Changeset 28017 for branches/czw_branch/20100427/pswarp/src/pswarpLoop.c

Timestamp:

May 18, 2010, 5:06:01 PM (16 years ago)

Author:

watersc1

Message:

Almost working version of the mask-stats/software revision code.

Almost.

Location:

branches/czw_branch/20100427

Files:

• . (modified) (1 prop)
• pswarp/src/pswarpLoop.c (modified) (5 diffs)

branches/czw_branch/20100427/pswarp/src/pswarpLoop.c

@@ -266 +266 @@
         goto DONE;
     }
-
-    // Correct image for change in the plate scale
-    {
-        psAssert(input && input->fpa && input->fpa->toSky, "Require astrometry for input");
-        psAssert(outFPA && outFPA && outFPA->toSky, "Require astrometry for output");
-
-        double inScale = input->fpa->toSky->Xs + input->fpa->toSky->Ys; // Plate scale for input
-        double outScale = outFPA->toSky->Xs + outFPA->toSky->Ys; // Plate scale for output
-        float correction = PS_SQR(outScale / inScale); // Correction factor to apply to image
-        psLogMsg("pswarp", PS_LOG_INFO, "Correcting flux by %f to account for pixel scales", correction);
-
-        psBinaryOp(output->image, output->image, "*", psScalarAlloc(correction, PS_TYPE_F32));
-        if (output->variance) {
-            psBinaryOp(output->variance, output->variance, "*",
-                       psScalarAlloc(PS_SQR(correction), PS_TYPE_F32));
-        }
-    }
-
-
+    bool doStats = psMetadataLookupBool(&mdok,recipe,"MASK.STATS");
+    if (doStats) {
+      if (!pswarpMaskStats(output, stats, config)) {
+        psError(psErrorCodeLast(), false, "Unable to calculate mask stats.");
+        psFree(cells);
+        psFree(view);
+        goto DONE;
+      }
+    }
     // Set covariance matrix for output
     {
@@ -291 +281 @@
         psAssert(covariances, "Should be there");
         psArray *covars = psListToArray(covariances); // Array of covariance matrices
-        output->covariance = psImageCovarianceAverage(covars);
+        psKernel *covar = psImageCovarianceAverage(covars);
         psFree(covars);
         psMetadataRemoveKey(output->analysis, PSWARP_ANALYSIS_COVARIANCES);
+
+        // Correct covariance matrix scale for the mean (square root of the) Jacobian
+        double jacobian = psMetadataLookupF64(NULL, output->analysis, PSWARP_ANALYSIS_JACOBIAN); // Jacobian
+        int goodPixels = psMetadataLookupS32(NULL, output->analysis, PSWARP_ANALYSIS_GOODPIX);   // Good pixels
+        jacobian /= goodPixels;
+        output->covariance = psImageCovarianceScale(covar, jacobian);
+        psFree(covar);
+
         psImageCovarianceTransfer(output->variance, output->covariance);
     }
@@ -340 +338 @@
 
     pswarpVersionHeader(hdu->header);
-
+    
     if (!pmAstromWriteWCS(hdu->header, outFPA, outChip, WCS_NONLIN_TOL)) {
         psError(psErrorCodeLast(), false, "Unable to generate WCS header.");
@@ -355 +353 @@
     pswarpFileActivation(config, detectorFiles, false);
     pswarpFileActivation(config, independentFiles, false);
+
 
     // We need a new PSF model for the warped frame.  It would be good to generate this analytically, but
@@ -429 +428 @@
         }
     }
+    
 
     // Add MD5 information for readout