; from a new GOGO-no-coda (1999/09)

;	Copyright (C) 1999 shigeo

;	special thanks to Keiichi SAKAI, URURI

; hacked and back-ported to LAME

;	 by Takehiro TOMINAGA Nov 2000

%include "nasm.h"

	globaldef fht_3DN

	segment_data

	align	16

costab	dd	0x80000000, 0

	dd	1.414213562,1.414213562

	dd	9.238795283293805e-01, 9.238795283293805e-01

	dd	3.826834424611044e-01, 3.826834424611044e-01

	dd	9.951847264044178e-01, 9.951847264044178e-01

	dd	9.801714304836734e-02, 9.801714304836734e-02

	dd	9.996988186794428e-01, 9.996988186794428e-01

	dd	2.454122920569705e-02, 2.454122920569705e-02

	dd	9.999811752815535e-01, 9.999811752815535e-01

	dd	6.135884819898878e-03, 6.135884819898878e-03

D_1_0_0_0	dd	0.0		, 1.0

	segment_code

PIC_OFFSETTABLE

;void fht_3DN(float *fz, int nn);

proc	fht_3DN

	pushd	ebp, ebx, esi, edi

	sub	esp, 20

	call	get_pc.bp

	add	ebp, PIC_BASE()

	mov	r0, [esp+40]		;fi

	mov	r1, [esp+44]		;r1 = nn

	lea	r3, [PIC_EBP_REL(costab)]		;tri = costab

	lea	r4, [r0+r1*8]		;r4 = fn = &fz[n]

	mov	[esp+16], r4

	mov	r4, 8			;kx = k1/2

	pmov	mm7, [r3]

	loopalign 16

.do1

	lea	r3, [r3+16]	;tri += 2;

	pmov	mm6, [PIC_EBP_REL(costab+8)]

	lea	r2, [r4+r4*2]		;k3*fsize/2

	mov	r5, 4		;i = 1*fsize

	loopalign 16

.do2:

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

	;f

	pmov	mm0, [r0]	;fi0

	pmov	mm1, [r0+r4*2]	;fi1

	pmov	mm2, [r0+r2*2]	;fi3

	pmov	mm3, [r0+r4*4]	;fi2

	pupldq	mm0, mm0	;fi0 | fi0

	pupldq	mm1, mm1	;fi1 | fi1

	pupldq	mm2, mm2	;fi2 | fi2

	pupldq	mm3, mm3	;fi3 | fi3

	pxor	mm1, mm7	;fi1 | -fi1

	pxor	mm3, mm7	;fi3 | -fi3

	pfsub	mm0, mm1	;f1 | f0

	pfsub	mm2, mm3	;f3 | f2

	pmov	mm4, mm0

	pfadd	mm0, mm2	;f1+f3|f0+f2 = fi1 | fi0

	pfsub	mm4, mm2	;f1-f3|f0-f2 = fi3 | fi2

	pmovd	[r0], mm0	;fi[0]

	puphdq	mm0, mm0

	pmovd	[r0+r4*4], mm4	;fi[k2]

	puphdq	mm4, mm4

	pmovd	[r0+r4*2], mm4	;fi[k1]

	pmovd	[r0+r2*2], mm0	;fi[k3]

	lea	r0, [r0+r4*8]

	;g

	pmov	mm0, [r1]	;gi0

	pmov	mm1, [r1+r4*2]	;gi1

	pmov	mm2, [r1+r4*4]	;gi2

	pmov	mm3, [r1+r2*2]	;gi3

	pupldq	mm1, mm1

	pupldq	mm0, mm0	;gi0 | gi0

	pupldq	mm2, mm3	;gi3 | gi2

	pxor	mm1, mm7	;gi1 | -gi1

	pfsub	mm0, mm1	;gi0-gi1|gi0+gi1 = g1 | g0

	pfmul	mm2, mm6	;gi3*SQRT2|gi2*SQRT2 = g3 | g2

	pmov	mm4, mm0

	pfadd	mm0, mm2	;g1+g3|g0+g2 = gi1 | gi0

	pfsub	mm4, mm2	;g1-g3|g0-g2 = gi3 | gi2

	pmovd	[r1], mm0	;gi[0]

	puphdq	mm0, mm0

	pmovd	[r1+r4*4], mm4	;gi[k2]

	puphdq	mm4, mm4

	cmp	r0, [esp + 16]

	pmovd	[r1+r4*2], mm0	;gi[k1]

	pmovd	[r1+r2*2], mm4	;gi[k3]

	jb near .do2

	pmov	mm6, [r3+r5]	; this is not aligned address!!

	loopalign 16

.for:

;

; mm6 = c1 | s1

; mm7 = 0x800000000 | 0

