motion.c
上传用户:ma_junhua
上传日期:2008-04-11
资源大小:2752k
文件大小:47k
开发平台:

C/C++

  1. /* motion.c, motion estimation                                              */
  2. /* Copyright (C) 1996, MPEG Software Simulation Group. All Rights Reserved. */
  3. /*
  4.  * Disclaimer of Warranty
  5.  *
  6.  * These software programs are available to the user without any license fee or
  7.  * royalty on an "as is" basis.  The MPEG Software Simulation Group disclaims
  8.  * any and all warranties, whether express, implied, or statuary, including any
  9.  * implied warranties or merchantability or of fitness for a particular
  10.  * purpose.  In no event shall the copyright-holder be liable for any
  11.  * incidental, punitive, or consequential damages of any kind whatsoever
  12.  * arising from the use of these programs.
  13.  *
  14.  * This disclaimer of warranty extends to the user of these programs and user's
  15.  * customers, employees, agents, transferees, successors, and assigns.
  16.  *
  17.  * The MPEG Software Simulation Group does not represent or warrant that the
  18.  * programs furnished hereunder are free of infringement of any third-party
  19.  * patents.
  20.  *
  21.  * Commercial implementations of MPEG-1 and MPEG-2 video, including shareware,
  22.  * are subject to royalty fees to patent holders.  Many of these patents are
  23.  * general enough such that they are unavoidable regardless of implementation
  24.  * design.
  25.  *
  26.  */
  27. #include <stdio.h>
  28. #include "config.h"
  29. #include "global.h"
  30. /* private prototypes */
  31. static void frame_ME _ANSI_ARGS_((unsigned char *oldorg, unsigned char *neworg,
  32.   unsigned char *oldref, unsigned char *newref, unsigned char *cur,
  33.   int i, int j, int sxf, int syf, int sxb, int syb, struct mbinfo *mbi));
  34. static void field_ME _ANSI_ARGS_((unsigned char *oldorg, unsigned char *neworg,
  35.   unsigned char *oldref, unsigned char *newref, unsigned char *cur,
  36.   unsigned char *curref, int i, int j, int sxf, int syf, int sxb, int syb,
  37.   struct mbinfo *mbi, int secondfield, int ipflag));
  38. static void frame_estimate _ANSI_ARGS_((unsigned char *org,
  39.   unsigned char *ref, unsigned char *mb,
  40.   int i, int j,
  41.   int sx, int sy, int *iminp, int *jminp, int *imintp, int *jmintp,
  42.   int *iminbp, int *jminbp, int *dframep, int *dfieldp,
  43.   int *tselp, int *bselp, int imins[2][2], int jmins[2][2]));
  44. static void field_estimate _ANSI_ARGS_((unsigned char *toporg,
  45.   unsigned char *topref, unsigned char *botorg, unsigned char *botref,
  46.   unsigned char *mb, int i, int j, int sx, int sy, int ipflag,
  47.   int *iminp, int *jminp, int *imin8up, int *jmin8up, int *imin8lp,
  48.   int *jmin8lp, int *dfieldp, int *d8p, int *selp, int *sel8up, int *sel8lp,
  49.   int *iminsp, int *jminsp, int *dsp));
  50. static void dpframe_estimate _ANSI_ARGS_((unsigned char *ref,
  51.   unsigned char *mb, int i, int j, int iminf[2][2], int jminf[2][2],
  52.   int *iminp, int *jminp, int *imindmvp, int *jmindmvp,
  53.   int *dmcp, int *vmcp));
  54. static void dpfield_estimate _ANSI_ARGS_((unsigned char *topref,
  55.   unsigned char *botref, unsigned char *mb,
  56.   int i, int j, int imins, int jmins, int *imindmvp, int *jmindmvp,
  57.   int *dmcp, int *vmcp));
  58. static int fullsearch _ANSI_ARGS_((unsigned char *org, unsigned char *ref,
  59.   unsigned char *blk,
  60.   int lx, int i0, int j0, int sx, int sy, int h, int xmax, int ymax,
  61.   int *iminp, int *jminp));
  62. static int dist1 _ANSI_ARGS_((unsigned char *blk1, unsigned char *blk2,
  63.   int lx, int hx, int hy, int h, int distlim));
  64. static int dist2 _ANSI_ARGS_((unsigned char *blk1, unsigned char *blk2,
  65.   int lx, int hx, int hy, int h));
  66. static int bdist1 _ANSI_ARGS_((unsigned char *pf, unsigned char *pb,
  67.   unsigned char *p2, int lx, int hxf, int hyf, int hxb, int hyb, int h));
  68. static int bdist2 _ANSI_ARGS_((unsigned char *pf, unsigned char *pb,
  69.   unsigned char *p2, int lx, int hxf, int hyf, int hxb, int hyb, int h));
  70. static int variance _ANSI_ARGS_((unsigned char *p, int lx));
  71. /*
  72.  * motion estimation for progressive and interlaced frame pictures
  73.  *
  74.  * oldorg: source frame for forward prediction (used for P and B frames)
  75.  * neworg: source frame for backward prediction (B frames only)
  76.  * oldref: reconstructed frame for forward prediction (P and B frames)
  77.  * newref: reconstructed frame for backward prediction (B frames only)
  78.  * cur:    current frame (the one for which the prediction is formed)
  79.  * sxf,syf: forward search window (frame coordinates)
  80.  * sxb,syb: backward search window (frame coordinates)
  81.  * mbi:    pointer to macroblock info structure
  82.  *
  83.  * results:
  84.  * mbi->
  85.  *  mb_type: 0, MB_INTRA, MB_FORWARD, MB_BACKWARD, MB_FORWARD|MB_BACKWARD
  86.  *  MV[][][]: motion vectors (frame format)
  87.  *  mv_field_sel: top/bottom field (for field prediction)
  88.  *  motion_type: MC_FRAME, MC_FIELD
  89.  *
  90.  * uses global vars: pict_type, frame_pred_dct
  91.  */
  92. void motion_estimation(oldorg,neworg,oldref,newref,cur,curref,
  93.   sxf,syf,sxb,syb,mbi,secondfield,ipflag)
  94. unsigned char *oldorg,*neworg,*oldref,*newref,*cur,*curref;
  95. int sxf,syf,sxb,syb;
  96. struct mbinfo *mbi;
  97. int secondfield,ipflag;
  98. {
  99.   int i, j;
  100.   /* loop through all macroblocks of the picture */
  101.   for (j=0; j<height2; j+=16)
  102.   {
  103.     for (i=0; i<width; i+=16)
  104.     {
  105.       if (pict_struct==FRAME_PICTURE)
  106.         frame_ME(oldorg,neworg,oldref,newref,cur,i,j,sxf,syf,sxb,syb,mbi);
  107.       else
  108.         field_ME(oldorg,neworg,oldref,newref,cur,curref,i,j,sxf,syf,sxb,syb,
  109.           mbi,secondfield,ipflag);
  110.       mbi++;
  111.     }
  112.     if (!quiet)
  113.     {
  114.       putc('.',stderr);
  115.       fflush(stderr);
  116.     }
  117.   }
  118.   if (!quiet)
  119.     putc('n',stderr);
  120. }
  121. static void frame_ME(oldorg,neworg,oldref,newref,cur,i,j,sxf,syf,sxb,syb,mbi)
  122. unsigned char *oldorg,*neworg,*oldref,*newref,*cur;
  123. int i,j,sxf,syf,sxb,syb;
  124. struct mbinfo *mbi;
  125. {
  126.   int imin,jmin,iminf,jminf,iminr,jminr;
  127.   int imint,jmint,iminb,jminb;
  128.   int imintf,jmintf,iminbf,jminbf;
  129.   int imintr,jmintr,iminbr,jminbr;
  130.   int var,v0;
  131.   int dmc,dmcf,dmcr,dmci,vmc,vmcf,vmcr,vmci;
  132.   int dmcfield,dmcfieldf,dmcfieldr,dmcfieldi;
  133.   int tsel,bsel,tself,bself,tselr,bselr;
  134.   unsigned char *mb;
  135.   int imins[2][2],jmins[2][2];
  136.   int imindp,jmindp,imindmv,jmindmv,dmc_dp,vmc_dp;
  137.   mb = cur + i + width*j;
  138.   var = variance(mb,width);
  139.   if (pict_type==I_TYPE)
  140.     mbi->mb_type = MB_INTRA;
  141.   else if (pict_type==P_TYPE)
  142.   {
  143.     if (frame_pred_dct)
  144.     {
  145.       dmc = fullsearch(oldorg,oldref,mb,
  146.                        width,i,j,sxf,syf,16,width,height,&imin,&jmin);
  147.       vmc = dist2(oldref+(imin>>1)+width*(jmin>>1),mb,
  148.                   width,imin&1,jmin&1,16);
  149.       mbi->motion_type = MC_FRAME;
  150.     }
  151.     else
  152.     {
  153.       frame_estimate(oldorg,oldref,mb,i,j,sxf,syf,
  154.         &imin,&jmin,&imint,&jmint,&iminb,&jminb,
  155.         &dmc,&dmcfield,&tsel,&bsel,imins,jmins);
  156.       if (M==1)
  157.         dpframe_estimate(oldref,mb,i,j>>1,imins,jmins,
  158.           &imindp,&jmindp,&imindmv,&jmindmv,&dmc_dp,&vmc_dp);
  159.       /* select between dual prime, frame and field prediction */
  160.       if (M==1 && dmc_dp<dmc && dmc_dp<dmcfield)
  161.       {
  162.         mbi->motion_type = MC_DMV;
  163.         dmc = dmc_dp;
  164.         vmc = vmc_dp;
  165.       }
  166.       else if (dmc<=dmcfield)
  167.       {
  168.         mbi->motion_type = MC_FRAME;
  169.         vmc = dist2(oldref+(imin>>1)+width*(jmin>>1),mb,
  170.                     width,imin&1,jmin&1,16);
  171.       }
  172.       else
  173.       {
  174.         mbi->motion_type = MC_FIELD;
  175.         dmc = dmcfield;
  176.         vmc = dist2(oldref+(tsel?width:0)+(imint>>1)+(width<<1)*(jmint>>1),
  177.                     mb,width<<1,imint&1,jmint&1,8);
  178.         vmc+= dist2(oldref+(bsel?width:0)+(iminb>>1)+(width<<1)*(jminb>>1),
  179.                     mb+width,width<<1,iminb&1,jminb&1,8);
  180.       }
  181.     }
  182.     /* select between intra or non-intra coding:
  183.      *
  184.      * selection is based on intra block variance (var) vs.
  185.      * prediction error variance (vmc)
  186.      *
  187.      * blocks with small prediction error are always coded non-intra
  188.      * even if variance is smaller (is this reasonable?)
  189.      */
  190.     if (vmc>var && vmc>=9*256)
  191.       mbi->mb_type = MB_INTRA;
  192.     else
  193.     {
  194.       /* select between MC / No-MC
  195.        *
  196.        * use No-MC if var(No-MC) <= 1.25*var(MC)
  197.        * (i.e slightly biased towards No-MC)
  198.        *
  199.        * blocks with small prediction error are always coded as No-MC
  200.        * (requires no motion vectors, allows skipping)
  201.        */
  202.       v0 = dist2(oldref+i+width*j,mb,width,0,0,16);
  203.       if (4*v0>5*vmc && v0>=9*256)
  204.       {
  205.         /* use MC */
  206.         var = vmc;
  207.         mbi->mb_type = MB_FORWARD;
  208.         if (mbi->motion_type==MC_FRAME)
  209.         {
  210.           mbi->MV[0][0][0] = imin - (i<<1);
  211.           mbi->MV[0][0][1] = jmin - (j<<1);
  212.         }
  213.         else if (mbi->motion_type==MC_DMV)
  214.         {
  215.           /* these are FRAME vectors */
  216.           /* same parity vector */
  217.           mbi->MV[0][0][0] = imindp - (i<<1);
  218.           mbi->MV[0][0][1] = (jmindp<<1) - (j<<1);
  219.           /* opposite parity vector */
  220.           mbi->dmvector[0] = imindmv;
  221.           mbi->dmvector[1] = jmindmv;
  222.         }
  223.         else
  224.         {
  225.           /* these are FRAME vectors */
  226.           mbi->MV[0][0][0] = imint - (i<<1);
  227.           mbi->MV[0][0][1] = (jmint<<1) - (j<<1);
  228.           mbi->MV[1][0][0] = iminb - (i<<1);
  229.           mbi->MV[1][0][1] = (jminb<<1) - (j<<1);
  230.           mbi->mv_field_sel[0][0] = tsel;
  231.           mbi->mv_field_sel[1][0] = bsel;
  232.         }
  233.       }
  234.       else
  235.       {
  236.         /* No-MC */
  237.         var = v0;
  238.         mbi->mb_type = 0;
  239.         mbi->motion_type = MC_FRAME;
  240.         mbi->MV[0][0][0] = 0;
  241.         mbi->MV[0][0][1] = 0;
  242.       }
  243.     }
  244.   }
  245.   else /* if (pict_type==B_TYPE) */
  246.   {
  247.     if (frame_pred_dct)
  248.     {
  249.       /* forward */
  250.       dmcf = fullsearch(oldorg,oldref,mb,
  251.                         width,i,j,sxf,syf,16,width,height,&iminf,&jminf);
  252.       vmcf = dist2(oldref+(iminf>>1)+width*(jminf>>1),mb,
  253.                    width,iminf&1,jminf&1,16);
  254.       /* backward */
  255.       dmcr = fullsearch(neworg,newref,mb,
  256.                         width,i,j,sxb,syb,16,width,height,&iminr,&jminr);
  257.       vmcr = dist2(newref+(iminr>>1)+width*(jminr>>1),mb,
  258.                    width,iminr&1,jminr&1,16);
  259.       /* interpolated (bidirectional) */
  260.       vmci = bdist2(oldref+(iminf>>1)+width*(jminf>>1),
  261.                     newref+(iminr>>1)+width*(jminr>>1),
  262.                     mb,width,iminf&1,jminf&1,iminr&1,jminr&1,16);
  263.       /* decisions */
  264.       /* select between forward/backward/interpolated prediction:
  265.        * use the one with smallest mean sqaured prediction error
  266.        */
  267.       if (vmcf<=vmcr && vmcf<=vmci)
  268.       {
  269.         vmc = vmcf;
  270.         mbi->mb_type = MB_FORWARD;
  271.       }
  272.       else if (vmcr<=vmci)
  273.       {
  274.         vmc = vmcr;
  275.         mbi->mb_type = MB_BACKWARD;
  276.       }
  277.       else
  278.       {
  279.         vmc = vmci;
  280.         mbi->mb_type = MB_FORWARD|MB_BACKWARD;
  281.       }
  282.       mbi->motion_type = MC_FRAME;
  283.     }
  284.     else
  285.     {
  286.       /* forward prediction */
  287.       frame_estimate(oldorg,oldref,mb,i,j,sxf,syf,
  288.         &iminf,&jminf,&imintf,&jmintf,&iminbf,&jminbf,
  289.         &dmcf,&dmcfieldf,&tself,&bself,imins,jmins);
  290.       /* backward prediction */
  291.       frame_estimate(neworg,newref,mb,i,j,sxb,syb,
  292.         &iminr,&jminr,&imintr,&jmintr,&iminbr,&jminbr,
  293.         &dmcr,&dmcfieldr,&tselr,&bselr,imins,jmins);
  294.       /* calculate interpolated distance */
  295.       /* frame */
  296.       dmci = bdist1(oldref+(iminf>>1)+width*(jminf>>1),
  297.                     newref+(iminr>>1)+width*(jminr>>1),
  298.                     mb,width,iminf&1,jminf&1,iminr&1,jminr&1,16);
  299.       /* top field */
  300.       dmcfieldi = bdist1(
  301.                     oldref+(imintf>>1)+(tself?width:0)+(width<<1)*(jmintf>>1),
  302.                     newref+(imintr>>1)+(tselr?width:0)+(width<<1)*(jmintr>>1),
  303.                     mb,width<<1,imintf&1,jmintf&1,imintr&1,jmintr&1,8);
  304.       /* bottom field */
  305.       dmcfieldi+= bdist1(
  306.                     oldref+(iminbf>>1)+(bself?width:0)+(width<<1)*(jminbf>>1),
  307.                     newref+(iminbr>>1)+(bselr?width:0)+(width<<1)*(jminbr>>1),
  308.                     mb+width,width<<1,iminbf&1,jminbf&1,iminbr&1,jminbr&1,8);
  309.       /* select prediction type of minimum distance from the
  310.        * six candidates (field/frame * forward/backward/interpolated)
  311.        */
  312.       if (dmci<dmcfieldi && dmci<dmcf && dmci<dmcfieldf
  313.           && dmci<dmcr && dmci<dmcfieldr)
  314.       {
  315.         /* frame, interpolated */
  316.         mbi->mb_type = MB_FORWARD|MB_BACKWARD;
  317.         mbi->motion_type = MC_FRAME;
  318.         vmc = bdist2(oldref+(iminf>>1)+width*(jminf>>1),
  319.                      newref+(iminr>>1)+width*(jminr>>1),
  320.                      mb,width,iminf&1,jminf&1,iminr&1,jminr&1,16);
  321.       }
  322.       else if (dmcfieldi<dmcf && dmcfieldi<dmcfieldf
  323.                && dmcfieldi<dmcr && dmcfieldi<dmcfieldr)
  324.       {
  325.         /* field, interpolated */
  326.         mbi->mb_type = MB_FORWARD|MB_BACKWARD;
  327.         mbi->motion_type = MC_FIELD;
  328.         vmc = bdist2(oldref+(imintf>>1)+(tself?width:0)+(width<<1)*(jmintf>>1),
  329.                      newref+(imintr>>1)+(tselr?width:0)+(width<<1)*(jmintr>>1),
  330.                      mb,width<<1,imintf&1,jmintf&1,imintr&1,jmintr&1,8);
  331.         vmc+= bdist2(oldref+(iminbf>>1)+(bself?width:0)+(width<<1)*(jminbf>>1),
  332.                      newref+(iminbr>>1)+(bselr?width:0)+(width<<1)*(jminbr>>1),
  333.                      mb+width,width<<1,iminbf&1,jminbf&1,iminbr&1,jminbr&1,8);
  334.       }
  335.       else if (dmcf<dmcfieldf && dmcf<dmcr && dmcf<dmcfieldr)
  336.       {
  337.         /* frame, forward */
  338.         mbi->mb_type = MB_FORWARD;
  339.         mbi->motion_type = MC_FRAME;
  340.         vmc = dist2(oldref+(iminf>>1)+width*(jminf>>1),mb,
  341.                     width,iminf&1,jminf&1,16);
  342.       }
  343.       else if (dmcfieldf<dmcr && dmcfieldf<dmcfieldr)
  344.       {
  345.         /* field, forward */
  346.         mbi->mb_type = MB_FORWARD;
  347.         mbi->motion_type = MC_FIELD;
  348.         vmc = dist2(oldref+(tself?width:0)+(imintf>>1)+(width<<1)*(jmintf>>1),
  349.                     mb,width<<1,imintf&1,jmintf&1,8);
  350.         vmc+= dist2(oldref+(bself?width:0)+(iminbf>>1)+(width<<1)*(jminbf>>1),
  351.                     mb+width,width<<1,iminbf&1,jminbf&1,8);
  352.       }
  353.       else if (dmcr<dmcfieldr)
  354.       {
  355.         /* frame, backward */
  356.         mbi->mb_type = MB_BACKWARD;
  357.         mbi->motion_type = MC_FRAME;
  358.         vmc = dist2(newref+(iminr>>1)+width*(jminr>>1),mb,
  359.                     width,iminr&1,jminr&1,16);
  360.       }
  361.       else
  362.       {
  363.         /* field, backward */
  364.         mbi->mb_type = MB_BACKWARD;
  365.         mbi->motion_type = MC_FIELD;
  366.         vmc = dist2(newref+(tselr?width:0)+(imintr>>1)+(width<<1)*(jmintr>>1),
  367.                     mb,width<<1,imintr&1,jmintr&1,8);
  368.         vmc+= dist2(newref+(bselr?width:0)+(iminbr>>1)+(width<<1)*(jminbr>>1),
  369.                     mb+width,width<<1,iminbr&1,jminbr&1,8);
  370.       }
  371.     }
  372.     /* select between intra or non-intra coding:
  373.      *
  374.      * selection is based on intra block variance (var) vs.
  375.      * prediction error variance (vmc)
  376.      *
  377.      * blocks with small prediction error are always coded non-intra
  378.      * even if variance is smaller (is this reasonable?)
  379.      */
  380.     if (vmc>var && vmc>=9*256)
  381.       mbi->mb_type = MB_INTRA;
  382.     else
  383.     {
  384.       var = vmc;
  385.       if (mbi->motion_type==MC_FRAME)
  386.       {
  387.         /* forward */
  388.         mbi->MV[0][0][0] = iminf - (i<<1);
  389.         mbi->MV[0][0][1] = jminf - (j<<1);
  390.         /* backward */
  391.         mbi->MV[0][1][0] = iminr - (i<<1);
  392.         mbi->MV[0][1][1] = jminr - (j<<1);
  393.       }
  394.       else
  395.       {
  396.         /* these are FRAME vectors */
  397.         /* forward */
  398.         mbi->MV[0][0][0] = imintf - (i<<1);
  399.         mbi->MV[0][0][1] = (jmintf<<1) - (j<<1);
  400.         mbi->MV[1][0][0] = iminbf - (i<<1);
  401.         mbi->MV[1][0][1] = (jminbf<<1) - (j<<1);
  402.         mbi->mv_field_sel[0][0] = tself;
  403.         mbi->mv_field_sel[1][0] = bself;
  404.         /* backward */
  405.         mbi->MV[0][1][0] = imintr - (i<<1);
  406.         mbi->MV[0][1][1] = (jmintr<<1) - (j<<1);
  407.         mbi->MV[1][1][0] = iminbr - (i<<1);
  408.         mbi->MV[1][1][1] = (jminbr<<1) - (j<<1);
  409.         mbi->mv_field_sel[0][1] = tselr;
  410.         mbi->mv_field_sel[1][1] = bselr;
  411.       }
  412.     }
  413.   }
  414.   mbi->var = var;
  415. }
  416. /*
  417.  * motion estimation for field pictures
  418.  *
  419.  * oldorg: original frame for forward prediction (P and B frames)
  420.  * neworg: original frame for backward prediction (B frames only)
  421.  * oldref: reconstructed frame for forward prediction (P and B frames)
  422.  * newref: reconstructed frame for backward prediction (B frames only)
  423.  * cur:    current original frame (the one for which the prediction is formed)
  424.  * curref: current reconstructed frame (to predict second field from first)
  425.  * sxf,syf: forward search window (frame coordinates)
  426.  * sxb,syb: backward search window (frame coordinates)
  427.  * mbi:    pointer to macroblock info structure
  428.  * secondfield: indicates second field of a frame (in P fields this means
  429.  *              that reference field of opposite parity is in curref instead
  430.  *              of oldref)
  431.  * ipflag: indicates a P type field which is the second field of a frame
  432.  *         in which the first field is I type (this restricts predictions
  433.  *         to be based only on the opposite parity (=I) field)
  434.  *
  435.  * results:
  436.  * mbi->
  437.  *  mb_type: 0, MB_INTRA, MB_FORWARD, MB_BACKWARD, MB_FORWARD|MB_BACKWARD
  438.  *  MV[][][]: motion vectors (field format)
  439.  *  mv_field_sel: top/bottom field
  440.  *  motion_type: MC_FIELD, MC_16X8
  441.  *
  442.  * uses global vars: pict_type, pict_struct
  443.  */
  444. static void field_ME(oldorg,neworg,oldref,newref,cur,curref,i,j,
  445.   sxf,syf,sxb,syb,mbi,secondfield,ipflag)
  446. unsigned char *oldorg,*neworg,*oldref,*newref,*cur,*curref;
  447. int i,j,sxf,syf,sxb,syb;
  448. struct mbinfo *mbi;
  449. int secondfield,ipflag;
  450. {
  451.   int w2;
  452.   unsigned char *mb, *toporg, *topref, *botorg, *botref;
  453.   int var,vmc,v0,dmc,dmcfieldi,dmc8i;
  454.   int imin,jmin,imin8u,jmin8u,imin8l,jmin8l,dmcfield,dmc8,sel,sel8u,sel8l;
  455.   int iminf,jminf,imin8uf,jmin8uf,imin8lf,jmin8lf,dmcfieldf,dmc8f,self,sel8uf,sel8lf;
  456.   int iminr,jminr,imin8ur,jmin8ur,imin8lr,jmin8lr,dmcfieldr,dmc8r,selr,sel8ur,sel8lr;
  457.   int imins,jmins,ds,imindmv,jmindmv,vmc_dp,dmc_dp;
  458.   w2 = width<<1;
  459.   mb = cur + i + w2*j;
  460.   if (pict_struct==BOTTOM_FIELD)
  461.     mb += width;
  462.   var = variance(mb,w2);
  463.   if (pict_type==I_TYPE)
  464.     mbi->mb_type = MB_INTRA;
  465.   else if (pict_type==P_TYPE)
  466.   {
  467.     toporg = oldorg;
  468.     topref = oldref;
  469.     botorg = oldorg + width;
  470.     botref = oldref + width;
  471.     if (secondfield)
  472.     {
  473.       /* opposite parity field is in same frame */
  474.       if (pict_struct==TOP_FIELD)
  475.       {
  476.         /* current is top field */
  477.         botorg = cur + width;
  478.         botref = curref + width;
  479.       }
  480.       else
  481.       {
  482.         /* current is bottom field */
  483.         toporg = cur;
  484.         topref = curref;
  485.       }
  486.     }
  487.     field_estimate(toporg,topref,botorg,botref,mb,i,j,sxf,syf,ipflag,
  488.                    &imin,&jmin,&imin8u,&jmin8u,&imin8l,&jmin8l,
  489.                    &dmcfield,&dmc8,&sel,&sel8u,&sel8l,&imins,&jmins,&ds);
  490.     if (M==1 && !ipflag)  /* generic condition which permits Dual Prime */
  491.       dpfield_estimate(topref,botref,mb,i,j,imins,jmins,&imindmv,&jmindmv,
  492.         &dmc_dp,&vmc_dp);
  493.     /* select between dual prime, field and 16x8 prediction */
  494.     if (M==1 && !ipflag && dmc_dp<dmc8 && dmc_dp<dmcfield)
  495.     {
  496.       /* Dual Prime prediction */
  497.       mbi->motion_type = MC_DMV;
  498.       dmc = dmc_dp;     /* L1 metric */
  499.       vmc = vmc_dp;     /* we already calculated L2 error for Dual */
  500.     }
  501.     else if (dmc8<dmcfield)
  502.     {
  503.       /* 16x8 prediction */
  504.       mbi->motion_type = MC_16X8;
  505.       /* upper half block */
  506.       vmc = dist2((sel8u?botref:topref) + (imin8u>>1) + w2*(jmin8u>>1),
  507.                   mb,w2,imin8u&1,jmin8u&1,8);
  508.       /* lower half block */
  509.       vmc+= dist2((sel8l?botref:topref) + (imin8l>>1) + w2*(jmin8l>>1),
  510.                   mb+8*w2,w2,imin8l&1,jmin8l&1,8);
  511.     }
  512.     else
  513.     {
  514.       /* field prediction */
  515.       mbi->motion_type = MC_FIELD;
  516.       vmc = dist2((sel?botref:topref) + (imin>>1) + w2*(jmin>>1),
  517.                   mb,w2,imin&1,jmin&1,16);
  518.     }
  519.     /* select between intra and non-intra coding */
  520.     if (vmc>var && vmc>=9*256)
  521.       mbi->mb_type = MB_INTRA;
  522.     else
  523.     {
  524.       /* zero MV field prediction from same parity ref. field
  525.        * (not allowed if ipflag is set)
  526.        */
  527.       if (!ipflag)
  528.         v0 = dist2(((pict_struct==BOTTOM_FIELD)?botref:topref) + i + w2*j,
  529.                    mb,w2,0,0,16);
  530.       if (ipflag || (4*v0>5*vmc && v0>=9*256))
  531.       {
  532.         var = vmc;
  533.         mbi->mb_type = MB_FORWARD;
  534.         if (mbi->motion_type==MC_FIELD)
  535.         {
  536.           mbi->MV[0][0][0] = imin - (i<<1);
  537.           mbi->MV[0][0][1] = jmin - (j<<1);
  538.           mbi->mv_field_sel[0][0] = sel;
  539.         }
  540.         else if (mbi->motion_type==MC_DMV)
  541.         {
  542.           /* same parity vector */
  543.           mbi->MV[0][0][0] = imins - (i<<1);
  544.           mbi->MV[0][0][1] = jmins - (j<<1);
  545.           /* opposite parity vector */
  546.           mbi->dmvector[0] = imindmv;
  547.           mbi->dmvector[1] = jmindmv;
  548.         }
  549.         else
  550.         {
  551.           mbi->MV[0][0][0] = imin8u - (i<<1);
  552.           mbi->MV[0][0][1] = jmin8u - (j<<1);
  553.           mbi->MV[1][0][0] = imin8l - (i<<1);
  554.           mbi->MV[1][0][1] = jmin8l - ((j+8)<<1);
  555.           mbi->mv_field_sel[0][0] = sel8u;
  556.           mbi->mv_field_sel[1][0] = sel8l;
  557.         }
  558.       }
  559.       else
  560.       {
  561.         /* No MC */
  562.         var = v0;
  563.         mbi->mb_type = 0;
  564.         mbi->motion_type = MC_FIELD;
  565.         mbi->MV[0][0][0] = 0;
  566.         mbi->MV[0][0][1] = 0;
  567.         mbi->mv_field_sel[0][0] = (pict_struct==BOTTOM_FIELD);
  568.       }
  569.     }
  570.   }
  571.   else /* if (pict_type==B_TYPE) */
  572.   {
  573.     /* forward prediction */
  574.     field_estimate(oldorg,oldref,oldorg+width,oldref+width,mb,
  575.                    i,j,sxf,syf,0,
  576.                    &iminf,&jminf,&imin8uf,&jmin8uf,&imin8lf,&jmin8lf,
  577.                    &dmcfieldf,&dmc8f,&self,&sel8uf,&sel8lf,&imins,&jmins,&ds);
  578.     /* backward prediction */
  579.     field_estimate(neworg,newref,neworg+width,newref+width,mb,
  580.                    i,j,sxb,syb,0,
  581.                    &iminr,&jminr,&imin8ur,&jmin8ur,&imin8lr,&jmin8lr,
  582.                    &dmcfieldr,&dmc8r,&selr,&sel8ur,&sel8lr,&imins,&jmins,&ds);
  583.     /* calculate distances for bidirectional prediction */
  584.     /* field */
  585.     dmcfieldi = bdist1(oldref + (self?width:0) + (iminf>>1) + w2*(jminf>>1),
  586.                        newref + (selr?width:0) + (iminr>>1) + w2*(jminr>>1),
  587.                        mb,w2,iminf&1,jminf&1,iminr&1,jminr&1,16);
  588.     /* 16x8 upper half block */
  589.     dmc8i = bdist1(oldref + (sel8uf?width:0) + (imin8uf>>1) + w2*(jmin8uf>>1),
  590.                    newref + (sel8ur?width:0) + (imin8ur>>1) + w2*(jmin8ur>>1),
  591.                    mb,w2,imin8uf&1,jmin8uf&1,imin8ur&1,jmin8ur&1,8);
  592.     /* 16x8 lower half block */
  593.     dmc8i+= bdist1(oldref + (sel8lf?width:0) + (imin8lf>>1) + w2*(jmin8lf>>1),
  594.                    newref + (sel8lr?width:0) + (imin8lr>>1) + w2*(jmin8lr>>1),
  595.                    mb+8*w2,w2,imin8lf&1,jmin8lf&1,imin8lr&1,jmin8lr&1,8);
  596.     /* select prediction type of minimum distance */
  597.     if (dmcfieldi<dmc8i && dmcfieldi<dmcfieldf && dmcfieldi<dmc8f
  598.         && dmcfieldi<dmcfieldr && dmcfieldi<dmc8r)
  599.     {
  600.       /* field, interpolated */
  601.       mbi->mb_type = MB_FORWARD|MB_BACKWARD;
  602.       mbi->motion_type = MC_FIELD;
  603.       vmc = bdist2(oldref + (self?width:0) + (iminf>>1) + w2*(jminf>>1),
  604.                    newref + (selr?width:0) + (iminr>>1) + w2*(jminr>>1),
  605.                    mb,w2,iminf&1,jminf&1,iminr&1,jminr&1,16);
  606.     }
  607.     else if (dmc8i<dmcfieldf && dmc8i<dmc8f
  608.              && dmc8i<dmcfieldr && dmc8i<dmc8r)
  609.     {
  610.       /* 16x8, interpolated */
  611.       mbi->mb_type = MB_FORWARD|MB_BACKWARD;
  612.       mbi->motion_type = MC_16X8;
  613.       /* upper half block */
  614.       vmc = bdist2(oldref + (sel8uf?width:0) + (imin8uf>>1) + w2*(jmin8uf>>1),
  615.                    newref + (sel8ur?width:0) + (imin8ur>>1) + w2*(jmin8ur>>1),
  616.                    mb,w2,imin8uf&1,jmin8uf&1,imin8ur&1,jmin8ur&1,8);
  617.       /* lower half block */
  618.       vmc+= bdist2(oldref + (sel8lf?width:0) + (imin8lf>>1) + w2*(jmin8lf>>1),
  619.                    newref + (sel8lr?width:0) + (imin8lr>>1) + w2*(jmin8lr>>1),
  620.                    mb+8*w2,w2,imin8lf&1,jmin8lf&1,imin8lr&1,jmin8lr&1,8);
  621.     }
  622.     else if (dmcfieldf<dmc8f && dmcfieldf<dmcfieldr && dmcfieldf<dmc8r)
  623.     {
  624.       /* field, forward */
  625.       mbi->mb_type = MB_FORWARD;
  626.       mbi->motion_type = MC_FIELD;
  627.       vmc = dist2(oldref + (self?width:0) + (iminf>>1) + w2*(jminf>>1),
  628.                   mb,w2,iminf&1,jminf&1,16);
  629.     }
  630.     else if (dmc8f<dmcfieldr && dmc8f<dmc8r)
  631.     {
  632.       /* 16x8, forward */
  633.       mbi->mb_type = MB_FORWARD;
  634.       mbi->motion_type = MC_16X8;
  635.       /* upper half block */
  636.       vmc = dist2(oldref + (sel8uf?width:0) + (imin8uf>>1) + w2*(jmin8uf>>1),
  637.                   mb,w2,imin8uf&1,jmin8uf&1,8);
  638.       /* lower half block */
  639.       vmc+= dist2(oldref + (sel8lf?width:0) + (imin8lf>>1) + w2*(jmin8lf>>1),
  640.                   mb+8*w2,w2,imin8lf&1,jmin8lf&1,8);
  641.     }
  642.     else if (dmcfieldr<dmc8r)
  643.     {
  644.       /* field, backward */
  645.       mbi->mb_type = MB_BACKWARD;
  646.       mbi->motion_type = MC_FIELD;
  647.       vmc = dist2(newref + (selr?width:0) + (iminr>>1) + w2*(jminr>>1),
  648.                   mb,w2,iminr&1,jminr&1,16);
  649.     }
  650.     else
  651.     {
  652.       /* 16x8, backward */
  653.       mbi->mb_type = MB_BACKWARD;
  654.       mbi->motion_type = MC_16X8;
  655.       /* upper half block */
  656.       vmc = dist2(newref + (sel8ur?width:0) + (imin8ur>>1) + w2*(jmin8ur>>1),
  657.                   mb,w2,imin8ur&1,jmin8ur&1,8);
  658.       /* lower half block */
  659.       vmc+= dist2(newref + (sel8lr?width:0) + (imin8lr>>1) + w2*(jmin8lr>>1),
  660.                   mb+8*w2,w2,imin8lr&1,jmin8lr&1,8);
  661.     }
  662.     /* select between intra and non-intra coding */
  663.     if (vmc>var && vmc>=9*256)
  664.       mbi->mb_type = MB_INTRA;
  665.     else
  666.     {
  667.       var = vmc;
  668.       if (mbi->motion_type==MC_FIELD)
  669.       {
  670.         /* forward */
  671.         mbi->MV[0][0][0] = iminf - (i<<1);
  672.         mbi->MV[0][0][1] = jminf - (j<<1);
  673.         mbi->mv_field_sel[0][0] = self;
  674.         /* backward */
  675.         mbi->MV[0][1][0] = iminr - (i<<1);
  676.         mbi->MV[0][1][1] = jminr - (j<<1);
  677.         mbi->mv_field_sel[0][1] = selr;
  678.       }
  679.       else /* MC_16X8 */
  680.       {
  681.         /* forward */
  682.         mbi->MV[0][0][0] = imin8uf - (i<<1);
  683.         mbi->MV[0][0][1] = jmin8uf - (j<<1);
  684.         mbi->mv_field_sel[0][0] = sel8uf;
  685.         mbi->MV[1][0][0] = imin8lf - (i<<1);
  686.         mbi->MV[1][0][1] = jmin8lf - ((j+8)<<1);
  687.         mbi->mv_field_sel[1][0] = sel8lf;
  688.         /* backward */
  689.         mbi->MV[0][1][0] = imin8ur - (i<<1);
  690.         mbi->MV[0][1][1] = jmin8ur - (j<<1);
  691.         mbi->mv_field_sel[0][1] = sel8ur;
  692.         mbi->MV[1][1][0] = imin8lr - (i<<1);
  693.         mbi->MV[1][1][1] = jmin8lr - ((j+8)<<1);
  694.         mbi->mv_field_sel[1][1] = sel8lr;
  695.       }
  696.     }
  697.   }
  698.   mbi->var = var;
  699. }
  700. /*
  701.  * frame picture motion estimation
  702.  *
  703.  * org: top left pel of source reference frame
  704.  * ref: top left pel of reconstructed reference frame
  705.  * mb:  macroblock to be matched
  706.  * i,j: location of mb relative to ref (=center of search window)
  707.  * sx,sy: half widths of search window
  708.  * iminp,jminp,dframep: location and value of best frame prediction
  709.  * imintp,jmintp,tselp: location of best field pred. for top field of mb
  710.  * iminbp,jminbp,bselp: location of best field pred. for bottom field of mb
  711.  * dfieldp: value of field prediction
  712.  */
  713. static void frame_estimate(org,ref,mb,i,j,sx,sy,
  714.   iminp,jminp,imintp,jmintp,iminbp,jminbp,dframep,dfieldp,tselp,bselp,
  715.   imins,jmins)
  716. unsigned char *org,*ref,*mb;
  717. int i,j,sx,sy;
  718. int *iminp,*jminp;
  719. int *imintp,*jmintp,*iminbp,*jminbp;
  720. int *dframep,*dfieldp;
  721. int *tselp,*bselp;
  722. int imins[2][2],jmins[2][2];
  723. {
  724.   int dt,db,dmint,dminb;
  725.   int imint,iminb,jmint,jminb;
  726.   /* frame prediction */
  727.   *dframep = fullsearch(org,ref,mb,width,i,j,sx,sy,16,width,height,
  728.                         iminp,jminp);
  729.   /* predict top field from top field */
  730.   dt = fullsearch(org,ref,mb,width<<1,i,j>>1,sx,sy>>1,8,width,height>>1,
  731.                   &imint,&jmint);
  732.   /* predict top field from bottom field */
  733.   db = fullsearch(org+width,ref+width,mb,width<<1,i,j>>1,sx,sy>>1,8,width,height>>1,
  734.                   &iminb,&jminb);
  735.   imins[0][0] = imint;
  736.   jmins[0][0] = jmint;
  737.   imins[1][0] = iminb;
  738.   jmins[1][0] = jminb;
  739.   /* select prediction for top field */
  740.   if (dt<=db)
  741.   {
  742.     dmint=dt; *imintp=imint; *jmintp=jmint; *tselp=0;
  743.   }
  744.   else
  745.   {
  746.     dmint=db; *imintp=iminb; *jmintp=jminb; *tselp=1;
  747.   }
  748.   /* predict bottom field from top field */
  749.   dt = fullsearch(org,ref,mb+width,width<<1,i,j>>1,sx,sy>>1,8,width,height>>1,
  750.                   &imint,&jmint);
  751.   /* predict bottom field from bottom field */
  752.   db = fullsearch(org+width,ref+width,mb+width,width<<1,i,j>>1,sx,sy>>1,8,width,height>>1,
  753.                   &iminb,&jminb);
  754.   imins[0][1] = imint;
  755.   jmins[0][1] = jmint;
  756.   imins[1][1] = iminb;
  757.   jmins[1][1] = jminb;
  758.   /* select prediction for bottom field */
  759.   if (db<=dt)
  760.   {
  761.     dminb=db; *iminbp=iminb; *jminbp=jminb; *bselp=1;
  762.   }
  763.   else
  764.   {
  765.     dminb=dt; *iminbp=imint; *jminbp=jmint; *bselp=0;
  766.   }
  767.   *dfieldp=dmint+dminb;
  768. }
  769. /*
  770.  * field picture motion estimation subroutine
  771.  *
  772.  * toporg: address of original top reference field
  773.  * topref: address of reconstructed top reference field
  774.  * botorg: address of original bottom reference field
  775.  * botref: address of reconstructed bottom reference field
  776.  * mb:  macroblock to be matched
  777.  * i,j: location of mb (=center of search window)
  778.  * sx,sy: half width/height of search window
  779.  *
  780.  * iminp,jminp,selp,dfieldp: location and distance of best field prediction
  781.  * imin8up,jmin8up,sel8up: location of best 16x8 pred. for upper half of mb
  782.  * imin8lp,jmin8lp,sel8lp: location of best 16x8 pred. for lower half of mb
  783.  * d8p: distance of best 16x8 prediction
  784.  * iminsp,jminsp,dsp: location and distance of best same parity field
  785.  *                    prediction (needed for dual prime, only valid if
  786.  *                    ipflag==0)
  787.  */
  788. static void field_estimate(toporg,topref,botorg,botref,mb,i,j,sx,sy,ipflag,
  789.   iminp,jminp,imin8up,jmin8up,imin8lp,jmin8lp,dfieldp,d8p,selp,sel8up,sel8lp,
  790.   iminsp,jminsp,dsp)
  791. unsigned char *toporg, *topref, *botorg, *botref, *mb;
  792. int i,j,sx,sy;
  793. int ipflag;
  794. int *iminp, *jminp;
  795. int *imin8up, *jmin8up, *imin8lp, *jmin8lp;
  796. int *dfieldp,*d8p;
  797. int *selp, *sel8up, *sel8lp;
  798. int *iminsp, *jminsp, *dsp;
  799. {
  800.   int dt, db, imint, jmint, iminb, jminb, notop, nobot;
  801.   /* if ipflag is set, predict from field of opposite parity only */
  802.   notop = ipflag && (pict_struct==TOP_FIELD);
  803.   nobot = ipflag && (pict_struct==BOTTOM_FIELD);
  804.   /* field prediction */
  805.   /* predict current field from top field */
  806.   if (notop)
  807.     dt = 65536; /* infinity */
  808.   else
  809.     dt = fullsearch(toporg,topref,mb,width<<1,
  810.                     i,j,sx,sy>>1,16,width,height>>1,
  811.                     &imint,&jmint);
  812.   /* predict current field from bottom field */
  813.   if (nobot)
  814.     db = 65536; /* infinity */
  815.   else
  816.     db = fullsearch(botorg,botref,mb,width<<1,
  817.                     i,j,sx,sy>>1,16,width,height>>1,
  818.                     &iminb,&jminb);
  819.   /* same parity prediction (only valid if ipflag==0) */
  820.   if (pict_struct==TOP_FIELD)
  821.   {
  822.     *iminsp = imint; *jminsp = jmint; *dsp = dt;
  823.   }
  824.   else
  825.   {
  826.     *iminsp = iminb; *jminsp = jminb; *dsp = db;
  827.   }
  828.   /* select field prediction */
  829.   if (dt<=db)
  830.   {
  831.     *dfieldp = dt; *iminp = imint; *jminp = jmint; *selp = 0;
  832.   }
  833.   else
  834.   {
  835.     *dfieldp = db; *iminp = iminb; *jminp = jminb; *selp = 1;
  836.   }
  837.   /* 16x8 motion compensation */
  838.   /* predict upper half field from top field */
  839.   if (notop)
  840.     dt = 65536;
  841.   else
  842.     dt = fullsearch(toporg,topref,mb,width<<1,
  843.                     i,j,sx,sy>>1,8,width,height>>1,
  844.                     &imint,&jmint);
  845.   /* predict upper half field from bottom field */
  846.   if (nobot)
  847.     db = 65536;
  848.   else
  849.     db = fullsearch(botorg,botref,mb,width<<1,
  850.                     i,j,sx,sy>>1,8,width,height>>1,
  851.                     &iminb,&jminb);
  852.   /* select prediction for upper half field */
  853.   if (dt<=db)
  854.   {
  855.     *d8p = dt; *imin8up = imint; *jmin8up = jmint; *sel8up = 0;
  856.   }
  857.   else
  858.   {
  859.     *d8p = db; *imin8up = iminb; *jmin8up = jminb; *sel8up = 1;
  860.   }
  861.   /* predict lower half field from top field */
  862.   if (notop)
  863.     dt = 65536;
  864.   else
  865.     dt = fullsearch(toporg,topref,mb+(width<<4),width<<1,
  866.                     i,j+8,sx,sy>>1,8,width,height>>1,
  867.                     &imint,&jmint);
  868.   /* predict lower half field from bottom field */
  869.   if (nobot)
  870.     db = 65536;
  871.   else
  872.     db = fullsearch(botorg,botref,mb+(width<<4),width<<1,
  873.                     i,j+8,sx,sy>>1,8,width,height>>1,
  874.                     &iminb,&jminb);
  875.   /* select prediction for lower half field */
  876.   if (dt<=db)
  877.   {
  878.     *d8p += dt; *imin8lp = imint; *jmin8lp = jmint; *sel8lp = 0;
  879.   }
  880.   else
  881.   {
  882.     *d8p += db; *imin8lp = iminb; *jmin8lp = jminb; *sel8lp = 1;
  883.   }
  884. }
  885. static void dpframe_estimate(ref,mb,i,j,iminf,jminf,
  886.   iminp,jminp,imindmvp, jmindmvp, dmcp, vmcp)
  887. unsigned char *ref, *mb;
  888. int i,j;
  889. int iminf[2][2], jminf[2][2];
  890. int *iminp, *jminp;
  891. int *imindmvp, *jmindmvp;
  892. int *dmcp,*vmcp;
  893. {
  894.   int pref,ppred,delta_x,delta_y;
  895.   int is,js,it,jt,ib,jb,it0,jt0,ib0,jb0;
  896.   int imins,jmins,imint,jmint,iminb,jminb,imindmv,jmindmv;
  897.   int vmc,local_dist;
  898.   /* Calculate Dual Prime distortions for 9 delta candidates
  899.    * for each of the four minimum field vectors
  900.    * Note: only for P pictures!
  901.    */
  902.   /* initialize minimum dual prime distortion to large value */
  903.   vmc = 1 << 30;
  904.   for (pref=0; pref<2; pref++)
  905.   {
  906.     for (ppred=0; ppred<2; ppred++)
  907.     {
  908.       /* convert Cartesian absolute to relative motion vector
  909.        * values (wrt current macroblock address (i,j)
  910.        */
  911.       is = iminf[pref][ppred] - (i<<1);
  912.       js = jminf[pref][ppred] - (j<<1);
  913.       if (pref!=ppred)
  914.       {
  915.         /* vertical field shift adjustment */
  916.         if (ppred==0)
  917.           js++;
  918.         else
  919.           js--;
  920.         /* mvxs and mvys scaling*/
  921.         is<<=1;
  922.         js<<=1;
  923.         if (topfirst == ppred)
  924.         {
  925.           /* second field: scale by 1/3 */
  926.           is = (is>=0) ? (is+1)/3 : -((-is+1)/3);
  927.           js = (js>=0) ? (js+1)/3 : -((-js+1)/3);
  928.         }
  929.         else
  930.           continue;
  931.       }
  932.       /* vector for prediction from field of opposite 'parity' */
  933.       if (topfirst)
  934.       {
  935.         /* vector for prediction of top field from bottom field */
  936.         it0 = ((is+(is>0))>>1);
  937.         jt0 = ((js+(js>0))>>1) - 1;
  938.         /* vector for prediction of bottom field from top field */
  939.         ib0 = ((3*is+(is>0))>>1);
  940.         jb0 = ((3*js+(js>0))>>1) + 1;
  941.       }
  942.       else
  943.       {
  944.         /* vector for prediction of top field from bottom field */
  945.         it0 = ((3*is+(is>0))>>1);
  946.         jt0 = ((3*js+(js>0))>>1) - 1;
  947.         /* vector for prediction of bottom field from top field */
  948.         ib0 = ((is+(is>0))>>1);
  949.         jb0 = ((js+(js>0))>>1) + 1;
  950.       }
  951.       /* convert back to absolute half-pel field picture coordinates */
  952.       is += i<<1;
  953.       js += j<<1;
  954.       it0 += i<<1;
  955.       jt0 += j<<1;
  956.       ib0 += i<<1;
  957.       jb0 += j<<1;
  958.       if (is >= 0 && is <= (width-16)<<1 &&
  959.           js >= 0 && js <= (height-16))
  960.       {
  961.         for (delta_y=-1; delta_y<=1; delta_y++)
  962.         {
  963.           for (delta_x=-1; delta_x<=1; delta_x++)
  964.           {
  965.             /* opposite field coordinates */
  966.             it = it0 + delta_x;
  967.             jt = jt0 + delta_y;
  968.             ib = ib0 + delta_x;
  969.             jb = jb0 + delta_y;
  970.             if (it >= 0 && it <= (width-16)<<1 &&
  971.                 jt >= 0 && jt <= (height-16) &&
  972.                 ib >= 0 && ib <= (width-16)<<1 &&
  973.                 jb >= 0 && jb <= (height-16))
  974.             {
  975.               /* compute prediction error */
  976.               local_dist = bdist2(
  977.                 ref + (is>>1) + (width<<1)*(js>>1),
  978.                 ref + width + (it>>1) + (width<<1)*(jt>>1),
  979.                 mb,             /* current mb location */
  980.                 width<<1,       /* adjacent line distance */
  981.                 is&1, js&1, it&1, jt&1, /* half-pel flags */
  982.                 8);             /* block height */
  983.               local_dist += bdist2(
  984.                 ref + width + (is>>1) + (width<<1)*(js>>1),
  985.                 ref + (ib>>1) + (width<<1)*(jb>>1),
  986.                 mb + width,     /* current mb location */
  987.                 width<<1,       /* adjacent line distance */
  988.                 is&1, js&1, ib&1, jb&1, /* half-pel flags */
  989.                 8);             /* block height */
  990.               /* update delta with least distortion vector */
  991.               if (local_dist < vmc)
  992.               {
  993.                 imins = is;
  994.                 jmins = js;
  995.                 imint = it;
  996.                 jmint = jt;
  997.                 iminb = ib;
  998.                 jminb = jb;
  999.                 imindmv = delta_x;
  1000.                 jmindmv = delta_y;
  1001.                 vmc = local_dist;
  1002.               }
  1003.             }
  1004.           }  /* end delta x loop */
  1005.         } /* end delta y loop */
  1006.       }
  1007.     }
  1008.   }
  1009.   /* Compute L1 error for decision purposes */
  1010.   local_dist = bdist1(
  1011.     ref + (imins>>1) + (width<<1)*(jmins>>1),
  1012.     ref + width + (imint>>1) + (width<<1)*(jmint>>1),
  1013.     mb,
  1014.     width<<1,
  1015.     imins&1, jmins&1, imint&1, jmint&1,
  1016.     8);
  1017.   local_dist += bdist1(
  1018.     ref + width + (imins>>1) + (width<<1)*(jmins>>1),
  1019.     ref + (iminb>>1) + (width<<1)*(jminb>>1),
  1020.     mb + width,
  1021.     width<<1,
  1022.     imins&1, jmins&1, iminb&1, jminb&1,
  1023.     8);
  1024.   *dmcp = local_dist;
  1025.   *iminp = imins;
  1026.   *jminp = jmins;
  1027.   *imindmvp = imindmv;
  1028.   *jmindmvp = jmindmv;
  1029.   *vmcp = vmc;
  1030. }
  1031. static void dpfield_estimate(topref,botref,mb,i,j,imins,jmins,
  1032.   imindmvp, jmindmvp, dmcp, vmcp)
  1033. unsigned char *topref, *botref, *mb;
  1034. int i,j;
  1035. int imins, jmins;
  1036. int *imindmvp, *jmindmvp;
  1037. int *dmcp,*vmcp;
  1038. {
  1039.   unsigned char *sameref, *oppref;
  1040.   int io0,jo0,io,jo,delta_x,delta_y,mvxs,mvys,mvxo0,mvyo0;
  1041.   int imino,jmino,imindmv,jmindmv,vmc_dp,local_dist;
  1042.   /* Calculate Dual Prime distortions for 9 delta candidates */
  1043.   /* Note: only for P pictures! */
  1044.   /* Assign opposite and same reference pointer */
  1045.   if (pict_struct==TOP_FIELD)
  1046.   {
  1047.     sameref = topref;    
  1048.     oppref = botref;
  1049.   }
  1050.   else 
  1051.   {
  1052.     sameref = botref;
  1053.     oppref = topref;
  1054.   }
  1055.   /* convert Cartesian absolute to relative motion vector
  1056.    * values (wrt current macroblock address (i,j)
  1057.    */
  1058.   mvxs = imins - (i<<1);
  1059.   mvys = jmins - (j<<1);
  1060.   /* vector for prediction from field of opposite 'parity' */
  1061.   mvxo0 = (mvxs+(mvxs>0)) >> 1;  /* mvxs // 2 */
  1062.   mvyo0 = (mvys+(mvys>0)) >> 1;  /* mvys // 2 */
  1063.   /* vertical field shift correction */
  1064.   if (pict_struct==TOP_FIELD)
  1065.     mvyo0--;
  1066.   else
  1067.     mvyo0++;
  1068.   /* convert back to absolute coordinates */
  1069.   io0 = mvxo0 + (i<<1);
  1070.   jo0 = mvyo0 + (j<<1);
  1071.   /* initialize minimum dual prime distortion to large value */
  1072.   vmc_dp = 1 << 30;
  1073.   for (delta_y = -1; delta_y <= 1; delta_y++)
  1074.   {
  1075.     for (delta_x = -1; delta_x <=1; delta_x++)
  1076.     {
  1077.       /* opposite field coordinates */
  1078.       io = io0 + delta_x;
  1079.       jo = jo0 + delta_y;
  1080.       if (io >= 0 && io <= (width-16)<<1 &&
  1081.           jo >= 0 && jo <= (height2-16)<<1)
  1082.       {
  1083.         /* compute prediction error */
  1084.         local_dist = bdist2(
  1085.           sameref + (imins>>1) + width2*(jmins>>1),
  1086.           oppref  + (io>>1)    + width2*(jo>>1),
  1087.           mb,             /* current mb location */
  1088.           width2,         /* adjacent line distance */
  1089.           imins&1, jmins&1, io&1, jo&1, /* half-pel flags */
  1090.           16);            /* block height */
  1091.         /* update delta with least distortion vector */
  1092.         if (local_dist < vmc_dp)
  1093.         {
  1094.           imino = io;
  1095.           jmino = jo;
  1096.           imindmv = delta_x;
  1097.           jmindmv = delta_y;
  1098.           vmc_dp = local_dist;
  1099.         }
  1100.       }
  1101.     }  /* end delta x loop */
  1102.   } /* end delta y loop */
  1103.   /* Compute L1 error for decision purposes */
  1104.   *dmcp = bdist1(
  1105.     sameref + (imins>>1) + width2*(jmins>>1),
  1106.     oppref  + (imino>>1) + width2*(jmino>>1),
  1107.     mb,             /* current mb location */
  1108.     width2,         /* adjacent line distance */
  1109.     imins&1, jmins&1, imino&1, jmino&1, /* half-pel flags */
  1110.     16);            /* block height */
  1111.   *imindmvp = imindmv;
  1112.   *jmindmvp = jmindmv;
  1113.   *vmcp = vmc_dp;
  1114. }
  1115. /*
  1116.  * full search block matching
  1117.  *
  1118.  * blk: top left pel of (16*h) block
  1119.  * h: height of block
  1120.  * lx: distance (in bytes) of vertically adjacent pels in ref,blk
  1121.  * org: top left pel of source reference picture
  1122.  * ref: top left pel of reconstructed reference picture
  1123.  * i0,j0: center of search window
  1124.  * sx,sy: half widths of search window
  1125.  * xmax,ymax: right/bottom limits of search area
  1126.  * iminp,jminp: pointers to where the result is stored
  1127.  *              result is given as half pel offset from ref(0,0)
  1128.  *              i.e. NOT relative to (i0,j0)
  1129.  */
  1130. static int fullsearch(org,ref,blk,lx,i0,j0,sx,sy,h,xmax,ymax,iminp,jminp)
  1131. unsigned char *org,*ref,*blk;
  1132. int lx,i0,j0,sx,sy,h,xmax,ymax;
  1133. int *iminp,*jminp;
  1134. {
  1135.   int i,j,imin,jmin,ilow,ihigh,jlow,jhigh;
  1136.   int d,dmin;
  1137.   int k,l,sxy;
  1138.   ilow = i0 - sx;
  1139.   ihigh = i0 + sx;
  1140.   if (ilow<0)
  1141.     ilow = 0;
  1142.   if (ihigh>xmax-16)
  1143.     ihigh = xmax-16;
  1144.   jlow = j0 - sy;
  1145.   jhigh = j0 + sy;
  1146.   if (jlow<0)
  1147.     jlow = 0;
  1148.   if (jhigh>ymax-h)
  1149.     jhigh = ymax-h;
  1150.   /* full pel search, spiraling outwards */
  1151.   imin = i0;
  1152.   jmin = j0;
  1153.   dmin = dist1(org+imin+lx*jmin,blk,lx,0,0,h,65536);
  1154.   sxy = (sx>sy) ? sx : sy;
  1155.   for (l=1; l<=sxy; l++)
  1156.   {
  1157.     i = i0 - l;
  1158.     j = j0 - l;
  1159.     for (k=0; k<8*l; k++)
  1160.     {
  1161.       if (i>=ilow && i<=ihigh && j>=jlow && j<=jhigh)
  1162.       {
  1163.         d = dist1(org+i+lx*j,blk,lx,0,0,h,dmin);
  1164.         if (d<dmin)
  1165.         {
  1166.           dmin = d;
  1167.           imin = i;
  1168.           jmin = j;
  1169.         }
  1170.       }
  1171.       if      (k<2*l) i++;
  1172.       else if (k<4*l) j++;
  1173.       else if (k<6*l) i--;
  1174.       else            j--;
  1175.     }
  1176.   }
  1177.   /* half pel */
  1178.   dmin = 65536;
  1179.   imin <<= 1;
  1180.   jmin <<= 1;
  1181.   ilow = imin - (imin>0);
  1182.   ihigh = imin + (imin<((xmax-16)<<1));
  1183.   jlow = jmin - (jmin>0);
  1184.   jhigh = jmin + (jmin<((ymax-h)<<1));
  1185.   for (j=jlow; j<=jhigh; j++)
  1186.     for (i=ilow; i<=ihigh; i++)
  1187.     {
  1188.       d = dist1(ref+(i>>1)+lx*(j>>1),blk,lx,i&1,j&1,h,dmin);
  1189.       if (d<dmin)
  1190.       {
  1191.         dmin = d;
  1192.         imin = i;
  1193.         jmin = j;
  1194.       }
  1195.     }
  1196.   *iminp = imin;
  1197.   *jminp = jmin;
  1198.   return dmin;
  1199. }
  1200. /*
  1201.  * total absolute difference between two (16*h) blocks
  1202.  * including optional half pel interpolation of blk1 (hx,hy)
  1203.  * blk1,blk2: addresses of top left pels of both blocks
  1204.  * lx:        distance (in bytes) of vertically adjacent pels
  1205.  * hx,hy:     flags for horizontal and/or vertical interpolation
  1206.  * h:         height of block (usually 8 or 16)
  1207.  * distlim:   bail out if sum exceeds this value
  1208.  */
  1209. static int dist1(blk1,blk2,lx,hx,hy,h,distlim)
  1210. unsigned char *blk1,*blk2;
  1211. int lx,hx,hy,h;
  1212. int distlim;
  1213. {
  1214.   unsigned char *p1,*p1a,*p2;
  1215.   int i,j;
  1216.   int s,v;
  1217.   s = 0;
  1218.   p1 = blk1;
  1219.   p2 = blk2;
  1220.   if (!hx && !hy)
  1221.     for (j=0; j<h; j++)
  1222.     {
  1223.       if ((v = p1[0]  - p2[0])<0)  v = -v; s+= v;
  1224.       if ((v = p1[1]  - p2[1])<0)  v = -v; s+= v;
  1225.       if ((v = p1[2]  - p2[2])<0)  v = -v; s+= v;
  1226.       if ((v = p1[3]  - p2[3])<0)  v = -v; s+= v;
  1227.       if ((v = p1[4]  - p2[4])<0)  v = -v; s+= v;
  1228.       if ((v = p1[5]  - p2[5])<0)  v = -v; s+= v;
  1229.       if ((v = p1[6]  - p2[6])<0)  v = -v; s+= v;
  1230.       if ((v = p1[7]  - p2[7])<0)  v = -v; s+= v;
  1231.       if ((v = p1[8]  - p2[8])<0)  v = -v; s+= v;
  1232.       if ((v = p1[9]  - p2[9])<0)  v = -v; s+= v;
  1233.       if ((v = p1[10] - p2[10])<0) v = -v; s+= v;
  1234.       if ((v = p1[11] - p2[11])<0) v = -v; s+= v;
  1235.       if ((v = p1[12] - p2[12])<0) v = -v; s+= v;
  1236.       if ((v = p1[13] - p2[13])<0) v = -v; s+= v;
  1237.       if ((v = p1[14] - p2[14])<0) v = -v; s+= v;
  1238.       if ((v = p1[15] - p2[15])<0) v = -v; s+= v;
  1239.       if (s >= distlim)
  1240.         break;
  1241.       p1+= lx;
  1242.       p2+= lx;
  1243.     }
  1244.   else if (hx && !hy)
  1245.     for (j=0; j<h; j++)
  1246.     {
  1247.       for (i=0; i<16; i++)
  1248.       {
  1249.         v = ((unsigned int)(p1[i]+p1[i+1]+1)>>1) - p2[i];
  1250.         if (v>=0)
  1251.           s+= v;
  1252.         else
  1253.           s-= v;
  1254.       }
  1255.       p1+= lx;
  1256.       p2+= lx;
  1257.     }
  1258.   else if (!hx && hy)
  1259.   {
  1260.     p1a = p1 + lx;
  1261.     for (j=0; j<h; j++)
  1262.     {
  1263.       for (i=0; i<16; i++)
  1264.       {
  1265.         v = ((unsigned int)(p1[i]+p1a[i]+1)>>1) - p2[i];
  1266.         if (v>=0)
  1267.           s+= v;
  1268.         else
  1269.           s-= v;
  1270.       }
  1271.       p1 = p1a;
  1272.       p1a+= lx;
  1273.       p2+= lx;
  1274.     }
  1275.   }
  1276.   else /* if (hx && hy) */
  1277.   {
  1278.     p1a = p1 + lx;
  1279.     for (j=0; j<h; j++)
  1280.     {
  1281.       for (i=0; i<16; i++)
  1282.       {
  1283.         v = ((unsigned int)(p1[i]+p1[i+1]+p1a[i]+p1a[i+1]+2)>>2) - p2[i];
  1284.         if (v>=0)
  1285.           s+= v;
  1286.         else
  1287.           s-= v;
  1288.       }
  1289.       p1 = p1a;
  1290.       p1a+= lx;
  1291.       p2+= lx;
  1292.     }
  1293.   }
  1294.   return s;
  1295. }
  1296. /*
  1297.  * total squared difference between two (16*h) blocks
  1298.  * including optional half pel interpolation of blk1 (hx,hy)
  1299.  * blk1,blk2: addresses of top left pels of both blocks
  1300.  * lx:        distance (in bytes) of vertically adjacent pels
  1301.  * hx,hy:     flags for horizontal and/or vertical interpolation
  1302.  * h:         height of block (usually 8 or 16)
  1303.  */
  1304. static int dist2(blk1,blk2,lx,hx,hy,h)
  1305. unsigned char *blk1,*blk2;
  1306. int lx,hx,hy,h;
  1307. {
  1308.   unsigned char *p1,*p1a,*p2;
  1309.   int i,j;
  1310.   int s,v;
  1311.   s = 0;
  1312.   p1 = blk1;
  1313.   p2 = blk2;
  1314.   if (!hx && !hy)
  1315.     for (j=0; j<h; j++)
  1316.     {
  1317.       for (i=0; i<16; i++)
  1318.       {
  1319.         v = p1[i] - p2[i];
  1320.         s+= v*v;
  1321.       }
  1322.       p1+= lx;
  1323.       p2+= lx;
  1324.     }
  1325.   else if (hx && !hy)
  1326.     for (j=0; j<h; j++)
  1327.     {
  1328.       for (i=0; i<16; i++)
  1329.       {
  1330.         v = ((unsigned int)(p1[i]+p1[i+1]+1)>>1) - p2[i];
  1331.         s+= v*v;
  1332.       }
  1333.       p1+= lx;
  1334.       p2+= lx;
  1335.     }
  1336.   else if (!hx && hy)
  1337.   {
  1338.     p1a = p1 + lx;
  1339.     for (j=0; j<h; j++)
  1340.     {
  1341.       for (i=0; i<16; i++)
  1342.       {
  1343.         v = ((unsigned int)(p1[i]+p1a[i]+1)>>1) - p2[i];
  1344.         s+= v*v;
  1345.       }
  1346.       p1 = p1a;
  1347.       p1a+= lx;
  1348.       p2+= lx;
  1349.     }
  1350.   }
  1351.   else /* if (hx && hy) */
  1352.   {
  1353.     p1a = p1 + lx;
  1354.     for (j=0; j<h; j++)
  1355.     {
  1356.       for (i=0; i<16; i++)
  1357.       {
  1358.         v = ((unsigned int)(p1[i]+p1[i+1]+p1a[i]+p1a[i+1]+2)>>2) - p2[i];
  1359.         s+= v*v;
  1360.       }
  1361.       p1 = p1a;
  1362.       p1a+= lx;
  1363.       p2+= lx;
  1364.     }
  1365.   }
  1366.   return s;
  1367. }
  1368. /*
  1369.  * absolute difference error between a (16*h) block and a bidirectional
  1370.  * prediction
  1371.  *
  1372.  * p2: address of top left pel of block
  1373.  * pf,hxf,hyf: address and half pel flags of forward ref. block
  1374.  * pb,hxb,hyb: address and half pel flags of backward ref. block
  1375.  * h: height of block
  1376.  * lx: distance (in bytes) of vertically adjacent pels in p2,pf,pb
  1377.  */
  1378. static int bdist1(pf,pb,p2,lx,hxf,hyf,hxb,hyb,h)
  1379. unsigned char *pf,*pb,*p2;
  1380. int lx,hxf,hyf,hxb,hyb,h;
  1381. {
  1382.   unsigned char *pfa,*pfb,*pfc,*pba,*pbb,*pbc;
  1383.   int i,j;
  1384.   int s,v;
  1385.   pfa = pf + hxf;
  1386.   pfb = pf + lx*hyf;
  1387.   pfc = pfb + hxf;
  1388.   pba = pb + hxb;
  1389.   pbb = pb + lx*hyb;
  1390.   pbc = pbb + hxb;
  1391.   s = 0;
  1392.   for (j=0; j<h; j++)
  1393.   {
  1394.     for (i=0; i<16; i++)
  1395.     {
  1396.       v = ((((unsigned int)(*pf++ + *pfa++ + *pfb++ + *pfc++ + 2)>>2) +
  1397.             ((unsigned int)(*pb++ + *pba++ + *pbb++ + *pbc++ + 2)>>2) + 1)>>1)
  1398.            - *p2++;
  1399.       if (v>=0)
  1400.         s+= v;
  1401.       else
  1402.         s-= v;
  1403.     }
  1404.     p2+= lx-16;
  1405.     pf+= lx-16;
  1406.     pfa+= lx-16;
  1407.     pfb+= lx-16;
  1408.     pfc+= lx-16;
  1409.     pb+= lx-16;
  1410.     pba+= lx-16;
  1411.     pbb+= lx-16;
  1412.     pbc+= lx-16;
  1413.   }
  1414.   return s;
  1415. }
  1416. /*
  1417.  * squared error between a (16*h) block and a bidirectional
  1418.  * prediction
  1419.  *
  1420.  * p2: address of top left pel of block
  1421.  * pf,hxf,hyf: address and half pel flags of forward ref. block
  1422.  * pb,hxb,hyb: address and half pel flags of backward ref. block
  1423.  * h: height of block
  1424.  * lx: distance (in bytes) of vertically adjacent pels in p2,pf,pb
  1425.  */
  1426. static int bdist2(pf,pb,p2,lx,hxf,hyf,hxb,hyb,h)
  1427. unsigned char *pf,*pb,*p2;
  1428. int lx,hxf,hyf,hxb,hyb,h;
  1429. {
  1430.   unsigned char *pfa,*pfb,*pfc,*pba,*pbb,*pbc;
  1431.   int i,j;
  1432.   int s,v;
  1433.   pfa = pf + hxf;
  1434.   pfb = pf + lx*hyf;
  1435.   pfc = pfb + hxf;
  1436.   pba = pb + hxb;
  1437.   pbb = pb + lx*hyb;
  1438.   pbc = pbb + hxb;
  1439.   s = 0;
  1440.   for (j=0; j<h; j++)
  1441.   {
  1442.     for (i=0; i<16; i++)
  1443.     {
  1444.       v = ((((unsigned int)(*pf++ + *pfa++ + *pfb++ + *pfc++ + 2)>>2) +
  1445.             ((unsigned int)(*pb++ + *pba++ + *pbb++ + *pbc++ + 2)>>2) + 1)>>1) - *p2++;
  1446.       s+=v*v;
  1447.     }
  1448.     p2+= lx-16;
  1449.     pf+= lx-16;
  1450.     pfa+= lx-16;
  1451.     pfb+= lx-16;
  1452.     pfc+= lx-16;
  1453.     pb+= lx-16;
  1454.     pba+= lx-16;
  1455.     pbb+= lx-16;
  1456.     pbc+= lx-16;
  1457.   }
  1458.   return s;
  1459. }
  1460. /*
  1461.  * variance of a (16*16) block, multiplied by 256
  1462.  * p:  address of top left pel of block
  1463.  * lx: distance (in bytes) of vertically adjacent pels
  1464.  */
  1465. static int variance(p,lx)
  1466. unsigned char *p;
  1467. int lx;
  1468. {
  1469.   int i,j;
  1470.   unsigned int v,s,s2;
  1471.   s = s2 = 0;
  1472.   for (j=0; j<16; j++)
  1473.   {
  1474.     for (i=0; i<16; i++)
  1475.     {
  1476.       v = *p++;
  1477.       s+= v;
  1478.       s2+= v*v;
  1479.     }
  1480.     p+= lx-16;
  1481.   }
  1482.   return s2 - (s*s)/256;
  1483. }