# include "psphotInternal.h"

bool psphotAddNoise (pmConfig *config, const pmFPAview *view, const char *filerule) {
    return psphotAddOrSubNoise (config, view, filerule, true);
}

bool psphotSubNoise (pmConfig *config, const pmFPAview *view, const char *filerule) {
    return psphotAddOrSubNoise (config, view, filerule, false);
}

// for now, let's store the detections on the readout->analysis for each readout
bool psphotAddOrSubNoise (pmConfig *config, const pmFPAview *view, const char *filerule, bool add)
{
    bool status = true;

    // select the appropriate recipe information
    psMetadata *recipe  = psMetadataLookupPtr (&status, config->recipes, PSPHOT_RECIPE);
    psAssert (recipe, "missing recipe?");

    int num = psphotFileruleCount(config, filerule);

    // loop over the available readouts
    for (int i = 0; i < num; i++) {
	if (!psphotAddOrSubNoiseReadout (config, view, filerule, i, recipe, add)) {
            psError (PSPHOT_ERR_CONFIG, false, "failed on to modify noise for %s entry %d", filerule, i);
	    return false;
	}
    }
    return true;
}

bool psphotAddOrSubNoiseReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe, bool add) {

    bool status = false;

    // find the currently selected readout
    pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest
    psAssert (file, "missing file?");

    pmReadout *readout = pmFPAviewThisReadout(view, file->fpa);
    psAssert (readout, "missing readout?");
    psAssert (readout->parent, "missing cell?");
    psAssert (readout->parent->concepts, "missing concepts?");

    pmDetections *detections = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.DETECTIONS");
    psAssert (detections, "missing detections?");

    psArray *sources = detections->allSources;
    psAssert (sources, "missing sources?");

    psTimerStart ("psphot.noise");

    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
    psImageMaskType maskVal = psMetadataLookupImageMask(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels
    psAssert (maskVal, "missing mask value?");

    // increase variance by factor*(object noise):
    // weight = flux/gain + rn^2/g^2
    // we are adding factor*flux/gain

    float FACTOR = psMetadataLookupF32 (&status, recipe, "NOISE.FACTOR");
    PS_ASSERT (status, false);
    float SIZE = psMetadataLookupF32 (&status, recipe, "NOISE.SIZE");
    PS_ASSERT (status, false);

    if (SIZE <= 0) {
       return true;
    }

    float GAIN = psMetadataLookupF32(&status, readout->parent->concepts, "CELL.GAIN"); // Cell gain
    PS_ASSERT (status, false);
    if (isfinite(GAIN)) {
        FACTOR /= GAIN;
    }

    // loop over all source
    for (int i = 0; i < sources->n; i++) {
        pmSource *source = sources->data[i];

        // skip sources which were not subtracted
	// NOTE: this bit is not modified when pmSourceOp applies to noise
        if (!(source->tmpFlags & PM_SOURCE_TMPF_SUBTRACTED)) continue;

	pmSourceNoiseOp (source, PM_MODEL_OP_FULL | PM_MODEL_OP_NOISE, FACTOR, SIZE, add, maskVal, 0, 0);
    }
    if (add) {
        psLogMsg ("psphot.noise", PS_LOG_WARN, "add noise for %ld objects: %f sec\n", sources->n, psTimerMark ("psphot.noise"));
    } else {
        psLogMsg ("psphot.noise", PS_LOG_WARN, "sub noise for %ld objects: %f sec\n", sources->n, psTimerMark ("psphot.noise"));
    }

    return true;
}