;

	pmov	mm1, mm6

	mov	r0, [esp+40]	; fz

	puphdq	mm1, mm1	; c1 | c1

	lea	r1, [r0+r4*2]

	pfadd	mm1, mm1	; c1+c1 | c1+c1

	pfmul	mm1, mm6	; 2*c1*c1 | 2*c1*s1

	pfsub	mm1, [PIC_EBP_REL(D_1_0_0_0)] ; 2*c1*c1-1.0 | 2*c1*s1 = -c2 | s2

	pmov	mm0, mm1

	pxor	mm7, mm6	; c1 | -s1

	pupldq	mm2, mm0

	pupldq	mm3, mm6	; ** | c1

	puphdq	mm0, mm2	; s2 | c2

	puphdq	mm6, mm3	;-s1 | c1

	pxor	mm0, [PIC_EBP_REL(costab)]	; c2 | -s2

; mm0 =  s2| c2

; mm1 = -c2| s2

; mm6 =  c1| s1

; mm7 =  s1|-c1 (we use the opposite sign. from GOGO here)

	pmov	[esp], mm0

	pmov	[esp+8], mm1

	sub	r1, r5		;r1 = gi

	add	r0, r5		;r0 = fi

	loopalign 16

.do3:

	pmov	mm2, [r0+r4*2] ; fi[k1]

	pmov	mm4, [r1+r4*2] ; gi[k1]

	pmov	mm3, [r0+r2*2] ; fi[k3]

	pmov	mm5, [r1+r2*2] ; gi[k3]

	pupldq	mm2, mm2	; fi1 | fi1

	pupldq	mm4, mm4	; gi1 | gi1

	pupldq	mm3, mm3	; fi3 | fi3

	pupldq	mm5, mm5	; gi3 | gi3

	pfmul	mm2, mm0	; s2 * fi1 | c2 * fi1

	pfmul	mm4, mm1	;-c2 * gi1 | s2 * gi1

	pfmul	mm3, mm0	; s2 * fi3 | c2 * fi3

	pfmul	mm5, mm1	;-c2 * gi3 | s2 * gi3

	pfadd	mm2, mm4		;b | a

	pfadd	mm3, mm5		;d | c

	pmov	mm0, [r0]

	pmov	mm4, [r1]

	pmov	mm1, [r0+r4*4]

	pmov	mm5, [r1+r4*4]

	pupldq	mm0, mm4		;gi0 | fi0

	pupldq	mm1, mm5		;gi2 | fi2

	pmov	mm4, mm2

	pmov	mm5, mm3

	pfadd	mm2, mm0		;g0 | f0

	pfadd	mm3, mm1		;g2 | f2

	pfsub	mm0, mm4		;g1 | f1

	pfsub	mm1, mm5		;g3 | f3

	pmov	mm4, mm3

	pmov	mm5, mm1

	pupldq	mm4, mm4		;f2 | f2

	puphdq	mm5, mm5		;g3 | g3

	puphdq	mm3, mm3		;g2 | g2

	pupldq	mm1, mm1		;f3 | f3

	pfmul	mm4, mm6		;f2 * c1 | f2 * s1

	pfmul	mm5, mm7		;g3 * s1 | g3 *-c1

	pfmul	mm3, mm6		;g2 * c1 | g2 * s1

	pfmul	mm1, mm7		;f3 * s1 | f3 *-c1

	pfadd	mm4, mm5		;a | b

	pfsub	mm3, mm1		;d | c

	pmov	mm5, mm2

	pmov	mm1, mm0

	pupldq	mm2, mm2		;f0 | f0

	pupldq	mm0, mm0		;f1 | f1

	puphdq	mm1, mm2		;f0 | g1

	puphdq	mm5, mm0		;f1 | g0

	pmov	mm2, mm4

	pmov	mm0, mm3

	pfadd	mm4, mm1		;fi0 | gi1

	pfadd	mm3, mm5		;fi1 | gi0

	pfsub	mm1, mm2		;fi2 | gi3

	pfsub	mm5, mm0		;fi3 | gi2

	pmovd	[r1+r4*2], mm4	;gi[k1]

	puphdq	mm4, mm4

	pmovd	[r1], mm3		;gi[0]

	puphdq	mm3, mm3

	pmovd	[r1+r2*2], mm1	;gi[k3]

	puphdq	mm1, mm1

	pmovd	[r1+r4*4], mm5	;gi[k2]

	puphdq	mm5, mm5

	pmovd	[r0], mm4	;fi[0]

	pmovd	[r0+r4*2], mm3	;fi[k1]

	pmovd	[r0+r4*4], mm1	;fi[k2]

	pmovd	[r0+r2*2], mm5	;fi[k3]

	lea	r0, [r0+r4*8]

	lea	r1, [r1+r4*8]

	cmp	r0, [esp + 16]

	pmov	mm0, [esp]

	pmov	mm1, [esp+8]

	jb near	.do3

	add	r5, 4

