Index: /trunk/psLib/src/dataManip/psFunctions.h
===================================================================
--- /trunk/psLib/src/dataManip/psFunctions.h	(revision 2268)
+++ /trunk/psLib/src/dataManip/psFunctions.h	(revision 2269)
@@ -12,6 +12,6 @@
 *  @author George Gusciora, MHPCC
 *
-*  @version $Revision: 1.31 $ $Name: not supported by cvs2svn $
-*  @date $Date: 2004-10-27 00:57:31 $
+*  @version $Revision: 1.32 $ $Name: not supported by cvs2svn $
+*  @date $Date: 2004-11-03 03:30:30 $
 *
 *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -409,12 +409,4 @@
                                const psSpline1D *spline);
 
-psS32 p_psVectorBinDisectF32(float *bins,
-                             psS32 numBins,
-                             float x);
-
-psS32 p_psVectorBinDisectS32(psS32 *bins,
-                             psS32 numBins,
-                             psS32 x);
-
 psS32 p_psVectorBinDisect(psVector *bins,
                           psScalar *x);
Index: /trunk/psLib/src/dataManip/psMinimize.c
===================================================================
--- /trunk/psLib/src/dataManip/psMinimize.c	(revision 2268)
+++ /trunk/psLib/src/dataManip/psMinimize.c	(revision 2269)
@@ -9,6 +9,6 @@
  *  @author GLG, MHPCC
  *
- *  @version $Revision: 1.83 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2004-11-02 19:08:33 $
+ *  @version $Revision: 1.84 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2004-11-03 03:30:30 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -61,10 +61,9 @@
 returned as a psVector sums.
  
-XXX: change name
 XXX: Use a static vector.
  *****************************************************************************/
-void p_psBuildSums1D(double x,
-                     psS32 polyOrder,
-                     psVector* sums)
+void psBuildSums1D(double x,
+                   psS32 polyOrder,
+                   psVector* sums)
 {
     psS32 i = 0;
@@ -74,5 +73,4 @@
         sums = psVectorAlloc(polyOrder, PS_TYPE_F64);
     }
-
     if (polyOrder > sums->n) {
         sums = psVectorRealloc(sums, polyOrder);
@@ -109,5 +107,4 @@
     psS32 i;
     psS32 k;
-    psBool mustFreeX = false;
     float sig;
     float p;
@@ -150,9 +147,5 @@
     }
 
-    if (mustFreeX == true) {
-        psFree(X);
-    }
     psFree(u);
-
     psTrace(".psLib.dataManip.CalculateSecondDerivs", 4,
             "---- CalculateSecondDerivs() end ----\n");
@@ -171,7 +164,7 @@
 polynomials which are stored in psSpline1D.
  
-XXX: check types/sizes
 XXX: This is F32 only
  *****************************************************************************/
+/*
 float p_psNRSpline1DEval(psSpline1D *spline,
                          const psVector* restrict x,
@@ -187,5 +180,5 @@
     PS_VECTOR_CHECK_NULL(y, NAN);
     PS_VECTOR_CHECK_TYPE(y, PS_TYPE_F32, NAN);
-
+ 
     psS32 n;
     psS32 klo;
@@ -197,7 +190,7 @@
     float D;
     float Y;
-
+ 
     n = spline->n;
-    klo = p_psVectorBinDisectF32(spline->domains, (spline->n)+1, X);
+    klo = p_psVectorBinDisect32(spline->domains, (spline->n)+1, X);
     khi = klo + 1;
     H = (spline->domains)[khi] - (spline->domains)[klo];
@@ -206,13 +199,13 @@
     C = ((A*A*A)-A) * (H*H/6.0);
     D = ((B*B*B)-B) * (H*H/6.0);
-
+ 
     Y = (A * y->data.F32[klo]) +
         (B * y->data.F32[khi]) +
         (C * (spline->p_psDeriv2)[klo]) +
         (D * (spline->p_psDeriv2)[khi]);
-
+ 
     return(Y);
 }
-
+*/
 /*****************************************************************************/
 /* FUNCTION IMPLEMENTATION - PUBLIC                                          */
@@ -239,4 +232,5 @@
  
 XXX: usage of yErr is not specified in IfA documentation.
+ 
 XXX: Is the x argument redundant?  What do we do if the x argument is
 supplied, but does not equal the domains specified in mySpline?
@@ -414,5 +408,4 @@
 /******************************************************************************
 XXX: We assume unnormalized gaussians.
-XXX: Currently, yErr is ignored.
  *****************************************************************************/
 psVector *psMinimizeLMChi2Gauss1D(psImage *deriv,
@@ -441,6 +434,4 @@
         x = ((psVector *) (coords->data[i]))->data.F32[0];
         out->data.F32[i] = psGaussian(x, mean, stdev, false);
-        //        psTrace(".psLib.dataManip.psMinimize", 6,
-        //            "out->data.F32[%d]: (x, m, s) is (%f, %f, %f) is %f\n", i, x, mean, stdev, out->data.F32[i]);
     }
 
@@ -449,5 +440,4 @@
         float tmp = (x - mean) * psGaussian(x, mean, stdev, false);
         deriv->data.F32[i][0] = tmp / (stdev * stdev);
-
         tmp = (x - mean) * (x - mean) *
               psGaussian(x, mean, stdev, 0);
@@ -482,5 +472,5 @@
 XXX: This must work for both F32 and F64.  F32 is currently implemented.
      Note: since the LUD routines are only implemented in F64, then we
->     will have to convert all F32 input vectors to F64 regardless.  So,
+     will have to convert all F32 input vectors to F64 regardless.  So,
      the F64 port might be.
  *****************************************************************************/
