Changeset 1188

Timestamp:

Jul 6, 2004, 6:05:29 PM (22 years ago)

Author:

gusciora

Message:

...

Location:

trunk/psLib/src

Files:

• dataManip/psMinimize.c (modified) (26 diffs)
• math/psMinimize.c (modified) (26 diffs)

trunk/psLib/src/dataManip/psMinimize.c

@@ -25 +25 @@
 #include <math.h>
 
-#define MAX_LMM_NUM_ITERATIONS 500
+#define MAX_LMM_ITERATIONS 100
+#define MAX_MINIMIZE_ITERATIONS 100
 typedef struct
 {
@@ -38 +39 @@
     float (*d_evalModel) (const psVector *, const psVector *, int);
 }
-psModelData;
-// The first argument to evalModel() and
-// d_evalModel() specifies the data point.
-// It must have the same size as the second
-// dimension of *domain.
-// The second argument must have the same size
-// as *initialGuess and *paramMask.
-
+psMinChi2Data;
 
 typedef struct
@@ -56 +50 @@
     float (*d_evalModel) (const psVector *, const psVector *, int);
 }
-psModelData2;
+psMinimizeData;
 
 
 /******************************************************************************
-p_psMinFunc(*v, *params): This routine calls the user supplied function to
-be minimized.  The "v" argument of this function corresponds to the
-parameters of the function to be minimized.
- 
- 
+p_psMinFunc(*params, *funcData): We use the GSL-supplied function
+gsl_multimin_fdfminimizer_iterate() to minimize an arbitary function supplied
+by the user.  The GSL function requires the user-supplied function to be in
+a different format than the psLib format.  The purpose of this procedure is
+to serve as a GSL-format wrapper for the psLib user-supplied function which
+is to be minimized.
+    *params: The parameters of the function to be minimized.  These will be
+ varied by GSL in order to minimize the function.
+    *funcData: a psLib struct which contains the data point to be minimized,
+ the function and derivative function pointers, an initial guess at
+ the parameters, an option parameter mask, etc.
  *****************************************************************************/
-double p_psMinFunc(const gsl_vector *v,
-                   void *params)
+double p_psMinFunc(const gsl_vector *params,
+                   void *funcData)
 {
     int i;    // Loop index variable.
     int j;    // Loop index variable.
     float tmpf;    // Temporary floating point variable.
-    const psVector *restrict coord   = ((psModelData2 *) params)->coord;
-    const psVector *restrict mask     = ((psModelData2 *) params)->paramMask;
-    psVector *restrict initialGuess = ((psModelData2 *)params)->initialGuess;
+    const psVector *restrict coord   = ((psMinimizeData *) funcData)->coord;
+    const psVector *restrict mask     = ((psMinimizeData *) funcData)->paramMask;
+    psVector *restrict initialGuess = ((psMinimizeData *)funcData)->initialGuess;
     float (*evalModel)(const psVector *, const psVector *) =
-        ((psModelData2 *) params)->evalModel;
+        ((psMinimizeData *) funcData)->evalModel;
     psVector *inputParameterList = NULL;
 
+    // The GSL routines will call this function with the masked parameters
+    // removed.  However, the user-supplied function (to be modified) does not
+    // have those parameters removed.  Here will create a new parameter list
+    // with the masked parameters added (we expand initialGuess).
     inputParameterList = psVectorAlloc(initialGuess->n, PS_TYPE_F32);
-
     if (mask != NULL) {
         j = 0;
@@ -87 +90 @@
                 inputParameterList->data.F32[i] = initialGuess->data.F32[i];
             } else {
-                inputParameterList->data.F32[i] = gsl_vector_get(v, j++);
-            }
-        }
-    } else {
-        for (i=0;i<initialGuess->n;i++) {
-            inputParameterList->data.F32[i] = gsl_vector_get(v, i);
-        }
-    }
-    //printf("HMMM: (%.1f %.1f %.1f %.1f)\n",
-    //       inputParameterList->data.F32[0],
-    //       inputParameterList->data.F32[1],
-    //       coord->data.F32[0],
-    //       coord->data.F32[1]);
-
+                inputParameterList->data.F32[i] = gsl_vector_get(params, j++);
+            }
+        }
+    } else {
+        for (i=0;i<initialGuess->n;i++) {
+            inputParameterList->data.F32[i] = gsl_vector_get(params, i);
+        }
+    }
+
+    // Call the user-supplied function.
     tmpf = evalModel(inputParameterList, coord);
+
+    // Free allocated memory and return the value of the function.
     psFree(inputParameterList);
-    //printf("Called p_psMinFunc()\n");
     return(tmpf);
 }
 
