Index: /branches/eam_branch_20080511/ppSim/src/Makefile.am =================================================================== --- /branches/eam_branch_20080511/ppSim/src/Makefile.am (revision 17957) +++ /branches/eam_branch_20080511/ppSim/src/Makefile.am (revision 17958) @@ -4,34 +4,35 @@ ppSim_LDFLAGS = $(PSLIB_LIBS) $(PSMODULE_LIBS) $(PSPHOT_LIBS) $(PSASTRO_LIBS) ppSim_SOURCES = \ - ppSim.c \ - ppSimArguments.c \ - ppSimCreate.c \ - ppSimLoadStars.c \ - ppSimMakeStars.c \ - ppSimMakeGalaxies.c \ - ppSimMakeBiassec.c \ - ppSimMakeBias.c \ - ppSimMakeDark.c \ - ppSimMakeSky.c \ - ppSimInsertStars.c \ - ppSimInsertGalaxies.c \ - ppSimAddOverscan.c \ - ppSimAddNoise.c \ - ppSimSaturate.c \ - ppSimBounds.c \ - ppSimStars.c \ - ppSimSetPSF.c \ - ppSimUtils.c \ - ppSimLoop.c \ - ppSimLoadSpots.c \ - ppSimPhotom.c \ - ppSimPhotomReadout.c \ - ppSimPhotomFiles.c \ - ppSimLoadForceSources.c \ - ppSimMergeReadouts.c \ - ppSimDetections.c \ - ppSimMergeSources.c \ - ppSimMosaicChip.c \ - ppSimRandomGaussian.c \ + ppSim.c \ + ppSimArguments.c \ + ppSimCreate.c \ + ppSimLoadStars.c \ + ppSimMakeStars.c \ + ppSimMakeGalaxies.c \ + ppSimMakeBiassec.c \ + ppSimMakeBias.c \ + ppSimMakeDark.c \ + ppSimMakeSky.c \ + ppSimInsertStars.c \ + ppSimInsertGalaxies.c \ + ppSimAddOverscan.c \ + ppSimAddNoise.c \ + ppSimSaturate.c \ + ppSimBounds.c \ + ppSimStars.c \ + ppSimSetPSF.c \ + ppSimUtils.c \ + ppSimLoop.c \ + ppSimLoadSpots.c \ + ppSimPhotom.c \ + ppSimPhotomReadoutFake.c \ + ppSimPhotomReadoutForce.c \ + ppSimPhotomFiles.c \ + ppSimLoadForceSources.c \ + ppSimMergeReadouts.c \ + ppSimDetections.c \ + ppSimMergeSources.c \ + ppSimMosaicChip.c \ + ppSimRandomGaussian.c \ ppSimBadPixels.c Index: /branches/eam_branch_20080511/ppSim/src/ppSim.h =================================================================== --- /branches/eam_branch_20080511/ppSim/src/ppSim.h (revision 17957) +++ /branches/eam_branch_20080511/ppSim/src/ppSim.h (revision 17958) @@ -187,9 +187,15 @@ void ppSimRandomGaussianFree(); -bool ppSimPhotomFiles (pmConfig *config, pmFPAfile *input); -bool ppSimPhotomReadout(pmConfig *config, const pmFPAview *view); +bool ppSimPhotomFiles (pmConfig *config, pmFPAfile *fakeFile, pmFPAfile *forceFile); + +bool ppSimPhotomReadoutFake(pmConfig *config, const pmFPAview *view); +bool ppSimPhotomReadoutForce(pmConfig *config, const pmFPAview *view); + psArray *ppSimLoadForceSources(pmConfig *config, const pmFPAview *view); bool ppSimDetections (psImage *significance, psMetadata *recipe, psArray *sources); psArray *ppSimMergeSources (psArray *in1, psArray *in2); +psArray *ppSimSelectSources (pmConfig *config, const pmFPAview *view, const char *filename); +bool ppSimDefinePixels (psArray *sources, pmReadout *readout, psMetadata *recipe); + #endif Index: /branches/eam_branch_20080511/ppSim/src/ppSimArguments.c =================================================================== --- /branches/eam_branch_20080511/ppSim/src/ppSimArguments.c (revision 17957) +++ /branches/eam_branch_20080511/ppSim/src/ppSimArguments.c (revision 17958) @@ -125,4 +125,10 @@ ppSimArgToRecipeF32(&status, options, "PIXEL.SCALE", arguments, "-scale"); // Plate scale + // XXX we need to be more flexible: get this from the input image header or WCS + float ra0 = psMetadataLookupF32(NULL, arguments, "-ra"); // Right Ascension of boresight + float dec0 = psMetadataLookupF32(NULL, arguments, "-dec"); // Declination of boresight + psMetadataAddF32(options, PS_LIST_TAIL, "RA", 0, "Boresight RA (radians)", ra0 * M_PI / 180.0); + psMetadataAddF32(options, PS_LIST_TAIL, "DEC", 0, "Boresight Declination (radians)", dec0 * M_PI / 180.0); + psFree (arguments); return true; Index: /branches/eam_branch_20080511/ppSim/src/ppSimCreate.c =================================================================== --- /branches/eam_branch_20080511/ppSim/src/ppSimCreate.c (revision 17957) +++ /branches/eam_branch_20080511/ppSim/src/ppSimCreate.c (revision 17958) @@ -3,188 +3,206 @@ pmFPAfile *ppSimCreate(pmConfig *config) { - bool status; - bool simImage = false; - PS_ASSERT_PTR_NON_NULL(config, NULL); - pmFPA *fpa = NULL; - - // the input image defines the camera. if it is not supplied, the user must have - // supplied a camera and other metadata on the command line - pmFPAfile *input = pmFPAfileDefineFromArgs (&status, config, "PPSIM.INPUT", "INPUT"); - if (!input) { - // if we have not specified the camera already, we need to interpolate the recipes associated with this camera, and read other command-line recipes - if (!pmConfigReadRecipes(config, PM_RECIPE_SOURCE_CL)) { - psError(PS_ERR_IO, false, "Error merging recipes from camera config for %s", config->cameraName); - return NULL; - } - } else { - // If an image is supplied, we still generate a fake image and merge them together downstream - // (otherwise, we get the variance wrong). - simImage = false; - if (input->type != PM_FPA_FILE_IMAGE) { - psError(PS_ERR_IO, true, "PPSIM.INPUT is not of type IMAGE"); - return NULL; - } - } - - // generate the fpa structure used by the output camera (determined from INPUT or specified) - assert (config->camera); - fpa = pmFPAConstruct(config->camera, config->cameraName); // FPA to contain the observation - if (!fpa) { - psError(PS_ERR_UNEXPECTED_NULL, false, "Unable to construct an FPA from camera configuration."); - return NULL; - } - - // define the output image file -- this is the basis for the ppSimLoop - pmFPAfile *output = pmFPAfileDefineOutput(config, fpa, "PPSIM.OUTPUT"); - if (!output) { - psError(PS_ERR_UNEXPECTED_NULL, false, "Unable to create output file from PPSIM.OUTPUT. Did you forget to specify the format?"); - return NULL; - } - if (output->type != PM_FPA_FILE_IMAGE) { - psError(PS_ERR_BAD_PARAMETER_TYPE, true, "PPSIM.OUTPUT type is not IMAGE"); - psFree(fpa); - return NULL; - } - // XXX we should not require the output image to be written - output->save = true; - - config->format = psMemIncrRefCounter (output->format); - config->formatName = psStringCopy (output->formatName); - - // the recipe is now fully realized for the desired camera - psMetadata *recipe = psMetadataLookupMetadata(&status, config->recipes, PPSIM_RECIPE); // Recipe - - char *typeStr = psMetadataLookupStr(NULL, recipe, "IMAGE.TYPE"); // Type of image to simulate - ppSimType type = ppSimTypeFromString (typeStr); // Type of image to simulate - - if (type == PPSIM_TYPE_OBJECT) { - // adjust the seeing by the scale - float seeing = psMetadataLookupF32(&status, recipe, "SEEING"); - if (isnan(seeing)) { - psError(PS_ERR_BAD_PARAMETER_TYPE, true, "seeing is not defined"); - psFree(fpa); - return NULL; - } - float scale = psMetadataLookupF32(&status, recipe, "PIXEL.SCALE"); - if (isnan(scale)) { - psError(PS_ERR_BAD_PARAMETER_TYPE, true, "pixel scale is not defined"); - psFree(fpa); - return NULL; - } - psMetadataAddF32(recipe, PS_LIST_TAIL, "SEEING", PS_META_REPLACE, "Seeing SIGMA (pixels)", seeing / 2.0 / sqrt(2.0 * log(2.0)) / scale); - } - - if ((type == PPSIM_TYPE_OBJECT) || (type == PPSIM_TYPE_FLAT)) { - // determine the zeropoint from the filter - float zp = psMetadataLookupF32(&status, recipe, "ZEROPOINT"); - if (isnan(zp)) { - char *filter = psMetadataLookupStr(&status, recipe, "FILTER"); - float zp = ppSimGetZeroPoint (recipe, filter); - psMetadataAddF32(recipe, PS_LIST_TAIL, "ZEROPOINT", PS_META_REPLACE, "Photometric zeropoint", zp); - } - } - - // For photometry, we operate on the chip-mosaicked image. we create a copy of the mosaicked - // image for psphot so we can write out a clean image - bool doPhotom = psMetadataLookupBool(&status, recipe, "PHOTOM"); // Density of fakes - if (doPhotom) { - - // we need a chip image if we perform photometry (is it necessary to build it if we don't use it?) - pmFPAfile *chipImage = pmFPAfileDefineChipMosaic(config, output->fpa, "PPSIM.CHIP"); - if (!chipImage) { - psError(PS_ERR_IO, false, _("Unable to generate new file from PPSIM.CHIP")); - psFree(fpa); - return NULL; - } - if (chipImage->type != PM_FPA_FILE_IMAGE) { - psError(PS_ERR_IO, true, "PPSIM.CHIP is not of type IMAGE"); - psFree(fpa); - return NULL; - } - - // define associated psphot input/output files - if (!ppSimPhotomFiles (config, chipImage)) { - psError(PSPHOT_ERR_CONFIG, false, "Trouble defining the additional input/output files for psphot"); - psFree(fpa); - return NULL; - } - } else { - // have we supplied a psf model? this happens in ppSimPhotomFiles if we request a photometry - // analysis. however, even if we do not, a psf model may be used to generate the fake - // sources. - if (psMetadataLookupPtr(NULL, config->arguments, "PSPHOT.PSF")) { - // tie the psf file to the chipMosaic - pmFPAfileBindFromArgs(&status, output, config, "PSPHOT.PSF.LOAD", "PSPHOT.PSF"); - if (!status) { - psError(PS_ERR_UNKNOWN, false, "Failed to find/build PSPHOT.PSF.LOAD"); + bool status; + bool simImage = false; + PS_ASSERT_PTR_NON_NULL(config, NULL); + pmFPA *fpa = NULL; + + // the input image defines the camera. if it is not supplied, the user must have + // supplied a camera and other metadata on the command line + pmFPAfile *input = pmFPAfileDefineFromArgs (&status, config, "PPSIM.INPUT", "INPUT"); + if (!input) { + // if we have not specified the camera already, we need to interpolate the recipes associated with this camera, and read other command-line recipes + if (!pmConfigReadRecipes(config, PM_RECIPE_SOURCE_CL)) { + psError(PS_ERR_IO, false, "Error merging recipes from camera config for %s", config->cameraName); + return NULL; + } + } else { + // If an image is supplied, we still generate a fake image and merge them together downstream + // (otherwise, we get the variance wrong). + simImage = false; + if (input->type != PM_FPA_FILE_IMAGE) { + psError(PS_ERR_IO, true, "PPSIM.INPUT is not of type IMAGE"); + return NULL; + } + } + + // generate the fpa structure used by the output camera (determined from INPUT or specified) + assert (config->camera); + fpa = pmFPAConstruct(config->camera, config->cameraName); // FPA to contain the observation + if (!fpa) { + psError(PS_ERR_UNEXPECTED_NULL, false, "Unable to construct an FPA from camera configuration."); + return NULL; + } + + // define the output image file -- this is the basis for the ppSimLoop + pmFPAfile *output = pmFPAfileDefineOutput(config, fpa, "PPSIM.OUTPUT"); + if (!output) { + psError(PS_ERR_UNEXPECTED_NULL, false, "Unable to create output file from PPSIM.OUTPUT. Did you forget to specify the format?"); + return NULL; + } + if (output->type != PM_FPA_FILE_IMAGE) { + psError(PS_ERR_BAD_PARAMETER_TYPE, true, "PPSIM.OUTPUT type is not IMAGE"); psFree(fpa); return NULL; - } - } - } - - // PPSIM.SOURCES carries the constructed, fake sources with their true parameters - // XXX only invoke this code for OBJECT types of images? - pmFPAfile *simSources = pmFPAfileDefineOutput (config, output->fpa, "PPSIM.SOURCES"); - if (!simSources) { - psError(PS_ERR_UNKNOWN, false, "Cannot find a rule for PPSIM.SOURCES"); - return false; - } - simSources->save = true; - - // if we have loaded an input image, we derive certain values from the image, if possible - if (input) { - // we need to extract certain metadata from the image and populate the recipe. - // or else we need to set the fpa concepts based on the recipe options... - + } + // XXX we should not require the output image to be written + output->save = true; + + config->format = psMemIncrRefCounter (output->format); + config->formatName = psStringCopy (output->formatName); + + // the recipe is now fully realized for the desired camera psMetadata *recipe = psMetadataLookupMetadata(&status, config->recipes, PPSIM_RECIPE); // Recipe - ppSimArgToRecipeF32(&status, recipe, "EXPTIME", input->fpa->concepts, "FPA.EXPOSURE"); - char *filter = ppSimArgToRecipeStr(&status, recipe, "FILTER", input->fpa->concepts, "FPA.FILTER"); - - float zp = ppSimGetZeroPoint (recipe, filter); - psMetadataAddF32(recipe, PS_LIST_TAIL, "ZEROPOINT", PS_META_REPLACE, "Photometric zeropoint", zp); - } - - pmFPALevel phuLevel = pmFPAPHULevel(output->format); // Level at which PHU goes - - pmFPAview *view = pmFPAviewAlloc(0);// View for current level - - if (phuLevel == PM_FPA_LEVEL_FPA) { - if (!pmFPAAddSourceFromView(fpa, "Simulation", view, output->format)) { - psError(PS_ERR_UNKNOWN, false, "Unable to add PHU to FPA."); - psFree(fpa); - psFree(view); - return NULL; - } - } - - pmChip *chip; // Chip from FPA - while ((chip = pmFPAviewNextChip(view, fpa, 1))) { - if (phuLevel == PM_FPA_LEVEL_CHIP) { - if (!pmFPAAddSourceFromView(fpa, "Simulation", view, output->format)) { - psError(PS_ERR_UNKNOWN, false, "Unable to add PHU to FPA."); - psFree(fpa); - psFree(view); - return NULL; - } - } - - pmCell *cell; // Cell from chip - while ((cell = pmFPAviewNextCell(view, fpa, 1))) { - if (phuLevel == PM_FPA_LEVEL_CELL) { + char *typeStr = psMetadataLookupStr(NULL, recipe, "IMAGE.TYPE"); // Type of image to simulate + ppSimType type = ppSimTypeFromString (typeStr); // Type of image to simulate + + if (type == PPSIM_TYPE_OBJECT) { + // adjust the seeing by the scale + float seeing = psMetadataLookupF32(&status, recipe, "SEEING"); + if (isnan(seeing)) { + psError(PS_ERR_BAD_PARAMETER_TYPE, true, "seeing is not defined"); + psFree(fpa); + return NULL; + } + float scale = psMetadataLookupF32(&status, recipe, "PIXEL.SCALE"); + if (isnan(scale)) { + psError(PS_ERR_BAD_PARAMETER_TYPE, true, "pixel scale is not defined"); + psFree(fpa); + return NULL; + } + psMetadataAddF32(recipe, PS_LIST_TAIL, "SEEING", PS_META_REPLACE, "Seeing SIGMA (pixels)", seeing / 2.0 / sqrt(2.0 * log(2.0)) / scale); + } + + if ((type == PPSIM_TYPE_OBJECT) || (type == PPSIM_TYPE_FLAT)) { + // determine the zeropoint from the filter + float zp = psMetadataLookupF32(&status, recipe, "ZEROPOINT"); + if (isnan(zp)) { + char *filter = psMetadataLookupStr(&status, recipe, "FILTER"); + float zp = ppSimGetZeroPoint (recipe, filter); + psMetadataAddF32(recipe, PS_LIST_TAIL, "ZEROPOINT", PS_META_REPLACE, "Photometric zeropoint", zp); + } + } + + // For photometry, we operate on the chip-mosaicked image. we create a copy of the mosaicked + // image for psphot so we can write out a clean image + bool doPhotom = psMetadataLookupBool(&status, recipe, "PHOTOM"); // Density of fakes + if (doPhotom) { + + // XXX optionally select fake and/or force + + // we need a chip image if we perform photometry (is it necessary to build it if we don't use it?) + pmFPAfile *fakeImage = pmFPAfileDefineChipMosaic(config, output->fpa, "PPSIM.FAKE.CHIP"); + if (!fakeImage) { + psError(PS_ERR_IO, false, _("Unable to generate new file from PPSIM.FAKE.CHIP")); + psFree(fpa); + return NULL; + } + if (fakeImage->type != PM_FPA_FILE_IMAGE) { + psError(PS_ERR_IO, true, "PPSIM.FAKE.CHIP is not of type IMAGE"); + psFree(fpa); + return NULL; + } + + // we need a chip image if we perform photometry (is it necessary to build it if we don't use it?) + pmFPAfile *forceImage = NULL; + if (input) { + forceImage = pmFPAfileDefineChipMosaic(config, input->fpa, "PPSIM.FORCE.CHIP"); + if (!forceImage) { + psError(PS_ERR_IO, false, _("Unable to generate new file from PPSIM.FORCE.CHIP")); + psFree(fpa); + return NULL; + } + if (forceImage->type != PM_FPA_FILE_IMAGE) { + psError(PS_ERR_IO, true, "PPSIM.FORCE.CHIP is not of type IMAGE"); + psFree(fpa); + return NULL; + } + } + + // define associated psphot input/output files + if (!ppSimPhotomFiles (config, fakeImage, forceImage)) { + psError(PSPHOT_ERR_CONFIG, false, "Trouble defining the additional input/output files for psphot"); + psFree(fpa); + return NULL; + } + } else { + // have we supplied a psf model? this happens in ppSimPhotomFiles if we request a photometry + // analysis. however, even if we do not, a psf model may be used to generate the fake + // sources. + if (psMetadataLookupPtr(NULL, config->arguments, "PSPHOT.PSF")) { + // tie the psf file to the chipMosaic + pmFPAfileBindFromArgs(&status, output, config, "PSPHOT.PSF.LOAD", "PSPHOT.PSF"); + if (!status) { + psError(PS_ERR_UNKNOWN, false, "Failed to find/build PSPHOT.PSF.LOAD"); + psFree(fpa); + return NULL; + } + } + } + + // PPSIM.SOURCES carries the constructed, fake sources with their true parameters + // XXX only invoke this code for OBJECT types of images? + pmFPAfile *simSources = pmFPAfileDefineOutput (config, output->fpa, "PPSIM.SOURCES"); + if (!simSources) { + psError(PS_ERR_UNKNOWN, false, "Cannot find a rule for PPSIM.SOURCES"); + return false; + } + simSources->save = true; + + // if we have loaded an input image, we derive certain values from the image, if possible + if (input) { + // we need to extract certain metadata from the image and populate the recipe. + // or else we need to set the fpa concepts based on the recipe options... + + psMetadata *recipe = psMetadataLookupMetadata(&status, config->recipes, PPSIM_RECIPE); // Recipe + + ppSimArgToRecipeF32(&status, recipe, "EXPTIME", input->fpa->concepts, "FPA.EXPOSURE"); + char *filter = ppSimArgToRecipeStr(&status, recipe, "FILTER", input->fpa->concepts, "FPA.FILTER"); + + float zp = ppSimGetZeroPoint (recipe, filter); + psMetadataAddF32(recipe, PS_LIST_TAIL, "ZEROPOINT", PS_META_REPLACE, "Photometric zeropoint", zp); + } + + pmFPALevel phuLevel = pmFPAPHULevel(output->format); // Level at which PHU goes + + pmFPAview *view = pmFPAviewAlloc(0);// View for current level + + if (phuLevel == PM_FPA_LEVEL_FPA) { if (!pmFPAAddSourceFromView(fpa, "Simulation", view, output->format)) { - psError(PS_ERR_UNKNOWN, false, "Unable to add PHU to FPA."); - psFree(fpa); - psFree(view); - return NULL; - } - } - } - } - - psFree(fpa); - psFree(view); - - return output; + psError(PS_ERR_UNKNOWN, false, "Unable to add PHU to FPA."); + psFree(fpa); + psFree(view); + return NULL; + } + } + + pmChip *chip; // Chip from FPA + while ((chip = pmFPAviewNextChip(view, fpa, 1))) { + if (phuLevel == PM_FPA_LEVEL_CHIP) { + if (!pmFPAAddSourceFromView(fpa, "Simulation", view, output->format)) { + psError(PS_ERR_UNKNOWN, false, "Unable to add PHU to FPA."); + psFree(fpa); + psFree(view); + return NULL; + } + } + + pmCell *cell; // Cell from chip + while ((cell = pmFPAviewNextCell(view, fpa, 1))) { + if (phuLevel == PM_FPA_LEVEL_CELL) { + if (!pmFPAAddSourceFromView(fpa, "Simulation", view, output->format)) { + psError(PS_ERR_UNKNOWN, false, "Unable to add PHU to FPA."); + psFree(fpa); + psFree(view); + return NULL; + } + } + } + } + + psFree(fpa); + psFree(view); + + return output; } Index: /branches/eam_branch_20080511/ppSim/src/ppSimDetectionLimits.c =================================================================== --- /branches/eam_branch_20080511/ppSim/src/ppSimDetectionLimits.c (revision 17957) +++ /branches/eam_branch_20080511/ppSim/src/ppSimDetectionLimits.c (revision 17958) @@ -17,5 +17,5 @@ // find the currently selected readout. // we always perform photometry on the mosaiced chip - pmReadout *readout = pmFPAfileThisReadout (config->files, view, "PPSIM.CHIP"); + pmReadout *readout = pmFPAfileThisReadout (config->files, view, "PPSIM.FAKE.CHIP"); PS_ASSERT_PTR_NON_NULL (readout, false); Index: /branches/eam_branch_20080511/ppSim/src/ppSimLoadForceSources.c =================================================================== --- /branches/eam_branch_20080511/ppSim/src/ppSimLoadForceSources.c (revision 17957) +++ /branches/eam_branch_20080511/ppSim/src/ppSimLoadForceSources.c (revision 17958) @@ -25,13 +25,15 @@ psAssert(srcFile, "missing PPSIM.REAL.SOURCES"); + // select the fpa and chip from PPSIM.REAL.SOURCES, which carries the astrometry in the headers pmFPA *fpa = srcFile->fpa; pmChip *chip = pmFPAviewThisChip (view, fpa); - // pmCell *cell = pmFPAviewThisCell (view, fpa); - pmReadout *readout = pmFPAviewThisReadout (view, fpa); + + // we use the readout for the PPSIM.FORCE.CHIP, which contains the chip-mosaicked image(s) + pmReadout *readout = pmFPAfileThisReadout (config->files, view, "PPSIM.FORCE.CHIP"); // XXX we should have only one cell and readout per chip, right? // read WCS from existing output from psphot & psastro - // XXX can we be more flexible? ie, use the header of PPSIM.CHIP? + // XXX can we be more flexible? ie, use the header of PPSIM.FORCE.CHIP? bool bilevelAstrometry = false; status = psastroAstromGuessSetFPA (fpa, &bilevelAstrometry); Index: /branches/eam_branch_20080511/ppSim/src/ppSimLoadSpots.c =================================================================== --- /branches/eam_branch_20080511/ppSim/src/ppSimLoadSpots.c (revision 17957) +++ /branches/eam_branch_20080511/ppSim/src/ppSimLoadSpots.c (revision 17958) @@ -9,12 +9,14 @@ bool status; - psArray *spots = psArrayAllocEmpty (1); - psMetadataAddPtr (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 psMetadata *recipe = psMetadataLookupMetadata(&status, config->recipes, PPSIM_RECIPE); // Recipe bool forcedPhot = psMetadataLookupBool(&status, recipe, "FORCED.PHOT"); // Density of fakes - if (!forcedPhot) return true; + 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 @@ -30,6 +32,4 @@ char *catdir = psMetadataLookupStr(NULL, recipe, "FORCED.CATDIR"); - char *photcode = psMetadataLookupStr(NULL, recipe, "FORCED.PHOTCODE"); - float MAGmax = 100.0; // Size of FPA @@ -47,35 +47,19 @@ psMetadataAddF32(astroRecipe, PS_LIST_TAIL, "DEC_MAX", PS_META_REPLACE, "", dec0 + radius); - // XXX these two values are somewhat bogus: can we have psastroLoadRefstars optionally ignore photcode and MAGmax? - psMetadataAddStr(astroRecipe, PS_LIST_TAIL, "DVO.GETSTAR.PHOTCODE", PS_META_REPLACE, "", photcode); - psMetadataAddF32(astroRecipe, PS_LIST_TAIL, "DVO.GETSTAR.MAG.MAX", PS_META_REPLACE, "", MAGmax); + // 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); - psArray *refStars = psastroLoadRefstars(config); - if (!refStars || refStars->n == 0) { + // 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(refStars); + psFree(spots); return false; } - psLogMsg("ppSim", PS_LOG_INFO, "Adding %ld reference stars", refStars->n); + psLogMsg("ppSim", PS_LOG_INFO, "Adding %ld reference stars", spots->n); - // convert the pmAstromObj sources to the storage format used by the CMF files - psArrayRealloc (spots, refStars->n); - for (long i = 0; i < refStars->n; i++) { - ppSimStar *spot = ppSimStarAlloc (); - - pmAstromObj *ref = refStars->data[i]; // Reference star - spot->ra = ref->sky->r; // RA of star - spot->dec = ref->sky->d; // Dec of star - spot->mag = ref->Mag; // Magnitude of star - - // Apply rotation, make FPA center of boresite - // convert the ra,dec to x,y using examples in psastro. - // spot->x = cos(pa) * xi - sin(pa) * eta + x0fpa; - // spot->y = sin(pa) * xi + cos(pa) * eta + y0fpa; - - spots->data[i] = spot; - - psTrace("ppSim", 10, "Adding catalogue star: %.1f,%.1f --> %.2f\n", spot->x, spot->y, ref->Mag); - } + 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; Index: /branches/eam_branch_20080511/ppSim/src/ppSimPhotom.c =================================================================== --- /branches/eam_branch_20080511/ppSim/src/ppSimPhotom.c (revision 17957) +++ /branches/eam_branch_20080511/ppSim/src/ppSimPhotom.c (revision 17958) @@ -16,10 +16,11 @@ int blankMask = 0; // XXX not sure what this should be set to... - ppSimMosaicChip(config, blankMask, view, "PPSIM.CHIP", "PPSIM.OUTPUT"); + ppSimMosaicChip(config, blankMask, view, "PPSIM.FORCE.CHIP", "PPSIM.INPUT"); + ppSimMosaicChip(config, blankMask, view, "PPSIM.FAKE.CHIP", "PPSIM.OUTPUT"); - // XXX have ppSimMosaicChip return the pointer? - pmFPAfile *input = psMetadataLookupPtr (&status, config->files, "PPSIM.CHIP"); + // use PPSIM.FAKE.CHIP as the base since it is guaranteed to exist (derived from PPSIM.OUTPUT) + pmFPAfile *input = psMetadataLookupPtr (&status, config->files, "PPSIM.FAKE.CHIP"); if (!status) { - psError(PSPHOT_ERR_CONFIG, false, "PPSIM.CHIP I/O file is not defined"); + psError(PSPHOT_ERR_CONFIG, false, "PPSIM.FAKE.CHIP I/O file is not defined"); return false; } @@ -46,5 +47,9 @@ // run the actual photometry analysis - if (!ppSimPhotomReadout (config, view)) { + if (!ppSimPhotomReadoutFake (config, view)) { + psError(psErrorCodeLast(), false, "failure in psphotReadout for chip %d, cell %d, readout %d\n", view->chip, view->cell, view->readout); + return false; + } + if (!ppSimPhotomReadoutForce (config, view)) { psError(psErrorCodeLast(), false, "failure in psphotReadout for chip %d, cell %d, readout %d\n", view->chip, view->cell, view->readout); return false; Index: /branches/eam_branch_20080511/ppSim/src/ppSimPhotomFiles.c =================================================================== --- /branches/eam_branch_20080511/ppSim/src/ppSimPhotomFiles.c (revision 17957) +++ /branches/eam_branch_20080511/ppSim/src/ppSimPhotomFiles.c (revision 17958) @@ -2,5 +2,5 @@ // define the needed / desired I/O files -bool ppSimPhotomFiles (pmConfig *config, pmFPAfile *input) { +bool ppSimPhotomFiles (pmConfig *config, pmFPAfile *fakeFile, pmFPAfile *forceFile) { bool status = false; @@ -13,5 +13,5 @@ int DX = psMetadataLookupS32 (&status, recipe, "BACKGROUND.XBIN"); int DY = psMetadataLookupS32 (&status, recipe, "BACKGROUND.YBIN"); - pmFPAfile *output = pmFPAfileDefineFromFile (config, input, DX, DY, "PSPHOT.BACKMDL"); + pmFPAfile *output = pmFPAfileDefineFromFile (config, fakeFile, DX, DY, "PSPHOT.BACKMDL"); if (!output) { psError(PSPHOT_ERR_CONFIG, false, "Cannot find a rule for PSPHOT.BACKMDL"); @@ -24,5 +24,5 @@ int DX = psMetadataLookupS32 (&status, recipe, "BACKGROUND.XBIN"); int DY = psMetadataLookupS32 (&status, recipe, "BACKGROUND.YBIN"); - pmFPAfile *output = pmFPAfileDefineFromFile (config, input, DX, DY, "PSPHOT.BACKMDL.STDEV"); + pmFPAfile *output = pmFPAfileDefineFromFile (config, fakeFile, DX, DY, "PSPHOT.BACKMDL.STDEV"); if (!output) { psError(PSPHOT_ERR_CONFIG, false, "Cannot find a rule for PSPHOT.BACKMDL.STDEV"); @@ -33,5 +33,5 @@ // optionally save the full background image if (psMetadataLookupBool(NULL, recipe, "SAVE.BACKGND")) { - pmFPAfile *output = pmFPAfileDefineFromFile (config, input, 1, 1, "PSPHOT.BACKGND"); + pmFPAfile *output = pmFPAfileDefineFromFile (config, fakeFile, 1, 1, "PSPHOT.BACKGND"); if (!output) { psError(PSPHOT_ERR_CONFIG, false, "Cannot find a rule for PSPHOT.BACKGND"); @@ -42,5 +42,5 @@ // optionally save the background-subtracted image if (psMetadataLookupBool(NULL, recipe, "SAVE.BACKSUB")) { - pmFPAfile *output = pmFPAfileDefineFromFile (config, input, 1, 1, "PSPHOT.BACKSUB"); + pmFPAfile *output = pmFPAfileDefineFromFile (config, fakeFile, 1, 1, "PSPHOT.BACKSUB"); if (!output) { psError(PSPHOT_ERR_CONFIG, false, "Cannot find a rule for PSPHOT.BACKSUB"); @@ -57,6 +57,6 @@ } - // the fake sources are carried on the input->fpa structures - pmFPAfile *outsources = pmFPAfileDefineOutputFromFile (config, input, "PPSIM.FAKE.SOURCES"); + // the fake sources are carried on the fakeFile->fpa structures (PPSIM.FAKE.CHIP) + pmFPAfile *outsources = pmFPAfileDefineOutput (config, fakeFile->fpa, "PPSIM.FAKE.SOURCES"); if (!outsources) { psError(PSPHOT_ERR_CONFIG, false, "Cannot find a rule for PPSIM.FAKE.SOURCES"); @@ -65,6 +65,14 @@ outsources->save = true; + // the forced sources are carried on the forceFile->fpa structures (PPSIM.FORCE.CHIP) + pmFPAfile *outforced = pmFPAfileDefineOutput (config, forceFile->fpa, "PPSIM.FORCE.SOURCES"); + if (!outforced) { + psError(PSPHOT_ERR_CONFIG, false, "Cannot find a rule for PPSIM.FORCE.SOURCES"); + return false; + } + outforced->save = true; + // we bind the psffile to the chip-mosaicked image so we have a reference to the image - pmFPAfile *psffile = pmFPAfileBindFromArgs(&status, input, config, "PSPHOT.PSF.LOAD", "PSPHOT.PSF"); + pmFPAfile *psffile = pmFPAfileBindFromArgs(&status, fakeFile, config, "PSPHOT.PSF.LOAD", "PSPHOT.PSF"); if (!psffile) { psError(PSPHOT_ERR_CONFIG, false, "Failed to find/build PSPHOT.PSF.LOAD"); Index: /branches/eam_branch_20080511/ppSim/src/ppSimPhotomReadout.c =================================================================== --- /branches/eam_branch_20080511/ppSim/src/ppSimPhotomReadout.c (revision 17957) +++ /branches/eam_branch_20080511/ppSim/src/ppSimPhotomReadout.c (revision 17958) @@ -158,9 +158,27 @@ // these are outputs: we need to generate the fpa structure + pmCell *forceCell = pmFPAfileThisCell (config->files, view, "PPSIM.FORCE.SOURCES"); psAssert (forceCell, "no cell?"); + pmChip *forceChip = pmFPAfileThisChip (config->files, view, "PPSIM.FORCE.SOURCES"); psAssert (forceChip, "no chip?"); + pmReadout *forceReadout = pmFPAfileThisReadout (config->files, view, "PPSIM.FORCE.SOURCES"); + if (!forceReadout) { + forceReadout = pmReadoutAlloc (forceCell); + psFree (forceReadout); // there is a copy on 'cell' as well + } + psAssert (forceReadout, "no forceReadout?"); + pmChipSetDataStatus (forceChip, true); + psMetadataAddArray (forceReadout->analysis, PS_LIST_TAIL, "PSPHOT.SOURCES", PS_META_REPLACE, "forced photometry ", forceSources); + psFree (forceSources); + + pmCell *fakeCell = pmFPAfileThisCell (config->files, view, "PPSIM.FAKE.SOURCES"); psAssert (fakeCell, "no cell?"); + pmChip *fakeChip = pmFPAfileThisChip (config->files, view, "PPSIM.FAKE.SOURCES"); psAssert (fakeChip, "no chip?"); pmReadout *fakeReadout = pmFPAfileThisReadout (config->files, view, "PPSIM.FAKE.SOURCES"); - psMetadataAdd (fakeReadout->analysis, PS_LIST_TAIL, "PSPHOT.SOURCES", PS_DATA_ARRAY, "fake photometry ", fakeSources); + if (!fakeReadout) { + fakeReadout = pmReadoutAlloc (fakeCell); + psFree (fakeReadout); // there is a copy on 'cell' as well + } + psAssert (fakeReadout, "no fakeReadout?"); + pmChipSetDataStatus (fakeChip, true); + psMetadataAddArray (fakeReadout->analysis, PS_LIST_TAIL, "PSPHOT.SOURCES", PS_META_REPLACE, "fake photometry ", fakeSources); - psMetadataAdd (readout->analysis, PS_LIST_TAIL, "PPSIM.FORCE.SOURCES", PS_DATA_ARRAY, "forced photometry ", forceSources); - psFree (forceSources); return true; } Index: /branches/eam_branch_20080511/ppSim/src/ppSimPhotomReadoutFake.c =================================================================== --- /branches/eam_branch_20080511/ppSim/src/ppSimPhotomReadoutFake.c (revision 17958) +++ /branches/eam_branch_20080511/ppSim/src/ppSimPhotomReadoutFake.c (revision 17958) @@ -0,0 +1,124 @@ +# include "ppSim.h" + +bool ppSimPhotomReadoutFake(pmConfig *config, const pmFPAview *view) { + + psTimerStart ("psphotReadout"); + + // select the current recipe + psMetadata *recipe = psMetadataLookupPtr (NULL, config->recipes, PSPHOT_RECIPE); + if (!recipe) { + psError(PSPHOT_ERR_CONFIG, false, "missing recipe %s", PSPHOT_RECIPE); + return false; + } + +# if 0 + // set the photcode for this image + if (!psphotAddPhotcode (recipe, config, view, "PPSIM.FAKE.CHIP")) { + psError (PSPHOT_ERR_CONFIG, false, "trouble defining the photcode"); + return false; + } +# endif + + // *** in this section, perform the photometry for real + fake sources on PPSIM.FAKE.CHIP *** + + // find the currently selected readout. + // we always perform photometry on the mosaiced chip + pmReadout *readout = pmFPAfileThisReadout (config->files, view, "PPSIM.FAKE.CHIP"); + PS_ASSERT_PTR_NON_NULL (readout, false); + + // the previously measured real sources are loaded into the PPSIM.REAL.SOURCES pmFPAfile + psArray *realMeasuredSources = ppSimSelectSources (config, view, "PPSIM.REAL.SOURCES"); + PS_ASSERT_PTR_NON_NULL (realMeasuredSources, false); + + // the newly injected fake sources are saved on the PPSIM.OUTPUT pmFPAfile + psArray *injectedSources = ppSimSelectSources (config, view, "PPSIM.OUTPUT"); + PS_ASSERT_PTR_NON_NULL (injectedSources, false); + + // Generate the mask and weight images, including the user-defined analysis region of interest + psphotSetMaskAndWeight (config, readout, recipe); + + // XXX need to define the source pixels + ppSimDefinePixels (realMeasuredSources, readout, recipe); + ppSimDefinePixels (injectedSources, readout, recipe); + + // make a copy of the injected fake-source parameters so these can be saved independently + // from the measured fake-source parameters + psArray *fakeSources = psArrayAlloc (injectedSources->n); + for (int i = 0; i < injectedSources->n; i++) { + fakeSources->data[i] = pmSourceCopy (injectedSources->data[i]); + } + + // make a copy of the measured real-source parameters so these can be fitted here and in + // ppSimPhotomReadoutForce, if desired. + psArray *realSources = psArrayAlloc (realMeasuredSources->n); + for (int i = 0; i < realMeasuredSources->n; i++) { + realSources->data[i] = pmSourceCopy (realMeasuredSources->data[i]); + } + + // load the psf model, if suppled. FWHM_X,FWHM_Y,etc are saved in the recipe + // this function uses PSPHOT.PSF.LOAD as the pmFPAfile + pmPSF *psf = psphotLoadPSF (config, view, recipe); + assert (psf); + + // remove all sources + psphotRemoveAllSources (realSources, recipe); + + // generate a background model (median, smoothed image) + if (!psphotModelBackground (config, view, "PPSIM.FAKE.CHIP")) { + return psphotReadoutCleanup (config, readout, recipe, NULL, NULL, NULL); + } + if (!psphotSubtractBackground (config, view, "PPSIM.FAKE.CHIP")) { + return psphotReadoutCleanup (config, readout, recipe, NULL, NULL, NULL); + } + + // add noise for real (subtracted) objects + psphotAddNoise (readout, realSources, recipe); + + // XXX fake sources should measure peak->x,y, force sources should not + psMaskType maskVal = 0xff; + psImage *significance = psphotSignificanceImage (readout, recipe, 1, maskVal); + ppSimDetections (significance, recipe, fakeSources); + psFree (significance); + + // remove noise for subtracted objects (ie, return to normal noise level) + psphotSubNoise (readout, realSources, recipe); + + // replace all sources + psphotReplaceAllSources (realSources, recipe); + + // construct an initial model for each object + psphotGuessModels (readout, realSources, recipe, psf); + psphotGuessModels (readout, fakeSources, recipe, psf); + + // linear fit to real + fake sources + psArray *sources = ppSimMergeSources (realSources, fakeSources); + psphotFitSourcesLinear (readout, sources, recipe, psf, FALSE); // XXX option to NOT subtract + psphotReplaceAllSources (sources, recipe); + psFree (sources); // only frees the merged references + + // calculate source magnitudes (for which set??) + pmReadout *background = psphotSelectBackground (config, view, false); + psphotMagnitudes(fakeSources, recipe, psf, background); + + // drop the references to the image pixels held by each source + psphotSourceFreePixels (realSources); + psphotSourceFreePixels (fakeSources); + + // create the exported-metadata and free local data + // XXX this places the sources on readout->analysis as PSPHOT.SOURCES. modify? + // (or don't supply the sources, and do this with a different function) + psphotReadoutCleanup(config, readout, recipe, NULL, psf, NULL); + + pmCell *fakeCell = pmFPAfileThisCell (config->files, view, "PPSIM.FAKE.SOURCES"); psAssert (fakeCell, "no cell?"); + pmChip *fakeChip = pmFPAfileThisChip (config->files, view, "PPSIM.FAKE.SOURCES"); psAssert (fakeChip, "no chip?"); + pmReadout *fakeReadout = pmFPAfileThisReadout (config->files, view, "PPSIM.FAKE.SOURCES"); + if (!fakeReadout) { + fakeReadout = pmReadoutAlloc (fakeCell); + psFree (fakeReadout); // there is a copy on 'cell' as well + } + psAssert (fakeReadout, "no fakeReadout?"); + pmChipSetDataStatus (fakeChip, true); + psMetadataAddArray (fakeReadout->analysis, PS_LIST_TAIL, "PSPHOT.SOURCES", PS_META_REPLACE, "fake photometry ", fakeSources); + + return true; +} Index: /branches/eam_branch_20080511/ppSim/src/ppSimPhotomReadoutForce.c =================================================================== --- /branches/eam_branch_20080511/ppSim/src/ppSimPhotomReadoutForce.c (revision 17958) +++ /branches/eam_branch_20080511/ppSim/src/ppSimPhotomReadoutForce.c (revision 17958) @@ -0,0 +1,117 @@ +# include "ppSim.h" + +bool ppSimPhotomReadoutForce(pmConfig *config, const pmFPAview *view) { + + psTimerStart ("psphotReadout"); + + // select the current recipe + psMetadata *recipe = psMetadataLookupPtr (NULL, config->recipes, PSPHOT_RECIPE); + if (!recipe) { + psError(PSPHOT_ERR_CONFIG, false, "missing recipe %s", PSPHOT_RECIPE); + return false; + } + +# if 0 + // set the photcode for this image + if (!psphotAddPhotcode (recipe, config, view, "PPSIM.FAKE.CHIP")) { + psError (PSPHOT_ERR_CONFIG, false, "trouble defining the photcode"); + return false; + } +# endif + + // *** in this section, perform the photometry for real + force sources on PPSIM.FORCE.CHIP *** + + // find the currently selected readout. + // we always perform photometry on the mosaiced chip + pmReadout *readout = pmFPAfileThisReadout (config->files, view, "PPSIM.FORCE.CHIP"); + PS_ASSERT_PTR_NON_NULL (readout, false); + + // the previously measured real sources are loaded into the PPSIM.REAL.SOURCES pmFPAfile + psArray *realMeasuredSources = ppSimSelectSources (config, view, "PPSIM.REAL.SOURCES"); + PS_ASSERT_PTR_NON_NULL (realMeasuredSources, false); + + // load the forced source lists + psArray *forceSources = ppSimLoadForceSources (config, view); + psAssert (forceSources, "failed to load force photometry sources"); + + // Generate the mask and weight images, including the user-defined analysis region of interest + psphotSetMaskAndWeight (config, readout, recipe); + + // XXX need to define the source pixels + ppSimDefinePixels (realMeasuredSources, readout, recipe); + ppSimDefinePixels (forceSources, readout, recipe); + + // make a copy of the measured real-source parameters so these can be fitted here and in + // ppSimPhotomReadoutForce, if desired. + psArray *realSources = psArrayAlloc (realMeasuredSources->n); + for (int i = 0; i < realMeasuredSources->n; i++) { + realSources->data[i] = pmSourceCopy (realMeasuredSources->data[i]); + } + + // load the psf model, if suppled. FWHM_X,FWHM_Y,etc are saved in the recipe + // this function uses PSPHOT.PSF.LOAD as the pmFPAfile + pmPSF *psf = psphotLoadPSF (config, view, recipe); + assert (psf); + + // remove all sources + psphotRemoveAllSources (realSources, recipe); + + // generate a background model (median, smoothed image) + if (!psphotModelBackground (config, view, "PPSIM.FORCE.CHIP")) { + return psphotReadoutCleanup (config, readout, recipe, NULL, NULL, NULL); + } + if (!psphotSubtractBackground (config, view, "PPSIM.FORCE.CHIP")) { + return psphotReadoutCleanup (config, readout, recipe, NULL, NULL, NULL); + } + + // add noise for real (subtracted) objects + psphotAddNoise (readout, realSources, recipe); + + // XXX fake sources should measure peak->x,y, force sources should not + psMaskType maskVal = 0xff; + psImage *significance = psphotSignificanceImage (readout, recipe, 1, maskVal); + ppSimDetections (significance, recipe, forceSources); + psFree (significance); + + // remove noise for subtracted objects (ie, return to normal noise level) + psphotSubNoise (readout, realSources, recipe); + + // replace all sources + psphotReplaceAllSources (realSources, recipe); + + // construct an initial model for each object + psphotGuessModels (readout, realSources, recipe, psf); + psphotGuessModels (readout, forceSources, recipe, psf); + + // linear fit to real + force sources + psArray *sources = ppSimMergeSources (realSources, forceSources); + psphotFitSourcesLinear (readout, sources, recipe, psf, FALSE); // XXX option to NOT subtract + psphotReplaceAllSources (sources, recipe); + psFree (sources); // only frees the merged references + + // calculate source magnitudes (for which set??) + pmReadout *background = psphotSelectBackground (config, view, false); + psphotMagnitudes(forceSources, recipe, psf, background); + + // drop the references to the image pixels held by each source + psphotSourceFreePixels (realSources); + psphotSourceFreePixels (forceSources); + + // create the exported-metadata and free local data + // XXX this places the sources on readout->analysis as PSPHOT.SOURCES. modify? + // (or don't supply the sources, and do this with a different function) + psphotReadoutCleanup(config, readout, recipe, NULL, psf, NULL); + + pmCell *forceCell = pmFPAfileThisCell (config->files, view, "PPSIM.FORCE.SOURCES"); psAssert (forceCell, "no cell?"); + pmChip *forceChip = pmFPAfileThisChip (config->files, view, "PPSIM.FORCE.SOURCES"); psAssert (forceChip, "no chip?"); + pmReadout *forceReadout = pmFPAfileThisReadout (config->files, view, "PPSIM.FORCE.SOURCES"); + if (!forceReadout) { + forceReadout = pmReadoutAlloc (forceCell); + psFree (forceReadout); // there is a copy on 'cell' as well + } + psAssert (forceReadout, "no forceReadout?"); + pmChipSetDataStatus (forceChip, true); + psMetadataAddArray (forceReadout->analysis, PS_LIST_TAIL, "PSPHOT.SOURCES", PS_META_REPLACE, "force photometry ", forceSources); + + return true; +} Index: /branches/eam_branch_20080511/ppSim/src/ppSimUtils.c =================================================================== --- /branches/eam_branch_20080511/ppSim/src/ppSimUtils.c (revision 17957) +++ /branches/eam_branch_20080511/ppSim/src/ppSimUtils.c (revision 17958) @@ -274,2 +274,28 @@ return NAN; } + +psArray *ppSimSelectSources (pmConfig *config, const pmFPAview *view, const char *filename) { + + pmReadout *readout = pmFPAfileThisReadout (config->files, view, filename); + PS_ASSERT_PTR_NON_NULL (readout, NULL); + + psArray *sources = psMetadataLookupPtr (NULL, readout->analysis, "PSPHOT.SOURCES"); + return sources; +} + +bool ppSimDefinePixels (psArray *sources, pmReadout *readout, psMetadata *recipe) { + + bool status; + + float OUTER = psMetadataLookupF32 (&status, recipe, "SKY_OUTER_RADIUS"); + if (!status) return NULL; + + for (int i = 0; i < sources->n; i++) { + pmSource *source = sources->data[i]; + + // allocate image, weight, mask arrays for each peak (square of radius OUTER) + pmSourceDefinePixels (source, readout, source->peak->x, source->peak->y, OUTER); + } + return true; +} +