;
;Schnelle Fouriertrafo mit
;genau 256 Stuetzstellen
;haftmann#software
;erstellt 15.10.1993
;geändert 21.10.1993
;
;Quelle: Taschenbuch mathematischer
;Formeln und moderner Verfahren,
;Verlag Harri Deutsch Thun
;Hrsg: Prof. Dr. Horst Stöcker
;Seiten 620 und 621 (Pseudocode)
;[maschinennah optimiert]
;=================================
os equ 0f003h
anfang equ 7000h
org anfang
;Sprungtabelle u. Übergabezelle
jr plotb ;anfang
what db tabhi ;anfang+2
jp fft ;anfang+3
jp coswin ;anfang+6
jp betrag ;anfang+9
;Schnittstellenbeschreibung
;(ShortInt = -128..127)
;(Integer = -32768..32767)
;
;*****************
;PLOT (Anfang +0)
;Zeichnet ein Oszillogramm in die
;linke "Bildschirmhälfte" unter
;Verwendung von CAOS-LINE
;
;PE: *(Anfang+2): High-Adresse
;der abzubildenden 100h-Tabelle
; *eine 100h lange Tabelle
;mit ShortInt-Werten belegt
;
;***************
;FFT (Anfang +3)
;Führt mit 256 komplexen ShortInt-
;Daten eine FFT aus. Das Ergebnis
;sind 256 komplexe Integer-Daten
;(=400h Bytes), also die A(i) und
;B(i).
;
;PE: *in (Anfang +400h) die
;Realteile als ShortInt
; *in (Anfang +500h) die
;Imaginaerteile als ShortInt
;
;PA: *in (Anfang +400h) die
;Low-Realteile als Integer
; *in (Anfang +500h) die
;Low-Imaginärteile als Integer
; *in (Anfang +600h) die
;High-Realteile als Integer
; *in (Anfang +700h) die
;High-Imaginärteile als Integer
;
;Hinweis: Die FFT-Routine wandelt
;als erstes die angebotenen ShortInts
;zu Integers und rechnet intern damit.
;
;******************
;COSWIN (Anfang +6)
;Legt über eine 100h lange ShortInt-
;Tabelle ein Cosinusquadratfenster
;(auch Hanning-Fenster genannt); d.h.
;alle Werte werden mit der COS^2-Funk-
;tion multipliziert.
;
;PE: *(Anfang +2): High-Adresse
;der zu fensternden 100h-Tabelle
; *die Tabelle selbst mit 100h
;ShortInt-Werten
;
;PA: *dieselbe Tabelle gefenstert
;
;******************
;BETRAG (Anfang +9)
;Rechnet die komplexen Ausgabewerte
;der FFT in Bytes um, bereitet die
;Rücktransformation vor und erzeugt
;eine Tabelle mit Beträgen (als Byte
;= unsigned ShortInt) und eine Ta-
;belle mit Phasen.
;Der Wert in (Anfang+2) wird als ge-
;meinsamer Exponent gewertet und ent-
;sprechend verschoben. Man beachte,
;daß bei Hin- und Rücktrafo sich die-
;ser Wert ganz normal um 8 erhöht, da
;die Normalisierung bei der Hintrafo
;nicht durchgeführt wird. Dadurch wird
;das Clipping der Vorgängerversion
;vermieden und trotzdem die Datenbrei-
;te von 8 Bit garantiert.
;
;PE: *Siehe PA der FFT!
;
;PA: *(Anfang +400h): 100h Imagi-
;närteile(!!!) als ShortInt
; *(Anfang +500h): 100h Real-
;teile als ShortInt
; *(Anfang +600h): 100h Beträge
;als Byte-80h (fuer PLOT günstig)
; *(Anfang +700h): 100h Phasen
;als ShortInt
;-128=-Pi/2, 0=0, 127= knapp Pi/2
;
;Kernel
;*** Unterfunktion PLOT fuer BASIC ***
plotb in a,88h
bit 2,a
jr nz,plot ;fuer MC-Progs
call 0f018h
call plot
call 0f01bh
ret
plot
ld hl,0b782h
ld b,8
xor a
plo1 ld m,a ;Args zu Null
inc l
djnz plo1
ld e,a
ld a,(what)
ld d,a
ld a,(de)
add 80h
ld l,88h
ld m,a
inc e
jr plot2
plot1
ld l,82h
inc m
plot2 ld l,86h
inc m
ld l,88h
ld a,m
ld l,84h
ld m,a
ld a,(de)
add 80h
ld l,88h
ld m,a
push hl
push de
call os
db 3eh
pop de
pop hl
inc e
jr nz,plot1
ret ;(plot)
;
;*** Die eigentliche FFT ***
fft
;Vorbereitung: High-Teile Vorzeichenset!
ld h,tabhi
call xcbw
inc h
call xcbw
dec h
;
ld a,80h ;hy=n2
ld b,1 ;b=argument
fftl1 ;Äußere Schleife
;8fach Durchlauf
ld hy,a
ld c,0 ;c=angle
ld l,0 ;l=j
fftl2 ;Mittlere Schleife
push hl
ld a,c
call cosinus
ld (cos),hl
ld a,c
call sinus
ld (sin),hl
pop hl
push hl ;l=i
ld a,l
add hy
ld e,a
ld d,h
push bc ;angle und arg
fftl3 ;Innere Schleife
;i=l
;hl=^x(i)
;de=^x(l)
ld c,m
inc h
inc h
ld b,m
push bc ;Stack:x(i)
ex de,hl
ld c,m
inc h
inc h
ld b,m ;bc:x(l)
ex (sp),hl ;Stack:^x(l)
or a
sbc hl,bc
ld (xdum),hl
add hl,bc
add hl,bc
ex de,hl ;hl:^x(i)
ld m,d
dec h
dec h
ld m,e
pop de
;
inc d
inc h
ld c,m
inc h
inc h
ld b,m
push bc ;Stack:y(i)
ex de,hl
ld b,m
dec h
dec h
ld c,m ;bc:y(l)
ex (sp),hl ;Stack:^y(l)
or a
sbc hl,bc
ld (ydum),hl
add hl,bc
add hl,bc
ex de,hl ;hl:^y(i)
ld m,d
dec h
dec h
ld m,e
pop de
dec h
dec d
;
push hl
push de
ld de,(xdum)
ld bc,(cos)
call mult
push hl
ld de,(ydum)
ld bc,(sin)
call mult
pop de
add hl,de
pop de
ex de,hl
ld m,e
inc h
inc h
ld m,d
push hl ;hl:^y(l)
ld de,(xdum)
; ld bc,(sin)
call mult
push hl
ld de,(ydum)
ld bc,(cos)
call mult
pop de
or a
sbc hl,de
pop de
ex de,hl
inc h
ld m,d
dec h
dec h
ld m,e
dec h
ex de,hl
pop hl
;
; call os
; db 23h,'i',0
ld a,hy
add a ;n2 2mal=n1
add e
ld e,a
ld a,l
add hy ;n2 2mal=n1
jr c,fftl3e
add hy ;Ueberlauf?
ld l,a
jp nc,fftl3
fftl3e ;Innere Schleife ENDE
; call os
; db 23h,'j',0
pop bc
pop hl
inc l ;j=Startindex
ld a,c
add b
ld c,a
rla
jp nc,fftl2
;Mittlere Schleife ENDE
; call os
; db 23h,'k',0
rlc b ;Argument
ld a,hy
rrca ;n2 halbieren
jp nc,fftl1
;Äußere Schleife ENDE
;*** fft 2.Teil ***
;z.Z.D=H=tabhi
;de=^x(j)
;hl=^x(i)
ld l,0 ;l=i=0
fftl4 ;Schleife 4
;Bitreversion (Ggs. Quelle)
ld a,l
ld b,8
brev rrca
rl e
djnz brev
;a ist noch l !
cp e ;i<j?
jr nc,fftv1
ld b,4
fftx ld c,m
ld a,(de)
ex de,hl
ld m,c
ld (de),a
inc d
inc h
djnz fftx
ld h,tabhi
ld d,h
fftv1
inc l ;i++
jr nz,fftl4
ret ;(END FFT)
;
;*** Cosinusquadratfenster ***
coswin
ld a,(what)
ld ly,0
ld hy,a
ld hl,0 ;1.Sinuswert
coswl1
push hl ;Sinuswert
call coswup
ld a,ly
inc a
srl a
call sinus
pop de
push hl
add hl,de
srl h
rr l ;Arit.Mittel
call coswup
pop hl ;Sinuswert
jr nz,coswl1
ret
;
coswup ld b,h
ld c,l
ld d,h
ld e,l
push iy
call mult
pop iy
ld c,l
ld b,h
ld e,(iy)
ld a,e
rla
sbc a
ld d,a ;cbw
push iy
call mult
pop iy
ld (iy),l
inc ly
ret
;
;*** Betragsbildung und ***
;*** Vorbereitung Rücktrafo ***
betrag
ld hl,tab
ld bc,0
call maxseek
inc h
call maxseek
dec h
ld a,0feh
maxsk scf
rl c
rl b
inc a
jr nc,maxsk
;A=Anz. mögliche Schiebungen
cp 0ffh
ccf
adc 0 ;ff=>0
call shift
inc h
call shift
dec h
ld c,a
ld a,(what)
add c ;Neuer Exp.
ld (what),a
inc h
ld d,h
ld e,l
inc h
call copy
inc h
dec d
call copy
inc d ;alter RealTeil
betrag1 ld a,(de)
call quad
ld a,(de)
ld c,a
dec d
ld a,(de) ;alter ImagTeil
ld b,a
push bc
push hl
call quad
pop bc
add hl,bc ;läuft nie über
call sqrt ;Pythagoras
inc d
inc d
ld h,a
sub 80h ;anzeigegünstig
ld (de),a ;Betrag
pop bc ;c=real b=imag
ld a,b
or c
jr z,betrag4 ;Null!
ld a,b
call abs ;a:=|a|
rl l ;VZ Imag
ld b,a
ld a,c
call abs ;a:=|a|
rl l ;VZ Real
cp b
rl l ;wenn Real<Imag
cp b ;CY noch mal
jr z,betrag2
jr c,betrag5
ld a,b
betrag5 call div ;a:=a/h, a<h
call arcsin ;a:=arcsin(a),
ld h,a ;a<=b5=√2/2
jr betr2a
betrag2 ld h,20h ;45°
betr2a ld bc,340h ;FU-Schleife
betrag6 rr l
jr c,betrag7
ld a,c
sub h
ld h,a
betrag7 sla c
djnz betrag6
betrag4 inc d
ld (de),a
dec d
dec d
; call os
; db 23h,'l',0
; call os
; db 1bh
inc e
jr nz,betrag1
ret
;
copy ;256 Bytes LDIR
push bc
ld bc,256
ldir
dec h
dec d
pop bc
ret
;
div ;a:=a/h, a<h (Nachkommastellen)
push bc
ld b,8
div1 add a
jr c,div2
cp h
jr nc,div2
sla c
jr div3
div2 sub h
sls c ;1 von rechts
div3 djnz div1
ld a,c
pop bc
ret
;
arcsin ;a:=arcsin(a)
;VR:AF
push hl
ld hl,sintab
arcs1 inc l
cp m
jr nc,arcs1
add a
sub m
dec l
sub m
rla
jr c,arcs2
inc l
arcs2 ld a,l
and 3fh
pop hl
ret
;
shift ;256 Bytes ShiftLeft um A Bits
ld e,m
inc h
inc h
ld d,m
or a
jr z,shift3;Nur runden
ex de,hl
ld b,a
shift1 add hl,hl
djnz shift1
ex de,hl
shift3 rl e ;runden
jr nc,shift2
inc d
shift2 ld m,d
dec h
dec h ;Kein Low-Teil
inc l
jr nz,shift
ret
;
maxseek ;Maximumsuche (Wert, NICHT
;Stelle! BC=Maximum (PE&PA)
ld e,m
inc h
inc h
ld d,m
dec h
dec h
bit 7,d ;negativ?
jr z,maxsk1
push hl ;"NEG DE"
ld hl,0
or a
sbc hl,de
ex de,hl
pop hl
maxsk1 ex de,hl
or a
sbc hl,bc
jr c,maxsk2;kleiner
add hl,bc ;Korrektur
ld c,l
ld b,h ;neues Max
maxsk2 ex de,hl
inc l
jr nz,maxseek
ret
;
;*** Unterprogramme ***
xcbw ;Byte-zu-Word-Trafo
ld l,0
xcbw1 ld a,m
rla
sbc a ;wie CBW
inc h
inc h
ld m,a
dec h
dec h
inc l
jr nz,xcbw1
ret
;
cosinus ;Cosinus
add 40h
sinus ;Sinuswert ermitteln (Tabelle)
;hl:=sin(a)
;hl-Aufbau: V000000Y YYYYYYYY
;(Vorzeichen - Betrag)
ld d,0
ld l,a
rla ;Bit 7
rr d ;=Vorzeichen
rla ;Bit 6
ld a,l
jr nc,sinu1
cpl
inc a
and 3fh ;Gipfel?
jr nz,sinu1
set 0,d ;Gipfel!!!
sinu1
or 0c0h
cp 0feh
jr c,sinu2
set 0,d ;100h
sinu2 ld l,a
ld h,sintab/256
ld e,m
ex de,hl
ret
;
quad ;Quadrat hl:=a^2 vzb=>vzl
;VR:A,B,HL
call abs
quad1 push de
ld b,8
ld d,0
ld e,a
ld h,d
ld l,d ;hl=0
quad2 add hl,hl
rla
jr nc,quad3
add hl,de
quad3 djnz quad2
pop de
ret
;
abs ;a:=abs(a); cy=1:negativ
or a
ret p ;positiv
neg
ret
;
sqrt ;Wurzel a:=sqrt(hl) vzl=>vzl
xor a
dec a
push de
ld e,a
ld d,a
sqr1 inc a
inc de
inc de
sbc hl,de
jr nc,sqr1
sra d ;halbieren
rr e ;cy=1
add hl,de
bit 7,h
jr nz,sqr2
inc a ;aufrunden
jr nz,sqr2 ;100h?
dec a ;ffh!
sqr2 pop de
ret
;
mult ;Multiplikation
;bc: 16bit Vorzeichen-Betrag
;de: 16bit ZWK
;hl: Ergebnis
;VR:hl,de,af,ly
ld hl,0
bit 0,b
jr z,mul1
add hl,de
jr mul2 ;praktisch!
mul1 push bc
xor a
ld b,a
ld ly,8
mul3 sra d
rr e
rr b
rlc c
jr nc,mul4
add b
adc hl,de
mul4 dec ly
jr nz,mul3
pop bc
rla ;aufrunden?
jr nc,mul2
inc hl
mul2 bit 7,b
ret z
ld e,l
ld d,h
ld hl,0
or a
sbc hl,de
ret
;
;*** Testteil ***
db 7fh,7fh,'SINUS',1
ld a,l
call sinus
outhl call os
db 1ah
ret
db 7fh,7fh,'MULT',1
ld b,h
ld c,l
call mult
jr outhl
db 7fh,7fh,'QUAD',1
ld a,l
call quad
jr outhl
db 7fh,7fh,'SQRT',1
call sqrt
outa ld l,a
ld h,0
jr outhl
db 7fh,7fh,'ARCSIN',1
ld a,l
call arcsin
jr outa
db 7fh,7fh,'JPHL',1
jp (hl)
;
;*** Variablen ***
xdum dw 0
ydum dw 0
sin dw 0
cos dw 0
;*** Tabellenteil ***
ds anfang+3c0h-$
sintab
DB 00H,06H,0DH,13H
DB 19H,1FH,26H,2CH
DB 32H,38H,3EH,44H
DB 4AH,50H,56H,5CH
DB 62H,68H,6DH,73H
DB 79H,7EH,84H,89H
DB 8EH,93H,98H,9DH
DB 0A2H,0A7H,0ACH,0B1H
DB 0B5H,0B9H,0BEH,0C2H
DB 0C6H,0CAH,0CEH,0D1H
DB 0D5H,0D8H,0DCH,0DFH
DB 0E2H,0E5H,0E7H,0EAH
DB 0EDH,0EFH,0F1H,0F3H
DB 0F5H,0F7H,0F8H,0FAH
DB 0FBH,0FCH,0FDH,0FEH
DB 0FFH,0FFH,00H,00H
;Datentabellen ab anfang+400h
tab
tabhi equ $/256
; Hintrafo Rücktrafo
;+400h: X Low B Low
;+500h: Y Low A Low
;+600h: X High Betraege
;+700h: Y High (Phasen)
END
Ç DEVEX KCC Ç TEXOREX KCC Ç FORTHEX KCC
| Detected encoding: OEM (CP437) | 1
|
|
|