/** @file  psImage.c
 *
 *  @brief Contains basic image definitions and operations.
 *
 *  This file defines the basic type for an image struct and functions useful
 *  in manupulating images.
 *
 *  @author Robert DeSonia, MHPCC
 *  @author Ross Harman, MHPCC
 *
 *  @version $Revision: 1.61 $ $Name: not supported by cvs2svn $
 *  @date $Date: 2005-03-18 02:35:14 $
 *
 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
 *
 *  That is the routine used to generate matrices.
 */

#include <string.h>
#include <math.h>

#include "psMemory.h"
#include "psError.h"
#include "psImage.h"
#include "psString.h"

#include "psImageErrors.h"

#define SQUARE(x) ((x)*(x))
#define MIN(x,y) (((x) > (y)) ? (y) : (x))
#define MAX(x,y) (((x) > (y)) ? (x) : (y))

static void imageFree(psImage* image)
{
    if (image == NULL) {
        return;
    }

    if (image->type.type == PS_TYPE_PTR) {
        // 2-D array of pointers -- must dereference elements
        psU32 oldNumRows = image->numRows;
        psU32 oldNumCols = image->numCols;
        psPtr* rowPtr;

        for (psU32 row = 0; row < oldNumRows; row++) {
            rowPtr = image->data.PTR[row];
            for (psU32 col = 0; col < oldNumCols; col++) {
                psMemDecrRefCounter(rowPtr[col]);
            }
        }
    }

    if (image->parent != NULL) {
        psArrayRemove(image->parent->children,image);
        image->parent = NULL;
    }

    psImageFreeChildren(image);

    psFree(image->rawDataBuffer);
    psFree(image->data.V);
}

psImage* psImageAlloc(psU32 numCols,
                      psU32 numRows,
                      const psElemType type)
{
    psS32 area = 0;
    psS32 elementSize = PSELEMTYPE_SIZEOF(type);  // element size in bytes
    psS32 rowSize = numCols * elementSize;        // row size in bytes.

    area = numCols * numRows;

    if (area < 1) {
        psError(PS_ERR_BAD_PARAMETER_VALUE, true,
                PS_ERRORTEXT_psImage_AREA_NEGATIVE,
                numRows, numCols);
        return NULL;
    }

    psImage* image = (psImage* ) psAlloc(sizeof(psImage));

    p_psMemSetDeallocator(image, (psFreeFcn) imageFree);

    image->data.V = psAlloc(sizeof(psPtr ) * numRows);

    image->rawDataBuffer = psAlloc(area * elementSize);

    // set the row pointers.
    image->data.V[0] = image->rawDataBuffer;
    for (psS32 i = 1; i < numRows; i++) {
        image->data.V[i] = (psPtr )((int8_t *) image->data.V[i - 1] + rowSize);
    }

    *(psS32 *)&image->col0 = 0;
    *(psS32 *)&image->row0 = 0;
    *(psU32 *)&image->numCols = numCols;
    *(psU32 *)&image->numRows = numRows;
    *(psDimen* ) & image->type.dimen = PS_DIMEN_IMAGE;
    *(psElemType* ) & image->type.type = type;
    image->parent = NULL;
    image->children = NULL;

    return image;
}

psRegion* psRegionAlloc(psF32 x0,
                        psF32 x1,
                        psF32 y0,
                        psF32 y1)
{
    psRegion* out = psAlloc(sizeof(psRegion));

    out->x0 = x0;
    out->y0 = y0;
    out->x1 = x1;
    out->y1 = y1;

    return out;
}