; mm6 =  c1| s1

; mm7 =  s1|-c1 (we use the opposite sign. from GOGO here)

	pfmul	mm6, [r3]	; c1*a | s1*a

	pfmul	mm7, [r3+8]	; s1*b |-c1*b

	cmp	r5, r4

	pfsub	mm6, mm7	; c1*a-s1*b | s1*a+c1*b

	pupldq	mm7,mm6

	puphdq	mm6,mm7

	pmov	mm7, [PIC_EBP_REL(costab)]

	jb near	.for

	mov	r0, [esp+40]	;fi

	cmp	r4, [esp+40+4]

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

	jb near	.do1

.exitttt

	femms

	add	esp,20

	popd	ebp, ebx, esi, edi

endproc

;void fht_E3DN(float *fz, int nn);

proc	fht_E3DN

	pushd	ebp, ebx, esi, edi

	sub	esp, 20

	call	get_pc.bp

	add	ebp, PIC_BASE()

	mov	r0, [esp+40]		;fi

	mov	r1, [esp+44]		;r1 = nn

	lea	r3, [PIC_EBP_REL(costab)]		;tri = costab

	lea	r4, [r0+r1*8]		;r4 = fn = &fz[n]

	mov	[esp+16], r4

	mov	r4, 8			;kx = k1/2

	pmov	mm7, [r3]

	loopalign 16

.do1

	lea	r3, [r3+16]	;tri += 2;

	pmov	mm6, [PIC_EBP_REL(costab+8)]

	lea	r2, [r4+r4*2]		;k3*fsize/2

	mov	r5, 4		;i = 1*fsize

	loopalign 16

.do2:

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

;f

	pmov	mm0, [r0]	; X  | fi0

	pmov	mm1, [r0+r4*4]	; X  | fi2

	pupldq	mm0, [r0+r4*2]	;fi1 | fi0

	pupldq	mm1, [r0+r2*2]	;fi3 | fi2

	pfpnacc	mm0, mm0	;fi0+fi1 | fi0-fi1 = f0|f1

	pfpnacc	mm1, mm1	;fi2+fi3 | fi2-fi3 = f2|f3

	pmov	mm2, mm0

	pfadd	mm0, mm1	;f0+f2|f1+f3 = fi0 | fi1

	pfsub	mm2, mm1	;f0-f2|f1-f3 = fi2 | fi3

	pmovd	[r0+r4*2], mm0	;fi[k1]

	pmovd	[r0+r2*2], mm2	;fi[k3]

	puphdq	mm0, mm0

	puphdq	mm2, mm2

	pmovd	[r0], mm0	;fi[0]

	pmovd	[r0+r4*4], mm2	;fi[k2]

	lea	r0, [r0+r4*8]

;g

	pmov	mm3, [r1]	;    gi0

	pmov	mm4, [r1+r2*2]	;    gi3

	pupldq	mm3, [r1+r4*2]	;gi1|gi0

	pupldq	mm4, [r1+r4*4]	;gi2|gi3

	pfpnacc	mm3, mm3	;gi0+gi1  |gi0-gi1   = f0|f1

	pfmul	mm4, mm6	;gi2*SQRT2|gi3*SQRT2 = f2|f3

	pmov	mm5, mm3

	pfadd	mm3, mm4	;f0+f2|f1+f3

	pfsub	mm5, mm4	;f0-f2|f1-f3

	cmp	r0, [esp + 16]

	pmovd	[r1+r4*2], mm3	;gi[k1]

	pmovd	[r1+r2*2], mm5	;gi[k3]

	puphdq	mm3, mm3

	puphdq	mm5, mm5

	pmovd	[r1], mm3	;gi[0]

	pmovd	[r1+r4*4], mm5	;gi[k2]

	jb near .do2

	pmov	mm6, [r3+r5]	; this is not aligned address!!

	loopalign 16

.for:

;

; mm6 = c1 | s1

; mm7 = 0x800000000 | 0

;

	pmov	mm5, mm6

	mov	r0, [esp+40]	; fz

	puphdq	mm5, mm5	; c1 | c1

	lea	r1, [r0+r4*2]

	pfadd	mm5, mm5	; c1+c1 | c1+c1

	pfmul	mm5, mm6	; 2*c1*c1 | 2*c1*s1

	pfsub	mm5, [PIC_EBP_REL(D_1_0_0_0)] ; 2*c1*c1-1.0 | 2*c1*s1 = -c2 | s2

	pswapd	mm4, mm5	; s2 |-c2

	pxor	mm4, mm7	; s2 | c2

	pxor	mm7, mm6	; c1 |-s1

	pswapd	mm6, mm6	; s1 | c1

