Changeset 38114 for trunk/psModules/src/objects

Timestamp:

Apr 13, 2015, 5:59:19 PM (11 years ago)

Author:

eugene

Message:

adding code to repair the petrosian mag errors

File:

• trunk/psModules/src/objects/pmSourceIO_CMF.c.in (modified) (4 diffs)

trunk/psModules/src/objects/pmSourceIO_CMF.c.in

@@ -552 +552 @@
     psMetadataAddStr (outhead, PS_LIST_TAIL, "EXTNAME", PS_META_REPLACE, "xsrc table extension", extname);
 
+    // we can either write out the petrosian mag errors correctly or inverted (by default, psphot writes them out inverted)
+    bool repairMagErrors = psMetadataLookupBool (&status, recipe, "REPAIR.PETROSIAN.MAG.ERRORS");
+    if (repairMagErrors) {
+      psMetadataAddBool (outhead, PS_LIST_TAIL, "REPAIR.PETROSIAN.MAG.ERRORS", PS_META_REPLACE, "", true);
+    }
+
     pmChip *chip = readout->parent->parent;
     pmFPA  *fpa  = chip->parent;
@@ -651 +657 @@
                 // XXX note that this mag is either calibrated or instrumental depending on existence of zero point 
                 float mag = (extpars->petrosianFlux > 0.0) ? -2.5*log10(extpars->petrosianFlux) + magOffset : NAN; // XXX zero point
-                // XXX NOTE EAM 20140806 : PETRO_MAG_ERR is inverted!! oops!!
+                // NOTE EAM 20140806 : PETRO_MAG_ERR was inverted!! this allows for it to be repaired
                 float magErr = (extpars->petrosianFlux > 0.0) ? extpars->petrosianFlux / extpars->petrosianFluxErr : NAN; // XXX zero point
+                if (repairMagErrors) {
+                  // I need to add the kron error in quadrature becasue pet_error ignores the object flux
+                  float Krf  = source->moments ? source->moments->KronFlux : NAN;
+                  float dKrf = source->moments ? source->moments->KronFluxErr : NAN;
+                  if (isfinite (Krf) && isfinite (dKrf)) {
+                    magErr = sqrt(PS_SQR(1.0 / magErr) + PS_SQR(dKrf / Krf));
+                  } else {
+                    magErr = 1.0 / magErr;
+                  }
+                }
                 psMetadataAdd (row, PS_LIST_TAIL, "PETRO_MAG",        PS_DATA_F32, "Petrosian Magnitude", mag);
                 psMetadataAdd (row, PS_LIST_TAIL, "PETRO_MAG_ERR",    PS_DATA_F32, "Petrosian Magnitude Error", magErr);
@@ -836 +852 @@
     }
 
+    bool repairMagErrors = psMetadataLookupBool (&status, tableHeader, "REPAIR.PETROSIAN.MAG.ERRORS");
+
     for (long i = 0; i < numSources; i++) {
         psMetadata *row = psFitsReadTableRow(fits, i); // Table row
@@ -867 +885 @@
         extpars->axes.theta = psMetadataLookupF32(&status, row, "F25_THETA");
 
+        // magErr may have been saved in inverted form
         float mag = psMetadataLookupF32(&status, row, "PETRO_MAG");
         float magErr = psMetadataLookupF32(&status, row, "PETRO_MAG_ERR");
         if (isfinite(mag)) {
-            extpars->petrosianFlux    = pow(10., (magOffset - mag) / 2.5);
-            if (isfinite(magErr)) {
-                extpars->petrosianFluxErr = extpars->petrosianFlux / magErr;
-            }
-        }
+          extpars->petrosianFlux    = pow(10., (magOffset - mag) / 2.5);
+          if (isfinite(magErr)) {
+            if (repairMagErrors) {
+              // I need to add the kron error in quadrature becasue pet_error ignores the object flux
+              float Krf  = source->moments ? source->moments->KronFlux : NAN;
+              float dKrf = source->moments ? source->moments->KronFluxErr : NAN;
+              if (isfinite (Krf) && isfinite (dKrf)) {
+                magErr = 1.0 / sqrt(PS_SQR(magErr) - PS_SQR(dKrf / Krf));
+              } else {
+                magErr = 1.0 / magErr;
+              }
+              extpars->petrosianFluxErr = extpars->petrosianFlux * magErr;
+            } else {
+              extpars->petrosianFluxErr = extpars->petrosianFlux / magErr;
+            }
+          }
+        }
 
         extpars->petrosianRadius   = psMetadataLookupF32(&status, row, "PETRO_RADIUS");