Index: /branches/eam_branches/ipp-20130509/psLib/src/math/psMinimizeLMM.c
===================================================================
--- /branches/eam_branches/ipp-20130509/psLib/src/math/psMinimizeLMM.c	(revision 35636)
+++ /branches/eam_branches/ipp-20130509/psLib/src/math/psMinimizeLMM.c	(revision 35637)
@@ -680,6 +680,7 @@
     }
 
-    psImage *alpha   = psImageAlloc(Alpha->numCols, Alpha->numRows, PS_TYPE_F32);
-    psVector *beta   = psVectorAlloc(Beta->n, PS_TYPE_F32);
+    int nFitParams = Beta->n;
+    psImage *alpha   = psImageAlloc(nFitParams, nFitParams, PS_TYPE_F32);
+    psVector *beta   = psVectorAlloc(nFitParams, PS_TYPE_F32);
     psVector *Params = psVectorAlloc(params->n, PS_TYPE_F32);
 
@@ -687,6 +688,6 @@
     psF32 Chisq = 0.0;
     psF32 lambda = 0.001;
+    psF32 dLinear = 0.0;
     psF32 nu = 2.0;
-    psF32 dLinear = 0.0;
 
     // the user provides the error or NULL.  we need to convert
@@ -700,5 +701,5 @@
 
     // number of degrees of freedom for this fit
-    int nDOF = dy->n - params->n;
+    int nDOF = dy->n - nFitParams;
 
     // calculate initial alpha and beta, set chisq (min->value)
@@ -722,9 +723,13 @@
     }
 
-    // iterate until: (a) nIter = min->iter or (b) (chisq / ndof) < maxChisq and deltaChisq < minTol (but don't stop unless Chisq is finite)
     bool done = (min->iter >= min->maxIter);
     while (!done) {
         psTrace("psLib.math", 5, "Iteration number %d.  (max iterations is %d).\n", min->iter, min->maxIter);
-        psTrace("psLib.math", 5, "Last delta is %f.  stop if < %f, accept if < %f\n", min->lastDelta, min->minTol, min->maxTol);
+
+	if (min->chisqConvergence) {
+	  psTrace("psLib.math", 5, "Last delta is %f.  stop if < %f, accept if < %f\n", min->lastDelta, min->minTol, min->maxTol);
+	} else {
+	  psTrace("psLib.math", 5, "Last delta is %f.  stop if < %f, accept if < %f\n", min->rParSigma, min->minTol*nFitParams, min->maxTol*nFitParams);
+	}
 
         // set a new guess for Alpha, Beta, Params
@@ -732,6 +737,6 @@
             min->iter ++;
 	    if (min->iter >=  min->maxIter) break;
-            // lambda *= 10.0;
-            lambda *= 2.0;
+            lambda *= 10.0;
+            // ALT? lambda *= 2.0;
             continue;
         }
@@ -754,17 +759,13 @@
 	// Beta is the actual parameter change for this pass
 
-	// note that Alpha & Beta only represent unmasked parameters
+	// note that Alpha & Beta only represent unmasked parameters, while params and Params have all 
 
 	// dParSigma = Alpha[i][i] : error (squared) on parameter i
-	// dParDelta = Beta[i]^2   : change (squared) in parameter i
+	// dParDelta = Params->data.F32[i] - params->data.F32[i]     : change on parameter i
 	float rParSigma = 0.0;
-        // for (int j = 0; j < Beta->n; j++) {
-	//     rParSigma += PS_SQR(Beta->data.F32[j]) / Alpha->data.F32[j][j];
-	// }
         for (int j = 0, J = 0; j < Params->n; j++) {
 	    if (paramMask && (paramMask->data.PS_TYPE_VECTOR_MASK_DATA[j])) {
 		continue;
 	    }
-	    // fprintf (stderr, "del: %f  - %f  vs %f  : %f\n", Params->data.F32[j], params->data.F32[j], sqrt(Alpha->data.F32[J][J]), PS_SQR(Params->data.F32[j] - params->data.F32[j]) / Alpha->data.F32[J][J]);
 	    rParSigma += PS_SQR(Params->data.F32[j] - params->data.F32[j]) / Alpha->data.F32[J][J];
 	    J++;
@@ -779,6 +780,6 @@
             min->iter ++;
 	    if (min->iter >= min->maxIter) break;
-            // lambda *= 10.0;
-            lambda *= 2.0;
+            lambda *= 10.0;
+            // ALT lambda *= 2.0;
             continue;
         }
@@ -807,16 +808,4 @@
 	// XXX the old version of lambda changes:
 	// XXX : Madsen gives suggestion for better use of rho
-# if (0)
-        if (rho >= -1e-6) {
-            min->value = Chisq;
-            alpha  = psImageCopy(alpha, Alpha, PS_TYPE_F32);
-            beta   = psVectorCopy(beta, Beta, PS_TYPE_F32);
-            params = psVectorCopy(params, Params, PS_TYPE_F32);
-            lambda *= 0.25;
-        } else {
-            lambda *= 10.0;
-        }
-# endif
-
         // rho is positive if the new chisq is smaller
         if (rho >= -1e-6) {
@@ -826,22 +815,33 @@
             params = psVectorCopy(params, Params, PS_TYPE_F32);
         } 
