# include "ppSim.h"

// Generate a header containing WCS keywords
// this function is called with only one of fpa, chip, cell not NULL
bool ppSimInitHeader(pmConfig *config,
                     pmFPA *fpa,
                     pmChip *chip,
                     pmCell *cell)
{

    float ra0 = psMetadataLookupF32(NULL, config->arguments, "RA"); // Boresight RA (radians)
    float dec0 = psMetadataLookupF32(NULL, config->arguments, "DEC"); // Boresight Dec (radians)
    float pa = psMetadataLookupF32(NULL, config->arguments, "PA"); // Position angle (radians)
    float scale = psMetadataLookupF32(NULL, config->arguments, "SCALE"); // plate scale in arcsec / pixel
    scale *= M_PI / 3600.0 / 180.0; // convert plate scale to radians/pixel

    int binning = psMetadataLookupS32(NULL, config->arguments, "BINNING"); // Binning in x and y

    float x0 = 0.0, y0 = 0.0;
    int xParity = 0, yParity = 0;
    if (cell) {
        int x0Chip = psMetadataLookupS32(NULL, cell->parent->concepts, "CHIP.X0");
        int y0Chip = psMetadataLookupS32(NULL, cell->parent->concepts, "CHIP.Y0");
        int xParityChip = psMetadataLookupS32(NULL, cell->parent->concepts, "CHIP.XPARITY");
        int yParityChip = psMetadataLookupS32(NULL, cell->parent->concepts, "CHIP.YPARITY");

        int x0Cell = psMetadataLookupS32(NULL, cell->concepts, "CELL.X0");
        int y0Cell = psMetadataLookupS32(NULL, cell->concepts, "CELL.Y0");
        int xParityCell = psMetadataLookupS32(NULL, cell->concepts, "CELL.XPARITY");
        int yParityCell = psMetadataLookupS32(NULL, cell->concepts, "CELL.YPARITY");

        x0 = PPSIM_FPA_TO_CELL(0.0, x0Cell, xParityCell, binning, x0Chip, xParityChip);
        y0 = PPSIM_FPA_TO_CELL(0.0, y0Cell, yParityCell, binning, y0Chip, yParityChip);
        xParity = xParityCell * xParityChip;
        yParity = yParityCell * yParityChip;
    }
    if (chip) {
        int x0Chip = psMetadataLookupS32(NULL, chip->concepts, "CHIP.X0");
        int y0Chip = psMetadataLookupS32(NULL, chip->concepts, "CHIP.Y0");
        int xParityChip = psMetadataLookupS32(NULL, chip->concepts, "CHIP.XPARITY");
        int yParityChip = psMetadataLookupS32(NULL, chip->concepts, "CHIP.YPARITY");

        x0 = PPSIM_FPA_TO_CELL(0.0, 0, 1, binning, x0Chip, xParityChip);
        y0 = PPSIM_FPA_TO_CELL(0.0, 0, 1, binning, y0Chip, yParityChip);
        xParity = xParityChip;
        yParity = yParityChip;
    }
    if (fpa) {
        psRegion *bounds = ppSimFPABounds (fpa);
        x0 = 0.5 * (bounds->x1 - bounds->x0);
        y0 = 0.5 * (bounds->y1 - bounds->y0);
        xParity = 1;
        yParity = 1;
        psFree (bounds);
    }
    assert(xParity != 0 && yParity != 0);

    psMetadata *header = psMetadataAlloc(); // Header, to return
    pmAstromWCS *wcs = pmAstromWCSAlloc(1, 1); // WCS structure
    wcs->toSky = psProjectionAlloc(ra0, dec0, scale * xParity, scale * yParity, PS_PROJ_TAN);
    wcs->cdelt1 = scale * PM_DEG_RAD * xParity;
    wcs->cdelt2 = scale * PM_DEG_RAD * yParity;
    wcs->crpix1 = x0;
    wcs->crpix2 = y0;
    wcs->trans->x->coeff[1][0] = cos(pa) * wcs->cdelt1;
    wcs->trans->x->coeff[0][1] = -sin(pa) * wcs->cdelt1;
    wcs->trans->y->coeff[1][0] = sin(pa) * wcs->cdelt2;
    wcs->trans->y->coeff[0][1] = cos(pa) * wcs->cdelt2;

    // These aren't used by pmAstromWCStoHeader, but set them anyway
    wcs->trans->x->coeff[0][0] = ra0;
    wcs->trans->y->coeff[0][0] = dec0;
    wcs->trans->x->coeff[1][1] = 0.0;
    wcs->trans->y->coeff[1][1] = 0.0;

    pmAstromWCStoHeader(header, wcs);
    psFree(wcs);

    if (cell) {
        cell->hdu->header = header;
        cell->data_exists = true;
    }
    if (chip) {
        chip->hdu->header = header;
        chip->data_exists = true;
    }
    if (fpa) {
        fpa->hdu->header = header;
    }

    return true;
}

