Changeset 28667 for trunk/psphot/src/psphotSignificanceImage.c

Timestamp:

Jul 13, 2010, 2:34:04 PM (16 years ago)

Author:

Paul Price

Message:

It's bothered me for a while that the noise model that we have been using is only kind of spatially variable. The variance map is intended to track the change in the noise properties as a function of position, and while it does keep track of different noise levels due to, e.g., different detectors, vignetting, etc., it has not tracked the different noise that comes from convolution kernels changing over an image. This had been an intentional oversight --- the prescription we were using for calculating the variance and covariance removed the spatial variation (by prescribing that the variance level under a convolution remain unchanged) and we were sweeping it under the carpet so we could calculate the covariance pseudo-matrix and use that for setting the absolute level of the noise.

I've now figured out how to have our cake and eat it too. By changing where the scalings are applied, the variance map can track the changing noise introduced by a varying convolution kernel, while the covariance pseudo-matrix tracks an absolute change. This fix is purely algorithmic --- it does not affect any of the variances/covariances stored on disk, nor is it affected by them. I hope that this means that any Magic recertification required is either minimal or completely unnecessary.

Again, all that has changed is where the scalings are applied.

File:

• trunk/psphot/src/psphotSignificanceImage.c (modified) (1 diff)

trunk/psphot/src/psphotSignificanceImage.c

@@ -76 +76 @@
     // Calculate correction factor for the covariance produced by the (potentially multiple) smoothing
     psKernel *kernel = psImageSmoothKernel(SIGMA_SMTH, NSIGMA_SMTH); // Kernel used for smoothing
-    float factor = 1.0 / psImageCovarianceCalculateFactor(kernel, readout->covariance);
+    double sum2 = 0.0;                                               // Sum of kernel squared
+    for (int y = kernel->yMin; y <= kernel->yMax; y++) {
+        for (int x = kernel->xMin; x <= kernel->xMax; x++) {
+            sum2 += PS_SQR(kernel->kernel[y][x]);
+        }
+    }
+    float factor = 1.0 / (sum2 * psImageCovarianceCalculateFactor(kernel, readout->covariance));
     psFree(kernel);