Changeset 9626 for trunk/psModules/src/astrom/pmAstrometryObjects.c

Timestamp:

Oct 18, 2006, 8:40:19 AM (20 years ago)

Author:

magnier

Message:

using a single mask for X and Y, iterate between X and Y for clip fitting in pmAstromMatchFit

File:

• trunk/psModules/src/astrom/pmAstrometryObjects.c (modified) (5 diffs)

trunk/psModules/src/astrom/pmAstrometryObjects.c

@@ -8 +8 @@
 *  @author EAM, IfA
 *
-*  @version $Revision: 1.15 $ $Name: not supported by cvs2svn $
-*  @date $Date: 2006-10-07 03:54:08 $
+*  @version $Revision: 1.16 $ $Name: not supported by cvs2svn $
+*  @date $Date: 2006-10-18 18:40:19 $
 *
 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
@@ -227 +227 @@
     }
 
-    psStats *stats = psStatsAlloc (PS_STAT_CLIPPED_MEAN | PS_STAT_CLIPPED_STDEV);
-    stats->clipSigma = psMetadataLookupF32 (&status, config, "PSASTRO.CHIP.NSIGMA");
-    stats->clipIter  = psMetadataLookupS32 (&status, config, "PSASTRO.CHIP.NITER");
-
     // XXX EAM : clip fit seems to only work for F32!
     // XXX EAM : clip fit fails to handle NULL error
@@ -260 +256 @@
 
     // constant errors
-    psVector *xMask = psVectorAlloc (match->n, PS_TYPE_U8);
-    psVector *yMask = psVectorAlloc (match->n, PS_TYPE_U8);
-    xMask->n = yMask->n = match->n;
-    psVectorInit (xMask, 0);
-    psVectorInit (yMask, 0);
+    psVector *mask = psVectorAlloc (match->n, PS_TYPE_U8);
+    mask->n = match->n;
+    psVectorInit (mask, 0);
 
     double dX, dY;
+
+    psStats *stats = psStatsAlloc (PS_STAT_CLIPPED_MEAN | PS_STAT_CLIPPED_STDEV);
+    stats->clipSigma = psMetadataLookupF32 (&status, config, "PSASTRO.CHIP.NSIGMA");
+    stats->clipIter = 1;
+    int clipIter  = psMetadataLookupS32 (&status, config, "PSASTRO.CHIP.NITER");
+    // we run 'clipIter' cycles clipping in each of x and y
 
     // fit chip-to-FPA transformation
     // XXX should force x->n and y->n to be the same
-    psVectorClipFitPolynomial2D (map->x, stats, xMask, 0, x, wt, X, Y);
-    psLogMsg ("psModules.astrom", 3, "x resid: %f +/- %f (%ld of %ld)\n", stats->clippedMean, stats->clippedStdev, stats->clippedNvalues, x->n);
-    dX = stats->clippedStdev;
-
-    psVectorClipFitPolynomial2D (map->y, stats, yMask, 0, y, wt, X, Y);
-    psLogMsg ("psModules.astrom", 3, "y resid: %f +/- %f (%ld of %ld)\n", stats->clippedMean, stats->clippedStdev, stats->clippedNvalues, y->n);
-    dY = stats->clippedStdev;
+    for (int i = 0; i < clipIter; i++) {
+        psVectorClipFitPolynomial2D (map->x, stats, mask, 0xff, x, wt, X, Y);
+        psLogMsg ("psModules.astrom", 3, "x resid: %f +/- %f (%ld of %ld)\n", stats->clippedMean, stats->clippedStdev, stats->clippedNvalues, x->n);
+        dX = stats->clippedStdev;
+
+        psVectorClipFitPolynomial2D (map->y, stats, mask, 0xff, y, wt, X, Y);
+        psLogMsg ("psModules.astrom", 3, "y resid: %f +/- %f (%ld of %ld)\n", stats->clippedMean, stats->clippedStdev, stats->clippedNvalues, y->n);
+        dY = stats->clippedStdev;
+    }
+
+    // XXXX test dump output code
+    # if (0)
+        psVector *xfit = psPolynomial2DEvalVector(map->x, X, Y);
+    psVector *yfit = psPolynomial2DEvalVector(map->y, X, Y);
+    FILE *f = fopen ("match.dat", "w");
+    for (int i = 0; i < x->n; i++) {
+        fprintf (f, "%3d %1d %1d %f %f  %f %f   %f %f\n",
+                 i, mask->data.U8[i], mask->data.U8[i],
+                 X->data.F32[i], Y->data.F32[i],
+                 xfit->data.F32[i], yfit->data.F32[i],
+                 x->data.F32[i], y->data.F32[i]);
+    }
+    fclose (f);
+    # endif
 
     // XXX this is somewhat ad-hoc: we need to be sure these values are in arcsec
@@ -294 +311 @@
     psFree (wt);
     psFree (stats);
-    psFree (xMask);
-    psFree (yMask);
+    psFree (mask);
 
     return (map);
@@ -402 +418 @@
     *obj->TP   = *old->TP;
     *obj->sky  = *old->sky;
+    obj->Mag   =  old->Mag;
+    obj->dMag  =  old->dMag;
 
     return(obj);