Index: /trunk/ppSub/src/ppSubArguments.c =================================================================== --- /trunk/ppSub/src/ppSubArguments.c (revision 19164) +++ /trunk/ppSub/src/ppSubArguments.c (revision 19165) @@ -200,5 +200,6 @@ psMetadataAddF32(arguments, PS_LIST_TAIL, "-rej", 0, "Rejection thresold (sigma)", NAN); psMetadataAddU8(arguments, PS_LIST_TAIL, "-mask-bad", 0, "Mask value for bad pixels", 0); - psMetadataAddU8(arguments, PS_LIST_TAIL, "-mask-blank", 0, "Mask value for blank region", 0); + psMetadataAddU8(arguments, PS_LIST_TAIL, "-mask-poor", 0, "Mask value for poor pixels", 0); + psMetadataAddF32(arguments, PS_LIST_TAIL, "-poor-frac", 0, "Fraction of weight 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); @@ -286,7 +287,9 @@ VALUE_ARG_RECIPE_FLOAT("-badfrac", "BADFRAC", F32); VALUE_ARG_RECIPE_FLOAT("-penalty", "PENALTY", F32); - - valueArgRecipeStr(arguments, recipe, "-mask-bad", "MASK.BAD", recipe); - valueArgRecipeStr(arguments, recipe, "-mask-blank", "MASK.BLANK", recipe); + VALUE_ARG_RECIPE_FLOAT("-poor-frac", "POOR.FRACTION", F32); + + 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); Index: /trunk/ppSub/src/ppSubLoop.c =================================================================== --- /trunk/ppSub/src/ppSubLoop.c (revision 19164) +++ /trunk/ppSub/src/ppSubLoop.c (revision 19165) @@ -142,5 +142,5 @@ return false; } - ppStatsFPA(stats, output->fpa, view, pmConfigMaskGet("BLANK", config), config); + ppStatsFPA(stats, output->fpa, view, pmConfigMaskGet("MASK.VALUE", config), config); } Index: /trunk/ppSub/src/ppSubReadout.c =================================================================== --- /trunk/ppSub/src/ppSubReadout.c (revision 19164) +++ /trunk/ppSub/src/ppSubReadout.c (revision 19165) @@ -98,6 +98,10 @@ int binning = psMetadataLookupS32(NULL, recipe, "SPAM.BINNING"); // Binning for SPAM kernel float penalty = psMetadataLookupF32(NULL, recipe, "PENALTY"); // Penalty for wideness - psMaskType maskIn = pmConfigMaskGet("MASK.VALUE", config); // Bits to mask going in to pmSubtractionMatch - psMaskType maskOut = pmConfigMaskGet("BLANK", config); // Bits to mask after pmSubtractionMatch + psString maskValStr = psMetadataLookupStr(NULL, recipe, "MASK.IN"); // Name of bits to mask going in + psMaskType maskVal = pmConfigMaskGet(maskValStr, config); // Bits to mask going in to pmSubtractionMatch + psString maskPoorStr = psMetadataLookupStr(NULL, recipe, "MASK.POOR"); // Name of bits to mask for poor + psMaskType maskPoor = pmConfigMaskGet(maskPoorStr, config); // Bits to mask for poor pixels + psString maskBadStr = psMetadataLookupStr(NULL, recipe, "MASK.BAD"); // Name of bits to mask for bad + psMaskType maskBad = pmConfigMaskGet(maskBadStr, config); // Bits to mask for bad pixels float badFrac = psMetadataLookupF32(NULL, recipe, "BADFRAC"); // Maximum bad fraction const char *stampsName = psMetadataLookupStr(&mdok, config->arguments, "STAMPS"); // Filename for stamps @@ -113,5 +117,13 @@ int renormWidth = psMetadataLookupS32(&mdok, recipe, "RENORM.WIDTH"); // Width for renormalise - pmSubtractionMode mode = dual ? PM_SUBTRACTION_MODE_DUAL : PM_SUBTRACTION_MODE_UNSURE; // Subtraction mode + psString interpModeStr = psMetadataLookupStr(&mdok, recipe, "INTERPOLATION"); // Interpolation mode + psImageInterpolateMode interpMode = psImageInterpolateModeFromString(interpModeStr); // Interp + if (interpMode == PS_INTERPOLATE_NONE) { + psError(PS_ERR_BAD_PARAMETER_VALUE, false, "Unknown interpolation mode: %s", interpModeStr); + return false; + } + float poorFrac = psMetadataLookupF32(&mdok, recipe, "POOR.FRACTION"); // Fraction for "poor" + + pmSubtractionMode subMode = dual ? PM_SUBTRACTION_MODE_DUAL : PM_SUBTRACTION_MODE_UNSURE; // Subtracn mode // Generate masks if they don't exist @@ -134,9 +146,9 @@ } - if (!pmReadoutMaskNonfinite(inRO, maskIn)) { + if (!pmReadoutMaskNonfinite(inRO, pmConfigMaskGet("SAT", config))) { psError(PS_ERR_UNKNOWN, false, "Unable to mask non-finite pixels in input."); return false; } - if (!pmReadoutMaskNonfinite(refRO, maskIn)) { + if (!pmReadoutMaskNonfinite(refRO, pmConfigMaskGet("SAT", config))) { psError(PS_ERR_UNKNOWN, false, "Unable to mask non-finite pixels in reference."); return false; @@ -153,8 +165,20 @@ } + // Interpolate over bad pixels, so the bad pixels don't explode + if (!pmReadoutInterpolateBadPixels(inRO, maskVal, interpMode, poorFrac, maskPoor, maskBad)) { + psError(PS_ERR_UNKNOWN, false, "Unable to interpolate bad pixels for input image."); + return false; + } + if (!pmReadoutInterpolateBadPixels(refRO, maskVal, interpMode, poorFrac, maskPoor, maskBad)) { + psError(PS_ERR_UNKNOWN, false, "Unable to interpolate bad pixels for reference image."); + return false; + } + maskVal |= maskBad; + + // Match the PSFs if (!pmSubtractionMatch(inConv, refConv, inRO, refRO, footprint, regionSize, spacing, threshold, sources, stampsName, type, size, order, widths, orders, inner, ringsOrder, binning, - penalty, optimum, optWidths, optOrder, optThresh, iter, rej, maskIn, - maskOut, badFrac, mode)) { + penalty, optimum, optWidths, optOrder, optThresh, iter, rej, maskVal, + maskBad, maskPoor, poorFrac, badFrac, subMode)) { psError(PS_ERR_UNKNOWN, false, "Unable to match images."); psFree(inConv); @@ -221,5 +245,5 @@ #ifdef TESTING { - pmReadoutMaskApply(minuend, maskOut); + pmReadoutMaskApply(minuend, maskVal); psFits *fits = psFitsOpen("minuend.fits", "w"); psFitsWriteImage(fits, NULL, minuend->image, 0, NULL); @@ -227,5 +251,5 @@ } { - pmReadoutMaskApply(subtrahend, maskOut); + pmReadoutMaskApply(subtrahend, maskVal); psFits *fits = psFitsOpen("subtrahend.fits", "w"); psFitsWriteImage(fits, NULL, subtrahend->image, 0, NULL); @@ -267,5 +291,5 @@ // set maskValue and markValue in the psphot recipe - psMaskType maskValue = pmConfigMaskGet("BLANK", config); // Bits to mask + psMaskType maskValue = maskVal; psMaskType markValue = pmConfigMaskGet("MARK.VALUE", config); // Bits to use for marking psMetadataAddU8(recipe, PS_LIST_TAIL, "MASK.PSPHOT", PS_META_REPLACE, "Bits to mask", maskValue); @@ -306,5 +330,5 @@ outRO->data_exists = outCell->data_exists = outCell->parent->data_exists = true; - pmReadoutMaskApply(outRO, maskOut); + pmReadoutMaskApply(outRO, maskBad); psFree(inConv); @@ -314,5 +338,5 @@ for (int y = 0; y < outRO->image->numRows; y++) { for (int x = 0; x < outRO->image->numCols; x++) { - if (isnan(outRO->image->data.F32[y][x]) && !(outRO->mask->data.U8[y][x] & maskOut)) { + if (isnan(outRO->image->data.F32[y][x]) && !(outRO->mask->data.U8[y][x] & maskVal)) { printf("Unmasked NAN at %d %d --> %d\n", x, y, outRO->mask->data.U8[y][x]); } @@ -329,5 +353,5 @@ if (renorm) { psRandom *rng = psRandomAlloc(PS_RANDOM_TAUS, 0); // Random number generator - if (!pmReadoutWeightRenorm(outRO, maskOut, PS_STAT_ROBUST_MEDIAN, PS_STAT_ROBUST_STDEV, + if (!pmReadoutWeightRenorm(outRO, maskBad, PS_STAT_ROBUST_MEDIAN, PS_STAT_ROBUST_STDEV, renormWidth, rng)) { psError(PS_ERR_UNKNOWN, false, "Unable to renormalise weights for photometry.");