-void p_psMinFuncDeriv(const gsl_vector *v,
-                      void *params,
+/******************************************************************************
+p_psMinFuncDeriv(*params, *funcData): We use the GSL-supplied function
+gsl_multimin_fdfminimizer_iterate() to minimize an arbitary function supplied
+by the user.  The GSL function requires the user-supplied function to be in
+a different format than the psLib format.  The purpose of this procedure is
+to serve as a GSL-format wrapper for the psLib user-supplied function which
+is to be minimized.
+    *params: The parameters of the function to be minimized.  These will be
+ varied by GSL in order to minimize the function.
+    *funcData: a psLib struct which contains the data point to be minimized,
+ the function and derivative function pointers, an initial guess at
+ the parameters, an option parameter mask, etc.
+    *df: we calculate the derivative of the function w.r.t. to each parameter
+ in "params" and return those derivatives in this psVector.
+ *****************************************************************************/
+void p_psMinFuncDeriv(const gsl_vector *params,
+                      void *funcData,
                       gsl_vector *df)
 {
@@ -114 +129 @@
     int j;    // Loop index variable.
     float tmpf;    // Temporary floating point variable.
-    const psVector *restrict coord   = ((psModelData2 *) params)->coord;
-    const psVector *restrict mask     = ((psModelData2 *) params)->paramMask;
-    psVector *restrict initialGuess = ((psModelData2 *)params)->initialGuess;
+    const psVector *restrict coord   = ((psMinimizeData *) funcData)->coord;
+    const psVector *restrict mask     = ((psMinimizeData *) funcData)->paramMask;
+    psVector *restrict initialGuess = ((psMinimizeData *)funcData)->initialGuess;
     float (*d_evalModel)(const psVector *, const psVector *, int) =
-        ((psModelData2 *) params)->d_evalModel;
+        ((psMinimizeData *) funcData)->d_evalModel;
     psVector *inputParameterList = NULL;
 
+    // The GSL routines will call this function with the masked parameters
+    // removed.  However, the user-supplied function (to be modified) does not
+    // have those parameters removed.  Here will create a new parameter list
+    // with the masked parameters added (we expand initialGuess).
     inputParameterList = psVectorAlloc(initialGuess->n, PS_TYPE_F32);
-
     if (mask != NULL) {
         j = 0;
@@ -129 +147 @@
                 inputParameterList->data.F32[i] = initialGuess->data.F32[i];
             } else {
-                inputParameterList->data.F32[i] = gsl_vector_get(v, j++);
-            }
-        }
-    } else {
-        for (i=0;i<initialGuess->n;i++) {
-            inputParameterList->data.F32[i] = gsl_vector_get(v, i);
-        }
-    }
-
+                inputParameterList->data.F32[i] = gsl_vector_get(params, j++);
+            }
+        }
+    } else {
+        for (i=0;i<initialGuess->n;i++) {
+            inputParameterList->data.F32[i] = gsl_vector_get(params, i);
+        }
+    }
+
+    // Evaluate the derivative w.r.t. each parameter.
+    // NOTE: we can probably remove the calls for masked parameters.
     for (i=0;i<initialGuess->n;i++) {
         tmpf = d_evalModel(inputParameterList, coord, i);
         gsl_vector_set(df, i, tmpf);
     }
+
+    // Free allocated memory.
     psFree(inputParameterList);
 }
 
-/* Compute both f and df together. */
-void p_psMinFuncFuncDeriv(const gsl_vector *x,
-                          void *params,
+/******************************************************************************
+    Compute both p_psMinFunc and p_psMinFuncDeriv together.
+ *****************************************************************************/
+void p_psMinFuncFuncDeriv(const gsl_vector *params,
+                          void *funcData,
                           double *f,
                           gsl_vector *df)
 {
-    *f = p_psMinFunc(x, params);
-    p_psMinFuncDeriv(x, params, df);
+    *f = p_psMinFunc(params, funcData);
+    p_psMinFuncDeriv(params, funcData, df);
 }
 
@@ -161 +185 @@
 This routine must minimize an arbitrary function; it must determine the set
 of parameters of that function such that the 
- 
  *****************************************************************************/
 psVector *
@@ -175 +198 @@
     int iter = 0;
     gsl_multimin_function_fdf f;
-    double *xInit = NULL;
     const gsl_multimin_fdfminimizer_type *T;
     gsl_multimin_fdfminimizer *s;
-    psModelData2 inputData;
+    psMinimizeData inputData;
     gsl_vector *x;
 
@@ -188 +210 @@
     inputData.paramCount = 0;
 
-    printf("psMinimize(): initialGuess->n is %d\n", initialGuess->n);
-    printf("psMinimize(): coord->n is %d\n", coord->n);
-    printf("psMinimize(): coord is (%.1f, %.1f)\n", coord->data.F32[0],
-           coord->data.F32[1]);
-
+    // If the user supplied a parameter mask, then count the number of
+    // non-masked elements.  This will be used later in allocating a vector
+    // for the parameters.
     if (paramMask != NULL) {
         for (i=0;i<paramMask->n;i++) {
@@ -204 +224 @@
 
     // The initial guess at the parameters for the function are written into
-    // the vector inputParameterList.  If the paramMask is not NULL, then those
-    // parameters are masked out.
-
-    xInit = (double *) psAlloc(inputData.paramCount * sizeof(double));
+    // the vector inputParameterList.  If the paramMask is not NULL, then
+    // masked parameters are masked out.
+    x = gsl_vector_alloc(inputData.paramCount);
     if (paramMask != NULL) {
         j = 0;
         for (i=0;i<initialGuess->n;i++) {
             if (paramMask->data.U8[i] == 0) {
-                xInit[j++] = initialGuess->data.F32[i];
-            }
-        }
-    } else {
-        for (i=0;i<initialGuess->n;i++) {
-            xInit[i] = initialGuess->data.F32[i];
+                gsl_vector_set(x, j++, initialGuess->data.F32[i]);
+            }
+        }
+    } else {
+        for (i=0;i<initialGuess->n;i++) {
+            gsl_vector_set(x, i, initialGuess->data.F32[i]);
         }
     }
@@ -224 +243 @@
     f.fdf = &p_psMinFuncFuncDeriv;
     f.n = inputData.paramCount;
-    // GUS: is this correct?
-    //f.params = xInit;
     f.params = &inputData;
 
-
-    printf("inputData.paramCount is %d\n", inputData.paramCount);
-    x = gsl_vector_alloc(inputData.paramCount);
-    for (i=0;i<inputData.paramCount;i++) {
-        gsl_vector_set(x, i, xInit[i]);
-    }
     T = gsl_multimin_fdfminimizer_conjugate_fr;
     s = gsl_multimin_fdfminimizer_alloc(T, inputData.paramCount);
-
     gsl_multimin_fdfminimizer_set(s, &f, x, 0.01, 1e-4);
-
     do {
         iter++;
@@ -251 +260 @@
             printf ("Minimum found at:\n");
 
-    } while (status == GSL_CONTINUE && iter < 100);
+    } while (status == GSL_CONTINUE && iter < MAX_MINIMIZE_ITERATIONS);
 
     // In the above steps we had removed the masked elements from the
@@ -270 +279 @@
         }
     }
