/******************************************************************************
 *
 * This software module was originally developed by
 *
 *   J. Ignacio Ronda (UPM-GTI / ACTS-MoMuSys)
 *
 * and edited by
 *
 *   Fernando Jaureguizar (UPM-GTI / ACTS-MoMuSys)
 *   Robert Danielsen (Telenor / ACTS-MoMuSys)
 *   Angel Pacheco (UPM-GTI / ACTS-MoMuSys)
 *   J. Ignacio Ronda (UPM-GTI / ACTS-MoMuSys)
 *   Ulrike Pestel (TUH / ACTS-MoMuSys)
 *   Minhua Zhou (HHI / ACTS-MoMuSys)
 *   Michael Wollborn (TUH / ACTS-MoMuSys)
 *   Ferran  Marques (UPC / ACTS-MoMuSys)
 *   Cecile Dufour (LEP / ACTS-MoMuSys)
 *   Bob Eifrig (NextLevel Systems)
 *
 * in the course of development of the MPEG-4 Video (ISO/IEC 14496-2) standard.
 * This software module is an implementation of a part of one or more MPEG-4
 * Video (ISO/IEC 14496-2) tools as specified by the MPEG-4 Video (ISO/IEC
 * 14496-2) standard.
 *
 * ISO/IEC gives users of the MPEG-4 Video (ISO/IEC 14496-2) standard free
 * license to this software module or modifications thereof for use in hardware
 * or software products claiming conformance to the MPEG-4 Video (ISO/IEC
 * 14496-2) standard.
 *
 * Those intending to use this software module in hardware or software products
 * are advised that its use may infringe existing patents. The original
 * developer of this software module and his/her company, the subsequent
 * editors and their companies, and ISO/IEC have no liability for use of this
 * software module or modifications thereof in an implementation. Copyright is
 * not released for non MPEG-4 Video (ISO/IEC 14496-2) standard conforming
 * products.
 *
 * ACTS-MoMuSys partners retain full right to use the code for his/her own
 * purpose, assign or donate the code to a third party and to inhibit third
 * parties from using the code for non MPEG-4 Video (ISO/IEC 14496-2) standard
 * conforming products. This copyright notice must be included in all copies or
 * derivative works.
 *
 * Copyright (c) 1997
 *
 *****************************************************************************/

/***********************************************************HeaderBegin*******
 *
 * File:        mot_est_mb.c
 *
 * Authors:     J.Ignacio Ronda / Angel Pacheco, UPM-GTI
 *
 * Created:     02.02.96
 *
 * Description: Functions to perform the motion estimation over one MB
 *
 * Flags:       -D_DEBUG_      : prints out of bound warnings
 *              -D_DEBUG_RANGE_: additional more verbose printing
 *              -D_SAD_        : prints the SAD for the selected MV
 *              -D_SAD_EXHAUS_ : additional more verbose printing (too
 *                               exhaustive with all the SAD's)
 *
 * Modified:
 *      21.04.96 Robert Danielsen: Reformatted. New headers.
 *      02.07.96 Fernando Jaureguizar: Reformatted.
 *      08.07.96 Michael Wollborn: In MB estimation for Intra the sad_min
 *               bug has been changed so that here the SAD is now
 *               computed without prior cancellation.
 *      16.07.96 Fernando Jaureguizar: Fixed two bugs for
 *              estimation flags in line 660 and 755.
 *      23.07.96 Fernando Jaureguizar: Fixed a bug in the MV
 *               range and in the zero vector for HalfPixel. A new
 *               "if" to allow shorter MV when the SAD is the same.
 *               Moved the DEFAULT_8_WIN define to mot_util.h
 *      18.09.96 Fernando Jaureguizar: Now new B_TRANSP and
 *               MB_TRANSP defines are in mot_util.h. Their calls
 *               are modified.
 *      22.10.96 Angel Pacheco: modifications in the range of
 *               the search area due to the new padding scheme (only
 *               the padded MB are valid for the MV search).
 *      28.11.96 Ferran Marques: Integration of changes from (22.10.96)
 *      05.02.97 Angel Pacheco: deletion of the
 *               *pel_pos arrays (VM5.1) having the valid positions
 *               (padded pixels) in which the search for the motion
 *               vectors must be done (as VM4.x said).
 *      14.05.97 Fernando Jaureguizar: new commentaries about new
 *               enlarged f_code range according to VM7.0 Modification
 *               in Obtain_Range() function to cope with enlarged
 *               f_code range.
 *      16.06.97 Angel Pacheco: unified the TRANSPARENT modes.
 *      17.07.97 Ulrike Pestel: ref_sel_code
 *      01.08.97 Fernando Jaureguizar: changed the wrong name
 *               advanced by the correct one enable_8x8_mv.
 *      23.10.97 Cecile Dufour: modifications for static
 *               sprite (low latency)
 *      13.11.97 Angel Pacheco: formating and new header.
 *      09.03.98 Fernando Jaureguizar: New formating. Deleted not used
 *               pmv_x and pmv_y parameters in RangeInSearchArea() and
 *               br_height in MBMotionEstimation()
 *
 ***********************************************************HeaderEnd*********/

/************************    INCLUDE FILES    ********************************/

#include "momusys.h"
#include "mot_util.h"
#include "mot_est.h"

/* ------------------------------------------------------------------------- */

/* Obtaining if two floating point values are equal*/
#define ABSDOUBLE(x)   (((x) > 0.0001) ? (x) : (((x) < -0.0001) ? -(x): 0.0 ))
#define ARE_EQUAL(x,y) ( (ABSDOUBLE((Float)(x)-(y))>0.1)?(0):(1) )

/* auxiliar define for indexing in MBMotionEstimation */
#define INDEX_BIG(x,y) ((x)+(y)*(vop_width+2*edge))
#define INDEX_NOR(x,y) ((x)+(y)*(MB_SIZE))

/* ------------------------------------------------------------------------- */

/***********************************************************CommentBegin******
 *
 * -- RangeInSearchArea -- computes the range of the search area
 *
 * Authors :
 *      UPM-GTI - J. Ignacio Ronda / Angel Pacheco
 *
 * Created :
 *      02.02.96
 *
 * Purpose :
 *      computes the range of the search area for the predicted MV's
 *      INSIDE the overlapped zone between reference and current vops
 *
 * Arguments in :
 *      Int     i,            horizontal MBcoordinate in pixels
 *      Int     j,            vertical MB coordinate in pixels
 *      Int     block,        block position (0 16x16; 1-2-3-4 8x8)
 *      Int     prev_x,       absolute horizontal position of the previous vop
 *      Int     prev_y,       absolute vertical position of the previous vop
 *      Int     vop_width,    horizontal vop dimension
 *      Int     vop_height,   vertical vop dimension
 *      Int     br_x,         absolute horizontal position of the current vop
 *      Int     br_y,         absolute vertical   position of the current vop
 *      Int     edge,         edge arround the reference vop
 *      Int     f_code,       MV search range 1/2 pel: 1=32,2=64,...,7=2048
 *
 * Arguments in/out :
 *
 *
 * Arguments out :
 *      Float   *mv_x_min,    min horizontal range
 *      Float   *mv_x_max,    max horizontal range
 *      Float   *mv_y_min,    min vertical range
 *      Float   *mv_y_max,    max vertical range
 *      Int     *out          the search area does not exist (the reference
 *                            and current BB does not overlap)
 *
 * Return values :
 *      none
 *
 * Side effects :
 *
 *
 * Description :
 *
 *
 * See also :
 *
 *
 * Modified :
 *
 *
 ***********************************************************CommentEnd********/

