[svn] / branches / dev-api-4 / xvidcore / src / motion / motion_comp.c Repository:
ViewVC logotype

Diff of /branches/dev-api-4/xvidcore/src/motion/motion_comp.c

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

trunk/xvidcore/src/motion/motion_comp.c revision 437, Sat Sep 7 09:12:22 2002 UTC branches/dev-api-4/xvidcore/src/motion/motion_comp.c revision 982, Thu Apr 10 13:05:54 2003 UTC
# Line 1  Line 1 
1  /*****************************************************************************  // 30.10.2002   corrected qpel chroma rounding
2   *  // 04.10.2002   added qpel support to MBMotionCompensation
3   *  XVID MPEG-4 VIDEO CODEC  // 01.05.2002   updated MBMotionCompensationBVOP
4   *  - Motion Compensation module -  // 14.04.2002   bframe compensation
5   *  
6   *  Copyright(C) 2002 Michael Militzer <michael@xvid.org>  #include <stdio.h>
  *  Copyright(C) 2002 Edouard Gomez <ed.gomez@wanadoo.fr>  
  *  Copyright(C) 2002 Christoph Lampert <gruel@web.de>  
  *  
  *  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.  
  *  
  *  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  
  *  
  *************************************************************************/  
7    
8  #include "../encoder.h"  #include "../encoder.h"
9  #include "../utils/mbfunctions.h"  #include "../utils/mbfunctions.h"
10  #include "../image/interpolate8x8.h"  #include "../image/interpolate8x8.h"
11    #include "../image/reduced.h"
12  #include "../utils/timer.h"  #include "../utils/timer.h"
13  #include "motion.h"  #include "motion.h"
14    
15  #define ABS(X) (((X)>0)?(X):-(X))  #ifndef RSHIFT
16  #define SIGN(X) (((X)>0)?1:-1)  #define RSHIFT(a,b) ((a) > 0 ? ((a) + (1<<((b)-1)))>>(b) : ((a) + (1<<((b)-1))-1)>>(b))
17    #endif
18    
19    /* assume b>0 */
20    #ifndef RDIV
21    #define RDIV(a,b) (((a)>0 ? (a) + ((b)>>1) : (a) - ((b)>>1))/(b))
22    #endif
23    
24    
25    /* This is borrowed from        decoder.c   */
26    static __inline int gmc_sanitize(int value, int quarterpel, int fcode)
27    {
28            int length = 1 << (fcode+4);
29    
30    //      if (quarterpel) value *= 2;
31    
32            if (value < -length)
33                    return -length;
34            else if (value >= length)
35                    return length-1;
36            else return value;
37    }
38    
39    /* And this is borrowed from   bitstream.c  until we find a common solution */
40    
41    static uint32_t __inline
42    log2bin(uint32_t value)
43    {
44    /* Changed by Chenm001 */
45    #if !defined(_MSC_VER)
46            int n = 0;
47    
48            while (value) {
49                    value >>= 1;
50                    n++;
51            }
52            return n;
53    #else
54            __asm {
55                    bsr eax, value
56                    inc eax
57            }
58    #endif
59    }
60    
61    
62  static __inline void  static __inline void
63  compensate8x8_halfpel(int16_t * const dct_codes,  compensate16x16_interpolate(int16_t * const dct_codes,
64                                            uint8_t * const cur,                                            uint8_t * const cur,
65                                            const uint8_t * const ref,                                            const uint8_t * const ref,
66                                            const uint8_t * const refh,                                            const uint8_t * const refh,
67                                            const uint8_t * const refv,                                            const uint8_t * const refv,
68                                            const uint8_t * const refhv,                                            const uint8_t * const refhv,
69                                            const uint32_t x,                                                          uint8_t * const tmp,
70                                            const uint32_t y,                                                          uint32_t x,
71                                                            uint32_t y,
72                                            const int32_t dx,                                            const int32_t dx,
73                                            const int dy,                                                          const int32_t dy,
74                                            const uint32_t stride)                                                          const int32_t stride,
75                                                            const int quarterpel,
76                                                            const int reduced_resolution,
77                                                            const int32_t rounding)
78  {  {
79          int32_t ddx, ddy;          const uint8_t * ptr;
80    
81          switch (((dx & 1) << 1) + (dy & 1))     // ((dx%2)?2:0)+((dy%2)?1:0)          if (!reduced_resolution) {
         {  
         case 0:  
                 ddx = dx / 2;  
                 ddy = dy / 2;  
                 transfer_8to16sub(dct_codes, cur + y * stride + x,  
                                                   ref + (int) ((y + ddy) * stride + x + ddx), stride);  
                 break;  
82    
83          case 1:                  if(quarterpel) {
84                  ddx = dx / 2;                          if ((dx&3) | (dy&3)) {
85                  ddy = (dy - 1) / 2;                                  interpolate16x16_quarterpel(tmp - y * stride - x,
86                  transfer_8to16sub(dct_codes, cur + y * stride + x,                                                                                          (uint8_t *) ref, tmp + 32,
87                                                    refv + (int) ((y + ddy) * stride + x + ddx), stride);                                                                                          tmp + 64, tmp + 96, x, y, dx, dy, stride, rounding);
88                  break;                                  ptr = tmp;
89                            } else ptr =  ref + (y + dy/4)*stride + x + dx/4; // fullpixel position
90    
91          case 2:                  } else ptr = get_ref(ref, refh, refv, refhv, x, y, 1, dx, dy, stride);
                 ddx = (dx - 1) / 2;  
                 ddy = dy / 2;  
                 transfer_8to16sub(dct_codes, cur + y * stride + x,  
                                                   refh + (int) ((y + ddy) * stride + x + ddx), stride);  
                 break;  
92    
         default:                                        // case 3:  
                 ddx = (dx - 1) / 2;  
                 ddy = (dy - 1) / 2;  
93                  transfer_8to16sub(dct_codes, cur + y * stride + x,                  transfer_8to16sub(dct_codes, cur + y * stride + x,
94                                                    refhv + (int) ((y + ddy) * stride + x + ddx), stride);                                                          ptr, stride);
95                  break;                  transfer_8to16sub(dct_codes+64, cur + y * stride + x + 8,
96                                                            ptr + 8, stride);
97                    transfer_8to16sub(dct_codes+128, cur + y * stride + x + 8*stride,
98                                                            ptr + 8*stride, stride);
99                    transfer_8to16sub(dct_codes+192, cur + y * stride + x + 8*stride+8,
100                                                            ptr + 8*stride + 8, stride);
101    
102            } else { //reduced_resolution
103    
104                    x *= 2; y *= 2;
105    
106                    ptr = get_ref(ref, refh, refv, refhv, x, y, 1, dx, dy, stride);
107    
108                    filter_18x18_to_8x8(dct_codes, cur+y*stride + x, stride);
109                    filter_diff_18x18_to_8x8(dct_codes, ptr, stride);
110    
111                    filter_18x18_to_8x8(dct_codes+64, cur+y*stride + x + 16, stride);
112                    filter_diff_18x18_to_8x8(dct_codes+64, ptr + 16, stride);
113    
114                    filter_18x18_to_8x8(dct_codes+128, cur+(y+16)*stride + x, stride);
115                    filter_diff_18x18_to_8x8(dct_codes+128, ptr + 16*stride, stride);
116    
117                    filter_18x18_to_8x8(dct_codes+192, cur+(y+16)*stride + x + 16, stride);
118                    filter_diff_18x18_to_8x8(dct_codes+192, ptr + 16*stride + 16, stride);
119    
120                    transfer32x32_copy(cur + y*stride + x, ptr, stride);
121          }          }
122  }  }
123    
124    static __inline void
125    compensate8x8_interpolate(      int16_t * const dct_codes,
126                                                            uint8_t * const cur,
127                                                            const uint8_t * const ref,
128                                                            const uint8_t * const refh,
129                                                            const uint8_t * const refv,
130                                                            const uint8_t * const refhv,
131                                                            uint8_t * const tmp,
132                                                            uint32_t x,
133                                                            uint32_t y,
134                                                            const int32_t dx,
135                                                            const int32_t dy,
136                                                            const int32_t stride,
137                                                            const int32_t quarterpel,
138                                                            const int reduced_resolution,
139                                                            const int32_t rounding)
140    {
141            const uint8_t * ptr;
142    
143            if (!reduced_resolution) {
144    
145                    if(quarterpel) {
146                            if ((dx&3) | (dy&3)) {
147                                    interpolate8x8_quarterpel(tmp - y*stride - x,
148                                                                                    (uint8_t *) ref, tmp + 32,
149                                                                                    tmp + 64, tmp + 96, x, y, dx, dy, stride, rounding);
150                                    ptr = tmp;
151                            } else ptr = ref + (y + dy/4)*stride + x + dx/4; // fullpixel position
152                    } else ptr = get_ref(ref, refh, refv, refhv, x, y, 1, dx, dy, stride);
153    
154                            transfer_8to16sub(dct_codes, cur + y * stride + x, ptr, stride);
155    
156            } else { //reduced_resolution
157    
158                    x *= 2; y *= 2;
159    
160                    ptr = get_ref(ref, refh, refv, refhv, x, y, 1, dx, dy, stride);
161    
162                    filter_18x18_to_8x8(dct_codes, cur+y*stride + x, stride);
163                    filter_diff_18x18_to_8x8(dct_codes, ptr, stride);
164    
165                    transfer16x16_copy(cur + y*stride + x, ptr, stride);
166            }
167    }
168    
169    /* XXX: slow, inelegant... */
170    static void
171    interpolate18x18_switch(uint8_t * const cur,
172                                                    const uint8_t * const refn,
173                                                    const uint32_t x,
174                                                    const uint32_t y,
175                                                    const int32_t dx,
176                                                    const int dy,
177                                                    const int32_t stride,
178                                                    const int32_t rounding)
179    {
180            interpolate8x8_switch(cur, refn, x-1, y-1, dx, dy, stride, rounding);
181            interpolate8x8_switch(cur, refn, x+7, y-1, dx, dy, stride, rounding);
182            interpolate8x8_switch(cur, refn, x+9, y-1, dx, dy, stride, rounding);
183    
184            interpolate8x8_switch(cur, refn, x-1, y+7, dx, dy, stride, rounding);
185            interpolate8x8_switch(cur, refn, x+7, y+7, dx, dy, stride, rounding);
186            interpolate8x8_switch(cur, refn, x+9, y+7, dx, dy, stride, rounding);
187    
188            interpolate8x8_switch(cur, refn, x-1, y+9, dx, dy, stride, rounding);
189            interpolate8x8_switch(cur, refn, x+7, y+9, dx, dy, stride, rounding);
190            interpolate8x8_switch(cur, refn, x+9, y+9, dx, dy, stride, rounding);
191    }
192    
193    static void
194    CompensateChroma(       int dx, int dy,
195                                            const int i, const int j,
196                                            IMAGE * const Cur,
197                                            const IMAGE * const Ref,
198                                            uint8_t * const temp,
199                                            int16_t * const coeff,
200                                            const int32_t stride,
201                                            const int rounding,
202                                            const int rrv)
203    { /* uv-block-based compensation */
204    
205            if (!rrv) {
206                    transfer_8to16sub(coeff, Cur->u + 8 * j * stride + 8 * i,
207                                                            interpolate8x8_switch2(temp, Ref->u, 8 * i, 8 * j,
208                                                                                                            dx, dy, stride, rounding),
209                                                            stride);
210                    transfer_8to16sub(coeff + 64, Cur->v + 8 * j * stride + 8 * i,
211                                                            interpolate8x8_switch2(temp, Ref->v, 8 * i, 8 * j,
212                                                                                                            dx, dy, stride, rounding),
213                                                            stride);
214            } else {
215                    uint8_t * current, * reference;
216    
217                    current = Cur->u + 16*j*stride + 16*i;
218                    reference = temp - 16*j*stride - 16*i;
219                    interpolate18x18_switch(reference, Ref->u, 16*i, 16*j, dx, dy, stride, rounding);
220                    filter_18x18_to_8x8(coeff, current, stride);
221                    filter_diff_18x18_to_8x8(coeff, temp, stride);
222                    transfer16x16_copy(current, temp, stride);
223    
224                    current = Cur->v + 16*j*stride + 16*i;
225                    interpolate18x18_switch(reference, Ref->v, 16*i, 16*j, dx, dy, stride, rounding);
226                    filter_18x18_to_8x8(coeff + 64, current, stride);
227                    filter_diff_18x18_to_8x8(coeff + 64, temp, stride);
228                    transfer16x16_copy(current, temp, stride);
229            }
230    }
231    
232  void  void
233  MBMotionCompensation(MACROBLOCK * const mb,  MBMotionCompensation(MACROBLOCK * const mb,
# Line 98  Line 237 
237                                           const IMAGE * const refh,                                           const IMAGE * const refh,
238                                           const IMAGE * const refv,                                           const IMAGE * const refv,
239                                           const IMAGE * const refhv,                                           const IMAGE * const refhv,
240                                            const IMAGE * const refGMC,
241                                           IMAGE * const cur,                                           IMAGE * const cur,
242                                           int16_t * dct_codes,                                           int16_t * dct_codes,
243                                           const uint32_t width,                                           const uint32_t width,
244                                           const uint32_t height,                                           const uint32_t height,
245                                           const uint32_t edged_width,                                           const uint32_t edged_width,
246                                           const uint32_t rounding)                                          const int32_t quarterpel,
247                                            const int reduced_resolution,
248                                            const int32_t rounding)
249  {  {
250          static const uint32_t roundtab[16] =          int32_t dx;
251                  { 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2 };          int32_t dy;
252    
253            uint8_t * const tmp = refv->u;
254    
255          if (mb->mode == MODE_INTER || mb->mode == MODE_INTER_Q) {          if ( (!reduced_resolution) && (mb->mode == MODE_NOT_CODED) ) {  /* quick copy for early SKIP */
256                  int32_t dx = mb->mvs[0].x;  /* early SKIP is only activated in P-VOPs, not in S-VOPs, so mcsel can never be 1 */
                 int32_t dy = mb->mvs[0].y;  
257    
258                  compensate8x8_halfpel(&dct_codes[0 * 64], cur->y, ref->y, refh->y,                  transfer16x16_copy(cur->y + 16 * (i + j * edged_width),
259                                                            refv->y, refhv->y, 16 * i, 16 * j, dx, dy,                                                    ref->y + 16 * (i + j * edged_width),
                                                           edged_width);  
                 compensate8x8_halfpel(&dct_codes[1 * 64], cur->y, ref->y, refh->y,  
                                                           refv->y, refhv->y, 16 * i + 8, 16 * j, dx, dy,  
260                                                            edged_width);                                                            edged_width);
                 compensate8x8_halfpel(&dct_codes[2 * 64], cur->y, ref->y, refh->y,  
                                                           refv->y, refhv->y, 16 * i, 16 * j + 8, dx, dy,  
                                                           edged_width);  
                 compensate8x8_halfpel(&dct_codes[3 * 64], cur->y, ref->y, refh->y,  
                                                           refv->y, refhv->y, 16 * i + 8, 16 * j + 8, dx,  
                                                           dy, edged_width);  
   
                 dx = (dx & 3) ? (dx >> 1) | 1 : dx / 2;  
                 dy = (dy & 3) ? (dy >> 1) | 1 : dy / 2;  
   
                 /* uv-image-based compensation */  
   
                 interpolate8x8_switch(refv->u, ref->u, 8 * i, 8 * j, dx, dy,  
                                                           edged_width / 2, rounding);  
                 transfer_8to16sub(&dct_codes[4 * 64],  
                                                   cur->u + 8 * j * edged_width / 2 + 8 * i,  
                                                   refv->u + 8 * j * edged_width / 2 + 8 * i,  
                                                   edged_width / 2);  
261    
262                  interpolate8x8_switch(refv->v, ref->v, 8 * i, 8 * j, dx, dy,                  transfer8x8_copy(cur->u + 8 * (i + j * edged_width/2),
263                                                            edged_width / 2, rounding);                                                          ref->u + 8 * (i + j * edged_width/2),
                 transfer_8to16sub(&dct_codes[5 * 64],  
                                                   cur->v + 8 * j * edged_width / 2 + 8 * i,  
                                                   refv->v + 8 * j * edged_width / 2 + 8 * i,  
264                                                    edged_width / 2);                                                    edged_width / 2);
265                    transfer8x8_copy(cur->v + 8 * (i + j * edged_width/2),
266                                                            ref->v + 8 * (i + j * edged_width/2),
267                                                            edged_width / 2);
268                    return;
269            }
270    
271            if ((mb->mode == MODE_NOT_CODED || mb->mode == MODE_INTER
272                                    || mb->mode == MODE_INTER_Q)) {
273    
274            /* reduced resolution + GMC:  not possible */
275    
276                    if (mb->mcsel) {
277    
278                            /* call normal routine once, easier than "if (mcsel)"ing all the time */
279    
280                            transfer_8to16sub(&dct_codes[0*64], cur->y + 16*j*edged_width + 16*i,
281                                                                                            refGMC->y + 16*j*edged_width + 16*i, edged_width);
282                            transfer_8to16sub(&dct_codes[1*64], cur->y + 16*j*edged_width + 16*i+8,
283                                                                                            refGMC->y + 16*j*edged_width + 16*i+8, edged_width);
284                            transfer_8to16sub(&dct_codes[2*64], cur->y + (16*j+8)*edged_width + 16*i,
285                                                                                            refGMC->y + (16*j+8)*edged_width + 16*i, edged_width);
286                            transfer_8to16sub(&dct_codes[3*64], cur->y + (16*j+8)*edged_width + 16*i+8,
287                                                                                            refGMC->y + (16*j+8)*edged_width + 16*i+8, edged_width);
288    
289    /* lumi is needed earlier for mode decision, but chroma should be done block-based, but it isn't, yet. */
290    
291                            transfer_8to16sub(&dct_codes[4 * 64], cur->u + 8 *j*edged_width/2 + 8*i,
292                                                                    refGMC->u + 8 *j*edged_width/2 + 8*i, edged_width/2);
293    
294                            transfer_8to16sub(&dct_codes[5 * 64], cur->v + 8*j* edged_width/2 + 8*i,
295                                                                    refGMC->v + 8*j* edged_width/2 + 8*i, edged_width/2);
296    
297                            return;
298                    }
299    
300                    /* ordinary compensation */
301    
302                    dx = (quarterpel ? mb->qmvs[0].x : mb->mvs[0].x);
303                    dy = (quarterpel ? mb->qmvs[0].y : mb->mvs[0].y);
304    
305                    if (reduced_resolution) {
306                            dx = RRV_MV_SCALEUP(dx);
307                            dy = RRV_MV_SCALEUP(dy);
308                    }
309    
310                    compensate16x16_interpolate(&dct_codes[0 * 64], cur->y, ref->y, refh->y,
311                                                            refv->y, refhv->y, tmp, 16 * i, 16 * j, dx, dy,
312                                                            edged_width, quarterpel, reduced_resolution, rounding);
313    
314                    if (quarterpel) { dx /= 2; dy /= 2; }
315    
316                    dx = (dx >> 1) + roundtab_79[dx & 0x3];
317                    dy = (dy >> 1) + roundtab_79[dy & 0x3];
318    
319            } else {                                        // mode == MODE_INTER4V
320                    int k, sumx = 0, sumy = 0;
321                    const VECTOR * const mvs = (quarterpel ? mb->qmvs : mb->mvs);
322    
323                    for (k = 0; k < 4; k++) {
324                            dx = mvs[k].x;
325                            dy = mvs[k].y;
326                            sumx += quarterpel ? dx/2 : dx;
327                            sumy += quarterpel ? dy/2 : dy;
328    
329                            if (reduced_resolution){
330                                    dx = RRV_MV_SCALEUP(dx);
331                                    dy = RRV_MV_SCALEUP(dy);
332                            }
333    
334          } else                                          // mode == MODE_INTER4V                          compensate8x8_interpolate(&dct_codes[k * 64], cur->y, ref->y, refh->y,
335                                                                            refv->y, refhv->y, tmp, 16 * i + 8*(k&1), 16 * j + 8*(k>>1), dx,
336                                                                            dy, edged_width, quarterpel, reduced_resolution, rounding);
337                    }
338                    dx = (sumx >> 3) + roundtab_76[sumx & 0xf];
339                    dy = (sumy >> 3) + roundtab_76[sumy & 0xf];
340            }
341    
342            CompensateChroma(dx, dy, i, j, cur, ref, tmp,
343                                            &dct_codes[4 * 64], edged_width / 2, rounding, reduced_resolution);
344    }
345    
346    
347    void
348    MBMotionCompensationBVOP(MBParam * pParam,
349                                                    MACROBLOCK * const mb,
350                                                    const uint32_t i,
351                                                    const uint32_t j,
352                                                    IMAGE * const cur,
353                                                    const IMAGE * const f_ref,
354                                                    const IMAGE * const f_refh,
355                                                    const IMAGE * const f_refv,
356                                                    const IMAGE * const f_refhv,
357                                                    const IMAGE * const b_ref,
358                                                    const IMAGE * const b_refh,
359                                                    const IMAGE * const b_refv,
360                                                    const IMAGE * const b_refhv,
361                                                    int16_t * dct_codes)
362          {          {
363                  int32_t sum, dx, dy;          const uint32_t edged_width = pParam->edged_width;
364            int32_t dx, dy, b_dx, b_dy, sumx, sumy, b_sumx, b_sumy;
365            int k;
366            const int quarterpel = pParam->vol_flags & XVID_VOL_QUARTERPEL;
367            const uint8_t * ptr1, * ptr2;
368            uint8_t * const tmp = f_refv->u;
369            const VECTOR * const fmvs = (quarterpel ? mb->qmvs : mb->mvs);
370            const VECTOR * const bmvs = (quarterpel ? mb->b_qmvs : mb->b_mvs);
371    
372            switch (mb->mode) {
373            case MODE_FORWARD:
374                    dx = fmvs->x; dy = fmvs->y;
375    
376                    compensate16x16_interpolate(&dct_codes[0 * 64], cur->y, f_ref->y, f_refh->y,
377                                                            f_refv->y, f_refhv->y, tmp, 16 * i, 16 * j, dx,
378                                                            dy, edged_width, quarterpel, 0, 0);
379    
380                    if (quarterpel) { dx /= 2; dy /= 2; }
381    
382                    CompensateChroma(       (dx >> 1) + roundtab_79[dx & 0x3],
383                                                            (dy >> 1) + roundtab_79[dy & 0x3],
384                                                            i, j, cur, f_ref, tmp,
385                                                            &dct_codes[4 * 64], edged_width / 2, 0, 0);
386    
387                    return;
388    
389            case MODE_BACKWARD:
390                    b_dx = bmvs->x; b_dy = bmvs->y;
391    
392                    compensate16x16_interpolate(&dct_codes[0 * 64], cur->y, b_ref->y, b_refh->y,
393                                                                                    b_refv->y, b_refhv->y, tmp, 16 * i, 16 * j, b_dx,
394                                                                                    b_dy, edged_width, quarterpel, 0, 0);
395    
396                    if (quarterpel) { b_dx /= 2; b_dy /= 2; }
397    
398                    CompensateChroma(       (b_dx >> 1) + roundtab_79[b_dx & 0x3],
399                                                            (b_dy >> 1) + roundtab_79[b_dy & 0x3],
400                                                            i, j, cur, b_ref, tmp,
401                                                            &dct_codes[4 * 64], edged_width / 2, 0, 0);
402    
403                    return;
404    
405            case MODE_INTERPOLATE: /* _could_ use DIRECT, but would be overkill (no 4MV there) */
406            case MODE_DIRECT_NO4V:
407                    dx = fmvs->x; dy = fmvs->y;
408                    b_dx = bmvs->x; b_dy = bmvs->y;
409    
410                    if (quarterpel) {
411    
412                            if ((dx&3) | (dy&3)) {
413                                    interpolate16x16_quarterpel(tmp - i * 16 - j * 16 * edged_width,
414                                            (uint8_t *) f_ref->y, tmp + 32,
415                                            tmp + 64, tmp + 96, 16*i, 16*j, dx, dy, edged_width, 0);
416                                    ptr1 = tmp;
417                            } else ptr1 = f_ref->y + (16*j + dy/4)*edged_width + 16*i + dx/4; // fullpixel position
418    
419                            if ((b_dx&3) | (b_dy&3)) {
420                                    interpolate16x16_quarterpel(tmp - i * 16 - j * 16 * edged_width + 16,
421                                            (uint8_t *) b_ref->y, tmp + 32,
422                                            tmp + 64, tmp + 96, 16*i, 16*j, b_dx, b_dy, edged_width, 0);
423                                    ptr2 = tmp + 16;
424                            } else ptr2 = b_ref->y + (16*j + b_dy/4)*edged_width + 16*i + b_dx/4; // fullpixel position
425    
426                            b_dx /= 2;
427                            b_dy /= 2;
428                            dx /= 2;
429                            dy /= 2;
430    
431                    } else {
432                            ptr1 = get_ref(f_ref->y, f_refh->y, f_refv->y, f_refhv->y,
433                                                            i, j, 16, dx, dy, edged_width);
434    
435                            ptr2 = get_ref(b_ref->y, b_refh->y, b_refv->y, b_refhv->y,
436                                                            i, j, 16, b_dx, b_dy, edged_width);
437                    }
438                    for (k = 0; k < 4; k++)
439                                    transfer_8to16sub2(&dct_codes[k * 64],
440                                                                            cur->y + (i * 16+(k&1)*8) + (j * 16+((k>>1)*8)) * edged_width,
441                                                                            ptr1 + (k&1)*8 + (k>>1)*8*edged_width,
442                                                                            ptr2 + (k&1)*8 + (k>>1)*8*edged_width, edged_width);
443    
444    
445                    dx = (dx >> 1) + roundtab_79[dx & 0x3];
446                    dy = (dy >> 1) + roundtab_79[dy & 0x3];
447    
448                  compensate8x8_halfpel(&dct_codes[0 * 64], cur->y, ref->y, refh->y,                  b_dx = (b_dx >> 1) + roundtab_79[b_dx & 0x3];
449                                                            refv->y, refhv->y, 16 * i, 16 * j, mb->mvs[0].x,                  b_dy = (b_dy >> 1) + roundtab_79[b_dy & 0x3];
450                                                            mb->mvs[0].y, edged_width);  
451                  compensate8x8_halfpel(&dct_codes[1 * 64], cur->y, ref->y, refh->y,                  break;
452                                                            refv->y, refhv->y, 16 * i + 8, 16 * j,  
453                                                            mb->mvs[1].x, mb->mvs[1].y, edged_width);          default: // MODE_DIRECT (or MODE_DIRECT_NONE_MV in case of bframes decoding)
454                  compensate8x8_halfpel(&dct_codes[2 * 64], cur->y, ref->y, refh->y,                  sumx = sumy = b_sumx = b_sumy = 0;
455                                                            refv->y, refhv->y, 16 * i, 16 * j + 8,  
456                                                            mb->mvs[2].x, mb->mvs[2].y, edged_width);                  for (k = 0; k < 4; k++) {
457                  compensate8x8_halfpel(&dct_codes[3 * 64], cur->y, ref->y, refh->y,  
458                                                            refv->y, refhv->y, 16 * i + 8, 16 * j + 8,                          dx = fmvs[k].x; dy = fmvs[k].y;
459                                                            mb->mvs[3].x, mb->mvs[3].y, edged_width);                          b_dx = bmvs[k].x; b_dy = bmvs[k].y;
460    
461                  sum = mb->mvs[0].x + mb->mvs[1].x + mb->mvs[2].x + mb->mvs[3].x;                          if (quarterpel) {
462                  dx = (sum ? SIGN(sum) *                                  sumx += dx/2; sumy += dy/2;
463                            (roundtab[ABS(sum) % 16] + (ABS(sum) / 16) * 2) : 0);                                  b_sumx += b_dx/2; b_sumy += b_dy/2;
464    
465                  sum = mb->mvs[0].y + mb->mvs[1].y + mb->mvs[2].y + mb->mvs[3].y;                                  if ((dx&3) | (dy&3)) {
466                  dy = (sum ? SIGN(sum) *                                          interpolate8x8_quarterpel(tmp - (i * 16+(k&1)*8) - (j * 16+((k>>1)*8)) * edged_width,
467                            (roundtab[ABS(sum) % 16] + (ABS(sum) / 16) * 2) : 0);                                                  (uint8_t *) f_ref->y,
468                                                    tmp + 32, tmp + 64, tmp + 96,
469                  /* uv-block-based compensation */                                                  16*i + (k&1)*8, 16*j + (k>>1)*8, dx, dy, edged_width, 0);
470                  interpolate8x8_switch(refv->u, ref->u, 8 * i, 8 * j, dx, dy,                                          ptr1 = tmp;
471                                                            edged_width / 2, rounding);                                  } else ptr1 = f_ref->y + (16*j + (k>>1)*8 + dy/4)*edged_width + 16*i + (k&1)*8 + dx/4;
472                  transfer_8to16sub(&dct_codes[4 * 64],  
473                                                    cur->u + 8 * j * edged_width / 2 + 8 * i,                                  if ((b_dx&3) | (b_dy&3)) {
474                                                    refv->u + 8 * j * edged_width / 2 + 8 * i,                                          interpolate8x8_quarterpel(tmp - (i * 16+(k&1)*8) - (j * 16+((k>>1)*8)) * edged_width + 16,
475                                                    (uint8_t *) b_ref->y,
476                                                    tmp + 16, tmp + 32, tmp + 48,
477                                                    16*i + (k&1)*8, 16*j + (k>>1)*8, b_dx, b_dy, edged_width, 0);
478                                            ptr2 = tmp + 16;
479                                    } else ptr2 = b_ref->y + (16*j + (k>>1)*8 + b_dy/4)*edged_width + 16*i + (k&1)*8 + b_dx/4;
480                            } else {
481                                    sumx += dx; sumy += dy;
482                                    b_sumx += b_dx; b_sumy += b_dy;
483    
484                                    ptr1 = get_ref(f_ref->y, f_refh->y, f_refv->y, f_refhv->y,
485                                                                    2*i + (k&1), 2*j + (k>>1), 8, dx, dy, edged_width);
486                                    ptr2 = get_ref(b_ref->y, b_refh->y, b_refv->y, b_refhv->y,
487                                                                    2*i + (k&1), 2*j + (k>>1), 8, b_dx, b_dy,  edged_width);
488                            }
489                            transfer_8to16sub2(&dct_codes[k * 64],
490                                                                    cur->y + (i * 16+(k&1)*8) + (j * 16+((k>>1)*8)) * edged_width,
491                                                                    ptr1, ptr2,     edged_width);
492    
493                    }
494    
495                    dx = (sumx >> 3) + roundtab_76[sumx & 0xf];
496                    dy = (sumy >> 3) + roundtab_76[sumy & 0xf];
497                    b_dx = (b_sumx >> 3) + roundtab_76[b_sumx & 0xf];
498                    b_dy = (b_sumy >> 3) + roundtab_76[b_sumy & 0xf];
499    
500                    break;
501            }
502    
503            // uv block-based chroma interpolation for direct and interpolate modes
504            transfer_8to16sub2(&dct_codes[4 * 64],
505                                                    cur->u + (j * 8) * edged_width / 2 + (i * 8),
506                                                    interpolate8x8_switch2(tmp, b_ref->u, 8 * i, 8 * j,
507                                                                                                    b_dx, b_dy, edged_width / 2, 0),
508                                                    interpolate8x8_switch2(tmp + 8, f_ref->u, 8 * i, 8 * j,
509                                                                                                    dx, dy, edged_width / 2, 0),
510                                                    edged_width / 2);                                                    edged_width / 2);
511    
512                  interpolate8x8_switch(refv->v, ref->v, 8 * i, 8 * j, dx, dy,          transfer_8to16sub2(&dct_codes[5 * 64],
513                                                            edged_width / 2, rounding);                                                  cur->v + (j * 8) * edged_width / 2 + (i * 8),
514                  transfer_8to16sub(&dct_codes[5 * 64],                                                  interpolate8x8_switch2(tmp, b_ref->v, 8 * i, 8 * j,
515                                                    cur->v + 8 * j * edged_width / 2 + 8 * i,                                                                                                  b_dx, b_dy, edged_width / 2, 0),
516                                                    refv->v + 8 * j * edged_width / 2 + 8 * i,                                                  interpolate8x8_switch2(tmp + 8, f_ref->v, 8 * i, 8 * j,
517                                                                                                    dx, dy, edged_width / 2, 0),
518                                                    edged_width / 2);                                                    edged_width / 2);
519          }          }
520    
521    
522    
523    void generate_GMCparameters( const int num_wp, const int res,
524                                                    const WARPPOINTS *const warp,
525                                                    const int width, const int height,
526                                                    GMC_DATA *const gmc)
527    {
528            const int du0 = warp->duv[0].x;
529            const int dv0 = warp->duv[0].y;
530            const int du1 = warp->duv[1].x;
531            const int dv1 = warp->duv[1].y;
532            const int du2 = warp->duv[2].x;
533            const int dv2 = warp->duv[2].y;
534    
535            gmc->W = width;
536            gmc->H = height;
537    
538            gmc->rho = 4 - log2bin(res-1);  // = {3,2,1,0} for res={2,4,8,16}
539    
540            gmc->alpha = log2bin(gmc->W-1);
541            gmc->Ws = (1 << gmc->alpha);
542    
543            gmc->dxF = 16*gmc->Ws + RDIV( 8*gmc->Ws*du1, gmc->W );
544            gmc->dxG = RDIV( 8*gmc->Ws*dv1, gmc->W );
545            gmc->Fo  = (res*du0 + 1) << (gmc->alpha+gmc->rho-1);
546            gmc->Go  = (res*dv0 + 1) << (gmc->alpha+gmc->rho-1);
547    
548            if (num_wp==2) {
549                    gmc->dyF = -gmc->dxG;
550                    gmc->dyG =  gmc->dxF;
551            } else if (num_wp==3) {
552                    gmc->beta = log2bin(gmc->H-1);
553                    gmc->Hs = (1 << gmc->beta);
554                    gmc->dyF =                       RDIV( 8*gmc->Hs*du2, gmc->H );
555                    gmc->dyG = 16*gmc->Hs + RDIV( 8*gmc->Hs*dv2, gmc->H );
556                    if (gmc->beta > gmc->alpha) {
557                            gmc->dxF <<= (gmc->beta - gmc->alpha);
558                            gmc->dxG <<= (gmc->beta - gmc->alpha);
559                            gmc->alpha = gmc->beta;
560                            gmc->Ws = 1<< gmc->beta;
561                    } else {
562                            gmc->dyF <<= gmc->alpha - gmc->beta;
563                            gmc->dyG <<= gmc->alpha - gmc->beta;
564                    }
565            }
566    
567            gmc->cFo = gmc->dxF + gmc->dyF + (1 << (gmc->alpha+gmc->rho+1));
568            gmc->cFo += 16*gmc->Ws*(du0-1);
569    
570            gmc->cGo = gmc->dxG + gmc->dyG + (1 << (gmc->alpha+gmc->rho+1));
571            gmc->cGo += 16*gmc->Ws*(dv0-1);
572    }
573    
574    void
575    generate_GMCimage(      const GMC_DATA *const gmc_data, // [input] precalculated data
576                                            const IMAGE *const pRef,                // [input]
577                                            const int mb_width,
578                                            const int mb_height,
579                                            const int stride,
580                                            const int stride2,
581                                            const int fcode,                                // [input] some parameters...
582                                            const int32_t quarterpel,               // [input] for rounding avgMV
583                                            const int reduced_resolution,   // [input] ignored
584                                            const int32_t rounding,                 // [input] for rounding image data
585                                            MACROBLOCK *const pMBs,                 // [output] average motion vectors
586                                            IMAGE *const pGMC)                              // [output] full warped image
587    {
588    
589            unsigned int mj,mi;
590            VECTOR avgMV;
591    
592            for (mj = 0; mj < (unsigned int)mb_height; mj++)
593                    for (mi = 0; mi < (unsigned int)mb_width; mi++) {
594    
595                            avgMV = generate_GMCimageMB(gmc_data, pRef, mi, mj,
596                                                    stride, stride2, quarterpel, rounding, pGMC);
597    
598                            pMBs[mj*mb_width+mi].amv.x = gmc_sanitize(avgMV.x, quarterpel, fcode);
599                            pMBs[mj*mb_width+mi].amv.y = gmc_sanitize(avgMV.y, quarterpel, fcode);
600                            pMBs[mj*mb_width+mi].mcsel = 0; /* until mode decision */
601            }
602    }
603    
604    
605    
606    #define MLT(i)  (((16-(i))<<16) + (i))
607    static const uint32_t MTab[16] = {
608      MLT( 0), MLT( 1), MLT( 2), MLT( 3), MLT( 4), MLT( 5), MLT( 6), MLT(7),
609      MLT( 8), MLT( 9), MLT(10), MLT(11), MLT(12), MLT(13), MLT(14), MLT(15)
610    };
611    #undef MLT
612    
613    VECTOR generate_GMCimageMB( const GMC_DATA *const gmc_data,
614                                                            const IMAGE *const pRef,
615                                                            const int mi, const int mj,
616                                                            const int stride,
617                                                            const int stride2,
618                                                            const int quarterpel,
619                                                            const int rounding,
620                                                            IMAGE *const pGMC)
621    {
622            const int W = gmc_data->W;
623            const int H = gmc_data->H;
624    
625            const int rho = gmc_data->rho;
626            const int alpha = gmc_data->alpha;
627    
628            const int rounder = ( 128 - (rounding<<(rho+rho)) ) << 16;
629    
630            const int dxF = gmc_data->dxF;
631            const int dyF = gmc_data->dyF;
632            const int dxG = gmc_data->dxG;
633            const int dyG = gmc_data->dyG;
634    
635            uint8_t *dstY, *dstU, *dstV;
636    
637            int I,J;
638            VECTOR avgMV = {0,0};
639    
640            int32_t Fj, Gj;
641    
642            dstY = &pGMC->y[(mj*16)*stride+mi*16] + 16;
643    
644            Fj = gmc_data->Fo + dyF*mj*16 + dxF*mi*16;
645            Gj = gmc_data->Go + dyG*mj*16 + dxG*mi*16;
646    
647            for (J = 16; J > 0; --J) {
648                    int32_t Fi, Gi;
649    
650                    Fi = Fj; Fj += dyF;
651                    Gi = Gj; Gj += dyG;
652                    for (I = -16; I < 0; ++I) {
653                            int32_t F, G;
654                            uint32_t ri, rj;
655    
656                            F = ( Fi >> (alpha+rho) ) << rho; Fi += dxF;
657                            G = ( Gi >> (alpha+rho) ) << rho; Gi += dxG;
658    
659                            avgMV.x += F;
660                            avgMV.y += G;
661    
662                            ri = MTab[F&15];
663                            rj = MTab[G&15];
664    
665                            F >>= 4;
666                            G >>= 4;
667    
668                            if (F < -1) F = -1;
669                            else if (F > W) F = W;
670                            if (G< -1) G=-1;
671                            else if (G>H) G=H;
672    
673                            {        // MMX-like bilinear...
674                                    const int offset = G*stride + F;
675                                    uint32_t f0, f1;
676                                    f0 = pRef->y[ offset +0 ];
677                                    f0 |= pRef->y[ offset +1 ] << 16;
678                                    f1 = pRef->y[ offset+stride +0 ];
679                                    f1 |= pRef->y[ offset+stride +1 ] << 16;
680                                    f0 = (ri*f0)>>16;
681                                    f1 = (ri*f1) & 0x0fff0000;
682                                    f0 |= f1;
683                                    f0 = ( rj*f0 + rounder ) >> 24;
684    
685                                    dstY[I] = (uint8_t)f0;
686                            }
687  }  }
688    
689                    dstY += stride;
690            }
691    
692            dstU = &pGMC->u[(mj*8)*stride2+mi*8] + 8;
693            dstV = &pGMC->v[(mj*8)*stride2+mi*8] + 8;
694    
695            Fj = gmc_data->cFo + dyF*4 *mj*8 + dxF*4 *mi*8;
696            Gj = gmc_data->cGo + dyG*4 *mj*8 + dxG*4 *mi*8;
697    
698            for (J = 8; J > 0; --J) {
699                    int32_t Fi, Gi;
700                    Fi = Fj; Fj += 4*dyF;
701                    Gi = Gj; Gj += 4*dyG;
702    
703                    for (I = -8; I < 0; ++I) {
704                            int32_t F, G;
705                            uint32_t ri, rj;
706    
707                            F = ( Fi >> (alpha+rho+2) ) << rho; Fi += 4*dxF;
708                            G = ( Gi >> (alpha+rho+2) ) << rho; Gi += 4*dxG;
709    
710                            ri = MTab[F&15];
711                            rj = MTab[G&15];
712    
713                            F >>= 4;
714                            G >>= 4;
715    
716                            if (F < -1) F=-1;
717                            else if (F >= W/2) F = W/2;
718                            if (G < -1) G = -1;
719                            else if (G >= H/2) G = H/2;
720    
721                            {
722                                    const int offset = G*stride2 + F;
723                                    uint32_t f0, f1;
724    
725                                    f0      = pRef->u[ offset                +0 ];
726                                    f0 |= pRef->u[ offset            +1 ] << 16;
727                                    f1      = pRef->u[ offset+stride2 +0 ];
728                                    f1 |= pRef->u[ offset+stride2 +1 ] << 16;
729                                    f0 = (ri*f0)>>16;
730                                    f1 = (ri*f1) & 0x0fff0000;
731                                    f0 |= f1;
732                                    f0 = ( rj*f0 + rounder ) >> 24;
733    
734                                    dstU[I] = (uint8_t)f0;
735    
736    
737                                    f0      = pRef->v[ offset                +0 ];
738                                    f0 |= pRef->v[ offset            +1 ] << 16;
739                                    f1      = pRef->v[ offset+stride2 +0 ];
740                                    f1 |= pRef->v[ offset+stride2 +1 ] << 16;
741                                    f0 = (ri*f0)>>16;
742                                    f1 = (ri*f1) & 0x0fff0000;
743                                    f0 |= f1;
744                                    f0 = ( rj*f0 + rounder ) >> 24;
745    
746                                    dstV[I] = (uint8_t)f0;
747                            }
748                    }
749                    dstU += stride2;
750                    dstV += stride2;
751            }
752    
753    
754            avgMV.x -= 16*((256*mi+120)<<4);        // 120 = 15*16/2
755            avgMV.y -= 16*((256*mj+120)<<4);
756    
757            avgMV.x = RSHIFT( avgMV.x, (4+7-quarterpel) );
758            avgMV.y = RSHIFT( avgMV.y, (4+7-quarterpel) );
759    
760            return avgMV;
761    }
762    
763    
764    
765    #ifdef OLD_GRUEL_GMC
766    void
767    generate_GMCparameters( const int num_wp,                       // [input]: number of warppoints
768                                                    const int res,                  // [input]: resolution
769                                                    const WARPPOINTS *const warp, // [input]: warp points
770                                                    const int width, const int height,
771                                                    GMC_DATA *const gmc)    // [output] precalculated parameters
772    {
773    
774    /* We follow mainly two sources: The original standard, which is ugly, and the
775       thesis from Andreas Dehnhardt, which is much nicer.
776    
777            Notation is: indices are written next to the variable,
778                                     primes in the standard are denoted by a suffix 'p'.
779            types are   "c"=constant, "i"=input parameter, "f"=calculated, then fixed,
780                                    "o"=output data, " "=other, "u" = unused, "p"=calc for every pixel
781    
782    type | variable name  |   ISO name (TeX-style) |  value or range  |  usage
783    -------------------------------------------------------------------------------------
784     c   | H                          |   H                                    |  [16 , ?]            |  image width (w/o edges)
785     c   | W                          |   W                                    |  [16 , ?]            |  image height (w/o edges)
786    
787     c   | i0                         |   i_0                                  |  0                           |  ref. point #1, X
788     c   | j0                         |   j_0                                  |  0                           |  ref. point #1, Y
789     c   | i1                         |   i_1                                  |  W                           |  ref. point #2, X
790     c   | j1                         |   j_1                                  |  0                           |  ref. point #2, Y
791     cu  | i2                         |   i_2                                  |  0                           |  ref. point #3, X
792     cu  | i2                         |   j_2                                  |  H                           |  ref. point #3, Y
793    
794     i   | du0                |   du[0]                        |  [-16863,16863]  |  warp vector #1, Y
795     i   | dv0                |   dv[0]                        |  [-16863,16863]  |  warp vector #1, Y
796     i   | du1                |   du[1]                        |  [-16863,16863]  |  warp vector #2, Y
797     i   | dv1                |   dv[1]                        |  [-16863,16863]  |  warp vector #2, Y
798     iu  | du2                |   du[2]                        |  [-16863,16863]  |  warp vector #3, Y
799     iu  | dv2                |   dv[2]                        |  [-16863,16863]  |  warp vector #3, Y
800    
801     i   | s                          |   s                                 |  {2,4,8,16}     |  interpol. resolution
802     f   | sigma              |             -                          |  log2(s)            |  X / s == X >> sigma
803     f   | r                          |   r                                 |  =16/s                   |  complementary res.
804     f   | rho                      |   \rho                                 |  log2(r)              |  X / r == X >> rho
805    
806     f   | i0s                      |   i'_0                                 |                                |
807     f   | j0s                      |   j'_0                                 |                                |
808     f       | i1s                  |   i'_1                                 |                                |
809     f       | j1s                  |   j'_1                                 |                                |
810     f       | i2s                  |   i'_2                                 |                                |
811     f       | j2s                  |   j'_2                                 |                                |
812    
813     f   | alpha              |   \alpha                       |                              |  2^{alpha-1} < W <= 2^alpha
814     f   | beta                |   \beta                            |                                 |  2^{beta-1} < H <= 2^beta
815    
816     f   | Ws                        |   W'                            | W = 2^{alpha}      |  scaled width
817     f   | Hs                        |   H'                            | W = 2^{beta}        |  scaled height
818    
819     f   | i1ss                |   i''_1                            |  "virtual sprite stuff"
820     f   | j1ss                |   j''_1                            |  "virtual sprite stuff"
821     f   | i2ss                |   i''_2                            |  "virtual sprite stuff"
822     f   | j2ss                |   j''_2                            |  "virtual sprite stuff"
823    */
824    
825    /* Some calculations are disabled because we only use 2 warppoints at the moment */
826    
827            int du0 = warp->duv[0].x;
828            int dv0 = warp->duv[0].y;
829            int du1 = warp->duv[1].x;
830            int dv1 = warp->duv[1].y;
831    //      int du2 = warp->duv[2].x;
832    //      int dv2 = warp->duv[2].y;
833    
834            gmc->num_wp = num_wp;
835    
836            gmc->s = res;                                           /* scaling parameters 2,4,8 or 16 */
837            gmc->sigma = log2bin(res-1);    /* log2bin(15)=4, log2bin(16)=5, log2bin(17)=5  */
838            gmc->r = 16/res;
839            gmc->rho = 4 - gmc->sigma;              /* = log2bin(r-1) */
840    
841            gmc->W = width;
842            gmc->H = height;                        /* fixed reference coordinates */
843    
844            gmc->alpha = log2bin(gmc->W-1);
845            gmc->Ws= 1<<gmc->alpha;
846    
847    //      gmc->beta = log2bin(gmc->H-1);
848    //      gmc->Hs= 1<<gmc->beta;
849    
850    //      printf("du0=%d dv0=%d du1=%d dv1=%d s=%d sigma=%d W=%d alpha=%d, Ws=%d, rho=%d\n",du0,dv0,du1,dv1,gmc->s,gmc->sigma,gmc->W,gmc->alpha,gmc->Ws,gmc->rho);
851    
852            /* i2s is only needed for num_wp >= 3, etc.  */
853            /* the 's' values are in 1/s pel resolution */
854            gmc->i0s = res/2 * ( du0 );
855            gmc->j0s = res/2 * ( dv0 );
856            gmc->i1s = res/2 * (2*width + du1 + du0 );
857            gmc->j1s = res/2 * ( dv1 + dv0 );
858    //      gmc->i2s = res/2 * ( du2 + du0 );
859    //      gmc->j2s = res/2 * (2*height + dv2 + dv0 );
860    
861            /* i2s and i2ss are only needed for num_wp == 3, etc.  */
862    
863            /* the 'ss' values are in 1/16 pel resolution */
864            gmc->i1ss = 16*gmc->Ws + ROUNDED_DIV(((gmc->W-gmc->Ws)*(gmc->r*gmc->i0s) + gmc->Ws*(gmc->r*gmc->i1s - 16*gmc->W)),gmc->W);
865            gmc->j1ss = ROUNDED_DIV( ((gmc->W - gmc->Ws)*(gmc->r*gmc->j0s) + gmc->Ws*gmc->r*gmc->j1s) ,gmc->W );
866    
867    //      gmc->i2ss = ROUNDED_DIV( ((gmc->H - gmc->Hs)*(gmc->r*gmc->i0s) + gmc->Hs*(gmc->r*gmc->i2s)), gmc->H);
868    //      gmc->j2ss = 16*gmc->Hs + ROUNDED_DIV( ((gmc->H-gmc->Hs)*(gmc->r*gmc->j0s) + gmc->Ws*(gmc->r*gmc->j2s - 16*gmc->H)), gmc->H);
869    
870            return;
871    }
872    
873    void
874    generate_GMCimage(      const GMC_DATA *const gmc_data,         // [input] precalculated data
875                                            const IMAGE *const pRef,                        // [input]
876                                            const int mb_width,
877                                            const int mb_height,
878                                            const int stride,
879                                            const int stride2,
880                                            const int fcode,                                        // [input] some parameters...
881                                            const int32_t quarterpel,                       // [input] for rounding avgMV
882                                            const int reduced_resolution,           // [input] ignored
883                                            const int32_t rounding,                 // [input] for rounding image data
884                                            MACROBLOCK *const pMBs,         // [output] average motion vectors
885                                            IMAGE *const pGMC)                      // [output] full warped image
886    {
887    
888            unsigned int mj,mi;
889            VECTOR avgMV;
890    
891            for (mj = 0;mj < mb_height; mj++)
892            for (mi = 0;mi < mb_width; mi++) {
893    
894                    avgMV = generate_GMCimageMB(gmc_data, pRef, mi, mj,
895                                            stride, stride2, quarterpel, rounding, pGMC);
896    
897                    pMBs[mj*mb_width+mi].amv.x = gmc_sanitize(avgMV.x, quarterpel, fcode);
898                    pMBs[mj*mb_width+mi].amv.y = gmc_sanitize(avgMV.y, quarterpel, fcode);
899                    pMBs[mj*mb_width+mi].mcsel = 0; /* until mode decision */
900            }
901    }
902    
903    
904    VECTOR generate_GMCimageMB(     const GMC_DATA *const gmc_data, /* [input] all precalc data */
905                                                            const IMAGE *const pRef,                /* [input] */
906                                                            const int mi, const int mj,     /* [input] MB position  */
907                                                            const int stride,                               /* [input] Lumi stride */
908                                                            const int stride2,                              /* [input] chroma stride */
909                                                            const int quarterpel,                   /* [input] for rounding of avgMV */
910                                                            const int rounding,                     /* [input] for rounding of imgae data */
911                                                            IMAGE *const pGMC)                              /* [outut] generate image */
912    
913    /*
914    type | variable name  |   ISO name (TeX-style) |  value or range  |  usage
915    -------------------------------------------------------------------------------------
916     p   | F                          |   F(i,j)                       |                              | pelwise motion vector X in s-th pel
917     p   | G                          |   G(i,j)                       |                              | pelwise motion vector Y in s-th pel
918     p   | Fc                        |   F_c(i,j)                    |                                |
919     p   | Gc                        |   G_c(i,j)                    |                                | same for chroma
920    
921     p   | Y00                      |   Y_{00}                         |  [0,255*s*s]        | first: 4 neighbouring Y-values
922     p   | Y01                      |   Y_{01}                         |  [0,255]            | at fullpel position, around the
923     p   | Y10                      |   Y_{10}                         |  [0,255*s]    | position where pelweise MV points to
924     p   | Y11                      |   Y_{11}                         |  [0,255]            | later: bilinear interpol Y-values in Y00
925    
926     p   | C00                      |   C_{00}                         |  [0,255*s*s]        | same for chroma Cb and Cr
927     p   | C01                      |   C_{01}                         |  [0,255]            |
928     p   | C10                      |   C_{10}                         |  [0,255*s]    |
929     p   | C11                      |   C_{11}                         |  [0,255]            |
930    
931    */
932    {
933            const int W = gmc_data->W;
934            const int H = gmc_data->H;
935    
936            const int s = gmc_data->s;
937            const int sigma = gmc_data->sigma;
938    
939            const int r = gmc_data->r;
940            const int rho = gmc_data->rho;
941    
942            const int i0s = gmc_data->i0s;
943            const int j0s = gmc_data->j0s;
944    
945            const int i1ss = gmc_data->i1ss;
946            const int j1ss = gmc_data->j1ss;
947    //      const int i2ss = gmc_data->i2ss;
948    //      const int j2ss = gmc_data->j2ss;
949    
950            const int alpha = gmc_data->alpha;
951            const int Ws    = gmc_data->Ws;
952    
953    //      const int beta  = gmc_data->beta;
954    //      const int Hs    = gmc_data->Hs;
955    
956            int I,J;
957            VECTOR avgMV = {0,0};
958    
959            for (J=16*mj;J<16*(mj+1);J++)
960            for (I=16*mi;I<16*(mi+1);I++)
961            {
962                    int F= i0s + ( ((-r*i0s+i1ss)*I + (r*j0s-j1ss)*J + (1<<(alpha+rho-1))) >>  (alpha+rho) );
963                    int G= j0s + ( ((-r*j0s+j1ss)*I + (-r*i0s+i1ss)*J + (1<<(alpha+rho-1))) >> (alpha+rho) );
964    
965    /* this naive implementation (with lots of multiplications) isn't slower (rather faster) than
966       working incremental. Don't ask me why... maybe the whole this is memory bound? */
967    
968                    const int ri= F & (s-1); // fractional part of pelwise MV X
969                    const int rj= G & (s-1); // fractional part of pelwise MV Y
970    
971                    int Y00,Y01,Y10,Y11;
972    
973    /* unclipped values are used for avgMV */
974                    avgMV.x += F-(I<<sigma);                /* shift position to 1/s-pel, as the MV is */
975                    avgMV.y += G-(J<<sigma);                /* TODO: don't do this (of course) */
976    
977                    F >>= sigma;
978                    G >>= sigma;
979    
980    /* clip values to be in range. Since we have edges, clip to 1 less than lower boundary
981       this way positions F+1/G+1 are still right */
982    
983                    if (F< -1)
984                            F=-1;
985                    else if (F>W)
986                            F=W;    /* W or W-1 doesn't matter, so save 1 subtract ;-) */
987                    if (G< -1)
988                            G=-1;
989                    else if (G>H)
990                            G=H;            /* dito */
991    
992                    Y00 = pRef->y[ G*stride + F ];                          // Lumi values
993                    Y01 = pRef->y[ G*stride + F+1 ];
994                    Y10 = pRef->y[ G*stride + F+stride ];
995                    Y11 = pRef->y[ G*stride + F+stride+1 ];
996    
997                    /* bilinear interpolation */
998                    Y00 = ((s-ri)*Y00 + ri*Y01);
999                    Y10 = ((s-ri)*Y10 + ri*Y11);
1000                    Y00 = ((s-rj)*Y00 + rj*Y10 + s*s/2 - rounding ) >> (sigma+sigma);
1001    
1002                    pGMC->y[J*stride+I] = (uint8_t)Y00;                                                                             /* output 1 Y-pixel */
1003            }
1004    
1005    
1006    /* doing chroma _here_ is even more stupid and slow, because won't be used until Compensation and
1007            most likely not even then (only if the block really _is_ GMC)
1008    */
1009    
1010            for (J=8*mj;J<8*(mj+1);J++)             /* this plays the role of j_c,i_c in the standard */
1011            for (I=8*mi;I<8*(mi+1);I++)             /* For I_c we have to use I_c = 4*i_c+1 ! */
1012            {
1013                    /* same positions for both chroma components, U=Cb and V=Cr */
1014                    int Fc=((-r*i0s+i1ss)*(4*I+1) + (r*j0s-j1ss)*(4*J+1) +2*Ws*r*i0s
1015                                                    -16*Ws +(1<<(alpha+rho+1)))>>(alpha+rho+2);
1016                    int Gc=((-r*j0s+j1ss)*(4*I+1) +(-r*i0s+i1ss)*(4*J+1) +2*Ws*r*j0s
1017                                                    -16*Ws +(1<<(alpha+rho+1))) >>(alpha+rho+2);
1018    
1019                    const int ri= Fc & (s-1); // fractional part of pelwise MV X
1020                    const int rj= Gc & (s-1); // fractional part of pelwise MV Y
1021    
1022                    int C00,C01,C10,C11;
1023    
1024                    Fc >>= sigma;
1025                    Gc >>= sigma;
1026    
1027                    if (Fc< -1)
1028                            Fc=-1;
1029                    else if (Fc>=W/2)
1030                            Fc=W/2;         /* W or W-1 doesn't matter, so save 1 subtraction ;-) */
1031                    if (Gc< -1)
1032                            Gc=-1;
1033                    else if (Gc>=H/2)
1034                            Gc=H/2;         /* dito */
1035    
1036    /* now calculate U data */
1037                    C00 = pRef->u[ Gc*stride2 + Fc ];                               // chroma-value Cb
1038                    C01 = pRef->u[ Gc*stride2 + Fc+1 ];
1039                    C10 = pRef->u[ (Gc+1)*stride2 + Fc ];
1040                    C11 = pRef->u[ (Gc+1)*stride2 + Fc+1 ];
1041    
1042                    /* bilinear interpolation */
1043                    C00 = ((s-ri)*C00 + ri*C01);
1044                    C10 = ((s-ri)*C10 + ri*C11);
1045                    C00 = ((s-rj)*C00 + rj*C10 + s*s/2 - rounding ) >> (sigma+sigma);
1046    
1047                    pGMC->u[J*stride2+I] = (uint8_t)C00;                                                                            /* output 1 U-pixel */
1048    
1049    /* now calculate V data */
1050                    C00 = pRef->v[ Gc*stride2 + Fc ];                               // chroma-value Cr
1051                    C01 = pRef->v[ Gc*stride2 + Fc+1 ];
1052                    C10 = pRef->v[ (Gc+1)*stride2 + Fc ];
1053                    C11 = pRef->v[ (Gc+1)*stride2 + Fc+1 ];
1054    
1055                    /* bilinear interpolation */
1056                    C00 = ((s-ri)*C00 + ri*C01);
1057                    C10 = ((s-ri)*C10 + ri*C11);
1058                    C00 = ((s-rj)*C00 + rj*C10 + s*s/2 - rounding ) >> (sigma+sigma);
1059    
1060                    pGMC->v[J*stride2+I] = (uint8_t)C00;                                                                            /* output 1 V-pixel */
1061            }
1062    
1063    /* The average vector is rounded from 1/s-pel to 1/2 or 1/4 using the '//' operator*/
1064    
1065            avgMV.x = RSHIFT( avgMV.x, (sigma+7-quarterpel) );
1066            avgMV.y = RSHIFT( avgMV.y, (sigma+7-quarterpel) );
1067    
1068            /* ^^^^ this is the way MS Reference Software does it */
1069    
1070            return avgMV;   /* clipping to fcode area is done outside! */
1071    }
1072    
1073    #endif

Legend:
Removed from v.437  
changed lines
  Added in v.982

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