Index: trunk/ppSub/src/Makefile.am =================================================================== --- trunk/ppSub/src/Makefile.am (revision 23688) +++ trunk/ppSub/src/Makefile.am (revision 23740) @@ -13,22 +13,25 @@ ppSub_LDFLAGS = $(PSLIB_LIBS) $(PSMODULE_LIBS) $(PPSTATS_LIBS) $(PSPHOT_LIBS) $(PPSUB_LIBS) -ppSub_SOURCES = \ - ppSub.c \ - ppSubArguments.c \ - ppSubVersion.c \ - ppSubBackground.c \ - ppSubCamera.c \ - ppSubData.c \ - ppSubLoop.c \ - ppSubReadout.c \ - ppSubDefineOutput.c \ - ppSubExtras.c \ - ppSubMakePSF.c \ - ppSubMatchPSFs.c \ - ppSubReadoutPhotometry.c \ - ppSubReadoutSubtract.c \ - ppSubReadoutUpdate.c \ - ppSubSetMasks.c \ - ppSubReadoutRenorm.c \ +ppSub_SOURCES = \ + ppSub.c \ + ppSubArguments.c \ + ppSubVersion.c \ + ppSubBackground.c \ + ppSubCamera.c \ + ppSubData.c \ + ppSubErrorCodes.c \ + ppSubFiles.c \ + ppSubLoop.c \ + ppSubDefineOutput.c \ + ppSubExtras.c \ + ppSubMakePSF.c \ + ppSubMatchPSFs.c \ + ppSubReadoutInverse.c \ + ppSubReadoutJpeg.c \ + ppSubReadoutPhotometry.c \ + ppSubReadoutStats.c \ + ppSubReadoutSubtract.c \ + ppSubSetMasks.c \ + ppSubReadoutRenorm.c \ ppSubVarianceFactors.c @@ -42,4 +45,16 @@ ppSub.h + +### Error codes. +BUILT_SOURCES = ppSubErrorCodes.h ppSubErrorCodes.c +CLEANFILES = ppSubErrorCodes.h ppSubErrorCodes.c + +ppSubErrorCodes.h : ppSubErrorCodes.dat ppSubErrorCodes.h.in + $(ERRORCODES) --data=ppSubErrorCodes.dat --outdir=. ppSubErrorCodes.h + +ppSubErrorCodes.c : ppSubErrorCodes.dat ppSubErrorCodes.c.in ppSubErrorCodes.h + $(ERRORCODES) --data=ppSubErrorCodes.dat --outdir=. ppSubErrorCodes.c + + clean-local: -rm -f TAGS @@ -48,3 +63,2 @@ tags: etags `find . -name \*.[ch] -print` - Index: trunk/ppSub/src/ppSub.c =================================================================== --- trunk/ppSub/src/ppSub.c (revision 23688) +++ trunk/ppSub/src/ppSub.c (revision 23740) @@ -28,4 +28,5 @@ psLibInit(NULL); + ppSubData *data = NULL; // Processing data pmConfig *config = pmConfigRead(&argc, argv, PPSUB_RECIPE); // Configuration if (!config) { @@ -49,5 +50,7 @@ } - if (!ppSubArgumentsSetup(argc, argv, config)) { + data = ppSubDataAlloc(); // Processing data + + if (!ppSubArguments(argc, argv, config, data)) { psErrorStackPrint(stderr, "Error reading arguments.\n"); exitValue = PS_EXIT_CONFIG_ERROR; @@ -55,5 +58,5 @@ } - if (!ppSubCamera(config)) { + if (!ppSubCamera(config, data)) { psErrorStackPrint(stderr, "Error setting up camera.\n"); exitValue = PS_EXIT_CONFIG_ERROR; @@ -61,11 +64,5 @@ } - if (!ppSubArgumentsParse(config)) { - psErrorStackPrint(stderr, "Error reading arguments.\n"); - exitValue = PS_EXIT_CONFIG_ERROR; - goto die; - } - - if (!ppSubLoop(config)) { + if (!ppSubLoop(config, data)) { psErrorStackPrint(stderr, "Error performing subtraction.\n"); exitValue = PS_EXIT_SYS_ERROR; @@ -77,6 +74,8 @@ psTimerStop(); + psFree(data); + psFree(config); + pmVisualClose(); //close plot windows, if -visual is set - psFree(config); pmModelClassCleanup(); pmConfigDone(); Index: trunk/ppSub/src/ppSub.h =================================================================== --- trunk/ppSub/src/ppSub.h (revision 23688) +++ trunk/ppSub/src/ppSub.h (revision 23740) @@ -14,6 +14,9 @@ #define PP_SUB_H +#include #include #include + +#include "ppSubErrorCodes.h" /// @addtogroup ppSub @@ -24,13 +27,24 @@ // Output files, for activation/deactivation typedef enum { - PPSUB_FILES_IMAGE = 0x01, // Image files - PPSUB_FILES_PHOT = 0x02, // Photometry files - PPSUB_FILES_ALL = 0xFF, // All files + PPSUB_FILES_INPUT = 0x01, // Input files + PPSUB_FILES_CONV = 0x02, // Convolved files (output) + PPSUB_FILES_SUB = 0x04, // Subtracted files (output) + PPSUB_FILES_INV = 0x08, // Inverse subtracted files (output) + PPSUB_FILES_PSF = 0x10, // PSF files (output) + PPSUB_FILES_PHOT_SUB = 0x20, // Subtraction photometry files (output) + PPSUB_FILES_PHOT_INV = 0x40, // Inverse subtraction photometry files (output) + PPSUB_FILES_PHOT = 0x80, // General photometry files (internal) + PPSUB_FILES_ALL = 0xFF, // All files } ppSubFiles; /// Data for processing typedef struct { - psErrorCode quality; /// Quality code; 0 for no problem - psMetadata *stats; /// Statistics + psErrorCode quality; // Quality code; 0 for no problem + bool photometry; // Perform photometry? + bool inverse; // Output inverse subtraction as well? + psString stamps; // Stamps file + pmPSF *psf; // Point Spread Function + FILE *statsFile; // Statistics file + psMetadata *stats; // Statistics } ppSubData; @@ -39,70 +53,62 @@ /// Setup the arguments parsing -bool ppSubArgumentsSetup(int argc, char *argv[], ///< Command-line arguments - pmConfig *config ///< Configuration - ); - -/// Parse the arguments -bool ppSubArgumentsParse(pmConfig *config ///< Configuration +bool ppSubArguments(int argc, char *argv[], ///< Command-line arguments + pmConfig *config, ///< Configuration + ppSubData *data ///< Processing data ); /// Parse the camera input -bool ppSubCamera(pmConfig *config ///< Configuration +bool ppSubCamera(pmConfig *config, ///< Configuration + ppSubData *data ///< Processing data ); /// Loop over the FPA hierarchy -bool ppSubLoop(pmConfig *config ///< Configuration +bool ppSubLoop(pmConfig *config, ///< Configuration + ppSubData *data ///< Processing data ); /// Perform PSF-matched image subtraction on the readout bool ppSubReadout(pmConfig *config, ///< Configuration - ppSubData *data, ///< Processing data - const pmFPAview *view ///< View of readout to subtract + ppSubData *data ///< Processing data ); /// Generate (if needed) and set or update the masks for input and reference images -bool ppSubSetMasks(pmConfig *config, ///< Configuration - const pmFPAview *view ///< View of active readout +bool ppSubSetMasks(pmConfig *config ///< Configuration ); /// Generate the PSF-matching kernel and convolve the images as needed. Most of this function involves /// looking up the parameters in the recipe and supplying them to the function pmSubtractionMatch() -bool ppSubMatchPSFs(pmConfig *config, ///< Configuration - ppSubData *data, ///< Processing data - const pmFPAview *view ///< View of active readout +bool ppSubMatchPSFs(pmConfig *config, ///< Configuration + ppSubData *data ///< Processing data ); /// Generate the output readout and pass the kernel info to the header -bool ppSubDefineOutput(pmConfig *config, ///< Configuration - const pmFPAview *view ///< View of active readout +bool ppSubDefineOutput(const char *name,///< Name of output to define + pmConfig *config ///< Configuration ); /// Photometry stage 1: measure the PSF from the minuend image -bool ppSubMakePSF(pmConfig *config, ///< Configuration - ppSubData *data, ///< Processing data - const pmFPAview *view ///< View of active readout +bool ppSubMakePSF(pmConfig *config, ///< Configuration + ppSubData *data ///< Processing data ); /// Perform the actual image subtraction, update output concepts -bool ppSubReadoutSubtract(pmConfig *config, ///< Configuration - const pmFPAview *view ///< View of active readout +bool ppSubReadoutSubtract(pmConfig *config ///< Configuration ); /// Photometry stage 2: find and measure sources on the subtracted image -bool ppSubReadoutPhotometry(pmConfig *config, ///< Configuration - ppSubData *data, ///< Processing data - const pmFPAview *view ///< View of active readout +bool ppSubReadoutPhotometry(const char *name, ///< Name of file to photometer + pmConfig *config, ///< Configuration + ppSubData *data ///< Processing data ); /// Renormalize, update headers and generate JPEGs bool ppSubReadoutUpdate(pmConfig *config, ///< Configuration - ppSubData *data, ///< Processing data - const pmFPAview *view ///< View of active readout + ppSubData *data ///< Processing data ); /// Higher-order background subtraction -bool ppSubBackground(pmConfig *config, ///< Configuration - const pmFPAview *view ///< View to readout +bool ppSubBackground(pmConfig *config ///< Configuration ); @@ -121,4 +127,40 @@ ); + +/// Activate or deactivate files +void ppSubFilesActivate(pmConfig *config, // Configuration + ppSubFiles files, // File to activate/deactivate + bool state // Activation state + ); + +/// Generate a view suitable for a readout +/// +/// Assumes we're working with skycells +pmFPAview *ppSubViewReadout(void); + +/// Iterate down the FPA hierarchy, opening files +bool ppSubFilesIterateDown(pmConfig *config, // Configuration + ppSubFiles files // Files to open + ); + +/// Iterate up the FPA hierarchy, closing files +bool ppSubFilesIterateUp(pmConfig *config, // Configuration + ppSubFiles files // Files to open + ); + +/// Collect statistics +bool ppSubReadoutStats(pmConfig *config,// Configuration + ppSubData *data // Processing data + ); + +/// Generate JPEG images +bool ppSubReadoutJpeg(pmConfig *config // Configuration + ); + +/// Generate inverse subtraction +bool ppSubReadoutInverse(pmConfig *config // Configuration + ); + + // Copy every instance of a single keyword from one metadata to another bool psMetadataCopySingle(psMetadata *target, psMetadata *source, const char *name); Index: trunk/ppSub/src/ppSubArguments.c =================================================================== --- trunk/ppSub/src/ppSubArguments.c (revision 23688) +++ trunk/ppSub/src/ppSubArguments.c (revision 23740) @@ -41,136 +41,4 @@ } -// Get a float-point value from the command-line or recipe, and add it to the arguments -#define VALUE_ARG_RECIPE_FLOAT(ARGNAME, RECIPENAME, TYPE) { \ - ps##TYPE value = psMetadataLookup##TYPE(NULL, config->arguments, ARGNAME); \ - if (isnan(value)) { \ - bool mdok; \ - value = psMetadataLookup##TYPE(&mdok, recipe, RECIPENAME); \ - if (!mdok) { \ - psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Unable to find %s in recipe %s", \ - RECIPENAME, PPSUB_RECIPE); \ - goto ERROR; \ - } \ - } \ - psMetadataAdd##TYPE(recipe, PS_LIST_TAIL, RECIPENAME, PS_META_REPLACE, NULL, value); \ -} - -// Get an integer value from the command-line or recipe, and add it to the arguments -#define VALUE_ARG_RECIPE_INT(ARGNAME, RECIPENAME, TYPE, UNSET) { \ - ps##TYPE value = psMetadataLookup##TYPE(NULL, config->arguments, ARGNAME); \ - if (value == UNSET) { \ - bool mdok; \ - value = psMetadataLookup##TYPE(&mdok, recipe, RECIPENAME); \ - if (!mdok) { \ - psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Unable to find %s in recipe %s", \ - RECIPENAME, PPSUB_RECIPE); \ - goto ERROR; \ - } \ - } \ - psMetadataAdd##TYPE(recipe, PS_LIST_TAIL, RECIPENAME, PS_META_REPLACE, NULL, value); \ -} - -/** - * Get a string value from the command-line and add it to the target - */ -static bool valueArgStr(psMetadata *arguments, // Command-line arguments - const char *argName, // Argument name in the command-line arguments - const char *mdName, // Name for value in the metadata - psMetadata *target // Target metadata to which to add value - ) -{ - psString value = psMetadataLookupStr(NULL, arguments, argName); // Value of interest - if (value && strlen(value) > 0) { - return psMetadataAddStr(target, PS_LIST_TAIL, mdName, PS_META_REPLACE, NULL, value); - } - return false; -} - -/** - * Get a string value from the command-line or recipe and add it to the target - */ -static bool valueArgRecipeStr(psMetadata *arguments, // Command-line arguments - psMetadata *recipe, // Recipe - const char *argName, // Argument name in the command-line arguments - const char *mdName, // Name for value in the metadata and recipe - psMetadata *target // Target metadata to which to add value - ) -{ - bool mdok; // Status of MD lookup - psString value = psMetadataLookupStr(&mdok, arguments, argName); // Value of interest - if (!value) { - value = psMetadataLookupStr(&mdok, recipe, mdName); - if (!value) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Unable to find %s in recipe %s", - mdName, PPSUB_RECIPE); - return false; - } - } - return psMetadataAddStr(target, PS_LIST_TAIL, mdName, PS_META_REPLACE, NULL, value); -} - -/** - * Get a vector from the command-line or recipe, and add it to the target - */ -static bool vectorArgRecipe(psMetadata *arguments, // Command-line arguments - const char *argName, // Argument name in the command-line arguments - const psMetadata *recipe, // Recipe - const char *recipeName, // Name for value in the recipe - psMetadata *target, // Target to which to add value - psElemType type // Type for vector - ) -{ - psVector *vector; // Vector - psString string = psMetadataLookupStr(NULL, arguments, argName); // String from arguments - if (string) { - psArray *array = psStringSplitArray(string, ", ", false); // Array of strings - vector = psVectorAlloc(array->n, type); - for (int i = 0; i < array->n; i++) { - const char *subString = array->data[i]; // String with a value - char *end; // Ptr to end of string parsed - - switch (type) { - case PS_TYPE_F32: - vector->data.F32[i] = strtof(subString, &end); - break; - case PS_TYPE_S32: { - long value = strtol(subString, &end, 10); - if (value > PS_MAX_S32 || value < PS_MIN_S32) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - "%s list includes value beyond S32 representation: %s", - argName, string); - psFree(vector); - return false; - } - vector->data.S32[i] = value; - break; - } - default: - psAbort("Unsupported type: %x\n", type); - } - if (end == subString) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Unable to decipher %s list: %s", - argName, string); - psFree(vector); - return false; - } - } - psFree(array); - } else { - vector = psMetadataLookupPtr(NULL, recipe, recipeName); - if (!psMemCheckVector(vector) || vector->type.type != type) { - psError(PS_ERR_BAD_PARAMETER_TYPE, false, "%s in recipe %s is not a vector of type F32.", - recipeName, PPSUB_RECIPE); - return false; - } - psMemIncrRefCounter(vector); - } - - psMetadataAddVector(target, PS_LIST_TAIL, recipeName, PS_META_REPLACE, NULL, vector); - psFree(vector); // Drop reference - - return true; -} - /** * Add a single filename to the arguments as an array, so that it can be used with pmFPAfileBindFromArgs, etc @@ -189,15 +57,7 @@ } -bool ppSubArgumentsSetup(int argc, char *argv[], pmConfig *config) +bool ppSubArguments(int argc, char *argv[], pmConfig *config, ppSubData *data) { - // int argnum; // Argument Number assert(config); - - psMetadata *recipe = psMetadataLookupMetadata(NULL, config->recipes, PPSUB_RECIPE); // Recipe for ppSim - if (!recipe) { - psError(PS_ERR_UNEXPECTED_NULL, false, "Unable to find recipe %s", PPSUB_RECIPE); - return false; - } - psMetadata *arguments = config->arguments; // Command-line arguments @@ -212,46 +72,9 @@ psMetadataAddStr(arguments, PS_LIST_TAIL, "-kernel", 0, "Precalculated kernel to apply", NULL); psMetadataAddStr(arguments, PS_LIST_TAIL, "-stats", 0, "Statistics file", NULL); - psMetadataAddF32(arguments, PS_LIST_TAIL, "-region", 0, "Size of iso-kernel region", NAN); - psMetadataAddS32(arguments, PS_LIST_TAIL, "-size", 0, "Kernel half-size (pixels)", 0); - psMetadataAddS32(arguments, PS_LIST_TAIL, "-order", 0, "Spatial polynomial order", -1); - psMetadataAddStr(arguments, PS_LIST_TAIL, "-type", 0, - "Kernel type (ISIS|POIS|SPAM|FRIES|GUNK|RINGS)", NULL); - psMetadataAddF32(arguments, PS_LIST_TAIL, "-penalty", 0, "Penalty for wideness", NAN); - psMetadataAddStr(arguments, PS_LIST_TAIL, "-isis-widths", 0, - "ISIS Gaussian FWHMs (comma-separated)", NULL); - psMetadataAddStr(arguments, PS_LIST_TAIL, "-isis-orders", 0, - "ISIS polynomial orders (comma-separated)", NULL); - psMetadataAddS32(arguments, PS_LIST_TAIL, "-rings-order", 0, "RINGS polynomial order", -1); - psMetadataAddS32(arguments, PS_LIST_TAIL, "-inner", 0, "SPAM and FRIES inner radius", -1); - psMetadataAddS32(arguments, PS_LIST_TAIL, "-spam-binning", 0, "SPAM kernel binning", 0); - psMetadataAddF32(arguments, PS_LIST_TAIL, "-spacing", 0, "Typical stamp spacing (pixels)", NAN); - psMetadataAddS32(arguments, PS_LIST_TAIL, "-footprint", 0, "Stamp footprint half-size (pixels)", -1); - psMetadataAddF32(arguments, PS_LIST_TAIL, "-source-radius", 0, "Source matching radius (pixels)", NAN); - psMetadataAddS32(arguments, PS_LIST_TAIL, "-stride", 0, "Size of convolution patches (pixels)", -1); - psMetadataAddF32(arguments, PS_LIST_TAIL, "-threshold", 0, "Minimum threshold for stamps (ADU)", NAN); - psMetadataAddS32(arguments, PS_LIST_TAIL, "-iter", 0, "Number of rejection iterations", -1); - psMetadataAddF32(arguments, PS_LIST_TAIL, "-rej", 0, "Rejection thresold (sigma)", NAN); - psMetadataAddF32(arguments, PS_LIST_TAIL, "-sys", 0, "Relative systematic error in kernel", NAN); - psMetadataAddImageMask(arguments, PS_LIST_TAIL, "-mask-bad", 0, "Mask value for bad pixels", 0); - psMetadataAddImageMask(arguments, PS_LIST_TAIL, "-mask-poor", 0, "Mask value for poor pixels", 0); - psMetadataAddF32(arguments, PS_LIST_TAIL, "-poor-frac", 0, "Fraction of variance for poor pixels", NAN); - psMetadataAddF32(arguments, PS_LIST_TAIL, "-badfrac", 0, "Maximum fraction of bad pixels to accept", 1.0); - psMetadataAddBool(arguments, PS_LIST_TAIL, "-reverse", 0, "Reverse sense of subtraction?", false); - psMetadataAddBool(arguments, PS_LIST_TAIL, "-generate-mask", 0, "Generate mask if not supplied?", false); - psMetadataAddStr(arguments, PS_LIST_TAIL, "-stamps", 0, - "Stamps filename; file has x,y on each line", NULL); - psMetadataAddBool(arguments, PS_LIST_TAIL, "-opt", 0, - "Derive optimum parameters for ISIS kernels?", false); - psMetadataAddF32(arguments, PS_LIST_TAIL, "-opt-min", 0, "Minimum value for optimum kernel search", NAN); - psMetadataAddF32(arguments, PS_LIST_TAIL, "-opt-max", 0, "Minimum value for optimum kernel search", NAN); - psMetadataAddF32(arguments, PS_LIST_TAIL, "-opt-step", 0, "Step value for optimum kernel search", NAN); - psMetadataAddF32(arguments, PS_LIST_TAIL, "-opt-tol", 0, "Tolerance for optimum kernel search", NAN); - psMetadataAddS32(arguments, PS_LIST_TAIL, "-opt-order", 0, "Maximum order for optimum kernel search", -1); - psMetadataAddBool(arguments, PS_LIST_TAIL, "-dual", 0, "Dual convolution", false); - psMetadataAddBool(arguments, PS_LIST_TAIL, "-photometry", 0, "Perform photometry?", false); + psMetadataAddStr(arguments, PS_LIST_TAIL, "-stamps", 0, "Stamps filename; x,y on each line", NULL); psMetadataAddS32(arguments, PS_LIST_TAIL, "-threads", 0, "Number of threads", 0); - psMetadataAddS32(arguments, PS_LIST_TAIL, "-bin1", 0, "Binning factor for first level", 0); - psMetadataAddS32(arguments, PS_LIST_TAIL, "-bin2", 0, "Binning factor for second level", 0); psMetadataAddStr(arguments, PS_LIST_TAIL, "-dumpconfig", 0, "file to dump configuration to", NULL); + psMetadataAddBool(arguments, PS_LIST_TAIL, "-photometry", 0, "Perform photometry?", NULL); + psMetadataAddBool(arguments, PS_LIST_TAIL, "-inverse", 0, "Generate inverse subtractions?", NULL); psMetadataAddBool(arguments, PS_LIST_TAIL, "-visual", 0, "Show diagnostic plots", NULL); @@ -302,80 +125,16 @@ } - if (psMetadataLookupBool(NULL, arguments, "-photometry")) { - psMetadataAddBool(recipe, PS_LIST_TAIL, "PHOTOMETRY", PS_META_REPLACE, "Perform photometry?", true); - } + data->stamps = psMemIncrRefCounter(psMetadataLookupStr(NULL, arguments, "-stamps")); - return true; -} - -bool ppSubArgumentsParse(pmConfig *config) -{ - assert(config); - - psMetadata *arguments = config->arguments; // Command-line arguments - - valueArgStr(arguments, "-stats", "STATS", arguments); - valueArgStr(arguments, "-stamps", "STAMPS", arguments); - - psMetadata *recipe = psMetadataLookupMetadata(NULL, config->recipes, PPSUB_RECIPE); // Recipe for ppSim - if (!recipe) { - psError(PS_ERR_UNEXPECTED_NULL, false, "Unable to find recipe %s", PPSUB_RECIPE); - goto ERROR; - } - - VALUE_ARG_RECIPE_FLOAT("-region", "REGION.SIZE", F32); - VALUE_ARG_RECIPE_INT("-size", "KERNEL.SIZE", S32, 0); - VALUE_ARG_RECIPE_INT("-order", "SPATIAL.ORDER", S32, -1); - VALUE_ARG_RECIPE_FLOAT("-spacing", "STAMP.SPACING", F32); - VALUE_ARG_RECIPE_INT("-rings-order", "RINGS.ORDER", S32, -1); - VALUE_ARG_RECIPE_INT("-inner", "INNER", S32, -1); - VALUE_ARG_RECIPE_INT("-spam-binning", "SPAM.BINNING", S32, 0); - VALUE_ARG_RECIPE_INT("-footprint", "STAMP.FOOTPRINT", S32, -1); - VALUE_ARG_RECIPE_FLOAT("-source-radius", "SOURCE.RADIUS", F32); - VALUE_ARG_RECIPE_INT("-stride", "STRIDE", S32, -1); - VALUE_ARG_RECIPE_FLOAT("-threshold", "STAMP.THRESHOLD", F32); - VALUE_ARG_RECIPE_INT("-iter", "ITER", S32, -1); - VALUE_ARG_RECIPE_FLOAT("-rej", "REJ", F32); - VALUE_ARG_RECIPE_FLOAT("-sys", "SYS", F32); - VALUE_ARG_RECIPE_FLOAT("-badfrac", "BADFRAC", F32); - VALUE_ARG_RECIPE_FLOAT("-penalty", "PENALTY", F32); - VALUE_ARG_RECIPE_FLOAT("-poor-frac", "POOR.FRACTION", F32); - VALUE_ARG_RECIPE_INT("-bin1", "BIN1", S32, 0); - VALUE_ARG_RECIPE_INT("-bin2", "BIN2", S32, 0); - - valueArgRecipeStr(arguments, recipe, "-mask-in", "MASK.IN", recipe); - valueArgRecipeStr(arguments, recipe, "-mask-bad", "MASK.BAD", recipe); - valueArgRecipeStr(arguments, recipe, "-mask-poor", "MASK.POOR", recipe); - - vectorArgRecipe(arguments, "-isis-widths", recipe, "ISIS.WIDTHS", recipe, PS_TYPE_F32); - vectorArgRecipe(arguments, "-isis-orders", recipe, "ISIS.ORDERS", recipe, PS_TYPE_S32); - - psVector *widths = psMetadataLookupPtr(NULL, recipe, "ISIS.WIDTHS"); // ISIS Gaussian widths - psVector *orders = psMetadataLookupPtr(NULL, recipe, "ISIS.ORDERS"); // ISIS Polynomial orders - if (widths->n != orders->n) { - psError(PS_ERR_BAD_PARAMETER_SIZE, true, "Size of vectors for ISIS widths and orders do not match."); - goto ERROR; - } - - if (psMetadataLookupBool(NULL, arguments, "-opt") || psMetadataLookupBool(NULL, recipe, "OPTIMUM")) { - psMetadataAddBool(recipe, PS_LIST_TAIL, "OPTIMUM", PS_META_REPLACE, - "Derive optimum parameters for ISIS kernels?", true); - VALUE_ARG_RECIPE_FLOAT("-opt-min", "OPTIMUM.MIN", F32); - VALUE_ARG_RECIPE_FLOAT("-opt-max", "OPTIMUM.MAX", F32); - VALUE_ARG_RECIPE_FLOAT("-opt-step","OPTIMUM.STEP", F32); - VALUE_ARG_RECIPE_FLOAT("-opt-tol", "OPTIMUM.TOL", F32); - VALUE_ARG_RECIPE_INT("-opt-order", "OPTIMUM.ORDER", S32, -1); - } - - psMetadataAddBool(arguments, PS_LIST_TAIL, "REVERSE", 0, "Reverse sense of subtraction", - psMetadataLookupBool(NULL, arguments, "-reverse")); - psMetadataAddBool(recipe, PS_LIST_TAIL, "MASK.GENERATE", PS_META_REPLACE, "Generate mask if not supplied", - psMetadataLookupBool(NULL, arguments, "-generate-mask")); - psMetadataAddBool(recipe, PS_LIST_TAIL, "DUAL", PS_META_REPLACE, "Dual convolution?", - psMetadataLookupBool(NULL, arguments, "-dual")); - - // Need to update this because it could have been overwritten by the camera's own recipe - if (psMetadataLookupBool(NULL, arguments, "-photometry")) { - psMetadataAddBool(recipe, PS_LIST_TAIL, "PHOTOMETRY", PS_META_REPLACE, "Perform photometry?", true); + const char *statsName = psMetadataLookupStr(NULL, arguments, "-stats"); // Filename for statistics + if (statsName && strlen(statsName) > 0) { + psString resolved = pmConfigConvertFilename(statsName, config, true, true); // Resolved filename + data->statsFile = fopen(resolved, "w"); + if (!data->statsFile) { + psError(PS_ERR_IO, true, "Unable to open statistics file %s for writing.\n", resolved); + psFree(resolved); + return false; + } + psFree(resolved); } @@ -384,19 +143,5 @@ } - // Translate the kernel type - psString type = psMetadataLookupStr(NULL, arguments, "-type"); // Name of kernel type - if (!type || strlen(type) == 0) { - type = psMetadataLookupStr(NULL, recipe, "KERNEL.TYPE"); - if (!type || strlen(type) == 0) { - psError(PS_ERR_BAD_PARAMETER_VALUE, false, "Unable to find KERNEL.TYPE specified."); - goto ERROR; - } - } - psMetadataAddStr(recipe, PS_LIST_TAIL, "KERNEL.TYPE", PS_META_REPLACE, "Type of kernel", type); - - psTrace("ppSub", 1, "Done reading command-line arguments\n"); - - // XXX move this to ppSubArguments - int threads = psMetadataLookupS32(NULL, config->arguments, "-threads"); // Number of threads + int threads = psMetadataLookupS32(NULL, arguments, "-threads"); // Number of threads if (threads > 0) { if (!psThreadPoolInit(threads)) { @@ -406,7 +151,6 @@ } + psTrace("ppSub", 1, "Done reading command-line arguments\n"); + return true; - -ERROR: - return false; } Index: trunk/ppSub/src/ppSubBackground.c =================================================================== --- trunk/ppSub/src/ppSubBackground.c (revision 23688) +++ trunk/ppSub/src/ppSubBackground.c (revision 23740) @@ -22,8 +22,7 @@ #include "ppSub.h" -bool ppSubBackground(pmConfig *config, const pmFPAview *view) +bool ppSubBackground(pmConfig *config) { psAssert(config, "Require configuration"); - psAssert(view, "Require view"); bool mdok; // Status of metadata lookups @@ -36,4 +35,5 @@ psImageMaskType maskBad = pmConfigMaskGet("BLANK", config); // Bits to mask + pmFPAview *view = ppSubViewReadout(); // View to readout pmReadout *outRO = pmFPAfileThisReadout(config->files, view, "PPSUB.OUTPUT"); // Output image pmReadout *modelRO = pmFPAfileThisReadout(config->files, view, "PSPHOT.BACKMDL"); // Background model @@ -44,4 +44,5 @@ if (!psphotModelBackground(config, view, "PPSUB.OUTPUT")) { psError(PS_ERR_UNKNOWN, false, "Unable to model background"); + psFree(view); return false; } @@ -50,7 +51,10 @@ if (!modelRO) { psError(PS_ERR_UNEXPECTED_NULL, true, "Unable to find background model"); + psFree(view); return false; } } + psFree(view); + psImageBinning *binning = psMetadataLookupPtr(&mdok, modelRO->analysis, "PSPHOT.BACKGROUND.BINNING"); // Binning for model Index: trunk/ppSub/src/ppSubCamera.c =================================================================== --- trunk/ppSub/src/ppSubCamera.c (revision 23688) +++ trunk/ppSub/src/ppSubCamera.c (revision 23740) @@ -134,5 +134,5 @@ -bool ppSubCamera(pmConfig *config) +bool ppSubCamera(pmConfig *config, ppSubData *data) { psAssert(config, "Require configuration"); @@ -147,5 +147,5 @@ pmFPAfile *inVar = defineInputFile(config, input, "PPSUB.INPUT.VARIANCE", "INPUT.VARIANCE", PM_FPA_FILE_VARIANCE); - defineInputFile(config, input, "PPSUB.INPUT.SOURCES", "INPUT.SOURCES", PM_FPA_FILE_CMF); + defineInputFile(config, NULL, "PPSUB.INPUT.SOURCES", "INPUT.SOURCES", PM_FPA_FILE_CMF); // Reference image @@ -158,25 +158,5 @@ pmFPAfile *refVar = defineInputFile(config, ref, "PPSUB.REF.VARIANCE", "REF.VARIANCE", PM_FPA_FILE_VARIANCE); - defineInputFile(config, ref, "PPSUB.REF.SOURCES", "REF.SOURCES", PM_FPA_FILE_CMF); - - - // Output image - pmFPAfile *output = defineOutputFile(config, input, true, "PPSUB.OUTPUT", PM_FPA_FILE_IMAGE); - pmFPAfile *outMask = defineOutputFile(config, output, false, "PPSUB.OUTPUT.MASK", PM_FPA_FILE_MASK); - if (!output || !outMask) { - psError(PS_ERR_UNKNOWN, false, "Unable to define output files"); - return false; - } - output->save = true; - outMask->save = true; - pmFPAfile *outVar = NULL; - if (inVar && refVar) { - outVar = defineOutputFile(config, output, false, "PPSUB.OUTPUT.VARIANCE", PM_FPA_FILE_VARIANCE); - if (!outVar) { - psError(PS_ERR_UNKNOWN, false, "Unable to define output files"); - return false; - } - outVar->save = true; - } + defineInputFile(config, NULL, "PPSUB.REF.SOURCES", "REF.SOURCES", PM_FPA_FILE_CMF); @@ -216,4 +196,66 @@ + // Now that the camera has been determined, we can read the recipe + psMetadata *recipe = psMetadataLookupMetadata(NULL, config->recipes, PPSUB_RECIPE); // Recipe for ppSim + if (!recipe) { + psError(PS_ERR_UNEXPECTED_NULL, false, "Unable to find recipe %s", PPSUB_RECIPE); + return false; + } + if (psMetadataLookupBool(NULL, config->arguments, "-photometry")) { + psMetadataAddBool(recipe, PS_LIST_TAIL, "PHOTOMETRY", PS_META_REPLACE, + "Perform photometry?", true); + } + if (psMetadataLookupBool(NULL, config->arguments, "-inverse")) { + psMetadataAddBool(recipe, PS_LIST_TAIL, "INVERSE", PS_META_REPLACE, + "Generate inverse subtractions?", true); + } + + data->inverse = psMetadataLookupBool(NULL, recipe, "INVERSE"); + data->photometry = psMetadataLookupBool(NULL, recipe, "PHOTOMETRY"); + + + // Output image + pmFPAfile *output = defineOutputFile(config, inConvImage, true, "PPSUB.OUTPUT", PM_FPA_FILE_IMAGE); + pmFPAfile *outMask = defineOutputFile(config, output, false, "PPSUB.OUTPUT.MASK", PM_FPA_FILE_MASK); + if (!output || !outMask) { + psError(PS_ERR_UNKNOWN, false, "Unable to define output files"); + return false; + } + output->save = true; + outMask->save = true; + if (inVar && refVar) { + pmFPAfile *outVar = defineOutputFile(config, output, false, "PPSUB.OUTPUT.VARIANCE", + PM_FPA_FILE_VARIANCE); + if (!outVar) { + psError(PS_ERR_UNKNOWN, false, "Unable to define output files"); + return false; + } + outVar->save = true; + } + + pmFPAfile *inverse = NULL; // Inverse output image + if (data->inverse) { + // Inverse output image + inverse = defineOutputFile(config, output, true, "PPSUB.INVERSE", PM_FPA_FILE_IMAGE); + pmFPAfile *invMask = defineOutputFile(config, inverse, false, "PPSUB.INVERSE.MASK", + PM_FPA_FILE_MASK); + if (!inverse || !invMask) { + psError(PS_ERR_UNKNOWN, false, "Unable to define output files"); + return false; + } + inverse->save = true; + invMask->save = true; + if (inVar && refVar) { + pmFPAfile *invVar = defineOutputFile(config, inverse, false, "PPSUB.INVERSE.VARIANCE", + PM_FPA_FILE_VARIANCE); + if (!invVar) { + psError(PS_ERR_UNKNOWN, false, "Unable to define output files"); + return false; + } + invVar->save = true; + } + } + + // Output JPEGs pmFPAfile *jpeg1 = pmFPAfileDefineOutput(config, NULL, "PPSUB.OUTPUT.JPEG1"); @@ -245,25 +287,7 @@ } - psMetadata *recipe = psMetadataLookupMetadata(NULL, config->recipes, PPSUB_RECIPE); // Recipe for ppSim - if (!recipe) { - psError(PS_ERR_UNEXPECTED_NULL, false, "Unable to find recipe %s", PPSUB_RECIPE); - return false; - } - // psPhot input - if (psMetadataLookupBool(NULL, recipe, "PHOTOMETRY")) { + if (data->photometry) { psphotModelClassInit(); // load implementation-specific models - - // Internal-ish file for getting the PSF from the minuend - pmFPAfile *psf = pmFPAfileDefineOutputFromFile(config, output, "PSPHOT.PSF.LOAD"); - if (!psf) { - psError(PS_ERR_IO, false, "Failed to build FPA from PSPHOT.PSF.LOAD"); - return false; - } - if (psf->type != PM_FPA_FILE_PSF) { - psError(PS_ERR_IO, true, "PSPHOT.PSF.LOAD is not of type PSF"); - return false; - } - pmFPAfileActivate(config->files, false, "PSPHOT.PSF.LOAD"); pmFPAfile *psphot = pmFPAfileDefineFromFPA(config, output->fpa, 1, 1, "PSPHOT.INPUT"); @@ -276,4 +300,17 @@ return false; } + pmFPAfileActivate(config->files, false, "PSPHOT.INPUT"); + + // Internal-ish file for getting the PSF from the minuend + pmFPAfile *psf = pmFPAfileDefineOutputFromFile(config, psphot, "PSPHOT.PSF.LOAD"); + if (!psf) { + psError(PS_ERR_IO, false, "Failed to build FPA from PSPHOT.PSF.LOAD"); + return false; + } + if (psf->type != PM_FPA_FILE_PSF) { + psError(PS_ERR_IO, true, "PSPHOT.PSF.LOAD is not of type PSF"); + return false; + } + pmFPAfileActivate(config->files, false, "PSPHOT.PSF.LOAD"); if (!psphotDefineFiles(config, psphot)) { @@ -281,4 +318,26 @@ return false; } + + // Deactivate psphot output sources --- we want to define output source files of our own + pmFPAfile *psphotOutput = pmFPAfileSelectSingle(config->files, "PSPHOT.OUTPUT", 0); + psphotOutput->save = false; + + pmFPAfile *outSources = defineOutputFile(config, output, false, "PPSUB.OUTPUT.SOURCES", + PM_FPA_FILE_CMF); + if (!outSources) { + psError(PS_ERR_UNKNOWN, false, "Unable to set up output source file."); + return false; + } + outSources->save = true; + + if (data->inverse) { + pmFPAfile *invSources = defineOutputFile(config, inverse, false, "PPSUB.INVERSE.SOURCES", + PM_FPA_FILE_CMF); + if (!invSources) { + psError(PS_ERR_UNKNOWN, false, "Unable to set up inverse source file."); + return false; + } + invSources->save = true; + } } Index: trunk/ppSub/src/ppSubData.c =================================================================== --- trunk/ppSub/src/ppSubData.c (revision 23688) +++ trunk/ppSub/src/ppSubData.c (revision 23740) @@ -11,16 +11,19 @@ -// Image files to activate/deactivate -static const char *imageFiles[] = { "PPSUB.OUTPUT", "PPSUB.OUTPUT.MASK", "PPSUB.OUTPUT.VARIANCE", - "PPSUB.OUTPUT.KERNELS", "PPSUB.OUTPUT.JPEG1", "PPSUB.OUTPUT.JPEG2", - "PPSUB.INPUT.CONV", "PPSUB.INPUT.CONV.MASK", "PPSUB.INPUT.CONV.VARIANCE", - "PPSUB.REF.CONV", "PPSUB.REF.CONV.MASK", "PPSUB.REF.CONV.VARIANCE", - NULL }; - - - -static void subOptionsFree(ppSubData *options) +static void subDataFree(ppSubData *data) { - psFree(options->stats); + if (data->statsFile) { + psString stats = psMetadataConfigFormat(data->stats); // Statistics to output + if (!stats || strlen(stats) == 0) { + psWarning("Unable to generate statistics file."); + } else { + fprintf(data->statsFile, "%s", stats); + } + psFree(stats); + fclose(data->statsFile); + } + psFree(data->stamps); + psFree(data->psf); + psFree(data->stats); return; } @@ -28,12 +31,17 @@ ppSubData *ppSubDataAlloc(void) { - ppSubData *options = psAlloc(sizeof(ppSubData)); // Processing data, to return - psMemSetDeallocator(options, (psFreeFunc)subOptionsFree); + ppSubData *data = psAlloc(sizeof(ppSubData)); // Processing data, to return + psMemSetDeallocator(data, (psFreeFunc)subDataFree); - options->quality = 0; - options->stats = psMetadataAlloc(); - psMetadataAddS32(options->stats, PS_LIST_TAIL, "QUALITY", 0, "Data quality", 0); + data->quality = 0; + data->photometry = false; + data->inverse = false; + data->stamps = NULL; + data->psf = NULL; + data->statsFile = NULL; + data->stats = psMetadataAlloc(); + psMetadataAddS32(data->stats, PS_LIST_TAIL, "QUALITY", 0, "Data quality", 0); - return options; + return data; } @@ -49,12 +57,5 @@ } - if (files & PPSUB_FILES_IMAGE) { - for (int i = 0; imageFiles[i]; i++) { - pmFPAfileActivate(config->files, imageFiles[i], false); - } - } - if (files & PPSUB_FILES_PHOT) { - psphotFilesActivate(config, false); - } + ppSubFilesActivate(config, files, false); psErrorClear(); Index: trunk/ppSub/src/ppSubDefineOutput.c =================================================================== --- trunk/ppSub/src/ppSubDefineOutput.c (revision 23688) +++ trunk/ppSub/src/ppSubDefineOutput.c (revision 23740) @@ -21,25 +21,29 @@ #include "ppSub.h" -bool ppSubDefineOutput(pmConfig *config, const pmFPAview *view) +bool ppSubDefineOutput(const char *name, pmConfig *config) { + psAssert(name, "Require name"); psAssert(config, "Require configuration"); - psAssert(view, "Require view"); - pmCell *outCell = pmFPAfileThisCell(config->files, view, "PPSUB.OUTPUT"); // Output cell + pmFPAview *view = ppSubViewReadout(); // View to readout + pmCell *outCell = pmFPAfileThisCell(config->files, view, name); // Output cell pmFPA *outFPA = outCell->parent->parent; // Output FPA - pmHDU *outHDU = outFPA->hdu; // Output HDU + pmHDU *outHDU = outFPA->hdu; // Output HDU if (!outHDU->header) { outHDU->header = psMetadataAlloc(); } - // generate an output readout (first check if it's already there by virtue of kernels) - pmReadout *outRO = pmFPAfileThisReadout(config->files, view, "PPSUB.OUTPUT.KERNELS"); // RO with kernel - if (!outRO) { - outRO = pmReadoutAlloc(outCell); // Output readout: subtraction + // The output readout may already be present if we read in the convolution kernel + pmReadout *outRO = NULL; // Output readout + if (outCell->readouts && outCell->readouts->n > 0 && outCell->readouts->data[0]) { + outRO = psMemIncrRefCounter(outCell->readouts->data[0]); + } else { + outRO = pmReadoutAlloc(outCell); } - // convolved input images + // Convolved input images pmReadout *inConv = pmFPAfileThisReadout(config->files, view, "PPSUB.INPUT.CONV"); // Input readout pmReadout *refConv = pmFPAfileThisReadout(config->files, view, "PPSUB.REF.CONV"); // Reference readout + psFree(view); // Add kernel descrption to header. @@ -55,6 +59,4 @@ if (!kernels) { psError(PS_ERR_UNEXPECTED_NULL, true, "Unable to find SUBTRACTION.KERNEL"); - psFree(inConv); - psFree(refConv); psFree(outRO); return false; @@ -63,15 +65,17 @@ kernels->description); + // Add additional data to the header + pmFPAfile *refFile = psMetadataLookupPtr(NULL, config->files, "PPSUB.REF"); // Reference file + pmFPAfile *inFile = psMetadataLookupPtr(NULL, config->files, "PPSUB.INPUT"); // Input file + psMetadataAddStr(outHDU->header, PS_LIST_TAIL, "PPSUB.REFERENCE", 0, + "Subtraction reference", refFile->filename); + psMetadataAddStr(outHDU->header, PS_LIST_TAIL, "PPSUB.INPUT", 0, + "Subtraction input", inFile->filename); + ppSubVersionHeader(outHDU->header); + + outRO->analysis = psMetadataCopy(outRO->analysis, analysis); -#ifdef TESTING - { - psImage *kernelImage = psMetadataLookupPtr(&mdok, analysis, - "SUBTRACTION.KERNEL.IMAGE"); // Image of kernel - psFits *fits = psFitsOpen("kernel.fits", "w"); - psFitsWriteImage(fits, NULL, kernelImage, 0, NULL); - psFitsClose(fits); - } -#endif + psFree(outRO); return true; Index: trunk/ppSub/src/ppSubErrorCodes.c.in =================================================================== --- trunk/ppSub/src/ppSubErrorCodes.c.in (revision 23740) +++ trunk/ppSub/src/ppSubErrorCodes.c.in (revision 23740) @@ -0,0 +1,37 @@ +/** @file ppSubErrorCodes.c.in + * + * @brief + * + * @ingroup ppSub + * + * @author IfA + * @version $Revision: 1.3 $ $Name: not supported by cvs2svn $ + * @date $Date: 2009-02-05 20:44:04 $ + * Copyright 2009 Institute for Astronomy, University of Hawaii + */ + +#include "pslib.h" +#include "ppSubErrorCodes.h" + +/* + * The line + { PPSUB_ERR_$X{ErrorCode}, "$X{ErrorDescription}"}, + * (without the Xs) + * will be replaced by values from errorCodes.dat + */ +void ppSubErrorRegister(void) +{ + static psErrorDescription errors[] = { + { PPSUB_ERR_BASE, "First value we use; lower values belong to psLib" }, + { PPSUB_ERR_${ErrorCode}, "${ErrorDescription}"}, + }; + static int nerror = PPSUB_ERR_NERROR - PPSUB_ERR_BASE; ///< number of values in enum + + for (int i = 0; i < nerror; i++) { + psErrorDescription *tmp = psAlloc(sizeof(psErrorDescription)); + *tmp = errors[i]; + psErrorRegister(tmp, 1); + psFree(tmp); /* it's on the internal list */ + } + nerror = 0; // don't register more than once +} Index: trunk/ppSub/src/ppSubErrorCodes.dat =================================================================== --- trunk/ppSub/src/ppSubErrorCodes.dat (revision 23740) +++ trunk/ppSub/src/ppSubErrorCodes.dat (revision 23740) @@ -0,0 +1,11 @@ +# +# This file is used to generate pswarpErrorClasses.h +# +BASE = 14000 First value we use; lower values belong to psLib +UNKNOWN Unknown ppSub error code +NOT_IMPLEMENTED Desired feature is not yet implemented +ARGUMENTS Incorrect arguments +CONFIG Problem in configure files +IO Problem in FITS I/O +DATA Problem in data values +NO_OVERLAP No overlap between input and skycell Index: trunk/ppSub/src/ppSubErrorCodes.h.in =================================================================== --- trunk/ppSub/src/ppSubErrorCodes.h.in (revision 23740) +++ trunk/ppSub/src/ppSubErrorCodes.h.in (revision 23740) @@ -0,0 +1,30 @@ +/** @file ppSubErrorCodes.h.in + * + * @brief + * + * @ingroup ppSub + * + * @author IfA + * @version $Revision: 1.2 $ $Name: not supported by cvs2svn $ + * @date $Date: 2009-02-05 20:44:04 $ + * Copyright 2009 Institute for Astronomy, University of Hawaii + */ + +#if !defined(PPSUB_ERROR_CODES_H) +#define PPSUB_ERROR_CODES_H +/* + * The line + * PPSUB_ERR_$X{ErrorCode}, + * (without the X) + * + * will be replaced by values from errorCodes.dat + */ +typedef enum { + PPSUB_ERR_BASE = 14000, + PPSUB_ERR_${ErrorCode}, + PPSUB_ERR_NERROR +} ppSubErrorCode; + +void ppSubErrorRegister(void); + +#endif Index: trunk/ppSub/src/ppSubFiles.c =================================================================== --- trunk/ppSub/src/ppSubFiles.c (revision 23740) +++ trunk/ppSub/src/ppSubFiles.c (revision 23740) @@ -0,0 +1,188 @@ +#ifdef HAVE_CONFIG_H +#include +#endif + +#include +#include +#include +#include + +#include "ppSub.h" + + +// Input files +static const char *inputFiles[] = { "PPSUB.INPUT", "PPSUB.INPUT.MASK", "PPSUB.INPUT.VARIANCE", + "PPSUB.INPUT.SOURCES", "PPSUB.REF", "PPSUB.REF.MASK", + "PPSUB.REF.VARIANCE", "PPSUB.REF.SOURCES", NULL }; + +// Convolved files +static const char *convFiles[] = { "PPSUB.INPUT.CONV", "PPSUB.INPUT.CONV.MASK", "PPSUB.INPUT.CONV.VARIANCE", + "PPSUB.REF.CONV", "PPSUB.REF.CONV.MASK", "PPSUB.REF.CONV.VARIANCE", + NULL }; + +// Subtraction files +static const char *subFiles[] = { "PPSUB.OUTPUT", "PPSUB.OUTPUT.MASK", "PPSUB.OUTPUT.VARIANCE", + "PPSUB.OUTPUT.JPEG1", "PPSUB.OUTPUT.JPEG2", + NULL }; + +// Subtraction photometry +static const char *subPhotFiles[] = { "PPSUB.OUTPUT.SOURCES", NULL }; + +// Inverse subtraction files +static const char *invFiles[] = { "PPSUB.INVERSE", "PPSUB.INVERSE.MASK", "PPSUB.INVERSE.VARIANCE", NULL }; + +// Inverse subtraction photometry +static const char *invPhotFiles[] = { "PPSUB.INVERSE.SOURCES", NULL }; + +// PSF files +static const char *psfFiles[] = { "PSPHOT.PSF.SAVE", NULL }; + +// Calculation (may be either input or output) files +static const char *calcFiles[] = { "PPSUB.OUTPUT.KERNELS", NULL }; + + +// Activate/deactivate a list of files +static void filesActivate(pmConfig *config, // Configuration + const char **files, // List of files + bool state // Activation status to set + ) +{ + for (int i = 0; files[i]; i++) { + pmFPAfileActivate(config->files, state, files[i]); + } + return; +} + +// Activate/deactivate a list of files depending on their 'save' boolean. +// This is so we can activate/deactivate the 'calculation' files, which may be either input or output, which +// is indicated by their 'save' boolean. +static void filesActivateSave(pmConfig *config, // Configuration + const char **files, // List of files + bool save, // Activate when this save state is set + bool state // Activation status to set + ) +{ + for (int i = 0; files[i]; i++) { + pmFPAfile *file = pmFPAfileSelectSingle(config->files, files[i], 0); + if (file && file->save == save) { + pmFPAfileActivate(config->files, state, files[i]); + } + } + return; +} + +void ppSubFilesActivate(pmConfig *config, ppSubFiles files, bool state) +{ + psAssert(config, "Require configuration"); + + if (files & PPSUB_FILES_INPUT) { + filesActivate(config, inputFiles, state); + filesActivateSave(config, calcFiles, false, state); + } + if (files & PPSUB_FILES_CONV) { + filesActivate(config, convFiles, state); + } + if (files & PPSUB_FILES_SUB) { + filesActivate(config, subFiles, state); + filesActivateSave(config, calcFiles, true, state); + } + if (files & PPSUB_FILES_INV) { + filesActivate(config, invFiles, state); + } + if (files & PPSUB_FILES_PHOT_SUB) { + filesActivate(config, subPhotFiles, state); + } + if (files & PPSUB_FILES_PHOT_INV) { + filesActivate(config, invPhotFiles, state); + } + if (files & PPSUB_FILES_PSF) { + filesActivate(config, psfFiles, state); + } + if (files & PPSUB_FILES_PHOT) { + psphotFilesActivate(config, state); + } + + return; +} + + +pmFPAview *ppSubViewReadout(void) +{ + pmFPAview *view = pmFPAviewAlloc(0); + view->chip = view->cell = view->readout = 0; + return view; +} + +bool ppSubFilesIterateDown(pmConfig *config, ppSubFiles files) +{ + psAssert(config, "Require configuration"); + + ppSubFilesActivate(config, PPSUB_FILES_ALL, false); + ppSubFilesActivate(config, files, true); + + pmFPAview *view = pmFPAviewAlloc(0);// View to FPA top + + // FPA + if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) { + return false; + } + + // Chip + view->chip = 0; + if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) { + return false; + } + + // Cell + view->cell = 0; + if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) { + return false; + } + + // Readout + view->readout = 0; + if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) { + return false; + } + + ppSubFilesActivate(config, PPSUB_FILES_ALL, false); + + return true; +} + +bool ppSubFilesIterateUp(pmConfig *config, ppSubFiles files) +{ + psAssert(config, "Require configuration"); + + ppSubFilesActivate(config, PPSUB_FILES_ALL, false); + ppSubFilesActivate(config, files, true); + + pmFPAview *view = ppSubViewReadout(); // View to readout + + // Readout + if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { + return false; + } + + // Cell + view->readout = -1; + if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { + return false; + } + + // Chip + view->cell = -1; + if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { + return false; + } + + // FPA + view->chip = -1; + if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { + return false; + } + + ppSubFilesActivate(config, PPSUB_FILES_ALL, false); + + return true; +} Index: trunk/ppSub/src/ppSubLoop.c =================================================================== --- trunk/ppSub/src/ppSubLoop.c (revision 23688) +++ trunk/ppSub/src/ppSubLoop.c (revision 23740) @@ -17,9 +17,8 @@ #include #include -#include #include "ppSub.h" -bool ppSubLoop(pmConfig *config) +bool ppSubLoop(pmConfig *config, ppSubData *data) { psAssert(config, "Require configuration."); @@ -28,151 +27,146 @@ pmConfigRecipesCull(config, "PPSUB,PPSTATS,PSPHOT,MASKS,JPEG"); - ppSubData *data = ppSubDataAlloc(); // Processing data - bool mdok; // Status of MD lookup - const char *statsName = psMetadataLookupStr(&mdok, config->arguments, "STATS"); // Filename for statistics - FILE *statsFile = NULL; // File stream for statistics - if (statsName && strlen(statsName) > 0) { - psString resolved = pmConfigConvertFilename(statsName, config, true, true); // Resolved filename - statsFile = fopen(resolved, "w"); - if (!statsFile) { - psError(PS_ERR_IO, true, "Unable to open statistics file %s for writing.\n", resolved); - psFree(resolved); - goto ERROR; - } - psFree(resolved); + pmFPAfile *input = psMetadataLookupPtr(NULL, config->files, "PPSUB.INPUT"); + pmFPAfile *reference = psMetadataLookupPtr(NULL, config->files, "PPSUB.REF"); + pmFPAfile *output = psMetadataLookupPtr(NULL, config->files, "PPSUB.OUTPUT"); + psAssert(input && reference && output, "Require files"); + + if (!ppSubFilesIterateDown(config, PPSUB_FILES_INPUT | PPSUB_FILES_CONV)) { + psError(PPSUB_ERR_IO, false, "Unable to load files."); + return false; } - pmFPAfile *input = psMetadataLookupPtr(NULL, config->files, "PPSUB.INPUT"); - if (!input) { - psError(PS_ERR_UNEXPECTED_NULL, false, "Can't find input data!\n"); - goto ERROR; + psTimerStart("PPSUB_MATCH"); + + if (!ppSubSetMasks(config)) { + psError(PS_ERR_UNKNOWN, false, "Unable to set masks."); + return false; } - pmFPAfile *reference = psMetadataLookupPtr(NULL, config->files, "PPSUB.REF"); - if (!reference) { - psError(PS_ERR_UNEXPECTED_NULL, false, "Can't find reference data!\n"); - goto ERROR; + if (!ppSubMatchPSFs(config, data)) { + psError(PS_ERR_UNKNOWN, false, "Unable to match PSFs."); + return false; + } + if (data->quality) { + // Can't do anything at all + return true; } - pmFPAfile *output = psMetadataLookupPtr(NULL, config->files, "PPSUB.OUTPUT"); - if (!output) { - psError(PS_ERR_UNEXPECTED_NULL, false, "Can't find output data!\n"); - goto ERROR; + psMetadataAddF32(data->stats, PS_LIST_TAIL, "TIME_MATCH", 0, "Time to match PSFs", + psTimerClear("PPSUB_MATCH")); + + // Close input files + if (!ppSubFilesIterateUp(config, PPSUB_FILES_INPUT)) { + psError(PPSUB_ERR_IO, false, "Unable to close input files."); + return false; } - pmFPAview *view = pmFPAviewAlloc(0); // Pointer into FPA hierarchy - - // Iterate over the FPA hierarchy - if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) { - goto ERROR; + // Set up subtraction files + if (!ppSubFilesIterateDown(config, PPSUB_FILES_SUB)) { + psError(PPSUB_ERR_IO, false, "Unable to set up subtraction files."); + return false; } - pmChip *inChip; // Input chip of interest - while ((inChip = pmFPAviewNextChip(view, input->fpa, 1)) != NULL) { - pmChip *refChip = pmFPAviewThisChip(view, reference->fpa); // Reference chip of interest - if ((!inChip->file_exists && refChip->file_exists) || - (inChip->file_exists && !refChip->file_exists)) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, "FPA format discrepency between input and reference"); - psFree(view); - goto ERROR; + if (!ppSubDefineOutput("PPSUB.OUTPUT", config)) { + psError(PS_ERR_UNKNOWN, false, "Unable to define output."); + return false; + } + + if (!data->quality && !ppSubMakePSF(config, data)) { + psError(PS_ERR_UNKNOWN, false, "Unable to generate PSF."); + return false; + } + + if (!ppSubReadoutSubtract(config)) { + psError(PS_ERR_UNKNOWN, false, "Unable to subtract images."); + return false; + } + + // Close convolved files + if (!ppSubFilesIterateUp(config, PPSUB_FILES_PSF | PPSUB_FILES_CONV)) { + psError(PPSUB_ERR_IO, false, "Unable to close input files."); + return false; + } + + // Higher order background subtraction using psphot + if (!ppSubBackground(config)) { + psError(PS_ERR_UNKNOWN, false, "Unable to subtract background."); + return false; + } + + if (!ppSubFilesIterateDown(config, PPSUB_FILES_PHOT_SUB)) { + psError(PPSUB_ERR_IO, false, "Unable to set up photometry files."); + return false; + } + + if (!data->quality && !ppSubReadoutPhotometry("PPSUB.OUTPUT", config, data)) { + psError(PS_ERR_UNKNOWN, false, "Unable to perform photometry."); + return false; + } + + if (!ppSubFilesIterateUp(config, PPSUB_FILES_PHOT_SUB)) { + psError(PPSUB_ERR_IO, false, "Unable to set up photometry files."); + return false; + } + + // Perform statistics on the cell + if (!ppSubReadoutStats(config, data)) { + psError(PS_ERR_UNKNOWN, false, "Unable to collect statistics"); + return false; + } + + if (!ppSubReadoutJpeg(config)) { + psError(PS_ERR_UNKNOWN, false, "Unable to update."); + return false; + } + + if (data->inverse) { + // Set up inverse subtraction files + if (!ppSubFilesIterateDown(config, PPSUB_FILES_INV)) { + psError(PPSUB_ERR_IO, false, "Unable to set up inverse files."); + return false; } - if (!inChip->file_exists) { - continue; + if (data->inverse && !ppSubDefineOutput("PPSUB.INVERSE", config)) { + psError(PS_ERR_UNKNOWN, false, "Unable to define inverse."); + return false; } - if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) { - goto ERROR; + if (!ppSubReadoutInverse(config)) { + psError(PS_ERR_UNKNOWN, false, "Unable to invert images."); + return false; } - pmCell *inCell; // Cell of interest - while ((inCell = pmFPAviewNextCell(view, input->fpa, 1)) != NULL) { - pmCell *refCell = pmFPAviewThisCell(view, reference->fpa); // Reference cell of interest - if ((!inCell->file_exists && refCell->file_exists) || - (inCell->file_exists && !refCell->file_exists)) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - "FPA format discrepency between input and reference"); - psFree(view); - goto ERROR; - } - if (!inCell->file_exists) { - continue; - } - if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) { - goto ERROR; - } - - pmReadout *inRO; // Readin of interest - while ((inRO = pmFPAviewNextReadout(view, input->fpa, 1))) { - if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) { - goto ERROR; - } - pmReadout *refRO = pmFPAviewThisReadout(view, reference->fpa);// Reference readout of interest - if (!refRO || (!inRO->data_exists && refRO->data_exists) || - (inRO->data_exists && !refRO->data_exists)) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - "FPA format discrepency between input and reference"); - psFree(view); - goto ERROR; - } - if (!inRO->data_exists) { - continue; - } - - // Perform the analysis - if (!ppSubReadout(config, data, view)) { - psError(PS_ERR_UNKNOWN, false, "Unable to subtract images.\n"); - goto ERROR; - } - - if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) { - goto ERROR; - } - } - - // Perform statistics on the cell - if (statsFile) { - pmFPAfile *output = psMetadataLookupPtr(NULL, config->files, "PPSUB.OUTPUT"); // Output file - if (!output) { - psError(PS_ERR_UNEXPECTED_NULL, true, "Unable to find file PPSUB.OUTPUT.\n"); - goto ERROR; - } - psImageMaskType maskValue = pmConfigMaskGet("MASK.VALUE", config); - ppStatsFPA(data->stats, output->fpa, view, maskValue, config); - } - - if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { - goto ERROR; - } + // Close subtraction files and open inverse photometry files + if (!ppSubFilesIterateUp(config, PPSUB_FILES_SUB)) { + psError(PPSUB_ERR_IO, false, "Unable to close subtraction files."); + return false; } - if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { - goto ERROR; + if (!ppSubFilesIterateDown(config, PPSUB_FILES_PHOT_INV)) { + psError(PPSUB_ERR_IO, false, "Unable to set up inverse files."); + return false; + } + + if (!data->quality && !ppSubReadoutPhotometry("PPSUB.INVERSE", config, data)) { + psError(PS_ERR_UNKNOWN, false, "Unable to perform photometry."); + return false; + } + + // Close inverse subtraction files + if (!ppSubFilesIterateUp(config, PPSUB_FILES_INV | PPSUB_FILES_PHOT_INV)) { + psError(PPSUB_ERR_IO, false, "Unable to close subtraction files."); + return false; + } + } else { + // Close subtraction files + if (!ppSubFilesIterateUp(config, PPSUB_FILES_SUB)) { + psError(PPSUB_ERR_IO, false, "Unable to close subtraction files."); + return false; } } - if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { - goto ERROR; - } - - psFree(view); - - // Write out summary statistics - if (statsFile) { - psMetadataAddF32(data->stats, PS_LIST_TAIL, "TIME_SUB", 0, "Time for subtraction completion", - psTimerMark("ppSub")); - - const char *statsMDC = psMetadataConfigFormat(data->stats); - if (!statsMDC || strlen(statsMDC) == 0) { - psWarning("Unable to generate statistics MDC file.\n"); - } else { - fprintf(statsFile, "%s", statsMDC); - } - psFree((void *)statsMDC); - fclose(statsFile); - } - - psString dump_file = psMetadataLookupStr(&mdok, config->arguments, "-dumpconfig"); + psString dump_file = psMetadataLookupStr(NULL, config->arguments, "-dumpconfig"); if (dump_file) { pmFPAfile *input = psMetadataLookupPtr(NULL, config->files, "PPSUB.INPUT"); // Input file @@ -180,9 +174,4 @@ } - psFree(data); return true; - -ERROR: - psFree(data); - return false; } Index: trunk/ppSub/src/ppSubMakePSF.c =================================================================== --- trunk/ppSub/src/ppSubMakePSF.c (revision 23688) +++ trunk/ppSub/src/ppSubMakePSF.c (revision 23740) @@ -22,8 +22,11 @@ #include "ppSub.h" -bool ppSubMakePSF(pmConfig *config, ppSubData *data, const pmFPAview *view) +bool ppSubMakePSF(pmConfig *config, ppSubData *data) { psAssert(config, "Require configuration"); - psAssert(view, "Require view"); + + if (!data->photometry) { + return true; + } psTimerStart("PPSUB_PHOT"); @@ -44,4 +47,5 @@ pmReadout *minuend = NULL; // Image that will be positive following subtraction pmFPAfile *minuendFile = NULL; // File for minuend image + pmFPAview *view = ppSubViewReadout(); // View to readout if (reverse) { minuend = pmFPAfileThisReadout(config->files, view, "PPSUB.REF.CONV"); @@ -53,7 +57,8 @@ #if 1 + pmFPAfile *photFile = psMetadataLookupPtr(&mdok, config->files, "PSPHOT.INPUT"); // Photometry file +#if 0 + pmReadout *photRO = pmFPAviewThisReadout(view, photFile->fpa); // Readout to photometer pmReadout *template = minuend; - pmFPAfile *photFile = psMetadataLookupPtr(&mdok, config->files, "PSPHOT.INPUT"); // Photometry file - pmReadout *photRO = pmFPAviewThisReadout(view, photFile->fpa); // Readout to photometer if (!photRO) { pmCell *cell = pmFPAviewThisCell(view, photFile->fpa); // Cell to photometer @@ -74,4 +79,17 @@ } #else + if (!pmFPACopy(photFile->fpa, minuendFile->fpa)) { + psError(PS_ERR_UNKNOWN, false, "Unable to copy FPA for photometry"); + psFree(view); + return false; + } + pmReadout *photRO = pmFPAviewThisReadout(view, photFile->fpa); // Readout to photometer + if (psMetadataLookup(photRO->analysis, "PSPHOT.SOURCES")) { + psMetadataRemoveKey(photRO->analysis, "PSPHOT.SOURCES"); + } +#endif + + +#else // Supply the minuend pmFPAfile to psphot as PSPHOT.INPUT: psMetadataAddPtr(config->files, PS_LIST_TAIL, "PSPHOT.INPUT", PS_DATA_UNKNOWN | PS_META_REPLACE, @@ -87,5 +105,6 @@ psErrorStackPrint(stderr, "Unable to determine PSF"); psWarning("Unable to determine PSF --- suspect bad data quality."); - ppSubDataQuality(config, data, PSPHOT_ERR_PSF, PPSUB_FILES_PHOT); + ppSubDataQuality(config, data, PSPHOT_ERR_PSF, PPSUB_FILES_PHOT_SUB | PPSUB_FILES_PHOT_INV); + psFree(view); return true; } @@ -93,4 +112,5 @@ // Record the FWHM in the output header pmReadout *outRO = pmFPAfileThisReadout(config->files, view, "PPSUB.OUTPUT"); // Output readout + psFree(view); pmHDU *hdu = pmHDUFromCell(outRO->parent); // HDU with header psMetadataItemSupplement(hdu->header, psphotRecipe, "FWHM_MAJ"); @@ -100,4 +120,7 @@ psMetadataRemoveKey(photRO->analysis, "PSPHOT.HEADER"); + data->psf = psMemIncrRefCounter(psMetadataLookupPtr(NULL, photRO->parent->parent->analysis, + "PSPHOT.PSF")); + return true; } Index: trunk/ppSub/src/ppSubMatchPSFs.c =================================================================== --- trunk/ppSub/src/ppSubMatchPSFs.c (revision 23688) +++ trunk/ppSub/src/ppSubMatchPSFs.c (revision 23740) @@ -22,12 +22,13 @@ #include "ppSub.h" -bool ppSubMatchPSFs(pmConfig *config, ppSubData *data, const pmFPAview *view) +bool ppSubMatchPSFs(pmConfig *config, ppSubData *data) { psAssert(config, "Require configuration"); - psAssert(view, "Require view"); // Look up recipe values psMetadata *recipe = psMetadataLookupMetadata(NULL, config->recipes, PPSUB_RECIPE); // Recipe for ppSim psAssert(recipe, "We checked this earlier, so it should be here."); + + pmFPAview *view = ppSubViewReadout(); // View to readout // Input images @@ -52,7 +53,11 @@ pmReadout *kernelRO = pmFPAfileThisReadout(config->files, view, "PPSUB.OUTPUT.KERNELS"); // RO with kernel + psFree(view); + // Sources in image, used for stamps: these must be loaded from previous analysis stages - psArray *inSources = psMetadataLookupPtr(&mdok, inRO->analysis, "PSPHOT.SOURCES"); // Input source list - psArray *refSources = psMetadataLookupPtr(&mdok, refRO->analysis, "PSPHOT.SOURCES"); // Ref source list + pmReadout *inSourceRO = pmFPAfileThisReadout(config->files, view, "PPSUB.INPUT.SOURCES"); + pmReadout *refSourceRO = pmFPAfileThisReadout(config->files, view, "PPSUB.REF.SOURCES"); + psArray *inSources = psMetadataLookupPtr(&mdok, inSourceRO->analysis, "PSPHOT.SOURCES"); // Source list + psArray *refSources = psMetadataLookupPtr(&mdok, refSourceRO->analysis, "PSPHOT.SOURCES"); // Source list psArray *sources = NULL; // Merged list of sources Index: trunk/ppSub/src/ppSubReadout.c =================================================================== --- trunk/ppSub/src/ppSubReadout.c (revision 23688) +++ trunk/ppSub/src/ppSubReadout.c (revision 23740) @@ -21,32 +21,22 @@ #include "ppSub.h" -bool ppSubReadout(pmConfig *config, ppSubData *data, const pmFPAview *view) +bool ppSubReadout(const char *name, bool reverse, pmConfig *config, ppSubData *data, const pmFPAview *view) { - psTimerStart("PPSUB_MATCH"); + psAssert(name, "Require name"); + psAssert(config, "Require configuration"); + psAssert(data, "Require data"); + psAssert(view, "Require view"); - if (!ppSubSetMasks(config, view)) { - psError(PS_ERR_UNKNOWN, false, "Unable to set masks."); - return false; - } - - if (!ppSubMatchPSFs(config, data, view)) { - psError(PS_ERR_UNKNOWN, false, "Unable to match PSFs."); - return false; - } else if (data->quality) { - // Can't do anything at all - return true; - } - - if (!ppSubDefineOutput(config, view)) { + if (!ppSubDefineOutput(name, config, data, view)) { psError(PS_ERR_UNKNOWN, false, "Unable to define output."); return false; } - if (!data->quality && !ppSubMakePSF(config, data, view)) { + if (!data->quality && !ppSubMakePSF(name, config, data, view)) { psError(PS_ERR_UNKNOWN, false, "Unable to generate PSF."); return false; } - if (!ppSubReadoutSubtract(config, view)) { + if (!ppSubReadoutSubtract(name, reverse, config, view)) { psError(PS_ERR_UNKNOWN, false, "Unable to subtract images."); return false; @@ -54,19 +44,22 @@ // Higher order background subtraction using psphot - if (!ppSubBackground(config, view)) { + if (!ppSubBackground(name, config, view)) { psError(PS_ERR_UNKNOWN, false, "Unable to subtract background."); return false; } - if (!data->quality && !ppSubReadoutPhotometry(config, data, view)) { + if (!data->quality && data->!ppSubReadoutPhotometry(name, config, data, view)) { psError(PS_ERR_UNKNOWN, false, "Unable to perform photometry."); return false; } - if (!ppSubReadoutUpdate(config, data, view)) { + if (!ppSubReadoutUpdate(name, config, data, view)) { psError(PS_ERR_UNKNOWN, false, "Unable to update."); return false; } + + + return true; } Index: trunk/ppSub/src/ppSubReadoutInverse.c =================================================================== --- trunk/ppSub/src/ppSubReadoutInverse.c (revision 23740) +++ trunk/ppSub/src/ppSubReadoutInverse.c (revision 23740) @@ -0,0 +1,27 @@ +#ifdef HAVE_CONFIG_H +#include +#endif + +#include +#include + +#include "ppSub.h" + + +bool ppSubReadoutInverse(pmConfig *config) +{ + psAssert(config, "Require configuration"); + + pmFPAview *view = ppSubViewReadout(); // View to readout + pmReadout *outRO = pmFPAfileThisReadout(config->files, view, "PPSUB.OUTPUT"); + pmReadout *invRO = pmFPAfileThisReadout(config->files, view, "PPSUB.INVERSE"); + + invRO->image = (psImage*)psBinaryOp(invRO->image, outRO->image, "*", psScalarAlloc(-1.0, PS_TYPE_F32)); + invRO->mask = psMemIncrRefCounter(outRO->mask); + invRO->variance = psMemIncrRefCounter(outRO->variance); + invRO->covariance = psMemIncrRefCounter(outRO->covariance); + + invRO->data_exists = invRO->parent->data_exists = invRO->parent->parent->data_exists = true; + + return true; +} Index: trunk/ppSub/src/ppSubReadoutJpeg.c =================================================================== --- trunk/ppSub/src/ppSubReadoutJpeg.c (revision 23740) +++ trunk/ppSub/src/ppSubReadoutJpeg.c (revision 23740) @@ -0,0 +1,50 @@ +#ifdef HAVE_CONFIG_H +#include +#endif + +#include +#include +#include + +#include "ppSub.h" + +bool ppSubReadoutJpeg(pmConfig *config) +{ + psAssert(config, "Require configuration"); + + pmFPAview *view = ppSubViewReadout(); // View to readout + pmReadout *outRO = pmFPAfileThisReadout(config->files, view, "PPSUB.OUTPUT"); // Readout to jpeg + + // Look up recipe values + psMetadata *recipe = psMetadataLookupMetadata(NULL, config->recipes, PPSUB_RECIPE); // Recipe for ppSim + psAssert(recipe, "We checked this earlier, so it should be here."); + + // Generate binned JPEGs + psImageMaskType maskBad = pmConfigMaskGet("BLANK", config); // Bits to mask for bad pixels + + int bin1 = psMetadataLookupS32(NULL, recipe, "BIN1"); // First binning level + int bin2 = psMetadataLookupS32(NULL, recipe, "BIN2"); // Second binning level + + // Target cells + pmCell *cell1 = pmFPAfileThisCell(config->files, view, "PPSUB.OUTPUT.JPEG1"); // Rebinned cell once + pmCell *cell2 = pmFPAfileThisCell(config->files, view, "PPSUB.OUTPUT.JPEG2"); // Rebinned cell twice + psFree(view); + + pmReadout *ro1 = pmReadoutAlloc(cell1), *ro2 = pmReadoutAlloc(cell2); // Binned readouts + if (!pmReadoutRebin(ro1, outRO, maskBad, bin1, bin1)) { + psError(PS_ERR_UNKNOWN, false, "Unable to bin output (1st binning)"); + psFree(ro1); + psFree(ro2); + return false; + } + if (!pmReadoutRebin(ro2, ro1, 0, bin2, bin2)) { + psError(PS_ERR_UNKNOWN, false, "Unable to bin output (2nd binning)"); + psFree(ro1); + psFree(ro2); + return false; + } + psFree(ro1); + psFree(ro2); + + return true; +} Index: trunk/ppSub/src/ppSubReadoutPhotometry.c =================================================================== --- trunk/ppSub/src/ppSubReadoutPhotometry.c (revision 23688) +++ trunk/ppSub/src/ppSubReadoutPhotometry.c (revision 23740) @@ -22,8 +22,11 @@ #include "ppSub.h" -bool ppSubReadoutPhotometry (pmConfig *config, ppSubData *data, const pmFPAview *view) +bool ppSubReadoutPhotometry(const char *name, pmConfig *config, ppSubData *data) { psAssert(config, "Require configuration"); - psAssert(view, "Require view"); + + if (!data->photometry) { + return false; + } // Look up recipe values @@ -39,4 +42,5 @@ // The PSF (measured in ppSubMakePSF) is stored on the chip->analysis of PSPHOT.INPUT // In order to use an incoming PSF, it must be stored on the chip->analysis of PSPHOT.PSF.LOAD + pmFPAview *view = ppSubViewReadout(); // View to readout pmChip *psfInputChip = pmFPAfileThisChip(config->files, view, "PSPHOT.INPUT"); // Chip with PSF psAssert (psfInputChip, "should have been generated for ppSubMakePSF"); @@ -47,5 +51,5 @@ psErrorStackPrint(stderr, "No PSF available"); psWarning("No PSF available --- suspect bad data quality."); - ppSubDataQuality(config, data, psErrorCodeLast(), PPSUB_FILES_PHOT); + ppSubDataQuality(config, data, psErrorCodeLast(), PPSUB_FILES_PHOT_SUB | PPSUB_FILES_PHOT_INV); return true; } @@ -58,8 +62,6 @@ // around the pointers so PSPHOT.INPUT corresponds to the output image; previously, it was // equivalent to the minuend image. - pmFPAfile *outputFile = psMetadataLookupPtr(&mdok, config->files, "PPSUB.OUTPUT"); // Output file - pmReadout *outRO = pmFPAviewThisReadout(view, outputFile->fpa); // Readout with the sources + pmReadout *inRO = pmFPAfileThisReadout(config->files, view, name); // Readout with image and sources - // XXX possibly rename this to PPSUB.RESID? pmFPAfile *photFile = psMetadataLookupPtr(&mdok, config->files, "PSPHOT.INPUT"); // Photometry file pmReadout *photRO = pmFPAviewThisReadout(view, photFile->fpa); // Readout to photometer @@ -68,17 +70,22 @@ photRO = pmReadoutAlloc(cell); // Output readout: subtraction } - photRO->image = psImageCopy(photRO->image, outRO->image, PS_TYPE_F32); - if (outRO->variance) { - photRO->variance = psImageCopy(photRO->variance, outRO->variance, PS_TYPE_F32); + photRO->image = psImageCopy(photRO->image, inRO->image, PS_TYPE_F32); + if (inRO->variance) { + photRO->variance = psImageCopy(photRO->variance, inRO->variance, PS_TYPE_F32); } else { psFree(photRO->variance); photRO->variance = NULL; } - if (outRO->mask) { - photRO->mask = psImageCopy(photRO->mask, outRO->mask, PS_TYPE_IMAGE_MASK); + if (inRO->mask) { + photRO->mask = psImageCopy(photRO->mask, inRO->mask, PS_TYPE_IMAGE_MASK); } else { psFree(photRO->mask); photRO->mask = NULL; } + psMetadataAddPtr(photRO->parent->parent->analysis, PS_LIST_TAIL, "PSPHOT.PSF", + PS_META_REPLACE | PS_DATA_UNKNOWN, "Point-spread function", data->psf); + + psFree(photRO->analysis); + photRO->analysis = psMetadataAlloc(); if (!psphotReadoutMinimal(config, view)) { @@ -87,17 +94,23 @@ psErrorStackPrint(stderr, "Unable to perform photometry on image"); psWarning("Unable to perform photometry on image --- suspect bad data quality."); - ppSubDataQuality(config, data, psErrorCodeLast(), PPSUB_FILES_PHOT); + ppSubDataQuality(config, data, psErrorCodeLast(), PPSUB_FILES_PHOT_SUB | PPSUB_FILES_PHOT_INV); } - photRO->data_exists = true; - photRO->parent->data_exists = true; - photRO->parent->parent->data_exists = true; + if (!data->quality && !psMetadataCopySingle(inRO->analysis, photRO->analysis, "PSPHOT.SOURCES")) { + psError(PS_ERR_UNEXPECTED_NULL, false, "Unable to copy PSPHOT.SOURCES"); + return false; + } if (data->stats) { psArray *sources = psMetadataLookupPtr(NULL, photRO->analysis, "PSPHOT.SOURCES"); // Sources - psMetadataAddS32(data->stats, PS_LIST_TAIL, "NUM_SOURCES", 0, "Number of sources detected", - sources ? sources->n : 0); - psMetadataAddF32(data->stats, PS_LIST_TAIL, "TIME_PHOT", 0, "Time to do photometry", - psTimerClear("PPSUB_PHOT")); + bool mdok; + int numSources = psMetadataLookupS32(&mdok, data->stats, "NUM_SOURCES"); // Number of sources + numSources += sources ? sources->n : 0; + psMetadataAddS32(data->stats, PS_LIST_TAIL, "NUM_SOURCES", PS_META_REPLACE, + "Total number of sources detected", numSources); + float newTime = psTimerClear("PPSUB_PHOT"); // Time for photometry + float oldTime = psMetadataLookupF32(&mdok, data->stats, "TIME_PHOT"); // Previous time for photometry + psMetadataAddF32(data->stats, PS_LIST_TAIL, "TIME_PHOT", PS_META_REPLACE, "Time to do photometry", + isfinite(oldTime) ? oldTime + newTime : newTime); } Index: trunk/ppSub/src/ppSubReadoutStats.c =================================================================== --- trunk/ppSub/src/ppSubReadoutStats.c (revision 23740) +++ trunk/ppSub/src/ppSubReadoutStats.c (revision 23740) @@ -0,0 +1,48 @@ +#ifdef HAVE_CONFIG_H +#include +#endif + +#include +#include +#include + +#include "ppSub.h" + + +bool ppSubReadoutStats(pmConfig *config, ppSubData *data) +{ + psAssert(config, "Require configuration"); + psAssert(data, "Require data"); + + if (!data->statsFile) { + // Nothing to do + return true; + } + + pmFPAfile *output = psMetadataLookupPtr(NULL, config->files, "PPSUB.OUTPUT"); // Output file + if (!output) { + psError(PS_ERR_UNEXPECTED_NULL, true, "Unable to find file PPSUB.OUTPUT.\n"); + return false; + } + psImageMaskType maskValue = pmConfigMaskGet("MASK.VALUE", config); + pmFPAview *view = ppSubViewReadout(); // View to readout + ppStatsFPA(data->stats, output->fpa, view, maskValue, config); + + pmReadout *outRO = pmFPAviewThisReadout(view, output->fpa); // Readout of interest + psFree(view); + + // Statistics on the matching + psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_MODE); + psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_STAMPS); + psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_DEV_MEAN); + psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_DEV_RMS); + psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_NORM); + psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_BGDIFF); + psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_MX); + psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_MY); + psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_MXX); + psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_MXY); + psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_MYY); + + return true; +} Index: trunk/ppSub/src/ppSubReadoutSubtract.c =================================================================== --- trunk/ppSub/src/ppSubReadoutSubtract.c (revision 23688) +++ trunk/ppSub/src/ppSubReadoutSubtract.c (revision 23740) @@ -23,9 +23,7 @@ #define WCS_TOLERANCE 0.001 // Tolerance for WCS -bool ppSubReadoutSubtract(pmConfig *config, const pmFPAview *view) +bool ppSubReadoutSubtract(pmConfig *config) { psAssert(config, "Require configuration"); - psAssert(view, "Require view"); - // Look up recipe values @@ -35,4 +33,6 @@ bool reverse = psMetadataLookupBool(&mdok, config->arguments, "REVERSE"); // Reverse sense of subtraction? + + pmFPAview *view = ppSubViewReadout(); // View to readout // Subtraction is: minuend - subtrahend @@ -91,4 +91,5 @@ psError(PS_ERR_UNKNOWN, false, "Unable to copy concepts from input to output."); psFree(outRO); + psFree(view); return false; } @@ -99,4 +100,5 @@ pmHDU *outHDU = outFPA->hdu; // Output HDU pmChip *outChip = pmFPAfileThisChip(config->files, view, "PPSUB.OUTPUT"); // Output chip + psFree(view); if (!outHDU || !inHDU) { psError(PS_ERR_UNKNOWN, false, "Unable to find HDU at FPA level to copy astrometry."); Index: trunk/ppSub/src/ppSubReadoutUpdate.c =================================================================== --- trunk/ppSub/src/ppSubReadoutUpdate.c (revision 23688) +++ trunk/ppSub/src/ppSubReadoutUpdate.c (revision 23740) @@ -20,88 +20,2 @@ #include "ppSub.h" - -bool ppSubReadoutUpdate(pmConfig *config, ppSubData *data, const pmFPAview *view) -{ - psAssert(config, "Require configuration"); - psAssert(view, "Require view"); - - bool mdok = false; // Status of MD lookup - - // Look up recipe values - psMetadata *recipe = psMetadataLookupMetadata(&mdok, config->recipes, PPSUB_RECIPE); // Recipe for ppSim - psAssert(recipe, "We checked this earlier, so it should be here."); - - pmFPAfile *outFile = psMetadataLookupPtr(&mdok, config->files, "PPSUB.OUTPUT"); // Output file - pmReadout *outRO = pmFPAfileThisReadout(config->files, view, "PPSUB.OUTPUT"); // Output image - pmFPA *outFPA = outFile->fpa; // Output FPA - pmHDU *outHDU = outFPA->hdu; // Output HDU - - // Add additional data to the header - pmFPAfile *refFile = psMetadataLookupPtr(NULL, config->files, "PPSUB.REF"); // Reference file - pmFPAfile *inFile = psMetadataLookupPtr(NULL, config->files, "PPSUB.INPUT"); // Input file - psMetadataAddStr(outHDU->header, PS_LIST_TAIL, "PPSUB.REFERENCE", 0, - "Subtraction reference", refFile->filename); - psMetadataAddStr(outHDU->header, PS_LIST_TAIL, "PPSUB.INPUT", 0, - "Subtraction input", inFile->filename); - ppSubVersionHeader(outHDU->header); - - // Statistics on the matching - if (data->stats) { - psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_MODE); - psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_STAMPS); - psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_DEV_MEAN); - psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_DEV_RMS); - psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_NORM); - psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_BGDIFF); - psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_MX); - psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_MY); - psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_MXX); - psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_MXY); - psMetadataCopySingle(data->stats, outRO->analysis, PM_SUBTRACTION_ANALYSIS_MYY); - - psMetadataAddF32(data->stats, PS_LIST_TAIL, "TIME_MATCH", 0, "Time to match PSFs", - psTimerClear("PPSUB_MATCH")); - } - - // Generate binned JPEGs - { - psImageMaskType maskBad = pmConfigMaskGet("BLANK", config); // Bits to mask for bad pixels - - int bin1 = psMetadataLookupS32(NULL, recipe, "BIN1"); // First binning level - int bin2 = psMetadataLookupS32(NULL, recipe, "BIN2"); // Second binning level - - // Target cells - pmCell *cell1 = pmFPAfileThisCell(config->files, view, "PPSUB.OUTPUT.JPEG1"); // Rebinned cell once - pmCell *cell2 = pmFPAfileThisCell(config->files, view, "PPSUB.OUTPUT.JPEG2"); // Rebinned cell twice - - pmReadout *ro1 = pmReadoutAlloc(cell1), *ro2 = pmReadoutAlloc(cell2); // Binned readouts - if (!pmReadoutRebin(ro1, outRO, maskBad, bin1, bin1)) { - psError(PS_ERR_UNKNOWN, false, "Unable to bin output (1st binning)"); - psFree(ro1); - psFree(ro2); - return false; - } - if (!pmReadoutRebin(ro2, ro1, 0, bin2, bin2)) { - psError(PS_ERR_UNKNOWN, false, "Unable to bin output (2nd binning)"); - psFree(ro1); - psFree(ro2); - return false; - } - psFree(ro1); - psFree(ro2); - } - -#ifdef TESTING - // Significance image - { - psImage *sig = (psImage*)psBinaryOp(NULL, outRO->image, "*", outRO->image); - psBinaryOp(sig, sig, "/", outRO->variance); - psFits *fits = psFitsOpen("significance.fits", "w"); - psFitsWriteImage(fits, NULL, sig, 0, NULL); - psFitsClose(fits); - psFree(sig); - } -#endif - - return true; -} Index: trunk/ppSub/src/ppSubSetMasks.c =================================================================== --- trunk/ppSub/src/ppSubSetMasks.c (revision 23688) +++ trunk/ppSub/src/ppSubSetMasks.c (revision 23740) @@ -22,8 +22,7 @@ #include "ppSub.h" -bool ppSubSetMasks(pmConfig *config, const pmFPAview *view) +bool ppSubSetMasks(pmConfig *config) { psAssert(config, "Require configuration"); - psAssert(view, "Require view"); psImageMaskType maskValue, markValue; // Mask values @@ -50,5 +49,6 @@ psAssert(lowValue, "LOW or BAD must be non-zero"); - // input images + // Input images + pmFPAview *view = ppSubViewReadout(); // View to readout pmReadout *inRO = pmFPAfileThisReadout(config->files, view, "PPSUB.INPUT"); // Input readout pmReadout *refRO = pmFPAfileThisReadout(config->files, view, "PPSUB.REF"); // Reference readout