найти сдвиг во времени между двумя одинаковыми формами волны - PullRequest
19 голосов
/ 14 января 2011

Я должен сравнить два сигнала времени-напряжения.Из-за особенностей источников этих сигналов одна из них может быть сдвинутой во времени версией другой.

Как я могу узнать, есть ли сдвиг во времени?и если да, то сколько это стоит.

Я делаю это на Python и хочу использовать библиотеки numpy / scipy.

Ответы [ 5 ]

35 голосов
/ 14 января 2011

scipy обеспечивает функцию корреляции, которая будет отлично работать для небольших входных данных, а также, если вы хотите некруговую корреляцию, означающую, что сигнал не будет распространяться. обратите внимание, что в mode='full' размер массива, возвращаемого signal.correlation, является суммой размеров сигнала минус один (т. е. len(a) + len(b) - 1), поэтому значение из argmax отключено (размер сигнала -1 = 20) из того, что вы, похоже, ожидаете .

from scipy import signal, fftpack
import numpy
a = numpy.array([0, 1, 2, 3, 4, 3, 2, 1, 0, 1, 2, 3, 4, 3, 2, 1, 0, 0, 0, 0, 0])
b = numpy.array([0, 0, 0, 0, 0, 1, 2, 3, 4, 3, 2, 1, 0, 1, 2, 3, 4, 3, 2, 1, 0])
numpy.argmax(signal.correlate(a,b)) -> 16
numpy.argmax(signal.correlate(b,a)) -> 24

Два разных значения соответствуют тому, находится ли сдвиг в a или b.

Если вам нужна круговая корреляция и для большого размера сигнала, вы можете использовать теорему свертки / преобразования Фурье с оговоркой, что корреляция очень похожа, но не идентична свертке.

A = fftpack.fft(a)
B = fftpack.fft(b)
Ar = -A.conjugate()
Br = -B.conjugate()
numpy.argmax(numpy.abs(fftpack.ifft(Ar*B))) -> 4
numpy.argmax(numpy.abs(fftpack.ifft(A*Br))) -> 17

снова эти два значения соответствуют тому, интерпретируете ли вы сдвиг в a или сдвиг в b.

Отрицательное сопряжение связано с переворачиванием одной из функций сверткой, но в корреляции переворот отсутствует. Вы можете отменить переключение, изменив один из сигналов и затем взяв БПФ, или взяв БПФ сигнала, а затем взяв отрицательное сопряжение. то есть верно следующее: Ar = -A.conjugate() = fft(a[::-1])

10 голосов
/ 14 января 2011

Если один сдвинут во времени другим, вы увидите пик в корреляции.Поскольку вычисление корреляции стоит дорого, лучше использовать FFT.Итак, что-то вроде этого должно работать:

af = scipy.fft(a)
bf = scipy.fft(b)
c = scipy.ifft(af * scipy.conj(bf))

time_shift = argmax(abs(c))
6 голосов
/ 15 января 2011

Эта функция, вероятно, более эффективна для реальных значений сигналов.Он использует rfft и нулевую площадку для входов до степени 2, достаточно большой, чтобы обеспечить линейную (т.е. некруговую) корреляцию:

def rfft_xcorr(x, y):
    M = len(x) + len(y) - 1
    N = 2 ** int(np.ceil(np.log2(M)))
    X = np.fft.rfft(x, N)
    Y = np.fft.rfft(y, N)
    cxy = np.fft.irfft(X * np.conj(Y))
    cxy = np.hstack((cxy[:len(x)], cxy[N-len(y)+1:]))
    return cxy

Возвращаемое значение равно длине M = len(x) + len(y) - 1 (взломано вместе с hstackубрать лишние нули от округления до степени 2).В неотрицательные лаги cxy[0], cxy[1], ..., cxy[len(x)-1], в то время как отрицательные лаги cxy[-1], cxy[-2], ..., cxy[-len(y)+1].

1009 * Для согласования опорного сигнала, я бы вычислить rfft_xcorr(x, ref) и искать пика.Например:
def match(x, ref):
    cxy = rfft_xcorr(x, ref)
    index = np.argmax(cxy)
    if index < len(x):
        return index
    else: # negative lag
        return index - len(cxy)   

In [1]: ref = np.array([1,2,3,4,5])
In [2]: x = np.hstack(([2,-3,9], 1.5 * ref, [0,3,8]))
In [3]: match(x, ref)
Out[3]: 3
In [4]: x = np.hstack((1.5 * ref, [0,3,8], [2,-3,-9]))
In [5]: match(x, ref)
Out[5]: 0
In [6]: x = np.hstack((1.5 * ref[1:], [0,3,8], [2,-3,-9,1]))
In [7]: match(x, ref)
Out[7]: -1

Это не надежный способ сопоставления сигналов, но он быстрый и простой.

2 голосов
/ 26 апреля 2017

Вот еще один вариант:

from scipy import signal, fftpack

def get_max_correlation(original, match):
    z = signal.fftconvolve(original, match[::-1])
    lags = np.arange(z.size) - (match.size - 1)
    return ( lags[np.argmax(np.abs(z))] )
2 голосов
/ 14 января 2011

Это зависит от типа вашего сигнала (периодического? ...), от того, имеют ли оба сигнала одинаковую амплитуду, и от того, какую точность вы ищете.

Корреляционная функция, упомянутая highBandWidth, может действительноработать для вас.Достаточно просто, чтобы вы попробовали.

Другой, более точный вариант, который я использую для высокоточной подгонки спектральной линии: вы моделируете свой "мастер" сигнал со сплайном и подгоняете время-перемещенный с ним сигнал (при этом возможно масштабирование сигнала, если это необходимо).Это дает очень точные временные сдвиги.Одним из преимуществ этого подхода является то, что вам не нужно изучать корреляционную функцию.Например, вы можете легко создать сплайн с помощью interpolate.UnivariateSpline() (из SciPy).SciPy возвращает функцию, которая затем легко устанавливается с помощью optimize.leastsq ().

...