Changeset 31154 for trunk/psphot/src/psphotExtendedSourceAnalysis.c

Timestamp:

Apr 4, 2011, 1:12:39 PM (15 years ago)

Author:

eugene

Message:

use the smoothed image for footprint culling; use a more stringent culling for saturated stars; more tweaks to get sat stars to be fitted; various updates to psphotStack; unify psphotImageLoop varients; be a bit careful about number of stars in a PSF clump

Location:

trunk/psphot

Files:

• . (modified) (2 props)
• src/psphotExtendedSourceAnalysis.c (modified) (9 diffs)

trunk/psphot

@@ -19 +19 @@
 psphot-config
 Doxyfile
+a.out.dSYM

@@ -19 +19 @@
 psphot-config
 Doxyfile
+a.out.dSYM

trunk/psphot/src/psphotExtendedSourceAnalysis.c

@@ -1 +1 @@
 # include "psphotInternal.h"
 
-// ?? these cannot happen here --> we would need to do this in psphotExtendedSourceAnalysis
-// XXX option to choose a consistent position
-// XXX option to choose a consistent elliptical contour
-// XXX SDSS uses the r-band petrosian radius to measure petrosian fluxes in all bands
-// XXX consistent choice of extendedness...
+// measure the elliptical radial profile and use this to measure the petrosian parameters for the sources
+// XXX this function needs to be threaded
 
 // for now, let's store the detections on the readout->analysis for each readout
@@ -40 +37 @@
     int Next = 0;
     int Npetro = 0;
-    int Nisophot = 0;
     int Nannuli = 0;
-    int Nkron = 0;
 
     psTimerStart ("psphot.extended");
@@ -68 +63 @@
     assert (maskVal);
 
-    // XXX require petrosian analysis for non-linear fits? 
-
-    // XXX temporary user-supplied systematic sky noise measurement (derive from background model)
-    float skynoise = psMetadataLookupF32 (&status, recipe, "SKY.NOISE");
+    // get the sky noise from the background analysis; if this is missing, get the user-supplied value
+    float skynoise = psMetadataLookupF32 (&status, readout->analysis, "SKY_STDEV");
+    if (!status) {
+        skynoise = psMetadataLookupF32 (&status, recipe, "SKY.NOISE");
+        psWarning ("failed to get sky noise level from background analysis; defaulting to user supplied value of %f\n", skynoise);
+    }
 
     // S/N limit to perform full non-linear fits
@@ -78 +75 @@
     // which extended source analyses should we perform?
     bool doPetrosian    = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_PETROSIAN");
-    bool doIsophotal    = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_ISOPHOTAL");
     bool doAnnuli       = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_ANNULI");
-    bool doKron         = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_KRON");
-
-# if (0)
-    // if backModel or backStdev are missing, the values of sky and/or skyErr will be set to NAN
-    // XXX use this to set skynoise
-    pmReadout *backModel = psphotSelectBackground (config, view);
-    pmReadout *backStdev = psphotSelectBackgroundStdev (config, view);
-# endif
+    bool doPetroStars   = psMetadataLookupBool (&status, recipe, "PETROSIAN_FOR_STARS");
 
     // source analysis is done in S/N order (brightest first)
-    sources = psArraySort (sources, pmSourceSortBySN);
+    sources = psArraySort (sources, pmSourceSortByFlux);
 
     // option to limit analysis to a specific region
@@ -103 +92 @@
 
         // skip PSF-like and non-astronomical objects
-        if (source->type == PM_SOURCE_TYPE_STAR) continue;
         if (source->type == PM_SOURCE_TYPE_DEFECT) continue;
         if (source->type == PM_SOURCE_TYPE_SATURATED) continue;
         if (source->mode & PM_SOURCE_MODE_DEFECT) continue;
         if (source->mode & PM_SOURCE_MODE_SATSTAR) continue;
-        if (!(source->mode & PM_SOURCE_MODE_EXT_LIMIT)) continue;
+
+        // optionally allow non-extended objects to get petrosians as well
+        if (!doPetroStars) {
+            if (!(source->mode & PM_SOURCE_MODE_EXT_LIMIT)) continue;
+            if (source->type == PM_SOURCE_TYPE_STAR) continue;
+        }
 
         // limit selection to some SN limit
         assert (source->peak); // how can a source not have a peak?
-        if (source->peak->SN < SN_LIM) continue;
+        if (sqrt(source->peak->detValue) < SN_LIM) continue;
 
         // limit selection by analysis region
@@ -119 +112 @@
         if (source->peak->x > AnalysisRegion.x1) continue;
         if (source->peak->y > AnalysisRegion.y1) continue;
-
-        // fprintf (stderr, "xsrc: %f, %f : %f\n", source->peak->xf, source->peak->yf, source->peak->SN); 
 
         // replace object in image
@@ -136 +127 @@
 
         // if we request any of these measurements, we require the radial profile
-        if (doPetrosian || doIsophotal || doAnnuli || doKron) {
+        if (doPetrosian || doAnnuli) {
             if (!psphotRadialProfile (source, recipe, skynoise, maskVal)) {
                 // all measurements below require the radial profile; skip them all
@@ -144 +135 @@
                 continue;
             }
+            Nannuli ++;
             source->mode |= PM_SOURCE_MODE_RADIAL_FLUX;
         }
@@ -169 +161 @@
     psLogMsg ("psphot", PS_LOG_INFO, "extended source analysis: %f sec for %d objects\n", psTimerMark ("psphot.extended"), Next);
     psLogMsg ("psphot", PS_LOG_INFO, "  %d petrosian\n", Npetro);
-    psLogMsg ("psphot", PS_LOG_INFO, "  %d isophotal\n", Nisophot);
     psLogMsg ("psphot", PS_LOG_INFO, "  %d annuli\n", Nannuli);
-    psLogMsg ("psphot", PS_LOG_INFO, "  %d kron\n", Nkron);
 
     psphotVisualShowResidualImage (readout, false);

@@ -19 +19 @@
 psphot-config
 Doxyfile
+a.out.dSYM