Changeset 41895 for trunk/psastro/src/psastroOneChipFit.c

Timestamp:

Nov 4, 2021, 6:08:21 PM (5 years ago)

Author:

eugene

Message:

merge changes from eam_branches/ipp-dev-20210817 (add CERSTD, NASTUSED to header, allow rejection no CERSTD, add psastroTIO)

Location:

trunk/psastro

Files:

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

trunk/psastro/src/psastroOneChipFit.c

@@ -40 +40 @@
     // allowed limits for valid solutions
     REQUIRED_RECIPE_VALUE (float maxError, "PSASTRO.MAX.ERROR", F32);
-    REQUIRED_RECIPE_VALUE (int minNstar, "PSASTRO.MIN.NSTAR", S32);
+    REQUIRED_RECIPE_VALUE (float maxStdev, "PSASTRO.MAX.STDEV", F32);
+    REQUIRED_RECIPE_VALUE (int   minNstar, "PSASTRO.MIN.NSTAR", S32);
 
     psArray *match = NULL;
@@ -121 +122 @@
         }
 
-        // create output toFPA; set masks appropriate to the Elixir DVO astrometry format
+        // Create output toFPA; set masks appropriate to the Elixir DVO astrometry format.
+        // Forward transformations (chip->fpa->tpa->sky) use Ordinary polynomials
         psFree (chip->toFPA);
-        chip->toFPA = psPlaneTransformAlloc (order, order);
+        chip->toFPA = psPlaneTransformAlloc (order, order, PS_POLYNOMIAL_ORD);
         for (int i = 0; i <= chip->toFPA->x->nX; i++) {
             for (int j = 0; j <= chip->toFPA->x->nY; j++) {
@@ -192 +194 @@
     float astError = 0.5*(rawXstdev + rawYstdev) * plateScale;
 
+    // astStdev is the average 1D stdev of median residuals in arcsec ('results' are in FPA units = microns)
+    // the median residuals are calculated in a grid of N x N bins 
+    float astStdev = 0.5*(results->dXstdev + results->dYstdev) * plateScale;
+
     // x and y are forced to use the same subset of values:
     int astNstar = results->yStats->clippedNvalues;
@@ -198 +204 @@
 
     // XXX should these result in errors or be handled another way?
-    psLogMsg ("psastro", PS_LOG_INFO, "astrometry solution: error: %f arcsec, Nstars: %d", astError, astNstar);
+    psLogMsg ("psastro", PS_LOG_INFO, "astrometry solution: error: %f arcsec, Nstars: %d, stdev: %f arcsec", astError, astNstar, astStdev);
     if (astError > maxError) {
         psLogMsg("psastro", PS_LOG_INFO, "residual error is too large, failed to find a solution: %f > %f", astError, maxError);
         validSolution = false;
     }
+    if (astStdev > maxStdev) {
+        psLogMsg("psastro", PS_LOG_INFO, "residual stdev is too large, failed to find a solution: %f > %f", astStdev, maxStdev);
+        validSolution = false;
+    }
     if (astNstar < minNstar) {
         psLogMsg("psastro", PS_LOG_INFO, "solution uses too few stars: %d < %d", astNstar, minNstar);
@@ -211 +221 @@
     psMetadataAddF32 (updates, PS_LIST_TAIL, "PERROR",   PS_META_REPLACE, "astrometry error (pixels)", pixError);
     psMetadataAddF32 (updates, PS_LIST_TAIL, "CERROR",   PS_META_REPLACE, "astrometry error (arcsec)", astError);
+    psMetadataAddF32 (updates, PS_LIST_TAIL, "CERSTD",   PS_META_REPLACE, "astrometry stdev (arcsec)", astStdev);
     if (validSolution) {
         psMetadataAddF32 (updates, PS_LIST_TAIL, "CPRECISE", PS_META_REPLACE, "astrometry precision (arcsec)", astError/sqrt(astNstar));
         psMetadataAddS32 (updates, PS_LIST_TAIL, "NASTRO",   PS_META_REPLACE, "number of astrometry stars", astNstar);
+        psMetadataAddS32 (updates, PS_LIST_TAIL, "NASTUSED", PS_META_REPLACE, "number of astrometry stars", astNstar);
     } else {
         psastroChipFailureHeader (updates);
+        psMetadataAddS32 (updates, PS_LIST_TAIL, "NASTUSED", PS_META_REPLACE, "number of astrometry stars", astNstar);
     }
     psMetadataAddF32 (updates, PS_LIST_TAIL, "EQUINOX",  PS_META_REPLACE, "equinox of ref catalog", 2000.0); // XXX this is bogus: should be defined based on equinox of refstars