-
-# if (0)
-        // adjust the gain ratio (lambda) based on rho
-	if (rho < 0.25) {
-            lambda *= 2.0;
-        } 
-	if (rho > 0.75) {
-            lambda *= 0.333;
-        }
-# else
-	if (rho > 0.0) {
-	  lambda *= PS_MAX(0.33, (1.0 - pow(2.0*rho - 1.0, 3.0)));
-	  nu = 2.0;
-	} else {
-	  lambda *= nu;
-	  nu *= 2.0;
+	switch (min->gainFactorMode) {
+	  case 0:
+	    if (rho >= -1e-6) {
+	      lambda *= 0.25;
+	    } else {
+	      lambda *= 10.0;
+	    }
+	    break;
+
+	  case 1:
+	    // adjust the gain ratio (lambda) based on rho
+	    if (rho < 0.25) {
+	      lambda *= 2.0;
+	    } 
+	    if (rho > 0.75) {
+	      lambda *= 0.333;
+	    }
+	    break;
+
+	  case 2:
+	    if (rho > 0.0) {
+	      lambda *= PS_MAX(0.33, (1.0 - pow(2.0*rho - 1.0, 3.0)));
+	      nu = 2.0;
+	    } else {
+	      lambda *= nu;
+	      nu *= 2.0;
+	    }
+	    break;
 	}
-# endif
         min->iter++;
 
@@ -850,5 +850,14 @@
 	done = (min->iter >= min->maxIter);
 	
+	// 2) require deltaChi > 1e-6 (ie, chisq is decreasing, but accept an insignificant change)
+	if (min->lastDelta < -1e-6) {
+	    continue;
+	}
+
+	// save this value in case we stop iterating
+	min->rParSigma = rParSigma;
+
 	// 2) require chisqDOF < maxChisqDOF (if maxChisqDOF is not NAN)
+	// keep iterating regardless of rParSigma in this case
 	float chisqDOF = Chisq / nDOF;
 	if (isfinite(min->maxChisqDOF) && (chisqDOF > min->maxChisqDOF)) {
@@ -856,13 +865,12 @@
 	}
 
-	// 3) require deltaChi > 1e-6 (ie, chisq is decreasing, but accept an insignificant change)
-	if (min->lastDelta < -1e-6) {
-	    continue;
+	// delta-chisq or rParSigma ?
+	if (min->chisqConvergence) {
+	  done |= (min->lastDelta < min->minTol);
+	} else {
+	  done |= (rParSigma < min->minTol*nFitParams);
 	}
-
-	// 4) require rParDelta < \eta * rParSigma
-	done |= (rParSigma < 0.01);
-    }
-    psTrace("psLib.math", 5, "chisq: %f, last delta: %f, Niter: %d\n", min->value, min->lastDelta, min->iter);
+    }
+    psTrace("psLib.math", 5, "chisq: %f, last delta: %f, rParSigma: %f, Niter: %d\n", min->value, min->lastDelta, min->rParSigma, min->iter);
 
     // construct & return the covariance matrix (if requested)
@@ -898,8 +906,14 @@
 
     // if the last improvement was at least as good as maxTol, accept the fit:
-    // XXX : maybe re-work this to use rParSigma & maxTol
-    if (min->lastDelta <= min->maxTol) {
+    if (min->chisqConvergence) {
+      if (min->lastDelta <= min->maxTol) {
 	psTrace("psLib.math", 6, "---- end (true) ----\n");
         return(true);
+      }
+    } else {
+      if (min->rParSigma <= min->maxTol*nFitParams) {
+	psTrace("psLib.math", 6, "---- end (true) ----\n");
+        return(true);
+      }
     }
     psTrace("psLib.math", 6, "---- end (false) ----\n");
@@ -952,6 +966,11 @@
     min->iter = 0;
     min->lastDelta = NAN;
+    min->rParSigma = NAN;
     min->maxChisqDOF = NAN;
 
+    // we default to the old algorithm for convergence
+    min->chisqConvergence = true;
+    min->gainFactorMode = 0;
+    
     return(min);
 }
Index: /branches/eam_branches/ipp-20130509/psLib/src/math/psMinimizeLMM.h
===================================================================
--- /branches/eam_branches/ipp-20130509/psLib/src/math/psMinimizeLMM.h	(revision 35636)
+++ /branches/eam_branches/ipp-20130509/psLib/src/math/psMinimizeLMM.h	(revision 35637)
@@ -30,4 +30,7 @@
 #include "psConstants.h"
 
+# define PS_MINIMIZE_LMM_GAIN_FACTOR_MODE 0
+# define PS_MINIMIZE_LMM_CHISQ_CONVERGENCE 1
+
 #define PS_DETERMINE_BRACKET_STEP_SIZE 0.10
 #define PS_MAX_LMM_ITERATIONS 100
@@ -82,5 +85,8 @@
     int iter;                          ///< Number of iterations to date
     float lastDelta;                   ///< The last difference for the fit
+    float rParSigma;		       ///< last fractional change in the parameters
     float maxChisqDOF;		       ///< for Chisq minimization, require that we reach here before checking tolerance
+    int gainFactorMode;
+    bool chisqConvergence;
 }
 psMinimization;