Index: trunk/psModules/test/objects/tap_pmSourceFitModel.c
===================================================================
--- trunk/psModules/test/objects/tap_pmSourceFitModel.c	(revision 10095)
+++ trunk/psModules/test/objects/tap_pmSourceFitModel.c	(revision 10098)
@@ -7,6 +7,8 @@
 #include "pstap.h"
 
-bool fitModels (psRandom *seed, float flux, float radius);
-bool fitModelFlux (psRandom *seed, float flux, float radius);
+bool fitModels (psRandom *seed, float flux, float radius, float sigma);
+
+// create a single source (sigma =
+bool fitModelFlux (psRandom *seed, float flux, float radius, float sigma);
 
 int main (void)
@@ -15,6 +17,4 @@
     pmSourceFitModelInit (15, 0.01, 1.0, true);
 
-    // psTraceSetLevel ("psLib.math.psMinimizeLMChi2", 6);
-
     plan_tests(240);
 
@@ -22,31 +22,15 @@
     psRandom *seed = psRandomAlloc (PS_RANDOM_TAUS, 0);
 
-    # if (0)
+    static float radius[] = {3.0, 5.0, 7.0, 10.0, 15.0, 25.0};
+    static float sigma[] = {1.0, 1.5, 2.0};
+    static float flux[] = {10000.0, 3000.0, 1000.0, 300.0, 100.0, 30.0, 10.0};
 
-        fitModelFlux (seed, 10000.0);
-    fitModelFlux (seed, 3000.0);
-    fitModelFlux (seed, 1000.0);
-    fitModelFlux (seed, 300.0);
-    fitModelFlux (seed, 100.0);
-    fitModelFlux (seed, 30.0);
-    fitModelFlux (seed, 30.0, 4.0);
-    fitModelFlux (seed, 30.0, 5.0);
-    fitModelFlux (seed, 30.0, 7.0);
-    fitModelFlux (seed, 30.0, 10.0);
-    fitModelFlux (seed, 30.0, 15.0);
-    fitModelFlux (seed, 30.0, 20.0);
-    fitModelFlux (seed, 30.0, 25.0);
-    fitModelFlux (seed, 30.0, 35.0);
-    # endif
-
-    fitModels (seed, 100.0,  7.0);
-    fitModels (seed,  80.0,  7.0);
-    fitModels (seed,  60.0,  7.0);
-    fitModels (seed,  40.0,  7.0);
-    fitModels (seed,  20.0,  7.0);
-
-    fitModels (seed,  30.0,  5.0);
-    fitModels (seed,  30.0, 10.0);
-    fitModels (seed,  30.0, 15.0);
+    for (int i = 0; i < sizeof(sigma)/sizeof(float); i++) {
+        for (int j = 0; j < sizeof(radius)/sizeof(float); j++) {
+            for (int k = 0; k < sizeof(flux)/sizeof(float); k++) {
+                fitModels (seed, flux[k], radius[j], sigma[i]);
+            }
+        }
+    }
 
     return exit_status();
@@ -59,36 +43,46 @@
 static psVector *par5 = NULL;
 