static Void
RangeInSearchArea(
  Int     i,          /* <-- horizontal MBcoordinate in pixels               */
  Int     j,          /* <-- vertical MB coordinate in pixels                */
  Int     block,      /* <-- block position (0 16x16; 1-2-3-4 8x8)           */
  Int     prev_x,     /* <-- absolute horizontal position of the previous vop*/
  Int     prev_y,     /* <-- absolute vertical position of the previous vop  */
  Int     vop_width,  /* <-- horizontal vop dimension                        */
  Int     vop_height, /* <-- vertical vop dimension                          */
  Int     br_x,       /* <-- absolute horizontal position of the current vop */
  Int     br_y,       /* <-- absolute vertical   position of the current vop */
  Int     edge,       /* <-- edge arround the reference vop                  */
  Int     f_code,     /* <-  MV search range 1/2 pel: 1=32,2=64,...,7=2048 */
  Float   *mv_x_min,  /* <-- min horizontal range                            */
  Float   *mv_x_max,  /* <-- max horizontal range                            */
  Float   *mv_y_min,  /* <-- min vertical range                              */
  Float   *mv_y_max,  /* <-- max vertical range                              */
  Int     *out        /* --> the search area does not exist (the reference   */
                      /*     and current BB does not overlap)                */
  )
{
  Int   dim_curr_x_max,
        dim_curr_y_max,
        dim_curr_x_min,
        dim_curr_y_min;
  Int   dim_prev_x_max,
        dim_prev_y_max,
        dim_prev_x_min,
        dim_prev_y_min;
  Int   mb_b_size,
        block_x,
        block_y;

  *out=0;

  switch (block)
    {
    case 0:                           /* 8x8 or 16x16 block search */
      block_x=0;                      /*****************************/
      block_y=0;                      /** 1 2 ********  0  *********/
      mb_b_size=MB_SIZE;              /** 3 4 ********     *********/
      break;                          /*****************************/
    case 1:
      block_x=0;
      block_y=0;
      mb_b_size=B_SIZE;
      break;
    case 2:
      block_x=B_SIZE;
      block_y=0;
      mb_b_size=B_SIZE;
      break;
    case 3:
      block_x=0;
      block_y=B_SIZE;
      mb_b_size=B_SIZE;
      break;
    case 4:
      block_x=B_SIZE;
      block_y=B_SIZE;
      mb_b_size=B_SIZE;
      break;
    default:
      fprintf(stderr,"RangeInSearchArea: error in block parameter (%d)\n",
        (int)block);
      return;
    }

  /* min x/y */
  dim_curr_x_min=(Int)(br_x+i*MB_SIZE+*mv_x_min+block_x);
  dim_curr_y_min=(Int)(br_y+j*MB_SIZE+*mv_y_min+block_y);
  dim_prev_x_min=prev_x-edge;
  dim_prev_y_min=prev_y-edge;

  /* max x/y */
  /*the MB right-pixels inside */
  dim_curr_x_max=(Int)(br_x+i*MB_SIZE+*mv_x_max+mb_b_size+block_x);
  /*the MB bottom-pixels inside */
  dim_curr_y_max=(Int)(br_y+j*MB_SIZE+*mv_y_max+mb_b_size+block_y);
  dim_prev_x_max=prev_x+vop_width +edge;
  dim_prev_y_max=prev_y+vop_height+edge;

#ifdef _DEBUG_RANGE_
  fprintf(stderr,"BEFORE: x's <%d,%d>  y's <%d,%d> \n",*mv_x_min,*mv_x_max,
          *mv_y_min,*mv_y_max);
#endif

  /* range x/y min */

  if (dim_curr_x_min > dim_prev_x_max)
    {
#ifdef _DEBUG_RANGE_
    fprintf(stderr,"mv_x_min out\n");
#endif
    *out=1;
    }
  else if(dim_curr_x_min < dim_prev_x_min)
    {
#ifdef _DEBUG_RANGE_
    fprintf(stderr,"mv_x_min %d -> ",*mv_x_min);
#endif
    *mv_x_min = *mv_x_min + ( dim_prev_x_min - dim_curr_x_min ) ;
#ifdef _DEBUG_RANGE_
    fprintf(stderr,"%d\n",*mv_x_min);
#endif
    }

  if(!(*out))
    {
    if (dim_curr_y_min > dim_prev_y_max)
      {
#ifdef _DEBUG_RANGE_
      fprintf(stderr,"mv_y_min out\n");
#endif
      *out=1;
      }
    else if(dim_curr_y_min < dim_prev_y_min)
      {
#ifdef _DEBUG_RANGE_
      fprintf(stderr,"mv_y_min %d -> ",*mv_y_min);
#endif
      *mv_y_min = *mv_y_min + ( dim_prev_y_min - dim_curr_y_min ) ;
#ifdef _DEBUG_RANGE_
      fprintf(stderr,"%d\n",*mv_y_min);
#endif
      }
    }

  /* range x/y max */
  if(!(*out))
    {
    if(dim_curr_x_max < dim_prev_x_min)
      {
#ifdef _DEBUG_RANGE_
      fprintf(stderr,"mv_x_max out\n");
#endif
      *out=1;
      }
    if ((!(*out))&&(dim_curr_x_max > dim_prev_x_max))
      {
#ifdef _DEBUG_RANGE_
      fprintf(stderr,"mv_x_max %d -> ",*mv_x_max);
#endif
      *mv_x_max = *mv_x_max - ( dim_curr_x_max - dim_prev_x_max) ;
#ifdef _DEBUG_RANGE_
      fprintf(stderr,"%d\n",*mv_x_max);
#endif
      }
    }

  if(!(*out))
    {
    if(dim_curr_y_max < dim_prev_y_min)
      {
#ifdef _DEBUG_RANGE_
      fprintf(stderr,"mv_x_max out\n");
#endif
      *out=1; /* already set */
      }
    if ((!(*out))&&(dim_curr_y_max > dim_prev_y_max))
      {
#ifdef _DEBUG_RANGE_
      fprintf(stderr,"mv_y_max %d -> ",*mv_y_max);
#endif
      *mv_y_max = *mv_y_max - ( dim_curr_y_max - dim_prev_y_max) ;
#ifdef _DEBUG_RANGE_
      fprintf(stderr,"%d\n",*mv_y_max);
#endif
      }
    }

  if(*mv_x_min>*mv_x_max)
    {
    *out=1;
    }

  if ( (!(*out)) && (*mv_y_min>*mv_y_max))
    {
    *out=1;
    }

#ifdef _DEBUG_
  if (f_code!=0&&(!(*out)) &&
    ((ABS(*mv_x_min)>=(32<<(f_code-1))) || (ABS(*mv_x_max)>=(32<<(f_code-1))) ||
    (ABS(*mv_y_min)>=(32<<(f_code-1))) || (ABS(*mv_y_max)>=(32<<(f_code-1))) ||
    (ABS(*mv_x_max-*mv_x_min)>=(32<<(f_code-1))) ||
    (ABS(*mv_y_max-*mv_y_min)>=(32<<(f_code-1))) ) )
    {
    fprintf(stderr,
    "error in RangeInSearchArea: this message will never be seen, hopefully\n");
    *out=1;
    }
#endif

#ifdef _DEBUG_RANGE_
  if (!(*out))
    fprintf(stderr,"AFTER: x's <%d,%d>  y's <%d,%d> \n",*mv_x_min,*mv_x_max,
            *mv_y_min,*mv_y_max);
  else
    fprintf(stderr,"AFTER out\n");
#endif
  return;
}