-    psFree(xInit);
     return(initialGuess);
 }
@@ -278 +286 @@
 
 
-// The first argument to evalModel() and
-// d_evalModel() specifies the data point.
-// It must have the same size as the second
-// dimension of *domain.
-// The second argument must have the same size
-// as *initialGuess and *paramMask.
+// The first argument to evalModel() and d_evalModel() specifies the data
+// point.  It must have the same size as the second dimension of *domain.
+// The second argument must have the same size as *initialGuess and
+// *paramMask.
 
 /******************************************************************************
- ******************************************************************************
-gsl_function_f(x, params, f): This function serves as a standardized wrapper
-for the user supplied function which is to be minimized.  The GSL
-minimization routines have no knowledge of the user-supplied function.  They
-only call this routine, which then calls the user-supplied function.  The
-arguments are:
- 
+p_psMinChi2Func(*x, *funcData, *outdata): We use the GSL-supplied function
+gsl_multifit_fdfsolver_iterate() to determine the function parameters that
+best fit the supllied set of data points.  That GSL function requires the
+user-supplied function to be in a different format than the psLib format.
+The purpose of this procedure is to serve as a GSL-format wrapper for the
+psLib user-supplied function which is to be minimized.
     x: These are the parameters which are to be varied by GSL in order to
  minimized chi2 over the data set.
-    params: this data structure contains the input values over which the
+    funcData: this data structure contains the input values over which the
  function will be evaluated, the expected value of the function at
  those points, the amount of error tolerable at those points, a mask
@@ -302 +307 @@
     outData: The function is evaluated at each point, then subtract the
  expected value and divide by the error.
- ******************************************************************************
  *****************************************************************************/
-int gsl_function_f(const gsl_vector *x,
-                   void *params,
-                   gsl_vector *outData)
+int p_psMinChi2Func(const gsl_vector *params,
+                    void *funcData,
+                    gsl_vector *outData)
 {
     int i;    // Loop index variable.
     int j;    // Loop index variable.
     float tmpf;    // Temporary floating point variable.
-    const psImage *restrict  domain   = ((psModelData *)params)->domain;
-    const psVector *restrict data     = ((psModelData *)params)->data;
-    const psVector *restrict errors   = ((psModelData *) params)->errors;
-    const psVector *restrict mask     = ((psModelData *) params)->paramMask;
-    psVector *restrict initialGuess = ((psModelData *)params)->initialGuess;
-    float (*evalModel)(const psVector *, const psVector *) = ((psModelData *) params)->evalModel;
+    const psImage *restrict  domain   = ((psMinChi2Data *)funcData)->domain;
+    const psVector *restrict data     = ((psMinChi2Data *)funcData)->data;
+    const psVector *restrict errors   = ((psMinChi2Data *) funcData)->errors;
+    const psVector *restrict mask     = ((psMinChi2Data *) funcData)->paramMask;
+    psVector *restrict initialGuess = ((psMinChi2Data *)funcData)->initialGuess;
+    float (*evalModel)(const psVector *, const psVector *) = ((psMinChi2Data *) funcData)->evalModel;
     psVector *inputParameterList = NULL;
     psVector *tmpVecPtr = NULL;
@@ -326 +330 @@
     // have those parameters removed.  Here will create a new parameter list
     // with the masked parameters added (we expand initialGuess).
+
     inputParameterList = psVectorAlloc(initialGuess->n, PS_TYPE_F32);
-
     if (mask != NULL) {
         j = 0;
@@ -334 +338 @@
                 inputParameterList->data.F32[i] = initialGuess->data.F32[i];
             } else {
-                inputParameterList->data.F32[i] = gsl_vector_get(x, j++);
-            }
-        }
-    } else {
-        for (i=0;i<initialGuess->n;i++) {
-            inputParameterList->data.F32[i] = gsl_vector_get(x, i);
-        }
-    }
-
+                inputParameterList->data.F32[i] = gsl_vector_get(params, j++);
+            }
+        }
+    } else {
+        for (i=0;i<initialGuess->n;i++) {
+            inputParameterList->data.F32[i] = gsl_vector_get(params, i);
+        }
+    }
+
+    // Evaluate the function at each data point.
     for (i=0;i<domain->numRows;i++) {
         for (j=0;j<domain->numCols;j++) {
@@ -353 +358 @@
     }
 
+    // Free allocated memory.
     psFree(inputParameterList);
     psFree(tmpVecPtr);
+
     return GSL_SUCCESS;
 }
 
