Changeset 18289

Timestamp:

Jun 23, 2008, 12:43:11 PM (18 years ago)

Author:

Paul Price

Message:

Merging in dual convolution development branch.

Location:

trunk/ppSub

Files:

• src/ppSubArguments.c (modified) (2 diffs)
• src/ppSubCamera.c (modified) (5 diffs)
• src/ppSubReadout.c (modified) (8 diffs)
• test/fake.c (added)

trunk/ppSub/src/ppSubArguments.c

@@ -186 +186 @@
     psMetadataAddS32(arguments, PS_LIST_TAIL, "-order", 0, "Spatial polynomial order", -1);
     psMetadataAddStr(arguments, PS_LIST_TAIL, "-type", 0, "Kernel type (ISIS|POIS|SPAM|FRIES|GUNK|RINGS)", NULL);
+    psMetadataAddF32(arguments, PS_LIST_TAIL, "-penalty", 0, "Penalty for wideness", NAN);
     psMetadataAddStr(arguments, PS_LIST_TAIL, "-isis-widths", 0, "ISIS Gaussian FWHMs (comma-separated)", NULL);
     psMetadataAddStr(arguments, PS_LIST_TAIL, "-isis-orders", 0, "ISIS polynomial orders (comma-separated)", NULL);
@@ -260 +261 @@
     VALUE_ARG_RECIPE_FLOAT("-rej",        "REJ",             F32);
     VALUE_ARG_RECIPE_FLOAT("-badfrac",    "BADFRAC",         F32);
+    VALUE_ARG_RECIPE_FLOAT("-penalty",    "PENALTY",         F32);
 
     valueArgRecipeStr(arguments, recipe, "-mask-bad",   "MASK.BAD",   config->arguments);

trunk/ppSub/src/ppSubCamera.c

@@ -82 +82 @@
 
     // Output image
