#include <FittingFunctionLorentzian.h>
Inherits FittingFunction.
Inheritance diagram for FittingFunctionLorentzian:


Public Member Functions | |
| FittingFunctionLorentzian () | |
| FittingFunctionLorentzian (long lNbIteration, double dLimitSearch) | |
| virtual | ~FittingFunctionLorentzian () |
| virtual PtrFunction | getFunction () const |
| virtual PtrFunctionDerivate | getFunctionDerivate () const |
| virtual PtrFunctionAndDerivate | getFunctionAndDerivate () const |
| virtual int | doFit (long points, gsl_vector *x, gsl_vector *y, gsl_vector *sigma, gsl_vector *guess, bool bEstimateGuess, int iSearchStoppingMethod=1, bool bUseScaled=true) |
| virtual int | generateFunctionFit (double dValMinX, double dResolutionX, long lNbPoints, gsl_vector *funcX, gsl_vector *funcY) |
| virtual int | generateFunctionFit (long lNbPoints, gsl_vector *x, gsl_vector *funcX, gsl_vector *funcY) |
| void | printState (size_t iter, gsl_multifit_fdfsolver *s) |
| void | printResults () |
| double | computeValue (double dX, double dPosition, double dWidth, double dHeight, double dBackground=0.0) |
| double | computeDerivateValue (double dX, double dPosition, double dWidth, double dHeight) |
| TO BE FINISHED. | |
| void | initializeInitialsParameters (gsl_vector *vExperimentalDataX, gsl_vector *vExperimentalDataY, gsl_vector *vInitialGuess, bool bEstimate) |
| void | estimateInitialGuess (gsl_vector *vExperimentalDataX, gsl_vector *vExperimentalDataY) |
| double | estimateInitialPosition (gsl_vector *vExperimentalDataX, gsl_vector *vExperimentalDataY) |
| double | estimateInitialSigma (gsl_vector *vExperimentalDataX, gsl_vector *vExperimentalDataY) |
| double | estimateInitialHeight (gsl_vector *vExperimentalDataX, gsl_vector *vExperimentalDataY) |
| double | estimateInitialBackground (gsl_vector *vExperimentalDataX, gsl_vector *vExperimentalDataY) |
| std::string | getEquation () const |
Static Public Member Functions | |
| static int | Function (const gsl_vector *peak_params, void *calc_data, gsl_vector *func) |
| Implements the Lorentzian function for the model :
The HWHM is equal to :
| |
| static int | FunctionDerivate (const gsl_vector *peak_params, void *calc_data, gsl_matrix *Jacobian) |
| Partial Derivates of the Lorentzian function according each of the variable parameters
The jacobian matrix is filled with the derivate of Yi with the parameters :
| |
| static int | FunctionAndDerivate (const gsl_vector *peak_params, void *calc_data, gsl_vector *func, gsl_matrix *Jacobian) |
Definition at line 18 of file FittingFunctionLorentzian.h.
|
|
Definition at line 22 of file FittingFunctionLorentzian.cpp. References FittingFunction::_dSearchLimit, FittingFunction::_iNbParameters, FittingFunction::_lNbMaxIteration, FittingFunction::_mFunctionMap, FittingFunction::_pParameters, FittingFunction::_pParametersErrors, FittingFunction::_vInitialGuess, NB_ITERATION_MAX, and SEARCH_LIMIT. 00023 { 00024 _iNbParameters = 3; 00025 _lNbMaxIteration = NB_ITERATION_MAX; 00026 _dSearchLimit = SEARCH_LIMIT; 00027 _pParameters = new double[_iNbParameters]; 00028 _pParametersErrors = new double[_iNbParameters]; 00029 _mFunctionMap["Position"] = 0; 00030 _mFunctionMap["Width"] = 1; 00031 _mFunctionMap["Height"] = 2; 00032 _vInitialGuess = gsl_vector_alloc(3); 00033 }
|
|
||||||||||||
|
Definition at line 36 of file FittingFunctionLorentzian.cpp. References FittingFunction::_dSearchLimit, FittingFunction::_iNbParameters, FittingFunction::_lNbMaxIteration, FittingFunction::_mFunctionMap, FittingFunction::_pParameters, FittingFunction::_pParametersErrors, and FittingFunction::_vInitialGuess. 00037 { 00038 _lNbMaxIteration = lNbIteration; 00039 _dSearchLimit = dLimitSearch; 00040 _iNbParameters = 3; 00041 _pParameters = new double[_iNbParameters]; 00042 _pParametersErrors = new double[_iNbParameters]; 00043 _mFunctionMap["Position"] = 0; 00044 _mFunctionMap["Width"] = 1; 00045 _mFunctionMap["Height"] = 2; 00046 _vInitialGuess = gsl_vector_alloc(3); 00047 }
|
|
|
Definition at line 50 of file FittingFunctionLorentzian.cpp. References FittingFunction::_pParameters, FittingFunction::_pParametersErrors, and FittingFunction::_vInitialGuess. 00051 { 00052 printf ("Destruction FittingFunctionLorentzian !!!\n"); 00053 delete[] _pParameters; 00054 delete[] _pParametersErrors; 00055 _pParameters = NULL; 00056 _pParametersErrors = NULL; 00057 gsl_vector_free(_vInitialGuess); 00058 }
|
|
||||||||||||||||||||
|
TO BE FINISHED.
Implements FittingFunction. Definition at line 447 of file FittingFunctionLorentzian.cpp.
|
|
||||||||||||||||||||||||
|
Implements FittingFunction. Definition at line 432 of file FittingFunctionLorentzian.cpp. 00433 { 00434 00435 double dXPos,dXPos2,dWidth2,dYi; 00436 00437 dXPos = dX - dPosition; 00438 dXPos2 = dXPos * dXPos; 00439 dWidth2 = dWidth * dWidth; 00440 00441 dYi = dHeight / (1 + dXPos2 / dWidth2); 00442 00443 return dYi; 00444 }
|
|
||||||||||||||||||||||||||||||||||||
|
Implements FittingFunction. Definition at line 182 of file FittingFunctionLorentzian.cpp. References FittingFunction::_dSearchLimit, FittingFunction::_iNbParameters, FittingFunction::_vInitialGuess, getFunction(), getFunctionAndDerivate(), getFunctionDerivate(), initializeInitialsParameters(), FittingData::n, printState(), FittingData::sigma, FittingData::x, and FittingData::y. 00183 { 00184 initializeInitialsParameters(x,y,guess,bEstimateGuess); 00185 long i; 00186 long iter; 00187 long status; 00188 00189 double* datax = new double[points]; 00190 double* datay = new double[points]; 00191 double* datas = new double[points]; 00192 00193 // fill the data arrays 00194 for (i=0; i<points; i++) 00195 { 00196 datax[i] = gsl_vector_get (x, i); // x = index 00197 datay[i] = gsl_vector_get (y, i); // y = data to be fit 00198 datas[i] = gsl_vector_get (sigma, i); // s = sigma (weight) of data values 00199 } 00200 00201 // gsl data structure 00202 struct FittingData gsl_data; 00203 gsl_data.n = points; 00204 gsl_data.x = datax; 00205 gsl_data.y = datay; 00206 gsl_data.sigma = datas; 00207 00208 // #ifdef DEBUG_FIT 00209 // printf ("Init = %g \t %g\t %g\t \n", gsl_vector_get(_vInitialGuess,0),gsl_vector_get(_vInitialGuess,1),gsl_vector_get(_vInitialGuess,2)); 00210 // #endif 00211 00212 cout << "\n---------------------------------------------" << endl; 00213 cout << "Initials Parameters For The Fit" << endl; 00214 // set-up the solver 00215 const gsl_multifit_fdfsolver_type *T; 00216 if (bUseScaled) 00217 { 00218 T = gsl_multifit_fdfsolver_lmsder; 00219 cout << "Jacobian Data Scaled" << endl; 00220 } 00221 else 00222 { 00223 T = gsl_multifit_fdfsolver_lmder; 00224 cout << "Jacobian Data Not Scaled" << endl; 00225 } 00226 gsl_multifit_fdfsolver *s = gsl_multifit_fdfsolver_alloc (T, points, _iNbParameters); 00227 gsl_multifit_function_fdf f; 00228 00229 f.f = getFunction(); // lorentz_f; 00230 f.df = getFunctionDerivate(); // lorentz_df; 00231 f.fdf = getFunctionAndDerivate(); // lorentz_fdf; 00232 f.n = points; 00233 f.p = _iNbParameters; 00234 f.params = &gsl_data; 00235 00236 switch(iSearchStoppingMethod) 00237 { 00238 case 1 : 00239 cout << "Test Delta Search Stopping Function Used" << endl; 00240 break; 00241 00242 case 2 : 00243 cout << "Test Gradient Search Stopping Function Used" << endl; 00244 break; 00245 } 00246 00247 cout << "\n----------------------------------------------" << endl; 00248 cout << "Initials Values For The Fit " << endl; 00249 cout << "Position = " << gsl_vector_get (_vInitialGuess, 0) << endl; 00250 cout << "Width = " << gsl_vector_get (_vInitialGuess, 1) << endl; 00251 cout << "Height = " << gsl_vector_get (_vInitialGuess, 2) << endl; 00252 00253 00254 gsl_multifit_fdfsolver_set ( s, &f, _vInitialGuess ); 00255 00256 cout << "\n---------------------------------------------" << endl; 00257 cout << "Iterations" << endl; 00258 iter = 0; 00259 00260 #ifdef DEBUG_FIT 00261 printState (iter, s); 00262 #endif 00263 00264 do 00265 { 00266 iter++; 00267 status = gsl_multifit_fdfsolver_iterate (s); 00268 00269 #ifdef DEBUG_FIT 00270 printf ("status = %s\n", gsl_strerror (status)); 00271 printState (iter, s); 00272 #endif 00273 00274 if (status) 00275 { 00276 break; 00277 } 00278 00279 00280 switch(iSearchStoppingMethod) 00281 { 00282 case 1 : 00283 status = gsl_multifit_test_delta (s->dx, s->x, /*_dSearchLimit*/0.0, _dSearchLimit); 00284 break; 00285 00286 case 2 : 00287 gsl_vector* gradient = gsl_vector_alloc(_iNbParameters); 00288 gsl_multifit_gradient(s->J,s->f,gradient); 00289 status = gsl_multifit_test_gradient (gradient, _dSearchLimit); 00290 gsl_vector_free(gradient); 00291 break; 00292 } 00293 } 00294 while (status == GSL_CONTINUE && iter < _lNbMaxIteration); 00295 00296 _lNbIterations = iter; 00297 // Check whether the fitting has converged or failed 00298 //if ( status == GSL_SUCCESS ) 00299 // { 00300 gsl_matrix *covar = gsl_matrix_alloc (_iNbParameters, _iNbParameters); 00301 gsl_multifit_covar (s->J, 0.0, covar); 00302 double chi = gsl_blas_dnrm2(s->f); 00303 00304 _pParameters[0] = FIT(0); 00305 _pParameters[1] = FIT(1); 00306 _pParameters[2] = FIT(2); 00307 00308 _pParametersErrors[0] = ERR(0); 00309 _pParametersErrors[1] = ERR(1); 00310 _pParametersErrors[2] = ERR(2); 00311 00312 //double dHeight = _pParameters[2]; 00313 double dEstimatedWidth = _pParameters[1]; 00314 _dQualityFactor = pow(chi, 2.0)/ (points - _iNbParameters); 00315 // _dSigma = 1.0/ (sqrt(2*M_PI)*dHeight); 00316 // _dFHWM = 2.355*_dSigma; 00317 00318 _dHWHM = dEstimatedWidth; 00319 _dFHWM = 2*_dHWHM; 00320 00321 #ifdef DEBUG_FIT 00322 cout << "\n---------------------------------------------" << endl; 00323 cout << "Results For Lorentzian" << endl; 00324 printResults(); 00325 #endif 00326 00327 gsl_matrix_free(covar); 00328 // } 00329 00330 #ifdef DEBUG_FIT 00331 printf ("status = %s\n", gsl_strerror (status)); 00332 #endif 00333 00334 // Free data 00335 gsl_multifit_fdfsolver_free (s); 00336 00337 delete[] datax; 00338 delete[] datay; 00339 delete[] datas; 00340 00341 datax = NULL; 00342 datax = NULL; 00343 datax = NULL; 00344 00345 return status; 00346 }
Here is the call graph for this function: ![]() |
|
||||||||||||
|
Definition at line 475 of file FittingFunctionLorentzian.cpp. Referenced by estimateInitialHeight(), and estimateInitialSigma().
|
|
||||||||||||
|
Implements FittingFunction. Definition at line 554 of file FittingFunctionLorentzian.cpp. Referenced by initializeInitialsParameters(). 00555 { 00556 // fill the initials parameters values estimation vector for the initialisation 00557 // [0] = position 00558 // [1] = width 00559 // [2] = height 00560 // Thanks to the max index of all the points, we can get the 00561 // position_value = x_value at index_max 00562 // width = standard deviation 00563 // height_value = y_value at index_max 00564 //max_index = gsl_vector_max_index (vExperimentalDataY); 00565 00566 double dInitialPosition = estimateInitialPosition(vExperimentalDataX,vExperimentalDataY); 00567 double dInitialSigma = estimateInitialSigma(vExperimentalDataX,vExperimentalDataY); 00568 double dInitialHeight = estimateInitialHeight(vExperimentalDataX,vExperimentalDataY); 00569 00570 gsl_vector_set(_vInitialGuess, 0, dInitialPosition); 00571 gsl_vector_set(_vInitialGuess, 1, dInitialSigma); 00572 gsl_vector_set(_vInitialGuess, 2, dInitialHeight); 00573 }
|
|
||||||||||||
|
Definition at line 464 of file FittingFunctionLorentzian.cpp. References estimateInitialBackground(). Referenced by estimateInitialSigma(). 00465 { 00466 size_t max_index = gsl_vector_max_index (vExperimentalDataY); 00467 00468 double background = estimateInitialBackground(vExperimentalDataX,vExperimentalDataY); 00469 double maxValue = gsl_vector_get(vExperimentalDataY, max_index); 00470 double height = maxValue - background; 00471 return height; 00472 }
Here is the call graph for this function: ![]() |
|
||||||||||||
|
Definition at line 453 of file FittingFunctionLorentzian.cpp. 00454 { 00455 size_t max_index = gsl_vector_max_index (vExperimentalDataY); 00456 double dPosition = gsl_vector_get(vExperimentalDataX,max_index); 00457 00458 return dPosition; 00459 }
|
|
||||||||||||
|
Definition at line 482 of file FittingFunctionLorentzian.cpp. References estimateInitialBackground(), and estimateInitialHeight(). 00483 { 00484 double dSigma; 00485 double dCovariance; 00486 double dHeight = estimateInitialHeight(vExperimentalDataX,vExperimentalDataY); 00487 double dBackground = estimateInitialBackground(vExperimentalDataX,vExperimentalDataY); 00488 00489 size_t iIndexMax = gsl_vector_max_index(vExperimentalDataY); 00490 double dYMaxValue = gsl_vector_get(vExperimentalDataY, iIndexMax); 00491 00492 double dYCurrentValue = dYMaxValue; 00493 00494 size_t iCurrentIndex = iIndexMax; 00495 size_t lVectorSize = vExperimentalDataX->size;//lSizeX; 00496 double lHalfVectorSize = 0.5 * lVectorSize; 00497 double dHalfMaximum = 0.5 * dHeight + dBackground; 00498 00499 size_t iIndexLow = iCurrentIndex; 00500 size_t iIndexHigh = iCurrentIndex; 00501 00502 if (iIndexMax < lHalfVectorSize) 00503 { 00504 // cout << "\n\t--> Look in increasing direction" << endl; 00505 //- look in increasing direction 00506 while(iCurrentIndex < lVectorSize - 1 && dYCurrentValue > dHalfMaximum) 00507 { 00508 dYCurrentValue = gsl_vector_get(vExperimentalDataY, ++iCurrentIndex); 00509 iIndexHigh = iCurrentIndex; 00510 iIndexLow = iIndexMax - iCurrentIndex; 00511 } 00512 } 00513 else 00514 { 00515 // cout << "\n\t--> Look in decreasing direction" << endl; 00516 //- look in decreasing direction 00517 while( iCurrentIndex > 0 && dYCurrentValue > dHalfMaximum) 00518 { 00519 dYCurrentValue = gsl_vector_get(vExperimentalDataY, --iCurrentIndex); 00520 iIndexLow = iCurrentIndex; 00521 iIndexHigh = iIndexMax + iCurrentIndex; 00522 } 00523 } 00524 00525 if (iCurrentIndex != iIndexMax) 00526 { 00527 double dPosition = gsl_vector_get(vExperimentalDataX, iIndexMax); 00528 double dX2 = gsl_vector_get(vExperimentalDataX, iCurrentIndex); 00529 dSigma = ::fabs(dX2 - dPosition); 00530 } 00531 else 00532 { 00533 //- try an arbitrary value that should be not so far from reality 00534 dSigma = lVectorSize / 6; 00535 } 00536 00537 /* cout << "--------------------------------" << endl; 00538 cout << "Peak : (x,y) --> (" << gsl_vector_get(vExperimentalDataX, iIndexMax) << "," 00539 << gsl_vector_get(vExperimentalDataY, iIndexMax) << ")" << endl; 00540 00541 00542 cout << "Low : index=" << iIndexLow << endl; 00543 cout << "\txlow=" << gsl_vector_get(vExperimentalDataX, iIndexLow) << " | ylow=" << gsl_vector_get(vExperimentalDataY, iIndexLow) << endl; 00544 cout << "\txhigh=" << gsl_vector_get(vExperimentalDataX, iIndexHigh) << " | yhigh=" << gsl_vector_get(vExperimentalDataY, iIndexHigh) << endl; 00545 00546 cout << "\txhigh - xlow= " << gsl_vector_get(vExperimentalDataX, iIndexHigh) - gsl_vector_get(vExperimentalDataX, iIndexLow) << endl; 00547 cout << "\tsigma=" << dSigma << endl; 00548 */ 00549 dCovariance = dSigma * dSigma; 00550 return dSigma; 00551 }
Here is the call graph for this function: ![]() |
|
||||||||||||||||
|
Implements the Lorentzian function for the model :
The HWHM is equal to :
ModelId=43849A02034A.
Definition at line 81 of file FittingFunctionLorentzian.cpp. References FittingData::sigma, FittingData::x, and FittingData::y. Referenced by FunctionAndDerivate(), and getFunction(). 00082 { 00083 size_t i; 00084 size_t iNbPoints = ((struct FittingData *)calc_data)->n; 00085 double *x = ((struct FittingData *)calc_data)->x; 00086 double *y = ((struct FittingData *)calc_data)->y; 00087 double *sigma = ((struct FittingData *) calc_data)->sigma; 00088 00089 double pos = gsl_vector_get (peak_params, 0); 00090 double width = gsl_vector_get (peak_params, 1); 00091 double height = gsl_vector_get (peak_params, 2); 00092 00093 //#ifdef DEBUG_FIT 00094 // printf("pos = %g, width = %g, height = %g\n", pos, width, height); 00095 //#endif 00096 00097 double dXPos,dXPos2,dWidth2,dYi; 00098 00099 for (i=0; i < iNbPoints; i++) 00100 { 00101 dXPos = x[i] - pos; 00102 dXPos2 = dXPos * dXPos; 00103 dWidth2 = width * width; 00104 00105 dYi = height / (1 + dXPos2 / dWidth2); 00106 00107 gsl_vector_set (func, i,( dYi - y[i]) / sigma[i]); 00108 } 00109 00110 return GSL_SUCCESS; 00111 }
|
|
||||||||||||||||||||
|
Definition at line 172 of file FittingFunctionLorentzian.cpp. References Function(), and FunctionDerivate(). Referenced by getFunctionAndDerivate(). 00173 { 00174 Function (peak_params, calc_data, func); 00175 FunctionDerivate (peak_params, calc_data, Jacobian); 00176 return GSL_SUCCESS; 00177 }
Here is the call graph for this function: ![]() |
|
||||||||||||||||
|
Partial Derivates of the Lorentzian function according each of the variable parameters
The jacobian matrix is filled with the derivate of Yi with the parameters :
with
Definition at line 126 of file FittingFunctionLorentzian.cpp. References FittingData::n, FittingData::sigma, and FittingData::x. Referenced by FunctionAndDerivate(), and getFunctionDerivate(). 00127 { 00128 00129 size_t i; 00130 00131 size_t n = ((struct FittingData *)calc_data)->n; 00132 double *x = ((struct FittingData *)calc_data)->x; 00133 double *sigma = ((struct FittingData *)calc_data)->sigma; 00134 00135 double pos = gsl_vector_get (peak_params, 0); 00136 double width = gsl_vector_get (peak_params, 1); 00137 double height = gsl_vector_get (peak_params, 2); 00138 00139 double d_pos, d_width, d_height; 00140 00141 // Fill the Jacobian Matrix 00142 // Jacobian matrix J(i,j) = dfi / dxj, 00143 // where fi = (Yi - yi)/sigma[i], 00144 // Yi = Yi = c / [ 1 + ((x-a)/b)^2] 00145 // and the xj are the parameters (a,b,c) 00146 00147 double s,alpha,beta,e,e2; 00148 for (i=0; i< n; i++) 00149 { 00150 00151 s = sigma[i]; 00152 alpha = (x[i] - pos)/width; 00153 beta = 2 * height / width; 00154 e = (1+alpha*alpha); 00155 e2 = e * e; 00156 00157 d_pos = beta * alpha / e2 / s; // derivate dYi/dpos 00158 d_width = beta * alpha * alpha / e2 / s; // derivate dYi/dwidth 00159 d_height = 1/e/s; // derivate dYi/dheight 00160 00161 // fill the jacobian matrix J(points, parameters) 00162 gsl_matrix_set (Jacobian, i, 0, d_pos); 00163 gsl_matrix_set (Jacobian, i, 1, d_width); 00164 gsl_matrix_set (Jacobian, i, 2, d_height); 00165 } 00166 00167 return GSL_SUCCESS; 00168 }
|
|
||||||||||||||||||||
|
Implements FittingFunction. Definition at line 377 of file FittingFunctionLorentzian.cpp. References FittingFunction::_pParameters. 00378 { 00379 00380 int i; 00381 00382 double dPos = _pParameters[0]; 00383 double dWidth = _pParameters[1]; 00384 double dHeight = _pParameters[2]; 00385 00386 double dXi,dXPos,dXPos2,dWidth2,dYi; 00387 00388 for (i=0; i < lNbPoints; i++) 00389 { 00390 dXi = gsl_vector_get(x,i); 00391 dXPos = dXi - dPos; 00392 dXPos2 = dXPos * dXPos; 00393 dWidth2 = dWidth * dWidth; 00394 00395 dYi = dHeight / (1.0 + dXPos2 / dWidth2); 00396 00397 gsl_vector_set (funcX, i, dXi); 00398 gsl_vector_set (funcY, i, dYi); 00399 } 00400 return 1; 00401 }
|
|
||||||||||||||||||||||||
|
Implements FittingFunction. Definition at line 350 of file FittingFunctionLorentzian.cpp. References FittingFunction::_pParameters. 00350 { 00351 00352 int i; 00353 00354 double dPos = _pParameters[0]; 00355 double dWidth = _pParameters[1]; 00356 double dHeight = _pParameters[2]; 00357 00358 double dXi,dXPos,dXPos2,dWidth2,dYi; 00359 00360 for (i=0; i < lNbPoints; i++) 00361 { 00362 dXi = dValMinX + i*dResolutionX; 00363 dXPos = dXi - dPos; 00364 dXPos2 = dXPos * dXPos; 00365 dWidth2 = dWidth * dWidth; 00366 00367 dYi = dHeight / (1.0 + dXPos2 / dWidth2); 00368 00369 gsl_vector_set (funcX, i, dXi); 00370 gsl_vector_set (funcY, i, dYi); 00371 } 00372 return 1; 00373 }
|
|
|
Implements FittingFunction. Definition at line 15 of file FittingFunctionLorentzian.cpp.
|
|
|
Implements FittingFunction. Definition at line 61 of file FittingFunctionLorentzian.cpp. References Function(). Referenced by doFit(). 00061 { 00062 return FittingFunctionLorentzian::Function; 00063 }
Here is the call graph for this function: ![]() |
|
|
Implements FittingFunction. Definition at line 71 of file FittingFunctionLorentzian.cpp. References FunctionAndDerivate(). Referenced by doFit(). 00071 { 00072 return FittingFunctionLorentzian::FunctionAndDerivate; 00073 }
Here is the call graph for this function: ![]() |
|
|
Implements FittingFunction. Definition at line 66 of file FittingFunctionLorentzian.cpp. References FunctionDerivate(). Referenced by doFit(). 00066 { 00067 return FittingFunctionLorentzian::FunctionDerivate; 00068 }
Here is the call graph for this function: ![]() |
|
||||||||||||||||||||
|
Implements FittingFunction. Definition at line 576 of file FittingFunctionLorentzian.cpp. References FittingFunction::_vInitialGuess, and estimateInitialGuess(). Referenced by doFit(). 00577 { 00578 if (bEstimate) 00579 { 00580 estimateInitialGuess(vExperimentalDataX,vExperimentalDataY); 00581 } 00582 else 00583 { 00584 gsl_vector_set(_vInitialGuess,0,gsl_vector_get(vInitialGuess,0)); 00585 gsl_vector_set(_vInitialGuess,1,gsl_vector_get(vInitialGuess,1)); 00586 gsl_vector_set(_vInitialGuess,2,gsl_vector_get(vInitialGuess,2)); 00587 } 00588 }
Here is the call graph for this function: ![]() |
|
|
Definition at line 420 of file FittingFunctionLorentzian.cpp. References FittingFunction::_dFHWM, FittingFunction::_dHWHM, FittingFunction::_dQualityFactor, FittingFunction::_dSigma, FittingFunction::_pParametersErrors, and FittingFunction::getParameter(). 00421 { 00422 printf ("pos = %g \t+/- %g\n", getParameter("Position"),_pParametersErrors[0]); 00423 printf ("width = %g \t+/- %g\n", getParameter("Width"), _pParametersErrors[1]); 00424 printf ("height = %g \t+/- %g\n", getParameter("Height"), _pParametersErrors[2]); 00425 printf ("chisq/dof = %g\n", _dQualityFactor); 00426 printf ("Sigma = %g\n", _dSigma); 00427 printf ("HWHM = %g\n", _dHWHM); 00428 printf ("FHWM = %g\n", _dFHWM); 00429 }
Here is the call graph for this function: ![]() |
|
||||||||||||
|
Definition at line 409 of file FittingFunctionLorentzian.cpp. Referenced by doFit(). 00410 { 00411 printf ("iter: %3u x : %g %g %g |f(x)|=%g\n", 00412 iter, 00413 gsl_vector_get (s->x, 0), 00414 gsl_vector_get (s->x, 1), 00415 gsl_vector_get (s->x, 2), 00416 gsl_blas_dnrm2 (s->f)); 00417 }
|
1.4.5