# include "ppSim.h"

// Load the relevant forced-photometry positions for this field.  This function does not determine
// the pixel coordinates, merely the celestial coordinates for sources in the general vicinity.  The
// sources are saved on the fpa->analysis metadata as FORCED.SPOTS.  We always create an entry; it
// will be empty if no sources were selected or forced photometry is not desired.
bool ppSimLoadSpots (pmFPA *fpa, pmConfig *config) {

    bool status;


    psMetadata *recipe = psMetadataLookupMetadata(&status, config->recipes, PPSIM_RECIPE); // Recipe

    bool forcedPhot = psMetadataLookupBool(&status, recipe, "FORCED.PHOT"); // Density of fakes
    if (!forcedPhot) {
	psArray *spots = psArrayAllocEmpty (1);
	psMetadataAddArray (fpa->analysis, PS_LIST_TAIL, "FORCED.SPOTS", PS_META_REPLACE, "forced photometry positions", spots);
	psFree (spots); // free the extra (local) reference; a copy remains on fpa->analysis
	return true;
    }

    // We read the catalogue stars using psastroLoadRefstars
    psMetadata *astroRecipe = psMetadataLookupPtr(NULL, config->recipes, PSASTRO_RECIPE);
    if (!astroRecipe) {
        psError(PSASTRO_ERR_CONFIG, false, "Can't find recipe %s", PSASTRO_RECIPE);
        return false;
    }

    float ra0     = psMetadataLookupF32(NULL, recipe, "RA");        // Boresight RA (radians)
    float dec0    = psMetadataLookupF32(NULL, recipe, "DEC");       // Boresight Dec (radians)
    float scale   = psMetadataLookupF32(NULL, recipe, "PIXEL.SCALE") * M_PI / 3600.0 / 180.0; // Plate scale (radians/pixel)

    if (isnan(ra0) || isnan(dec0)) {
        psError(PS_ERR_UNKNOWN, false, "image boresite coords not defined.");
	return false;
    }

    char *catdir = psMetadataLookupStr(NULL, recipe, "FORCED.CATDIR");

    // Size of FPA
    psRegion *bounds = ppSimFPABounds (fpa);
    float radius = 0.5 * PS_MAX(bounds->x1 - bounds->x0, bounds->y1 - bounds->y0) * scale;
    psFree(bounds);

    // modify the PSASTRO recipe to use the values desired for FORCED.PHOT.  we can use a view on
    // the PSASTRO recipe because we are not calling psastro elsewhere in this program.  XXX need to
    // use the WCS to define the overlap region
    psMetadataAddStr(astroRecipe, PS_LIST_TAIL, "PSASTRO.CATDIR",  PS_META_REPLACE, "", catdir);
    psMetadataAddF32(astroRecipe, PS_LIST_TAIL, "RA_MIN",  PS_META_REPLACE, "", ra0 - radius);
    psMetadataAddF32(astroRecipe, PS_LIST_TAIL, "RA_MAX",  PS_META_REPLACE, "", ra0 + radius);
    psMetadataAddF32(astroRecipe, PS_LIST_TAIL, "DEC_MIN", PS_META_REPLACE, "", dec0 - radius);
    psMetadataAddF32(astroRecipe, PS_LIST_TAIL, "DEC_MAX", PS_META_REPLACE, "", dec0 + radius);
    
    // tell psastroLoadRefstars to ignore photcode and maglim
    psMetadataAddStr(astroRecipe, PS_LIST_TAIL, "DVO.GETSTAR.PHOTCODE",  PS_META_REPLACE, "", "NONE");
    psMetadataAddF32(astroRecipe, PS_LIST_TAIL, "DVO.GETSTAR.MAG.MAX",  PS_META_REPLACE, "", NAN);

    // load refstars from the catalog
    psArray *spots = psastroLoadRefstars(config);
    if (!spots || spots->n == 0) {
        psError(PS_ERR_UNKNOWN, false, "Unable to find reference stars.");
        psFree(spots);
        return false;
    }
    psLogMsg("ppSim", PS_LOG_INFO, "Adding %ld reference stars", spots->n);

    psMetadataAddArray (fpa->analysis, PS_LIST_TAIL, "FORCED.SPOTS", PS_META_REPLACE, "forced photometry positions", spots);
    psFree (spots); // free the extra (local) reference; a copy remains on fpa->analysis

    return true;
}
