]> git.uio.no Git - u/mrichter/AliRoot.git/blob - PWG3/vertexingHF/AliHFPtSpectrum.h
37ed7302846254116fe910beb2d5eaa0b7384bcb
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliHFPtSpectrum.h
1 #ifndef ALIHFPTSPECTRUM_H
2 #define ALIHFPTSPECTRUM_H
3
4 /* Copyright(c) 1998-2010, ALICE Experiment at CERN, All rights reserved. *
5  * See cxx source for full Copyright notice                               */
6
7 //***********************************************************************
8 // Class AliHFPtSpectrum
9 // Base class for feed-down corrections on heavy-flavour decays
10 // computes the cross-section via one of the three implemented methods:
11 //   0) Consider no feed-down prediction 
12 //   1) Subtract the feed-down with the "fc" method 
13 //       Yield = Reco * fc;  where fc = 1 / ( 1 + (eff_b/eff_c)*(N_b/N_c) ) ;
14 //   2) Subtract the feed-down with the "Nb" method
15 //       Yield = Reco - Feed-down (exact formula on the function implementation)
16 //
17 //  (the corrected yields per bin are divided by the bin-width)
18 //
19 // Author: Z.Conesa, zconesa@in2p3.fr
20 //***********************************************************************
21
22 #include "TNamed.h"
23 #include "TH1.h"
24 #include "TMath.h"
25 #include "TGraphAsymmErrors.h"
26
27 class AliHFPtSpectrum: public TNamed
28 {
29
30  public:
31   
32   // Constructor
33   AliHFPtSpectrum(const char* name="AliHFPtSpectrum", const char* title="HF feed down correction class", Int_t option=1);
34   // Copy constructor
35   AliHFPtSpectrum(const AliHFPtSpectrum &rhs);
36   // Assignment operator
37   AliHFPtSpectrum& operator=(const AliHFPtSpectrum &source);
38   // Destructor
39   virtual ~AliHFPtSpectrum();
40
41   //
42   // Setters
43   //
44   // Set the theoretical direct & feeddown pt spectrum
45   void SetMCptSpectra(TH1D *hDirect, TH1D *hFeedDown);
46   // Set the theoretical feeddown pt spectrum
47   void SetFeedDownMCptSpectra(TH1D *hFeedDown);
48   // Set the theoretical direct & feeddown pt spectrum upper and lower bounds
49   void SetMCptDistributionsBounds(TH1D *hDirectMax, TH1D *hDirectMin, TH1D *hFeedDownMax, TH1D *hFeedDownMin);
50   // Set the theoretical feeddown pt spectrum upper and lower bounds
51   void SetFeedDownMCptDistributionsBounds(TH1D *hFeedDownMax, TH1D *hFeedDownMin);
52   // Set the acceptance and efficiency corrections for direct  
53   void SetDirectAccEffCorrection(TH1D *hDirectEff);
54   // Set the acceptance and efficiency corrections for direct & feeddown 
55   void SetAccEffCorrection(TH1D *hDirectEff, TH1D *hFeedDownEff);
56   // Set the reconstructed spectrum
57   void SetReconstructedSpectrum(TH1D *hRec);
58   void SetReconstructedSpectrumSystematics(TGraphAsymmErrors *gRec); 
59   // Set the calculation option flag for feed-down correction: 0=none, 1=fc , 2=Nb 
60   void SetFeedDownCalculationOption(Int_t option){ fFeedDownOption = option; }
61   // Set if the calculation has to consider asymmetric uncertaInt_ties or not
62   void SetComputeAsymmetricUncertainties(Bool_t flag){ fAsymUncertainties = flag; }
63   // Set the luminosity and its uncertainty
64   void SetLuminosity(Double_t luminosity, Double_t unc){
65     fLuminosity[0]=luminosity;  fLuminosity[1]=unc;
66   }
67   // Set the trigger efficiency and its uncertainty
68   void SetTriggerEfficiency(Double_t efficiency, Double_t unc){
69     fTrigEfficiency[0]=efficiency; fTrigEfficiency[1]=unc;
70   }
71   // Set global acceptance x efficiency correction uncertainty (in percentages)
72   void SetAccEffPercentageUncertainty(Double_t globalEffUnc, Double_t globalBCEffRatioUnc){
73     fGlobalEfficiencyUncertainties[0] = globalEffUnc;
74     fGlobalEfficiencyUncertainties[1] = globalBCEffRatioUnc;
75   }
76   // Set the normalization factors
77   void SetNormalization(Double_t normalization){
78     fLuminosity[0]=normalization; fTrigEfficiency[0]=1.0;
79   }
80   void SetNormalization(Double_t nevents, Double_t sigma){
81     fLuminosity[0]=nevents/sigma; fTrigEfficiency[0]=1.0;
82   }
83   void SetNormalization(Double_t nevents, Double_t sigma, Double_t sigmaunc){
84     fLuminosity[0] = nevents/sigma; 
85     fTrigEfficiency[0] = 1.0;
86     fLuminosity[1] = fLuminosity[0] * TMath::Sqrt( (1/nevents) + (sigmaunc/sigma)*(sigmaunc/sigma) );
87   }
88
89   //
90   // Getters
91   //
92   // Return the theoretical predictions used for the calculation (rebinned if needed)
93   TH1D * GetDirectTheoreticalSpectrum() { return (fhDirectMCpt ? (TH1D*)fhDirectMCpt : NULL); }
94   TH1D * GetDirectTheoreticalUpperLimitSpectrum() { return (fhDirectMCptMax ? (TH1D*)fhDirectMCptMax : NULL); }
95   TH1D * GetDirectTheoreticalLowerLimitSpectrum() { return (fhDirectMCptMin ? (TH1D*)fhDirectMCptMin : NULL); }
96   TH1D * GetFeedDownTheoreticalSpectrum() { return (fhFeedDownMCpt ? (TH1D*)fhFeedDownMCpt : NULL); }
97   TH1D * GetFeedDownTheoreticalUpperLimitSpectrum() { return (fhFeedDownMCptMax ? (TH1D*)fhFeedDownMCptMax : NULL); }
98   TH1D * GetFeedDownTheoreticalLowerLimitSpectrum() { return (fhFeedDownMCptMin ? (TH1D*)fhFeedDownMCptMin : NULL); }
99   // Return the acceptance and efficiency corrections (rebinned if needed)
100   TH1D * GetDirectAccEffCorrection() { return (fhDirectEffpt ? (TH1D*)fhDirectEffpt : NULL); }
101   TH1D * GetFeedDownAccEffCorrection() { return (fhFeedDownEffpt ? (TH1D*)fhFeedDownEffpt : NULL); }
102   // Return the TGraphAsymmErrors of the feed-down correction (extreme systematics)
103   TGraphAsymmErrors * GetFeedDownCorrectionFcExtreme() { return (fgFcExtreme ?  fgFcExtreme : NULL); }
104   // Return the TGraphAsymmErrors of the feed-down correction (conservative systematics)
105   TGraphAsymmErrors * GetFeedDownCorrectionFcConservative() { return (fgFcConservative ?  fgFcConservative : NULL); }
106   // Return the histogram of the feed-down correction
107   TH1D * GetHistoFeedDownCorrectionFc() { return (fhFc ?  (TH1D*)fhFc : NULL); }
108   // Return the histograms of the feed-down correction bounds
109   TH1D * GetHistoUpperLimitFeedDownCorrectionFc() { return (fhFcMax ? (TH1D*)fhFcMax : NULL); }
110   TH1D * GetHistoLowerLimitFeedDownCorrectionFc() { return (fhFcMin ? (TH1D*)fhFcMin : NULL); }
111   // Return the TGraphAsymmErrors of the yield after feed-down correction (systematics but feed-down) 
112   TGraphAsymmErrors * GetFeedDownCorrectedSpectrum() { return (fgYieldCorr ? fgYieldCorr : NULL); }
113   // Return the TGraphAsymmErrors of the yield after feed-down correction (feed-down extreme systematics)
114   TGraphAsymmErrors * GetFeedDownCorrectedSpectrumExtreme() { return (fgYieldCorrExtreme ? fgYieldCorrExtreme : NULL); }
115   // Return the TGraphAsymmErrors of the yield after feed-down correction (feed-down conservative systematics)
116   TGraphAsymmErrors * GetFeedDownCorrectedSpectrumConservative() { return (fgYieldCorrConservative ? fgYieldCorrConservative : NULL); }
117   // Return the histogram of the yield after feed-down correction 
118   TH1D * GetHistoFeedDownCorrectedSpectrum() { return (fhYieldCorr ? (TH1D*)fhYieldCorr : NULL); }
119   // Return the histogram of the yield after feed-down correction bounds
120   TH1D * GetHistoUpperLimitFeedDownCorrectedSpectrum() { return (fhYieldCorrMax ? (TH1D*)fhYieldCorrMax : NULL); }
121   TH1D * GetHistoLowerLimitFeedDownCorrectedSpectrum() { return (fhYieldCorrMin ? (TH1D*)fhYieldCorrMin : NULL); }
122   // Return the equivalent invariant cross-section TGraphAsymmErrors (systematics but feed-down) 
123   TGraphAsymmErrors * GetCrossSectionFromYieldSpectrum() { return (fgSigmaCorr ? fgSigmaCorr : NULL); }
124   // Return the equivalent invariant cross-section TGraphAsymmErrors (feed-down extreme systematics)
125   TGraphAsymmErrors * GetCrossSectionFromYieldSpectrumExtreme() { return (fgSigmaCorrExtreme ? fgSigmaCorrExtreme : NULL); }
126   // Return the equivalent invariant cross-section TGraphAsymmErrors (feed-down conservative systematics)
127   TGraphAsymmErrors * GetCrossSectionFromYieldSpectrumConservative() { return (fgSigmaCorrConservative ? fgSigmaCorrConservative : NULL); }
128   // Return the equivalent invariant cross-section histogram
129   TH1D * GetHistoCrossSectionFromYieldSpectrum() { return (fhSigmaCorr ? (TH1D*)fhSigmaCorr : NULL); }
130   // Return the equivalent invariant cross-section histogram bounds
131   TH1D * GetHistoUpperLimitCrossSectionFromYieldSpectrum() { return (fhSigmaCorrMax ? (TH1D*)fhSigmaCorrMax : NULL); }
132   TH1D * GetHistoLowerLimitCrossSectionFromYieldSpectrum() { return (fhSigmaCorrMin ? (TH1D*)fhSigmaCorrMin : NULL); }
133
134   //
135   // Main function:
136   //    Compute the invariant cross-section from the yield (correct it)
137   // variables : analysed delta_y, BR for the final correction, BR b --> decay (relative to the input theoretical prediction)
138   void ComputeHFPtSpectrum(Double_t deltaY=1.0, Double_t branchingRatioC=1.0, Double_t branchingRatioBintoFinalDecay=1.0);
139
140   // Compute the systematic uncertainties
141   //   taking as input the AliHFSystErr uncertainties
142   void ComputeSystUncertainties(Int_t decay, Bool_t combineFeedDown);
143   //
144   // Drawing the corrected spectrum comparing to theoretical prediction
145   void DrawSpectrum(TGraphAsymmErrors *gPrediction);
146
147   //
148   // Basic functions
149   // 
150   void EstimateAndSetDirectEfficiencyRecoBin(TH1D *hSimu, TH1D *hReco);
151   void EstimateAndSetFeedDownEfficiencyRecoBin(TH1D *hSimu, TH1D *hReco);
152   
153   //
154   // Functions to  reweight histograms for testing purposes: 
155   //   to reweight the simulation: hToReweight is reweighted as hReference/hToReweight
156   TH1D * ReweightHisto(TH1D *hToReweight, TH1D *hReference);
157   //   to reweight the reco-histos: hRecToReweight is reweighted as hReference/hMCToReweight
158   TH1D * ReweightRecHisto(TH1D *hRecToReweight, TH1D *hMCToReweight, TH1D *hMCReference);
159
160
161  protected:
162
163   // Initialization 
164   Bool_t Initialize();
165   
166   // Basic functions
167   //
168   // Compute the feed-down correction via fc-method
169   void CalculateFeedDownCorrectionFc(); 
170   // Correct the yield for feed-down correction via fc-method
171   void CalculateFeedDownCorrectedSpectrumFc(); 
172   // Correct the yield for feed-down correction via Nb-method
173   void CalculateFeedDownCorrectedSpectrumNb(Double_t deltaY, Double_t branchingRatioBintoFinalDecay); 
174
175   // Check histograms consistency function
176   Bool_t CheckHistosConsistency(TH1D *h1, TH1D *h2);
177   // Function to rebin the theoretical spectra in the data-reconstructed spectra binning
178   TH1D * RebinTheoreticalSpectra(TH1D *hTheory, const char *name);
179   // Function to estimate the efficiency in the data-reconstructed spectra binning
180   TH1D * EstimateEfficiencyRecoBin(TH1D *hSimu, TH1D *hReco, const char *name);
181
182
183   //
184   // Input spectra
185   //
186   TH1D *fhDirectMCpt;            // Input MC c-->D spectra
187   TH1D *fhFeedDownMCpt;          // Input MC b-->D spectra
188   TH1D *fhDirectMCptMax;         // Input MC maximum c-->D spectra
189   TH1D *fhDirectMCptMin;         // Input MC minimum c-->D spectra
190   TH1D *fhFeedDownMCptMax;       // Input MC maximum b-->D spectra
191   TH1D *fhFeedDownMCptMin;       // Input MC minimum b-->D spectra
192   TH1D *fhDirectEffpt;           // c-->D Acceptance and efficiency correction
193   TH1D *fhFeedDownEffpt;         // b-->D Acceptance and efficiency correction
194   TH1D *fhRECpt;                 // all reconstructed D
195   //
196   TGraphAsymmErrors *fgRECSystematics; // all reconstructed D Systematic uncertainties
197   //
198   // Normalization factors
199   Double_t fLuminosity[2];           // analyzed luminosity & uncertainty
200   Double_t fTrigEfficiency[2];       // trigger efficiency & uncertainty
201   Double_t fGlobalEfficiencyUncertainties[2]; // uncertainties on the efficiency [0]=c, b, [1]=b/c
202
203   //
204   // Output spectra
205   //
206   TH1D *fhFc;                            // Correction histo fc = 1 / ( 1 + (eff_b/eff_c)*(N_b/N_c) ) 
207   TH1D *fhFcMax;                         // Maximum fc histo
208   TH1D *fhFcMin;                         // Minimum fc histo
209   TGraphAsymmErrors * fgFcExtreme;       // Extreme correction as TGraphAsymmErrors
210   TGraphAsymmErrors * fgFcConservative;  // Extreme correction as TGraphAsymmErrors
211   TH1D *fhYieldCorr;                     // Corrected yield (stat unc. only)
212   TH1D *fhYieldCorrMax;                  // Maximum corrected yield  
213   TH1D *fhYieldCorrMin;                  // Minimum corrected yield  
214   TGraphAsymmErrors * fgYieldCorr;              // Corrected yield as TGraphAsymmErrors  (syst but feed-down)
215   TGraphAsymmErrors * fgYieldCorrExtreme;       // Extreme corrected yield as TGraphAsymmErrors  (syst from feed-down)
216   TGraphAsymmErrors * fgYieldCorrConservative;  // Conservative corrected yield as TGraphAsymmErrors  (syst from feed-down) 
217   TH1D *fhSigmaCorr;                     // Corrected cross-section (stat unc. only)
218   TH1D *fhSigmaCorrMax;                  // Maximum corrected cross-section  
219   TH1D *fhSigmaCorrMin;                  // Minimum corrected cross-section
220   TGraphAsymmErrors * fgSigmaCorr;              // Corrected cross-section as TGraphAsymmErrors (syst but feed-down)
221   TGraphAsymmErrors * fgSigmaCorrExtreme;       // Extreme corrected cross-section as TGraphAsymmErrors (syst from feed-down)
222   TGraphAsymmErrors * fgSigmaCorrConservative;  // Conservative corrected cross-section as TGraphAsymmErrors  (syst from feed-down)
223
224   //
225   Int_t fFeedDownOption;            // feed-down correction flag: 0=none, 1=fc, 2=Nb 
226   Bool_t fAsymUncertainties;        // flag: asymmetric uncertainties are (1) or not (0) considered
227
228
229   ClassDef(AliHFPtSpectrum,1) // Class for Heavy Flavor spectra corrections
230 };
231
232 #endif