Changeset 23594 for branches/cnb_branches/cnb_branch_20090301/psLib/src/math/psMinimizeLMM.c

Timestamp:

Mar 29, 2009, 6:15:31 PM (17 years ago)

Author:

beaumont

Message:

merged with head

Location:

branches/cnb_branches/cnb_branch_20090301

Files:

• . (modified) (1 prop)
• psLib/src/math/psMinimizeLMM.c (modified) (14 diffs)

branches/cnb_branches/cnb_branch_20090301/psLib/src/math/psMinimizeLMM.c

@@ -95 +95 @@
     // check for non-finite entries in result
     for (int i = 0; i < Beta->n; i++) {
-        if (!isfinite(Beta->data.F32[i])) {
-            // psError(PS_ERR_BAD_PARAMETER_VALUE, 3, "Fit value diverges: vector[%d] is %.2f\n", i, Beta->data.F32[i]);
-            return false;
-        }
+        if (!isfinite(Beta->data.F32[i])) {
+            // psError(PS_ERR_BAD_PARAMETER_VALUE, 3, "Fit value diverges: vector[%d] is %.2f\n", i, Beta->data.F32[i]);
+            return false;
+        }
     }
 
@@ -104 +104 @@
     // (we must do this before truncating Beta below)
     if (dLinear) {
-        *dLinear = psMinLM_dLinear(Beta, beta, lambda);
-    }
-
-    // full-length Beta for checkLimits functions 
+        *dLinear = psMinLM_dLinear(Beta, beta, lambda);
+    }
+
+    // full-length Beta for checkLimits functions
     psVector *tmpBeta = psVectorAlloc(params->n, PS_TYPE_F32);
     psVectorInit (tmpBeta, 0.0);
@@ -113 +113 @@
     // set tmpBeta values which are not masked
     for (int j = 0, n = 0; j < params->n; j++) {
-        if (paramMask && (paramMask->data.PS_TYPE_VECTOR_MASK_DATA[j])) continue;
-        tmpBeta->data.F32[j] = Beta->data.F32[n];
-        n++;
+        if (paramMask && (paramMask->data.PS_TYPE_VECTOR_MASK_DATA[j])) continue;
+        tmpBeta->data.F32[j] = Beta->data.F32[n];
+        n++;
     }
 
@@ -127 +127 @@
         if (checkLimits) {
             checkLimits (PS_MINIMIZE_BETA_LIMIT, j, Params->data.F32, tmpBeta->data.F32);
-        }
+        }
 
         Params->data.F32[j] = params->data.F32[j] - tmpBeta->data.F32[j];
@@ -135 +135 @@
             checkLimits (PS_MINIMIZE_PARAM_MIN,  j, Params->data.F32, tmpBeta->data.F32);
             checkLimits (PS_MINIMIZE_PARAM_MAX,  j, Params->data.F32, tmpBeta->data.F32);
-        }
+        }
     }
 
     // apply tmpBeta after limits have been checked
     for (int j = 0, n = 0; j < params->n; j++) {
-        if (paramMask && (paramMask->data.PS_TYPE_VECTOR_MASK_DATA[j])) continue;
-        Beta->data.F32[n] = tmpBeta->data.F32[j];
-        n++;
+        if (paramMask && (paramMask->data.PS_TYPE_VECTOR_MASK_DATA[j])) continue;
+        Beta->data.F32[n] = tmpBeta->data.F32[j];
+        n++;
     }
 
@@ -169 +169 @@
     // if none are available, return false
     if (!psMinLM_AllocAB (&Alpha, &Beta, params, paramMask)) {
-        return false;
+        return false;
     }
 
@@ -241 +241 @@
     // beta only counts unmasked parameters
     for (int i = 0; i < beta->n; i++) {
-        dh = lambda*B[i] + b[i];
-        sh = 0.5*B[i]*dh;
-        Sh += sh;
+        dh = lambda*B[i] + b[i];
+        sh = 0.5*B[i]*dh;
+        Sh += sh;
         dLinear += lambda*PS_SQR(B[i]) + B[i]*b[i];
     }
