# include "psphotInternal.h"
# define RADIUS_TYPE int

static float PSF_FIT_NSIGMA;
static float PSF_FIT_PADDING;
static float PSF_FIT_RADIUS = 0;	// radius to use in fitting (ignored if <= 0,
					// and a per-object radius is calculated)

bool psphotInitRadiusPSF(const psMetadata *recipe, const pmModelType type) {

    bool status = true;

    PSF_FIT_NSIGMA  = psMetadataLookupF32(&status, recipe, "PSF_FIT_NSIGMA");
    PSF_FIT_PADDING = psMetadataLookupF32(&status, recipe, "PSF_FIT_PADDING");
    PSF_FIT_RADIUS =  psMetadataLookupF32(&status, recipe, "PSF_FIT_RADIUS");

    return true;
}

// call this function whenever you (re)-define the PSF model
bool psphotCheckRadiusPSF (pmReadout *readout, pmSource *source, pmModel *model, psImageMaskType markVal)
{
    psF32 *PAR = model->params->data.F32;

    // XXX do we have a better value for the sky noise level?  not really...
    pmMoments *moments = source->moments;

    // set the fit radius based on the object flux limit and the model
    float radiusFit = PSF_FIT_RADIUS;
    if (radiusFit <= 0) {		// use fixed radius
	if (moments == NULL) {
	    radiusFit = model->modelRadius(model->params, PSF_FIT_NSIGMA*moments->dSky);
	} else {
	    radiusFit = model->modelRadius(model->params, 1.0);
	}
    }
    model->radiusFit = (RADIUS_TYPE)(radiusFit + PSF_FIT_PADDING);

    if (isnan(model->radiusFit)) psAbort("error in radius");
	
    if (source->mode & PM_SOURCE_MODE_SATSTAR) {
	model->radiusFit *= 2;
    }

    bool status = pmSourceRedefinePixels (source, readout, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit);

    // set the mask to flag the excluded pixels
    psImageKeepCircle (source->maskObj, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit, "OR", markVal);
    return status;
}

bool psphotCheckRadiusPSFBlend (pmReadout *readout, pmSource *source, pmModel *model, psImageMaskType markVal, float dR) {

    psF32 *PAR = model->params->data.F32;

    pmMoments *moments = source->moments;
    if (moments == NULL) return false;

    // set the fit radius based on the object flux limit and the model
    model->radiusFit = (RADIUS_TYPE) (model->modelRadius (model->params, PSF_FIT_NSIGMA*moments->dSky) + dR + PSF_FIT_PADDING);
    if (isnan(model->radiusFit)) psAbort("error in radius");
	
    if (source->mode &  PM_SOURCE_MODE_SATSTAR) {
	model->radiusFit *= 2;
    }

    bool status = pmSourceRedefinePixels (source, readout, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit);

    // set the mask to flag the excluded pixels
    psImageKeepCircle (source->maskObj, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit, "OR", markVal);
    return status;
}

static float EXT_FIT_NSIGMA;
static float EXT_FIT_PADDING;

bool psphotInitRadiusEXT (psMetadata *recipe, pmModelType type) {

    bool status;

    EXT_FIT_NSIGMA   = psMetadataLookupF32 (&status, recipe, "EXT_FIT_NSIGMA");
    EXT_FIT_PADDING  = psMetadataLookupF32 (&status, recipe, "EXT_FIT_PADDING");

    return true;
}

// call this function whenever you (re)-define the EXT model
bool psphotCheckRadiusEXT (pmReadout *readout, pmSource *source, pmModel *model, psImageMaskType markVal) {

    psF32 *PAR = model->params->data.F32;

    pmMoments *moments = source->moments;
    if (moments == NULL) return false;

    // set the fit radius based on the object flux limit and the model
    float rawRadius = model->modelRadius (model->params, EXT_FIT_NSIGMA*moments->dSky);

    model->radiusFit = rawRadius + EXT_FIT_PADDING;
    if (isnan(model->radiusFit)) psAbort("error in radius");

    // redefine the pixels if needed
    bool status = pmSourceRedefinePixels (source, readout, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit);

    // set the mask to flag the excluded pixels
    psImageKeepCircle (source->maskObj, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit, "OR", markVal);
    return status;
}
