; back port from GOGO-no coda 2.24b by Takehiro TOMINAGA

; GOGO-no-coda

;	Copyright (C) 1999 shigeo

;	special thanks to URURI

%include "nasm.h"

	externdef costab_fft

	externdef sintab_fft

	segment_data

	align 32

D_1_41421	dd	1.41421356

D_1_0	dd	1.0

D_0_5	dd	0.5

D_0_25	dd	0.25

D_0_0005	dd	0.0005

D_0_0	dd	0.0

	segment_code

;void fht(float *fz, int n);

proc	fht_FPU

%$fz	arg	4

%$n	arg	4

%$k	local	4

%$f0	local	4

%$f1	local	4

%$f2	local	4

%$f3	local	4

%$g0	local	4

%$g1	local	4

%$g2	local	4

%$g3	local	4

%$s1	local	4

%$c1	local	4

%$s2	local	4

%$c2	local	4

%$t_s	local	4

%$t_c	local	4

	alloc

	pushd	ebp, ebx, esi, edi

fht_FPU_1st_part:

fht_FPU_2nd_part:

fht_FPU_3rd_part:

.do_init:

	mov	r3, 16		;k1*fsize = 4*fsize = k4

	mov	r4, 8		;kx = k1/2

	mov	r2, 48		;k3*fsize

	mov	dword [sp(%$k)], 2	;k = 2

	mov	r0, [sp(%$fz)]	;fi

	lea	r1, [r0+8]		;gi = fi + kx

.do:

.do2:

	;f

	fld	dword [r0]

	fsub	dword [r0+r3]

	fld	dword [r0]

	fadd	dword [r0+r3]

	fld	dword [r0+r3*2]

	fsub	dword [r0+r2]

	fld	dword [r0+r3*2]

	fadd	dword [r0+r2]		;f2 f3 f0 f1

	fld	st2			;f0 f2 f3 f0 f1

	fadd	st0, st1

	fstp	dword [r0]		;fi[0]

	fld	st3			;f1 f2 f3 f0 f1

	fadd	st0, st2

	fstp	dword [r0+r3]		;fi[k1]

	fsubr	st0, st2		;f0-f2 f3 f0 f1

	fstp	dword [r0+r3*2]		;fi[k2]

	fsubr	st0, st2		;f1-f3 f0 f1

	fstp	dword [r0+r2]		;fi[k3]

	fcompp

	;g

	fld	dword [r1]

	fsub	dword [r1+r3]

	fld	dword [r1]

	fadd	dword [r1+r3]

	fld	dword [D_1_41421]

	fmul	dword [r1+r2]

	fld	dword [D_1_41421]

	fmul	dword [r1+r3*2]		;g2 g3 g0 g1

	fld	st2			;g0 g2 g3 g0 g1

	fadd	st0, st1

	fstp	dword [r1]		;gi[0]

	fld	st3			;g1 g2 g3 g0 g1

	fadd	st0, st2

	fstp	dword [r1+r3]		;gi[k1]

	fsubr	st0, st2		;g0-g2 g3 g0 g1

	fstp	dword [r1+r3*2]		;gi[k2]

	fsubr	st0, st2		;g1-g3 g0 g1

	fstp	dword [r1+r2]		;gi[k3]

	fcompp

	lea	r0, [r0+r3*4]

	lea	r1, [r1+r3*4]

	cmp	r0, r6

	jb	.do2

	mov	r0, [sp(%$k)]

	fld	dword [costab_fft +r0*4]

	fstp	dword [sp(%$t_c)]

	fld	dword [sintab_fft +r0*4]

	fstp	dword [sp(%$t_s)]

	fld	dword [D_1_0]

	fstp	dword [sp(%$c1)]

	fld	dword [D_0_0]

	fstp	dword [sp(%$s1)]

.for_init:

	mov	r5, 4		;i = 1*fsize

