Index: trunk/psModules/src/detrend/pmPattern.c
===================================================================
--- trunk/psModules/src/detrend/pmPattern.c	(revision 33326)
+++ trunk/psModules/src/detrend/pmPattern.c	(revision 33340)
@@ -680,13 +680,15 @@
 	j = 8 * (x + 1) + y; // Determine that neighbor's index
 	if (fabs(C->data.F64[y][x]) > fabs(D->data.F64[y][x+1])) {
-	  critical_value = 2.0 * fabs(D->data.F64[y][x+1] / C->data.F64[y][x]);
+	  critical_value = 2.0 * fabs(D->data.F64[y][x+1]);
 	}
 	else {
-	  critical_value = 2.0 * fabs(C->data.F64[y][x] / D->data.F64[y][x+1]);
-	}
-	psTrace("psModules.detrend.cont",5,"CmD %d %d %d %d %g %g", // diagnostic
+	  critical_value = 2.0 * fabs(C->data.F64[y][x]);
+	}
+	if (critical_value < 25) { critical_value = 25; }
+	psTrace("psModules.detrend.cont",5,"CmD %d %d %d %d %g %g %g", // diagnostic
 		i,x,y,j,
 		C->data.F64[y][x],
-		D->data.F64[y][x+1]
+		D->data.F64[y][x+1],
+		critical_value
 		);
 	if (!(meanMask->data.PS_TYPE_VECTOR_MASK_DATA[j] & PM_PATTERN_IGNORE)&&  // If there are no errors with the neighbor,
@@ -701,15 +703,18 @@
       if (x - 1 > -1) { // etc.
 	j = 8 * (x - 1) + y;
-	psTrace("psModules.detrend.cont",5,"DmC %d %d %d %d %g %g",
+	if (fabs(C->data.F64[y][x-1]) > fabs(D->data.F64[y][x])) {
+	  critical_value = 2.0 * fabs(D->data.F64[y][x]);
+	}
+	else {
+	  critical_value = 2.0 * fabs(C->data.F64[y][x-1]);
+	}
+	if (critical_value < 25) { critical_value = 25; }
+	psTrace("psModules.detrend.cont",5,"DmC %d %d %d %d %g %g %g",
 		i,x,y,j,
 		D->data.F64[y][x],
-		C->data.F64[y][x-1]
+		C->data.F64[y][x-1],
+		critical_value
 		);
-	if (fabs(C->data.F64[y][x-1]) > fabs(D->data.F64[y][x])) {
-	  critical_value = 2.0 * fabs(D->data.F64[y][x] / C->data.F64[y][x-1]);
-	}
-	else {
-	  critical_value = 2.0 * fabs(C->data.F64[y][x-1] / D->data.F64[y][x]);
-	}
+
 	if (!(meanMask->data.PS_TYPE_VECTOR_MASK_DATA[j] & PM_PATTERN_IGNORE)&&
 	    (isfinite(D->data.F64[y][x]))&&(isfinite(C->data.F64[y][x-1]))&&
@@ -729,9 +734,10 @@
 		);
 	if (fabs(A->data.F64[y][x]) > fabs(B->data.F64[y+1][x])) {
-	  critical_value = 2.0 * fabs(B->data.F64[y+1][x] / A->data.F64[y][x]);
+	  critical_value = 2.0 * fabs(B->data.F64[y+1][x]);
 	}
 	else {
-	  critical_value = 2.0 * fabs(A->data.F64[y][x] / B->data.F64[y+1][x]);
-	}
+	  critical_value = 2.0 * fabs(A->data.F64[y][x]);
+	}
+	if (critical_value < 25) { critical_value = 25; }
 	if (!(meanMask->data.PS_TYPE_VECTOR_MASK_DATA[j] & PM_PATTERN_IGNORE)&&
 	    (isfinite(A->data.F64[y][x]))&&(isfinite(B->data.F64[y+1][x]))&&
@@ -751,9 +757,10 @@
 		);
 	if (fabs(A->data.F64[y-1][x]) > fabs(B->data.F64[y][x])) {
-	  critical_value = 2.0 * fabs(B->data.F64[y][x] / A->data.F64[y-1][x]);
+	  critical_value = 2.0 * fabs(B->data.F64[y][x]);
 	}
 	else {
-	  critical_value = 2.0 * fabs(A->data.F64[y-1][x] / B->data.F64[y][x]);
-	}
+	  critical_value = 2.0 * fabs(A->data.F64[y-1][x]);
+	}
+	if (critical_value < 25) { critical_value = 25; }
 	if (!(meanMask->data.PS_TYPE_VECTOR_MASK_DATA[j] & PM_PATTERN_IGNORE)&&
 	    (isfinite(B->data.F64[y][x]))&&(isfinite(A->data.F64[y-1][x]))&&
@@ -766,10 +773,38 @@
       }
     }
+    double max_XX = 0;
+    double solution_V = 0;
+    int i_peak = -1;
     for (int i = 0; i < numCells; i++) { // If any cells have no value of themself, set the matrix to 1.0.
       if (XX->data.F64[i][i] == 0.0) {
 	XX->data.F64[i][i] = 1.0;
       }
-    }
-#if (0)
+      if (XX->data.F64[i][i] > max_XX) {
+	max_XX = XX->data.F64[i][i];
+	solution_V = solution->data.F64[i];
+	i_peak = i;
+      }
+    }
+    psTrace("psModules.detrend.cont",5,"fixed point: %d %g\n",
+	    i_peak,solution_V);
+
+    for (int i = 0; i < numCells; i++) {
+/*        if (!((XX->data.F64[i][i] == 1.0)&& */
+/*  	    (solution->data.F64[i] == 0.0))) { */
+	solution->data.F64[i] -= solution_V;
+	if (i != i_peak) {
+	  for (int j = 0; j < numCells; j++) {
+	    XX->data.F64[i][j] -= XX->data.F64[i_peak][j];
+	  }
+	}
+/*        } */
+    }
+    for (int i = 0; i < numCells; i++) {
+      XX->data.F64[i_peak][i] = 0.0;
+    }
+    XX->data.F64[i_peak][i_peak] = 1.0;
+    
+    
+#if (1)
     for (int i = 0; i < numCells; i++) { // print matrix A
       psTrace("psModules.detrend.cont",5,"A: %3d % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f % 2.0f",
@@ -805,7 +840,7 @@
     //    psMatrixLUSolve(XX,solution);
     psMatrixGJSolve(XX,solution);
-#if (0)
+#if (1)
     for (int i = 0; i < numCells; i++) { // print vector b
-      psTrace("psModules.detrend.cont",5,"b: %d %f",
+      psTrace("psModules.detrend.cont",5,"x: %d %f",
 	      i,
 	      solution->data.F64[i]