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

Diff of /branches/dev-api-3/xvidcore/src/motion/motion_est.c

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

revision 768, Thu Jan 9 11:36:33 2003 UTC revision 818, Fri Feb 7 15:02:56 2003 UTC
# Line 32  Line 32 
32  #include <stdio.h>  #include <stdio.h>
33  #include <stdlib.h>  #include <stdlib.h>
34  #include <string.h>     // memcpy  #include <string.h>     // memcpy
35    #include <math.h>       // lrint
36    
37  #include "../encoder.h"  #include "../encoder.h"
38  #include "../utils/mbfunctions.h"  #include "../utils/mbfunctions.h"
# Line 56  Line 57 
57  d_mv_bits(int x, int y, const VECTOR pred, const uint32_t iFcode, const int qpel, const int rrv)  d_mv_bits(int x, int y, const VECTOR pred, const uint32_t iFcode, const int qpel, const int rrv)
58  {  {
59          int xb, yb;          int xb, yb;
60          x += x * qpel; y += y * qpel;          x = qpel ? x<<1 : x;
61            y = qpel ? y<<1 : y;
62          if (rrv) { x = RRV_MV_SCALEDOWN(x); y = RRV_MV_SCALEDOWN(y); }          if (rrv) { x = RRV_MV_SCALEDOWN(x); y = RRV_MV_SCALEDOWN(y); }
63    
64          x -= pred.x;          x -= pred.x;
65          y -= pred.y;          y -= pred.y;
66    
# Line 79  Line 82 
82          return xb + yb;          return xb + yb;
83  }  }
84    
85    static int32_t ChromaSAD2(int fx, int fy, int bx, int by, const SearchData * const data)
86    {
87            int sad;
88            const uint32_t stride = data->iEdgedWidth/2;
89            uint8_t * f_refu = data->RefQ,
90                    * f_refv = data->RefQ + 8,
91                    * b_refu = data->RefQ + 16,
92                    * b_refv = data->RefQ + 24;
93    
94            switch (((fx & 1) << 1) | (fy & 1))     {
95                    case 0:
96                            fx = fx / 2; fy = fy / 2;
97                            f_refu = (uint8_t*)data->RefCU + fy * stride + fx, stride;
98                            f_refv = (uint8_t*)data->RefCV + fy * stride + fx, stride;
99                            break;
100                    case 1:
101                            fx = fx / 2; fy = (fy - 1) / 2;
102                            interpolate8x8_halfpel_v(f_refu, data->RefCU + fy * stride + fx, stride, data->rounding);
103                            interpolate8x8_halfpel_v(f_refv, data->RefCV + fy * stride + fx, stride, data->rounding);
104                            break;
105                    case 2:
106                            fx = (fx - 1) / 2; fy = fy / 2;
107                            interpolate8x8_halfpel_h(f_refu, data->RefCU + fy * stride + fx, stride, data->rounding);
108                            interpolate8x8_halfpel_h(f_refv, data->RefCV + fy * stride + fx, stride, data->rounding);
109                            break;
110                    default:
111                            fx = (fx - 1) / 2; fy = (fy - 1) / 2;
112                            interpolate8x8_halfpel_hv(f_refu, data->RefCU + fy * stride + fx, stride, data->rounding);
113                            interpolate8x8_halfpel_hv(f_refv, data->RefCV + fy * stride + fx, stride, data->rounding);
114                            break;
115            }
116    
117            switch (((bx & 1) << 1) | (by & 1))     {
118                    case 0:
119                            bx = bx / 2; by = by / 2;
120                            b_refu = (uint8_t*)data->b_RefCU + by * stride + bx, stride;
121                            b_refv = (uint8_t*)data->b_RefCV + by * stride + bx, stride;
122                            break;
123                    case 1:
124                            bx = bx / 2; by = (by - 1) / 2;
125                            interpolate8x8_halfpel_v(b_refu, data->b_RefCU + by * stride + bx, stride, data->rounding);
126                            interpolate8x8_halfpel_v(b_refv, data->b_RefCV + by * stride + bx, stride, data->rounding);
127                            break;
128                    case 2:
129                            bx = (bx - 1) / 2; by = by / 2;
130                            interpolate8x8_halfpel_h(b_refu, data->b_RefCU + by * stride + bx, stride, data->rounding);
131                            interpolate8x8_halfpel_h(b_refv, data->b_RefCV + by * stride + bx, stride, data->rounding);
132                            break;
133                    default:
134                            bx = (bx - 1) / 2; by = (by - 1) / 2;
135                            interpolate8x8_halfpel_hv(b_refu, data->b_RefCU + by * stride + bx, stride, data->rounding);
136                            interpolate8x8_halfpel_hv(b_refv, data->b_RefCV + by * stride + bx, stride, data->rounding);
137                            break;
138            }
139    
140            sad = sad8bi(data->CurU, b_refu, f_refu, stride);
141            sad += sad8bi(data->CurV, b_refv, f_refv, stride);
142    
143            return sad;
144    }
145    
146    
147  static int32_t  static int32_t
148  ChromaSAD(int dx, int dy, const SearchData * const data)  ChromaSAD(int dx, int dy, const SearchData * const data)
149  {  {
# Line 139  Line 204 
204  {  {
205          switch ( ((x&1)<<1) | (y&1) ) {          switch ( ((x&1)<<1) | (y&1) ) {
206                  case 0 : return data->Ref + x/2 + (y/2)*(data->iEdgedWidth);                  case 0 : return data->Ref + x/2 + (y/2)*(data->iEdgedWidth);
207                    case 3 : return data->RefHV + (x-1)/2 + ((y-1)/2)*(data->iEdgedWidth);
208                  case 1 : return data->RefV + x/2 + ((y-1)/2)*(data->iEdgedWidth);                  case 1 : return data->RefV + x/2 + ((y-1)/2)*(data->iEdgedWidth);
209                  case 2 : return data->RefH + (x-1)/2 + (y/2)*(data->iEdgedWidth);                  default : return data->RefH + (x-1)/2 + (y/2)*(data->iEdgedWidth);      //case 2
                 default : return data->RefHV + (x-1)/2 + ((y-1)/2)*(data->iEdgedWidth);  
210          }          }
211  }  }
212    
# Line 202  Line 267 
267    
268          ref1 = GetReferenceB(halfpel_x, halfpel_y, dir, data);          ref1 = GetReferenceB(halfpel_x, halfpel_y, dir, data);
269          switch( ((x&1)<<1) + (y&1) ) {          switch( ((x&1)<<1) + (y&1) ) {
270          case 0: // pure halfpel position          case 3: // x and y in qpel resolution - the "corners" (top left/right and
271                  return (uint8_t *) ref1;                           // bottom left/right) during qpel refinement
272                    ref2 = GetReferenceB(halfpel_x, y - halfpel_y, dir, data);
273                    ref3 = GetReferenceB(x - halfpel_x, halfpel_y, dir, data);
274                    ref4 = GetReferenceB(x - halfpel_x, y - halfpel_y, dir, data);
275                    interpolate8x8_avg4(Reference, ref1, ref2, ref3, ref4, iEdgedWidth, rounding);
276                    interpolate8x8_avg4(Reference+8, ref1+8, ref2+8, ref3+8, ref4+8, iEdgedWidth, rounding);
277                    interpolate8x8_avg4(Reference+8*iEdgedWidth, ref1+8*iEdgedWidth, ref2+8*iEdgedWidth, ref3+8*iEdgedWidth, ref4+8*iEdgedWidth, iEdgedWidth, rounding);
278                    interpolate8x8_avg4(Reference+8*iEdgedWidth+8, ref1+8*iEdgedWidth+8, ref2+8*iEdgedWidth+8, ref3+8*iEdgedWidth+8, ref4+8*iEdgedWidth+8, iEdgedWidth, rounding);
279                    break;
280    
281          case 1: // x halfpel, y qpel - top or bottom during qpel refinement          case 1: // x halfpel, y qpel - top or bottom during qpel refinement
282                  ref2 = GetReferenceB(halfpel_x, y - halfpel_y, dir, data);                  ref2 = GetReferenceB(halfpel_x, y - halfpel_y, dir, data);
283                  interpolate8x8_avg2(Reference, ref1, ref2, iEdgedWidth, rounding, 8);                  interpolate8x8_avg2(Reference, ref1, ref2, iEdgedWidth, rounding, 8);
# Line 220  Line 294 
294                  interpolate8x8_avg2(Reference+8*iEdgedWidth+8, ref1+8*iEdgedWidth+8, ref2+8*iEdgedWidth+8, iEdgedWidth, rounding, 8);                  interpolate8x8_avg2(Reference+8*iEdgedWidth+8, ref1+8*iEdgedWidth+8, ref2+8*iEdgedWidth+8, iEdgedWidth, rounding, 8);
295                  break;                  break;
296    
297          default: // x and y in qpel resolution - the "corners" (top left/right and          case 0: // pure halfpel position
298                           // bottom left/right) during qpel refinement                  return (uint8_t *) ref1;
                 ref2 = GetReferenceB(halfpel_x, y - halfpel_y, dir, data);  
                 ref3 = GetReferenceB(x - halfpel_x, halfpel_y, dir, data);  
                 ref4 = GetReferenceB(x - halfpel_x, y - halfpel_y, dir, data);  
                 interpolate8x8_avg4(Reference, ref1, ref2, ref3, ref4, iEdgedWidth, rounding);  
                 interpolate8x8_avg4(Reference+8, ref1+8, ref2+8, ref3+8, ref4+8, iEdgedWidth, rounding);  
                 interpolate8x8_avg4(Reference+8*iEdgedWidth, ref1+8*iEdgedWidth, ref2+8*iEdgedWidth, ref3+8*iEdgedWidth, ref4+8*iEdgedWidth, iEdgedWidth, rounding);  
                 interpolate8x8_avg4(Reference+8*iEdgedWidth+8, ref1+8*iEdgedWidth+8, ref2+8*iEdgedWidth+8, ref3+8*iEdgedWidth+8, ref4+8*iEdgedWidth+8, iEdgedWidth, rounding);  
                 break;  
299          }          }
300          return Reference;          return Reference;
301  }  }
# Line 239  Line 305 
305  static void  static void
306  CheckCandidate16(const int x, const int y, const int Direction, int * const dir, const SearchData * const data)  CheckCandidate16(const int x, const int y, const int Direction, int * const dir, const SearchData * const data)
307  {  {
308          int t, xc, yc;          int xc, yc;
309          const uint8_t * Reference;          const uint8_t * Reference;
310          VECTOR * current;          VECTOR * current;
311            int32_t sad; uint32_t t;
312    
313          if ( (x > data->max_dx) | (x < data->min_dx)          if ( (x > data->max_dx) || (x < data->min_dx)
314                  | (y > data->max_dy) | (y < data->min_dy) ) return;                  || (y > data->max_dy) || (y < data->min_dy) ) return;
315    
316          if (data->qpel_precision) { // x and y are in 1/4 precision          if (!data->qpel_precision) {
                 Reference = Interpolate16x16qpel(x, y, 0, data);  
                 xc = x/2; yc = y/2; //for chroma sad  
                 current = data->currentQMV;  
         } else {  
317                  Reference = GetReference(x, y, data);                  Reference = GetReference(x, y, data);
318                  current = data->currentMV;                  current = data->currentMV;
319                  xc = x; yc = y;                  xc = x; yc = y;
320            } else { // x and y are in 1/4 precision
321                    Reference = Interpolate16x16qpel(x, y, 0, data);
322                    xc = x/2; yc = y/2; //for chroma sad
323                    current = data->currentQMV;
324          }          }
         t = d_mv_bits(x, y, data->predMV, data->iFcode, data->qpel^data->qpel_precision, 0);  
325    
326          data->temp[0] = sad16v(data->Cur, Reference, data->iEdgedWidth, data->temp + 1);          sad = sad16v(data->Cur, Reference, data->iEdgedWidth, data->temp + 1);
327            t = d_mv_bits(x, y, data->predMV, data->iFcode, data->qpel^data->qpel_precision, 0);
328    
329          data->temp[0] += (data->lambda16 * t * data->temp[0])>>10;          sad += (data->lambda16 * t * sad)>>10;
330          data->temp[1] += (data->lambda8 * t * (data->temp[1] + NEIGH_8X8_BIAS))>>10;          data->temp[1] += (data->lambda8 * t * (data->temp[1] + NEIGH_8X8_BIAS))>>10;
331    
332          if (data->chroma) data->temp[0] += ChromaSAD((xc >> 1) + roundtab_79[xc & 0x3],          if (data->chroma) sad += ChromaSAD((xc >> 1) + roundtab_79[xc & 0x3],
333                                                                                                          (yc >> 1) + roundtab_79[yc & 0x3], data);                                                                                                          (yc >> 1) + roundtab_79[yc & 0x3], data);
334    
335          if (data->temp[0] < data->iMinSAD[0]) {          if (sad < data->iMinSAD[0]) {
336                  data->iMinSAD[0] = data->temp[0];                  data->iMinSAD[0] = sad;
337                  current[0].x = x; current[0].y = y;                  current[0].x = x; current[0].y = y;
338                  *dir = Direction; }                  *dir = Direction;
339            }
340    
341          if (data->temp[1] < data->iMinSAD[1]) {          if (data->temp[1] < data->iMinSAD[1]) {
342                  data->iMinSAD[1] = data->temp[1]; current[1].x = x; current[1].y= y; }                  data->iMinSAD[1] = data->temp[1]; current[1].x = x; current[1].y= y; }
# Line 282  Line 350 
350  }  }
351    
352  static void  static void
353    CheckCandidate8(const int x, const int y, const int Direction, int * const dir, const SearchData * const data)
354    {
355            int32_t sad; uint32_t t;
356            const uint8_t * Reference;
357    
358            if ( (x > data->max_dx) || (x < data->min_dx)
359                    || (y > data->max_dy) || (y < data->min_dy) ) return;
360    
361            if (!data->qpel_precision) Reference = GetReference(x, y, data);
362            else Reference = Interpolate8x8qpel(x, y, 0, 0, data);
363    
364            sad = sad8(data->Cur, Reference, data->iEdgedWidth);
365            t = d_mv_bits(x, y, data->predMV, data->iFcode, data->qpel^data->qpel_precision, 0);
366    
367            sad += (data->lambda8 * t * (sad+NEIGH_8X8_BIAS))>>10;
368    
369            if (sad < *(data->iMinSAD)) {
370                    *(data->iMinSAD) = sad;
371                    data->currentMV->x = x; data->currentMV->y = y;
372                    *dir = Direction;
373            }
374    }
375    
376    
377    static void
378  CheckCandidate32(const int x, const int y, const int Direction, int * const dir, const SearchData * const data)  CheckCandidate32(const int x, const int y, const int Direction, int * const dir, const SearchData * const data)
379  {  {
380          uint32_t t;          uint32_t t;
381          const uint8_t * Reference;          const uint8_t * Reference;
382    
383          if ( (!(x&1) && x !=0) | (!(y&1) && y !=0) || //non-zero integer value          if ( (!(x&1) && x !=0) || (!(y&1) && y !=0) || //non-zero integer value
384                  (x > data->max_dx) | (x < data->min_dx)                  (x > data->max_dx) || (x < data->min_dx)
385                  | (y > data->max_dy) | (y < data->min_dy) ) return;                  || (y > data->max_dy) || (y < data->min_dy) ) return;
386    
387          Reference = GetReference(x, y, data);          Reference = GetReference(x, y, data);
388          t = d_mv_bits(x, y, data->predMV, data->iFcode, 0, 1);          t = d_mv_bits(x, y, data->predMV, data->iFcode, 0, 1);
# Line 317  Line 410 
410  static void  static void
411  CheckCandidate16no4v(const int x, const int y, const int Direction, int * const dir, const SearchData * const data)  CheckCandidate16no4v(const int x, const int y, const int Direction, int * const dir, const SearchData * const data)
412  {  {
413          int32_t sad;          int32_t sad, xc, yc;
414          const uint8_t * Reference;          const uint8_t * Reference;
415          uint32_t t;          uint32_t t;
416          VECTOR * current;          VECTOR * current;
# Line 330  Line 423 
423          if (data->qpel_precision) { // x and y are in 1/4 precision          if (data->qpel_precision) { // x and y are in 1/4 precision
424                  Reference = Interpolate16x16qpel(x, y, 0, data);                  Reference = Interpolate16x16qpel(x, y, 0, data);
425                  current = data->currentQMV;                  current = data->currentQMV;
426                    xc = x/2; yc = y/2;
427          } else {          } else {
428                  Reference = GetReference(x, y, data);                  Reference = GetReference(x, y, data);
429                  current = data->currentMV;                  current = data->currentMV;
430                    xc = x; yc = y;
431          }          }
432          t = d_mv_bits(x, y, data->predMV, data->iFcode,          t = d_mv_bits(x, y, data->predMV, data->iFcode,
433                                          data->qpel^data->qpel_precision, data->rrv);                                          data->qpel^data->qpel_precision, data->rrv);
# Line 340  Line 435 
435          sad = sad16(data->Cur, Reference, data->iEdgedWidth, 256*4096);          sad = sad16(data->Cur, Reference, data->iEdgedWidth, 256*4096);
436          sad += (data->lambda16 * t * sad)>>10;          sad += (data->lambda16 * t * sad)>>10;
437    
438            if (data->chroma) sad += ChromaSAD((xc >> 1) + roundtab_79[xc & 0x3],
439                                                                                    (yc >> 1) + roundtab_79[yc & 0x3], data);
440    
441    
442          if (sad < *(data->iMinSAD)) {          if (sad < *(data->iMinSAD)) {
443                  *(data->iMinSAD) = sad;                  *(data->iMinSAD) = sad;
444                  current->x = x; current->y = y;                  current->x = x; current->y = y;
445                  *dir = Direction; }                  *dir = Direction;
446            }
447  }  }
448    
449  static void  static void
450  CheckCandidate32I(const int x, const int y, const int Direction, int * const dir, const SearchData * const data)  CheckCandidate32I(const int x, const int y, const int Direction, int * const dir, const SearchData * const data)
451  {  {
452  // maximum speed - for P/B/I decision  // maximum speed - for P/B/I decision
453            int32_t sad;
454    
455          if ( (x > data->max_dx) | (x < data->min_dx)          if ( (x > data->max_dx) || (x < data->min_dx)
456                  | (y > data->max_dy) | (y < data->min_dy) ) return;                  || (y > data->max_dy) || (y < data->min_dy) ) return;
457    
458          data->temp[0] = sad32v_c(data->Cur, data->Ref + x/2 + (y/2)*(data->iEdgedWidth),          sad = sad32v_c(data->Cur, data->Ref + x/2 + (y/2)*(data->iEdgedWidth),
459                                                          data->iEdgedWidth, data->temp+1);                                                          data->iEdgedWidth, data->temp+1);
460    
461          if (data->temp[0] < *(data->iMinSAD)) {          if (sad < *(data->iMinSAD)) {
462                  *(data->iMinSAD) = data->temp[0];                  *(data->iMinSAD) = sad;
463                  data->currentMV[0].x = x; data->currentMV[0].y = y;                  data->currentMV[0].x = x; data->currentMV[0].y = y;
464                  *dir = Direction; }                  *dir = Direction;
465            }
466          if (data->temp[1] < data->iMinSAD[1]) {          if (data->temp[1] < data->iMinSAD[1]) {
467                  data->iMinSAD[1] = data->temp[1]; data->currentMV[1].x = x; data->currentMV[1].y = y; }                  data->iMinSAD[1] = data->temp[1]; data->currentMV[1].x = x; data->currentMV[1].y = y; }
468          if (data->temp[2] < data->iMinSAD[2]) {          if (data->temp[2] < data->iMinSAD[2]) {
# Line 375  Line 477 
477  static void  static void
478  CheckCandidateInt(const int xf, const int yf, const int Direction, int * const dir, const SearchData * const data)  CheckCandidateInt(const int xf, const int yf, const int Direction, int * const dir, const SearchData * const data)
479  {  {
480          int32_t sad, xb, yb;          int32_t sad, xb, yb, xcf, ycf, xcb, ycb;
481          uint32_t t;          uint32_t t;
482          const uint8_t *ReferenceF, *ReferenceB;          const uint8_t *ReferenceF, *ReferenceB;
483          VECTOR *current;          VECTOR *current;
# Line 388  Line 490 
490                  xb = data->currentMV[1].x; yb = data->currentMV[1].y;                  xb = data->currentMV[1].x; yb = data->currentMV[1].y;
491                  ReferenceB = GetReferenceB(xb, yb, 1, data);                  ReferenceB = GetReferenceB(xb, yb, 1, data);
492                  current = data->currentMV;                  current = data->currentMV;
493                    xcf = xf; ycf = yf;
494                    xcb = xb; ycb = yb;
495          } else {          } else {
496                  ReferenceF = Interpolate16x16qpel(xf, yf, 0, data);                  ReferenceF = Interpolate16x16qpel(xf, yf, 0, data);
497                  xb = data->currentQMV[1].x; yb = data->currentQMV[1].y;                  xb = data->currentQMV[1].x; yb = data->currentQMV[1].y;
498                  current = data->currentQMV;                  current = data->currentQMV;
499                  ReferenceB = Interpolate16x16qpel(xb, yb, 1, data);                  ReferenceB = Interpolate16x16qpel(xb, yb, 1, data);
500                    xcf = xf/2; ycf = yf/2;
501                    xcb = xb/2; ycb = yb/2;
502          }          }
503    
504          t = d_mv_bits(xf, yf, data->predMV, data->iFcode, data->qpel^data->qpel_precision, 0)          t = d_mv_bits(xf, yf, data->predMV, data->iFcode, data->qpel^data->qpel_precision, 0)
# Line 401  Line 507 
507          sad = sad16bi(data->Cur, ReferenceF, ReferenceB, data->iEdgedWidth);          sad = sad16bi(data->Cur, ReferenceF, ReferenceB, data->iEdgedWidth);
508          sad += (data->lambda16 * t * sad)>>10;          sad += (data->lambda16 * t * sad)>>10;
509    
510            if (data->chroma) sad += ChromaSAD2((xcf >> 1) + roundtab_79[xcf & 0x3],
511                                                                                    (ycf >> 1) + roundtab_79[ycf & 0x3],
512                                                                                    (xcb >> 1) + roundtab_79[xcb & 0x3],
513                                                                                    (ycb >> 1) + roundtab_79[ycb & 0x3], data);
514    
515          if (sad < *(data->iMinSAD)) {          if (sad < *(data->iMinSAD)) {
516                  *(data->iMinSAD) = sad;                  *(data->iMinSAD) = sad;
517                  current->x = xf; current->y = yf;                  current->x = xf; current->y = yf;
518                  *dir = Direction; }                  *dir = Direction;
519            }
520  }  }
521    
522  static void  static void
523  CheckCandidateDirect(const int x, const int y, const int Direction, int * const dir, const SearchData * const data)  CheckCandidateDirect(const int x, const int y, const int Direction, int * const dir, const SearchData * const data)
524  {  {
525          int32_t sad = 0;          int32_t sad = 0, xcf = 0, ycf = 0, xcb = 0, ycb = 0;
526          uint32_t k;          uint32_t k;
527          const uint8_t *ReferenceF;          const uint8_t *ReferenceF;
528          const uint8_t *ReferenceB;          const uint8_t *ReferenceB;
# Line 429  Line 541 
541                          data->directmvB[k].y                          data->directmvB[k].y
542                          : mvs.y - data->referencemv[k].y);                          : mvs.y - data->referencemv[k].y);
543    
544                  if (( mvs.x > data->max_dx ) || ( mvs.x < data->min_dx )                  if ( (mvs.x > data->max_dx) | (mvs.x < data->min_dx)
545                          || ( mvs.y > data->max_dy ) || ( mvs.y < data->min_dy )                          | (mvs.y > data->max_dy) | (mvs.y < data->min_dy)
546                          || ( b_mvs.x > data->max_dx ) || ( b_mvs.x < data->min_dx )                          | (b_mvs.x > data->max_dx) | (b_mvs.x < data->min_dx)
547                          || ( b_mvs.y > data->max_dy ) || ( b_mvs.y < data->min_dy )) return;                          | (b_mvs.y > data->max_dy) | (b_mvs.y < data->min_dy) ) return;
548    
549                    if (data->qpel) {
550                  mvs.x *= 2 - data->qpel; mvs.y *= 2 - data->qpel;                          xcf += mvs.x/2; ycf += mvs.y/2;
551                  b_mvs.x *= 2 - data->qpel; b_mvs.y *= 2 - data->qpel; //we move to qpel precision anyway                          xcb += b_mvs.x/2; ycb += b_mvs.y/2;
552                    } else {
553                            xcf += mvs.x; ycf += mvs.y;
554                            xcb += b_mvs.x; ycb += b_mvs.y;
555                            mvs.x *= 2; mvs.y *= 2; //we move to qpel precision anyway
556                            b_mvs.x *= 2; b_mvs.y *= 2;
557                    }
558    
559                  ReferenceF = Interpolate8x8qpel(mvs.x, mvs.y, k, 0, data);                  ReferenceF = Interpolate8x8qpel(mvs.x, mvs.y, k, 0, data);
560                  ReferenceB = Interpolate8x8qpel(b_mvs.x, b_mvs.y, k, 1, data);                  ReferenceB = Interpolate8x8qpel(b_mvs.x, b_mvs.y, k, 1, data);
# Line 448  Line 566 
566    
567          sad += (data->lambda16 * d_mv_bits(x, y, zeroMV, 1, 0, 0) * sad)>>10;          sad += (data->lambda16 * d_mv_bits(x, y, zeroMV, 1, 0, 0) * sad)>>10;
568    
569            if (data->chroma) sad += ChromaSAD2((xcf >> 3) + roundtab_76[xcf & 0xf],
570                                                                                    (ycf >> 3) + roundtab_76[ycf & 0xf],
571                                                                                    (xcb >> 3) + roundtab_76[xcb & 0xf],
572                                                                                    (ycb >> 3) + roundtab_76[ycb & 0xf], data);
573    
574          if (sad < *(data->iMinSAD)) {          if (sad < *(data->iMinSAD)) {
575                  *(data->iMinSAD) = sad;                  *(data->iMinSAD) = sad;
576                  data->currentMV->x = x; data->currentMV->y = y;                  data->currentMV->x = x; data->currentMV->y = y;
577                  *dir = Direction; }                  *dir = Direction;
578            }
579  }  }
580    
581  static void  static void
582  CheckCandidateDirectno4v(const int x, const int y, const int Direction, int * const dir, const SearchData * const data)  CheckCandidateDirectno4v(const int x, const int y, const int Direction, int * const dir, const SearchData * const data)
583  {  {
584          int32_t sad;          int32_t sad, xcf, ycf, xcb, ycb;
585          const uint8_t *ReferenceF;          const uint8_t *ReferenceF;
586          const uint8_t *ReferenceB;          const uint8_t *ReferenceB;
587          VECTOR mvs, b_mvs;          VECTOR mvs, b_mvs;
# Line 474  Line 598 
598                  data->directmvB[0].y                  data->directmvB[0].y
599                  : mvs.y - data->referencemv[0].y);                  : mvs.y - data->referencemv[0].y);
600    
601          if (( mvs.x > data->max_dx ) || ( mvs.x < data->min_dx )          if ( (mvs.x > data->max_dx) | (mvs.x < data->min_dx)
602                  || ( mvs.y > data->max_dy ) || ( mvs.y < data->min_dy )                  | (mvs.y > data->max_dy) | (mvs.y < data->min_dy)
603                  || ( b_mvs.x > data->max_dx ) || ( b_mvs.x < data->min_dx )                  | (b_mvs.x > data->max_dx) | (b_mvs.x < data->min_dx)
604                  || ( b_mvs.y > data->max_dy ) || ( b_mvs.y < data->min_dy )) return;                  | (b_mvs.y > data->max_dy) | (b_mvs.y < data->min_dy) ) return;
605    
606          mvs.x *= 2 - data->qpel; mvs.y *= 2 - data->qpel;          if (data->qpel) {
607          b_mvs.x *= 2 - data->qpel; b_mvs.y *= 2 - data->qpel; //we move to qpel precision anyway                  xcf = 4*(mvs.x/2); ycf = 4*(mvs.y/2);
608                    xcb = 4*(b_mvs.x/2); ycb = 4*(b_mvs.y/2);
609          ReferenceF = Interpolate16x16qpel(mvs.x, mvs.y, 0, data);          ReferenceF = Interpolate16x16qpel(mvs.x, mvs.y, 0, data);
610          ReferenceB = Interpolate16x16qpel(b_mvs.x, b_mvs.y, 1, data);          ReferenceB = Interpolate16x16qpel(b_mvs.x, b_mvs.y, 1, data);
611            } else {
612                    xcf = 4*mvs.x; ycf = 4*mvs.y;
613                    xcb = 4*b_mvs.x; ycb = 4*b_mvs.y;
614                    ReferenceF = GetReference(mvs.x, mvs.y, data);
615                    ReferenceB = GetReferenceB(b_mvs.x, b_mvs.y, 1, data);
616            }
617    
618          sad = sad16bi(data->Cur, ReferenceF, ReferenceB, data->iEdgedWidth);          sad = sad16bi(data->Cur, ReferenceF, ReferenceB, data->iEdgedWidth);
619          sad += (data->lambda16 * d_mv_bits(x, y, zeroMV, 1, 0, 0) * sad)>>10;          sad += (data->lambda16 * d_mv_bits(x, y, zeroMV, 1, 0, 0) * sad)>>10;
620    
621          if (sad < *(data->iMinSAD)) {          if (data->chroma) sad += ChromaSAD2((xcf >> 3) + roundtab_76[xcf & 0xf],
622                  *(data->iMinSAD) = sad;                                                                                  (ycf >> 3) + roundtab_76[ycf & 0xf],
623                  data->currentMV->x = x; data->currentMV->y = y;                                                                                  (xcb >> 3) + roundtab_76[xcb & 0xf],
624                  *dir = Direction; }                                                                                  (ycb >> 3) + roundtab_76[ycb & 0xf], data);
 }  
   
 static void  
 CheckCandidate8(const int x, const int y, const int Direction, int * const dir, const SearchData * const data)  
 {  
         int32_t sad; uint32_t t;  
         const uint8_t * Reference;  
   
         if ( (x > data->max_dx) | (x < data->min_dx)  
                 | (y > data->max_dy) | (y < data->min_dy) ) return;  
   
         if (data->qpel) Reference = Interpolate16x16qpel(x, y, 0, data);  
         else Reference = GetReference(x, y, data);  
   
         sad = sad8(data->Cur, Reference, data->iEdgedWidth);  
         t = d_mv_bits(x, y, data->predMV, data->iFcode, data->qpel && !data->qpel_precision, 0);  
   
         sad += (data->lambda8 * t * (sad+NEIGH_8X8_BIAS))>>10;  
625    
626          if (sad < *(data->iMinSAD)) {          if (sad < *(data->iMinSAD)) {
627                  *(data->iMinSAD) = sad;                  *(data->iMinSAD) = sad;
628                  data->currentMV->x = x; data->currentMV->y = y;                  data->currentMV->x = x; data->currentMV->y = y;
629                  *dir = Direction; }                  *dir = Direction;
630            }
631  }  }
632    
633  /* CHECK_CANDIATE FUNCTIONS END */  /* CHECK_CANDIATE FUNCTIONS END */
# Line 752  Line 865 
865          Data.iFcode = current->fcode;          Data.iFcode = current->fcode;
866          Data.rounding = pParam->m_rounding_type;          Data.rounding = pParam->m_rounding_type;
867          Data.qpel = pParam->m_quarterpel;          Data.qpel = pParam->m_quarterpel;
868          Data.chroma = current->global_flags & XVID_ME_COLOUR;          Data.chroma = current->motion_flags & PMV_CHROMA16;
869          Data.rrv = current->global_flags & XVID_REDUCED;          Data.rrv = current->global_flags & XVID_REDUCED;
870    
871          if ((current->global_flags & XVID_REDUCED)) {          if ((current->global_flags & XVID_REDUCED)) {
# Line 801  Line 914 
914    
915  //initial skip decision  //initial skip decision
916  /* no early skip for GMC (global vector = skip vector is unknown!)      */  /* no early skip for GMC (global vector = skip vector is unknown!)      */
917                          if (current->coding_type == P_VOP)      { /* no fast SKIP for S(GMC)-VOPs */                          if (!(current->global_flags & XVID_GMC))        { /* no fast SKIP for S(GMC)-VOPs */
918                                  if (pMB->dquant == NO_CHANGE && sad00 < pMB->quant * INITIAL_SKIP_THRESH * (Data.rrv ? 4:1) )                                  if (pMB->dquant == NO_CHANGE && sad00 < pMB->quant * INITIAL_SKIP_THRESH * (Data.rrv ? 4:1) )
919                                          if (Data.chroma || SkipDecisionP(pCurrent, pRef, x, y, iEdgedWidth/2, pMB->quant, Data.rrv)) {                                          if (Data.chroma || SkipDecisionP(pCurrent, pRef, x, y, iEdgedWidth/2, pMB->quant, Data.rrv)) {
920                                                  SkipMacroblockP(pMB, sad00);                                                  SkipMacroblockP(pMB, sad00);
# Line 815  Line 928 
928                                                  current->global_flags & XVID_INTER4V, pMB);                                                  current->global_flags & XVID_INTER4V, pMB);
929    
930  /* final skip decision, a.k.a. "the vector you found, really that good?" */  /* final skip decision, a.k.a. "the vector you found, really that good?" */
931                          if (current->coding_type == P_VOP)      {                          if (!(current->global_flags & XVID_GMC))        {
932                                  if ( (pMB->dquant == NO_CHANGE) && (sad00 < pMB->quant * MAX_SAD00_FOR_SKIP)                                  if ( (pMB->dquant == NO_CHANGE) && (sad00 < pMB->quant * MAX_SAD00_FOR_SKIP)
933                                          && ((100*pMB->sad16)/(sad00+1) > FINAL_SKIP_THRESH * (Data.rrv ? 4:1)) )                                          && ((100*pMB->sad16)/(sad00+1) > FINAL_SKIP_THRESH * (Data.rrv ? 4:1)) )
934                                          if (Data.chroma || SkipDecisionP(pCurrent, pRef, x, y, iEdgedWidth/2, pMB->quant, Data.rrv)) {                                          if (Data.chroma || SkipDecisionP(pCurrent, pRef, x, y, iEdgedWidth/2, pMB->quant, Data.rrv)) {
# Line 855  Line 968 
968                  }                  }
969          }          }
970    
971          if (current->coding_type == S_VOP)      /* first GMC step only for S(GMC)-VOPs */          if (current->global_flags & XVID_GMC )  /* GMC only for S(GMC)-VOPs */
972                  current->GMC_MV = GlobalMotionEst( pMBs, pParam, current->fcode );          {
973          else                  current->warp = GlobalMotionEst( pMBs, pParam, current, reference, pRefH, pRefV, pRefHV);
974                  current->GMC_MV = zeroMV;          }
975    
976          return 0;          return 0;
977  }  }
# Line 1269  Line 1382 
1382    
1383  /* search backward or forward */  /* search backward or forward */
1384  static void  static void
1385  SearchBF(       const uint8_t * const pRef,  SearchBF(       const IMAGE * const pRef,
1386                          const uint8_t * const pRefH,                          const uint8_t * const pRefH,
1387                          const uint8_t * const pRefV,                          const uint8_t * const pRefV,
1388                          const uint8_t * const pRefHV,                          const uint8_t * const pRefHV,
# Line 1291  Line 1404 
1404          *Data->iMinSAD = MV_MAX_ERROR;          *Data->iMinSAD = MV_MAX_ERROR;
1405          Data->iFcode = iFcode;          Data->iFcode = iFcode;
1406          Data->qpel_precision = 0;          Data->qpel_precision = 0;
1407            Data->temp[5] = Data->temp[6] = Data->temp[7] = 256*4096; // reset chroma-sad cache
1408    
1409          Data->Ref = pRef + (x + y * Data->iEdgedWidth) * 16;          Data->Ref = pRef->y + (x + y * Data->iEdgedWidth) * 16;
1410          Data->RefH = pRefH + (x + y * Data->iEdgedWidth) * 16;          Data->RefH = pRefH + (x + y * Data->iEdgedWidth) * 16;
1411          Data->RefV = pRefV + (x + y * Data->iEdgedWidth) * 16;          Data->RefV = pRefV + (x + y * Data->iEdgedWidth) * 16;
1412          Data->RefHV = pRefHV + (x + y * Data->iEdgedWidth) * 16;          Data->RefHV = pRefHV + (x + y * Data->iEdgedWidth) * 16;
1413            Data->RefCU = pRef->u + (x + y * Data->iEdgedWidth/2) * 8;
1414            Data->RefCV = pRef->v + (x + y * Data->iEdgedWidth/2) * 8;
1415    
1416          Data->predMV = *predMV;          Data->predMV = *predMV;
1417    
# Line 1435  Line 1551 
1551          Data->bRefH = b_RefH + k;          Data->bRefH = b_RefH + k;
1552          Data->bRefV = b_RefV + k;          Data->bRefV = b_RefV + k;
1553          Data->bRefHV = b_RefHV + k;          Data->bRefHV = b_RefHV + k;
1554            Data->RefCU = f_Ref->u + (x + (Data->iEdgedWidth/2) * y) * 8;
1555            Data->RefCV = f_Ref->v + (x + (Data->iEdgedWidth/2) * y) * 8;
1556            Data->b_RefCU = b_Ref->u + (x + (Data->iEdgedWidth/2) * y) * 8;
1557            Data->b_RefCV = b_Ref->v + (x + (Data->iEdgedWidth/2) * y) * 8;
1558    
1559          k = Data->qpel ? 4 : 2;          k = Data->qpel ? 4 : 2;
1560          Data->max_dx = k * (pParam->width - x * 16);          Data->max_dx = k * (pParam->width - x * 16);
# Line 1473  Line 1593 
1593          (*CheckCandidate)(0, 0, 255, &k, Data);          (*CheckCandidate)(0, 0, 255, &k, Data);
1594    
1595  // initial (fast) skip decision  // initial (fast) skip decision
1596          if (*Data->iMinSAD < pMB->quant * INITIAL_SKIP_THRESH * 2) {          if (*Data->iMinSAD < pMB->quant * INITIAL_SKIP_THRESH * (2 + Data->chroma?1:0)) {
1597                  //possible skip - checking chroma                  //possible skip
1598                    if (Data->chroma) {
1599                            pMB->mode = MODE_DIRECT_NONE_MV;
1600                            return *Data->iMinSAD; // skip.
1601                    } else {
1602                  SkipDecisionB(pCur, f_Ref, b_Ref, pMB, x, y, Data);                  SkipDecisionB(pCur, f_Ref, b_Ref, pMB, x, y, Data);
1603                  if (pMB->mode == MODE_DIRECT_NONE_MV) return *Data->iMinSAD; // skip.                  if (pMB->mode == MODE_DIRECT_NONE_MV) return *Data->iMinSAD; // skip.
1604          }          }
1605            }
1606    
1607          skip_sad = *Data->iMinSAD;          skip_sad = *Data->iMinSAD;
1608    
# Line 1527  Line 1652 
1652  }  }
1653    
1654  static void  static void
1655  SearchInterpolate(const uint8_t * const f_Ref,  SearchInterpolate(const IMAGE * const f_Ref,
1656                                  const uint8_t * const f_RefH,                                  const uint8_t * const f_RefH,
1657                                  const uint8_t * const f_RefV,                                  const uint8_t * const f_RefV,
1658                                  const uint8_t * const f_RefHV,                                  const uint8_t * const f_RefHV,
1659                                  const uint8_t * const b_Ref,                                  const IMAGE * const b_Ref,
1660                                  const uint8_t * const b_RefH,                                  const uint8_t * const b_RefH,
1661                                  const uint8_t * const b_RefV,                                  const uint8_t * const b_RefV,
1662                                  const uint8_t * const b_RefHV,                                  const uint8_t * const b_RefHV,
# Line 1559  Line 1684 
1684          fData->iFcode = bData.bFcode = fcode; fData->bFcode = bData.iFcode = bcode;          fData->iFcode = bData.bFcode = fcode; fData->bFcode = bData.iFcode = bcode;
1685    
1686          i = (x + y * fData->iEdgedWidth) * 16;          i = (x + y * fData->iEdgedWidth) * 16;
1687          bData.bRef = fData->Ref = f_Ref + i;          bData.bRef = fData->Ref = f_Ref->y + i;
1688          bData.bRefH = fData->RefH = f_RefH + i;          bData.bRefH = fData->RefH = f_RefH + i;
1689          bData.bRefV = fData->RefV = f_RefV + i;          bData.bRefV = fData->RefV = f_RefV + i;
1690          bData.bRefHV = fData->RefHV = f_RefHV + i;          bData.bRefHV = fData->RefHV = f_RefHV + i;
1691          bData.Ref = fData->bRef = b_Ref + i;          bData.Ref = fData->bRef = b_Ref->y + i;
1692          bData.RefH = fData->bRefH = b_RefH + i;          bData.RefH = fData->bRefH = b_RefH + i;
1693          bData.RefV = fData->bRefV = b_RefV + i;          bData.RefV = fData->bRefV = b_RefV + i;
1694          bData.RefHV = fData->bRefHV = b_RefHV + i;          bData.RefHV = fData->bRefHV = b_RefHV + i;
1695            bData.b_RefCU = fData->RefCU = f_Ref->u + (x + (fData->iEdgedWidth/2) * y) * 8;
1696            bData.b_RefCV = fData->RefCV = f_Ref->v + (x + (fData->iEdgedWidth/2) * y) * 8;
1697            bData.RefCU = fData->b_RefCU = b_Ref->u + (x + (fData->iEdgedWidth/2) * y) * 8;
1698            bData.RefCV = fData->b_RefCV = b_Ref->v + (x + (fData->iEdgedWidth/2) * y) * 8;
1699    
1700    
1701          bData.bpredMV = fData->predMV = *f_predMV;          bData.bpredMV = fData->predMV = *f_predMV;
1702          fData->bpredMV = bData.predMV = *b_predMV;          fData->bpredMV = bData.predMV = *b_predMV;
# Line 1683  Line 1813 
1813          int32_t iMinSAD;          int32_t iMinSAD;
1814          VECTOR currentMV[3];          VECTOR currentMV[3];
1815          VECTOR currentQMV[3];          VECTOR currentQMV[3];
1816            int32_t temp[8];
1817          memset(&Data, 0, sizeof(SearchData));          memset(&Data, 0, sizeof(SearchData));
1818          Data.iEdgedWidth = pParam->edged_width;          Data.iEdgedWidth = pParam->edged_width;
1819          Data.currentMV = currentMV; Data.currentQMV = currentQMV;          Data.currentMV = currentMV; Data.currentQMV = currentQMV;
1820          Data.iMinSAD = &iMinSAD;          Data.iMinSAD = &iMinSAD;
1821          Data.lambda16 = lambda_vec16[frame->quant];          Data.lambda16 = lambda_vec16[frame->quant];
         Data.chroma = frame->quant;  
1822          Data.qpel = pParam->m_quarterpel;          Data.qpel = pParam->m_quarterpel;
1823          Data.rounding = 0;          Data.rounding = 0;
1824            Data.chroma = frame->motion_flags & PMV_CHROMA8;
1825            Data.temp = temp;
1826    
1827          Data.RefQ = f_refV->u; // a good place, also used in MC (for similar purpose)          Data.RefQ = f_refV->u; // a good place, also used in MC (for similar purpose)
1828          // note: i==horizontal, j==vertical          // note: i==horizontal, j==vertical
# Line 1710  Line 1842 
1842                                  }                                  }
1843    
1844                          Data.Cur = frame->image.y + (j * Data.iEdgedWidth + i) * 16;                          Data.Cur = frame->image.y + (j * Data.iEdgedWidth + i) * 16;
1845                            Data.CurU = frame->image.u + (j * Data.iEdgedWidth/2 + i) * 8;
1846                            Data.CurV = frame->image.v + (j * Data.iEdgedWidth/2 + i) * 8;
1847                          pMB->quant = frame->quant;                          pMB->quant = frame->quant;
1848    
1849  /* direct search comes first, because it (1) checks for SKIP-mode  /* direct search comes first, because it (1) checks for SKIP-mode
# Line 1728  Line 1862 
1862                          if (pMB->mode == MODE_DIRECT_NONE_MV) { n_count++; continue; }                          if (pMB->mode == MODE_DIRECT_NONE_MV) { n_count++; continue; }
1863    
1864                          // forward search                          // forward search
1865                          SearchBF(f_ref->y, f_refH->y, f_refV->y, f_refHV->y,                          SearchBF(f_ref, f_refH->y, f_refV->y, f_refHV->y,
1866                                                  &frame->image, i, j,                                                  &frame->image, i, j,
1867                                                  frame->motion_flags,                                                  frame->motion_flags,
1868                                                  frame->fcode, pParam,                                                  frame->fcode, pParam,
# Line 1736  Line 1870 
1870                                                  MODE_FORWARD, &Data);                                                  MODE_FORWARD, &Data);
1871    
1872                          // backward search                          // backward search
1873                          SearchBF(b_ref->y, b_refH->y, b_refV->y, b_refHV->y,                          SearchBF(b_ref, b_refH->y, b_refV->y, b_refHV->y,
1874                                                  &frame->image, i, j,                                                  &frame->image, i, j,
1875                                                  frame->motion_flags,                                                  frame->motion_flags,
1876                                                  frame->bcode, pParam,                                                  frame->bcode, pParam,
# Line 1744  Line 1878 
1878                                                  MODE_BACKWARD, &Data);                                                  MODE_BACKWARD, &Data);
1879    
1880                          // interpolate search comes last, because it uses data from forward and backward as prediction                          // interpolate search comes last, because it uses data from forward and backward as prediction
1881                          SearchInterpolate(f_ref->y, f_refH->y, f_refV->y, f_refHV->y,                          SearchInterpolate(f_ref, f_refH->y, f_refV->y, f_refHV->y,
1882                                                  b_ref->y, b_refH->y, b_refV->y, b_refHV->y,                                                  b_ref, b_refH->y, b_refV->y, b_refHV->y,
1883                                                  &frame->image,                                                  &frame->image,
1884                                                  i, j,                                                  i, j,
1885                                                  frame->fcode, frame->bcode,                                                  frame->fcode, frame->bcode,
# Line 1906  Line 2040 
2040                  }                  }
2041          }          }
2042          sSAD /= (pParam->mb_height-2)*(pParam->mb_width-2);          sSAD /= (pParam->mb_height-2)*(pParam->mb_width-2);
2043          if (sSAD > IntraThresh + INTRA_BIAS) return I_VOP;  //      if (sSAD > IntraThresh + INTRA_BIAS) return I_VOP;
2044          if (sSAD > InterThresh ) return P_VOP;          if (sSAD > InterThresh ) return P_VOP;
2045          emms();          emms();
2046          return B_VOP;          return B_VOP;
2047    
2048  }  }
2049    
2050  static void  
2051  CheckGMC(int x, int y, const int dir, int * iDirection,  static WARPPOINTS
2052                  const MACROBLOCK * const pMBs, uint32_t * bestcount, VECTOR * GMC,  GlobalMotionEst(const MACROBLOCK * const pMBs,
2053                  const MBParam * const pParam)                                  const MBParam * const pParam,
2054                                    const FRAMEINFO * const current,
2055                                    const FRAMEINFO * const reference,
2056                                    const IMAGE * const pRefH,
2057                                    const IMAGE * const pRefV,
2058                                    const IMAGE * const pRefHV      )
2059  {  {
         uint32_t mx, my, a, count = 0;  
2060    
2061          for (my = 1; my < pParam->mb_height-1; my++)          const int deltax=8;             // upper bound for difference between a MV and it's neighbour MVs
2062                  for (mx = 1; mx < pParam->mb_width-1; mx++) {          const int deltay=8;
2063                          VECTOR mv;          const int grad=512;             // lower bound for deviation in MB
2064                          const MACROBLOCK *pMB = &pMBs[mx + my * pParam->mb_width];  
2065                          if (pMB->mode == MODE_INTRA || pMB->mode == MODE_NOT_CODED) continue;          WARPPOINTS gmc;
2066                          mv = pMB->mvs[0];  
2067                          a = ABS(mv.x - x) + ABS(mv.y - y);          uint32_t mx, my;
2068                          if (a < 6) count += 6 - a;  
2069                  }          int MBh = pParam->mb_height;
2070            int MBw = pParam->mb_width;
2071    
2072            int *MBmask= calloc(MBh*MBw,sizeof(int));
2073            double DtimesF[4] = { 0.,0., 0., 0. };
2074            double sol[4] = { 0., 0., 0., 0. };
2075            double a,b,c,n,denom;
2076            double meanx,meany;
2077            int num,oldnum;
2078    
2079            if (!MBmask) { fprintf(stderr,"Mem error\n"); return gmc;}
2080    
2081    // filter mask of all blocks
2082    
2083            for (my = 1; my < MBh-1; my++)
2084            for (mx = 1; mx < MBw-1; mx++)
2085            {
2086                    const int mbnum = mx + my * MBw;
2087                    const MACROBLOCK *pMB = &pMBs[mbnum];
2088                    const VECTOR mv = pMB->mvs[0];
2089    
2090          if (count > *bestcount) {                  if (pMB->mode == MODE_INTRA || pMB->mode == MODE_NOT_CODED)
2091                  *bestcount = count;                          continue;
2092                  *iDirection = dir;  
2093                  GMC->x = x; GMC->y = y;                  if ( ( (ABS(mv.x -   (pMB-1)->mvs[0].x) < deltax) && (ABS(mv.y -   (pMB-1)->mvs[0].y) < deltay) )
2094                    &&   ( (ABS(mv.x -   (pMB+1)->mvs[0].x) < deltax) && (ABS(mv.y -   (pMB+1)->mvs[0].y) < deltay) )
2095                    &&   ( (ABS(mv.x - (pMB-MBw)->mvs[0].x) < deltax) && (ABS(mv.y - (pMB-MBw)->mvs[0].y) < deltay) )
2096                    &&   ( (ABS(mv.x - (pMB+MBw)->mvs[0].x) < deltax) && (ABS(mv.y - (pMB+MBw)->mvs[0].y) < deltay) ) )
2097                            MBmask[mbnum]=1;
2098          }          }
2099    
2100            for (my = 1; my < MBh-1; my++)
2101            for (mx = 1; mx < MBw-1; mx++)
2102            {
2103                    const uint8_t *const pCur = current->image.y + 16*my*pParam->edged_width + 16*mx;
2104    
2105                    const int mbnum = mx + my * MBw;
2106                    if (!MBmask[mbnum])
2107                            continue;
2108    
2109                    if (sad16 ( pCur, pCur+1 , pParam->edged_width, 65536) <= grad )
2110                            MBmask[mbnum] = 0;
2111                    if (sad16 ( pCur, pCur+pParam->edged_width, pParam->edged_width, 65536) <= grad )
2112                            MBmask[mbnum] = 0;
2113    
2114  }  }
2115    
2116            emms();
2117    
2118  static VECTOR          do {            /* until convergence */
 GlobalMotionEst(const MACROBLOCK * const pMBs, const MBParam * const pParam, const uint32_t iFcode)  
 {  
2119    
2120          uint32_t count, bestcount = 0;          a = b = c = n = 0;
2121          int x, y;          DtimesF[0] = DtimesF[1] = DtimesF[2] = DtimesF[3] = 0.;
2122          VECTOR gmc = {0,0};          for (my = 0; my < MBh; my++)
2123          int step, min_x, max_x, min_y, max_y;                  for (mx = 0; mx < MBw; mx++)
2124          uint32_t mx, my;                  {
2125          int iDirection, bDirection;                          const int mbnum = mx + my * MBw;
2126                            const MACROBLOCK *pMB = &pMBs[mbnum];
2127                            const VECTOR mv = pMB->mvs[0];
2128    
2129          min_x = min_y = -32<<iFcode;                          if (!MBmask[mbnum])
2130          max_x = max_y = 32<<iFcode;                                  continue;
2131    
2132  //step1: let's find a rough camera panning                          n++;
2133          for (step = 32; step >= 2; step /= 2) {                          a += 16*mx+8;
2134                  bestcount = 0;                          b += 16*my+8;
2135                  for (y = min_y; y <= max_y; y += step)                          c += (16*mx+8)*(16*mx+8)+(16*my+8)*(16*my+8);
2136                          for (x = min_x ; x <= max_x; x += step) {  
2137                                  count = 0;                          DtimesF[0] += (double)mv.x;
2138                                  //for all macroblocks                          DtimesF[1] += (double)mv.x*(16*mx+8) + (double)mv.y*(16*my+8);
2139                                  for (my = 1; my < pParam->mb_height-1; my++)                          DtimesF[2] += (double)mv.x*(16*my+8) - (double)mv.y*(16*mx+8);
2140                                          for (mx = 1; mx < pParam->mb_width-1; mx++) {                          DtimesF[3] += (double)mv.y;
2141                                                  const MACROBLOCK *pMB = &pMBs[mx + my * pParam->mb_width];                  }
2142                                                  VECTOR mv;  
2143            denom = a*a+b*b-c*n;
2144    
2145    /* Solve the system:     sol = (D'*E*D)^{-1} D'*E*F   */
2146    /* D'*E*F has been calculated in the same loop as matrix */
2147    
2148            sol[0] = -c*DtimesF[0] + a*DtimesF[1] + b*DtimesF[2];
2149            sol[1] =  a*DtimesF[0] - n*DtimesF[1]                + b*DtimesF[3];
2150            sol[2] =  b*DtimesF[0]                - n*DtimesF[2] - a*DtimesF[3];
2151            sol[3] =                 b*DtimesF[1] - a*DtimesF[2] - c*DtimesF[3];
2152    
2153            sol[0] /= denom;
2154            sol[1] /= denom;
2155            sol[2] /= denom;
2156            sol[3] /= denom;
2157    
2158            meanx = meany = 0.;
2159            oldnum = 0;
2160            for (my = 0; my < MBh; my++)
2161                    for (mx = 0; mx < MBw; mx++)
2162                    {
2163                            const int mbnum = mx + my * MBw;
2164                            const MACROBLOCK *pMB = &pMBs[mbnum];
2165                            const VECTOR mv = pMB->mvs[0];
2166    
2167                                                  if (pMB->mode == MODE_INTRA || pMB->mode == MODE_NOT_CODED)                          if (!MBmask[mbnum])
2168                                                          continue;                                                          continue;
2169    
2170                                                  mv = pMB->mvs[0];                          oldnum++;
2171                                                  if ( ABS(mv.x - x) <= step && ABS(mv.y - y) <= step )   /* GMC translation is always halfpel-res */                          meanx += ABS(( sol[0] + (16*mx+8)*sol[1] + (16*my+8)*sol[2] ) - mv.x );
2172                                                          count++;                          meany += ABS(( sol[3] - (16*mx+8)*sol[2] + (16*my+8)*sol[1] ) - mv.y );
                                         }  
                                 if (count >= bestcount) { bestcount = count; gmc.x = x; gmc.y = y; }  
2173                          }                          }
2174                  min_x = gmc.x - step;  
2175                  max_x = gmc.x + step;          if (4*meanx > oldnum)   /* better fit than 0.25 is useless */
2176                  min_y = gmc.y - step;                  meanx /= oldnum;
2177                  max_y = gmc.y + step;          else
2178                    meanx = 0.25;
2179    
2180            if (4*meany > oldnum)
2181                    meany /= oldnum;
2182            else
2183                    meany = 0.25;
2184    
2185    /*      fprintf(stderr,"sol = (%8.5f, %8.5f, %8.5f, %8.5f)\n",sol[0],sol[1],sol[2],sol[3]);
2186            fprintf(stderr,"meanx = %8.5f  meany = %8.5f   %d\n",meanx,meany, oldnum);
2187    */
2188            num = 0;
2189            for (my = 0; my < MBh; my++)
2190                    for (mx = 0; mx < MBw; mx++)
2191                    {
2192                            const int mbnum = mx + my * MBw;
2193                            const MACROBLOCK *pMB = &pMBs[mbnum];
2194                            const VECTOR mv = pMB->mvs[0];
2195    
2196                            if (!MBmask[mbnum])
2197                                    continue;
2198    
2199                            if  ( ( ABS(( sol[0] + (16*mx+8)*sol[1] + (16*my+8)*sol[2] ) - mv.x ) > meanx )
2200                               || ( ABS(( sol[3] - (16*mx+8)*sol[2] + (16*my+8)*sol[1] ) - mv.y ) > meany ) )
2201                                    MBmask[mbnum]=0;
2202                            else
2203                                    num++;
2204          }          }
2205    
2206          if (bestcount < (pParam->mb_height-2)*(pParam->mb_width-2)/10)          } while ( (oldnum != num) && (num>=4) );
                 gmc.x = gmc.y = 0; //no camara pan, no GMC  
2207    
2208  // step2: let's refine camera panning using gradiend-descent approach          if (num < 4)
2209  // TODO: more warping points may be evaluated here (like in interpolate mode search - two vectors in one diamond)          {
2210          bestcount = 0;                  gmc.duv[0].x= gmc.duv[0].y= gmc.duv[1].x= gmc.duv[1].y= gmc.duv[2].x= gmc.duv[2].y=0;
2211          CheckGMC(gmc.x, gmc.y, 255, &iDirection, pMBs, &bestcount, &gmc, pParam);          } else {
         do {  
                 x = gmc.x; y = gmc.y;  
                 bDirection = iDirection; iDirection = 0;  
                 if (bDirection & 1) CheckGMC(x - 1, y, 1+4+8, &iDirection, pMBs, &bestcount, &gmc, pParam);  
                 if (bDirection & 2) CheckGMC(x + 1, y, 2+4+8, &iDirection, pMBs, &bestcount, &gmc, pParam);  
                 if (bDirection & 4) CheckGMC(x, y - 1, 1+2+4, &iDirection, pMBs, &bestcount, &gmc, pParam);  
                 if (bDirection & 8) CheckGMC(x, y + 1, 1+2+8, &iDirection, pMBs, &bestcount, &gmc, pParam);  
2212    
2213          } while (iDirection);                  gmc.duv[0].x=(int)(sol[0]+0.5);
2214                    gmc.duv[0].y=(int)(sol[3]+0.5);
2215    
2216                    gmc.duv[1].x=(int)(sol[1]*pParam->width+0.5);
2217                    gmc.duv[1].y=(int)(-sol[2]*pParam->width+0.5);
2218    
2219          if (pParam->m_quarterpel) {                  gmc.duv[2].x=0;
2220                  gmc.x *= 2;                  gmc.duv[2].y=0;
                 gmc.y *= 2;     /* we store the halfpel value as pseudo-qpel to make comparison easier */  
2221          }          }
2222    //      fprintf(stderr,"wp1 = ( %4d, %4d)  wp2 = ( %4d, %4d) \n", gmc.duv[0].x, gmc.duv[0].y, gmc.duv[1].x, gmc.duv[1].y);
2223    
2224            free(MBmask);
2225    
2226          return gmc;          return gmc;
2227  }  }
Legend:
Removed from v.768  
changed lines
  Added in v.818
No admin address has been configured
 ViewVC Help
Powered by ViewVC 1.0.4