Index: /branches/eam_branches/ipp-20130419/pswarp/src/Makefile.am =================================================================== --- /branches/eam_branches/ipp-20130419/pswarp/src/Makefile.am (revision 35526) +++ /branches/eam_branches/ipp-20130419/pswarp/src/Makefile.am (revision 35527) @@ -51,5 +51,6 @@ pswarpTransformSources.c \ pswarpTransformTile.c \ - pswarpUpdateStatistics.c \ + pswarpUpdateStatistics.c \ + pswarpUpdateMetadata.c \ pswarpVersion.c \ pswarpFiles.c Index: /branches/eam_branches/ipp-20130419/pswarp/src/pswarp.h =================================================================== --- /branches/eam_branches/ipp-20130419/pswarp/src/pswarp.h (revision 35526) +++ /branches/eam_branches/ipp-20130419/pswarp/src/pswarp.h (revision 35527) @@ -201,2 +201,17 @@ bool pswarpMakePSF (pmConfig *config, pmFPAfile *output, psMetadata *stats); bool pswarpUpdateStatistics (pmFPA *output, psMetadata *stats, pmFPA *input, pmFPA *astrom, pmConfig *config); +bool pswarpUpdateMetadata (pmFPA *output, pmFPA *input, pmFPA *astrom, pmConfig *config); + +// XXX functions in pswarpParseCamera + +bool AddStringAsArray (psMetadata *md, char *string, char *name); + +pmFPAfile *pswarpDefineInputFile(pmConfig *config,// Configuration + pmFPAfile *bind, // File to which to bind, or NULL + char *filerule, // Name of file rule + char *argname, // Argument name + pmFPAfileType fileType // Type of file + ); + +bool pswarpParseSingleInput (pmConfig *config); +bool pswarpParseMultiInput (pmConfig *config, psMetadata *fileListMD); Index: /branches/eam_branches/ipp-20130419/pswarp/src/pswarpArguments.c =================================================================== --- /branches/eam_branches/ipp-20130419/pswarp/src/pswarpArguments.c (revision 35526) +++ /branches/eam_branches/ipp-20130419/pswarp/src/pswarpArguments.c (revision 35527) @@ -83,14 +83,50 @@ pswarpSetThreads(); - pmConfigFileSetsMD (config->arguments, &argc, argv, "INPUT", "-file", "-list"); - pmConfigFileSetsMD (config->arguments, &argc, argv, "MASK", "-mask", "-masklist"); - pmConfigFileSetsMD (config->arguments, &argc, argv, "VARIANCE", "-variance", "-variancelist"); - pmConfigFileSetsMD (config->arguments, &argc, argv, "BACKGROUND", "-background", "-bkglist"); + // there are three mutually exclusive ways of providing the input + // 1) supply -file (filename) [-mask .. -variance ..] on the command line + // 2) supply -input (input.mdc) on the command line + // 3) load inputs from RUN config info + + // below, we check first for -file then for -input. failure to find either implies use of + // the configuration metadata file. + + bool singleInput = pmConfigFileSetsMD (config->arguments, &argc, argv, "INPUT", "-file", "-list"); + if (singleInput) { + if (psArgumentGet(argc, argv, "-input")) { + psErrorStackPrint(stderr, "error in arguments : -input and -file / -list are mutually exclusive"); + exit(PS_EXIT_CONFIG_ERROR); + } + pmConfigFileSetsMD (config->arguments, &argc, argv, "MASK", "-mask", "-masklist"); + pmConfigFileSetsMD (config->arguments, &argc, argv, "VARIANCE", "-variance", "-variancelist"); + pmConfigFileSetsMD (config->arguments, &argc, argv, "BACKGROUND", "-background", "-bkglist"); + } else { + // find the input data file (an mdc file) + if ((N = psArgumentGet(argc, argv, "-input"))) { + if (argc <= N+1) { + psErrorStackPrint(stderr, "Expected to see 1 more argument; saw %d", argc - 1); + exit(PS_EXIT_CONFIG_ERROR); + } + psArgumentRemove(N, &argc, argv); + + unsigned int numBad = 0; // Number of bad lines + psMetadata *inputs = psMetadataConfigRead(NULL, &numBad, argv[N], false); // Input file info + if (!inputs || numBad > 0) { + psErrorStackPrint(stderr, "Unable to cleanly read MDC file with inputs."); + exit(PS_EXIT_CONFIG_ERROR); + } + psMetadataAddMetadata(config->arguments, PS_LIST_TAIL, "INPUTS", 0, "Metadata with input details", inputs); + psFree(inputs); + + psArgumentRemove(N, &argc, argv); + } + } + if (argc != 3) { + usage(); + } + if (psErrorCodeLast() != PS_ERR_NONE) { + psErrorStackPrint(stderr, "error in arguments"); + exit(PS_EXIT_CONFIG_ERROR); + } - if (argc != 3) { - usage(); - } - - psArray *array; Index: /branches/eam_branches/ipp-20130419/pswarp/src/pswarpDefineLayout.c =================================================================== --- /branches/eam_branches/ipp-20130419/pswarp/src/pswarpDefineLayout.c (revision 35526) +++ /branches/eam_branches/ipp-20130419/pswarp/src/pswarpDefineLayout.c (revision 35527) @@ -14,21 +14,7 @@ bool pswarpDefineLayout (pmConfig *config) { - // place input astrometry transformations in 'input' - pmFPAfile *input = psMetadataLookupPtr(NULL, config->files, "PSWARP.INPUT"); - if (!input) { - psError(PSWARP_ERR_CONFIG, false, "Can't find input data!\n"); - return false; - } - // input astrometry may be embedded in 'input' or supplied separately - pmFPAfile *astrom = psMetadataLookupPtr(NULL, config->files, "PSWARP.ASTROM"); - if (!astrom) { - astrom = input; - } - if (!pswarpLoadAstrometry (input, astrom, config)) { - psError(PSWARP_ERR_CONFIG, false, "problem loading input astrometry\n"); - return false; - } - - pmFPAfile *output = psMetadataLookupPtr(NULL, config->files, "PSWARP.OUTPUT"); + bool status = false; + + pmFPAfile *output = psMetadataLookupPtr(&status, config->files, "PSWARP.OUTPUT"); if (!output) { psError(PSWARP_ERR_CONFIG, false, "Can't find output data!\n"); @@ -36,5 +22,5 @@ } // select the input data sources - pmFPAfile *skycell = psMetadataLookupPtr (NULL, config->files, "PSWARP.SKYCELL"); + pmFPAfile *skycell = psMetadataLookupPtr (&status, config->files, "PSWARP.SKYCELL"); if (!skycell) { psError(PSWARP_ERR_CONFIG, false, "Can't find skycell data!\n"); @@ -46,16 +32,46 @@ } - // given the input data, determine which output elements should be generated - // and generate the appropriate images - - // find the R,D center for the skycell, use for common projection - psProjection *frame = pswarpLocalFrame (skycell->fpa); - - // generate Lmin,max, Mmin,max for both datasets - pswarpBounds *srcBounds = pswarpMakeBounds (input->fpa, frame); - pswarpBounds *tgtBounds = pswarpMakeBounds (skycell->fpa, frame); - - // find the output (tgt) chips which overlap the input (src) chips - pswarpFindOverlap (input->fpa, output->fpa, srcBounds, tgtBounds); + // chips are not processed unless we have determined they overlap the inputs + pmFPAExcludeChips (output->fpa); + + int nInputs = psMetadataLookupS32(&status, config->arguments, "NUM_INPUTS"); + if (!status) { + psError(PSWARP_ERR_DATA, true, "number of inputs is not defined (programming error)"); + return false; + } + + for (int i = 0; i < nInputs; i++) { + // place input astrometry transformations in 'input' + pmFPAfile *input = pmFPAfileSelectSingle(config->files, "PSWARP.INPUT", i); + if (!input) { + psError(PSWARP_ERR_CONFIG, false, "Can't find input data!\n"); + return false; + } + + // input astrometry may be embedded in 'input' or supplied separately + pmFPAfile *astrom = pmFPAfileSelectSingle(config->files, "PSWARP.ASTROM", i); + if (!astrom) { + astrom = input; + } + if (!pswarpLoadAstrometry (input, astrom, config)) { + psError(PSWARP_ERR_CONFIG, false, "problem loading input astrometry\n"); + return false; + } + + // given the input data, determine which output elements should be generated + // and generate the appropriate images + + // find the R,D center for the skycell, use for common projection + psProjection *frame = pswarpLocalFrame (skycell->fpa); + + // generate Lmin,max, Mmin,max for both datasets + pswarpBounds *srcBounds = pswarpMakeBounds (input->fpa, frame); + pswarpBounds *tgtBounds = pswarpMakeBounds (skycell->fpa, frame); + + // find the output (tgt) chips which overlap the input (src) chips + pswarpFindOverlap (input->fpa, output->fpa, srcBounds, tgtBounds); + psFree (srcBounds); + psFree (tgtBounds); + } // The loop below generates the output pixels. XXX Should this be deferred until we @@ -86,6 +102,6 @@ return false; } - int numCols = psMetadataLookupS32(NULL, hdu->header, "NAXIS1"); ///< Number of columns - int numRows = psMetadataLookupS32(NULL, hdu->header, "NAXIS2"); ///< Number of rows + int numCols = psMetadataLookupS32(&status, hdu->header, "NAXIS1"); ///< Number of columns + int numRows = psMetadataLookupS32(&status, hdu->header, "NAXIS2"); ///< Number of rows if ((numCols == 0) || (numRows == 0)) { psError(PSWARP_ERR_DATA, false, "skycell has invalid dimensions %d x %d", numCols, numRows); @@ -97,4 +113,5 @@ pmReadout *readout = pmReadoutAlloc(cell); ///< output readout readout->image = psImageAlloc(numCols, numRows, PS_TYPE_F32); + psImageInit(readout->image, NAN); psFree(readout); // Drop reference @@ -129,4 +146,6 @@ # include "pswarpFileNames.h" +static int Nout = 0; + // XX function for tests bool pswarpDumpOutput (pmConfig *config) { @@ -144,25 +163,28 @@ pmFPAfileActivate(config->files, true, "PSWARP.OUTPUT"); + pmFPAview *view = pmFPAviewAlloc(0); + pmChip *chip; - pmFPAview *view = pmFPAviewAlloc(0); - if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) { - psError(psErrorCodeLast(), false, "Unable to read files."); - goto DONE; - } while ((chip = pmFPAviewNextChip (view, output->fpa, 1)) != NULL) { psTrace ("pswarp", 4, "Chip %d: %x %x\n", view->chip, chip->file_exists, chip->process); if (!chip->process || !chip->file_exists) { continue; } - if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) { - psError(psErrorCodeLast(), false, "Unable to read files."); - goto DONE; - } + pmCell *cell; while ((cell = pmFPAviewNextCell (view, output->fpa, 1)) != NULL) { psTrace ("pswarp", 4, "Cell %d: %x %x\n", view->cell, cell->file_exists, cell->process); if (!cell->process || !cell->file_exists) { continue; } - if (!pmFPAfileIOChecks (config, view, PM_FPA_BEFORE) || - !pmFPAfileIOChecks (config, view, PM_FPA_AFTER)) { - psError(psErrorCodeLast(), false, "Unable to read files."); - goto DONE; + + // process each of the readouts + pmReadout *readout; + while ((readout = pmFPAviewNextReadout(view, output->fpa, 1)) != NULL) { + if (!readout->data_exists) { + continue; + } + + char name[64]; + snprintf (name, 64, "dumpwarp.%02d.fits", Nout); + Nout ++; + + psphotSaveImage (NULL, readout->image, name); } } Index: /branches/eam_branches/ipp-20130419/pswarp/src/pswarpLoop.c =================================================================== --- /branches/eam_branches/ipp-20130419/pswarp/src/pswarpLoop.c (revision 35526) +++ /branches/eam_branches/ipp-20130419/pswarp/src/pswarpLoop.c (revision 35527) @@ -36,22 +36,4 @@ // select the input data sources - pmFPAfile *input = psMetadataLookupPtr(&status, config->files, "PSWARP.INPUT"); - if (!input) { - psError(PSWARP_ERR_CONFIG, true, "Can't find input data!\n"); - return false; - } - - // use the external astrometry source if supplied - pmFPAfile *astrom = psMetadataLookupPtr(&status, config->files, "PSWARP.ASTROM"); - if (!astrom) { - astrom = input; - } - - if (astrom->camera != input->camera) { - psError(PSWARP_ERR_DATA, true, "Input camera and astrometry camera do not match."); - return false; - } - - // select the input data sources pmFPAfile *output = psMetadataLookupPtr(&status, config->files, "PSWARP.OUTPUT"); if (!output) { @@ -67,97 +49,133 @@ } - // ensure everyone is off except what we need below - pswarpFileActivation(config, detectorFiles, false); - pswarpFileActivation(config, photFiles, false); - pswarpFileActivation(config, independentFiles, false); - pswarpFileActivation(config, skycellFiles, false); + pmFPAview *view = pmFPAviewAlloc(0); + + int nInputs = psMetadataLookupS32(&status, config->arguments, "NUM_INPUTS"); + if (!status) { + psError(PSWARP_ERR_DATA, true, "number of inputs is not defined (programming error)"); + return false; + } + + // loop over this section once per input group + for (int i = 0; i < nInputs; i++) { + // select the input data sources + pmFPAfile *input = pmFPAfileSelectSingle(config->files, "PSWARP.INPUT", i); + if (!input) { + psError(PSWARP_ERR_CONFIG, true, "Can't find input data!\n"); + return false; + } + + fprintf (stderr, "loading file %s\n", input->filename); + + // use the external astrometry source if supplied + pmFPAfile *astrom = pmFPAfileSelectSingle(config->files, "PSWARP.ASTROM", i); + if (!astrom) { + astrom = input; + } + + if (astrom->camera != input->camera) { + psError(PSWARP_ERR_DATA, true, "Input camera and astrometry camera do not match."); + return false; + } + + // ensure everyone is off except what we need below + pswarpFileActivation(config, detectorFiles, false); + pswarpFileActivation(config, photFiles, false); + pswarpFileActivation(config, independentFiles, false); + pswarpFileActivation(config, skycellFiles, false); - // the loops below load the input data - pswarpFileActivation(config, detectorFiles, true); - - // need to read in these (but only if CMF, right?) - pmFPAfileActivate(config->files, true, "PSWARP.ASTROM"); - - // pswarpFileActivation(config, independentFiles, true); - - pmFPAview *view = pmFPAviewAlloc(0); - - // files associated with the science image - if (!pmFPAfileIOChecks (config, view, PM_FPA_BEFORE)) { - psError(psErrorCodeLast(), false, "Unable to read files."); - goto FAIL; - } - - // *** main transformation block - // *** this section loops over the input chips/cells and reads them one at a time - // *** the output chips/cells are filled where appropriate, but not yet written to disk - pmChip *chip; - while ((chip = pmFPAviewNextChip (view, input->fpa, 1)) != NULL) { - psTrace ("pswarp", 4, "Chip %d: %x %x\n", view->chip, chip->file_exists, chip->process); - if (!chip->process || !chip->file_exists) { continue; } - if (!pmFPAfileIOChecks (config, view, PM_FPA_BEFORE)) { - psError(psErrorCodeLast(), false, "Unable to read files."); - goto FAIL; - } - - pmCell *cell; - while ((cell = pmFPAviewNextCell (view, input->fpa, 1)) != NULL) { - psTrace ("pswarp", 4, "Cell %d: %x %x\n", view->cell, cell->file_exists, cell->process); - if (!cell->process || !cell->file_exists) { continue; } - if (!pmFPAfileIOChecks (config, view, PM_FPA_BEFORE)) { - psError(psErrorCodeLast(), false, "Unable to read files."); - goto FAIL; - } - - // process each of the readouts - pmReadout *readout; - while ((readout = pmFPAviewNextReadout(view, input->fpa, 1)) != NULL) { - if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) { - psError(psErrorCodeLast(), false, "Unable to read files."); - goto FAIL; - } - if (!readout->data_exists) { - continue; - } - - // Copy the detections from the astrometry carrier to the input, so they can be accessed by - // pswarpTransformReadout - if (astrom != input) { - pmReadout *astromRO = pmFPAviewThisReadout(view, astrom->fpa); // Readout for astrometry - pmDetections *detections = psMetadataLookupPtr(&status, astromRO->analysis, "PSPHOT.DETECTIONS"); // Sources from astrometry - if (detections) { - psMetadataAddPtr(readout->analysis, PS_LIST_TAIL, "PSPHOT.DETECTIONS", PS_DATA_ARRAY, "Sources from input astrometry", detections); - } + // the loops below load the input data + pswarpFileActivation(config, detectorFiles, true); + + // need to read in these (but only if CMF, right?) + pmFPAfileActivate(config->files, true, "PSWARP.ASTROM"); + + // pswarpFileActivation(config, independentFiles, true); + + pmFPAviewReset (view); + + // files associated with the science image + if (!pmFPAfileIOChecks (config, view, PM_FPA_BEFORE)) { + psError(psErrorCodeLast(), false, "Unable to read files."); + goto FAIL; + } + + // *** main transformation block + // *** this section loops over the input chips/cells and reads them one at a time + // *** the output chips/cells are filled where appropriate, but not yet written to disk + pmChip *chip; + while ((chip = pmFPAviewNextChip (view, input->fpa, 1)) != NULL) { + psTrace ("pswarp", 4, "Chip %d: %x %x\n", view->chip, chip->file_exists, chip->process); + if (!chip->process || !chip->file_exists) { continue; } + if (!pmFPAfileIOChecks (config, view, PM_FPA_BEFORE)) { + psError(psErrorCodeLast(), false, "Unable to read files."); + goto FAIL; + } + + pmCell *cell; + while ((cell = pmFPAviewNextCell (view, input->fpa, 1)) != NULL) { + psTrace ("pswarp", 4, "Cell %d: %x %x\n", view->cell, cell->file_exists, cell->process); + if (!cell->process || !cell->file_exists) { continue; } + if (!pmFPAfileIOChecks (config, view, PM_FPA_BEFORE)) { + psError(psErrorCodeLast(), false, "Unable to read files."); + goto FAIL; } - pswarpTransformToTarget(output->fpa, readout, config); + // process each of the readouts + pmReadout *readout; + while ((readout = pmFPAviewNextReadout(view, input->fpa, 1)) != NULL) { + if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) { + psError(psErrorCodeLast(), false, "Unable to read files."); + goto FAIL; + } + if (!readout->data_exists) { + continue; + } + + // Copy the detections from the astrometry carrier to the input, so they can be accessed by + // pswarpTransformReadout + if (astrom != input) { + pmReadout *astromRO = pmFPAviewThisReadout(view, astrom->fpa); // Readout for astrometry + pmDetections *detections = psMetadataLookupPtr(&status, astromRO->analysis, "PSPHOT.DETECTIONS"); // Sources from astrometry + if (detections) { + psMetadataAddPtr(readout->analysis, PS_LIST_TAIL, "PSPHOT.DETECTIONS", PS_DATA_ARRAY, "Sources from input astrometry", detections); + } + } + + // psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "SET.VALUE", PS_META_REPLACE, "test", (i + 1)*1.0); + pswarpTransformToTarget(output->fpa, readout, config); - if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { - psError(psErrorCodeLast(), false, "Unable to write files."); - goto FAIL; - } - } - if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { - psError(psErrorCodeLast(), false, "Unable to write files."); - goto FAIL; - } - } - if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { - psError(psErrorCodeLast(), false, "Unable to write files."); - goto FAIL; - } - } - if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { - psError(psErrorCodeLast(), false, "Unable to write files."); - goto FAIL; - } - - // Done with the detector side of things - pswarpFileActivation(config, detectorFiles, false); - pswarpFileActivation(config, independentFiles, false); - - if (!pswarpUpdateStatistics (output->fpa, stats, input->fpa, astrom->fpa, config)) { - psError(psErrorCodeLast(), false, "problem generating statistics."); - goto FAIL; + if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { + psError(psErrorCodeLast(), false, "Unable to write files."); + goto FAIL; + } + } + if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { + psError(psErrorCodeLast(), false, "Unable to write files."); + goto FAIL; + } + } + if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { + psError(psErrorCodeLast(), false, "Unable to write files."); + goto FAIL; + } + } + if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { + psError(psErrorCodeLast(), false, "Unable to write files."); + goto FAIL; + } + + // Done with the detector side of things + pswarpFileActivation(config, detectorFiles, false); + pswarpFileActivation(config, independentFiles, false); + + if (!pswarpUpdateStatistics (output->fpa, stats, input->fpa, astrom->fpa, config)) { + psError(psErrorCodeLast(), false, "problem generating statistics."); + goto FAIL; + } + if (!pswarpUpdateMetadata (output->fpa, input->fpa, astrom->fpa, config)) { + psError(psErrorCodeLast(), false, "problem generating statistics."); + goto FAIL; + } } @@ -166,4 +184,8 @@ goto FAIL; } + + // XXX + // pswarpDumpOutput (config); + // exit (0); psFree(view); @@ -174,4 +196,6 @@ return false; } + +// static int Nout = 0; // once the output fpa elements have been built, loop over the fpa and generate stats @@ -195,4 +219,9 @@ while ((readout = pmFPAviewNextReadout(view, output, 1)) != NULL) { pswarpTransformReadout (readout, input, config); + // char name[64]; + // snprintf (name, 64, "dumpwarp.%02d.fits", Nout); + // Nout ++; + // + // psphotSaveImage (NULL, readout->image, name); } } Index: /branches/eam_branches/ipp-20130419/pswarp/src/pswarpParseCamera.c =================================================================== --- /branches/eam_branches/ipp-20130419/pswarp/src/pswarpParseCamera.c (revision 35526) +++ /branches/eam_branches/ipp-20130419/pswarp/src/pswarpParseCamera.c (revision 35527) @@ -13,53 +13,151 @@ #include "pswarp.h" -// Define an input file -static pmFPAfile *defineInputFile(pmConfig *config,// Configuration - pmFPAfile *bind, // File to which to bind, or NULL - char *filerule, // Name of file rule - char *argname, // Argument name - pmFPAfileType fileType // Type of file - ) -{ - bool status; - - pmFPAfile *file = NULL; - // look for the file on the argument list - if (bind) { - file = pmFPAfileBindFromArgs(&status, bind, config, filerule, argname); - } else { - file = pmFPAfileDefineFromArgs(&status, config, filerule, argname); - } - if (!status) { - psError(psErrorCodeLast(), false, "Failed to load file definition for %s", filerule); - return false; - } - if (!file) { - // look for the file on the RUN metadata - file = pmFPAfileDefineFromRun(&status, bind, config, filerule); // File to return - if (!status) { - psError(psErrorCodeLast(), false, "Failed to load file definition for %s", filerule); - return NULL; - } - } - - if (!file) { - return NULL; - } - - if (file->type != fileType) { - psError(PSWARP_ERR_CONFIG, true, "%s is not of type %s", filerule, pmFPAfileStringFromType(fileType)); - return NULL; - } - - return file; -} +static bool foundAstrom = false; +static bool foundVariance = false; +static bool foundBackground = false; bool pswarpParseCamera(pmConfig *config) { psAssert(config, "Require configuration"); - bool mdok; // Status of MD lookup + bool status; // Status of MD lookup + + // *** parse the input information (either from -file or from -input) + + // if INPUTS exists, we have a metadata file to provide input, variance, mask, etc + psMetadata *inputFile = psMetadataLookupPtr(&status, config->arguments, "INPUTS"); + if (inputFile) { + if (!pswarpParseMultiInput (config, inputFile)) { + psError(PSWARP_ERR_CONFIG, false, "failed to parse multi-input file"); + return false; + } + } else { + if (!pswarpParseSingleInput (config)) { + psError(PSWARP_ERR_CONFIG, false, "failed to parse multi-input file"); + return false; + } + } + + // *** parse the output information (output target, output astrometry [skycell]) + + // The input skycell is a required argument: it defines the output image + pmConfig *skyConfig = NULL; + pmFPAfile *skycell = NULL; + + skycell = pmFPAfileDefineNewConfig(&status, &skyConfig, config, "PSWARP.SKYCELL", "SKYCELL"); + if (!status || !skycell) { + psError(psErrorCodeLast(), false, "Failed to build FPA from PSWARP.SKYCELL"); + return false; + } + psMetadataAddStr(config->arguments, PS_LIST_TAIL, "SKYCELL.CAMERA", 0, "Name of camera for skycell", skyConfig->cameraName); + psFree(skyConfig); + + psMetadata *recipe = psMetadataLookupPtr (&status, config->recipes, PSWARP_RECIPE); + if (!recipe) { + psError(PSWARP_ERR_CONFIG, false, "missing recipe %s", PSWARP_RECIPE); + return false; + } + + // The output skycell + pmFPAfile *output = pmFPAfileDefineOutputForFormat(config, NULL, "PSWARP.OUTPUT", skycell->cameraName, skycell->formatName); + if (!output) { + psError(psErrorCodeLast(), false, "Failed to build FPA from PSWARP.OUTPUT"); + return false; + } + output->save = true; + + pmFPAfile *outMask = pmFPAfileDefineOutputForFormat(config, output->fpa, "PSWARP.OUTPUT.MASK", skycell->cameraName, skycell->formatName); + if (!outMask) { + psError(psErrorCodeLast(), false, "Failed to build FPA from PSWARP.OUTPUT.MASK"); + return false; + } + outMask->save = true; + + // only create an output variance in we supply an input variance + if (foundVariance) { + pmFPAfile *outVariance = pmFPAfileDefineOutputForFormat(config, output->fpa, "PSWARP.OUTPUT.VARIANCE", skycell->cameraName, skycell->formatName); + if (!outVariance) { + psError(psErrorCodeLast(), false, "Failed to build FPA from PSWARP.OUTPUT.VARIANCE"); + return false; + } + outVariance->save = true; + } + + // XXX we assume input sources come from the input astrom description, but this need not be true (we could define an input sources file) + if (foundAstrom && psMetadataLookupBool(&status, recipe, "SOURCES")) { + pmFPAfile *outSources = pmFPAfileDefineOutputForFormat(config, output->fpa, "PSWARP.OUTPUT.SOURCES", skycell->cameraName, skycell->formatName); + if (!outSources) { + psError(psErrorCodeLast(), false, "Failed to build FPA from PSWARP.OUTPUT.SOURCES"); + return false; + } + outSources->save = true; + } + + // only create an output background if an input background is supplied + if (foundBackground && psMetadataLookupBool(&status, recipe, "BACKGROUND.MODEL")) { + pmFPAfile *outBackground = pmFPAfileDefineSkycell(config, NULL, "PSWARP.OUTPUT.BKGMODEL"); + outBackground->xBin = psMetadataLookupS32(&status, recipe, "BKG.XGRID"); + outBackground->yBin = psMetadataLookupS32(&status, recipe, "BKG.YGRID"); + + if (!outBackground) { + psError(psErrorCodeLast(), false, "Failed to build FPA from PSWARP.OUTPUT.BKGMODEL"); + return false; + } + outBackground->save = true; + } + + if (psMetadataLookupBool(&status, recipe, "PSF")) { + // This file, PSPHOT.INPUT, is just used as a carrier; output files (eg, PSPHOT.RESID) are defined by psphotDefineFiles + pmFPAfile *psphotInput = pmFPAfileDefineOutputForFormat(config, NULL, "PSPHOT.INPUT", skycell->cameraName, skycell->formatName); + if (!psphotInput) { + psError(psErrorCodeLast(), false, _("Unable to generate output file from PSPHOT.INPUT")); + return false; + } + psphotInput->src = psMemIncrRefCounter(output->fpa); + + // specify the number of psphot input images (psphotReadout loops over all input images) + psMetadataAddS32 (config->arguments, PS_LIST_TAIL, "PSPHOT.INPUT.NUM", PS_META_REPLACE, "number of inputs", 1); + + // the input sources (read from the input astrometry file) are transformed (in pswarpLoop) to the readout->analysis + // entries of the output file PSWARP.OUTPUT.SOURCES, associated with the main output pmFPAfile PSWARP.OUTPUT + + // the PSPHOT.INPUT.CMF is used to supply those sources to psphot (specifically psphotLoadPSFSources). + + // pmFPAfile *psphotInSources = pmFPAfileDefineSkycell(config, output->fpa, "PSPHOT.INPUT.CMF"); + pmFPAfile *psphotInSources = pmFPAfileDefineOutputForFormat(config, output->fpa, "PSPHOT.INPUT.CMF", skycell->cameraName, skycell->formatName); + if (!psphotInSources) { + psError(psErrorCodeLast(), false, _("Unable to generate output file from PSPHOT.INPUT.CMF")); + return false; + } + + // Ensure psphot defines an output file for the PSF. + psMetadata *psphotRecipe = psMetadataLookupPtr(NULL, config->recipes, PSPHOT_RECIPE); // Recipe + if (!psphotRecipe) { + psError(PSWARP_ERR_CONFIG, false, "Unable to find %s recipe.", PSPHOT_RECIPE); + return false; + } + psMetadataAddBool(psphotRecipe, PS_LIST_TAIL, "SAVE.PSF", PS_META_REPLACE, "Save PSF?", true); + + // Define associated psphot input/output files + if (!psphotDefineFiles(config, psphotInput)) { + psError(PSWARP_ERR_CONFIG, false, + "Unable to define the additional input/output files for psphot"); + return false; + } + + // Turn off writing the sources we find as part of deriving the PSF --- they might clobber the ones we + // get from transforming the input coordinates. + pmFPAfile *psphotSources = psMetadataLookupPtr(NULL, config->files, "PSPHOT.OUTPUT"); + assert(psphotSources); + psphotSources->save = false; + } + + psTrace("pswarp", 1, "Done with pswarpParseCamera...\n"); + return true; +} + +bool pswarpParseSingleInput (pmConfig *config) { // The input image(s) is required: it defines the camera - pmFPAfile *input = defineInputFile(config, NULL, "PSWARP.INPUT", "INPUT", PM_FPA_FILE_IMAGE); + pmFPAfile *input = pswarpDefineInputFile(config, NULL, "PSWARP.INPUT", "INPUT", PM_FPA_FILE_IMAGE); if (!input) { psError(psErrorCodeLast(), false, "Failed to build FPA from PSWARP.INPUT"); @@ -77,11 +175,16 @@ } } - if (astrom && (astrom->type != PM_FPA_FILE_CMF) && (astrom->type != PM_FPA_FILE_WCS)) { - psLogMsg("pswarp", PS_LOG_INFO, "%s is neither CMF nor WCS", "PSWARP.ASTROM"); + if (astrom) { + if ((astrom->type != PM_FPA_FILE_CMF) && (astrom->type != PM_FPA_FILE_WCS)) { + psLogMsg("pswarp", PS_LOG_INFO, "%s is neither CMF nor WCS", "PSWARP.ASTROM"); + // XXX raise an error here?? + } else { + foundAstrom = true; + } } psLogMsg("pswarp", PS_LOG_INFO, "Astrometry source: %s", astrom ? "supplied" : "header"); // Define the input mask file(s) - pmFPAfile *inMask = defineInputFile(config, input, "PSWARP.MASK", "MASK", PM_FPA_FILE_MASK); + pmFPAfile *inMask = pswarpDefineInputFile(config, input, "PSWARP.MASK", "MASK", PM_FPA_FILE_MASK); if (!inMask) { psLogMsg("pswarp", PS_LOG_INFO, "No mask supplied"); @@ -89,145 +192,21 @@ // Define the input variance file(s) - pmFPAfile *inVariance = defineInputFile(config, input, "PSWARP.VARIANCE", "VARIANCE", - PM_FPA_FILE_VARIANCE); - if (!inVariance) { + pmFPAfile *inVariance = pswarpDefineInputFile(config, input, "PSWARP.VARIANCE", "VARIANCE", PM_FPA_FILE_VARIANCE); + if (inVariance) { + foundVariance = true; + } else { psLogMsg("pswarp", PS_LOG_INFO, "No variance supplied"); } - pmFPAfile *inBackground = defineInputFile(config, NULL, "PSWARP.BKGMODEL", "BACKGROUND", - PM_FPA_FILE_IMAGE); - if (!inBackground) { + pmFPAfile *inBackground = pswarpDefineInputFile(config, NULL, "PSWARP.BKGMODEL", "BACKGROUND", PM_FPA_FILE_IMAGE); + if (inBackground) { + foundBackground = true; + } else { psLogMsg("pswarp", PS_LOG_INFO, "No background models supplied"); } - - // The input skycell is a required argument: it defines the output image - pmConfig *skyConfig = NULL; - pmFPAfile *skycell = NULL; - -# define USE_PSWARP_DEFINE_SKYCELL 0 -# if (USE_PSWARP_DEFINE_SKYCELL) - // XXX this function is no longer needed; drop - status = pswarpDefineSkycell(&skycell, &skyConfig, config, "PSWARP.SKYCELL", "SKYCELL"); -# else - skycell = pmFPAfileDefineNewConfig(&status, &skyConfig, config, "PSWARP.SKYCELL", "SKYCELL"); -# endif - - if (!status || !skycell) { - psError(psErrorCodeLast(), false, "Failed to build FPA from PSWARP.SKYCELL"); - return false; - } - psMetadataAddStr(config->arguments, PS_LIST_TAIL, "SKYCELL.CAMERA", 0, "Name of camera for skycell", skyConfig->cameraName); - psFree(skyConfig); - - psMetadata *recipe = psMetadataLookupPtr (&status, config->recipes, PSWARP_RECIPE); - if (!recipe) { - psError(PSWARP_ERR_CONFIG, false, "missing recipe %s", PSWARP_RECIPE); - return false; - } - - // The output skycell - // XXX The output needs to be in a more generic format, not just skycell. - // I need to understand how badly things depend on SKYCELL being a special format - // (is it just pswarp, or is it everything downstream as well??) - // XXX for a test, generate a generic output (based on output target astrometry) - // XXX pmFPAfile *output = pmFPAfileDefineSkycell(config, NULL, "PSWARP.OUTPUT"); - pmFPAfile *output = pmFPAfileDefineOutputForFormat(config, NULL, "PSWARP.OUTPUT", skycell->cameraName, skycell->formatName); - if (!output) { - psError(psErrorCodeLast(), false, "Failed to build FPA from PSWARP.OUTPUT"); - return false; - } - output->save = true; - - // pmFPAfile *outMask = pmFPAfileDefineSkycell(config, output->fpa, "PSWARP.OUTPUT.MASK"); - pmFPAfile *outMask = pmFPAfileDefineOutputForFormat(config, output->fpa, "PSWARP.OUTPUT.MASK", skycell->cameraName, skycell->formatName); - if (!outMask) { - psError(psErrorCodeLast(), false, "Failed to build FPA from PSWARP.OUTPUT.MASK"); - return false; - } - outMask->save = true; - - if (inVariance) { - // pmFPAfile *outVariance = pmFPAfileDefineSkycell(config, output->fpa, "PSWARP.OUTPUT.VARIANCE"); - pmFPAfile *outVariance = pmFPAfileDefineOutputForFormat(config, output->fpa, "PSWARP.OUTPUT.VARIANCE", skycell->cameraName, skycell->formatName); - if (!outVariance) { - psError(psErrorCodeLast(), false, "Failed to build FPA from PSWARP.OUTPUT.VARIANCE"); - return false; - } - outVariance->save = true; - } - - if (astrom && psMetadataLookupBool(&mdok, recipe, "SOURCES")) { - // pmFPAfile *outSources = pmFPAfileDefineSkycell(config, output->fpa, "PSWARP.OUTPUT.SOURCES"); - pmFPAfile *outSources = pmFPAfileDefineOutputForFormat(config, output->fpa, "PSWARP.OUTPUT.SOURCES", skycell->cameraName, skycell->formatName); - if (!outSources) { - psError(psErrorCodeLast(), false, "Failed to build FPA from PSWARP.OUTPUT.SOURCES"); - return false; - } - outSources->save = true; - } - - if (inBackground && psMetadataLookupBool(&mdok, recipe, "BACKGROUND.MODEL")) { - pmFPAfile *outBackground = pmFPAfileDefineSkycell(config, NULL, "PSWARP.OUTPUT.BKGMODEL"); - outBackground->xBin = psMetadataLookupS32(&mdok, recipe, "BKG.XGRID"); - outBackground->yBin = psMetadataLookupS32(&mdok, recipe, "BKG.YGRID"); - - if (!outBackground) { - psError(psErrorCodeLast(), false, "Failed to build FPA from PSWARP.OUTPUT.BKGMODEL"); - return false; - } - outBackground->save = true; - } - - if (psMetadataLookupBool(&mdok, recipe, "PSF")) { - // This file, PSPHOT.INPUT, is just used as a carrier; output files (eg, PSPHOT.RESID) are defined by - // psphotDefineFiles - // XXX old : pmFPAfile *psphotInput = pmFPAfileDefineSkycell(config, NULL, "PSPHOT.INPUT"); - pmFPAfile *psphotInput = pmFPAfileDefineOutputForFormat(config, NULL, "PSPHOT.INPUT", skycell->cameraName, skycell->formatName); - if (!psphotInput) { - psError(psErrorCodeLast(), false, _("Unable to generate output file from PSPHOT.INPUT")); - return false; - } - psphotInput->src = psMemIncrRefCounter(output->fpa); - - // specify the number of psphot input images (psphotReadout loops over all input images) - psMetadataAddS32 (config->arguments, PS_LIST_TAIL, "PSPHOT.INPUT.NUM", PS_META_REPLACE, "number of inputs", 1); - - // the input sources (read from the input astrometry file) are transformed (in pswarpLoop) to the readout->analysis - // entries of the output file PSWARP.OUTPUT.SOURCES, associated with the main output pmFPAfile PSWARP.OUTPUT - - // the PSPHOT.INPUT.CMF is used to supply those sources to psphot (specifically psphotLoadPSFSources). - - // pmFPAfile *psphotInSources = pmFPAfileDefineSkycell(config, output->fpa, "PSPHOT.INPUT.CMF"); - pmFPAfile *psphotInSources = pmFPAfileDefineOutputForFormat(config, output->fpa, "PSPHOT.INPUT.CMF", skycell->cameraName, skycell->formatName); - if (!psphotInSources) { - psError(psErrorCodeLast(), false, _("Unable to generate output file from PSPHOT.INPUT.CMF")); - return false; - } - - // Ensure psphot defines an output file for the PSF. - psMetadata *psphotRecipe = psMetadataLookupPtr(NULL, config->recipes, PSPHOT_RECIPE); // Recipe - if (!psphotRecipe) { - psError(PSWARP_ERR_CONFIG, false, "Unable to find %s recipe.", PSPHOT_RECIPE); - return false; - } - psMetadataAddBool(psphotRecipe, PS_LIST_TAIL, "SAVE.PSF", PS_META_REPLACE, "Save PSF?", true); - - // Define associated psphot input/output files - if (!psphotDefineFiles(config, psphotInput)) { - psError(PSWARP_ERR_CONFIG, false, - "Unable to define the additional input/output files for psphot"); - return false; - } - - // Turn off writing the sources we find as part of deriving the PSF --- they might clobber the ones we - // get from transforming the input coordinates. - pmFPAfile *psphotSources = psMetadataLookupPtr(NULL, config->files, "PSPHOT.OUTPUT"); - assert(psphotSources); - psphotSources->save = false; - } // Chip selection: turn on only the chips specified - char *chipLine = psMetadataLookupStr(&mdok, config->arguments, "CHIP_SELECTIONS"); - if (mdok) { + char *chipLine = psMetadataLookupStr(&status, config->arguments, "CHIP_SELECTIONS"); + if (status) { psArray *chips = psStringSplitArray (chipLine, ",", false); if (chips->n > 0) { @@ -244,5 +223,137 @@ } - psTrace("pswarp", 1, "Done with pswarpParseCamera...\n"); + psMetadataAddS32(config->arguments, PS_LIST_TAIL, "NUM_INPUTS", PS_META_REPLACE, "input file sets", 1); return true; } + +// read input file information from a metadata file +// XXX for now, in the multi-input format, we disable PSF, BACKGROUND, and CHIP_SELECTIONS +bool pswarpParseMultiInput (pmConfig *config, psMetadata *fileListMD) { + + // the multi-input file consists of a set of metadata blocks, each with entries + // for INPUT and (possibly) ASTROM, MASK, VARIANCE, ??? + + bool status = false; + + for (int i = 0; i < fileListMD->list->n; i++) { + psMetadataItem *item = psMetadataGet(fileListMD, i); + if (item->type != PS_DATA_METADATA) { + psError(PS_ERR_BAD_PARAMETER_TYPE, true, "Component %s of the input metadata is not of type METADATA", item->name); + return false; + } + + // next input metadata block of interest + psMetadata *fileBlockMD = item->data.md; + + psString inputName = psMetadataLookupStr(&status, fileBlockMD, "INPUT"); // Name of image + if (!inputName || !strlen(inputName)) { + psError(PS_ERR_UNKNOWN, false, "input file not found in metadata block %d (%s)", i, item->name); + return false; + } + AddStringAsArray (config->arguments, inputName, "FILENAMES"); + + pmFPAfile *input = pswarpDefineInputFile (config, NULL, "PSWARP.INPUT", "FILENAMES", PM_FPA_FILE_IMAGE); + if (!input) { + psError(psErrorCodeLast(), false, "Failed to build FPA from PSWARP.INPUT"); + return false; + } + + psString astromName = psMetadataLookupStr(&status, fileBlockMD, "ASTROM"); // Name of astrom file + if (astromName && strlen(astromName)) { + AddStringAsArray (config->arguments, astromName, "FILENAMES"); + + pmFPAfile *astrom = pmFPAfileDefineFromArgs (&status, config, "PSWARP.ASTROM", "FILENAMES"); + if (!status) { + psLogMsg("pswarp", PS_LOG_INFO, "Failed to load file definition for %s", "PSWARP.ASTROM"); + } + if (astrom) { + if ((astrom->type != PM_FPA_FILE_CMF) && (astrom->type != PM_FPA_FILE_WCS)) { + psLogMsg("pswarp", PS_LOG_INFO, "%s is neither CMF nor WCS", "PSWARP.ASTROM"); + } else { + foundAstrom = true; + } + } + psLogMsg("pswarp", PS_LOG_INFO, "Astrometry source: %s", astrom ? "supplied" : "header"); + } + + pmFPAfile *mask = NULL; + psString maskName = psMetadataLookupStr(&status, fileBlockMD, "MASK"); // Name of mask + if (maskName && strlen(maskName)) { + AddStringAsArray (config->arguments, maskName, "FILENAMES"); + + mask = pswarpDefineInputFile (config, input, "PSWARP.MASK", "FILENAMES", PM_FPA_FILE_MASK); + } + if (!mask) { + psLogMsg("pswarp", PS_LOG_INFO, "No mask supplied"); + } + + pmFPAfile *variance = NULL; + psString varianceName = psMetadataLookupStr(&status, fileBlockMD, "VARIANCE"); // Name of variance + if (varianceName && strlen(varianceName)) { + AddStringAsArray (config->arguments, varianceName, "FILENAMES"); + + variance = pswarpDefineInputFile (config, input, "PSWARP.VARIANCE", "FILENAMES", PM_FPA_FILE_VARIANCE); + } + if (variance) { + foundVariance = true; + } else { + psLogMsg("pswarp", PS_LOG_INFO, "No variance supplied"); + } + } + + psMetadataAddS32(config->arguments, PS_LIST_TAIL, "NUM_INPUTS", PS_META_REPLACE, "input file sets", fileListMD->list->n); + return true; +} + +// Define an input file based on name in config->arguments or config RUN block +pmFPAfile *pswarpDefineInputFile(pmConfig *config,// Configuration + pmFPAfile *bind, // File to which to bind, or NULL + char *filerule, // Name of file rule + char *argname, // Argument name + pmFPAfileType fileType // Type of file + ) +{ + bool status; + + pmFPAfile *file = NULL; + // look for the file on the argument list + if (bind) { + file = pmFPAfileBindFromArgs(&status, bind, config, filerule, argname); + } else { + file = pmFPAfileDefineFromArgs(&status, config, filerule, argname); + } + if (!status) { + psError(psErrorCodeLast(), false, "Failed to load file definition for %s", filerule); + return false; + } + if (!file) { + // look for the file on the RUN metadata + file = pmFPAfileDefineFromRun(&status, bind, config, filerule); // File to return + if (!status) { + psError(psErrorCodeLast(), false, "Failed to load file definition for %s", filerule); + return NULL; + } + } + + if (!file) { + return NULL; + } + + if (file->type != fileType) { + psError(PSWARP_ERR_CONFIG, true, "%s is not of type %s", filerule, pmFPAfileStringFromType(fileType)); + return NULL; + } + + return file; +} + +bool AddStringAsArray (psMetadata *md, char *string, char *name) { + + psArray *dummy = psArrayAlloc(1); // Dummy array of filenames for this FPA + dummy->data[0] = psStringCopy(string); + + psMetadataAddArray(md, PS_LIST_TAIL, name, PS_META_REPLACE, "string added as array", dummy); + psFree(dummy); + return true; +} + Index: /branches/eam_branches/ipp-20130419/pswarp/src/pswarpTransformTile.c =================================================================== --- /branches/eam_branches/ipp-20130419/pswarp/src/pswarpTransformTile.c (revision 35526) +++ /branches/eam_branches/ipp-20130419/pswarp/src/pswarpTransformTile.c (revision 35527) @@ -89,4 +89,7 @@ for (int y = yMin; y < yMax; y++) { + + int yOut = y - outRow0; ///< Position on output image + for (int x = xMin; x < xMax; x++) { // Only transform those pixels requested @@ -94,4 +97,11 @@ continue; } + + int xOut = x - outCol0; + + // XXX the existing image may already have valid data -- probably should keep + // the best, but for the moment, just keep the pixel which is not NAN + if (isfinite(outImageData[yOut][xOut])) continue; + // pswarpMapApply converts the output coordinate (x,y) to the input coordinate. // both are in the parent frames of the input and output images. @@ -127,8 +137,7 @@ } - int xOut = x - outCol0, yOut = y - outRow0; ///< Position on output image - if (outImageData) { - outImageData[yOut][xOut] = imageValue * jacobian; + // XXX TEST outImageData[yOut][xOut] = value; + outImageData[yOut][xOut] = imageValue * jacobian; } if (outVarData) { Index: /branches/eam_branches/ipp-20130419/pswarp/src/pswarpUpdateMetadata.c =================================================================== --- /branches/eam_branches/ipp-20130419/pswarp/src/pswarpUpdateMetadata.c (revision 35527) +++ /branches/eam_branches/ipp-20130419/pswarp/src/pswarpUpdateMetadata.c (revision 35527) @@ -0,0 +1,113 @@ +/** @file pswarpUpdateMetadata.c + * + * @brief update generic metadata info for warp + * @ingroup pswarp + * + * @author IfA + * @date $Date: 2009-02-13 21:54:32 $ + * Copyright 2009 Institute for Astronomy, University of Hawaii + */ + +#include "pswarp.h" + +// once the output fpa elements have been built, loop over the fpa and generate stats +// for each readout +bool pswarpUpdateMetadata (pmFPA *output, pmFPA *input, pmFPA *astrom, pmConfig *config) { + + pmFPAview *view = pmFPAviewAlloc(0); + + pmChip *chip; + while ((chip = pmFPAviewNextChip (view, output, 1)) != NULL) { + psTrace ("pswarp", 4, "Chip %d: %x %x\n", view->chip, chip->file_exists, chip->process); + if (!chip->process || !chip->file_exists) { continue; } + + // use this in output metadata info (MD5 sums) + const char *chipName = psMetadataLookupStr(NULL, chip->concepts, "CHIP.NAME"); + + pmCell *cell; + while ((cell = pmFPAviewNextCell (view, output, 1)) != NULL) { + psTrace ("pswarp", 4, "Cell %d: %x %x\n", view->cell, cell->file_exists, cell->process); + if (!cell->process || !cell->file_exists) { continue; } + + // use this in output metadata info (MD5 sums) + const char *cellName = psMetadataLookupStr(NULL, cell->concepts, "CELL.NAME"); + + // process each of the readouts + pmReadout *readout; + while ((readout = pmFPAviewNextReadout(view, output, 1)) != NULL) { + // skip empty output readouts + if (!readout->data_exists) continue; + + // Set covariance matrix for output +# if (0) + psList *covariances = psMetadataLookupPtr(&mdok, output->analysis, PSWARP_ANALYSIS_COVARIANCES); // Covariance matrices + psAssert(covariances, "Should be there"); + psArray *covars = psListToArray(covariances); // Array of covariance matrices + psKernel *covar = psImageCovarianceAverage(covars); + psFree(covars); + psMetadataRemoveKey(output->analysis, PSWARP_ANALYSIS_COVARIANCES); + + // Correct covariance matrix scale for the mean (square root of the) Jacobian + double jacobian = psMetadataLookupF64(NULL, output->analysis, PSWARP_ANALYSIS_JACOBIAN); // Jacobian + int goodPixels = psMetadataLookupS32(NULL, output->analysis, PSWARP_ANALYSIS_GOODPIX); // Good pixels + jacobian /= goodPixels; + output->covariance = psImageCovarianceScale(covar, jacobian); + psFree(covar); + + if (output->variance) { + psImageCovarianceTransfer(output->variance, output->covariance); + } +# endif + + { + // Add MD5 information for readout + psString headerName = NULL; ///< Header name for MD5 + psVector *md5 = psImageMD5(readout->image); ///< md5 hash + psString md5string = psMD5toString(md5); ///< String + + pmHDU *hdu = pmHDUFromReadout(readout); + psStringAppend(&headerName, "MD5_%s_%s_%d", chipName, cellName, view->readout); + psMetadataAddStr(hdu->header, PS_LIST_TAIL, headerName, PS_META_REPLACE, "Image MD5", md5string); + psFree(md5); + psFree(md5string); + psFree(headerName); + } + } + + psList *inputCells = psMetadataLookupPtr (NULL, cell->analysis, "INPUT.CELLS"); + if (!pmConceptsAverageCells(cell, inputCells, NULL, NULL, false)) { + psError(psErrorCodeLast(), false, "Unable to average cell concepts."); + psFree(view); + return false; + } + + // XXX Is this too ad-hoc? + psRegion *trimsec = psMetadataLookupPtr(NULL, cell->concepts, "CELL.TRIMSEC"); ///< Trim section + trimsec->x0 = trimsec->x1 = trimsec->y0 = trimsec->y1 = 0; ///< All pixels + + } + } + + if (!psMetadataCopy(output->concepts, input->concepts)) { + psError(psErrorCodeLast(), false, "Unable to copy FPA concepts from input to output."); + return false; + } + + // Update ZP from the astrometry + { + psMetadataItem *item = psMetadataLookup(output->concepts, "FPA.ZP"); + item->data.F32 = psMetadataLookupF32(NULL, astrom->concepts, "FPA.ZP"); + } + + // need to update the chip and fpa level astrometry appropriately, see + // code in psastro + +# if (0) + if (!pmAstromWriteWCS(hdu->header, outFPA, outChip, WCS_NONLIN_TOL)) { + psError(psErrorCodeLast(), false, "Unable to generate WCS header."); + return false; + } +# endif + + return true; +} Index: /branches/eam_branches/ipp-20130419/pswarp/src/pswarpUpdateStatistics.c =================================================================== --- /branches/eam_branches/ipp-20130419/pswarp/src/pswarpUpdateStatistics.c (revision 35526) +++ /branches/eam_branches/ipp-20130419/pswarp/src/pswarpUpdateStatistics.c (revision 35527) @@ -15,4 +15,9 @@ // for each readout bool pswarpUpdateStatistics (pmFPA *output, psMetadata *fpaStats, pmFPA *input, pmFPA *astrom, pmConfig *config) { + + if (!fpaStats) { + psLogMsg("pswarp", PS_LOG_INFO, "stats not requested, skipping"); + return true; + } // load the recipe @@ -95,69 +100,15 @@ } - // Set covariance matrix for output - // XXX is this the right place for the covariance calculation? -# if (0) - psList *covariances = psMetadataLookupPtr(&mdok, output->analysis, PSWARP_ANALYSIS_COVARIANCES); // Covariance matrices - psAssert(covariances, "Should be there"); - psArray *covars = psListToArray(covariances); // Array of covariance matrices - psKernel *covar = psImageCovarianceAverage(covars); - psFree(covars); - psMetadataRemoveKey(output->analysis, PSWARP_ANALYSIS_COVARIANCES); - - // Correct covariance matrix scale for the mean (square root of the) Jacobian - double jacobian = psMetadataLookupF64(NULL, output->analysis, PSWARP_ANALYSIS_JACOBIAN); // Jacobian - int goodPixels = psMetadataLookupS32(NULL, output->analysis, PSWARP_ANALYSIS_GOODPIX); // Good pixels - jacobian /= goodPixels; - output->covariance = psImageCovarianceScale(covar, jacobian); - psFree(covar); - - if (output->variance) { - psImageCovarianceTransfer(output->variance, output->covariance); - } -# endif - - { - // Add MD5 information for readout - psString headerName = NULL; ///< Header name for MD5 - psVector *md5 = psImageMD5(readout->image); ///< md5 hash - psString md5string = psMD5toString(md5); ///< String - - pmHDU *hdu = pmHDUFromReadout(readout); - psStringAppend(&headerName, "MD5_%s_%s_%d", chipName, cellName, view->readout); - psMetadataAddStr(hdu->header, PS_LIST_TAIL, headerName, PS_META_REPLACE, "Image MD5", md5string); - psFree(md5); - psFree(md5string); - psFree(headerName); - } psFree (readoutName); } - - psList *inputCells = psMetadataLookupPtr (NULL, cell->analysis, "INPUT.CELLS"); - if (!pmConceptsAverageCells(cell, inputCells, NULL, NULL, false)) { - psError(psErrorCodeLast(), false, "Unable to average cell concepts."); - psFree(view); - return false; - } - - // XXX Is this too ad-hoc? - psRegion *trimsec = psMetadataLookupPtr(NULL, cell->concepts, "CELL.TRIMSEC"); ///< Trim section - trimsec->x0 = trimsec->x1 = trimsec->y0 = trimsec->y1 = 0; ///< All pixels - } } // Perform statistics on the output image (ppStatsFPA loops down the hierarchy) - if (fpaStats) { - pmFPAviewReset (view); - if (!ppStatsFPA(fpaStats, output, view, maskValue, config)) { + pmFPAviewReset (view); + if (!ppStatsFPA(fpaStats, output, view, maskValue, config)) { psWarning("Unable to perform statistics on warped image."); - } } psFree(view); - - if (!psMetadataCopy(output->concepts, input->concepts)) { - psError(psErrorCodeLast(), false, "Unable to copy FPA concepts from input to output."); - return false; - } if (Nreadout == 0) { @@ -168,20 +119,4 @@ } - // Update ZP from the astrometry - { - psMetadataItem *item = psMetadataLookup(output->concepts, "FPA.ZP"); - item->data.F32 = psMetadataLookupF32(NULL, astrom->concepts, "FPA.ZP"); - } - - // need to update the chip and fpa level astrometry appropriately, see - // code in psastro - -# if (0) - if (!pmAstromWriteWCS(hdu->header, outFPA, outChip, WCS_NONLIN_TOL)) { - psError(psErrorCodeLast(), false, "Unable to generate WCS header."); - return false; - } -# endif - return true; }