Changeset 30784 for branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourcePhotometry.c

Timestamp:

Mar 3, 2011, 3:11:09 PM (15 years ago)

Author:

eugene

Message:

fix the calculation of chisq values for constant errors (re-calculated, do not attempt to correct raw value)

File:

• branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourcePhotometry.c (modified) (3 diffs)

branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourcePhotometry.c

@@ -741 +741 @@
 # endif
 
-// determine chisq, etc for linear normalization-only fit
-bool pmSourceChisq (pmModel *model, psImage *image, psImage *mask, psImage *variance, psImageMaskType maskVal, const float covarFactor, int nParams)
+// determine chisq, nPix, nDOF, chisqNorm : model->nPar must be set
+bool pmSourceChisq (pmModel *model, psImage *image, psImage *mask, psImage *variance, psImageMaskType maskVal)
 {
     PS_ASSERT_PTR_NON_NULL(model, false);
@@ -757 +757 @@
             if (variance->data.F32[j][i] <= 0)
                 continue;
-            // dC += PS_SQR (image->data.F32[j][i]) / (covarFactor * variance->data.F32[j][i]);
             dC += PS_SQR (image->data.F32[j][i]) / variance->data.F32[j][i];
             Npix ++;
         }
     }
-
     model->nPix = Npix;
-    model->nDOF = Npix - nParams - 1;
+    model->nDOF = Npix - model->nPar;
     model->chisq = dC;
     model->chisqNorm = dC / model->nDOF;
@@ -771 +769 @@
 }
 
+
+// return source aperture magnitude
+bool pmSourceChisqUnsubtracted (pmSource *source, pmModel *model, psImageMaskType maskVal)
+{
+    PS_ASSERT_PTR_NON_NULL(source, false);
+    PS_ASSERT_PTR_NON_NULL(model, false);
+
+    float dC = 0.0;
+    int Npix = 0;
+
+    // the model function returns the source flux at a position
+    psVector *coord = psVectorAlloc(2, PS_TYPE_F32);
+
+    psVector *params = model->params;
+    psImage  *image = source->pixels;
+    psImage  *mask = source->maskObj;
+    psImage  *variance = source->variance;
+
+    int dX = image->col0;
+    int dY = image->row0;
+
+    for (int iy = 0; iy < image->numRows; iy++) {
+        for (int ix = 0; ix < image->numCols; ix++) {
+
+            // skip pixels which are masked
+            if (mask->data.PS_TYPE_IMAGE_MASK_DATA[iy][ix] & maskVal) continue;
+
+            if (variance->data.F32[iy][ix] <= 0) continue;
+
+            coord->data.F32[0] = (psF32) ix + dX + 0.5;
+            coord->data.F32[1] = (psF32) iy + dY + 0.5;
+
+            // for the full model, add all points
+            float value = model->modelFunc (NULL, params, coord);
+
+            // fprintf (stderr, "%d, %d : %f, %f : %f - %f : %f\n", 
+            // ix, iy, coord->data.F32[0], coord->data.F32[1], image->data.F32[iy][ix], value, dC);
+
+            dC += PS_SQR (image->data.F32[iy][ix] - value) / variance->data.F32[iy][ix];
+            Npix ++;
+        }
+    }
+    model->nPix = Npix;
+    model->nDOF = Npix - model->nPar;
+    model->chisq = dC;
+    model->chisqNorm = dC / model->nDOF;
+
+    psFree (coord);
+    return (true);
+}
 
 double pmSourceModelWeight(const pmSource *Mi, int term, const bool unweighted_sum, const float covarFactor, psImageMaskType maskVal)