-bool fitModels (psRandom *seed, float flux, float radius)
+bool fitModels (psRandom *seed, float flux, float radius, float sigma)
 {
 
     psMemId id = psMemGetId();
 
-    diag("test model fit");
+    diag("test model fit - flux: %f, radius: %f, sigma: %f", flux, radius, sigma);
 
-    par1 = psVectorAllocEmpty (100, PS_TYPE_F32);
-    par2 = psVectorAllocEmpty (100, PS_TYPE_F32);
-    par3 = psVectorAllocEmpty (100, PS_TYPE_F32);
-    par4 = psVectorAllocEmpty (100, PS_TYPE_F32);
-    par5 = psVectorAllocEmpty (100, PS_TYPE_F32);
+    par1 = psVectorAllocEmpty (200, PS_TYPE_F32);
+    par2 = psVectorAllocEmpty (200, PS_TYPE_F32);
+    par3 = psVectorAllocEmpty (200, PS_TYPE_F32);
+    par4 = psVectorAllocEmpty (200, PS_TYPE_F32);
+    par5 = psVectorAllocEmpty (200, PS_TYPE_F32);
 
-    for (int i = 0; i < 100; i++) {
-        fitModelFlux (seed, flux, radius);
+    for (int i = 0; i < 200; i++) {
+        fitModelFlux (seed, flux, radius, sigma);
     }
+
+    float signal = 2*M_PI*sigma*sigma*flux;
+    float noise = sqrt(signal + 4*M_PI*sigma*sigma*(100 + PS_SQR(5)));
+    float dMag = noise / signal;
+    float dPos = sigma * dMag;
+    diag ("signal: %f, noise: %f, dMag: %f, dPos: %f", signal, noise, dMag, dPos);
 
     psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEAN | PS_STAT_SAMPLE_STDEV);
     psVectorStats (stats, par1, NULL, NULL, 0);
-    fprintf (stderr, "Io: %f +/- %f\n", stats->sampleMean, stats->sampleStdev);
-    // ok_float_tol (stats->sampleStdev, 0.0, 0.0, "Io ref/fit stdev: %f", stats->sampleStdev);
+    // fprintf (stderr, "Io: %f +/- %f : %f sig\n", stats->sampleMean, stats->sampleStdev, stats->sampleStdev/dMag);
+    ok ((stats->sampleStdev/dMag < 2.0), "Io ref/fit stdev: %f : %f sigma", stats->sampleStdev, stats->sampleStdev/dMag);
     psVectorStats (stats, par2, NULL, NULL, 0);
-    fprintf (stderr, "Xo: %f +/- %f\n", stats->sampleMean, stats->sampleStdev);
+    // fprintf (stderr, "Xo: %f +/- %f : %f sig\n", stats->sampleMean, stats->sampleStdev, stats->sampleStdev/dPos);
+    ok ((stats->sampleStdev/dPos < 2.0), "Xo ref/fit stdev: %f : %f sigma", stats->sampleStdev, stats->sampleStdev/dPos);
     // ok_float_tol (stats->sampleStdev, 0.0, 0.0, "Xo ref-fit stdev: %f", stats->sampleStdev);
     psVectorStats (stats, par3, NULL, NULL, 0);
-    fprintf (stderr, "Yo: %f +/- %f\n", stats->sampleMean, stats->sampleStdev);
+    // fprintf (stderr, "Yo: %f +/- %f : %f sig\n", stats->sampleMean, stats->sampleStdev, stats->sampleStdev/dPos);
+    ok ((stats->sampleStdev/dPos < 2.0), "Yo ref/fit stdev: %f : %f sigma", stats->sampleStdev, stats->sampleStdev/dPos);
     // ok_float_tol (stats->sampleStdev, 0.0, 0.0, "Yo ref-fit stdev: %f", stats->sampleStdev);
     psVectorStats (stats, par4, NULL, NULL, 0);
-    fprintf (stderr, "Sx: %f +/- %f\n", stats->sampleMean, stats->sampleStdev);
+    // fprintf (stderr, "Sx: %f +/- %f : %f sig\n", stats->sampleMean, stats->sampleStdev, stats->sampleStdev/dMag);
+    ok ((stats->sampleStdev/dMag < 2.0), "Sx ref/fit stdev: %f : %f sigma", stats->sampleStdev, stats->sampleStdev/dMag);
     // ok_float_tol (stats->sampleStdev, 0.0, 0.0, "Sx ref/fit stdev: %f", stats->sampleStdev);
     psVectorStats (stats, par5, NULL, NULL, 0);
-    fprintf (stderr, "Sy: %f +/- %f\n", stats->sampleMean, stats->sampleStdev);
+    // fprintf (stderr, "Sy: %f +/- %f : %f sig\n", stats->sampleMean, stats->sampleStdev, stats->sampleStdev/dMag);
+    ok ((stats->sampleStdev/dMag < 2.0), "Sy ref/fit stdev: %f : %f sigma", stats->sampleStdev, stats->sampleStdev/dMag);
     // ok_float_tol (stats->sampleStdev, 0.0, 0.0, "Sy ref/fit stdev: %f", stats->sampleStdev);
 
@@ -98,4 +92,5 @@
     psFree (par4);
     psFree (par5);
+    psFree (stats);
 
     ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
@@ -103,5 +98,5 @@
 }
 
-bool fitModelFlux (psRandom *seed, float flux, float radius)
+bool fitModelFlux (psRandom *seed, float flux, float radius, float sigma)
 {
 
@@ -122,6 +117,6 @@
     source->modelEXT->params->data.F32[2] = 50;
     source->modelEXT->params->data.F32[3] = 50;
-    source->modelEXT->params->data.F32[4] = 2.0*sqrt(2.0);
-    source->modelEXT->params->data.F32[5] = 2.0*sqrt(2.0);
+    source->modelEXT->params->data.F32[4] = 2.0*sqrt(sigma);
+    source->modelEXT->params->data.F32[5] = 2.0*sqrt(sigma);
     source->modelEXT->params->data.F32[6] = 0;
 
@@ -176,4 +171,8 @@
     psVectorExtend (par5, 100, 1);
 
+    psFree (rnd);
+    psFree (source);
+    psFree (guess);
+
     return true;
 }