ROOT logo
/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

/* $Id$ */

#include <TCanvas.h>
#include <TF1.h>
#include <TFormula.h>
#include <TH1F.h>
#include <TH2F.h>
#include <TMath.h>
#include <TMatrix.h>
#include <TRandom.h>
#include <TString.h>
#include <TVector.h>
#include <TVectorD.h>
#include <TVirtualFitter.h>

#include "AliTMinuitToolkit.h"

//--------------------------------------------------------------------------------------
//
// The AliTMinuitToolkit serves as an easy to use interface for the TMinuit 
// package: 
// 
// - It allows to fit a curve to one and two dimensional histograms 
//   (TH2F::Fit() only allows to fit a hyperplane).
// - Or n points can be specified directly via a n x 2 matrix.
// - An option for robust fitting with non-linear functions is implemented.
//
// A small example illustrating the usage of AliTMinuitToolkit is given in the function 
// "AliTMinuitToolkit::Test()".
// 
// 
// 1. Setting the formula:
//
//  The formula is simply set via "void SetFitFunction(TFormula * formula)".
//
//
// 2. Adding the data points
//
//  - In order to fit a histogram, use "void FitHistogram(TH1F * his)" or 
//    "void FitHistogram(TH2F * his)". The fitter is started automatically
//  - Alternatively, the direct specification of the points is possible via
//    "void SetPoints(TMatrixD * points)". Note, that the each point 
//    corresponds to one row in the matrix. The fitter is then started with 
//    the command "void Fit()". The weight of each point can be specified
//    with an n-dimensional vector using "void SetWeights(TVectorD * weights)".
//
//
// 3. Accessing the fit results
//
//  The N parameters of the formula are stored in a N-dimensional vector which
//  is returned by "TVectorD * GetParameters()". In a similar way the covariance 
//  matrix of the fit is returned via "TMatrixD * GetCovarianceMatrix()" which
//  is of the type N x N.
//
//
// 4. Non-linear robust fitting:
//
//  Even a few outliers can lead to wrong results of a least-squares fitting 
//  procedure. In this case the use of robust(resistant) methods can be 
//  helpful, but a stronger dependence on starting values or convergence to
//  local minima can occur.
//
//  The robust option becomes active if EnableRobust(true, sigma) is called. It is
//  very much recommended that a normalization value (scale variable) corresponding 
//  to an expected deviation (sigma) is specified via 
//  "EnableRobust(Bool_t b, Double_t sigma)".
//
//  Performing the fit without knowledge of sigma is also possible if only
//  "EnableRobust(true)" is activated, but this is NOT RECOMMENDED.
//
//  The method is based on another estimator instead of chi^2. For small deviations 
//  the function behaves like x^2 and for larger deviations like |x| - the so 
//  called Huber estimator:
//
//   h(x) = x^2                              , for x < 2.5*sigma
//   h(x) = 2*(2.5*sigma)*x - (2.5*sigma)^2  , for x > 2.5*sigma
//
//  If a weighting function is specified in addition, a second run with the ordinary 
//  metric is started, but before entering the iteration every point is weighted 
//  according to its distance to the outcoming function of the first run. The weighting
//  function w(x) must be defined on the intervall x in [0,1]. w(0) then 
//  corresponds to the weight of the closest point and w(1) to the point with the
//  largest distance.
//
//  Some standard weighting functions are predefined in 
//  "SetWeightFunction(Char_t * name, Float_t param1, Float_t param2 = 0)":
//   - "BOX" equals to 1 if x < param1 and to 0 if x > param1.
//   - "EXPONENTIAL" corresponds to "Math::Exp(-TMath::Log(param1)*x)"
//   - "ERRORFUNCTION" corresponds to "TMath::Erfc((x-param1)/param2)"
//
//
//  REFERENCE for non-linear robust fitting:
//  Ekblom H. and Madsen K. (1988), Alogrithms for non-linear Huber estimation,
//  BIT Numerical Mathematics 29 (1989) 60-76.
//  internet: http://www.springerlink.com/content/m277218542988344/
//
//
// 5. examples:
//
//  A small example illustrating the working principles of AliTMinuitToolkit is given
//  in the function "AliTMinuitToolkit::Test()".
//
//
//
// Comments and questions are always welcome: A.Kalweit@gsi.de
//--------------------------------------------------------------------------------------


ClassImp(AliTMinuitToolkit)

AliTMinuitToolkit::AliTMinuitToolkit() : 
   TNamed(),
   fFormula(0),
   fWeightFunction(0),
   fFitAlgorithm(""),
   fPoints(0),
   fWeights(0),
   fParam(0),
   fParamLimits(0),
   fCovar(0),
   fChi2(0),
   fMaxCalls(0),
   fPrecision(0),
   fUseRobust(0),
   fExpectedSigma(0)
{
 //
 // standard constructor
 //
 fMaxCalls = 500;
 fPrecision = 1;
 fUseRobust = false;
 fExpectedSigma = 0;
}


