fdct_sse2_skal.asm

wince下的xvidcore开发库,可用于MP4等视频播放开发

 ;/****************************************************************************
 ; *
 ; *  XVID MPEG-4 VIDEO CODEC
 ; *  - SSE2 forward discrete cosine transform -
 ; *
 ; *  Copyright(C) 2003 Pascal Massimino <[e-mail]>
 ; *
 ; *  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: fdct_sse2_skal.asm,v 1.1 2005/07/21 09:08:25 klschoef Exp $
 ; *
 ; ***************************************************************************/
 
 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_COFF
 SECTION .rodata data
 %else
 SECTION .rodata data align=16
 %endif
 
 ALIGN 16
 tan1:    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)-1
 sqrt2:   times 8 dw 0x5a82    ; 0.5/sqrt(2)
 
 ;-----------------------------------------------------------------------------
 ; Inverse DCT tables
 ;-----------------------------------------------------------------------------
 
 ALIGN 16
 iTab1:
   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, 0xa73b
 
 iTab2:
   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, 0x84df
 
 iTab3:
   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, 0x8c04
 
 iTab4:
   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!) <[e-mail]>
 
 ALIGN 16
 Idct_Rnd0: dd  65535, 65535, 65535, 65535
 Idct_Rnd1: dd   3612,  3612,  3612,  3612
 Idct_Rnd2: dd   2271,  2271,  2271,  2271
 Idct_Rnd3: dd   1203,  1203,  1203,  1203
 Idct_Rnd4: dd   1023,  1023,  1023,  1023
 Idct_Rnd5: dd    102,   102,   102,   102
 Idct_Rnd6: dd    398,   398,   398,   398
 Idct_Rnd7: dd    469,   469,   469,   469
 
 Idct_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 16
 fTab1:
   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, 0xcdb7
 
 fTab2:
   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, 0xba41
 
 fTab3:
   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, 0xbe4d
 
 fTab4:
   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, 0xc4df
 
 
 ALIGN 16
 Fdct_Rnd0: dw  6,8,8,8, 6,8,8,8
 Fdct_Rnd1: dw  8,8,8,8, 8,8,8,8
 Fdct_Rnd2: dw 10,8,8,8, 8,8,8,8
 Rounder1:  dw  1,1,1,1, 1,1,1,1
 
 ;=============================================================================
 ; Code
 ;=============================================================================
 
 SECTION .text
 
 cglobal idct_sse2_skal
 cglobal idct_sse2_sparse_skal
 cglobal 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