-int gsl_function_df(const gsl_vector *x,
-                    void *params,
-                    gsl_matrix *J)
-{
-    const psImage *restrict domain   = ((psModelData *)params)->domain;
-    const psVector *restrict errors   = ((psModelData *) params)->errors;
-    const psVector *restrict mask     = ((psModelData *) params)->paramMask;
-    psVector *restrict initialGuess = ((psModelData *)params)->initialGuess;
+int p_psMinChi2FuncDeriv(const gsl_vector *params,
+                         void *funcData,
+                         gsl_matrix *J)
+{
+    const psImage *restrict domain   = ((psMinChi2Data *)funcData)->domain;
+    const psVector *restrict errors   = ((psMinChi2Data *) funcData)->errors;
+    const psVector *restrict mask     = ((psMinChi2Data *) funcData)->paramMask;
+    psVector *restrict initialGuess = ((psMinChi2Data *)funcData)->initialGuess;
     psVector *inputParameterList = NULL;
     psVector *tmpVecPtr = NULL;
-    float (*d_evalModel)(const psVector *, const psVector *, int) = ((psModelData *) params)->d_evalModel;
+    float (*d_evalModel)(const psVector *, const psVector *, int) = ((psMinChi2Data *) funcData)->d_evalModel;
 
     size_t i;
@@ -380 +387 @@
     // have those parameters removed.  Here will create a new parameter list
     // with the masked parameters added (we expand initialGuess).
+
     inputParameterList = psVectorAlloc(initialGuess->n, PS_TYPE_F32);
-
     if (mask != NULL) {
         j = 0;
@@ -388 +395 @@
                 inputParameterList->data.F32[i] = initialGuess->data.F32[i];
             } else {
-                inputParameterList->data.F32[i] = gsl_vector_get(x, j++);
-            }
-        }
-    } else {
-        for (i=0;i<initialGuess->n;i++) {
-            inputParameterList->data.F32[i] = gsl_vector_get(x, i);
-        }
-    }
-
+                inputParameterList->data.F32[i] = gsl_vector_get(params, j++);
+            }
+        }
+    } else {
+        for (i=0;i<initialGuess->n;i++) {
+            inputParameterList->data.F32[i] = gsl_vector_get(params, i);
+        }
+    }
+
+    // Evaluate the derivtaive at each data point, and w.r.t. each parameter.
     for (i=0;i<domain->numRows;i++) {
         for (j=0;j<tmpVecPtr->n;j++) {
@@ -413 +421 @@
 }
 
-int gsl_my_function_fdf(const gsl_vector *x,
-                        void *params,
-                        gsl_vector *f,
-                        gsl_matrix *J)
-{
-    gsl_function_f(x, params, f);
-    gsl_function_df(x, params, J);
+
+int p_psMinChi2FuncFuncDeriv(const gsl_vector *params,
+                             void *funcData,
+                             gsl_vector *f,
+                             gsl_matrix *J)
+{
+    p_psMinChi2Func(params, funcData, f);
+    p_psMinChi2FuncDeriv(params, funcData, J);
 
     return GSL_SUCCESS;
@@ -451 +460 @@
     // minimization.
     gsl_multifit_fdfsolver *s; // GSL data structure.
-    psModelData inputData;
+    psMinChi2Data inputData;
     float chiSqOld = 0.0;
 
@@ -464 +473 @@
     inputData.d_evalModel = DevalModel;
 
+    // If the user supplied a parameter mask, then count the number of
+    // non-masked elements.  This will be used later in allocating a vector
+    // for the parameters.
     if (paramMask != NULL) {
         for (i=0;i<paramMask->n;i++) {
@@ -477 +489 @@
     // the vector inputParameterList.  If the paramMask is not NULL, then those
     // parameters are masked out.
-
     xInit = (double *) psAlloc(inputData.paramCount * sizeof(double));
     if (paramMask != NULL) {
@@ -504 +515 @@
     // and the data structures those functions use.
 
-    f.f = &gsl_function_f;
-    f.df = &gsl_function_df;
-    f.fdf = &gsl_my_function_fdf;
+    f.f = &p_psMinChi2Func;
+    f.df = &p_psMinChi2FuncDeriv;
+    f.fdf = &p_psMinChi2FuncFuncDeriv;
     f.n = numData;
     f.p = inputData.paramCount;
@@ -542 +553 @@
         chiSqOld = *chiSq;
 
-    } while (status == GSL_CONTINUE && iter < MAX_LMM_NUM_ITERATIONS);
+    } while (status == GSL_CONTINUE && iter < MAX_LMM_ITERATIONS);
 
 

trunk/psLib/src/math/psMinimize.c

@@ -25 +25 @@
 #include <math.h>
 
-#define MAX_LMM_NUM_ITERATIONS 500
+#define MAX_LMM_ITERATIONS 100
+#define MAX_MINIMIZE_ITERATIONS 100
 typedef struct
 {
@@ -38 +39 @@
     float (*d_evalModel) (const psVector *, const psVector *, int);
 }
-psModelData;
-// The first argument to evalModel() and
-// d_evalModel() specifies the data point.
-// It must have the same size as the second
-// dimension of *domain.
-// The second argument must have the same size
-// as *initialGuess and *paramMask.
-
+psMinChi2Data;
 
 typedef struct
@@ -56 +50 @@
     float (*d_evalModel) (const psVector *, const psVector *, int);
 }
-psModelData2;
+psMinimizeData;
 
 
 /******************************************************************************
-p_psMinFunc(*v, *params): This routine calls the user supplied function to
-be minimized.  The "v" argument of this function corresponds to the
-parameters of the function to be minimized.
- 
- 
+p_psMinFunc(*params, *funcData): We use the GSL-supplied function
+gsl_multimin_fdfminimizer_iterate() to minimize an arbitary function supplied
+by the user.  The GSL function requires the user-supplied function to be in
+a different format than the psLib format.  The purpose of this procedure is
+to serve as a GSL-format wrapper for the psLib user-supplied function which
+is to be minimized.
+    *params: The parameters of the function to be minimized.  These will be
+ varied by GSL in order to minimize the function.
+    *funcData: a psLib struct which contains the data point to be minimized,
+ the function and derivative function pointers, an initial guess at
+ the parameters, an option parameter mask, etc.
  *****************************************************************************/
-double p_psMinFunc(const gsl_vector *v,
-                   void *params)
+double p_psMinFunc(const gsl_vector *params,
+                   void *funcData)
 {
     int i;    // Loop index variable.
     int j;    // Loop index variable.
     float tmpf;    // Temporary floating point variable.
-    const psVector *restrict coord   = ((psModelData2 *) params)->coord;
-    const psVector *restrict mask     = ((psModelData2 *) params)->paramMask;
-    psVector *restrict initialGuess = ((psModelData2 *)params)->initialGuess;
+    const psVector *restrict coord   = ((psMinimizeData *) funcData)->coord;
+    const psVector *restrict mask     = ((psMinimizeData *) funcData)->paramMask;
+    psVector *restrict initialGuess = ((psMinimizeData *)funcData)->initialGuess;
     float (*evalModel)(const psVector *, const psVector *) =
-        ((psModelData2 *) params)->evalModel;
+        ((psMinimizeData *) funcData)->evalModel;
     psVector *inputParameterList = NULL;
 
+    // The GSL routines will call this function with the masked parameters
+    // removed.  However, the user-supplied function (to be modified) does not
+    // have those parameters removed.  Here will create a new parameter list
+    // with the masked parameters added (we expand initialGuess).
     inputParameterList = psVectorAlloc(initialGuess->n, PS_TYPE_F32);
-
     if (mask != NULL) {
         j = 0;
@@ -87 +90 @@
                 inputParameterList->data.F32[i] = initialGuess->data.F32[i];
             } else {
-                inputParameterList->data.F32[i] = gsl_vector_get(v, j++);
-            }
-        }
-    } else {
-        for (i=0;i<initialGuess->n;i++) {
-            inputParameterList->data.F32[i] = gsl_vector_get(v, i);
-        }
-    }
-    //printf("HMMM: (%.1f %.1f %.1f %.1f)\n",
-    //       inputParameterList->data.F32[0],
-    //       inputParameterList->data.F32[1],
-    //       coord->data.F32[0],
-    //       coord->data.F32[1]);
-
+                inputParameterList->data.F32[i] = gsl_vector_get(params, j++);
+            }
+        }
+    } else {
+        for (i=0;i<initialGuess->n;i++) {
+            inputParameterList->data.F32[i] = gsl_vector_get(params, i);
+        }
+    }
+
+    // Call the user-supplied function.
     tmpf = evalModel(inputParameterList, coord);
