Changeset 14669

Timestamp:

Aug 27, 2007, 9:04:17 AM (19 years ago)

Author:

magnier

Message:

calculate radius for sersic correctly

File:

• trunk/psModules/src/objects/models/pmModel_SERSIC.c (modified) (2 diffs)

trunk/psModules/src/objects/models/pmModel_SERSIC.c

@@ -299 +299 @@
 psF64 PM_MODEL_RADIUS (const psVector *params, psF64 flux)
 {
-    psF64 z, f;
-    int Nstep = 0;
     psEllipseShape shape;
 
@@ -321 +319 @@
     psF64 limit = flux / PAR[PM_PAR_I0];
 
-    // use the fact that f is monotonically decreasing
-    z = 0;
-    Nstep = 0;
-
-    // choose a z value guaranteed to be beyond our limit
-    float z0 = pow((1.0 / limit), (1.0 / 2.25));
-    float z1 = (1.0 / limit) / PAR[PM_PAR_7];
-    z1 = PS_MAX (z0, z1);
-    z0 = 0.0;
-
-    // perform a type of bisection to find the value
-    float f0 = 1.0 / (1 + PAR[PM_PAR_7]*z0 + pow(z0, 2.25));
-    float f1 = 1.0 / (1 + PAR[PM_PAR_7]*z1 + pow(z1, 2.25));
-    while ((Nstep < 10) && (fabs(z1 - z0) > 0.5)) {
-        z = 0.5*(z0 + z1);
-        f = 1.0 / (1 + PAR[PM_PAR_7]*z + pow(z, 2.25));
-        if (f > limit) {
-            z0 = z;
-            f0 = f;
-        } else {
-            z1 = z;
-            f1 = f;
-        }
-        Nstep ++;
-    }
+    psF64 z = pow (-log(limit), (1.0 / PAR[PM_PAR_7]));
+
     psF64 radius = sigma * sqrt (2.0 * z);