.for:

	fld	dword [sp(%$c1)]

	fmul	dword [sp(%$t_c)]

	fld	dword [sp(%$s1)]

	fmul	dword [sp(%$t_s)]

	fsubp	st1, st0		;c1

	fld	dword [sp(%$c1)]

	fmul	dword [sp(%$t_s)]

	fld	dword [sp(%$s1)]

	fmul	dword [sp(%$t_c)]

	faddp	st1, st0		;s1 c1

	fld	st1

	fmul	st0, st0		;c1c1 s1 c1

	fld	st1

	fmul	st0, st0		;s1s1 c1c1 s1 c1

	fsubp	st1, st0		;c2 s1 c1

	fstp	dword [sp(%$c2)]	;s1 c1

	fld	st1			;c1 s1 c1

	fmul	st0, st1		;c1s1 s1 c1

	fadd	st0, st0		;s2 s1 c1

	fstp	dword [sp(%$s2)]	;s1 c1

	fstp	dword [sp(%$s1)]	;c1

	fstp	dword [sp(%$c1)]	;

	mov	r0, [sp(%$fz)]

	add	r0, r5		;r0 = fi

	mov	r1, [sp(%$fz)]

	add	r1, r3

	sub	r1, r5		;r1 = gi

.do3:

	fld	dword [sp(%$s2)]

	fmul	dword [r0+r3]

	fld	dword [sp(%$c2)]

	fmul	dword [r1+r3]

	fsubp	st1, st0		;b = s2*fi[k1] - c2*gi[k1]

	fld	dword [sp(%$c2)]

	fmul	dword [r0+r3]

	fld	dword [sp(%$s2)]

	fmul	dword [r1+r3]

	faddp	st1, st0		;a = c2*fi[k1] + s2*gi[k1]  b

	fld	dword [r0]

	fsub	st0, st1		;f1 a b

	fstp	dword [sp(%$f1)]	;a b

	fadd	dword [r0]		;f0 b

	fstp	dword [sp(%$f0)]	;b

	fld	dword [r1]

	fsub	st0, st1		;g1 b

	fstp	dword [sp(%$g1)]	;b

	fadd	dword [r1]		;g0

	fstp	dword [sp(%$g0)]	;

	fld	dword [sp(%$s2)]

	fmul	dword [r0+r2]

	fld	dword [sp(%$c2)]

	fmul	dword [r1+r2]

	fsubp	st1, st0		;b = s2*fi[k3] - c2*gi[k3]

	fld	dword [sp(%$c2)]

	fmul	dword [r0+r2]

	fld	dword [sp(%$s2)]

	fmul	dword [r1+r2]

	faddp	st1, st0		;a = c2*fi[k3] + s2*gi[k3]  b

	fld	dword [r0+r3*2]

	fsub	st0, st1		;f3 a b

	fstp	dword [sp(%$f3)]	;a b

	fadd	dword [r0+r3*2]	;f2 b

	fstp	dword [sp(%$f2)]	;b

	fld	dword [r1+r3*2]

	fsub	st0, st1		;g3 b

	fstp	dword [sp(%$g3)]	;b

	fadd	dword [r1+r3*2]	;g2

	fstp	dword [sp(%$g2)]	;

	fld	dword [sp(%$s1)]

	fmul	dword [sp(%$f2)]

	fld	dword [sp(%$c1)]

	fmul	dword [sp(%$g3)]

	fsubp	st1, st0		;b = s1*f2 - c1*g3

	fld	dword [sp(%$c1)]

	fmul	dword [sp(%$f2)]

	fld	dword [sp(%$s1)]

	fmul	dword [sp(%$g3)]

	faddp	st1, st0		;a = c1*f2 + s1*g3  b

	fld	dword [sp(%$f0)]

	fsub	st0, st1		;fi[k2] a b

	fstp	dword [r0+r3*2]

	fadd	dword [sp(%$f0)]	;fi[0] b

	fstp	dword [r0]

	fld	dword [sp(%$g1)]

	fsub	st0, st1		;gi[k3] b

	fstp	dword [r1+r2]

	fadd	dword [sp(%$g1)]	;gi[k1]

	fstp	dword [r1+r3]

	fld	dword [sp(%$c1)]

	fmul	dword [sp(%$g2)]

	fld	dword [sp(%$s1)]

	fmul	dword [sp(%$f3)]

	fsubp	st1, st0		;b = c1*g2 - s1*f3

	fld	dword [sp(%$s1)]

	fmul	dword [sp(%$g2)]

	fld	dword [sp(%$c1)]

	fmul	dword [sp(%$f3)]

	faddp	st1, st0		;a = s1*g2 + c1*f3  b

	fld	dword [sp(%$g0)]

	fsub	st0, st1		;gi[k2] a b

	fstp	dword [r1+r3*2]

	fadd	dword [sp(%$g0)]	;gi[0] b

	fstp	dword [r1]

	fld	dword [sp(%$f1)]

	fsub	st0, st1		;fi[k3] b

	fstp	dword [r0+r2]

	fadd	dword [sp(%$f1)]	;fi[k1]

	fstp	dword [r0+r3]

	lea	r0, [r0+r3*4]

	lea	r1, [r1+r3*4]

	cmp	r0, r6

	jb near	.do3

	add	r5, 4

	cmp	r5, r4

	jb near	.for

	cmp	r3, [sp(%$n)]

	jae	.exit

	add	dword [sp(%$k)], 2	;k  += 2;

	lea	r3, [r3*4]		;k1 *= 4

	lea	r2, [r2*4]		;k3 *= 4

	lea	r4, [r4*4]		;kx *= 4

	mov	r0, [sp(%$fz)]	;fi

	lea	r1, [r0+r4]		;gi = fi + kx

	jmp	.do