+
+    // Free allocated memory and return the value of the function.
     psFree(inputParameterList);
-    //printf("Called p_psMinFunc()\n");
     return(tmpf);
 }
 
-void p_psMinFuncDeriv(const gsl_vector *v,
-                      void *params,
+/******************************************************************************
+p_psMinFuncDeriv(*params, *funcData): We use the GSL-supplied function
+gsl_multimin_fdfminimizer_iterate() to minimize an arbitary function supplied
+by the user.  The GSL function requires the user-supplied function to be in
+a different format than the psLib format.  The purpose of this procedure is
+to serve as a GSL-format wrapper for the psLib user-supplied function which
+is to be minimized.
+    *params: The parameters of the function to be minimized.  These will be
+ varied by GSL in order to minimize the function.
+    *funcData: a psLib struct which contains the data point to be minimized,
+ the function and derivative function pointers, an initial guess at
+ the parameters, an option parameter mask, etc.
+    *df: we calculate the derivative of the function w.r.t. to each parameter
+ in "params" and return those derivatives in this psVector.
+ *****************************************************************************/
+void p_psMinFuncDeriv(const gsl_vector *params,
+                      void *funcData,
                       gsl_vector *df)
 {
@@ -114 +129 @@
     int j;    // Loop index variable.
     float tmpf;    // Temporary floating point variable.
-    const psVector *restrict coord   = ((psModelData2 *) params)->coord;
-    const psVector *restrict mask     = ((psModelData2 *) params)->paramMask;
-    psVector *restrict initialGuess = ((psModelData2 *)params)->initialGuess;
+    const psVector *restrict coord   = ((psMinimizeData *) funcData)->coord;
+    const psVector *restrict mask     = ((psMinimizeData *) funcData)->paramMask;
+    psVector *restrict initialGuess = ((psMinimizeData *)funcData)->initialGuess;
     float (*d_evalModel)(const psVector *, const psVector *, int) =
-        ((psModelData2 *) params)->d_evalModel;
+        ((psMinimizeData *) funcData)->d_evalModel;
     psVector *inputParameterList = NULL;
 
+    // The GSL routines will call this function with the masked parameters
+    // removed.  However, the user-supplied function (to be modified) does not
+    // have those parameters removed.  Here will create a new parameter list
+    // with the masked parameters added (we expand initialGuess).
     inputParameterList = psVectorAlloc(initialGuess->n, PS_TYPE_F32);
-
     if (mask != NULL) {
         j = 0;
@@ -129 +147 @@
                 inputParameterList->data.F32[i] = initialGuess->data.F32[i];
             } else {
-                inputParameterList->data.F32[i] = gsl_vector_get(v, j++);
-            }
-        }
-    } else {
-        for (i=0;i<initialGuess->n;i++) {
-            inputParameterList->data.F32[i] = gsl_vector_get(v, i);
-        }
-    }
-
+                inputParameterList->data.F32[i] = gsl_vector_get(params, j++);
+            }
+        }
+    } else {
+        for (i=0;i<initialGuess->n;i++) {
+            inputParameterList->data.F32[i] = gsl_vector_get(params, i);
+        }
+    }
+
+    // Evaluate the derivative w.r.t. each parameter.
+    // NOTE: we can probably remove the calls for masked parameters.
     for (i=0;i<initialGuess->n;i++) {
         tmpf = d_evalModel(inputParameterList, coord, i);
         gsl_vector_set(df, i, tmpf);
     }
