# include "psphotInternal.h"

/*
 * Cull a set of peaks contained in a psArray of pmFootprints
 */
psErrorCode
psphotCullPeaks(const pmReadout *readout,   // the image wherein lives the footprint
                const psMetadata *recipe,
                psArray *footprints) {  // array of pmFootprints

    bool status = false;

    float nsigma_delta = psMetadataLookupF32(&status, recipe, "FOOTPRINT_CULL_NSIGMA_DELTA");
    if (!status) {
        nsigma_delta = 0; // min.
    }
    float nsigma_min = psMetadataLookupF32(&status, recipe, "FOOTPRINT_CULL_NSIGMA_MIN");
    if (!status) {
        nsigma_min = 0;
    }
    float fPadding = psMetadataLookupF32(&status, recipe, "FOOTPRINT_CULL_NSIGMA_PAD");
    if (!status) {
        fPadding = 0;
    }
    const float skyStdev = psMetadataLookupF32(&status, readout->analysis, "SKY_STDEV");
    psAssert (status, "cannot find SKY_STDEV in readout->analysis");

    const float min_threshold = nsigma_min*skyStdev;
    
    for (int i = 0; i < footprints->n; i++) {
	// if (i % 50 == 0) fprintf (stderr, "cull %d\n", i);
	pmFootprint *fp = footprints->data[i];
	if (fp->npix > 30000) {
	    fprintf (stderr, "big footprint: %f %f to %f %f (%d pix)\n", fp->bbox.x0, fp->bbox.y0, fp->bbox.x1, fp->bbox.y1, fp->npix);
	}
	psTimerStart ("psphot.cull.footprints");
	if (pmFootprintCullPeaks(readout->image, readout->variance, fp, nsigma_delta, fPadding, min_threshold) != PS_ERR_NONE) {
	    return psError(PS_ERR_UNKNOWN, false, "Culling pmFootprint %d", fp->id);
	}
	float dtime = psTimerMark ("psphot.cull.footprints");
	if (dtime > 1.0) {
	    fprintf (stderr, "slow cull for %d (%f sec)\n", i, dtime);
	}
    }
    return PS_ERR_NONE;
}

// psLogMsg ("psphot", PS_LOG_INFO, "%ld peaks, %ld total footprints: %f sec\n", detections->peaks->n, detections->footprints->n, 