psRegion* psRegionFromString(char* region)
{
    psS32 col0;
    psS32 col1;
    psS32 row0;
    psS32 row1;

    // section should be of the form '[col0:col1,row0:row1]'
    if (region == NULL) {
        psError(PS_ERR_BAD_PARAMETER_NULL, true,
                PS_ERRORTEXT_psImage_SUBSECTION_NULL);
        return NULL;
    }

    if (sscanf(region,"[%d:%d,%d:%d]",&col0,&col1,&row0,&row1) < 4) {
        psError(PS_ERR_BAD_PARAMETER_NULL, true,
                PS_ERRORTEXT_psImage_SUBSECTION_INVALID,
                region);
        return NULL;
    }

    if (col0 > col1 || row0 > row1) {
        psError(PS_ERR_BAD_PARAMETER_VALUE, true,
                PS_ERRORTEXT_psImage_SUBSET_RANGE_MALFORMED,
                col0,col1,row0,row1);
        return NULL;
    }

    return psRegionAlloc(col0,col1,row0,row1);
}

char* psRegionToString(psRegion* region)
{
    char tmpText[256]; // big enough to store any region as text

    if (region == NULL) {
        psError(PS_ERR_BAD_PARAMETER_NULL, true,
                PS_ERRORTEXT_psImage_REGION_NULL);
        return NULL;
    }

    snprintf(tmpText,256,"[%g:%g,%g:%g]",
             region->x0, region->x1,
             region->y0, region->y1);

    return psStringCopy(tmpText);
}

psImage* psImageRecycle(psImage* old,
                        psU32 numCols,
                        psU32 numRows,
                        const psElemType type)
{
    psS32 elementSize = PSELEMTYPE_SIZEOF(type);  // element size in bytes
    psS32 rowSize = numCols * elementSize;        // row size in bytes.

    if (old == NULL) {
        old = psImageAlloc(numCols, numRows, type);
        return old;
    }

    if (old->type.dimen != PS_DIMEN_IMAGE) {
        psFree(old);
        psError(PS_ERR_BAD_PARAMETER_TYPE, true,
                PS_ERRORTEXT_psImage_NOT_AN_IMAGE);
        return NULL;
    }

    if (old->type.type == PS_TYPE_PTR) {
        // 2-D array of pointers -- must
        // dereference
        psU32 oldNumRows = old->numRows;
        psU32 oldNumCols = old->numCols;
        psPtr* rowPtr;

        for (psU32 row = 0; row < oldNumRows; row++) {
            rowPtr = old->data.PTR[row];
            for (psU32 col = 0; col < oldNumCols; col++) {
                psMemDecrRefCounter(rowPtr[col]);
                rowPtr[col] = NULL;
            }
        }
    }

    /* image already the right size/type? */
    if (numCols == old->numCols && numRows == old->numRows &&
            type == old->type.type) {
        return old;
    }
    // Resize the image buffer
    old->rawDataBuffer = psRealloc(old->data.V[0],
                                   numCols * numRows * elementSize);
    old->data.V = (psPtr *)psRealloc(old->data.V, numRows * sizeof(psPtr ));

    // recreate the row pointers
    old->data.V[0] = old->rawDataBuffer;
    for (psS32 i = 1; i < numRows; i++) {
        old->data.V[i] = (psPtr )((int8_t *) old->data.V[i - 1] + rowSize);
    }

    *(psU32 *)&old->numCols = numCols;
    *(psU32 *)&old->numRows = numRows;
    *(psElemType* ) & old->type.type = type;

    return old;
}

