home *** CD-ROM | disk | FTP | other *** search
/ Club Amiga de Montreal - CAM / CAM_CD_1.iso / files / 326.lha / KFFT_v1.1 / fft.asm < prev    next >
Assembly Source File  |  1989-12-23  |  9KB  |  362 lines
  1. \ FFT.ASM - Fast Fourier Transform assembly support words.
  2. \
  3. \ KFFT V1.1 (C)Copyright 1989, Jerry Kallaus.  All rights reserved. 
  4. \ May be freely redistributed for non-commercial use (FREEWARE).
  5.  
  6. \ See file fft.asm.doc.
  7. \ These words assume that registors A0-A1,D0-D3 are trashable.
  8.  
  9. ANEW TASK-fft.asm
  10.  
  11. anew task-fft.asm
  12.  
  13. auto_scale_fft? NOT CONSTANT nonauto?
  14.  
  15. variable save-here
  16.  
  17. : bpl$+4  $ 6a02 w, ; immediate \ used to keep assembler from scribling
  18.                                 \ over erased code
  19.  
  20. : markhere ( -- ) here save-here ! ;
  21. : backhere ( flag -- , backup to markhere if flag )
  22.      if save-here @ here - allot then ;
  23.  
  24. asm 2**                ( n -- 2**n )
  25.     moveq.l #1,d0
  26.     asl.l   tos,d0
  27.     move.l  d0,tos
  28.     forth{ both }
  29. end-code
  30.  
  31. asm 2CELL+             ( n -- n+8 )
  32.     addq.l  #8,tos
  33.     forth{ both }
  34. end-code
  35.  
  36. asm 2CELL-             ( n -- n-8 )
  37.     subq.l  #8,tos
  38.     forth{ both }
  39. end-code
  40.  
  41. asm 2CELLS             ( n -- n*8 )
  42.     asl.l   #3,tos
  43.     forth{ both }
  44. end-code
  45.  
  46. asm 4DUP               ( 1 2 3 4 -- 1 2 3 4 1 2 3 4 )
  47.     move.l  dsp,a0
  48.     movem.l (a0)+,d1-d3
  49.     move.l  tos,-(dsp)
  50.     movem.l d1-d3,-(dsp)
  51.     forth{ both }
  52. end-code
  53.  
  54. asm Z@                 ( addr -- real imag )
  55.     move.l  $0(org,tos.l),-(dsp)
  56.     move.l  $4(org,tos.l),tos
  57.     forth{ both }
  58. end-code
  59.  
  60. asm Z!                 ( real imag addr -- )
  61.     move.l  (dsp)+,$4(org,tos.l)
  62.     move.l  (dsp)+,$0(org,tos.l)
  63.     move.l  (dsp)+,tos
  64.     forth{ both }
  65. end-code
  66.  
  67. asm Z+                    ( z1 z2 -- z1+z2 )
  68.     movem.l (dsp)+,d0-d2
  69.     add.l   d0,d2
  70.     add.l   d1,tos
  71.     move.l  d2,-(dsp)
  72.     forth{ both }
  73. end-code
  74.  
  75. asm Z-                    ( z1 z2 -- z1-z2 )
  76.     movem.l (dsp)+,d0-d2
  77.     sub.l   d0,d2
  78.     sub.l   d1,tos
  79.     neg.l   tos
  80.     move.l  d2,-(dsp)
  81.     forth{ both }
  82. end-code
  83.  
  84.  
  85. \ Z* - complex multiply -  ( a b c d -- ac-bd ad+bc )
  86. \      scaled 2**14
  87. \      registers  3 2 1 7 = a b c d
  88.  
  89. asm Z*                    ( z1 z2 -- z1*z2 )
  90.     movem.l (dsp)+,d1-d3
  91.     move.l  d3,d0         a
  92.     muls    d1,d0         ac
  93.     muls    d2,d1         bc
  94.     muls    tos,d2        bd
  95.     sub.l   d2,d0         ac-bd
  96.     asl.l   #2,d0
  97.     swap    d0
  98.     ext.l   d0
  99.     move.l  d0,-(dsp)
  100.     muls    d3,tos        ad
  101.     add.l   d1,tos        ad+bc
  102.     asl.l   #2,tos
  103.     swap    tos
  104.     ext.l   tos
  105.     forth{ both }
  106. end-code
  107.  
  108. asm ZNEGATE              ( z -- -z )
  109.     neg.l   tos
  110.     move.l  (dsp),d0
  111.     neg.l   d0
  112.     move.l  d0,(dsp)
  113.     forth{ both }
  114. end-code
  115.  
  116. asm Z2/                   ( z -- z/2 )
  117.     asr.l   #1,tos
  118.     move.l  (dsp),d0
  119.     asr.l   #1,d0
  120.     move.l  d0,(dsp)
  121.     forth{ both }
  122. end-code
  123.  
  124. asm Z/2**N                ( z n -- z/2**n )
  125.     movem.l (dsp)+,d0-d1
  126.     asr.l   tos,d0
  127.     asr.l   tos,d1
  128.     move.l  d0,tos
  129.     move.l  d1,-(dsp)
  130.     forth{ both }
  131. end-code
  132.  
  133. asm NSBITS    ( value --  number of significant abs bits plus 1 sign bit )
  134.     tst.l   tos
  135.     beq.s   9$          exception is return zero if value is zero
  136.     bgt.s   1$
  137.     not.l   tos
  138.     bgt.s   1$
  139.     moveq.l #1,tos
  140.     bra.s   9$
  141. 1$: moveq.l #33,d0
  142. 2$: asl.l   #1,tos
  143.     dblt.w  d0,2$
  144.     move.l  d0,tos
  145. 9$:
  146.     forth{ both }
  147. end-code
  148.  
  149.  
  150. asm OR.ABS.ARRAY     ( addr ncells -- or'd-magnitude-bits-of-array )
  151.     move.l  (dsp)+,a0
  152.     adda.l  org,a0
  153.     moveq.l #0,d1
  154.     move.l  tos,d2        64k counter
  155.     swap    d2
  156.     bra.s   3$
  157. 1$: move.l  (a0)+,d0
  158.     bpl.s   2$
  159.     not.l   d0
  160. 2$: or.l    d0,d1
  161. 3$: dbra.w  tos,1$
  162.     dbra.w  d2,1$
  163.     move.l  d1,tos
  164.     forth{ both }
  165. end-code
  166.  
  167.                       ( Arithmetic shift array of n-cells by n-bits. )
  168.                       ( Left for n-bits positive, right for n-bits neg.)
  169. asm ASHIFT.ARRAY      ( array-addr n-cells n-bits -- )
  170.     movem.l  d4-d5,-(rp)   ( Limited to arrays up to 256k cells )
  171.     move.l  (dsp)+,d4          n
  172.     ble.l   8$
  173.     move.l  d4,a1
  174.     lsr.l   #2,d4
  175.     move.l  (dsp)+,d0
  176.     lea     $0(org,d0.l),a0    addr
  177.     moveq.l #16,d5
  178.     tst.l   tos
  179.     beq.l   9$
  180.     bgt.l   4$
  181.     neg.l   tos
  182.     subq.l  #1,d4           start of right shift code
  183.     bmi.s   2$
  184. 1$: movem.l (a0)+,d0-d3
  185.     asr.l   tos,d0
  186.     asr.l   tos,d1
  187.     asr.l   tos,d2
  188.     asr.l   tos,d3
  189.     movem.l d0-d3,-(a0)
  190.     adda.l  d5,a0
  191.     dbra.w  d4,1$
  192. 2$: move.l  a1,d4
  193.     moveq.l #3,d0
  194.     and.l   d0,d4
  195.     beq.s   9$
  196.     subq.l  #1,d4
  197. 3$: move.l  (a0),d0
  198.     asr.l   tos,d0
  199.     move.l  d0,(a0)+
  200.     dbra.w  d4,3$
  201.     bra.s   9$
  202. 4$: subq.l  #1,d4           start of left shift code
  203.     bmi.s   6$
  204. 5$: movem.l (a0)+,d0-d3
  205.     asl.l   tos,d0
  206.     asl.l   tos,d1
  207.     asl.l   tos,d2
  208.     asl.l   tos,d3
  209.     movem.l d0-d3,-(a0)
  210.     adda.l  d5,a0
  211.     dbra.w  d4,5$
  212. 6$: move.l  a1,d4
  213.     moveq.l #3,d0
  214.     and.l   d0,d4
  215.     beq.s   9$
  216.     subq.l  #1,d4
  217. 7$: move.l  (a0),d0
  218.     asl.l   tos,d0
  219.     move.l  d0,(a0)+
  220.     dbra.w  d4,7$
  221.     bra.s   9$
  222. 8$: adda.w  #4,dsp         pop data addr off stack
  223. 9$: movem.l (rp)+,d4-d5
  224.     move.l  (dsp)+,tos
  225. end-code
  226.  
  227.  
  228. asm STATS.ARRAY       ( array-addr n -- max min sumlo sumhi )
  229.     movem.l d4-d5,-(rp)
  230.     move.l  (dsp)+,a0
  231.     adda.l  org,a0
  232.     move.l  tos,d0
  233.     move.l  tos,d5           64k counter
  234.     swap    d5
  235.     moveq.l #0,tos           init sum to 0
  236.     moveq.l #0,d1
  237.     move.l  #$-80000000,d3   init max to -inf
  238.     move.l  #$3FFFFFFF,d2    init min to +inf
  239.     bra.s   4$
  240. 1$: move.l  (a0)+,d4
  241.     cmp.l   d4,d3
  242.     bge.s   2$
  243.     move.l  d4,d3            max
  244. 2$: cmp.l   d4,d2
  245.     ble.s   3$
  246.     move.l  d4,d2            min
  247. 3$: add.l   d4,d1            sum
  248.     tst.l   d4
  249.     smi.b   d4
  250.     ext.w   d4
  251.     ext.l   d4
  252.     addx.l  d4,tos
  253. 4$: dbra.w  d0,1$
  254.     dbra.w  d5,1$
  255.     movem.l  d1-d3,-(dsp)
  256. \    move.l  d3,-(dsp)
  257. \    move.l  d1,-(dsp)
  258.     movem.l (rp)+,d4-d5
  259. end-code
  260.  
  261.  
  262. asm QUICK.REVERSAL    ( array-data  reversal-map-of-swap-pairs -- )
  263.     move.l  a3,-(rp)           save regs on return stack
  264.     move.l  a5,-(rp)
  265.     lea     $0(org,tos.l),a0   r
  266.     move.l  (dsp)+,a1
  267.     adda.l  org,a1             a
  268.     move.l  (a0)+,tos          i
  269. 1$: move.l  (a0)+,d0           next j
  270.     lea     $0(a1,tos.l),a3    abs i
  271.     lea     $0(a1,d0.l),a5     abs j
  272.     move.l  (a3),d1            swap cmplx a[i] with a[j]
  273.     move.l  (a5),d2
  274.     move.l  d1,(a5)+
  275.     move.l  d2,(a3)+
  276.     move.l  (a3),d1
  277.     move.l  (a5),d2
  278.     move.l  d1,(a5)
  279.     move.l  d2,(a3)
  280.     move.l  (a0)+,tos          next i
  281.     bne.l   1$                 zero terminator in swap map
  282.     move.l  (dsp)+,tos         cache tos
  283.     move.l  (rp)+,a5           restore regs
  284.     move.l  (rp)+,a3
  285. end-code
  286.  
  287. \  inner.loop register usage
  288. \  i   a   n           ss  le le1                ui      ur
  289. \  0dr 1   2   3   4   5   6   7     0ar 1   2   3   4   5   6   7
  290. \  i  air aii ur-i ur+i ss ur le1    ai aip  an  ur      ui
  291. \                             high               le
  292.  
  293. asm INNER.LOOP    ( u le ss n a i le1 -- hi )
  294.     movem.l d4-d6/a2-a3/a5,-(rp)       save regs on return stack
  295.     movem.l (dsp)+,d0-d2/d5-d6/a3/a5
  296.     lea     $0(a4,d1.l),a0       a
  297.     lea     $0(a0,d2.l),a2       an
  298.     adda.l  d0,a0                ai
  299.     lea     $0(a0,tos.l),a1      aip
  300.     moveq.l #0,tos               hi - abs all output or'd in 7dr
  301.     move.l  a5,d3                ur
  302.     move.l  a5,d4
  303.     sub.l   a3,d3                ur-ui
  304.     add.l   a3,d4                ur+ui
  305.     subq.l  #4,d6
  306.     move.l  d6,a3                le
  307. 1$:
  308.     move.l  (a1),d1              a[ip]
  309.     move.l  $4(a1),d2
  310.     asr.l   d5,d1                scale-down
  311.     asr.l   d5,d2
  312.     move.l  d1,d0
  313.     sub.l   d2,d0                c-d
  314.     move.l  a5,d6                ur
  315.     muls    d6,d0                z
  316.     muls    d3,d2                fd
  317.     add.l   d0,d2                fd+z
  318.     muls    d4,d1                gc
  319.     sub.l   d0,d1                gc-z
  320.     asl.l   #2,d1                scale-down cmplx * result
  321.     swap    d1
  322.     ext.l   d1
  323.     asl.l   #2,d2
  324.     swap    d2
  325.     ext.l   d2
  326.     move.l  (a0),d0              a[i] real
  327.     asr.l   d5,d0                scale-down
  328.     move.l  d0,d6
  329.     sub.l   d2,d6                a[i]-t
  330.     move.l  d6,(a1)+             a[ip]
  331.     forth{ markhere bpl$+4 }
  332.     neg.l   d6
  333.     or.l    d6,tos
  334.     forth{ nonauto? backhere }
  335.     add.l   d2,d0                a[i]+t
  336.     move.l  d0,(a0)+             a[i]
  337.     forth{ markhere bpl$+4 }
  338.     neg.l   d0
  339.     or.l    d0,tos
  340.     forth{ nonauto? backhere }
  341.     move.l  (a0),d0              a[i] imag
  342.     asr.l   d5,d0                scale-down
  343.     move.l  d0,d6
  344.     sub.l   d1,d6                a[i]-t
  345.     move.l  d6,(a1)              a[ip]
  346.     forth{ markhere bpl$+4 }
  347.     neg.l   d6
  348.     or.l    d6,tos
  349.     forth{ nonauto? backhere }
  350.     add.l   d1,d0                a[i]+t
  351.     move.l  d0,(a0)              a[i]
  352.     forth{ markhere bpl$+4 }
  353.     neg.l   d0
  354.     or.l    d0,tos
  355.     forth{ nonauto? backhere }
  356.     adda.l  a3,a1
  357.     adda.l  a3,a0
  358.     cmp.l   a2,a0
  359.     blt.l   1$
  360.     movem.l (rp)+,d4-d6/a2-a3/a5
  361. end-code
  362.