; DISCLAIMER
; 
; Copyright (C) 1997 by Adrian Philip Whichello. All rights reserved.
; 
; I, Adrian Philip Whichello, MAKE NO WARRANTY ON THIS SOFTWARE, EITHER EXPRESS
; OR IMPLIED, WITH RESPECT TO QUALITY, ACCURACY, MERCHANTABILITY, OR FITNESS FOR
; A PARTICULAR PURPOSE.
; 
; IN NO EVENT WILL I BE HELD RESPONSIBLE FOR DIRECT, INDIRECT, SPECIAL,
; INCIDENTAL, OR CONSEQUENTIAL DAMAGES RESULTING FROM ANY DEFECT OF THE
; SOFTWARE.
; 
; I RELEASE THIS PROGRAM AS FREEWARE. YOU MAY COPY, GIVE AWAY AND USE THIS
; PROGRAM, PROVIDED THIS NOTICE REMAINS ATTACHED AT ALL TIMES.
; 

return  .equ 0
real    .equ return+2
imag    .equ real+2
N       .equ imag+2
sign    .equ N+2
scalfac .equ sign+2
i       .equ scalfac+2
j       .equ i+1
costab  .equ j+1
mplier  .equ costab+2
micand  .equ mplier+2
result  .equ micand+2
icos    .equ result+2
isin    .equ icos+1
tmp_re  .equ isin+1
tmp_im  .equ tmp_re+1
k       .equ tmp_im+1
prcount .equ k+1
celnum  .equ prcount+1
delta   .equ celnum+1
celdist .equ delta+1
celcount.equ celdist+1
pairnum .equ celcount+1
log2_N  .equ pairnum+1
N_1     .equ log2_N+1

cout    .equ 0C305

        .macro pop
        pla
        sta %1
        pla
        sta %1+1
        .endm
        
        .macro push
        lda %1+1
        pha
        lda %1
        pha
        .endm
        
        .proc byte_fft, 5
        pop return
        pop scalfac
        pop sign
        pop N
        pop imag
        pop real
        
        ; pre-transform bit swap
        
        lda #0FF        ; -1
        sta log2_N
        lda N
        sta i
        sta N_1         ; useful later
        dec N_1
        lda N+1
        sta j
$1      lsr j           ; find power of 2
        ror i
        inc log2_N
        bcc $1
        
        ldy #1
bitrev  lda #0
        sta j
        tya
        clc
        ldx log2_N
$1      ror a
        rol j
        dex
        bne $1
        cpy j
        bmi $2
        sty i
        lda (real),y
        pha
        lda (imag),y
        pha
        ldy j
        lda (real),y
        pha
        lda (imag),y
        ldy i
        sta (imag),y
        pla
        sta (real),y
        ldy j
        pla
        sta (imag),y
        pla
        sta (real),y
        ldy i
$2      iny
        cpy N_1
        bne bitrev
        
        lda cosadr
        sta costab
        lda cosadr+1
        sta costab+1
        lda #0
        tay
        sta (scalfac),y
        iny
        sta (scalfac),y
        
fft     lda N+1
        lsr a
        lda N
        ror a
        sta celnum
        lda #080
        sta delta
        ldx #1
        stx celdist
        stx pairnum
npass   jsr scale
        lda celnum
        sta celcount
        lda #0
        sta i
ncell   ldx #0
        stx k
        ldx pairnum
nc1     lda i
        clc
        adc celdist
        sta j
        stx prcount
        
        ldy k
        lda (costab),y
        sta icos
        lda sign
        bpl $2
        tya
        clc
        adc #040  ; go forward 64 places in cos table to get -sine
        jmp $3
$2      tya
        sec
        sbc #040  ; go back 64 places in cos table to get +sine
$3      tay
        lda (costab),y
        sta isin
        
        ldy j
        lda (real),y
        pha
        sta mplier
        lda icos
        sta micand
        jsr multply
        sta tmp_re      ;  tmp_re = real(j) * cos
        pla
        sta mplier
        lda isin
        sta micand
        jsr multply
        sta tmp_im      ;  tmp_im = real(j) * sign * -sin
        
        lda (imag),y
        pha
        sta mplier
        lda isin
        sta micand
        jsr multply
        clc
        adc tmp_re
        sta tmp_re      ;  tmp_re = real(j) * cos + imag(j) * sign * -sin
        pla
        sta mplier
        lda icos
        sta micand
        jsr multply
        sec
        sbc tmp_im
        sta tmp_im      ;  tmp_im = imag(j) * cos - real(j) * sign * -sin
        
        ldy i
        lda (imag),y
        pha
        clc
        adc tmp_im
        sta (imag),y    ; imag(i) = imag(i) + tmp_im
        lda (real),y
        pha
        clc
        adc tmp_re
        sta (real),y    ; real(i) = real(i) + tmp_re
        pla
        sec
        sbc tmp_re
        ldy j
        sta (real),y    ; real(j) = real(i) - tmp_re
        pla
        sec
        sbc tmp_im
        sta (imag),y    ; imag(j) = imag(i) - tmp_im
        
        lda k
        clc
        adc delta
        sta k
        inc i
        ldx prcount
        dex
        beq $1
        jmp nc1
