7.1.8. FFT valós vektor esetén
A (3.34) szerinti általános rész-transzformáltra áll, hogy
(7.40)
Ebből az fi valós esetben egyenesen következik, hogy
(7.41)
Ahol * a konjugált képzést jelenti. Más szóval:
(7.42)
számsorozat redundáns, hiszen az r=N/2p fölötti rész konjugálással képezhető a kisebb indexűekből.
Nézzük, mennyi a transzformált nem redundáns része. Ez a következőkből áll:
valós érték (7.34) szerint egyértelmű,
valós érték,
komplex értékek.
Ha összeszámláljuk, kiderül, hogy ez éppen N darab adatot jelent (a komplex számot itt két adatnak vettük). Ez megfelel a várakozásnak, hiszen a transzformáció N darab adatból indul.
A valós bemenő vektorra egyszerűsített algoritmusnál tehát nem tároljuk a teljes rész-transzformáltat, csak a nem redundáns részt. Ennek megfelelően át kell rendeznünk a (7.39) szerinti rekurziós formulát. Ez a rendezés az alábbi eredményt adja, amelyet nem részletezek:
(7.43)
(7.44)
(7.45)
Egy N/2 elemű komplex vektorban a 0 indexű elem nem tényleges komplex szám, hanem a valós illetőleg a képzetes rész helyén a transzformált két tisztán valós elemét, D0 és DN/2 –t tároljuk. A transzformáció eredményét ezen rend szerint kapjuk. Ez a tárolási séma a rész-transzformáltakra is vonatkozik.
A transzformációt végző C függvény [6] alapján a következő:
//A megfelelő paraméterátadási konvenció beállítása: WINAPI,
//és a függvény láthatóvá tétele a külvilág számára: __declspec(dllexport)
extern “C” __declspec(dllexport)
void WINAPI refft(double *poi) //transzformálandó vektor, mint paraméter
{ unsigned long int p,p2,q,r,i,nper2p,nper4p;
double s,*poi1,*poi2;
complex *p_exp;
complex *p_from,*p_to,*p1_from,*p1_to,*p2_from,*p2_to;
complex cx,cxexp;
//Első rekurziós lépés megvalósítása
p=n/2;
poi1=poi; //Munkaterületek kezdőcímeinek kijelölése
poi2=poi+p;
p1_to=fftbuff;
for(q=0;q<p;q++)
{ p1_to->re=*poi1+(*poi2);
p1_to->im=*poi1-(*poi2);
poi1++;
poi2++;
p1_to++;
}
//A további rekurziós lépések ciklusa
p_from=fftbuff;
p_to=(complex *)poi;
for(i=1;i<logn;i++)
{ p2=p;
p>>=1;
nper2p=n/p2;
//Az r=0 és r=n/2p eset
p2_from=p+(p1_from=p_from);
p1_to=p_to;
for(q=0;q<p;q++) { p1_to->re=p1_from->re+p2_from->re;
p1_to->im=p1_from->re-p2_from->re;
p1_from++;
p2_from++;
p1_to++;
}
p1_from+=p;
p2_from+=p;
p2_to=p_to+n/2-p;
//Általános eset
if(n>=4*p)
{ p_exp=exp_tabl+p-1;
nper4p=nper2p/2;
for(r=1;r<nper4p;r++)
{ cx_copy(p_exp,&cxexp);
for(q=0;q<p;q++)
{ cx_mul(p2_from,&cxexp,&cx);
cx_add(p1_from,&cx,p1_to);
cx_sub(p1_from,&cx,p2_to);
p2_to->im=-p2_to->im;
p1_from++;
p2_from++;
p1_to++;
p2_to++;
}
p_exp+=p;
p1_from+=p;
p2_from+=p;
p2_to-=p2;
}
}
//Az r=n/4p eset
p2_from=p+(p1_from=p_from);
for(q=0;q<p;q++)
{ p1_to->re=p1_from->im;
p1_to->im=-p2_from->im;
p1_from++;
p2_from++;
p1_to++;
}
p1_from=p_from;
p_from=p_to;
p_to=p1_from;
}
if(logn%2) //Visszamásolás a poi által mutatott adaterületre
{ poi1=(double *)fftbuff;
poi2=poi;
for(i=0;i<n;i++)
*(poi2++)=*(poi1++);
}
s=1.0/(double)n; //Az eredmény végigosztása n-nel
for(i=0;i<n;i++) *(poi++)*=s;
return;
}
Az inverz transzformáció hasonló meggondolások szerint kódolható, listája megtalálható a [6] alapján készült F.27-ben.
A fenti függvényt Delphi-ben a következőképpen kell deklarálni:
procedure ReFft(poi:pdouble); //pdouble=^double
stdcall; external ‘prjconvoldll.dll’ name ‘refft’;
ahol prjconvoldll.dll az állomány, illetve a dinamikusan szerkesztett könyvtár, ahonnét a kérdéses függvényt behívjuk, refft pedig a C függvény eredeti neve.
[6] Székely Vladimír: Képkorrekció, hanganalízis, térszámítás PC-n (Gyors Fourier transzformációs módszerek), Computer Books, Budapest, 1994