+
+    // Free allocated memory.
     psFree(inputParameterList);
 }
 
-/* Compute both f and df together. */
-void p_psMinFuncFuncDeriv(const gsl_vector *x,
-                          void *params,
+/******************************************************************************
+    Compute both p_psMinFunc and p_psMinFuncDeriv together.
+ *****************************************************************************/
+void p_psMinFuncFuncDeriv(const gsl_vector *params,
+                          void *funcData,
                           double *f,
                           gsl_vector *df)
 {
-    *f = p_psMinFunc(x, params);
-    p_psMinFuncDeriv(x, params, df);
+    *f = p_psMinFunc(params, funcData);
+    p_psMinFuncDeriv(params, funcData, df);
 }
 
@@ -161 +185 @@
 This routine must minimize an arbitrary function; it must determine the set
 of parameters of that function such that the 
- 
  *****************************************************************************/
 psVector *
@@ -175 +198 @@
     int iter = 0;
     gsl_multimin_function_fdf f;
-    double *xInit = NULL;
     const gsl_multimin_fdfminimizer_type *T;
     gsl_multimin_fdfminimizer *s;
-    psModelData2 inputData;
+    psMinimizeData inputData;
     gsl_vector *x;
 
@@ -188 +210 @@
     inputData.paramCount = 0;
 
-    printf("psMinimize(): initialGuess->n is %d\n", initialGuess->n);
-    printf("psMinimize(): coord->n is %d\n", coord->n);
-    printf("psMinimize(): coord is (%.1f, %.1f)\n", coord->data.F32[0],
-           coord->data.F32[1]);
-
+    // If the user supplied a parameter mask, then count the number of
+    // non-masked elements.  This will be used later in allocating a vector
+    // for the parameters.
     if (paramMask != NULL) {
         for (i=0;i<paramMask->n;i++) {
@@ -204 +224 @@
 
     // The initial guess at the parameters for the function are written into
-    // the vector inputParameterList.  If the paramMask is not NULL, then those
-    // parameters are masked out.
-
-    xInit = (double *) psAlloc(inputData.paramCount * sizeof(double));
+    // the vector inputParameterList.  If the paramMask is not NULL, then
+    // masked parameters are masked out.
+    x = gsl_vector_alloc(inputData.paramCount);
     if (paramMask != NULL) {
         j = 0;
         for (i=0;i<initialGuess->n;i++) {
             if (paramMask->data.U8[i] == 0) {
-                xInit[j++] = initialGuess->data.F32[i];
-            }
-        }
-    } else {
-        for (i=0;i<initialGuess->n;i++) {
-            xInit[i] = initialGuess->data.F32[i];
+                gsl_vector_set(x, j++, initialGuess->data.F32[i]);
+            }
+        }
+    } else {
+        for (i=0;i<initialGuess->n;i++) {
+            gsl_vector_set(x, i, initialGuess->data.F32[i]);
         }
     }
@@ -224 +243 @@
     f.fdf = &p_psMinFuncFuncDeriv;
     f.n = inputData.paramCount;
-    // GUS: is this correct?
-    //f.params = xInit;
     f.params = &inputData;
 
-
-    printf("inputData.paramCount is %d\n", inputData.paramCount);
-    x = gsl_vector_alloc(inputData.paramCount);
-    for (i=0;i<inputData.paramCount;i++) {
-        gsl_vector_set(x, i, xInit[i]);
-    }
     T = gsl_multimin_fdfminimizer_conjugate_fr;
     s = gsl_multimin_fdfminimizer_alloc(T, inputData.paramCount);
-
     gsl_multimin_fdfminimizer_set(s, &f, x, 0.01, 1e-4);
-
     do {
         iter++;
@@ -251 +260 @@
             printf ("Minimum found at:\n");
 
-    } while (status == GSL_CONTINUE && iter < 100);
+    } while (status == GSL_CONTINUE && iter < MAX_MINIMIZE_ITERATIONS);
 
     // In the above steps we had removed the masked elements from the
@@ -270 +279 @@
         }
     }
