diff options
author | Erich Eckner <git@eckner.net> | 2017-07-21 10:39:43 +0200 |
---|---|---|
committer | Erich Eckner <git@eckner.net> | 2017-07-21 10:39:43 +0200 |
commit | df7061ed77b9c7e412e4bb068b2eda81025cff70 (patch) | |
tree | 80b72b0ccf980aade4b52783b3dd97968c570d61 | |
parent | 2798b395bc067293d3a7d338a792224da17854e3 (diff) | |
download | units-df7061ed77b9c7e412e4bb068b2eda81025cff70.tar.xz |
fftunit: getrenntes Einlesen von Real- und Imaginaerteilen jetzt auch von Pointer und fuer doAlleResIms und doAlleResSmi (braucht man fuer voll-komplexe fft2d)
-rw-r--r-- | fftunit.inc | 64 | ||||
-rw-r--r-- | fftunit.pas | 82 |
2 files changed, 85 insertions, 61 deletions
diff --git a/fftunit.inc b/fftunit.inc index 7f5cc34..8d465f6 100644 --- a/fftunit.inc +++ b/fftunit.inc @@ -10,7 +10,7 @@ begin {$IFDEF fftNormierung} lenFak:=1/sqrt(len); {$ENDIF} - case ino of + case inO of doRes: for i:=0 to len-1 do begin {$IFDEF tFFTCooleyTukey} @@ -32,7 +32,7 @@ begin {$ENDIF} end; doAlleResSmi: - for i:=0 to hlen + (len and $1) - 1 do begin + for i:=0 to hLen + (len and $1) - 1 do begin {$IFDEF tFFTCooleyTukey} res[perm[i]]:=(q+i*schritt)^ {$IFDEF fftNormierung} * lenFak{$ENDIF}; ims[perm[i]]:=(q+(2*len-1-i)*schritt)^ * inv {$IFDEF fftNormierung} * lenFak{$ENDIF}; @@ -44,33 +44,33 @@ begin doResIms: begin res[0]:=q^ {$IFDEF fftNormierung} * lenFak{$ENDIF}; // perm[0]=0 ims[0]:=0; - for i:=1 to hlen + (len and $1) - 1 do begin + for i:=1 to hLen + (len and $1) - 1 do begin {$IFDEF tFFTCooleyTukey} res[perm[i]]:=(q+i*schritt)^ {$IFDEF fftNormierung} * lenFak{$ENDIF}; - ims[perm[i]]:=(q+(i+hlen)*schritt)^ * inv {$IFDEF fftNormierung} * lenFak{$ENDIF}; + ims[perm[i]]:=(q+(i+hLen)*schritt)^ * inv {$IFDEF fftNormierung} * lenFak{$ENDIF}; res[perm[len-i]]:=res[perm[i]]; ims[perm[len-i]]:=-ims[perm[i]]; {$ELSE} res[i]:=(q+i*schritt)^ {$IFDEF fftNormierung} * lenFak{$ENDIF}; - ims[i]:=(q+(i+hlen)*schritt)^ * inv {$IFDEF fftNormierung} * lenFak{$ENDIF}; + ims[i]:=(q+(i+hLen)*schritt)^ * inv {$IFDEF fftNormierung} * lenFak{$ENDIF}; res[len-i]:=res[i]; ims[len-i]:=-ims[i]; {$ENDIF} end; if not odd(len) then begin {$IFDEF tFFTCooleyTukey} - res[1]:=(q+hlen*schritt)^ {$IFDEF fftNormierung} * lenFak{$ENDIF}; // perm[hlen]=1 + res[1]:=(q+hLen*schritt)^ {$IFDEF fftNormierung} * lenFak{$ENDIF}; // perm[hLen]=1 ims[1]:=0; {$ELSE} - res[hlen]:=(q+hlen*schritt)^ {$IFDEF fftNormierung} * lenFak{$ENDIF}; - ims[hlen]:=0; + res[hLen]:=(q+hLen*schritt)^ {$IFDEF fftNormierung} * lenFak{$ENDIF}; + ims[hLen]:=0; {$ENDIF} end; end; doResSmi: begin res[0]:=q^ {$IFDEF fftNormierung} * lenFak{$ENDIF}; // perm[0]=0 ims[0]:=0; - for i:=1 to hlen + (len and $1) - 1 do begin + for i:=1 to hLen + (len and $1) - 1 do begin {$IFDEF tFFTCooleyTukey} res[perm[i]]:=(q+i*schritt)^ {$IFDEF fftNormierung} * lenFak{$ENDIF}; ims[perm[i]]:=(q+(len-i)*schritt)^ * inv {$IFDEF fftNormierung} * lenFak{$ENDIF}; @@ -85,11 +85,11 @@ begin end; if not odd(len) then begin {$IFDEF tFFTCooleyTukey} - res[1]:=(q+hlen*schritt)^ {$IFDEF fftNormierung} * lenFak{$ENDIF}; // perm[hlen]=1 + res[1]:=(q+hLen*schritt)^ {$IFDEF fftNormierung} * lenFak{$ENDIF}; // perm[hLen]=1 ims[1]:=0; {$ELSE} - res[hlen]:=(q+hlen*schritt)^ {$IFDEF fftNormierung} * lenFak{$ENDIF}; - ims[hlen]:=0; + res[hLen]:=(q+hLen*schritt)^ {$IFDEF fftNormierung} * lenFak{$ENDIF}; + ims[hLen]:=0; {$ENDIF} end; end; @@ -102,7 +102,7 @@ end; {$ENDIF} {$IFDEF fftunitGetrenntLaden} -//procedure tFFTAlgorithmus.laden(invers: boolean; qRe,qIm: tExtendedArray); +//procedure tFFTAlgorithmus.laden(invers: boolean; qRe,qIm: pSingle; schritt: longint); var i,inv: longint; {$IFDEF fftNormierung} @@ -113,15 +113,17 @@ begin {$IFDEF fftNormierung} lenFak:=1/sqrt(len); {$ENDIF} - case ino of - doGetrennt: + case inO of + doGetrennt,doAlleResIms,doAlleResSmi: + // Real- und Imaginärteile sind _immer_ in gleicher Reihenfolge; + // die ggf. verdrehte Dimension ist senkrecht dazu! (benötigt für fft2d) for i:=0 to len-1 do begin {$IFDEF tFFTCooleyTukey} - res[perm[i]]:=qRe[i] {$IFDEF fftNormierung} * lenFak{$ENDIF}; - ims[perm[i]]:=qIm[i] * inv {$IFDEF fftNormierung} * lenFak{$ENDIF}; + res[perm[i]]:=(qRe+i*schritt)^ {$IFDEF fftNormierung} * lenFak{$ENDIF}; + ims[perm[i]]:=(qIm+i*schritt)^ * inv {$IFDEF fftNormierung} * lenFak{$ENDIF}; {$ELSE} - res[i]:=qRe[i] {$IFDEF fftNormierung} * lenFak{$ENDIF}; - ims[i]:=qIm[i] * inv {$IFDEF fftNormierung} * lenFak{$ENDIF}; + res[i]:=(qRe+i*schritt)^ {$IFDEF fftNormierung} * lenFak{$ENDIF}; + ims[i]:=(qIm+i*schritt)^ * inv {$IFDEF fftNormierung} * lenFak{$ENDIF}; {$ENDIF} end; else @@ -136,7 +138,7 @@ var i,inv: longint; begin inv:=1-2*byte(invers); - case outo of + case outO of doRes: for i:=0 to len-1 do (z+i*schritt)^:=res[i]; @@ -152,21 +154,21 @@ begin end; doResIms: begin z^:=res[0]; - for i:=1 to hlen + (len and $1) - 1 do begin + for i:=1 to hLen + (len and $1) - 1 do begin (z+i*schritt)^:=res[i]; - (z+(i+hlen)*schritt)^:=ims[i] * inv; + (z+(i+hLen)*schritt)^:=ims[i] * inv; end; if not odd(len) then - (z+hlen*schritt)^:=res[hlen]; + (z+hLen*schritt)^:=res[hLen]; end; doResSmi: begin z^:=res[0]; - for i:=1 to hlen + (len and $1) - 1 do begin + for i:=1 to hLen + (len and $1) - 1 do begin (z+i*schritt)^:=res[i]; (z+(len-i)*schritt)^:=ims[i] * inv; end; if not odd(len) then - (z+hlen*schritt)^:=res[hlen]; + (z+hLen*schritt)^:=res[hLen]; end; doBetr: for i:=0 to len-1 do @@ -181,16 +183,18 @@ end; {$ENDIF} {$IFDEF fftunitGetrenntSpeichern} -//procedure tFFTAlgorithmus.speichern(invers: boolean; zRe,zIm: tExtendedArray); +//procedure tFFTAlgorithmus.speichern(invers: boolean; zRe,zIm: pSingle; schritt: longint); var i,inv: longint; begin inv:=1-2*byte(invers); - case outo of - doGetrennt: + case outO of + doGetrennt,doAlleResIms,doAlleResSmi: + // Real- und Imaginärteile sind _immer_ in gleicher Reihenfolge; + // die ggf. verdrehte Dimension ist senkrecht dazu! (benötigt für fft2d) for i:=0 to len-1 do begin - zRe[i]:=res[i]; - zIm[i]:=ims[i] * inv; + (zRe+i*schritt)^:=res[i]; + (zIm+i*schritt)^:=ims[i] * inv; end; else fehler('Diese Outputordnung ist in dieser Methode nicht implementiert!');; diff --git a/fftunit.pas b/fftunit.pas index b4b03bf..7b5995f 100644 --- a/fftunit.pas +++ b/fftunit.pas @@ -14,9 +14,9 @@ type tFFTAlgorithmus = class private - len,hlen: longint; + len,hLen: longint; res,ims: tExtendedArray; - ino,outo: tFFTDatenordnung; + inO,outO: tFFTDatenordnung; procedure dumpResIms; public constructor create(laenge: longint; inputOrdnung,outputOrdnung: tFFTDatenordnung); overload; @@ -25,14 +25,20 @@ type procedure laden(invers: boolean; q: pSingle; schritt: longint); overload; virtual; procedure laden(invers: boolean; q: pDouble; schritt: longint); overload; virtual; procedure laden(invers: boolean; q: pExtended; schritt: longint); overload; virtual; - procedure laden(invers: boolean; qRe,qIm: tExtendedArray); overload; virtual; + procedure laden(invers: boolean; qRe,qIm: pSingle; schritt: longint); overload; virtual; + procedure laden(invers: boolean; qRe,qIm: pDouble; schritt: longint); overload; virtual; + procedure laden(invers: boolean; qRe,qIm: pExtended; schritt: longint); overload; virtual; procedure ausfuehren; virtual; abstract; procedure speichern(invers: boolean; z: pSingle; schritt: longint); overload; virtual; procedure speichern(invers: boolean; z: pDouble; schritt: longint); overload; virtual; procedure speichern(invers: boolean; z: pExtended; schritt: longint); overload; virtual; - procedure speichern(invers: boolean; zRe,zIm: tExtendedArray); overload; virtual; + procedure speichern(invers: boolean; zRe,zIm: pSingle; schritt: longint); overload; virtual; + procedure speichern(invers: boolean; zRe,zIm: pDouble; schritt: longint); overload; virtual; + procedure speichern(invers: boolean; zRe,zIm: pExtended; schritt: longint); overload; virtual; procedure summen(amAnfang: boolean; var energie: extended; var nurNullen: boolean); function dumpParams: string; + property inOrdnung: tFFTDatenordnung read inO; + property outOrdnung: tFFTDatenordnung read outO; end; tFFTCooleyTukey = class(tFFTAlgorithmus) @@ -46,7 +52,9 @@ type procedure laden(invers: boolean; q: pSingle; schritt: longint); overload; override; procedure laden(invers: boolean; q: pDouble; schritt: longint); overload; override; procedure laden(invers: boolean; q: pExtended; schritt: longint); overload; override; - procedure laden(invers: boolean; qRe,qIm: tExtendedArray); overload; override; + procedure laden(invers: boolean; qRe,qIm: pSingle; schritt: longint); overload; override; + procedure laden(invers: boolean; qRe,qIm: pDouble; schritt: longint); overload; override; + procedure laden(invers: boolean; qRe,qIm: pExtended; schritt: longint); overload; override; procedure ausfuehren; override; end; @@ -157,10 +165,10 @@ end; constructor tFFTAlgorithmus.create(laenge: longint; inputOrdnung,outputOrdnung: tFFTDatenordnung); begin inherited create; - ino:=inputOrdnung; - outo:=outputOrdnung; + inO:=inputOrdnung; + outO:=outputOrdnung; len:=laenge; - hlen:=laenge div 2; + hLen:=laenge div 2; fillchar(res,sizeOf(res),#0); setLength(res,0); fillchar(ims,sizeOf(ims),#0); @@ -170,10 +178,10 @@ end; constructor tFFTAlgorithmus.create(original: tFFTAlgorithmus); begin inherited create; - ino:=original.ino; - outo:=original.outo; + inO:=original.inO; + outO:=original.outO; len:=original.len; - hlen:=original.hlen; + hLen:=original.hLen; fillchar(res,sizeOf(res),#0); setLength(res,0); fillchar(ims,sizeOf(ims),#0); @@ -207,7 +215,11 @@ procedure tFFTAlgorithmus.laden(invers: boolean; q: pExtended; schritt: longint) {$UNDEF fftunitLaden} {$DEFINE fftunitGetrenntLaden} -procedure tFFTAlgorithmus.laden(invers: boolean; qRe,qIm: tExtendedArray); +procedure tFFTAlgorithmus.laden(invers: boolean; qRe,qIm: pSingle; schritt: longint); +{$INCLUDE fftunit.inc} +procedure tFFTAlgorithmus.laden(invers: boolean; qRe,qIm: pDouble; schritt: longint); +{$INCLUDE fftunit.inc} +procedure tFFTAlgorithmus.laden(invers: boolean; qRe,qIm: pExtended; schritt: longint); {$INCLUDE fftunit.inc} {$UNDEF fftunitGetrenntLaden} @@ -221,7 +233,11 @@ procedure tFFTAlgorithmus.speichern(invers: boolean; z: pExtended; schritt: long {$UNDEF fftunitSpeichern} {$DEFINE fftunitGetrenntSpeichern} -procedure tFFTAlgorithmus.speichern(invers: boolean; zRe,zIm: tExtendedArray); +procedure tFFTAlgorithmus.speichern(invers: boolean; zRe,zIm: pSingle; schritt: longint); +{$INCLUDE fftunit.inc} +procedure tFFTAlgorithmus.speichern(invers: boolean; zRe,zIm: pDouble; schritt: longint); +{$INCLUDE fftunit.inc} +procedure tFFTAlgorithmus.speichern(invers: boolean; zRe,zIm: pExtended; schritt: longint); {$INCLUDE fftunit.inc} {$UNDEF fftunitGetrenntSpeichern} @@ -245,7 +261,7 @@ end; function tFFTAlgorithmus.dumpParams: string; begin - result:=className+': '+intToStr(len)+' ('+intToStr(hlen)+')'; + result:=className+': '+intToStr(len)+' ('+intToStr(hLen)+')'; end; {$UNDEF tFFTAlgorithmus} @@ -261,8 +277,8 @@ begin setLength(res,len); setLength(ims,len); loglen:=round(ln(len)/ln(2)); - setLength(wRe,hlen); - setLength(wIm,hlen); + setLength(wRe,hLen); + setLength(wIm,hLen); setLength(perm,len); for i:=0 to len-1 do begin res[i]:=nan; @@ -271,7 +287,7 @@ begin for j:=0 to loglen-1 do perm[i]:=perm[i] or (byte(odd(i shr j)) shl (loglen-1-j)); end; - for i:=0 to hlen-1 do begin + for i:=0 to hLen-1 do begin wRe[i]:=cos(2*pi*i/len); wIm[i]:=sin(2*pi*i/len); end; @@ -308,36 +324,40 @@ procedure tFFTCooleyTukey.laden(invers: boolean; q: pExtended; schritt: longint) {$UNDEF fftunitLaden} {$DEFINE fftunitGetrenntLaden} -procedure tFFTCooleyTukey.laden(invers: boolean; qRe,qIm: tExtendedArray); +procedure tFFTCooleyTukey.laden(invers: boolean; qRe,qIm: pSingle; schritt: longint); +{$INCLUDE fftunit.inc} +procedure tFFTCooleyTukey.laden(invers: boolean; qRe,qIm: pDouble; schritt: longint); +{$INCLUDE fftunit.inc} +procedure tFFTCooleyTukey.laden(invers: boolean; qRe,qIm: pExtended; schritt: longint); {$INCLUDE fftunit.inc} {$UNDEF fftunitGetrenntLaden} procedure tFFTCooleyTukey.ausfuehren; var - i,dist,absch,wnum,wstep: longint; + i,dist,absch,wNum,wStep: longint; t1,t2: extended; begin dist:=1; - wstep:=hlen; + wStep:=hLen; while dist<len do begin absch:=0; while absch<len do begin - wnum:=0; + wNum:=0; for i:=0 to dist-1 do begin // x_links: [absch+j] // x_rechts: [absch+j+dist] - t1:=wRe[wnum]*res[absch+i+dist] - wIm[wnum]*ims[absch+i+dist]; - t2:=wRe[wnum]*ims[absch+i+dist] + wIm[wnum]*res[absch+i+dist]; + t1:=wRe[wNum]*res[absch+i+dist] - wIm[wNum]*ims[absch+i+dist]; + t2:=wRe[wNum]*ims[absch+i+dist] + wIm[wNum]*res[absch+i+dist]; res[absch+i+dist]:=res[absch+i]-t1; ims[absch+i+dist]:=ims[absch+i]-t2; res[absch+i]:=res[absch+i]+t1; ims[absch+i]:=ims[absch+i]+t2; - wnum:=wnum+wstep; + wNum:=wNum+wStep; end; absch:=absch+2*dist; end; - wstep:=wstep div 2; + wStep:=wStep div 2; dist:=dist*2; end; end; @@ -368,9 +388,9 @@ begin end; preFFT:=tFFTCooleyTukey.create(length(rRe),doGetrennt,doGetrennt); - preFFT.laden(false,rRe,rIm); + preFFT.laden(false,pExtended(@(rRe[0])),pExtended(@(rIm[0])),1); preFFT.ausfuehren; - preFFT.speichern(false,sRe,sIm); + preFFT.speichern(false,pExtended(@(sRe[0])),pExtended(@(sIm[0])),1); preFFT.free; {$IFDEF fftNormierung} @@ -427,9 +447,9 @@ begin end; // a[j] := iFFT(a[j]) - subFFT.laden(true,res,ims); + subFFT.laden(true,pExtended(@(res[0])),pExtended(@(ims[0])),1); subFFT.ausfuehren; - subFFT.speichern(true,res,ims); + subFFT.speichern(true,pExtended(@(res[0])),pExtended(@(ims[0])),1); for i:=0 to length(res)-1 do begin res[i]:=res[i]; @@ -444,9 +464,9 @@ begin end; // a[j] := iFFT(a[j]) - subFFT.laden(false,res,ims); + subFFT.laden(false,pExtended(@(res[0])),pExtended(@(ims[0])),1); subFFT.ausfuehren; - subFFT.speichern(false,res,ims); + subFFT.speichern(false,pExtended(@(res[0])),pExtended(@(ims[0])),1); for i:=0 to length(res)-1 do begin res[i]:=res[i]; |