Calculer l'autocorrélation en utilisant FFT dans matlab
j'ai lu quelques explications sur la façon dont l'autocorrélation peut être calculée plus efficacement en utilisant le fft d'un signal, en multipliant la partie réelle par le conjugué complexe (domaine de fourier), puis en utilisant le FFT inverse, mais j'ai du mal à m'en rendre compte dans matlab parce qu'à un niveau détaillé, Je ne sais pas vraiment ce que je fais. : o) Toutes les âmes aimables là-bas se soucient de partager un peu de code et de sagesse?
Merci!
2 réponses
comme vous l'avez dit, Prenez le fft et multipliez le pointwise par son conjugué complexe, puis utilisez l'inverse fft (ou dans le cas de la corrélation croisée de deux signaux: Corr(x,y) <=> FFT(x)FFT(y)*)
x = rand(100,1);
len = length(x);
%# autocorrelation
nfft = 2^nextpow2(2*len-1);
r = ifft( fft(x,nfft) .* conj(fft(x,nfft)) );
%# rearrange and keep values corresponding to lags: -(len-1):+(len-1)
r = [r(end-len+2:end) ; r(1:len)];
%# compare with MATLAB's XCORR output
all( (xcorr(x)-r) < 1e-10 )
En fait, si vous regardez le code de xcorr.m, c'est exactement ce qu'il fait (seulement il doit traiter tous les cas de capitonnage, de normalisation, d'entrée vectorielle/matricielle, etc...)
Par théorème de Wiener–Khinchin, la densité spectrale de puissance (DSP) d'une fonction est la Transformée de Fourier de l'autocorrélation. Pour les signaux déterministes, le PSD est simplement la magnitude au carré de la Transformée de Fourier. Voir aussi le théorème de convolution.
quand il s'agit de transformations de Fourier discrètes (c.-à-d. en utilisant FFTs), vous obtenez en fait l'autocorrélation cyclique. Afin d'obtenir l'autocorrélation (linéaire) appropriée, vous devez zéroter le les données originales à deux fois sa longueur originale avant de prendre la Transformée de Fourier. Alors quelque chose comme:
x = [ ... ];
x_pad = [x zeros(size(x))];
X = fft(x_pad);
X_psd = abs(X).^2;
r_xx = ifft(X_psd);