Flag for In/Out of Plane Analysis
authorzconesa <zconesa@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 13 Jul 2012 19:39:46 +0000 (19:39 +0000)
committerzconesa <zconesa@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 13 Jul 2012 19:39:46 +0000 (19:39 +0000)
PWGHF/vertexingHF/AliHFPtSpectrum.cxx
PWGHF/vertexingHF/AliHFPtSpectrum.h
PWGHF/vertexingHF/macros/HFPtSpectrum.C

index bf84b74..c8f4ba9 100644 (file)
@@ -100,6 +100,7 @@ AliHFPtSpectrum::AliHFPtSpectrum(const char* name, const char* title, Int_t opti
   fPbPbElossHypothesis(kFALSE),
   fIsStatUncEff(kTRUE),
   fParticleAntiParticle(2),
+  fIsEventPlane(kFALSE),
   fhStatUncEffcSigma(NULL),
   fhStatUncEffbSigma(NULL),
   fhStatUncEffcFD(NULL),
@@ -162,6 +163,7 @@ AliHFPtSpectrum::AliHFPtSpectrum(const AliHFPtSpectrum &rhs):
   fPbPbElossHypothesis(rhs.fPbPbElossHypothesis),
   fIsStatUncEff(rhs.fIsStatUncEff),
   fParticleAntiParticle(rhs.fParticleAntiParticle),
+  fIsEventPlane(rhs.fIsEventPlane),
   fhStatUncEffcSigma(NULL),
   fhStatUncEffbSigma(NULL),
   fhStatUncEffcFD(NULL),
@@ -227,6 +229,7 @@ AliHFPtSpectrum &AliHFPtSpectrum::operator=(const AliHFPtSpectrum &source){
   fPbPbElossHypothesis = source.fPbPbElossHypothesis;
   fIsStatUncEff = source.fIsStatUncEff;
   fParticleAntiParticle = source.fParticleAntiParticle;
+  fIsEventPlane = source.fIsEventPlane;
   
   for(Int_t i=0; i<2; i++){
     fLuminosity[i] = source.fLuminosity[i];
@@ -1431,7 +1434,8 @@ void AliHFPtSpectrum::CalculateFeedDownCorrectedSpectrumNb(Double_t deltaY, Doub
     //
     Double_t frac = 1.0, errfrac =0.;
     if(fPbPbElossHypothesis) {
-      frac = fTab[0]*fNevts; 
+      frac = fTab[0]*fNevts;
+      if(fIsEventPlane) frac = frac/2.0;
       errfrac = frac * TMath::Sqrt( (fTab[1]/fTab[0])*(fTab[1]/fTab[0]) + (1/fNevts) );
     } else {
       frac = fLuminosity[0]; 
@@ -1440,7 +1444,7 @@ void AliHFPtSpectrum::CalculateFeedDownCorrectedSpectrumNb(Double_t deltaY, Doub
     
     value = ( fhRECpt->GetBinContent(ibin)>0. && fhRECpt->GetBinContent(ibin)!=0. && 
              fhFeedDownMCpt->GetBinContent(ibin)>0. && fhFeedDownEffpt->GetBinContent(ibin)>0. ) ?
-      fhRECpt->GetBinContent(ibin) - frac*(deltaY*branchingRatioBintoFinalDecay*fParticleAntiParticle*fTrigEfficiency[0]*fhFeedDownEffpt->GetBinContent(ibin)*fhFeedDownMCpt->GetBinContent(ibin) * fhRECpt->GetBinWidth(ibin) ) 
+      fhRECpt->GetBinContent(ibin) - frac*(deltaY*branchingRatioBintoFinalDecay*fParticleAntiParticle*fTrigEfficiency[0]*fhFeedDownEffpt->GetBinContent(ibin)*fhFeedDownMCpt->GetBinContent(ibin) * fhRECpt->GetBinWidth(ibin) )
       : 0. ;
     value /= fhRECpt->GetBinWidth(ibin);
     if (value<0.) value =0.;
@@ -1457,6 +1461,8 @@ void AliHFPtSpectrum::CalculateFeedDownCorrectedSpectrumNb(Double_t deltaY, Doub
       1 - (frac*deltaY*branchingRatioBintoFinalDecay*fParticleAntiParticle*fTrigEfficiency[0]*fhFeedDownEffpt->GetBinContent(ibin)*fhFeedDownMCpt->GetBinContent(ibin) * fhRECpt->GetBinWidth(ibin) ) / fhRECpt->GetBinContent(ibin)  : 0. ;
     if (correction<0.) correction = 0.;
 
+    //    cout << " pt="<< fhRECpt->GetBinCenter(ibin) << " rec="<< fhRECpt->GetBinContent(ibin) << ", corr="<<correction<<" = 1 - "<<  (frac*deltaY*branchingRatioBintoFinalDecay*fParticleAntiParticle*fTrigEfficiency[0]*fhFeedDownEffpt->GetBinContent(ibin)*fhFeedDownMCpt->GetBinContent(ibin) * fhRECpt->GetBinWidth(ibin) ) / fhRECpt->GetBinContent(ibin) << endl;
+
     // Systematic uncertainties
     //     (syst but feed-down)  delta_physics = sqrt ( (delta_reco_syst)^2 )  / bin-width
     //         (feed-down syst)  delta_physics = sqrt ( (k*delta_lumi/lumi)^2 + (k*delta_eff_trig/eff_trig)^2 
@@ -1512,13 +1518,14 @@ void AliHFPtSpectrum::CalculateFeedDownCorrectedSpectrumNb(Double_t deltaY, Doub
     //
     // Estimate how the result varies vs charm/beauty Eloss hypothesis
     //
-    if ( correction>0.0001 && fPbPbElossHypothesis){
+    if ( correction>1.0e-16 && fPbPbElossHypothesis){
       // Loop over the Eloss hypothesis
       //      Int_t rbin=0;
       for (Float_t rval=0.0025; rval<4.0; rval+=0.005){
        // correction =  [ 1  - (Tab *Nevt * delta_y * BR_b * eff_trig * eff_b * Nb_th *binwidth )* (rval) /reco ] 
        Double_t fcRcbvalue = 1 - (fTab[0]*fNevts*deltaY*branchingRatioBintoFinalDecay*fParticleAntiParticle*fTrigEfficiency[0]*fhFeedDownEffpt->GetBinContent(ibin)*fhFeedDownMCpt->GetBinContent(ibin)*fhRECpt->GetBinWidth(ibin) * rval ) / fhRECpt->GetBinContent(ibin) ;
-       if(fcRcbvalue<0.) fcRcbvalue=0.;
+       //      cout << " rval="<<rval <<" , fc="<<fcRcbvalue<<"    ";
+       if(fcRcbvalue<1.e-16) fcRcbvalue=0.;
        fhFcRcb->Fill( fhRECpt->GetBinCenter(ibin) , rval, fcRcbvalue );
        //    physics = reco * fcRcb / bin-width
        Double_t Rcbvalue = (fhRECpt->GetBinContent(ibin) && fcRcbvalue) ? 
index 4c1c591..7601de8 100644 (file)
@@ -78,6 +78,7 @@ class AliHFPtSpectrum: public TNamed
     if (flag) { fParticleAntiParticle = 2; AliInfo(" Setting for particle + anti-particle yields"); }
     else { fParticleAntiParticle = 1; AliInfo(" Setting for only (anti)particle yields, not the sum of both"); }
   }
+  void SetIsEventPlaneAnalysis(Bool_t flag){ fIsEventPlane = flag; }
   //
   void SetfIsStatUncEff(Bool_t flag){ fIsStatUncEff = flag; }
   // Set if the calculation has to consider Ratio(c/b eloss) hypothesis 
@@ -293,6 +294,7 @@ class AliHFPtSpectrum: public TNamed
   Bool_t fPbPbElossHypothesis;      // flag: whether to do estimates vs Ratio(c/b eloss) hypothesis
   Bool_t fIsStatUncEff;             // flag : consider (1) or not (0) the stat unc on the efficiencies
   Int_t fParticleAntiParticle;      // 1: only one sign, 2: yield is for particle+anti-particle
+  Bool_t fIsEventPlane;             // flag : when the analysis is done for In/Out of plane, divide the B-prediction by two
 
   //
   TH1D *fhStatUncEffcSigma;       // Uncertainty on the cross-section due to the prompt efficiency statistical uncertainty
@@ -300,7 +302,7 @@ class AliHFPtSpectrum: public TNamed
   TH1D *fhStatUncEffcFD;          // Uncertainty on the feed-down correction due to the prompt efficiency statistical uncertainty
   TH1D *fhStatUncEffbFD;          // Uncertainty on the feed-down correction due to the feed-down efficiency statistical uncertainty
 
-  ClassDef(AliHFPtSpectrum,3) // Class for Heavy Flavor spectra corrections
+  ClassDef(AliHFPtSpectrum,4) // Class for Heavy Flavor spectra corrections
 };
 
 #endif
index c56db05..d2b5f1a 100644 (file)
@@ -46,7 +46,7 @@ void HFPtSpectrum ( const char *mcfilename="FeedDownCorrectionMC.root",
                    const char *recofilename="Reconstructed.root", const char *recohistoname="hRawSpectrumD0",
                    const char *outfilename="HFPtSpectrum.root",
                    Int_t fdMethod=kNb, Double_t nevents=1.0, Double_t sigma=1.0, // sigma[nb]
-                   Bool_t isParticlePlusAntiParticleYield=true, Int_t cc=kpp7, Bool_t PbPbEloss=false ) {
+                   Bool_t isParticlePlusAntiParticleYield=true, Int_t cc=kpp7, Bool_t PbPbEloss=false, Bool_t kRaavsEP=kFALSE) {
 
 
   gROOT->Macro("$ALICE_ROOT/PWGHF/vertexingHF/macros/LoadLibraries.C");
@@ -257,6 +257,7 @@ void HFPtSpectrum ( const char *mcfilename="FeedDownCorrectionMC.root",
   spectra->SetLuminosity(lumi,lumiUnc);
   Double_t effTrig = 1.0;
   spectra->SetTriggerEfficiency(effTrig,0.);
+  if(kRaavsEP) spectra->SetIsEventPlaneAnalysis(kTRUE);
 
   // Set the global uncertainties on the efficiencies (in percent)
   Double_t globalEffUnc = 0.15; 
@@ -286,10 +287,7 @@ void HFPtSpectrum ( const char *mcfilename="FeedDownCorrectionMC.root",
     systematics->SetIsLowEnergy(true);
   } else if( cc!=kpp7 )  {
     systematics->SetCollisionType(1);
-    if ( cc == k07half ) { 
-      cout<<endl<<" Beware, you're using the systematics of the 0-10% of 2010 data !! FIX ME !!"<<endl<<endl; 
-      systematics->SetCentrality("010"); 
-    }
+    if ( cc == k07half ) systematics->SetCentrality("07half");
     else if ( cc == k010 )  systematics->SetCentrality("010");
     else if ( cc == k1020 )  systematics->SetCentrality("1020");
     else if ( cc == k020 )  systematics->SetCentrality("020");