Changeset 27840 for branches/simtest_nebulous_branches/psphot/src/psphotRadiusChecks.c

Timestamp:

May 3, 2010, 8:50:52 AM (16 years ago)

Author:

eugene

Message:

updates from trunk

Location:

branches/simtest_nebulous_branches

Files:

• . (modified) (1 prop)
• psphot (modified) (1 prop)
• psphot/src (modified) (1 prop)
• psphot/src/psphotRadiusChecks.c (modified) (7 diffs)

branches/simtest_nebulous_branches/psphot/src

@@ -18 +18 @@
 psphotVersionDefinitions.h
 psphotMomentsStudy
+psphotPetrosianStudy
+psphotForced
+psphotMakePSF
+psphotStack

@@ -18 +18 @@
 psphotVersionDefinitions.h
 psphotMomentsStudy
+psphotPetrosianStudy
+psphotForced
+psphotMakePSF
+psphotStack

branches/simtest_nebulous_branches/psphot/src/psphotRadiusChecks.c

@@ -4 +4 @@
 static float PSF_FIT_NSIGMA;
 static float PSF_FIT_PADDING;
+static float PSF_APERTURE = 0;  // radius to use in PSF aperture mags
 static float PSF_FIT_RADIUS = 0;        // radius to use in fitting (ignored if <= 0,
                                         // and a per-object radius is calculated)
 
-bool psphotInitRadiusPSF(const psMetadata *recipe, const pmModelType type) {
+bool psphotInitRadiusPSF(const psMetadata *recipe, const psMetadata *analysis, const pmModelType type) {
 
     bool status = true;
 
-    PSF_FIT_NSIGMA  = psMetadataLookupF32(&status, recipe, "PSF_FIT_NSIGMA");
+    PSF_FIT_NSIGMA = psMetadataLookupF32(&status, recipe, "PSF_FIT_NSIGMA");
     PSF_FIT_PADDING = psMetadataLookupF32(&status, recipe, "PSF_FIT_PADDING");
-    PSF_FIT_RADIUS =  psMetadataLookupF32(&status, recipe, "PSF_FIT_RADIUS");
+
+    PSF_FIT_RADIUS =  psMetadataLookupF32(&status, analysis, "PSF_FIT_RADIUS");
+    if (!status) {
+        PSF_FIT_RADIUS = psMetadataLookupF32(&status, recipe, "PSF_FIT_RADIUS");
+    }
+
+    PSF_APERTURE =  psMetadataLookupF32(&status, analysis, "PSF_APERTURE");
+    if (!status) {
+        PSF_APERTURE =  psMetadataLookupF32(&status, recipe, "PSF_APERTURE");
+    }
 
     return true;
@@ -34 +44 @@
             radiusFit = model->modelRadius(model->params, 1.0);
         }
+        model->fitRadius = (RADIUS_TYPE)(radiusFit + PSF_FIT_PADDING);
+    } else {
+        model->fitRadius = radiusFit;
     }
-    model->radiusFit = (RADIUS_TYPE)(radiusFit + PSF_FIT_PADDING);
-
-    if (isnan(model->radiusFit)) psAbort("error in radius");
+    if (isnan(model->fitRadius)) psAbort("error in radius");
         
     if (source->mode & PM_SOURCE_MODE_SATSTAR) {
-        model->radiusFit *= 2;
+        model->fitRadius *= 2;
     }
 
-    bool status = pmSourceRedefinePixels (source, readout, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit);
+    // radius used to measure aperture photometry
+    source->apRadius = PSF_APERTURE;
+
+    bool status = pmSourceRedefinePixels (source, readout, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->fitRadius);
 
     // set the mask to flag the excluded pixels
-    psImageKeepCircle (source->maskObj, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit, "OR", markVal);
+    psImageKeepCircle (source->maskObj, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->fitRadius, "OR", markVal);
     return status;
 }
@@ -58 +72 @@
 
     // set the fit radius based on the object flux limit and the model
-    model->radiusFit = (RADIUS_TYPE) (model->modelRadius (model->params, PSF_FIT_NSIGMA*moments->dSky) + dR + PSF_FIT_PADDING);
-    if (isnan(model->radiusFit)) psAbort("error in radius");
-        
+    float radiusFit = PSF_FIT_RADIUS;
+    if (radiusFit <= 0) {               // use fixed radius
+        if (moments == NULL) {
+            radiusFit = model->modelRadius(model->params, PSF_FIT_NSIGMA*moments->dSky);
+        } else {
+            radiusFit = model->modelRadius(model->params, 1.0);
+        }
+        model->fitRadius = (RADIUS_TYPE)(radiusFit + PSF_FIT_PADDING);
+    } else {
+        model->fitRadius = radiusFit;
+    }
+    if (isnan(model->fitRadius)) psAbort("error in radius");
+
+    // above sets a radius for a single star, bump by blend separation
+    model->fitRadius += dR;
+
     if (source->mode &  PM_SOURCE_MODE_SATSTAR) {
-        model->radiusFit *= 2;
+        model->fitRadius *= 2;
     }
 
