From a8730c5330d1748f4e60b213a5d2786f4892a057 Mon Sep 17 00:00:00 2001 From: "jiyun.park" Date: Sat, 3 Oct 2026 11:31:41 +0200 Subject: [PATCH 01/11] Audio features: take the modulus of the complex CQT and STFT CQTProcessor, CQTSpectralFluxProcessor (OPHMM) and RawSpectrumProcessor (SKF) passed dtype=np.float32 to librosa.cqt / librosa.stft. librosa 0.11 then returns a real array, dropping the imaginary part, so np.abs gave |Re X| instead of |X| (about 36% mean magnitude error). Compute the complex transform and cast its modulus to float32. OPHMM's score templates already used the complex modulus. --- matchmaker/features/audio.py | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/matchmaker/features/audio.py b/matchmaker/features/audio.py index 704efb7..00ddd1b 100644 --- a/matchmaker/features/audio.py +++ b/matchmaker/features/audio.py @@ -116,11 +116,11 @@ def __call__( sr=self.sample_rate, hop_length=self.hop_length, norm=self.norm, - dtype=np.float32, fmin=librosa.note_to_hz("A0"), n_bins=88, ) - return np.abs(cqt).T[1:-1], f_time + # the modulus of the complex CQT (a real dtype would make librosa drop the imaginary part) + return np.abs(cqt).T.astype(np.float32)[1:-1], f_time class CQTSpectralFluxProcessor(Processor): @@ -162,9 +162,9 @@ def __call__( n_bins=self.n_bins, bins_per_octave=self.bins_per_octave, norm=self.norm, - dtype=np.float32, ) - cqt_features = np.abs(cqt).T + # the modulus of the complex CQT (a real dtype would make librosa drop the imaginary part) + cqt_features = np.abs(cqt).T.astype(np.float32) if self.include_spectral_flux: if self.prev_magnitude is None: @@ -323,9 +323,9 @@ def __call__(self, data: InputAudioFrame): win_length=self.n_fft, hop_length=self.hop_length, center=False, - dtype=np.float32, ) - return np.abs(stft).T, f_time # (n_frames, n_bins) + # the modulus of the complex STFT (a real dtype would make librosa drop the imaginary part) + return np.abs(stft).T.astype(np.float32), f_time # (n_frames, n_bins) def compute_features_from_audio( From 011b476a2887db1a66da1136cd8428a9098b580d Mon Sep 17 00:00:00 2001 From: "jiyun.park" Date: Sat, 3 Oct 2026 11:31:41 +0200 Subject: [PATCH 02/11] outerhmm: neighbourhood transitions and s_j of Nakamura et al. (2016) The audio OPHMM used the symmetric +-3 transition table of the MIDI model (Nakamura et al. 2014, Table 3, with the +1 entry mistyped as 0.94531 for 0.84531), filled the whole matrix with 1e-6 outside the band that the forward step uses, boosted the exit probability by spectral flux, and defaulted to s_j = 1e-5. Follow the audio paper instead (Sec. III-B, IV-B2): nbh(i) = {j; 0 <= i - j <= 2} with a_{i,i+2} = 1e-50 and no backward moves, no mass outside the band, no flux boost, and s_j = 1e-100 (Sec. IV-C). --- matchmaker/prob/outer_product_hmm_audio.py | 37 ++++++++-------------- 1 file changed, 14 insertions(+), 23 deletions(-) diff --git a/matchmaker/prob/outer_product_hmm_audio.py b/matchmaker/prob/outer_product_hmm_audio.py index f325e7c..12fffab 100644 --- a/matchmaker/prob/outer_product_hmm_audio.py +++ b/matchmaker/prob/outer_product_hmm_audio.py @@ -14,25 +14,21 @@ NDArrayFloat = NDArray[np.float32] NDArrayInt = NDArray[np.int32] -# Score-position jump probs (d = i - j) from Nakamura 2014 Table 3, |d| <= 3. -# d=0 (staying on the same note) is omitted: it is zeroed by fill_diagonal and -# instead derived per-chord from note duration and tempo (a00, set in __init__). +# Neighbourhood transitions of Nakamura et al. (2016), Sec. III-B and IV-B2: +# nbh(i) = {j; 0 <= i - j <= 2}, with a_{i,i+2} = 1e-50 for deletion errors and +# a_{i,i} = 0 for insertion errors (staying on a note is the bottom-level a00). DEFAULT_TRANSITIONS = [ - (-3, 0.00509), - (-2, 0.00516), - (-1, 0.00886), - (1, 0.94531), # normal forward progression - (2, 0.00610), # deletion (skip one note) - (3, 0.00073), + (1, 1.0 - 1e-50), # normal forward progression + (2, 1e-50), # deletion (skip one note) ] -DEFAULT_D1 = 3 -DEFAULT_D2 = 3 +# neighbourhood reach in the forward step: D1 states back, D2 states ahead +DEFAULT_D1 = 0 +DEFAULT_D2 = 2 IOI_THRESHOLD = 0.035 # seconds SUSTAINED_DECAY: float = 0.3 # exponential decay rate for sustained notes -_FLUX_EXIT_BOOST: float = 1.0 _OTHER_PROB: float = 1e-6 _EMISSION_BETA: float = 5.0 _UNIFORM_FLOOR: float = 0.2 @@ -209,16 +205,16 @@ def compute_transition_matrix( if transitions is None: transitions = DEFAULT_TRANSITIONS - # Initialize transition matrix with epsilons - alpha = np.full((N, N), 1e-6, dtype=float) + # only the band is used in the forward step; moves outside it are the s_j r_i term + alpha = np.zeros((N, N), dtype=float) for delta, prob in transitions: for i in range(N): j = i + delta if 0 <= j < N: alpha[i, j] = prob - alpha += np.finfo(float).eps - alpha /= alpha.sum(axis=1, keepdims=True) + row_sums = alpha.sum(axis=1, keepdims=True) + alpha = np.divide(alpha, row_sums, out=np.zeros_like(alpha), where=row_sums > 0) return alpha, D1, D2 @@ -231,7 +227,7 @@ def __init__( tempo: float = 120.0, sample_rate: int = 16000, hop_length: int = 320, - s_j: float = 1e-5, + s_j: float = 1e-100, r_i: Optional[np.ndarray] = None, **kwargs, ) -> None: @@ -511,12 +507,7 @@ def forward_step( ) emit_pause = np.full(N, emit_pause_scalar, dtype=float) - # Spectral-flux-driven exit boost - flux = float(obs[88]) if obs.size > 88 else 0.0 - f = flux / (flux + 1.0) # [0,1) - boost = 1.0 + _FLUX_EXIT_BOOST * f - e0 = np.clip(self.e0 * boost, 1e-10, 1.0 - self.a01 - 1e-10) - a00 = np.clip(1.0 - self.a01 - e0, 1e-10, 1.0 - 1e-10) + e0, a00 = self.e0, self.a00 # Exit masses from each top state j (Eq.(6)) exit_mass = prev_sound * e0 + prev_pause * self.e1 # (N,) From 7688c45292db6612b5688857f01d921405760ba4 Mon Sep 17 00:00:00 2001 From: "jiyun.park" Date: Sat, 3 Oct 2026 11:31:41 +0200 Subject: [PATCH 03/11] SKF chord lengths in quarter notes; score tempo in quarter-note BPM build_chord_sequence measured note lengths in score beats divided by 4, which is wrong whenever the beat is not a quarter note; use onset_quarter. The tempo read from a score marking now honours the marking's unit (e.g. half = 100 is 200 quarter-note BPM). --- matchmaker/prob/skf.py | 12 ++++++------ matchmaker/utils/misc.py | 4 +++- 2 files changed, 9 insertions(+), 7 deletions(-) diff --git a/matchmaker/prob/skf.py b/matchmaker/prob/skf.py index 4e104f0..4d0fe76 100644 --- a/matchmaker/prob/skf.py +++ b/matchmaker/prob/skf.py @@ -62,9 +62,10 @@ def build_chord_sequence(note_array: np.ndarray): """ onsets = note_array["onset_beat"] pitches = note_array["pitch"] - durations = note_array["duration_beat"] + # score time in quarters is independent of the time signature's beat unit + quarters = note_array["onset_quarter"] if "onset_quarter" in note_array.dtype.names else onsets - unique_onsets = np.unique(onsets) + unique_onsets, first = np.unique(onsets, return_index=True) K = len(unique_onsets) chords: List[List[int]] = [] @@ -75,11 +76,10 @@ def build_chord_sequence(note_array: np.ndarray): chords.append(pitches[mask].tolist()) onset_beats[k] = ub - # Nominal length l_k = onset_{k+1} - onset_k (in beat units, then /4 - # to convert to whole-note units since 1 whole note = 4 beats) + # Nominal length l_k = onset_{k+1} - onset_k in whole notes (4 quarters) + onset_quarters = quarters[first] lengths = np.zeros(K) - for k in range(K - 1): - lengths[k] = (onset_beats[k + 1] - onset_beats[k]) / 4.0 + lengths[:-1] = np.diff(onset_quarters) / 4.0 # Last chord: use same length as previous, or 1 beat lengths[K - 1] = lengths[K - 2] if K > 1 else 0.25 diff --git a/matchmaker/utils/misc.py b/matchmaker/utils/misc.py index 884442d..73640a7 100644 --- a/matchmaker/utils/misc.py +++ b/matchmaker/utils/misc.py @@ -15,6 +15,7 @@ import numpy as np import partitura from partitura.score import ScoreLike +from partitura.utils.music import to_quarter_tempo # Tempo marking to BPM mapping # Reference: https://en.wikipedia.org/wiki/Tempo#Basic_tempo_markings @@ -262,7 +263,8 @@ def get_tempo_from_score( try: for tempo_obj in score_part.iter_all(partitura.score.Tempo): if hasattr(tempo_obj, "bpm") and tempo_obj.bpm is not None: - return float(tempo_obj.bpm) + # the marking's unit ("h" = 100 is 200 quarters per minute) + return float(to_quarter_tempo(tempo_obj.unit or "q", tempo_obj.bpm)) except Exception: pass From c7b09a489e6f9123de05fdd8f104fcba9af01df0 Mon Sep 17 00:00:00 2001 From: "jiyun.park" Date: Sat, 3 Oct 2026 11:31:41 +0200 Subject: [PATCH 04/11] Add arzt_tempo and arzt_multi_fold baselines arzt_tempo: frame-level Arzt OLTW with the simple tempo model of Arzt & Widmer (SMC 2010), computed from a backward path (dtw_loop records backpointers). arzt_multi_fold: for folded scores, one Arzt matcher per unfolded version of the score, outputting the matcher with the lowest recent cost (adapted from Arzt, Widmer & Dixon, ECAI 2008, and Widmer & Arzt, Vienna Talk 2010). --- README.md | 2 + matchmaker/dp/__init__.py | 2 + matchmaker/dp/dtw_loop.pyx | 17 +- matchmaker/dp/oltw_arzt.py | 304 +++++++++++++++++++++- matchmaker/dp/oltw_arzt_multi_fold.py | 355 ++++++++++++++++++++++++++ matchmaker/methods.yaml | 31 +++ tests/test_registry.py | 4 +- tests/test_tempo_models.py | 201 +++++++++++++++ 8 files changed, 911 insertions(+), 5 deletions(-) create mode 100644 matchmaker/dp/oltw_arzt_multi_fold.py create mode 100644 tests/test_tempo_models.py diff --git a/README.md b/README.md index a181751..0822c1d 100644 --- a/README.md +++ b/README.md @@ -349,6 +349,8 @@ Default method: `"arzt"` | `"dixon"` | On-line time warping by Dixon (2005) | | `"outerhmm"` | Outer-product HMM score follower by Nakamura (2014) | | `"skf"` | Switching Kalman Filter with hidden tempo by Jiang and Raphael (2020) | +| `"arzt_tempo"` | On-line time warping by Arzt with a tempo model (Arzt and Widmer, 2010) | +| `"arzt_multi_fold"` | Arzt OLTW tracking several unfolded score variants in parallel | ### MIDI (`input_type="midi"`) diff --git a/matchmaker/dp/__init__.py b/matchmaker/dp/__init__.py index 4a063d3..7d7e8ed 100644 --- a/matchmaker/dp/__init__.py +++ b/matchmaker/dp/__init__.py @@ -12,9 +12,11 @@ OnlineTimeWarpingArzt, OnlineTimeWarpingArztEvent, OnlineTimeWarpingArztFrame, + OnlineTimeWarpingArztTempoFrame, ) from .oltw_dixon import ( OnlineTimeWarpingDixon, OnlineTimeWarpingDixonEvent, OnlineTimeWarpingDixonFrame, ) +from .oltw_arzt_multi_fold import OnlineTimeWarpingArztMultiFold diff --git a/matchmaker/dp/dtw_loop.pyx b/matchmaker/dp/dtw_loop.pyx index 4ec4379..b1f8242 100644 --- a/matchmaker/dp/dtw_loop.pyx +++ b/matchmaker/dp/dtw_loop.pyx @@ -23,7 +23,8 @@ def oltw_arzt_loop( int window_end, int input_index, float min_costs, - int min_index + int min_index, + signed char[:] backpointers = None ): """ Cythonized main loop for online time warping. @@ -46,6 +47,8 @@ def oltw_arzt_loop( Minimum costs of current window. min_index: int Index of minimum costs of current window. + backpointers: numpy.ndarray, optional + 1D array of shape (window_end - window_start,) to record predecessor steps. Returns ------- @@ -56,7 +59,7 @@ def oltw_arzt_loop( min_costs: float Minimum costs of current window. """ - cdef float dist1, dist2, dist3, local_dist, norm_cost + cdef float dist1, dist2, dist3, min_dist, local_dist, norm_cost cdef int idx = 0, score_index = window_start cdef Py_ssize_t k, wk = window_cost.shape[0] cdef Py_ssize_t N = global_cost_matrix.shape[0] @@ -84,6 +87,14 @@ def oltw_arzt_loop( min_dist = min(dist1, dist2, dist3) global_cost_matrix[score_index + 1, 1] = min_dist + if backpointers is not None: + if dist1 <= dist2 and dist1 <= dist3: + backpointers[idx] = 1 + elif dist2 <= dist3: + backpointers[idx] = 2 + else: + backpointers[idx] = 3 + norm_cost = min_dist / (input_index + score_index + 1.0) # check if new cell has lower costs and might be current position @@ -96,7 +107,7 @@ def oltw_arzt_loop( update_cost_matrix(global_cost_matrix, N) - return global_cost_matrix, min_index, min_costs + return np.asarray(global_cost_matrix), min_index, min_costs @cython.boundscheck(False) @cython.wraparound(False) diff --git a/matchmaker/dp/oltw_arzt.py b/matchmaker/dp/oltw_arzt.py index c74c9bd..1dc789d 100644 --- a/matchmaker/dp/oltw_arzt.py +++ b/matchmaker/dp/oltw_arzt.py @@ -84,7 +84,6 @@ def __init__( ) self.N_ref: int = self.reference_features.shape[0] self.step_size: int = step_size - self.reset() def reset(self) -> None: self.current_index = 0 @@ -165,6 +164,7 @@ def __init__( "min_latency": float("inf"), } self._init_distance_func(distance_func) + self.reset() def __call__(self, observation: Any, perf_time: float) -> float: t0 = time.time() @@ -278,6 +278,307 @@ def step(self, input_features: NDArray[np.float32]) -> None: self.input_index += 1 +class OnlineTimeWarpingArztTempoFrame(OnlineTimeWarpingArztFrame): + """Frame-level OLTW with Arzt & Widmer (2010) Simple Tempo Model. + + Dynamically tracks relative tempo via rectified backward path + and stretches/compresses reference score features on the fly. + """ + + def __init__( + self, + reference_features: NDArray[np.float32], + score_positions: NDArray[np.float32], + window_size: int = WINDOW_SIZE, + step_size: int = STEP_SIZE, + distance_func: Union[ + str, Callable, Tuple[str, Dict[str, Any]] + ] = OnlineTimeWarpingArztFrame.DEFAULT_DISTANCE_FUNC, + start_window_size: Union[float, int] = START_WINDOW_SIZE, + frame_rate: int = FRAME_RATE, + ref_frame_to_beat: NDArray = None, + queue: Optional[RECVQueue] = None, + tempo_window_size: float = 3.0, + tempo_min_past: float = 1.0, + tempo_n_onsets: int = 20, + random_seed: Optional[int] = 1984, + use_tempo_model: bool = True, + **kwargs, + ) -> None: + if ref_frame_to_beat is None: + raise ValueError( + "Frame-level Arzt requires `ref_frame_to_beat` (per-frame beat mapping)." + ) + self.use_tempo_model = use_tempo_model + self.tempo_window_size = tempo_window_size + self.tempo_min_past = tempo_min_past + self.tempo_n_onsets = tempo_n_onsets + self.random_seed = random_seed + self._original_reference_features = reference_features.copy() + self._original_ref_frame_to_beat = ref_frame_to_beat.copy() + super().__init__( + reference_features=reference_features, + score_positions=score_positions, + window_size=window_size, + step_size=step_size, + distance_func=distance_func, + start_window_size=start_window_size, + frame_rate=frame_rate, + ref_frame_to_beat=ref_frame_to_beat, + queue=queue, + **kwargs, + ) + self.reset() + + def reset(self) -> None: + super().reset() + self.current_relative_tempo: float = 1.0 + self._consecutive_alterations: int = 0 + self._backpointers: list = [] + self._rng = np.random.RandomState(self.random_seed) + # The tempo model stretches the reference on the fly; restore the + # pristine copy so each run starts from the notated score. + self.reference_features = self._original_reference_features.copy() + self._ref_frame_to_beat = self._original_ref_frame_to_beat.copy() + self.N_ref = self.reference_features.shape[0] + self.global_cost_matrix = np.full( + (self.N_ref + 1, 2), np.inf, dtype=np.float32 + ) + self._init_onset_mask() + + def _init_onset_mask(self) -> None: + self.is_onset_frame = np.zeros(self.N_ref, dtype=bool) + onset_frames = np.searchsorted(self._ref_frame_to_beat, self.score_positions) + valid = onset_frames < self.N_ref + self.is_onset_frame[onset_frames[valid]] = True + + def _get_backward_path(self, max_history_frames: int = 500) -> list: + """Trace backward path from current position using recorded backpointers.""" + if not self._backpointers or self.input_index <= 0: + return [(self._current_frame, 0)] + + path = [] + s = self._current_frame + i = self.input_index - 1 + min_i = max(0, i - max_history_frames) + while i >= min_i and s >= 0: + path.append((s, i)) + if i >= len(self._backpointers): + i -= 1 + s = max(0, s - 1) + continue + w_start, w_end, bp = self._backpointers[i] + if bp is not None and w_start <= s < w_end: + code = bp[s - w_start] + if code == 1: + s -= 1 + elif code == 2: + i -= 1 + elif code == 3: + s -= 1 + i -= 1 + else: + break + else: + s = max(0, s - 1) + i -= 1 + path.reverse() + return path + + def _compute_relative_tempo(self) -> float: + """Estimate current relative tempo via rectified backward path (Arzt 2010 Sec. 4.1).""" + if self.input_index < int(self.frame_rate * 1.5): + return 1.0 + + cur_perf_time = self.input_index / self.frame_rate + max_history = int(self.frame_rate * 10) + path = self._get_backward_path(max_history_frames=max_history) + if len(path) < 2: + return 1.0 + + score_frames = np.array([p[0] for p in path]) + perf_frames = np.array([p[1] for p in path]) + + s_min, s_max = int(score_frames[0]), int(score_frames[-1]) + if s_max <= s_min: + return 1.0 + + local_onsets = np.where(self.is_onset_frame[s_min : s_max + 1])[0] + s_min + if len(local_onsets) < 2: + return 1.0 + + t_s_list = [] + t_p_list = [] + for sf in local_onsets: + idx = int(np.searchsorted(score_frames, sf)) + if idx < len(score_frames): + t_s_list.append(sf / self.frame_rate) + t_p_list.append(perf_frames[idx] / self.frame_rate) + + if len(t_s_list) < 2: + return 1.0 + + t_s_arr = np.array(t_s_list) + t_p_arr = np.array(t_p_list) + + # Onsets at least tempo_min_past seconds in the past + valid_mask = t_p_arr <= (cur_perf_time - self.tempo_min_past) + valid_indices = np.where(valid_mask)[0] + if len(valid_indices) == 0: + return 1.0 + + if len(valid_indices) > self.tempo_n_onsets: + valid_indices = valid_indices[-self.tempo_n_onsets:] + + half_win = self.tempo_window_size / 2.0 + tempi = [] + for idx in valid_indices: + ts = t_s_arr[idx] + w_start = ts - half_win + w_end = ts + half_win + if w_start >= t_s_arr[0] and w_end <= t_s_arr[-1]: + tp_start = float(np.interp(w_start, t_s_arr, t_p_arr)) + tp_end = float(np.interp(w_end, t_s_arr, t_p_arr)) + delta_p = tp_end - tp_start + delta_s = w_end - w_start + if delta_p > 1e-4: + tempi.append(delta_s / delta_p) + else: + tempi.append(1.0) + + if not tempi: + return 1.0 + + m = len(tempi) + weights = np.arange(1, m + 1) + t_est = float(np.sum(np.array(tempi) * weights) / np.sum(weights)) + return float(np.clip(t_est, 0.33, 3.0)) + + def _alter_score_representation(self, t: float) -> None: + """Alter score representation on the fly according to relative tempo (Arzt 2010 Sec. 4.2).""" + ps = self._current_frame + if ps + 2 >= self.N_ref: + return + + if self._consecutive_alterations >= 3: + self._consecutive_alterations = 0 + return + + r = self._rng.uniform(0.0, 1.0) + + if t > 1.0: + threshold = 1.0 / t + if r > threshold: + if self.is_onset_frame[ps + 1] or self.is_onset_frame[ps + 2]: + self._consecutive_alterations = 0 + return + + mean_feat = 0.5 * ( + self.reference_features[ps + 1] + self.reference_features[ps + 2] + ) + self.reference_features[ps + 1] = mean_feat + self.reference_features = np.delete( + self.reference_features, ps + 2, axis=0 + ) + + mean_beat = 0.5 * ( + self._ref_frame_to_beat[ps + 1] + self._ref_frame_to_beat[ps + 2] + ) + self._ref_frame_to_beat[ps + 1] = mean_beat + self._ref_frame_to_beat = np.delete(self._ref_frame_to_beat, ps + 2) + + self.is_onset_frame = np.delete(self.is_onset_frame, ps + 2) + self.global_cost_matrix = np.delete( + self.global_cost_matrix, ps + 3, axis=0 + ) + self.N_ref = self.reference_features.shape[0] + self._consecutive_alterations += 1 + else: + self._consecutive_alterations = 0 + + elif t < 1.0: + threshold = t + if r > threshold: + if self.is_onset_frame[ps] or self.is_onset_frame[ps + 1]: + self._consecutive_alterations = 0 + return + + mean_feat = 0.5 * ( + self.reference_features[ps] + self.reference_features[ps + 1] + ) + self.reference_features = np.insert( + self.reference_features, ps + 1, mean_feat, axis=0 + ) + + mean_beat = 0.5 * ( + self._ref_frame_to_beat[ps] + self._ref_frame_to_beat[ps + 1] + ) + self._ref_frame_to_beat = np.insert( + self._ref_frame_to_beat, ps + 1, mean_beat + ) + + self.is_onset_frame = np.insert(self.is_onset_frame, ps + 1, False) + row1 = np.asarray(self.global_cost_matrix[ps + 1]) + row2 = np.asarray(self.global_cost_matrix[ps + 2]) + new_cost_row = 0.5 * (row1 + row2) + new_cost_row[1] = np.inf + self.global_cost_matrix = np.insert( + self.global_cost_matrix, ps + 2, new_cost_row, axis=0 + ) + self.N_ref = self.reference_features.shape[0] + self._consecutive_alterations += 1 + else: + self._consecutive_alterations = 0 + else: + self._consecutive_alterations = 0 + + def step(self, input_features: NDArray[np.float32]) -> None: + if self.use_tempo_model and self.input_index > 0: + self.current_relative_tempo = self._compute_relative_tempo() + self._alter_score_representation(self.current_relative_tempo) + + min_costs = np.inf + min_index = max(self.window_index - self.step_size, 0) + + window_start, window_end = self.get_window() + window_cost = self.vdist( + self.reference_features[window_start:window_end], + input_features.squeeze(), + self.distance_func, + ) + + bp = None + if self.use_tempo_model: + bp = np.zeros(window_end - window_start, dtype=np.int8) + + self.global_cost_matrix, min_index, min_costs = oltw_arzt_loop( + global_cost_matrix=self.global_cost_matrix, + window_cost=window_cost, + window_start=window_start, + window_end=window_end, + input_index=self.input_index, + min_costs=min_costs, + min_index=min_index, + backpointers=bp, + ) + + if self.use_tempo_model: + self._backpointers.append((window_start, window_end, bp)) + + if self.input_index > 0: + self._current_frame = int( + np.clip( + min_index, + self._current_frame, + self._current_frame + self.step_size, + ) + ) + else: + self._current_frame = min_index + self.current_index = self._frame_to_score_idx(self._current_frame) + self.input_index += 1 + + # --------------------------------------------------------------------------- # Event-level subclass (MIDI) # --------------------------------------------------------------------------- @@ -313,6 +614,7 @@ def __init__( self._window_size = window_size self._start_window_size = start_window_size self._scipy_metric = distance_func + self.reset() def get_window(self) -> Tuple[int, int]: w = self._window_size diff --git a/matchmaker/dp/oltw_arzt_multi_fold.py b/matchmaker/dp/oltw_arzt_multi_fold.py new file mode 100644 index 0000000..9dfabf7 --- /dev/null +++ b/matchmaker/dp/oltw_arzt_multi_fold.py @@ -0,0 +1,355 @@ +#!/usr/bin/python +# -*- coding: utf-8 -*- +""" +Parallel Multi-Fold Online Time Warping (OLTWArztMultiFold) + +A baseline for folded scores with repeats. The score is unfolded into a few versions +(every repeat taken, the folded score itself, and the first four versions partitura's +``iter_unfolded_parts`` yields, duplicates removed); one frame-level Arzt OLTW matcher +(``OnlineTimeWarpingArztFrame``) runs on each version from the start, and the output, in +notated (folded) beats, is the matcher whose local matching cost averaged over the last +``selection_window`` frames is lowest, switching only when another is lower by +``switch_margin``. + +The idea of running several OLTW matchers in parallel and comparing their alignment costs +is adapted from + A. Arzt, G. Widmer and S. Dixon, "Automatic page turning for musicians via real-time + machine listening", Proc. ECAI, 2008 (Strategy 3) + G. Widmer and A. Arzt, "Robust real-time music tracking", Proc. Vienna Talk, 2010. +The unfolded versions, the cost window and the switching margin are this implementation's +own choices, not those papers'. +""" + +from __future__ import annotations + +import collections +import copy +import re +import time +from typing import Any, Callable, Dict, Generator, List, Optional, Tuple, Union + +import numpy as np +import partitura as pt +from numpy.typing import NDArray +from partitura.score import Part, merge_parts + +from matchmaker.base import OnlineAlignment +from matchmaker.dp.oltw_arzt import ( + OnlineTimeWarpingArztFrame, + START_WINDOW_SIZE, + STEP_SIZE, + WINDOW_SIZE, +) +from matchmaker.features.audio import FRAME_RATE +from matchmaker.features.processor import Processor +from matchmaker.io.audio import QUEUE_TIMEOUT +from matchmaker.io.queue import RECVQueue +from matchmaker.utils.misc import generate_score_audio, set_latency_stats +from partitura.io.exportmidi import get_ppq + + +class OnlineTimeWarpingArztMultiFold(OnlineTimeWarpingArztFrame): + """Multi-fold parallel OLTW tracker for folded scores with repeats. + + Parameters + ---------- + score_part : Part + The original (folded) partitura Part. + tempo : float + Nominal score tempo in BPM. + sample_rate : int + Audio sample rate. + processor : Processor + Audio feature processor (e.g. ChromaProcessor). + frame_rate : int + Audio frame rate (default: FRAME_RATE, 30). + window_size : int + Search window size in seconds (default: 10). + step_size : int + Maximum forward advance per input frame (default: 3). + start_window_size : float or int + Window size during warmup (default: 0.1). + selection_window : int + Number of trailing frames over which running matching cost is averaged + for hypothesis selection (default: 50 frames, about 1.7 s at 30 fps). + switch_margin : float + Hysteresis factor requiring an alternative hypothesis to have cost + strictly lower by this fraction before switching (default: 0.05). + queue : Optional[RECVQueue] + Input queue for streaming observations. + """ + + def __init__( + self, + reference_features: Optional[NDArray[np.float32]] = None, + score_positions: Optional[NDArray[np.float32]] = None, + score_part: Optional[Part] = None, + tempo: Optional[float] = None, + sample_rate: int = 44100, + processor: Optional[Processor] = None, + frame_rate: int = FRAME_RATE, + window_size: int = WINDOW_SIZE, + step_size: int = STEP_SIZE, + start_window_size: Union[float, int] = START_WINDOW_SIZE, + distance_func: str = "Manhattan", + selection_window: int = 50, + switch_margin: float = 0.05, + queue: Optional[RECVQueue] = None, + ref_frame_to_beat: Optional[NDArray] = None, + **kwargs, + ) -> None: + if score_part is None: + raise ValueError("OnlineTimeWarpingArztMultiFold requires `score_part`.") + if tempo is None: + tempo = 120.0 + + self.score_part = score_part + self.tempo = float(tempo) + self.sample_rate = sample_rate + self.processor = processor + self.frame_rate = frame_rate + self.window_size = window_size + self.step_size = step_size + self.start_window_size = start_window_size + self.distance_func_name = distance_func + self.selection_window = selection_window + self.switch_margin = switch_margin + self.queue_timeout = QUEUE_TIMEOUT + + # Folded note lookup: original note ID -> notated beat + folded_notes = self.score_part.note_array() + self.folded_note_dict = { + str(n["id"]): float(n["onset_beat"]) for n in folded_notes + } + notated_positions = np.unique(folded_notes["onset_beat"]).astype(np.float32) + + ref_feat = ( + reference_features + if reference_features is not None + else np.empty((0, 12), dtype=np.float32) + ) + dummy_ref_beat = ( + ref_frame_to_beat + if ref_frame_to_beat is not None + else np.zeros(max(1, len(ref_feat)), dtype=np.float32) + ) + + # Build parallel track hypotheses from the score repeat structure + # before the parent __init__ so that reset() always finds them. + self._build_hypotheses() + + super().__init__( + reference_features=ref_feat, + score_positions=notated_positions, + window_size=window_size, + step_size=step_size, + start_window_size=start_window_size, + distance_func=distance_func, + frame_rate=frame_rate, + ref_frame_to_beat=dummy_ref_beat, + queue=queue, + **kwargs, + ) + self.reset() + + def _build_hypotheses(self) -> None: + """Enumerate distinct structural fold clones safely and initialize trackers.""" + import itertools + import sys + + _prev_limit = sys.getrecursionlimit() + clones = [] + try: + sys.setrecursionlimit(max(_prev_limit, 25_000)) + # 1. Maximal unfold (canonical repeat-performed performance) + try: + max_unfolded = pt.score.unfold_part_maximal( + self.score_part, ignore_leaps=False + ) + clones.append( + merge_parts(max_unfolded.parts) + if hasattr(max_unfolded, "parts") + else max_unfolded + ) + except Exception: + pass + + # 2. Minimal / folded score (repeat-skipped performance) + clones.append(self.score_part) + + # 3. Intermediate unfold variants (bounded to prevent combinatorial explosion) + try: + gen = pt.score.iter_unfolded_parts(self.score_part, update_ids=True) + for c in itertools.islice(gen, 4): + clones.append(c) + except Exception: + pass + finally: + sys.setrecursionlimit(_prev_limit) + + # Deduplicate clones by note count and duration + unique_clones = [] + seen = set() + for c in clones: + na = c.note_array() + key = (len(na), round(float(na["onset_beat"].max()), 2)) + if key not in seen: + seen.add(key) + unique_clones.append(c) + + if not unique_clones: + unique_clones = [self.score_part] + + self.hypotheses: List[Dict[str, Any]] = [] + tick = get_ppq(self.score_part) + + for i, clone in enumerate(unique_clones): + # 1. Synthesize audio and extract reference features + clone_audio = generate_score_audio(clone, self.tempo, self.sample_rate) + if self.processor is not None: + features, _ = self.processor((clone_audio.astype(np.float32), 0.0)) + self.processor.reset() + else: + raise ValueError( + "OnlineTimeWarpingArztMultiFold requires an audio processor." + ) + + n_ref_frames = features.shape[0] + + # 2. Frame to unfolded beat + times = ( + (np.arange(n_ref_frames) / self.frame_rate) + * tick + * (self.tempo / 60.0) + ) + ref_frame_to_unfolded_beat = np.array( + [float(clone.beat_map(t)) for t in times], dtype=np.float32 + ) + score_positions_k = np.unique(clone.note_array()["onset_beat"]).astype( + np.float32 + ) + + # 3. Unfolded beat -> Notated beat mapping via note ID correspondence + clone_notes = clone.note_array() + u_beats = [] + n_beats = [] + for n in clone_notes: + orig_id = re.sub(r"-\d+$", "", str(n["id"])) + if orig_id in self.folded_note_dict: + u_beats.append(float(n["onset_beat"])) + n_beats.append(self.folded_note_dict[orig_id]) + + u_beats = np.asarray(u_beats, dtype=float) + n_beats = np.asarray(n_beats, dtype=float) + sort_idx = np.argsort(u_beats) + u_sorted = u_beats[sort_idx] + n_sorted = n_beats[sort_idx] + u_uniq, u_first = np.unique(u_sorted, return_index=True) + n_uniq = n_sorted[u_first] + + # Precompute frame to notated beat directly + if len(n_uniq) > 0: + ref_frame_to_notated_beat = np.interp( + ref_frame_to_unfolded_beat, + u_uniq, + n_uniq, + left=n_uniq[0], + right=n_uniq[-1], + ).astype(np.float32) + else: + ref_frame_to_notated_beat = ref_frame_to_unfolded_beat.astype(np.float32) + + # 4. Instantiate follower + tracker = OnlineTimeWarpingArztFrame( + reference_features=features, + score_positions=score_positions_k, + window_size=self.window_size, + step_size=self.step_size, + start_window_size=self.start_window_size, + distance_func=self.distance_func_name, + frame_rate=self.frame_rate, + ref_frame_to_beat=ref_frame_to_unfolded_beat, + queue=None, + ) + + self.hypotheses.append( + { + "index": i, + "clone": clone, + "tracker": tracker, + "ref_frame_to_notated_beat": ref_frame_to_notated_beat, + "recent_costs": collections.deque( + maxlen=self.selection_window + ), + "mean_recent_cost": float("inf"), + } + ) + + def reset(self) -> None: + super().reset() + self.active_hypothesis_idx = 0 + self.latency_stats = { + "total_latency": 0, + "total_frames": 0, + "max_latency": 0, + "min_latency": float("inf"), + } + for hyp in self.hypotheses: + hyp["tracker"].reset() + hyp["recent_costs"].clear() + hyp["mean_recent_cost"] = float("inf") + + def __call__(self, observation: Any, perf_time: float) -> float: + t0 = time.time() + self.step(observation) + current_beat = self.get_current_position() + self._alignment_path.append((float(perf_time), float(current_beat))) + self.latency_stats = set_latency_stats( + time.time() - t0, self.latency_stats, self.input_index + ) + return current_beat + + def step(self, input_features: NDArray[np.float32]) -> None: + obs_squeezed = input_features.squeeze() + + # 1. Advance every hypothesis tracker + for i, hyp in enumerate(self.hypotheses): + tracker: OnlineTimeWarpingArztFrame = hyp["tracker"] + tracker.step(input_features) + + # Measure local matching error at current hypothesis frame + curr_frame = tracker._current_frame + ref_feat = tracker.reference_features[curr_frame] + loc_cost = float(tracker.distance_func(ref_feat, obs_squeezed)) + + hyp["recent_costs"].append(loc_cost) + mean_cost = float(np.mean(hyp["recent_costs"])) + hyp["mean_recent_cost"] = mean_cost + + # 2. Decision Maker: select winning hypothesis with hysteresis + current_hyp_cost = self.hypotheses[self.active_hypothesis_idx][ + "mean_recent_cost" + ] + min_hyp_idx = min( + range(len(self.hypotheses)), + key=lambda i: self.hypotheses[i]["mean_recent_cost"], + ) + min_hyp_cost = self.hypotheses[min_hyp_idx]["mean_recent_cost"] + + if min_hyp_cost < current_hyp_cost * (1.0 - self.switch_margin): + self.active_hypothesis_idx = min_hyp_idx + + self.input_index += 1 + + def get_current_position(self) -> float: + """Return the current position projected into notated (folded) score beats.""" + winning_hyp = self.hypotheses[self.active_hypothesis_idx] + curr_frame = winning_hyp["tracker"]._current_frame + notated_beat = winning_hyp["ref_frame_to_notated_beat"][ + min(curr_frame, len(winning_hyp["ref_frame_to_notated_beat"]) - 1) + ] + return float(notated_beat) + + def run(self, verbose: bool = True) -> Generator[float, None, NDArray]: + self.reset() + return (yield from super().run(verbose=verbose)) diff --git a/matchmaker/methods.yaml b/matchmaker/methods.yaml index 9b96d1d..cabe986 100644 --- a/matchmaker/methods.yaml +++ b/matchmaker/methods.yaml @@ -169,6 +169,37 @@ methods: start_window_size: 0.1 step_size: 3 + arzt_tempo: + family: oltw + class: matchmaker.dp:OnlineTimeWarpingArztTempoFrame + reference: score_audio + args: *oltw_frame_args + config_passthrough: true + default_kwargs: + window_size: 10 + start_window_size: 0.1 + step_size: 3 + use_tempo_model: true + random_seed: 1984 + + arzt_multi_fold: + family: oltw + class: matchmaker.dp:OnlineTimeWarpingArztMultiFold + reference: score_audio + args: + <<: *oltw_frame_args + score_part: {from: score_part} + tempo: {from: tempo} + sample_rate: {from: sample_rate} + processor: {from: processor} + config_passthrough: true + default_kwargs: + window_size: 10 + start_window_size: 0.1 + step_size: 3 + selection_window: 50 + switch_margin: 0.05 + dixon: family: oltw class: matchmaker.dp:OnlineTimeWarpingDixonFrame diff --git a/tests/test_registry.py b/tests/test_registry.py index 629f627..ef88f08 100644 --- a/tests/test_registry.py +++ b/tests/test_registry.py @@ -97,7 +97,9 @@ def test_default_method_and_processor_exist(self): ) def test_families_are_populated(self): - self.assertEqual(OLTW_METHODS, {"arzt", "dixon"}) + self.assertEqual( + OLTW_METHODS, {"arzt", "dixon", "arzt_tempo", "arzt_multi_fold"} + ) self.assertEqual(PARANGONAR_METHODS, {"OPTM", "OTM", "SLT_OLTW", "SL_OLTW"}) def test_obj_tag_resolves_to_the_real_object(self): diff --git a/tests/test_tempo_models.py b/tests/test_tempo_models.py new file mode 100644 index 0000000..131683a --- /dev/null +++ b/tests/test_tempo_models.py @@ -0,0 +1,201 @@ +#!/usr/bin/env python +# -*- coding: utf-8 -*- +""" +Tests for Arzt (2010) Simple Tempo Model and Raphael/Jiang (2020) SKF. +""" + +import unittest +import numpy as np + +from matchmaker.dp.oltw_arzt import ( + OnlineTimeWarpingArztFrame, + OnlineTimeWarpingArztTempoFrame, +) +from matchmaker.prob.skf import ( + SwitchingKalmanFilterFollower, + build_chord_sequence, + build_spectral_templates, +) +from matchmaker.registry import REGISTRY +from tests.utils import generate_example_sequences + +RNG = np.random.RandomState(1984) + + +class TestArztTempoModel(unittest.TestCase): + """Tests for Arzt & Widmer (2010 SMC) Simple Tempo Model.""" + + def setUp(self): + # Generate synthetic reference and performance sequences + self.frame_rate = 50 + # 10 seconds of score audio -> 500 frames + n_frames = 500 + n_features = 12 + self.X = RNG.randn(n_frames, n_features).astype(np.float32) + # 1 beat per second -> 10 beats total + # Every 50 frames is 1 beat + self.ref_frame_to_beat = (np.arange(n_frames, dtype=np.float32) / self.frame_rate) + # Score positions at beats 0.0, 1.0, 2.0, ..., 9.0 + self.score_positions = np.arange(10, dtype=np.float32) + + def test_initialization(self): + tracker = OnlineTimeWarpingArztTempoFrame( + reference_features=self.X, + score_positions=self.score_positions, + ref_frame_to_beat=self.ref_frame_to_beat, + frame_rate=self.frame_rate, + window_size=5, + step_size=3, + ) + self.assertTrue(tracker.use_tempo_model) + self.assertEqual(tracker.current_relative_tempo, 1.0) + self.assertEqual(tracker._consecutive_alterations, 0) + self.assertEqual(len(tracker._backpointers), 0) + self.assertTrue(tracker.is_onset_frame[0]) + self.assertTrue(tracker.is_onset_frame[50]) + + def test_run_steps(self): + tracker = OnlineTimeWarpingArztTempoFrame( + reference_features=self.X, + score_positions=self.score_positions, + ref_frame_to_beat=self.ref_frame_to_beat, + frame_rate=self.frame_rate, + window_size=5, + step_size=3, + ) + + # Feed 100 frames of observations + for i in range(100): + obs = self.X[min(i, len(self.X) - 1)] + RNG.randn(12).astype(np.float32) * 0.01 + perf_time = i / self.frame_rate + beat = tracker(obs, perf_time) + self.assertIsInstance(beat, float) + self.assertGreaterEqual(beat, 0.0) + + # Check that backpointers were recorded + self.assertEqual(len(tracker._backpointers), 100) + path = tracker.alignment_path + self.assertEqual(path.shape[0], 2) + self.assertEqual(path.shape[1], 100) + + def test_reset(self): + tracker = OnlineTimeWarpingArztTempoFrame( + reference_features=self.X, + score_positions=self.score_positions, + ref_frame_to_beat=self.ref_frame_to_beat, + frame_rate=self.frame_rate, + ) + orig_len = tracker.N_ref + + # Run 50 frames + for i in range(50): + obs = self.X[i] + tracker(obs, i / self.frame_rate) + + # Reset + tracker.reset() + self.assertEqual(tracker.current_index, 0) + self.assertEqual(tracker._current_frame, 0) + self.assertEqual(tracker.current_relative_tempo, 1.0) + self.assertEqual(len(tracker._backpointers), 0) + self.assertEqual(tracker.N_ref, orig_len) + + def test_registry_lookup(self): + spec = REGISTRY.method("audio", "arzt_tempo") + self.assertEqual(spec.cls_path, "matchmaker.dp:OnlineTimeWarpingArztTempoFrame") + self.assertTrue(spec.default_kwargs.get("use_tempo_model")) + + +class TestRaphaelSwitchingStateSpace(unittest.TestCase): + """Tests for Raphael & Jiang (2020 ISMIR) Switching State-Space Model.""" + + def setUp(self): + # Create a simple synthetic note_array + self.note_array = np.zeros( + 4, + dtype=[ + ("pitch", "i4"), + ("onset_beat", "f4"), + ("duration_beat", "f4"), + ], + ) + self.note_array["pitch"] = [60, 64, 67, 72] # C major chord arpeggio + self.note_array["onset_beat"] = [0.0, 1.0, 2.0, 3.0] + self.note_array["duration_beat"] = [1.0, 1.0, 1.0, 1.0] + + def test_build_chord_sequence(self): + chords, lengths, onset_beats = build_chord_sequence(self.note_array) + self.assertEqual(len(chords), 4) + self.assertEqual(chords[0], [60]) + self.assertEqual(chords[1], [64]) + np.testing.assert_array_equal(onset_beats, [0.0, 1.0, 2.0, 3.0]) + # 1 beat = 0.25 whole note + np.testing.assert_allclose(lengths, [0.25, 0.25, 0.25, 0.25]) + + def test_build_spectral_templates(self): + chords, _, _ = build_chord_sequence(self.note_array) + templates = build_spectral_templates(chords, n_fft=512, sample_rate=8000) + self.assertEqual(templates.shape[0], 4) + self.assertEqual(templates.shape[1], 257) + # Each template row should sum to 1 + np.testing.assert_allclose(templates.sum(axis=1), np.ones(4), atol=1e-5) + + def test_initialization_and_step(self): + follower = SwitchingKalmanFilterFollower( + reference_features=self.note_array, + tempo=120.0, + sample_rate=8000, + n_fft=512, + hop_length=128, + max_hypotheses=50, + ) + self.assertEqual(follower.K, 4) + self.assertEqual(follower.current_index, 0) + # Initial tempo: 240 / 120 = 2.0 s/whole-note + self.assertAlmostEqual(follower.init_tempo, 2.0) + + # Process 20 frames + delta = 128 / 8000.0 + for i in range(20): + dummy_spectrum = np.random.rand(257).astype(np.float32) + pos = follower(dummy_spectrum, perf_time=i * delta) + self.assertIsInstance(pos, float) + self.assertGreaterEqual(pos, 0.0) + self.assertLessEqual(pos, 4.0) + + # Verify hypotheses count is bounded by max_hypotheses + self.assertLessEqual(len(follower.hypotheses), 50) + path = follower.alignment_path + self.assertEqual(path.shape, (2, 20)) + + def test_kalman_tempo_adaptation(self): + follower = SwitchingKalmanFilterFollower( + reference_features=self.note_array, + tempo=120.0, + sample_rate=8000, + n_fft=512, + hop_length=128, + ) + # Test transition Kalman update directly: + # Expected duration = length * tempo = 0.25 * 2.0 = 0.5s. + # Suppose performer played chord 0 in 0.25s (twice as fast, tempo should decrease s/wn): + fast_age = int(0.25 / follower.delta) + mu_new, var_new = follower._kalman_update( + mu_t_prev=2.0, var_t_prev=0.04, k_prev=0, age_prev=fast_age + ) + self.assertLess(mu_new, 2.0) + + # Suppose performer played chord 0 in 1.0s (twice as slow, tempo should increase s/wn): + slow_age = int(1.0 / follower.delta) + mu_new_slow, var_new_slow = follower._kalman_update( + mu_t_prev=2.0, var_t_prev=0.04, k_prev=0, age_prev=slow_age + ) + self.assertGreater(mu_new_slow, 2.0) + + def test_registry_lookup(self): + spec = REGISTRY.method("audio", "skf") + self.assertEqual(spec.cls_path, "matchmaker.prob.skf:SwitchingKalmanFilterFollower") + + +if __name__ == "__main__": + unittest.main() From 770b6aa62d1b398005b6891c57187dbf5f46355f Mon Sep 17 00:00:00 2001 From: "jiyun.park" Date: Sat, 3 Oct 2026 11:31:41 +0200 Subject: [PATCH 05/11] run_examples: choose the example piece, custom files, unfolding and method kwargs --- run_examples.py | 125 ++++++++++++++++++++++++++++++++++++++++-------- 1 file changed, 105 insertions(+), 20 deletions(-) diff --git a/run_examples.py b/run_examples.py index 5665621..23cf9ad 100644 --- a/run_examples.py +++ b/run_examples.py @@ -87,16 +87,31 @@ def plot_alignment_path(perf_sec, score_beat, gt_perf, gt_score, save_path, run_ plt.close(fig) -def save_results(save_dir, run_name, results, perf_sec, score_beat, gt_perf, gt_score): - """Save wp/gt TSVs, results JSON, and the alignment path plot.""" +def save_results( + save_dir, + run_name, + results, + perf_sec, + score_beat, + gt_perf, + gt_score, + save_plots=True, +): + """Save predicted warping path (wp) and ground truth (gt) TSVs, results JSON, and alignment plot.""" save_dir.mkdir(parents=True, exist_ok=True) save_tsv(np.column_stack([perf_sec, score_beat]), save_dir / f"wp_{run_name}.tsv") save_tsv(np.column_stack([gt_perf, gt_score]), save_dir / f"gt_{run_name}.tsv") with open(save_dir / f"{run_name}.json", "w") as f: json.dump(results, f, indent=4) - plot_alignment_path( - perf_sec, score_beat, gt_perf, gt_score, save_dir / f"{run_name}.png", run_name - ) + if save_plots: + plot_alignment_path( + perf_sec, + score_beat, + gt_perf, + gt_score, + save_dir / f"{run_name}.png", + run_name, + ) def main(): @@ -105,18 +120,73 @@ def main(): group = parser.add_mutually_exclusive_group() group.add_argument("--audio", action="store_true", help="Use audio input mode") group.add_argument("--midi", action="store_true", help="Use MIDI input mode") + parser.add_argument( + "--piece", + type=str, + default="simple_mozart", + choices=list(EXAMPLE_PIECES.keys()), + help="Built-in example piece to run (default: simple_mozart)", + ) parser.add_argument( "--method", type=str, default=None, - help="Score following method (e.g., arzt, dixon, outerhmm, audio_outerhmm)", + help="Score following method (e.g., arzt, arzt_tempo, dixon, outerhmm)", + ) + parser.add_argument( + "--score", type=str, default=None, help="Path to custom score XML/MusicXML" + ) + parser.add_argument( + "--audio-file", + type=str, + default=None, + help="Path to custom performance audio file", + ) + parser.add_argument( + "--midi-file", + type=str, + default=None, + help="Path to custom performance MIDI file", + ) + parser.add_argument( + "--match", + type=str, + default=None, + help="Path to custom .match ground truth file", + ) + parser.add_argument( + "--unfold", action="store_true", help="Unfold score repetitions during load" + ) + parser.add_argument("--no-plots", action="store_true", help="Skip alignment plots") + parser.add_argument( + "--output-dir", type=Path, default=ROOT_DIR / "results", help="Result directory" + ) + parser.add_argument( + "--kwargs", type=json.loads, default={}, help="Method option overrides" ) args = parser.parse_args() input_mode = "midi" if args.midi else "audio" - performance_file = select_performance_file(input_mode) - print(f"Running matchmaker with the score file ({SCORE_FILE.name})...") + # Resolve piece files + if args.score: + score_file = Path(args.score) + performance_file = Path( + args.midi_file if input_mode == "midi" else args.audio_file + ) + match_file = Path(args.match) if args.match else None + run_name = score_file.stem + else: + piece_cfg = EXAMPLE_PIECES[args.piece] + score_file = Path(piece_cfg["score"]) + performance_file = Path( + piece_cfg["midi"] if input_mode == "midi" else piece_cfg["audio"] + ) + match_file = Path(piece_cfg["match"]) if "match" in piece_cfg else None + run_name = args.piece + + print(f"Performance file: {performance_file.name}, with input mode: {input_mode}") + print(f"Running matchmaker with the score file ({score_file.name})...") print("-" * 50) if args.method is not None: @@ -127,10 +197,12 @@ def main(): # Initialize matchmaker (simulation mode) try: mm = Matchmaker( - score_file=SCORE_FILE, + score_file=score_file, performance_file=performance_file, input_type=input_mode, method=method, + unfold_score=args.unfold, + kwargs=args.kwargs or None, ) except Empty as e: print(f"Error initializing Matchmaker: {e}") @@ -141,17 +213,30 @@ def main(): timestamp = datetime.datetime.now().strftime("%H:%M:%S.%f")[:-3] print(f"[{timestamp}] Current beat position: {current_position}") - print("-" * 50) - print(f"Running evaluation using the match file ({MATCH_FILE.name})...") - - results, perf_sec, score_beat, gt_perf, gt_score = evaluate(mm, MATCH_FILE) - print(f"Evaluation Result: {json.dumps(results, indent=4)}") - - results_dir = ROOT_DIR / "results" - save_results( - results_dir, "simple_example", results, perf_sec, score_beat, gt_perf, gt_score - ) - print(f"Detailed evaluation results saved in {results_dir}") + if match_file and match_file.is_file(): + print("-" * 50) + print(f"Running evaluation using the match file ({match_file.name})...") + + results, perf_sec, score_beat, gt_perf, gt_score = evaluate(mm, match_file) + print(f"Evaluation Result: {json.dumps(results, indent=4)}") + + results_dir = args.output_dir + save_results( + results_dir, + f"{run_name}_{method}", + results, + perf_sec, + score_beat, + gt_perf, + gt_score, + save_plots=not args.no_plots, + ) + print(f"Detailed evaluation results saved in {results_dir}") + else: + print( + "\nAlignment finished (no match file provided for quantitative evaluation)." + ) + return if __name__ == "__main__": From 82bb65cd171991c485f74d1bb4ba6a86fbee0fd7 Mon Sep 17 00:00:00 2001 From: "jiyun.park" Date: Sat, 3 Oct 2026 22:22:20 +0200 Subject: [PATCH 06/11] GT onset times from note_on, not the float32 onset_sec gt_from_match took performance onsets from note_array()'s onset_sec, which partitura stores as float32. A MIDI follower stamps its updates with the float64 onset it receives, so the GT time of the same onset could fall a few tens of nanoseconds earlier, and the perf -> score lookup ("last position at or before t") then missed the update made at that onset and scored the previous event. In the eval fold about a third of MIDI onsets were scored one event late this way; 6/8 scores, where one event is up to a beat, failed the tracking criterion on paths that are correct. Use the PerformedPart notes' note_on (float64) instead: GT and follower times now coincide exactly. --- matchmaker/utils/eval.py | 9 ++++++++- 1 file changed, 8 insertions(+), 1 deletion(-) diff --git a/matchmaker/utils/eval.py b/matchmaker/utils/eval.py index 5e02c0a..57623fd 100644 --- a/matchmaker/utils/eval.py +++ b/matchmaker/utils/eval.py @@ -234,7 +234,14 @@ def gt_from_match(match_path, score_notes): ) for n in match_notes } - onset_at = {str(n["id"]): float(n["onset_sec"]) for n in perf.note_array()} + # note_on, not note_array()'s float32 onset_sec: a MIDI follower stamps its + # updates with the float64 onset, and a GT time rounded below it would hide + # the update made at that very onset from the perf -> score lookup. + onset_at = { + str(n["id"]): float(n["note_on"]) + for part in perf.performedparts + for n in part.notes + } beats, secs = [], [] for a in alignment: From c471813dd744b9252327b8e39b5590583683290b Mon Sep 17 00:00:00 2001 From: "jiyun.park" Date: Mon, 5 Oct 2026 14:55:57 +0200 Subject: [PATCH 07/11] Audio OPHMM: chord durations in quarter notes The bottom-level self-transitions took the chord durations from score beats times 60 / tempo. The tempo is quarter-note BPM, but score beats follow the time signature, so in 6/8 (eighth-note beats) every chord was expected to last twice as long, and in 2/2 half as long. Use onset_quarter; the reported score positions stay in beats. --- matchmaker/prob/outer_product_hmm_audio.py | 21 +++++++++++++++------ 1 file changed, 15 insertions(+), 6 deletions(-) diff --git a/matchmaker/prob/outer_product_hmm_audio.py b/matchmaker/prob/outer_product_hmm_audio.py index 12fffab..2e345c9 100644 --- a/matchmaker/prob/outer_product_hmm_audio.py +++ b/matchmaker/prob/outer_product_hmm_audio.py @@ -306,8 +306,16 @@ def __init__( } # Bottom transitions a_{l',l}^{(i)} and exit probs e_l^{(i)} (Eq.(5)) frame_rate = float(self.sample_rate) / float(self.hop_length) + # chord durations in quarter notes: `tempo` is quarter-note BPM, while + # score beats follow the time signature (an eighth in 6/8, a half in 2/2) + notes = ( + self.reference_features + if isinstance(self.reference_features, np.ndarray) + else self.reference_features.note_array() + ) + quarter_at = dict(zip(notes["onset_beat"].tolist(), notes["onset_quarter"].tolist())) self.a00 = self._compute_chord_self_transition_probs( - unique_onsets=unique_onsets, + onset_quarters=np.array([quarter_at[b] for b in unique_onsets.tolist()]), tempo=tempo, frame_rate=frame_rate, ) @@ -346,20 +354,21 @@ def _pause_self_transition_prob( @staticmethod def _compute_chord_self_transition_probs( - unique_onsets: np.ndarray, + onset_quarters: np.ndarray, tempo: float, frame_rate: float, ) -> np.ndarray: """ Compute self-transition probabilities from chord durations (Eq.5). - a_i = 1 - 1/d_i, where d_i = duration_sec / frame_time. + a_i = 1 - 1/d_i, where d_i = duration_sec / frame_time. ``onset_quarters`` + are the chord onsets in quarter notes and ``tempo`` is quarter-note BPM. """ - N = len(unique_onsets) + N = len(onset_quarters) frame_time = 1.0 / max(frame_rate, 1e-6) - # Convert onset beats to seconds, then compute inter-onset durations - onset_sec = unique_onsets * (60.0 / tempo) + # Convert onset quarters to seconds, then compute inter-onset durations + onset_sec = np.asarray(onset_quarters, dtype=float) * (60.0 / tempo) dur_sec = np.zeros(N, dtype=float) if N >= 2: dur_sec[:-1] = np.diff(onset_sec) From 948aca722f915853f89dedf572de8f601ea04148 Mon Sep 17 00:00:00 2001 From: "jiyun.park" Date: Wed, 7 Oct 2026 15:21:40 +0200 Subject: [PATCH 08/11] GT onset times from the streamed performance when it holds the match notes A .match file and the MIDI file a follower streams can store the same notes at different tick resolutions ((n)ASAP: 480 vs 384 ppq), so their onsets differ by up to half a tick. A GT time a fraction of a millisecond before the follower's update was scored against the position before that update. resolve_gt/gt_from_match take the streamed performance and use its onsets when every matched note has the same id and pitch in it. --- matchmaker/utils/eval.py | 51 ++++++++++++++++++++++++++++++++++--- tests/test_utils_eval_gt.py | 41 +++++++++++++++++++++++++++++ 2 files changed, 88 insertions(+), 4 deletions(-) create mode 100644 tests/test_utils_eval_gt.py diff --git a/matchmaker/utils/eval.py b/matchmaker/utils/eval.py index 57623fd..1c2e864 100644 --- a/matchmaker/utils/eval.py +++ b/matchmaker/utils/eval.py @@ -214,7 +214,39 @@ def evaluate_alignment( return {"beat": beat_results, "ms": ms_results} -def gt_from_match(match_path, score_notes): +def _streamed_onsets(performance, match_perf, performance_ids): + """Onset times of ``performance`` by note id, if it holds the match's notes. + + A .match file and the MIDI file a follower streams can store the same notes + with different tick resolutions (e.g. (n)ASAP: 480 vs 384 ppq), so their + times differ by up to half a tick. The follower stamps its updates with the + times it receives, so GT must use them: a GT time a fraction of a + millisecond early would be scored against the position before that update. + Returns None unless every matched note has the same id and pitch in both. + """ + if performance is None: + return None + if isinstance(performance, (str, Path)): + performance = pt.load_performance(str(performance)) + parts = getattr(performance, "performedparts", [performance]) + streamed = {str(n["id"]): n for part in parts for n in part.notes} + pitch_in_match = { + str(n["id"]): int(n["pitch"]) for part in match_perf.performedparts for n in part.notes + } + for pid in performance_ids: + note = streamed.get(pid) + if note is None or int(note["pitch"]) != pitch_in_match[pid]: + return None + return {pid: float(note["note_on"]) for pid, note in streamed.items()} + + +def gt_from_match(match_path, score_notes, performance=None): + """GT (perf_sec, score_beat) from a .match file. + + ``performance`` is the performance the follower receives (a MIDI path or a + partitura performance); when its notes are those of the match file, their + onset times are used instead of the match file's own. + """ perf, alignment, match_score = pt.load_match(str(match_path), create_score=True) match_notes = match_score.note_array() shared_ids = {str(n["id"]) for n in score_notes} & { @@ -242,6 +274,14 @@ def gt_from_match(match_path, score_notes): for part in perf.performedparts for n in part.notes } + matched_ids = [ + str(a["performance_id"]) + for a in alignment + if a.get("label") == "match" and str(a["performance_id"]) in onset_at + ] + streamed = _streamed_onsets(performance, perf, matched_ids) + if streamed is not None: + onset_at = streamed beats, secs = [], [] for a in alignment: @@ -261,14 +301,17 @@ def gt_from_match(match_path, score_notes): return secs[keep], beats[keep] -def resolve_gt(gt, score_notes): - """Return GT as (perf_sec, score_beat) arrays, from an ndarray, .match, or .tsv.""" +def resolve_gt(gt, score_notes, performance=None): + """Return GT as (perf_sec, score_beat) arrays, from an ndarray, .match, or .tsv. + + ``performance``: the performance the follower receives; see gt_from_match. + """ if isinstance(gt, np.ndarray): # gt array columns: [perf_sec, score_beat] arr = np.asarray(gt, dtype=float) return arr[:, 0], arr[:, 1] if Path(gt).suffix.lower() == ".match": - return gt_from_match(gt, score_notes) + return gt_from_match(gt, score_notes, performance) with open(gt) as f: header = f.readline().rstrip("\n").split("\t") data = np.loadtxt(gt, delimiter="\t", skiprows=1, ndmin=2) diff --git a/tests/test_utils_eval_gt.py b/tests/test_utils_eval_gt.py new file mode 100644 index 0000000..9983d1c --- /dev/null +++ b/tests/test_utils_eval_gt.py @@ -0,0 +1,41 @@ +import copy +import unittest + +import numpy as np +import partitura as pt + +from matchmaker.utils.eval import gt_from_match, resolve_gt + +MATCH = "matchmaker/assets/Bach-fugue_bwv_858.match" + + +class TestGTFromStreamedPerformance(unittest.TestCase): + """GT times come from the performance the follower streams when it holds the same notes.""" + + @classmethod + def setUpClass(cls): + cls.perf, _, score = pt.load_match(MATCH, create_score=True) + cls.score_notes = score.note_array() + cls.secs, cls.beats = gt_from_match(MATCH, cls.score_notes) + + def shifted(self, offset): + perf = copy.deepcopy(self.perf) + for note in perf.performedparts[0].notes: + note["note_on"] = float(note["note_on"]) + offset + return perf + + def test_streamed_onsets_replace_the_match_times(self): + secs, beats = resolve_gt(MATCH, self.score_notes, performance=self.shifted(4e-4)) + np.testing.assert_allclose(secs, self.secs + 4e-4) + np.testing.assert_array_equal(beats, self.beats) + + def test_other_notes_fall_back_to_the_match_times(self): + perf = self.shifted(4e-4) + perf.performedparts[0].notes[0]["pitch"] += 1 + secs, beats = gt_from_match(MATCH, self.score_notes, performance=perf) + np.testing.assert_array_equal(secs, self.secs) + np.testing.assert_array_equal(beats, self.beats) + + +if __name__ == "__main__": + unittest.main() From 28847937bfd351bd89fcf7a79979d2568cc1e19d Mon Sep 17 00:00:00 2001 From: "jiyun.park" Date: Wed, 7 Oct 2026 15:21:40 +0200 Subject: [PATCH 09/11] PF MIDI processor: stamp an observation with its first onset ParticleFilterMidiProcessor returned the middle of the polling window, which precedes the onsets in the window's second half and follows those in its first half. Use the first onset of the frame, as the other MIDI processors do. --- matchmaker/features/midi.py | 11 +++++++++-- 1 file changed, 9 insertions(+), 2 deletions(-) diff --git a/matchmaker/features/midi.py b/matchmaker/features/midi.py index d56d721..555b46e 100644 --- a/matchmaker/features/midi.py +++ b/matchmaker/features/midi.py @@ -525,13 +525,19 @@ def __call__( Returns ------- - (KorzeniowskiObservation, frame_time) + (KorzeniowskiObservation, time) + ``time`` is the first onset in the frame, or the frame time when + the frame has none. The frame time is the middle of the polling + window, which lies before the onsets in its second half and after + those in its first half; the other MIDI processors also stamp an + observation with its first note. """ data, f_time = frame # Reset onset list for this frame. self.current_onsets = [] + onset_time = None for msg, m_time in data: @@ -540,6 +546,7 @@ def __call__( # new onset in this frame self.current_onsets.append(msg.note) self.velocities[msg.note] = msg.velocity + onset_time = m_time if onset_time is None else min(onset_time, m_time) self.active_notes[msg.note] = ( msg.velocity, @@ -590,7 +597,7 @@ def __call__( loudness=self.compute_loudness(), ) - return observation, f_time + return observation, (f_time if onset_time is None else onset_time) def compute_loudness(self) -> float: """ From 6c56d3c4affb10e226fff511ce66c453b66be66f Mon Sep 17 00:00:00 2001 From: "jiyun.park" Date: Thu, 15 Oct 2026 10:17:05 +0200 Subject: [PATCH 10/11] feat(dp): implement and optimize Arzt tempo model integration in OnlineTimeWarpingArztFrame - Implement ArztTempoModel in matchmaker/utils/tempo_models.py (Arzt & Widmer 2010 simple tempo model with dynamic score alterations and standalone onset support) - Integrate ArztTempoModel into OnlineTimeWarpingArztFrame with bounded backward step and tempo projection - Support tempo model kwargs and benchmark evaluation in matchmaker/methods.yaml and utils/eval.py - Add unit tests for ArztTempoModel and OnlineTimeWarpingArztFrame --- .pre-commit-config.yaml | 1 + matchmaker/dp/__init__.py | 3 +- matchmaker/dp/dtw_loop.pyx | 7 +- matchmaker/dp/oltw_arzt.py | 620 +++++++++++++++---------------- matchmaker/features/audio.py | 23 +- matchmaker/methods.yaml | 6 +- matchmaker/utils/eval.py | 78 +++- matchmaker/utils/tempo_models.py | 313 +++++++++++++++- tests/test_oltw_arzt.py | 49 +++ tests/test_tempo_models.py | 165 +++++--- 10 files changed, 857 insertions(+), 408 deletions(-) diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index 460d66d..1fd957e 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -3,6 +3,7 @@ repos: rev: v4.4.0 hooks: - id: check-yaml + args: [--unsafe] - id: end-of-file-fixer - id: trailing-whitespace - repo: https://github.com/astral-sh/ruff-pre-commit diff --git a/matchmaker/dp/__init__.py b/matchmaker/dp/__init__.py index 7d7e8ed..c5ff815 100644 --- a/matchmaker/dp/__init__.py +++ b/matchmaker/dp/__init__.py @@ -12,11 +12,10 @@ OnlineTimeWarpingArzt, OnlineTimeWarpingArztEvent, OnlineTimeWarpingArztFrame, - OnlineTimeWarpingArztTempoFrame, ) +from .oltw_arzt_multi_fold import OnlineTimeWarpingArztMultiFold from .oltw_dixon import ( OnlineTimeWarpingDixon, OnlineTimeWarpingDixonEvent, OnlineTimeWarpingDixonFrame, ) -from .oltw_arzt_multi_fold import OnlineTimeWarpingArztMultiFold diff --git a/matchmaker/dp/dtw_loop.pyx b/matchmaker/dp/dtw_loop.pyx index b1f8242..94b918d 100644 --- a/matchmaker/dp/dtw_loop.pyx +++ b/matchmaker/dp/dtw_loop.pyx @@ -47,8 +47,9 @@ def oltw_arzt_loop( Minimum costs of current window. min_index: int Index of minimum costs of current window. - backpointers: numpy.ndarray, optional - 1D array of shape (window_end - window_start,) to record predecessor steps. + backpointers: numpy.ndarray or None + Optional array of shape (window_end - window_start,) to record backpointers + (1=vertical, 2=horizontal, 3=diagonal). Returns ------- @@ -59,7 +60,7 @@ def oltw_arzt_loop( min_costs: float Minimum costs of current window. """ - cdef float dist1, dist2, dist3, min_dist, local_dist, norm_cost + cdef float dist1, dist2, dist3, local_dist, norm_cost cdef int idx = 0, score_index = window_start cdef Py_ssize_t k, wk = window_cost.shape[0] cdef Py_ssize_t N = global_cost_matrix.shape[0] diff --git a/matchmaker/dp/oltw_arzt.py b/matchmaker/dp/oltw_arzt.py index 1dc789d..63558af 100644 --- a/matchmaker/dp/oltw_arzt.py +++ b/matchmaker/dp/oltw_arzt.py @@ -37,11 +37,18 @@ MatchmakerInvalidParameterTypeError, ) from matchmaker.utils.misc import set_latency_stats +from matchmaker.utils.tempo_models import ArztTempoModel, make_tempo_model STEP_SIZE: int = 3 START_WINDOW_SIZE: Union[float, int] = 0.1 WINDOW_SIZE: int = 10 +_ADV_IN, _ADV_SC, _BOTH = ( + 0, + 1, + 2, +) # Dixon (2005) GetInc: Row (input), Column (score), Both + # --------------------------------------------------------------------------- # Base class @@ -84,6 +91,7 @@ def __init__( ) self.N_ref: int = self.reference_features.shape[0] self.step_size: int = step_size + self.reset() def reset(self) -> None: self.current_index = 0 @@ -136,18 +144,12 @@ def __init__( start_window_size: Union[float, int] = START_WINDOW_SIZE, frame_rate: int = FRAME_RATE, ref_frame_to_beat: NDArray = None, + tempo_model=None, + backward_forward: Optional[bool] = None, + max_run_count: Optional[int] = None, queue: Optional[RECVQueue] = None, **kwargs, ) -> None: - super().__init__( - reference_features=reference_features, - score_positions=score_positions, - window_size=window_size, - step_size=step_size, - start_window_size=start_window_size, - queue=queue, - **kwargs, - ) if ref_frame_to_beat is None: raise ValueError( "Frame-level Arzt requires `ref_frame_to_beat` (per-frame beat mapping)." @@ -157,6 +159,24 @@ def __init__( self._window_size = int(np.round(window_size * self.frame_rate)) self._start_window_size = int(np.round(start_window_size * frame_rate)) self.queue_timeout = QUEUE_TIMEOUT + + self.tempo_model = make_tempo_model(tempo_model) + self.backward_forward = bool(backward_forward) + self.max_run_count = max_run_count + + if self.tempo_model is not None: + self._orig_features = np.array(reference_features, copy=True) + self._orig_f2b = np.asarray(ref_frame_to_beat, dtype=np.float64).copy() + + super().__init__( + reference_features=reference_features, + score_positions=score_positions, + window_size=window_size, + step_size=step_size, + start_window_size=start_window_size, + queue=queue, + **kwargs, + ) self.latency_stats: Dict[str, float] = { "total_latency": 0, "total_frames": 0, @@ -164,7 +184,6 @@ def __init__( "min_latency": float("inf"), } self._init_distance_func(distance_func) - self.reset() def __call__(self, observation: Any, perf_time: float) -> float: t0 = time.time() @@ -176,9 +195,238 @@ def __call__(self, observation: Any, perf_time: float) -> float: return beat def reset(self) -> None: + if self.tempo_model is not None: + self.reference_features = self._orig_features.copy() + self._ref_frame_to_beat = self._orig_f2b.copy() + self.N_ref = len(self.reference_features) + on = np.searchsorted(self._ref_frame_to_beat, self.score_positions) + self.is_onset = np.zeros(self.N_ref, dtype=bool) + self.is_onset[on[on < self.N_ref]] = True + if self.tempo_model is not None or self.backward_forward: + self._backpointers = {} + self._last_window_end = 0 + self._min_x_needed = 0 + self._run_count, self._run_prev = 0, None + if self.backward_forward: + K = 4 * self._window_size + self._cost_band = np.full( + (K, 2 * self._window_size + 64), np.inf, np.float32 + ) + self._cost_start = np.zeros(K, dtype=np.int64) + self._cost_len = np.zeros(K, dtype=np.int64) + self._n_bf = 0 super().reset() self.input_features: list = [] self._current_frame = 0 + if self.tempo_model is not None: + self.tempo_model.reset(self) + + @property + def current_relative_tempo(self) -> float: + return 1.0 if self.tempo_model is None else self.tempo_model.current_tempo + + @property + def n_deleted(self) -> int: + return self.tempo_model.n_deleted if self.tempo_model else 0 + + @property + def n_inserted(self) -> int: + return self.tempo_model.n_inserted if self.tempo_model else 0 + + @property + def use_tempo_model(self) -> bool: + return self.tempo_model is not None + + @property + def N(self) -> int: + return len(self.reference_features) + + @property + def c(self) -> int: + return self._window_size + + @property + def V(self) -> np.ndarray: + return self.reference_features + + # Interface for the tempo model (TempoModel). + + def tempo_last_input(self) -> int: + return self.input_index - 1 + + def tempo_current_cell(self) -> Tuple[int, int]: + return self.input_index - 1, self._current_frame + + def tempo_backward_path(self, start_x: int, start_y: int): + x, y = start_x, start_y + while x > 0 and y > 0: + entry = self._backpointers.get(x) + if entry is None: + break + ws, codes = entry + offset = y - ws + if offset < 0 or offset >= len(codes): + break + code = codes[offset] + if code == 1: # vertical: (x, y-1) + y -= 1 + elif code == 2: # horizontal: (x-1, y) + x -= 1 + elif code == 3: # diagonal: (x-1, y-1) + x -= 1 + y -= 1 + else: + break + yield x, y + + def tempo_note_path_end(self, x: int, y: int) -> None: + self._min_x_needed = x if not self._min_x_needed else min(self._min_x_needed, x) + + def tempo_last_score_frame(self) -> int: + return self._current_frame + + def tempo_score(self): + return self.reference_features, self._ref_frame_to_beat, self.is_onset + + def tempo_set_score(self, features, frame_to_beat, is_onset) -> None: + n = len(features) + ps = self._current_frame + G = self.global_cost_matrix # row k + 1 = score frame k + if n < self.N_ref: + r = ps + 2 + self.global_cost_matrix = np.ascontiguousarray(np.delete(G, r + 1, axis=0)) + for x, (s, c) in list(self._backpointers.items()): + if s > r: + self._backpointers[x] = (s - 1, c) + elif r - s < len(c): + self._backpointers[x] = (s, np.delete(c, r - s)) + if self.backward_forward: + self._band_alter(r, -1) + elif n > self.N_ref: + r = ps + 1 + if r + 1 < len(G): + interp = 0.5 * (G[r] + G[r + 1]) + else: + interp = np.full(2, np.inf, dtype=np.float32) + self.global_cost_matrix = np.ascontiguousarray( + np.insert(G, r + 1, interp, axis=0) + ) + for x, (s, c) in list(self._backpointers.items()): + if s >= r: + self._backpointers[x] = (s + 1, c) + elif r - s < len(c): + self._backpointers[x] = (s, np.insert(c, r - s, 1).astype(np.int8)) + if self.backward_forward: + self._band_alter(r, 1) + self.reference_features = features + self._ref_frame_to_beat = frame_to_beat + self.is_onset = is_onset + self.N_ref = n + + # Strategy 1 (Arzt, Widmer & Dixon 2008, Sec. 3.2). + + def _band_alter(self, r, sign): + st, ln, band = self._cost_start, self._cost_len, self._cost_band + for k in range(len(band)): + s, n = st[k], ln[k] + if sign < 0: + if s > r: + st[k] = s - 1 + elif r - s < n: + band[k, r - s : n - 1] = band[k, r - s + 1 : n].copy() + band[k, n - 1] = np.inf + ln[k] = n - 1 + else: + if s >= r: + st[k] = s + 1 + elif r - s < n: + m = min(n + 1, band.shape[1]) + band[k, r - s + 1 : m] = band[k, r - s : m - 1].copy() + band[k, r - s] = np.inf + ln[k] = m + + def _cost_row(self, y, x0, x1): + xs = np.arange(x0, x1 + 1) + K = len(self._cost_band) + sl = xs % K + k = y - self._cost_start[sl] + ok = ( + (xs > self.input_index - K) + & (xs <= self.input_index) + & (k >= 0) + & (k < self._cost_len[sl]) + ) + out = np.full(len(xs), np.inf) + out[ok] = self._cost_band[sl[ok], k[ok]] + return out + + def _cost_col(self, x, y0, y1): + out = np.full(y1 - y0 + 1, np.inf) + K = len(self._cost_band) + if x > self.input_index or x <= self.input_index - K or x < 0: + return out + sl = x % K + s, n = self._cost_start[sl], self._cost_len[sl] + a, b = max(y0, s), min(y1, s + n - 1) + if a <= b: + out[a - y0 : b - y0 + 1] = self._cost_band[sl, a - s : b - s + 1] + return out + + def _cost(self, x, y): + return float(self._cost_col(x, y, y)[0]) + + def _update_runs(self, move: int) -> None: + if move == self._run_prev: + self._run_count += 1 + else: + self._run_count, self._run_prev = 1, move + + def _get_inc(self, x, y): + M = self.max_run_count + if M is not None and self._run_count > M: + return _ADV_SC if self._run_prev == _ADV_IN else _ADV_IN + c = self._window_size + best = self._cost(x, y) / (x + y + 1.0) + bx, by = x, y + x0, y0 = max(0, x - c + 1), max(0, y - c + 1) + r = self._cost_row(y, x0, x) / (np.arange(x0, x + 1) + y + 1.0) + k = int(np.argmin(r)) + if r[k] < best: + best, bx, by = r[k], x0 + k, y + cl = self._cost_col(x, y0, y) / (x + np.arange(y0, y + 1) + 1.0) + k = int(np.argmin(cl)) + if cl[k] < best: + best, bx, by = cl[k], x, y0 + k + if bx < x: + return _ADV_SC + if by < y: + return _ADV_IN + return _BOTH + + def _strategy1(self, T, y): + b = 50 if self._n_bf % 5 == 4 else 10 + self._n_bf += 1 + x, target = T, y - b + for x, y in self.tempo_backward_path(T, y): + if y <= target: + break + self.tempo_note_path_end(x, y) + J = self._last_window_end - 1 + self._run_count, self._run_prev = 0, None + while x < T and y < J: + if self._get_inc(x, y) != _ADV_SC: + if not np.isfinite(self._cost(x + 1, y)) and not np.isfinite( + self._cost(x + 1, y + 1) + ): + break + x += 1 + if x <= T and y < J and self._get_inc(x, y) != _ADV_IN: + y += 1 + if not np.isfinite(self._cost(x, y)): + break + if self.max_run_count is not None: + self._update_runs(self._get_inc(x, y)) + return y @property def window_index(self) -> int: @@ -197,6 +445,23 @@ def _frame_to_score_idx(self, frame: int) -> int: return max(0, min(idx, len(self.score_positions) - 1)) def get_current_position(self) -> float: + if self.tempo_model is not None and self.input_index > 5: + x, y = self.input_index - 1, self._current_frame + path = [(x, y)] + for bx, by in self.tempo_backward_path(x, y): + path.append((bx, by)) + if len(path) >= 3: + break + if len(path) >= 2: + bx, by = path[-1] + dx = x - bx + tempo = self.tempo_model.current_tempo + pred_y = by + dx * tempo + return float( + self._frame_to_beat( + int(round(np.clip(pred_y, 0, len(self._ref_frame_to_beat) - 1))) + ) + ) return self._frame_to_beat(self._current_frame) def _init_distance_func(self, distance_func): @@ -244,298 +509,8 @@ def get_window(self) -> Tuple[int, int]: return start, end def step(self, input_features: NDArray[np.float32]) -> None: - min_costs = np.inf - min_index = max(self.window_index - self.step_size, 0) - - window_start, window_end = self.get_window() - window_cost = self.vdist( - self.reference_features[window_start:window_end], - input_features.squeeze(), - self.distance_func, - ) - - self.global_cost_matrix, min_index, min_costs = oltw_arzt_loop( - global_cost_matrix=self.global_cost_matrix, - window_cost=window_cost, - window_start=window_start, - window_end=window_end, - input_index=self.input_index, - min_costs=min_costs, - min_index=min_index, - ) - - if self.input_index > 0: - self._current_frame = int( - np.clip( - min_index, - self._current_frame, - self._current_frame + self.step_size, - ) - ) - else: - self._current_frame = min_index - self.current_index = self._frame_to_score_idx(self._current_frame) - self.input_index += 1 - - -class OnlineTimeWarpingArztTempoFrame(OnlineTimeWarpingArztFrame): - """Frame-level OLTW with Arzt & Widmer (2010) Simple Tempo Model. - - Dynamically tracks relative tempo via rectified backward path - and stretches/compresses reference score features on the fly. - """ - - def __init__( - self, - reference_features: NDArray[np.float32], - score_positions: NDArray[np.float32], - window_size: int = WINDOW_SIZE, - step_size: int = STEP_SIZE, - distance_func: Union[ - str, Callable, Tuple[str, Dict[str, Any]] - ] = OnlineTimeWarpingArztFrame.DEFAULT_DISTANCE_FUNC, - start_window_size: Union[float, int] = START_WINDOW_SIZE, - frame_rate: int = FRAME_RATE, - ref_frame_to_beat: NDArray = None, - queue: Optional[RECVQueue] = None, - tempo_window_size: float = 3.0, - tempo_min_past: float = 1.0, - tempo_n_onsets: int = 20, - random_seed: Optional[int] = 1984, - use_tempo_model: bool = True, - **kwargs, - ) -> None: - if ref_frame_to_beat is None: - raise ValueError( - "Frame-level Arzt requires `ref_frame_to_beat` (per-frame beat mapping)." - ) - self.use_tempo_model = use_tempo_model - self.tempo_window_size = tempo_window_size - self.tempo_min_past = tempo_min_past - self.tempo_n_onsets = tempo_n_onsets - self.random_seed = random_seed - self._original_reference_features = reference_features.copy() - self._original_ref_frame_to_beat = ref_frame_to_beat.copy() - super().__init__( - reference_features=reference_features, - score_positions=score_positions, - window_size=window_size, - step_size=step_size, - distance_func=distance_func, - start_window_size=start_window_size, - frame_rate=frame_rate, - ref_frame_to_beat=ref_frame_to_beat, - queue=queue, - **kwargs, - ) - self.reset() - - def reset(self) -> None: - super().reset() - self.current_relative_tempo: float = 1.0 - self._consecutive_alterations: int = 0 - self._backpointers: list = [] - self._rng = np.random.RandomState(self.random_seed) - # The tempo model stretches the reference on the fly; restore the - # pristine copy so each run starts from the notated score. - self.reference_features = self._original_reference_features.copy() - self._ref_frame_to_beat = self._original_ref_frame_to_beat.copy() - self.N_ref = self.reference_features.shape[0] - self.global_cost_matrix = np.full( - (self.N_ref + 1, 2), np.inf, dtype=np.float32 - ) - self._init_onset_mask() - - def _init_onset_mask(self) -> None: - self.is_onset_frame = np.zeros(self.N_ref, dtype=bool) - onset_frames = np.searchsorted(self._ref_frame_to_beat, self.score_positions) - valid = onset_frames < self.N_ref - self.is_onset_frame[onset_frames[valid]] = True - - def _get_backward_path(self, max_history_frames: int = 500) -> list: - """Trace backward path from current position using recorded backpointers.""" - if not self._backpointers or self.input_index <= 0: - return [(self._current_frame, 0)] - - path = [] - s = self._current_frame - i = self.input_index - 1 - min_i = max(0, i - max_history_frames) - while i >= min_i and s >= 0: - path.append((s, i)) - if i >= len(self._backpointers): - i -= 1 - s = max(0, s - 1) - continue - w_start, w_end, bp = self._backpointers[i] - if bp is not None and w_start <= s < w_end: - code = bp[s - w_start] - if code == 1: - s -= 1 - elif code == 2: - i -= 1 - elif code == 3: - s -= 1 - i -= 1 - else: - break - else: - s = max(0, s - 1) - i -= 1 - path.reverse() - return path - - def _compute_relative_tempo(self) -> float: - """Estimate current relative tempo via rectified backward path (Arzt 2010 Sec. 4.1).""" - if self.input_index < int(self.frame_rate * 1.5): - return 1.0 - - cur_perf_time = self.input_index / self.frame_rate - max_history = int(self.frame_rate * 10) - path = self._get_backward_path(max_history_frames=max_history) - if len(path) < 2: - return 1.0 - - score_frames = np.array([p[0] for p in path]) - perf_frames = np.array([p[1] for p in path]) - - s_min, s_max = int(score_frames[0]), int(score_frames[-1]) - if s_max <= s_min: - return 1.0 - - local_onsets = np.where(self.is_onset_frame[s_min : s_max + 1])[0] + s_min - if len(local_onsets) < 2: - return 1.0 - - t_s_list = [] - t_p_list = [] - for sf in local_onsets: - idx = int(np.searchsorted(score_frames, sf)) - if idx < len(score_frames): - t_s_list.append(sf / self.frame_rate) - t_p_list.append(perf_frames[idx] / self.frame_rate) - - if len(t_s_list) < 2: - return 1.0 - - t_s_arr = np.array(t_s_list) - t_p_arr = np.array(t_p_list) - - # Onsets at least tempo_min_past seconds in the past - valid_mask = t_p_arr <= (cur_perf_time - self.tempo_min_past) - valid_indices = np.where(valid_mask)[0] - if len(valid_indices) == 0: - return 1.0 - - if len(valid_indices) > self.tempo_n_onsets: - valid_indices = valid_indices[-self.tempo_n_onsets:] - - half_win = self.tempo_window_size / 2.0 - tempi = [] - for idx in valid_indices: - ts = t_s_arr[idx] - w_start = ts - half_win - w_end = ts + half_win - if w_start >= t_s_arr[0] and w_end <= t_s_arr[-1]: - tp_start = float(np.interp(w_start, t_s_arr, t_p_arr)) - tp_end = float(np.interp(w_end, t_s_arr, t_p_arr)) - delta_p = tp_end - tp_start - delta_s = w_end - w_start - if delta_p > 1e-4: - tempi.append(delta_s / delta_p) - else: - tempi.append(1.0) - - if not tempi: - return 1.0 - - m = len(tempi) - weights = np.arange(1, m + 1) - t_est = float(np.sum(np.array(tempi) * weights) / np.sum(weights)) - return float(np.clip(t_est, 0.33, 3.0)) - - def _alter_score_representation(self, t: float) -> None: - """Alter score representation on the fly according to relative tempo (Arzt 2010 Sec. 4.2).""" - ps = self._current_frame - if ps + 2 >= self.N_ref: - return - - if self._consecutive_alterations >= 3: - self._consecutive_alterations = 0 - return - - r = self._rng.uniform(0.0, 1.0) - - if t > 1.0: - threshold = 1.0 / t - if r > threshold: - if self.is_onset_frame[ps + 1] or self.is_onset_frame[ps + 2]: - self._consecutive_alterations = 0 - return - - mean_feat = 0.5 * ( - self.reference_features[ps + 1] + self.reference_features[ps + 2] - ) - self.reference_features[ps + 1] = mean_feat - self.reference_features = np.delete( - self.reference_features, ps + 2, axis=0 - ) - - mean_beat = 0.5 * ( - self._ref_frame_to_beat[ps + 1] + self._ref_frame_to_beat[ps + 2] - ) - self._ref_frame_to_beat[ps + 1] = mean_beat - self._ref_frame_to_beat = np.delete(self._ref_frame_to_beat, ps + 2) - - self.is_onset_frame = np.delete(self.is_onset_frame, ps + 2) - self.global_cost_matrix = np.delete( - self.global_cost_matrix, ps + 3, axis=0 - ) - self.N_ref = self.reference_features.shape[0] - self._consecutive_alterations += 1 - else: - self._consecutive_alterations = 0 - - elif t < 1.0: - threshold = t - if r > threshold: - if self.is_onset_frame[ps] or self.is_onset_frame[ps + 1]: - self._consecutive_alterations = 0 - return - - mean_feat = 0.5 * ( - self.reference_features[ps] + self.reference_features[ps + 1] - ) - self.reference_features = np.insert( - self.reference_features, ps + 1, mean_feat, axis=0 - ) - - mean_beat = 0.5 * ( - self._ref_frame_to_beat[ps] + self._ref_frame_to_beat[ps + 1] - ) - self._ref_frame_to_beat = np.insert( - self._ref_frame_to_beat, ps + 1, mean_beat - ) - - self.is_onset_frame = np.insert(self.is_onset_frame, ps + 1, False) - row1 = np.asarray(self.global_cost_matrix[ps + 1]) - row2 = np.asarray(self.global_cost_matrix[ps + 2]) - new_cost_row = 0.5 * (row1 + row2) - new_cost_row[1] = np.inf - self.global_cost_matrix = np.insert( - self.global_cost_matrix, ps + 2, new_cost_row, axis=0 - ) - self.N_ref = self.reference_features.shape[0] - self._consecutive_alterations += 1 - else: - self._consecutive_alterations = 0 - else: - self._consecutive_alterations = 0 - - def step(self, input_features: NDArray[np.float32]) -> None: - if self.use_tempo_model and self.input_index > 0: - self.current_relative_tempo = self._compute_relative_tempo() - self._alter_score_representation(self.current_relative_tempo) + if self.tempo_model is not None and self.input_index > 0: + self.tempo_model.before_step(self) min_costs = np.inf min_index = max(self.window_index - self.step_size, 0) @@ -547,9 +522,10 @@ def step(self, input_features: NDArray[np.float32]) -> None: self.distance_func, ) - bp = None - if self.use_tempo_model: - bp = np.zeros(window_end - window_start, dtype=np.int8) + store = (self.tempo_model is not None) or self.backward_forward + backpointers = None + if store: + backpointers = np.zeros(window_end - window_start, dtype=np.int8) self.global_cost_matrix, min_index, min_costs = oltw_arzt_loop( global_cost_matrix=self.global_cost_matrix, @@ -559,22 +535,43 @@ def step(self, input_features: NDArray[np.float32]) -> None: input_index=self.input_index, min_costs=min_costs, min_index=min_index, - backpointers=bp, + backpointers=backpointers, ) - if self.use_tempo_model: - self._backpointers.append((window_start, window_end, bp)) + if self.backward_forward: + k = self.input_index % len(self._cost_band) + n = min(window_end - window_start, self._cost_band.shape[1]) + self._cost_band[k, :n] = self.global_cost_matrix[ + window_start + 1 : window_start + 1 + n, 0 + ] + self._cost_start[k], self._cost_len[k] = window_start, n + + if store: + self._backpointers[self.input_index] = (window_start, backpointers) + self._last_window_end = window_end + if (self.input_index + 1) % self._window_size == 0: + keep = min(self.input_index - 2 * self._window_size, self._min_x_needed) + for x in [x for x in list(self._backpointers.keys()) if x < keep - 1]: + del self._backpointers[x] + self._min_x_needed = 0 if self.input_index > 0: - self._current_frame = int( + cf = self._current_frame + min_pos = max(0, cf - 1) if self.tempo_model is not None else cf + pos = int( np.clip( min_index, - self._current_frame, - self._current_frame + self.step_size, + min_pos, + cf + self.step_size, ) ) + if self.backward_forward and (self.input_index + 1) % 2 == 0: + pos = self._strategy1(self.input_index, pos) + pos = int(np.clip(pos, min_pos, cf + self.step_size)) + self._current_frame = int(pos) else: self._current_frame = min_index + self.current_index = self._frame_to_score_idx(self._current_frame) self.input_index += 1 @@ -614,7 +611,6 @@ def __init__( self._window_size = window_size self._start_window_size = start_window_size self._scipy_metric = distance_func - self.reset() def get_window(self) -> Tuple[int, int]: w = self._window_size diff --git a/matchmaker/features/audio.py b/matchmaker/features/audio.py index 00ddd1b..1232fab 100644 --- a/matchmaker/features/audio.py +++ b/matchmaker/features/audio.py @@ -9,7 +9,7 @@ import librosa import numpy as np -from matchmaker.features.processor import Processor, KorzeniowskiObservation +from matchmaker.features.processor import KorzeniowskiObservation, Processor SAMPLE_RATE = 44100 FRAME_RATE = 30 @@ -382,7 +382,6 @@ def __init__( win_length: int = 2048, n_fft: int = 4096, ): - super().__init__() self.sample_rate = sample_rate @@ -401,7 +400,6 @@ def __init__( self.frame_index = 0 def reset(self): - self.previous_spectrum = None self.frame_index = 0 @@ -410,17 +408,12 @@ def __call__( self, data: InputAudioFrame, ): - frame, f_time = data spectrum = self.compute_spectrum(frame) - onset = self.compute_onset( - frame - ) + onset = self.compute_onset(frame) - loudness = self.compute_loudness( - frame - ) + loudness = self.compute_loudness(frame) observation = KorzeniowskiObservation( spectrum=spectrum, @@ -429,7 +422,7 @@ def __call__( ) return observation, f_time - + def compute_spectrum( self, frame: np.ndarray, @@ -438,7 +431,7 @@ def compute_spectrum( Compute the normalized magnitude spectrum. """ - frame = frame[:self.win_length] * self.window + frame = frame[: self.win_length] * self.window magnitude = np.abs( np.fft.rfft( @@ -454,7 +447,6 @@ def compute_spectrum( return magnitude - def compute_onset(self, frame: np.ndarray) -> float: """ Compute normalized causal spectral-flux onset activation. @@ -469,7 +461,7 @@ def compute_onset(self, frame: np.ndarray) -> float: float Non-negative normalized onset activation. """ - windowed = frame[:self.win_length] * self.window + windowed = frame[: self.win_length] * self.window spectrum = np.abs(np.fft.rfft(windowed, n=self.n_fft)) @@ -490,7 +482,6 @@ def compute_onset(self, frame: np.ndarray) -> float: return float(flux / (energy + 1e-8)) - def compute_loudness( self, frame: np.ndarray, @@ -516,4 +507,4 @@ def compute_loudness( ref=1.0, )[0, 0] - return float(loudness) \ No newline at end of file + return float(loudness) diff --git a/matchmaker/methods.yaml b/matchmaker/methods.yaml index cabe986..cb32522 100644 --- a/matchmaker/methods.yaml +++ b/matchmaker/methods.yaml @@ -171,16 +171,16 @@ methods: arzt_tempo: family: oltw - class: matchmaker.dp:OnlineTimeWarpingArztTempoFrame + class: matchmaker.dp:OnlineTimeWarpingArztFrame reference: score_audio args: *oltw_frame_args config_passthrough: true default_kwargs: + processor: chroma window_size: 10 start_window_size: 0.1 step_size: 3 - use_tempo_model: true - random_seed: 1984 + tempo_model: !obj matchmaker.utils.tempo_models:ArztTempoModel arzt_multi_fold: family: oltw diff --git a/matchmaker/utils/eval.py b/matchmaker/utils/eval.py index 1c2e864..ccedbe8 100644 --- a/matchmaker/utils/eval.py +++ b/matchmaker/utils/eval.py @@ -145,6 +145,69 @@ def get_evaluation_results( return results +def compute_sparc( + score_beat: np.ndarray, + perf_sec: np.ndarray, + fs: float = 50.0, + fc: float = 10.0, + padlevel: int = 4, +) -> float: + """Compute Spectral Arc Length (SPARC) for an alignment path. + + Parameters + ---------- + score_beat : np.ndarray + Score positions (in beats). + perf_sec : np.ndarray + Performance timestamps (in seconds). + fs : float + Uniform resampling frequency in Hz (default: 50.0 Hz). + fc : float + Cutoff frequency in Hz for arc length integration (default: 10.0 Hz). + padlevel : int + Zero-padding level for FFT (default: 4). + + Returns + ------- + float + SPARC value (strictly negative; values closer to 0 indicate smoother tracking). + """ + score_beat = np.asarray(score_beat, dtype=float) + perf_sec = np.asarray(perf_sec, dtype=float) + + valid = np.isfinite(score_beat) & np.isfinite(perf_sec) + score_beat = score_beat[valid] + perf_sec = perf_sec[valid] + + if len(perf_sec) < 2: + return 0.0 + + t_start, t_end = float(perf_sec[0]), float(perf_sec[-1]) + duration = t_end - t_start + if duration <= 0.1: + return 0.0 + + n_pts = max(int(np.round(duration * fs)), 2) + t_uniform = np.linspace(t_start, t_end, n_pts) + s_interp = np.interp(t_uniform, perf_sec, score_beat) + vel = np.gradient(s_interp, 1.0 / fs) + + nfft = int(2 ** (np.ceil(np.log2(len(vel))) + padlevel)) + freq = np.fft.rfftfreq(nfft, d=1.0 / fs) + mask = freq <= fc + freq_filtered = freq[mask] + + Mf = np.abs(np.fft.rfft(vel, n=nfft))[mask] + max_mf = Mf.max() + if max_mf > 0: + Mf = Mf / max_mf + + d_freq = np.diff(freq_filtered) / fc + d_mf = np.diff(Mf) + arc = np.sum(np.sqrt(d_freq**2 + d_mf**2)) + return float(-arc) + + def evaluate_alignment( aligned_score_beats, aligned_perf_sec, @@ -173,7 +236,7 @@ def evaluate_alignment( Returns ------- dict - {"beat": {...}, "ms": {...}} + {"beat": {...}, "ms": {...}, "sparc": float} """ # GT interpolator: score beat → perf time valid = np.isfinite(gt_perf_sec) & np.isfinite(gt_score_beats) @@ -211,7 +274,9 @@ def evaluate_alignment( tolerances=ms_tolerances, ) - return {"beat": beat_results, "ms": ms_results} + sparc_val = compute_sparc(aligned_score_beats, aligned_perf_sec) + + return {"beat": beat_results, "ms": ms_results, "sparc": sparc_val} def _streamed_onsets(performance, match_perf, performance_ids): @@ -231,7 +296,9 @@ def _streamed_onsets(performance, match_perf, performance_ids): parts = getattr(performance, "performedparts", [performance]) streamed = {str(n["id"]): n for part in parts for n in part.notes} pitch_in_match = { - str(n["id"]): int(n["pitch"]) for part in match_perf.performedparts for n in part.notes + str(n["id"]): int(n["pitch"]) + for part in match_perf.performedparts + for n in part.notes } for pid in performance_ids: note = streamed.get(pid) @@ -266,9 +333,6 @@ def gt_from_match(match_path, score_notes, performance=None): ) for n in match_notes } - # note_on, not note_array()'s float32 onset_sec: a MIDI follower stamps its - # updates with the float64 onset, and a GT time rounded below it would hide - # the update made at that very onset from the perf -> score lookup. onset_at = { str(n["id"]): float(n["note_on"]) for part in perf.performedparts @@ -311,7 +375,7 @@ def resolve_gt(gt, score_notes, performance=None): arr = np.asarray(gt, dtype=float) return arr[:, 0], arr[:, 1] if Path(gt).suffix.lower() == ".match": - return gt_from_match(gt, score_notes, performance) + return gt_from_match(gt, score_notes, performance=performance) with open(gt) as f: header = f.readline().rstrip("\n").split("\t") data = np.loadtxt(gt, delimiter="\t", skiprows=1, ndmin=2) diff --git a/matchmaker/utils/tempo_models.py b/matchmaker/utils/tempo_models.py index 897bd50..0a4478e 100644 --- a/matchmaker/utils/tempo_models.py +++ b/matchmaker/utils/tempo_models.py @@ -3,12 +3,17 @@ """ Tempo Models -TODO ----- -* Adapt models from ACCompanion +Includes: +* ReactiveTempoModel: Baseline raw IOI ratio model. +* MovingAverageTempoModel: Smoothed moving average IOI model. +* KalmanTempoModel: Linear Kalman filter tempo tracker (ACCompanion). +* LinearTempoModel: Sensorimotor synchronization linear model. +* JointAdaptationAnticipationModel: Joint adaptation & anticipation model. +* ArztTempoModel: Simple tempo model of Arzt & Widmer (SMC 2010) with backward + path estimation and dynamic score alterations for OLTW, also usable on onsets. """ -from typing import Optional, Tuple +from typing import List, Optional, Tuple, Union import numpy as np @@ -296,7 +301,7 @@ def update_beat_period( self.prev_perf_onset = performed_onset self.prev_score_onset = score_onset self.counter += 1 - return # skip update on first observation + return # skip update on first observation else: performed_ioi = abs(performed_onset - self.prev_perf_onset) @@ -527,5 +532,303 @@ def update_beat_period( self.mns.append(self.motor_noise.normal(0, 0.0025)) +class ArztTempoModel(TempoModel): + """ + Simple tempo model of Arzt & Widmer (2010). + + References + ---------- + A. Arzt & G. Widmer (2010), "Simple Tempo Models for Real-Time Music + Tracking", Sound and Music Computing Conference (SMC), Sec. 4. + + When used with an on-line time warping tracker (e.g., + ``OnlineTimeWarpingArztFrame``), after every input frame, the relative + tempo t is the Eq. 1-weighted mean of the local tempi at the ``n_onsets`` + most recent score onsets at least ``min_past`` seconds in the past, each + over a ``window_size``-second window of the rectified backward path (as in + Mueller et al. (2009), ISMIR, Sec. 3.1 / 3.3: path rectified at onsets, + extended with slope 1 past its ends). The score representation is then + altered at the last processed score frame ps: with probability 1 - 1/t + (t > 1) frames ps+1 and ps+2 are merged into their mean, with probability + 1 - t (t < 1) the mean of ps and ps+1 is inserted; onset frames are not + duplicated (the insertion is postponed), and at most 3 alterations are + made in a row. + + When used as a standalone tempo model (e.g. with discrete performed and score + onsets), ``update_beat_period`` calculates local tempi over recent onset intervals + and computes the Eq. 1 weighted average. + + Parameters + ---------- + window_size : float + Window around an onset for its local tempo, in seconds (default: 3.0). + min_past : float + Onsets used must be at least this many seconds in the past (default: 1.0). + n_onsets : int + Number of recent onsets averaged (default: 20). + random_seed : int or None + Seed for the random score alterations (default: 1984). + tempo_min : float + Lower bound on relative tempo factor (default: 0.85). + tempo_max : float + Upper bound on relative tempo factor (default: 1.20). + min_advance : int + Minimum score frames between successive alterations (default: 10). + init_beat_period : float + Initial tempo in seconds per beat (default: 0.5). + init_score_onset : float + Initial score onset time in beats (default: 0.0). + score_onsets : Optional[np.ndarray] + Score onsets array in beats. + """ + + window_size: float + min_past: float + n_onsets: int + random_seed: Optional[int] + tempo_min: float + tempo_max: float + min_advance: int + n_inserted: int + n_deleted: int + n_postponed: int + + def __init__( + self, + window_size: float = 3.0, + min_past: float = 1.0, + n_onsets: int = 20, + random_seed: Optional[int] = 1984, + tempo_min: float = 0.85, + tempo_max: float = 1.20, + min_advance: int = 10, + init_beat_period: float = 0.5, + init_score_onset: float = 0.0, + score_onsets: Optional[np.ndarray] = None, + ) -> None: + super().__init__( + init_beat_period=init_beat_period, + init_score_onset=init_score_onset, + score_onsets=score_onsets, + ) + self.init_beat_period = init_beat_period + self.window_size = window_size + self.min_past = min_past + self.n_onsets = n_onsets + self.random_seed = random_seed + self.tempo_min = tempo_min + self.tempo_max = tempo_max + self.min_advance = min_advance + self._current_tempo = 1.0 + self.n_inserted = self.n_deleted = self.n_postponed = 0 + self._last_altered_ps: int = -1 + self._consecutive = 0 + self._insert_pending = False + self._w = 0 + self._min_past = 0 + self._rng = ( + np.random.RandomState(random_seed) + if random_seed is not None + else np.random.RandomState() + ) + self._history_perf: List[float] = [] + self._history_score: List[float] = [] + + def reset(self, tracker=None) -> None: + """Reset the tempo model state.""" + self.counter = 0 + self.prev_perf_onset = None + self.est_onset = None + self.asynchrony = 0.0 + self.beat_period = getattr(self, "init_beat_period", 0.5) + self._consecutive = 0 + self._insert_pending = False + self._current_tempo = 1.0 + self.n_inserted = self.n_deleted = self.n_postponed = 0 + self._last_altered_ps = -1 + self._history_perf = [] + self._history_score = [] + if tracker is not None: + fr = tracker.frame_rate + self._w = int(round(self.window_size * fr)) + self._min_past = int(round(self.min_past * fr)) + self._rng = ( + np.random.RandomState(self.random_seed) + if self.random_seed is not None + else np.random.RandomState() + ) + elif self.random_seed is not None: + self._rng = np.random.RandomState(self.random_seed) + + @property + def current_tempo(self) -> float: + """Current relative tempo (performance / score representation).""" + return self._current_tempo + + # --- Tracker integration for OLTW --- + + def before_step(self, tracker) -> None: + """Estimate the tempo and alter the tracker's score representation.""" + t = self.relative_tempo(tracker) + if t is not None: + self._current_tempo = t + base_bp = getattr(self, "init_beat_period", 0.5) + self.beat_period = base_bp / max(1e-4, t) + self.alter(tracker, t) + + def relative_tempo(self, tracker) -> Optional[float]: + """Sec. 4.1: relative tempo from the rectified backward path.""" + half = (self._w - 1) // 2 + limit = tracker.tempo_last_input() - self._min_past + is_onset = tracker.tempo_score()[2] + x, y = tracker.tempo_current_cell() + y_start = y + phi = {} # row -> min input index on the backward path + onsets = [] # onset rows on the path, most recent first + last_row = y + phi[y] = x + need_below = None + for x, y in tracker.tempo_backward_path(x, y): + phi[y] = x + if y != last_row: + if is_onset[last_row] and phi[last_row] <= limit: + onsets.append(last_row) + last_row = y + if need_below is None and len(onsets) >= self.n_onsets: + need_below = float(onsets[self.n_onsets - 1]) - half + if need_below is not None and float(y) < need_below and is_onset[y]: + break + if is_onset[last_row] and phi[last_row] <= limit and last_row not in onsets: + onsets.append(last_row) + tracker.tempo_note_path_end(x, y) + if not onsets: + return None + # anchors: onsets on the path plus both ends (Mueller et al. 2009, Sec. 3.3) + rows = sorted(set([y, y_start] + [r for r in phi if is_onset[r]])) + g = np.array(rows, dtype=np.float64) + p = np.array([phi[r] for r in rows], dtype=np.float64) + + def phi_r(n): + if n <= g[0]: + return p[0] - (g[0] - n) # extended with slope 1 + if n >= g[-1]: + return p[-1] + (n - g[-1]) + return float(np.interp(n, g, p)) + + sel = onsets[: self.n_onsets][::-1] # oldest first (Eq. 1 weights 1..n) + tempi = [] + for o in sel: + go = float(o) + n1, n2 = go - half, go + half + tempi.append((n2 - n1 + 1) / (phi_r(n2) - phi_r(n1) + 1)) + w = np.arange(1, len(tempi) + 1) + return float(np.sum(np.array(tempi) * w) / np.sum(w)) + + def alter(self, tracker, t: Optional[float]) -> None: + """Sec. 4.2: alter the score representation at the last processed frame.""" + if t is None: + self._consecutive = 0 + return + t = max(self.tempo_min, min(self.tempo_max, t)) + if t >= 1.0: + self._insert_pending = False + ps = tracker.tempo_last_score_frame() + if ps <= self._last_altered_ps + self.min_advance: + return + V, f2b, on = tracker.tempo_score() + if ps < 0 or ps + 2 >= len(V): + return + if self._consecutive >= 3: # at most 3 alterations in a row + self._consecutive = 0 + return + r = self._rng.uniform(0.0, 1.0) + if t > 1.0 and r > 1.0 / t: + V[ps + 1] = 0.5 * (V[ps + 1] + V[ps + 2]) + if on[ps + 2] and not on[ps + 1]: + f2b[ps + 1] = f2b[ps + 2] + elif not on[ps + 1]: + f2b[ps + 1] = 0.5 * (f2b[ps + 1] + f2b[ps + 2]) + on[ps + 1] |= on[ps + 2] + tracker.tempo_set_score( + np.delete(V, ps + 2, axis=0), + np.delete(f2b, ps + 2), + np.delete(on, ps + 2), + ) + self.n_deleted += 1 + self._consecutive += 1 + self._last_altered_ps = ps + elif t < 1.0 and (self._insert_pending or r > t): + if on[ps] or on[ps + 1]: # onset vectors are not duplicated + self._insert_pending = True + self.n_postponed += 1 + self._consecutive = 0 + return + self._insert_pending = False + tracker.tempo_set_score( + np.insert(V, ps + 1, 0.5 * (V[ps] + V[ps + 1]), axis=0), + np.insert(f2b, ps + 1, 0.5 * (f2b[ps] + f2b[ps + 1])), + np.insert(on, ps + 1, False), + ) + self.n_inserted += 1 + self._consecutive += 1 + self._last_altered_ps = ps + else: + self._consecutive = 0 + + # --- Standalone onset update (Eq. 1) --- + + def update_beat_period( + self, + performed_onset: float, + score_onset: float, + *args, + **kwargs, + ) -> None: + """Update beat period using the Eq. 1 weighted average of recent local tempi.""" + self.est_onset = performed_onset + self._history_perf.append(performed_onset) + self._history_score.append(score_onset) + self.prev_perf_onset = performed_onset + self.prev_score_onset = score_onset + + if len(self._history_perf) < 2: + return + + perf = np.array(self._history_perf) + score = np.array(self._history_score) + n = min(self.n_onsets, len(perf) - 1) + tempi = [] + for i in range(1, n + 1): + dp = perf[-1] - perf[-1 - i] + ds = score[-1] - score[-1 - i] + if dp > 0 and ds > 0: + tempi.append(dp / ds) + + if tempi: + tempi = tempi[::-1] + w = np.arange(1, len(tempi) + 1) + est_bp = float(np.sum(np.array(tempi) * w) / np.sum(w)) + base_bp = getattr(self, "init_beat_period", 0.5) + rel_t = base_bp / max(1e-4, est_bp) + rel_t = max(self.tempo_min, min(self.tempo_max, rel_t)) + self._current_tempo = rel_t + self.beat_period = base_bp / rel_t + + +# Aliases for backward compatibility +ArztWidmer2010TempoModel = ArztTempoModel + + +def make_tempo_model(tempo_model, **kwargs) -> Optional[TempoModel]: + """None, a ``TempoModel`` instance, or a class (built with ``kwargs``).""" + if tempo_model is None or isinstance(tempo_model, TempoModel): + return tempo_model + if isinstance(tempo_model, type) and issubclass(tempo_model, TempoModel): + return tempo_model(**kwargs) + raise TypeError( + f"tempo_model must be None, a TempoModel or its class, got {tempo_model!r}" + ) + + if __name__ == "__main__": pass diff --git a/tests/test_oltw_arzt.py b/tests/test_oltw_arzt.py index c22b21d..4c1aaf5 100644 --- a/tests/test_oltw_arzt.py +++ b/tests/test_oltw_arzt.py @@ -139,6 +139,55 @@ def test_distance_func(self): # __call__ now returns the float beat position self.assertTrue(isinstance(current_position, float)) + def test_tempo_model(self): + """Test OnlineTimeWarpingArztFrame with ArztTempoModel.""" + from matchmaker.utils.tempo_models import ArztTempoModel + + X, Y, path = generate_example_sequences( + lenX=50, + centers=3, + n_features=3, + maxreps=4, + minreps=1, + noise_scale=0.00, + random_state=RNG, + dtype=np.float32, + ) + score_positions = np.arange(X.shape[0], dtype=np.float32) + + # Baseline tracker without tempo model + oltw_no_tm = OnlineTimeWarpingArzt( + reference_features=X, + score_positions=score_positions, + ref_frame_to_beat=score_positions, + window_size=5, + step_size=3, + frame_rate=50, + ) + self.assertEqual(oltw_no_tm.current_relative_tempo, 1.0) + self.assertFalse(oltw_no_tm.backward_forward) + + # Tracker with ArztTempoModel class passed + oltw_tm = OnlineTimeWarpingArzt( + reference_features=X, + score_positions=score_positions, + ref_frame_to_beat=score_positions, + window_size=5, + step_size=3, + frame_rate=50, + tempo_model=ArztTempoModel, + ) + self.assertFalse(oltw_tm.backward_forward) + self.assertIsInstance(oltw_tm.tempo_model, ArztTempoModel) + + positions = [] + for obs in Y: + pos = oltw_tm(obs, 0.0) + self.assertIsInstance(pos, float) + positions.append(pos) + self.assertGreater(len(positions), 0) + self.assertIsInstance(oltw_tm.current_relative_tempo, float) + # Test that error is raised if incorrect name self.assertRaises( MatchmakerInvalidOptionError, diff --git a/tests/test_tempo_models.py b/tests/test_tempo_models.py index 131683a..578213e 100644 --- a/tests/test_tempo_models.py +++ b/tests/test_tempo_models.py @@ -5,105 +5,148 @@ """ import unittest + import numpy as np -from matchmaker.dp.oltw_arzt import ( - OnlineTimeWarpingArztFrame, - OnlineTimeWarpingArztTempoFrame, -) +from matchmaker.dp.oltw_arzt import OnlineTimeWarpingArztFrame from matchmaker.prob.skf import ( SwitchingKalmanFilterFollower, build_chord_sequence, build_spectral_templates, ) from matchmaker.registry import REGISTRY +from matchmaker.utils.tempo_models import ArztTempoModel from tests.utils import generate_example_sequences RNG = np.random.RandomState(1984) class TestArztTempoModel(unittest.TestCase): - """Tests for Arzt & Widmer (2010 SMC) Simple Tempo Model.""" + """Tests for Arzt, Widmer & Dixon (2008) OLTW with the Arzt & Widmer (2010) tempo model.""" def setUp(self): # Generate synthetic reference and performance sequences self.frame_rate = 50 - # 10 seconds of score audio -> 500 frames - n_frames = 500 + # 20 seconds of score audio -> 1000 frames + n_frames = 1000 n_features = 12 self.X = RNG.randn(n_frames, n_features).astype(np.float32) - # 1 beat per second -> 10 beats total - # Every 50 frames is 1 beat - self.ref_frame_to_beat = (np.arange(n_frames, dtype=np.float32) / self.frame_rate) - # Score positions at beats 0.0, 1.0, 2.0, ..., 9.0 - self.score_positions = np.arange(10, dtype=np.float32) + # 1 beat per second; every 50 frames is 1 beat + self.ref_frame_to_beat = np.arange(n_frames, dtype=np.float32) / self.frame_rate + # Score positions at beats 0.0, 1.0, 2.0, ..., 19.0 + self.score_positions = np.arange(20, dtype=np.float32) - def test_initialization(self): - tracker = OnlineTimeWarpingArztTempoFrame( + def _tracker(self, tempo_model=ArztTempoModel, **kwargs): + return OnlineTimeWarpingArztFrame( reference_features=self.X, score_positions=self.score_positions, ref_frame_to_beat=self.ref_frame_to_beat, frame_rate=self.frame_rate, - window_size=5, - step_size=3, + tempo_model=tempo_model, + **kwargs, ) + + def _perform(self, tracker, speed, n_steps): + """Feed the reference played at ``speed`` times the score tempo.""" + beats = [] + for i in range(n_steps): + idx = min(int(i * speed), len(self.X) - 1) + obs = self.X[idx] + RNG.randn(12).astype(np.float32) * 0.01 + beats.append(tracker(obs, i / self.frame_rate)) + return np.array(beats) + + def test_initialization(self): + tracker = self._tracker(window_size=5) self.assertTrue(tracker.use_tempo_model) + self.assertEqual(tracker.c, 250) self.assertEqual(tracker.current_relative_tempo, 1.0) - self.assertEqual(tracker._consecutive_alterations, 0) - self.assertEqual(len(tracker._backpointers), 0) - self.assertTrue(tracker.is_onset_frame[0]) - self.assertTrue(tracker.is_onset_frame[50]) + self.assertEqual(tracker.N, len(self.X)) + self.assertTrue(tracker.is_onset[0]) + self.assertTrue(tracker.is_onset[50]) + self.assertFalse(tracker.is_onset[25]) def test_run_steps(self): - tracker = OnlineTimeWarpingArztTempoFrame( - reference_features=self.X, - score_positions=self.score_positions, - ref_frame_to_beat=self.ref_frame_to_beat, - frame_rate=self.frame_rate, - window_size=5, - step_size=3, - ) - - # Feed 100 frames of observations - for i in range(100): - obs = self.X[min(i, len(self.X) - 1)] + RNG.randn(12).astype(np.float32) * 0.01 - perf_time = i / self.frame_rate - beat = tracker(obs, perf_time) - self.assertIsInstance(beat, float) - self.assertGreaterEqual(beat, 0.0) - - # Check that backpointers were recorded - self.assertEqual(len(tracker._backpointers), 100) + tracker = self._tracker() + beats = self._perform(tracker, 1.0, 200) + self.assertTrue(all(isinstance(b, float) for b in beats)) + # Played at the score tempo: the position follows the input closely + expected = np.arange(200) / self.frame_rate + self.assertLess(np.abs(beats[100:] - expected[100:]).max(), 0.2) path = tracker.alignment_path - self.assertEqual(path.shape[0], 2) - self.assertEqual(path.shape[1], 100) + self.assertEqual(path.shape, (2, 200)) + + def test_tempo_model_follows_faster_performance(self): + tracker = self._tracker() + beats = self._perform(tracker, 1.5, 400) + self.assertGreater(tracker.current_relative_tempo, 1.2) + self.assertGreater(tracker.n_deleted, 0) + self.assertEqual(tracker.n_inserted, 0) + self.assertLess(tracker.N, len(self.X)) + expected = 1.5 * np.arange(400) / self.frame_rate + self.assertLess(np.abs(beats[200:] - expected[200:]).max(), 0.5) + + def test_without_tempo_model(self): + tracker = self._tracker(tempo_model=None) + self._perform(tracker, 1.5, 200) + self.assertEqual(tracker.n_deleted + tracker.n_inserted, 0) + self.assertEqual(tracker.N, len(self.X)) + self.assertEqual(tracker.current_relative_tempo, 1.0) def test_reset(self): - tracker = OnlineTimeWarpingArztTempoFrame( - reference_features=self.X, - score_positions=self.score_positions, - ref_frame_to_beat=self.ref_frame_to_beat, - frame_rate=self.frame_rate, - ) - orig_len = tracker.N_ref - - # Run 50 frames - for i in range(50): - obs = self.X[i] - tracker(obs, i / self.frame_rate) + tracker = self._tracker() + self._perform(tracker, 1.5, 300) + self.assertLess(tracker.N, len(self.X)) - # Reset tracker.reset() self.assertEqual(tracker.current_index, 0) - self.assertEqual(tracker._current_frame, 0) self.assertEqual(tracker.current_relative_tempo, 1.0) - self.assertEqual(len(tracker._backpointers), 0) - self.assertEqual(tracker.N_ref, orig_len) + self.assertEqual(tracker.N, len(self.X)) + self.assertEqual(tracker.alignment_path.shape, (2, 0)) + np.testing.assert_array_equal(tracker.V, self.X) def test_registry_lookup(self): + from matchmaker.utils.tempo_models import ArztTempoModel, TempoModel + + self.assertTrue(issubclass(ArztTempoModel, TempoModel)) + spec = REGISTRY.method("audio", "arzt_tempo") - self.assertEqual(spec.cls_path, "matchmaker.dp:OnlineTimeWarpingArztTempoFrame") - self.assertTrue(spec.default_kwargs.get("use_tempo_model")) + self.assertEqual(spec.cls_path, "matchmaker.dp:OnlineTimeWarpingArztFrame") + self.assertIs(spec.default_kwargs.get("tempo_model"), ArztTempoModel) + self.assertEqual(spec.default_kwargs.get("processor"), "chroma") + + def test_arzt_tempo_model_standalone_onsets(self): + from matchmaker.utils.tempo_models import ( + ArztTempoModel, + ArztWidmer2010TempoModel, + make_tempo_model, + ) + + tm = ArztTempoModel(init_beat_period=0.5, n_onsets=5) + self.assertIs(ArztWidmer2010TempoModel, ArztTempoModel) + self.assertIs(make_tempo_model(tm), tm) + built = make_tempo_model(ArztTempoModel, init_beat_period=0.6) + self.assertIsInstance(built, ArztTempoModel) + self.assertEqual(built.init_beat_period, 0.6) + + # First onset + bp, est = tm(performed_onset=0.0, score_onset=0.0) + self.assertEqual(bp, 0.5) + self.assertEqual(est, 0.0) + + # Performance running 20% slower: performed IOI = 0.6s per score beat + for i in range(1, 10): + perf_t = i * 0.6 + score_b = float(i) + bp, est = tm(performed_onset=perf_t, score_onset=score_b) + + # Relative tempo should reflect slower tempo (t < 1.0, beat_period > 0.5) + self.assertLess(tm.current_tempo, 1.0) + self.assertGreater(tm.beat_period, 0.5) + + # Reset + tm.reset() + self.assertEqual(tm.counter, 0) + self.assertEqual(tm.current_tempo, 1.0) class TestRaphaelSwitchingStateSpace(unittest.TestCase): @@ -194,7 +237,9 @@ def test_kalman_tempo_adaptation(self): def test_registry_lookup(self): spec = REGISTRY.method("audio", "skf") - self.assertEqual(spec.cls_path, "matchmaker.prob.skf:SwitchingKalmanFilterFollower") + self.assertEqual( + spec.cls_path, "matchmaker.prob.skf:SwitchingKalmanFilterFollower" + ) if __name__ == "__main__": From 62ed7be68329c1c58460a1cf1a16ee2e01cac7f9 Mon Sep 17 00:00:00 2001 From: "jiyun.park" Date: Thu, 22 Oct 2026 11:47:38 +0200 Subject: [PATCH 11/11] Address PR review: validate tempo and unify OLTW directions Validate tempo values and score alteration boundaries. Share direction enums with Dixon and count actual refinement advances. Clarify backpointer fields, separate SKF regression tests, and fix audio test discovery on headless CI. --- matchmaker/dp/direction.py | 12 + matchmaker/dp/dtw_loop.pyx | 19 +- matchmaker/dp/oltw_arzt.py | 367 ++++++++++++++++++------------- matchmaker/dp/oltw_dixon.py | 9 +- matchmaker/prob/skf.py | 9 +- matchmaker/utils/tempo_models.py | 38 +++- tests/test_oltw_arzt.py | 191 ++++++++++++++++ tests/test_prob_skf.py | 129 +++++++++++ tests/test_tempo_models.py | 146 ++++-------- tests/test_utils_audio.py | 45 ++-- 10 files changed, 681 insertions(+), 284 deletions(-) create mode 100644 matchmaker/dp/direction.py create mode 100644 tests/test_prob_skf.py diff --git a/matchmaker/dp/direction.py b/matchmaker/dp/direction.py new file mode 100644 index 0000000..47b8da7 --- /dev/null +++ b/matchmaker/dp/direction.py @@ -0,0 +1,12 @@ +"""Shared reference/input directions for online time warping.""" + +from enum import IntEnum + + +class Direction(IntEnum): + """Axes advanced by a step; NONE marks a missing predecessor.""" + + NONE = -1 + REF = 0 + INPUT = 1 + BOTH = 2 diff --git a/matchmaker/dp/dtw_loop.pyx b/matchmaker/dp/dtw_loop.pyx index 94b918d..0a5ec59 100644 --- a/matchmaker/dp/dtw_loop.pyx +++ b/matchmaker/dp/dtw_loop.pyx @@ -8,6 +8,13 @@ cimport numpy as np import numpy as np +from matchmaker.dp.direction import Direction + +cdef signed char NO_DIRECTION = Direction.NONE +cdef signed char REF_DIRECTION = Direction.REF +cdef signed char INPUT_DIRECTION = Direction.INPUT +cdef signed char BOTH_DIRECTION = Direction.BOTH + # from numpy.math cimport INFINITY from libc.math cimport INFINITY @@ -49,7 +56,7 @@ def oltw_arzt_loop( Index of minimum costs of current window. backpointers: numpy.ndarray or None Optional array of shape (window_end - window_start,) to record backpointers - (1=vertical, 2=horizontal, 3=diagonal). + Values follow Direction: NONE=-1, REF=0, INPUT=1, BOTH=2. Returns ------- @@ -72,6 +79,8 @@ def oltw_arzt_loop( # global_cost_matrix[1, 1] = sum(window_cost) min_costs = global_cost_matrix[1, 1] min_index = 0 + if backpointers is not None: + backpointers[0] = NO_DIRECTION while score_index < window_end: @@ -90,11 +99,11 @@ def oltw_arzt_loop( if backpointers is not None: if dist1 <= dist2 and dist1 <= dist3: - backpointers[idx] = 1 + backpointers[idx] = REF_DIRECTION elif dist2 <= dist3: - backpointers[idx] = 2 + backpointers[idx] = INPUT_DIRECTION else: - backpointers[idx] = 3 + backpointers[idx] = BOTH_DIRECTION norm_cost = min_dist / (input_index + score_index + 1.0) @@ -108,6 +117,8 @@ def oltw_arzt_loop( update_cost_matrix(global_cost_matrix, N) + # Row arithmetic in the tempo model needs an ndarray, not a Cython memoryview. + # np.asarray shares the existing buffer without copying its data. return np.asarray(global_cost_matrix), min_index, min_costs @cython.boundscheck(False) diff --git a/matchmaker/dp/oltw_arzt.py b/matchmaker/dp/oltw_arzt.py index 63558af..d309d21 100644 --- a/matchmaker/dp/oltw_arzt.py +++ b/matchmaker/dp/oltw_arzt.py @@ -15,13 +15,14 @@ """ import time -from typing import Any, Callable, Dict, Generator, Optional, Tuple, Union +from typing import Any, Callable, Dict, Generator, NamedTuple, Optional, Tuple, Union import numpy as np import scipy.spatial.distance from numpy.typing import NDArray from matchmaker.base import OnlineAlignment +from matchmaker.dp.direction import Direction from matchmaker.dp.dtw_loop import oltw_arzt_loop from matchmaker.features.audio import FRAME_RATE from matchmaker.io.audio import QUEUE_TIMEOUT @@ -43,11 +44,12 @@ START_WINDOW_SIZE: Union[float, int] = 0.1 WINDOW_SIZE: int = 10 -_ADV_IN, _ADV_SC, _BOTH = ( - 0, - 1, - 2, -) # Dixon (2005) GetInc: Row (input), Column (score), Both + +class BackpointerColumn(NamedTuple): + """Predecessor directions in one input frame's reference window.""" + + window_start: int + predecessor_directions: NDArray[np.int8] # --------------------------------------------------------------------------- @@ -206,7 +208,8 @@ def reset(self) -> None: self._backpointers = {} self._last_window_end = 0 self._min_x_needed = 0 - self._run_count, self._run_prev = 0, None + self.run_count_ref = 0 + self.run_count_input = 0 if self.backward_forward: K = 4 * self._window_size self._cost_band = np.full( @@ -257,27 +260,26 @@ def tempo_last_input(self) -> int: def tempo_current_cell(self) -> Tuple[int, int]: return self.input_index - 1, self._current_frame - def tempo_backward_path(self, start_x: int, start_y: int): - x, y = start_x, start_y - while x > 0 and y > 0: - entry = self._backpointers.get(x) - if entry is None: + def tempo_backward_path(self, start_input_index: int, start_ref_index: int): + input_index, ref_index = start_input_index, start_ref_index + while input_index > 0 and ref_index > 0: + column = self._backpointers.get(input_index) + if column is None: break - ws, codes = entry - offset = y - ws - if offset < 0 or offset >= len(codes): + offset = ref_index - column.window_start + if offset < 0 or offset >= len(column.predecessor_directions): break - code = codes[offset] - if code == 1: # vertical: (x, y-1) - y -= 1 - elif code == 2: # horizontal: (x-1, y) - x -= 1 - elif code == 3: # diagonal: (x-1, y-1) - x -= 1 - y -= 1 + direction = column.predecessor_directions[offset] + if direction == Direction.REF: + ref_index -= 1 + elif direction == Direction.INPUT: + input_index -= 1 + elif direction == Direction.BOTH: + input_index -= 1 + ref_index -= 1 else: break - yield x, y + yield input_index, ref_index def tempo_note_path_end(self, x: int, y: int) -> None: self._min_x_needed = x if not self._min_x_needed else min(self._min_x_needed, x) @@ -289,144 +291,211 @@ def tempo_score(self): return self.reference_features, self._ref_frame_to_beat, self.is_onset def tempo_set_score(self, features, frame_to_beat, is_onset) -> None: - n = len(features) - ps = self._current_frame - G = self.global_cost_matrix # row k + 1 = score frame k - if n < self.N_ref: - r = ps + 2 - self.global_cost_matrix = np.ascontiguousarray(np.delete(G, r + 1, axis=0)) - for x, (s, c) in list(self._backpointers.items()): - if s > r: - self._backpointers[x] = (s - 1, c) - elif r - s < len(c): - self._backpointers[x] = (s, np.delete(c, r - s)) + n_frames = len(features) + if len(frame_to_beat) != n_frames or len(is_onset) != n_frames: + raise ValueError("Score features, beat mapping and onset flags must agree.") + if abs(n_frames - self.N_ref) > 1: + raise ValueError("A tempo alteration must insert or delete a single frame.") + + cost_matrix = self.global_cost_matrix # row k + 1 = reference frame k + if n_frames < self.N_ref: + deleted_ref_index = self._current_frame + 2 + if not 0 <= deleted_ref_index < self.N_ref: + raise ValueError("The reference frame to delete is out of bounds.") + self.global_cost_matrix = np.ascontiguousarray( + np.delete(cost_matrix, deleted_ref_index + 1, axis=0) + ) + for input_index, column in list(self._backpointers.items()): + if column.window_start > deleted_ref_index: + self._backpointers[input_index] = column._replace( + window_start=column.window_start - 1 + ) + else: + offset = deleted_ref_index - column.window_start + if offset < len(column.predecessor_directions): + self._backpointers[input_index] = column._replace( + predecessor_directions=np.delete( + column.predecessor_directions, offset + ) + ) if self.backward_forward: - self._band_alter(r, -1) - elif n > self.N_ref: - r = ps + 1 - if r + 1 < len(G): - interp = 0.5 * (G[r] + G[r + 1]) + self._band_alter(deleted_ref_index, -1) + elif n_frames > self.N_ref: + inserted_ref_index = self._current_frame + 1 + if not 0 <= inserted_ref_index <= self.N_ref: + raise ValueError("The reference insertion position is out of bounds.") + if inserted_ref_index + 1 < len(cost_matrix): + interpolated_cost = 0.5 * ( + cost_matrix[inserted_ref_index] + + cost_matrix[inserted_ref_index + 1] + ) else: - interp = np.full(2, np.inf, dtype=np.float32) + interpolated_cost = np.full(2, np.inf, dtype=np.float32) self.global_cost_matrix = np.ascontiguousarray( - np.insert(G, r + 1, interp, axis=0) + np.insert( + cost_matrix, inserted_ref_index + 1, interpolated_cost, axis=0 + ) ) - for x, (s, c) in list(self._backpointers.items()): - if s >= r: - self._backpointers[x] = (s + 1, c) - elif r - s < len(c): - self._backpointers[x] = (s, np.insert(c, r - s, 1).astype(np.int8)) + for input_index, column in list(self._backpointers.items()): + if column.window_start >= inserted_ref_index: + self._backpointers[input_index] = column._replace( + window_start=column.window_start + 1 + ) + else: + offset = inserted_ref_index - column.window_start + if offset < len(column.predecessor_directions): + self._backpointers[input_index] = column._replace( + predecessor_directions=np.insert( + column.predecessor_directions, offset, Direction.REF + ).astype(np.int8) + ) if self.backward_forward: - self._band_alter(r, 1) + self._band_alter(inserted_ref_index, 1) self.reference_features = features self._ref_frame_to_beat = frame_to_beat self.is_onset = is_onset - self.N_ref = n + self.N_ref = n_frames # Strategy 1 (Arzt, Widmer & Dixon 2008, Sec. 3.2). - def _band_alter(self, r, sign): - st, ln, band = self._cost_start, self._cost_len, self._cost_band - for k in range(len(band)): - s, n = st[k], ln[k] - if sign < 0: - if s > r: - st[k] = s - 1 - elif r - s < n: - band[k, r - s : n - 1] = band[k, r - s + 1 : n].copy() - band[k, n - 1] = np.inf - ln[k] = n - 1 + def _band_alter(self, ref_index, change): + for slot in range(len(self._cost_band)): + window_start = self._cost_start[slot] + window_length = self._cost_len[slot] + offset = ref_index - window_start + if change < 0: + if window_start > ref_index: + self._cost_start[slot] = window_start - 1 + elif offset < window_length: + self._cost_band[slot, offset : window_length - 1] = self._cost_band[ + slot, offset + 1 : window_length + ].copy() + self._cost_band[slot, window_length - 1] = np.inf + self._cost_len[slot] = window_length - 1 else: - if s >= r: - st[k] = s + 1 - elif r - s < n: - m = min(n + 1, band.shape[1]) - band[k, r - s + 1 : m] = band[k, r - s : m - 1].copy() - band[k, r - s] = np.inf - ln[k] = m - - def _cost_row(self, y, x0, x1): - xs = np.arange(x0, x1 + 1) - K = len(self._cost_band) - sl = xs % K - k = y - self._cost_start[sl] - ok = ( - (xs > self.input_index - K) - & (xs <= self.input_index) - & (k >= 0) - & (k < self._cost_len[sl]) + if window_start >= ref_index: + self._cost_start[slot] = window_start + 1 + elif offset < window_length: + new_length = min(window_length + 1, self._cost_band.shape[1]) + self._cost_band[slot, offset + 1 : new_length] = self._cost_band[ + slot, offset : new_length - 1 + ].copy() + self._cost_band[slot, offset] = np.inf + self._cost_len[slot] = new_length + + def _cost_row(self, ref_index, input_start, input_end): + input_indices = np.arange(input_start, input_end + 1) + capacity = len(self._cost_band) + slots = input_indices % capacity + offsets = ref_index - self._cost_start[slots] + valid = ( + (input_indices > self.input_index - capacity) + & (input_indices <= self.input_index) + & (offsets >= 0) + & (offsets < self._cost_len[slots]) ) - out = np.full(len(xs), np.inf) - out[ok] = self._cost_band[sl[ok], k[ok]] - return out - - def _cost_col(self, x, y0, y1): - out = np.full(y1 - y0 + 1, np.inf) - K = len(self._cost_band) - if x > self.input_index or x <= self.input_index - K or x < 0: - return out - sl = x % K - s, n = self._cost_start[sl], self._cost_len[sl] - a, b = max(y0, s), min(y1, s + n - 1) - if a <= b: - out[a - y0 : b - y0 + 1] = self._cost_band[sl, a - s : b - s + 1] - return out - - def _cost(self, x, y): - return float(self._cost_col(x, y, y)[0]) - - def _update_runs(self, move: int) -> None: - if move == self._run_prev: - self._run_count += 1 - else: - self._run_count, self._run_prev = 1, move - - def _get_inc(self, x, y): - M = self.max_run_count - if M is not None and self._run_count > M: - return _ADV_SC if self._run_prev == _ADV_IN else _ADV_IN - c = self._window_size - best = self._cost(x, y) / (x + y + 1.0) - bx, by = x, y - x0, y0 = max(0, x - c + 1), max(0, y - c + 1) - r = self._cost_row(y, x0, x) / (np.arange(x0, x + 1) + y + 1.0) - k = int(np.argmin(r)) - if r[k] < best: - best, bx, by = r[k], x0 + k, y - cl = self._cost_col(x, y0, y) / (x + np.arange(y0, y + 1) + 1.0) - k = int(np.argmin(cl)) - if cl[k] < best: - best, bx, by = cl[k], x, y0 + k - if bx < x: - return _ADV_SC - if by < y: - return _ADV_IN - return _BOTH - - def _strategy1(self, T, y): - b = 50 if self._n_bf % 5 == 4 else 10 + costs = np.full(len(input_indices), np.inf) + costs[valid] = self._cost_band[slots[valid], offsets[valid]] + return costs + + def _cost_col(self, input_index, ref_start, ref_end): + costs = np.full(ref_end - ref_start + 1, np.inf) + capacity = len(self._cost_band) + if ( + input_index > self.input_index + or input_index <= self.input_index - capacity + or input_index < 0 + ): + return costs + slot = input_index % capacity + window_start = self._cost_start[slot] + window_length = self._cost_len[slot] + overlap_start = max(ref_start, window_start) + overlap_end = min(ref_end, window_start + window_length - 1) + if overlap_start <= overlap_end: + costs[overlap_start - ref_start : overlap_end - ref_start + 1] = ( + self._cost_band[ + slot, overlap_start - window_start : overlap_end - window_start + 1 + ] + ) + return costs + + def _cost(self, input_index, ref_index): + return float(self._cost_col(input_index, ref_index, ref_index)[0]) + + def update_run_counts(self, adv_ref, adv_input): + """Count actual single-axis advances, following Dixon's convention.""" + if adv_ref and adv_input: + self.run_count_ref = 0 + self.run_count_input = 0 + elif adv_ref: + self.run_count_ref += 1 + self.run_count_input = 0 + elif adv_input: + self.run_count_input += 1 + self.run_count_ref = 0 + + def get_inc(self, input_index, ref_index) -> Direction: + if self.max_run_count is not None: + if self.run_count_ref >= self.max_run_count: + return Direction.INPUT + if self.run_count_input >= self.max_run_count: + return Direction.REF + window_size = self._window_size + best_cost = self._cost(input_index, ref_index) / (input_index + ref_index + 1.0) + best_input, best_ref = input_index, ref_index + input_start = max(0, input_index - window_size + 1) + ref_start = max(0, ref_index - window_size + 1) + row_costs = self._cost_row(ref_index, input_start, input_index) / ( + np.arange(input_start, input_index + 1) + ref_index + 1.0 + ) + best_offset = int(np.argmin(row_costs)) + if row_costs[best_offset] < best_cost: + best_cost = row_costs[best_offset] + best_input, best_ref = input_start + best_offset, ref_index + column_costs = self._cost_col(input_index, ref_start, ref_index) / ( + input_index + np.arange(ref_start, ref_index + 1) + 1.0 + ) + best_offset = int(np.argmin(column_costs)) + if column_costs[best_offset] < best_cost: + best_input, best_ref = input_index, ref_start + best_offset + if best_input < input_index: + return Direction.REF + if best_ref < ref_index: + return Direction.INPUT + return Direction.BOTH + + def _strategy1(self, input_end, ref_index): + backward_distance = 50 if self._n_bf % 5 == 4 else 10 self._n_bf += 1 - x, target = T, y - b - for x, y in self.tempo_backward_path(T, y): - if y <= target: + input_index, target_ref = input_end, ref_index - backward_distance + for input_index, ref_index in self.tempo_backward_path(input_end, ref_index): + if ref_index <= target_ref: break - self.tempo_note_path_end(x, y) - J = self._last_window_end - 1 - self._run_count, self._run_prev = 0, None - while x < T and y < J: - if self._get_inc(x, y) != _ADV_SC: - if not np.isfinite(self._cost(x + 1, y)) and not np.isfinite( - self._cost(x + 1, y + 1) - ): + self.tempo_note_path_end(input_index, ref_index) + ref_end = self._last_window_end - 1 + self.run_count_ref = 0 + self.run_count_input = 0 + while input_index < input_end and ref_index < ref_end: + previous_input, previous_ref = input_index, ref_index + if self.get_inc(input_index, ref_index) != Direction.REF: + if not np.isfinite( + self._cost(input_index + 1, ref_index) + ) and not np.isfinite(self._cost(input_index + 1, ref_index + 1)): break - x += 1 - if x <= T and y < J and self._get_inc(x, y) != _ADV_IN: - y += 1 - if not np.isfinite(self._cost(x, y)): + input_index += 1 + if ( + input_index <= input_end + and ref_index < ref_end + and self.get_inc(input_index, ref_index) != Direction.INPUT + ): + ref_index += 1 + if not np.isfinite(self._cost(input_index, ref_index)): break - if self.max_run_count is not None: - self._update_runs(self._get_inc(x, y)) - return y + self.update_run_counts( + ref_index > previous_ref, input_index > previous_input + ) + return ref_index @property def window_index(self) -> int: @@ -525,7 +594,9 @@ def step(self, input_features: NDArray[np.float32]) -> None: store = (self.tempo_model is not None) or self.backward_forward backpointers = None if store: - backpointers = np.zeros(window_end - window_start, dtype=np.int8) + backpointers = np.full( + window_end - window_start, Direction.NONE, dtype=np.int8 + ) self.global_cost_matrix, min_index, min_costs = oltw_arzt_loop( global_cost_matrix=self.global_cost_matrix, @@ -547,7 +618,9 @@ def step(self, input_features: NDArray[np.float32]) -> None: self._cost_start[k], self._cost_len[k] = window_start, n if store: - self._backpointers[self.input_index] = (window_start, backpointers) + self._backpointers[self.input_index] = BackpointerColumn( + window_start=window_start, predecessor_directions=backpointers + ) self._last_window_end = window_end if (self.input_index + 1) % self._window_size == 0: keep = min(self.input_index - 2 * self._window_size, self._min_x_needed) diff --git a/matchmaker/dp/oltw_dixon.py b/matchmaker/dp/oltw_dixon.py index b93affa..e7efe24 100644 --- a/matchmaker/dp/oltw_dixon.py +++ b/matchmaker/dp/oltw_dixon.py @@ -16,7 +16,6 @@ """ import time -from enum import IntEnum from typing import Any, Dict, Generator, Optional, Tuple import numpy as np @@ -24,18 +23,12 @@ from numpy.typing import NDArray from matchmaker.base import OnlineAlignment +from matchmaker.dp.direction import Direction from matchmaker.features.audio import FRAME_RATE from matchmaker.io.audio import QUEUE_TIMEOUT from matchmaker.io.queue import RECVQueue from matchmaker.utils.misc import set_latency_stats - -class Direction(IntEnum): - REF = 0 - INPUT = 1 - BOTH = 2 - - MAX_RUN_COUNT: int = 3 WINDOW_SIZE = 10 # seconds for frame, events for event diff --git a/matchmaker/prob/skf.py b/matchmaker/prob/skf.py index 4d0fe76..c9ef791 100644 --- a/matchmaker/prob/skf.py +++ b/matchmaker/prob/skf.py @@ -62,8 +62,13 @@ def build_chord_sequence(note_array: np.ndarray): """ onsets = note_array["onset_beat"] pitches = note_array["pitch"] - # score time in quarters is independent of the time signature's beat unit - quarters = note_array["onset_quarter"] if "onset_quarter" in note_array.dtype.names else onsets + # Quarter-note positions are independent of the time signature's beat unit. + # Legacy arrays without onset_quarter are assumed to use quarter-note beats. + quarters = ( + note_array["onset_quarter"] + if "onset_quarter" in note_array.dtype.names + else onsets + ) unique_onsets, first = np.unique(onsets, return_index=True) K = len(unique_onsets) diff --git a/matchmaker/utils/tempo_models.py b/matchmaker/utils/tempo_models.py index 0a4478e..a64fc0b 100644 --- a/matchmaker/utils/tempo_models.py +++ b/matchmaker/utils/tempo_models.py @@ -611,6 +611,16 @@ def __init__( init_score_onset=init_score_onset, score_onsets=score_onsets, ) + if not np.isfinite(init_beat_period) or init_beat_period <= 0: + raise ValueError("init_beat_period must be finite and positive.") + if ( + not np.isfinite(tempo_min) + or not np.isfinite(tempo_max) + or not 0 < tempo_min <= tempo_max + ): + raise ValueError( + "Tempo bounds must be finite and satisfy 0 < tempo_min <= tempo_max." + ) self.init_beat_period = init_beat_period self.window_size = window_size self.min_past = min_past @@ -669,12 +679,18 @@ def current_tempo(self) -> float: def before_step(self, tracker) -> None: """Estimate the tempo and alter the tracker's score representation.""" - t = self.relative_tempo(tracker) - if t is not None: - self._current_tempo = t - base_bp = getattr(self, "init_beat_period", 0.5) - self.beat_period = base_bp / max(1e-4, t) - self.alter(tracker, t) + relative_tempo = self.relative_tempo(tracker) + if relative_tempo is not None: + if not np.isfinite(relative_tempo) or relative_tempo <= 0: + relative_tempo = None + else: + beat_period = self.init_beat_period / max(1e-4, relative_tempo) + if not np.isfinite(beat_period) or beat_period <= 0: + relative_tempo = None + else: + self._current_tempo = relative_tempo + self.beat_period = beat_period + self.alter(tracker, relative_tempo) def relative_tempo(self, tracker) -> Optional[float]: """Sec. 4.1: relative tempo from the rectified backward path.""" @@ -726,7 +742,7 @@ def phi_r(n): def alter(self, tracker, t: Optional[float]) -> None: """Sec. 4.2: alter the score representation at the last processed frame.""" - if t is None: + if t is None or not np.isfinite(t) or t <= 0: self._consecutive = 0 return t = max(self.tempo_min, min(self.tempo_max, t)) @@ -785,6 +801,8 @@ def update_beat_period( **kwargs, ) -> None: """Update beat period using the Eq. 1 weighted average of recent local tempi.""" + if not np.isfinite(performed_onset) or not np.isfinite(score_onset): + return self.est_onset = performed_onset self._history_perf.append(performed_onset) self._history_score.append(score_onset) @@ -811,8 +829,10 @@ def update_beat_period( base_bp = getattr(self, "init_beat_period", 0.5) rel_t = base_bp / max(1e-4, est_bp) rel_t = max(self.tempo_min, min(self.tempo_max, rel_t)) - self._current_tempo = rel_t - self.beat_period = base_bp / rel_t + beat_period = base_bp / rel_t + if np.isfinite(beat_period) and beat_period > 0: + self._current_tempo = rel_t + self.beat_period = beat_period # Aliases for backward compatibility diff --git a/tests/test_oltw_arzt.py b/tests/test_oltw_arzt.py index 4c1aaf5..be93408 100644 --- a/tests/test_oltw_arzt.py +++ b/tests/test_oltw_arzt.py @@ -5,11 +5,17 @@ """ import unittest +from unittest.mock import patch import numpy as np from scipy.spatial import distance as sp_distance +from matchmaker.dp.direction import Direction +from matchmaker.dp.dtw_loop import oltw_arzt_loop +from matchmaker.dp.oltw_arzt import BackpointerColumn from matchmaker.dp.oltw_arzt import OnlineTimeWarpingArztFrame as OnlineTimeWarpingArzt +from matchmaker.dp.oltw_dixon import Direction as DixonDirection +from matchmaker.dp.oltw_dixon import OnlineTimeWarpingDixonFrame from matchmaker.utils import ( CYTHONIZED_METRICS_W_ARGUMENTS, CYTHONIZED_METRICS_WO_ARGUMENTS, @@ -18,6 +24,7 @@ MatchmakerInvalidOptionError, MatchmakerInvalidParameterTypeError, ) +from matchmaker.utils.tempo_models import ArztTempoModel from tests.utils import generate_example_sequences RNG = np.random.RandomState(1984) @@ -216,3 +223,187 @@ def test_tempo_model(self): # scipy metrics may not guarantee correct alignments; type check only current_position = oltw(Y[0], 0.0) self.assertTrue(isinstance(current_position, float)) + + +class TestArztTempoIntegration(unittest.TestCase): + def make_follower(self, **kwargs): + return OnlineTimeWarpingArzt( + reference_features=np.eye(8, dtype=np.float32), + score_positions=np.arange(8, dtype=np.float32), + ref_frame_to_beat=np.arange(8, dtype=np.float32), + frame_rate=1, + tempo_model=ArztTempoModel, + **kwargs, + ) + + def test_cython_predecessors_use_shared_directions(self): + self.assertIs(DixonDirection, Direction) + for direction, predecessor in ( + (Direction.REF, (1, 1)), + (Direction.INPUT, (2, 0)), + (Direction.BOTH, (1, 0)), + ): + with self.subTest(direction=direction): + costs = np.full((4, 2), 10.0, dtype=np.float32) + costs[predecessor] = 0.0 + directions = np.full(1, Direction.NONE, dtype=np.int8) + result, _, _ = oltw_arzt_loop( + costs, + np.ones(1, dtype=np.float32), + 1, + 2, + 1, + np.inf, + 0, + directions, + ) + self.assertEqual(directions[0], direction) + self.assertIsInstance(result, np.ndarray) + self.assertTrue(np.shares_memory(costs, result)) + + def test_origin_has_no_predecessor(self): + directions = np.zeros(2, dtype=np.int8) + oltw_arzt_loop( + np.full((3, 2), np.inf, dtype=np.float32), + np.ones(2, dtype=np.float32), + 0, + 2, + 0, + np.inf, + 0, + directions, + ) + self.assertEqual(directions[0], Direction.NONE) + + def test_backward_path_follows_each_direction_and_stops_at_none(self): + follower = self.make_follower() + follower._backpointers = { + 3: BackpointerColumn( + 0, + np.array( + [ + Direction.NONE, + Direction.NONE, + Direction.INPUT, + Direction.REF, + ], + dtype=np.int8, + ), + ), + 2: BackpointerColumn( + 0, + np.array( + [ + Direction.NONE, + Direction.NONE, + Direction.BOTH, + ], + dtype=np.int8, + ), + ), + 1: BackpointerColumn(0, np.full(2, Direction.NONE, dtype=np.int8)), + } + self.assertEqual( + list(follower.tempo_backward_path(3, 3)), [(3, 2), (2, 2), (1, 1)] + ) + self.assertEqual(list(follower.tempo_backward_path(4, 3)), []) + self.assertEqual(list(follower.tempo_backward_path(3, 5)), []) + + def test_last_deletable_frame_updates_costs_and_predecessors(self): + follower = self.make_follower() + follower._current_frame = 5 + follower.global_cost_matrix[:] = np.arange(18).reshape(9, 2) + expected_costs = follower.global_cost_matrix[:-1].copy() + follower._backpointers[0] = BackpointerColumn( + 0, np.full(8, Direction.BOTH, dtype=np.int8) + ) + follower.tempo_set_score( + follower.V[:-1], follower._ref_frame_to_beat[:-1], follower.is_onset[:-1] + ) + self.assertEqual(follower.N_ref, 7) + np.testing.assert_array_equal(follower.global_cost_matrix, expected_costs) + self.assertEqual(len(follower._backpointers[0].predecessor_directions), 7) + + def test_insertion_interpolates_costs_and_adds_ref_predecessor(self): + follower = self.make_follower() + follower._current_frame = 2 + follower.global_cost_matrix[:] = np.arange(18).reshape(9, 2) + expected_row = 0.5 * ( + follower.global_cost_matrix[3] + follower.global_cost_matrix[4] + ) + follower._backpointers[0] = BackpointerColumn( + 0, np.full(8, Direction.BOTH, dtype=np.int8) + ) + follower.tempo_set_score( + np.insert(follower.V, 3, follower.V[2], axis=0), + np.insert(follower._ref_frame_to_beat, 3, 2.5), + np.insert(follower.is_onset, 3, False), + ) + np.testing.assert_array_equal(follower.global_cost_matrix[4], expected_row) + self.assertEqual(follower.N_ref, 9) + self.assertEqual( + follower._backpointers[0].predecessor_directions[3], Direction.REF + ) + + def test_invalid_alterations_leave_score_and_costs_unchanged(self): + for current_frame, size_change in ((6, -1), (8, 1), (2, -2)): + with self.subTest(current_frame=current_frame, size_change=size_change): + follower = self.make_follower() + follower._current_frame = current_frame + costs = follower.global_cost_matrix.copy() + if size_change < 0: + features = follower.V[:size_change] + beats = follower._ref_frame_to_beat[:size_change] + onsets = follower.is_onset[:size_change] + else: + features = np.insert(follower.V, 0, follower.V[0], axis=0) + beats = np.insert(follower._ref_frame_to_beat, 0, 0.0) + onsets = np.insert(follower.is_onset, 0, False) + with self.assertRaises(ValueError): + follower.tempo_set_score(features, beats, onsets) + self.assertEqual(follower.N_ref, 8) + np.testing.assert_array_equal(follower.V, np.eye(8, dtype=np.float32)) + np.testing.assert_array_equal(follower.global_cost_matrix, costs) + + def test_run_counts_and_limit_follow_dixon(self): + follower = self.make_follower(backward_forward=True, max_run_count=3) + dixon = OnlineTimeWarpingDixonFrame( + reference_features=follower.V, + score_positions=follower.score_positions, + ref_frame_to_beat=follower._ref_frame_to_beat, + frame_rate=1, + window_size=1, + ) + for adv_ref, adv_input in ( + (True, False), + (True, False), + (True, False), + (True, True), + (False, True), + (False, True), + (False, True), + ): + follower.update_run_counts(adv_ref, adv_input) + dixon.update_run_counts(adv_ref, adv_input) + self.assertEqual(follower.run_count_ref, dixon.run_count_ref) + self.assertEqual(follower.run_count_input, dixon.run_count_input) + if follower.run_count_ref == 3: + self.assertEqual(follower.get_inc(3, 3), Direction.INPUT) + if follower.run_count_input == 3: + self.assertEqual(follower.get_inc(3, 3), Direction.REF) + follower.reset() + self.assertEqual((follower.run_count_ref, follower.run_count_input), (0, 0)) + + def test_refinement_counts_actual_moves_without_extra_direction_lookup(self): + follower = self.make_follower(backward_forward=True, max_run_count=3) + follower._last_window_end = 8 + with ( + patch.object(follower, "tempo_backward_path", return_value=iter([(1, 1)])), + patch.object( + follower, "get_inc", side_effect=[Direction.BOTH, Direction.INPUT] + ) as get_inc, + patch.object(follower, "_cost", return_value=1.0), + ): + follower._strategy1(2, 2) + self.assertEqual(get_inc.call_count, 2) + self.assertEqual((follower.run_count_ref, follower.run_count_input), (0, 1)) diff --git a/tests/test_prob_skf.py b/tests/test_prob_skf.py new file mode 100644 index 0000000..70f0293 --- /dev/null +++ b/tests/test_prob_skf.py @@ -0,0 +1,129 @@ +"""Tests for the switching Kalman filter audio follower.""" + +import unittest + +import numpy as np +import partitura as pt + +from matchmaker.prob.skf import ( + SwitchingKalmanFilterFollower, + build_chord_sequence, + build_spectral_templates, +) +from matchmaker.registry import REGISTRY + + +class TestSwitchingKalmanFilterFollower(unittest.TestCase): + """Tests for Raphael & Jiang (2020 ISMIR) Switching State-Space Model.""" + + def setUp(self): + # Create a simple synthetic note_array + self.note_array = np.zeros( + 4, + dtype=[ + ("pitch", "i4"), + ("onset_beat", "f4"), + ("duration_beat", "f4"), + ], + ) + self.note_array["pitch"] = [60, 64, 67, 72] # C major chord arpeggio + self.note_array["onset_beat"] = [0.0, 1.0, 2.0, 3.0] + self.note_array["duration_beat"] = [1.0, 1.0, 1.0, 1.0] + + def test_build_chord_sequence(self): + chords, lengths, onset_beats = build_chord_sequence(self.note_array) + self.assertEqual(len(chords), 4) + self.assertEqual(chords[0], [60]) + self.assertEqual(chords[1], [64]) + np.testing.assert_array_equal(onset_beats, [0.0, 1.0, 2.0, 3.0]) + # 1 beat = 0.25 whole note + np.testing.assert_allclose(lengths, [0.25, 0.25, 0.25, 0.25]) + + def test_chord_lengths_are_independent_of_time_signature(self): + for beats, beat_type, expected_beats in ( + (4, 4, [0.0, 1.0, 2.0, 3.0]), + (6, 8, [0.0, 2.0, 4.0, 6.0]), + (2, 2, [0.0, 0.5, 1.0, 1.5]), + ): + with self.subTest(time_signature=(beats, beat_type)): + part = pt.score.Part("P1", quarter_duration=480) + part.add(pt.score.TimeSignature(beats, beat_type), start=0) + for note_index in range(4): + part.add( + pt.score.Note(id=f"n{note_index}", step="C", octave=4), + start=note_index * 480, + end=(note_index + 1) * 480, + ) + _, lengths, onset_beats = build_chord_sequence(part.note_array()) + np.testing.assert_allclose(lengths, [0.25] * 4) + np.testing.assert_allclose(onset_beats, expected_beats) + + def test_build_spectral_templates(self): + chords, _, _ = build_chord_sequence(self.note_array) + templates = build_spectral_templates(chords, n_fft=512, sample_rate=8000) + self.assertEqual(templates.shape[0], 4) + self.assertEqual(templates.shape[1], 257) + # Each template row should sum to 1 + np.testing.assert_allclose(templates.sum(axis=1), np.ones(4), atol=1e-5) + + def test_initialization_and_step(self): + follower = SwitchingKalmanFilterFollower( + reference_features=self.note_array, + tempo=120.0, + sample_rate=8000, + n_fft=512, + hop_length=128, + max_hypotheses=50, + ) + self.assertEqual(follower.K, 4) + self.assertEqual(follower.current_index, 0) + # Initial tempo: 240 / 120 = 2.0 s/whole-note + self.assertAlmostEqual(follower.init_tempo, 2.0) + + # Process 20 frames + delta = 128 / 8000.0 + for i in range(20): + dummy_spectrum = np.random.rand(257).astype(np.float32) + pos = follower(dummy_spectrum, perf_time=i * delta) + self.assertIsInstance(pos, float) + self.assertGreaterEqual(pos, 0.0) + self.assertLessEqual(pos, 4.0) + + # Verify hypotheses count is bounded by max_hypotheses + self.assertLessEqual(len(follower.hypotheses), 50) + path = follower.alignment_path + self.assertEqual(path.shape, (2, 20)) + + def test_kalman_tempo_adaptation(self): + follower = SwitchingKalmanFilterFollower( + reference_features=self.note_array, + tempo=120.0, + sample_rate=8000, + n_fft=512, + hop_length=128, + ) + # Test transition Kalman update directly: + # Expected duration = length * tempo = 0.25 * 2.0 = 0.5s. + # Suppose performer played chord 0 in 0.25s (twice as fast, tempo should decrease s/wn): + fast_age = int(0.25 / follower.delta) + mu_new, var_new = follower._kalman_update( + mu_t_prev=2.0, var_t_prev=0.04, k_prev=0, age_prev=fast_age + ) + self.assertLess(mu_new, 2.0) + + # Suppose performer played chord 0 in 1.0s (twice as slow, tempo should increase s/wn): + slow_age = int(1.0 / follower.delta) + mu_new_slow, var_new_slow = follower._kalman_update( + mu_t_prev=2.0, var_t_prev=0.04, k_prev=0, age_prev=slow_age + ) + self.assertGreater(mu_new_slow, 2.0) + + def test_registry_lookup(self): + spec = REGISTRY.method("audio", "skf") + self.assertEqual( + spec.cls_path, "matchmaker.prob.skf:SwitchingKalmanFilterFollower" + ) + + +if __name__ == "__main__": + unittest.main() diff --git a/tests/test_tempo_models.py b/tests/test_tempo_models.py index 578213e..9d308d3 100644 --- a/tests/test_tempo_models.py +++ b/tests/test_tempo_models.py @@ -1,22 +1,17 @@ #!/usr/bin/env python # -*- coding: utf-8 -*- """ -Tests for Arzt (2010) Simple Tempo Model and Raphael/Jiang (2020) SKF. +Tests for the Arzt (2010) simple tempo model. """ import unittest +from unittest.mock import patch import numpy as np from matchmaker.dp.oltw_arzt import OnlineTimeWarpingArztFrame -from matchmaker.prob.skf import ( - SwitchingKalmanFilterFollower, - build_chord_sequence, - build_spectral_templates, -) from matchmaker.registry import REGISTRY from matchmaker.utils.tempo_models import ArztTempoModel -from tests.utils import generate_example_sequences RNG = np.random.RandomState(1984) @@ -65,6 +60,41 @@ def test_initialization(self): self.assertTrue(tracker.is_onset[50]) self.assertFalse(tracker.is_onset[25]) + def test_invalid_tempo_configuration(self): + for value in (0.0, -1.0, np.nan, np.inf, -np.inf): + for parameter in ("init_beat_period", "tempo_min", "tempo_max"): + with self.subTest(parameter=parameter, value=value): + with self.assertRaises(ValueError): + ArztTempoModel(**{parameter: value}) + with self.assertRaises(ValueError): + ArztTempoModel(tempo_min=1.2, tempo_max=0.8) + + def test_invalid_estimates_preserve_last_valid_tempo(self): + tracker = self._tracker() + model = tracker.tempo_model + with patch.object(model, "relative_tempo", return_value=1.1): + model.before_step(tracker) + for value in (None, 0.0, -1.0, np.nan, np.inf, -np.inf): + with self.subTest(value=value): + with patch.object(model, "relative_tempo", return_value=value): + model.before_step(tracker) + model.alter(tracker, value) + self.assertEqual(model.current_tempo, 1.1) + self.assertEqual(model.beat_period, 0.5 / 1.1) + self.assertEqual(model.n_inserted + model.n_deleted, 0) + np.testing.assert_array_equal(tracker.V, self.X) + + def test_nonfinite_onsets_do_not_poison_tempo_history(self): + model = ArztTempoModel() + model(0.0, 0.0) + for value in (np.nan, np.inf, -np.inf): + model(value, 1.0) + model(1.0, value) + model(0.5, 1.0) + self.assertEqual(model._history_perf, [0.0, 0.5]) + self.assertEqual(model._history_score, [0.0, 1.0]) + self.assertEqual(model.beat_period, 0.5) + def test_run_steps(self): tracker = self._tracker() beats = self._perform(tracker, 1.0, 200) @@ -85,6 +115,15 @@ def test_tempo_model_follows_faster_performance(self): expected = 1.5 * np.arange(400) / self.frame_rate self.assertLess(np.abs(beats[200:] - expected[200:]).max(), 0.5) + def test_backward_forward_refinement_tracks_with_shared_directions(self): + for speed in (0.8, 1.0, 1.5): + with self.subTest(speed=speed): + tracker = self._tracker(backward_forward=True, max_run_count=3) + beats = self._perform(tracker, speed, 400) + expected = speed * np.arange(400) / self.frame_rate + self.assertTrue(np.all(np.isfinite(beats))) + self.assertLess(np.abs(beats[200:] - expected[200:]).max(), 0.5) + def test_without_tempo_model(self): tracker = self._tracker(tempo_model=None) self._perform(tracker, 1.5, 200) @@ -149,98 +188,5 @@ def test_arzt_tempo_model_standalone_onsets(self): self.assertEqual(tm.current_tempo, 1.0) -class TestRaphaelSwitchingStateSpace(unittest.TestCase): - """Tests for Raphael & Jiang (2020 ISMIR) Switching State-Space Model.""" - - def setUp(self): - # Create a simple synthetic note_array - self.note_array = np.zeros( - 4, - dtype=[ - ("pitch", "i4"), - ("onset_beat", "f4"), - ("duration_beat", "f4"), - ], - ) - self.note_array["pitch"] = [60, 64, 67, 72] # C major chord arpeggio - self.note_array["onset_beat"] = [0.0, 1.0, 2.0, 3.0] - self.note_array["duration_beat"] = [1.0, 1.0, 1.0, 1.0] - - def test_build_chord_sequence(self): - chords, lengths, onset_beats = build_chord_sequence(self.note_array) - self.assertEqual(len(chords), 4) - self.assertEqual(chords[0], [60]) - self.assertEqual(chords[1], [64]) - np.testing.assert_array_equal(onset_beats, [0.0, 1.0, 2.0, 3.0]) - # 1 beat = 0.25 whole note - np.testing.assert_allclose(lengths, [0.25, 0.25, 0.25, 0.25]) - - def test_build_spectral_templates(self): - chords, _, _ = build_chord_sequence(self.note_array) - templates = build_spectral_templates(chords, n_fft=512, sample_rate=8000) - self.assertEqual(templates.shape[0], 4) - self.assertEqual(templates.shape[1], 257) - # Each template row should sum to 1 - np.testing.assert_allclose(templates.sum(axis=1), np.ones(4), atol=1e-5) - - def test_initialization_and_step(self): - follower = SwitchingKalmanFilterFollower( - reference_features=self.note_array, - tempo=120.0, - sample_rate=8000, - n_fft=512, - hop_length=128, - max_hypotheses=50, - ) - self.assertEqual(follower.K, 4) - self.assertEqual(follower.current_index, 0) - # Initial tempo: 240 / 120 = 2.0 s/whole-note - self.assertAlmostEqual(follower.init_tempo, 2.0) - - # Process 20 frames - delta = 128 / 8000.0 - for i in range(20): - dummy_spectrum = np.random.rand(257).astype(np.float32) - pos = follower(dummy_spectrum, perf_time=i * delta) - self.assertIsInstance(pos, float) - self.assertGreaterEqual(pos, 0.0) - self.assertLessEqual(pos, 4.0) - - # Verify hypotheses count is bounded by max_hypotheses - self.assertLessEqual(len(follower.hypotheses), 50) - path = follower.alignment_path - self.assertEqual(path.shape, (2, 20)) - - def test_kalman_tempo_adaptation(self): - follower = SwitchingKalmanFilterFollower( - reference_features=self.note_array, - tempo=120.0, - sample_rate=8000, - n_fft=512, - hop_length=128, - ) - # Test transition Kalman update directly: - # Expected duration = length * tempo = 0.25 * 2.0 = 0.5s. - # Suppose performer played chord 0 in 0.25s (twice as fast, tempo should decrease s/wn): - fast_age = int(0.25 / follower.delta) - mu_new, var_new = follower._kalman_update( - mu_t_prev=2.0, var_t_prev=0.04, k_prev=0, age_prev=fast_age - ) - self.assertLess(mu_new, 2.0) - - # Suppose performer played chord 0 in 1.0s (twice as slow, tempo should increase s/wn): - slow_age = int(1.0 / follower.delta) - mu_new_slow, var_new_slow = follower._kalman_update( - mu_t_prev=2.0, var_t_prev=0.04, k_prev=0, age_prev=slow_age - ) - self.assertGreater(mu_new_slow, 2.0) - - def test_registry_lookup(self): - spec = REGISTRY.method("audio", "skf") - self.assertEqual( - spec.cls_path, "matchmaker.prob.skf:SwitchingKalmanFilterFollower" - ) - - if __name__ == "__main__": unittest.main() diff --git a/tests/test_utils_audio.py b/tests/test_utils_audio.py index 22c3077..f2e64b1 100644 --- a/tests/test_utils_audio.py +++ b/tests/test_utils_audio.py @@ -4,6 +4,8 @@ This module contains tests for the utils/audio.py module. """ +import os +import runpy import unittest from io import StringIO from unittest.mock import patch @@ -18,19 +20,8 @@ list_audio_devices, ) - -def has_audio_devices(): - """Check if system has audio devices""" - p = pyaudio.PyAudio() - - DEVICE_COUNT = p.get_device_count() - - has_audio = DEVICE_COUNT > 0 - - return has_audio - - -HAS_AUDIO_DEVICES = has_audio_devices() +# Reuse the CI guard and initialization-error handling in the audio utilities. +HAS_AUDIO_DEVICES = bool(get_audio_devices()) SKIP_REASON = (not HAS_AUDIO_DEVICES, "No audio device detected") @@ -64,7 +55,6 @@ def test_check_input_audio_devices(self) -> None: self.assertTrue(isinstance(has_audio_inputs, bool)) - @unittest.skipIf(*SKIP_REASON) @patch( "matchmaker.utils.audio.get_audio_devices", return_value=["Device 1", "Device 2"], @@ -96,3 +86,30 @@ def test_get_device_index_from_name(self, mock_stdout) -> None: # Test raising error with self.assertRaises(ValueError): get_device_index_from_name("NoT a ReAl AuDiO dEvIcE") + + +class TestAudioDeviceAvailability(unittest.TestCase): + def test_module_import_in_ci_does_not_initialize_audio(self): + with ( + patch.dict(os.environ, {"CI": "true"}), + patch( + "pyaudio.PyAudio", + side_effect=OSError(-9999, "Unanticipated host error"), + ) as initialize_audio, + ): + module = runpy.run_path(__file__) + self.assertFalse(module["HAS_AUDIO_DEVICES"]) + initialize_audio.assert_not_called() + + def test_module_import_handles_unavailable_audio_host(self): + with ( + patch.dict(os.environ, {"CI": "", "GITHUB_ACTIONS": ""}), + patch( + "pyaudio.PyAudio", + side_effect=OSError(-9999, "Unanticipated host error"), + ) as initialize_audio, + patch("builtins.print"), + ): + module = runpy.run_path(__file__) + self.assertFalse(module["HAS_AUDIO_DEVICES"]) + initialize_audio.assert_called_once()