[PD] Cross correlation using FFT

[Jamie Bullock]

Hi,

On Fri, 7 Apr 2006 12:38:53 -0500
"Charles Henry" <czhenry at gmail.com> wrote:

> Hi,
>   I encountered this problem before when I wanted to do a 'windowed'
> cross-correlation function in pd too.  Your algorithm is simpler than
> mine, I wrote an external to do it on its own.
> 
> (x + x'i)(y-y'i) = (xy + x'y') + (x'y - xy')*i
> 

And if you weren't time reversing y, it would be:

(xy + x'y') + (xy' + x'y)i 

.. as Matju suggests?


> Now then, there's the problem of finding what the output means:
> the cross-correlation (in this example) is being done as a convolution
> between signal x, and time-reversed signal y.
> 
> So, the summation notation can shed some light on it:
> 
> 1st: simple convolution: x * y (i) = sum( j=1, N-1: x(j) y(i-j) )
> In the case of ffts , we do circular convolution, so if it's y(N) it
> means y(0), and y(-5)=y(N-5)
> 
> The real thing:
> x * rev(y) (i) = sum ( j=1, N :  x(j) rev(y) (i-j)  )
> = sum ( j=1, N :  x(j) y (j-i)  )
> 

I'm afraid I don't understand your notation here. Are we still dealing with FFTs? ...and is i the subscript for the input vector array y?

Something like:
	
	for(i = 0; i < ??? ; i++){
		for(j = 1; j < N; j++)
			correlation_coefficient[i] = x[j] * y[j - i];
	}

.. or are x and y fft(x) and fft(y), making the above code somewhat more complex ?

> So, the output should reflect that for each positive i , we get the
> signal y, shifted i samples forward in time, in product with x.  For
> negative i (the samples N-i), we get the signal y delayed by i samples
> in product with x.
> 
> In my implementation, I kept the middle half of y the same, and let x
> be completely the same, like this:
> 
> x:|/\/\/\/\/\/\/\/\/\/\|
> y:|-----/\/\/\/\/\-----|
> 
> Then, the cross correlation has a maximum value in either (1, N/4)
> where y lags behind x, or in (3*N/4, N-1) where y leads x.  Then, the
> samples from (N/4, 3*N/4) were just garbage, left over from the fft.
> 

Did you devise that approach yourself, or is there a paper you can point me towards?

> There are some problems with this, to be sure.  It all depends on what
> you want to use it for?
> 

I would like to try a few things like detecting proximity between two microphones, and detection of acoustic spill. Maybe also fundamental estimation by auto correlation. It just seems a useful tool.


> Chuck
> 
> P.S.  I don't have the external available at the moment....but I can
> send it later, if you'd like.
> 

That would be really good. I think I probably understand code better than I understand maths ;-)

We can do a swap - you can find my time-domain  circular cross correlation external as part of my flib library in the PD CVS under postlude. It's just a brute force approach, not especially pretty or efficient. An optional argument gives the maximum delay (in samples) between the two inputs for each correlation coefficient. Conceptually the summations are performed for -maxdelay < 0 < maxdelay, but in reality the values are 'wrapped' so that the index to the delayed signal stays within the array bounds.

best,
Jamie