AliTMinuitToolkit::AliTMinuitToolkit(const AliTMinuitToolkit&) :
   TNamed(),
   fFormula(0),
   fWeightFunction(0),
   fFitAlgorithm(""),
   fPoints(0),
   fWeights(0),
   fParam(0),
   fParamLimits(0),
   fCovar(0),
   fChi2(0),
   fMaxCalls(0),
   fPrecision(0),
   fUseRobust(0),
   fExpectedSigma(0)
{


}


AliTMinuitToolkit& AliTMinuitToolkit::operator=(const AliTMinuitToolkit&) {

 return *this;
}



AliTMinuitToolkit::~AliTMinuitToolkit(){
  //
  // destructor
  //
  delete fPoints;
  delete fWeights;
  delete fWeightFunction;
  delete fParamLimits;
  delete fFormula;
  delete fParam;
  delete fCovar;
  delete fChi2;
}

void AliTMinuitToolkit::FitHistogram(TH1F *const his) {
 //
 // Fit a one dimensional histogram
 //
 fPoints = new TMatrixD(his->GetNbinsX(), 2);
 
 for(Int_t ibin=0; ibin < his->GetNbinsX(); ibin++) {
  Double_t x = his->GetXaxis()->GetBinCenter(ibin+1);
  Double_t y = his->GetBinContent(ibin+1);
  
  (*fPoints)(ibin, 0) = x;
  (*fPoints)(ibin, 1) = y;
 }
 
 Fit();
}


void AliTMinuitToolkit::FitHistogram(TH2F *const his) {
 //
 // Fit a curve to a two dimensional histogram
 //
 fPoints = new TMatrixD((Long64_t)his->GetEntries(), 2);
 Long64_t entry = 0;
 
 for(Int_t ibin=0; ibin < his->GetNbinsX(); ibin++) {
  Double_t x = his->GetXaxis()->GetBinCenter(ibin);
  for(Int_t jbin=0; jbin < his->GetNbinsY(); jbin++) {   
   Long64_t n = his->GetBin(ibin, jbin);
   Double_t y = his->GetYaxis()->GetBinCenter(jbin);
   for(Int_t ientries=0; ientries < his->GetBinContent(n); ientries++) {
    (*fPoints)(entry,0) = x;
    (*fPoints)(entry,1) = y;
    entry++;
   }
   
  }
 }

 Fit();
}


void AliTMinuitToolkit::SetWeightFunction(const Char_t *name, Float_t param1, Float_t param2) {
 //
 // Set the weight function which must be defined on the interval [0,1].
 //
 TString FuncType(name);
 FuncType.ToUpper();
 
 if (FuncType == "EXPONENTIAL") fWeightFunction = new TFormula("exp", Form("TMath::Exp(-TMath::Log(%f)*x)", param1));
 if (FuncType == "BOX") fWeightFunction = new TFormula("box", Form("TMath::Erfc((x-%f)/0.0001)", param1));
 if (FuncType == "ERRORFUNCTION") fWeightFunction = new TFormula("err", Form("TMath::Erfc((x-%f)/%f)", param1, param2));
 
}


void AliTMinuitToolkit::FitterFCN(int &/*npar*/, double */*dummy*/, double &fchisq, double *gin, int /*iflag*/){
  //
  // internal function which gives the specified function to the TMinuit function
  //  

  //
  AliTMinuitToolkit * fitter = (AliTMinuitToolkit*)TVirtualFitter::GetFitter()->GetObjectFit();
  fchisq = 0;
  Int_t nvar       = fitter->GetPoints()->GetNcols()-1;
  Int_t npoints    = fitter->GetPoints()->GetNrows();
  
  // calculate mean deviation for normalization or use user-defined sigma
  Double_t dev = 0.;
  if (fitter->GetExpectedSigma() == 0 && fitter->GetStatus() == true) {
   for (Int_t ipoint=0; ipoint<npoints; ipoint++){
    Double_t x[100];
     for (Int_t ivar=0; ivar<nvar; ivar++){
      x[ivar] = (*fitter->GetPoints())(ipoint, ivar);      
     }
    Float_t funx = fitter->GetFormula()->EvalPar(x,gin);
    Double_t delta = (*fitter->GetPoints())(ipoint, nvar) - funx;
    dev += TMath::Sqrt(TMath::Abs(delta));
   }
   dev = dev/npoints; 
  } else {
   dev = fitter->GetExpectedSigma();
  }
  // calculate chisquare  
  for (Int_t ipoint=0; ipoint<npoints; ipoint++){
    Double_t x[100];
    for (Int_t ivar=0; ivar<nvar; ivar++){
      x[ivar] = (*fitter->GetPoints())(ipoint, ivar);      
    }
    Float_t funx = fitter->GetFormula()->EvalPar(x,gin);   
    Double_t delta = TMath::Abs((*fitter->GetPoints())(ipoint, nvar) - funx);
    if (fitter->GetStatus() == true) {
     delta = delta/dev; // normalization
     if (delta <= 2.5) fchisq+= delta*delta; // new metric: Huber-k-estimator
     if (delta > 2.5) fchisq+= 2*(2.5)*delta - (2.5*2.5);
    } else {
     Double_t weight = (*fitter->GetWeights())(ipoint);
     fchisq+= delta*delta*weight; //old metric
    }
  }
 }
 

