Changeset 26703 for branches/eam_branches/20091201/ppSub/src/ppSubMatchPSFs.c

Timestamp:

Jan 28, 2010, 4:16:19 AM (16 years ago)

Author:

Paul Price

Message:

Solve for the normalisation of the convolution kernel completely separately from that of the other kernel terms. We use a 'normalisation window', basically an aperture in which we measure the flux of the stamps to determine the normalisation value. This introduces a new recipe parameter, NORM.FRAC, which is the fraction of the total in the window at which we set the normalisation window/aperture. The quality of the result seems to be rather sensitive to this value. Also found a bug in the imposition of the penalty function in dual convolution mode: it wasn't being set for the second solution part. This seems to produce decent dual convolution subtractions.

File:

• branches/eam_branches/20091201/ppSub/src/ppSubMatchPSFs.c (modified) (9 diffs)

branches/eam_branches/20091201/ppSub/src/ppSubMatchPSFs.c

@@ -22 +22 @@
 #include "ppSub.h"
 
+#define COVAR_FRAC 0.01                 // Truncation fraction for covariance matrix
+
 // Normalise a region on an image
 static void normaliseRegion(psImage *image, // Image to normalise
@@ -38 +40 @@
 }
 
+#if 0
 // Measure the PSF for an image
 static float subImagePSF(ppSubData *data, // Processing data
@@ -85 +88 @@
     return fwhm;
 }
+#endif
 
 // Scale the kernel parameters according to the PSFs
@@ -114 +118 @@
     }
 
+#if 0
     // Input images
     pmReadout *inRO = pmFPAfileThisReadout(config->files, view, "PPSUB.INPUT"); // Input readout
@@ -133 +138 @@
     float inFWHM = subImagePSF(data, inRO, inSources); // FWHM for input
     float refFWHM = subImagePSF(data, refRO, refSources); // FWHM for reference
+#else
+    float inFWHM = 7.631371;
+    float refFWHM = 10.005879;
+#endif
+
     psLogMsg("ppSub", PS_LOG_INFO, "Input FWHM: %f\nReference FWHM: %f\n", inFWHM, refFWHM);
     if (!isfinite(inFWHM) || !isfinite(refFWHM)) {
@@ -252 +262 @@
     float rej = psMetadataLookupF32(NULL, recipe, "REJ"); // Rejection threshold
     float kernelErr = psMetadataLookupF32(NULL, recipe, "KERNEL.ERR"); // Relative systematic error in kernel
+    float normFrac = psMetadataLookupF32(NULL, recipe, "NORM.FRAC"); // Fraction of window for normalisn windw
     float sysErr = psMetadataLookupF32(NULL, recipe, "SYS.ERR"); // Relative systematic error in images
     float skyErr = psMetadataLookupF32(NULL, recipe, "SKY.ERR"); // Additional error in sky
@@ -306 +317 @@
     }
 
+    if (inRO->covariance) {
+        psImageCovarianceTruncate(inRO->covariance, COVAR_FRAC);
+    }
+    if (refRO->covariance) {
+        psImageCovarianceTruncate(refRO->covariance, COVAR_FRAC);
+    }
+
     // Match the PSFs
     bool success = false;               // Operation was successful?
@@ -315 +333 @@
                                      spacing, threshold, sources, data->stamps, type, size, order,
                                      widths, orders, inner, ringsOrder, binning, penalty, optimum,
-                                     optWidths, optOrder, optThresh, iter, rej, sysErr, skyErr, kernelErr, maskVal,
-                                     maskBad, maskPoor, poorFrac, badFrac, subMode);
+                                     optWidths, optOrder, optThresh, iter, rej, normFrac,
+                                     sysErr, skyErr, kernelErr, maskVal, maskBad, maskPoor,
+                                     poorFrac, badFrac, subMode);
     }
 
@@ -407 +426 @@
     pmConceptsCopyFPA(refConv->parent->parent->parent, refRO->parent->parent->parent, true, true);
 
+    if (inConv->covariance) {
+        psImageCovarianceTruncate(inConv->covariance, COVAR_FRAC);
+    }
+    if (refConv->covariance) {
+        psImageCovarianceTruncate(refConv->covariance, COVAR_FRAC);
+    }
+
     if (inConv->variance) {
         psImageCovarianceTransfer(inConv->variance, inConv->covariance);