python - 对时间戳数组使用互相关有意义吗？

Question

我想找到两个时间戳数组之间的偏移量。比方说，它们可以表示两个音轨中开始发出哔哔声。

注意:任一音轨中可能有额外或缺失的开始。

我发现了一些看起来很有希望的有关互相关的信息(例如 https://dsp.stackexchange.com/questions/736/how-do-i-implement-cross-correlation-to-prove-two-audio-files-are-similar)。

我假设每个音轨的持续时间为 10 秒，并将蜂鸣声表示为采样率为 44.1 kHz 的“方波”的峰值:

import numpy as np

rfft = np.fft.rfft
irfft = np.fft.irfft

track_1 = np.array([..., 5.2, 5.5, 7.0, ...])
# The onset in track_2 at 8.0 is "extra," it has no
# corresponding onset in track_1
track_2 = np.array([..., 7.2, 7.45, 8.0, 9.0, ...])
frequency = 44100
num_samples = 10 * frequency
wave_1 = np.zeros(num_samples)
wave_1[(track_1 * frequency).astype(int)] = 1
wave_2 = np.zeros(num_samples)
wave_2[(track_2 * frequency).astype(int)] = 1
xcor = irfft(rfft(wave_1) * np.conj(rfft(wave_2)))
offset = xcor.argmax()

这种方法不是特别快，但即使频率非常低，我也能够获得相当一致的结果。但是……我不知道这是不是个好主意!有没有比互相关更好的方法来找到这个偏移量？

编辑:添加了关于缺失和额外发作的注释。

Answer 1

如果 track_1 和 track_2 是蜂鸣声的时间戳，并且都捕捉到所有的蜂鸣声，那么就没有必要构建波形并进行互相关。只需找到两个蜂鸣时间戳数组之间的平均延迟即可:

import numpy as np

frequency = 44100
num_samples = 10 * frequency
actual_time_delay = 0.020
timestamp_noise = 0.001

# Assumes track_1 and track_2 are timestamps of beeps, with a delay and noise added
track_1 = np.array([1,2,10000])
track_2 = np.array([1,2,10000]) + actual_time_delay*frequency + 
          frequency*np.random.uniform(-timestamp_noise,timestamp_noise,len(track_1))

calculated_time_delay = np.mean(track_2 - track_1) / frequency
print('Calculated time delay (s):',calculated_time_delay)

这会产生大约 0.020 秒，具体取决于引入的随机错误，并且会尽可能快。

编辑:如果需要处理额外或缺失的时间戳，则可以按如下方式修改代码。本质上，如果时间戳随机性的误差是有界的，则对所有值之间的延迟执行统计“模式”函数。时间戳随机性范围内的任何事物都聚集在一起并被识别，然后对原始识别的延迟进行平均。

import numpy as np
from scipy.stats import mode

frequency = 44100
num_samples = 10 * frequency
actual_time_delay = 0.020

# Assumes track_1 and track_2 are timestamps of beeps, with a delay, noise added,
#   and extra/missing beeps
timestamp_noise = 0.001
timestamp_noise_threshold = 0.002
track_1 = np.array([1,2,5,10000,100000,200000])
track_2 = np.array([1,2,44,10000,30000]) + actual_time_delay*frequency 
track_2 = track_2 + frequency*np.random.uniform(-timestamp_noise,timestamp_noise,len(track_2))

deltas = []
for t1 in track_1:
    for t2 in track_2:
        deltas.append( t2 - t1)
sorted_deltas = np.sort(deltas)/frequency
truncated_deltas = np.array(sorted_deltas/(timestamp_noise_threshold)+0.5, 
dtype=int)*timestamp_noise_threshold

truncated_time_delay = min(mode(truncated_deltas).mode,key=abs)
calculated_time_delay = np.mean(sorted_deltas[truncated_deltas == truncated_time_delay])

print('Calculated time delay (s):',calculated_time_delay)

显然，如果缺少太多蜂鸣声或存在额外的蜂鸣声，那么代码会在某个时候开始失败，但通常它应该运行良好，并且比尝试生成整个波形并执行关联要快得多。