#include "VsipTypes.h"

/****************************************************************************
 *
 *  Function    : VsipConvolution
 *
 *  Description : Computes the linear convolution with a nxn kernel
 *                Kernel sizes allowed are 3x3, 5x5, 7x7
 *                The difference between the output (dest) size and
 *                input (src) size must be 0 or +/- (n - 1)
 *
 *  Parameters  : src - input image object
 *                kernel - nxn kernel object
 *                rsz - row size of the output image
 *                csz - column size of the output image
 *
 *  Returns     : convolved image
 *
 ****************************************************************************/
$generic

  $VsipConvolution = 'VsipConvolutionBit',    $InType = 'uint1',  $OutType = 'int32',
       $Per = 'AuxPerimeterBit';
  $VsipConvolution = 'VsipConvolutionByte',   $InType = 'uint8',  $OutType = 'uint8',
       $Per = 'AuxPerimeterByte';
  $VsipConvolution = 'VsipConvolutionInt',    $InType = 'int32',  $OutType = 'int32',
       $Per = 'AuxPerimeterInt';
  $VsipConvolution = 'VsipConvolutionFloat',  $InType = 'float',  $OutType = 'float',
       $Per = 'AuxPerimeterFloat';
  $VsipConvolution = 'VsipConvolutionDouble', $InType = 'double', $OutType = 'double',
       $Per = 'AuxPerimeterDouble';

$in

$OutType[:,:] $VsipConvolution($InType src[:,:], $OutType kernel[:,:],
                             VsipIndexesType rsz, VsipIndexesType csz)
{

  // ***** recover the src and kernel sizes *****
  VsipIndexesType rs, VsipIndexesType cs = extents(src);
  VsipIndexesAuxType rk, VsipIndexesAuxType ck = extents(kernel);

  // ***** test the sizes of the image, kernel and dest (see header)  *****
  assert(((rs > 0) && (cs > 0)),
         "ERROR: Image can not have zero dimension. (",rs,"x",cs,")\n");
  assert((((rk == 3) && (ck == 3))  || ((rk == 5) && (ck == 5))  ||
          ((rk == 7) && (ck == 7))),
         "ERROR: Invalid Kernel size (", rk, "'", ck, ").\n");
  assert((((rsz == rs) && (csz == cs)) ||
          ((rsz == rs - (rk - 1)) && (csz == cs - (ck - 1))) ||
          ((rsz == rs + (rk - 1)) && (csz == cs + (ck - 1)))),
         "ERROR: Invalid destination size (",rsz,",",csz,").\n");

  // ***** reflects the kernel over the origin *****
  $OutType invkernel[:,:] =
       for elem in kernel at (VsipIndexesType i, VsipIndexesType j)
       return(array(kernel[rk - 1 - i, ck - 1 - j]));

  // ***** creates the necessary perimeter around src *****
  $OutType persrc[:,:] =
       if (rsz == rs - (rk - 1))
          return(src)
       elif (rsz == rs)
          return($Per(src, rk / 2, rk / 2, rk / 2, rk / 2, 0))
       else
          return($Per(src, rk - 1, rk - 1, rk - 1, rk - 1, 0));
 

  // ***** computes the convolution for the modified src *****
  $OutType result[:,:] =
       forwindow win[rk,ck] in persrc
            { $OutType conv =
                 for elem1 in win dot elem2 in invkernel
                 return(sum(elem1*elem2));
            }
       return(array(conv));

} return (result);

$end_generic

/****************************************************************************
 *
 *  Function    : VsipConvolutionCx
 *
 *  Description : Computes the linear convolution with a nxn kernel
 *                Kernel sizes allowed are 3x3, 5x5, 7x7
 *                The difference between the output (dest) size and
 *                input (src) size must be 0 or +/- (n - 1)
 *
 *  Parameters  : src - input image object
 *                kernel - nxn kernel object
 *                rsz - row size of the output image
 *                csz - column size of the output image
 *
 *  Returns     : convolved image
 *
 ****************************************************************************/
$generic

  $VsipConvolutionCx = 'VsipConvolutionCxInt',  $InType = 'complex int32',
              $OutType = 'complex int32', $Per = 'AuxPerimeterCxInt';
  $VsipConvolutionCx = 'VsipConvolutionCxFloat',$InType = 'complexfloat',
              $OutType = 'complexfloat', $Per = 'AuxPerimeterCxFloat';

$in

$OutType[:,:] $VsipConvolutionCx($InType src[:,:], $OutType kernel[:,:],
                               VsipIndexesType rsz, VsipIndexesType csz)
{

  // ***** recover the src and kernel sizes *****
  VsipIndexesType rs, VsipIndexesType cs = extents(src);
  VsipIndexesAuxType rk, VsipIndexesAuxType ck = extents(kernel);

  // ***** test the sizes of the image, kernel and dest (see header)  *****
  assert(((rs > 0) && (cs > 0)),
         "ERROR: Image can not have zero dimension. (",rs,"x",cs,")\n");
  assert((((rk == 3) && (ck == 3))  || ((rk == 5) && (ck == 5))  ||
          ((rk == 7) && (ck == 7))),
         "ERROR: Invalid Kernel size (", rk, "'", ck, ").\n");
  assert((((rsz == rs) && (csz == cs)) ||
          ((rsz == rs - (rk - 1)) && (csz == cs - (ck - 1))) ||
          ((rsz == rs + (rk - 1)) && (csz == cs + (ck - 1)))),
         "ERROR: Invalid destination size (",rsz,",",csz,").\n");

  // ***** reflects the kernel over the origin *****
  $OutType invkernel[:,:] =
       for elem in kernel at (VsipIndexesType i, VsipIndexesType j)
       return(array(kernel[rk - 1 - i, ck - 1 - j]));

  // ***** creates the necessary perimeter around src *****
  $OutType persrc[:,:] =
       if (rsz == rs - (rk - 1))
          return(src)
       elif (rsz == rs)
          return($Per(src, rk / 2, rk / 2, rk / 2, rk / 2, (0,0)))
       else
          return($Per(src, rk - 1, rk - 1, rk - 1, rk - 1, (0,0)));
 

  // ***** computes the convolution for the modified src *****
  $OutType result[:,:] =
       forwindow win[rk,ck] in persrc
            { $OutType conv =
                 for elem1 in win dot elem2 in invkernel
                 return(sum(elem1*elem2));
            }
       return(array(conv));

} return (result);

$end_generic