-    bool status = pmSourceRedefinePixels (source, readout, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit);
+    // radius used to measure aperture photometry
+    source->apRadius = PSF_APERTURE;
+
+    bool status = pmSourceRedefinePixels (source, readout, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->fitRadius);
 
     // set the mask to flag the excluded pixels
-    psImageKeepCircle (source->maskObj, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit, "OR", markVal);
+    psImageKeepCircle (source->maskObj, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->fitRadius, "OR", markVal);
     return status;
 }
@@ -74 +104 @@
 static float EXT_FIT_NSIGMA;
 static float EXT_FIT_PADDING;
+static float EXT_FIT_MAX_RADIUS;
 
 bool psphotInitRadiusEXT (psMetadata *recipe, pmModelType type) {
@@ -79 +110 @@
     bool status;
 
-    EXT_FIT_NSIGMA   = psMetadataLookupF32 (&status, recipe, "EXT_FIT_NSIGMA");
-    EXT_FIT_PADDING  = psMetadataLookupF32 (&status, recipe, "EXT_FIT_PADDING");
+    EXT_FIT_NSIGMA     = psMetadataLookupF32 (&status, recipe, "EXT_FIT_NSIGMA");
+    EXT_FIT_PADDING    = psMetadataLookupF32 (&status, recipe, "EXT_FIT_PADDING");
+    EXT_FIT_MAX_RADIUS = psMetadataLookupF32 (&status, recipe, "EXT_FIT_MAX_RADIUS");
 
     return true;
@@ -86 +118 @@
 
 // call this function whenever you (re)-define the EXT model
+float psphotSetRadiusEXT (pmReadout *readout, pmSource *source, psImageMaskType markVal) {
+
+    psAssert (source, "source not defined??");
+    psAssert (source->peak, "peak not defined??");
+
+    pmPeak *peak = source->peak;
+
+    // set the radius based on the footprint:
+    if (!peak->footprint) goto escape;
+    pmFootprint *footprint = peak->footprint;
+    if (!footprint->spans) goto escape;
+    if (footprint->spans->n < 1) goto escape;
+
+    // find the max radius
+    float radius = 0.0;
+    for (int j = 0; j < footprint->spans->n; j++) {
+        pmSpan *span = footprint->spans->data[j];
+
+        float dY  = span->y  - peak->yf;
+        float dX0 = span->x0 - peak->xf;
+        float dX1 = span->x1 - peak->xf;
+
+        radius = PS_MAX (radius, hypot(dY, dX0));
+        radius = PS_MAX (radius, hypot(dY, dX1));
+    }
+
+    radius += EXT_FIT_PADDING;
+    if (isnan(radius)) psAbort("error in radius");
+
+    radius = PS_MIN (radius, EXT_FIT_MAX_RADIUS);
+
+    // redefine the pixels if needed
+    pmSourceRedefinePixels (source, readout, peak->xf, peak->yf, radius);
+
+    // set the mask to flag the excluded pixels
+    psImageKeepCircle (source->maskObj, peak->xf, peak->yf, radius, "OR", markVal);
+    return radius;
+
+escape:
+    return NAN;
+    // bool result = psphotCheckRadiusEXT (readout, source, model, markVal);
+    // return result;
+}
+
+// alternative EXT radius based on model guess (for use without footprints)
 bool psphotCheckRadiusEXT (pmReadout *readout, pmSource *source, pmModel *model, psImageMaskType markVal) {
+
+    psAbort ("do not use this function");
 
     psF32 *PAR = model->params->data.F32;
@@ -96 +175 @@
     float rawRadius = model->modelRadius (model->params, EXT_FIT_NSIGMA*moments->dSky);
 
-    model->radiusFit = rawRadius + EXT_FIT_PADDING;
-    if (isnan(model->radiusFit)) psAbort("error in radius");
+    model->fitRadius = rawRadius + EXT_FIT_PADDING;
+    if (isnan(model->fitRadius)) psAbort("error in radius");
 
     // redefine the pixels if needed
-    bool status = pmSourceRedefinePixels (source, readout, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit);
+    bool status = pmSourceRedefinePixels (source, readout, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->fitRadius);
 
     // set the mask to flag the excluded pixels
-    psImageKeepCircle (source->maskObj, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit, "OR", markVal);
+    psImageKeepCircle (source->maskObj, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->fitRadius, "OR", markVal);
     return status;
 }

@@ -18 +18 @@
 psphotVersionDefinitions.h
 psphotMomentsStudy
+psphotPetrosianStudy
+psphotForced
+psphotMakePSF
+psphotStack