; mm4 =  s2| c2

; mm5 = -c2| s2

; mm6 =  c1| s1

; mm7 =  s1|-c1 (we use the opposite sign. from GOGO here)

	pmov	[esp], mm4

	pmov	[esp+8], mm5

	sub	r1, r5		;r1 = gi

	add	r0, r5		;r0 = fi

	loopalign 16

.do3:

	pmov	mm0, [r0+r2*2] ; fi[k1]

	pmov	mm2, [r1+r2*2] ; gi[k1]

	pmov	mm1, [r0+r4*2] ; fi[k3]

	pmov	mm3, [r1+r4*2] ; gi[k3]

	pupldq	mm0, mm0

	pupldq	mm2, mm2

	pupldq	mm1, mm1

	pupldq	mm3, mm3

	pfmul	mm0, mm4

	pfmul	mm2, mm5

	pfmul	mm1, mm4

	pfmul	mm3, mm5

	pfadd	mm0, mm2		;d | c

	pfadd	mm1, mm3		;b | a

	pmov	mm2, [r0+r4*4]		;fi2

	pupldq	mm3, [r1+r4*4]		;gi2 | -

	pmov	mm4, [r0]		;fi0

	pupldq	mm5, [r1]		;gi0 | -

	pupldq	mm2, mm0		;c | fi2

	puphdq	mm3, mm0		;d | gi2

	pupldq	mm4, mm1		;a | fi0

	puphdq	mm5, mm1		;b | gi0

	pfpnacc	mm2, mm2		;f2 | f3

	pfpnacc	mm3, mm3		;g2 | g3

	pfpnacc	mm4, mm4		;f0 | f1

	pfpnacc	mm5, mm5		;g0 | g1

	pmov	mm0, mm2

	pmov	mm1, mm3

	pupldq	mm2, mm2		;f3 | f3

	pupldq	mm3, mm3		;g3 | g3

	puphdq	mm0, mm0		;f2 | f2

	puphdq	mm1, mm1		;g2 | g2

	pswapd	mm4, mm4		;f1 | f0

	pswapd	mm5, mm5		;g1 | g0

	pfmul	mm0, mm7		;f2 * s1 | f2 *-c1

	pfmul	mm3, mm6		;g3 * c1 | g3 * s1

	pfmul	mm1, mm6		;g2 * c1 | g2 * s1

	pfmul	mm2, mm7		;f3 * s1 | f3 *-c1

	pfsub	mm0, mm3		; b |-a

	pfsub	mm1, mm2		; d | c

	pmov	mm2, mm5

	pmov	mm3, mm4

	pupldq	mm4, mm0		;-a | f0

	pupldq	mm5, mm1		; c | g0

	puphdq	mm2, mm0		; b | g1

	puphdq	mm3, mm1		; d | f1

	pfpnacc	mm4, mm4		;fi2 | fi0

	pfpnacc	mm5, mm5		;gi0 | gi2

	pfpnacc	mm2, mm2		;gi1 | gi3

	pfpnacc	mm3, mm3		;fi1 | fi3

	pmovd	[r0], mm4		;fi[0]

	pmovd	[r1+r4*4], mm5		;gi[k2]

	pmovd	[r1+r2*2], mm2		;gi[k3]

	pmovd	[r0+r2*2], mm3		;fi[k3]

	puphdq	mm4, mm4

	puphdq	mm5, mm5

	puphdq	mm2, mm2

	puphdq	mm3, mm3

	pmovd	[r0+r4*4], mm4		;fi[k2]

	pmovd	[r1], mm5		;gi[0]

	pmovd	[r1+r4*2], mm2		;gi[k1]

	pmovd	[r0+r4*2], mm3		;fi[k1]

	lea	r0, [r0+r4*8]

	lea	r1, [r1+r4*8]

	cmp	r0, [esp + 16]

	pmov	mm4, [esp]

	pmov	mm5, [esp+8]

	jb near	.do3

	add	r5, 4

; mm6 =  c1| s1

; mm7 =  s1|-c1 (we use the opposite sign. from GOGO here)

	pfmul	mm6, [r3]	; c1*a | s1*a

	pfmul	mm7, [r3+8]	; s1*b |-c1*b

	cmp	r5, r4

	pfsub	mm6, mm7	; c1*a-s1*b | s1*a+c1*b

	pswapd	mm6, mm6 ; ???	; s1*a+c1*b | c1*a-s1*b

	pmov	mm7, [PIC_EBP_REL(costab)]

	jb near	.for

	mov	r0, [esp+40]	;fi

	cmp	r4, [esp+40+4]

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

	jb near	.do1

.exitttt

	femms

	add	esp,20

	popd	ebp, ebx, esi, edi

endproc