/***********************************************************CommentBegin******
 *
 * -- Obtain_Range -- computes the range of the search area
 *
 * Authors :
 *      UPM-GTI - J. Ignacio Ronda / Angel Pacheco
 *
 * Created :
 *      02.02.96
 *
 * Purpose :
 *      computes the range of the search area for the predicted MV's
 *      (it can be outside the overlapped zone between the reference
 *      and current vops)
 *
 * Arguments in :
 *      Int     f_code,      MV search range 1/2 pel: 1=32,2=64,...,7=2048
 *      Int     type,        MBM_INTER16==16x16 search;
 *                           MBM_INTER8==8x8 search
 *      Float   pmv_x,       predicted horizontal motion vector (0.0 in 16x16)
 *      Float   pmv_y,       predicted horizontal motion vector (0.0 in 16x16)
 *
 * Arguments in/out :
 *
 *
 * Arguments out :
 *      Float   *mv_x_min,   min horizontal range
 *      Float   *mv_x_max,   max horizontal range
 *      Float   *mv_y_min,   min vertical range
 *      Float   *mv_y_max    max vertical range
 *
 * Return values :
 *      1 on success, 0 on error
 *
 * Side effects :
 *
 *
 * Description :
 *
 *
 * See also :
 *
 *
 * Modified :
 *      14.05.97 Fernando Jaureguizar: Modification to cope with enlarged
 *               f_code range (VM7.0).
 *      23.10.97 Cecile Dufour: f_code=0 included for static sprite (low
 *               latency)
 *      12.11.97 Bob Eifrig Use explicit search radius not calculated from
 *               f_code.
 *
 ***********************************************************CommentEnd********/

static Int
Obtain_Range(
  Int     f_code,    /* <-- MV search range 1/2 pel: 1=32,2=64,...,7=2048   */
  Int     sr,        /* <-- Serach range (radius)                           */
  Int     type,      /* <-- MBM_INTER16==16x16 search;
                            MBM_INTER8==8x8 search                          */
  Float   pmv_x,     /* <-- predicted horizontal motion vector              */
  Float   pmv_y,     /* <-- predicted horizontal motion vector              */
  Float   *mv_x_min, /* --> min horizontal range                            */
  Float   *mv_x_max, /* --> max horizontal range                            */
  Float   *mv_y_min, /* --> min vertical range                              */
  Float   *mv_y_max  /* --> max vertical range                              */
  )
{
  Int    error;

  Float  aux_x_min, aux_y_min,
         aux_x_max, aux_y_max;

  Float  range;

  error=0;

  if (f_code < 0 || f_code > 8)
    fprintf(stderr,"Obtain_Range: f_code (%d) not supported\n",(int)f_code);
  else
    if (f_code==0) /* for Low Latency Sprite */
    {
    *mv_x_min=0;
    *mv_x_max=0;
    *mv_y_min=0;
    *mv_y_max=0;
    }
  else
    {
    range = sr;

    if ((pmv_x<-range)||(pmv_x>=range))
      {
      fprintf(stderr, "Obtain_Range: Error (f_code %d): pmv_x %f\n",
              (int)f_code,pmv_x);
      error=1;
      }
    *mv_x_min=-range; *mv_x_max= range - 0.5;

    if ((pmv_y<-range)||(pmv_y>=range))
      {
      fprintf(stderr, "Obtain_Range: Error (f_code %d): pmv_y %f\n",
              (int)f_code,pmv_y);
      error=1;
      }
    *mv_y_min=-range; *mv_y_max= range - 0.5;
    }

  if (type==MBM_INTER8)
    {
    aux_x_min=pmv_x - DEFAULT_8_WIN;
    aux_y_min=pmv_y - DEFAULT_8_WIN;
    aux_x_max=pmv_x + DEFAULT_8_WIN;
    aux_y_max=pmv_y + DEFAULT_8_WIN;

    if(*mv_x_min < aux_x_min)
       *mv_x_min = aux_x_min;
    if(*mv_y_min < aux_y_min)
       *mv_y_min = aux_y_min;
    if(*mv_x_max > aux_x_max)
       *mv_x_max = aux_x_max;
    if(*mv_y_max > aux_y_max)
       *mv_y_max = aux_y_max;
    }
  else if (type!=MBM_INTER16)
    {
    fprintf(stderr, "Obtain_Range: erroneous type (%d)\n",(int)type);
    return(0);
    }

  if (error==1)
    return(0);
  else
    return(1);
}

/***********************************************************CommentBegin******
 *
 * -- MBMotionEstimation -- Computes INTEGER PRECISION MV's for a MB
 *
 * Authors :
 *      UPM-GTI - J. Ignacio Ronda / Angel Pacheco
 *
 * Created :
 *      02.02.96
 *
 * Purpose :
 *      Computes INTEGER PRECISION MV's for a MB. Also returns
 *      prediction errors. Requires the whole MVs data images to
 *      obtain prediction for the current MV, which determines search
 *      range
 *
 * Arguments in :
 *      Vop     *curr_vop,     current VOP
 *      SInt    *alpha,        alpha plane of current picture
 *      SInt    *prev,         extended reference picture
 *      SInt    *alpha_sub,    subsampled alpha plane
 *      Int     br_x,          horizontal bounding rectangle coordinate
 *      Int     br_y,          vertical bounding rectangle coordinate
 *      Int     br_width,      bounding rectangule width
 *      Int     i,             horizontal MBcoordinate in pixels
 *      Int     j,             vertical MB coordinate in pixels
 *      Int     prev_x,        absolute horizontal position of the previous
 *                             vop
 *      Int     prev_y,        absolute vertical position of the previous
 *                             vop
 *      Int     vop_width,     horizontal vop dimension
 *      Int     vop_height,    vertical vop dimension
 *      Int     enable_8x8_mv, 8x8 MV (=1) or only 16x16 MV (=0)
 *      Int     edge,          edge
 *      Int     f_code,        MV search range 1/2 pel: 1=32,2=64,...,7=2048
 *      Int     ref_sel_code   for scalability (UPS)
 *
 * Arguments in/out :
 *
 *
 * Arguments out :
 *      Float   *mv16_w,       predicted horizontal 16x16 MV (if approp.)
 *      Float   *mv16_h,       predicted vertical 16x16 MV   (if approp.)
 *      Float   *mv8_w,        predicted horizontal 8x8 MV   (if approp.)
 *      Float   *mv8_h,        predicted vertical 8x8 MV     (if approp.)
 *      Float   *mvF_w,        predicted horizontal field MVs (if approp.)
 *      Float   *mvF_h,        predicted vertical field MVs   (if approp.)
 *      Int     *min_error,    minimum prediction error
 *      SInt    *flags
 *
 * Return values :
 *      none
 *
 * Side effects :
 *
 *
 * Description :
 *
 *
 * See also :
 *
 *
 * Modified :
 *      23.10.96: Angel Pacheco: new input parameter, pel_pos,
 *                to restrict the search area to the reference padded MB's.
 *      05.02.97: Angel Pacheco: deletion of the pel_pos array (VM5.1)
 *      12.11.97  Bob Eifrig: add field motion vector search code.
 *
 ***********************************************************CommentEnd********/