.exit:

	popd	ebp, ebx, esi, edi

endproc

;*************************************************************

;void fht_FPU_FXCH(float *fz, int n);

proc	fht_FPU_FXCH

%$fz	arg	4

%$n	arg	4

%$k	local	4

%$f0	local	4

%$f1	local	4

%$f2	local	4

%$f3	local	4

%$g0	local	4

%$g1	local	4

%$g2	local	4

%$g3	local	4

%$s1	local	4

%$c1	local	4

%$s2	local	4

%$c2	local	4

%$t_s	local	4

%$t_c	local	4

	alloc

	pushd	ebp, ebx, esi, edi

fht_FPU_FXCH_1st_part:

fht_FPU_FXCH_2nd_part:

fht_FPU_FXCH_3rd_part:

.do_init:

	mov	r3, 16		;k1*fsize = 4*fsize = k4

	mov	r4, 8		;kx = k1/2

	mov	r2, 48		;k3*fsize

	mov	dword [sp(%$k)], 2	;k = 2

	mov	r0, [sp(%$fz)]	;fi

	lea	r1, [r0+8]		;gi = fi + kx

.do:

.do2:

	;f

	fld	dword [r0]

	fsub	dword [r0+r3]

	fld	dword [r0]

	fadd	dword [r0+r3]

	fld	dword [r0+r3*2]

	fsub	dword [r0+r2]

	fld	dword [r0+r3*2]

	fadd	dword [r0+r2]		;f2 f3 f0 f1

	fld	st3

	fld	st3

	fxch	st5

	fadd	st0, st3

	fxch	st4

	fadd	st0, st2

	fxch	st3

	fsubp	st1, st0

	fxch	st1

	fsubp	st4, st0

	fxch	st2

	fstp	dword [r0+r3]		;fi[k1]

	fstp	dword [r0]		;fi[0]

	fstp	dword [r0+r2]		;fi[k3]

	fstp	dword [r0+r3*2]		;fi[k2]

	;g

	fld	dword [r1]

	fsub	dword [r1+r3]

	fld	dword [r1]

	fadd	dword [r1+r3]

	fld	dword [D_1_41421]

	fmul	dword [r1+r2]

	fld	dword [D_1_41421]

	fmul	dword [r1+r3*2]		;g2 g3 g0 g1

	fld	st3

	fld	st3

	fxch	st5

	fadd	st0, st3

	fxch	st4

	fadd	st0, st2

	fxch	st3

	fsubp	st1, st0

	fxch	st1

	fsubp	st4, st0

	fxch	st2

	fstp	dword [r1+r3]		;gi[k1]

	fstp	dword [r1]		;gi[0]

	fstp	dword [r1+r2]		;gi[k3]

	fstp	dword [r1+r3*2]		;gi[k2]

	lea	r0, [r0+r3*4]

	lea	r1, [r1+r3*4]

	cmp	r0, r6

	jb	.do2

	mov	r0, [sp(%$k)]

	fld	dword [costab_fft +r0*4]

	fld	dword [sintab_fft +r0*4]

	fld	dword [D_1_0]

	fld	dword [D_0_0]

	fxch	st3

	fstp	dword [sp(%$t_c)]

	fxch	st1

	fstp	dword [sp(%$t_s)]

	fstp	dword [sp(%$c1)]

	fstp	dword [sp(%$s1)]

