/** @file ppSubReadoutPhotometry.c
 *
 *  @brief
 *
 *  @ingroup ppSub
 *
 *  @author IfA
 *  @version $Revision: 1.4 $ $Name: not supported by cvs2svn $
 *  @date $Date: 2009-02-18 00:31:20 $
 *  Copyright 2009 Institute for Astronomy, University of Hawaii
 */

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

#include <stdio.h>
#include <pslib.h>
#include <psmodules.h>
#include <psphot.h>

#include "ppSub.h"

bool ppSubReadoutPhotometry (pmConfig *config, psMetadata *stats, const pmFPAview *view)
{
    psAssert(config, "Require configuration");
    psAssert(view, "Require view");

    // Look up recipe values
    psMetadata *recipe = psMetadataLookupMetadata(NULL, config->recipes, PPSUB_RECIPE); // Recipe for ppSim
    psAssert(recipe, "We checked this earlier, so it should be here.");

    if (!psMetadataLookupBool(NULL, recipe, "PHOTOMETRY")) {
        pmFPAfileDropInternal(config->files, "PSPHOT.BACKMDL");
        pmFPAfileDropInternal(config->files, "PSPHOT.BACKMDL.STDEV");
        return true;
    }

    // The PSF (measured in ppSubMakePSF) is stored on the chip->analysis of PSPHOT.INPUT
    // In order to use an incoming PSF, it must be stored on the chip->analysis of PSPHOT.PSF.LOAD
    pmChip *psfInputChip = pmFPAfileThisChip(config->files, view, "PSPHOT.INPUT"); // Chip with PSF
    psAssert (psfInputChip, "should have been generated for ppSubMakePSF");
    pmChip *psfLoadChip = pmFPAfileThisChip(config->files, view, "PSPHOT.PSF.LOAD"); // Chip to have PSF
    psAssert (psfLoadChip, "PSPHOT.PSF.LOAD should have been defined in ppSubCamera");
    pmPSF *psf = psMetadataLookupPtr(NULL, psfInputChip->analysis, "PSPHOT.PSF"); // PSF for photometry
    if (!psf) {
        psError(PS_ERR_UNEXPECTED_NULL, true, "Unable to find PSF from psphot");
        return false;
    }
    psMetadataAddPtr(psfLoadChip->analysis, PS_LIST_TAIL, "PSPHOT.PSF", PS_DATA_UNKNOWN | PS_META_REPLACE,
                     "PSF from ppSubMakePSF", psf);

    bool mdok = false;

    // psphotReadoutMinimal performs the photometry analysis on PSPHOT.INPUT; we need to move
    // around the pointers so PSPHOT.INPUT corresponds to the output image; previously, it was
    // equivalent to the minuend image.
    pmFPAfile *outputFile = psMetadataLookupPtr(&mdok, config->files, "PPSUB.OUTPUT"); // Output file
    pmReadout *outRO = pmFPAviewThisReadout(view, outputFile->fpa); // Readout with the sources

    // XXX possibly rename this to PPSUB.RESID?
    pmFPAfile *photFile = psMetadataLookupPtr(&mdok, config->files, "PSPHOT.INPUT"); // Photometry file
    pmReadout *photRO = pmFPAviewThisReadout(view, photFile->fpa); // Readout to photometer
    if (!photRO) {
        pmCell *cell = pmFPAviewThisCell(view, photFile->fpa); // Cell to photometer
        photRO = pmReadoutAlloc(cell); // Output readout: subtraction
    }
    photRO->image = psImageCopy(photRO->image, outRO->image, PS_TYPE_F32);
    if (outRO->variance) {
        photRO->variance = psImageCopy(photRO->variance, outRO->variance, PS_TYPE_F32);
    } else {
        psFree(photRO->variance);
        photRO->variance = NULL;
    }
    if (outRO->mask) {
        photRO->mask = psImageCopy(photRO->mask, outRO->mask, PS_TYPE_IMAGE_MASK);
    } else {
        psFree(photRO->mask);
        photRO->mask = NULL;
    }

#if 0
    psMetadataAddPtr(config->files, PS_LIST_TAIL, "PSPHOT.INPUT", PS_DATA_UNKNOWN | PS_META_REPLACE,
                     "psphot input: view on another pmFPAfile", photFile);
#endif

    if (!psphotReadoutMinimal(config, view)) {
        psWarning("Unable to perform photometry on subtracted image.");
        psErrorStackPrint(stderr, "Error stack from photometry:");
        psErrorClear();
    }
#if 1
    photRO->data_exists = true;
    photRO->parent->data_exists = true;
    photRO->parent->parent->data_exists = true;
#endif

    if (stats) {
        psArray *sources = psMetadataLookupPtr(NULL, photRO->analysis, "PSPHOT.SOURCES"); // Sources
        psMetadataAddS32(stats, PS_LIST_TAIL, "NUM_SOURCES", 0, "Number of sources detected",
                         sources ? sources->n : 0);
        psMetadataAddF32(stats, PS_LIST_TAIL, "TIME_PHOT", 0, "Time to do photometry",
                         psTimerClear("PPSUB_PHOT"));
    }

    pmFPAfileDropInternal(config->files, "PSPHOT.BACKMDL");
    pmFPAfileDropInternal(config->files, "PSPHOT.BACKMDL.STDEV");

    return true;
}





#ifdef TESTING
    // Record data about sources: not everything gets into the output CMF files
    {
        pmReadout *photRO = pmFPAviewThisReadout(view, photFile->fpa); // Readout with the sources
        psArray *sources = psMetadataLookupPtr(NULL, photRO->analysis, "PSPHOT.SOURCES"); // Sources
        FILE *sourceFile = fopen("sources.dat", "w"); // File for sources
        fprintf(sourceFile,
                "# x y mag mag_err psf_chisq cr_nsigma ext_nsigma psf_qf flags m_x m_y m_xx m_xy m_yy\n");
        for (int i = 0; i < sources->n; i++) {
            pmSource *source = sources->data[i];
            if (!source) {
                continue;
            }

            float x, y;             // Position of source
            float chi2;             // chi^2 for source
            if (source->modelPSF) {
                x = source->modelPSF->params->data.F32[PM_PAR_XPOS];
                y = source->modelPSF->params->data.F32[PM_PAR_YPOS];
                chi2 = source->modelPSF->chisq;
            } else if (source->peak) {
                x = source->peak->xf;
                y = source->peak->yf;
                chi2 = NAN;
            } else {
                psWarning("No position available for source.");
                continue;
            }

            float xMoment = NAN, yMoment = NAN, xxMoment = NAN, xyMoment = NAN, yyMoment = NAN;
            if (source->moments) {
                xMoment = source->moments->Mx;
                yMoment = source->moments->My;
                xxMoment = source->moments->Mxx;
                xyMoment = source->moments->Mxy;
                yyMoment = source->moments->Myy;
            }

            fprintf(sourceFile, "%f %f %f %f %f %f %f %f %d %f %f %f %f %f\n",
                    x, y, source->psfMag, source->errMag, chi2, source->crNsigma, source->extNsigma,
                    source->pixWeight, source->mode, xMoment, yMoment, xxMoment, xyMoment, yyMoment);
        }
        fclose(sourceFile);
    }
#endif