void AliTMinuitToolkit::Fit() {
  //
  // internal function that calls the fitter
  //
  Int_t nparam  = fParam->GetNrows();
  Int_t npoints = fPoints->GetNrows();
  Int_t nvar    = fPoints->GetNcols()-1;
  
  // set all paramter limits to infinity as default
  if (fParamLimits == 0) {
   fParamLimits = new TMatrixD(nparam ,2);
   for (Int_t iparam=0; iparam<nparam; iparam++){
    (*fParamLimits)(iparam, 0) = 0;
    (*fParamLimits)(iparam, 1) = 0;
   }
  }
  
  // set all weights to 1 as default
  Bool_t weightFlag = false;
  if (fWeightFunction == 0) {
   fWeightFunction = new TFormula("constant", "1");
  } else {
   weightFlag = true;
  }
  
  // migrad fit algorithm as default
  if (fFitAlgorithm == "") {
   fFitAlgorithm = "migrad";
  }
  
  // assign weights
  if (fWeights == 0) {
   fWeights = new TVectorD(npoints);
   for (Int_t ipoint=0; ipoint<npoints; ipoint++) (*fWeights)(ipoint) = 1;
  }
  
  // set up the fitter
  TVirtualFitter *minuit = TVirtualFitter::Fitter(0, nparam);
  minuit->SetObjectFit(this); 
  minuit->SetFCN(AliTMinuitToolkit::FitterFCN);
  
  // initialize paramters (step size???)
  for (Int_t iparam=0; iparam<nparam; iparam++){
   minuit->SetParameter(iparam, Form("p[%d]",iparam), (*fParam)(iparam), (*fParam)(iparam)/10, (*fParamLimits)(iparam, 0), (*fParamLimits)(iparam, 1));
  }
   
  //
  Double_t argList[2];
  argList[0] = fMaxCalls; //maximal number of calls 
  argList[1] = fPrecision; //tolerance normalized to 0.001 
  if (fMaxCalls == 500 && fPrecision == 1) minuit->ExecuteCommand(fFitAlgorithm, 0, 0); 
  if (fMaxCalls != 500 || fPrecision != 1) minuit->ExecuteCommand(fFitAlgorithm, argList, 2);
  // two additional arguments can be specified ExecuteCommand("migrad", argList, 2) - use 0,0 for default
  
  // fill parameter vector
  for (Int_t ivar=0; ivar<nparam; ivar++){
   (*fParam)(ivar) = minuit->GetParameter(ivar);
   fFormula->SetParameter(ivar, minuit->GetParameter(ivar));
  }

  // if a weight function is specified -> enter 2nd run with weights
  if (weightFlag == true && fUseRobust == true) {
   // sort points for weighting
   Double_t *sortList = new Double_t[npoints];
   Int_t  *indexList = new Int_t[npoints];   
   for (Int_t ipoint=0; ipoint<npoints; ipoint++){
    Double_t funx = fFormula->Eval((*fPoints)(ipoint, 0));
    Double_t delta = TMath::Abs((*fPoints)[ipoint][nvar] - funx);
    sortList[ipoint] = delta;
   } 
   TMath::Sort(npoints, sortList, indexList, false);
   for (Int_t ip=0; ip<npoints; ip++){
    Double_t t = ip/(Double_t)npoints;
    (*fWeights)(indexList[ip]) = fWeightFunction->Eval(t);
   }
   
   // set up the fitter
   fUseRobust = false;
   for (Int_t iparam=0; iparam<nparam; iparam++){
    minuit->SetParameter(iparam, Form("p[%d]",iparam), (*fParam)(iparam), (*fParam)(iparam)/10, (*fParamLimits)(iparam, 0), (*fParamLimits)(iparam, 1));
   }
   // start fitting
   if (fMaxCalls == 500 && fPrecision == 1) minuit->ExecuteCommand(fFitAlgorithm, 0, 0); 
   if (fMaxCalls != 500 || fPrecision != 1) minuit->ExecuteCommand(fFitAlgorithm, argList, 2);
   fUseRobust = true;
   
   delete [] sortList; 
   delete [] indexList;    
  }
  
  // fill parameter vector
  for (Int_t ivar=0; ivar<nparam; ivar++){
   (*fParam)(ivar) = minuit->GetParameter(ivar);
   fFormula->SetParameter(ivar, minuit->GetParameter(ivar));
  }
  
  // fill covariance matrix
  fCovar = new TMatrixD(nparam, nparam);
  for(Int_t i=0; i < nparam; i++) {
   for(Int_t j=0; j < nparam; j++) {
    (*fCovar)(i,j) = minuit->GetCovarianceMatrixElement(i,j);
   }
  }
  
  if (weightFlag == false) fWeightFunction = 0;
}



