Index: trunk/archive/pslib/include/psFunctions.h
===================================================================
--- trunk/archive/pslib/include/psFunctions.h	(revision 217)
+++ trunk/archive/pslib/include/psFunctions.h	(revision 220)
@@ -31,4 +31,75 @@
 } psPolynomial2D;
 
+/** Three-dimensional polynomial */
+typedef struct {
+    int nX, nY, nZ;			//!< Number of terms in x, y and z
+    float *restrict *restrict *restrict coeff; //!< Coefficients
+    float *restrict *restrict *restrict coeffErr; //!< Error in coefficients
+    char *restrict *restrict *restrict mask; //!< Coefficients mask
+} psPolynomial3D;
+
+/** Four-dimensional polynomial */
+typedef struct {
+    int nW, nX, nY, nZ;			//!< Number of terms in w, x, y and z
+    float *restrict *restrict *restrict *restrict coeff; //!< Coefficients
+    float *restrict *restrict *restrict *restrict coeffErr; //!< Error in coefficients
+    char *restrict *restrict *restrict *restrict mask; //!< Coefficients mask
+} psPolynomial4D;
+
+
+/** Constructors */
+psPolynomial1D *psPolynomial1DAlloc(int n //!< Number of terms
+				    );
+psPolynomial2D *psPolynomial2DAlloc(int nX, int nY //!< Number of terms in x and y
+				    );
+psPolynomial3D *psPolynomial3DAlloc(int nX, int nY, int nZ //!< Number of terms in x, y and z
+				    );
+psPolynomial4D *psPolynomial4DAlloc(int nW, int nX, int nY, int nZ //!< Number of terms in w, x, y and z
+				    );
+
+/** Destructors */
+void psPolynomial1DFree(psPolynomial1D *restrict myPoly //!< Polynomial to destroy
+			);
+void psPolynomial2DFree(psPolynomial2D *restrict myPoly //!< Polynomial to destroy
+			);
+void psPolynomial3DFree(psPolynomial3D *restrict myPoly //!< Polynomial to destroy
+			);
+void psPolynomial4DFree(psPolynomial4D *restrict myPoly //!< Polynomial to destroy
+			);
+
+
+/** Evaluate 1D polynomial */
+float
+psEvalPolynomial1D(float x,		//!< Value at which to evaluate
+		   const psPolynomial1D *restrict myPoly //!< Coefficients for the polynomial
+		   );
+
+/** Evaluate 2D polynomial */
+float
+psEvalPolynomial2D(float x,		//!< Value x at which to evaluate
+		   float y,		//!< Value y at which to evaluate
+		   const psPolynomial2D *restrict myPoly //!< Coefficients for the polynomial
+		   );
+
+/** Evaluate 3D polynomial */
+float
+psEvalPolynomial3D(float x,		//!< Value x at which to evaluate
+		   float y,		//!< Value y at which to evaluate
+		   float z,		//!< Value z at which to evaluate
+		   const psPolynomial3D *restrict myPoly //!< Coefficients for the polynomial
+		   );
+
+/** Evaluate 4D polynomial */
+float
+psEvalPolynomial4D(float w,		//!< Value w at which to evaluate
+		   float x,		//!< Value x at which to evaluate
+		   float y,		//!< Value y at which to evaluate
+		   float z,		//!< Value z at which to evaluate
+		   const psPolynomial4D *restrict myPoly //!< Coefficients for the polynomial
+		   );
+
+/************************************************************************************************************/
+
+/* Double-precision polynomials, mainly for use in astrometry */
 
 /** Double-precision one-dimensional polynomial */
@@ -48,38 +119,41 @@
 } psDPolynomial2D;
 
+/** Double-precision three-dimensional polynomial */
+typedef struct {
+    int nX, nY, nZ;			//!< Number of terms in x, y and z
+    double *restrict *restrict *restrict coeff;	//!< Coefficients
+    double *restrict *restrict *restrict coeffErr; //!< Error in coefficients
+    char *restrict *restrict *restrict mask; //!< Coefficient mask
+} psDPolynomial3D;
+
+/** Double-precision four-dimensional polynomial */
+typedef struct {
+    int nW, nX, nY, nZ;			//!< Number of terms in w, x, y and z
+    double *restrict *restrict *restrict *restrict coeff; //!< Coefficients
+    double *restrict *restrict *restrict *restrict coeffErr; //!< Error in coefficients
+    char *restrict *restrict *restrict *restrict mask; //!< Coefficients mask
+} psDPolynomial4D;
 
 /** Constructors */
-psPolynomial1D *psPolynomial1DAlloc(int n //!< Number of terms
-				    );
-psPolynomial2D *psPolynomial2DAlloc(int nX, int nY //!< Number of terms in x and y
-				    );
 psPolynomial1D *psDPolynomial1DAlloc(int n //!< Number of terms
 				     );
 psPolynomial2D *psDPolynomial2DAlloc(int nX, int nY //!< Number of terms in x and y
 				     );
+psPolynomial3D *psDPolynomial3DAlloc(int nX, int nY, int nZ //!< Number of terms in x, y and z
+				     );
+psPolynomial4D *psDPolynomial4DAlloc(int nW, int nX, int nY, int nZ //!< Number of terms in w, x, y and z
+				     );
+
 
 /** Destructors */
-void psPolynomial1DFree(psPolynomial1D *restrict myPoly //!< Polynomial to destroy
-			);
-void psPolynomial2DFree(psPolynomial2D *restrict myPoly //!< Polynomial to destroy
-			);
 void psDPolynomial1DFree(psPolynomial1D *restrict myPoly //!< Polynomial to destroy
 			 );
-void psDPolynomial2DFree(psPolynomial1D *restrict myPoly //!< Polynomial to destroy
+void psDPolynomial2DFree(psPolynomial2D *restrict myPoly //!< Polynomial to destroy
+			 );
+void psDPolynomial3DFree(psPolynomial3D *restrict myPoly //!< Polynomial to destroy
+			 );
+void psDPolynomial4DFree(psPolynomial4D *restrict myPoly //!< Polynomial to destroy
 			 );
 
-
-/** Evaluate 1D polynomial */
-float
-psEvalPolynomial1D(float x,		//!< Value at which to evaluate
-		   const psPolynomial1D *restrict myPoly //!< Coefficients for the polynomial
-		   );
-
-/** Evaluate 2D polynomial */
-float
-psEvalPolynomial2D(float x,		//!< Value x at which to evaluate
-		   float y,		//!< Value y at which to evaluate
-		   const psPolynomial2D *restrict myPoly //!< Coefficients for the polynomial
-		   );
 
 /** Evaluate 1D polynomial (double precision) */
@@ -96,3 +170,21 @@
 		    );
 
+/** Evaluate 3D polynomial (double precision) */
+double
+psEvalDPolynomial3D(double x,		//!< Value x at which to evaluate
+		    double y,		//!< Value y at which to evaluate
+		    double z,		//!< Value z at which to evaluate
+		    const psDPolynomial3D *restrict myPoly //!< Coefficients for the polynomial
+		    );
+
+/** Evaluate 4D polynomial (double precision) */
+double
+psEvalDPolynomial4D(double w,		//!< Value w at which to evaluate
+		    double x,		//!< Value x at which to evaluate
+		    double y,		//!< Value y at which to evaluate
+		    double z,		//!< Value z at which to evaluate
+		    const psDPolynomial4D *restrict myPoly //!< Coefficients for the polynomial
+		    );
+
+
 #endif