Index: /trunk/pswarp/src/pswarp.h =================================================================== --- /trunk/pswarp/src/pswarp.h (revision 21367) +++ /trunk/pswarp/src/pswarp.h (revision 21368) @@ -6,6 +6,6 @@ * * @author IfA - * @version $Revision: 1.17 $ $Name: not supported by cvs2svn $ - * @date $Date: 2009-02-05 20:44:04 $ + * @version $Revision: 1.18 $ $Name: not supported by cvs2svn $ + * @date $Date: 2009-02-06 03:10:36 $ * Copyright 2009 Institute for Astronomy, University of Hawaii */ @@ -31,4 +31,5 @@ #define PSWARP_ANALYSIS_VARFACTOR "PSWARP.VARFACTOR" ///< Name for variance factor in analysis metadata #define PSWARP_ANALYSIS_GOODPIX "PSWARP.GOODPIX" ///< Name for number of good pixels in analysis metadata +#define PSWARP_ANALYSIS_COVARIANCES "PSWARP.COVARIANCES" ///< Name for covariance matrices on analysis MD /** @@ -43,5 +44,5 @@ } pswarpMap; -/* the pswarpMapGrid carries a collection of pswarpMag structures representing +/* the pswarpMapGrid carries a collection of pswarpMag structures representing * the local value of the pswarpMap at different locations in the image. */ @@ -66,6 +67,7 @@ /** output values for this tile */ - long goodPixels; //< Number of good pixels - int xMin, xMax, yMin, yMax; //< Bounds of tile + long goodPixels; ///< Number of good pixels + int xMin, xMax, yMin, yMax; ///< Bounds of tile + psKernel *covariance; ///< Covariance matrix } pswarpTransformTileArgs; Index: /trunk/pswarp/src/pswarpArguments.c =================================================================== --- /trunk/pswarp/src/pswarpArguments.c (revision 21367) +++ /trunk/pswarp/src/pswarpArguments.c (revision 21368) @@ -6,6 +6,6 @@ * * @author IfA - * @version $Revision: 1.23 $ $Name: not supported by cvs2svn $ - * @date $Date: 2009-02-05 20:44:04 $ + * @version $Revision: 1.24 $ $Name: not supported by cvs2svn $ + * @date $Date: 2009-02-06 03:10:36 $ * Copyright 2009 Institute for Astronomy, University of Hawaii */ @@ -92,7 +92,7 @@ } - // the mask and weight entries are optional (build from gain?) + // the mask and variance entries are optional (build from gain?) pmConfigFileSetsMD (config->arguments, &argc, argv, "MASK", "-mask", "-masklist"); - pmConfigFileSetsMD (config->arguments, &argc, argv, "WEIGHT", "-weight", "-weightlist"); + pmConfigFileSetsMD (config->arguments, &argc, argv, "VARIANCE", "-variance", "-variancelist"); if (argc != 3) { Index: /trunk/pswarp/src/pswarpFileNames.h =================================================================== --- /trunk/pswarp/src/pswarpFileNames.h (revision 21367) +++ /trunk/pswarp/src/pswarpFileNames.h (revision 21368) @@ -8,22 +8,22 @@ // Lists of file rules for the detectors -static char *detectorFiles[] = { - "PSWARP.INPUT", - "PSWARP.MASK", - "PSWARP.WEIGHT", - NULL +static char *detectorFiles[] = { + "PSWARP.INPUT", + "PSWARP.MASK", + "PSWARP.VARIANCE", + NULL }; // Lists of file rules for the skycells -static char *skycellFiles[] = { - "PSWARP.OUTPUT", - "PSWARP.OUTPUT.MASK", - "PSWARP.OUTPUT.WEIGHT", - NULL +static char *skycellFiles[] = { + "PSWARP.OUTPUT", + "PSWARP.OUTPUT.MASK", + "PSWARP.OUTPUT.VARIANCE", + NULL }; // Lists of file rules for photometry -static char *photFiles[] = { - "PSPHOT.INPUT", +static char *photFiles[] = { + "PSPHOT.INPUT", "PSPHOT.OUTPUT", #if PSPHOT_FIND_PSF @@ -31,21 +31,21 @@ #endif "PSPHOT.RESID", - "PSPHOT.BACKMDL", - "PSPHOT.BACKMDL.STDEV", - "PSPHOT.BACKGND", + "PSPHOT.BACKMDL", + "PSPHOT.BACKMDL.STDEV", + "PSPHOT.BACKGND", "PSPHOT.BACKSUB", - "PSPHOT.PSF.SAVE", - "SOURCE.PLOT.MOMENTS", + "PSPHOT.PSF.SAVE", + "SOURCE.PLOT.MOMENTS", "SOURCE.PLOT.PSFMODEL", - "SOURCE.PLOT.APRESID", - NULL + "SOURCE.PLOT.APRESID", + NULL }; // Lists of file rules for the detectors and skycells, and an independent list of everything else -static char *independentFiles[] = { +static char *independentFiles[] = { "PSWARP.ASTROM", // Read independently from the input pixels "PSWARP.SKYCELL", // Don't care about the skycell once we have its WCS "PSWARP.OUTPUT.SOURCES", // Save these independently so we can do the PSF - NULL + NULL }; Index: /trunk/pswarp/src/pswarpLoop.c =================================================================== --- /trunk/pswarp/src/pswarpLoop.c (revision 21367) +++ /trunk/pswarp/src/pswarpLoop.c (revision 21368) @@ -6,6 +6,6 @@ * * @author IfA - * @version $Revision: 1.36 $ $Name: not supported by cvs2svn $ - * @date $Date: 2009-02-05 20:44:04 $ + * @version $Revision: 1.37 $ $Name: not supported by cvs2svn $ + * @date $Date: 2009-02-06 03:10:36 $ * Copyright 2009 Institute for Astronomy, University of Hawaii */ @@ -13,18 +13,16 @@ #include "pswarp.h" #include - -#define WCS_NONLIN_TOL 0.001 ///< Non-linear tolerance for header WCS -#define PSPHOT_FIND_PSF 1 ///< Use psphot's findPSF function? -#define TESTING 0 ///< Testing output? - -#include "pswarpFileNames.h" ///< Lists of file rules used at different stages - -/** - * XXX these are generic functions which should be moved to psModules - * Activate a list of files - */ -static void fileActivation(pmConfig *config, ///< Configuration - char **files, ///< Files to turn on/off - bool state ///< Activation state +#include "pswarpFileNames.h" // Lists of file rules used at different stages + +#define WCS_NONLIN_TOL 0.001 // Non-linear tolerance for header WCS +#define PSPHOT_FIND_PSF 1 // Use psphot's findPSF function? +#define TESTING 0 // Testing output? + + +// XXX these are generic functions which should be moved to psModules +// Activate a list of files +static void fileActivation(pmConfig *config, // Configuration + char **files, // Files to turn on/off + bool state // Activation state ) { @@ -293,4 +291,15 @@ } + // Set covariance matrix for output + { + psList *covariances = psMetadataLookupPtr(&mdok, output->analysis, + PSWARP_ANALYSIS_COVARIANCES); // Covariance matrices + psAssert(covariances, "Should be there"); + psArray *covars = psListToArray(covariances); // Array of covariance matrices + output->covariance = psImageCovarianceAverage(covars); + psFree(covars); + psMetadataRemoveKey(output->analysis, PSWARP_ANALYSIS_COVARIANCES); + } + if (!pmConceptsAverageCells(outCell, cells, NULL, NULL, false)) { psError(PS_ERR_UNKNOWN, false, "Unable to average cell concepts."); @@ -348,5 +357,5 @@ ioChecksBefore(config); - // supply the readout and fpa of interest to psphot + // supply the readout and fpa of interest to psphot pmFPAfile *photFile = psMetadataLookupPtr(NULL, config->files, "PSPHOT.INPUT"); pmFPACopy(photFile->fpa, outFPA); @@ -359,15 +368,15 @@ // We can adjust the weight directly since this is a deep copy pmReadout *photRO = pmFPAviewThisReadout(view, photFile->fpa); - float vf = psMetadataLookupF32(NULL, photRO->parent->concepts, "CELL.VARFACTOR"); ///< Variance factor - psBinaryOp(photRO->weight, photRO->weight, "*", psScalarAlloc(vf, PS_TYPE_F32)); - - // grab the sources of interest from the storage location (pmFPAfile PSPHOT.INPUT.CMF) - psArray *sources = psphotLoadPSFSources (config, view); - if (!sources) { + float vf = psMetadataLookupF32(NULL, photRO->parent->concepts, "CELL.VARFACTOR"); // Variance factor + psBinaryOp(photRO->variance, photRO->variance, "*", psScalarAlloc(vf, PS_TYPE_F32)); + + // grab the sources of interest from the storage location (pmFPAfile PSPHOT.INPUT.CMF) + psArray *sources = psphotLoadPSFSources (config, view); + if (!sources) { psError(PS_ERR_UNKNOWN, false, "No sources supplied to measure PSF"); return false; } - // measure the PSF using these sources + // measure the PSF using these sources if (!psphotReadoutFindPSF(config, view, sources)) { psError(PS_ERR_UNKNOWN, false, "Unable to determine PSF for warped image."); Index: /trunk/pswarp/src/pswarpParseCamera.c =================================================================== --- /trunk/pswarp/src/pswarpParseCamera.c (revision 21367) +++ /trunk/pswarp/src/pswarpParseCamera.c (revision 21368) @@ -6,6 +6,6 @@ * * @author IfA - * @version $Revision: 1.23 $ $Name: not supported by cvs2svn $ - * @date $Date: 2009-02-05 20:44:04 $ + * @version $Revision: 1.24 $ $Name: not supported by cvs2svn $ + * @date $Date: 2009-02-06 03:10:36 $ * Copyright 2009 Institute for Astronomy, University of Hawaii */ @@ -57,11 +57,11 @@ } - pmFPAfile *inWeight = pmFPAfileBindFromArgs(&status, input, config, "PSWARP.WEIGHT", "WEIGHT"); + pmFPAfile *inVariance = pmFPAfileBindFromArgs(&status, input, config, "PSWARP.VARIANCE", "VARIANCE"); if (!status) { psError (PS_ERR_UNKNOWN, false, "failed to load find definition"); return NULL; } - if (!inWeight) { - psLogMsg ("pswarp", 3, "no weight supplied\n"); + if (!inVariance) { + psLogMsg ("pswarp", 3, "no variance supplied\n"); } @@ -90,11 +90,11 @@ outMask->save = true; - if (inWeight) { - pmFPAfile *outWeight = pmFPAfileDefineSkycell(config, output->fpa, "PSWARP.OUTPUT.WEIGHT"); - if (!outWeight) { - psError(PSWARP_ERR_CONFIG, false, "Failed to build FPA from PSWARP.OUTPUT.WEIGHT"); + if (inVariance) { + pmFPAfile *outVariance = pmFPAfileDefineSkycell(config, output->fpa, "PSWARP.OUTPUT.VARIANCE"); + if (!outVariance) { + psError(PSWARP_ERR_CONFIG, false, "Failed to build FPA from PSWARP.OUTPUT.VARIANCE"); return false; } - outWeight->save = true; + outVariance->save = true; } Index: /trunk/pswarp/src/pswarpSetMaskBits.c =================================================================== --- /trunk/pswarp/src/pswarpSetMaskBits.c (revision 21367) +++ /trunk/pswarp/src/pswarpSetMaskBits.c (revision 21368) @@ -6,6 +6,6 @@ * * @author IfA - * @version $Revision: 1.4 $ $Name: not supported by cvs2svn $ - * @date $Date: 2009-02-05 20:44:04 $ + * @version $Revision: 1.5 $ $Name: not supported by cvs2svn $ + * @date $Date: 2009-02-06 03:10:36 $ * Copyright 2009 Institute for Astronomy, University of Hawaii */ @@ -20,29 +20,29 @@ * updated in the config metadata. */ -bool pswarpSetMaskBits (pmConfig *config) { - - psImageMaskType maskIn = 0x00; ///< mask for the input image - psImageMaskType markIn = 0x00; ///< mark for the input image - psImageMaskType maskOut = 0x00; ///< mask for the output image +bool pswarpSetMaskBits (pmConfig *config) +{ + psImageMaskType maskIn = 0x00; // mask for the input image + psImageMaskType markIn = 0x00; // mark for the input image + psImageMaskType maskOut = 0x00; // mask for the output image // this function sets the required single-image mask bits if (!pmConfigMaskSetBits (&maskIn, &markIn, config)) { - psError (PS_ERR_UNKNOWN, true, "Unable to define the mask bit values"); - return false; + psError (PS_ERR_UNKNOWN, true, "Unable to define the mask bit values"); + return false; } // mask for non-linear flat regions (default to DETECTOR if not defined) - psImageMaskType badMask = pmConfigMaskGet("BAD.WARP", config); + psImageMaskType badMask = pmConfigMaskGet("BAD.WARP", config); if (!badMask) { - badMask = 0x01; - pmConfigMaskSet (config, "BAD.WARP", badMask); + badMask = 0x01; + pmConfigMaskSet (config, "BAD.WARP", badMask); } maskOut |= badMask; // mask for non-linear flat regions (default to DETECTOR if not defined) - psImageMaskType poorMask = pmConfigMaskGet("POOR.WARP", config); + psImageMaskType poorMask = pmConfigMaskGet("POOR.WARP", config); if (!poorMask) { - poorMask = 0x02; - pmConfigMaskSet (config, "POOR.WARP", poorMask); + poorMask = 0x02; + pmConfigMaskSet (config, "POOR.WARP", poorMask); } maskOut |= poorMask; Index: /trunk/pswarp/src/pswarpTransformReadout.c =================================================================== --- /trunk/pswarp/src/pswarpTransformReadout.c (revision 21367) +++ /trunk/pswarp/src/pswarpTransformReadout.c (revision 21368) @@ -6,6 +6,6 @@ * * @author IfA - * @version $Revision: 1.15 $ $Name: not supported by cvs2svn $ - * @date $Date: 2009-02-05 20:44:04 $ + * @version $Revision: 1.16 $ $Name: not supported by cvs2svn $ + * @date $Date: 2009-02-06 03:10:36 $ * Copyright 2009 Institute for Astronomy, University of Hawaii */ @@ -41,5 +41,9 @@ psAssert(mdok, "MASK.INPUT was not defined"); - int nThreads = psMetadataLookupS32(&mdok, config->arguments, "NTHREADS"); ///< Number of threads + if (!input->covariance) { + input->covariance = psImageCovarianceNone(); + } + + int nThreads = psMetadataLookupS32(&mdok, config->arguments, "NTHREADS"); // Number of threads if (!mdok) { nThreads = 0; @@ -85,11 +89,11 @@ // Interpolation options : move these from the arguments to explicit assignments psImageInterpolation *interp = psImageInterpolationAlloc(interpolationMode, input->image, - input->weight, input->mask, maskIn, + input->variance, input->mask, maskIn, NAN, NAN, maskBad, maskPoor, poorFrac, numKernels); - if (input->weight && !output->weight) { - output->weight = psImageAlloc(output->image->numCols, output->image->numRows, PS_TYPE_F32); - psImageInit(output->weight, NAN); + if (input->variance && !output->variance) { + output->variance = psImageAlloc(output->image->numCols, output->image->numRows, PS_TYPE_F32); + psImageInit(output->variance, NAN); } if ((input->mask || maskPoor || maskBad) && !output->mask) { @@ -135,5 +139,12 @@ psThreadJob *job = NULL; int xMin = output->image->numCols, xMax = 0, yMin = output->image->numRows, yMax = 0; // Bounds - int goodPixels = 0; ///< total number of good pixels across all tiles + int goodPixels = 0; // total number of good pixels across all tiles + psList *covariances = psMetadataLookupPtr(&mdok, output->analysis, + PSWARP_ANALYSIS_COVARIANCES); // Collection of covar. matrices + if (!covariances) { + covariances = psListAlloc(NULL); + psMetadataAddList(output->analysis, PS_LIST_TAIL, PSWARP_ANALYSIS_COVARIANCES, 0, + "Collection of covariance matrices", covariances); + } while ((job = psThreadJobGetDone()) != NULL) { if (job->args->n < 1) { @@ -147,4 +158,7 @@ yMin = PS_MIN(args->yMin, yMin); yMax = PS_MAX(args->yMax, yMax); + if (args->covariance) { + psListAdd(covariances, PS_LIST_TAIL, args->covariance); + } } psFree(job); Index: /trunk/pswarp/src/pswarpTransformTile.c =================================================================== --- /trunk/pswarp/src/pswarpTransformTile.c (revision 21367) +++ /trunk/pswarp/src/pswarpTransformTile.c (revision 21368) @@ -6,6 +6,6 @@ * * @author IfA - * @version $Revision: 1.8 $ $Name: not supported by cvs2svn $ - * @date $Date: 2009-02-05 20:44:04 $ + * @version $Revision: 1.9 $ $Name: not supported by cvs2svn $ + * @date $Date: 2009-02-06 03:10:36 $ * Copyright 2009 Institute for Astronomy, University of Hawaii */ @@ -20,4 +20,5 @@ psFree(args->interp); psFree(args->region); + psFree(args->covariance); return; } @@ -42,4 +43,5 @@ args->yMin = PS_MAX_S32; args->yMax = PS_MIN_S32; + args->covariance = NULL; return args; @@ -60,5 +62,5 @@ // Dereference images for convenience psF32 **outImageData = args->output->image->data.F32; - psF32 **outVarData = (args->output->weight) ? args->output->weight->data.F32 : NULL; + psF32 **outVarData = (args->output->variance) ? args->output->variance->data.F32 : NULL; psImageMaskType **outMaskData = (args->output->mask) ? args->output->mask->data.PS_TYPE_IMAGE_MASK_DATA : NULL; psImageMaskType **inMaskData = (args->input->mask) ? args->input->mask->data.PS_TYPE_IMAGE_MASK_DATA : NULL; @@ -92,9 +94,7 @@ } - double imageValue, varValue; // Value of image and variance map - // psImagePixelInterpolate determines the value at pixel coordinate (x,y) in child coordinates + double imageValue, varValue; // Value of image and variance map psImageMaskType maskValue = inMaskData ? inMaskData[(int)yIn][(int)xIn] : 0; // Value of mask - if (!psImageInterpolate(&imageValue, &varValue, &maskValue, xIn, yIn, args->interp)) { psError(PS_ERR_UNKNOWN, false, "Unable to interpolate image."); @@ -118,4 +118,17 @@ } + if (goodPixels > 0) { + float xOut = 0.5 * (xMin + xMax), yOut = 0.5 * (yMin + yMax); // Position of interest on output + double xIn, yIn; // Position of interest on input + pswarpMapApply(&xIn, &yIn, map, xOut + 0.5, yOut + 0.5); + // XXX Why are we subtracting the *output* col0,row0 from the *input* coordinates? + // I expect these are zero, so that it makes no difference, but it's distracting. + xIn -= outCol0; + yIn -= outRow0; + psKernel *kernel = psImageInterpolationKernel(xIn, yIn, args->interp->mode); // Interpolation kernel + args->covariance = psImageCovarianceCalculate(kernel, args->input->covariance); + psFree(kernel); + } + args->goodPixels = goodPixels; args->xMin = xMin;