void AliTMinuitToolkit::Test() {
 //
 // This test function shows the basic working principles of this class 
 // and illustrates how a robust fit can improve the results
 //
 
 // 1. provide some example histogram
 TH1F * hist = new TH1F("test", "with (red) and without (black) robust option", 20,0,4);
 TRandom * rand = new TRandom();
 for (Int_t i = 0; i < 10000; i++) {
  hist->Fill(rand->Exp(1));
  if (i < 1000) hist->Fill(3); //"outliers"
  if (i < 1070) hist->Fill(3.5);
  if (i < 670) hist->Fill(2);
  if (i < 770) hist->Fill(1.5);//"outliers"
  if (i < 740) hist->Fill(1);
 }
 TCanvas * canv = new TCanvas();
 canv->cd(1);
 hist->Draw();
 
 // 2. example fit without robust option
 AliTMinuitToolkit * tool = new AliTMinuitToolkit();
 TFormula *aFormExp = new TFormula("formExp", "[0]*TMath::Exp(-[1]*x)");
 tool->SetFitFunction(aFormExp);
 TVectorD *vec1 = new TVectorD(2); // Set initial values
 (*vec1)(0) = 1800;
 (*vec1)(1) = 1;
 tool->SetInitialParam(vec1);
 tool->FitHistogram(hist);
 
 // draw fit function
 TF1 *func = new TF1("test", "[0]*TMath::Exp(-[1]*x)", 0, 6);
 func->SetParameters((*tool->GetParameters())(0), (*tool->GetParameters())(1));
 func->Draw("same");
 
 // 3 . robust fit 
 TVectorD *vec2 = new TVectorD(2);
 (*vec2)(0) = 1800;
 (*vec2)(1) = 1;
 tool->SetInitialParam(vec2);
 tool->EnableRobust(true, 10);
 tool->SetWeightFunction("box", 0.75);
 tool->FitHistogram(hist);
 TF1 *func2 = new TF1("test2", "[0]*TMath::Exp(-[1]*x)", 0, 6);
 func2->SetParameter(0, (*tool->GetParameters())(0));
 func2->SetParameter(1, (*tool->GetParameters())(1));
 func2->SetLineColor(kRed);
 func2->Draw("same");
 
}

 AliTMinuitToolkit.cxx:1
 AliTMinuitToolkit.cxx:2
 AliTMinuitToolkit.cxx:3
 AliTMinuitToolkit.cxx:4
 AliTMinuitToolkit.cxx:5
 AliTMinuitToolkit.cxx:6
 AliTMinuitToolkit.cxx:7
 AliTMinuitToolkit.cxx:8
 AliTMinuitToolkit.cxx:9
 AliTMinuitToolkit.cxx:10
 AliTMinuitToolkit.cxx:11
 AliTMinuitToolkit.cxx:12
 AliTMinuitToolkit.cxx:13
 AliTMinuitToolkit.cxx:14
 AliTMinuitToolkit.cxx:15
 AliTMinuitToolkit.cxx:16
 AliTMinuitToolkit.cxx:17
 AliTMinuitToolkit.cxx:18
 AliTMinuitToolkit.cxx:19
 AliTMinuitToolkit.cxx:20
 AliTMinuitToolkit.cxx:21
 AliTMinuitToolkit.cxx:22
 AliTMinuitToolkit.cxx:23
 AliTMinuitToolkit.cxx:24
 AliTMinuitToolkit.cxx:25
 AliTMinuitToolkit.cxx:26
 AliTMinuitToolkit.cxx:27
 AliTMinuitToolkit.cxx:28
 AliTMinuitToolkit.cxx:29
 AliTMinuitToolkit.cxx:30
 AliTMinuitToolkit.cxx:31
 AliTMinuitToolkit.cxx:32
 AliTMinuitToolkit.cxx:33
 AliTMinuitToolkit.cxx:34
 AliTMinuitToolkit.cxx:35
 AliTMinuitToolkit.cxx:36
 AliTMinuitToolkit.cxx:37
 AliTMinuitToolkit.cxx:38
 AliTMinuitToolkit.cxx:39
 AliTMinuitToolkit.cxx:40
 AliTMinuitToolkit.cxx:41
 AliTMinuitToolkit.cxx:42
 AliTMinuitToolkit.cxx:43
 AliTMinuitToolkit.cxx:44
 AliTMinuitToolkit.cxx:45
 AliTMinuitToolkit.cxx:46
 AliTMinuitToolkit.cxx:47
 AliTMinuitToolkit.cxx:48
 AliTMinuitToolkit.cxx:49
 AliTMinuitToolkit.cxx:50
 AliTMinuitToolkit.cxx:51
 AliTMinuitToolkit.cxx:52
 AliTMinuitToolkit.cxx:53
 AliTMinuitToolkit.cxx:54
 AliTMinuitToolkit.cxx:55
 AliTMinuitToolkit.cxx:56
 AliTMinuitToolkit.cxx:57
 AliTMinuitToolkit.cxx:58
 AliTMinuitToolkit.cxx:59
 AliTMinuitToolkit.cxx:60
 AliTMinuitToolkit.cxx:61
 AliTMinuitToolkit.cxx:62
 AliTMinuitToolkit.cxx:63
 AliTMinuitToolkit.cxx:64
 AliTMinuitToolkit.cxx:65
 AliTMinuitToolkit.cxx:66
 AliTMinuitToolkit.cxx:67
 AliTMinuitToolkit.cxx:68
 AliTMinuitToolkit.cxx:69
 AliTMinuitToolkit.cxx:70
 AliTMinuitToolkit.cxx:71
 AliTMinuitToolkit.cxx:72
 AliTMinuitToolkit.cxx:73
 AliTMinuitToolkit.cxx:74
 AliTMinuitToolkit.cxx:75
 AliTMinuitToolkit.cxx:76
 AliTMinuitToolkit.cxx:77
 AliTMinuitToolkit.cxx:78
 AliTMinuitToolkit.cxx:79
 AliTMinuitToolkit.cxx:80
 AliTMinuitToolkit.cxx:81
 AliTMinuitToolkit.cxx:82
 AliTMinuitToolkit.cxx:83
 AliTMinuitToolkit.cxx:84
 AliTMinuitToolkit.cxx:85
 AliTMinuitToolkit.cxx:86
 AliTMinuitToolkit.cxx:87
 AliTMinuitToolkit.cxx:88
 AliTMinuitToolkit.cxx:89
 AliTMinuitToolkit.cxx:90
 AliTMinuitToolkit.cxx:91
 AliTMinuitToolkit.cxx:92
 AliTMinuitToolkit.cxx:93
 AliTMinuitToolkit.cxx:94
 AliTMinuitToolkit.cxx:95
 AliTMinuitToolkit.cxx:96
 AliTMinuitToolkit.cxx:97
 AliTMinuitToolkit.cxx:98
 AliTMinuitToolkit.cxx:99
 AliTMinuitToolkit.cxx:100
 AliTMinuitToolkit.cxx:101
 AliTMinuitToolkit.cxx:102
 AliTMinuitToolkit.cxx:103
 AliTMinuitToolkit.cxx:104
 AliTMinuitToolkit.cxx:105
 AliTMinuitToolkit.cxx:106
 AliTMinuitToolkit.cxx:107
 AliTMinuitToolkit.cxx:108
 AliTMinuitToolkit.cxx:109
 AliTMinuitToolkit.cxx:110
 AliTMinuitToolkit.cxx:111
 AliTMinuitToolkit.cxx:112
 AliTMinuitToolkit.cxx:113
 AliTMinuitToolkit.cxx:114
 AliTMinuitToolkit.cxx:115
 AliTMinuitToolkit.cxx:116
 AliTMinuitToolkit.cxx:117
 AliTMinuitToolkit.cxx:118
 AliTMinuitToolkit.cxx:119
 AliTMinuitToolkit.cxx:120
 AliTMinuitToolkit.cxx:121
 AliTMinuitToolkit.cxx:122
 AliTMinuitToolkit.cxx:123
 AliTMinuitToolkit.cxx:124
 AliTMinuitToolkit.cxx:125
 AliTMinuitToolkit.cxx:126
 AliTMinuitToolkit.cxx:127
 AliTMinuitToolkit.cxx:128
 AliTMinuitToolkit.cxx:129
 AliTMinuitToolkit.cxx:130
 AliTMinuitToolkit.cxx:131
 AliTMinuitToolkit.cxx:132
 AliTMinuitToolkit.cxx:133
 AliTMinuitToolkit.cxx:134
 AliTMinuitToolkit.cxx:135
 AliTMinuitToolkit.cxx:136
 AliTMinuitToolkit.cxx:137
 AliTMinuitToolkit.cxx:138
 AliTMinuitToolkit.cxx:139
 AliTMinuitToolkit.cxx:140
 AliTMinuitToolkit.cxx:141
 AliTMinuitToolkit.cxx:142
 AliTMinuitToolkit.cxx:143
 AliTMinuitToolkit.cxx:144
 AliTMinuitToolkit.cxx:145
 AliTMinuitToolkit.cxx:146
 AliTMinuitToolkit.cxx:147
 AliTMinuitToolkit.cxx:148
 AliTMinuitToolkit.cxx:149
 AliTMinuitToolkit.cxx:150
 AliTMinuitToolkit.cxx:151
 AliTMinuitToolkit.cxx:152
 AliTMinuitToolkit.cxx:153
 AliTMinuitToolkit.cxx:154
 AliTMinuitToolkit.cxx:155
 AliTMinuitToolkit.cxx:156
 AliTMinuitToolkit.cxx:157
 AliTMinuitToolkit.cxx:158
 AliTMinuitToolkit.cxx:159
 AliTMinuitToolkit.cxx:160
 AliTMinuitToolkit.cxx:161
 AliTMinuitToolkit.cxx:162
 AliTMinuitToolkit.cxx:163
 AliTMinuitToolkit.cxx:164
 AliTMinuitToolkit.cxx:165
 AliTMinuitToolkit.cxx:166
 AliTMinuitToolkit.cxx:167
 AliTMinuitToolkit.cxx:168
 AliTMinuitToolkit.cxx:169
 AliTMinuitToolkit.cxx:170
 AliTMinuitToolkit.cxx:171
 AliTMinuitToolkit.cxx:172
 AliTMinuitToolkit.cxx:173
 AliTMinuitToolkit.cxx:174
 AliTMinuitToolkit.cxx:175
 AliTMinuitToolkit.cxx:176
 AliTMinuitToolkit.cxx:177
 AliTMinuitToolkit.cxx:178
 AliTMinuitToolkit.cxx:179
 AliTMinuitToolkit.cxx:180
 AliTMinuitToolkit.cxx:181
 AliTMinuitToolkit.cxx:182
 AliTMinuitToolkit.cxx:183
 AliTMinuitToolkit.cxx:184
 AliTMinuitToolkit.cxx:185
 AliTMinuitToolkit.cxx:186
 AliTMinuitToolkit.cxx:187
 AliTMinuitToolkit.cxx:188
 AliTMinuitToolkit.cxx:189
 AliTMinuitToolkit.cxx:190
 AliTMinuitToolkit.cxx:191
 AliTMinuitToolkit.cxx:192
 AliTMinuitToolkit.cxx:193
 AliTMinuitToolkit.cxx:194
 AliTMinuitToolkit.cxx:195
 AliTMinuitToolkit.cxx:196
 AliTMinuitToolkit.cxx:197
 AliTMinuitToolkit.cxx:198
 AliTMinuitToolkit.cxx:199
 AliTMinuitToolkit.cxx:200
 AliTMinuitToolkit.cxx:201
 AliTMinuitToolkit.cxx:202
 AliTMinuitToolkit.cxx:203
 AliTMinuitToolkit.cxx:204
 AliTMinuitToolkit.cxx:205
 AliTMinuitToolkit.cxx:206
 AliTMinuitToolkit.cxx:207
 AliTMinuitToolkit.cxx:208
 AliTMinuitToolkit.cxx:209
 AliTMinuitToolkit.cxx:210
 AliTMinuitToolkit.cxx:211
 AliTMinuitToolkit.cxx:212
 AliTMinuitToolkit.cxx:213
 AliTMinuitToolkit.cxx:214
 AliTMinuitToolkit.cxx:215
 AliTMinuitToolkit.cxx:216
 AliTMinuitToolkit.cxx:217
 AliTMinuitToolkit.cxx:218
 AliTMinuitToolkit.cxx:219
 AliTMinuitToolkit.cxx:220
 AliTMinuitToolkit.cxx:221
 AliTMinuitToolkit.cxx:222
 AliTMinuitToolkit.cxx:223
 AliTMinuitToolkit.cxx:224
 AliTMinuitToolkit.cxx:225
 AliTMinuitToolkit.cxx:226
 AliTMinuitToolkit.cxx:227
 AliTMinuitToolkit.cxx:228
 AliTMinuitToolkit.cxx:229
 AliTMinuitToolkit.cxx:230
 AliTMinuitToolkit.cxx:231
 AliTMinuitToolkit.cxx:232
 AliTMinuitToolkit.cxx:233
 AliTMinuitToolkit.cxx:234
 AliTMinuitToolkit.cxx:235
 AliTMinuitToolkit.cxx:236
 AliTMinuitToolkit.cxx:237
 AliTMinuitToolkit.cxx:238
 AliTMinuitToolkit.cxx:239
 AliTMinuitToolkit.cxx:240
 AliTMinuitToolkit.cxx:241
 AliTMinuitToolkit.cxx:242
 AliTMinuitToolkit.cxx:243
 AliTMinuitToolkit.cxx:244
 AliTMinuitToolkit.cxx:245
 AliTMinuitToolkit.cxx:246
 AliTMinuitToolkit.cxx:247
 AliTMinuitToolkit.cxx:248
 AliTMinuitToolkit.cxx:249
 AliTMinuitToolkit.cxx:250
 AliTMinuitToolkit.cxx:251
 AliTMinuitToolkit.cxx:252
 AliTMinuitToolkit.cxx:253
 AliTMinuitToolkit.cxx:254
 AliTMinuitToolkit.cxx:255
 AliTMinuitToolkit.cxx:256
 AliTMinuitToolkit.cxx:257
 AliTMinuitToolkit.cxx:258
 AliTMinuitToolkit.cxx:259
 AliTMinuitToolkit.cxx:260
 AliTMinuitToolkit.cxx:261
 AliTMinuitToolkit.cxx:262
 AliTMinuitToolkit.cxx:263
 AliTMinuitToolkit.cxx:264
 AliTMinuitToolkit.cxx:265
 AliTMinuitToolkit.cxx:266
 AliTMinuitToolkit.cxx:267
 AliTMinuitToolkit.cxx:268
 AliTMinuitToolkit.cxx:269
 AliTMinuitToolkit.cxx:270
 AliTMinuitToolkit.cxx:271
 AliTMinuitToolkit.cxx:272
 AliTMinuitToolkit.cxx:273
 AliTMinuitToolkit.cxx:274
 AliTMinuitToolkit.cxx:275
 AliTMinuitToolkit.cxx:276
 AliTMinuitToolkit.cxx:277
 AliTMinuitToolkit.cxx:278
 AliTMinuitToolkit.cxx:279
 AliTMinuitToolkit.cxx:280
 AliTMinuitToolkit.cxx:281
 AliTMinuitToolkit.cxx:282
 AliTMinuitToolkit.cxx:283
 AliTMinuitToolkit.cxx:284
 AliTMinuitToolkit.cxx:285
 AliTMinuitToolkit.cxx:286
 AliTMinuitToolkit.cxx:287
 AliTMinuitToolkit.cxx:288
 AliTMinuitToolkit.cxx:289
 AliTMinuitToolkit.cxx:290
 AliTMinuitToolkit.cxx:291
 AliTMinuitToolkit.cxx:292
 AliTMinuitToolkit.cxx:293
 AliTMinuitToolkit.cxx:294
 AliTMinuitToolkit.cxx:295
 AliTMinuitToolkit.cxx:296
 AliTMinuitToolkit.cxx:297
 AliTMinuitToolkit.cxx:298
 AliTMinuitToolkit.cxx:299
 AliTMinuitToolkit.cxx:300
 AliTMinuitToolkit.cxx:301
 AliTMinuitToolkit.cxx:302
 AliTMinuitToolkit.cxx:303
 AliTMinuitToolkit.cxx:304
 AliTMinuitToolkit.cxx:305
 AliTMinuitToolkit.cxx:306
 AliTMinuitToolkit.cxx:307
 AliTMinuitToolkit.cxx:308
 AliTMinuitToolkit.cxx:309
 AliTMinuitToolkit.cxx:310
 AliTMinuitToolkit.cxx:311
 AliTMinuitToolkit.cxx:312
 AliTMinuitToolkit.cxx:313
 AliTMinuitToolkit.cxx:314
 AliTMinuitToolkit.cxx:315
 AliTMinuitToolkit.cxx:316
 AliTMinuitToolkit.cxx:317
 AliTMinuitToolkit.cxx:318
 AliTMinuitToolkit.cxx:319
 AliTMinuitToolkit.cxx:320
 AliTMinuitToolkit.cxx:321
 AliTMinuitToolkit.cxx:322
 AliTMinuitToolkit.cxx:323
 AliTMinuitToolkit.cxx:324
 AliTMinuitToolkit.cxx:325
 AliTMinuitToolkit.cxx:326
 AliTMinuitToolkit.cxx:327
 AliTMinuitToolkit.cxx:328
 AliTMinuitToolkit.cxx:329
 AliTMinuitToolkit.cxx:330
 AliTMinuitToolkit.cxx:331
 AliTMinuitToolkit.cxx:332
 AliTMinuitToolkit.cxx:333
 AliTMinuitToolkit.cxx:334
 AliTMinuitToolkit.cxx:335
 AliTMinuitToolkit.cxx:336
 AliTMinuitToolkit.cxx:337
 AliTMinuitToolkit.cxx:338
 AliTMinuitToolkit.cxx:339
 AliTMinuitToolkit.cxx:340
 AliTMinuitToolkit.cxx:341
 AliTMinuitToolkit.cxx:342
 AliTMinuitToolkit.cxx:343
 AliTMinuitToolkit.cxx:344
 AliTMinuitToolkit.cxx:345
 AliTMinuitToolkit.cxx:346
 AliTMinuitToolkit.cxx:347
 AliTMinuitToolkit.cxx:348
 AliTMinuitToolkit.cxx:349
 AliTMinuitToolkit.cxx:350
 AliTMinuitToolkit.cxx:351
 AliTMinuitToolkit.cxx:352
 AliTMinuitToolkit.cxx:353
 AliTMinuitToolkit.cxx:354
 AliTMinuitToolkit.cxx:355
 AliTMinuitToolkit.cxx:356
 AliTMinuitToolkit.cxx:357
 AliTMinuitToolkit.cxx:358
 AliTMinuitToolkit.cxx:359
 AliTMinuitToolkit.cxx:360
 AliTMinuitToolkit.cxx:361
 AliTMinuitToolkit.cxx:362
 AliTMinuitToolkit.cxx:363
 AliTMinuitToolkit.cxx:364
 AliTMinuitToolkit.cxx:365
 AliTMinuitToolkit.cxx:366
 AliTMinuitToolkit.cxx:367
 AliTMinuitToolkit.cxx:368
 AliTMinuitToolkit.cxx:369
 AliTMinuitToolkit.cxx:370
 AliTMinuitToolkit.cxx:371
 AliTMinuitToolkit.cxx:372
 AliTMinuitToolkit.cxx:373
 AliTMinuitToolkit.cxx:374
 AliTMinuitToolkit.cxx:375
 AliTMinuitToolkit.cxx:376
 AliTMinuitToolkit.cxx:377
 AliTMinuitToolkit.cxx:378
 AliTMinuitToolkit.cxx:379
 AliTMinuitToolkit.cxx:380
 AliTMinuitToolkit.cxx:381
 AliTMinuitToolkit.cxx:382
 AliTMinuitToolkit.cxx:383
 AliTMinuitToolkit.cxx:384
 AliTMinuitToolkit.cxx:385
 AliTMinuitToolkit.cxx:386
 AliTMinuitToolkit.cxx:387
 AliTMinuitToolkit.cxx:388
 AliTMinuitToolkit.cxx:389
 AliTMinuitToolkit.cxx:390
 AliTMinuitToolkit.cxx:391
 AliTMinuitToolkit.cxx:392
 AliTMinuitToolkit.cxx:393
 AliTMinuitToolkit.cxx:394
 AliTMinuitToolkit.cxx:395
 AliTMinuitToolkit.cxx:396
 AliTMinuitToolkit.cxx:397
 AliTMinuitToolkit.cxx:398
 AliTMinuitToolkit.cxx:399
 AliTMinuitToolkit.cxx:400
 AliTMinuitToolkit.cxx:401
 AliTMinuitToolkit.cxx:402
 AliTMinuitToolkit.cxx:403
 AliTMinuitToolkit.cxx:404
 AliTMinuitToolkit.cxx:405
 AliTMinuitToolkit.cxx:406
 AliTMinuitToolkit.cxx:407
 AliTMinuitToolkit.cxx:408
 AliTMinuitToolkit.cxx:409
 AliTMinuitToolkit.cxx:410
 AliTMinuitToolkit.cxx:411
 AliTMinuitToolkit.cxx:412
 AliTMinuitToolkit.cxx:413
 AliTMinuitToolkit.cxx:414
 AliTMinuitToolkit.cxx:415
 AliTMinuitToolkit.cxx:416
 AliTMinuitToolkit.cxx:417
 AliTMinuitToolkit.cxx:418
 AliTMinuitToolkit.cxx:419
 AliTMinuitToolkit.cxx:420
 AliTMinuitToolkit.cxx:421
 AliTMinuitToolkit.cxx:422
 AliTMinuitToolkit.cxx:423
 AliTMinuitToolkit.cxx:424
 AliTMinuitToolkit.cxx:425
 AliTMinuitToolkit.cxx:426
 AliTMinuitToolkit.cxx:427
 AliTMinuitToolkit.cxx:428
 AliTMinuitToolkit.cxx:429
 AliTMinuitToolkit.cxx:430
 AliTMinuitToolkit.cxx:431
 AliTMinuitToolkit.cxx:432
 AliTMinuitToolkit.cxx:433
 AliTMinuitToolkit.cxx:434
 AliTMinuitToolkit.cxx:435
 AliTMinuitToolkit.cxx:436
 AliTMinuitToolkit.cxx:437
 AliTMinuitToolkit.cxx:438
 AliTMinuitToolkit.cxx:439
 AliTMinuitToolkit.cxx:440
 AliTMinuitToolkit.cxx:441
 AliTMinuitToolkit.cxx:442
 AliTMinuitToolkit.cxx:443
 AliTMinuitToolkit.cxx:444
 AliTMinuitToolkit.cxx:445
 AliTMinuitToolkit.cxx:446
 AliTMinuitToolkit.cxx:447
 AliTMinuitToolkit.cxx:448
 AliTMinuitToolkit.cxx:449
 AliTMinuitToolkit.cxx:450
 AliTMinuitToolkit.cxx:451
 AliTMinuitToolkit.cxx:452
 AliTMinuitToolkit.cxx:453
 AliTMinuitToolkit.cxx:454
 AliTMinuitToolkit.cxx:455
 AliTMinuitToolkit.cxx:456
 AliTMinuitToolkit.cxx:457
 AliTMinuitToolkit.cxx:458