;/****************************************************************************; *; *  XVID MPEG-4 VIDEO CODEC; *  - SSE2 forward discrete cosine transform -; *; *  Copyright(C) 2003 Pascal Massimino <skal@planet-d.net>; *; *  This program is free software; you can redistribute it and/or modify it; *  under the terms of the GNU General Public License as published by; *  the Free Software Foundation; either version 2 of the License, or; *  (at your option) any later version.; *; *  This program is distributed in the hope that it will be useful,; *  but WITHOUT ANY WARRANTY; without even the implied warranty of; *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the; *  GNU General Public License for more details.; *; *  You should have received a copy of the GNU General Public License; *  along with this program; if not, write to the Free Software; *  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA; *; * $Id$; *; ***************************************************************************/BITS 32%macro cglobal 1  %ifdef PREFIX    global _%1    %define %1 _%1  %else    global %1  %endif%endmacro;-----------------------------------------------------------------------------;;                          -=FDCT=-;; Vertical pass is an implementation of the scheme:;  Loeffler C., Ligtenberg A., and Moschytz C.S.:;  Practical Fast 1D DCT Algorithm with Eleven Multiplications,;  Proc. ICASSP 1989, 988-991.;; Horizontal pass is a double 4x4 vector/matrix multiplication,; (see also Intel's Application Note 922:;  http://developer.intel.com/vtune/cbts/strmsimd/922down.htm;  Copyright (C) 1999 Intel Corporation);  ; Notes:;  * tan(3pi/16) is greater than 0.5, and would use the;    sign bit when turned into 16b fixed-point precision. So,;    we use the trick: x*tan3 = x*(tan3-1)+x; ;  * There's only one SSE-specific instruction (pshufw).;;  * There's still 1 or 2 ticks to save in fLLM_PASS, but;    I prefer having a readable code, instead of a tightly;    scheduled one...;;  * Quantization stage (as well as pre-transposition for the;    idct way back) can be included in the fTab* constants;    (with induced loss of precision, somehow);;  * Some more details at: http://skal.planet-d.net/coding/dct.html;;;//////////////////////////////////////////////////////////////////////;;  == Mean square errors ==;   0.000 0.001 0.001 0.002 0.000 0.002 0.001 0.000    [0.001];   0.035 0.029 0.032 0.032 0.031 0.032 0.034 0.035    [0.032];   0.026 0.028 0.027 0.027 0.025 0.028 0.028 0.025    [0.027];   0.037 0.032 0.031 0.030 0.028 0.029 0.026 0.031    [0.030];   0.000 0.001 0.001 0.002 0.000 0.002 0.001 0.001    [0.001];   0.025 0.024 0.022 0.022 0.022 0.022 0.023 0.023    [0.023];   0.026 0.028 0.025 0.028 0.030 0.025 0.026 0.027    [0.027];   0.021 0.020 0.020 0.022 0.020 0.022 0.017 0.019    [0.020];  ;  == Abs Mean errors ==;   0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000    [0.000];   0.020 0.001 0.003 0.003 0.000 0.004 0.002 0.003    [0.002];   0.000 0.001 0.001 0.001 0.001 0.004 0.000 0.000    [0.000];   0.027 0.001 0.000 0.002 0.002 0.002 0.001 0.000    [0.003];   0.000 0.000 0.000 0.000 0.000 0.001 0.000 0.001    [-0.000];   0.001 0.003 0.001 0.001 0.002 0.001 0.000 0.000    [-0.000];   0.000 0.002 0.002 0.001 0.001 0.002 0.001 0.000    [-0.000];   0.000 0.002 0.001 0.002 0.001 0.002 0.001 0.001    [-0.000];;  =========================;  Peak error:   1.0000;  Peak MSE:     0.0365;  Overall MSE:  0.0201;  Peak ME:      0.0265;  Overall ME:   0.0006;;-----------------------------------------------------------------------------;;                          -=IDCT=-;; A little slower than fdct, because the final stages (butterflies and; descaling) require some unpairable shifting and packing, all on; the same CPU unit.;;   THIS IDCT IS NOT IEEE-COMPLIANT: IT WILL FAIL THE [-300,300];   INPUT RANGE TEST (because of overflow). But the [-256,255] one;   is OK, and I'm fine with it (for now;);;  == Mean square errors ==;   0.007 0.006 0.005 0.007 0.006 0.007 0.005 0.007    [0.006];   0.006 0.008 0.007 0.007 0.007 0.008 0.008 0.008    [0.007];   0.008 0.008 0.008 0.008 0.007 0.009 0.010 0.007    [0.008];   0.007 0.007 0.006 0.007 0.008 0.007 0.006 0.008    [0.007];   0.007 0.006 0.006 0.006 0.006 0.005 0.006 0.006    [0.006];   0.008 0.007 0.006 0.008 0.007 0.008 0.009 0.009    [0.008];   0.008 0.006 0.010 0.008 0.008 0.008 0.007 0.007    [0.008];   0.007 0.006 0.006 0.007 0.007 0.006 0.006 0.007    [0.006];  ;  == Abs Mean errors ==;   0.001 0.000 0.000 0.001 0.001 0.000 0.000 0.000    [0.000];   0.000 0.002 0.002 0.000 0.001 0.001 0.000 0.002    [0.000];   0.001 0.002 0.001 0.001 0.001 0.001 0.000 0.001    [-0.001];   0.000 0.002 0.000 0.000 0.001 0.000 0.000 0.001    [-0.000];   0.000 0.001 0.001 0.001 0.000 0.001 0.000 0.001    [0.000];   0.000 0.001 0.001 0.001 0.001 0.000 0.001 0.000    [0.000];   0.001 0.001 0.002 0.001 0.001 0.002 0.001 0.001    [0.001];   0.000 0.000 0.001 0.000 0.000 0.000 0.000 0.000    [0.000];;  =========================;;  Peak error:   1.0000;  Peak MSE:     0.0096;  Overall MSE:  0.0070;  Peak ME:      0.0024;  Overall ME:   0.0001;;-----------------------------------------------------------------------------;=============================================================================; Read only data;=============================================================================%ifdef FORMAT_COFFSECTION .rodata data%elseSECTION .rodata data align=16%endifALIGN 16tan1:    times 8 dw 0x32ec    ; tan( pi/16)tan2:    times 8 dw 0x6a0a    ; tan(2pi/16)  (=sqrt(2)-1)tan3:    times 8 dw 0xab0e    ; tan(3pi/16)-1sqrt2:   times 8 dw 0x5a82    ; 0.5/sqrt(2);-----------------------------------------------------------------------------; Inverse DCT tables;-----------------------------------------------------------------------------ALIGN 16iTab1:  dw 0x4000, 0x539f, 0x4000, 0x22a3  dw 0x4000, 0xdd5d, 0x4000, 0xac61  dw 0x4000, 0x22a3, 0xc000, 0xac61  dw 0xc000, 0x539f, 0x4000, 0xdd5d  dw 0x58c5, 0x4b42, 0x4b42, 0xee58  dw 0x3249, 0xa73b, 0x11a8, 0xcdb7  dw 0x3249, 0x11a8, 0xa73b, 0xcdb7  dw 0x11a8, 0x4b42, 0x4b42, 0xa73biTab2:  dw 0x58c5, 0x73fc, 0x58c5, 0x300b  dw 0x58c5, 0xcff5, 0x58c5, 0x8c04  dw 0x58c5, 0x300b, 0xa73b, 0x8c04  dw 0xa73b, 0x73fc, 0x58c5, 0xcff5  dw 0x7b21, 0x6862, 0x6862, 0xe782  dw 0x45bf, 0x84df, 0x187e, 0xba41  dw 0x45bf, 0x187e, 0x84df, 0xba41  dw 0x187e, 0x6862, 0x6862, 0x84dfiTab3:  dw 0x539f, 0x6d41, 0x539f, 0x2d41  dw 0x539f, 0xd2bf, 0x539f, 0x92bf  dw 0x539f, 0x2d41, 0xac61, 0x92bf  dw 0xac61, 0x6d41, 0x539f, 0xd2bf  dw 0x73fc, 0x6254, 0x6254, 0xe8ee  dw 0x41b3, 0x8c04, 0x1712, 0xbe4d  dw 0x41b3, 0x1712, 0x8c04, 0xbe4d  dw 0x1712, 0x6254, 0x6254, 0x8c04iTab4:  dw 0x4b42, 0x6254, 0x4b42, 0x28ba  dw 0x4b42, 0xd746, 0x4b42, 0x9dac  dw 0x4b42, 0x28ba, 0xb4be, 0x9dac  dw 0xb4be, 0x6254, 0x4b42, 0xd746  dw 0x6862, 0x587e, 0x587e, 0xeb3d  dw 0x3b21, 0x979e, 0x14c3, 0xc4df  dw 0x3b21, 0x14c3, 0x979e, 0xc4df  dw 0x14c3, 0x587e, 0x587e, 0x979e  ; the original rounding trick is by  ; Michel Lespinasse (hi Walken!) <walken@zoy.org>ALIGN 16Idct_Rnd0: dd  65535, 65535, 65535, 65535Idct_Rnd1: dd   3612,  3612,  3612,  3612Idct_Rnd2: dd   2271,  2271,  2271,  2271Idct_Rnd3: dd   1203,  1203,  1203,  1203Idct_Rnd4: dd   1023,  1023,  1023,  1023Idct_Rnd5: dd    102,   102,   102,   102Idct_Rnd6: dd    398,   398,   398,   398Idct_Rnd7: dd    469,   469,   469,   469Idct_Sparse_Rnd0: times 4 dw  (65535>>11)Idct_Sparse_Rnd1: times 4 dw  ( 3612>>11)Idct_Sparse_Rnd2: times 4 dw  ( 2271>>11)  ; other rounders are zero...;-----------------------------------------------------------------------------; Forward DCT tables;-----------------------------------------------------------------------------ALIGN 16fTab1:  dw 0x4000, 0x4000, 0x58c5, 0x4b42,  dw 0xdd5d, 0xac61, 0xa73b, 0xcdb7,  dw 0x4000, 0x4000, 0x3249, 0x11a8,  dw 0x539f, 0x22a3, 0x4b42, 0xee58,  dw 0x4000, 0xc000, 0x3249, 0xa73b,  dw 0x539f, 0xdd5d, 0x4b42, 0xa73b,  dw 0xc000, 0x4000, 0x11a8, 0x4b42,  dw 0x22a3, 0xac61, 0x11a8, 0xcdb7fTab2:  dw 0x58c5, 0x58c5, 0x7b21, 0x6862,  dw 0xcff5, 0x8c04, 0x84df, 0xba41,  dw 0x58c5, 0x58c5, 0x45bf, 0x187e,  dw 0x73fc, 0x300b, 0x6862, 0xe782,  dw 0x58c5, 0xa73b, 0x45bf, 0x84df,  dw 0x73fc, 0xcff5, 0x6862, 0x84df,  dw 0xa73b, 0x58c5, 0x187e, 0x6862,  dw 0x300b, 0x8c04, 0x187e, 0xba41fTab3:  dw 0x539f, 0x539f, 0x73fc, 0x6254,  dw 0xd2bf, 0x92bf, 0x8c04, 0xbe4d,  dw 0x539f, 0x539f, 0x41b3, 0x1712,  dw 0x6d41, 0x2d41, 0x6254, 0xe8ee,  dw 0x539f, 0xac61, 0x41b3, 0x8c04,  dw 0x6d41, 0xd2bf, 0x6254, 0x8c04,  dw 0xac61, 0x539f, 0x1712, 0x6254,  dw 0x2d41, 0x92bf, 0x1712, 0xbe4dfTab4:  dw 0x4b42, 0x4b42, 0x6862, 0x587e,  dw 0xd746, 0x9dac, 0x979e, 0xc4df,  dw 0x4b42, 0x4b42, 0x3b21, 0x14c3,  dw 0x6254, 0x28ba, 0x587e, 0xeb3d,  dw 0x4b42, 0xb4be, 0x3b21, 0x979e,  dw 0x6254, 0xd746, 0x587e, 0x979e,  dw 0xb4be, 0x4b42, 0x14c3, 0x587e,  dw 0x28ba, 0x9dac, 0x14c3, 0xc4dfALIGN 16Fdct_Rnd0: dw  6,8,8,8, 6,8,8,8Fdct_Rnd1: dw  8,8,8,8, 8,8,8,8Fdct_Rnd2: dw 10,8,8,8, 8,8,8,8Rounder1:  dw  1,1,1,1, 1,1,1,1;=============================================================================; Code;=============================================================================SECTION .textcglobal idct_sse2_skalcglobal idct_sse2_sparse_skalcglobal fdct_sse2_skal;-----------------------------------------------------------------------------; Helper macro iMTX_MULT;-----------------------------------------------------------------------------%macro iMTX_MULT 4   ; %1=src, %2 = Table to use, %3=rounder, %4=Shift  movdqa  xmm0, [ecx+%1*16]     ; xmm0 = [01234567]  pshuflw xmm0, xmm0, 11011000b ; [0213]  pshufhw xmm0, xmm0, 11011000b ; [02134657]  pshufd  xmm4, xmm0, 00000000b ; [02020202]  pshufd  xmm5, xmm0, 10101010b ; [46464646]  pshufd  xmm6, xmm0, 01010101b ; [13131313]  pshufd  xmm7, xmm0, 11111111b ; [57575757]  pmaddwd xmm4, [%2+ 0]   ; dot [M00,M01][M04,M05][M08,M09][M12,M13]  pmaddwd xmm5, [%2+16]   ; dot [M02,M03][M06,M07][M10,M11][M14,M15]  pmaddwd xmm6, [%2+32]   ; dot [M16,M17][M20,M21][M24,M25][M28,M29]  pmaddwd xmm7, [%2+48]   ; dot [M18,M19][M22,M23][M26,M27][M30,M31]  paddd   xmm4, [%3]      ; Round  paddd   xmm6, xmm7      ; [b0|b1|b2|b3]  paddd   xmm4, xmm5      ; [a0|a1|a2|a3]  movdqa  xmm7, xmm6  paddd   xmm6, xmm4      ; mm6=a+b  psubd   xmm4, xmm7      ; mm4=a-b  psrad   xmm6, %4        ; => out [0123]  psrad   xmm4, %4        ; => out [7654]  packssdw xmm6, xmm4     ; [01237654]  pshufhw xmm6, xmm6, 00011011b ; [01234567]  movdqa  [ecx+%1*16], xmm6%endmacro;-----------------------------------------------------------------------------; Helper macro iLLM_PASS;-----------------------------------------------------------------------------%macro iLLM_PASS 1  ; %1: src/dst  movdqa xmm0, [tan3]    ; t3-1  movdqa xmm3, [%1+16*3] ; x3  movdqa xmm1, xmm0       ; t3-1  movdqa xmm5, [%1+16*5] ; x5  movdqa xmm4, [tan1]    ; t1  movdqa xmm6, [%1+16*1] ; x1  movdqa xmm7, [%1+16*7] ; x7  movdqa xmm2, xmm4       ; t1  pmulhw xmm0, xmm3       ; x3*(t3-1)  pmulhw xmm1, xmm5       ; x5*(t3-1)  paddsw xmm0, xmm3       ; x3*t3  paddsw xmm1, xmm5       ; x5*t3  psubsw xmm0, xmm5       ; x3*t3-x5 = tm35  paddsw xmm1, xmm3       ; x3+x5*t3 = tp35  pmulhw xmm4, xmm7       ; x7*t1  pmulhw xmm2, xmm6       ; x1*t1  paddsw xmm4, xmm6       ; x1+t1*x7 = tp17  psubsw xmm2, xmm7       ; x1*t1-x7 = tm17