Void
MBMotionEstimation(
  Vop     *curr_vop,     /* <-- Current vop pointer                          */
  SInt    *alpha,        /* <-- alpha plane of current picture               */
  SInt    *prev,         /* <-- extended reference picture                   */
  SInt    *alpha_sub,    /* <-- subsampled alpha plane                       */
  Int     br_x,          /* <-- horizontal bounding rectangle coordinate     */
  Int     br_y,          /* <-- vertical bounding rectangle coordinate       */
  Int     br_width,      /* <-- bounding rectangule width                    */
  Int     i,             /* <-- horizontal MBcoordinate in pixels            */
  Int     j,             /* <-- vertical MB coordinate in pixels             */
  Int     prev_x,        /* <-- absolute horiz. position of previous vop     */
  Int     prev_y,        /* <-- absolute verti. position of previous vop     */
  Int     vop_width,     /* <-- horizontal vop dimension                     */
  Int     vop_height,    /* <-- vertical vop dimension                       */
  Int     enable_8x8_mv, /* <-- 8x8 MV (=1) or only 16x16 MV (=0)            */
  Int     edge,          /* <-- edge                                         */
  Int     f_code,        /* <-- MV search range 1/2 pel: 1=32,2=64,...,7=2048*/
  SInt    *MB_Nb,        /* <-- array with the Nb per MB                     */
  Int     ref_sel_code,  /* <-- ref_sel_code for scalability                 */
  Float   *mv16_w,       /* --> predicted horizontal 16x16 MV(if approp.)    */
  Float   *mv16_h,       /* --> predicted vertical 16x16 MV  (if approp.)    */
  Float   *mv8_w,        /* --> predicted horizontal 8x8 MV  (if approp.)    */
  Float   *mv8_h,        /* --> predicted vertical 8x8 MV    (if approp.)    */
  Float   *mvF_w,        /* --> predicted horizontal field MV(if approp.)    */
  Float   *mvF_h,        /* --> predicted vertical field MV  (if approp.)    */
  Int     *min_error,    /* --> Minimum prediction error                     */
  SInt    *flags
  )
{
  Int     x, y;
  Int     sad, sad_min=MV_MAX_ERROR;
  Int     mv_x;
  Int     mv_y;
  Int     block;
  Float   pmv_x,
          pmv_y;
  SInt    *act_block,
          *alpha_block;

  Float   mv_x_min, mv_x_max,
          mv_y_min, mv_y_max;
  Int     int_mv_x_min=0, int_mv_x_max=0,
          int_mv_y_min=0, int_mv_y_max=0;

  Int     pos16, pos8;
  Int     mvm_width   = br_width/MB_SIZE;

  Int     x_curr      = i*MB_SIZE,
          y_curr      = j*MB_SIZE;
  Int     hb,vb;
  Int     out;

  Int     rel_ori_x;
  Int     rel_ori_y;

#ifdef _SAD_
  Int     sad0;
#endif
  Int     min_error16, min_errorF[4], min_error8 = 0;
  SInt    *curr = (SInt *)GetImageData(GetVopY(curr_vop));
  Int     ewidth = vop_width + 2*edge;

  rel_ori_x=br_x-prev_x;
  rel_ori_y=br_y-prev_y;

  /* Load act_block and alpha_block */

  act_block=  (SInt*)LoadArea(curr, x_curr,y_curr,MB_SIZE,MB_SIZE,br_width);
  alpha_block=(SInt*)LoadArea(alpha,x_curr,y_curr,MB_SIZE,MB_SIZE,br_width);

  /* Compute 16x16 integer MVs */


  /* Obtain range */

  Obtain_Range (f_code, GetVopSearchRangeFor(curr_vop), MBM_INTER16,
                0.0, 0.0, &mv_x_min, &mv_x_max,
                &mv_y_min, &mv_y_max);

  RangeInSearchArea (i,j,0, prev_x, prev_y, vop_width, vop_height,
                     br_x, br_y, edge,f_code, &mv_x_min, &mv_x_max,
                     &mv_y_min, &mv_y_max,&out);

  /* Compute */

  if(!out)
    {
    int_mv_x_min=(int)ceil((double)mv_x_min);
    int_mv_y_min=(int)ceil((double)mv_y_min);
    int_mv_x_max=(int)floor((double)mv_x_max);
    int_mv_y_max=(int)floor((double)mv_y_max);
    flags[0]=ARE_EQUAL(int_mv_x_min,mv_x_min);
    flags[1]=ARE_EQUAL(int_mv_x_max,mv_x_max);
    flags[2]=ARE_EQUAL(int_mv_y_min,mv_y_min);
    flags[3]=ARE_EQUAL(int_mv_y_max,mv_y_max);

    sad_min=MV_MAX_ERROR;
    mv_x = mv_y = 2000; /* A very large MV */

    if(ref_sel_code !=3 )
      {
      for (y=int_mv_y_min; y<=int_mv_y_max; y++)
        for (x=int_mv_x_min; x<=int_mv_x_max; x++)
          {
          if (x==0 && y==0)
#ifdef _SAD_
            sad0 =
#endif
            sad=SAD_Macroblock(prev+INDEX_BIG(x_curr+rel_ori_x,
                               y_curr+rel_ori_y), act_block, alpha_block,
                               (vop_width+2*edge), MV_MAX_ERROR)
                    - (MB_Nb[j*mvm_width +i]/2 + 1);
          else
            sad=SAD_Macroblock(prev+INDEX_BIG(x_curr+x+rel_ori_x,
                               y_curr+y+rel_ori_y), act_block, alpha_block,
                               (vop_width+2*edge), sad_min);
#ifdef _SAD_EXHAUS_
          fprintf(stdout,"+o+ [%2d,%2d] sad16(%3d,%3d)=%4d\n",i,j,x,y,sad);
#endif

          if (sad<sad_min)
              {
              sad_min=sad;
              mv_x=x;
              mv_y=y;
              }
          else if (sad==sad_min)
            if((ABS(x)+ABS(y)) < (ABS(mv_x)+ABS(mv_y)))
              {
              sad_min=sad;
              mv_x=x;
              mv_y=y;
              }
          }/*for*/
      }
    else /* scalability P-VOP ->ref_sel_code==3: MV=0  */
      {
      sad=SAD_Macroblock(prev+INDEX_BIG(x_curr+rel_ori_x,
                         y_curr+rel_ori_y), act_block, alpha_block,
                         (vop_width+2*edge), sad_min);
      sad_min=sad;
      mv_x=0;
      mv_y=0;
      }
    flags[0]=flags[0] && (mv_x==int_mv_x_min);
    flags[1]=flags[1] && (mv_x==int_mv_x_max);
    flags[2]=flags[2] && (mv_y==int_mv_y_min);
    flags[3]=flags[3] && (mv_y==int_mv_y_max);
    }
  else
    {
    mv_x=mv_y=0;
    sad_min=MV_MAX_ERROR;
    }

  /* Store result */

  out=((mv_x==2000)||(mv_y==2000));
  if(out)
    {
    mv_x=mv_y=0;
    sad_min=MV_MAX_ERROR;
    }

  pos16=2*j*2*mvm_width + 2*i;
  mv16_w[pos16]=   mv_x; mv16_h[pos16]=   mv_y;
  mv16_w[pos16+1]= mv_x; mv16_h[pos16+1]= mv_y;
  pos16+=2*mvm_width;
  mv16_w[pos16]=   mv_x; mv16_h[pos16]=   mv_y;
  mv16_w[pos16+1]= mv_x; mv16_h[pos16+1]= mv_y;
  min_error16 = sad_min;

  if (curr_vop->interlaced && !out)
    {
    if ((int_mv_y_min + y_curr) & 1)
      int_mv_y_min++;
    pos16 -= 2*mvm_width;

    /* cur=Top, ref=Top */
    sad_min = MV_MAX_ERROR;
    for (y=int_mv_y_min; y<=int_mv_y_max; y += 2)
      {
      for (x=int_mv_x_min; x<=int_mv_x_max; x++)
        {
        sad = SAD_Field16x8(prev + INDEX_BIG(x_curr + x, y_curr + y),
                            act_block, alpha_block, ewidth, sad_min);

        if ((sad<sad_min) ||
            ((sad==sad_min) && (MVLEN(x,y) < MVLEN(mv_x,mv_y))))
          {
          sad_min=sad;
          mv_x=x;
          mv_y=y;
          }
        }
      }
    min_errorF[0] = sad_min;
    mvF_w[pos16] = (Float)mv_x; mvF_h[pos16] = (Float)mv_y;

    /* cur=Top, ref=Bot */
    sad_min = MV_MAX_ERROR;
    for (y=int_mv_y_min; y<=int_mv_y_max; y += 2)
      {
      for (x=int_mv_x_min; x<=int_mv_x_max; x++)
        {
        sad = SAD_Field16x8(prev + INDEX_BIG(x_curr + x, y_curr + y + 1),
                            act_block, alpha_block, ewidth, sad_min);

        if ((sad<sad_min) ||
            ((sad==sad_min) && (MVLEN(x,y) < MVLEN(mv_x,mv_y))))
          {
          sad_min=sad;
          mv_x=x;
          mv_y=y;
          }
        }
      }
    min_errorF[1] = sad_min;
    mvF_w[pos16+1] = (Float)mv_x; mvF_h[pos16+1] = (Float)mv_y;

    /* cur=Bot, ref=Top */
    sad_min = MV_MAX_ERROR;
    for (y=int_mv_y_min; y<=int_mv_y_max; y += 2)
      {
      for (x=int_mv_x_min; x<=int_mv_x_max; x++)
        {
        sad = SAD_Field16x8(prev + INDEX_BIG(x_curr + x, y_curr + y),
                            act_block + MB_SIZE, alpha_block + MB_SIZE,
                            ewidth, sad_min);


        if ((sad<sad_min) ||
            ((sad==sad_min) && (MVLEN(x,y) < MVLEN(mv_x,mv_y))))
          {
          sad_min=sad;
          mv_x=x;
          mv_y=y;
          }
        }
      }
    min_errorF[2] = sad_min;
    pos16 += 2*mvm_width;
    mvF_w[pos16] = (Float)mv_x; mvF_h[pos16] = (Float)mv_y;

    /* cur=Bot, ref=Bot */
    sad_min = MV_MAX_ERROR;
    for (y=int_mv_y_min; y<=int_mv_y_max; y += 2)
      {
      for (x=int_mv_x_min; x<=int_mv_x_max; x++)
        {
        sad = SAD_Field16x8(prev + INDEX_BIG(x_curr + x, y_curr + y + 1),
                            act_block + MB_SIZE, alpha_block + MB_SIZE,
                            ewidth, sad_min);

        if ((sad<sad_min) ||
            ((sad==sad_min) && (MVLEN(x,y) < MVLEN(mv_x,mv_y))))
          {
          sad_min=sad;
          mv_x=x;
          mv_y=y;
          }
        }
      }
    min_errorF[3] = sad_min;
    mvF_w[pos16+1] = (Float)mv_x; mvF_h[pos16+1] = (Float)mv_y;

    *min_error = MIN(min_errorF[0], min_errorF[1]) + /* Best F-top predictor */
                 MIN(min_errorF[2], min_errorF[3]);  /* Best F-bot predictor */
    if (*min_error > min_error16)
      *min_error = min_error16;
    }
  else
    {
    pos16 -= 2*mvm_width;
    mvF_w[pos16] = mvF_w[pos16+1] = 0.0;
    mvF_h[pos16] = mvF_h[pos16+1] = 0.0;
    pos16 += 2*mvm_width;
    mvF_w[pos16] = mvF_w[pos16+1] = 0.0;
    mvF_h[pos16] = mvF_h[pos16+1] = 0.0;
    *min_error = min_error16;
    }

  /* Compute 8x8 MVs */

  if(enable_8x8_mv==1)
    {
    if(!out)
      {
      for (block=0;block<4;block++)
        {
        /* Obtain range */
        if(block==0)
          {
          hb=vb=0;
          }
        else if (block==1)
          {
          hb=1;vb=0;
          }
        else if (block==2)
          {
          hb=0;vb=1;
          }
        else
          {
          hb=vb=1;
          }
        if(!B_TRANSP(alpha_sub,br_width/B_SIZE,i,j,hb,vb,MBM_TRANSPARENT))
          {
          /* VM 2.*: restrict the search range based on the current 16x16 MV
           *         inside a window around it
           */

          pmv_x=mv16_w[pos16]; pmv_y=mv16_h[pos16];

          Obtain_Range(f_code, GetVopSearchRangeFor(curr_vop), MBM_INTER8,
                       pmv_x, pmv_y, &mv_x_min,
                       &mv_x_max, &mv_y_min, &mv_y_max);

          RangeInSearchArea(i,j, block+1, prev_x, prev_y,
                     vop_width, vop_height, br_x, br_y, edge,f_code,
                     &mv_x_min, &mv_x_max, &mv_y_min,&mv_y_max,&out);

          if(!out)
            {
            int_mv_x_min=(int)ceil((double)mv_x_min);
            int_mv_y_min=(int)ceil((double)mv_y_min);
            int_mv_x_max=(int)floor((double)mv_x_max);
            int_mv_y_max=(int)floor((double)mv_y_max);
            flags[4+block*4]=ARE_EQUAL(int_mv_x_min,mv_x_min);
            flags[4+block*4+1]=ARE_EQUAL(int_mv_x_max,mv_x_max);
            flags[4+block*4+2]=ARE_EQUAL(int_mv_y_min,mv_y_min);
            flags[4+block*4+3]=ARE_EQUAL(int_mv_y_max,mv_y_max);

            sad_min=MV_MAX_ERROR;
            mv_x = mv_y = 2000; /* A very large MV */

            for (y=int_mv_y_min; y<=int_mv_y_max; y++)
              for (x=int_mv_x_min; x<=int_mv_x_max; x++)
                {
                sad=SAD_Block(prev+
                      INDEX_BIG(x_curr + x + 8*(block==1||block==3)+rel_ori_x,
                                y_curr + y + 8*(block==2||block==3)+rel_ori_y),
                       act_block+INDEX_NOR(8*(block==1||block==3),
                                           8*(block==2||block==3)),
                       alpha_block+INDEX_NOR(8*(block==1||block==3),
                                             8*(block==2||block==3)),
                       (vop_width +2*edge), sad_min);
#ifdef _SAD_EXHAUS_
                fprintf(stdout,"-o- [%2d,%2d] sad8_%d(%3d,%3d)=%4d\n",
                        i,j,block,x,y,sad);
#endif

                if (sad<sad_min)
                  {
                  sad_min=sad;
                  mv_x=x;
                  mv_y=y;
                  }
                else if (sad==sad_min)
                  if((ABS(x)+ABS(y)) < (ABS(mv_x)+ABS(mv_y)))
                    {
                    sad_min=sad;
                    mv_x=x;
                    mv_y=y;
                    }
                }/*for*/

            flags[4+block*4]   = flags[4+block*4]   && (mv_x==int_mv_x_min);
            flags[4+block*4+1] = flags[4+block*4+1] && (mv_x==int_mv_x_max);
            flags[4+block*4+2] = flags[4+block*4+2] && (mv_y==int_mv_y_min);
            flags[4+block*4+3] = flags[4+block*4+3] && (mv_y==int_mv_y_max);
            }
          else
            {
            mv_x=mv_y=0;
            sad_min=MV_MAX_ERROR;
            /* punish the OUTER blocks with MV_MAX_ERROR in order to
               be INTRA CODED */
            }
          }
        else /* if B_TRANS */
          {
          mv_x=mv_y=0;
          sad_min=0; /* don't punish the TRANSP blocks with MV_MAX_ERROR; */
          }

        /* Store result */

        if(block==0)
          {
          pos8=2*j*2*mvm_width + 2*i;
          min_error8 += sad_min;
          }
        else if (block==1)
          {
          pos8=2*j*2*mvm_width + 2*i+1;
          min_error8 += sad_min;
          }
        else if (block==2)
          {
          pos8=(2*j+1)*2*mvm_width + 2*i;
          min_error8 += sad_min;
          }
        else
          {
          pos8=(2*j+1)*2*mvm_width + 2*i+1;
          min_error8 += sad_min;
          }

        /* Store result: absolute coordinates (note that in restricted mode
           pmv is (0,0)) */
        mv8_w[pos8]=mv_x;
        mv8_h[pos8]=mv_y;
        }/*for block*/

      if (min_error8 < *min_error)
        *min_error = min_error8;
      }
    else
      { /* all the four blocks are out */
      pos8=2*j*2*mvm_width + 2*i;      mv8_w[pos8]=mv8_h[pos8]=0.0;
      pos8=2*j*2*mvm_width + 2*i+1;    mv8_w[pos8]=mv8_h[pos8]=0.0;
      pos8=(2*j+1)*2*mvm_width + 2*i;  mv8_w[pos8]=mv8_h[pos8]=0.0;
      pos8=(2*j+1)*2*mvm_width + 2*i+1;mv8_w[pos8]=mv8_h[pos8]=0.0;
      min_error8 = MV_MAX_ERROR;
      }
    }/* enable_8x8_mv*/
#ifdef _SAD_
  pos16 = 2*j*2*mvm_width + 2*i;
  fprintf(stdout,"### [%2d,%2d] sad0=%4d         sad16(%5d,%5d)=%4d", i, j,
                  sad0, (int)mv16_w[pos16], (int)mv16_h[pos16], *min_error16);
  if (enable_8x8_mv)
    {
    fprintf(stdout," sad8=%4d\n",
              *min_error8_0+*min_error8_1+*min_error8_2+*min_error8_3);
    pos8=2*j*2*mvm_width + 2*i;
    fprintf(stdout,"### [%2d,%2d]      sad8_0(%5d,%5d)=%4d\n", i, j,
              (int)mv8_w[pos8],(int)mv8_h[pos8],*min_error8_0);
    pos8=2*j*2*mvm_width + 2*i+1;
    fprintf(stdout,"### [%2d,%2d]      sad8_1(%5d,%5d)=%4d\n", i, j,
              (int)mv8_w[pos8],(int)mv8_h[pos8],*min_error8_1);
    pos8=(2*j+1)*2*mvm_width + 2*i;
    fprintf(stdout,"### [%2d,%2d]      sad8_2(%5d,%5d)=%4d\n", i, j,
              (int)mv8_w[pos8],(int)mv8_h[pos8],*min_error8_2);
    pos8=(2*j+1)*2*mvm_width + 2*i+1;
    fprintf(stdout,"### [%2d,%2d]      sad8_3(%5d,%5d)=%4d", i, j,
              (int)mv8_w[pos8],(int)mv8_h[pos8],*min_error8_3);
    }
  fprintf(stdout,"\n");
#endif

  /* Free dynamic memory */
  free((Char*)act_block);
  free((Char*)alpha_block);
}

