Interpolation of cross-correlation from zero padded idft of power spectrum shows odd behavior
13 visualizzazioni (ultimi 30 giorni)
Mostra commenti meno recenti
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))
0 Commenti
Risposta accettata
Paul
il 1 Dic 2022
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
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
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))
0 Commenti
Più risposte (0)
Vedere anche
Categorie
Scopri di più su Fourier Analysis and Filtering in Help Center e File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!