PWGPP/EVCHAR/makeCentralityFit.C

   1 void makeCentralityFit(const char * run="188362",const char * system="TKL", int Rebin=1,int Nevt=1e6)
   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/GlauberNtuplePA/GlauberMC_pPb_ntuple_sigma70_mind4_r662_a546_Rpro4.root";
  15   TString finname = Form("/home/atoia/analysis/data2012/multRef/EventStat_temp_%s.root",run);
  16   TString foutname = Form("/home/atoia/analysis/data2012/fit/%s_fit_%s.root",system,run);
  17   TString foutnameGlau = Form("/home/atoia/analysis/data2012/fit/%s_ntuple_%s.root",system,run);
  18   const char *histname=Form("fHOutMult%s",system);
  19
  20   AliCentralityGlauberFit *mPM = new AliCentralityGlauberFit(finnameGlau);
  21   mPM->SetInputFile(finname);
  22   mPM->SetInputNtuple(finnameGlau);
  23   mPM->SetOutputFile(foutname);
  24   mPM->SetOutputNtuple(foutnameGlau);
  25   mPM->AddHisto(histname);
  26
  27   mPM->SetRebin(Rebin);
  28   mPM->SetNevents(Nevt);
  29   mPM->SetAncestorMode(1); // 1: Npart**alpha, 2: alpha*Npart + (1-alpha)*Ncoll
  30   mPM->SetFastFitMode(0);  // 1:NBD, 2:Gauss
  31   mPM->UseChi2(kTRUE);     // If TRUE minimize Chi2
  32   mPM->UseAverage(kFALSE); // If TRUE use Average
  33   mPM->SetNtrials(1);
  34
  35   // ----------range to fit--------------
  36   if (strncmp (system,"V0A",1) == 0) {
  37     mPM->SetRangeToFit(15., 600.);   // range to fit
  38     mPM->SetRangeToScale(15.); // range to scale
  39     mPM->SetGlauberParam(1,12.8,13, 20,0.5,2.5, 1,1,1); // fit parameters
  40   }
  41   else if (strncmp (system,"TKL",1) == 0) {
  42     mPM->SetRangeToFit(15., 200.);   // range to fit
  43     mPM->SetRangeToScale(15.); // range to scale
  44     mPM->SetGlauberParam(1,4.9,5.5, 1,0.61,0.62, 1,1,1); // fit parameters
  45   }
  46   else if (strncmp (system,"CL1",1) == 0) {
  47     mPM->SetRangeToFit(40., 400.);
  48     mPM->SetRangeToScale(40.);
  49     mPM->SetGlauberParam(1,7.9,8, 1,0.43,0.46, 1,1,1);
  50   }
  51   else if (strncmp (system,"TRK",1) == 0) {
  52     mPM->SetRangeToFit(10., 2600.);
  53     mPM->SetRangeToScale(10.);
  54     mPM->SetGlauberParam(1,3.9,4.2, 1,1.3,2.5, 1,0.801,0.81);
  55   }
  56
  57   mPM->MakeFits();
  58
  59   // ----------for Minuit--------------
  60   //mPM->MakeFitsMinuitNBD(1,8,0.5);          // initial parameters
  61
  62
  63   TFile * f = new TFile (foutname);
  64   TH1 * hd = (TH1*) gDirectory->Get(Form("fHOutMult%s",system));
  65   TH1 * hg = (TH1*) gDirectory->Get(Form("fHOutMult%s_GLAU",system));
  66   hg->SetLineColor(kRed);
  67   hd->Draw("e");
  68   hg->Draw("same");
  69
  70
  71
  72 }
  73