Changeset 41808 for branches/eam_branches/ipp-dev-20210817/psModules

Timestamp:

Sep 13, 2021, 6:44:06 AM (5 years ago)

Author:

eugene

Message:

some tests of pattern.row implementation

File:

• branches/eam_branches/ipp-dev-20210817/psModules/src/detrend/pmPattern.c (modified) (6 diffs)

branches/eam_branches/ipp-dev-20210817/psModules/src/detrend/pmPattern.c

@@ -139 +139 @@
 // Measurement and application
 //////////////////////////////////////////////////////////////////////////////////////////////////////////////
+
+// USE_BACKGROUND_STDEV: if true, the analysis will use the measured robust stdev to clip the out-of-range pixles
+// if false, the stdev will be estimated based on the Poisson statistics of the median (sky level).
+# define USE_BACKGROUND_STDEV 0
 
 bool pmPatternRow(pmReadout *ro, int order, int iter, float rej, float thresh,
@@ -172 +176 @@
     }
 
-# if (0) 
+# if (USE_BACKGROUND_STDEV) 
     float lower = stats->robustMedian - thresh * stats->robustStdev; // Lower bound for data
     float upper = stats->robustMedian + thresh * stats->robustStdev; // Upper bound for data
     float background = stats->robustMedian;
 # else
-
     // the signal we are looking for is a small variation on top of the background.  if
     // the background is uniform with only read noise + sky noise, then the pixel-to-pixel
@@ -201 +204 @@
 # endif    
     
-    // the expected amplitude of the 2nd order term is 
-    // if the predicted sky poisson signal is larger than 
-
-    // I want to add a constraint to the fit so the amplitude or coeffs are |value| < x
+    pmCell *cell = ro->parent;
+    const char *cellName = psMetadataLookupStr(NULL, cell->concepts, "CELL.NAME"); // Name of cell
+    fprintf (stderr, "pattern row background %s: %f - %f - %f\n", cellName,lower, background, upper);
+
+
+    // XXX add code to skip fit if background stdev is large compared to expected signal
+    // XXX add a constraint to the fit so the |amplitude| < x (x ~ 20 - 25)
 
     psFree(stats);
@@ -217 +223 @@
 
     psStats *clip = psStatsAlloc(clipMean | clipStdev); // Clipping statistics
-    clip->clipIter = iter;
+    // XXX clip->clipIter = iter;
+    clip->clipIter = 1; // XXX skip iteration for a test
     clip->clipSigma = rej;
     psVector *clipMask = psVectorAlloc(numCols, PS_TYPE_VECTOR_MASK); // Mask for clipping
@@ -259 +266 @@
         float validXmax = -1;
 
+        // XXX can we do just as well fitting 1/3 of the pixels? (NOT REALLY)
+        // (x % 3) ||
         for (int x = 0; x < numCols; x++) {
-            if ((mask && mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & maskVal) ||
-                data->data.F32[x] < lower || data->data.F32[x] > upper) {
-                clipMask->data.PS_TYPE_VECTOR_MASK_DATA[x] = 0xFF;
+            if ((mask && mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & maskVal) ||
+                data->data.F32[x] < lower || data->data.F32[x] > upper) {
+                clipMask->data.PS_TYPE_VECTOR_MASK_DATA[x] = 0xFF;
             } else {
                 clipMask->data.PS_TYPE_VECTOR_MASK_DATA[x] = 0;
@@ -271 +280 @@
         }
 
+        // XXX how much time is spent in the fitting
         if (num < validNmin) {
             // Not enough points to fit