-    pmFPAfile *output = pmFPAfileDefineSkycell(config, NULL, "PPSUB.OUTPUT");
+    pmFPAfile *output = pmFPAfileDefineFromFile(config, input, 1, 1, "PPSUB.OUTPUT");
     if (!output) {
         psError(PS_ERR_IO, false, _("Unable to generate output file from PPSUB.OUTPUT"));
@@ -94 +94 @@
 
     // Output mask
-    pmFPAfile *outMask = pmFPAfileDefineSkycell(config, output->fpa, "PPSUB.OUTPUT.MASK");
+    pmFPAfile *outMask = pmFPAfileDefineOutput(config, output->fpa, "PPSUB.OUTPUT.MASK");
     if (!outMask) {
         psError(PS_ERR_IO, false, _("Unable to generate output file from PPSUB.OUTPUT.MASK"));
@@ -107 +107 @@
     // Output weight
     if (inputWeight && refWeight) {
-        pmFPAfile *outWeight = pmFPAfileDefineSkycell(config, output->fpa, "PPSUB.OUTPUT.WEIGHT");
+        pmFPAfile *outWeight = pmFPAfileDefineOutput(config, output->fpa, "PPSUB.OUTPUT.WEIGHT");
         if (!outWeight) {
             psError(PS_ERR_IO, false, _("Unable to generate output file from PPSUB.OUTPUT.WEIGHT"));
@@ -119 +119 @@
     }
 
+#if 0
     if (!pmFPAAddSourceFromFormat(output->fpa, "Subtraction", output->format)) {
         psError(PS_ERR_UNKNOWN, false, "Unable to generate output FPA.");
         return false;
     }
+#endif
 
     pmFPAfile *sources = pmFPAfileDefineFromArgs(&status, config, "PPSUB.SOURCES", "PPSUB.SOURCES");
@@ -153 +155 @@
 
         // Internal-ish file for getting the PSF from the matched addition
-        pmFPAfile *psf = pmFPAfileDefineSkycell(config, psphot->fpa, "PSPHOT.PSF.LOAD");
+        pmFPAfile *psf = pmFPAfileDefineOutputFromFile(config, output, "PSPHOT.PSF.LOAD");
         if (!psf) {
             psError(PS_ERR_IO, false, "Failed to build FPA from PSPHOT.PSF.LOAD");

trunk/ppSub/src/ppSubReadout.c

@@ -12 +12 @@
 
 #define WCS_TOLERANCE 0.001             // Tolerance for WCS
-//#define TESTING                         // For test output
+#define TESTING                         // For test output
 
 
@@ -73 +73 @@
     int ringsOrder = psMetadataLookupS32(NULL, config->arguments, "RINGS.ORDER"); // RINGS polynomial order
     int binning = psMetadataLookupS32(NULL, config->arguments, "SPAM.BINNING"); // Binning for SPAM kernel
+    float penalty = psMetadataLookupF32(NULL, config->arguments, "PENALTY"); // Penalty for wideness
     psMaskType maskBad = pmConfigMask(psMetadataLookupStr(NULL, config->arguments, "MASK.BAD"),
                                     config); // Value to mask
@@ -131 +132 @@
 
     if (!pmSubtractionMatch(inConv, refConv, inRO, refRO, footprint, regionSize, spacing, threshold, sources,
-                            stampsName, type, size, order, widths, orders, inner, ringsOrder,
-                            binning, optimum, optWidths, optOrder, optThresh, iter, rej, maskBad,
+                            stampsName, type, size, order, widths, orders, inner, ringsOrder, binning,
+                            penalty, optimum, optWidths, optOrder, optThresh, iter, rej, maskBad,
                             maskBlank, badFrac, mode)) {
         psError(PS_ERR_UNKNOWN, false, "Unable to match images.");
@@ -184 +185 @@
 #endif
 
-    // "Subtract" the mask and weight map
+    // Subtraction is: minuend - subtrahend
+    pmReadout *minuend = inConv;
+    pmReadout *subtrahend = refConv;
+
+    if (reverse) {
+        pmReadout *temp = subtrahend;
+        subtrahend = minuend;
+        minuend = temp;
+    }
+
+#ifdef TESTING
+    {
+        pmReadoutMaskApply(minuend, maskBlank);
+        psFits *fits = psFitsOpen("minuend.fits", "w");
+        psFitsWriteImage(fits, NULL, minuend->image, 0, NULL);
+        psFitsClose(fits);
+    }
+    {
+        pmReadoutMaskApply(subtrahend, maskBlank);
+        psFits *fits = psFitsOpen("subtrahend.fits", "w");
+        psFitsWriteImage(fits, NULL, subtrahend->image, 0, NULL);
+            psFitsClose(fits);
+    }
+#endif
+
     outRO->mask = (psImage*)psBinaryOp(outRO->mask, inConv->mask, "|", refConv->mask);
-    if (inConv->weight && refConv->weight) {
-        outRO->weight = (psImage*)psBinaryOp(outRO->weight, inConv->weight, "+", refConv->weight);
-    }
-    outRO->data_exists = outCell->data_exists = outCell->parent->data_exists = true;
+    outRO->data_exists = outCell->data_exists = outCell->parent->data_exists = true; // It'll be there soon
 
     // Photometry is to be performed in two stages:
@@ -198 +220 @@
     pmPSF *psf = NULL;                  // PSF for photometry
     if (psMetadataLookupBool(NULL, config->arguments, "PHOTOMETRY")) {
-        // We use a summed image as the basis for the PSF: this will have the maximum S/N.
-        outRO->image = (psImage*)psBinaryOp(outRO->image, inConv->image, "+", refConv->image);
+        outRO->image = psImageCopy(outRO->image, minuend->image, PS_TYPE_F32);
+        if (minuend->weight) {
+            outRO->weight = psImageCopy(outRO->weight, minuend->weight, PS_TYPE_F32);
+        }
 
         pmFPAfile *photFile = psMetadataLookupPtr(NULL, config->files, "PSPHOT.INPUT");
@@ -244 +268 @@
     }
 
-    // Do the subtraction
-    {
-        // Subtraction is: minuend - subtrahend
-        pmReadout *minuend = inConv;
-        pmReadout *subtrahend = refConv;
-
-        if (reverse) {
-            pmReadout *temp = subtrahend;
-            subtrahend = minuend;
-            minuend = temp;
-        }
-
-        outRO->image = (psImage*)psBinaryOp(outRO->image, minuend->image, "-", subtrahend->image);
-
-#ifdef TESTING
-        {
-            psFits *fits = psFitsOpen("minuend.fits", "w");
-            psFitsWriteImage(fits, NULL, minuend->image, 0, NULL);
-            psFitsClose(fits);
-        }
-        {
-            psFits *fits = psFitsOpen("subtrahend.fits", "w");
-            psFitsWriteImage(fits, NULL, subtrahend->image, 0, NULL);
-            psFitsClose(fits);
-        }
-#endif
-
-        pmReadoutMaskApply(outRO, maskBlank);
-    }
+    // Do the actual subtraction
+    outRO->image = (psImage*)psBinaryOp(outRO->image, minuend->image, "-", subtrahend->image);
+    if (inConv->weight && refConv->weight) {
+        outRO->weight = (psImage*)psBinaryOp(outRO->weight, inConv->weight, "+", refConv->weight);
+    }
+    outRO->data_exists = outCell->data_exists = outCell->parent->data_exists = true;
+
+    pmReadoutMaskApply(outRO, maskBlank);
 
     psFree(inConv);
@@ -304 +307 @@
     }
 
-    psTraceSetLevel("psphot", 6);
-
     // Photometry stage 2: find and measure sources on the subtracted image
     if (psMetadataLookupBool(NULL, config->arguments, "PHOTOMETRY")) {
@@ -314 +315 @@
         pmReadout *psfRO = pmFPAfileThisReadout(config->files, view, "PSPHOT.PSF.LOAD");
         if (!psfRO) {
-            psError(PS_ERR_UNEXPECTED_NULL, true, "Unable to find file PSPHOT.PSF.LOAD");
+            psError(PS_ERR_UNEXPECTED_NULL, true, "Unable to find readout for file PSPHOT.PSF.LOAD");
             return false;
         }