Re: Fwd: Potential serious problem with xcorr

[Mike Miller]

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