Index: trunk/archive/pslib/include/psFunctions.h
===================================================================
--- trunk/archive/pslib/include/psFunctions.h	(revision 1089)
+++ trunk/archive/pslib/include/psFunctions.h	(revision 1423)
@@ -26,6 +26,15 @@
 /************************************************************************************************************/
 
+/** Type of polynomial */
+typedef enum {
+    PS_POLYNOMIAL_ORD,			///< Ordinary polynomial
+    PS_POLYNOMIAL_CHEB			///< Chebyshev polynomial
+} psPolynomialType;
+
+/************************************************************************************************************/
+
 /** one-dimensional polynomial */
 typedef struct {
+    psPolynomialType type;		///< Polynomial type
     int n;				///< Number of terms
     float *restrict coeff;		///< Coefficients
@@ -36,4 +45,5 @@
 /** Two-dimensional polynomial */
 typedef struct {
+    psPolynomialType type;		///< Polynomial type
     int nX, nY;				///< Number of terms in x and y
     float *restrict *restrict coeff;	///< Coefficients
@@ -44,4 +54,5 @@
 /** Three-dimensional polynomial */
 typedef struct {
+    psPolynomialType type;		///< Polynomial type
     int nX, nY, nZ;			///< Number of terms in x, y and z
     float *restrict *restrict *restrict coeff; ///< Coefficients
@@ -52,4 +63,5 @@
 /** Four-dimensional polynomial */
 typedef struct {
+    psPolynomialType type;		///< Polynomial type
     int nW, nX, nY, nZ;			///< Number of terms in w, x, y and z
     float *restrict *restrict *restrict *restrict coeff; ///< Coefficients
@@ -66,13 +78,17 @@
 /** Constructor */
 psPolynomial1D *psPolynomial1DAlloc(int n ///< Number of terms
+				    psPolynomialType type ///< Polynomial type
 				    );
 /** Constructor */
 psPolynomial2D *psPolynomial2DAlloc(int nX, int nY ///< Number of terms in x and y
+				    psPolynomialType type ///< Polynomial type
 				    );
 /** Constructor */
 psPolynomial3D *psPolynomial3DAlloc(int nX, int nY, int nZ ///< Number of terms in x, y and z
+				    psPolynomialType type ///< Polynomial type
 				    );
 /** Constructor */
 psPolynomial4D *psPolynomial4DAlloc(int nW, int nX, int nY, int nZ ///< Number of terms in w, x, y and z
+				    psPolynomialType type ///< Polynomial type
 				    );
 
@@ -121,4 +137,35 @@
 		   );
 
+/** Evaluate 1D polynomial: vector input */
+psVector *
+psPolynomial1DEvalVector(psVector *x,	///< Value at which to evaluate
+			 const psPolynomial1D *restrict myPoly ///< Coefficients for the polynomial
+			 );
+
+/** Evaluate 2D polynomial: vector input */
+psVector *
+psPolynomial2DEvalVector(psVector *x,	///< Value x at which to evaluate
+			 psVector *y,	///< Value y at which to evaluate
+			 const psPolynomial2D *restrict myPoly ///< Coefficients for the polynomial
+			 );
+
+/** Evaluate 3D polynomial: vector input */
+psVector *
+psPolynomial3DEvalVector(psVector *x,	///< Value x at which to evaluate
+			 psVector *y,	///< Value y at which to evaluate
+			 psVector *z,	///< Value z at which to evaluate
+			 const psPolynomial3D *restrict myPoly ///< Coefficients for the polynomial
+			 );
+
+/** Evaluate 4D polynomial: vector input */
+psVector *
+psPolynomial4DEvalVector(psVector *w,	///< Value w at which to evaluate
+			 psVector *x,	///< Value x at which to evaluate
+			 psVector *y,	///< Value y at which to evaluate
+			 psVector *z,	///< Value z at which to evaluate
+			 const psPolynomial4D *restrict myPoly ///< Coefficients for the polynomial
+			 );
+
+
 /* \} */ // End of MathGroup Functions
 
@@ -129,4 +176,5 @@
 /** Double-precision one-dimensional polynomial */
 typedef struct {
+    psPolynomialType type;		///< Polynomial type
     int n;				///< Number of terms
     double *restrict coeff;		///< Coefficients
@@ -137,4 +185,5 @@
 /** Double-precision two-dimensional polynomial */
 typedef struct {
+    psPolynomialType type;		///< Polynomial type
     int nX, nY;				///< Number of terms in x and y
     double *restrict *restrict coeff;	///< Coefficients
@@ -145,4 +194,5 @@
 /** Double-precision three-dimensional polynomial */
 typedef struct {
+    psPolynomialType type;		///< Polynomial type
     int nX, nY, nZ;			///< Number of terms in x, y and z
     double *restrict *restrict *restrict coeff;	///< Coefficients
@@ -153,4 +203,5 @@
 /** Double-precision four-dimensional polynomial */
 typedef struct {
+    psPolynomialType type;		///< Polynomial type
     int nW, nX, nY, nZ;			///< Number of terms in w, x, y and z
     double *restrict *restrict *restrict *restrict coeff; ///< Coefficients
@@ -165,15 +216,19 @@
 
 /** Constructor */
-psDPolynomial1D *psDPolynomial1DAlloc(int n ///< Number of terms
-    );
+psDPolynomial1D *psDPolynomial1DAlloc(int n, ///< Number of terms
+				      psPolynomialType type ///< Polynomial type
+				      );
 /** Constructor */
 psDPolynomial2D *psDPolynomial2DAlloc(int nX, int nY ///< Number of terms in x and y
-    );
+				      psPolynomialType type ///< Polynomial type
+				      );
 /** Constructor */
 psDPolynomial3D *psDPolynomial3DAlloc(int nX, int nY, int nZ ///< Number of terms in x, y and z
-    );
+				      psPolynomialType type ///< Polynomial type
+				      );
 /** Constructor */
 psDPolynomial4D *psDPolynomial4DAlloc(int nW, int nX, int nY, int nZ ///< Number of terms in w, x, y and z
-    );
+				      psPolynomialType type ///< Polynomial type
+				      );
 
 
@@ -222,5 +277,74 @@
 		    );
 
+/** Evaluate 1D polynomial: vector input */
+psVector *
+psDPolynomial1DEvalVector(psVector *x,	///< Value at which to evaluate
+			  const psDPolynomial1D *restrict myPoly ///< Coefficients for the polynomial
+			  );
+
+/** Evaluate 2D polynomial: vector input */
+psVector *
+psDPolynomial2DEvalVector(psVector *x,	///< Value x at which to evaluate
+			  psVector *y,	///< Value y at which to evaluate
+			  const psDPolynomial2D *restrict myPoly ///< Coefficients for the polynomial
+			  );
+
+/** Evaluate 3D polynomial: vector input */
+psVector *
+psDPolynomial3DEvalVector(psVector *x,	///< Value x at which to evaluate
+			  psVector *y,	///< Value y at which to evaluate
+			  psVector *z,	///< Value z at which to evaluate
+			  const psDPolynomial3D *restrict myPoly ///< Coefficients for the polynomial
+			  );
+
+/** Evaluate 4D polynomial: vector input */
+psVector *
+psDPolynomial4DEvalVector(psVector *w,	///< Value w at which to evaluate
+			  psVector *x,	///< Value x at which to evaluate
+			  psVector *y,	///< Value y at which to evaluate
+			  psVector *z,	///< Value z at which to evaluate
+			  const psDPolynomial4D *restrict myPoly ///< Coefficients for the polynomial
+			  );
+
 /* \} */ // End of MathGroup Functions
 
+/************************************************************************************************************/
+
+/* Splines */
+
+/** A 1D cubic-spline */
+typedef struct {
+    int n;				///< Number of spline pieces
+    psPolynomial1D *spline;		///< Array of n splines
+    float *domains;			///< The boundaries between each spline piece.  Size is n+1.
+} psSpline1D;
+
+/** Functions **************************************************************/
+/** \addtogroup MathGroup Math Utilities
+ *  \{
+ */
+
+/** Constructors */
+psSpline1D *psSpline1DAlloc(int n,	///< Number of spline pieces
+			    int order,	///< Order of spline pieces
+			    float min,	///< Minimum data value
+			    float max	///< Maximum data value
+			    );
+psSpline1D *psSpline1DAllocGeneric(const psVector *bounds ///< Boundaries between each spline piece.  Number
+							  ///< of spline pieces can hence be inferred.
+				   int order	///< Order of spline pieces
+				   );
+/** Destructor */
+void psSpline1DFree(psSpline1D *mySpline ///< Spline to destroy
+		    );
+/** Evaluator */
+float psSpline1DEval(float x,		///< Coordinate at which to evaluate
+		     const psSpline1D *spline ///< Spline to evaluate
+		     );
+
+/** Evaluator: vector version */
+psVector *psSpline1DEvalVector(psVector *x, ///< Coordinates at which to evaluate
+			       const psSpline1D *spline ///< Spline to evaluate
+			       );
+
 #endif
