#include <stdio.h>
#include <strings.h>
#include <assert.h>
#include "pslib.h"

#include "pmFPA.h"
#include "pmFPARead.h"
#include "pmHDU.h"
#include "pmHDUUtils.h"
#include "pmConcepts.h"
#include "psRegionIsBad.h"
#include "pmHDUUtils.h"

//////////////////////////////////////////////////////////////////////////////////////////////////////////////
// File-static functions
//////////////////////////////////////////////////////////////////////////////////////////////////////////////

// Carve a readout from the image pixels
static pmReadout *readoutCarve(pmCell *cell, // Cell into which the new readout will be placed
                               psImage *image, // Image that will be carved
                               const psRegion *trimsec, // Trim section
                               const psList *biassecs // Bias sections
                              )
{
    pmReadout *readout = pmReadoutAlloc(cell);

    // The image corresponding to the trim region
    if (psRegionIsBad(*trimsec)) {
        psString regionString = psRegionToString(*trimsec);
        psError(PS_ERR_IO, true, "Invalid trim section: %s\n", regionString);
        psFree(regionString);
        psFree(readout);
        return NULL;
    }
    readout->image = psMemIncrRefCounter(psImageSubset(image, *trimsec));

    // Get the list of overscans
    psListIterator *iter = psListIteratorAlloc((psList*)biassecs, PS_LIST_HEAD, false); // Iterator
    psRegion *biassec = NULL;       // A BIASSEC region from the list
    while ((biassec = psListGetAndIncrement(iter))) {
        if (psRegionIsBad(*biassec)) {
            psString regionString = psRegionToString(*biassec);
            psError(PS_ERR_IO, true, "Invalid bias section: %s\n", regionString);
            psFree(regionString);
            psFree(readout);
            return NULL;
        }
        psImage *overscan = psMemIncrRefCounter(psImageSubset(image, *biassec));
        psListAdd(readout->bias, PS_LIST_TAIL, overscan);
        psFree(overscan);
    }
    psFree(iter);

    return readout;
}

//////////////////////////////////////////////////////////////////////////////////////////////////////////////
// Public functions
//////////////////////////////////////////////////////////////////////////////////////////////////////////////

// Read the next readout; return true if we read pixels in.
bool pmReadoutReadNext(pmReadout *readout, // Readout into which to read
                       psFits *fits,    // FITS file from which to read
                       int z,           // Readout number/plane; zero-offset indexing
                       int numRows      // The number of rows to read
                      )
{
    // Get the HDU and read the header
    pmCell *cell = readout->parent;     // The parent cell
    pmHDU *hdu = pmHDUFromCell(cell);   // The HDU
    if (!hdu) {
        return false;
    }

    if (!pmHDUReadHeader(hdu, fits)) {
        psError(PS_ERR_IO, false, "Unable to read HDU header!\n");
        return false;
    }

    // Check the third dimension
    bool mdok = true;                   // Status of MD lookup
    int naxis = psMetadataLookupS32(&mdok, hdu->header, "NAXIS"); // The number of axes
    if (!mdok) {
        psError(PS_ERR_IO, true, "Unable to find NAXIS in header for extension %s\n", hdu->extname);
        return false;
    }
    if (naxis == 0) {
        // No pixels to read, as for a PHU.
        return false;
    }
    if (naxis < 2 || naxis > 3) {
        psError(PS_ERR_IO, true, "NAXIS in header of extension %s (= %d) is not valid.\n",
                hdu->extname, naxis);
        return false;
    }
    int naxis3 = 1;                     // The number of image planes
    if (naxis == 3) {
        naxis3 = psMetadataLookupS32(&mdok, hdu->header, "NAXIS3");
        if (!mdok) {
            psError(PS_ERR_IO, true, "Unable to find NAXIS3 in header for extension %s\n", hdu->extname);
            return false;
        }
    }
    if (z >= naxis3) {
        // Nothing to see here.  Move along.
        return false;
    }

    // Get the size of the image plane
    int naxis1 = psMetadataLookupS32(&mdok, hdu->header, "NAXIS1"); // The number of columns
    if (!mdok) {
        psError(PS_ERR_IO, true, "Unable to find NAXIS1 in header for extension %s\n", hdu->extname);
        return false;
    }
    int naxis2 = psMetadataLookupS32(&mdok, hdu->header, "NAXIS2"); // The number of columns
    if (!mdok) {
        psError(PS_ERR_IO, true, "Unable to find NAXIS2 in header for extension %s\n", hdu->extname);
        return false;
    }

    // Read the FITS image
    int offset = readout->row0;         // The row number to read
    if (readout->image) {
        offset += readout->image->numRows;
    }
    if (offset >= naxis2) {
        // We've read everything there is
        return false;
    }
    int upper = offset + numRows;       // The upper limit for the pixel read
    if (upper > naxis2) {
        upper = naxis2;
    }
    psRegion region = psRegionSet(0, naxis1, offset, upper); // Region to be read
    psImage *image = psFitsReadImage(fits, region, z); // The image
    if (readout->image) {
        psFree(readout->image);
    }

    // Stick the image into the HDU and the readout
    if (!hdu->images) {
        hdu->images = psArrayAlloc(naxis3);
    }
    if (hdu->images->n != naxis3) {
        hdu->images = psArrayRealloc(hdu->images, naxis3);
    }
    if (hdu->images->data[z]) {
        psFree(hdu->images->data[z]);
    }
    hdu->images->data[z] = image;
    readout->image = psImageSubset(hdu->images->data[z], psRegionSet(0,0,0,0));
    readout->row0 = region.y0;

    return true;
}


