# Fwd: Potential serious problem with xcorr

5 messages
Open this post in threaded view
|

## Fwd: Potential serious problem with xcorr

 Hi Mike, Just wanted to ask you something regarding xcorr in signal. First of all, I get the same apparent mistake with matlab online, so nothing to worry on that level. However, have a look at the output of this code: N    = 256; t    = linspace (0, 1, N).' * 2 * pi; dt   = t(2) - t(1); freq = 3; x    = sin (freq * t); y    = x.^2; [C lag] = xcorr (x, y); Cw      = arrayfun (@(i)sum (x .* circshift (y, i)), lag.'); plot(lag, C, lag, Cw) indeed the equivalence is obtained with plot(lag, C, lag, cumtrapz(lag, Cw.*sign(lag.'))) So the C obtained with xcorr is sort of the integral of the actual result for periodic signals. I find this gross difference whenever I use periodic signals (both x and y). I think the issue comes from the fact that xcorr is padding with zeros, but here the signals are actually periodic. Don't you think it would be nice to offer a warning about this in the doc (which needs some fixing, will send patch), and maybe offer xcorr_periodic or something like that for the case in which the signals are periodic? Regards, -- JuanPi Carbajal https://goo.gl/ayiJzi----- “An article about computational result is advertising, not scholarship. The actual scholarship is the full software environment, code and data, that produced  the  result.” - Buckheit and Donoho -- JuanPi Carbajal https://goo.gl/ayiJzi----- “An article about computational result is advertising, not scholarship. The actual scholarship is the full software environment, code and data, that produced  the  result.” - Buckheit and Donoho
Open this post in threaded view
|

## Re: Fwd: Potential serious problem with xcorr

 On Tue, Apr 02, 2019 at 00:42:56 +0200, JuanPi wrote: > Just wanted to ask you something regarding xcorr in signal. > First of all, I get the same apparent mistake with matlab online, so > nothing to worry on that level. > However, have a look at the output of this code: > N    = 256; > t    = linspace (0, 1, N).' * 2 * pi; > dt   = t(2) - t(1); > freq = 3; > x    = sin (freq * t); > y    = x.^2; > > [C lag] = xcorr (x, y); > Cw      = arrayfun (@(i)sum (x .* circshift (y, i)), lag.'); > > plot(lag, C, lag, Cw) > > indeed the equivalence is obtained with > > plot(lag, C, lag, cumtrapz(lag, Cw.*sign(lag.'))) > > So the C obtained with xcorr is sort of the integral of the actual > result for periodic signals. > > I find this gross difference whenever I use periodic signals (both x and y). > I think the issue comes from the fact that xcorr is padding with > zeros, but here the signals are actually periodic. Right, the 'xcorr' function does not compute a circular correlation. I wouldn't call this a serious problem or a mistake. Try comparing the above against     Cw = arrayfun (@(i)sum (x .* zeroshift (y, i)), lag.'); instead, where 'zeroshift' is defined as     function y = zeroshift (x, i)       if (i < 0)         y = [x(1-i:end); zeros(-i, 1)];       else         y = [zeros(i, 1); x(1:end-i)];       endif     endfunction > Don't you think it would be nice to offer a warning about this in the > doc (which needs some fixing, will send patch), and maybe offer > xcorr_periodic or something like that for the case in which the > signals are periodic? Sure, the docs can always be improved. It is only very briefly mentioned currently with     where data not provided (for example x(-1), y(N+1)) is zero. Thanks, -- mike
Open this post in threaded view
|

## Re: Fwd: Potential serious problem with xcorr

 In reply to this post by JuanPi On Tue, Apr 02, 2019 at 00:42:56 +0200, JuanPi wrote: > maybe offer > xcorr_periodic or something like that for the case in which the > signals are periodic? I forgot about this last part of your message. The simplest and fastest way to compute the circular cross-correlation today is with 'fft'. Continuing with your example     corr_fft = ifft (fft (x) .* conj (fft (y)));     corr_circ = real ([corr_fft(2:end); corr_fft]); This agrees with the circular cross-correlation you computed using 'arrayfun' and 'circshift'. You can also use the function 'cconv' to compute the same thing     cc = cconv (x, conj (flip (y)), length (x));     cc = circshift (cc, 1);     corr_circ2 = real ([cc(2:end); cc]); I would be in favor of documenting ways to compute a circular cross-correlation using 'fft' and/or 'cconv'. -- mike