-    psFree(xInit);
     return(initialGuess);
 }
@@ -278 +286 @@
 
 
-// The first argument to evalModel() and
-// d_evalModel() specifies the data point.
-// It must have the same size as the second
-// dimension of *domain.
-// The second argument must have the same size
-// as *initialGuess and *paramMask.
+// The first argument to evalModel() and d_evalModel() specifies the data
+// point.  It must have the same size as the second dimension of *domain.
+// The second argument must have the same size as *initialGuess and
+// *paramMask.
 
 /******************************************************************************
- ******************************************************************************
-gsl_function_f(x, params, f): This function serves as a standardized wrapper
-for the user supplied function which is to be minimized.  The GSL
-minimization routines have no knowledge of the user-supplied function.  They
-only call this routine, which then calls the user-supplied function.  The
-arguments are:
- 
+p_psMinChi2Func(*x, *funcData, *outdata): We use the GSL-supplied function
+gsl_multifit_fdfsolver_iterate() to determine the function parameters that
+best fit the supllied set of data points.  That GSL function requires the
+user-supplied function to be in a different format than the psLib format.
+The purpose of this procedure is to serve as a GSL-format wrapper for the
+psLib user-supplied function which is to be minimized.
     x: These are the parameters which are to be varied by GSL in order to
  minimized chi2 over the data set.
-    params: this data structure contains the input values over which the
+    funcData: this data structure contains the input values over which the
  function will be evaluated, the expected value of the function at
  those points, the amount of error tolerable at those points, a mask
@@ -302 +307 @@
     outData: The function is evaluated at each point, then subtract the
  expected value and divide by the error.
- ******************************************************************************
  *****************************************************************************/
-int gsl_function_f(const gsl_vector *x,
-                   void *params,
-                   gsl_vector *outData)
+int p_psMinChi2Func(const gsl_vector *params,
+                    void *funcData,
+                    gsl_vector *outData)
 {
     int i;    // Loop index variable.
     int j;    // Loop index variable.
     float tmpf;    // Temporary floating point variable.
-    const psImage *restrict  domain   = ((psModelData *)params)->domain;
-    const psVector *restrict data     = ((psModelData *)params)->data;
-    const psVector *restrict errors   = ((psModelData *) params)->errors;
-    const psVector *restrict mask     = ((psModelData *) params)->paramMask;
-    psVector *restrict initialGuess = ((psModelData *)params)->initialGuess;
-    float (*evalModel)(const psVector *, const psVector *) = ((psModelData *) params)->evalModel;
+    const psImage *restrict  domain   = ((psMinChi2Data *)funcData)->domain;
+    const psVector *restrict data     = ((psMinChi2Data *)funcData)->data;
+    const psVector *restrict errors   = ((psMinChi2Data *) funcData)->errors;
+    const psVector *restrict mask     = ((psMinChi2Data *) funcData)->paramMask;
+    psVector *restrict initialGuess = ((psMinChi2Data *)funcData)->initialGuess;
+    float (*evalModel)(const psVector *, const psVector *) = ((psMinChi2Data *) funcData)->evalModel;
     psVector *inputParameterList = NULL;
     psVector *tmpVecPtr = NULL;
@@ -326 +330 @@
     // have those parameters removed.  Here will create a new parameter list
     // with the masked parameters added (we expand initialGuess).
+
     inputParameterList = psVectorAlloc(initialGuess->n, PS_TYPE_F32);
-
     if (mask != NULL) {
         j = 0;
@@ -334 +338 @@
                 inputParameterList->data.F32[i] = initialGuess->data.F32[i];
             } else {
-                inputParameterList->data.F32[i] = gsl_vector_get(x, j++);
-            }
-        }
-    } else {
-        for (i=0;i<initialGuess->n;i++) {
-            inputParameterList->data.F32[i] = gsl_vector_get(x, i);
-        }
-    }
-
+                inputParameterList->data.F32[i] = gsl_vector_get(params, j++);
+            }
+        }
+    } else {
+        for (i=0;i<initialGuess->n;i++) {
+            inputParameterList->data.F32[i] = gsl_vector_get(params, i);
+        }
+    }
+
+    // Evaluate the function at each data point.
     for (i=0;i<domain->numRows;i++) {
         for (j=0;j<domain->numCols;j++) {
@@ -353 +358 @@
     }
 
+    // Free allocated memory.
     psFree(inputParameterList);
     psFree(tmpVecPtr);
+
     return GSL_SUCCESS;
 }
 