psImage* psImageCopy(psImage* restrict output,
                     const psImage* input,
                     psElemType type)
{
    psElemType inDatatype;
    psS32 elementSize;
    psS32 elements;
    psS32 numRows;
    psS32 numCols;

    if (input == NULL || input->data.V == NULL) {
        psError(PS_ERR_BAD_PARAMETER_NULL, true,
                PS_ERRORTEXT_psImage_IMAGE_NULL);
        psFree(output);
        return NULL;
    }

    if (input == output) {
        psError(PS_ERR_BAD_PARAMETER_NULL, true,
                PS_ERRORTEXT_psImage_INPLACE_NOTSUPPORTED);
        psFree(output);
        return NULL;
    }

    if (input->type.dimen != PS_DIMEN_IMAGE) {
        psError(PS_ERR_BAD_PARAMETER_TYPE, true,
                PS_ERRORTEXT_psImage_NOT_AN_IMAGE);
        psFree(output);
        return NULL;
    }

    inDatatype = input->type.type;
    numRows = input->numRows;
    numCols = input->numCols;
    elements = numRows * numCols;
    elementSize = PSELEMTYPE_SIZEOF(inDatatype);

    if (inDatatype == PS_TYPE_PTR || type == PS_TYPE_PTR) {
        psError(PS_ERR_BAD_PARAMETER_TYPE, true,
                PS_ERRORTEXT_psImage_IMAGE_TYPE_UNSUPPORTED,
                PS_TYPE_PTR_NAME);
        psFree(output);
        return NULL;
    }

    output = psImageRecycle(output, numCols, numRows, type);

    // cover the trival case of copy of the same
    // datatype.
    if (type == inDatatype) {
        for (psS32 row=0;row<numRows;row++) {
            memcpy(output->data.V[row], input->data.V[row], elementSize * numCols);
        }
        return output;
    }

    #define PSIMAGE_ELEMENT_COPY(IN,INTYPE,OUT,OUTTYPE,ELEMENTS) { \
        ps##INTYPE *in; \
        ps##OUTTYPE *out; \
        for(psS32 row=0;row<numRows;row++) { \
            in = IN->data.INTYPE[row]; \
            out = OUT->data.OUTTYPE[row]; \
            for (psS32 col=0;col<numCols;col++) { \
                *(out++) = *(in++); \
            } \
        } \
    }

    #define PSIMAGE_COPY_CASE(OUT,OUTTYPE) { \
        switch (inDatatype) { \
        case PS_TYPE_S8: \
            PSIMAGE_ELEMENT_COPY(input,S8,OUT,OUTTYPE,elements); \
            break; \
        case PS_TYPE_S16: \
            PSIMAGE_ELEMENT_COPY(input,S16,OUT,OUTTYPE,elements); \
            break; \
        case PS_TYPE_S32: \
            PSIMAGE_ELEMENT_COPY(input,S32,OUT,OUTTYPE,elements); \
            break; \
        case PS_TYPE_S64: \
            PSIMAGE_ELEMENT_COPY(input,S64,OUT,OUTTYPE,elements); \
            break; \
        case PS_TYPE_U8: \
            PSIMAGE_ELEMENT_COPY(input,U8,OUT,OUTTYPE,elements); \
            break; \
        case PS_TYPE_U16: \
            PSIMAGE_ELEMENT_COPY(input,U16,OUT,OUTTYPE,elements); \
            break; \
        case PS_TYPE_U32: \
            PSIMAGE_ELEMENT_COPY(input,U32,OUT,OUTTYPE,elements); \
            break; \
        case PS_TYPE_U64: \
            PSIMAGE_ELEMENT_COPY(input,U64,OUT,OUTTYPE,elements); \
            break; \
        case PS_TYPE_F32: \
            PSIMAGE_ELEMENT_COPY(input,F32,OUT,OUTTYPE,elements); \
            break; \
        case PS_TYPE_F64: \
            PSIMAGE_ELEMENT_COPY(input,F64,OUT,OUTTYPE,elements); \
            break; \
        case PS_TYPE_C32: \
            PSIMAGE_ELEMENT_COPY(input,C32,OUT,OUTTYPE,elements); \
            break; \
        case PS_TYPE_C64: \
            PSIMAGE_ELEMENT_COPY(input,C64,OUT,OUTTYPE,elements); \
            break; \
        default: \
            break; \
        } \
    }

    switch (type) {
    case PS_TYPE_S8:
        PSIMAGE_COPY_CASE(output, S8);
        break;
    case PS_TYPE_S16:
        PSIMAGE_COPY_CASE(output, S16);
        break;
    case PS_TYPE_S32:
        PSIMAGE_COPY_CASE(output, S32);
        break;
    case PS_TYPE_S64:
        PSIMAGE_COPY_CASE(output, S64);
        break;
    case PS_TYPE_U8:
        PSIMAGE_COPY_CASE(output, U8);
        break;
    case PS_TYPE_U16:
        PSIMAGE_COPY_CASE(output, U16);
        break;
    case PS_TYPE_U32:
        PSIMAGE_COPY_CASE(output, U32);
        break;
    case PS_TYPE_U64:
        PSIMAGE_COPY_CASE(output, U64);
        break;
    case PS_TYPE_F32:
        PSIMAGE_COPY_CASE(output, F32);
        break;
    case PS_TYPE_F64:
        PSIMAGE_COPY_CASE(output, F64);
        break;
    case PS_TYPE_C32:
        PSIMAGE_COPY_CASE(output, C32);
        break;
    case PS_TYPE_C64:
        PSIMAGE_COPY_CASE(output, C64);
        break;
    default: {
            char* typeStr;
            PS_TYPE_NAME(typeStr,type);
            psError(PS_ERR_BAD_PARAMETER_TYPE, true,
                    PS_ERRORTEXT_psImage_IMAGE_TYPE_UNSUPPORTED,
                    typeStr);
            psFree(output);

            break;
        }
    }
    return output;
}

