.text
	.even
.globl _ASM_FFT_main_loop
_ASM_FFT_main_loop:
	link a6,#-8
		  movem.l  d2-d7/a2-a6,-(sp)
		  move.l   a2,a4          | a4 = cos
		  move.l   a3,a5          | a5 = sin
		  moveq.l  #1,d1          | d1 = d = 1
		  clr.l    d2
		  bset     d0,d2
		  asr.l    #1,d2          | d2 = li = 1<<( power-1 )

		  subq.w   #1,d0          | d0 = i
FFT_main_loop_1:
		  movem.l  a0-a1,sp@-
		  lea.l    (a0,d1.w*4),a2 | a2 = &x_real[ d ]
		  lea.l    (a1,d1.w*4),a3 | a3 = &x_imag[ d ]

		  move.w   d2,d3          |
		  subq.w   #1,d3          | d3 = j
FFT_main_loop_2:
		  move.w   d1,d4
		  subq.w   #1,d4          | d4 = l
FFT_main_loop_3:
		  fmove.s  a0@,fp0       | fp0 = *real_ptr
		  fmove.s  a1@,fp1       | fp1 = *imag_ptr
		  fmove.s  a2@,fp2       | fp2 = *t_real_ptr
		  fmove.s  a3@,fp3       | fp3 = *t_imag_ptr
		  fmove.x  fp2,fp4
		  fmove.x  fp3,fp5
		  fadd.x   fp0,fp2        | fp2 = *real_ptr + *t_real_ptr
		  fadd.x   fp1,fp3        | fp3 = *imag_ptr + *t_imag_ptr
		  fsub.x   fp4,fp0        | fp0 = *real_ptr - *t_real_ptr
		  fsub.x   fp5,fp1        | fp1 = *imag_ptr - *t_imag_ptr
		  fmove.s  fp2,a0@+      | fp2 -> *real_ptr++
		  fmove.s  fp3,a1@+      | fp3 -> *imag_ptr++
		  fmove.s  fp0,a2@+      | fp4 -> *t_real_ptr++
		  fmove.s  fp1,a3@+      | fp5 -> *t_imag_ptr++
		  dbra     d4,FFT_main_loop_3
		  lea.l    (a0,d1.w*4),a0 | real_ptr += d|
		  lea.l    (a1,d1.w*4),a1 | imag_ptr += d|
		  lea.l    (a2,d1.w*4),a2 | t_real_ptr += d|
		  lea.l    (a3,d1.w*4),a3 | t_imag_ptr += d|
		  dbra     d3,FFT_main_loop_2

		  movem.l  sp@+,a0-a1
		  movem.l  a0-a1,sp@-

		  asl.l    #1,d1          | d = d<<1|
		  asr.l    #1,d2          | li = li>>1|
		  beq      FFT_main_loop_9 | li = 0 -> end of FFT
		  lea.l    (a0,d1.w*4),a0 | a0 = real_ptr = &x_real[ d ]
		  lea.l    (a1,d1.w*4),a1 | a1 = imag_ptr = &x_imag[ d ]
		  move.w   d2,d3
		  subq.w   #1,d3          | d3 = j
FFT_main_loop_5:
		  move.l   a4,a2          | a2 = cos_ptr = cos
		  move.l   a5,a3          | a3 = sin_ptr = sin
		  move.w   d1,d4
		  subq.w   #1,d4          | d4 = l
FFT_main_loop_6:
		  fmove.s  a0@,fp0       | fp0 = *real_ptr
		  fmove.s  a1@,fp1       | fp1 = *imag_ptr
		  fmove.s  a2@,fp2       | fp2 = *cos_ptr
		  fmove.s  a3@,fp3       | fp3 = *sin_ptr
		  fmove.x  fp2,fp4
		  fmove.x  fp3,fp5
		  fmul.x   fp0,fp2        | fp2 = (*real_ptr)*(*cos_ptr) (#2)
		  fmul.x   fp1,fp3        | fp3 = (*imag_ptr)*(*sin_ptr) (#2)
		  fadd.x   fp2,fp3        | fp3 = real
		  fmul.x   fp1,fp4        | fp4 = (*imag_ptr)*(*cos_ptr) (#2)
		  fmul.x   fp0,fp5        | fp5 = (*real_ptr)*(*sin_ptr) (#2)
		  fsub.x   fp5,fp4        | fp4 = imag
		  fmove.s  fp3,a0@+      | *real_ptr++ = real
		  fmove.s  fp4,a1@+      | *imag_ptr++ = imag
		  lea.l    (a2,d2.w*4),a2 | cos_ptr += li
		  lea.l    (a3,d2.w*4),a3 | sin_ptr += li
		  dbra     d4,FFT_main_loop_6
		  lea.l    (a0,d1.w*4),a0 | real_ptr += d
		  lea.l    (a1,d1.w*4),a1 | imag_ptr += d
		  dbra     d3,FFT_main_loop_5
FFT_main_loop_9:
		  movem.l  sp@+,a0-a1

		  dbra     d0,FFT_main_loop_1

		  movem.l  sp@+,d2-d7/a2-a6
		  unlk a6
		  rts