Can multiple cross-correlation be expressed in terms of multiple Fourier transforms?

A post explains that the cross correlation of two signals can be given by $$\mathcal{F}^{-1} \left\{\mathcal{F}\left\{x(t) \right\} \mathcal{F}\left\{y(t) \right\}^* \right\}$$

Does this generalize to a cross-correlation of multiple signals in the form in Equation $$(1)$$? $$\mathcal{F}^{-1} \left\{\mathcal{F}\left\{x(t) \right\} \prod_{j}^{n-1} \mathcal{F}\left\{y_j(t) \right\}^* \right\}\tag{1}$$

Or perhaps this in Equation $$(2)$$? $$\mathcal{F}^{-1} \left\{\mathcal{F}\left\{x(t) \right\} \left(\prod_{j}^{n-1} \mathcal{F}\left\{y_j(t) \right\}\right)^* \right\}\tag{2}$$

• What's the definition of cross-correlation of multiple signals? Feb 4, 2022 at 9:32

Assuming we interpret the "Cross Correlation of Multiple Signals" to be a cross correlation result for each of a group of signals with a reference signal, this would result in a correlation versus two parameters, $$\rho(\tau,n)$$ where $$\tau$$ represents the typical delay parameter we would use in a cross-correlation function, and $$n$$ represents an index associated with each of the multiple signals. When we would like the cross correlation between a refence signal $$x(t)$$ and multiple other signals given as $$y_n(t)$$, the result for this would be given by:

$$\rho(\tau,n) = \mathcal{F}^{-1} \left\{\mathcal{F}\left\{x(t) \right\} \mathcal{F}\left\{y_n(t) \right\}^* \right\}$$

To note, the FFT approach given in general results in a circular cross-correlation function as opposed to a linear cross-correlation (unless we zero pad to twice the length). The above form can be solved as a single compact matrix equation by placing each result for the internal FFT's in column form (given the default inverse Fourier transform operates across columns in MATLAB, Octave and Python), where the Fourier transform for $$x(t)$$ would just repeat in each column. The product of the two matrices would be element by element rather than vector dot products ($$c_{rc} = a_{rc}b_{rc}$$ for row $$r$$ and column $$c$$).

The resulting matrix can be further processed if the interest was in something else, such as the average of each cross correlation, or whatever parameter the OP is ultimately interested in related to the cross correlation of the multiple waveforms.

To demonstrate, I used this technique for a condensed cross correlation of multiple signals as a joint solution for the Doppler offset and code delay from a single capture of a GPS signal, as demonstrated in this graphic below. This technique can be used for the fast acquisition of any signal where we would use correlation techniques to resolve frequency and time offsets. The result is a "one-shot" coarse frequency and time correction that would place the signal right in the tracking range of typical carrier recovery and timing recovery loops in a receiver.

The complete MATLAB code is available here with shortened version of the main operations below.

Given two time domain signals, $$a$$ and $$b$$, multiple cross correlation functions (correlation vs delay) are computed between the two signals using

$$\text{ifft}\{ \text{fft}\{b_m\} (\text{fft}\{a_m\})^* \}$$

Where $$a_m$$ is a matrix where each column is the FFT result circularly rotated by one bin and thus representing Doppler frequency offsets spaced by $$f_s/N$$ where $$f_s$$ is the sampling rate and $$N$$ is the number of bins (or number of samples in the time domain), and $$b_m$$ is a matrix where each column is the FFT of $$b$$ (repeated).

%============================================================
% create all doppler versions of a by shifting the fft:
%============================================================
a=a(:);    % force column vectors
b=b(:);
N=max(length(a),length(b));
fft_a= fft(a,N);
% frequency axis for all doppler bins -pi to pi
freq_axis=[-floor((N-1)/2):floor(N/2)]*2*pi/N;

% create Hankel matrix for all Dopplers (each column is the fft
% of vector a at a different Doppler offset)
fft_idx=(1:N)';
fft_dopplers = fft_idx(:,ones(N,1))+freq_index(ones(N,1),:);
fft_dopplers= mod(fft_dopplers-1,N)+1;    %Hankel subscripts
fft_dopplers(:)=  fft_a(fft_dopplers);      %fill data
% alternate method if solving for all Doppler bins (processing intensive):
% fft_all_dopplers= hankel(fft_a,[fft_a(end); fft_a(2:end)]);
%==================================================================
% Circular Correlation using circ_corr= ifft(fft(a).*conj(fft(b))
%==================================================================
fft_code_mtx= fft(b,N)*ones(1,N);
fdout=(ifft(fft_code_mtx.*conj(fft_dopplers)));