]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWGHF/vertexingHF/macros/HFPtSpectrumRaa.C
Merge branch 'feature-movesplit'
[u/mrichter/AliRoot.git] / PWGHF / vertexingHF / macros / HFPtSpectrumRaa.C
index 18d2dd7fd05e6e75e3dcf4ad7f5fb3d24faaf578..7cb960cc208379353d7ad22d26f837a3636c0fd9 100644 (file)
@@ -50,7 +50,8 @@
 //
 ///////////////////////////////////////////////////////////////////////////////
 
-enum centrality{ kpp, k07half, kpPb0100, k010, k1020, k020, k2040, k2030, k3040, k4050, k3050, k5060, k4060, k6080, k4080, k5080, k80100 };
+enum centrality{ kpp, k07half, kpPb0100, k010, k1020, k020, k2040, k2030, k3040, k4050, k3050, k5060, k4060, k6080, k4080, k5080, k80100, kpPb020, kpPb2040, kpPb4060, kpPb60100 };
+enum centestimator{ kV0M, kV0A, kZNA };
 enum energy{ k276, k5dot023, k55 };
 enum BFDSubtrMethod { kfc, kNb };
 enum RaavsEP {kPhiIntegrated, kInPlane, kOutOfPlane};
@@ -102,6 +103,7 @@ void HFPtSpectrumRaa(const char *ppfile="HFPtSpectrum_D0Kpi_method2_rebinnedth_2
                     Int_t fdMethod = kNb, Int_t cc=kpp, Int_t Energy=k276,
                     Double_t MinHypo=1./3., Double_t MaxHypo=3.0, Double_t MaxRb=6.0,
                     Bool_t isRbHypo=false, Double_t CentralHypo = 1.0,
+                    Int_t ccestimator = kV0M,
                     Bool_t isUseTaaForRaa=true, const char *shadRbcFile="", Int_t nSigmaShad=3.0,
                     Int_t isRaavsEP=kPhiIntegrated, Bool_t isScaledAndExtrapRef=kFALSE)
 {
@@ -121,44 +123,68 @@ void HFPtSpectrumRaa(const char *ppfile="HFPtSpectrum_D0Kpi_method2_rebinnedth_2
   }
   // Values from Alberica's twiki:
   //   https://twiki.cern.ch/twiki/bin/viewauth/ALICE/CentStudies
-  if ( cc == k07half ) {
-    Tab = 24.81; TabSyst = 0.8037;
-  } else if ( cc == k010 ) {
-    Tab = 23.48; TabSyst = 0.97;
-  } else if ( cc == k1020 ) {
-    Tab = 14.4318; TabSyst = 0.5733;
-  } else if ( cc == k020 ) {
-    Tab = 18.93; TabSyst = 0.74;
-  } else if ( cc == k2040 ) {
-    Tab = 6.86; TabSyst = 0.28;
-  } else if ( cc == k2030 ) {
-    Tab = 8.73769; TabSyst = 0.370219;
-  } else if ( cc == k3040 ) {
-    Tab = 5.02755; TabSyst = 0.22099;
-  } else if ( cc == k4050 ) {
-    Tab = 2.68327; TabSyst = 0.137073;
-  } else if ( cc == k3050 ) {
-    Tab = 3.87011; TabSyst = 0.183847;
-  } else if ( cc == k4060 ) {
-    Tab = 2.00;  TabSyst= 0.11;
-  } else if ( cc == k4080 ) {
-    Tab = 1.20451; TabSyst = 0.071843;
-  } else if ( cc == k5060 ) {
-    Tab = 1.32884; TabSyst = 0.0929536;
-  } else if ( cc == k6080 ) {
-    Tab = 0.419; TabSyst = 0.033;
-  } else if ( cc == k5080 ) {
-    Tab = 0.719; TabSyst = 0.054;
-  } else if ( cc == k80100 ){
-    Tab = 0.0690; TabSyst = 0.0062;
+  if( ccestimator == kV0M ) {
+    if ( cc == k07half ) {
+      Tab = 24.81; TabSyst = 0.8037;
+    } else if ( cc == k010 ) {
+      Tab = 23.48; TabSyst = 0.97;
+    } else if ( cc == k1020 ) {
+      Tab = 14.4318; TabSyst = 0.5733;
+    } else if ( cc == k020 ) {
+      Tab = 18.93; TabSyst = 0.74;
+    } else if ( cc == k2040 ) {
+      Tab = 6.86; TabSyst = 0.28;
+    } else if ( cc == k2030 ) {
+      Tab = 8.73769; TabSyst = 0.370219;
+    } else if ( cc == k3040 ) {
+      Tab = 5.02755; TabSyst = 0.22099;
+    } else if ( cc == k4050 ) {
+      Tab = 2.68327; TabSyst = 0.137073;
+    } else if ( cc == k3050 ) {
+      Tab = 3.87011; TabSyst = 0.183847;
+    } else if ( cc == k4060 ) {
+      Tab = 2.00;  TabSyst= 0.11;
+    } else if ( cc == k4080 ) {
+      Tab = 1.20451; TabSyst = 0.071843;
+    } else if ( cc == k5060 ) {
+      Tab = 1.32884; TabSyst = 0.0929536;
+    } else if ( cc == k6080 ) {
+      Tab = 0.419; TabSyst = 0.033;
+    } else if ( cc == k5080 ) {
+      Tab = 0.719; TabSyst = 0.054;
+    } else if ( cc == k80100 ){
+      Tab = 0.0690; TabSyst = 0.0062;
+    }
   }
 
   // pPb Glauber (A. Toia)
   // https://twiki.cern.ch/twiki/bin/viewauth/ALICE/PACentStudies#Glauber_Calculations_with_sigma
-  else if( cc == kpPb0100 ){
+  if( cc == kpPb0100 ){
     Tab = 0.098334; TabSyst = 0.0070679;
     A=207.2; B=1.;
   }
+  else if( ccestimator == kV0A ){
+    if ( cc == kpPb020 ) {
+      Tab = 0.183; TabSyst = 0.006245;
+    } else if ( cc == kpPb2040 ) {
+      Tab = 0.134; TabSyst = 0.004899;
+    } else if ( cc == kpPb4060 ) {
+      Tab = 0.092; TabSyst = 0.004796;
+    } else if ( cc == kpPb60100 ) {
+      Tab = 0.041; TabSyst = 0.008832;
+    }
+  }
+  else if( ccestimator == kZNA ){
+    if ( cc == kpPb020 ) {
+      Tab = 0.164; TabSyst = 0.010724;
+    } else if ( cc == kpPb2040 ) {
+      Tab = 0.137; TabSyst = 0.005099;
+    } else if ( cc == kpPb4060 ) {
+      Tab = 0.1011; TabSyst = 0.006;
+    } else if ( cc == kpPb60100 ) {
+      Tab = 0.0459; TabSyst = 0.003162;
+    }
+  }
 
   //
   // Reading the pp file 
@@ -190,8 +216,8 @@ void HFPtSpectrumRaa(const char *ppfile="HFPtSpectrum_D0Kpi_method2_rebinnedth_2
   // 
   TFile * ABf = new TFile(ABfile,"read");
   TH1D *hSigmaAB = (TH1D*)ABf->Get("histoSigmaCorr");
-  TH2D *hSigmaABRcb = (TH2D*)ABf->Get("histoSigmaCorrRcb");
-  TGraphAsymmErrors * gSigmaABSyst = (TGraphAsymmErrors*)ABf->Get("gSigmaCorr");
+  //  TH2D *hSigmaABRcb = (TH2D*)ABf->Get("histoSigmaCorrRcb");
+  //  TGraphAsymmErrors * gSigmaABSyst = (TGraphAsymmErrors*)ABf->Get("gSigmaCorr");
   TGraphAsymmErrors * gSigmaABSystFeedDown = (TGraphAsymmErrors*)ABf->Get("gSigmaCorrConservative");
   TNtuple * nSigmaAB = (TNtuple*)ABf->Get("fnSigma");
   //
@@ -239,8 +265,24 @@ void HFPtSpectrumRaa(const char *ppfile="HFPtSpectrum_D0Kpi_method2_rebinnedth_2
     else if (isRaavsEP == kOutOfPlane) systematicsAB->SetCentrality("3050OutOfPlane");
   }
   //
-  else if ( cc == kpPb0100 ){ 
-    systematicsAB->SetCollisionType(2); 
+  else if ( cc == kpPb0100 || cc == kpPb020 || cc == kpPb2040 || cc == kpPb4060 || cc == kpPb60100 ) {
+    systematicsAB->SetCollisionType(2);
+    if(ccestimator==kV0A) {
+      if(cc == kpPb020) systematicsAB->SetCentrality("020V0A");
+      else if(cc == kpPb2040) systematicsAB->SetCentrality("2040V0A");
+      else if(cc == kpPb4060) systematicsAB->SetCentrality("4060V0A");
+      else if(cc == kpPb60100) systematicsAB->SetCentrality("60100V0A");
+    } else if (ccestimator==kZNA) {
+      if(cc == kpPb020) systematicsAB->SetCentrality("020ZNA");
+      else if(cc == kpPb2040) systematicsAB->SetCentrality("2040ZNA");
+      else if(cc == kpPb4060) systematicsAB->SetCentrality("4060ZNA");
+      else if(cc == kpPb60100) systematicsAB->SetCentrality("60100ZNA");
+    } else {
+      if(!(cc == kpPb0100)) {
+       cout <<" Error on the pPb options"<<endl;
+       return;
+      }
+    }
   }
   else { 
     cout << " Systematics not yet implemented " << endl;
@@ -296,6 +338,10 @@ void HFPtSpectrumRaa(const char *ppfile="HFPtSpectrum_D0Kpi_method2_rebinnedth_2
     TH1D *hRbcShadCentral = (TH1D*)fshad->Get("hDfromBoverDfromc_L0");
     TH1D *hRbcShadMin = (TH1D*)fshad->Get("hDfromBoverDfromc_L0");
     TH1D *hRbcShadMax = (TH1D*)fshad->Get("hDfromBoverDfromc_L1");
+    if(!hRbcShadCentral || !hRbcShadMin || !hRbcShadMax) {
+      cout<< endl <<">> Shadowing input histograms are not ok !! "<<endl<<endl;
+      return;
+    }
     //       nSigmaShad
     //    nSigmaShad
     for(Int_t i=1; i<=nbins; i++) {
@@ -322,7 +368,7 @@ void HFPtSpectrumRaa(const char *ppfile="HFPtSpectrum_D0Kpi_method2_rebinnedth_2
   //
   TH2D * hRABvsRcb = new TH2D("hRABvsRcb"," R_{AB}(c) vs Rcb Eloss hypothesis; p_{T} [GeV/c] ;  R_{AB}(c) ; Rcb Eloss hypothesis  ",nbins,limits,800,0.,4.);
   TH2D * hRABvsRb = new TH2D("hRABvsRb"," R_{AB}(c) vs Rb Eloss hypothesis; p_{T} [GeV/c] ;  R_{AB}(c) ; Rb Eloss hypothesis ",nbins,limits,800,0.,4.);
-  TH2D * hRABBeautyvsRCharm = new TH2D("hRABBeautyvsRCharm"," R_{AB}(c) vs Rb Eloss hypothesis; p_{T} [GeV/c] ;  R_{AB}(b) ;  R_{AB}(c) ",nbins,limits,800,0.,4.);
+  //  TH2D * hRABBeautyvsRCharm = new TH2D("hRABBeautyvsRCharm"," R_{AB}(c) vs Rb Eloss hypothesis; p_{T} [GeV/c] ;  R_{AB}(b) ;  R_{AB}(c) ",nbins,limits,800,0.,4.);
   Int_t nbinsHypo=800;//200;
   Double_t *limitsHypo = new Double_t[nbinsHypo+1];
   for(Int_t i=1; i<=nbinsHypo+1; i++) limitsHypo[i-1]= i*4./800.;
@@ -609,8 +655,6 @@ void HFPtSpectrumRaa(const char *ppfile="HFPtSpectrum_D0Kpi_method2_rebinnedth_2
     Double_t minFdSyst = 0., maxFdSyst = 0.;
     if ( ElossHypo == ElossCentral[ hABbin ] ) {
 
-      cout <<endl<<">>>>>>>>>>> tete >>>>>>>>>>>>>>"<<endl<<endl;
-
       //
       // Data stat uncertainty
       //
@@ -770,9 +814,11 @@ void HFPtSpectrumRaa(const char *ppfile="HFPtSpectrum_D0Kpi_method2_rebinnedth_2
     //
     // Filling Eloss scenarii information
     //
+    //    trick in case not fine enough Rb hypothesis to cope with the min/max range
+    //    if( RaaCharm>0 && ( (ElossHypo >= MinHypo && ElossHypo <=MaxHypo) || ElossHypo == ElossCentral[ hABbin ] ) && RaaBeauty<=MaxRb ) {
+    //      by default better not use it, to monitor when this happens (could affect results)
     if( RaaCharm>0 && ElossHypo >= MinHypo && ElossHypo <=MaxHypo && RaaBeauty<=MaxRb ) {
 
-      cout <<endl<<">>>>>>>>>>> tetuu >>>>>>>>>>>>>>"<<endl<<endl;
       Double_t Ehigh =  ElossMax[ hABbin ] ;
       Double_t Elow =  ElossMin[ hABbin ] ;
       if ( RaaCharm > Ehigh ) ElossMax[ hABbin ] = RaaCharm ;
@@ -1045,7 +1091,7 @@ void HFPtSpectrumRaa(const char *ppfile="HFPtSpectrum_D0Kpi_method2_rebinnedth_2
   legrcb->Draw();
   TLatex* tc;
   TString system = "Pb-Pb   #sqrt{s_{NN}}=2.76 TeV";
-  if(cc==kpPb0100) system = "p-Pb   #sqrt{s_{NN}}=5.023 TeV";
+  if( cc==kpPb0100 || cc==kpPb020 || cc==kpPb2040 || cc==kpPb4060 || cc==kpPb60100 ) system = "p-Pb   #sqrt{s_{NN}}=5.023 TeV";
   if(decay==1) tc =new TLatex(0.18,0.82,Form("D^{0},  %s ",system.Data()));
   else if(decay==2) tc =new TLatex(0.18,0.82,Form("D^{+},  %s ",system.Data()));
   else if(decay==3) tc =new TLatex(0.18,0.82,Form("D^{*+},  %s ",system.Data()));