Changeset 25433 for branches/eam_branches/20090715/psphot/src/psphotRadiiFromProfiles.c

Timestamp:

Sep 17, 2009, 2:26:32 PM (17 years ago)

Author:

eugene

Message:

change radius for extended sources to use footprint; clean up some of the visualizations; plug some leaks

File:

• branches/eam_branches/20090715/psphot/src/psphotRadiiFromProfiles.c (modified) (6 diffs)

branches/eam_branches/20090715/psphot/src/psphotRadiiFromProfiles.c

@@ -3 +3 @@
 // Given the Petrosian data (radii, fluxes) determine the radius for each profile at the desisred isophote
 
-bool psphotRadiiFromProfiles (pmPetrosian *petrosian, float fluxMin, float fluxMax) {
+bool psphotRadiiFromProfiles (pmSource *source, pmPetrosian *petrosian, float fluxMin, float fluxMax) {
 
   petrosian->isophotalRadii = psVectorAlloc(petrosian->theta->n, PS_TYPE_F32);
@@ -10 +10 @@
       psVector *radii = petrosian->radii->data[i];
       psVector *fluxes = petrosian->fluxes->data[i];
-      float radius = psphotRadiusFromProfile (radii, fluxes, fluxMin, fluxMax);
+      float radius = psphotRadiusFromProfile (source, radii, fluxes, fluxMin, fluxMax);
 
       // psphotPetrosianVisualProfileByAngle (radii, fluxes, radius);
@@ -20 +20 @@
 }
 
-float psphotRadiusFromProfile (psVector *radius, psVector *flux, float fluxMin, float fluxMax) {
+float psphotRadiusFromProfile (pmSource *source, psVector *radius, psVector *flux, float fluxMin, float fluxMax) {
 
     // 'flux' is a noisy sample of the galaxy radial profile at points 'radius'
@@ -49 +49 @@
             psVectorAppend (values, radius->data.F32[i]);
         }
-        psVectorStats (stats, values, NULL, NULL, 0);
+        if (values->n > 1) {
+            psVectorStats (stats, values, NULL, NULL, 0);
 
-        // if we have a valid range, rebin with bin size 1/2 of median radius
-        if (isfinite(stats->sampleMedian)) {
-            Rbin = MAX(1, 0.5*stats->sampleMedian);
+            // if we have a valid range, rebin with bin size 1/2 of median radius
+            if (isfinite(stats->sampleMedian)) {
+                Rbin = MAX(1, 0.5*stats->sampleMedian);
+            }
         }
+        psFree (values);
+        psFree (stats);
     }
     Rbin = 3;
@@ -90 +94 @@
                 // calculate the value for the nOut bin
                 // XXX need to fix this as well psStats (stats, values);
-                psVectorStats (stats, values, NULL, NULL, 0);
-                fluxBinned->data.F32[nOut] = stats->sampleMedian;
+                float value;
+                if (values->n > 0) {
+                    psVectorStats (stats, values, NULL, NULL, 0);
+                    value = stats->sampleMedian;
+                } else {
+                    value = NAN;
+                }
+                fluxBinned->data.F32[nOut] = value;
                 nOut ++;
                 radiusBinned->data.F32[nOut] = (nOut + 0.5)*Rbin;
@@ -118 +128 @@
             // XXX is there a macro in psLib that does this interpolation?
             if (i == 0) { 
-                psWarning ("impossible condition f[0] < Fo");
-                continue;
-            } else {
-                Ro = radiusBinned->data.F32[i-1] + (radiusBinned->data.F32[i] - radiusBinned->data.F32[i-1]) * (Fo - fluxBinned->data.F32[i-1]) / (fluxBinned->data.F32[i] - fluxBinned->data.F32[i-1]);
-            }
+                psLogMsg ("psphot", 3, "bogus radial profile for ..., skipping");
+                psFree (fluxBinned);
+                psFree (radiusBinned);
+                return NAN;
+            } 
+            Ro = radiusBinned->data.F32[i-1] + (radiusBinned->data.F32[i] - radiusBinned->data.F32[i-1]) * (Fo - fluxBinned->data.F32[i-1]) / (fluxBinned->data.F32[i] - fluxBinned->data.F32[i-1]);
             above = FALSE;
         }