$1      lda i
        clc
        adc celdist
        sta i
        dec celcount
        beq np1
        jmp ncell
        
np1     lsr celnum
        beq exit
        asl pairnum
        asl celdist
        lsr delta
        jmp npass
        
exit    push return
        rts
        
cosadr  .word cosine

cosine  .byte  127.,  127.,  127.,  127.,  127.,  127.,  127.,  126.
        .byte  126.,  125.,  124.,  123.,  122.,  122.,  121.,  119.
        .byte  118.,  117.,  116.,  114.,  113.,  111.,  110.,  108.
        .byte  106.,  105.,  103.,  101.,   99.,   97.,   95.,   93.
        .byte   91.,   88.,   86.,   84.,   81.,   79.,   76.,   74.
        .byte   71.,   68.,   66.,   63.,   60.,   58.,   55.,   52.
        .byte   49.,   46.,   43.,   40.,   37.,   34.,   31.,   28.
        .byte   25.,   22.,   19.,   16.,   13.,    9.,    6.,    3.
        .byte    0.,   -3.,   -6.,   -9.,  -13.,  -16.,  -19.,  -22.
        .byte  -25.,  -28.,  -31.,  -34.,  -37.,  -40.,  -43.,  -46.
        .byte  -49.,  -52.,  -55.,  -58.,  -60.,  -63.,  -66.,  -68.
        .byte  -71.,  -74.,  -76.,  -79.,  -81.,  -84.,  -86.,  -88.
        .byte  -91.,  -93.,  -95.,  -97.,  -99., -101., -103., -105.
        .byte -106., -108., -110., -111., -113., -114., -116., -117.
        .byte -118., -119., -121., -122., -122., -123., -124., -125.
        .byte -126., -126., -127., -127., -127., -128., -128., -128.
        .byte -128., -128., -128., -128., -127., -127., -127., -126.
        .byte -126., -125., -124., -123., -122., -122., -121., -119.
        .byte -118., -117., -116., -114., -113., -111., -110., -108.
        .byte -106., -105., -103., -101.,  -99.,  -97.,  -95.,  -93.
        .byte  -91.,  -88.,  -86.,  -84.,  -81.,  -79.,  -76.,  -74.
        .byte  -71.,  -68.,  -66.,  -63.,  -60.,  -58.,  -55.,  -52.
        .byte  -49.,  -46.,  -43.,  -40.,  -37.,  -34.,  -31.,  -28.
        .byte  -25.,  -22.,  -19.,  -16.,  -13.,   -9.,   -6.,   -3.
        .byte    0.,    3.,    6.,    9.,   13.,   16.,   19.,   22.
        .byte   25.,   28.,   31.,   34.,   37.,   40.,   43.,   46.
        .byte   49.,   52.,   55.,   58.,   60.,   63.,   66.,   68.
        .byte   71.,   74.,   76.,   79.,   81.,   84.,   86.,   88.
        .byte   91.,   93.,   95.,   97.,   99.,  101.,  103.,  105.
        .byte  106.,  108.,  110.,  111.,  113.,  114.,  116.,  117.
        .byte  118.,  119.,  121.,  122.,  122.,  123.,  124.,  125.
        .byte  126.,  126.,  127.,  127.,  127.,  127.,  127.,  127.
        
multply lda #0
        sta mplier+1
        sta micand+1    ; clear MSB's
        sta result
        sta result+1
        lda micand
        beq $9
        bpl $1
        lda #0FF
        sta micand+1
$1      lda mplier
        beq $9
        bpl $2
        lda #0FF
        sta mplier+1
$2      ldx #0F
$3      lsr mplier+1
        ror mplier
        bcc $4
        lda micand
        clc
        adc result
        sta result
        lda micand+1
        adc result+1
        sta result+1
$4      asl micand
        rol micand+1
        dex
        bne $3
        
$9      lda result+1
        rol result
        rol a
        rts

scale   ldy #0
$2      lda (real),y
        beq $1
        cmp #0C0
        bpl $1
        cmp #040
        bmi $1
        jmp shift
        lda (imag),y
        beq $1
        cmp #0C0
        bpl $1
        cmp #040
        bmi $1
        jmp shift
$1      iny
        cpy N
        bne $2
        rts

shift   ldy #0
        lda (scalfac),y
        clc
        adc #1
        sta (scalfac),y
$3      lda (real),y
        clc
        adc #080
        lsr a
        sec
        sbc #040
        sta (real),y
        lda (imag),y
        clc
        adc #080
        lsr a
        sec
        sbc #040
        sta (imag),y
        iny
        cpy N
        bne $3
        rts

        .end