/***********************************************************CommentBegin******
 *
 * -- FindHalfPel -- Computes MV with half-pixel accurary
 *
 * Authors :
 *      UPM-GTI - J. Ignacio Ronda / Angel Pacheco
 *
 * Created :
 *      02.02.96
 *
 * Purpose :
 *      Computes MV with half-pixel accurary for a 16x16 or 8x8 block
 *
 * Arguments in :
 *      Int     x,             horizontal block coordinate in pixels
 *      Int     y,             vertical blocl coordinate in pixels
 *      SInt    *prev,         interpolated reference vop
 *      Int     *curr,         current block
 *      Int     *curr_alpha,   current block, alpha
 *      Int     bs,            block size (8/16)
 *      Int     comp,          block position in MB (0,1,2,3)
 *      Int     rel_x,         relative horizontal position between
 *                             curr & prev vops
 *      Int     rel_y,         relative vertical position between
 *                             curr & prev vops
 *      Int     pels,          vop width
 *      Int     lines,         vop height
 *      Int     edge,          edge
 *      SInt    flags,
 *
 * Arguments in/out :
 *      Float   *mvx,          horizontal motion vector
 *      Float   *mvy,          vertical motion vector
 *
 * Arguments out :
 *      Int     *min_error     prediction error
 *
 * Return values :
 *      none
 *
 * Side effects :
 *
 *
 * Description :
 *    1) The preference for the 16x16 null motion vector must be applied
 *       again (now, the SAD is recomputed for the recontructed interpolated
 *       reference)
 *
 * See also :
 *
 *
 * Modified :
 *      23.10.96: Angel Pacheco: new input parameter, half_pel_pos,
 *                to restrict the search area to the reference padded MB's.
 *      05.02.97: Angel Pacheco: deleting of the half_pel_pos array
 *                (VM5.1)
 *      23.10.97: Cecile Dufour: added flag_static_sprite for low latency
 *                static sprite
 *
 ***********************************************************CommentEnd********/