// Read in the cell, and allocate the readouts
bool pmCellRead(pmCell *cell,           // Cell to read into
                psFits *fits,           // FITS file from which to read
                psDB *db                // Database handle, for "concepts" ingest
               )
{
    pmHDU *hdu = pmHDUFromCell(cell);   // The HDU
    if (!hdu) {
        return false;                    // Nothing to see here; move along
    }
    if (!hdu->images && !pmHDURead(hdu, fits)) {
        psError(PS_ERR_IO, false, "Unable to read HDU for cell.\n");
        return false;
    }

    pmConceptsReadCell(cell, PM_CONCEPT_SOURCE_HEADER, false, NULL);

    // Having read the cell, we now have to cut it up
    psRegion *trimsec = psMetadataLookupPtr(NULL, cell->concepts, "CELL.TRIMSEC");
    psList *biassecs = psMetadataLookupPtr(NULL, cell->concepts, "CELL.BIASSEC");

    // Iterate over each of the image planes
    for (int i = 0; i < hdu->images->n; i++) {
        psImage *image = hdu->images->data[i]; // The i-th plane

        // XXX: Type conversion here to support the modules, which don't have multiple type support yet
        if (image->type.type != PS_TYPE_F32) {
            psImage *temp = psImageCopy(NULL, image, PS_TYPE_F32); // Temporary image
            psFree(hdu->images->data[i]);
            hdu->images->data[i] = temp;
            image = temp;
        }

        pmReadout *readout = readoutCarve(cell, image, trimsec, biassecs);
        readout->mask = NULL;
        readout->weight = NULL;
        psFree(readout);                // Drop reference
    }

    pmCellSetDataStatus (cell, true);
    return true;
}

// Read in the component cells
bool pmChipRead(pmChip *chip,           // Chip to read into
                psFits *fits,           // FITS file from which to read
                psDB *db                // Database handle, for "concepts" ingest
               )
{
    bool success = false;               // Were we able to read at least one HDU?
    psArray *cells = chip->cells;       // Array of cells
    for (int i = 0; i < cells->n; i++) {
        pmCell *cell = cells->data[i];  // The cell of interest
        success |= pmCellRead(cell, fits, db);
    }
    if (success) {
        pmConceptsReadChip(chip, PM_CONCEPT_SOURCE_HEADER, false, NULL);
        // XXX probably could just use chip->data_exists
        pmChipSetDataStatus (chip, true);
    }

    return success;
}

// Read in the component chips
bool pmFPARead(pmFPA *fpa,              // FPA to read into
               psFits *fits,            // FITS file from which to read
               psDB *db                 // Database handle, for "concepts" ingest
              )
{
    bool success = false;               // Were we able to read at least one HDU?
    psArray *chips = fpa->chips;        // Array of chips
    for (int i = 0; i < chips->n; i++) {
        pmChip *chip = chips->data[i];  // The cell of interest
        success |= pmChipRead(chip, fits, db);
    }
    if (success) {
        pmConceptsReadFPA(fpa, PM_CONCEPT_SOURCE_HEADER, NULL);
    } else {
        psLogMsg(__func__, PS_LOG_WARN, "Unable to read any chips in FPA.\n");
    }

    return success;
}

