PWG1/TPC/AliPerfAnalyzeInvPt.cxx

   1
   2 /**************************************************************************
   3  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   4  *                                                                        *
   5  * Author: The ALICE Off-line Project.                                    *
   6  * Contributors are mentioned in the code where appropriate.              *
   7  *                                                                        *
   8  * Permission to use, copy, modify and distribute this software and its   *
   9  * documentation strictly for non-commercial purposes is hereby granted   *
  10  * without fee, provided that the above copyright notice appears in all   *
  11  * copies and that both the copyright notice and this permission notice   *
  12  * appear in the supporting documentation. The authors make no claims     *
  13  * about the suitability of this software for any purpose. It is          *
  14  * provided "as is" without express or implied warranty.                  *
  15  **************************************************************************/
  16
  17 //------------------------------------------------------------------------------
  18 // Implementation of AliPerfAnalyzeInvPt class. It analyzes the output of
  19 // AliPerformancePtCalib.cxx  and AliPerformancePtCalibMC.cxx: via
  20 // AliPerformancePtCalib::Analyse() or AliPerformancePtCalibMC::Analyse()
  21 //
  22 // Projection of charge/pt vs theta and vs phi resp. histoprams will be fitted with either
  23 // polynomial or gaussian fit function to extract minimum position of charge/pt.
  24 // Fit options and theta, phi bins can be set by user.
  25 // Attention: use the Set* functions of AliPerformancePtCalib.h and AliPerformancePtCalibMC.h
  26 // when running AliPerformancePtCalib*::Analyse()
  27 // The result of the analysis (histograms/graphs) are stored in the folder which is
  28 // a data member of AliPerformancePtCalib*.
  29 //
  30 // Author: S.Schuchmann 11/13/2009
  31 //------------------------------------------------------------------------------
  32
  33 /*
  34
  35 // after running the performance task of AliPerformancePtCalib* or AliPerformancePtCalibMC*
  36 // read the file, and get component:
  37
  38 TFile f("Output.root");
  39 AliPerformancePtCalib *compObj = (AliPerformancePtCalib*)coutput->FindObject("AliPerformancePtCalib");
  40 //or
  41 AliPerformancePtCalibMC *compObj = (AliPerformancePtCalibMC*)coutput->FindObject("AliPerformancePtCalibMC");
  42
  43 // set phi and theta bins for fitting and analyse comparison data
  44 compObj->SetProjBinsTheta(thetaBins,nThetaBins,minPhi, maxPhi);
  45 compObj->SetProjBinsPhi(phiBins,nPhiBins,minTheta,maxTheta);
  46 compObj->SetMakeFitOption(kFALSE,exclRange,fitRange);
  47 compObj->SetDoRebin(rebin);
  48 compObj->Analyse();
  49 //for details see functions of class AliPerformancePtCalib*
  50
  51 // the output histograms/graphs will be stored in the folder "folderRes"
  52 compObj->GetAnalysisFolder()->ls("*");
  53
  54 // user can save whole comparison object (or only folder with anlysed histograms)
  55 // in the seperate output file (e.g.)
  56 TFile fout("Analysed_InvPt.root","recreate");
  57 compObj->Write(); // compObj->GetAnalysisFolder()->Write();
  58 fout.Close();
  59
  60 */
  61
  62 #include "TNamed.h"
  63 #include "TF1.h"
  64 #include "TH1F.h"
  65 #include "TH2F.h"
  66 #include "TH1D.h"
  67 #include "TMath.h"
  68 #include "TGraphErrors.h"
  69 #include "TCanvas.h"
  70 #include "TObjArray.h"
  71
  72 #include "AliPerfAnalyzeInvPt.h"
  73
  74 ClassImp(AliPerfAnalyzeInvPt)
  75
  76
  77 // fit functions
  78 //____________________________________________________________________________________________________________________________________________
  79    Double_t AliPerfAnalyzeInvPt::Polynomial(Double_t *x, const Double_t *par)
  80 {
  81    // fit function for fitting of 1/pt with a polynomial of 4th order, rejecting points between  +/-par[4]
  82
  83
  84    if (x[0] > -par[4] && x[0] < par[4]) {
  85       TF1::RejectPoint();
  86       return 0;
  87    }
  88
  89    return  par[2]+par[0]*pow((x[0]-par[1]),2)+par[3]*pow((x[0]-par[1]),4);
  90
  91 }
  92 //____________________________________________________________________________________________________________________________________________
  93 Double_t AliPerfAnalyzeInvPt::PolynomialRejP(Double_t *x, const Double_t *par)
  94 {
  95    // fit function for fitting of 1/pt with a polynomial of 4th order to improve result (fit range and rejection zone is adjusted to first guess of minimum position)
  96
  97    Double_t pos  = par[5];
  98    Double_t neg = - par[5];
  99    pos += par[4];
 100    neg += par[4];
 101
 102    if (x[0] > neg && x[0] < pos) {
 103       TF1::RejectPoint();
 104       return 0;
 105    }
 106
 107    return   par[2]+par[0]*pow((x[0]-par[1]),2)+par[3]*pow((x[0]-par[1]),4);
 108
 109 }
 110 //____________________________________________________________________________________________________________________________________________
 111 Double_t AliPerfAnalyzeInvPt::InvGauss(Double_t *x, const Double_t *par)
 112 {
 113    // fit function for fitting of 1/pt with gaussian, rejecting points between  +/-par[4]
 114    if (x[0] > -par[6] && x[0] < par[6]) {
 115       TF1::RejectPoint();
 116       return 0;
 117    }
 118
 119    return par[3]+par[0]*(TMath::Exp(-0.5*(TMath::Power((x[0]-par[1])/par[2], 2.0)))+par[4]*pow((x[0]-par[1]),2)+par[5]*pow((x[0]-par[1]),4)) ;
 120 }
 121 //____________________________________________________________________________________________________________________________________________
 122 Double_t AliPerfAnalyzeInvPt::InvGaussRejP(Double_t *x, const Double_t *par)
 123 {
 124    // fit function for fitting of 1/pt with gaussian to improve result (fit range and rejection zone is adjusted to first guess of minimum position)
 125    Double_t pos  = par[7];//0.12;
 126    Double_t neg = - par[7];//0.12;
 127    pos += par[6];
 128    neg += par[6];
 129
 130    if (x[0] > neg && x[0] < pos) {
 131       TF1::RejectPoint();
 132       return 0;
 133    }
 134
 135
 136    return par[3]+par[0]*TMath::Exp(-0.5*(TMath::Power((x[0]-par[1])/par[2], 2.0)))+par[4]*pow((x[0]-par[1]),2)+par[5]*pow((x[0]-par[1]),4) ;
 137 }
 138
 139
 140 //_____________________________________________________________________________________________________________________________________________
 141 AliPerfAnalyzeInvPt::AliPerfAnalyzeInvPt():
 142    TNamed("AliPerfAnalyzeInvPt","AliPerfAnalyzeInvPt"),
 143    fNThetaBins(0),
 144    fNPhiBins(0),
 145    fRange(0),
 146    fExclRange(0),
 147    fFitGaus(0) ,
 148    fDoRebin(0),
 149    fRebin(0),
 150    //histograms,graphs and functions
 151    fHistH2InvPtTheta(0),
 152    fHistH2InvPtPhi(0),
 153    fGrMinPosTheta(0),
 154    fGrMinPosPhi(0),
 155    fFitMinPos(0),
 156    fFitMinPosRejP(0),
 157    fFitInvGauss(0),
 158    fFitInvGaussRejP(0)
 159 {
 160    // Default constructor
 161
 162    fFitGaus = kFALSE;// flag for gaussian fit
 163    fNThetaBins = 0; // theta bins for projections
 164    fNPhiBins = 0; //phi bins for projections
 165    fRange = 0; //fit range
 166    fExclRange = 0; //fit exclusion range
 167    fDoRebin = kFALSE; // flag for rebin
 168    fRebin = 0; // bins for rebin
 169
 170    for(Int_t i=0;i<100;i++){
 171
 172       fHistFitTheta[i] = NULL;
 173       fHistFitPhi[i] = NULL;
 174    }
 175
 176
 177 }
 178 //_____________________________________________________________________________________________________________________________________________
 179 AliPerfAnalyzeInvPt::AliPerfAnalyzeInvPt(const Char_t* name="AliAnalyzeInvPt",const Char_t* title="AliAnalyzeInvPt"):
 180    TNamed(name, title),
 181    fNThetaBins(0),
 182    fNPhiBins(0),
 183    fRange(0),
 184    fExclRange(0),
 185    fFitGaus(0) ,
 186    fDoRebin(0),
 187    fRebin(0),
 188    //histograms,graphs and functions
 189    fHistH2InvPtTheta(0),
 190    fHistH2InvPtPhi(0),
 191    fGrMinPosTheta(0),
 192    fGrMinPosPhi(0),
 193    fFitMinPos(0),
 194    fFitMinPosRejP(0),
 195    fFitInvGauss(0),
 196    fFitInvGaussRejP(0)
 197 {
 198    fFitGaus = kFALSE;// flag for gaussian fit
 199    fNThetaBins = 0; // theta bins for projections
 200    fNPhiBins = 0; //phi bins for projections
 201    fRange = 0; //fit range
 202    fExclRange = 0; //fit exclusion range
 203    fDoRebin = kFALSE; // flag for rebin
 204    fRebin = 0; // bins for rebin
 205
 206    for(Int_t i=0;i<100;i++){
 207
 208       fHistFitTheta[i] = NULL;
 209       fHistFitPhi[i] = NULL;
 210    }
 211
 212 }
 213
 214
 215 //______________________________________________________________________________________________________________________________________
 216 void AliPerfAnalyzeInvPt::InitGraphs(Double_t *binsXTheta,Double_t *fitParamTheta,Double_t *errFitParamTheta,Double_t *binsXPhi,Double_t *fitParamPhi,Double_t *errFitParamPhi){
 217
 218    // initialize graphs for storing fit parameter
 219
 220    fGrMinPosTheta = new TGraphErrors(fNThetaBins,binsXTheta,fitParamTheta,0, errFitParamTheta);
 221    fGrMinPosTheta->SetMarkerStyle(20);
 222    fGrMinPosTheta->SetMarkerColor(2);
 223    fGrMinPosTheta->SetLineColor(2);
 224    fGrMinPosTheta->GetYaxis()->SetTitle("min pos (GeV/c)^{-1}");
 225    fGrMinPosTheta->GetXaxis()->SetTitle("#theta bin no.");
 226    fGrMinPosTheta->GetYaxis()->SetTitleOffset(1.2);
 227    fGrMinPosTheta->SetTitle("#theta bins");
 228
 229    fGrMinPosPhi = new TGraphErrors(fNPhiBins,binsXPhi,fitParamPhi,0,errFitParamPhi);
 230    fGrMinPosPhi->SetMarkerStyle(20);
 231    fGrMinPosPhi->SetMarkerColor(4);
 232    fGrMinPosPhi->SetLineColor(4);
 233    fGrMinPosPhi->GetYaxis()->SetTitle("min pos (GeV/c)^{-1}");
 234    fGrMinPosPhi->GetXaxis()->SetTitle("#phi bin no.");
 235    fGrMinPosPhi->GetYaxis()->SetTitleOffset(1.2);
 236    fGrMinPosPhi->SetTitle("#phi bins");
 237 }
 238
 239 //______________________________________________________________________________________________________________________________________
 240
 241 void AliPerfAnalyzeInvPt::InitFitFcn(){
 242    // fit functions
 243    fFitMinPos = new TF1("fFitMinPos", Polynomial,-4.0,4.0,5);
 244    fFitMinPos->SetLineColor(4);
 245    fFitMinPos->SetLineWidth(1);
 246    fFitMinPos->SetParameter(0,1.0);
 247    fFitMinPos->SetParameter(1,0.0);
 248    fFitMinPos->SetParameter(2,0.0);
 249
 250    fFitMinPosRejP = new TF1("fFitMinPosRejP",PolynomialRejP,-4.0,4.0,6);
 251    fFitMinPosRejP->SetLineColor(2);
 252    fFitMinPosRejP->SetLineWidth(1);
 253    fFitMinPosRejP->SetParameter(0,1.0);
 254    fFitMinPosRejP->SetParameter(1,0.0);
 255    fFitMinPosRejP->SetParameter(2,0.0);
 256
 257
 258    fFitInvGauss = new TF1("fFitInvGauss", InvGauss,-4.0,4.0,6);
 259    fFitInvGauss->SetLineColor(4);
 260    fFitInvGauss->SetLineWidth(1);
 261    fFitInvGauss->SetParameter(2,1.0);
 262
 263    fFitInvGaussRejP = new TF1("fFitInvGaussRejP", InvGaussRejP,-4.0,4.0,7);
 264    fFitInvGaussRejP->SetLineColor(2);
 265    fFitInvGaussRejP->SetLineWidth(1);
 266    fFitInvGaussRejP->SetParameter(2,1.0);
 267
 268 }
 269 //______________________________________________________________________________________________________________________________________
 270 void AliPerfAnalyzeInvPt::StartAnalysis(const TH2F *histThetaInvPt, const TH2F *histPhiInvPt, TObjArray* aFolderObj){
 271    //start analysis: fitting charge/pt spectra
 272
 273
 274
 275    if(!histThetaInvPt) {
 276       Printf("warning: no 1/pt histogram to analyse in theta bins!");
 277       return;
 278    }
 279    if(!histPhiInvPt) {
 280       Printf("warning: no 1/pt histogram to analyse in phit bins!");
 281       return;
 282    }
 283
 284    Double_t thetaBins[9] = {0.77,0.97,1.17,1.37,1.57,1.77,1.97,2.17,2.37};                 // default theta bins
 285    Double_t phiBins[13] = {0.0,0.5,1.0,1.5,2.0,2.5,3.0,3.5,4.0,4.5,5.0,5.5,6.5};           // default phi bins
 286
 287    Int_t nThetaBins = 9;
 288    Int_t nPhiBins = 13;
 289    Double_t range = 1.0;        // fit range
 290    Double_t  exclRange  = 0.13; // range of point rejection in fit
 291    Bool_t fitGaus=kFALSE;       // fit with Gaussian or with polynomial (kFALSE)
 292
 293    if(!fNThetaBins){
 294       fNThetaBins = nThetaBins;
 295       for(Int_t k = 0;k<fNThetaBins;k++){
 296          fThetaBins[k] = thetaBins[k];
 297       }
 298       Printf("warning: for theta bins no user intput! bins are set to default values.");
 299    }
 300
 301    if(!fNPhiBins){
 302       fNPhiBins = nPhiBins;
 303       for(Int_t k = 0;k<fNPhiBins;k++){
 304          fPhiBins[k] = phiBins[k];
 305       }
 306
 307       Printf("warning: for phi bins no user intput! bins are set to default values.");
 308
 309    }
 310
 311    if(!fRange){
 312       fRange = range;
 313       fExclRange = exclRange;
 314       fFitGaus = fitGaus;
 315       Printf("warning: no fit range is set. default fitting conditions are used.");
 316    }
 317
 318
 319    //fit parameter arrays
 320    Double_t fitParamTheta[100],fitParamPhi[100];
 321    Double_t errFitParamTheta[100], errFitParamPhi[100];
 322    Double_t binsXTheta[100],binsXPhi[100];
 323
 324    //copy histograms
 325    fHistH2InvPtTheta  = (TH2F*)histThetaInvPt->Clone("invPtVsTheta");
 326    fHistH2InvPtPhi    = (TH2F*)histPhiInvPt->Clone("invPtVsPhi");
 327
 328    // initialize fit functions
 329    InitFitFcn();
 330
 331    // canvas for 2D histograms (input)
 332    TCanvas *thephiCan =new TCanvas("thephiCan","theta and phi vs invPt",800,400);
 333    thephiCan->Divide(2,1);
 334    thephiCan->cd(1);
 335    fHistH2InvPtTheta->SetTitle("#theta vs charge/pt for selected #phi range");
 336    fHistH2InvPtTheta->GetXaxis()->SetTitle("charge/p_{t} (GeV/c)^{-1}");
 337    fHistH2InvPtTheta->GetYaxis()->SetTitle("#theta (rad)");
 338    fHistH2InvPtTheta->DrawCopy("col1z");
 339
 340    thephiCan->cd(2);
 341    fHistH2InvPtPhi->SetTitle("#phi vs charge/pt for selected #theta range");
 342    fHistH2InvPtPhi->GetXaxis()->SetTitle("charge/p_{t} (GeV/c)^{-1}");
 343    fHistH2InvPtPhi->GetYaxis()->SetTitle("#phi (rad)");
 344    fHistH2InvPtPhi->DrawCopy("col1z");
 345
 346    // canvas for 1D histograms (projections)
 347    TCanvas *theCan =new TCanvas("theCan","invPt theta bins",1200,900);
 348    theCan->Divide((fNThetaBins+2)/2,2);
 349    TCanvas *phiCan =new TCanvas("phiCan","invPt phi bins",1200,900);
 350    phiCan->Divide((fNPhiBins+2)/2,2);
 351
 352    Int_t countPad = 1;
 353
 354    // analyse charge/pt in bins of theta
 355    for(Int_t i=0;i<fNThetaBins;i++){
 356
 357       //set name for each projection
 358       TString name = "fit_theta_";
 359       name +=i;
 360       fHistH2InvPtTheta->SetName(name.Data());
 361
 362       //find bins for projection
 363       Int_t firstBin= fHistH2InvPtTheta->GetYaxis()->FindBin(fThetaBins[i]);
 364       if(i>0) firstBin +=1;
 365       Int_t lastBin = fHistH2InvPtTheta->GetYaxis()->FindBin(fThetaBins[i+1]);
 366       if( i == fNThetaBins-1) {
 367          firstBin= fHistH2InvPtTheta->GetYaxis()->FindBin(fThetaBins[0]);
 368          lastBin = fHistH2InvPtTheta->GetYaxis()->FindBin(fThetaBins[i]);
 369       }
 370
 371       //make projection
 372       fHistFitTheta[i] =  (TH1D*)fHistH2InvPtTheta->ProjectionX("_px",firstBin,lastBin,"e");
 373
 374       Char_t titleTheta[50];
 375       if(TMath::Abs(i-(fNThetaBins-1))<0.5) snprintf(titleTheta,50,"charge/pt (GeV/c) integrated over #theta");
 376       else  snprintf(titleTheta,50,"charge/pt (GeV/c) for #theta range: %1.3f - %1.3f",fThetaBins[i],fThetaBins[i+1]);
 377
 378       fHistFitTheta[i]->SetTitle(titleTheta);
 379
 380       Double_t invPtMinPos  = 0;
 381       Double_t invPtMinPosErr = 0;
 382       Double_t invPtMinPosImpr  = 0;
 383       Double_t invPtMinPosErrImpr = 0;
 384
 385       if(fDoRebin) fHistFitTheta[i]->Rebin(fRebin);
 386
 387       //start fitting
 388       if(!fFitGaus){
 389          Printf("making polynomial fit in charge/pt in theta bins");
 390          theCan->cd(countPad);
 391          MakeFit(fHistFitTheta[i],fFitMinPos, invPtMinPos,invPtMinPosErr, fExclRange,fRange);
 392          MakeFitBetter(fHistFitTheta[i],fFitMinPosRejP, invPtMinPosImpr, invPtMinPosErrImpr, invPtMinPos,fExclRange,fRange);
 393
 394          //plot projection and fit
 395          fHistFitTheta[i]->DrawCopy();
 396          fFitMinPos->DrawCopy("L,same");
 397          fFitMinPosRejP->DrawCopy("L,same");
 398       }
 399       else{
 400          Printf("making gauss fit in charge/pt in theta bins");
 401          theCan->cd(countPad);
 402          MakeFitInvGauss(fHistFitTheta[i],fFitInvGauss, invPtMinPos, invPtMinPosErr,fExclRange,fRange);
 403          MakeFitInvGaussBetter(fHistFitTheta[i],fFitInvGaussRejP, invPtMinPosImpr, invPtMinPosErrImpr,invPtMinPos,fExclRange,fRange);
 404
 405          //plot projection and fit
 406          fHistFitTheta[i]->DrawCopy();
 407          fFitInvGauss->DrawCopy("L,same");
 408          fFitInvGaussRejP->DrawCopy("L,same");
 409       }
 410       //add objects for analysis folder
 411       aFolderObj->Add(fHistFitTheta[i]);
 412
 413       //store fit parameter
 414       fitParamTheta[i] = invPtMinPosImpr;
 415       errFitParamTheta[i] = invPtMinPosErrImpr;
 416
 417       //count bins and pad number for displaying
 418       binsXTheta[i] = i+1.0;
 419       countPad++;
 420    }
 421
 422
 423    countPad = 1;
 424
 425
 426    // analyse charge/pt in bins of phi
 427
 428    for(Int_t i=0;i<fNPhiBins;i++){
 429
 430       //set name for each projection
 431       TString name = "fit_phi_";
 432       name +=i;
 433       fHistH2InvPtPhi->SetName(name.Data());
 434
 435       //find bins for projection
 436       Int_t  firstBin = fHistH2InvPtPhi->GetYaxis()->FindBin(fPhiBins[i]);
 437       if(i>0) firstBin +=1;
 438       Int_t   lastBin =  fHistH2InvPtPhi->GetYaxis()->FindBin(fPhiBins[i+1]);
 439       if(TMath::Abs(i-(fNPhiBins-1))<0.5){
 440          firstBin = fHistH2InvPtPhi->GetYaxis()->FindBin(fPhiBins[0]);
 441          lastBin =  fHistH2InvPtPhi->GetYaxis()->FindBin(fPhiBins[i]);
 442       }
 443
 444       //make projection
 445       fHistFitPhi[i] =  (TH1D*) fHistH2InvPtPhi->ProjectionX("_px",firstBin,lastBin,"e");
 446
 447       Char_t titlePhi[50];
 448       if(TMath::Abs(i-(fNPhiBins-1))<0.5) snprintf(titlePhi,50,"charge/pt (GeV/c) integrated over #phi");
 449       else  snprintf(titlePhi,50,"charge/pt (GeV/c) for #phi range: %1.3f - %1.3f",fPhiBins[i],fPhiBins[i+1]);
 450
 451       fHistFitPhi[i]->SetTitle(titlePhi);
 452
 453       Double_t invPtMinPos  = 0;
 454       Double_t invPtMinPosErr = 0;
 455       Double_t invPtMinPosImpr  = 0;
 456       Double_t invPtMinPosErrImpr = 0;
 457       if(fDoRebin) fHistFitPhi[i]->Rebin(fRebin);
 458
 459       //start fitting
 460       if(!fFitGaus){
 461          Printf("making polynomial fit in charge/pt in phi bins");
 462          phiCan->cd(countPad);
 463          MakeFit(fHistFitPhi[i],fFitMinPos, invPtMinPos, invPtMinPosErr,fExclRange,fRange);
 464          MakeFitBetter(fHistFitPhi[i],fFitMinPosRejP, invPtMinPosImpr, invPtMinPosErrImpr,invPtMinPos,fExclRange,fRange);
 465
 466          //plot projection and fit
 467          fHistFitPhi[i]->DrawCopy();
 468          fFitMinPos->DrawCopy("L,same");
 469          fFitMinPosRejP->DrawCopy("L,same");
 470
 471       }
 472       else {
 473          Printf("making gauss fit in charge/pt in phi bins");
 474          phiCan->cd(countPad);
 475          MakeFitInvGauss(fHistFitPhi[i],fFitInvGauss, invPtMinPos, invPtMinPosErr, exclRange,fRange);
 476          MakeFitInvGaussBetter(fHistFitPhi[i],fFitInvGaussRejP, invPtMinPosImpr, invPtMinPosErrImpr, invPtMinPos, fExclRange,fRange);
 477
 478          //plot projection and fit
 479          fHistFitPhi[i]->DrawCopy();
 480          fFitInvGauss->DrawCopy("L,same");
 481          fFitInvGaussRejP->DrawCopy("L,same");
 482       }
 483       //add objects for analysis folder
 484       aFolderObj->Add(fHistFitPhi[i]);
 485
 486       //store fit parameter
 487       fitParamPhi[i] = invPtMinPosImpr;
 488       errFitParamPhi[i] = invPtMinPosErrImpr;
 489
 490       //count bins and pad number for displaying
 491       binsXPhi[i] = i+1.0;
 492       countPad++;
 493    }
 494
 495    //initialize graphs for displaying the fit parameter
 496    InitGraphs(binsXTheta,fitParamTheta,errFitParamTheta,binsXPhi,fitParamPhi,errFitParamPhi);
 497
 498    //plot fit values = minimum positions of charge/pt
 499    TCanvas *canFitVal = new TCanvas("canFitVal","min pos histos",800,400);
 500    canFitVal->Divide(2,1);
 501
 502    canFitVal->cd(1);
 503    fGrMinPosTheta->Draw("AP");
 504    canFitVal->cd(2);
 505    fGrMinPosPhi->Draw("AP");
 506
 507    Printf("AliPerfAnalyzeInvPt::StartAnalysis: ***** NOTE: last bin is always fit result of integral over all angle BINS which have been set by user! angle ranges of 2D projections are kept.*****");
 508
 509    //add objects to folder
 510    aFolderObj->Add(fGrMinPosTheta);
 511    aFolderObj->Add(fGrMinPosPhi);
 512    aFolderObj->Add(fHistH2InvPtTheta);
 513    aFolderObj->Add(fHistH2InvPtPhi);
 514
 515
 516 }
 517
 518
 519 //____________________________________________________________________________________________________________________________________________
 520
 521 void AliPerfAnalyzeInvPt::MakeFit(TH1D *hproy, TF1 * fitpb, Double_t &mean, Double_t &errMean, Double_t &excl,Double_t &range)
 522 {
 523    // fit charge/pt and extract minimum position
 524
 525    fitpb->SetRange(-range,range);
 526    fitpb->SetParLimits(1,-0.05,0.05);
 527    fitpb->SetParameter(0,1.0);
 528    fitpb->FixParameter(4,excl);
 529    fitpb->FixParameter(2,0.0);
 530
 531    hproy->Fit(fitpb,"RM");
 532    mean = fitpb->GetParameter(1);
 533    errMean = fitpb->GetParError(1);
 534    if(!mean)  errMean = 0.0;
 535 }
 536 //____________________________________________________________________________________________________________________________________________
 537 void AliPerfAnalyzeInvPt::MakeFitBetter(TH1D *hproy, TF1 * fitpb2, Double_t &mean, Double_t &errMean, Double_t &f, Double_t &excl, Double_t &range)
 538 {
 539    // adjust fit range to minimum position of AliPerfAnalyzeInvPt::MakeFit and fit 1/pt and extract new minimum position
 540
 541    fitpb2->FixParameter(5,excl);
 542    fitpb2->FixParameter(4,f);
 543    // fitpb2->FixParameter(2,0.0);
 544    fitpb2->SetRange(-range+f,range+f);
 545    fitpb2->SetParLimits(1,-0.05,0.05);
 546    fitpb2->SetParameter(0,1.0);
 547    hproy->Fit(fitpb2,"RM");
 548    mean = fitpb2->GetParameter(1);
 549    errMean = fitpb2->GetParError(1);
 550    if(!mean)   errMean = 0.0;
 551
 552 }
 553 //____________________________________________________________________________________________________________________________________________
 554 void AliPerfAnalyzeInvPt::MakeFitInvGauss(TH1D *hproy, TF1 * fitpb, Double_t &mean, Double_t &errMean, Double_t &excl,Double_t &range)
 555 {
 556    // fit charge/pt and extract minimum position
 557
 558    fitpb->FixParameter(6,excl);
 559    fitpb->SetRange(-range,range);
 560    fitpb->SetParameter(0,-1.0);
 561    fitpb->SetParameter(2,1.0);
 562    fitpb->SetParameter(3,25000.0);
 563    fitpb->SetParameter(4,-1.0);
 564    fitpb->SetParLimits(1,-0.02,0.02);
 565    hproy->Fit(fitpb,"RM");
 566    mean = fitpb->GetParameter(1);
 567    errMean = fitpb->GetParError(1);
 568    if(!mean)   errMean = 0.0;
 569
 570 }
 571 //____________________________________________________________________________________________________________________________________________
 572 void AliPerfAnalyzeInvPt::MakeFitInvGaussBetter(TH1D *hproy, TF1 * fitpb2, Double_t &mean, Double_t &errMean, Double_t &f, Double_t &excl, Double_t &range)
 573 {
 574    // adjust fit range to minimum position of AliPerfAnalyzeInvPt::MakeFitInvGauss and fit 1/pt and extract new minimum position
 575
 576    fitpb2->FixParameter(7,excl);
 577    fitpb2->FixParameter(6,f);
 578    fitpb2->SetRange(-range+f,range+f);
 579    fitpb2->SetParameter(0,-1.0);
 580    fitpb2->SetParameter(2,1.0);
 581    fitpb2->SetParameter(3,25000.0);
 582    fitpb2->SetParameter(4,-1.0);
 583    fitpb2->SetParLimits(1,-0.02,0.02);
 584    hproy->Fit(fitpb2,"RM");
 585    mean = fitpb2->GetParameter(1);
 586    errMean = fitpb2->GetParError(1);
 587    if(!mean)   errMean = 0.0;
 588
 589 }
 590
 591
 592 //____________________________________________________________________________________________________________________________________________
 593
 594 // set variables for analysis
 595 void AliPerfAnalyzeInvPt::SetProjBinsPhi(const Double_t *phiBinArray,const Int_t nphBins){
 596
 597    // set theta bins for Analyse()
 598    //set phi bins as array and set number of this array which is equal to number of bins analysed
 599    //the last analysed bin will always be the projection from first to last bin in the array
 600    if(nphBins){
 601       fNPhiBins = nphBins;
 602
 603       for(Int_t k = 0;k<fNPhiBins;k++){
 604          fPhiBins[k] = phiBinArray[k];
 605       }
 606       Printf("AliPerfAnalyzeInvPt::SetProjBinsPhi:number of bins in phi set to %i",fNPhiBins);
 607    }
 608    else Printf("AliPerfAnalyzeInvPt::SetProjBinsPhi: number of bins in theta NOT set. Default values are used.");
 609 }
 610 //____________________________________________________________________________________________________________________________________________
 611 void AliPerfAnalyzeInvPt::SetProjBinsTheta(const Double_t *thetaBinArray, const Int_t nthBins){
 612    // set theta bins for Analyse()
 613    //set theta bins as array and set number of this array which is equal to number of bins analysed
 614    //the last analysed bin will always be the projection from first to last bin in the array
 615    if(nthBins){
 616       fNThetaBins = nthBins;
 617       for(Int_t k = 0;k<fNThetaBins;k++){
 618          fThetaBins[k] = thetaBinArray[k];
 619       }
 620       Printf("AliPerfAnalyzeInvPt::SetProjBinsTheta: number of bins in theta set to %i",fNThetaBins);
 621    }
 622    else Printf("AliPerfAnalyzeInvPt::SetProjBinsTheta: number of bins in theta NOT set. Default values are used.");
 623 }
 624 //____________________________________________________________________________________________________________________________________________
 625 void AliPerfAnalyzeInvPt::SetMakeFitOption(const Bool_t setGausFit, const Double_t exclusionR,const Double_t fitR ){
 626
 627    //set the fit options:
 628    //for usage of gaussian function instead of polynomial (default) set setGausFit=kTRUE
 629    //set the range of rejection of points around 0 via exclusionR
 630    //set the fit range around 0 with fitR
 631    if(fitR){
 632       fFitGaus = setGausFit;
 633       fExclRange  = exclusionR;
 634       fRange = fitR;
 635
 636       if(fFitGaus) Printf("AliPerfAnalyzeInvPt::SetMakeFitOption: set MakeGausFit with fit range %2.3f and exclusion range in fabs(1/pt): %2.3f",fRange,fExclRange);
 637       else  Printf("AliPerfAnalyzeInvPt::SetMakeFitOption: set standard polynomial fit with fit range %2.3f and exclusion range in fabs(1/pt): %2.3f",fRange,fExclRange);
 638    }
 639    else Printf(" AliPerfAnalyzeInvPt::SetMakeFitOption: no user input. Set standard polynomial fit with fit range 1.0 and exclusion range 0.13.");
 640 }