#ifndef STAC_H
#define STAC_H

#include "pslib.h"
#define abs(x) ((x) >= 0 ? (x) : (-(x)))
#define MAXCHAR 80

/************************************************************************************************************/

// time.c

// Get the current time
double getTime(void);

/************************************************************************************************************/

// stacRead.c

// Read the input files and return an array of images
psArray *stacReadImages(psArray **headers, // The image headers, to be returned
                        psArray *filenames // The file names of the images
    );

// Read a file of coordinates (x,y)
psArray *stacReadCoords(const char *filename);

// Read a single map file and return the transformation
psPlaneTransform *stacReadMap(const char *filename);

// Read the map files and return an array of transformations
psArray *stacReadMaps(psArray *filenames // The file names of the images whose maps are to be read
    );

/************************************************************************************************************/

// stacErrorImages.c

// Calculate the error images
psArray *stacErrorImages(psArray *inputs, // Array of input images
                         float gain,    // Gain, in e/ADU
                         float rn       // Read noise, in e
    );

/************************************************************************************************************/

// stacTransform.c

// Transform input images, return success
bool stacTransform(psArray **outputs,   // Transformed images for output
                   psArray **outErrors, // Transformed error images for output
                   const psArray *images, // Array of images to be transformed
                   const psArray *maps, // Array of polynomials that do the transformation
                   const psArray *errors, // Array of error images to be transformed
                   const psArray *masks, // Masks of input images
                   const psImage *region, // Region of interest for transformation
                   const psVector *scales, // Relative scales
                   const psVector *offsets, // Relative offsets
                   int outnx, int outny // Size of output images
    );

/************************************************************************************************************/

// stacCheckMemory.c

// Check memory
void stacCheckMemory(void);

// Print out the problem
void stacMemoryProblem(const psMemBlock* ptr, ///< the pointer to the problematic memory block.
                       const char *file, ///< the file in which the problem originated
                       psS32 lineno     ///< the line number in which the problem originated
    );

// Print out a memblock when it's allocated --- this function used as a callback
psMemId stacMemPrint(const psMemBlock *ptr);


/************************************************************************************************************/

// stacCombine.c

// Get the mean for a bunch of values
float stacCombineMean(psVector *values, // Values for which to take the mean
                      psVector *errors, // Errors in the values
                      psVector *masks   // Masks for the values, 0 = don't use, 1 = use
    );

// Get the median for a bunch of values
float stacCombineMedian(psVector *values, // Values for which to take the median
                        psVector *errors, // Errors in the values
                        psVector *masks // Masks for the values, 0 = don't use, 1 = use
    );

// Combine the transformed images
bool stacCombine(psImage **combined,    // The combined image for output
                 psArray **rejected,    // Array of rejection masks
                 psArray *images,       // Array of transformed images
                 psArray *errors,       // Array of transformed error images
                 int nReject,           // Number of rejection iterations
                 psImage *region,       // Region to combine
                 psVector *saturated,   // Saturation limits for each image
                 psVector *bad,         // Bad pixel limits for each image
                 float reject           // Rejection (k-sigma)
    );

/************************************************************************************************************/

// stacInvertMaps.c

// Invert an array of maps
psArray *stacInvertMaps(const psArray *maps, // Array of maps to invert
                        int outnx, int outny // Size of output image
    );

/************************************************************************************************************/

// stacRejection.c

// Return the relative gradient for a given pixel
float stacGradient(psImage *image,      // Input for which to measure the gradient
                   int x, int y         // Coordinates at which to measure the gradient
    );

// Transform the rejection masks back to the source frame
psArray *stacRejection(psArray *inputs, // Input images
                       psArray *rejected, // Rejected images
                       psArray *regions, // Regions of interest
                       psArray *maps,   // Maps from input to transformed image
                       psArray *inverseMaps, // Maps from transformed to input image
                       float frac,      // Fraction of pixel rejected before looking more carefully
                       float grad,      // Gradient limit for rejection
                       psArray *names   // Names of original images (only for writing out when TESTING)
    );

/************************************************************************************************************/

// stacAreaOfInterest.c

// Returns the change in x' and y' for a change in x and y of 1
bool stacPlaneTransformDeriv(psPlane *out, // Output derivative, assumed already allocated
                             psPlaneTransform *transform, // The transform for which to obtain the derivative
                             psPlane *in // The position at which to obtain the derivative
    );

// Return a mask image designating the region of interest in the source frame
// Basically, uses derivatives of the transformation to work out what pixels in the source might be of interest
psImage *stacAreaOfInterest(psImage *mask, // Mask that is to be used to determine the area of interest
                            psPlaneTransform *map, // Map from the mask to the area of interest
                            int nxOut, int nyOut // Size of the area of interest
    );

/************************************************************************************************************/

// stacSize.c

// Calculate the size of the output image
bool stacSize(int *outnx, int *outny,   // Size of output image (returned)
              float *xMapDiff, float *yMapDiff, // Difference applied to maps
              const psVector *xSizes,   // Sizes of images in x
              const psVector *ySizes,   // Sizes of images in y
              psArray *maps             // Transformation maps
    );

/************************************************************************************************************/

// stacScale.c

// Calculate the background in an image
float stacBackground(const psImage *image, // Image for which to get the background
                     int sample         // Sample in increments of this value
    );


// Calculate the relative scales and offsets between the images
bool stacScales(psVector **scalesPtr,   // Scales to return
                psVector **offsetsPtr,  // Offsets to return
                const psArray *images,  // Images on which to measure the scales and offsets
                const char *starFile,   // File containing coordinates to photometer
                const char *starMapFile, // Map for coodinates to the common output frame
                float xMapDiff, float yMapDiff, // Difference from the map to apply (due to size difference)
                float aper              // Aperture to use for photometry (radius)
    );

// Rescale images
bool stacRescale(psArray *images,       // Images to rescale
                 psArray *errImages,    // Variance images to rescale
                 const psImage *mask,   // Mask indicating which pixels to scale
                 const psVector *scales,// Scales for images
                 const psVector *offsets // Offsets for images
    );



/************************************************************************************************************/

// stacWrite.c

// Write map out
bool stacWriteMap(const char *mapName,  // Filename to write to
                  psPlaneTransform *map // Map to write
    );

// Write multiple maps
bool stacWriteMaps(const psArray *names, // Filenames of the input images (will add ".map")
                   const psArray *maps  // Maps to write
    );

/************************************************************************************************************/

#endif
