Changeset 35798

Timestamp:

Jul 11, 2013, 8:00:59 AM (13 years ago)

Author:

watersc1

Message:

Fix to psMM that avoids memory corruption when Ncensored = 0.

Location:

trunk/psLib/src/math

Files:

• psMinimizeLMM.c (modified) (3 diffs)
• psMixtureModels.c (modified) (4 diffs, 1 prop)

trunk/psLib/src/math/psMinimizeLMM.c

@@ -501 +501 @@
         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);
+        psTrace("psLib.math.dLinear", 5, "Iteration number %d.  (max iterations is %d).\n", min->iter, min->maxIter);
+        psTrace("psLib.math.dLinear", 5, "Last delta is %f.  stop if < %f, accept if < %f\n", min->lastDelta, min->minTol, min->maxTol);
 
         // set a new guess for Alpha, Beta, Params
@@ -519 +521 @@
             p_psVectorPrint(psTraceGetDestination(), Params, "params guess (1)");
         }
+        if (psTraceGetLevel("psLib.math.dLinear") >= 6) {
+          p_psImagePrint(psTraceGetDestination(), Alpha, "alpha guess (1)");
+          p_psVectorPrint(psTraceGetDestination(), Beta, "beta guess (1)");
+          p_psVectorPrint(psTraceGetDestination(), Params, "params guess (1)");
+          p_psVectorPrint(psTraceGetDestination(), params, "params guess (1)");
+        }         
 
         // calculate Chisq for new guess, update Alpha & Beta
@@ -538 +546 @@
 
         psTrace("psLib.math.dLinear", 5, "last chisq: %f, new chisq %f, delta: %f, dLinear: %f, rho: %f, lambda: %g\n", min->value, Chisq, min->lastDelta, dLinear, rho, lambda);
-
+        if (psTraceGetLevel("psLib.math.dLinear") >= 6) {
+          p_psImagePrint(psTraceGetDestination(), Alpha, "alpha guess (2)");
+          p_psVectorPrint(psTraceGetDestination(), Beta, "beta guess (2)");
+          p_psVectorPrint(psTraceGetDestination(), Params, "params guess (2)");
+        }         
         // dump some useful info if trace is defined
         if (psTraceGetLevel("psLib.math") >= 6) {

trunk/psLib/src/math/psMixtureModels.c

@@ -116 +116 @@
                         psVector *offsets) {
   // Allocate offset vector
-  offsets = psVectorRecycle(offsets,in->numRows,PS_TYPE_S32);
+  //  offsets = psVectorRecycle(offsets,in->numRows,PS_TYPE_S32);
 
   *Ncensored = 0;
@@ -169 +169 @@
   psImage  *Dcensored = psMMCensorData(D,&Ncensored,offsets);
 
-  psVector *tempModes = psVectorAlloc(N - Ncensored,PS_TYPE_F32);
-  bool success = psMMkmeansUncensored(Dcensored,
-                                      dim,
-                                      N - Ncensored,
-                                      tempModes,
-                                      means,
-                                      m,
-                                      iterations,
-                                      V);
+
+  bool success;
   if (Ncensored == 0) { // We can directly return the modes and means calculated
-    if (modes) {
-      psFree(modes);
-    }
-    modes = psMemIncrRefCounter(tempModes);
+    success = psMMkmeansUncensored(Dcensored,
+                                   dim,
+                                   N,
+                                   modes,
+                                   means,
+                                   m,
+                                   iterations,
+                                   V);
   }
   else {
+    psVector *tempModes = psVectorAlloc(N - Ncensored,PS_TYPE_F32);
+    success = psMMkmeansUncensored(Dcensored,
+                                   dim,
+                                   N - Ncensored,
+                                   tempModes,
+                                   means,
+                                   m,
+                                   iterations,
+                                   V);
+    
     if (!modes) {
       modes = psVectorAlloc(N,PS_TYPE_F32);
     }
-
+    
     psVectorInit(modes,NAN);
     for (int i = 0; i < N - Ncensored; i++) {
       modes->data.F32[i - offsets->data.S32[i]] = tempModes->data.F32[i];
     }
+    psFree(tempModes);
   }
   psFree(offsets);
   psFree(Dcensored);
-  psFree(tempModes);
+
   return(success);
 }
@@ -373 +381 @@
   // We need to scan the input for invalid d
   long Ncensored;
-  psVector *offsets = psVectorAlloc(N,PS_TYPE_F32);
+  psVector *offsets = psVectorAlloc(N,PS_TYPE_S32);
   psImage  *Dcensored = psMMCensorData(D,&Ncensored,offsets);
-
-  psVector *tempModes = psVectorAlloc(N - Ncensored,PS_TYPE_F32);
-  psImage  *tempP     = psImageAlloc(m,N - Ncensored,PS_TYPE_F32);
-  bool success = psMMGMMUncensored(Dcensored,
-                                   dim,
-                                   N - Ncensored,
-                                   tempModes,
-                                   means,
-                                   sigma,
-                                   pi,
-                                   tempP,
-                                   m,
-                                   iterations,
-                                   V);
+  bool success;
   if (Ncensored == 0) { // We can directly return the modes and means calculated
-    if (modes) {
-      psFree(modes);
-    }
-    if (P) {
-      psFree(P);
-    }
-    modes = psMemIncrRefCounter(tempModes);
-    P     = psMemIncrRefCounter(tempP);
+    success = psMMGMMUncensored(Dcensored,
+                                dim,
+                                N,
+                                modes,
+                                means,
+                                sigma,
+                                pi,
+                                P,
+                                m,
+                                iterations,
+                                V);
   }
   else {
+    psVector *tempModes = psVectorAlloc(N - Ncensored,PS_TYPE_F32);
+    psImage  *tempP     = psImageAlloc(m,N - Ncensored,PS_TYPE_F32);
+
+    success = psMMGMMUncensored(Dcensored,
+                                dim,
+                                N,
+                                tempModes,
+                                means,
+                                sigma,
+                                pi,
+                                tempP,
+                                m,
+                                iterations,
+                                V);
+      
+    
     if (!modes) {
       modes = psVectorAlloc(N,PS_TYPE_F32);
@@ -415 +429 @@
       }
     }
+    psFree(tempModes);
+    psFree(tempP);
   }
   psFree(offsets);
   psFree(Dcensored);
-  psFree(tempModes);
-  psFree(tempP);
   return(success);
 }