# include "psphotInternal.h"

// XXX we can probably move this into pmGrowthCurve.c
// XXX need to change the way we grab an image, or else use
//     the 'center' option

// XXX add a option to turn off the curve-of-growth (ie, make the apMag = fitMag everywhere);

bool psphotGrowthCurve (pmReadout *readout, pmPSF *psf, bool ignore, psMaskType maskVal) {

    // bool status;
    float xc, yc, dx, dy;
    float fitMag, apMag;
    float radius;

    // create template model
    pmModel *modelRef = pmModelAlloc(psf->type);

    // use the center of the center pixel of the image
    xc = 0.5*readout->image->numCols + readout->image->col0 + 0.5;
    yc = 0.5*readout->image->numRows + readout->image->row0 + 0.5;
    dx = psf->growth->maxRadius + 1;
    dy = psf->growth->maxRadius + 1;

    // assign the x and y coords to the image center
    // create an object with center intensity of 1000
    modelRef->params->data.F32[PM_PAR_SKY] = 0;
    modelRef->params->data.F32[PM_PAR_I0] = 1000;
    modelRef->params->data.F32[PM_PAR_XPOS] = xc;
    modelRef->params->data.F32[PM_PAR_YPOS] = yc;

    // create modelPSF from this model
    pmModel *model = pmModelFromPSF (modelRef, psf);

    // measure the fitMag for this model
    pmSourcePhotometryModel (&fitMag, model);
    psf->growth->fitMag = fitMag;

    // generate working image for this source
    psRegion region = {xc - dx, xc + dx, yc - dy, yc + dy};
    psImage *view = psImageSubset (readout->image, region);
    psImage *image = psImageCopy (NULL, view, PS_TYPE_F32);
    psImage *mask = psImageCopy (NULL, view, PS_TYPE_U8);
    psImageInit (mask, 0);

    // loop over a range of source fluxes
    // no need to interpolate since we have force the object center
    // to 0.5, 0.5 above
    for (int i = 0; i < psf->growth->radius->n; i++) {

        psImageInit (image, 0.0);

        radius = psf->growth->radius->data.F32[i];

        // NOTE: we use pmModelAdd not pmSourceAdd because we are not working with a normal source
        psImageKeepCircle (mask, xc, yc, radius, "OR", PM_MASK_MARK);
        pmModelAdd (image, mask, model, PM_MODEL_OP_FULL, maskVal);
        pmSourcePhotometryAper (&apMag, model, image, mask, maskVal);
        psImageKeepCircle (mask, xc, yc, radius, "AND", PS_NOT_U8(PM_MASK_MARK));

        // the 'ignore' mode is for testing
        if (ignore) {
            psf->growth->apMag->data.F32[i] = fitMag;
        } else {
            psf->growth->apMag->data.F32[i] = apMag;
        }
    }

    psf->growth->apRef = psVectorInterpolate (psf->growth->radius, psf->growth->apMag, psf->growth->refRadius);
    psf->growth->apLoss = psf->growth->fitMag - psf->growth->apRef;

    psLogMsg ("psphot.growth", 4, "GrowthCurve : apLoss : %f\n", psf->growth->apLoss);

    psFree (view);
    psFree (image);
    psFree (mask);
    psFree (model);
    psFree (modelRef);

    return true;
}
