Changeset 26007 for branches/pap/ppStack/src/ppStackMatch.c

Timestamp:

Nov 2, 2009, 5:08:31 PM (17 years ago)

Author:

Paul Price

Message:

Serious rework of pmStackCombine for clarity and some additional functionality. Split the combinePixels function into two operating modes (test and combine) to make clear what we're doing; pmStackCombine could be split into two separate functions (e.g., pmStackTest and pmStackCombine) in the future, since the testing part doesn't need to produce any output image, but this involves moving a lot of stuff around in ppStack, so not doing it yet. Verified that variance scaling is correct, and softening of variance is correct. Big improvement in minimising loss of good pixels comes from throwing out only the most variant pixel on each rejection pass (consequently changed iterations parameter to iterations per input); this involved another rework of the core combinePixels function (moving the switch on the number of valid inputs inside the iteration). Much happier with the outlier rejection now.

File:

• branches/pap/ppStack/src/ppStackMatch.c (modified) (7 diffs)

branches/pap/ppStack/src/ppStackMatch.c

@@ -18 +18 @@
 #define COVAR_FRAC 0.01                 // Truncation fraction for covariance matrix
 
-#define TESTING                         // Enable debugging output
+//#define TESTING                         // Enable debugging output
 
 #ifdef TESTING
@@ -167 +167 @@
     )
 {
+#if 1
     bool mdok; // Status of metadata lookups
 
@@ -192 +193 @@
     psImageMaskType maskBad = pmConfigMaskGet("BLANK", config); // Bits to mask
 
+    psImageCovarianceTransfer(readout->variance, readout->covariance);
     return pmReadoutVarianceRenormalise(readout, maskBad, num, minValid, maxValid);
+#else
+    return true;
+#endif
 }
 
@@ -212 +217 @@
     int size = psMetadataLookupS32(NULL, ppsub, "KERNEL.SIZE"); // Kernel half-size
 
-    psString maskValStr = psMetadataLookupStr(NULL, recipe, "MASK.VAL"); // Name of bits to mask going in
+    psString maskValStr = psMetadataLookupStr(NULL, ppsub, "MASK.VAL"); // Name of bits to mask going in
     psImageMaskType 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
@@ -377 +382 @@
             }
 #endif
+
+            fprintf(stderr, "vf = %f\n", psImageCovarianceFactor(readout->covariance));
+
 
             if (threads > 0) {
@@ -516 +524 @@
             psFree(iter);
             options->matchChi2->data.F32[index] = sum / (psImageCovarianceFactor(readout->covariance) * num);
+            fprintf(stderr, "chi2 = %f ; vf = %f\n", sum/num, psImageCovarianceFactor(readout->covariance));
         }
 
@@ -542 +551 @@
     psRandom *rng = psRandomAlloc(PS_RANDOM_TAUS); // Random number generator
     if (!psImageBackground(bg, NULL, readout->image, readout->mask, maskVal | maskBad, rng)) {
-      psWarning("Can't measure background for image.");
-      psErrorClear();
+        psWarning("Can't measure background for image.");
+        psErrorClear();
     } else {
-      if (!psMetadataLookupBool(NULL, config->arguments, "PPSTACK.SKIP.BG.SUB")) {
-        psLogMsg("ppStack", PS_LOG_INFO, "Correcting convolved image background by %lf (+/- %lf)",
-                 psStatsGetValue(bg, PS_STAT_ROBUST_MEDIAN), psStatsGetValue(bg, PS_STAT_ROBUST_STDEV));
-        (void)psBinaryOp(readout->image, readout->image, "-",
-                         psScalarAlloc(psStatsGetValue(bg, PS_STAT_ROBUST_MEDIAN), PS_TYPE_F32));
-      }
+        if (!psMetadataLookupBool(NULL, config->arguments, "PPSTACK.SKIP.BG.SUB")) {
+            psLogMsg("ppStack", PS_LOG_INFO, "Correcting convolved image background by %lf (+/- %lf)",
+                     psStatsGetValue(bg, PS_STAT_ROBUST_MEDIAN), psStatsGetValue(bg, PS_STAT_ROBUST_STDEV));
+            (void)psBinaryOp(readout->image, readout->image, "-",
+                             psScalarAlloc(psStatsGetValue(bg, PS_STAT_ROBUST_MEDIAN), PS_TYPE_F32));
+        }
     }