Index: /trunk/psLib/src/dataManip/psMinimize.c =================================================================== --- /trunk/psLib/src/dataManip/psMinimize.c (revision 2201) +++ /trunk/psLib/src/dataManip/psMinimize.c (revision 2202) @@ -9,6 +9,6 @@ * @author GLF, MHPCC * - * @version $Revision: 1.59 $ $Name: not supported by cvs2svn $ - * @date $Date: 2004-10-26 21:24:42 $ + * @version $Revision: 1.60 $ $Name: not supported by cvs2svn $ + * @date $Date: 2004-10-26 23:14:04 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -393,6 +393,6 @@ /****************************************************************************** -psMinimizeLMChi2(): This routine will take an procedure which calculates an -arbitrary function and it's derivative and minimize the chi-squared match +psMinimizeLMChi2(): This routine will take an procedure which calculates +an arbitrary function and it's derivative and minimize the chi-squared match between that function at the specified coords and the specified value at those coords. @@ -417,266 +417,11 @@ *****************************************************************************/ bool psMinimizeLMChi2(psMinimization *min, + psImage *covar, psVector *params, const psVector *paramMask, - psImage *covar, - const psArray *coords, - const psVector *value, + const psArray *x, + const psVector *y, + const psVector *yErr, psMinimizeLMChi2Func func) -{ - PS_CHECK_NULL_PTR_RETURN_NULL(min); - PS_CHECK_NULL_VECTOR_RETURN_NULL(params); - PS_CHECK_EMPTY_VECTOR_RETURN_NULL(params); - PS_CHECK_NULL_PTR_RETURN_NULL(coords); - PS_CHECK_NULL_VECTOR_RETURN_NULL(value); - PS_CHECK_EMPTY_VECTOR_RETURN_NULL(value); - if (paramMask != NULL) { - PS_CHECK_VECTOR_SIZE_EQUAL_RETURN_NULL(params, paramMask); - } - psTrace(".psLib.dataManip.psMinimizeLMChi2", 4, - "---- psMinimizeLMChi2() begin ----\n"); - int numData = value->n; - int numParams = params->n; - int i; - int j; - int k; - int l; - int n; - int p; - psVector *beta = psVectorAlloc(numParams, PS_TYPE_F64); - psVector *perm = psVectorAlloc(numParams, PS_TYPE_F64); - - psVector *paramDeltasF64 = psVectorAlloc(numParams, PS_TYPE_F64); - psVector *origParams = psVectorAlloc(numParams, PS_TYPE_F32); - psVector *newParams = psVectorAlloc(numParams, PS_TYPE_F32); - - psImage *alpha = psImageAlloc(numParams, numParams, PS_TYPE_F32); - psImage *A = psImageAlloc(numParams, numParams, PS_TYPE_F64); - psImage *aOut = psImageAlloc(numParams, numParams, PS_TYPE_F64); - - psVector **deriv = (psVector **) psAlloc(numData * sizeof(psVector *)); - for (i=0;in, PS_TYPE_F32); - psVector *newValueVec = psVectorAlloc(value->n, PS_TYPE_F32); - - float currChi2 = 0.0; - float newChi2 = 0.0; - float lamda = 0.00005; - - psTrace(".psLib.dataManip.psMinimizeLMChi2", 6, - "min->maxIter is %d\n", min->maxIter); - psTrace(".psLib.dataManip.psMinimizeLMChi2", 6, - "min->tol is %f\n", min->tol); - - for (p=0;pdata.F32[p] = params->data.F32[p]; - } - - min->lastDelta = HUGE; - min->iter = 0; - - while ((min->lastDelta > min->tol) && (min->iter < min->maxIter)) { - psTrace(".psLib.dataManip.psMinimizeLMChi2", 4, - "------------------------------------------------------\n"); - psTrace(".psLib.dataManip.psMinimizeLMChi2", 4, - "Iteration %d. Delta is %f\n", min->iter, min->lastDelta); - - // - // Calculate the current values and chi-squared of the function. - // - currChi2 = 0.0; - for (n=0;ndata.F32[n] = func(deriv[n], - params, - (psVector *) coords->data[n]); - currChi2+= (currValueVec->data.F32[n] * currValueVec->data.F32[n]); - } - - for (p=0;pdata.F32[%d] is %f.\n", p, params->data.F32[p]); - } - psTrace(".psLib.dataManip.psMinimizeLMChi2", 6, - "Current chi-squared is (%f)\n", currChi2); - - // - // Mask elements of the derivative for each data point. - // - for (p=0;pdata.U8[p] != 0)) { - for (n=0;ndata.F32[p] = 0.0; - } - } - } - - // - // Calculate the BETA vector. - // - for (p=0;pdata.F64[p] = 0.0; - for (n=0;ndata.F64[p])+= - (value->data.F32[n] - currValueVec->data.F32[n]) * - (deriv[n])->data.F32[p]; - } - // XXX: multiple by -1 here? - (beta->data.F64[p])*= -1.0; - psTrace(".psLib.dataManip.psMinimizeLMChi2", 6, - "beta->data.F64[%d] is %f.\n", p, beta->data.F64[p]); - } - - // - // Calculate the ALPHA matrix. - // - for (k=0;kdata.F32[k][l] = 0.0; - for (n=0;ndata.F32[k][l]+= (deriv[n])->data.F32[k] * - (deriv[n])->data.F32[l]; - } - } - } - - // - // Calculate the matrix A. - // - for (j=0;jdata.F64[j][k] = - (double) ((1.0 + lamda) * alpha->data.F32[j][k]); - } else { - A->data.F64[j][k] = (double) alpha->data.F32[j][k]; - } - } - } - for (j=0;jdata.F64[j][k]); - } - psTrace(".psLib.dataManip.psMinimizeLMChi2", 6, "Matrix A[][]:\n"); - } - - // - // Solve A * alpha = Beta - // - aOut = psMatrixLUD(aOut, perm, A); - paramDeltasF64 = psMatrixLUSolve(paramDeltasF64, aOut, beta, perm); - - // - // Mask any masked parameters. - // - for (i=0;idata.F64[%d] is %f.\n", i, paramDeltasF64->data.F64[i]); - if ((paramMask != NULL) && (paramMask->data.U8[i] != 0)) { - newParams->data.F32[i] = origParams->data.F32[i]; - } else { - newParams->data.F32[i] = params->data.F32[i] - - (float) paramDeltasF64->data.F64[i]; - } - } - - psTrace(".psLib.dataManip.psMinimizeLMChi2", 6, - "Calling func() with new parameters:\n"); - for (i=0;idata.F32[%d] is %f.\n", i, newParams->data.F32[i]); - } - - - // - // Calculate new function values. - // - newChi2 = 0.0; - for (n=0;ndata.F32[n] = func(deriv[n], newParams, - (psVector *) coords->data[n]); - newChi2+= (newValueVec->data.F32[n] * newValueVec->data.F32[n]); - } - psTrace(".psLib.dataManip.psMinimizeLMChi2", 4, - "old/new chi-squareds are (%f, %f)\n", currChi2, newChi2); - - // - // If the new chi-squared is lower, then keep it. - // - if (currChi2 > newChi2) { - min->lastDelta = currChi2 - newChi2; - min->value = newChi2; - - // We already masked params. - for (i=0;idata.F32[i] = (float) newParams->data.F32[i]; - } - lamda*= 0.1; - } else { - lamda*= 10.0; - } - psTrace(".psLib.dataManip.psMinimizeLMChi2", 4, - "lamda is %f\n", lamda); - min->iter++; - } - psFree(beta); - psFree(perm); - psFree(paramDeltasF64); - psFree(origParams); - psFree(newParams); - psFree(alpha); - psFree(A); - psFree(aOut); - for (i=0;iiter < min->maxIter) || - (min->lastDelta <= min->tol)) { - return(true); - } - - psTrace(".psLib.dataManip.psMinimizeLMChi2", 4, - "---- psMinimizeLMChi2() end (false) ----\n"); - return(false); -} - -/****************************************************************************** -psMyMinimizeLMChi2(): This routine will take an procedure which calculates -an arbitrary function and it's derivative and minimize the chi-squared match -between that function at the specified coords and the specified value at -those coords. - -XXX: Do this: - After checking that all entries in the paramMask are 1 or 0, when - forming the A matrix from alpha, try this: - - A[i][i] = (1 + lambda*paramask[i]) * alpha[i][i]; - -XXX: This is very different from what is specified in the SDR. Must -coordinate with IfA on new SDR. - -XXX: Do vector/image recycles. - -XXX: probably yErr will be part of the SDR. - -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, - the F64 port might be. - *****************************************************************************/ -bool psMyMinimizeLMChi2(psMinimization *min, - psImage *covar, - psVector *params, - const psVector *paramMask, - const psArray *x, - const psVector *y, - const psVector *yErr, - psMyMinimizeLMChi2Func func) { PS_CHECK_NULL_PTR_RETURN_NULL(min); Index: /trunk/psLib/src/dataManip/psMinimize.h =================================================================== --- /trunk/psLib/src/dataManip/psMinimize.h (revision 2201) +++ /trunk/psLib/src/dataManip/psMinimize.h (revision 2202) @@ -8,6 +8,6 @@ * @author GLG, MHPCC * - * @version $Revision: 1.30 $ $Name: not supported by cvs2svn $ - * @date $Date: 2004-10-26 21:24:43 $ + * @version $Revision: 1.31 $ $Name: not supported by cvs2svn $ + * @date $Date: 2004-10-26 23:14:04 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -83,110 +83,17 @@ ); - -/* - Bug 203: - -typedef float (*psMinimizeLMChi2Func)(psVector *deriv, const psVector *params, -const psVector *x); - -This function takes the current guess for the parameters for which we are -trying to get the best values (params), and a single vector (x) of -conditions from the array of x values fed into the minimiser. - - bool psMinimizeLMChi2(psMinimization *min, - psImage *covar, - psVector *params, - const psVector *paramMask, - const psArray *x, - const psVector *y, - const psVector *yErr, - psMinimizeLMChi2Func func); - -This takes the minimization specs (min), returns the covariance matrix -(covar), takes the best guess of initial parameters (params), the parameter -mask (paramMask), and takes multiple vectors of conditions in an array (x), -the corresponding measured values (y) and errors (yErr) and the function to -fit (func). - -For example, for GRB afterglows, I have flux as a function of time and -frequency. So I stuff into the "psArray *x" all my time and frequency values -(so a whole heap of vectors of size 2), I have the measured values in "y" and -errors in "yErr". Then each of the time-frequency pairs are passed to my model -function with the current parameters, and the model function returns the flux -for that time-frequency pair, and the derivative with respect to each of the -parameters. - -This seems reasonable. The only thing we could change would be to have the -function be defined: - -typedef psVector* (*psMinimizeLMChi2Func)(psImage *deriv, - const psVector *params, - const psArray *x); - -So it would return the model value for each of the measurements at once, and -return for each the derivatives (so that it returns a matrix). - -What do you think? - - - -I'm not sure I understand how LM chi-squared minimization will work with the -following function that you define: - - typedef psVector* (*psMinimizeLMChi2Func)(psImage *deriv, - const psVector *params, - const psArray *x); - -The chi-squared minimization algorithm, as defined in NR, requires that -function to be minimized be evaluated at each data point, and that all -derivatives, with respect to each parameter, be calculated at each data -point. In the above, can I assume that - - x is an array of psVectors, with each vector corresponding to a single - data point. - - The returned value has the same length as x. It contains the value of - the function at each data point in x. - - deriv: an n-by-p matrix where "n" is the number of data points, and "p" - is the number of parameters. The [i][j] element of this matrix - holds the derivative of the function at the i-th data point with - respect to the j-th parameter. - - - -*/ - - -// XXX: What if any of these arguments are NULL? - - typedef -float (*psMinimizeLMChi2Func)(psVector *deriv, - const psVector *params, - const psVector *coords); +psVector* (*psMinimizeLMChi2Func)(psImage *deriv, + const psVector *params, + const psArray *x); bool psMinimizeLMChi2(psMinimization *min, + psImage *covar, psVector *params, const psVector *paramMask, - psImage *covar, - const psArray *coords, - const psVector *value, + const psArray *x, + const psVector *y, + const psVector *yErr, psMinimizeLMChi2Func func); - - -typedef -psVector* (*psMyMinimizeLMChi2Func)(psImage *deriv, - const psVector *params, - const psArray *x); - -bool psMyMinimizeLMChi2(psMinimization *min, - psImage *covar, - psVector *params, - const psVector *paramMask, - const psArray *x, - const psVector *y, - const psVector *yErr, - psMyMinimizeLMChi2Func func); typedef Index: /trunk/psLib/src/math/psMinimize.c =================================================================== --- /trunk/psLib/src/math/psMinimize.c (revision 2201) +++ /trunk/psLib/src/math/psMinimize.c (revision 2202) @@ -9,6 +9,6 @@ * @author GLF, MHPCC * - * @version $Revision: 1.59 $ $Name: not supported by cvs2svn $ - * @date $Date: 2004-10-26 21:24:42 $ + * @version $Revision: 1.60 $ $Name: not supported by cvs2svn $ + * @date $Date: 2004-10-26 23:14:04 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -393,6 +393,6 @@ /****************************************************************************** -psMinimizeLMChi2(): This routine will take an procedure which calculates an -arbitrary function and it's derivative and minimize the chi-squared match +psMinimizeLMChi2(): This routine will take an procedure which calculates +an arbitrary function and it's derivative and minimize the chi-squared match between that function at the specified coords and the specified value at those coords. @@ -417,266 +417,11 @@ *****************************************************************************/ bool psMinimizeLMChi2(psMinimization *min, + psImage *covar, psVector *params, const psVector *paramMask, - psImage *covar, - const psArray *coords, - const psVector *value, + const psArray *x, + const psVector *y, + const psVector *yErr, psMinimizeLMChi2Func func) -{ - PS_CHECK_NULL_PTR_RETURN_NULL(min); - PS_CHECK_NULL_VECTOR_RETURN_NULL(params); - PS_CHECK_EMPTY_VECTOR_RETURN_NULL(params); - PS_CHECK_NULL_PTR_RETURN_NULL(coords); - PS_CHECK_NULL_VECTOR_RETURN_NULL(value); - PS_CHECK_EMPTY_VECTOR_RETURN_NULL(value); - if (paramMask != NULL) { - PS_CHECK_VECTOR_SIZE_EQUAL_RETURN_NULL(params, paramMask); - } - psTrace(".psLib.dataManip.psMinimizeLMChi2", 4, - "---- psMinimizeLMChi2() begin ----\n"); - int numData = value->n; - int numParams = params->n; - int i; - int j; - int k; - int l; - int n; - int p; - psVector *beta = psVectorAlloc(numParams, PS_TYPE_F64); - psVector *perm = psVectorAlloc(numParams, PS_TYPE_F64); - - psVector *paramDeltasF64 = psVectorAlloc(numParams, PS_TYPE_F64); - psVector *origParams = psVectorAlloc(numParams, PS_TYPE_F32); - psVector *newParams = psVectorAlloc(numParams, PS_TYPE_F32); - - psImage *alpha = psImageAlloc(numParams, numParams, PS_TYPE_F32); - psImage *A = psImageAlloc(numParams, numParams, PS_TYPE_F64); - psImage *aOut = psImageAlloc(numParams, numParams, PS_TYPE_F64); - - psVector **deriv = (psVector **) psAlloc(numData * sizeof(psVector *)); - for (i=0;in, PS_TYPE_F32); - psVector *newValueVec = psVectorAlloc(value->n, PS_TYPE_F32); - - float currChi2 = 0.0; - float newChi2 = 0.0; - float lamda = 0.00005; - - psTrace(".psLib.dataManip.psMinimizeLMChi2", 6, - "min->maxIter is %d\n", min->maxIter); - psTrace(".psLib.dataManip.psMinimizeLMChi2", 6, - "min->tol is %f\n", min->tol); - - for (p=0;pdata.F32[p] = params->data.F32[p]; - } - - min->lastDelta = HUGE; - min->iter = 0; - - while ((min->lastDelta > min->tol) && (min->iter < min->maxIter)) { - psTrace(".psLib.dataManip.psMinimizeLMChi2", 4, - "------------------------------------------------------\n"); - psTrace(".psLib.dataManip.psMinimizeLMChi2", 4, - "Iteration %d. Delta is %f\n", min->iter, min->lastDelta); - - // - // Calculate the current values and chi-squared of the function. - // - currChi2 = 0.0; - for (n=0;ndata.F32[n] = func(deriv[n], - params, - (psVector *) coords->data[n]); - currChi2+= (currValueVec->data.F32[n] * currValueVec->data.F32[n]); - } - - for (p=0;pdata.F32[%d] is %f.\n", p, params->data.F32[p]); - } - psTrace(".psLib.dataManip.psMinimizeLMChi2", 6, - "Current chi-squared is (%f)\n", currChi2); - - // - // Mask elements of the derivative for each data point. - // - for (p=0;pdata.U8[p] != 0)) { - for (n=0;ndata.F32[p] = 0.0; - } - } - } - - // - // Calculate the BETA vector. - // - for (p=0;pdata.F64[p] = 0.0; - for (n=0;ndata.F64[p])+= - (value->data.F32[n] - currValueVec->data.F32[n]) * - (deriv[n])->data.F32[p]; - } - // XXX: multiple by -1 here? - (beta->data.F64[p])*= -1.0; - psTrace(".psLib.dataManip.psMinimizeLMChi2", 6, - "beta->data.F64[%d] is %f.\n", p, beta->data.F64[p]); - } - - // - // Calculate the ALPHA matrix. - // - for (k=0;kdata.F32[k][l] = 0.0; - for (n=0;ndata.F32[k][l]+= (deriv[n])->data.F32[k] * - (deriv[n])->data.F32[l]; - } - } - } - - // - // Calculate the matrix A. - // - for (j=0;jdata.F64[j][k] = - (double) ((1.0 + lamda) * alpha->data.F32[j][k]); - } else { - A->data.F64[j][k] = (double) alpha->data.F32[j][k]; - } - } - } - for (j=0;jdata.F64[j][k]); - } - psTrace(".psLib.dataManip.psMinimizeLMChi2", 6, "Matrix A[][]:\n"); - } - - // - // Solve A * alpha = Beta - // - aOut = psMatrixLUD(aOut, perm, A); - paramDeltasF64 = psMatrixLUSolve(paramDeltasF64, aOut, beta, perm); - - // - // Mask any masked parameters. - // - for (i=0;idata.F64[%d] is %f.\n", i, paramDeltasF64->data.F64[i]); - if ((paramMask != NULL) && (paramMask->data.U8[i] != 0)) { - newParams->data.F32[i] = origParams->data.F32[i]; - } else { - newParams->data.F32[i] = params->data.F32[i] - - (float) paramDeltasF64->data.F64[i]; - } - } - - psTrace(".psLib.dataManip.psMinimizeLMChi2", 6, - "Calling func() with new parameters:\n"); - for (i=0;idata.F32[%d] is %f.\n", i, newParams->data.F32[i]); - } - - - // - // Calculate new function values. - // - newChi2 = 0.0; - for (n=0;ndata.F32[n] = func(deriv[n], newParams, - (psVector *) coords->data[n]); - newChi2+= (newValueVec->data.F32[n] * newValueVec->data.F32[n]); - } - psTrace(".psLib.dataManip.psMinimizeLMChi2", 4, - "old/new chi-squareds are (%f, %f)\n", currChi2, newChi2); - - // - // If the new chi-squared is lower, then keep it. - // - if (currChi2 > newChi2) { - min->lastDelta = currChi2 - newChi2; - min->value = newChi2; - - // We already masked params. - for (i=0;idata.F32[i] = (float) newParams->data.F32[i]; - } - lamda*= 0.1; - } else { - lamda*= 10.0; - } - psTrace(".psLib.dataManip.psMinimizeLMChi2", 4, - "lamda is %f\n", lamda); - min->iter++; - } - psFree(beta); - psFree(perm); - psFree(paramDeltasF64); - psFree(origParams); - psFree(newParams); - psFree(alpha); - psFree(A); - psFree(aOut); - for (i=0;iiter < min->maxIter) || - (min->lastDelta <= min->tol)) { - return(true); - } - - psTrace(".psLib.dataManip.psMinimizeLMChi2", 4, - "---- psMinimizeLMChi2() end (false) ----\n"); - return(false); -} - -/****************************************************************************** -psMyMinimizeLMChi2(): This routine will take an procedure which calculates -an arbitrary function and it's derivative and minimize the chi-squared match -between that function at the specified coords and the specified value at -those coords. - -XXX: Do this: - After checking that all entries in the paramMask are 1 or 0, when - forming the A matrix from alpha, try this: - - A[i][i] = (1 + lambda*paramask[i]) * alpha[i][i]; - -XXX: This is very different from what is specified in the SDR. Must -coordinate with IfA on new SDR. - -XXX: Do vector/image recycles. - -XXX: probably yErr will be part of the SDR. - -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, - the F64 port might be. - *****************************************************************************/ -bool psMyMinimizeLMChi2(psMinimization *min, - psImage *covar, - psVector *params, - const psVector *paramMask, - const psArray *x, - const psVector *y, - const psVector *yErr, - psMyMinimizeLMChi2Func func) { PS_CHECK_NULL_PTR_RETURN_NULL(min); Index: /trunk/psLib/src/math/psMinimize.h =================================================================== --- /trunk/psLib/src/math/psMinimize.h (revision 2201) +++ /trunk/psLib/src/math/psMinimize.h (revision 2202) @@ -8,6 +8,6 @@ * @author GLG, MHPCC * - * @version $Revision: 1.30 $ $Name: not supported by cvs2svn $ - * @date $Date: 2004-10-26 21:24:43 $ + * @version $Revision: 1.31 $ $Name: not supported by cvs2svn $ + * @date $Date: 2004-10-26 23:14:04 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -83,110 +83,17 @@ ); - -/* - Bug 203: - -typedef float (*psMinimizeLMChi2Func)(psVector *deriv, const psVector *params, -const psVector *x); - -This function takes the current guess for the parameters for which we are -trying to get the best values (params), and a single vector (x) of -conditions from the array of x values fed into the minimiser. - - bool psMinimizeLMChi2(psMinimization *min, - psImage *covar, - psVector *params, - const psVector *paramMask, - const psArray *x, - const psVector *y, - const psVector *yErr, - psMinimizeLMChi2Func func); - -This takes the minimization specs (min), returns the covariance matrix -(covar), takes the best guess of initial parameters (params), the parameter -mask (paramMask), and takes multiple vectors of conditions in an array (x), -the corresponding measured values (y) and errors (yErr) and the function to -fit (func). - -For example, for GRB afterglows, I have flux as a function of time and -frequency. So I stuff into the "psArray *x" all my time and frequency values -(so a whole heap of vectors of size 2), I have the measured values in "y" and -errors in "yErr". Then each of the time-frequency pairs are passed to my model -function with the current parameters, and the model function returns the flux -for that time-frequency pair, and the derivative with respect to each of the -parameters. - -This seems reasonable. The only thing we could change would be to have the -function be defined: - -typedef psVector* (*psMinimizeLMChi2Func)(psImage *deriv, - const psVector *params, - const psArray *x); - -So it would return the model value for each of the measurements at once, and -return for each the derivatives (so that it returns a matrix). - -What do you think? - - - -I'm not sure I understand how LM chi-squared minimization will work with the -following function that you define: - - typedef psVector* (*psMinimizeLMChi2Func)(psImage *deriv, - const psVector *params, - const psArray *x); - -The chi-squared minimization algorithm, as defined in NR, requires that -function to be minimized be evaluated at each data point, and that all -derivatives, with respect to each parameter, be calculated at each data -point. In the above, can I assume that - - x is an array of psVectors, with each vector corresponding to a single - data point. - - The returned value has the same length as x. It contains the value of - the function at each data point in x. - - deriv: an n-by-p matrix where "n" is the number of data points, and "p" - is the number of parameters. The [i][j] element of this matrix - holds the derivative of the function at the i-th data point with - respect to the j-th parameter. - - - -*/ - - -// XXX: What if any of these arguments are NULL? - - typedef -float (*psMinimizeLMChi2Func)(psVector *deriv, - const psVector *params, - const psVector *coords); +psVector* (*psMinimizeLMChi2Func)(psImage *deriv, + const psVector *params, + const psArray *x); bool psMinimizeLMChi2(psMinimization *min, + psImage *covar, psVector *params, const psVector *paramMask, - psImage *covar, - const psArray *coords, - const psVector *value, + const psArray *x, + const psVector *y, + const psVector *yErr, psMinimizeLMChi2Func func); - - -typedef -psVector* (*psMyMinimizeLMChi2Func)(psImage *deriv, - const psVector *params, - const psArray *x); - -bool psMyMinimizeLMChi2(psMinimization *min, - psImage *covar, - psVector *params, - const psVector *paramMask, - const psArray *x, - const psVector *y, - const psVector *yErr, - psMyMinimizeLMChi2Func func); typedef Index: /trunk/psLib/src/sys/psTrace.h =================================================================== --- /trunk/psLib/src/sys/psTrace.h (revision 2201) +++ /trunk/psLib/src/sys/psTrace.h (revision 2202) @@ -7,8 +7,8 @@ * * @author Robert Lupton, Princeton University - * @author George Gusciora, MHPCC + * @author GLG, MHPCC * - * @version $Revision: 1.23 $ $Name: not supported by cvs2svn $ - * @date $Date: 2004-09-20 22:03:35 $ + * @version $Revision: 1.24 $ $Name: not supported by cvs2svn $ + * @date $Date: 2004-10-26 23:14:04 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii Index: /trunk/psLib/src/sysUtils/psTrace.h =================================================================== --- /trunk/psLib/src/sysUtils/psTrace.h (revision 2201) +++ /trunk/psLib/src/sysUtils/psTrace.h (revision 2202) @@ -7,8 +7,8 @@ * * @author Robert Lupton, Princeton University - * @author George Gusciora, MHPCC + * @author GLG, MHPCC * - * @version $Revision: 1.23 $ $Name: not supported by cvs2svn $ - * @date $Date: 2004-09-20 22:03:35 $ + * @version $Revision: 1.24 $ $Name: not supported by cvs2svn $ + * @date $Date: 2004-10-26 23:14:04 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii Index: /trunk/psLib/test/dataManip/tst_psMinimize06.c =================================================================== --- /trunk/psLib/test/dataManip/tst_psMinimize06.c (revision 2201) +++ /trunk/psLib/test/dataManip/tst_psMinimize06.c (revision 2202) @@ -8,5 +8,5 @@ #define NUM_ITERATIONS 10000 #define ERR_TOL 0.0 -#define N 5 +#define N 20 #define MIN_VALUE 5.0 #define NUM_PARAMS 3 @@ -24,57 +24,7 @@ *****************************************************************************/ -float myFunc(psVector *myDeriv, - psVector *myParams, - psVector *myCoords) -{ - float sum = 0.0; - // float coordData = 0.0; - // float expData = 0.0; - int i; - - if (myDeriv == NULL) { - myDeriv = psVectorAlloc(myParams->n, PS_TYPE_F32); - psError(__func__, "myDeriv is NULL.\n"); - } - - // Simply test that coords were passed in correctly. - /* - for (i=0;idata.F32[0]; - expData = (float) (i+10); - if (fabs(coordData - expData) > FLT_EPSILON) { - printf("ERROR(1): coordinate data was incorrectly passed to myFunc()\n"); - printf("ERROR(1): was (%f) should be (%f)\n", coordData, expData); - testStatus = false; - } - coordData = myCoords->data.F32[1]; - expData = (float) (i+3); - if (fabs(coordData - expData) > FLT_EPSILON) { - printf("ERROR(2): coordinate data was incorrectly passed to myFunc()\n"); - printf("ERROR(2): was (%f) should be (%f)\n", coordData, expData); - testStatus = false; - } - } - */ - - sum = 0.0; - for (i=0;idata.F32[i] - expectedParm[i]) * - (myParams->data.F32[i] - expectedParm[i]); - myDeriv->data.F32[i] = (2.0 * myParams->data.F32[i]) - - (2.0 * expectedParm[i]); - } - // for (i=0;idata.F32[%d] is %f\n", i, myParams->data.F32[i]); - // for (i=0;idata.F32[%d] is %f\n", i, myDeriv->data.F32[i]); - - sum+= MIN_VALUE; - return(sum); -} - -psVector *myFunc2(psImage *myDeriv, - psVector *myParams, - psArray *myCoords) +psVector *myFunc(psImage *myDeriv, + psVector *myParams, + psArray *myCoords) { psVector *sum = psVectorAlloc(myCoords->n, PS_TYPE_F32); @@ -136,12 +86,12 @@ } - psMyMinimizeLMChi2(min, - myCovar, - myParams, - NULL, - myCoords, - y, - NULL, - (psMyMinimizeLMChi2Func) myFunc2); + psMinimizeLMChi2(min, + myCovar, + myParams, + NULL, + myCoords, + y, + NULL, + (psMinimizeLMChi2Func) myFunc); printf("\nThe chi-squared is %f\n", min->value); @@ -165,73 +115,6 @@ } -int t02() -{ - int currentId = psMemGetId(); - int memLeaks = 0; - int i = 0; - psArray *myCoords; - psVector *myParams; - psVector *myParamMask; - psImage *myCovar; - psMinimization *min; - psVector *y = psVectorAlloc(N, PS_TYPE_F32); - - psTraceSetLevel(".psLib", 0); - t02(); - /************************************************************************** - *************************************************************************/ - myParams = psVectorAlloc(NUM_PARAMS, PS_TYPE_F32); - myParamMask = psVectorAlloc(NUM_PARAMS, PS_TYPE_F32); - min = psMinimizationAlloc(NUM_ITERATIONS, ERR_TOL); - myCovar = psImageAlloc(NUM_PARAMS, NUM_PARAMS, PS_TYPE_F32); - - myCoords = psArrayAlloc(N); - for (i=0;idata[i] = (psPTR *) psVectorAlloc(2, PS_TYPE_F32); - ((psVector *) (myCoords->data[i]))->data.F32[0] = (float) (i+10); - ((psVector *) (myCoords->data[i]))->data.F32[1] = (float) (i+3); - y->data.F32[i] = (float) i; - } - for (i=0;idata.F32[i] = (float) i; - myParams->data.F32[i] = expectedParm[i] * 1.3; - myParams->data.F32[i] = (float) (5 + i); - myParams->data.F32[i] = 0.0; - myParamMask->data.U8[i] = 0; - } - - psMinimizeLMChi2(min, - myParams, - NULL, - myCovar, - myCoords, - y, - (psMinimizeLMChi2Func) myFunc); - - printf("\nThe chi-squared is %f\n", min->value); - for (i=0;idata.F32[i], expectedParm[i]); - } - - psFree(myCoords); - psFree(myParams); - psFree(myParamMask); - psFree(min); - psFree(y); - psFree(myCovar); - - psMemCheckCorruption(1); - memLeaks = psMemCheckLeaks(currentId,NULL,NULL); - if (0 != memLeaks) { - psAbort(__func__,"Memory Leaks! (%d leaks)", memLeaks); - } - return (!testStatus); -} - int main() { t01(); - // t02(); } Index: /trunk/psLib/test/sysUtils/tst_psTrace.c =================================================================== --- /trunk/psLib/test/sysUtils/tst_psTrace.c (revision 2201) +++ /trunk/psLib/test/sysUtils/tst_psTrace.c (revision 2202) @@ -23,5 +23,5 @@ {testTrace00, 0, "psTraceSetLevel() and psTraceGetLevel()", 0, false}, {testTrace01, 1, "psTraceSetLevel(): set multiple components in one call", 0, false}, - {testTrace02, 2, "psTraceSetLevel(): test static inheritance", 0, false}, + {testTrace02, 2, "psTraceSetLevel(): test static/dynamic inheritance", 0, false}, {testTrace03, 3, "psTraceReset()", 0, false}, {testTrace04, 4, "psTrace()", 0, false},