PWGPP/EVCHAR/makeCentralityFit.C

   1 void makeCentralityFit(const char * run="167693",const char * system="V0M", int Rebin=100,int Nevt=1e5, int Ncell=0)
   2 {
   3  //load libraries
   4   gSystem->SetBuildDir("/tmp/");
   5   gSystem->Load("libCore.so");
   6   gSystem->Load("libTree.so");
   7   gSystem->Load("libGeom.so");
   8   gSystem->Load("libVMC.so");
   9   gSystem->Load("libPhysics.so");
  10   gSystem->Load("libSTEERBase.so");
  11   gROOT->ProcessLine(".include $ALICE_ROOT/include");
  12   gROOT->LoadMacro("AliCentralityGlauberFit.cxx+");
  13
  14   const char *finnameGlau ="/home/atoia/GlauberNtuple/GlauberMC_PbPb_ntuple_sigma64_mind4_r662_a546.root";
  15
  16
  17   TString finname = Form("/home/atoia/analysis/data2011/multRef/EventStat_temp_%s.root",run);
  18   TString foutname = Form("/home/atoia/analysis/data2011/fit/%s_fitTEST_%s.root",system,run);
  19   TString foutnameGlau = Form("/home/atoia/analysis/data2011/fit/%s_ntupleTEST_%s.root",system,run);
  20   const char *histname=Form("fHOutMult%s",system);
  21
  22
  23   /*
  24   TString finname = Form("/home/atoia/analysis/EPVzero/VZEROEquaFactorStat.root");
  25   TString foutname = Form("/home/atoia/analysis/data2011/fit/%sCell%d_fit_%s.root",system,Ncell,run);
  26   TString foutnameGlau = Form("/home/atoia/analysis/data2011/fit/%sCell%d_ntuple_%s.root",system,Ncell,run);
  27   const char *histname=Form("fMultCell_%d",Ncell);
  28   */
  29   AliCentralityGlauberFit *mPM = new AliCentralityGlauberFit(finnameGlau);
  30   mPM->SetInputFile(finname);
  31   mPM->SetInputNtuple(finnameGlau);
  32   mPM->SetOutputFile(foutname);
  33   mPM->SetOutputNtuple(foutnameGlau);
  34   mPM->AddHisto(histname);
  35
  36   mPM->SetRebin(Rebin);
  37   mPM->SetNevents(Nevt);
  38   mPM->SetAncestorMode(2); // 1: Npart**alpha, 2: alpha*Npart + (1-alpha)*Ncoll
  39   mPM->SetFastFitMode(0);  // 1:NBD, 2:Gauss
  40   mPM->UseChi2(kTRUE);     // If TRUE minimize Chi2
  41   mPM->UseAverage(kFALSE); // If TRUE use Average
  42   mPM->SetNtrials(1);
  43
  44   // ----------range to fit--------------
  45   if (strncmp (system,"V0M",1) == 0) {
  46     mPM->SetRangeToFit(100., 22000.);   // range to fit
  47     mPM->SetRangeToScale(100.); // range to scale
  48     mPM->SetGlauberParam(1,28.7,29., 1,1.601,1.5, 1,0.80,0.805); // fit parameters
  49   }
  50   else if (strncmp (system,"CL1",1) == 0) {
  51     mPM->SetRangeToFit(40., 5400.);
  52     mPM->SetRangeToScale(40.);
  53     mPM->SetGlauberParam(1,7.13,7.3, 1,1.217,1.9, 1,0.802,0.815);
  54   }
  55   else if (strncmp (system,"TRK",1) == 0) {
  56     mPM->SetRangeToFit(10., 2600.);
  57     mPM->SetRangeToScale(10.);
  58     mPM->SetGlauberParam(1,3.9,4.2, 1,1.3,2.5, 1,0.801,0.81);
  59   }
  60
  61   mPM->MakeFits();
  62
  63   // ----------for Minuit--------------
  64   //mPM->MakeFitsMinuitNBD(0.8,28.,1.29);          // initial parameters
  65
  66
  67   TFile * f = new TFile (foutname);
  68   TH1 * hd = (TH1*) gDirectory->Get(Form("fHOutMult%s",system));
  69   TH1 * hg = (TH1*) gDirectory->Get(Form("fHOutMult%s_GLAU",system));
  70   //TH1 * hd = (TH1*) gDirectory->Get(Form("fMultCell_%d",Ncell));
  71   //TH1 * hg = (TH1*) gDirectory->Get(Form("fMultCell_%d_GLAU",Ncell));
  72   hg->SetLineColor(kRed);
  73   hd->Draw("e");
  74   hg->Draw("same");
  75
  76
  77
  78 }
  79