在时间戳数组上使用互相关是否有意义？

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

是的，它确实有道理。这通常在matlab中完成。这是一个类似应用程序的链接：

http://www.mathworks.com/help/signal/ug/cross-correlation-of-delayed-signal-in-noise.html

几个考虑

互相关通常用于相关信号噪声过多的情况。但是，如果您没有任何噪音可担心，我会使用不同的方法。

Answer 2

如果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 s，具体取决于引入的随机错误，并且速度差不多。

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

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)

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