psS32 psImageFreeChildren(psImage* image)
{
    psS32 numFreed = 0;

    if (image == NULL) {
        return numFreed;
    }

    if (image->children != NULL) {
        psImage** children = (psImage**)image->children->data;
        numFreed = image->children->n;

        // orphan the children first
        // (so psFree doesn't try to modify the parent's children array while I'm using it)
        for (psS32 i=0;i<numFreed;i++) {
            children[i]->parent = NULL;
        }

        psFree(image->children);
        image->children = NULL;
    }

    return numFreed;
}

psC64 psImagePixelInterpolate(const psImage* input,
                              float x,
                              float y,
                              const psImage* mask,
                              psU32 maskVal,
                              psC64 unexposedValue,
                              psImageInterpolateMode mode)
{

    if (input == NULL) {
        psError(PS_ERR_BAD_PARAMETER_NULL,true,
                PS_ERRORTEXT_psImage_IMAGE_NULL);
        return unexposedValue;
    }

    #define PSIMAGE_PIXEL_INTERPOLATE_CASE(TYPE)                             \
case PS_TYPE_##TYPE:                                                 \
    switch (mode) {                                                  \
    case PS_INTERPOLATE_FLAT:                                        \
        return p_psImagePixelInterpolateFLAT_##TYPE(                 \
                input,                                               \
                x,                                                   \
                y,                                                   \
                mask,                                                \
                maskVal,                                             \
                unexposedValue);                                     \
        break;                                                       \
    case PS_INTERPOLATE_BILINEAR:                                    \
        return p_psImagePixelInterpolateBILINEAR_##TYPE(             \
                input,                                               \
                x,                                                   \
                y,                                                   \
                mask,                                                \
                maskVal,                                             \
                unexposedValue);                                     \
        break;                                                       \
    case PS_INTERPOLATE_BILINEAR_VARIANCE:                           \
        return p_psImagePixelInterpolateBILINEAR_VARIANCE_##TYPE(    \
                input,                                               \
                x,                                                   \
                y,                                                   \
                mask,                                                \
                maskVal,                                             \
                unexposedValue);                                     \
        break;                                                       \
    default:                                                         \
        psError(PS_ERR_BAD_PARAMETER_VALUE,true,                     \
                PS_ERRORTEXT_psImage_INTERPOLATE_METHOD_INVALID,     \
                mode);                                               \
    }                                                                \
    break;

    switch (input->type.type) {
        PSIMAGE_PIXEL_INTERPOLATE_CASE(U8);
        PSIMAGE_PIXEL_INTERPOLATE_CASE(U16);
        PSIMAGE_PIXEL_INTERPOLATE_CASE(U32);
        PSIMAGE_PIXEL_INTERPOLATE_CASE(U64);
        PSIMAGE_PIXEL_INTERPOLATE_CASE(S8);
        PSIMAGE_PIXEL_INTERPOLATE_CASE(S16);
        PSIMAGE_PIXEL_INTERPOLATE_CASE(S32);
        PSIMAGE_PIXEL_INTERPOLATE_CASE(S64);
        PSIMAGE_PIXEL_INTERPOLATE_CASE(F32);
        PSIMAGE_PIXEL_INTERPOLATE_CASE(F64);
        PSIMAGE_PIXEL_INTERPOLATE_CASE(C32);
        PSIMAGE_PIXEL_INTERPOLATE_CASE(C64);
    default: {
            char* typeStr;
            PS_TYPE_NAME(typeStr,input->type.type);
            psError(PS_ERR_BAD_PARAMETER_TYPE,true,
                    PS_ERRORTEXT_psImage_IMAGE_TYPE_UNSUPPORTED,
                    typeStr);
        }
    }

    return unexposedValue;
}

