Logo ROOT   6.30.04
Reference Guide
 All Namespaces Files Pages
TFitter.cxx
Go to the documentation of this file.
1 // @(#)root/minuit:$Id$
2 // Author: Rene Brun 31/08/99
3 /*************************************************************************
4  * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
5  * All rights reserved. *
6  * *
7  * For the licensing terms see $ROOTSYS/LICENSE. *
8  * For the list of contributors see $ROOTSYS/README/CREDITS. *
9  *************************************************************************/
10 
11 #include "TMinuit.h"
12 #include "TFitter.h"
13 #include "TH1.h"
14 #include "TF1.h"
15 #include "TF2.h"
16 #include "TF3.h"
17 #include "TList.h"
18 #include "TGraph.h"
19 #include "TGraph2D.h"
20 #include "TMultiGraph.h"
21 #include "TMath.h"
22 
23 ////////////////////////////////////////////////////////////////////////////////
24 /// \class TFitter
25 ///
26 /// The ROOT standard fitter based on TMinuit
27 ///
28 ////////////////////////////////////////////////////////////////////////////////
29 
30 // extern void H1FitChisquare(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag);
31 // extern void H1FitLikelihood(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag);
32 // extern void GraphFitChisquare(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag);
33 // extern void Graph2DFitChisquare(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag);
34 // extern void MultiGraphFitChisquare(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag);
35 extern void F2Fit(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag);
36 extern void F3Fit(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag);
37 
38 ClassImp(TFitter);
39 
40 ////////////////////////////////////////////////////////////////////////////////
41 /// Default constructor
42 
43 TFitter::TFitter(Int_t maxpar)
44 {
45  fMinuit = new TMinuit(maxpar);
46  fNlog = 0;
47  fSumLog = 0;
48  fCovar = 0;
49  SetName("MinuitFitter");
50 }
51 
52 ////////////////////////////////////////////////////////////////////////////////
53 /// Default destructor
54 
55 TFitter::~TFitter()
56 {
57  if (fCovar) delete [] fCovar;
58  if (fSumLog) delete [] fSumLog;
59  delete fMinuit;
60 }
61 
62 // ////////////////////////////////////////////////////////////////////////////////
63 // /// return a chisquare equivalent
64 
65 Double_t TFitter::Chisquare(Int_t , Double_t * ) const
66 {
67  Error("Chisquare","This function is deprecated - use ROOT::Fit::Chisquare class");
68  //Double_t amin = 0;
69  //H1FitChisquare(npar,params,amin,params,1);
70  return TMath::QuietNaN();
71 }
72 
73 ////////////////////////////////////////////////////////////////////////////////
74 /// reset the fitter environment
75 
76 void TFitter::Clear(Option_t *)
77 {
78  if (fCovar) {delete [] fCovar; fCovar = 0;}
79  fMinuit->mncler();
80 
81  //reset the internal Minuit random generator to its initial state
82  Double_t val = 3;
83  Int_t inseed = 12345;
84  fMinuit->mnrn15(val,inseed);
85 }
86 
87 ////////////////////////////////////////////////////////////////////////////////
88 /// Execute a fitter command;
89 /// command : command string
90 /// args : list of nargs command arguments
91 
92 Int_t TFitter::ExecuteCommand(const char *command, Double_t *args, Int_t nargs)
93 {
94  if (fCovar) {delete [] fCovar; fCovar = 0;}
95  Int_t ierr = 0;
96  fMinuit->mnexcm(command,args,nargs,ierr);
97  return ierr;
98 }
99 
100 ////////////////////////////////////////////////////////////////////////////////
101 /// Fix parameter ipar.
102 
103 void TFitter::FixParameter(Int_t ipar)
104 {
105  if (fCovar) {delete [] fCovar; fCovar = 0;}
106  fMinuit->FixParameter(ipar);
107 }
108 
109 ////////////////////////////////////////////////////////////////////////////////
110 ///Computes point-by-point confidence intervals for the fitted function
111 ///
112 ///Parameters:
113 /// - n - number of points
114 /// - ndim - dimensions of points
115 /// - x - points, at which to compute the intervals, for ndim > 1
116 /// should be in order: (x0,y0, x1, y1, ... xn, yn)
117 /// - ci - computed intervals are returned in this array
118 /// - cl - confidence level, default=0.95
119 ///
120 ///NOTE, that the intervals are approximate for nonlinear(in parameters) models
121 
122 void TFitter::GetConfidenceIntervals(Int_t n, Int_t ndim, const Double_t *x, Double_t *ci, Double_t cl)
123 {
124  TF1 *f = (TF1*)fUserFunc;
125  Int_t npar = f->GetNumberFreeParameters();
126  Int_t npar_real = f->GetNpar();
127  Double_t *grad = new Double_t[npar_real];
128  Double_t *sum_vector = new Double_t[npar];
129  Bool_t *fixed=0;
130  Double_t al, bl;
131  if (npar_real != npar){
132  fixed = new Bool_t[npar_real];
133  memset(fixed,0,npar_real*sizeof(Bool_t));
134 
135  for (Int_t ipar=0; ipar<npar_real; ipar++){
136  fixed[ipar]=0;
137  f->GetParLimits(ipar,al,bl);
138  if (al*bl != 0 && al >= bl) {
139  //this parameter is fixed
140  fixed[ipar]=1;
141  }
142  }
143  }
144  Double_t c=0;
145 
146  Double_t *matr = GetCovarianceMatrix();
147  if (!matr){
148  delete [] grad;
149  delete [] sum_vector;
150  if (fixed)
151  delete [] fixed;
152  return;
153  }
154 
155  Double_t t = TMath::StudentQuantile(0.5 + cl/2, f->GetNDF());
156  Double_t chidf = TMath::Sqrt(f->GetChisquare()/f->GetNDF());
157  Int_t igrad, ifree=0;
158  for (Int_t ipoint=0; ipoint<n; ipoint++){
159  c=0;
160  f->GradientPar(x+ndim*ipoint, grad);
161  //multiply the covariance matrix by gradient
162  for (Int_t irow=0; irow<npar; irow++){
163  sum_vector[irow]=0;
164  igrad = 0;
165  for (Int_t icol=0; icol<npar; icol++){
166  igrad = 0;
167  ifree=0;
168  if (fixed) {
169  //find the free parameter #icol
170  while (ifree<icol+1){
171  if (fixed[igrad]==0) ifree++;
172  igrad++;
173  }
174  igrad--;
175  //now the [igrad] element of gradient corresponds to [icol] element of cov.matrix
176  } else {
177  igrad = icol;
178  }
179  sum_vector[irow]+=matr[irow*npar_real+icol]*grad[igrad];
180  }
181  }
182  igrad = 0;
183  for (Int_t i=0; i<npar; i++){
184  igrad = 0; ifree=0;
185  if (fixed) {
186  //find the free parameter #icol
187  while (ifree<i+1){
188  if (fixed[igrad]==0) ifree++;
189  igrad++;
190  }
191  igrad--;
192  } else {
193  igrad = i;
194  }
195  c+=grad[igrad]*sum_vector[i];
196  }
197 
198  c=TMath::Sqrt(c);
199  ci[ipoint]=c*t*chidf;
200  }
201 
202  delete [] grad;
203  delete [] sum_vector;
204  if (fixed)
205  delete [] fixed;
206 }
207 
208 ////////////////////////////////////////////////////////////////////////////////
209 ///Computes confidence intervals at level cl. Default is 0.95
210 ///
211 ///The TObject parameter can be a TGraphErrors, a TGraph2DErrors or a TH1,2,3.
212 ///For Graphs, confidence intervals are computed for each point,
213 ///the value of the graph at that point is set to the function value at that
214 ///point, and the graph y-errors (or z-errors) are set to the value of
215 ///the confidence interval at that point.
216 ///For Histograms, confidence intervals are computed for each bin center
217 ///The bin content of this bin is then set to the function value at the bin
218 ///center, and the bin error is set to the confidence interval value.
219 ///NOTE: confidence intervals are approximate for nonlinear models!
220 ///
221 ///Allowed combinations:
222 ///
223 /// - Fitted object Passed object
224 /// - TGraph TGraphErrors, TH1
225 /// - TGraphErrors, AsymmErrors TGraphErrors, TH1
226 /// - TH1 TGraphErrors, TH1
227 /// - TGraph2D TGraph2DErrors, TH2
228 /// - TGraph2DErrors TGraph2DErrors, TH2
229 /// - TH2 TGraph2DErrors, TH2
230 /// - TH3 TH3
231 
232 void TFitter::GetConfidenceIntervals(TObject *obj, Double_t cl)
233 {
234  if (obj->InheritsFrom(TGraph::Class())) {
235  TGraph *gr = (TGraph*)obj;
236  if (!gr->GetEY()){
237  Error("GetConfidenceIntervals", "A TGraphErrors should be passed instead of a graph");
238  return;
239  }
240  if (fObjectFit->InheritsFrom(TGraph2D::Class())){
241  Error("GetConfidenceIntervals", "A TGraph2DErrors should be passed instead of a graph");
242  return;
243  }
244  if (fObjectFit->InheritsFrom(TH1::Class())){
245  if (((TH1*)(fObjectFit))->GetDimension()>1){
246  Error("GetConfidenceIntervals", "A TGraph2DErrors or a TH23 should be passed instead of a graph");
247  return;
248  }
249  }
250  GetConfidenceIntervals(gr->GetN(),1,gr->GetX(), gr->GetEY(), cl);
251  for (Int_t i=0; i<gr->GetN(); i++)
252  gr->SetPoint(i, gr->GetX()[i], ((TF1*)(fUserFunc))->Eval(gr->GetX()[i]));
253  }
254 
255  //TGraph2D/////////////////
256  else if (obj->InheritsFrom(TGraph2D::Class())) {
257  TGraph2D *gr2 = (TGraph2D*)obj;
258  if (!gr2->GetEZ()){
259  Error("GetConfidenceIntervals", "A TGraph2DErrors should be passed instead of a TGraph2D");
260  return;
261  }
262  if (fObjectFit->InheritsFrom(TGraph::Class())){
263  Error("GetConfidenceIntervals", "A TGraphErrors should be passed instead of a TGraph2D");
264  return;
265  }
266  if (fObjectFit->InheritsFrom(TH1::Class())){
267  if (((TH1*)(fObjectFit))->GetDimension()==1){
268  Error("GetConfidenceIntervals", "A TGraphErrors or a TH1 should be passed instead of a graph");
269  return;
270  }
271  }
272  TF2 *f=(TF2*)fUserFunc;
273  Double_t xy[2];
274  Int_t np = gr2->GetN();
275  Int_t npar = f->GetNpar();
276  Double_t *grad = new Double_t[npar];
277  Double_t *sum_vector = new Double_t[npar];
278  Double_t *x = gr2->GetX();
279  Double_t *y = gr2->GetY();
280  Double_t t = TMath::StudentQuantile(0.5 + cl/2, f->GetNDF());
281  Double_t chidf = TMath::Sqrt(f->GetChisquare()/f->GetNDF());
282  Double_t *matr=GetCovarianceMatrix();
283  Double_t c = 0;
284  for (Int_t ipoint=0; ipoint<np; ipoint++){
285  xy[0]=x[ipoint];
286  xy[1]=y[ipoint];
287  f->GradientPar(xy, grad);
288  for (Int_t irow=0; irow<f->GetNpar(); irow++){
289  sum_vector[irow]=0;
290  for (Int_t icol=0; icol<npar; icol++)
291  sum_vector[irow]+=matr[irow*npar+icol]*grad[icol];
292  }
293  c = 0;
294  for (Int_t i=0; i<npar; i++)
295  c+=grad[i]*sum_vector[i];
296  c=TMath::Sqrt(c);
297  gr2->SetPoint(ipoint, xy[0], xy[1], f->EvalPar(xy));
298  gr2->GetEZ()[ipoint]=c*t*chidf;
299 
300  }
301  delete [] grad;
302  delete [] sum_vector;
303  }
304 
305  //TH1/////////////////////////
306  else if (obj->InheritsFrom(TH1::Class())) {
307  if (fObjectFit->InheritsFrom(TGraph::Class())){
308  if (((TH1*)obj)->GetDimension()>1){
309  Error("GetConfidenceIntervals", "Fitted graph and passed histogram have different number of dimensions");
310  return;
311  }
312  }
313  if (fObjectFit->InheritsFrom(TGraph2D::Class())){
314  if (((TH1*)obj)->GetDimension()!=2){
315  Error("GetConfidenceIntervals", "Fitted graph and passed histogram have different number of dimensions");
316  return;
317  }
318  }
319  if (fObjectFit->InheritsFrom(TH1::Class())){
320  if (((TH1*)(fObjectFit))->GetDimension()!=((TH1*)(obj))->GetDimension()){
321  Error("GetConfidenceIntervals", "Fitted and passed histograms have different number of dimensions");
322  return;
323  }
324  }
325 
326 
327  TH1 *hfit = (TH1*)obj;
328  TF1 *f = (TF1*)GetUserFunc();
329  Int_t npar = f->GetNpar();
330  Double_t *grad = new Double_t[npar];
331  Double_t *sum_vector = new Double_t[npar];
332  Double_t x[3];
333 
334  Int_t hxfirst = hfit->GetXaxis()->GetFirst();
335  Int_t hxlast = hfit->GetXaxis()->GetLast();
336  Int_t hyfirst = hfit->GetYaxis()->GetFirst();
337  Int_t hylast = hfit->GetYaxis()->GetLast();
338  Int_t hzfirst = hfit->GetZaxis()->GetFirst();
339  Int_t hzlast = hfit->GetZaxis()->GetLast();
340 
341  TAxis *xaxis = hfit->GetXaxis();
342  TAxis *yaxis = hfit->GetYaxis();
343  TAxis *zaxis = hfit->GetZaxis();
344  Double_t t = TMath::StudentQuantile(0.5 + cl/2, f->GetNDF());
345  Double_t chidf = TMath::Sqrt(f->GetChisquare()/f->GetNDF());
346  Double_t *matr=GetCovarianceMatrix();
347  Double_t c=0;
348  for (Int_t binz=hzfirst; binz<=hzlast; binz++){
349  x[2]=zaxis->GetBinCenter(binz);
350  for (Int_t biny=hyfirst; biny<=hylast; biny++) {
351  x[1]=yaxis->GetBinCenter(biny);
352  for (Int_t binx=hxfirst; binx<=hxlast; binx++) {
353  x[0]=xaxis->GetBinCenter(binx);
354  f->GradientPar(x, grad);
355  for (Int_t irow=0; irow<npar; irow++){
356  sum_vector[irow]=0;
357  for (Int_t icol=0; icol<npar; icol++)
358  sum_vector[irow]+=matr[irow*npar+icol]*grad[icol];
359  }
360  c = 0;
361  for (Int_t i=0; i<npar; i++)
362  c+=grad[i]*sum_vector[i];
363  c=TMath::Sqrt(c);
364  hfit->SetBinContent(binx, biny, binz, f->EvalPar(x));
365  hfit->SetBinError(binx, biny, binz, c*t*chidf);
366  }
367  }
368  }
369  delete [] grad;
370  delete [] sum_vector;
371  }
372  else {
373  Error("GetConfidenceIntervals", "This object type is not supported");
374  return;
375  }
376 
377 }
378 
379 ////////////////////////////////////////////////////////////////////////////////
380 /// return a pointer to the covariance matrix
381 
382 Double_t *TFitter::GetCovarianceMatrix() const
383 {
384  if (fCovar) return fCovar;
385  Int_t npars = fMinuit->GetNumPars();
386  ((TFitter*)this)->fCovar = new Double_t[npars*npars];
387  fMinuit->mnemat(fCovar,npars);
388  return fCovar;
389 }
390 
391 ////////////////////////////////////////////////////////////////////////////////
392 /// return element i,j from the covariance matrix
393 
394 Double_t TFitter::GetCovarianceMatrixElement(Int_t i, Int_t j) const
395 {
396  GetCovarianceMatrix();
397  Int_t npars = fMinuit->GetNumPars();
398  if (i < 0 || i >= npars || j < 0 || j >= npars) {
399  Error("GetCovarianceMatrixElement","Illegal arguments i=%d, j=%d",i,j);
400  return 0;
401  }
402  return fCovar[j+npars*i];
403 }
404 
405 ////////////////////////////////////////////////////////////////////////////////
406 /// return current errors for a parameter
407 /// ipar : parameter number
408 /// eplus : upper error
409 /// eminus : lower error
410 /// eparab : parabolic error
411 /// globcc : global correlation coefficient
412 
413 Int_t TFitter::GetErrors(Int_t ipar,Double_t &eplus, Double_t &eminus, Double_t &eparab, Double_t &globcc) const
414 {
415 
416  Int_t ierr = 0;
417  fMinuit->mnerrs(ipar, eplus,eminus,eparab,globcc);
418  return ierr;
419 }
420 
421 
422 ////////////////////////////////////////////////////////////////////////////////
423 /// return the total number of parameters (free + fixed)
424 
425 Int_t TFitter::GetNumberTotalParameters() const
426 {
427  return fMinuit->fNpar + fMinuit->fNpfix;
428 }
429 
430 ////////////////////////////////////////////////////////////////////////////////
431 /// return the number of free parameters
432 
433 Int_t TFitter::GetNumberFreeParameters() const
434 {
435  return fMinuit->fNpar;
436 }
437 
438 
439 ////////////////////////////////////////////////////////////////////////////////
440 /// return error of parameter ipar
441 
442 Double_t TFitter::GetParError(Int_t ipar) const
443 {
444  Int_t ierr = 0;
445  TString pname;
446  Double_t value,verr,vlow,vhigh;
447 
448  fMinuit->mnpout(ipar, pname,value,verr,vlow,vhigh,ierr);
449  return verr;
450 }
451 
452 
453 ////////////////////////////////////////////////////////////////////////////////
454 /// return current value of parameter ipar
455 
456 Double_t TFitter::GetParameter(Int_t ipar) const
457 {
458  Int_t ierr = 0;
459  TString pname;
460  Double_t value,verr,vlow,vhigh;
461 
462  fMinuit->mnpout(ipar, pname,value,verr,vlow,vhigh,ierr);
463  return value;
464 }
465 
466 ////////////////////////////////////////////////////////////////////////////////
467 /// return current values for a parameter
468 ///
469 /// - ipar : parameter number
470 /// - parname : parameter name
471 /// - value : initial parameter value
472 /// - verr : initial error for this parameter
473 /// - vlow : lower value for the parameter
474 /// - vhigh : upper value for the parameter
475 ///
476 /// WARNING! parname must be suitably dimensionned in the calling function.
477 
478 Int_t TFitter::GetParameter(Int_t ipar, char *parname,Double_t &value,Double_t &verr,Double_t &vlow, Double_t &vhigh) const
479 {
480  Int_t ierr = 0;
481  TString pname;
482  fMinuit->mnpout(ipar, pname,value,verr,vlow,vhigh,ierr);
483  strcpy(parname,pname.Data());
484  return ierr;
485 }
486 
487 ////////////////////////////////////////////////////////////////////////////////
488 /// return name of parameter ipar
489 
490 const char *TFitter::GetParName(Int_t ipar) const
491 {
492  if (!fMinuit || ipar < 0 || ipar > fMinuit->fNu) return "";
493  return fMinuit->fCpnam[ipar];
494 }
495 
496 ////////////////////////////////////////////////////////////////////////////////
497 /// return global fit parameters
498 ///
499 /// - amin : chisquare
500 /// - edm : estimated distance to minimum
501 /// - errdef
502 /// - nvpar : number of variable parameters
503 /// - nparx : total number of parameters
504 
505 Int_t TFitter::GetStats(Double_t &amin, Double_t &edm, Double_t &errdef, Int_t &nvpar, Int_t &nparx) const
506 {
507  Int_t ierr = 0;
508  fMinuit->mnstat(amin,edm,errdef,nvpar,nparx,ierr);
509  return ierr;
510 }
511 
512 ////////////////////////////////////////////////////////////////////////////////
513 /// return Sum(log(i) i=0,n
514 /// used by log likelihood fits
515 
516 Double_t TFitter::GetSumLog(Int_t n)
517 {
518  if (n < 0) return 0;
519  if (n > fNlog) {
520  if (fSumLog) delete [] fSumLog;
521  fNlog = 2*n+1000;
522  fSumLog = new Double_t[fNlog+1];
523  Double_t fobs = 0;
524  for (Int_t j=0;j<=fNlog;j++) {
525  if (j > 1) fobs += TMath::Log(j);
526  fSumLog[j] = fobs;
527  }
528  }
529  if (fSumLog) return fSumLog[n];
530  return 0;
531 }
532 
533 
534 ////////////////////////////////////////////////////////////////////////////////
535 ///return kTRUE if parameter ipar is fixed, kFALSE othersise)
536 
537 Bool_t TFitter::IsFixed(Int_t ipar) const
538 {
539  if (fMinuit->fNiofex[ipar] == 0 ) return kTRUE;
540  return kFALSE;
541 }
542 
543 
544 ////////////////////////////////////////////////////////////////////////////////
545 /// Print fit results
546 
547 void TFitter::PrintResults(Int_t level, Double_t amin) const
548 {
549  fMinuit->mnprin(level,amin);
550 }
551 
552 ////////////////////////////////////////////////////////////////////////////////
553 /// Release parameter ipar.
554 
555 void TFitter::ReleaseParameter(Int_t ipar)
556 {
557  if (fCovar) {delete [] fCovar; fCovar = 0;}
558  fMinuit->Release(ipar);
559 }
560 
561 ////////////////////////////////////////////////////////////////////////////////
562 /// Specify the address of the fitting algorithm
563 
564 void TFitter::SetFCN(void (*fcn)(Int_t &, Double_t *, Double_t &f, Double_t *, Int_t))
565 {
566  if (fCovar) {delete [] fCovar; fCovar = 0;}
567  TVirtualFitter::SetFCN(fcn);
568  fMinuit->SetFCN(fcn);
569 }
570 
571 ////////////////////////////////////////////////////////////////////////////////
572 /// ret fit method (chisquare or loglikelihood)
573 
574 void TFitter::SetFitMethod(const char *name)
575 {
576  if (fCovar) {delete [] fCovar; fCovar = 0;}
577  // if (!strcmp(name,"H1FitChisquare")) SetFCN(H1FitChisquare);
578  // if (!strcmp(name,"H1FitLikelihood")) SetFCN(H1FitLikelihood);
579  // if (!strcmp(name,"GraphFitChisquare")) SetFCN(GraphFitChisquare);
580  // if (!strcmp(name,"Graph2DFitChisquare")) SetFCN(Graph2DFitChisquare);
581  // if (!strcmp(name,"MultiGraphFitChisquare")) SetFCN(MultiGraphFitChisquare);
582  if (!strcmp(name,"F2Minimizer")) SetFCN(F2Fit);
583  if (!strcmp(name,"F3Minimizer")) SetFCN(F3Fit);
584 }
585 
586 ////////////////////////////////////////////////////////////////////////////////
587 /// set initial values for a parameter
588 ///
589 /// - ipar : parameter number
590 /// - parname : parameter name
591 /// - value : initial parameter value
592 /// - verr : initial error for this parameter
593 /// - vlow : lower value for the parameter
594 /// - vhigh : upper value for the parameter
595 
596 Int_t TFitter::SetParameter(Int_t ipar,const char *parname,Double_t value,Double_t verr,Double_t vlow, Double_t vhigh)
597 {
598  if (fCovar) {delete [] fCovar; fCovar = 0;}
599  Int_t ierr = 0;
600  fMinuit->mnparm(ipar,parname,value,verr,vlow,vhigh,ierr);
601  return ierr;
602 }
603 
604 
605 
606 ////////////////////////////////////////////////////////////////////////////////
607 
608 void F2Fit(Int_t &/*npar*/, Double_t * /*gin*/, Double_t &f,Double_t *u, Int_t /*flag*/)
609 {
610  TVirtualFitter *fitter = TVirtualFitter::GetFitter();
611  TF2 *f2 = (TF2*)fitter->GetObjectFit();
612  f2->InitArgs(u, f2->GetParameters() );
613  f = f2->EvalPar(u);
614 }
615 
616 ////////////////////////////////////////////////////////////////////////////////
617 
618 void F3Fit(Int_t &/*npar*/, Double_t * /*gin*/, Double_t &f,Double_t *u, Int_t /*flag*/)
619 {
620  TVirtualFitter *fitter = TVirtualFitter::GetFitter();
621  TF3 *f3 = (TF3*)fitter->GetObjectFit();
622  f3->InitArgs(u, f3->GetParameters() );
623  f = f3->EvalPar(u);
624 }