Void
FindHalfPel(
  Int    x,           /* <-- horizontal block coordinate in pixels            */
  Int    y,           /* <-- vertical blocl coordinate in pixels              */
  SInt   *prev,       /* <-- interpolated reference vop                       */
  Int    *curr,       /* <-- current block                                    */
  Int    *curr_alpha, /* <-- current block, alpha                             */
  Int    bs,          /* <-- block size (8/16)                                */
  Int    comp,        /* <-- block position in MB (0,1,2,3)                   */
  Int    rel_x,       /* <-- relative horiz. position between curr & prev vops*/
  Int    rel_y,       /* <-- relative verti. position between curr & prev vops*/
  Int    pels,        /* <-- vop width                                        */
  Int    lines,       /* <-- vop height                                       */
  Int    edge,        /* <-- edge                                             */
  SInt   *flags,      /* <-- flags                                            */
  SInt   Nb,          /* <-- the Nb value for (0,0)                           */
  Float  *mvx,        /* <-> horizontal motion vector                         */
  Float  *mvy,        /* <-> vertical motion vector                           */
  Int    *min_error,  /* --> prediction error                                 */
  Int    flag_static_sprite /* for low latency sprite */
  )
{
  Int        i, m, n;
  Int        half_pel=0;
  Int        start_x, start_y,
             stop_x, stop_y,
             new_x, new_y,
             lx;
  Int        min_pos;
  Int        AE, AE_min;
  PixPoint   search[9];
  Int        flag_pos;
  Int        max_i_loop;

#ifdef _SAD_
  Int        sad0;
#endif

  if (flag_static_sprite==1) /* Low latency STATIC_SPRITE */
     max_i_loop=1;
  else
     max_i_loop=9;

  start_x = -1;
  stop_x = 1;
  start_y = -1;
  stop_y = 1;

  new_x = x + (Int)(*mvx) + rel_x; /* referred to prev_ipol image origin */
  new_y = y + (Int)(*mvy) + rel_y; /* referred to prev_ipol image origin */

  new_x += ((comp&1)<<3);
  new_y += ((comp&2)<<2);

  lx = pels + 2*edge;

  /** in 1-pixel search we check if we are inside the range; so don't check
      it again
  **/

  /* avoids trying outside the reference image */
  /* we really move in relative coordinates    */
  /* edge==0 if restricted mode, so it works   */
  /* we also check with flags if we are inside */
  /* search range (1==don't make 1/2 search by */
  /* this side                                 */

  if (bs==8)
    flag_pos=4+comp*4;
  else /* ==16*/
    flag_pos=0;

  if ((new_x <= -edge) || *(flags+flag_pos))
    start_x=-2;
  else if  ((new_x >= (pels-bs+edge)) || *(flags+flag_pos+1))
    stop_x=-2;

  if ((new_y <= -edge) || *(flags+flag_pos+2))
    start_y=-2;
  else if  ((new_y >= (lines-bs+edge)) || *(flags+flag_pos+3))
    stop_y=-2;

  search[0].x = 0;              search[0].y = 0;
  search[1].x = start_x;        search[1].y = start_y; /*   1 2 3   */
  search[2].x = 0;              search[2].y = start_y; /*   4 0 5   */
  search[3].x = stop_x;         search[3].y = start_y; /*   6 7 8   */
  search[4].x = start_x;        search[4].y = 0;
  search[5].x = stop_x;         search[5].y = 0;
  search[6].x = start_x;        search[6].y = stop_y;
  search[7].x = 0;              search[7].y = stop_y;
  search[8].x = stop_x;         search[8].y = stop_y;

  AE_min = MV_MAX_ERROR;
  min_pos = 0;
  for (i = 0; i < max_i_loop; i++)
    {
    if ((search[i].x!=-2)&&(search[i].y!=-2)) /* inside the reference */
      {
      AE = 0;
      for (n = 0; n < bs; n++)
        for (m = 0; m < bs; m++)
          {
          /* Find absolute error */
          half_pel = *(prev + 2*new_x + 2*m + search[i].x +
                      (2*new_y + 2*n + search[i].y)*lx*2);
          AE+=(*(curr_alpha + m + n*16)!=0)*ABS(half_pel - *(curr + m + n*16));
          }
#ifdef _SAD_EXHAUS_
      if (bs == 16)
        fprintf(stdout,"  +o+ half  sad16(%3d,%3d)=%4d\n",search[i].x,
                search[i].y,AE);
      else
        fprintf(stdout,"  -o- half  sad8_%d(%3d,%3d)=%4d\n",comp,search[i].x,
                search[i].y,AE);
#endif

#ifdef _SAD_
      if (i == 0 && *mvx==0 && *mvy==0)
        sad0=AE;
#endif

      if (i==0 && bs==16 && *mvx==0 && *mvy==0)
        AE -= (Nb/2 + 1);

      if (AE < AE_min)
        {
        AE_min = AE;
        min_pos = i;
        }
      }
    /* else don't look outside the reference */
    }
#ifdef _SAD_
    if (bs==16 && *mvx==0 && *mvy==0)
      fprintf(stdout,"  ### half  sad0=%4d(-%3d)",sad0,Nb/2 + 1);
    else
      fprintf(stdout,"  ### half  sad0=    (    )");
#endif

  /* Store optimal values */

  *min_error = AE_min;
  *mvx += ((Float)search[min_pos].x)/2.0;
  *mvy += ((Float)search[min_pos].y)/2.0;

#ifdef _SAD_
  if (bs==16)
    fprintf(stdout,"   sad16(%5.1f,%5.1f)=%4d\n",*mvx,*mvy,*min_error);
  else
    fprintf(stdout,"   sad8_%d(%5.1f,%5.1f)=%4d\n",comp,*mvx,*mvy,*min_error);
#endif

  return;
}