.for_init:

	mov	r5, 4		;i = 1*fsize

.for:

	fld	dword [sp(%$c1)]

	fmul	dword [sp(%$t_c)]

	fld	dword [sp(%$s1)]

	fmul	dword [sp(%$t_s)]

	fld	dword [sp(%$c1)]

	fmul	dword [sp(%$t_s)]

	fld	dword [sp(%$s1)]

	fmul	dword [sp(%$t_c)]

	fxch	st2

	fsubp	st3, st0		;c1

	faddp	st1, st0		;s1 c1

	fld	st1

	fxch	st2

	fmul	st0, st0		;c1c1 s1 c1

	fld	st1

	fxch	st2

	fmul	st0, st0		;s1s1 c1c1 s1 c1

	fxch	st3

	fst	dword [sp(%$c1)]	;c1

	fxch	st2

	fst	dword [sp(%$s1)]	;s1 c1c1 c1 s1s1

	fmulp	st2, st0

	fsubrp	st2, st0

	fadd	st0, st0		;s2 c2

	fxch	st1

	fstp	dword [sp(%$c2)]

	fstp	dword [sp(%$s2)]

	mov	r0, [sp(%$fz)]

	mov	r1, [sp(%$fz)]

	add	r0, r5		;r0 = fi

	add	r1, r3

	sub	r1, r5		;r1 = gi

