Changeset 21149

Timestamp:

Jan 22, 2009, 9:03:46 AM (17 years ago)

Author:

Paul Price

Message:

Adding test for deconvolution: measures the fraction of the kernel flux within a range of apertures about the centre and takes the maximum; if significantly in excess of 1 (e.g., 1.1), the convolved images are generally bad.

Location:

trunk/psModules/src/imcombine

Files:

• pmSubtractionAnalysis.c (modified) (1 diff)
• pmSubtractionAnalysis.h (modified) (1 diff)
• pmSubtractionMatch.c (modified) (2 diffs)

trunk/psModules/src/imcombine/pmSubtractionAnalysis.c

@@ -165 +165 @@
     UPDATE_VARFACTOR(vf2, PM_SUBTRACTION_ANALYSIS_VARFACTOR_2);
 
+    // Kernel shape
+    {
+        psImage *image = pmSubtractionKernelImage(kernels, 0.5, 0.5, false); // Image of the kernel
+        if (!image) {
+            psError(PS_ERR_UNKNOWN, false, "Unable to generate image of kernel.");
+            return false;
+        }
+        int size = kernels->size;       // Half-size of kernel
+        int fullSize = 2 * size + 1;    // Full size of kernel
+        float norm = 0.0;               // Normalisation (kernel sum)
+        for (int y = 0; y < fullSize; y++) {
+            for (int x = 0; x < fullSize; x++) {
+                norm += image->data.F32[y][x];
+            }
+        }
+        float max = -INFINITY;          // Maximum fraction
+        for (int r = 1; r < size; r++) {
+            unsigned long r2 = PS_SQR(r); // r^2
+            float sum = 0.0;            // Sum within circle
+            for (int y = 0, v = -size; y < fullSize; y++, v++) {
+                unsigned long v2 = PS_SQR(v); // y^2
+                for (int x = 0, u = -size; x < fullSize; x++, u++) {
+                    unsigned long u2 = PS_SQR(u); // u^2
+                    if (u2 + v2 <= r2) {
+                        sum += image->data.F32[y][x];
+                    }
+                }
+            }
+            float frac = sum / norm;    // Fraction of flux moving towards centre
+            psTrace("psModules.imcombine", 5, "Deconvolution fraction at %d: %f\n", r, frac);
+            max = PS_MAX(max, frac);
+        }
+        psFree(image);
+
+        psMetadataAddF32(analysis, PS_LIST_TAIL, PM_SUBTRACTION_ANALYSIS_DECONV_MAX,
+                         PS_META_DUPLICATE_OK, "Maximum deconvolution fraction", max);
+    }
+
     // Kernel moments
     {

trunk/psModules/src/imcombine/pmSubtractionAnalysis.h

@@ -22 +22 @@
 #define PM_SUBTRACTION_ANALYSIS_MXY          "SUBTRACTION.MXY"          // Kernel moment in xy
 #define PM_SUBTRACTION_ANALYSIS_MYY          "SUBTRACTION.MYY"          // Kernel moment in yy
+#define PM_SUBTRACTION_ANALYSIS_DECONV_MAX   "SUBTRACTION.DECONV.MAX"   // Maximum deconvolution fraction
 
 // Derive QA information about the subtraction

trunk/psModules/src/imcombine/pmSubtractionMatch.c

@@ -569 +569 @@
         ratios->data.F32[index] = (float)width1 / (float)width2;
         mask->data.PS_TYPE_MASK_DATA[index] = 0;
+        psTrace("psModules.imcombine", 3, "Stamp %d (%.1f,%.1f) widths: %d, %d --> %f\n",
+                index, stamp->x, stamp->y, width1, width2, ratios->data.F32[index]);
     }
 
@@ -680 +682 @@
     psFree(stats);
 
+
+    return PM_SUBTRACTION_MODE_2;
+
+
+    exit(1);
+
     return mode;
 }