Index: trunk/stac/src/stacInvertMaps.c
===================================================================
--- trunk/stac/src/stacInvertMaps.c	(revision 2897)
+++ trunk/stac/src/stacInvertMaps.c	(revision 3375)
@@ -1,3 +1,4 @@
 #include <stdio.h>
+#include <assert.h>
 #include "pslib.h"
 #include "stac.h"
@@ -24,16 +25,6 @@
 	psPlaneTransform *oldMap = (psPlaneTransform*)maps->data[mapNum]; // Uninverted map
 	// Check input
-	if (oldMap->x->nX != oldMap->x->nY) {
-	    psError("stac.invertMaps", "Polynomial order in x and y don't match for map %d\n", mapNum);
-	    return NULL;
-	}
-	if (oldMap->y->nX != oldMap->y->nY) {
-	    psError("stac.invertMaps", "Polynomial order in x and y don't match for map %d\n", mapNum);
-	    return NULL;
-	}
-	if (oldMap->x->nX != oldMap->y->nX) {
-	    psError("stac.invertMaps", "Polynomial order in x and y don't match for map %d\n", mapNum);
-	    return NULL;
-	}
+	assert(oldMap->x->nX == oldMap->x->nY && oldMap->y->nX == oldMap->y->nY &&
+	       oldMap->x->nX == oldMap->y->nX);
 	int order = oldMap->x->nX;	// Polynomial order
 	psTrace("stac.invertMaps", 4, "Generating order %d polynomial inverse transformation.\n", order);
@@ -91,5 +82,5 @@
 
 		    fakePoly->mask[i][j] = 0;
-		    double ijPoly = psDPolynomial2DEval(xIn->data.F32[g], yIn->data.F32[g], fakePoly);
+		    double ijPoly = psDPolynomial2DEval(fakePoly, xIn->data.F32[g], yIn->data.F32[g]);
 		    fakePoly->mask[i][j] = 1;
 
@@ -98,5 +89,5 @@
 			    
 			    fakePoly->mask[m][n] = 0;
-			    double mnPoly = psDPolynomial2DEval(xIn->data.F32[g], yIn->data.F32[g], fakePoly);
+			    double mnPoly = psDPolynomial2DEval(fakePoly, xIn->data.F32[g], yIn->data.F32[g]);
 			    fakePoly->mask[m][n] = 1;
 
@@ -113,5 +104,5 @@
 	// Solution via LU Decomposition
 	psVector *permutation = psVectorAlloc(nCoeff, PS_TYPE_F64); // Permutation vector for LU Decomposition
-	psImage *luMatrix = psMatrixLUD(NULL, permutation, matrix); // LU decomposed matrix
+	psImage *luMatrix = psMatrixLUD(NULL, &permutation, matrix); // LU decomposed matrix
 	psVector *xSolution = psMatrixLUSolve(NULL, luMatrix, xVector, permutation); // Solution in x
 	psVector *ySolution = psMatrixLUSolve(NULL, luMatrix, yVector, permutation); // Solution in y
@@ -178,14 +169,8 @@
 #if 0
 	// Can't handle higher order than linear yet
-	if (((psPlaneTransform*)maps->data[i])->x->nX != 2 ||
-	    ((psPlaneTransform*)maps->data[i])->x->nY != 2 ||
-	    ((psPlaneTransform*)maps->data[i])->y->nX != 2 ||
-	    ((psPlaneTransform*)maps->data[i])->y->nY != 2) {
-	    psError("stac.invertMaps",
-		    "STAC cannot currently support orders other than linear.\n");
-	    psFree(inverted);
-	    return NULL;
-	}
-
+	assert(((psPlaneTransform*)maps->data[i])->x->nX == 2 &&
+	       ((psPlaneTransform*)maps->data[i])->x->nY == 2 &&
+	       ((psPlaneTransform*)maps->data[i])->y->nX == 2 &&
+	       ((psPlaneTransform*)maps->data[i])->y->nY == 2);
 	psPlaneTransform *newMap = psPlaneTransformAlloc(2, 2); // Inverted map
 	psPlaneTransform *oldMap = (psPlaneTransform*)maps->data[i]; // Uninverted map
