pythonnumpysignal-processingcorrelation

find time shift between two similar waveforms


I have to compare two time-vs-voltage waveforms. Because of the peculiarity of the sources of these waveforms, one of them can be a time shifted version of the other.

How can i find whether there is a time shift? and if yes, how much is it.

I am doing this in Python and wish to use numpy/scipy libraries.


Solution

  • scipy provides a correlation function which will work fine for small input and also if you want non-circular correlation meaning that the signal will not wrap around. note that in mode='full' , the size of the array returned by signal.correlation is sum of the signal sizes minus one (i.e. len(a) + len(b) - 1), so the value from argmax is off by (signal size -1 = 20) from what you seem to expect.

    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
    

    The two different values correspond to whether the shift is in a or b.

    If you want circular correlation and for big signal size, you can use the convolution/Fourier transform theorem with the caveat that correlation is very similar to but not identical to convolution.

    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
    

    again the two values correspond to whether your interpreting a shift in a or a shift in b.

    The negative conjugation is due to convolution flipping one of the functions, but in correlation there is no flipping. You can undo the flipping by either reversing one of the signals and then taking the FFT, or taking the FFT of the signal and then taking the negative conjugate. i.e. the following is true: Ar = -A.conjugate() = fft(a[::-1])