#define PSIMAGE_PIXEL_INTERPOLATE_FLAT(TYPE,RETURNTYPE) \
inline RETURNTYPE p_psImagePixelInterpolateFLAT_##TYPE( \
        const psImage* input, \
        float x, \
        float y, \
        const psImage* mask, \
        psU32 maskVal, \
        RETURNTYPE unexposedValue) \
{ \
    psS32 intX = (psS32) round((psF64)(x) - 0.5 + FLT_EPSILON); \
    psS32 intY = (psS32) round((psF64)(y) - 0.5 + FLT_EPSILON); \
    psS32 lastX = input->numCols - 1; \
    psS32 lastY = input->numRows - 1; \
    \
    if ((intX < 0) || \
            (intX > lastX) || \
            (intY < 0) || \
            (intY > lastY) || \
            ( (mask!=NULL) && \
              ((mask->data.PS_TYPE_MASK_DATA[intY][intX] & maskVal) != 0) ) ) { \
        return unexposedValue; \
    } \
    \
    return input->data.TYPE[intY][intX]; \
}

PSIMAGE_PIXEL_INTERPOLATE_FLAT(U8,psF64)
PSIMAGE_PIXEL_INTERPOLATE_FLAT(U16,psF64)
PSIMAGE_PIXEL_INTERPOLATE_FLAT(U32,psF64)
PSIMAGE_PIXEL_INTERPOLATE_FLAT(U64,psF64)
PSIMAGE_PIXEL_INTERPOLATE_FLAT(S8,psF64)
PSIMAGE_PIXEL_INTERPOLATE_FLAT(S16,psF64)
PSIMAGE_PIXEL_INTERPOLATE_FLAT(S32,psF64)
PSIMAGE_PIXEL_INTERPOLATE_FLAT(S64,psF64)
PSIMAGE_PIXEL_INTERPOLATE_FLAT(F32,psF64)
PSIMAGE_PIXEL_INTERPOLATE_FLAT(F64,psF64)
PSIMAGE_PIXEL_INTERPOLATE_FLAT(C32,psC64)
PSIMAGE_PIXEL_INTERPOLATE_FLAT(C64,psC64)