/***********************************************************CommentBegin******
 *
 * -- FindHalfPelField -- Computes field MV with half-pixel accurary
 *
 * Author :
 *      Bob Eifrig
 *
 * Created :
 *      20.01.97
 *
 * Purpose :
 *      Computes field MV with half-pixel accurary for a 16x8 block
 *
 * Arguments in :
 *      Int     x,             horizontal block coordinate in pixels
 *      Int     y,             vertical blocl coordinate in pixels
 *      SInt    *prev,         interpolated reference vop
 *      Int     *curr,         current block
 *      Int     *curr_alpha,   current block, alpha
 *      Int     rel_x,         relative horizontal position between
 *                             curr & prev vops
 *      Int     rel_y,         relative vertical position between
 *                             curr & prev vops
 *      Int     pels,          vop width
 *      Int     lines,         vop height
 *      Int     edge,          edge
 *
 * Arguments in/out :
 *      Float   *mvx,          horizontal motion vector
 *      Float   *mvy,          vertical motion vector
 *
 * Arguments out :
 *      Int     *min_error     prediction error
 *
 * Return values :
 *      none
 *
 * Side effects :
 *
 *
 * See also :
 *
 *
 * Modified :
 *       05.02.97 X. C.: deletion of the *pel_pos arrays (VM5.1)
 *                having the valid positions (padded pixels) in which the
 *                search for the motion vectors must be done (as VM4.x said).
 *
 ***********************************************************CommentEnd********/