bool ppSimUpdateConceptsFPA (pmFPA *fpa, pmConfig *config) {

    float expTime = psMetadataLookupF32(NULL, config->arguments, "EXPTIME"); // Exposure time

    const char *filter = psMetadataLookupStr(NULL, config->arguments, "FILTER"); // Filter name
    if (!filter) {
        filter = "NONE";
    }

    ppSimType type = psMetadataLookupS32(NULL, config->arguments, "TYPE"); // Type of image to simulate

    // Update FPA concepts
    const char *typeStr;                // Exposure type String
    switch (type) {
      case PPSIM_TYPE_BIAS:   typeStr = "BIAS";   break;
      case PPSIM_TYPE_DARK:   typeStr = "DARK";   break;
      case PPSIM_TYPE_FLAT:   typeStr = "FLAT";   break;
      case PPSIM_TYPE_OBJECT: typeStr = "OBJECT"; break;
      default:
        psAbort("Should never get here.");
    }

    psMetadataAddStr(fpa->concepts, PS_LIST_TAIL, "FPA.OBSTYPE", PS_META_REPLACE,
                     "Observation type", typeStr);
    psMetadataAddStr(fpa->concepts, PS_LIST_TAIL, "FPA.OBJECT", PS_META_REPLACE,
                     "Observation name", typeStr);
    psMetadataAddF32(fpa->concepts, PS_LIST_TAIL, "FPA.EXPTIME", PS_META_REPLACE,
                     "Exposure time (sec)", expTime);
    psMetadataAddStr(fpa->concepts, PS_LIST_TAIL, "FPA.FILTERID", PS_META_REPLACE, "Filter name", filter);
    psMetadataAddStr(fpa->concepts, PS_LIST_TAIL, "FPA.FILTER", PS_META_REPLACE, "Filter name", filter);

    return true;
}

bool ppSimUpdateConceptsCell (pmCell *cell, pmConfig *config) {

    int binning = psMetadataLookupS32(NULL, config->arguments, "BINNING"); // Binning in x and y
    float expTime = psMetadataLookupF32(NULL, config->arguments, "EXPTIME"); // Exposure time

    psMetadataAddF32(cell->concepts, PS_LIST_TAIL, "CELL.EXPOSURE", PS_META_REPLACE,
                     "Exposure time (sec)", expTime);
    psMetadataAddF32(cell->concepts, PS_LIST_TAIL, "CELL.DARKTIME", PS_META_REPLACE,
                     "Dark time (sec)", expTime);
    psMetadataAddS32(cell->concepts, PS_LIST_TAIL, "CELL.XBIN", PS_META_REPLACE,
                     "Binning in x", binning);
    psMetadataAddS32(cell->concepts, PS_LIST_TAIL, "CELL.YBIN", PS_META_REPLACE,
                     "Binning in y", binning);

    return true;
}
