Interpolation of cross-correlation from zero padded idft of power spectrum shows odd behavior

Hello,

I'm trying to interpolate the cross correlation between two singals by zero padding the ifft of the signals' PSD. Can anyone explain the odd spikiness of the interpolated correlation?

When the dft is evalauted directly via goertzel's method at non-integer indices, the same spikiness shows up.

Thanks for taking a look at this!

Signals being correlated

d1 = designfilt("lowpassiir",FilterOrder=12, ...
HalfPowerFrequency=0.02,DesignMethod="butter");


N = 2^10;
y1 = rand(N,1);

y1 = filtfilt(d1,y1)';

N2 = 2^2;

shiftImposed = 1*N2;

y2 = interp1((1:N),y1,(1+shiftImposed:N+shiftImposed));


y3  = imresize(y1,1/N2);
y4  = imresize(y2,1/N2);


y5 = y3(length(y3)/8:length(y3)-length(y3)/8+1);
y6 = y4(length(y3)/8:length(y3)-length(y3)/8+1);

figure(54);plot(y5);hold on;plot(y6)

Correlation without interpolation

g1 = y5-mean(y5);
g2 = y6-mean(y6);

G1 = fft(g1);
G2 = fft(g2);

R = G1.*conj(G2);

N = length(R);
p=0;
N2 = 2^p*N;
x = (1:N2);

R2 = padarray(R,[0 (N2-N)/2] ,0,'both');

rr2 = fftshift(ifft(R2,'symmetric'));

figure(43);plot(abs(rr2))

Correlation with 2N interpolation via zero padding ifft of cross power spectrum

g1 = y5-mean(y5);
g2 = y6-mean(y6);

G1 = fft(g1);
G2 = fft(g2);

R = G1.*conj(G2);

N = length(R);
p=1;
N2 = 2^p*N;
x = (1:N2);

R2 = padarray(R,[0 (N2-N)/2] ,0,'both');

rr2 = fftshift(ifft(R2,'symmetric'));

figure(43);plot(abs(rr2))

Goertzel's Method to calculate dft bins of arbitrary location on unit circle shows the same thing as figure 3

x = (1:.5:N);

DD = fftshift(goertzel(R,x));

figure(457);plot(abs(DD))

Answer 1

[Paul]

Hi matt,

I think the issue is with how the padding is being done. See code below for that and other questions. Full disclosure: I'm only looking at how the fft/ifft is being used; I'm not familiar with this type of problem.

rng(100);

d1 = designfilt("lowpassiir",FilterOrder=12, ...
HalfPowerFrequency=0.02,DesignMethod="butter");
N = 2^10;
y1 = rand(N,1);
y1 = filtfilt(d1,y1)';
N2 = 2^2;
shiftImposed = 1*N2;
y2 = interp1((1:N),y1,(1+shiftImposed:N+shiftImposed));
y3  = imresize(y1,1/N2);
y4  = imresize(y2,1/N2);
y5 = y3(length(y3)/8:length(y3)-length(y3)/8+1);
y6 = y4(length(y3)/8:length(y3)-length(y3)/8+1);
figure(54);plot(y5);hold on;plot(y6)

Correlation without interpolation

g1 = y5-mean(y5);
g2 = y6-mean(y6);
G1 = fft(g1);
G2 = fft(g2);
R = G1.*conj(G2);
N = length(R);
p=0;
N2 = 2^p*N;
x = (1:N2);
% no padding here
N2 - N

ans = 0

R2 = padarray(R,[0 (N2-N)/2] ,0,'both');

Why is fftshift applied here?

rr2 = fftshift(ifft(R2,'symmetric'));
figure(43);plot(abs(rr2))

Correlation with 2N interpolation via zero padding ifft of cross power spectrum

g1 = y5-mean(y5);
g2 = y6-mean(y6);
G1 = fft(g1);
G2 = fft(g2);
R = G1.*conj(G2);
N = length(R);
p=1;
N2 = 2^p*N;
x = (1:N2);
% pad with this many 0's in both directions
N2-N

ans = 194

I think the issue is here. Need to fftshift R, and then pad the left and the right for negative and positive frequencies, and then ifftshift back for use with ifft in the next step.

%R2 = padarray(R,[0 (N2-N)/2] ,0,'both');
R2 = ifftshift(padarray(fftshift(R),[0 (N2-N)/2] ,0,'both'));

Why is fftshift applied here?

rr2 = fftshift(ifft(R2,'symmetric'));
figure(43);plot(abs(rr2))