Changeset 2250 for trunk/psLib/test/dataManip/tst_psMinimize07.c

Timestamp:

Nov 1, 2004, 1:57:08 PM (22 years ago)

Author:

gusciora

Message:

Fixed the robust stats routines for stdev/mean. Fixed minimization routines
for fitting gaussians. CVS:

---

Enter Log. Lines beginning with `CVS:' are removed automatically CVS: CVS:
Committing in . CVS: CVS: Modified Files: CVS:
src/dataManip/psMinimize.c src/dataManip/psMinimize.h CVS:
src/dataManip/psStats.c src/dataManip/makedir/psMinimize.d CVS:
src/dataManip/makedir/psStats.d CVS: test/dataManip/tst_psMinimize07.c

---

File:

• trunk/psLib/test/dataManip/tst_psMinimize07.c (modified) (6 diffs)

trunk/psLib/test/dataManip/tst_psMinimize07.c

@@ -10 +10 @@
 #define N 24
 #define NUM_PARAMS 2
-#define MEAN 8.5
-#define STDEV 1.5
+#define MEAN 14.5
+#define STDEV 1.2
 
 float expectedParm[NUM_PARAMS];
 psS32 testStatus = true;
 
-psVector *genGaussianVector()
+// This routine generates a vector of length "n" with the specified mean and
+// stdev.
+psVector *genGaussianVector(int n, float mean, float stdev)
 {
-    psVector *y = psVectorAlloc(N, PS_TYPE_F32);
+    psVector *y = psVectorAlloc(n, PS_TYPE_F32);
 
-    for (psS32 i=0;i<N;i++) {
-        y->data.F32[i] = psGaussian((float) i, MEAN, STDEV, false);
+    for (psS32 i=0;i<n;i++) {
+        y->data.F32[i] = psGaussian((float) i, mean, stdev, false);
     }
     return(y);
 }
 
-psVector *myFunc(psVector *params,
-                 psArray *coords)
-{
-    PS_PTR_CHECK_NULL(coords, NULL);
-    PS_PTR_CHECK_NULL(params, NULL);
-
-    float x;
-    int i;
-    float mean = params->data.F32[0];
-    float stdev = params->data.F32[1];
-    psVector *out = psVectorAlloc(coords->n, PS_TYPE_F32);
-
-    for (i=0;i<coords->n;i++) {
-        x = ((psVector *) (coords->data[i]))->data.F32[0];
-        out->data.F32[i] = psGaussian(x, mean, stdev, false);
-    }
-
-    return(out);
-}
-
+// This routine tries to fit a Gaussian to a set of data points via the
+// psMinimizeLMChi2() and psMinimizeLMChi2Gauss1D() functions.
 psS32 t00()
 {
@@ -56 +40 @@
     psMinimization *min;
 
-    psTraceSetLevel(".psLib.dataManip.psMinimize", 6);
-    psTraceSetLevel(".psLib.dataManip.p_psDetermineBracket", 0);
+    psTraceSetLevel(".psLib.dataManip.psMinimize", 0);
     /**************************************************************************
      *************************************************************************/
@@ -63 +46 @@
     myCovar = psImageAlloc(NUM_PARAMS, NUM_PARAMS, PS_TYPE_F32);
     myParams = psVectorAlloc(NUM_PARAMS, PS_TYPE_F32);
-    myParams->data.F32[0] = (float) (N/2);
-    myParams->data.F32[1] = 0.5;
+    myParams->data.F32[0] = MEAN * 0.7;
+    myParams->data.F32[1] = STDEV * 0.7;
     myCoords = psArrayAlloc(N);
-    psVector *y = genGaussianVector();
+    psVector *y = genGaussianVector(N, MEAN, STDEV);
 
     for (i=0;i<N;i++) {
         myCoords->data[i] = (psPtr *) psVectorAlloc(2, PS_TYPE_F32);
         ((psVector *) (myCoords->data[i]))->data.F32[0] = (float) i;
+        printf("Initial data %d: (%f, %f)\n", i, ((psVector *) (myCoords->data[i]))->data.F32[0], y->data.F32[i]);
     }
 
@@ -105 +89 @@
 }
 
+// This routine tries to fit a Gaussian to a set of data points via the
+// psMinimizeChi2Powell() function.
 psS32 t01()
 {
@@ -115 +101 @@
     psMinimization *min;
 
-    psTraceSetLevel(".psLib.dataManip.psMinimizePowell", 10);
-    psTraceSetLevel(".psLib.dataManip.myPowellChi2Func", 10);
+    psTraceSetLevel(".psLib.dataManip.psMinimize", 0);
+    psTraceSetLevel(".psLib.dataManip.psMinimizeLMChi2Gauss1D", 0);
+
+    psTraceSetLevel(".psLib.dataManip.psMinimizePowell", 0);
+    psTraceSetLevel(".psLib.dataManip.myPowellChi2Func", 0);
+    psTraceSetLevel(".psLib.dataManip.p_psLineMin", 0);
+    psTraceSetLevel(".psLib.dataManip.p_psDetermineBracket", 0);
     psTraceSetLevel(".psLib.dataManip.psFunctions.psGaussian", 0);
     psTraceSetLevel(".psLib.dataManip.psFunctions", 0);
-    psTraceSetLevel(".", 10);
-    psTracePrintLevels();
+    psTraceSetLevel(".", 0);
+    //    psTracePrintLevels();
     /**************************************************************************
      *************************************************************************/
@@ -126 +117 @@
     myCovar = psImageAlloc(NUM_PARAMS, NUM_PARAMS, PS_TYPE_F32);
     myParams = psVectorAlloc(NUM_PARAMS, PS_TYPE_F32);
-    myParams->data.F32[0] = (float) (N/2);
-    myParams->data.F32[0] = 2.0 * MEAN;
-    myParams->data.F32[1] = 0.5;
-    myParams->data.F32[1] = 1.5 * STDEV;
+    myParams->data.F32[0] = MEAN * 0.7;
+    myParams->data.F32[1] = STDEV * 0.7;
     myCoords = psArrayAlloc(N);
-    psVector *y = genGaussianVector();
+    psVector *y = genGaussianVector(N, MEAN, STDEV);
 
     for (i=0;i<N;i++) {
         myCoords->data[i] = (psPtr *) psVectorAlloc(2, PS_TYPE_F32);
         ((psVector *) (myCoords->data[i]))->data.F32[0] = (float) i;
+        printf("Initial data %d: (%f, %f)\n", i, ((psVector *) (myCoords->data[i]))->data.F32[0], y->data.F32[i]);
     }
 
+    //    psMinimizeChi2Powell(min, myParams, NULL, myCoords, y, NULL,
+    //                         (psMinimizeChi2PowellFunc) myFunc);
     psMinimizeChi2Powell(min, myParams, NULL, myCoords, y, NULL,
-                         (psMinimizeChi2PowellFunc) myFunc);
+                         (psMinimizeChi2PowellFunc) psMinimizePowellChi2Gauss1D);
 
     printf("\nThe chi-squared is %f\n", min->value);