Index: branches/eam_branches/ipp-20100823/psModules/src/imcombine/pmSubtractionEquation.c
===================================================================
--- branches/eam_branches/ipp-20100823/psModules/src/imcombine/pmSubtractionEquation.c	(revision 29489)
+++ branches/eam_branches/ipp-20100823/psModules/src/imcombine/pmSubtractionEquation.c	(revision 29490)
@@ -1104,5 +1104,14 @@
 	// XXX TEST: try some constraint on the svd solution
 	// solution = psMatrixSolveSVD(solution, sumMatrix, sumVector, NAN);
-	solution = psMatrixSolveSVD(solution, sumMatrix, sumVector, NAN);
+	// SINGLE solution
+	if (1) {
+	    solution = psMatrixSolveSVD(solution, sumMatrix, sumVector, NAN);
+	} else {
+	    psVector *PERM = NULL;
+	    psImage *LU = psMatrixLUDecomposition(NULL, &PERM, sumMatrix);
+	    solution = psMatrixLUSolution(solution, LU, sumVector, PERM);
+	    psFree (LU);
+	    psFree (PERM);
+	}
 
         for (int i = 0; i < solution->n; i++) {
@@ -1175,4 +1184,5 @@
         psVectorInit(solution, 0);
 
+	// DUAL solution
 	solution = psMatrixSolveSVD(solution, sumMatrix, sumVector, NAN);
 
Index: branches/eam_branches/ipp-20100823/psModules/src/imcombine/pmSubtractionKernels.c
===================================================================
--- branches/eam_branches/ipp-20100823/psModules/src/imcombine/pmSubtractionKernels.c	(revision 29489)
+++ branches/eam_branches/ipp-20100823/psModules/src/imcombine/pmSubtractionKernels.c	(revision 29490)
@@ -182,4 +182,7 @@
     return true;
 }
+
+// XXX *** this code used the central pixel to force zero net flux,
+// Alard actually uses kernel(0) for this purpose (for even order, kernel[i] = kernel'[i] - kernel[0])
 
 static bool pmSubtractionKernelPreCalcNormalize(pmSubtractionKernels *kernels, pmSubtractionKernelPreCalc *preCalc,