|
| Spørgsmål: autokorrelation, fft, power spe~ Fra : Joe |
Dato : 27-02-05 15:59 |
|
Hej
Jeg har lavet et matlab-program som udregner autokorrelationen for sekvensen
y på 2 måder. Den ene måde tager udgangspunkt i y's power spektrum mens den
anden måde udregner autokorrelationen direkte vha. xcorr-kommandoen. De 2
resultater afviger dog fra hinanden. Hvordan kan det være? Tak på forhånd
Her er koden:
clc
close all
clear
N=1000;
% Metode 1
y=[randn(1,N)];
Syy=(abs(fft(y)).^2)./N;
plot(real(ifft(Syy)),'--');
% Metode 2
Ryy=xcorr(y,'coeff')
[maxRyy,lagzero] = max(Ryy);
hold
plot(Ryy(lagzero:2*N-1),'r--')
| |
TA (07-03-2005)
| Kommentar Fra : TA |
Dato : 07-03-05 21:03 |
|
"Joe" <JoeNOSPAM@yahoo.com> wrote in message
news:cvsn7f$10kp$1@news.cybercity.dk...
>
> Jeg har lavet et matlab-program som udregner autokorrelationen for
sekvensen
> y på 2 måder. Den ene måde tager udgangspunkt i y's power spektrum mens
den
> anden måde udregner autokorrelationen direkte vha. xcorr-kommandoen. De 2
> resultater afviger dog fra hinanden. Hvordan kan det være? Tak på forhånd
Hey Joe,
Du får formentlig forskellige resultater fordi du i frekvensdomænet udfører
en cirkulær foldning af signalerne. At korrelere to signaler med hinanden
svarer til at folde det ene signal med en tidsreverseret version af det
andet. I frekvensdomænet kan et signal på simpel vis reverseres ved at
kompleks konjugere det. Desuden kan to signaler foldes ved at multiplicere
dem med hinanden. Korrelationsfunktionen kan derfor beregnes som:
FFT(x).*conj(FFT(x)). Dette er også effektspektret FFT(x).^2.
Effekterne af en cirkulær foldning kan undgåes ved at zero-pad'de
tidssignalet med zero'es (duh) således at signalet ikke overlapper sig selv
på enhedscirklen. Hvis tidssignalet zero-pad'des til dobbelt længde vil
dette aldrig ske. Har man til gengæld et meget langt tidssignal, og ønsker
kun at beregne korrelationsfunktionen for ganske få tau-værdier ('lags'), er
det tilstrækkeligt kun at zero-pad'de med et tilsvarende få antal zero'es.
Se endvidere følgende kode.
clc
close all
clear
N = 1024;
% tidssignal
x = randn(1, N);
% Metode 1 - korrelationsfunktionen beregnes vha fft
% for alle lags og tidssignalet zero-pad'des til
% dobbelt længde
%
tmp = [x zeros(1, N)];
Syy1 = (abs( fft( tmp)).^2)./N;
R1 = real( ifft( Syy1));
% Metode 2 - korrelationsfunktionen beregnes vha fft
% for de første 16 lags (to-sidet korrelationsfunktion)
% De øvrige korrelationsværdier indeholder fejl pga
% cirkulær foldning
%
tmp = [x zeros(1, 16)];
TMP = fft( tmp);
Syy2 = TMP.*conj(TMP)./length(tmp);
R2 = real( ifft( Syy2));
% Metode 3 - korrelationsfunktionen beregnes direkte fra
% tidssignalet vha xcorr
%
R3 = xcorr( x, x)./N;
% Plot korrelationsfunktionerne for alle positive lags
%
plot( R1(1:1024), 'y'), hold,
plot( R2(1:1024), 'r')
plot( R3(1024:1024+1023), 'g')
% R1 og R3 er identiske for alle lags,
% hvorimod R2 har højere korrelationer ved lange
% lags pga cirkulær korrelation
% Effekten af cirkulær korrelation er mere
% tydelig hvis man blot ser på de første 50 lags
%
% De to funktioner er ens de første 16 lags, men
% for højere lags er de forskellige
%
hold off
plot( R1(1:50), 'y'), hold
plot( R2(1:50), 'r')
--
TA
| |
|
|