Weźmy do analizy dokładnie ten sam przykład okienkowany w sposób symetryczny. Poprzednio przedział wynikający z okienkowania był niesymetryczny, gdyż dotyczył (0,T). Nowy przedział, czyli funkcja okienkująca niech będzie wyglądała jak na ilustracji 1, czyli w obszarze zainteresowania znajduje się przedział (-T/2 .. T/2).
Nawet zadbałem, by próbkowana funkcja była dokładnie ta sama jak poprzednio. Została przesunięte w dziedzinie czasu w taki sposób, by pokryła się dla dodatnich argumentów czasu w poprzednim przebiegu. Pokazuje to ilustracja 2.
Spowodowało to pewną ingerencję w sposób generacji ciągu próbek. Fragment kodu jest następujący (argument czasu został przesunięty o połowę długości czasu okienkowania):
Kod: Zaznacz cały
procedure TDFTForm.InputSampleTable ( SampleInterv : double ) ;
var
Loop : word ;
Time : double ;
begin (* TDFTForm.InputSampleTable *)
for Loop := 0 to SampleNum - 1 do
begin (* 1 *)
Time := SampleInterv * ( Loop - SampleNumDiv2 ) ;
TimeTable [ Loop ] := Time ;
InpSampleTable [ Loop ] := 10.0 + 5.0 * sin ( DoublePi * 15.0 * Time ) +
7.0 * cos ( DoublePi * 25.0 * Time ) ;
end (* 1 *) ;
TotalTime := SampleNum * SampleInterv ; (* Okres calego sprobkowanego sygnalu *)
end (* TDFTForm.InputSampleTable *) ;
Oczywiście to ma spore implikacje w algorytmie. Transformata Fouriera obecnie jest obrabiana w innym przedziale czasowym a to ma swoje odbicie na egzystencję w innej przestrzeni. Zmienia się początek nawijania.
Kod: Zaznacz cały
procedure DFT ( var DFTTable : array of ComplexT ;
var SampleTable : array of double ;
SampleTabSize : word ) ;
var
HarmInx : word ;
SampleInx : word ;
DAngle : double ;
Angle : double ;
begin (* DFT *)
SampleTabSize := SampleTabSize - 1 ;
for HarmInx := 0 to SampleTabSize do
begin (* 1 *)
DFTTable [ HarmInx ] := CAssign ( 0.0 , 0.0 ) ;
end (* 1 *) ;
for HarmInx := 0 to SampleTabSize do
begin (* 1 *)
DAngle := DoublePi * HarmInx / SampleNum ;
for SampleInx := 0 to SampleTabSize do
begin (* 2 *)
Angle := DAngle * SampleInx - Pi ;
DFTTable [ HarmInx ] . Re := DFTTable [ HarmInx ] . Re +
SampleTable [ SampleInx ] * cos ( Angle ) ;
DFTTable [ HarmInx ] . Im := DFTTable [ HarmInx ] . Im -
SampleTable [ SampleInx ] * sin ( Angle ) ;
end (* 2 *) ;
end (* 1 *) ;
end (* DFT *) ;
Nawijanie „tam” musie mieć ekwiwalentne rozwiązanie w nawijaniu „z powrotem”, inaczej nie trafimy do domu.
Kod: Zaznacz cały
procedure IDFT ( var SampleTable : array of double ;
var DFTTable : array of ComplexT ;
SampleTabSize : word ) ;
var
HarmInx : word ;
SampleInx : word ;
DAngle : double ;
Angle : double ;
begin (* IDFT *)
SampleTabSize := SampleTabSize - 1 ;
for SampleInx := 0 to SampleTabSize do
begin (* 1 *)
SampleTable [ SampleInx ] := 0.0 ;
end (* 1 *) ;
for SampleInx := 0 to SampleTabSize do
begin (* 1 *)
DAngle := - DoublePi * SampleInx / SampleNum ;
for HarmInx := 0 to SampleTabSize do
begin (* 2 *)
Angle := DAngle * HarmInx + Pi ;
SampleTable [ SampleInx ] := SampleTable [ SampleInx ] +
DFTTable [ HarmInx ] . Re * cos ( Angle ) +
DFTTable [ HarmInx ] . Im * sin ( Angle ) ;
end (* 2 *) ;
SampleTable [ SampleInx ] := SampleTable [ SampleInx ] / SampleTabSize ;
end (* 1 *) ;
end (* IDFT *) ;
Rezultat działania programu powinien spowodować pewną refleksję. Dane wejściowe dotyczą tej samej funkcji natomiast wyniki w znacznym stopniu się różnią. Przede wszystkim ma to swoje odbicie w zawartości harmonicznych. Część rzeczywistą pokazuje poniższa tabelka i ilustracja 3.
Kod: Zaznacz cały
Okres calego sygnalu: 1024,000 ms
*****************************************************************************
Okres probkowania: 2,000000 ms
Czestotliwosc probkowania: 500,000 Hz
Czestotliwosc podstawowa (pierwsza harmoniczna): 0,977 Hz
Pulsacja podstawowa (pierwsza harmoniczna): 6,135923 rad/s
Harm. 0 .re=-5082,4909935847 .im=-0,0000020849 Czest. =0,000 Hz
Harm. 1 .re=37,5737730002 .im=-3,1477393408 Czest. =0,977 Hz
( . . . )
Harm. 24 .re=345,1209218922 .im=51,8125053015 Czest. =23,438 Hz
Harm. 25 .re=910,0019285673 .im=47,1448851844 Czest. =24,414 Hz
Harm. 26 .re=-1350,5951458007 .im=43,3204757202 Czest. =25,391 Hz
Harm. 27 .re=-382,0569497582 .im=40,1240654503 Czest. =26,367 Hz
( . . . )
Harm. 483 .re=-154,4966348165 .im=-35,0690825601 Czest. =471,680 Hz
Harm. 484 .re=-220,7931480689 .im=-37,4083567146 Czest. =472,656 Hz
Harm. 485 .re=-382,0569525196 .im=-40,1240822263 Czest. =473,633 Hz
Harm. 486 .re=-1350,5951174269 .im=-43,3205377508 Czest. =474,609 Hz
Harm. 487 .re=910,0019632004 .im=-47,1448413477 Czest. =475,586 Hz
Harm. 488 .re=345,1209292216 .im=-51,8124877931 Czest. =476,563 Hz
Harm. 489 .re=214,9390322480 .im=-57,6509126488 Czest. =477,539 Hz
Widmo to wykazuje jedynie harmoniczną o częstotliwości 25Hz (podobnie jak poprzednio energia rozlała się i zasiliła sąsiednie prążki widma → próbkowanie nie trafiło precyzyjnie w częstotliwość sygnału). Patrząc na formułę sygnału, z której był generowany ciąg próbek:10+5 sin ( 2 π f1 t ) + 7 cos ( 2 π f2 t), gdzie f1=15Hz i F2=25Hz, nie powinno dziwić, że widmo funkcji parzystej (cos) przebiło się przez te wszystkie „zawirowania na wałku” natomiast widmo funkcji nieparzystej (sin) nie było w stanie dotrzeć.
Oczywiście wystąpiła harmoniczna symetryczna w stosunku do częstotliwości Nyquista.
Wybrane elementy widma części urojonej pokazuje poniższa tabelka i całe widmo pokazuje ilustracja 4.
Kod: Zaznacz cały
Harm. 13 .re=52,2355272581 .im=-143,3330134370 Czest. =12,695 Hz
Harm. 14 .re=55,5921916630 .im=-258,6466245991 Czest. =13,672 Hz
Harm. 15 .re=59,6666444390 .im=-1012,0509247227 Czest. =14,648 Hz
Harm. 16 .re=64,6764661718 .im=587,6381113817 Czest. =15,625 Hz
Harm. 17 .re=70,9411580662 .im=236,0272549218 Czest. =16,602 Hz
( . . . )
Harm. 494 .re=78,9490741259 .im=-150,5277796696 Czest. =482,422 Hz
Harm. 495 .re=70,9411655781 .im=-236,0272444121 Czest. =483,398 Hz
Harm. 496 .re=64,6764897364 .im=-587,6380824086 Czest. =484,375 Hz
Harm. 497 .re=59,6665927333 .im=1012,0509472941 Czest. =485,352 Hz
Harm. 498 .re=55,5921747719 .im=258,6466245501 Czest. =486,328 Hz
Harm. 499 .re=52,2355152801 .im=143,3330137584 Czest. =487,305 Hz
Wskazuje ono jedynie harmoniczną o częstotliwości 15Hz (podobnie jak poprzednio energia rozlała się i zasiliła sąsiednie prążki widma). Analogicznie jak w części rzeczywistej, tu na wierzch wydobyły się jedynie echa funkcji nieparzystych (o częstotliwości 15Hz). Oczywiście wystąpiła harmoniczna po drugiej stronie częstotliwości Nyquista.
Odwrotna transformata zrobiła to co do niej należało, odtworzony sygnał pokazuje ilustracja 5.
Błąd odtworzenia funkcji (wynik transformaty odwrotnej) jest na podobnym poziomie (ilustracja 6).
Patrząc na historię zaistniałych zdarzeń, można zadać sobie pytanie: co się stało, że wyszło jak wyszło?
Generalnie zagadnienie wykazuje dużą wrażliwość na zaistniałe warunki. Wybranie jakiejś funkcji okienkującej wpływa na wynik końcowy: wystąpienie określonej przyczyny generuje ściśle określone skutki. Nastąpił przeciek energii nawet pomiędzy częścią rzeczywistą a częścią urojoną widma. Tu również występuje swoista symetria: określone skutki są następstwem ściśle określonych przyczyn.
Program dla tropicieli: