[svn] / branches / release-1_3-branch / xvidcore / src / utils / mbtransquant.c Repository:
ViewVC logotype

Diff of /branches/release-1_3-branch/xvidcore/src/utils/mbtransquant.c

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 605, Sat Oct 19 12:20:33 2002 UTC revision 1713, Mon Jul 10 08:09:59 2006 UTC
# Line 1  Line 1 
1  /*****************************************************************************  /*****************************************************************************
2   *   *
3   *  XVID MPEG-4 VIDEO CODEC   *  XVID MPEG-4 VIDEO CODEC
4   *  - MacroBlock transfer and quantization -   *  - MB Transfert/Quantization functions -
5   *   *
6   *  Copyright(C) 2002-2001 Christoph Lampert <gruel@web.de>   *  Copyright(C) 2001-2003  Peter Ross <pross@xvid.org>
7   *               2002-2001 Michael Militzer <isibaar@xvid.org>   *               2001-2003  Michael Militzer <isibaar@xvid.org>
8   *               2002-2001 Peter Ross <pross@xvid.org>   *               2003       Edouard Gomez <ed.gomez@free.fr>
  *               2002      Daniel Smith <danielsmith@astroboymail.com>  
  *  
  *  This program is an implementation of a part of one or more MPEG-4  
  *  Video tools as specified in ISO/IEC 14496-2 standard.  Those intending  
  *  to use this software module in hardware or software products are  
  *  advised that its use may infringe existing patents or copyrights, and  
  *  any such use would be at such party's own risk.  The original  
  *  developer of this software module and his/her company, and subsequent  
  *  editors and their companies, will have no liability for use of this  
  *  software or modifications or derivatives thereof.  
9   *   *
10   *  This program is free software; you can redistribute it and/or modify   *  This program is free software; you can redistribute it and/or modify
11   *  it under the terms of the GNU General Public License as published by   *  it under the terms of the GNU General Public License as published by
# Line 31  Line 21 
21   *  along with this program; if not, write to the Free Software   *  along with this program; if not, write to the Free Software
22   *  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA   *  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
23   *   *
24   * $Id: mbtransquant.c,v 1.18 2002-10-19 12:20:33 edgomez Exp $   * $Id: mbtransquant.c,v 1.32 2006-07-10 08:09:59 syskin Exp $
25   *   *
26   ****************************************************************************/   ****************************************************************************/
27    
28    #include <stdio.h>
29    #include <stdlib.h>
30  #include <string.h>  #include <string.h>
31    
32  #include "../portab.h"  #include "../portab.h"
# Line 43  Line 35 
35  #include "../global.h"  #include "../global.h"
36  #include "mem_transfer.h"  #include "mem_transfer.h"
37  #include "timer.h"  #include "timer.h"
38    #include "../bitstream/mbcoding.h"
39    #include "../bitstream/zigzag.h"
40  #include "../dct/fdct.h"  #include "../dct/fdct.h"
41  #include "../dct/idct.h"  #include "../dct/idct.h"
42  #include "../quant/quant_mpeg4.h"  #include "../quant/quant.h"
 #include "../quant/quant_h263.h"  
43  #include "../encoder.h"  #include "../encoder.h"
44    
45  #define MIN(X, Y) ((X)<(Y)?(X):(Y))  #include  "../quant/quant_matrix.h"
 #define MAX(X, Y) ((X)>(Y)?(X):(Y))  
46    
47  #define TOOSMALL_LIMIT 3                /* skip blocks having a coefficient sum below this value */  MBFIELDTEST_PTR MBFieldTest;
48    
49  /* this isnt pretty, but its better than 20 ifdefs */  /*
50     * Skip blocks having a coefficient sum below this value. This value will be
51  void   * corrected according to the MB quantizer to avoid artifacts for quant==1
52  MBTransQuantIntra(const MBParam * pParam,   */
53                                    FRAMEINFO * frame,  #define PVOP_TOOSMALL_LIMIT 1
54                                    MACROBLOCK * pMB,  #define BVOP_TOOSMALL_LIMIT 3
                                   const uint32_t x_pos,  
                                   const uint32_t y_pos,  
                                   int16_t data[6 * 64],  
                                   int16_t qcoeff[6 * 64])  
 {  
   
         uint32_t stride = pParam->edged_width;  
         uint32_t stride2 = stride / 2;  
         uint32_t next_block = stride * 8;  
         uint32_t i;  
         uint32_t iQuant = frame->quant;  
         uint8_t *pY_Cur, *pU_Cur, *pV_Cur;  
         IMAGE *pCurrent = &frame->image;  
   
         pY_Cur = pCurrent->y + (y_pos << 4) * stride + (x_pos << 4);  
         pU_Cur = pCurrent->u + (y_pos << 3) * stride2 + (x_pos << 3);  
         pV_Cur = pCurrent->v + (y_pos << 3) * stride2 + (x_pos << 3);  
   
         start_timer();  
         transfer_8to16copy(&data[0 * 64], pY_Cur, stride);  
         transfer_8to16copy(&data[1 * 64], pY_Cur + 8, stride);  
         transfer_8to16copy(&data[2 * 64], pY_Cur + next_block, stride);  
         transfer_8to16copy(&data[3 * 64], pY_Cur + next_block + 8, stride);  
         transfer_8to16copy(&data[4 * 64], pU_Cur, stride2);  
         transfer_8to16copy(&data[5 * 64], pV_Cur, stride2);  
         stop_transfer_timer();  
   
         start_timer();  
         pMB->field_dct = 0;  
         if ((frame->global_flags & XVID_INTERLACING) &&  
                 (x_pos>0) && (x_pos<pParam->mb_width-1) &&  
                 (y_pos>0) && (y_pos<pParam->mb_height-1)) {  
                 pMB->field_dct = MBDecideFieldDCT(data);  
         }  
         stop_interlacing_timer();  
   
         for (i = 0; i < 6; i++) {  
                 uint32_t iDcScaler = get_dc_scaler(iQuant, i < 4);  
   
                 start_timer();  
                 fdct(&data[i * 64]);  
                 stop_dct_timer();  
   
                 if (pParam->m_quant_type == H263_QUANT) {  
                         start_timer();  
                         quant_intra(&qcoeff[i * 64], &data[i * 64], iQuant, iDcScaler);  
                         stop_quant_timer();  
   
                         start_timer();  
                         dequant_intra(&data[i * 64], &qcoeff[i * 64], iQuant, iDcScaler);  
                         stop_iquant_timer();  
                 } else {  
                         start_timer();  
                         quant4_intra(&qcoeff[i * 64], &data[i * 64], iQuant, iDcScaler);  
                         stop_quant_timer();  
   
                         start_timer();  
                         dequant4_intra(&data[i * 64], &qcoeff[i * 64], iQuant, iDcScaler);  
                         stop_iquant_timer();  
                 }  
55    
56                  start_timer();  /*****************************************************************************
57                  idct(&data[i * 64]);   * Local functions
58                  stop_idct_timer();   ****************************************************************************/
         }  
59    
60          if (pMB->field_dct) {  /* permute block and return field dct choice */
61                  next_block = stride;  static __inline uint32_t
62                  stride *= 2;  MBDecideFieldDCT(int16_t data[6 * 64])
63          }  {
64            uint32_t field = MBFieldTest(data);
65    
66          start_timer();          if (field)
67          transfer_16to8copy(pY_Cur, &data[0 * 64], stride);                  MBFrameToField(data);
         transfer_16to8copy(pY_Cur + 8, &data[1 * 64], stride);  
         transfer_16to8copy(pY_Cur + next_block, &data[2 * 64], stride);  
         transfer_16to8copy(pY_Cur + next_block + 8, &data[3 * 64], stride);  
         transfer_16to8copy(pU_Cur, &data[4 * 64], stride2);  
         transfer_16to8copy(pV_Cur, &data[5 * 64], stride2);  
         stop_transfer_timer();  
68    
69            return field;
70  }  }
71    
72    /* Performs Forward DCT on all blocks */
73  uint8_t  static __inline void
74  MBTransQuantInter(const MBParam * pParam,  MBfDCT(const MBParam * const pParam,
75                                    FRAMEINFO * frame,             const FRAMEINFO * const frame,
76                                    MACROBLOCK * pMB,             MACROBLOCK * const pMB,
77                                    const uint32_t x_pos,             uint32_t x_pos,
78                                    const uint32_t y_pos,             uint32_t y_pos,
79                                    int16_t data[6 * 64],             int16_t data[6 * 64])
                                   int16_t qcoeff[6 * 64])  
80  {  {
81            /* Handles interlacing */
         uint32_t stride = pParam->edged_width;  
         uint32_t stride2 = stride / 2;  
         uint32_t next_block = stride * 8;  
         uint32_t i;  
         uint32_t iQuant = frame->quant;  
         uint8_t *pY_Cur, *pU_Cur, *pV_Cur;  
         uint8_t cbp = 0;  
         uint32_t sum;  
         IMAGE *pCurrent = &frame->image;  
   
         pY_Cur = pCurrent->y + (y_pos << 4) * stride + (x_pos << 4);  
         pU_Cur = pCurrent->u + (y_pos << 3) * stride2 + (x_pos << 3);  
         pV_Cur = pCurrent->v + (y_pos << 3) * stride2 + (x_pos << 3);  
   
82          start_timer();          start_timer();
83          pMB->field_dct = 0;          pMB->field_dct = 0;
84          if ((frame->global_flags & XVID_INTERLACING) &&          if ((frame->vol_flags & XVID_VOL_INTERLACING) &&
85                  (x_pos>0) && (x_pos<pParam->mb_width-1) &&                  (x_pos>0) && (x_pos<pParam->mb_width-1) &&
86                  (y_pos>0) && (y_pos<pParam->mb_height-1)) {                  (y_pos>0) && (y_pos<pParam->mb_height-1)) {
87                  pMB->field_dct = MBDecideFieldDCT(data);                  pMB->field_dct = MBDecideFieldDCT(data);
88          }          }
89          stop_interlacing_timer();          stop_interlacing_timer();
90    
91          for (i = 0; i < 6; i++) {          /* Perform DCT */
                 /*  
                  *  no need to transfer 8->16-bit  
                  * (this is performed already in motion compensation)  
                  */  
92                  start_timer();                  start_timer();
93                  fdct(&data[i * 64]);          fdct((short * const)&data[0 * 64]);
94            fdct((short * const)&data[1 * 64]);
95            fdct((short * const)&data[2 * 64]);
96            fdct((short * const)&data[3 * 64]);
97            fdct((short * const)&data[4 * 64]);
98            fdct((short * const)&data[5 * 64]);
99                  stop_dct_timer();                  stop_dct_timer();
   
                 if (pParam->m_quant_type == 0) {  
                         start_timer();  
                         sum = quant_inter(&qcoeff[i * 64], &data[i * 64], iQuant);  
                         stop_quant_timer();  
                 } else {  
                         start_timer();  
                         sum = quant4_inter(&qcoeff[i * 64], &data[i * 64], iQuant);  
                         stop_quant_timer();  
                 }  
   
                 if ((sum >= TOOSMALL_LIMIT) || (qcoeff[i*64] != 0) ||  
                         (qcoeff[i*64+1] != 0) || (qcoeff[i*64+8] != 0)) {  
   
                         if (pParam->m_quant_type == H263_QUANT) {  
                                 start_timer();  
                                 dequant_inter(&data[i * 64], &qcoeff[i * 64], iQuant);  
                                 stop_iquant_timer();  
                         } else {  
                                 start_timer();  
                                 dequant4_inter(&data[i * 64], &qcoeff[i * 64], iQuant);  
                                 stop_iquant_timer();  
100                          }                          }
101    
102                          cbp |= 1 << (5 - i);  /* Performs Inverse DCT on all blocks */
103    static __inline void
104                          start_timer();  MBiDCT(int16_t data[6 * 64],
105                          idct(&data[i * 64]);             const uint8_t cbp)
                         stop_idct_timer();  
                 }  
         }  
   
         if (pMB->field_dct) {  
                 next_block = stride;  
                 stride *= 2;  
         }  
   
         start_timer();  
         if (cbp & 32)  
                 transfer_16to8add(pY_Cur, &data[0 * 64], stride);  
         if (cbp & 16)  
                 transfer_16to8add(pY_Cur + 8, &data[1 * 64], stride);  
         if (cbp & 8)  
                 transfer_16to8add(pY_Cur + next_block, &data[2 * 64], stride);  
         if (cbp & 4)  
                 transfer_16to8add(pY_Cur + next_block + 8, &data[3 * 64], stride);  
         if (cbp & 2)  
                 transfer_16to8add(pU_Cur, &data[4 * 64], stride2);  
         if (cbp & 1)  
                 transfer_16to8add(pV_Cur, &data[5 * 64], stride2);  
         stop_transfer_timer();  
   
         return cbp;  
   
 }  
   
 void  
 MBTransQuantIntra2(const MBParam * pParam,  
                                   FRAMEINFO * frame,  
                                   MACROBLOCK * pMB,  
                                   const uint32_t x_pos,  
                                   const uint32_t y_pos,  
                                   int16_t data[6 * 64],  
                                   int16_t qcoeff[6 * 64])  
 {  
         MBTrans(pParam,frame,pMB,x_pos,y_pos,data);  
         MBfDCT(pParam,frame,pMB,data);  
         MBQuantIntra(pParam,frame,pMB,data,qcoeff);  
         MBDeQuantIntra(pParam,frame->quant,data,qcoeff);  
         MBiDCT(data,0x3F);  
         MBTransAdd(pParam,frame,pMB,x_pos,y_pos,data,0x3F);  
 }  
   
   
 uint8_t  
 MBTransQuantInter2(const MBParam * pParam,  
                                   FRAMEINFO * frame,  
                                   MACROBLOCK * pMB,  
                                   const uint32_t x_pos,  
                                   const uint32_t y_pos,  
                                   int16_t data[6 * 64],  
                                   int16_t qcoeff[6 * 64])  
 {  
         uint8_t cbp;  
   
 /* there is no MBTrans for Inter block, that's done in motion compensation already */  
   
         MBfDCT(pParam,frame,pMB,data);  
         cbp = MBQuantInter(pParam,frame->quant,data,qcoeff);  
         MBDeQuantInter(pParam,frame->quant,data,qcoeff,cbp);  
         MBiDCT(data,cbp);  
         MBTransAdd(pParam,frame,pMB,x_pos,y_pos,data,cbp);  
   
         return cbp;  
 }  
   
 uint8_t  
 MBTransQuantInterBVOP(const MBParam * pParam,  
                                   FRAMEINFO * frame,  
                                   MACROBLOCK * pMB,  
                                   int16_t data[6 * 64],  
                                   int16_t qcoeff[6 * 64])  
 {  
         uint8_t cbp;  
   
 /* there is no MBTrans for Inter block, that's done in motion compensation already */  
   
         MBfDCT(pParam,frame,pMB,data);  
         cbp = MBQuantInter(pParam,frame->quant,data,qcoeff);  
   
 /* we don't have to DeQuant, iDCT and Transfer back data for B-frames */  
   
         return cbp;  
 }  
   
   
 void  
 MBfDCT(const MBParam * pParam,  
                                   FRAMEINFO * frame,  
                                   MACROBLOCK * pMB,  
                                   int16_t data[6 * 64])  
106  {  {
         int i;  
   
107          start_timer();          start_timer();
108          pMB->field_dct = 0;          if(cbp & (1 << (5 - 0))) idct((short * const)&data[0 * 64]);
109          if ((frame->global_flags & XVID_INTERLACING)) {          if(cbp & (1 << (5 - 1))) idct((short * const)&data[1 * 64]);
110                  pMB->field_dct = MBDecideFieldDCT(data);          if(cbp & (1 << (5 - 2))) idct((short * const)&data[2 * 64]);
111          }          if(cbp & (1 << (5 - 3))) idct((short * const)&data[3 * 64]);
112          stop_interlacing_timer();          if(cbp & (1 << (5 - 4))) idct((short * const)&data[4 * 64]);
113            if(cbp & (1 << (5 - 5))) idct((short * const)&data[5 * 64]);
114          for (i = 0; i < 6; i++) {          stop_idct_timer();
                 start_timer();  
                 fdct(&data[i * 64]);  
                 stop_dct_timer();  
         }  
115  }  }
116    
117  void  /* Quantize all blocks -- Intra mode */
118  MBQuantDeQuantIntra(const MBParam * pParam,  static __inline void
119                                          FRAMEINFO * frame,  MBQuantIntra(const MBParam * pParam,
120                                          MACROBLOCK * pMB,                           const FRAMEINFO * const frame,
121                             const MACROBLOCK * pMB,
122                                          int16_t qcoeff[6 * 64],                                          int16_t qcoeff[6 * 64],
123                                          int16_t data[6*64])                                          int16_t data[6*64])
124  {  {
125          int i;          int scaler_lum, scaler_chr;
126          int iQuant = frame->quant;          quant_intraFuncPtr quant;
127    
128          start_timer();          /* check if quant matrices need to be re-initialized with new quant */
129          pMB->field_dct = 0;          if (pParam->vol_flags & XVID_VOL_MPEGQUANT) {
130          if ((frame->global_flags & XVID_INTERLACING)) {                  if (pParam->last_quant_initialized_intra != pMB->quant) {
131                  pMB->field_dct = MBDecideFieldDCT(data);                          init_intra_matrix(pParam->mpeg_quant_matrices, pMB->quant);
132          }          }
133          stop_interlacing_timer();                  quant = quant_mpeg_intra;
   
         for (i = 0; i < 6; i++) {  
                 uint32_t iDcScaler = get_dc_scaler(iQuant, i < 4);  
   
                 if (pParam->m_quant_type == H263_QUANT) {  
                         start_timer();  
                         quant_intra(&qcoeff[i * 64], &data[i * 64], iQuant, iDcScaler);  
                         stop_quant_timer();  
   
                         start_timer();  
                         dequant_intra(&data[i * 64], &qcoeff[i * 64], iQuant, iDcScaler);  
                         stop_iquant_timer();  
134                  } else {                  } else {
135                          start_timer();                  quant = quant_h263_intra;
                         quant4_intra(&qcoeff[i * 64], &data[i * 64], iQuant, iDcScaler);  
                         stop_quant_timer();  
   
                         start_timer();  
                         dequant4_intra(&data[i * 64], &qcoeff[i * 64], iQuant, iDcScaler);  
                         stop_iquant_timer();  
                 }  
         }  
 }  
   
 void  
 MBQuantIntra(const MBParam * pParam,  
                          FRAMEINFO * frame,  
                          MACROBLOCK *pMB,  
                          int16_t data[6 * 64],  
                      int16_t qcoeff[6 * 64])  
 {  
         int i;  
         int iQuant = frame->quant;  
   
         start_timer();  
         pMB->field_dct = 0;  
         if ((frame->global_flags & XVID_INTERLACING)) {  
                 pMB->field_dct = MBDecideFieldDCT(data);  
136          }          }
         stop_interlacing_timer();  
137    
138          for (i = 0; i < 6; i++) {          scaler_lum = get_dc_scaler(pMB->quant, 1);
139                  uint32_t iDcScaler = get_dc_scaler(iQuant, i < 4);          scaler_chr = get_dc_scaler(pMB->quant, 0);
140    
141                  if (pParam->m_quant_type == H263_QUANT) {          /* Quantize the block */
                         start_timer();  
                         quant_intra(&qcoeff[i * 64], &data[i * 64], iQuant, iDcScaler);  
                         stop_quant_timer();  
                 } else {  
142                          start_timer();                          start_timer();
143                          quant4_intra(&qcoeff[i * 64], &data[i * 64], iQuant, iDcScaler);          quant(&data[0 * 64], &qcoeff[0 * 64], pMB->quant, scaler_lum, pParam->mpeg_quant_matrices);
144            quant(&data[1 * 64], &qcoeff[1 * 64], pMB->quant, scaler_lum, pParam->mpeg_quant_matrices);
145            quant(&data[2 * 64], &qcoeff[2 * 64], pMB->quant, scaler_lum, pParam->mpeg_quant_matrices);
146            quant(&data[3 * 64], &qcoeff[3 * 64], pMB->quant, scaler_lum, pParam->mpeg_quant_matrices);
147            quant(&data[4 * 64], &qcoeff[4 * 64], pMB->quant, scaler_chr, pParam->mpeg_quant_matrices);
148            quant(&data[5 * 64], &qcoeff[5 * 64], pMB->quant, scaler_chr, pParam->mpeg_quant_matrices);
149                          stop_quant_timer();                          stop_quant_timer();
150                  }                  }
         }  
 }  
151    
152  void  /* DeQuantize all blocks -- Intra mode */
153    static __inline void
154  MBDeQuantIntra(const MBParam * pParam,  MBDeQuantIntra(const MBParam * pParam,
155                             const int iQuant,                             const int iQuant,
156                                    int16_t qcoeff[6 * 64],                                    int16_t qcoeff[6 * 64],
157                                    int16_t data[6*64])                                    int16_t data[6*64])
158  {  {
159          int i;          int mpeg;
160            int scaler_lum, scaler_chr;
161    
162          for (i = 0; i < 6; i++) {          quant_intraFuncPtr const dequant[2] =
163                  uint32_t iDcScaler = get_dc_scaler(iQuant, i < 4);                  {
164                            dequant_h263_intra,
165                            dequant_mpeg_intra
166                    };
167    
168            mpeg = !!(pParam->vol_flags & XVID_VOL_MPEGQUANT);
169            scaler_lum = get_dc_scaler(iQuant, 1);
170            scaler_chr = get_dc_scaler(iQuant, 0);
171    
                 if (pParam->m_quant_type == H263_QUANT) {  
172                          start_timer();                          start_timer();
173                          dequant_intra(&data[i * 64], &qcoeff[i * 64], iQuant, iDcScaler);          dequant[mpeg](&qcoeff[0 * 64], &data[0 * 64], iQuant, scaler_lum, pParam->mpeg_quant_matrices);
174            dequant[mpeg](&qcoeff[1 * 64], &data[1 * 64], iQuant, scaler_lum, pParam->mpeg_quant_matrices);
175            dequant[mpeg](&qcoeff[2 * 64], &data[2 * 64], iQuant, scaler_lum, pParam->mpeg_quant_matrices);
176            dequant[mpeg](&qcoeff[3 * 64], &data[3 * 64], iQuant, scaler_lum, pParam->mpeg_quant_matrices);
177            dequant[mpeg](&qcoeff[4 * 64], &data[4 * 64], iQuant, scaler_chr, pParam->mpeg_quant_matrices);
178            dequant[mpeg](&qcoeff[5 * 64], &data[5 * 64], iQuant, scaler_chr, pParam->mpeg_quant_matrices);
179                          stop_iquant_timer();                          stop_iquant_timer();
                 } else {  
                         start_timer();  
                         dequant4_intra(&data[i * 64], &qcoeff[i * 64], iQuant, iDcScaler);  
                         stop_iquant_timer();  
                 }  
         }  
180  }  }
181    
182  uint8_t  static int
183    dct_quantize_trellis_c(int16_t *const Out,
184                                               const int16_t *const In,
185                                               int Q,
186                                               const uint16_t * const Zigzag,
187                                               const uint16_t * const QuantMatrix,
188                                               int Non_Zero,
189                                               int Sum,
190                                               int Lambda_Mod);
191    
192    /* Quantize all blocks -- Inter mode */
193    static __inline uint8_t
194  MBQuantInter(const MBParam * pParam,  MBQuantInter(const MBParam * pParam,
195                           const int iQuant,                           const FRAMEINFO * const frame,
196                             const MACROBLOCK * pMB,
197                                    int16_t data[6 * 64],                                    int16_t data[6 * 64],
198                                    int16_t qcoeff[6 * 64])                           int16_t qcoeff[6 * 64],
199                             int bvop,
200                             int limit)
201  {  {
202    
203          int i;          int i;
204          uint8_t cbp = 0;          uint8_t cbp = 0;
205          int sum;          int sum;
206            int code_block, mpeg;
207    
208            quant_interFuncPtr const quant[2] =
209                    {
210                            quant_h263_inter,
211                            quant_mpeg_inter
212                    };
213    
214            mpeg = !!(pParam->vol_flags & XVID_VOL_MPEGQUANT);
215    
216          for (i = 0; i < 6; i++) {          for (i = 0; i < 6; i++) {
217    
218                  if (pParam->m_quant_type == 0) {                  /* Quantize the block */
219                          start_timer();                          start_timer();
220                          sum = quant_inter(&qcoeff[i * 64], &data[i * 64], iQuant);  
221                    sum = quant[mpeg](&qcoeff[i*64], &data[i*64], pMB->quant, pParam->mpeg_quant_matrices);
222    
223                    if(sum && (pMB->quant > 2) && (frame->vop_flags & XVID_VOP_TRELLISQUANT)) {
224                            const uint16_t *matrix;
225                            const static uint16_t h263matrix[] =
226                                    {
227                                            16, 16, 16, 16, 16, 16, 16, 16,
228                                            16, 16, 16, 16, 16, 16, 16, 16,
229                                            16, 16, 16, 16, 16, 16, 16, 16,
230                                            16, 16, 16, 16, 16, 16, 16, 16,
231                                            16, 16, 16, 16, 16, 16, 16, 16,
232                                            16, 16, 16, 16, 16, 16, 16, 16,
233                                            16, 16, 16, 16, 16, 16, 16, 16,
234                                            16, 16, 16, 16, 16, 16, 16, 16
235                                    };
236    
237                            matrix = (mpeg)?get_inter_matrix(pParam->mpeg_quant_matrices):h263matrix;
238                            sum = dct_quantize_trellis_c(&qcoeff[i*64], &data[i*64],
239                                                                                     pMB->quant, &scan_tables[0][0],
240                                                                                     matrix,
241                                                                                     63,
242                                                                                     sum,
243                                                                                     pMB->lambda[i]);
244                    }
245                          stop_quant_timer();                          stop_quant_timer();
246    
247                    /*
248                     * We code the block if the sum is higher than the limit and if the first
249                     * two AC coefficients in zig zag order are not zero.
250                     */
251                    code_block = 0;
252                    if ((sum >= limit) || (qcoeff[i*64+1] != 0) || (qcoeff[i*64+8] != 0)) {
253                            code_block = 1;
254                  } else {                  } else {
                         start_timer();  
                         sum = quant4_inter(&qcoeff[i * 64], &data[i * 64], iQuant);  
                         stop_quant_timer();  
                 }  
255    
256                  if (sum >= TOOSMALL_LIMIT) {    // skip block ?                          if (bvop && (pMB->mode == MODE_DIRECT || pMB->mode == MODE_DIRECT_NO4V)) {
257                          cbp |= 1 << (5 - i);                                  /* dark blocks prevention for direct mode */
258                                    if ((qcoeff[i*64] < -1) || (qcoeff[i*64] > 0))
259                                            code_block = 1;
260                            } else {
261                                    /* not direct mode */
262                                    if (qcoeff[i*64] != 0)
263                                            code_block = 1;
264                  }                  }
265          }          }
266          return cbp;  
267                    /* Set the corresponding cbp bit */
268                    cbp |= code_block << (5 - i);
269  }  }
270    
271  void          return(cbp);
272    }
273    
274    /* DeQuantize all blocks -- Inter mode */
275    static __inline void
276  MBDeQuantInter( const MBParam * pParam,  MBDeQuantInter( const MBParam * pParam,
277                                  const int iQuant,                                  const int iQuant,
278                                    int16_t data[6 * 64],                                    int16_t data[6 * 64],
279                                    int16_t qcoeff[6 * 64],                                    int16_t qcoeff[6 * 64],
280                                    const uint8_t cbp)                                    const uint8_t cbp)
281  {  {
282          int i;          int mpeg;
283    
284          for (i = 0; i < 6; i++) {          quant_interFuncPtr const dequant[2] =
                 if (cbp & (1 << (5 - i)))  
285                  {                  {
286                          if (pParam->m_quant_type == H263_QUANT) {                          dequant_h263_inter,
287                                  start_timer();                          dequant_mpeg_inter
288                                  dequant_inter(&data[i * 64], &qcoeff[i * 64], iQuant);                  };
289                                  stop_iquant_timer();  
290                          } else {          mpeg = !!(pParam->vol_flags & XVID_VOL_MPEGQUANT);
291    
292                                  start_timer();                                  start_timer();
293                                  dequant4_inter(&data[i * 64], &qcoeff[i * 64], iQuant);          if(cbp & (1 << (5 - 0))) dequant[mpeg](&data[0 * 64], &qcoeff[0 * 64], iQuant, pParam->mpeg_quant_matrices);
294            if(cbp & (1 << (5 - 1))) dequant[mpeg](&data[1 * 64], &qcoeff[1 * 64], iQuant, pParam->mpeg_quant_matrices);
295            if(cbp & (1 << (5 - 2))) dequant[mpeg](&data[2 * 64], &qcoeff[2 * 64], iQuant, pParam->mpeg_quant_matrices);
296            if(cbp & (1 << (5 - 3))) dequant[mpeg](&data[3 * 64], &qcoeff[3 * 64], iQuant, pParam->mpeg_quant_matrices);
297            if(cbp & (1 << (5 - 4))) dequant[mpeg](&data[4 * 64], &qcoeff[4 * 64], iQuant, pParam->mpeg_quant_matrices);
298            if(cbp & (1 << (5 - 5))) dequant[mpeg](&data[5 * 64], &qcoeff[5 * 64], iQuant, pParam->mpeg_quant_matrices);
299                                  stop_iquant_timer();                                  stop_iquant_timer();
300                          }                          }
                 }  
         }  
 }  
   
 void  
 MBiDCT( int16_t data[6 * 64],  
                 const uint8_t cbp)  
 {  
         int i;  
   
         for (i = 0; i < 6; i++) {  
                 if (cbp & (1 << (5 - i)))  
                 {  
                         start_timer();  
                         idct(&data[i * 64]);  
                         stop_idct_timer();  
301    
302                  }  typedef void (transfer_operation_8to16_t) (int16_t *Dst, const uint8_t *Src, int BpS);
303          }  typedef void (transfer_operation_16to8_t) (uint8_t *Dst, const int16_t *Src, int BpS);
 }  
304    
305    
306  void  static __inline void
307  MBTrans(const MBParam * pParam,  MBTrans8to16(const MBParam * const pParam,
308                                    FRAMEINFO * frame,                           const FRAMEINFO * const frame,
309                                    MACROBLOCK * pMB,                           const MACROBLOCK * const pMB,
310                                    const uint32_t x_pos,                                    const uint32_t x_pos,
311                                    const uint32_t y_pos,                                    const uint32_t y_pos,
312                                    int16_t data[6 * 64])                                    int16_t data[6 * 64])
# Line 502  Line 315 
315          uint32_t stride2 = stride / 2;          uint32_t stride2 = stride / 2;
316          uint32_t next_block = stride * 8;          uint32_t next_block = stride * 8;
317          uint8_t *pY_Cur, *pU_Cur, *pV_Cur;          uint8_t *pY_Cur, *pU_Cur, *pV_Cur;
318          IMAGE *pCurrent = &frame->image;          const IMAGE * const pCurrent = &frame->image;
319    
320            /* Image pointers */
321          pY_Cur = pCurrent->y + (y_pos << 4) * stride + (x_pos << 4);          pY_Cur = pCurrent->y + (y_pos << 4) * stride + (x_pos << 4);
322          pU_Cur = pCurrent->u + (y_pos << 3) * stride2 + (x_pos << 3);          pU_Cur = pCurrent->u + (y_pos << 3) * stride2 + (x_pos << 3);
323          pV_Cur = pCurrent->v + (y_pos << 3) * stride2 + (x_pos << 3);          pV_Cur = pCurrent->v + (y_pos << 3) * stride2 + (x_pos << 3);
324    
325            /* Do the transfer */
326          start_timer();          start_timer();
327          transfer_8to16copy(&data[0 * 64], pY_Cur, stride);          transfer_8to16copy(&data[0 * 64], pY_Cur, stride);
328          transfer_8to16copy(&data[1 * 64], pY_Cur + 8, stride);          transfer_8to16copy(&data[1 * 64], pY_Cur + 8, stride);
# Line 518  Line 333 
333          stop_transfer_timer();          stop_transfer_timer();
334  }  }
335    
336  void  static __inline void
337  MBTransAdd(const MBParam * pParam,  MBTrans16to8(const MBParam * const pParam,
338                                    FRAMEINFO * frame,                           const FRAMEINFO * const frame,
339                                    MACROBLOCK * pMB,                           const MACROBLOCK * const pMB,
340                                    const uint32_t x_pos,                                    const uint32_t x_pos,
341                                    const uint32_t y_pos,                                    const uint32_t y_pos,
342                                    int16_t data[6 * 64],                                    int16_t data[6 * 64],
343                             const uint32_t add, /* Must be 1 or 0 */
344                                    const uint8_t cbp)                                    const uint8_t cbp)
345  {  {
346          uint8_t *pY_Cur, *pU_Cur, *pV_Cur;          uint8_t *pY_Cur, *pU_Cur, *pV_Cur;
347          uint32_t stride = pParam->edged_width;          uint32_t stride = pParam->edged_width;
348          uint32_t stride2 = stride / 2;          uint32_t stride2 = stride / 2;
349          uint32_t next_block = stride * 8;          uint32_t next_block = stride * 8;
350          IMAGE *pCurrent = &frame->image;          const IMAGE * const pCurrent = &frame->image;
351    
352            /* Array of function pointers, indexed by [add] */
353            transfer_operation_16to8_t  * const functions[2] =
354                    {
355                            (transfer_operation_16to8_t*)transfer_16to8copy,
356                            (transfer_operation_16to8_t*)transfer_16to8add,
357                    };
358    
359            transfer_operation_16to8_t *transfer_op = NULL;
360    
361            /* Image pointers */
362          pY_Cur = pCurrent->y + (y_pos << 4) * stride + (x_pos << 4);          pY_Cur = pCurrent->y + (y_pos << 4) * stride + (x_pos << 4);
363          pU_Cur = pCurrent->u + (y_pos << 3) * stride2 + (x_pos << 3);          pU_Cur = pCurrent->u + (y_pos << 3) * stride2 + (x_pos << 3);
364          pV_Cur = pCurrent->v + (y_pos << 3) * stride2 + (x_pos << 3);          pV_Cur = pCurrent->v + (y_pos << 3) * stride2 + (x_pos << 3);
# Line 542  Line 368 
368                  stride *= 2;                  stride *= 2;
369          }          }
370    
371            /* Operation function */
372            transfer_op = functions[add];
373    
374            /* Do the operation */
375          start_timer();          start_timer();
376          if (cbp & 32)          if (cbp&32) transfer_op(pY_Cur,                                         &data[0 * 64], stride);
377                  transfer_16to8add(pY_Cur, &data[0 * 64], stride);          if (cbp&16) transfer_op(pY_Cur + 8,                                     &data[1 * 64], stride);
378          if (cbp & 16)          if (cbp& 8) transfer_op(pY_Cur + next_block,            &data[2 * 64], stride);
379                  transfer_16to8add(pY_Cur + 8, &data[1 * 64], stride);          if (cbp& 4) transfer_op(pY_Cur + next_block + 8,        &data[3 * 64], stride);
380          if (cbp & 8)          if (cbp& 2) transfer_op(pU_Cur,                                         &data[4 * 64], stride2);
381                  transfer_16to8add(pY_Cur + next_block, &data[2 * 64], stride);          if (cbp& 1) transfer_op(pV_Cur,                                         &data[5 * 64], stride2);
         if (cbp & 4)  
                 transfer_16to8add(pY_Cur + next_block + 8, &data[3 * 64], stride);  
         if (cbp & 2)  
                 transfer_16to8add(pU_Cur, &data[4 * 64], stride2);  
         if (cbp & 1)  
                 transfer_16to8add(pV_Cur, &data[5 * 64], stride2);  
382          stop_transfer_timer();          stop_transfer_timer();
383  }  }
384    
385    /*****************************************************************************
386     * Module functions
387     ****************************************************************************/
388    
389    void
390    MBTransQuantIntra(const MBParam * const pParam,
391                                      const FRAMEINFO * const frame,
392                                      MACROBLOCK * const pMB,
393                                      const uint32_t x_pos,
394                                      const uint32_t y_pos,
395                                      int16_t data[6 * 64],
396                                      int16_t qcoeff[6 * 64])
397    {
398    
399  /* if sum(diff between field lines) < sum(diff between frame lines), use field dct */          /* Transfer data */
400            MBTrans8to16(pParam, frame, pMB, x_pos, y_pos, data);
401    
402            /* Perform DCT (and field decision) */
403            MBfDCT(pParam, frame, pMB, x_pos, y_pos, data);
404    
405  uint32_t          /* Quantize the block */
406  MBDecideFieldDCT(int16_t data[6 * 64])          MBQuantIntra(pParam, frame, pMB, data, qcoeff);
407    
408            /* DeQuantize the block */
409            MBDeQuantIntra(pParam, pMB->quant, data, qcoeff);
410    
411            /* Perform inverse DCT*/
412            MBiDCT(data, 0x3F);
413    
414            /* Transfer back the data -- Don't add data */
415            MBTrans16to8(pParam, frame, pMB, x_pos, y_pos, data, 0, 0x3F);
416    }
417    
418    
419    uint8_t
420    MBTransQuantInter(const MBParam * const pParam,
421                                      const FRAMEINFO * const frame,
422                                      MACROBLOCK * const pMB,
423                                      const uint32_t x_pos,
424                                      const uint32_t y_pos,
425                                      int16_t data[6 * 64],
426                                      int16_t qcoeff[6 * 64])
427  {  {
428            uint8_t cbp;
429            uint32_t limit;
430    
431            /* There is no MBTrans8to16 for Inter block, that's done in motion compensation
432             * already */
433    
434            /* Perform DCT (and field decision) */
435            MBfDCT(pParam, frame, pMB, x_pos, y_pos, data);
436    
437            /* Set the limit threshold */
438            limit = PVOP_TOOSMALL_LIMIT + ((pMB->quant == 1)? 1 : 0);
439    
440            if (frame->vop_flags & XVID_VOP_CARTOON)
441                    limit *= 3;
442    
443            /* Quantize the block */
444            cbp = MBQuantInter(pParam, frame, pMB, data, qcoeff, 0, limit);
445    
446            /* DeQuantize the block */
447            MBDeQuantInter(pParam, pMB->quant, data, qcoeff, cbp);
448    
449            /* Perform inverse DCT*/
450            MBiDCT(data, cbp);
451    
452            /* Transfer back the data -- Add the data */
453            MBTrans16to8(pParam, frame, pMB, x_pos, y_pos, data, 1, cbp);
454    
455            return(cbp);
456    }
457    
458    uint8_t
459    MBTransQuantInterBVOP(const MBParam * pParam,
460                                              FRAMEINFO * frame,
461                                              MACROBLOCK * pMB,
462                                              const uint32_t x_pos,
463                                              const uint32_t y_pos,
464                                              int16_t data[6 * 64],
465                                              int16_t qcoeff[6 * 64])
466    {
467            uint8_t cbp;
468            uint32_t limit;
469    
470            /* There is no MBTrans8to16 for Inter block, that's done in motion compensation
471             * already */
472    
473            /* Perform DCT (and field decision) */
474            MBfDCT(pParam, frame, pMB, x_pos, y_pos, data);
475    
476            /* Set the limit threshold */
477            limit = BVOP_TOOSMALL_LIMIT;
478    
479            if (frame->vop_flags & XVID_VOP_CARTOON)
480                    limit *= 2;
481    
482            /* Quantize the block */
483            cbp = MBQuantInter(pParam, frame, pMB, data, qcoeff, 1, limit);
484    
485            /*
486             * History comment:
487             * We don't have to DeQuant, iDCT and Transfer back data for B-frames.
488             *
489             * BUT some plugins require the rebuilt original frame to be passed so we
490             * have to take care of that here
491             */
492            if((pParam->plugin_flags & XVID_REQORIGINAL)) {
493    
494                    /* DeQuantize the block */
495                    MBDeQuantInter(pParam, pMB->quant, data, qcoeff, cbp);
496    
497                    /* Perform inverse DCT*/
498                    MBiDCT(data, cbp);
499    
500                    /* Transfer back the data -- Add the data */
501                    MBTrans16to8(pParam, frame, pMB, x_pos, y_pos, data, 1, cbp);
502            }
503    
504            return(cbp);
505    }
506    
507    /* if sum(diff between field lines) < sum(diff between frame lines), use field dct */
508    uint32_t
509    MBFieldTest_c(int16_t data[6 * 64])
510    {
511          const uint8_t blocks[] =          const uint8_t blocks[] =
512                  { 0 * 64, 0 * 64, 0 * 64, 0 * 64, 2 * 64, 2 * 64, 2 * 64, 2 * 64 };                  { 0 * 64, 0 * 64, 0 * 64, 0 * 64, 2 * 64, 2 * 64, 2 * 64, 2 * 64 };
513          const uint8_t lines[] = { 0, 16, 32, 48, 0, 16, 32, 48 };          const uint8_t lines[] = { 0, 16, 32, 48, 0, 16, 32, 48 };
# Line 577  Line 518 
518          for (i = 0; i < 7; ++i) {          for (i = 0; i < 7; ++i) {
519                  for (j = 0; j < 8; ++j) {                  for (j = 0; j < 8; ++j) {
520                          frame +=                          frame +=
521                                  ABS(data[0 * 64 + (i + 1) * 8 + j] - data[0 * 64 + i * 8 + j]);                                  abs(data[0 * 64 + (i + 1) * 8 + j] - data[0 * 64 + i * 8 + j]);
522                          frame +=                          frame +=
523                                  ABS(data[1 * 64 + (i + 1) * 8 + j] - data[1 * 64 + i * 8 + j]);                                  abs(data[1 * 64 + (i + 1) * 8 + j] - data[1 * 64 + i * 8 + j]);
524                          frame +=                          frame +=
525                                  ABS(data[2 * 64 + (i + 1) * 8 + j] - data[2 * 64 + i * 8 + j]);                                  abs(data[2 * 64 + (i + 1) * 8 + j] - data[2 * 64 + i * 8 + j]);
526                          frame +=                          frame +=
527                                  ABS(data[3 * 64 + (i + 1) * 8 + j] - data[3 * 64 + i * 8 + j]);                                  abs(data[3 * 64 + (i + 1) * 8 + j] - data[3 * 64 + i * 8 + j]);
528    
529                          field +=                          field +=
530                                  ABS(data[blocks[i + 1] + lines[i + 1] + j] -                                  abs(data[blocks[i + 1] + lines[i + 1] + j] -
531                                          data[blocks[i] + lines[i] + j]);                                          data[blocks[i] + lines[i] + j]);
532                          field +=                          field +=
533                                  ABS(data[blocks[i + 1] + lines[i + 1] + 8 + j] -                                  abs(data[blocks[i + 1] + lines[i + 1] + 8 + j] -
534                                          data[blocks[i] + lines[i] + 8 + j]);                                          data[blocks[i] + lines[i] + 8 + j]);
535                          field +=                          field +=
536                                  ABS(data[blocks[i + 1] + 64 + lines[i + 1] + j] -                                  abs(data[blocks[i + 1] + 64 + lines[i + 1] + j] -
537                                          data[blocks[i] + 64 + lines[i] + j]);                                          data[blocks[i] + 64 + lines[i] + j]);
538                          field +=                          field +=
539                                  ABS(data[blocks[i + 1] + 64 + lines[i + 1] + 8 + j] -                                  abs(data[blocks[i + 1] + 64 + lines[i + 1] + 8 + j] -
540                                          data[blocks[i] + 64 + lines[i] + 8 + j]);                                          data[blocks[i] + 64 + lines[i] + 8 + j]);
541                  }                  }
542          }          }
543    
544          if (frame > (field + 350)) {          return (frame >= (field + 350));
                 MBFrameToField(data);  
         }  
   
         return (frame > (field + 350));  
545  }  }
546    
547    
# Line 620  Line 557 
557    
558          /* left blocks */          /* left blocks */
559    
560          // 1=2, 2=4, 4=8, 8=1          /* 1=2, 2=4, 4=8, 8=1 */
561          MOVLINE(tmp, LINE(0, 1));          MOVLINE(tmp, LINE(0, 1));
562          MOVLINE(LINE(0, 1), LINE(0, 2));          MOVLINE(LINE(0, 1), LINE(0, 2));
563          MOVLINE(LINE(0, 2), LINE(0, 4));          MOVLINE(LINE(0, 2), LINE(0, 4));
564          MOVLINE(LINE(0, 4), LINE(2, 0));          MOVLINE(LINE(0, 4), LINE(2, 0));
565          MOVLINE(LINE(2, 0), tmp);          MOVLINE(LINE(2, 0), tmp);
566    
567          // 3=6, 6=12, 12=9, 9=3          /* 3=6, 6=12, 12=9, 9=3 */
568          MOVLINE(tmp, LINE(0, 3));          MOVLINE(tmp, LINE(0, 3));
569          MOVLINE(LINE(0, 3), LINE(0, 6));          MOVLINE(LINE(0, 3), LINE(0, 6));
570          MOVLINE(LINE(0, 6), LINE(2, 4));          MOVLINE(LINE(0, 6), LINE(2, 4));
571          MOVLINE(LINE(2, 4), LINE(2, 1));          MOVLINE(LINE(2, 4), LINE(2, 1));
572          MOVLINE(LINE(2, 1), tmp);          MOVLINE(LINE(2, 1), tmp);
573    
574          // 5=10, 10=5          /* 5=10, 10=5 */
575          MOVLINE(tmp, LINE(0, 5));          MOVLINE(tmp, LINE(0, 5));
576          MOVLINE(LINE(0, 5), LINE(2, 2));          MOVLINE(LINE(0, 5), LINE(2, 2));
577          MOVLINE(LINE(2, 2), tmp);          MOVLINE(LINE(2, 2), tmp);
578    
579          // 7=14, 14=13, 13=11, 11=7          /* 7=14, 14=13, 13=11, 11=7 */
580          MOVLINE(tmp, LINE(0, 7));          MOVLINE(tmp, LINE(0, 7));
581          MOVLINE(LINE(0, 7), LINE(2, 6));          MOVLINE(LINE(0, 7), LINE(2, 6));
582          MOVLINE(LINE(2, 6), LINE(2, 5));          MOVLINE(LINE(2, 6), LINE(2, 5));
# Line 648  Line 585 
585    
586          /* right blocks */          /* right blocks */
587    
588          // 1=2, 2=4, 4=8, 8=1          /* 1=2, 2=4, 4=8, 8=1 */
589          MOVLINE(tmp, LINE(1, 1));          MOVLINE(tmp, LINE(1, 1));
590          MOVLINE(LINE(1, 1), LINE(1, 2));          MOVLINE(LINE(1, 1), LINE(1, 2));
591          MOVLINE(LINE(1, 2), LINE(1, 4));          MOVLINE(LINE(1, 2), LINE(1, 4));
592          MOVLINE(LINE(1, 4), LINE(3, 0));          MOVLINE(LINE(1, 4), LINE(3, 0));
593          MOVLINE(LINE(3, 0), tmp);          MOVLINE(LINE(3, 0), tmp);
594    
595          // 3=6, 6=12, 12=9, 9=3          /* 3=6, 6=12, 12=9, 9=3 */
596          MOVLINE(tmp, LINE(1, 3));          MOVLINE(tmp, LINE(1, 3));
597          MOVLINE(LINE(1, 3), LINE(1, 6));          MOVLINE(LINE(1, 3), LINE(1, 6));
598          MOVLINE(LINE(1, 6), LINE(3, 4));          MOVLINE(LINE(1, 6), LINE(3, 4));
599          MOVLINE(LINE(3, 4), LINE(3, 1));          MOVLINE(LINE(3, 4), LINE(3, 1));
600          MOVLINE(LINE(3, 1), tmp);          MOVLINE(LINE(3, 1), tmp);
601    
602          // 5=10, 10=5          /* 5=10, 10=5 */
603          MOVLINE(tmp, LINE(1, 5));          MOVLINE(tmp, LINE(1, 5));
604          MOVLINE(LINE(1, 5), LINE(3, 2));          MOVLINE(LINE(1, 5), LINE(3, 2));
605          MOVLINE(LINE(3, 2), tmp);          MOVLINE(LINE(3, 2), tmp);
606    
607          // 7=14, 14=13, 13=11, 11=7          /* 7=14, 14=13, 13=11, 11=7 */
608          MOVLINE(tmp, LINE(1, 7));          MOVLINE(tmp, LINE(1, 7));
609          MOVLINE(LINE(1, 7), LINE(3, 6));          MOVLINE(LINE(1, 7), LINE(3, 6));
610          MOVLINE(LINE(3, 6), LINE(3, 5));          MOVLINE(LINE(3, 6), LINE(3, 5));
611          MOVLINE(LINE(3, 5), LINE(3, 3));          MOVLINE(LINE(3, 5), LINE(3, 3));
612          MOVLINE(LINE(3, 3), tmp);          MOVLINE(LINE(3, 3), tmp);
613  }  }
614    
615    /*****************************************************************************
616     *               Trellis based R-D optimal quantization
617     *
618     *   Trellis Quant code (C) 2003 Pascal Massimino skal(at)planet-d.net
619     *
620     ****************************************************************************/
621    
622    /*----------------------------------------------------------------------------
623     *
624     *        Trellis-Based quantization
625     *
626     * So far I understand this paper:
627     *
628     *  "Trellis-Based R-D Optimal Quantization in H.263+"
629     *    J.Wen, M.Luttrell, J.Villasenor
630     *    IEEE Transactions on Image Processing, Vol.9, No.8, Aug. 2000.
631     *
632     * we are at stake with a simplified Bellmand-Ford / Dijkstra Single
633     * Source Shortest Path algo. But due to the underlying graph structure
634     * ("Trellis"), it can be turned into a dynamic programming algo,
635     * partially saving the explicit graph's nodes representation. And
636     * without using a heap, since the open frontier of the DAG is always
637     * known, and of fixed size.
638     *--------------------------------------------------------------------------*/
639    
640    
641    
642    /* Codes lengths for relevant levels. */
643    
644    /* let's factorize: */
645    static const uint8_t Code_Len0[64] = {
646            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
647            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
648    static const uint8_t Code_Len1[64] = {
649            20,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
650            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
651    static const uint8_t Code_Len2[64] = {
652            19,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
653            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
654    static const uint8_t Code_Len3[64] = {
655            18,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
656            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
657    static const uint8_t Code_Len4[64] = {
658            17,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
659            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
660    static const uint8_t Code_Len5[64] = {
661            16,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
662            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
663    static const uint8_t Code_Len6[64] = {
664            15,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
665            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
666    static const uint8_t Code_Len7[64] = {
667            13,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
668            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
669    static const uint8_t Code_Len8[64] = {
670            11,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
671            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
672    static const uint8_t Code_Len9[64] = {
673            12,21,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
674            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
675    static const uint8_t Code_Len10[64] = {
676            12,20,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
677            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
678    static const uint8_t Code_Len11[64] = {
679            12,19,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
680            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
681    static const uint8_t Code_Len12[64] = {
682            11,17,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
683            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
684    static const uint8_t Code_Len13[64] = {
685            11,15,21,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
686            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
687    static const uint8_t Code_Len14[64] = {
688            10,12,19,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
689            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
690    static const uint8_t Code_Len15[64] = {
691            10,13,17,19,21,21,21,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
692            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
693    static const uint8_t Code_Len16[64] = {
694            9,12,13,18,18,19,19,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
695            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30};
696    static const uint8_t Code_Len17[64] = {
697            8,11,13,14,14,14,15,19,19,19,21,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
698            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
699    static const uint8_t Code_Len18[64] = {
700            7, 9,11,11,13,13,13,15,15,15,16,22,22,22,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
701            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
702    static const uint8_t Code_Len19[64] = {
703            5, 7, 9,10,10,11,11,11,11,11,13,14,16,17,17,18,18,18,18,18,18,18,18,20,20,21,21,30,30,30,30,30,
704            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30 };
705    static const uint8_t Code_Len20[64] = {
706            3, 4, 5, 6, 6, 6, 7, 7, 7, 7, 8, 8, 8, 9, 9,10,10,10,10,10,10,10,10,12,12,13,13,12,13,14,15,15,
707            15,16,16,16,16,17,17,17,18,18,19,19,19,19,19,19,19,19,21,21,22,22,30,30,30,30,30,30,30,30,30,30 };
708    
709    /* a few more table for LAST table: */
710    static const uint8_t Code_Len21[64] = {
711            13,20,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
712            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30};
713    static const uint8_t Code_Len22[64] = {
714            12,15,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,
715            30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30};
716    static const uint8_t Code_Len23[64] = {
717            10,12,15,15,15,16,16,16,16,17,17,17,17,17,17,17,17,18,18,18,18,18,18,18,18,19,19,19,19,20,20,20,
718            20,21,21,21,21,21,21,21,21,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30};
719    static const uint8_t Code_Len24[64] = {
720            5, 7, 7, 7, 7, 8, 8, 8, 8, 9, 9, 9, 9, 9, 9, 9, 9,10,10,10,10,10,10,10,10,11,11,11,11,12,12,12,
721            12,13,13,13,13,13,13,13,13,14,16,16,16,16,17,17,17,17,18,18,18,18,18,18,18,18,19,19,19,19,19,19};
722    
723    
724    static const uint8_t * const B16_17_Code_Len[24] = { /* levels [1..24] */
725            Code_Len20,Code_Len19,Code_Len18,Code_Len17,
726            Code_Len16,Code_Len15,Code_Len14,Code_Len13,
727            Code_Len12,Code_Len11,Code_Len10,Code_Len9,
728            Code_Len8, Code_Len7 ,Code_Len6 ,Code_Len5,
729            Code_Len4, Code_Len3, Code_Len3 ,Code_Len2,
730            Code_Len2, Code_Len1, Code_Len1, Code_Len1,
731    };
732    
733    static const uint8_t * const B16_17_Code_Len_Last[6] = { /* levels [1..6] */
734            Code_Len24,Code_Len23,Code_Len22,Code_Len21, Code_Len3, Code_Len1,
735    };
736    
737    /* TL_SHIFT controls the precision of the RD optimizations in trellis
738     * valid range is [10..16]. The bigger, the more trellis is vulnerable
739     * to overflows in cost formulas.
740     *  - 10 allows ac values up to 2^11 == 2048
741     *  - 16 allows ac values up to 2^8 == 256
742     */
743    #define TL_SHIFT 11
744    #define TL(q) ((0xfe00>>(16-TL_SHIFT))/(q*q))
745    
746    static const int Trellis_Lambda_Tabs[31] = {
747            TL( 1),TL( 2),TL( 3),TL( 4),TL( 5),TL( 6), TL( 7),
748            TL( 8),TL( 9),TL(10),TL(11),TL(12),TL(13),TL(14), TL(15),
749            TL(16),TL(17),TL(18),TL(19),TL(20),TL(21),TL(22), TL(23),
750            TL(24),TL(25),TL(26),TL(27),TL(28),TL(29),TL(30), TL(31)
751    };
752    #undef TL
753    
754    static int __inline
755    Find_Last(const int16_t *C, const uint16_t *Zigzag, int i)
756    {
757            while(i>=0)
758                    if (C[Zigzag[i]])
759                            return i;
760                    else i--;
761            return -1;
762    }
763    
764    #define TRELLIS_MIN_EFFORT      3
765    
766    /* this routine has been strippen of all debug code */
767    static int
768    dct_quantize_trellis_c(int16_t *const Out,
769                                               const int16_t *const In,
770                                               int Q,
771                                               const uint16_t * const Zigzag,
772                                               const uint16_t * const QuantMatrix,
773                                               int Non_Zero,
774                                               int Sum,
775                                               int Lambda_Mod)
776    {
777    
778            /* Note: We should search last non-zero coeffs on *real* DCT input coeffs
779             * (In[]), not quantized one (Out[]). However, it only improves the result
780             * *very* slightly (~0.01dB), whereas speed drops to crawling level :)
781             * Well, actually, taking 1 more coeff past Non_Zero into account sometimes
782             * helps. */
783            typedef struct { int16_t Run, Level; } NODE;
784    
785            NODE Nodes[65], Last = { 0, 0};
786            uint32_t Run_Costs0[64+1];
787            uint32_t * const Run_Costs = Run_Costs0 + 1;
788    
789            /* it's 1/lambda, actually */
790            const int Lambda = (Lambda_Mod*Trellis_Lambda_Tabs[Q-1])>>LAMBDA_EXP;
791    
792            int Run_Start = -1;
793            uint32_t Min_Cost = 2<<TL_SHIFT;
794    
795            int Last_Node = -1;
796            uint32_t Last_Cost = 0;
797    
798            int i, j;
799    
800            /* source (w/ CBP penalty) */
801            Run_Costs[-1] = 2<<TL_SHIFT;
802    
803            Non_Zero = Find_Last(Out, Zigzag, Non_Zero);
804            if (Non_Zero < TRELLIS_MIN_EFFORT)
805                    Non_Zero = TRELLIS_MIN_EFFORT;
806    
807            for(i=0; i<=Non_Zero; i++) {
808                    const int q = ((Q*QuantMatrix[Zigzag[i]])>>4);
809                    const int Mult = 2*q;
810                    const int Bias = (q-1) | 1;
811                    const int Lev0 = Mult + Bias;
812    
813                    const int AC = In[Zigzag[i]];
814                    const int Level1 = Out[Zigzag[i]];
815                    const unsigned int Dist0 = Lambda* AC*AC;
816                    uint32_t Best_Cost = 0xf0000000;
817                    Last_Cost += Dist0;
818    
819                    /* very specialized loop for -1,0,+1 */
820                    if ((uint32_t)(Level1+1)<3) {
821                            int dQ;
822                            int Run;
823                            uint32_t Cost0;
824    
825                            if (AC<0) {
826                                    Nodes[i].Level = -1;
827                                    dQ = Lev0 + AC;
828                            } else {
829                                    Nodes[i].Level = 1;
830                                    dQ = Lev0 - AC;
831                            }
832                            Cost0 = Lambda*dQ*dQ;
833    
834                            Nodes[i].Run = 1;
835                            Best_Cost = (Code_Len20[0]<<TL_SHIFT) + Run_Costs[i-1]+Cost0;
836                            for(Run=i-Run_Start; Run>0; --Run) {
837                                    const uint32_t Cost_Base = Cost0 + Run_Costs[i-Run];
838                                    const uint32_t Cost = Cost_Base + (Code_Len20[Run-1]<<TL_SHIFT);
839                                    const uint32_t lCost = Cost_Base + (Code_Len24[Run-1]<<TL_SHIFT);
840    
841                                    /* TODO: what about tie-breaks? Should we favor short runs or
842                                     * long runs? Although the error is the same, it would not be
843                                     * spread the same way along high and low frequencies... */
844    
845                                    /* Gruel: I'd say, favour short runs => hifreq errors (HVS) */
846    
847                                    if (Cost<Best_Cost) {
848                                            Best_Cost        = Cost;
849                                            Nodes[i].Run = Run;
850                                    }
851    
852                                    if (lCost<Last_Cost) {
853                                            Last_Cost  = lCost;
854                                            Last.Run   = Run;
855                                            Last_Node  = i;
856                                    }
857                            }
858                            if (Last_Node==i)
859                                    Last.Level = Nodes[i].Level;
860                    } else if (51U>(uint32_t)(Level1+25)) {
861                            /* "big" levels (not less than ESC3, though) */
862                            const uint8_t *Tbl_L1, *Tbl_L2, *Tbl_L1_Last, *Tbl_L2_Last;
863                            int Level2;
864                            int dQ1, dQ2;
865                            int Run;
866                            uint32_t Dist1,Dist2;
867                            int dDist21;
868    
869                            if (Level1>1) {
870                                    dQ1 = Level1*Mult-AC + Bias;
871                                    dQ2 = dQ1 - Mult;
872                                    Level2 = Level1-1;
873                                    Tbl_L1          = (Level1<=24) ? B16_17_Code_Len[Level1-1]         : Code_Len0;
874                                    Tbl_L2          = (Level2<=24) ? B16_17_Code_Len[Level2-1]         : Code_Len0;
875                                    Tbl_L1_Last = (Level1<=6) ? B16_17_Code_Len_Last[Level1-1] : Code_Len0;
876                                    Tbl_L2_Last = (Level2<=6) ? B16_17_Code_Len_Last[Level2-1] : Code_Len0;
877                            } else { /* Level1<-1 */
878                                    dQ1 = Level1*Mult-AC - Bias;
879                                    dQ2 = dQ1 + Mult;
880                                    Level2 = Level1 + 1;
881                                    Tbl_L1          = (Level1>=-24) ? B16_17_Code_Len[Level1^-1]      : Code_Len0;
882                                    Tbl_L2          = (Level2>=-24) ? B16_17_Code_Len[Level2^-1]      : Code_Len0;
883                                    Tbl_L1_Last = (Level1>=- 6) ? B16_17_Code_Len_Last[Level1^-1] : Code_Len0;
884                                    Tbl_L2_Last = (Level2>=- 6) ? B16_17_Code_Len_Last[Level2^-1] : Code_Len0;
885                            }
886    
887                            Dist1 = Lambda*dQ1*dQ1;
888                            Dist2 = Lambda*dQ2*dQ2;
889                            dDist21 = Dist2-Dist1;
890    
891                            for(Run=i-Run_Start; Run>0; --Run)
892                            {
893                                    const uint32_t Cost_Base = Dist1 + Run_Costs[i-Run];
894                                    uint32_t Cost1, Cost2;
895                                    int bLevel;
896    
897                                    /* for sub-optimal (but slightly worth it, speed-wise) search,
898                                     * uncomment the following:
899                                     *              if (Cost_Base>=Best_Cost) continue;
900                                     * (? doesn't seem to have any effect -- gruel ) */
901    
902                                    Cost1 = Cost_Base + (Tbl_L1[Run-1]<<TL_SHIFT);
903                                    Cost2 = Cost_Base + (Tbl_L2[Run-1]<<TL_SHIFT) + dDist21;
904    
905                                    if (Cost2<Cost1) {
906                                            Cost1 = Cost2;
907                                            bLevel = Level2;
908                                    } else {
909                                            bLevel = Level1;
910                                    }
911    
912                                    if (Cost1<Best_Cost) {
913                                            Best_Cost = Cost1;
914                                            Nodes[i].Run   = Run;
915                                            Nodes[i].Level = bLevel;
916                                    }
917    
918                                    Cost1 = Cost_Base + (Tbl_L1_Last[Run-1]<<TL_SHIFT);
919                                    Cost2 = Cost_Base + (Tbl_L2_Last[Run-1]<<TL_SHIFT) + dDist21;
920    
921                                    if (Cost2<Cost1) {
922                                            Cost1 = Cost2;
923                                            bLevel = Level2;
924                                    } else {
925                                            bLevel = Level1;
926                                    }
927    
928                                    if (Cost1<Last_Cost) {
929                                            Last_Cost  = Cost1;
930                                            Last.Run   = Run;
931                                            Last.Level = bLevel;
932                                            Last_Node  = i;
933                                    }
934                            } /* end of "for Run" */
935                    } else {
936                            /* Very very high levels, with no chance of being optimizable
937                             * => Simply pick best Run. */
938                            int Run;
939                            for(Run=i-Run_Start; Run>0; --Run) {
940                                    /* 30 bits + no distortion */
941                                    const uint32_t Cost = (30<<TL_SHIFT) + Run_Costs[i-Run];
942                                    if (Cost<Best_Cost) {
943                                            Best_Cost = Cost;
944                                            Nodes[i].Run   = Run;
945                                            Nodes[i].Level = Level1;
946                                    }
947    
948                                    if (Cost<Last_Cost) {
949                                            Last_Cost  = Cost;
950                                            Last.Run   = Run;
951                                            Last.Level = Level1;
952                                            Last_Node  = i;
953                                    }
954                            }
955                    }
956    
957    
958                    Run_Costs[i] = Best_Cost;
959    
960                    if (Best_Cost < Min_Cost + Dist0) {
961                            Min_Cost = Best_Cost;
962                            Run_Start = i;
963                    } else {
964                            /* as noticed by Michael Niedermayer (michaelni at gmx.at),
965                             * there's a code shorter by 1 bit for a larger run (!), same
966                             * level. We give it a chance by not moving the left barrier too
967                             * much. */
968                            while( Run_Costs[Run_Start]>Min_Cost+(1<<TL_SHIFT) )
969                                    Run_Start++;
970    
971                            /* spread on preceding coeffs the cost incurred by skipping this
972                             * one */
973                            for(j=Run_Start; j<i; ++j) Run_Costs[j] += Dist0;
974                            Min_Cost += Dist0;
975                    }
976            }
977    
978            /* It seems trellis doesn't give good results... just leave the block untouched
979             * and return the original sum value */
980            if (Last_Node<0)
981                    return Sum;
982    
983            /* reconstruct optimal sequence backward with surviving paths */
984            memset(Out, 0x00, 64*sizeof(*Out));
985            Out[Zigzag[Last_Node]] = Last.Level;
986            i = Last_Node - Last.Run;
987            Sum = abs(Last.Level);
988            while(i>=0) {
989                    Out[Zigzag[i]] = Nodes[i].Level;
990                    Sum += abs(Nodes[i].Level);
991                    i -= Nodes[i].Run;
992            }
993    
994            return Sum;
995    }
996    
997    /* original version including heavy debugging info */
998    
999    #ifdef DBGTRELL
1000    
1001    #define DBG 0
1002    
1003    static __inline uint32_t Evaluate_Cost(const int16_t *C, int Mult, int Bias,
1004                                                                               const uint16_t * Zigzag, int Max, int Lambda)
1005    {
1006    #if (DBG>0)
1007            const int16_t * const Ref = C + 6*64;
1008            int Last = Max;
1009            int Bits = 0;
1010            int Dist = 0;
1011            int i;
1012            uint32_t Cost;
1013    
1014            while(Last>=0 && C[Zigzag[Last]]==0)
1015                    Last--;
1016    
1017            if (Last>=0) {
1018                    int j=0, j0=0;
1019                    int Run, Level;
1020    
1021                    Bits = 2;   /* CBP */
1022                    while(j<Last) {
1023                            while(!C[Zigzag[j]])
1024                                    j++;
1025                            if (j==Last)
1026                                    break;
1027                            Level=C[Zigzag[j]];
1028                            Run = j - j0;
1029                            j0 = ++j;
1030                            if (Level>=-24 && Level<=24)
1031                                    Bits += B16_17_Code_Len[(Level<0) ? -Level-1 : Level-1][Run];
1032                            else
1033                                    Bits += 30;
1034                    }
1035                    Level = C[Zigzag[Last]];
1036                    Run = j - j0;
1037                    if (Level>=-6 && Level<=6)
1038                            Bits += B16_17_Code_Len_Last[(Level<0) ? -Level-1 : Level-1][Run];
1039                    else
1040                            Bits += 30;
1041            }
1042    
1043            for(i=0; i<=Last; ++i) {
1044                    int V = C[Zigzag[i]]*Mult;
1045                    if (V>0)
1046                            V += Bias;
1047                    else
1048                            if (V<0)
1049                                    V -= Bias;
1050                    V -= Ref[Zigzag[i]];
1051                    Dist += V*V;
1052            }
1053            Cost = Lambda*Dist + (Bits<<TL_SHIFT);
1054            if (DBG==1)
1055                    printf( " Last:%2d/%2d Cost = [(Bits=%5.0d) + Lambda*(Dist=%6.0d) = %d ] >>12= %d ", Last,Max, Bits, Dist, Cost, Cost>>12 );
1056            return Cost;
1057    
1058    #else
1059            return 0;
1060    #endif
1061    }
1062    
1063    
1064    static int
1065    dct_quantize_trellis_h263_c(int16_t *const Out, const int16_t *const In, int Q, const uint16_t * const Zigzag, int Non_Zero)
1066    {
1067    
1068        /*
1069             * Note: We should search last non-zero coeffs on *real* DCT input coeffs (In[]),
1070             * not quantized one (Out[]). However, it only improves the result *very*
1071             * slightly (~0.01dB), whereas speed drops to crawling level :)
1072             * Well, actually, taking 1 more coeff past Non_Zero into account sometimes helps.
1073             */
1074            typedef struct { int16_t Run, Level; } NODE;
1075    
1076            NODE Nodes[65], Last;
1077            uint32_t Run_Costs0[64+1];
1078            uint32_t * const Run_Costs = Run_Costs0 + 1;
1079            const int Mult = 2*Q;
1080            const int Bias = (Q-1) | 1;
1081            const int Lev0 = Mult + Bias;
1082            const int Lambda = Trellis_Lambda_Tabs[Q-1];    /* it's 1/lambda, actually */
1083    
1084            int Run_Start = -1;
1085            Run_Costs[-1] = 2<<TL_SHIFT;                          /* source (w/ CBP penalty) */
1086            uint32_t Min_Cost = 2<<TL_SHIFT;
1087    
1088            int Last_Node = -1;
1089            uint32_t Last_Cost = 0;
1090    
1091            int i, j;
1092    
1093    #if (DBG>0)
1094            Last.Level = 0; Last.Run = -1; /* just initialize to smthg */
1095    #endif
1096    
1097            Non_Zero = Find_Last(Out, Zigzag, Non_Zero);
1098            if (Non_Zero<0)
1099                    return -1;
1100    
1101            for(i=0; i<=Non_Zero; i++)
1102            {
1103                    const int AC = In[Zigzag[i]];
1104                    const int Level1 = Out[Zigzag[i]];
1105                    const int Dist0 = Lambda* AC*AC;
1106                    uint32_t Best_Cost = 0xf0000000;
1107                    Last_Cost += Dist0;
1108    
1109                    if ((uint32_t)(Level1+1)<3)                 /* very specialized loop for -1,0,+1 */
1110                    {
1111                            int dQ;
1112                            int Run;
1113                            uint32_t Cost0;
1114    
1115                            if (AC<0) {
1116                                    Nodes[i].Level = -1;
1117                                    dQ = Lev0 + AC;
1118                            } else {
1119                                    Nodes[i].Level = 1;
1120                                    dQ = Lev0 - AC;
1121                            }
1122                            Cost0 = Lambda*dQ*dQ;
1123    
1124                            Nodes[i].Run = 1;
1125                            Best_Cost = (Code_Len20[0]<<TL_SHIFT) + Run_Costs[i-1]+Cost0;
1126                            for(Run=i-Run_Start; Run>0; --Run)
1127                            {
1128                                    const uint32_t Cost_Base = Cost0 + Run_Costs[i-Run];
1129                                    const uint32_t Cost = Cost_Base + (Code_Len20[Run-1]<<TL_SHIFT);
1130                                    const uint32_t lCost = Cost_Base + (Code_Len24[Run-1]<<TL_SHIFT);
1131    
1132                                    /*
1133                                     * TODO: what about tie-breaks? Should we favor short runs or
1134                                     * long runs? Although the error is the same, it would not be
1135                                     * spread the same way along high and low frequencies...
1136                                     */
1137                                    if (Cost<Best_Cost) {
1138                                            Best_Cost    = Cost;
1139                                            Nodes[i].Run = Run;
1140                                    }
1141    
1142                                    if (lCost<Last_Cost) {
1143                                            Last_Cost  = lCost;
1144                                            Last.Run   = Run;
1145                                            Last_Node  = i;
1146                                    }
1147                            }
1148                            if (Last_Node==i)
1149                                    Last.Level = Nodes[i].Level;
1150    
1151                            if (DBG==1) {
1152                                    Run_Costs[i] = Best_Cost;
1153                                    printf( "Costs #%2d: ", i);
1154                                    for(j=-1;j<=Non_Zero;++j) {
1155                                            if (j==Run_Start)            printf( " %3.0d|", Run_Costs[j]>>12 );
1156                                            else if (j>Run_Start && j<i) printf( " %3.0d|", Run_Costs[j]>>12 );
1157                                            else if (j==i)               printf( "(%3.0d)", Run_Costs[j]>>12 );
1158                                            else                         printf( "  - |" );
1159                                    }
1160                                    printf( "<%3.0d %2d %d>", Min_Cost>>12, Nodes[i].Level, Nodes[i].Run );
1161                                    printf( "  Last:#%2d {%3.0d %2d %d}", Last_Node, Last_Cost>>12, Last.Level, Last.Run );
1162                                    printf( " AC:%3.0d Dist0:%3d Dist(%d)=%d", AC, Dist0>>12, Nodes[i].Level, Cost0>>12 );
1163                                    printf( "\n" );
1164                            }
1165                    }
1166                    else                      /* "big" levels */
1167                    {
1168                            const uint8_t *Tbl_L1, *Tbl_L2, *Tbl_L1_Last, *Tbl_L2_Last;
1169                            int Level2;
1170                            int dQ1, dQ2;
1171                            int Run;
1172                            uint32_t Dist1,Dist2;
1173                            int dDist21;
1174    
1175                            if (Level1>1) {
1176                                    dQ1 = Level1*Mult-AC + Bias;
1177                                    dQ2 = dQ1 - Mult;
1178                                    Level2 = Level1-1;
1179                                    Tbl_L1      = (Level1<=24) ? B16_17_Code_Len[Level1-1]     : Code_Len0;
1180                                    Tbl_L2      = (Level2<=24) ? B16_17_Code_Len[Level2-1]     : Code_Len0;
1181                                    Tbl_L1_Last = (Level1<=6) ? B16_17_Code_Len_Last[Level1-1] : Code_Len0;
1182                                    Tbl_L2_Last = (Level2<=6) ? B16_17_Code_Len_Last[Level2-1] : Code_Len0;
1183                            } else { /* Level1<-1 */
1184                                    dQ1 = Level1*Mult-AC - Bias;
1185                                    dQ2 = dQ1 + Mult;
1186                                    Level2 = Level1 + 1;
1187                                    Tbl_L1      = (Level1>=-24) ? B16_17_Code_Len[Level1^-1]      : Code_Len0;
1188                                    Tbl_L2      = (Level2>=-24) ? B16_17_Code_Len[Level2^-1]      : Code_Len0;
1189                                    Tbl_L1_Last = (Level1>=- 6) ? B16_17_Code_Len_Last[Level1^-1] : Code_Len0;
1190                                    Tbl_L2_Last = (Level2>=- 6) ? B16_17_Code_Len_Last[Level2^-1] : Code_Len0;
1191                            }
1192                            Dist1 = Lambda*dQ1*dQ1;
1193                            Dist2 = Lambda*dQ2*dQ2;
1194                            dDist21 = Dist2-Dist1;
1195    
1196                            for(Run=i-Run_Start; Run>0; --Run)
1197                            {
1198                                    const uint32_t Cost_Base = Dist1 + Run_Costs[i-Run];
1199                                    uint32_t Cost1, Cost2;
1200                                    int bLevel;
1201    
1202    /*
1203     * for sub-optimal (but slightly worth it, speed-wise) search, uncomment the following:
1204     *        if (Cost_Base>=Best_Cost) continue;
1205     */
1206                                    Cost1 = Cost_Base + (Tbl_L1[Run-1]<<TL_SHIFT);
1207                                    Cost2 = Cost_Base + (Tbl_L2[Run-1]<<TL_SHIFT) + dDist21;
1208    
1209                                    if (Cost2<Cost1) {
1210                                            Cost1 = Cost2;
1211                                            bLevel = Level2;
1212                                    } else
1213                                            bLevel = Level1;
1214    
1215                                    if (Cost1<Best_Cost) {
1216                                            Best_Cost = Cost1;
1217                                            Nodes[i].Run   = Run;
1218                                            Nodes[i].Level = bLevel;
1219                                    }
1220    
1221                                    Cost1 = Cost_Base + (Tbl_L1_Last[Run-1]<<TL_SHIFT);
1222                                    Cost2 = Cost_Base + (Tbl_L2_Last[Run-1]<<TL_SHIFT) + dDist21;
1223    
1224                                    if (Cost2<Cost1) {
1225                                            Cost1 = Cost2;
1226                                            bLevel = Level2;
1227                                    } else
1228                                            bLevel = Level1;
1229    
1230                                    if (Cost1<Last_Cost) {
1231                                            Last_Cost  = Cost1;
1232                                            Last.Run   = Run;
1233                                            Last.Level = bLevel;
1234                                            Last_Node  = i;
1235                                    }
1236                            } /* end of "for Run" */
1237    
1238                            if (DBG==1) {
1239                                    Run_Costs[i] = Best_Cost;
1240                                    printf( "Costs #%2d: ", i);
1241                                    for(j=-1;j<=Non_Zero;++j) {
1242                                            if (j==Run_Start)            printf( " %3.0d|", Run_Costs[j]>>12 );
1243                                            else if (j>Run_Start && j<i) printf( " %3.0d|", Run_Costs[j]>>12 );
1244                                            else if (j==i)               printf( "(%3.0d)", Run_Costs[j]>>12 );
1245                                            else                         printf( "  - |" );
1246                                    }
1247                                    printf( "<%3.0d %2d %d>", Min_Cost>>12, Nodes[i].Level, Nodes[i].Run );
1248                                    printf( "  Last:#%2d {%3.0d %2d %d}", Last_Node, Last_Cost>>12, Last.Level, Last.Run );
1249                                    printf( " AC:%3.0d Dist0:%3d Dist(%2d):%3d Dist(%2d):%3d", AC, Dist0>>12, Level1, Dist1>>12, Level2, Dist2>>12 );
1250                                    printf( "\n" );
1251                            }
1252                    }
1253    
1254                    Run_Costs[i] = Best_Cost;
1255    
1256                    if (Best_Cost < Min_Cost + Dist0) {
1257                            Min_Cost = Best_Cost;
1258                            Run_Start = i;
1259                    }
1260                    else
1261                    {
1262                            /*
1263                             * as noticed by Michael Niedermayer (michaelni at gmx.at), there's
1264                             * a code shorter by 1 bit for a larger run (!), same level. We give
1265                             * it a chance by not moving the left barrier too much.
1266                             */
1267    
1268                            while( Run_Costs[Run_Start]>Min_Cost+(1<<TL_SHIFT) )
1269                                    Run_Start++;
1270    
1271                            /* spread on preceding coeffs the cost incurred by skipping this one */
1272                            for(j=Run_Start; j<i; ++j) Run_Costs[j] += Dist0;
1273                            Min_Cost += Dist0;
1274                    }
1275            }
1276    
1277            if (DBG) {
1278                    Last_Cost = Evaluate_Cost(Out,Mult,Bias, Zigzag,Non_Zero, Lambda);
1279                    if (DBG==1) {
1280                            printf( "=> " );
1281                            for(i=0; i<=Non_Zero; ++i) printf( "[%3.0d] ", Out[Zigzag[i]] );
1282                            printf( "\n" );
1283                    }
1284            }
1285    
1286            if (Last_Node<0)
1287                    return -1;
1288    
1289            /* reconstruct optimal sequence backward with surviving paths */
1290            memset(Out, 0x00, 64*sizeof(*Out));
1291            Out[Zigzag[Last_Node]] = Last.Level;
1292            i = Last_Node - Last.Run;
1293            while(i>=0) {
1294                    Out[Zigzag[i]] = Nodes[i].Level;
1295                    i -= Nodes[i].Run;
1296            }
1297    
1298            if (DBG) {
1299                    uint32_t Cost = Evaluate_Cost(Out,Mult,Bias, Zigzag,Non_Zero, Lambda);
1300                    if (DBG==1) {
1301                            printf( "<= " );
1302                            for(i=0; i<=Last_Node; ++i) printf( "[%3.0d] ", Out[Zigzag[i]] );
1303                            printf( "\n--------------------------------\n" );
1304                    }
1305                    if (Cost>Last_Cost) printf( "!!! %u > %u\n", Cost, Last_Cost );
1306            }
1307            return Last_Node;
1308    }
1309    
1310    #undef DBG
1311    
1312    #endif

Legend:
Removed from v.605  
changed lines
  Added in v.1713

No admin address has been configured
ViewVC Help
Powered by ViewVC 1.0.4