Кросс-корреляция Python

У меня есть пара одномерных массивов (разной длины), как показано ниже:

data1 = [0,0,0,1,1,1,0,1,0,0,1]
data2 = [0,1,1,0,1,0,0,1]

Я хотел бы получить максимальную взаимную корреляцию 2 серии в Python. В Matlab, xcorr() функция вернет это ОК

Я попробовал следующие 2 метода:

  1. numpy.correlate(data1, data2)
  2. signal.fftconvolve(data2, data1[::-1], mode='full')

Оба метода дают мне одинаковые значения, но значения, которые я получаю от python, отличаются от значений, которые получаются из matlab. Python дает мне целые значения> 1, тогда как Matlab дает фактические значения корреляции между 0 и 1.

Сначала я попытался нормализовать 2 массива (среднее значение /SD), но полученные значения взаимной корреляции исчисляются тысячами, что кажется неправильным.

Matlab также даст вам значение запаздывания, при котором взаимная корреляция является наибольшей. Я предполагаю, что это легко сделать с помощью индексов, но каков наиболее подходящий способ сделать это, если мои массивы содержат десятки тысяч значений?

Я хотел бы подражать xcorr() функция, которая есть у Matlab, есть мысли о том, как бы я это сделал в Python?

3 ответа

Решение
numpy.correlate(arr1,arr2,"full")

дал мне тот же вывод, что и

xcorr(arr1,arr2)

дает в Matlab

Реализация MATLAB xcorr(x,y) и сравнение результата с примером.

      import scipy.signal as signal
def xcorr(x,y):
    """
    Perform Cross-Correlation on x and y
    x    : 1st signal
    y    : 2nd signal

    returns
    lags : lags of correlation
    corr : coefficients of correlation
    """
    corr = signal.correlate(x, y, mode="full")
    lags = signal.correlation_lags(len(x), len(y), mode="full")
    return lags, corr

n = np.array([i for i in range(0,15)])
x = 0.84**n
y = np.roll(x,5);
lags,c = xcorr(x,y);
plt.figure()
plt.stem(lags,c)
plt.show()

Этот код поможет найти задержку между двумя каналами в аудиофайле.

xin, fs = sf.read('recording1.wav')
frame_len = int(fs*5*1e-3)
dim_x =xin.shape
M = dim_x[0] # No. of rows
N= dim_x[1] # No. of col
sample_lim = frame_len*100
tau = [0]
M_lim = 20000 # for testing as processing takes time
for i in range(1,N):
    c = np.correlate(xin[0:M_lim,0],xin[0:M_lim,i],"full")
    maxlags = M_lim-1
    c = c[M_lim -1 -maxlags: M_lim + maxlags]
    Rmax_pos = np.argmax(c)
    pos = Rmax_pos-M_lim+1
    tau.append(pos)
print(tau)
Другие вопросы по тегам