.do3:

	fld	dword [sp(%$s2)]

	fmul	dword [r0+r3]

	fld	dword [sp(%$c2)]

	fmul	dword [r1+r3]

	fld	dword [sp(%$c2)]

	fmul	dword [r0+r3]

	fld	dword [sp(%$s2)]

	fmul	dword [r1+r3]

	fxch	st2

	fsubp	st3, st0		;b = s2*fi[k1] - c2*gi[k1]

	faddp	st1, st0		;a = c2*fi[k1] + s2*gi[k1]  b

	fld	dword [r1]

	fsub	st0, st2		;g1 a b

	fxch	st2

	fadd	dword [r1]		;g0 a g1

	fld	dword [r0]

	fsub	st0, st2		;f1 g0 a g1

	fxch	st2

	fadd	dword [r0]		;f0 g0 f1 g1

	fxch	st3

	fstp	dword [sp(%$g1)]

	fstp	dword [sp(%$g0)]

	fstp	dword [sp(%$f1)]

	fstp	dword [sp(%$f0)]

	fld	dword [sp(%$s2)]

	fmul	dword [r0+r2]

	fld	dword [sp(%$c2)]

	fmul	dword [r1+r2]

	fld	dword [sp(%$c2)]

	fmul	dword [r0+r2]

	fld	dword [sp(%$s2)]

	fmul	dword [r1+r2]

	fxch	st2

	fsubp	st3, st0		;b = s2*fi[k3] - c2*gi[k3]

	faddp	st1, st0		;a = c2*fi[k3] + s2*gi[k3]  b

	fld	dword [r1+r3*2]

	fsub	st0, st2		;g3 a b

	fxch	st2

	fadd	dword [r1+r3*2]	;g2 a g3

	fld	dword [r0+r3*2]

	fsub	st0, st2		;f3 g2 a g3

	fxch	st2

	fadd	dword [r0+r3*2]	;f2 g2 f3 g3

	fxch	st3

	fstp	dword [sp(%$g3)]

	fstp	dword [sp(%$g2)]

	fstp	dword [sp(%$f3)]

	fstp	dword [sp(%$f2)]

	fld	dword [sp(%$s1)]

	fmul	dword [sp(%$f2)]

	fld	dword [sp(%$c1)]

	fmul	dword [sp(%$g3)]

	fld	dword [sp(%$c1)]

	fmul	dword [sp(%$f2)]

	fld	dword [sp(%$s1)]

	fmul	dword [sp(%$g3)]

	fxch	st2

	fsubp	st3, st0		;b = s1*f2 - c1*g3

	faddp	st1, st0		;a = c1*f2 + s1*g3  b

	fld	dword [sp(%$g1)]

	fsub	st0, st2		;gi[k3] a b

	fxch	st2

	fadd	dword [sp(%$g1)]	;gi[k1] a gi[k3]

	fld	dword [sp(%$f0)]

	fsub	st0, st2		;fi[k2] gi[k1] a gi[k3]

	fxch	st2

	fadd	dword [sp(%$f0)]	;fi[0] gi[k1] fi[k2] gi[k3]

	fxch	st3

	fstp	dword [r1+r2]

	fstp	dword [r1+r3]

	fstp	dword [r0+r3*2]

	fstp	dword [r0]

	fld	dword [sp(%$c1)]

	fmul	dword [sp(%$g2)]

	fld	dword [sp(%$s1)]

	fmul	dword [sp(%$f3)]

	fld	dword [sp(%$s1)]

	fmul	dword [sp(%$g2)]

	fld	dword [sp(%$c1)]

	fmul	dword [sp(%$f3)]

	fxch	st2

	fsubp	st3, st0		;b = c1*g2 - s1*f3

	faddp	st1, st0		;a = s1*g2 + c1*f3  b

	fld	dword [sp(%$f1)]

	fsub	st0, st2		;fi[k3] a b

	fxch	st2

	fadd	dword [sp(%$f1)]	;fi[k1] a fi[k3]

	fld	dword [sp(%$g0)]

	fsub	st0, st2		;gi[k2] fi[k1] a fi[k3]

	fxch	st2

	fadd	dword [sp(%$g0)]	;gi[0] fi[k1] gi[k2] fi[k3]

	fxch	st3

	fstp	dword [r0+r2]

	fstp	dword [r0+r3]

	fstp	dword [r1+r3*2]

	fstp	dword [r1]

	lea	r0, [r0+r3*4]

	lea	r1, [r1+r3*4]

	cmp	r0, r6

	jb near	.do3

	add	r5, 4

	cmp	r5, r4

	jb near	.for

	cmp	r3, [sp(%$n)]

	jae	.exit

	add	dword [sp(%$k)], 2	;k  += 2;

	lea	r3, [r3*4]		;k1 *= 4

	lea	r2, [r2*4]		;k3 *= 4

	lea	r4, [r4*4]		;kx *= 4

	mov	r0, [sp(%$fz)]	;fi

	lea	r1, [r0+r4]		;gi = fi + kx

	jmp	.do

.exit:

	popd	ebp, ebx, esi, edi

endproc

	end