Index: trunk/psModules/src/detrend/pmPattern.c
===================================================================
--- trunk/psModules/src/detrend/pmPattern.c	(revision 33299)
+++ trunk/psModules/src/detrend/pmPattern.c	(revision 33304)
@@ -565,4 +565,8 @@
     psImage *C = psImageAlloc(8,8,PS_TYPE_F64); // Right edge
     psImage *D = psImageAlloc(8,8,PS_TYPE_F64); // Left edge
+    psImageInit(A,0.0);
+    psImageInit(B,0.0);
+    psImageInit(C,0.0);
+    psImageInit(D,0.0);
     
     for (int i = 0; i < numCells; i++) {
@@ -611,5 +615,5 @@
             continue; // Move on to next edge, as only part of this cell may be a problem
 	  }
-
+ 
 	  // If the returned value is zero, assume something is wrong.  Do I still need this?
 	  if (psStatsGetValue(bgStats,bgStat) < 1e-6) {
@@ -632,4 +636,12 @@
 	    else if (j == 3) { D->data.F64[y][x] = psStatsGetValue(bgStats,bgStat); }
 	  }
+
+	  for (int u = 0; u < subset->numCols; u++) {
+	    for (int v = 0; v < subset->numRows; v++) {
+	      psTrace("psModules.detrend.cont",10,"BOX: %d %d (%d %d) (%d %d) %f %d",
+		      i,j,x,y,u,v,subset->data.F32[v][u],submask->data.PS_TYPE_IMAGE_MASK_DATA[v][u]);
+	    }
+	  }	  
+	  
 	  psFree(subset);
 	  psFree(submask);
@@ -638,5 +650,5 @@
 	psTrace("psModules.detrend.cont",5, "OTA: %d (%d %d) A: %f B: %f C: %f D: %f",
 		i,x,y,
-		A->data.F32[y][x],B->data.F32[y][x],C->data.F32[y][x],D->data.F32[y][x]);		
+		A->data.F64[y][x],B->data.F64[y][x],C->data.F64[y][x],D->data.F64[y][x]);		
     }
     psFree(bgStats);
@@ -661,4 +673,5 @@
       int x = (i - y) / 8;
       int j;
+      double critical_value = 0.0;
       if (meanMask->data.PS_TYPE_VECTOR_MASK_DATA[i] & PM_PATTERN_IGNORE) {
 	continue;
@@ -666,4 +679,10 @@
       if (x + 1 < 8) {  // We have a neighbor adjacent in the +x direction
 	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]);
+	}
+	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
 		i,x,y,j,
@@ -672,5 +691,7 @@
 		);
 	if (!(meanMask->data.PS_TYPE_VECTOR_MASK_DATA[j] & PM_PATTERN_IGNORE)&&  // If there are no errors with the neighbor,
-	    (isfinite(C->data.F64[y][x]))&&(isfinite(D->data.F64[y][x+1]))) {    // and all edges have valid values:
+	    (isfinite(C->data.F64[y][x]))&&(isfinite(D->data.F64[y][x+1]))&&     // and all edges have valid values,
+	    (fabs(C->data.F64[y][x] - D->data.F64[y][x+1]) < critical_value)     // and there are no large discontinuities,
+	    ) {    
 	  solution->data.F64[i] += C->data.F64[y][x] - D->data.F64[y][x+1];     // Take the difference
 	  XX->data.F64[i][i] += 1;                                              // increment our relation with ourself
@@ -685,7 +706,14 @@
 		C->data.F64[y][x-1]
 		);
-
+	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]))) {
+	    (isfinite(D->data.F64[y][x]))&&(isfinite(C->data.F64[y][x-1]))&&
+	    (fabs(D->data.F64[y][x] - C->data.F64[y][x-1]) < critical_value)
+	    ) {
 	  solution->data.F64[i] += D->data.F64[y][x] - C->data.F64[y][x-1];
 	  XX->data.F64[i][i] += 1;
@@ -700,7 +728,14 @@
 		B->data.F64[y+1][x]
 		);
-
+	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]);
+	}
+	else {
+	  critical_value = 2.0 * fabs(A->data.F64[y][x] / B->data.F64[y+1][x]);
+	}
 	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]))) {
+	    (isfinite(A->data.F64[y][x]))&&(isfinite(B->data.F64[y+1][x]))&&
+	    (fabs(A->data.F64[y][x] - B->data.F64[y+1][x]) < critical_value)
+	    ) {
 	  solution->data.F64[i] += A->data.F64[y][x] - B->data.F64[y+1][x];
 	  XX->data.F64[i][i] += 1;
@@ -715,7 +750,14 @@
 		A->data.F64[y-1][x]
 		);
-
+	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]);
+	}
+	else {
+	  critical_value = 2.0 * fabs(A->data.F64[y-1][x] / B->data.F64[y][x]);
+	}
 	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]))) {
+	    (isfinite(B->data.F64[y][x]))&&(isfinite(A->data.F64[y-1][x]))&&
+	    (fabs(B->data.F64[y][x] - A->data.F64[y-1][x]) < critical_value)
+	    ) {
 	  solution->data.F64[i] += B->data.F64[y][x] - A->data.F64[y-1][x];
 	  XX->data.F64[i][i] += 1;