Index: /branches/eam_branches/ipp-20100823/psModules/src/camera/pmFPAMaskWeight.c =================================================================== --- /branches/eam_branches/ipp-20100823/psModules/src/camera/pmFPAMaskWeight.c (revision 29214) +++ /branches/eam_branches/ipp-20100823/psModules/src/camera/pmFPAMaskWeight.c (revision 29215) @@ -476,6 +476,29 @@ } - - +// find any pixels which are not already masked (with maskTest) which are not valid and raise maskSet bits +bool pmReadoutMaskInvalid (const pmReadout *readout, psImageMaskType maskTest, psImageMaskType maskSet) { + + if (!readout) return true; + + psImage *image = readout->image; + psImage *mask = readout->mask; + psImage *variance = readout->variance; + for (int y = 0; y < image->numRows; y++) { + for (int x = 0; x < image->numCols; x++) { + if (mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & maskTest) continue; + bool valid = false; + valid = isfinite(image->data.F32[y][x]); + if (valid && variance) { + valid &= isfinite(variance->data.F32[y][x]); + } + if (valid) continue; + mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] |= maskSet; + } + } + + return true; +} + +// raise maskVal for any invalid pixels bool pmReadoutMaskNonfinite(pmReadout *readout, psImageMaskType maskVal) { Index: /branches/eam_branches/ipp-20100823/psModules/src/camera/pmFPAMaskWeight.h =================================================================== --- /branches/eam_branches/ipp-20100823/psModules/src/camera/pmFPAMaskWeight.h (revision 29214) +++ /branches/eam_branches/ipp-20100823/psModules/src/camera/pmFPAMaskWeight.h (revision 29215) @@ -99,4 +99,7 @@ ); +// find any pixels which are not already masked (with maskTest) which are not valid and raise maskSet bits +bool pmReadoutMaskInvalid (const pmReadout *readout, psImageMaskType maskTest, psImageMaskType maskSet); + /// Apply a mask to the image and variance map /// Index: /branches/eam_branches/ipp-20100823/psModules/src/imcombine/pmStack.c =================================================================== --- /branches/eam_branches/ipp-20100823/psModules/src/imcombine/pmStack.c (revision 29214) +++ /branches/eam_branches/ipp-20100823/psModules/src/imcombine/pmStack.c (revision 29215) @@ -19,4 +19,5 @@ #include +#include #include @@ -34,9 +35,16 @@ -//#define TESTING // Enable test output -//#define TEST_X 843-1 // x coordinate to examine -//#define TEST_Y 813-1 // y coordinate to examine -//#define TEST_RADIUS 0 // Radius to examine - +#define TESTING // Enable test output +// #define TEST_X 4963 // x coordinate to examine +// #define TEST_Y 5353 // y coordinate to examine +#define TEST_X 40 // x coordinate to examine +#define TEST_Y 40 // y coordinate to examine +#define TEST_RADIUS 0.5 // Radius to examine + +# ifdef TESTING +# define CHECKPIX(XPIX,YPIX,MSG,...) { if (PS_SQR(XPIX - TEST_X) + PS_SQR(YPIX - TEST_Y) <= PS_SQR(TEST_RADIUS)) { fprintf(stderr,MSG,__VA_ARGS__); } } +# else +# define CHECKPIX(XPIX,YPIX,MSG,...) { } +# endif // Data structure for use as a buffer when combining pixels @@ -250,9 +258,5 @@ ) { -#ifdef TESTING - if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { - fprintf(stderr, "Marking image %d, pixel %d,%d for inspection\n", source, x, y); - } -#endif + CHECKPIX(x, y, "Marking image %d, pixel %d,%d for inspection\n", source, x, y); pmStackData *data = inputs->data[source]; // Stack data of interest if (!data) { @@ -272,9 +276,5 @@ ) { -#ifdef TESTING - if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { - fprintf(stderr, "Marking pixel image %d, pixel %d,%d for rejection\n", source, x, y); - } -#endif + CHECKPIX(x, y, "Marking pixel image %d, pixel %d,%d for rejection\n", source, x, y); pmStackData *data = inputs->data[source]; // Stack data of interest if (!data) { @@ -290,4 +290,6 @@ static void combineExtract(int *num, // Number of good pixels bool *suspect, // Any suspect pixels? + psImageMaskType *badMask, // OR of all bad (masked) pixels + psImageMaskType *goodMask, // OR of all good (unmasked) pixels combineBuffer *buffer, // Buffer with vectors psImage *image, // Combined image, for output @@ -300,6 +302,6 @@ const psVector *reject, // Indices of pixels to reject, or NULL int x, int y, // Coordinates of interest; frame of output image - psImageMaskType maskVal, // Value to mask - psImageMaskType maskSuspect // Value to suspect + psImageMaskType badMaskBits, // Value to mask as 'bad' + psImageMaskType suspectMaskBits // Value to mask as 'suspect' ) { @@ -322,4 +324,8 @@ } + // mask values to store possible mask combinations + *badMask = 0; + *goodMask = 0xffff; + // Extract the pixel and mask data int numGood = 0; // Number of good pixels @@ -328,5 +334,22 @@ // should be because of how pixelMapGenerate works if (reject && reject->data.U16[j] == i) { + // pixels can be rejected because: + // 1) only 1 input pixel (and 'safe' is true) + // 2) only 2 input pixels were available and they were mutually inconsistent (or variance info was missing) + // 3) NOTE : ifdef'ed out code for 3 inputs case + // 4) outlier from sample of N pixels + // 5) some of these may have been suspect. + // XXX raise a specific mask bit for these (currently results in BLANK) j++; + + pmStackData *data = inputs->data[i]; // Stack data of interest + if (data) { + int xIn = x - data->readout->col0, yIn = y - data->readout->row0; // Coordinates on input readout + psImage *mask = data->readout->mask; // Mask of interest + *badMask |= mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn]; // save the bad bits used + CHECKPIX(x, y, "reject: adding bit to mask: %d : %x (badMask = %x)\n", i, mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn], *badMask); + } else { + CHECKPIX(x, y, "reject: no item in data (badMask = %x)\n", *badMask); + } continue; } @@ -334,4 +357,5 @@ pmStackData *data = inputs->data[i]; // Stack data of interest if (!data) { + CHECKPIX(x, y, "skip: no item in data (badMask = %x)\n", *badMask); continue; } @@ -339,10 +363,13 @@ int xIn = x - data->readout->col0, yIn = y - data->readout->row0; // Coordinates on input readout psImage *mask = data->readout->mask; // Mask of interest - if (mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn] & maskVal) { + if (mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn] & badMaskBits) { + *badMask |= (mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn] & badMaskBits); // save the bad bits used + CHECKPIX(x, y, "skip: adding bit to mask: %d : %x (badMask = %x)\n", i, mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn], *badMask); continue; } - pixelSuspects->data.U8[numGood] = mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn] & maskSuspect ? - true : false; + pixelSuspects->data.U8[numGood] = mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn] & suspectMaskBits ? true : false; + + *goodMask &= mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn]; // save the mask bits still used psImage *image = data->readout->image; // Image of interest @@ -356,4 +383,6 @@ pixelSources->data.U16[numGood] = i; numGood++; + + CHECKPIX(x, y, "keep: %d : %x (badMask = %x)\n", i, mask->data.PS_TYPE_IMAGE_MASK_DATA[yIn][xIn], *badMask); } pixelData->n = numGood; @@ -370,8 +399,8 @@ if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { for (int i = 0; i < numGood; i++) { - fprintf(stderr, "Input %d, pixel %d,%d (%" PRIu16 "): %f %f (%f) %f %f %d\n", + fprintf(stderr, "Input %d, pixel %d,%d (%" PRIu16 "): %f %f (%f) %f %f %d %x %x -> %x %x\n", i, x, y, pixelSources->data.U16[i], pixelData->data.F32[i], pixelVariances->data.F32[i], addVariance->data.F32[i], pixelWeights->data.F32[i], pixelExps->data.F32[i], - pixelSuspects->data.U8[i]); + pixelSuspects->data.U8[i], badMaskBits, suspectMaskBits, *badMask, *goodMask); } } @@ -380,5 +409,4 @@ return; } - // Combine pixels @@ -392,5 +420,7 @@ combineBuffer *buffer, // Buffer with vectors int x, int y, // Coordinates of interest; frame of output image - psImageMaskType bad, // Value for bad pixels + psImageMaskType blankMask, // Value for empty pixels + psImageMaskType badMask, // Value for bad pixels + psImageMaskType goodMask, // Value for good pixels bool safe, // Safe combination? float invTotalWeight // Inverse of total weight for all inputs @@ -404,15 +434,17 @@ // Default option is that the pixel is bad float imageValue = NAN, varianceValue = NAN; // Value for combined image and variance map - psImageMaskType maskValue = bad; // Value for combined mask float expValue = 0.0, expWeightValue = NAN; // Exposure value (straight, and weighted) + + // default output mask value. badMask is the OR of all unused input pixels. + // if there are no input pixels, this value will be 0, in which case we want to set the output pixel to BLANK. + // if there are only good input pixels, and they do not result in a valid pixel, we still want to set this to BLANK. + psImageMaskType maskValue = badMask ? badMask : blankMask; // Value for combined mask + + CHECKPIX(x, y, "bad vs good : %x %x %x\n", maskValue, badMask, blankMask); switch (num) { case 0: { // Nothing to combine: it's bad -#ifdef TESTING - if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { - fprintf(stderr, "No inputs to combine, pixel %d,%d is bad.\n", x, y); - } -#endif + CHECKPIX(x, y, "No inputs to combine, pixel %d,%d is bad.\n", x, y); break; } @@ -428,17 +460,9 @@ expWeightValue = pixelExps->data.F32[0]; } - maskValue = 0; -#ifdef TESTING - if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { - fprintf(stderr, "Single input to combine, safety off, pixel %d,%d --> %f\n", - x, y, imageValue); - } -#endif - } -#ifdef TESTING - else if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { - fprintf(stderr, "Single input to combine, safety on, pixel %d,%d is bad.\n", x, y); - } -#endif + maskValue = goodMask; + CHECKPIX(x, y, "Single input to combine, safety off, pixel %d,%d --> %f\n", x, y, imageValue); + } else { + CHECKPIX(x, y, "Single input to combine, safety on, pixel %d,%d is bad.\n", x, y); + } break; } @@ -446,22 +470,11 @@ // Automatically accept the mean of the pixels only if we're not playing safe if (!safe) { - if (combinationMeanVariance(&imageValue, &varianceValue, &expValue, &expWeightValue, - pixelData, pixelVariances, pixelExps, pixelWeights)) { - maskValue = 0; -#ifdef TESTING - if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { - fprintf(stderr, "Two inputs to combine using unsafe, pixel %d,%d --> %f %f\n", - x, y, imageValue, varianceValue); - } -#endif + if (combinationMeanVariance(&imageValue, &varianceValue, &expValue, &expWeightValue, pixelData, pixelVariances, pixelExps, pixelWeights)) { + maskValue = goodMask; + CHECKPIX(x, y, "Two inputs to combine using unsafe, pixel %d,%d --> %f %f\n", x, y, imageValue, varianceValue); } + } else { + CHECKPIX(x, y, "Two inputs to combine, safety on, pixel %d,%d is bad\n", x, y); } -#ifdef TESTING - else { - if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { - fprintf(stderr, "Two inputs to combine, safety on, pixel %d,%d is bad\n", x, y); - } - } -#endif break; } @@ -472,10 +485,6 @@ break; } - maskValue = 0; -#ifdef TESTING - if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { - fprintf(stderr, "Combined inputs, pixel %d,%d --> %f %f\n", x, y, imageValue, varianceValue); - } -#endif + maskValue = goodMask; + CHECKPIX(x, y, "Combined inputs, pixel %d,%d --> %f %f\n", x, y, imageValue, varianceValue); break; } @@ -522,10 +531,5 @@ int numIter = PS_MAX(iter * num, 1); // Number of iterations -#ifdef TESTING - if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { - fprintf(stderr, "Testing pixel %d,%d: %d %f %f %f %d %d\n", - x, y, numIter, rej, sys, olympic, useVariance, safe); - } -#endif + CHECKPIX(x, y, "Testing pixel %d,%d: %d %f %f %f %d %d\n", x, y, numIter, rej, sys, olympic, useVariance, safe); psVector *pixelData = buffer->pixels; // Values for the pixel of interest @@ -548,12 +552,7 @@ // Systematic error contributes to the rejection level float sysVar = PS_SQR(sys * pixelData->data.F32[i]); -#ifdef TESTING - // Correct variance for comparison against weighted mean including itself - float compare = 1.0 - pixelWeights->data.F32[i] / sumWeights; - if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { - fprintf(stderr, "Variance %d (%d), pixel %d,%d: %f %f %f\n", i, pixelSources->data.U16[i], - x, y, pixelVariances->data.F32[i], sysVar, compare); - } -#endif + CHECKPIX(x, y, "Variance %d (%d), pixel %d,%d: %f %f %f\n", + i, pixelSources->data.U16[i], x, y, + pixelVariances->data.F32[i], sysVar, 1.0 - pixelWeights->data.F32[i] / sumWeights); pixelLimits->data.F32[i] = rej2 * (pixelVariances->data.F32[i] + sysVar); } @@ -593,10 +592,5 @@ combineMarkInspect(inputs, x, y, pixelSources->data.U16[1]); } -#ifdef TESTING - if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { - fprintf(stderr, "Flagged both inputs for pixel %d,%d (%f > %f x %f\n)", - x, y, diff, rej, sqrtf(sigma2)); - } -#endif + CHECKPIX(x, y, "Flagged both inputs for pixel %d,%d (%f > %f x %f\n)", x, y, diff, rej, sqrtf(sigma2)); } } else if (i == 0 && safe) { @@ -708,17 +702,9 @@ float median = combinationWeightedOlympic(pixelData, pixelWeights, olympic, buffer->sort); // Median for stack -#ifdef TESTING - if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { - fprintf(stderr, "Flag with variance pixel %d,%d: median = %f\n", x, y, median); - } -#endif + CHECKPIX(x, y, "Flag with variance pixel %d,%d: median = %f\n", x, y, median); float worst = -INFINITY; // Largest deviation for (int j = 0; j < num; j++) { float diff = pixelData->data.F32[j] - median; // Difference from expected -#ifdef TESTING - if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { - fprintf(stderr, "Testing input %d for pixel %d,%d: %f\n", j, x, y, diff); - } -#endif + CHECKPIX(x, y, "Testing input %d for pixel %d,%d: %f\n", j, x, y, diff); // Comparing squares --- cheaper than lots of sqrts @@ -737,11 +723,5 @@ combinationMedianStdev(&median, &stdev, pixelData, buffer->sort); float limit = rej * stdev; // Rejection limit -#ifdef TESTING - if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { - fprintf(stderr, - "Flag without variance pixel %d,%d; median = %f, stdev = %f, limit = %f\n", - x, y, median, stdev, limit); - } -#endif + CHECKPIX(x, y, "Flag without variance pixel %d,%d; median = %f, stdev = %f, limit = %f\n", x, y, median, stdev, limit); float worst = -INFINITY; // Largest deviation for (int j = 0; j < num; j++) { @@ -763,9 +743,5 @@ if (maskIndex >= 0) { if (suspect) { -#ifdef TESTING - if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { - fprintf(stderr, "Throwing out all suspect pixels for %d,%d\n", x, y); - } -#endif + CHECKPIX(x, y, "Throwing out all suspect pixels for %d,%d\n", x, y); // Throw out all suspect pixels int numGood = 0; // Number of good pixels @@ -796,9 +772,5 @@ } else { // Throw out masked pixel -#ifdef TESTING - if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { - fprintf(stderr, "Throwing out input %d for pixel %d,%d\n", maskIndex, x, y); - } -#endif + CHECKPIX(x, y, "Throwing out input %d for pixel %d,%d\n", maskIndex, x, y); combineMarkInspect(inputs, x, y, pixelSources->data.U16[maskIndex]); int numGood = 0; // Number of good pixels @@ -999,9 +971,19 @@ /// Stack input images -bool pmStackCombine(pmReadout *combined, pmReadout *expmaps, psArray *input, - psImageMaskType maskVal, psImageMaskType maskSuspect, - psImageMaskType bad, int kernelSize, - float iter, float rej, float sys, float olympic, - bool useVariance, bool safe, bool rejection) +bool pmStackCombine( + pmReadout *combined, // output stacked readout + pmReadout *expmaps, // output exposure map information + psArray *input, // input exposures + psImageMaskType badMaskBits, // treat these bits as 'bad' + psImageMaskType suspectMaskBits, // treat these bits as 'suspect' + psImageMaskType blankMaskBits, // use this mask value for pixels missing input data (distinguish between Ninput = 0 and Ngood = 0?) + int kernelSize, + float iter, + float rej, + float sys, + float olympic, + bool useVariance, + bool safe, + bool rejection) { bool haveVariances; // Do we have the variance maps? @@ -1012,5 +994,5 @@ } PS_ASSERT_INT_NONNEGATIVE(kernelSize, false); - PS_ASSERT_INT_POSITIVE(bad, false); + PS_ASSERT_INT_POSITIVE(blankMaskBits, false); if (isnan(rej)) { PS_ASSERT_FLOAT_EQUAL(iter, 0, false); @@ -1070,5 +1052,5 @@ } psFree(stack); - pmReadoutUpdateSize(combined, minInputCols, minInputRows, xSize, ySize, true, true, bad); + pmReadoutUpdateSize(combined, minInputCols, minInputRows, xSize, ySize, true, true, blankMaskBits); psTrace("psModules.imcombine", 1, "Have for combination [%d:%d,%d:%d] (%dx%d)\n", minInputCols, maxInputCols, minInputRows, maxInputRows, xSize, ySize); @@ -1131,10 +1113,12 @@ for (int y = minInputRows; y < maxInputRows; y++) { for (int x = minInputCols; x < maxInputCols; x++) { + CHECKPIX(x, y, "Combining pixel %d,%d: %x %x %f %f %f %f %d %d %d\n", x, y, badMaskBits, blankMaskBits, iter, rej, sys, olympic, useVariance, safe, rejection); + #ifdef TESTING - if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { - fprintf(stderr, "Combining pixel %d,%d: %x %x %f %f %f %f %d %d %d\n", - x, y, maskVal, bad, iter, rej, sys, olympic, useVariance, safe, rejection); - } -#endif + if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) { + fprintf (stderr, "combine\n"); + } +# endif + psVector *reject = NULL; // Images to reject for this pixel if (rejection) { @@ -1154,11 +1138,11 @@ #endif } - - int num; // Number of good pixels - bool suspect; // Suspect pixels in stack? - combineExtract(&num, &suspect, buffer, combinedImage, combinedMask, combinedVariance, - input, weights, exps, addVariance, reject, x, y, maskVal, maskSuspect); - combinePixels(combinedImage, combinedMask, combinedVariance, exp, expnum, expweight, - num, buffer, x, y, bad, safe, totalExpWeight); + + int num; // Number of good pixels + bool suspect; // Suspect pixels in stack? + psImageMaskType badMask = 0; // OR of mask bits in all bad input pixels + psImageMaskType goodMask = 0; // OR of mask bits in all good input pixels + combineExtract(&num, &suspect, &badMask, &goodMask, buffer, combinedImage, combinedMask, combinedVariance, input, weights, exps, addVariance, reject, x, y, badMaskBits, suspectMaskBits); + combinePixels(combinedImage, combinedMask, combinedVariance, exp, expnum, expweight, num, buffer, x, y, blankMaskBits, badMask, goodMask, safe, totalExpWeight); if (iter > 0) { Index: /branches/eam_branches/ipp-20100823/psModules/src/imcombine/pmSubtractionMatch.c =================================================================== --- /branches/eam_branches/ipp-20100823/psModules/src/imcombine/pmSubtractionMatch.c (revision 29214) +++ /branches/eam_branches/ipp-20100823/psModules/src/imcombine/pmSubtractionMatch.c (revision 29215) @@ -275,26 +275,26 @@ } -bool pmSubtractionMaskInvalid (const pmReadout *readout, psImageMaskType maskVal) { - - if (!readout) return true; - - psImage *image = readout->image; - psImage *mask = readout->mask; - psImage *variance = readout->variance; - for (int y = 0; y < image->numRows; y++) { - for (int x = 0; x < image->numCols; x++) { - if (mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & maskVal) continue; - bool valid = false; - valid = isfinite(image->data.F32[y][x]); - if (variance) { - valid &= isfinite(variance->data.F32[y][x]); - } - if (valid) continue; - mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] = maskVal; - } - } - - return true; -} +// bool pmSubtractionMaskInvalid (const pmReadout *readout, psImageMaskType maskVal) { +// +// if (!readout) return true; +// +// psImage *image = readout->image; +// psImage *mask = readout->mask; +// psImage *variance = readout->variance; +// for (int y = 0; y < image->numRows; y++) { +// for (int x = 0; x < image->numCols; x++) { +// if (mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & maskVal) continue; +// bool valid = false; +// valid = isfinite(image->data.F32[y][x]); +// if (variance) { +// valid &= isfinite(variance->data.F32[y][x]); +// } +// if (valid) continue; +// mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] = maskVal; +// } +// } +// +// return true; +// } bool pmSubtractionMatchPrecalc(pmReadout *conv1, pmReadout *conv2, const pmReadout *ro1, const pmReadout *ro2, @@ -387,6 +387,7 @@ } - pmSubtractionMaskInvalid(ro1, maskVal); - pmSubtractionMaskInvalid(ro2, maskVal); + // XXX this is done before calling this function + // pmSubtractionMaskInvalid(ro1, maskVal); + // pmSubtractionMaskInvalid(ro2, maskVal); // General background subtraction, since this is done in pmSubtractionMatch @@ -560,6 +561,6 @@ memCheck("start"); - pmSubtractionMaskInvalid(ro1, maskVal); - pmSubtractionMaskInvalid(ro2, maskVal); + // pmSubtractionMaskInvalid(ro1, maskVal); + // pmSubtractionMaskInvalid(ro2, maskVal); psRegion bounds = psRegionSet(NAN, NAN, NAN, NAN); // Bounds of valid pixels Index: /branches/eam_branches/ipp-20100823/psModules/src/imcombine/pmSubtractionStamps.c =================================================================== --- /branches/eam_branches/ipp-20100823/psModules/src/imcombine/pmSubtractionStamps.c (revision 29214) +++ /branches/eam_branches/ipp-20100823/psModules/src/imcombine/pmSubtractionStamps.c (revision 29215) @@ -106,10 +106,10 @@ if (xMax < xMin) { - fprintf (stderr, "%f,%f : x-border\n", xRaw, yRaw); - return false; + // fprintf (stderr, "%f,%f : x-border\n", xRaw, yRaw); + return false; } if (yMax < yMin) { - fprintf (stderr, "%f,%f : y-border\n", xRaw, yRaw); - return false; + // fprintf (stderr, "%f,%f : y-border\n", xRaw, yRaw); + return false; } @@ -458,5 +458,5 @@ // Take stamps off the top of the (sorted) list for (int j = xList->n - 1; j >= 0 && !goodStamp; j--) { - fprintf (stderr, "%d : xList: %ld elements\n", i, xList->n); + // fprintf (stderr, "%d : xList: %ld elements\n", i, xList->n); // Chop off the top of the list xList->n = j; @@ -706,5 +706,5 @@ // storage vector for flux data - psStats *stats = psStatsAlloc (PS_STAT_ROBUST_MEDIAN | PS_STAT_ROBUST_STDEV); + psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEDIAN); psVector *flux1 = psVectorAllocEmpty(2*stamps->num, PS_TYPE_F32); psVector *flux2 = psVectorAllocEmpty(2*stamps->num, PS_TYPE_F32); @@ -734,5 +734,5 @@ psAbort ("failed to generate stats"); } - float f1 = stats->robustMedian; + float f1 = stats->sampleMedian; psStatsInit (stats); @@ -740,5 +740,5 @@ psAbort ("failed to generate stats"); } - float f2 = stats->robustMedian; + float f2 = stats->sampleMedian; stamps->window1->kernel[y][x] = f1; @@ -755,4 +755,16 @@ } } + +#if 0 + { + psFits *fits = NULL; + fits = psFitsOpen ("window1.fits", "w"); + psFitsWriteImage (fits, NULL, stamps->window1->image, 0, NULL); + psFitsClose (fits); + fits = psFitsOpen ("window2.fits", "w"); + psFitsWriteImage (fits, NULL, stamps->window2->image, 0, NULL); + psFitsClose (fits); + } +#endif psTrace("psModules.imcombine", 3, "Window total (1): %f, threshold: %f\n", sum1, (1.0 - stamps->normFrac) * sum1);