#define PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(TYPE, RETURNTYPE, SUFFIX, FRACFUNC) \
inline RETURNTYPE p_psImagePixelInterpolateBILINEAR_##SUFFIX( \
        const psImage* input, \
        float x, \
        float y, \
        const psImage* mask, \
        psU32 maskVal, \
        RETURNTYPE unexposedValue) \
{ \
    double floorX = floor((psF64)(x) - 0.5); \
    double floorY = floor((psF64)(y) - 0.5); \
    psF64 fracX = x - 0.5 - floorX; \
    psF64 fracY = y - 0.5 - floorY; \
    psS32 intFloorX = (psS32) floorX; \
    psS32 intFloorY = (psS32) floorY; \
    psS32 lastX = input->numCols - 1; \
    psS32 lastY = input->numRows - 1; \
    ps##TYPE V00 = 0; \
    ps##TYPE V01 = 0; \
    ps##TYPE V10 = 0; \
    ps##TYPE V11 = 0; \
    psBool valid00 = false; \
    psBool valid01 = false; \
    psBool valid10 = false; \
    psBool valid11 = false; \
    \
    if (intFloorY >= 0 && intFloorY <= lastY) { \
        if (intFloorX >= 0 && intFloorX <= lastX) { \
            V00 = input->data.TYPE[intFloorY][intFloorX]; \
            valid00 = (mask == NULL) || \
                      ((mask->data.PS_TYPE_MASK_DATA[intFloorY][intFloorX] & maskVal) == 0); \
        } \
        if (intFloorX >= -1 && intFloorX < lastX) { \
            V10 = input->data.TYPE[intFloorY][intFloorX+1]; \
            valid10 = (mask == NULL) || \
                      ((mask->data.PS_TYPE_MASK_DATA[intFloorY][intFloorX+1] & maskVal) == 0); \
        } \
    } \
    if (intFloorY >= -1 && intFloorY < lastY) { \
        if (intFloorX >= 0 && intFloorX <= lastX) { \
            V01 = input->data.TYPE[intFloorY+1][intFloorX]; \
            valid01 = (mask == NULL) || \
                      ((mask->data.PS_TYPE_MASK_DATA[intFloorY+1][intFloorX] & maskVal) == 0); \
        } \
        if (intFloorX >= -1 && intFloorX < lastX) { \
            V11 = input->data.TYPE[intFloorY+1][intFloorX+1]; \
            valid11 = (mask == NULL) || \
                      ((mask->data.PS_TYPE_MASK_DATA[intFloorY+1][intFloorX+1] & maskVal) == 0); \
        } \
    } \
    \
    /* cover likely case of all pixels being valid more efficiently */  \
    if (valid00 && valid10 && valid01 && valid11) { \
        /* formula from the ADD */ \
        return V00*FRACFUNC((1.0-fracX)*(1.0-fracY)) + V10*FRACFUNC(fracX*(1.0-fracY)) + \
               V01*FRACFUNC(fracY*(1.0-fracX)) + V11*FRACFUNC(fracX*fracY); \
    } \
    \
    /* OK, at least one pixel is not valid - need to do it piecemeal */ \
    \
    RETURNTYPE V0 = 0.0; \
    psBool valid0 = true; \
    if (valid00 && valid10) { \
        V0 = V00*FRACFUNC(1-fracX)+V10*FRACFUNC(fracX); \
    } else if (valid00) { \
        V0 = V00; \
    } else if (valid10) { \
        V0 = V10; \
    } else { \
        valid0 = false; \
    } \
    \
    RETURNTYPE V1 = 0.0; \
    psBool valid1 = true; \
    if (valid01 && valid11) { \
        V1 = V01*FRACFUNC(1-fracX)+V11*FRACFUNC(fracX); \
    } else if (valid01) { \
        V1 = V01; \
    } else if (valid11) { \
        V1 = V11; \
    } else { \
        valid1 = false; \
    } \
    \
    if (valid0 && valid1) { \
        return V0*FRACFUNC(1-fracY) + V1*FRACFUNC(fracY); \
    } else if (valid0) { \
        return V0; \
    } else if (valid1) { \
        return V1; \
    } \
    \
    return unexposedValue; \
}

PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(U8,psF64,U8,)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(U16,psF64,U16,)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(U32,psF64,U32,)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(U64,psF64,U64,)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(S8,psF64,S8,)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(S16,psF64,S16,)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(S32,psF64,S32,)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(S64,psF64,S64,)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(F32,psF64,F32,)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(F64,psF64,F64,)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(C32,psC64,C32,)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(C64,psC64,C64,)

// Variance Version
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(U8,psF64,VARIANCE_U8,SQUARE)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(U16,psF64,VARIANCE_U16,SQUARE)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(U32,psF64,VARIANCE_U32,SQUARE)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(U64,psF64,VARIANCE_U64,SQUARE)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(S8,psF64,VARIANCE_S8,SQUARE)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(S16,psF64,VARIANCE_S16,SQUARE)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(S32,psF64,VARIANCE_S32,SQUARE)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(S64,psF64,VARIANCE_S64,SQUARE)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(F32,psF64,VARIANCE_F32,SQUARE)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(F64,psF64,VARIANCE_F64,SQUARE)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(C32,psC64,VARIANCE_C32,SQUARE)
PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(C64,psC64,VARIANCE_C64,SQUARE)
