/** @file  ImageSubtract.h
 *
 *  This file will ...
 *
 *  @author GLG, MHPCC
 *
 *  @version $Revision: 1.1 $ $Name: not supported by cvs2svn $
 *  @date $Date: 2005-05-25 23:00:50 $
 *
 *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
 *
 */

#if !defined(PM_IMAGE_COMBINE_H)
#define PM_IMAGE_COMBINE_H

#if HAVE_CONFIG_H
#include <config.h>
#endif

#include<stdio.h>
#include<math.h>
#include "pslib.h"
#include "psConstants.h"

typedef enum {
    PM_SUBTRACTION_KERNEL_POIS,    ///< POIS kernel --- delta functions
    PM_SUBTRACTION_KERNEL_ISIS    ///< ISIS kernel --- gaussians modified by polynomials
} pmSubtractionKernelsType;

typedef struct
{
    pmSubtractionKernelsType type;   ///< Type ofKernels --- allowing the use of multiple kernels
    psVector *u, *v;     ///< Offset (for POIS) or polynomial order (for ISIS)
    psVector *sigma;     ///< Width of Gaussian (for ISIS)
    psVector *xOrder, *yOrder;    ///< Spatial Polynomial order (for all)
    int subIndex;     ///< Index of kernel to be subtracted to maintain flux conservation
    psArray *preCalc;     ///< Array of images containing pre-calculated kernel (to
    ///< accelerate ISIS; don't use for POIS)
    psS32 p_size;                                       ///< The halfsize of the kernel
    psS32 p_spatialOrder;                               ///< The spatial order of the kernels
}
psSubtractionKernels;

psSubtractionKernels *pmSubtractionKernelsAllocPOIS(int size,
        int SpatialOrder);

psSubtractionKernels *pmSubtractionKernelsAllocISIS(const psVector *sigmas,
        const psVector *orders,
        int size,
        int SpatialOrder);

typedef enum {
    PM_STAMP_INIT,     ///< Initial state
    PM_STAMP_USED,     ///< Use this stamp
    PM_STAMP_REJECTED,     ///< This stamp has been rejected
    PM_STAMP_RECALC,     ///< Having been reset, this stamp is to be recalculated
    PM_STAMP_NONE     ///< No stamp in this region
} pmStampStatus;

typedef struct
{
    int x, y;      ///< Position
    psImage *matrix;     ///< Associated matrix
    psVector *vector;     ///< Assoicated vector
    pmStampStatus status;    ///< Status ofstamp
}
pmStamp;

psArray *pmSubtractionFindStamps(psArray *stamps, ///< Output stamps, or NULL
                                 const psImage *image, ///< Image for which to find stamps
                                 const psImage *mask, ///< Mask
                                 psU32 maskVal,  ///< Value for mask
                                 psF32 threshold, ///< Threshold for stamps in the image
                                 psS32 xNum,  ///< Number of stamps in x
                                 psS32 yNum,  ///< Number of stamps in y
                                 psS32 border  ///< Border around image to ignore (should be size of kernel)
                                );

bool pmSubtractionCalculateEquation(psArray *stamps, ///< The stamps for which to calculate the equation,
                                    const psImage *reference, ///< Reference image
                                    const psImage *input, ///< Input image
                                    const psSubtractionKernels *kernels, ///< The kernel basis functions
                                    psS32 footprint ///< Half-size of region over which to calculate equation
                                   );


psVector *pmSubtractionSolveEquation(psVector *solution, ///< Solution vector, or NULL
                                     const psArray *stamps ///< Array of stamps
                                    );

bool pmSubtractionRejectStamps(psArray *stamps,  ///< Array of stamps to check for rejection
                               psImage *mask,  ///< Mask image
                               psU32 badStampMaskVal, ///< Value to use in mask for bad stamp
                               psS32 footprint,  ///< Region to mask if stamp is bad
                               psF32 sigmaRej,  ///< Number of RMS deviations above zero at which to reject
                               const psImage *refImage, ///< Reference image
                               const psImage *inImage, ///< Input image
                               const psVector *solution, ///< Solution vector
                               const psSubtractionKernels *kernels ///< Array of kernel parameters
                              );

psImage *pmSubtractionKernelImage(psImage *out,
                                  const psVector *solution,
                                  const psSubtractionKernels *kernels,
                                  psF32 x,
                                  psF32 y
                                 );

#endif
