FittingFunctionLorentzianWithBackground.cpp

Go to the documentation of this file.
00001 // FittingFunctionLorentzianWithBackground.cpp: implementation of the FittingFunctionLorentzianWithBackground class.
00002 //
00004 #ifdef WIN32 
00005 #pragma warning(disable:4786)
00006 #endif
00007 
00008 #include "FittingFunctionLorentzianWithBackground.h"
00009 
00011 // Construction/Destruction
00013 
00014 std::string     FittingFunctionLorentzianWithBackground::getEquation() const
00015 {
00016         return "f(x) = Height / [1 + (x - Position)˛ / Width˛] + Background";
00017 }
00018 
00019 
00020 //##ModelId=43849A0202D6
00021 FittingFunctionLorentzianWithBackground::FittingFunctionLorentzianWithBackground()
00022 {       
00023         _iNbParameters                          =       4;
00024         _lNbMaxIteration                        =       NB_ITERATION_MAX;
00025         _dSearchLimit                           =       SEARCH_LIMIT;
00026         _pParameters                            =       new double[_iNbParameters];
00027         _pParametersErrors                      =       new double[_iNbParameters];
00028         _mFunctionMap["Position"]       =       0;
00029         _mFunctionMap["Width"]          =       1;
00030         _mFunctionMap["Height"]         =       2;
00031         _mFunctionMap["Background"]     =       3;
00032         _vInitialGuess = gsl_vector_alloc(4);
00033 }
00034 
00035 //##ModelId=43849A0202D7
00036 FittingFunctionLorentzianWithBackground::FittingFunctionLorentzianWithBackground(long lNbIteration,double dLimitSearch) 
00037 {
00038         _lNbMaxIteration                        =       lNbIteration;
00039         _dSearchLimit                           =       dLimitSearch;
00040         _iNbParameters                          =       4;
00041         _pParameters                            =       new double[_iNbParameters];
00042         _pParametersErrors                      =       new double[_iNbParameters];
00043         _mFunctionMap["Position"]       =       0;
00044         _mFunctionMap["Width"]          =       1;
00045         _mFunctionMap["Height"]         =       2;
00046         _mFunctionMap["Background"]     =       3;
00047         _vInitialGuess = gsl_vector_alloc(4);
00048 }
00049 
00050 //##ModelId=43849A0202DA
00051 FittingFunctionLorentzianWithBackground::~FittingFunctionLorentzianWithBackground()
00052 {
00053         printf ("Destruction FittingFunctionLorentzianWithBackground !!!\n");
00054         delete[] _pParameters;
00055         delete[] _pParametersErrors;
00056         _pParameters = NULL;
00057         _pParametersErrors = NULL;
00058         gsl_vector_free(_vInitialGuess);
00059 }
00060 
00061 //##ModelId=43849A0202DC
00062 PtrFunction FittingFunctionLorentzianWithBackground::getFunction() const{
00063         return FittingFunctionLorentzianWithBackground::Function;
00064 }
00065 
00066 //##ModelId=43849A0202E5
00067 PtrFunctionDerivate FittingFunctionLorentzianWithBackground::getFunctionDerivate() const{
00068         return FittingFunctionLorentzianWithBackground::FunctionDerivate;
00069 }
00070 
00071 //##ModelId=43849A0202E7
00072 PtrFunctionAndDerivate FittingFunctionLorentzianWithBackground::getFunctionAndDerivate() const {
00073         return FittingFunctionLorentzianWithBackground::FunctionAndDerivate;
00074 }
00075 
00076 
00081 //##ModelId=43849A0202E9
00082 int FittingFunctionLorentzianWithBackground::Function(const gsl_vector *peak_params, void *calc_data, gsl_vector *func)
00083 {
00084         size_t i;
00085         size_t iNbPoints        = ((struct FittingData *)calc_data)->n;
00086         double *x                       = ((struct FittingData *)calc_data)->x;
00087         double *y                       = ((struct FittingData *)calc_data)->y;
00088         double *sigma           = ((struct FittingData *)calc_data)->sigma;
00089 
00090         double pos                      = gsl_vector_get (peak_params, 0);
00091         double width            = gsl_vector_get (peak_params, 1);
00092         double height           = gsl_vector_get (peak_params, 2);
00093         double background       = gsl_vector_get (peak_params, 3);
00094 
00095 //#ifdef DEBUG_FIT      
00096  // printf("pos = %g, width = %g, height = %g, background = %g\n", pos, width, height, background);
00097 //#endif        
00098         
00099         double dXPos,dXPos2,dWidth2,dYi;
00100 
00101         for (i=0; i < iNbPoints; i++)
00102     {           
00103                 dXPos   = x[i] - pos;
00104                 dXPos2  = dXPos * dXPos;
00105                 dWidth2 = width * width;
00106                         
00107                 dYi = background + height / (1 + dXPos2 / dWidth2);
00108 
00109                 gsl_vector_set (func, i,( dYi - y[i]) / sigma[i]);              
00110     }
00111 
00112   return GSL_SUCCESS;
00113 }
00114 
00115 
00128 //##ModelId=43849A0202EE
00129 int FittingFunctionLorentzianWithBackground::FunctionDerivate(const gsl_vector *peak_params, void *calc_data, gsl_matrix *Jacobian)
00130 {
00131 
00132         size_t i;
00133         
00134         size_t  n = ((struct FittingData *)calc_data)->n;
00135         double *x = ((struct FittingData *)calc_data)->x;
00136         double *sigma = ((struct FittingData *)calc_data)->sigma;
00137 
00138         double pos                      = gsl_vector_get (peak_params, 0);
00139         double width            = gsl_vector_get (peak_params, 1);
00140         double height           = gsl_vector_get (peak_params, 2);
00141         //double background     = gsl_vector_get (peak_params, 3);
00142   
00143         double d_pos, d_width, d_height, d_background;
00144 
00145         // Fill the Jacobian Matrix     
00146         // Jacobian matrix J(i,j) = dfi / dxj, 
00147         // where fi = (Yi - yi)/sigma[i],  
00148         // Yi = d + c / [ 1 + ((x-a)/b)^2] 
00149         // and the xj are the parameters (a,b,c,d) 
00150 
00151         double s,alpha,beta,e,e2;
00152         
00153         for (i=0; i< n; i++)
00154     {
00155                 s                               = sigma[i]; 
00156                 alpha                   = (x[i] - pos)/width;
00157                 beta                    = 2 * height / width;
00158                 e                               = (1+alpha*alpha);
00159                 e2                              = e * e;
00160 
00161                 d_pos                   = beta * alpha / e2 / s;                        // derivate dYi/dPos            
00162                 d_width                 = beta * alpha * alpha / e2 / s;        // derivate dYi/dWidth          
00163                 d_height                = 1/e/s;                                                        // derivate dYi/dHeight
00164                 d_background    = 1/s;                                                          // derivate dYi/dBackground
00165         
00166                 // fill the jacobion matrix J(points, parameters)
00167                 gsl_matrix_set (Jacobian, i, 0, d_pos); 
00168                 gsl_matrix_set (Jacobian, i, 1, d_width);
00169                 gsl_matrix_set (Jacobian, i, 2, d_height);
00170                 gsl_matrix_set (Jacobian, i, 3, d_background);
00171     }
00172          
00173         return GSL_SUCCESS;
00174 }
00175 
00176 // Interface function for the gsl_multifit_fdfsolver
00177 //##ModelId=43849A0202F6
00178 int FittingFunctionLorentzianWithBackground::FunctionAndDerivate(const gsl_vector *peak_params, void *calc_data,gsl_vector *func, gsl_matrix *Jacobian) 
00179 {       
00180         Function  (peak_params, calc_data, func);
00181         FunctionDerivate (peak_params, calc_data, Jacobian);
00182         return GSL_SUCCESS;
00183 }
00184 
00185 
00186 
00187 //##ModelId=43849A0202FC
00188 int FittingFunctionLorentzianWithBackground::doFit(long points, gsl_vector *x, gsl_vector *y,gsl_vector *sigma,gsl_vector *guess,  bool bEstimateGuess,int iSearchStoppingMethod,bool bUseScaled)
00189 {
00190 initializeInitialsParameters(x,y,guess,bEstimateGuess);
00191         long            i;
00192         long            iter;
00193         long            status;
00194         
00195         double* datax = new double[points];
00196         double* datay = new double[points];
00197         double* datas = new double[points];
00198                 
00199         // fill the data arrays 
00200         for (i=0; i<points; i++)
00201         {               
00202                 datax[i] = gsl_vector_get (x, i);               // x = index            
00203                 datay[i] = gsl_vector_get (y, i);               // y = data to be fit    
00204                 datas[i] = gsl_vector_get (sigma, i);   // s = sigma (weight) of data values
00205         }
00206         
00207         // gsl data structure   
00208         struct FittingData gsl_data;
00209         gsl_data.n      = points;
00210         gsl_data.x      = datax;
00211         gsl_data.y              = datay;
00212         gsl_data.sigma  = datas;
00213         
00214 //      #ifdef DEBUG_FIT                
00215 //              printf ("Init = pos=%g \t width=%g\t height=%g\t background=%g\t\n", gsl_vector_get(_vInitialGuess,0),gsl_vector_get(_vInitialGuess,1),gsl_vector_get(_vInitialGuess,2),gsl_vector_get(_vInitialGuess,3));
00216 //      #endif
00217           cout << "\n---------------------------------------------" << endl;
00218         cout << "Initials Parameters For The Fit" << endl;      
00219         // set-up the solver 
00220         const gsl_multifit_fdfsolver_type *T;
00221         if (bUseScaled)
00222         {
00223                 T       = gsl_multifit_fdfsolver_lmsder;
00224                 cout << "Jacobian Data Scaled" << endl;
00225         }
00226         else
00227         {
00228                 T       = gsl_multifit_fdfsolver_lmder;
00229                 cout << "Jacobian Data Not Scaled" << endl;
00230         }
00231         gsl_multifit_fdfsolver *s                               = gsl_multifit_fdfsolver_alloc (T, points, _iNbParameters);
00232         gsl_multifit_function_fdf f;
00233 
00234         f.f             = getFunction();            
00235         f.df            = getFunctionDerivate();    
00236         f.fdf           = getFunctionAndDerivate(); 
00237         f.n             = points;
00238         f.p             = _iNbParameters;
00239         f.params        = &gsl_data;
00240   
00241         gsl_multifit_fdfsolver_set ( s,  &f,  _vInitialGuess );
00242         switch(iSearchStoppingMethod)
00243         {
00244                 case 1 :
00245                         cout << "Test Delta Search Stopping Function Used" << endl;
00246                         break;
00247 
00248                 case 2 :
00249                         cout << "Test Gradient Search Stopping Function Used" << endl;
00250                         break;
00251         }
00252         cout << "\n-------------------------------------------" << endl;
00253         cout << "Initials Values For The Fit " << endl;
00254         cout << "Position   = "         << gsl_vector_get (_vInitialGuess, 0) << endl;
00255         cout << "Width      = "         << gsl_vector_get (_vInitialGuess, 1) << endl;
00256         cout << "Height     = "         << gsl_vector_get (_vInitialGuess, 2) << endl;
00257         cout << "Background = "         << gsl_vector_get (_vInitialGuess, 3) << endl;
00258 
00259  
00260         cout << "\n---------------------------------------------" << endl;
00261         cout << "Iterations" << endl;
00262         iter = 0;
00263         #ifdef DEBUG_FIT
00264           printState (iter, s);
00265         #endif
00266 
00267         do
00268     {
00269                 iter++;
00270                 status = gsl_multifit_fdfsolver_iterate (s);
00271 
00272                 #ifdef DEBUG_FIT
00273                           printf ("status = %s\n", gsl_strerror (status));
00274                           printState (iter, s);
00275                 #endif
00276                 
00277                 if (status)
00278                 {
00279                         break;
00280                 }
00281 
00282                 switch(iSearchStoppingMethod)
00283                 {
00284                 case 1 :
00285                         status = gsl_multifit_test_delta (s->dx, s->x,  /*_dSearchLimit*/0.0, _dSearchLimit);
00286                         break;
00287 
00288                 case 2 :
00289                         gsl_vector* gradient = gsl_vector_alloc(_iNbParameters);
00290                         gsl_multifit_gradient(s->J,s->f,gradient);
00291                         status = gsl_multifit_test_gradient (gradient, _dSearchLimit);
00292                         gsl_vector_free(gradient);
00293                         break;
00294                 }
00295     }
00296         while (status == GSL_CONTINUE && iter < _lNbMaxIteration);
00297 
00298         _lNbIterations = iter;
00299         // Check whether the fitting has converged or failed  
00300         //if ( status == GSL_SUCCESS )
00301         //{
00302         gsl_matrix *covar       = gsl_matrix_alloc (_iNbParameters, _iNbParameters);
00303         gsl_multifit_covar (s->J, 0.0, covar);
00304         double chi                      = gsl_blas_dnrm2(s->f);
00305 
00306         _pParameters[0]                 =       FIT(0);
00307         _pParameters[1]                 =       FIT(1);
00308         _pParameters[2]                 =       FIT(2);
00309         _pParameters[3]                 =       FIT(3);
00310 
00311         _pParametersErrors[0]   =       ERR(0);
00312         _pParametersErrors[1]   =       ERR(1);
00313         _pParametersErrors[2]   =       ERR(2);
00314         _pParametersErrors[3]   =       ERR(3);
00315 
00316         //double dHeight = _pParameters[2];
00317         double dEstimatedWidth = _pParameters[1];
00318         _dQualityFactor =       pow(chi, 2.0)/ (points - _iNbParameters);
00319 //      _dSigma                 =       1.0/ (sqrt(2*M_PI)*dHeight);
00320 //      _dFHWM                  =       2.355*_dSigma;
00321                 
00322         _dHWHM                  = dEstimatedWidth;
00323         _dFHWM                  = 2*_dHWHM;
00324         #ifdef DEBUG_FIT
00325                         cout << "\n---------------------------------------------" << endl;
00326                         cout << "Results For Lorentzian With Background" << endl;
00327                         printResults(); 
00328         #endif  
00329         
00330         gsl_matrix_free(covar);         
00331         //      }
00332 
00333         #ifdef DEBUG_FIT  
00334                 printf ("status = %s\n", gsl_strerror (status));
00335         #endif
00336         
00337         // Free data 
00338         gsl_multifit_fdfsolver_free (s);
00339 
00340         delete[] datax;
00341         delete[] datay;
00342         delete[] datas;
00343         
00344         datax = NULL;
00345         datax = NULL;
00346         datax = NULL;
00347 
00348         return status;
00349 }
00350 
00351 
00352 //##ModelId=43849A02030A
00353 int FittingFunctionLorentzianWithBackground::generateFunctionFit(double dValMinX,double dResolutionX,long lNbPoints,gsl_vector* funcX,gsl_vector* funcY){
00354 
00355         int i;
00356 
00357         double dPos                     =       _pParameters[0];
00358         double dWidth           =       _pParameters[1];
00359         double dHeight          =       _pParameters[2];
00360         double dBackground      =       _pParameters[3];
00361 
00362         double dXi,dXPos,dXPos2,dWidth2,dYi;
00363 
00364         for (i=0; i < lNbPoints; i++)
00365     {           
00366                 dXi             = dValMinX + i*dResolutionX;
00367                 dXPos   = dXi  - dPos;
00368                 dXPos2  = dXPos * dXPos;
00369                 dWidth2 = dWidth * dWidth;
00370                                 
00371                 dYi             = dBackground + dHeight / (1.0 + dXPos2 / dWidth2);
00372 
00373                 gsl_vector_set (funcX, i, dXi);
00374                 gsl_vector_set (funcY, i, dYi);         
00375     }   
00376         return 1;
00377 }
00378 
00379 
00380 
00381 
00382 int FittingFunctionLorentzianWithBackground::generateFunctionFit(long lNbPoints, gsl_vector *x, gsl_vector* funcX,gsl_vector* funcY)
00383 {
00384         int i;
00385 
00386         double dPos                     =       _pParameters[0];
00387         double dWidth           =       _pParameters[1];
00388         double dHeight          =       _pParameters[2];
00389         double dBackground      =       _pParameters[3];
00390 
00391         double dXi,dXPos,dXPos2,dWidth2,dYi;
00392 
00393         for (i=0; i < lNbPoints; i++)
00394     {           
00395                 dXi             = gsl_vector_get(x,i);
00396                 dXPos   = dXi  - dPos;
00397                 dXPos2  = dXPos * dXPos;
00398                 dWidth2 = dWidth * dWidth;
00399                                 
00400                 dYi             = dBackground + dHeight / (1.0 + dXPos2 / dWidth2);
00401 
00402                 gsl_vector_set (funcX, i, dXi);
00403                 gsl_vector_set (funcY, i, dYi);         
00404     }   
00405         return 1;
00406 }
00407 
00408 
00409 
00410 
00411 // print the state of the fitting during iterations 
00412 //##ModelId=43849A020315
00413 void FittingFunctionLorentzianWithBackground::printState (size_t iter, gsl_multifit_fdfsolver * s)
00414 {
00415   printf ("iter: %3u x : %g %g %g %g |f(x)|=%g\n",
00416           iter,
00417           gsl_vector_get (s->x, 0), 
00418           gsl_vector_get (s->x, 1),
00419           gsl_vector_get (s->x, 2), 
00420                   gsl_vector_get (s->x, 3), 
00421           gsl_blas_dnrm2 (s->f));
00422 }
00423 
00424 //##ModelId=43849A020318
00425 void FittingFunctionLorentzianWithBackground::printResults()
00426 {
00427                 printf ("pos            = %g \t+/- %g\n",       getParameter("Position")        ,       _pParametersErrors[0]);
00428                 printf ("width          = %g \t+/- %g\n",       getParameter("Width")           ,       _pParametersErrors[1]);
00429                 printf ("height         = %g \t+/- %g\n",       getParameter("Height")          ,       _pParametersErrors[2]);
00430                 printf ("background     = %g \t+/- %g\n",       getParameter("Background")      ,       _pParametersErrors[3]);
00431         printf ("chisq/dof      = %g\n",                                _dQualityFactor);
00432                 printf ("Sigma          = %g\n",                                _dSigma);
00433                 printf ("HWHM           = %g\n",                                _dHWHM);
00434                 printf ("FHWM           = %g\n",                                _dFHWM);
00435 }
00436 
00437 double  FittingFunctionLorentzianWithBackground::computeValue(double dX,double dPosition, double dWidth, double dHeight, double dBackground)
00438 {
00439 
00440         double dXPos,dXPos2,dWidth2,dYi;
00441 
00442         dXPos   = dX - dPosition;
00443         dXPos2  = dXPos * dXPos;
00444         dWidth2 = dWidth * dWidth;
00445                         
00446         dYi = dBackground + dHeight / (1 + dXPos2 / dWidth2);
00447 
00448         return dYi;
00449 }
00450 
00451 
00453 double  FittingFunctionLorentzianWithBackground::computeDerivateValue(double dX,double dPosition, double dWidth, double dHeight)
00454 {
00455         return 0.0;
00456 }
00457 
00458 
00459 
00460 
00461 double FittingFunctionLorentzianWithBackground::estimateInitialPosition(gsl_vector *vExperimentalDataX, gsl_vector *vExperimentalDataY)
00462 {
00463         size_t max_index        = gsl_vector_max_index (vExperimentalDataY);
00464         double dPosition        = gsl_vector_get(vExperimentalDataX,max_index);
00465 
00466         return dPosition;
00467 }
00468 
00469 
00470 
00471 
00472 double FittingFunctionLorentzianWithBackground::estimateInitialHeight(gsl_vector *vExperimentalDataX, gsl_vector *vExperimentalDataY)
00473 {
00474         size_t max_index        = gsl_vector_max_index (vExperimentalDataY);
00475 
00476         double background       = estimateInitialBackground(vExperimentalDataX,vExperimentalDataY);
00477         double maxValue         = gsl_vector_get(vExperimentalDataY, max_index);
00478         double height           = maxValue - background;
00479         return height;
00480 }
00481 
00482 
00483 double FittingFunctionLorentzianWithBackground::estimateInitialBackground(gsl_vector *vExperimentalDataX, gsl_vector *vExperimentalDataY)
00484 {
00485         double dBackground = 0.0;
00486         size_t min_index = gsl_vector_min_index (vExperimentalDataY);
00487         //TODO do a mean with more points to better estimate background
00488         dBackground =  gsl_vector_get(vExperimentalDataY, min_index);
00489         
00490         return dBackground;
00491 }
00492 
00493 
00494 double FittingFunctionLorentzianWithBackground::estimateInitialSigma(gsl_vector *vExperimentalDataX, gsl_vector *vExperimentalDataY)
00495 {
00496         double dSigma;
00497         double dCovariance;
00498         double dHeight                          = estimateInitialHeight(vExperimentalDataX,vExperimentalDataY);
00499         double dBackground                      = estimateInitialBackground(vExperimentalDataX,vExperimentalDataY);
00500         
00501         size_t iIndexMax                        = gsl_vector_max_index(vExperimentalDataY);
00502         double dYMaxValue                       = gsl_vector_get(vExperimentalDataY, iIndexMax);
00503   
00504         double  dYCurrentValue          = dYMaxValue;
00505 
00506         size_t iCurrentIndex            = iIndexMax;
00507         size_t lVectorSize                      = vExperimentalDataX->size;//lSizeX;
00508         double lHalfVectorSize          = 0.5 * lVectorSize;
00509         double dHalfMaximum                     = 0.5 * dHeight + dBackground;
00510   
00511         size_t iIndexLow                        = iCurrentIndex;
00512         size_t iIndexHigh                       = iCurrentIndex;
00513 
00514     if (iIndexMax < lHalfVectorSize)
00515     {
00516                 //cout << "\n\t--> Look in increasing direction" << endl;
00517                 //- look in increasing direction
00518                 while(iCurrentIndex < lVectorSize - 1 && dYCurrentValue > dHalfMaximum)
00519                 {
00520                         dYCurrentValue  = gsl_vector_get(vExperimentalDataY, ++iCurrentIndex);
00521                         iIndexHigh              = iCurrentIndex;
00522                         iIndexLow               = iIndexMax - iCurrentIndex;
00523                 }
00524     }
00525     else
00526     {
00527                 //cout << "\n\t--> Look in decreasing direction" << endl;
00528                 //- look in decreasing direction
00529                 while( iCurrentIndex > 0 && dYCurrentValue > dHalfMaximum)
00530                 {
00531                         dYCurrentValue =  gsl_vector_get(vExperimentalDataY, --iCurrentIndex);
00532                         iIndexLow               = iCurrentIndex;
00533                         iIndexHigh              = iIndexMax + iCurrentIndex;
00534                 }
00535     }
00536 
00537     if (iCurrentIndex != iIndexMax)
00538     {
00539                 double dPosition        = gsl_vector_get(vExperimentalDataX, iIndexMax);
00540                 double dX2                      = gsl_vector_get(vExperimentalDataX, iCurrentIndex);
00541                 dSigma = ::fabs(dX2 - dPosition);
00542     }
00543     else
00544     {
00545                 //- try an arbitrary value that should be not so far from reality
00546                 dSigma = lVectorSize / 6;
00547     }
00548 
00549 /*      cout << "----------------------------------" << endl;
00550         cout << "Peak : (x,y) --> ("    << gsl_vector_get(vExperimentalDataX, iIndexMax) << "," 
00551                                                                         << gsl_vector_get(vExperimentalDataY, iIndexMax) << ")" << endl;
00552                 
00553         
00554         cout << "Low : index=" << iIndexLow << endl;
00555         cout << "\txlow=" << gsl_vector_get(vExperimentalDataX, iIndexLow) << " | ylow=" << gsl_vector_get(vExperimentalDataY, iIndexLow) << endl;
00556         cout << "\txhigh=" << gsl_vector_get(vExperimentalDataX, iIndexHigh) << " | yhigh=" << gsl_vector_get(vExperimentalDataY, iIndexHigh) << endl;
00557         
00558         cout << "\txhigh - xlow= " << gsl_vector_get(vExperimentalDataX, iIndexHigh) - gsl_vector_get(vExperimentalDataX, iIndexLow) << endl;
00559         cout << "\tsigma=" << dSigma << endl;
00560 */
00561     dCovariance = dSigma * dSigma;
00562         return dSigma;
00563 }
00564 
00565 
00566 
00567 void    FittingFunctionLorentzianWithBackground::estimateInitialGuess(gsl_vector* vExperimentalDataX,gsl_vector* vExperimentalDataY)
00568 {               
00569         // fill the initials parameters values estimation vector for the initialisation
00570         // [0] = position 
00571         // [1] = width  
00572         // [2] = height
00573         // [3] = background
00574         // Thanks to the max index of all the points, we can get the
00575         // position_value       = x_value at index_max
00576         // width                        = standard deviation
00577         // height_value         = y_value - min_value (=background)     at index_max
00578         // background_value = minimal value
00579         //min_index = gsl_vector_min_index (vExperimentalDataY);
00580         //max_index = gsl_vector_max_index (vExperimentalDataY);
00581         double dInitialPosition = estimateInitialPosition(vExperimentalDataX,vExperimentalDataY);
00582         double dInitialSigma    = estimateInitialSigma(vExperimentalDataX,vExperimentalDataY);
00583         double dInitialHeight   = estimateInitialHeight(vExperimentalDataX,vExperimentalDataY);
00584         double dInitialBackground = estimateInitialBackground(vExperimentalDataX,vExperimentalDataY);
00585 
00586         gsl_vector_set(_vInitialGuess, 0, dInitialPosition);
00587         gsl_vector_set(_vInitialGuess, 1, dInitialSigma);
00588         gsl_vector_set(_vInitialGuess, 2, dInitialHeight);
00589         gsl_vector_set(_vInitialGuess, 3, dInitialBackground);  
00590 }
00591 
00592 void    FittingFunctionLorentzianWithBackground::initializeInitialsParameters(gsl_vector* vExperimentalDataX,gsl_vector* vExperimentalDataY,gsl_vector* vInitialGuess,bool bEstimate)
00593 {
00594         if (bEstimate) 
00595         {
00596                 estimateInitialGuess(vExperimentalDataX,vExperimentalDataY);
00597         }
00598         else
00599         {
00600                 gsl_vector_set(_vInitialGuess,0,gsl_vector_get(vInitialGuess,0));
00601                 gsl_vector_set(_vInitialGuess,1,gsl_vector_get(vInitialGuess,1));
00602                 gsl_vector_set(_vInitialGuess,2,gsl_vector_get(vInitialGuess,2));
00603                 gsl_vector_set(_vInitialGuess,3,gsl_vector_get(vInitialGuess,3));
00604         }
00605 }
00606 

Generated on Tue Apr 14 09:34:28 2009 for Data Fitting Library by  doxygen 1.4.5