#include "VsipTypes.h"

/****************************************************************************
 *
 *  Function    : VsipExpand
 *
 *  Description : Expands the image by using bilinear interpolation
 *                in bilinear interpolation the pixel receives an value
 *                that is the compbination of its four neighbours values
 *                weighted by the corresponding area of influence of each
 *                neighbour.
 *
 *  Parameters  : src - input image object
 *                rsz - row size of the output image
 *                csz - column size of the output image
 *
 *  Returns     : Expanded image
 *
 ****************************************************************************/
$generic

   $VsipExpand = 'VsipExpandBit',    $IOType = 'uint1', $CompType = 'int16';
   $VsipExpand = 'VsipExpandByte',   $IOType = 'uint8', $CompType = 'int32';
   $VsipExpand = 'VsipExpandInt',    $IOType = 'int32', $CompType = 'int32';
   $VsipExpand = 'VsipExpandFloat',  $IOType = 'float', $CompType = 'float';
 

$in

$IOType[:,:] $VsipExpand($IOType src[:,:], VsipIndexesType rsz, VsipIndexesType csz)
{

  // ***** Recover the src size *****
  VsipIndexesType rs, VsipIndexesType cs = extents(src);

  // ***** Test the image and dest(expanded)  size *****
  assert(((rs > 1) && (cs > 1)),
         "ERROR: Image can not have zero dimention. (",rs,"x",cs,")\n");
  assert((rsz > rs) && (csz > cs),
         "ERROR: Destination must be bigger than source (",rsz,",",csz,").\n");
 

  VsipIndexesType rsrc = rs - 1;
  VsipIndexesType csrc = cs - 1;
  VsipIndexesType rres = rsz - 1;
  VsipIndexesType cres = csz - 1;
 

  // ***** Computes the expansion *****
  $IOType result[:,:] =
       // ***** for all pixels in the expanded image *****
       forVsipIndexesType i in [0~rres] cross VsipIndexesType j in [0~cres]
            { // ***** Computes the floor and ceiling of row and column *****
             // ***** in the original image                            *****
             VsipIndexesType floorr = i * rsrc / rres;
             VsipIndexesType floorc = j * csrc / cres;
             VsipIndexesType remainderr = i * rsrc % rres;
             VsipIndexesType remainderc = j * csrc % cres;
             VsipIndexesType ceilr =
                  if (remainderr == 0)
                       return(floorr)
                  else
                       return(floorr + 1);
             VsipIndexesType ceilc =
                  if (remainderc == 0)
                       return(floorc)
                  else
                       return(floorc + 1);

              // ***** Computes the new value. There are 4 different   *****
              // ***** possibilities. 1- The new pixel matches an old  *****
              // ***** pixel, 2- The new pixel matches the row of an   *****
              // ***** old pixel (linear interpolation), 3- The new    *****
              // ***** pixel matchs the column of an old pixel (linear *****
              // ***** interpolation), and 4- The new pixel does not   *****
              // ***** match row or column of an old pixel (bilinear   *****
              // ***** interpolation)                                  *****
              $IOType newval =
                  if ((remainderr == 0) && (remainderc == 0))
                       return(src[floorr,floorc])
                  elif (remainderr == 0)
                       return ((($CompType)src[floorr,floorc]*(ceilc*cres - j*csrc) +
                                ($CompType)src[floorr,ceilc]*(j*csrc - floorc*cres)) /
                                 cres)
                  elif (remainderc == 0)
                       return ((($CompType)src[floorr,floorc]*(ceilr*rres - i*rsrc) +
                                ($CompType)src[ceilr,floorc]*(i*rsrc - floorr*rres)) /
                                 rres)
                  else
                       return ((($CompType)src[floorr,floorc] * (ceilr*rres-i*rsrc) *
                                                                 (ceilc*cres-j*csrc) +
                                ($CompType)src[ceilr,floorc] * (i*rsrc-floorr*rres) *
                                                               (ceilc*cres-j*csrc) +
                                ($CompType)src[ceilr,ceilc] * (i*rsrc-floorr*rres) *
                                                              (j*csrc-floorc*cres)+
                                ($CompType)src[floorr,ceilc] * (ceilr*rres-i*rsrc) *
                                                               (j*csrc-floorc*cres)) /
                                (rres*cres));
             } return(array(newval));
 

} return (result);

$end_generic