Changeset 37778

Timestamp:

Jan 2, 2015, 11:46:31 AM (12 years ago)

Author:

eugene

Message:

repair nan-pixels in psf image maps

Location:

tags/ipp-20141224/psModules/src/objects

Files:

• pmPSFtryMakePSF.c (modified) (2 diffs)
• pmTrend2D.c (modified) (2 diffs)

tags/ipp-20141224/psModules/src/objects/pmPSFtryMakePSF.c

@@ -45 +45 @@
 #include "pmSourceVisual.h"
 
+bool pmPSF_DataDump (char *filename, psVector *x, psVector *y, psVector *e0, psVector *e1, psVector *e2, psVector *mask);
+
 /*****************************************************************************
 pmPSFFromPSFtry (psfTry): build a PSF model from a collection of source->modelEXT entries
@@ -356 +358 @@
 }
 
+bool pmPSF_DataDump (char *filename, psVector *x, psVector *y, psVector *e0, psVector *e1, psVector *e2, psVector *mask) {
+
+
+  FILE *f = fopen (filename, "w");
+  if (!f) return false;
+
+  for (int i = 0; i < x->n; i++) {
+
+    fprintf (f, "%6.1f %6.1f : %5.2f %5.2f %5.2f : %2d\n", 
+             x->data.F32[i], 
+             y->data.F32[i], 
+             e0->data.F32[i], 
+             e1->data.F32[i], 
+             e2->data.F32[i], 
+             mask->data.U8[i]);
+  }
+  fclose (f);
+  return true;
+}

tags/ipp-20141224/psModules/src/objects/pmTrend2D.c

@@ -211 +211 @@
         // XXX need to add the API which adjusts the scale
         status = psImageMapClipFit(pGoodFit, trend->map, trend->stats, mask, maskVal, x, y, f, df);
+        if (!status) {
+          psError(PS_ERR_UNKNOWN, true, "failed to build PSF model map");
+          return false;
+        }
+
+        // the psf model map can have nan pixels: repair this by extrapolation / interpolation
+        // XXX TEST: p_psImagePrint(0, trend->map->map, "before");
+        status = psImageMapRepair (trend->map->map);
+        // XXX TEST: p_psImagePrint(0, trend->map->map, "after");
+        if (!status) {
+          psError(PS_ERR_UNKNOWN, true, "failed to repair PSF model map");
+          return false;
+        }
         break;
 
@@ -217 +230 @@
     }
     return status;
+}
+
+bool psImageMapRepair (psImage *map) {
+
+
+  // XXX why is my repair not working??
+
+    // patch over bad regions (use average of 8 possible neighbor pixels)
+    // XXX consider testing all pixels against the 8 neighbors and replacing outliers...
+    double Count = 0;                   // number of good pixels
+    double Value = 0;                   // sum of good pixel's value
+    for (int iy = 0; iy < map->numRows; iy++) {
+        for (int ix = 0; ix < map->numCols; ix++) {
+            if (isfinite(map->data.F32[iy][ix])) {
+                Value += map->data.F32[iy][ix];
+                Count++;
+                continue;
+            }
+
+            double value = 0;
+            double count = 0;
+            for (int jy = iy - 1; jy <= iy + 1; jy++) {
+                if (jy <   0) continue;
+                if (jy >= map->numRows) continue;
+                for (int jx = ix - 1; jx <= ix + 1; jx++) {
+                    if (!jx && !jy) continue;
+                    if (jx   <   0) continue;
+                    if (jx   >= map->numCols) continue;
+                    if (!isfinite(map->data.F32[jy][jx])) continue;
+                    value += map->data.F32[jy][jx];
+                    count += 1.0;
+                }
+            }
+            if (count > 0) {
+              // psLogMsg ("psphot", PS_LOG_DETAIL, "patching image map %d, %d: %f (%d pts)\n", ix, iy, (value / count), (int) count);
+                map->data.F32[iy][ix] = value / count;
+            }
+        }
+    }
+    if (Count == 0) {
+        psError(PS_ERR_UNKNOWN, true, "failed to repair PSF model map");
+        return false;
+    }
+    Value /= Count;
+
+    // patch over remaining bad regions (use global average)
+    for (int iy = 0; iy < map->numRows; iy++) {
+        for (int ix = 0; ix < map->numCols; ix++) {
+            if (!isnan(map->data.F32[iy][ix])) continue;
+            map->data.F32[iy][ix] = Value;
+        }
+    }
+    return true;
 }