Void
FindHalfPelField(
  Int    x,           /* <-- X coord of pel in prev[] of MV=(0,0)             */
  Int    y,           /* <-- Y coord of pel in prev[] of MV=(0,0)             */
  SInt   *prev,       /* <-- ptr to 1st pel in top/bottom fld in reference    */
  Int    *curr,       /* <-- current block                                    */
  Int    *curr_alpha, /* <-- current block, alpha                             */
  Int    pels,        /* <-- vop width                                        */
  Int    lines,       /* <-- vop height                                       */
  Int    edge,        /* <-- edge                                             */
  Float  *mvx,        /* <-> horizontal motion vector                         */
  Float  *mvy,        /* <-> vertical motion vector                           */
  Int    *min_error   /* --> prediction error                                 */
  )
{
  Int        i, j, m, n, v;
  Int        start_x, start_y,
             stop_x, stop_y,
             xx, yy, lx;
  Int        min_pos;
  Int        AE, AE_min;
  PixPoint   search[9];

  start_x = -1;
  stop_x  =  1;
  start_y = -2;
  stop_y  =  2;

  if (y & 1)
    fprintf(stderr,"FindHalfPelField() with odd Y offset between ref & cur\n");
  if ((Int)(*mvy) & 1)
    fprintf(stderr,"FindHalfPelField() with odd Y fullpel MV\n");

  x += (Int)*mvx;   /* Coordinates in prev[] of center of search */
  y += (Int)*mvy;

  lx = pels + 2*edge;

  /** in 1-pixel search we check if we are inside the range; so don't check
      it again
  **/
  if (x <= -edge)
    start_x = -4;
  else if  (x >= (pels - MB_SIZE + edge))
    stop_x = -4;

  if (y <= -edge)
    start_y = -4;
  else if  (y >= (lines - MB_SIZE + edge))
    stop_y = -4;

  search[0].x = 0;              search[0].y = 0;
  search[1].x = start_x;        search[1].y = start_y; /*   1 2 3   */
  search[2].x = 0;              search[2].y = start_y; /*   4 0 5   */
  search[3].x = stop_x;         search[3].y = start_y; /*   6 7 8   */
  search[4].x = start_x;        search[4].y = 0;
  search[5].x = stop_x;         search[5].y = 0;
  search[6].x = start_x;        search[6].y = stop_y;
  search[7].x = 0;              search[7].y = stop_y;
  search[8].x = stop_x;         search[8].y = stop_y;

  AE_min = MV_MAX_ERROR;
  min_pos = 0;

  x <<= 1;          /* Convert to half-pel units */
  y <<= 1;

  for (i = 0; i < 9; i++)
    {
    if ((search[i].x == -4) || (search[i].y == -4))
      continue;
    xx =  (x + search[i].x) >> 1;       /* Truncate to left integral pel */
    yy = ((y + search[i].y) >> 1) & ~1; /* Truncate to upper integral pel
                                           in same field (yy must be even) */
    AE = 0;
    j = xx + yy*lx;
    switch ((0x3B5DC >> (i+i)) & 3)
      {
      case 0:                         /* no interp */
        for (n = 0; n < MB_SIZE; n += 2)
          for (m = 0; m < MB_SIZE; m++)
            {
            if (!curr_alpha[m + n*MB_SIZE])
              continue;
            v = prev[xx + m + lx*(yy + n)] - curr[m + n*MB_SIZE];
            AE += ABS(v);
            }
        break;
      case 1:                         /* Horiz interp */
        for (n = 0; n < MB_SIZE; n += 2)
          for (m = 0; m < MB_SIZE; m++)
            {
            if (!curr_alpha[m + n*MB_SIZE])
              continue;
            j = xx + m + lx*(yy + n);
            v = ((prev[j] + prev[j + 1] + 1) >> 1) - curr[m + n*MB_SIZE];
            AE += ABS(v);
            }
        break;
      case 2:                         /* Vert interp */
        for (n = 0; n < MB_SIZE; n += 2)
          for (m = 0; m < MB_SIZE; m++)
            {
            if (!curr_alpha[m + n*MB_SIZE])
              continue;
            j = xx + m + lx*(yy + n);
            v = ((prev[j] + prev[j + 2*lx] + 1) >> 1) - curr[m + n*MB_SIZE];
            AE += ABS(v);
            }
        break;
      case 3:                         /* Horiz & vert interp */
        for (n = 0; n < MB_SIZE; n += 2)
          for (m = 0; m < MB_SIZE; m++)
            {
            if (!curr_alpha[m + n*MB_SIZE])
              continue;
            j = xx + m + lx*(yy + n);
            v = ((prev[j] + prev[j + 1] + prev[j + 2*lx] +
                  prev[j + 2*lx + 1] + 2) >> 2) - curr[m + n*MB_SIZE];
            AE += ABS(v);
            }
        break;
      }

#ifdef _SAD_EXHAUS_
    fprintf(stdout,"  +o+ half  sad16x8(%3d,%3d)=%4d\n",search[i].x,
            search[i].y,AE);
#endif

    if (AE < AE_min)
      {
      AE_min = AE;
      min_pos = i;
      }
    }

  /* Store optimal values */

  *min_error = AE_min;
  *mvx += (Float)(search[min_pos].x/2.0);
  *mvy += (Float)(search[min_pos].y/2.0); /* Vector is in "frame" coordinates */

#ifdef _SAD_
  fprintf(stdout, "  ### half  sad0=    (    )   sad16x8(%5.1f,%5.1f)=%4d\n",
      *mvx,*mvy,*min_error);
#endif
}