@@ -295 +295 @@
         assert (!isnan(chisq));
 
-        // we track alpha,beta and params,deriv separately
+        // we track alpha,beta and params,deriv separately
         for (int j = 0, J = 0; j < params->n; j++) {
             if (paramMask && (paramMask->data.PS_TYPE_VECTOR_MASK_DATA[j])) continue;
@@ -304 +304 @@
                 if (paramMask && (paramMask->data.PS_TYPE_VECTOR_MASK_DATA[k])) continue;
                 alpha->data.F32[J][K] += weight * deriv->data.F32[k];
-                K++;
+                K++;
             }
             beta->data.F32[J] += weight * delta;
-            J++;
+            J++;
         }
     }
@@ -328 +328 @@
 between that function at the specified coords and the specified value at those
 coords.
- 
+
 This requires F32 input data; all internal calls use F32.
 XXX Make an F64 version?
@@ -386 +386 @@
     // Alpha & Beta only contain elements to represent the unmasked parameters
     if (!psMinLM_AllocAB (&Alpha, &Beta, params, paramMask)) {
-        psAbort ("programming error: no unmasked parameters to be fit\n");
-    }
-    
+        psAbort ("programming error: no unmasked parameters to be fit\n");
+    }
+
     psImage *alpha   = psImageAlloc(Alpha->numCols, Alpha->numRows, PS_TYPE_F32);
     psVector *beta   = psVectorAlloc(Beta->n, PS_TYPE_F32);
@@ -490 +490 @@
             psTrace ("psLib.math", 5, "failure to calculate covariance matrix\n");
         }
-        // set covar values which are not masked
-        psImageInit (covar, 0.0);
-        for (int j = 0, J = 0; j < params->n; j++) {
-            if (paramMask && (paramMask->data.PS_TYPE_VECTOR_MASK_DATA[j])) {
-                covar->data.F32[j][j] = 1.0;
-                continue;
-            }
-            for (int k = 0, K = 0; k < params->n; k++) {
-                if (paramMask && (paramMask->data.PS_TYPE_VECTOR_MASK_DATA[k])) continue;
-                covar->data.F32[j][k] = Alpha->data.F32[J][K];
-                K++;
-            }
-            J++;
-        }
+        // set covar values which are not masked
+        psImageInit (covar, 0.0);
+        for (int j = 0, J = 0; j < params->n; j++) {
+            if (paramMask && (paramMask->data.PS_TYPE_VECTOR_MASK_DATA[j])) {
+                covar->data.F32[j][j] = 1.0;
+                continue;
+            }
+            for (int k = 0, K = 0; k < params->n; k++) {
+                if (paramMask && (paramMask->data.PS_TYPE_VECTOR_MASK_DATA[k])) continue;
+                covar->data.F32[j][k] = Alpha->data.F32[J][K];
+                K++;
+            }
+            J++;
+        }
     }
 
@@ -533 +533 @@
     // count unmasked parameters
     if (paramMask) {
-        nParams = 0;
-        for (int i = 0; i < paramMask->n; i++) {
-            if (paramMask->data.PS_TYPE_VECTOR_MASK_DATA[i]) continue;
-            nParams ++;
-        }
-    }
-
-    if (nParams == 0) { 
-        return false;
+        nParams = 0;
+        for (int i = 0; i < paramMask->n; i++) {
+            if (paramMask->data.PS_TYPE_VECTOR_MASK_DATA[i]) continue;
+            nParams ++;
+        }
+    }
+
+    if (nParams == 0) {
+        return false;
     }
 
@@ -585 +585 @@
 }
 
-psMinConstraint* psMinConstraintAlloc()
+psMinConstraint* psMinConstraintAlloc(void)
 {
     psMinConstraint *tmp = psAlloc(sizeof(psMinConstraint));