Suppose there are two functions $f(x)$ and $g(x)$. Each function is nonzero only in the specified interval, $[m_f, M_f]$ and $[m_g, M_g]$. What I want to calculate is the cross-correlation between two functions: $$ (f\star g)(t) = \int_{-\infty}^{\infty} f^\ast(x)g(x+t)dx $$ which can be obtained by applying convolution theorem: $$ (f\star g)(t) = \mathcal{F}^{-1}\left\{{\mathcal{F}\left\{f\right\}}^\ast\cdot\mathcal{F}\left\{g\right\}\right\}=\mathcal{F}^{-1}\left\{\hat{f}^\ast\cdot\hat{g}\right\} $$ where $\hat{.}$ denotes forward transformation.

Now, to write this in FFT, I prepared four sets of arrays: x_f, x_g, y_f, y_g. x_f and x_g denote the $x$ domain of respective functions, and they are 1. zero-pad and 2. sampled in the same sampling. Basically, I took whichever was greater between $M_f-m_f$ and $M_g-m_g$, defined it to be $L/2$, made it as the interval axes. So, x_f[0] is $m_f$, and x_f[-1] is $m_f+L\geq M_f$. The same goes for x_g. y_f and y_g are defined to be $f(x)$ and $g(x)$ at each x_f and x_g.

If we take y_c = ifft(fft(y_f).conjugate()*fft(y_g)), we get a somewhat valid result. However, it is not easy to map y_c corresponds to what x_c so that y_c is equal to $(f\star g)(x)$ at x_c. How do you calculate this?

Also, what happens if there are three terms? i.e., $(f\star (g\star h))(t)$?


1 Answer 1


Let's say we have two discrete functions $f[n]$ and $g[n]$ that have finite support on $[m_f, M_f]\ m_f,M_f \in \mathbb{Z}$ and $[m_g, M_g] \ m_g,M_g \in \mathbb{Z}$. Than the lengths of the functions are $N_f = M_f-m_f+1$ and $N_g = M_g-m_g+1$

The length of the cross correlation $h[n]$ (or convolution) is simply

$$N_h = N_f+N_g-1$$

We treat the cross correlation as "convolve with the time flipped version of g[n]", i.e.

$$h[n] = f[n]*g[-n] $$

where $*$ is the convolution operator.

Since we are time flipping around $n=0$ The support of the time-flipped version is $[-M_g,-m_g]$, i.e. the limits change signs and flip positions. For example $[3,8]$ turns into $[-8,-3]$

The FFT implements the DFT (Discrete Fourier Transform). Multiplication in the DFT domain implements circular (not linear) convolution. Hence we need to zero pad both sequences to the length $N_h$ to avoid time domain aliasing. Note, that you need to zero-pad $g[n]$ AFTER time flipping it.

For convolution the support intervals just add, i.e. $h[n]$ will be non-zero on $[m_h,M_h] = [m_f-M_g,M_f-m_g]$

Programming languages like Matlab don't have a built-in way for managing the x-interval of an array so you have to keep track of it manually. Below is a code example that demonstrate the cross correlation between two different triangular waves.

%% create two tri-angular waves and cross correlate them
% wave 1: width 9, starting at n = 3;
f = conv(ones(5,1),ones(5,1)); f = f/max(f);
mf = [3,11];
% wave 2: width 13, starting at n = -5
g = conv(ones(7,1),ones(7,1)); g = g/max(g);
mg = [-5,7];

%% do the crosscorrelation
Nh = mf(2)+mg(2)-mf(1)-mg(1) + 1;       % length of cross correlation
h2 = ifft(fft(f,Nh).*fft(flip(g),Nh));  % F-domain convolution with padding
mh = mf + flip(-mg);                    % support interval of result

%% now plot it all
hold on;

grid on
set(gca,'xlim',[-5 18]);
xlabel('Time in samples');

enter image description here

Since the peak of $g[n]$ precedes the peak of $f[n]$ by 6 samples, the maximum of the cross correlation will be at $n=6$.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.