-int gsl_function_df(const gsl_vector *x,
-                    void *params,
-                    gsl_matrix *J)
-{
-    const psImage *restrict domain   = ((psModelData *)params)->domain;
-    const psVector *restrict errors   = ((psModelData *) params)->errors;
-    const psVector *restrict mask     = ((psModelData *) params)->paramMask;
-    psVector *restrict initialGuess = ((psModelData *)params)->initialGuess;
+int p_psMinChi2FuncDeriv(const gsl_vector *params,
+                         void *funcData,
+                         gsl_matrix *J)
+{
+    const psImage *restrict domain   = ((psMinChi2Data *)funcData)->domain;
+    const psVector *restrict errors   = ((psMinChi2Data *) funcData)->errors;
+    const psVector *restrict mask     = ((psMinChi2Data *) funcData)->paramMask;
+    psVector *restrict initialGuess = ((psMinChi2Data *)funcData)->initialGuess;
     psVector *inputParameterList = NULL;
     psVector *tmpVecPtr = NULL;
-    float (*d_evalModel)(const psVector *, const psVector *, int) = ((psModelData *) params)->d_evalModel;
+    float (*d_evalModel)(const psVector *, const psVector *, int) = ((psMinChi2Data *) funcData)->d_evalModel;
 
     size_t i;
@@ -380 +387 @@
     // have those parameters removed.  Here will create a new parameter list
     // with the masked parameters added (we expand initialGuess).
+
     inputParameterList = psVectorAlloc(initialGuess->n, PS_TYPE_F32);
-
     if (mask != NULL) {
         j = 0;
@@ -388 +395 @@
                 inputParameterList->data.F32[i] = initialGuess->data.F32[i];
             } else {
-                inputParameterList->data.F32[i] = gsl_vector_get(x, j++);
-            }
-        }
-    } else {
-        for (i=0;i<initialGuess->n;i++) {
-            inputParameterList->data.F32[i] = gsl_vector_get(x, i);
-        }
-    }
-
+                inputParameterList->data.F32[i] = gsl_vector_get(params, j++);
+            }
+        }
+    } else {
+        for (i=0;i<initialGuess->n;i++) {
+            inputParameterList->data.F32[i] = gsl_vector_get(params, i);
+        }
+    }
+
+    // Evaluate the derivtaive at each data point, and w.r.t. each parameter.
     for (i=0;i<domain->numRows;i++) {
         for (j=0;j<tmpVecPtr->n;j++) {
@@ -413 +421 @@
 }
 
-int gsl_my_function_fdf(const gsl_vector *x,
-                        void *params,
-                        gsl_vector *f,
-                        gsl_matrix *J)
-{
-    gsl_function_f(x, params, f);
-    gsl_function_df(x, params, J);
+
+int p_psMinChi2FuncFuncDeriv(const gsl_vector *params,
+                             void *funcData,
+                             gsl_vector *f,
+                             gsl_matrix *J)
+{
+    p_psMinChi2Func(params, funcData, f);
+    p_psMinChi2FuncDeriv(params, funcData, J);
 
     return GSL_SUCCESS;
@@ -451 +460 @@
     // minimization.
     gsl_multifit_fdfsolver *s; // GSL data structure.
-    psModelData inputData;
+    psMinChi2Data inputData;
     float chiSqOld = 0.0;
 
@@ -464 +473 @@
     inputData.d_evalModel = DevalModel;
 
+    // If the user supplied a parameter mask, then count the number of
+    // non-masked elements.  This will be used later in allocating a vector
+    // for the parameters.
     if (paramMask != NULL) {
         for (i=0;i<paramMask->n;i++) {
@@ -477 +489 @@
     // the vector inputParameterList.  If the paramMask is not NULL, then those
     // parameters are masked out.
-
     xInit = (double *) psAlloc(inputData.paramCount * sizeof(double));
     if (paramMask != NULL) {
@@ -504 +515 @@
     // and the data structures those functions use.
 
-    f.f = &gsl_function_f;
-    f.df = &gsl_function_df;
-    f.fdf = &gsl_my_function_fdf;
+    f.f = &p_psMinChi2Func;
+    f.df = &p_psMinChi2FuncDeriv;
+    f.fdf = &p_psMinChi2FuncFuncDeriv;
     f.n = numData;
     f.p = inputData.paramCount;
@@ -542 +553 @@
         chiSqOld = *chiSq;
 
-    } while (status == GSL_CONTINUE && iter < MAX_LMM_NUM_ITERATIONS);
+    } while (status == GSL_CONTINUE && iter < MAX_LMM_ITERATIONS);