@@ -730,7 +720,7 @@
 
 /******************************************************************************
-p_psVectorFitPolynomial1DCheb():  This routine will fit a Chebyshev
-polynomial of degree myPoly to the data points (x, y) and return the
-coefficients of that polynomial.
+VectorFitPolynomial1DCheb():  This routine will fit a Chebyshev polynomial of
+degree myPoly to the data points (x, y) and return the coefficients of that
+polynomial.
  
 XXX: yErr is currently ignored.
@@ -738,5 +728,5 @@
 XXX: Use private name?
 *****************************************************************************/
-psPolynomial1D *p_psVectorFitPolynomial1DCheby(psPolynomial1D* myPoly,
+psPolynomial1D *VectorFitPolynomial1DCheby(psPolynomial1D* myPoly,
         const psVector* restrict x,
         const psVector* restrict y,
@@ -785,5 +775,5 @@
         psFree(fScalar);
 
-        psTrace(".psLib.dataManip.p_psVectorFitPolynomial1DCheby", 6,
+        psTrace(".psLib.dataManip.VectorFitPolynomial1DCheby", 6,
                 "(x, X, y, f(X)) is (%f, %f, %f, %f)\n",
                 x->data.F64[i], X, y->data.F64[i], f->data.F64[i]);
@@ -809,12 +799,12 @@
 
 /******************************************************************************
-p_psVectorFitPolynomial1DOrd():  This routine will fit an ordinary
-polynomial of degree myPoly to the data points (x, y) and return the
-coefficients of that polynomial.
+VectorFitPolynomial1DOrd():  This routine will fit an ordinary polynomial of
+degree myPoly to the data points (x, y) and return the coefficients of that
+polynomial.
  
 XXX: Use private name?
 XXX: Use recycled vectors.
  *****************************************************************************/
-psPolynomial1D* p_psVectorFitPolynomial1DOrd(psPolynomial1D* myPoly,
+psPolynomial1D* VectorFitPolynomial1DOrd(psPolynomial1D* myPoly,
         const psVector* restrict x,
         const psVector* restrict y,
@@ -833,7 +823,7 @@
     psVector* xSums = NULL;
 
-    psTrace(".psLib.dataManip.psVectorFitPolynomial1D", 4,
-            "---- psVectorFitPolynomial1D() begin ----\n");
-    // printf("psVectorFitPolynomial1D()\n");
+    psTrace(".psLib.dataManip.VectorFitPolynomial1DOrd", 4,
+            "---- VectorFitPolynomial1DOrd() begin ----\n");
+    // printf("VectorFitPolynomial1D()\n");
     // for (i=0;i<x->n;i++) {
     // printf("(x, y, yErr) is (%f, %f, %f)\n", x->data.F64[i], y->data.F64[i], yErr->data.F64[i]);
@@ -865,5 +855,5 @@
     if (yErr == NULL) {
         for (i = 0; i < X->n; i++) {
-            p_psBuildSums1D(X->data.F64[i], 2 * polyOrder, xSums);
+            psBuildSums1D(X->data.F64[i], 2 * polyOrder, xSums);
 
             for (k = 0; k < polyOrder; k++) {
@@ -879,5 +869,5 @@
     } else {
         for (i = 0; i < X->n; i++) {
-            p_psBuildSums1D(X->data.F64[i], 2 * polyOrder, xSums);
+            psBuildSums1D(X->data.F64[i], 2 * polyOrder, xSums);
 
             for (k = 0; k < polyOrder; k++) {
@@ -916,6 +906,6 @@
     psFree(xSums);
 
-    psTrace(".psLib.dataManip.psVectorFitPolynomial1D", 4,
-            "---- psVectorFitPolynomial1D() begin ----\n");
+    psTrace(".psLib.dataManip.VectorFitPolynomial1DOrd", 4,
+            "---- VectorFitPolynomial1DOrd() begin ----\n");
     return (myPoly);
 }
@@ -962,5 +952,5 @@
         PS_VECTOR_GEN_X_INDEX_STATIC_F64(x64Static, y->n);
         if (myPoly->type == PS_POLYNOMIAL_CHEB) {
-            p_psNormalizeVector(x64Static);
+            p_psNormalizeVectorRange(x64Static, -1.0, 1.0);
         }
         x64 = x64Static;
@@ -974,7 +964,7 @@
     // Call the appropriate vector fitting routine.
     if (myPoly->type == PS_POLYNOMIAL_CHEB) {
-        p_psVectorFitPolynomial1DCheby(myPoly, x64, y64, yErr64);
+        VectorFitPolynomial1DCheby(myPoly, x64, y64, yErr64);
     } else if (myPoly->type == PS_POLYNOMIAL_ORD) {
-        tmpPoly = p_psVectorFitPolynomial1DOrd(myPoly, x64, y64, yErr64);
+        tmpPoly = VectorFitPolynomial1DOrd(myPoly, x64, y64, yErr64);
     } else {
         psError(__func__, "unknown polynomial type.\n");
Index: /trunk/psLib/src/dataManip/psStats.c
===================================================================
--- /trunk/psLib/src/dataManip/psStats.c	(revision 2268)
+++ /trunk/psLib/src/dataManip/psStats.c	(revision 2269)
@@ -9,6 +9,6 @@
  *  @author GLG, MHPCC
  *
- *  @version $Revision: 1.83 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2004-11-02 19:13:03 $
+ *  @version $Revision: 1.84 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2004-11-03 03:30:30 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -403,16 +403,18 @@
  *****************************************************************************/
 void p_psVectorSampleMedian(const psVector* restrict myVector,
-                            const psVector* restrict maskVector, psU32 maskVal, psStats* stats)
+                            const psVector* restrict maskVector,
+                            psU32 maskVal,
+                            psStats* stats)
 {
     psVector* unsortedVector = NULL;    // Temporary vector
     psVector* sortedVector = NULL;      // Temporary vector
-    psS32 i = 0;                  // Loop index variable
-    psS32 count = 0;              // # of points in this mean?
-    psS32 nValues = 0;            // # of points in vector
-    float rangeMin = 0.0;       // Exclude data below this
-    float rangeMax = 0.0;       // Exclude date above this
+    psS32 i = 0;                        // Loop index variable
+    psS32 count = 0;                    // # of points in this mean?
+    psS32 nValues = 0;                  // # of points in vector
+    float rangeMin = 0.0;               // Exclude data below this
+    float rangeMax = 0.0;               // Exclude date above this
 
     // Determine how many data points fit inside this min/max range
-    // and are not masked, IF the maskVector is not NULL>
+    // and are not masked, if the maskVector is not NULL>
     nValues = p_psVectorNValues(myVector, maskVector, maskVal, stats);
 
@@ -551,10 +553,12 @@
  *****************************************************************************/
 void p_psVectorSampleQuartiles(const psVector* restrict myVector,
-                               const psVector* restrict maskVector, psU32 maskVal, psStats* stats)
+                               const psVector* restrict maskVector,
+                               psU32 maskVal,
+                               psStats* stats)
 {
     psVector* unsortedVector = NULL;    // Temporary vector
     psVector* sortedVector = NULL;      // Temporary vector
     psS32 i = 0;                  // Loop index variable
-    psS32 count = 0;              // # of points in this mean?
+    psS32 count = 0;              // # of points in this mean.
     psS32 nValues = 0;            // # data points
     float rangeMin = 0.0;       // Exclude data below this
@@ -631,5 +635,7 @@
  *****************************************************************************/
 void p_psVectorSampleStdev(const psVector* restrict myVector,
-                           const psVector* restrict maskVector, psU32 maskVal, psStats* stats)
+                           const psVector* restrict maskVector,
+                           psU32 maskVal,
+                           psStats* stats)
 {
     psS32 i = 0;                  // Loop index variable
@@ -713,8 +719,10 @@
     stats
 Returns
-    NULL
- *****************************************************************************/
-void p_psVectorClippedStats(const psVector* restrict myVector,
-                            const psVector* restrict maskVector, psU32 maskVal, psStats* stats)
+    0 for success.
+ *****************************************************************************/
+int p_psVectorClippedStats(const psVector* restrict myVector,
+                           const psVector* restrict maskVector,
+                           psU32 maskVal,
+                           psStats* stats)
 {
     psS32 i = 0;                  // Loop index variable
@@ -726,11 +734,13 @@
     psVector* tmpMask = NULL;   // Temporary vector
 
-    // Endure that stats->clipIter is within the proper range.
+    // Ensure that stats->clipIter is within the proper range.
     if (!((PS_CLIPPED_NUM_ITER_LB <= stats->clipIter) && (stats->clipIter <= PS_CLIPPED_NUM_ITER_UB))) {
         psError(__func__, "Unallowed value for clipIter (%d).\n", stats->clipIter);
-    }
-    // Endure that stats->clipSigma is within the proper range.
+        return(-1);
+    }
+    // Ensure that stats->clipSigma is within the proper range.
     if (!((PS_CLIPPED_SIGMA_LB <= stats->clipSigma) && (stats->clipSigma <= PS_CLIPPED_SIGMA_UB))) {
         psError(__func__, "Unallowed value for clipSigma (%f).\n", stats->clipSigma);
+        return(-1);
     }
     // We allocate a temporary mask vector since during the iterative
@@ -798,7 +808,12 @@
 
     psFree(tmpMask);
-}
-
-void p_psNormalizeVectorF32_0(psVector* myData)
+    return(0);
+}
+
+/*****************************************************************************
+ *****************************************************************************/
+void p_psNormalizeVectorRangeF32(psVector* myData,
+                                 float outLow,
+                                 float outHigh)
 {
     float min = (float)HUGE;
@@ -815,22 +830,11 @@
     }
 
-    //  float range = max - min;
-    //    for (i = 0; i < myData->n; i++) {
-    //        myData->data.F32[i] = (myData->data.F32[i] - min) / range;
-    //    }
     for (i = 0; i < myData->n; i++) {
-        myData->data.F32[i] = (myData->data.F32[i] - min) / (max - min);
-    }
-}
-
-
-
-/*****************************************************************************
-p_psNormalizeVectorF32(myData): this is a private function which normalizes the
-elements of a vector to a range between 0.0 and 1.0.
- 
-XXX: how about an arbitrary p_psNormalizeVectorF32(myData, min, max)?
- *****************************************************************************/
-void p_psNormalizeVectorF32(psVector* myData)
+        myData->data.F32[i] = outLow + (myData->data.F32[i] - min) * (outHigh - outLow) / (max - min);
+    }
+}
+void p_psNormalizeVectorRangeF64(psVector* myData,
+                                 float outLow,
+                                 float outHigh)
 {
     float min = (float)HUGE;
@@ -840,63 +844,28 @@
     for (i = 0; i < myData->n; i++) {
         if (myData->data.F32[i] < min) {
-            min = myData->data.F32[i];
+            min = myData->data.F64[i];
         }
         if (myData->data.F32[i] > max) {
-            max = myData->data.F32[i];
-        }
-    }
-
-    //  float range = max - min;
-    //    for (i = 0; i < myData->n; i++) {
-    //        myData->data.F32[i] = (myData->data.F32[i] - min) / range;
-    //    }
+            max = myData->data.F64[i];
+        }
+    }
+
     for (i = 0; i < myData->n; i++) {
-        myData->data.F32[i] = (myData->data.F32[i] - 0.5 * (min + max)) /
-                              (0.5 * (max - min));
+        myData->data.F64[i] = outLow + (myData->data.F64[i] - min) * (outHigh - outLow) / (max - min);
     }
 }
 
 /*****************************************************************************
-p_psNormalizeVectorF64(myData): this is a private function which normalizes the
-elements of a vector to a range between -1.0 and 1.0.
-XXX: 0-1 or -1:1?
- 
-XXX: how about an arbitrary p_psNormalizeVectorF64(myData, min, max)?
- *****************************************************************************/
-void p_psNormalizeVectorF64(psVector* myData)
-{
-    float min = (double)HUGE;
-    float max = (double)-HUGE;
-    double range = 0.0;
-    psS32 i = 0;
-
-    for (i = 0; i < myData->n; i++) {
-        if (myData->data.F64[i] < min) {
-            min = myData->data.F64[i];
-        }
-        if (myData->data.F64[i] > max) {
-            max = myData->data.F64[i];
-        }
-    }
-
-    range = max - min;
-    for (i = 0; i < myData->n; i++) {
-        myData->data.F64[i] = -1.0 + 2.0 *
-                              ((myData->data.F64[i] - min) / range);
-    }
-
-    //    for (i = 0; i < myData->n; i++) {
-    //        myData->data.F64[i] = (myData->data.F64[i] - 0.5 * (min + max)) /
-    //                              (0.5 * (max - min));
-    //    }
-}
-
-// XXX: how about an arbitrary p_psNormalizeVector(myData, min, max)?
-void p_psNormalizeVector(psVector* myData)
+psNormalizeVectorRange(myData, low, high): this is a private function which
+normalizes the elements of a vector to a range between low and high.
+ *****************************************************************************/
+void psNormalizeVectorRange(psVector* myData,
+                            float low,
+                            float high)
 {
     if (myData->type.type == PS_TYPE_F32) {
-        p_psNormalizeVectorF32(myData);
+        p_psNormalizeVectorRangeF32(myData, low, high);
     } else if (myData->type.type == PS_TYPE_F64) {
-        p_psNormalizeVectorF64(myData);
+        p_psNormalizeVectorRangeF64(myData, low, high);
     } else {
         psError(__func__, "Unalowable data type.\n");
@@ -914,4 +883,6 @@
  
 XXX: Terminate when f(x)-getThisValue is within some error tolerance.
+ 
+XXX: Solve for X analytically.
  *****************************************************************************/
 float p_ps1DPolyMedian(psPolynomial1D* myPoly,
@@ -924,6 +895,4 @@
     float oldMidpoint = 1.0;
     float f = 0.0;
-
-    // printf("p_ps1DPolyMedian(%f, %f, %f) \n", rangeLow, rangeHigh, getThisValue);
 
     while (numIterations < PS_POLY_MEDIAN_MAX_ITERATIONS) {
@@ -956,4 +925,6 @@
 and i+1 is used for x[i]).  It then determines for what value x does that
 quadratic f(x) = yVal (the input parameter).
+ 
+XXX: After you fit the polynomial, solve for X analytically.
 *****************************************************************************/
 float fitQuadraticSearchForYThenReturnX(psVector *xVec,
@@ -1034,5 +1005,5 @@
     stats
 Returns
-    NULL
+    0 on success.
 *****************************************************************************/
 int p_psVectorRobustStats(const psVector* restrict myVector,
@@ -1119,8 +1090,4 @@
                             tmpStats->clippedStdev / 4.0f);
 
-    // The following was necessary to fit a gaussian to the data, since
-    // gaussian functions produce data between 0.0 and 1.0.
-    // p_psNormalizeVectorF32(robustHistogramVector);
-
     /**************************************************************************
     Determine the median/lower/upper quartile bin numbers.
@@ -1207,5 +1174,5 @@
     psVector *y = psVectorAlloc(robustHistogramVector->n, PS_TYPE_F32);
 
-    p_psNormalizeVectorF32_0(robustHistogramVector);
+    p_psNormalizeVectorRange(robustHistogramVector, 0.0, 1.0);
     for (i=0;i<robustHistogramVector->n;i++) {
         myCoords->data[i] = (psPtr *) psVectorAlloc(2, PS_TYPE_F32);
@@ -1620,6 +1587,4 @@
     }
     // ************************************************************************
-    // XXX: The Stdev calculation requires the mean.  Should we assume the
-    // mean has already been calculated? Or should we always calculate it?
     if (stats->options & PS_STAT_SAMPLE_STDEV) {
         p_psVectorSampleMean(inF32, mask, maskVal, stats);
@@ -1643,5 +1608,7 @@
 
     if ((stats->options & PS_STAT_CLIPPED_MEAN) || (stats->options & PS_STAT_CLIPPED_STDEV)) {
-        p_psVectorClippedStats(inF32, mask, maskVal, stats);
+        if (0 != p_psVectorClippedStats(inF32, mask, maskVal, stats)) {
+            psError(__func__, "p_psVectorClippedStats() failed.\n");
+        }
     }
     // ************************************************************************
Index: /trunk/psLib/src/dataManip/psStats.h
===================================================================
--- /trunk/psLib/src/dataManip/psStats.h	(revision 2268)
+++ /trunk/psLib/src/dataManip/psStats.h	(revision 2269)
@@ -10,6 +10,6 @@
  *  @author George Gusciora, MHPCC
  *
- *  @version $Revision: 1.32 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2004-10-27 00:57:31 $
+ *  @version $Revision: 1.33 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2004-11-03 03:30:30 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -174,8 +174,7 @@
 );
 
-
-void p_psNormalizeVector(psVector* myData);
-
-void p_psNormalizeVectorF64(psVector* myData);
+void psNormalizeVectorRange(psVector* myData,
+                            float low,
+                            float high);
 
 /// @}
Index: /trunk/psLib/src/math/psMinimize.c
===================================================================
--- /trunk/psLib/src/math/psMinimize.c	(revision 2268)
+++ /trunk/psLib/src/math/psMinimize.c	(revision 2269)
@@ -9,6 +9,6 @@
  *  @author GLG, MHPCC
  *
- *  @version $Revision: 1.83 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2004-11-02 19:08:33 $
+ *  @version $Revision: 1.84 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2004-11-03 03:30:30 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -61,10 +61,9 @@
 returned as a psVector sums.
  
-XXX: change name
 XXX: Use a static vector.
  *****************************************************************************/
-void p_psBuildSums1D(double x,
-                     psS32 polyOrder,
-                     psVector* sums)
+void psBuildSums1D(double x,
+                   psS32 polyOrder,
+                   psVector* sums)
 {
     psS32 i = 0;
@@ -74,5 +73,4 @@
         sums = psVectorAlloc(polyOrder, PS_TYPE_F64);
     }
-
     if (polyOrder > sums->n) {
         sums = psVectorRealloc(sums, polyOrder);
@@ -109,5 +107,4 @@
     psS32 i;
     psS32 k;
-    psBool mustFreeX = false;
     float sig;
     float p;
@@ -150,9 +147,5 @@
     }
 
-    if (mustFreeX == true) {
-        psFree(X);
-    }
     psFree(u);
-
     psTrace(".psLib.dataManip.CalculateSecondDerivs", 4,
             "---- CalculateSecondDerivs() end ----\n");
@@ -171,7 +164,7 @@
 polynomials which are stored in psSpline1D.
  
-XXX: check types/sizes
 XXX: This is F32 only
  *****************************************************************************/
+/*
 float p_psNRSpline1DEval(psSpline1D *spline,
                          const psVector* restrict x,
@@ -187,5 +180,5 @@
     PS_VECTOR_CHECK_NULL(y, NAN);
     PS_VECTOR_CHECK_TYPE(y, PS_TYPE_F32, NAN);
-
+ 
     psS32 n;
     psS32 klo;
@@ -197,7 +190,7 @@
     float D;
     float Y;
-
+ 
     n = spline->n;
-    klo = p_psVectorBinDisectF32(spline->domains, (spline->n)+1, X);
+    klo = p_psVectorBinDisect32(spline->domains, (spline->n)+1, X);
     khi = klo + 1;
     H = (spline->domains)[khi] - (spline->domains)[klo];
@@ -206,13 +199,13 @@
     C = ((A*A*A)-A) * (H*H/6.0);
     D = ((B*B*B)-B) * (H*H/6.0);
-
+ 
     Y = (A * y->data.F32[klo]) +
         (B * y->data.F32[khi]) +
         (C * (spline->p_psDeriv2)[klo]) +
         (D * (spline->p_psDeriv2)[khi]);
-
+ 
     return(Y);
 }
-
+*/
 /*****************************************************************************/
 /* FUNCTION IMPLEMENTATION - PUBLIC                                          */
@@ -239,4 +232,5 @@
  
 XXX: usage of yErr is not specified in IfA documentation.
+ 
 XXX: Is the x argument redundant?  What do we do if the x argument is
 supplied, but does not equal the domains specified in mySpline?
@@ -414,5 +408,4 @@
 /******************************************************************************
 XXX: We assume unnormalized gaussians.
-XXX: Currently, yErr is ignored.
  *****************************************************************************/
 psVector *psMinimizeLMChi2Gauss1D(psImage *deriv,
@@ -441,6 +434,4 @@
         x = ((psVector *) (coords->data[i]))->data.F32[0];
         out->data.F32[i] = psGaussian(x, mean, stdev, false);
-        //        psTrace(".psLib.dataManip.psMinimize", 6,
-        //            "out->data.F32[%d]: (x, m, s) is (%f, %f, %f) is %f\n", i, x, mean, stdev, out->data.F32[i]);
     }
 
@@ -449,5 +440,4 @@
         float tmp = (x - mean) * psGaussian(x, mean, stdev, false);
         deriv->data.F32[i][0] = tmp / (stdev * stdev);
-
         tmp = (x - mean) * (x - mean) *
               psGaussian(x, mean, stdev, 0);
@@ -482,5 +472,5 @@
 XXX: This must work for both F32 and F64.  F32 is currently implemented.
      Note: since the LUD routines are only implemented in F64, then we
->     will have to convert all F32 input vectors to F64 regardless.  So,
+     will have to convert all F32 input vectors to F64 regardless.  So,
      the F64 port might be.
  *****************************************************************************/
@@ -730,7 +720,7 @@
 
 /******************************************************************************
-p_psVectorFitPolynomial1DCheb():  This routine will fit a Chebyshev
-polynomial of degree myPoly to the data points (x, y) and return the
-coefficients of that polynomial.
+VectorFitPolynomial1DCheb():  This routine will fit a Chebyshev polynomial of
+degree myPoly to the data points (x, y) and return the coefficients of that
+polynomial.
  
 XXX: yErr is currently ignored.
@@ -738,5 +728,5 @@
 XXX: Use private name?
 *****************************************************************************/
-psPolynomial1D *p_psVectorFitPolynomial1DCheby(psPolynomial1D* myPoly,
+psPolynomial1D *VectorFitPolynomial1DCheby(psPolynomial1D* myPoly,
         const psVector* restrict x,
         const psVector* restrict y,
@@ -785,5 +775,5 @@
         psFree(fScalar);
 
-        psTrace(".psLib.dataManip.p_psVectorFitPolynomial1DCheby", 6,
+        psTrace(".psLib.dataManip.VectorFitPolynomial1DCheby", 6,
                 "(x, X, y, f(X)) is (%f, %f, %f, %f)\n",
                 x->data.F64[i], X, y->data.F64[i], f->data.F64[i]);
@@ -809,12 +799,12 @@
 
 /******************************************************************************
-p_psVectorFitPolynomial1DOrd():  This routine will fit an ordinary
-polynomial of degree myPoly to the data points (x, y) and return the
-coefficients of that polynomial.
+VectorFitPolynomial1DOrd():  This routine will fit an ordinary polynomial of
+degree myPoly to the data points (x, y) and return the coefficients of that
+polynomial.
  
 XXX: Use private name?
 XXX: Use recycled vectors.
  *****************************************************************************/
-psPolynomial1D* p_psVectorFitPolynomial1DOrd(psPolynomial1D* myPoly,
+psPolynomial1D* VectorFitPolynomial1DOrd(psPolynomial1D* myPoly,
         const psVector* restrict x,
         const psVector* restrict y,
@@ -833,7 +823,7 @@
     psVector* xSums = NULL;
 
-    psTrace(".psLib.dataManip.psVectorFitPolynomial1D", 4,
-            "---- psVectorFitPolynomial1D() begin ----\n");
-    // printf("psVectorFitPolynomial1D()\n");
+    psTrace(".psLib.dataManip.VectorFitPolynomial1DOrd", 4,
+            "---- VectorFitPolynomial1DOrd() begin ----\n");
+    // printf("VectorFitPolynomial1D()\n");
     // for (i=0;i<x->n;i++) {
     // printf("(x, y, yErr) is (%f, %f, %f)\n", x->data.F64[i], y->data.F64[i], yErr->data.F64[i]);
@@ -865,5 +855,5 @@
     if (yErr == NULL) {
         for (i = 0; i < X->n; i++) {
-            p_psBuildSums1D(X->data.F64[i], 2 * polyOrder, xSums);
+            psBuildSums1D(X->data.F64[i], 2 * polyOrder, xSums);
 
             for (k = 0; k < polyOrder; k++) {
@@ -879,5 +869,5 @@
     } else {
         for (i = 0; i < X->n; i++) {
-            p_psBuildSums1D(X->data.F64[i], 2 * polyOrder, xSums);
+            psBuildSums1D(X->data.F64[i], 2 * polyOrder, xSums);
 
             for (k = 0; k < polyOrder; k++) {
@@ -916,6 +906,6 @@
     psFree(xSums);
 
-    psTrace(".psLib.dataManip.psVectorFitPolynomial1D", 4,
-            "---- psVectorFitPolynomial1D() begin ----\n");
+    psTrace(".psLib.dataManip.VectorFitPolynomial1DOrd", 4,
+            "---- VectorFitPolynomial1DOrd() begin ----\n");
     return (myPoly);
 }
@@ -962,5 +952,5 @@
         PS_VECTOR_GEN_X_INDEX_STATIC_F64(x64Static, y->n);
         if (myPoly->type == PS_POLYNOMIAL_CHEB) {
-            p_psNormalizeVector(x64Static);
+            p_psNormalizeVectorRange(x64Static, -1.0, 1.0);
         }
         x64 = x64Static;
@@ -974,7 +964,7 @@
     // Call the appropriate vector fitting routine.
     if (myPoly->type == PS_POLYNOMIAL_CHEB) {
-        p_psVectorFitPolynomial1DCheby(myPoly, x64, y64, yErr64);
+        VectorFitPolynomial1DCheby(myPoly, x64, y64, yErr64);
     } else if (myPoly->type == PS_POLYNOMIAL_ORD) {
-        tmpPoly = p_psVectorFitPolynomial1DOrd(myPoly, x64, y64, yErr64);
+        tmpPoly = VectorFitPolynomial1DOrd(myPoly, x64, y64, yErr64);
     } else {
         psError(__func__, "unknown polynomial type.\n");
Index: /trunk/psLib/src/math/psPolynomial.h
===================================================================
--- /trunk/psLib/src/math/psPolynomial.h	(revision 2268)
+++ /trunk/psLib/src/math/psPolynomial.h	(revision 2269)
@@ -12,6 +12,6 @@
 *  @author George Gusciora, MHPCC
 *
-*  @version $Revision: 1.31 $ $Name: not supported by cvs2svn $
-*  @date $Date: 2004-10-27 00:57:31 $
+*  @version $Revision: 1.32 $ $Name: not supported by cvs2svn $
+*  @date $Date: 2004-11-03 03:30:30 $
 *
 *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -409,12 +409,4 @@
                                const psSpline1D *spline);
 
-psS32 p_psVectorBinDisectF32(float *bins,
-                             psS32 numBins,
-                             float x);
-
-psS32 p_psVectorBinDisectS32(psS32 *bins,
-                             psS32 numBins,
-                             psS32 x);
-
 psS32 p_psVectorBinDisect(psVector *bins,
                           psScalar *x);
Index: /trunk/psLib/src/math/psSpline.h
===================================================================
--- /trunk/psLib/src/math/psSpline.h	(revision 2268)
+++ /trunk/psLib/src/math/psSpline.h	(revision 2269)
@@ -12,6 +12,6 @@
 *  @author George Gusciora, MHPCC
 *
-*  @version $Revision: 1.31 $ $Name: not supported by cvs2svn $
-*  @date $Date: 2004-10-27 00:57:31 $
+*  @version $Revision: 1.32 $ $Name: not supported by cvs2svn $
+*  @date $Date: 2004-11-03 03:30:30 $
 *
 *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -409,12 +409,4 @@
                                const psSpline1D *spline);
 
-psS32 p_psVectorBinDisectF32(float *bins,
-                             psS32 numBins,
-                             float x);
-
-psS32 p_psVectorBinDisectS32(psS32 *bins,
-                             psS32 numBins,
-                             psS32 x);
-
 psS32 p_psVectorBinDisect(psVector *bins,
                           psScalar *x);
Index: /trunk/psLib/src/math/psStats.c
===================================================================
--- /trunk/psLib/src/math/psStats.c	(revision 2268)
+++ /trunk/psLib/src/math/psStats.c	(revision 2269)
@@ -9,6 +9,6 @@
  *  @author GLG, MHPCC
  *
- *  @version $Revision: 1.83 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2004-11-02 19:13:03 $
+ *  @version $Revision: 1.84 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2004-11-03 03:30:30 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -403,16 +403,18 @@
  *****************************************************************************/
 void p_psVectorSampleMedian(const psVector* restrict myVector,
-                            const psVector* restrict maskVector, psU32 maskVal, psStats* stats)
+                            const psVector* restrict maskVector,
+                            psU32 maskVal,
+                            psStats* stats)
 {
     psVector* unsortedVector = NULL;    // Temporary vector
     psVector* sortedVector = NULL;      // Temporary vector
-    psS32 i = 0;                  // Loop index variable
-    psS32 count = 0;              // # of points in this mean?
-    psS32 nValues = 0;            // # of points in vector
-    float rangeMin = 0.0;       // Exclude data below this
-    float rangeMax = 0.0;       // Exclude date above this
+    psS32 i = 0;                        // Loop index variable
+    psS32 count = 0;                    // # of points in this mean?
+    psS32 nValues = 0;                  // # of points in vector
+    float rangeMin = 0.0;               // Exclude data below this
+    float rangeMax = 0.0;               // Exclude date above this
 
     // Determine how many data points fit inside this min/max range
-    // and are not masked, IF the maskVector is not NULL>
+    // and are not masked, if the maskVector is not NULL>
     nValues = p_psVectorNValues(myVector, maskVector, maskVal, stats);
 
@@ -551,10 +553,12 @@
  *****************************************************************************/
 void p_psVectorSampleQuartiles(const psVector* restrict myVector,
-                               const psVector* restrict maskVector, psU32 maskVal, psStats* stats)
+                               const psVector* restrict maskVector,
+                               psU32 maskVal,
+                               psStats* stats)
 {
     psVector* unsortedVector = NULL;    // Temporary vector
     psVector* sortedVector = NULL;      // Temporary vector
     psS32 i = 0;                  // Loop index variable
-    psS32 count = 0;              // # of points in this mean?
+    psS32 count = 0;              // # of points in this mean.
     psS32 nValues = 0;            // # data points
     float rangeMin = 0.0;       // Exclude data below this
@@ -631,5 +635,7 @@
  *****************************************************************************/
 void p_psVectorSampleStdev(const psVector* restrict myVector,
-                           const psVector* restrict maskVector, psU32 maskVal, psStats* stats)
+                           const psVector* restrict maskVector,
+                           psU32 maskVal,
+                           psStats* stats)
 {
     psS32 i = 0;                  // Loop index variable
@@ -713,8 +719,10 @@
     stats
 Returns
-    NULL
- *****************************************************************************/
-void p_psVectorClippedStats(const psVector* restrict myVector,
-                            const psVector* restrict maskVector, psU32 maskVal, psStats* stats)
+    0 for success.
+ *****************************************************************************/
+int p_psVectorClippedStats(const psVector* restrict myVector,
+                           const psVector* restrict maskVector,
+                           psU32 maskVal,
+                           psStats* stats)
 {
     psS32 i = 0;                  // Loop index variable
@@ -726,11 +734,13 @@
     psVector* tmpMask = NULL;   // Temporary vector
 
-    // Endure that stats->clipIter is within the proper range.
+    // Ensure that stats->clipIter is within the proper range.
     if (!((PS_CLIPPED_NUM_ITER_LB <= stats->clipIter) && (stats->clipIter <= PS_CLIPPED_NUM_ITER_UB))) {
         psError(__func__, "Unallowed value for clipIter (%d).\n", stats->clipIter);
-    }
-    // Endure that stats->clipSigma is within the proper range.
+        return(-1);
+    }
+    // Ensure that stats->clipSigma is within the proper range.
     if (!((PS_CLIPPED_SIGMA_LB <= stats->clipSigma) && (stats->clipSigma <= PS_CLIPPED_SIGMA_UB))) {
         psError(__func__, "Unallowed value for clipSigma (%f).\n", stats->clipSigma);
+        return(-1);
     }
     // We allocate a temporary mask vector since during the iterative
@@ -798,7 +808,12 @@
 
     psFree(tmpMask);
-}
-
-void p_psNormalizeVectorF32_0(psVector* myData)
+    return(0);
+}
+
+/*****************************************************************************
+ *****************************************************************************/
+void p_psNormalizeVectorRangeF32(psVector* myData,
+                                 float outLow,
+                                 float outHigh)
 {
     float min = (float)HUGE;
@@ -815,22 +830,11 @@
     }
 
-    //  float range = max - min;
-    //    for (i = 0; i < myData->n; i++) {
-    //        myData->data.F32[i] = (myData->data.F32[i] - min) / range;
-    //    }
     for (i = 0; i < myData->n; i++) {
-        myData->data.F32[i] = (myData->data.F32[i] - min) / (max - min);
-    }
-}
-
-
-
-/*****************************************************************************
-p_psNormalizeVectorF32(myData): this is a private function which normalizes the
-elements of a vector to a range between 0.0 and 1.0.
- 
-XXX: how about an arbitrary p_psNormalizeVectorF32(myData, min, max)?
- *****************************************************************************/
-void p_psNormalizeVectorF32(psVector* myData)
+        myData->data.F32[i] = outLow + (myData->data.F32[i] - min) * (outHigh - outLow) / (max - min);
+    }
+}
+void p_psNormalizeVectorRangeF64(psVector* myData,
+                                 float outLow,
+                                 float outHigh)
 {
     float min = (float)HUGE;
@@ -840,63 +844,28 @@
     for (i = 0; i < myData->n; i++) {
         if (myData->data.F32[i] < min) {
-            min = myData->data.F32[i];
+            min = myData->data.F64[i];
         }
         if (myData->data.F32[i] > max) {
-            max = myData->data.F32[i];
-        }
-    }
-
-    //  float range = max - min;
-    //    for (i = 0; i < myData->n; i++) {
-    //        myData->data.F32[i] = (myData->data.F32[i] - min) / range;
-    //    }
+            max = myData->data.F64[i];
+        }
+    }
+
     for (i = 0; i < myData->n; i++) {
-        myData->data.F32[i] = (myData->data.F32[i] - 0.5 * (min + max)) /
-                              (0.5 * (max - min));
+        myData->data.F64[i] = outLow + (myData->data.F64[i] - min) * (outHigh - outLow) / (max - min);
     }
 }
 
 /*****************************************************************************
-p_psNormalizeVectorF64(myData): this is a private function which normalizes the
-elements of a vector to a range between -1.0 and 1.0.
-XXX: 0-1 or -1:1?
- 
-XXX: how about an arbitrary p_psNormalizeVectorF64(myData, min, max)?
- *****************************************************************************/
-void p_psNormalizeVectorF64(psVector* myData)
-{
-    float min = (double)HUGE;
-    float max = (double)-HUGE;
-    double range = 0.0;
-    psS32 i = 0;
-
-    for (i = 0; i < myData->n; i++) {
-        if (myData->data.F64[i] < min) {
-            min = myData->data.F64[i];
-        }
-        if (myData->data.F64[i] > max) {
-            max = myData->data.F64[i];
-        }
-    }
-
-    range = max - min;
-    for (i = 0; i < myData->n; i++) {
-        myData->data.F64[i] = -1.0 + 2.0 *
-                              ((myData->data.F64[i] - min) / range);
-    }
-
-    //    for (i = 0; i < myData->n; i++) {
-    //        myData->data.F64[i] = (myData->data.F64[i] - 0.5 * (min + max)) /
-    //                              (0.5 * (max - min));
-    //    }
-}
-
-// XXX: how about an arbitrary p_psNormalizeVector(myData, min, max)?
-void p_psNormalizeVector(psVector* myData)
+psNormalizeVectorRange(myData, low, high): this is a private function which
+normalizes the elements of a vector to a range between low and high.
+ *****************************************************************************/
+void psNormalizeVectorRange(psVector* myData,
+                            float low,
+                            float high)
 {
     if (myData->type.type == PS_TYPE_F32) {
-        p_psNormalizeVectorF32(myData);
+        p_psNormalizeVectorRangeF32(myData, low, high);
     } else if (myData->type.type == PS_TYPE_F64) {
-        p_psNormalizeVectorF64(myData);
+        p_psNormalizeVectorRangeF64(myData, low, high);
     } else {
         psError(__func__, "Unalowable data type.\n");
@@ -914,4 +883,6 @@
  
 XXX: Terminate when f(x)-getThisValue is within some error tolerance.
+ 
+XXX: Solve for X analytically.
  *****************************************************************************/
 float p_ps1DPolyMedian(psPolynomial1D* myPoly,
@@ -924,6 +895,4 @@
     float oldMidpoint = 1.0;
     float f = 0.0;
-
-    // printf("p_ps1DPolyMedian(%f, %f, %f) \n", rangeLow, rangeHigh, getThisValue);
 
     while (numIterations < PS_POLY_MEDIAN_MAX_ITERATIONS) {
@@ -956,4 +925,6 @@
 and i+1 is used for x[i]).  It then determines for what value x does that
 quadratic f(x) = yVal (the input parameter).
+ 
+XXX: After you fit the polynomial, solve for X analytically.
 *****************************************************************************/
 float fitQuadraticSearchForYThenReturnX(psVector *xVec,
@@ -1034,5 +1005,5 @@
     stats
 Returns
-    NULL
+    0 on success.
 *****************************************************************************/
 int p_psVectorRobustStats(const psVector* restrict myVector,
@@ -1119,8 +1090,4 @@
                             tmpStats->clippedStdev / 4.0f);
 
-    // The following was necessary to fit a gaussian to the data, since
-    // gaussian functions produce data between 0.0 and 1.0.
-    // p_psNormalizeVectorF32(robustHistogramVector);
-
     /**************************************************************************
     Determine the median/lower/upper quartile bin numbers.
@@ -1207,5 +1174,5 @@
     psVector *y = psVectorAlloc(robustHistogramVector->n, PS_TYPE_F32);
 
-    p_psNormalizeVectorF32_0(robustHistogramVector);
+    p_psNormalizeVectorRange(robustHistogramVector, 0.0, 1.0);
     for (i=0;i<robustHistogramVector->n;i++) {
         myCoords->data[i] = (psPtr *) psVectorAlloc(2, PS_TYPE_F32);
@@ -1620,6 +1587,4 @@
     }
     // ************************************************************************
-    // XXX: The Stdev calculation requires the mean.  Should we assume the
-    // mean has already been calculated? Or should we always calculate it?
     if (stats->options & PS_STAT_SAMPLE_STDEV) {
         p_psVectorSampleMean(inF32, mask, maskVal, stats);
@@ -1643,5 +1608,7 @@
 
     if ((stats->options & PS_STAT_CLIPPED_MEAN) || (stats->options & PS_STAT_CLIPPED_STDEV)) {
-        p_psVectorClippedStats(inF32, mask, maskVal, stats);
+        if (0 != p_psVectorClippedStats(inF32, mask, maskVal, stats)) {
+            psError(__func__, "p_psVectorClippedStats() failed.\n");
+        }
     }
     // ************************************************************************
Index: /trunk/psLib/src/math/psStats.h
===================================================================
--- /trunk/psLib/src/math/psStats.h	(revision 2268)
+++ /trunk/psLib/src/math/psStats.h	(revision 2269)
@@ -10,6 +10,6 @@
  *  @author George Gusciora, MHPCC
  *
- *  @version $Revision: 1.32 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2004-10-27 00:57:31 $
+ *  @version $Revision: 1.33 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2004-11-03 03:30:30 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -174,8 +174,7 @@
 );
 
-
-void p_psNormalizeVector(psVector* myData);
-
-void p_psNormalizeVectorF64(psVector* myData);
+void psNormalizeVectorRange(psVector* myData,
+                            float low,
+                            float high);
 
 /// @}