]> git.uio.no Git - u/mrichter/AliRoot.git/blob - PWGJE/UserTasks/AliAnalysisTaskPID.h
- Added threshold for spline systematics and different sys. scale factors below and...
[u/mrichter/AliRoot.git] / PWGJE / UserTasks / AliAnalysisTaskPID.h
1 #ifndef ALI_ANALYSIS_TASK_PID_H
2 #define ALI_ANALYSIS_TASK_PID_H
3
4 /*
5 This task collects PID output from different detectors.
6 Only tracks fulfilling some standard quality cuts are taken into account.
7 At the moment, only data from TPC and TOF is collected. But in future,
8 data from e.g. HMPID is also foreseen.
9
10 Class written by Benjamin Hess.
11 Contact: bhess@cern.ch
12 */
13
14 class TF1;
15 class TRandom3;
16 class AliAnalysisFilter;
17 class AliCFContainer;
18 class AliESDEvent;
19 class AliMCEvent;
20 class AliMCParticle;
21 class AliPID;
22 class AliPIDCombined;
23 class AliPIDResponse;
24 class AliTOFPIDResponse;
25 class AliVEvent;
26 class AliVTrack;
27
28 #include "TH1D.h"
29 #include "TH2D.h"
30 #include "TH3D.h"
31 #include "THnSparse.h"
32 #include "TProfile.h"
33 #include "TString.h"
34
35 #include "AliCentrality.h"
36 #include "AliCFContainer.h"
37
38 #include "AliPID.h"
39 #include "AliAnalysisTaskPIDV0base.h"
40
41 class AliAnalysisTaskPID : public AliAnalysisTaskPIDV0base {
42  public:
43   AliAnalysisTaskPID();
44   AliAnalysisTaskPID(const char *name);
45   virtual ~AliAnalysisTaskPID();
46   
47   virtual void   UserCreateOutputObjects();
48   virtual void   UserExec(Option_t *option);
49   virtual void   Terminate(const Option_t*);
50   
51   enum ErrorCode { kNoErrors = 1, kWarning = 0, kError = -1};
52   
53   enum dataAxes { kDataMCID = 0, kDataSelectSpecies = 1, kDataPt = 2, kDataDeltaPrimeSpecies = 3, kDataCentrality = 4,
54                   kDataJetPt = 5, kDataZ = 6, kDataXi = 7, kDataCharge = 8, kDataTOFpidInfo = 9, kDataNumAxes = 10 };
55
56   enum genAxes { kGenMCID = 0, kGenSelectSpecies = 1, kGenPt = 2, kGenDeltaPrimeSpecies = 3, kGenCentrality = 4,
57                  kGenJetPt = 5, kGenZ = 6, kGenXi = 7, kGenCharge = 8, kGenTOFpidInfo = 9, kGenNumAxes = 10 };
58   
59   enum genYieldAxes { kGenYieldMCID = 0, kGenYieldPt = 1, kGenYieldCentrality = 2, kGenYieldJetPt = 3, kGenYieldZ = 4, kGenYieldXi = 5,
60                       kGenYieldCharge = 6, kGenYieldNumAxes = 7 };
61   
62   enum ptResolutionAxes { kPtResJetPt = 0, kPtResGenPt = 1, kPtResRecPt = 2, kPtResCharge = 3, kPtResCentrality = 4, kPtResNumAxes = 5 };
63   
64   enum efficiencyAxes { kEffMCID = 0, kEffTrackPt = 1, kEffTrackEta = 2, kEffTrackCharge = 3, kEffCentrality = 4, kEffJetPt = 5,
65                         kEffZ = 6, kEffXi = 7, kEffNumAxes = 8 };
66   
67   enum EffSteps { kStepGenWithGenCuts = 0, kStepRecWithGenCuts = 1, kStepRecWithGenCutsMeasuredObs = 2,
68                   kStepRecWithRecCutsMeasuredObs = 3, kStepRecWithRecCutsMeasuredObsPrimaries = 4,
69                   kStepRecWithRecCutsMeasuredObsStrangenessScaled = 5, kStepRecWithRecCutsPrimaries = 6, kNumSteps = 7};
70   
71   enum TOFpidInfo { kNoTOFinfo = -2, kNoTOFpid = -1, kTOFpion = 0, kTOFkaon = 1, kTOFproton = 2, kNumTOFspecies = 3,
72                     kNumTOFpidInfoBins = 5 };
73   
74   enum EventCounterType { kTriggerSel = 0, kTriggerSelAndVtxCut = 1, kTriggerSelAndVtxCutAndZvtxCut = 2 };
75   
76   static Int_t PDGtoMCID(Int_t pdg);
77   
78   static void GetJetTrackObservables(Double_t trackPt, Double_t jetPt, Double_t& z, Double_t& xi);
79   
80   static Double_t GetMCStrangenessFactorCMS(Int_t motherPDG, Double_t motherGenPt);
81   static Double_t GetMCStrangenessFactorCMS(AliMCEvent* mcEvent, AliMCParticle* daughter);
82   
83   static Bool_t IsSecondaryWithStrangeMotherMC(AliMCEvent* mcEvent, Int_t partLabel);
84   
85   Int_t GetIndexOfChargeAxisData() const
86     { return fStoreAdditionalJetInformation ? kDataCharge : kDataCharge - fgkNumJetAxes; };
87   Int_t GetIndexOfChargeAxisGen() const
88     { return fStoreAdditionalJetInformation ? kGenCharge : kGenCharge - fgkNumJetAxes; };
89   Int_t GetIndexOfChargeAxisGenYield() const
90     { return fStoreAdditionalJetInformation ? kGenYieldCharge : kGenYieldCharge - fgkNumJetAxes; };
91   
92   Int_t GetIndexOfTOFpidInfoAxisData() const
93     { return fStoreAdditionalJetInformation ? kDataTOFpidInfo : kDataTOFpidInfo - fgkNumJetAxes; };
94   Int_t GetIndexOfTOFpidInfoAxisGen() const
95     { return fStoreAdditionalJetInformation ? kGenTOFpidInfo : kGenTOFpidInfo - fgkNumJetAxes; };
96   
97   Bool_t FillXsec(Double_t xsection)
98     { if (!fh1Xsec) return kFALSE; fh1Xsec->Fill("<#sigma>", xsection); return kTRUE; };
99   Bool_t FillPythiaTrials(Double_t avgTrials)
100     { if (!fh1Trials) return kFALSE; fh1Trials->Fill("#sum{ntrials}", avgTrials); return kTRUE; };
101     
102   Bool_t FillEfficiencyContainer(const Double_t* values, EffSteps step, Double_t weight = 1.0);
103   
104   Bool_t FillGeneratedYield(const Double_t* values, Double_t weight = 1.0);
105   Bool_t FillPtResolution(Int_t mcID, const Double_t* values);
106   Bool_t FillGenJets(Double_t centralityPercentile, Double_t jetPt, Double_t norm = -1.);
107   Bool_t FillRecJets(Double_t centralityPercentile, Double_t jetPt, Double_t norm = -1.);
108   
109   Bool_t IncrementEventCounter(Double_t centralityPercentile, EventCounterType type);
110   
111   void PostOutputData();
112   
113   void PrintSettings(Bool_t printSystematicsSettings = kFALSE) const;
114   
115   void PrintSystematicsSettings() const;
116   
117   Bool_t ProcessTrack(const AliVTrack* track, Int_t particlePDGcode, Double_t centralityPercentile, Double_t jetPt) ;
118   
119   ErrorCode GenerateDetectorResponse(ErrorCode errCode, Double_t mean, Double_t sigma, Double_t* responses,
120                                      Int_t nResponses,
121                                      Bool_t usePureGaus = kFALSE);
122   ErrorCode SetParamsForConvolutedGaus(Double_t gausMean, Double_t gausSigma);
123   
124   const TString GetCentralityEstimator() const { return fCentralityEstimator; };
125   void SetCentralityEstimator(TString estimator) { fCentralityEstimator = estimator; };
126   
127   Double_t GetCentralityPercentile(AliVEvent* evt) const;
128   
129   inline Double_t GetConvolutedGaussTransitionPar(Int_t index) const;
130   
131   Bool_t SetConvolutedGaussLambdaParameter(Double_t lambda);
132
133   Bool_t GetInputFromOtherTask() const { return fInputFromOtherTask; };
134   void SetInputFromOtherTask(Bool_t flag) { fInputFromOtherTask = flag; };
135   
136   Bool_t GetDoPID() const { return fDoPID; };
137   void SetDoPID(Bool_t flag) { fDoPID = flag; };
138   
139   Bool_t GetDoEfficiency() const { return fDoEfficiency; };
140   void SetDoEfficiency(Bool_t flag) { fDoEfficiency = flag; };
141   
142   Bool_t GetDoPtResolution() const { return fDoPtResolution; };
143   void SetDoPtResolution(Bool_t flag) { fDoPtResolution = flag; };
144   
145   Bool_t GetStoreCentralityPercentile() const { return fStoreCentralityPercentile; };
146   void SetStoreCentralityPercentile(Bool_t flag) { fStoreCentralityPercentile = flag; };
147   
148   Bool_t GetStoreAdditionalJetInformation() const { return fStoreAdditionalJetInformation; };
149   void SetStoreAdditionalJetInformation(Bool_t flag) { fStoreAdditionalJetInformation = flag; };
150   
151   Bool_t GetUseMCidForGeneration() const { return fUseMCidForGeneration; };
152   void SetUseMCidForGeneration(Bool_t flag) { fUseMCidForGeneration = flag; };
153   
154   Bool_t GetUseConvolutedGaus() const { return fUseConvolutedGaus; };
155   void SetUseConvolutedGaus(Bool_t flag) { fUseConvolutedGaus = flag; };
156   
157   Double_t GetAccuracyNonGaussianTail() const { return fAccuracyNonGaussianTail; };
158   void SetAccuracyNonGaussianTail(Double_t value) { fAccuracyNonGaussianTail = value; };
159   
160   Bool_t GetTakeIntoAccountMuons() const { return fTakeIntoAccountMuons; };
161   void SetTakeIntoAccountMuons(Bool_t flag) { fTakeIntoAccountMuons = flag; };
162   
163   Int_t GetTOFmode() const { return fTOFmode; };
164   void SetTOFmode(Int_t tofMode) { fTOFmode = tofMode; };
165   
166   Bool_t GetUseTPCDefaultPriors() const { return fTPCDefaultPriors; };
167   void SetUseTPCDefaultPriors(Bool_t flag) { fTPCDefaultPriors = flag; };
168   
169   Bool_t GetUsePriors() const { return fUsePriors; };
170   void SetUsePriors(Bool_t flag) { fUsePriors = flag; };
171   
172   Bool_t GetUseITS() const { return fUseITS; };
173   void SetUseITS(Bool_t flag) { fUseITS = flag; };
174   
175   Bool_t GetUseTOF() const { return fUseTOF; };
176   void SetUseTOF(Bool_t flag) { fUseTOF = flag; };
177   
178   Double_t GetEtaAbsCutLow() const { return fEtaAbsCutLow; };
179   Double_t GetEtaAbsCutUp() const { return fEtaAbsCutUp; };
180   Bool_t SetEtaAbsCutRange(Double_t lowerLimit, Double_t upperLimit);
181   
182   Bool_t IsInAcceptedEtaRange(Double_t etaAbs) const { return (etaAbs >= fEtaAbsCutLow && etaAbs <= fEtaAbsCutUp); };
183   
184   Double_t GetSystematicScalingSplinesThreshold() const { return fSystematicScalingSplinesThreshold; };
185   void SetSystematicScalingSplinesThreshold(Double_t threshold) { fSystematicScalingSplinesThreshold = threshold; };
186   
187   Double_t GetSystematicScalingSplinesBelowThreshold() const { return fSystematicScalingSplinesBelowThreshold; };
188   void SetSystematicScalingSplinesBelowThreshold(Double_t scaleFactor) 
189     { fSystematicScalingSplinesBelowThreshold = scaleFactor; CheckDoAnyStematicStudiesOnTheExpectedSignal(); };
190     
191   Double_t GetSystematicScalingSplinesAboveThreshold() const { return fSystematicScalingSplinesAboveThreshold; };
192   void SetSystematicScalingSplinesAboveThreshold(Double_t scaleFactor) 
193     { fSystematicScalingSplinesAboveThreshold = scaleFactor; CheckDoAnyStematicStudiesOnTheExpectedSignal(); };
194   
195   Double_t GetSystematicScalingEtaCorrectionMomentumThr() const { return fSystematicScalingEtaCorrectionMomentumThr; };
196   void SetSystematicScalingEtaCorrectionMomentumThr(Double_t threshold) { fSystematicScalingEtaCorrectionMomentumThr = threshold; };
197   
198   Double_t GetSystematicScalingEtaCorrectionLowMomenta() const { return fSystematicScalingEtaCorrectionLowMomenta; };
199   void SetSystematicScalingEtaCorrectionLowMomenta(Double_t scaleFactor) 
200     { fSystematicScalingEtaCorrectionLowMomenta = scaleFactor; CheckDoAnyStematicStudiesOnTheExpectedSignal(); };
201   
202   Double_t GetSystematicScalingEtaCorrectionHighMomenta() const { return fSystematicScalingEtaCorrectionHighMomenta; };
203   void SetSystematicScalingEtaCorrectionHighMomenta(Double_t scaleFactor) 
204     { fSystematicScalingEtaCorrectionHighMomenta = scaleFactor; CheckDoAnyStematicStudiesOnTheExpectedSignal(); };
205   
206   Double_t GetSystematicScalingEtaSigmaPara() const { return fSystematicScalingEtaSigmaPara; };
207   void SetSystematicScalingEtaSigmaPara(Double_t scaleFactor)
208     { fSystematicScalingEtaSigmaPara = scaleFactor; CheckDoAnyStematicStudiesOnTheExpectedSignal(); };
209   
210   Double_t GetSystematicScalingMultCorrection() const { return fSystematicScalingMultCorrection; };
211   void SetSystematicScalingMultCorrection(Double_t scaleFactor) 
212     { fSystematicScalingMultCorrection = scaleFactor; CheckDoAnyStematicStudiesOnTheExpectedSignal(); };
213   
214   void CleanupParticleFractionHistos();
215   Bool_t GetParticleFraction(Double_t trackPt, Double_t jetPt, Double_t multiplicity,
216                              AliPID::EParticleType species, Double_t& fraction, Double_t& fractionErrorStat,
217                              Double_t& fractionErrorSys) const;
218   Bool_t GetParticleFractions(Double_t trackPt, Double_t jetPt, Double_t centralityPercentile,
219                               Double_t* prob, Int_t smearSpeciesByError, Int_t takeIntoAccountSpeciesSysError,
220                               Bool_t uniformSystematicError = kFALSE) const;
221   const TH3D* GetParticleFractionHisto(Int_t species, Bool_t sysError = kFALSE) const;
222   Bool_t SetParticleFractionHisto(const TH3D* hist, Int_t species, Bool_t sysError = kFALSE);
223   Int_t GetParticleFractionHistoNbinsTrackPt() const;
224   Int_t GetParticleFractionHistoNbinsJetPt() const;
225   Int_t GetParticleFractionHistoNbinsCentrality() const;
226   Bool_t SetParticleFractionHistosFromFile(const TString filePathName, Bool_t sysError = kFALSE);
227   Int_t GetRandomParticleTypeAccordingToParticleFractions(Double_t trackPt, Double_t jetPt, 
228                                                           Double_t centralityPercentile,
229                                                           Bool_t smearByError,
230                                                           Bool_t takeIntoAccountSysError = kFALSE) const;
231   
232   TOFpidInfo GetTOFType(const AliVTrack* track, Int_t tofMode) const;
233   
234  protected:
235   void CheckDoAnyStematicStudiesOnTheExpectedSignal();
236   Double_t ConvolutedGaus(const Double_t* xx, const Double_t* par) const;
237   inline Double_t FastGaus(Double_t x, Double_t mean, Double_t sigma) const;
238   inline Double_t FastNormalisedGaus(Double_t x, Double_t mean, Double_t sigma) const;
239   Int_t FindBinWithinRange(TAxis* axis, Double_t value) const;
240   Int_t FindFirstBinAboveIn3dSubset(const TH3* hist, Double_t threshold, Int_t yValue, Int_t zValue) const;
241   Int_t FindLastBinAboveIn3dSubset(const TH3* hist, Double_t threshold, Int_t yValue, Int_t zValue) const;
242   virtual void SetUpGenHist(THnSparse* hist, Double_t* binsPt, Double_t* binsDeltaPrime, Double_t* binsCent, Double_t* binsJetPt) const;
243   virtual void SetUpGenYieldHist(THnSparse* hist, Double_t* binsPt, Double_t* binsCent, Double_t* binsJetPt) const;
244   virtual void SetUpHist(THnSparse* hist, Double_t* binsPt, Double_t* binsDeltaPrime, Double_t* binsCent, Double_t* binsJetPt) const;
245   virtual void SetUpPtResHist(THnSparse* hist, Double_t* binsPt, Double_t* binsJetPt, Double_t* binsCent) const;
246   virtual void SetUpPIDcombined();
247   
248   static const Int_t fgkNumJetAxes; // Number of additional axes for jets
249   static const Double_t fgkEpsilon; // Double_t threshold above zero
250   static const Int_t fgkMaxNumGenEntries; // Maximum number of generated detector responses per track and delta(Prime) and associated species
251
252  private:
253   static const Double_t fgkOneOverSqrt2; // = 1. / TMath::Sqrt2();
254   
255   AliPIDCombined* fPIDcombined; //! PID combined object
256   
257   Bool_t fInputFromOtherTask; // If set to kTRUE, no events are processed and the input must be fed in from another task. If set to kFALSE, normal event processing
258   
259   Bool_t fDoPID; // Only do PID processing (and post the output), if flag is set to kTRUE
260   Bool_t fDoEfficiency; // Only do efficiency processing (and post the output), if flag is set to kTRUE
261   Bool_t fDoPtResolution; // Only do pT resolution processing (and post the output), if flag is set to kTRUE
262   
263   Bool_t fStoreCentralityPercentile; // If set to kTRUE, store centrality percentile for each event. In case of kFALSE (appropriate for pp), centrality percentile will be set to -1 for every event
264   Bool_t fStoreAdditionalJetInformation; // If set to kTRUE, additional jet information like jetPt, z, xi will be stored in the THnSparses
265
266   Bool_t fTakeIntoAccountMuons; // Also take into account muons for the generation of the expected response and the most probable PID
267   Bool_t fUseITS; // Use ITS for PID combined probabilities
268   Bool_t fUseTOF; // Use TOF for PID combined probabilities
269   Bool_t fUsePriors; // Use priors for PID combined probabilities
270   Bool_t fTPCDefaultPriors; // Use TPC default priors for PID combined probabilities, if priors are enabled
271     
272   Bool_t fUseMCidForGeneration; // If MC, use MCid instead of PIDcombined to generate the signals
273   
274   Bool_t fUseConvolutedGaus; // Use convoluted gaus to generate detector response instead of pure gaus  
275   const Int_t fkConvolutedGausNPar; // Number of parameters for convoluted gaus
276   Double_t fAccuracyNonGaussianTail; // Accuracy of the non-gaussian tail (fraction of the maximum)
277   const Double_t fkDeltaPrimeLowLimit; // Lower deltaPrime limit
278   const Double_t fkDeltaPrimeUpLimit; // Upper deltaPrime limit
279   TF1* fConvolutedGausDeltaPrime; // Gaus convoluted with exponential tail to generate detector response (deltaPrime)
280   
281   Int_t fTOFmode; // TOF mode used for TOF PID info (affects num sigma inclusion/exclusion)
282   Double_t fConvolutedGaussTransitionPars[3]; // Parameter for transition from gaussian parameters to asymmetric shape
283   static const Double_t fgkSigmaReferenceForTransitionPars; // Reference sigma chosen to calculate transition parameters
284   
285   Double_t fEtaAbsCutLow; // Lower cut value on |eta|
286   Double_t fEtaAbsCutUp;  // Upper cut value on |eta|
287   
288   // For systematic studies
289   Bool_t   fDoAnySystematicStudiesOnTheExpectedSignal; // Internal flag indicating whether any systematic studies are going to be performed
290   Double_t fSystematicScalingSplinesThreshold;         // beta-gamma threshold for the systematic spline scale factor
291   Double_t fSystematicScalingSplinesBelowThreshold;        // Systematic scale factor for the splines (1. = no systematics) below threshold
292   Double_t fSystematicScalingSplinesAboveThreshold;        // Systematic scale factor for the splines (1. = no systematics) above threshold
293   Double_t fSystematicScalingEtaCorrectionMomentumThr;  // Momentum threshold for the systematic scale factor for the eta correction (separates low-p from high-p
294   Double_t fSystematicScalingEtaCorrectionLowMomenta;   // Systematic scale factor for the eta correction (1. = no systematics) at low momenta
295   Double_t fSystematicScalingEtaCorrectionHighMomenta;  // Systematic scale factor for the eta correction (1. = no systematics) at high momenta
296   Double_t fSystematicScalingEtaSigmaPara;   // Systematic scale factor for the parametrisation of the relative signal width (1. = no systematics) 
297   Double_t fSystematicScalingMultCorrection; // Systematic scale factor for the multiplicity correction (1. = no systematics) 
298   
299   TH3D* fFractionHists[AliPID::kSPECIES]; // 3D histos of particle fraction as function  with trackPt, jetPt (-1 for inclusive spectra), centralityPercentile (-1 for pp)
300   TH3D* fFractionSysErrorHists[AliPID::kSPECIES]; // 3D histos of sys. error of particle fraction as function  with trackPt, jetPt (-1 for inclusive spectra), centralityPercentile (-1 for pp)
301   
302   TString fCentralityEstimator; // Estimator for the centrality (e.g. V0A, V0M)
303   
304   THnSparseD* fhPIDdataAll; //! Data histo
305   
306   // Generated response information
307   THnSparseD* fhGenEl; //! Generated response for el
308   THnSparseD* fhGenKa; //! Generated response for ka
309   THnSparseD* fhGenPi; //! Generated response for pi
310   THnSparseD* fhGenMu; //! Generated response for mu
311   THnSparseD* fhGenPr; //! Generated response for pr
312   
313   // Generated responses for a single track
314   Double_t* fGenRespElDeltaPrimeEl; //! Generated responses for a single track
315   Double_t* fGenRespElDeltaPrimeKa; //! Generated responses for a single track
316   Double_t* fGenRespElDeltaPrimePi; //! Generated responses for a single track
317   Double_t* fGenRespElDeltaPrimePr; //! Generated responses for a single track
318   Double_t* fGenRespKaDeltaPrimeEl; //! Generated responses for a single track
319   Double_t* fGenRespKaDeltaPrimeKa; //! Generated responses for a single track
320   Double_t* fGenRespKaDeltaPrimePi; //! Generated responses for a single track
321   Double_t* fGenRespKaDeltaPrimePr; //! Generated responses for a single track
322   Double_t* fGenRespPiDeltaPrimeEl; //! Generated responses for a single track
323   Double_t* fGenRespPiDeltaPrimeKa; //! Generated responses for a single track
324   Double_t* fGenRespPiDeltaPrimePi; //! Generated responses for a single track
325   Double_t* fGenRespPiDeltaPrimePr; //! Generated responses for a single track
326   Double_t* fGenRespMuDeltaPrimeEl; //! Generated responses for a single track
327   Double_t* fGenRespMuDeltaPrimeKa; //! Generated responses for a single track
328   Double_t* fGenRespMuDeltaPrimePi; //! Generated responses for a single track
329   Double_t* fGenRespMuDeltaPrimePr; //! Generated responses for a single track
330   Double_t* fGenRespPrDeltaPrimeEl; //! Generated responses for a single track
331   Double_t* fGenRespPrDeltaPrimeKa; //! Generated responses for a single track
332   Double_t* fGenRespPrDeltaPrimePi; //! Generated responses for a single track
333   Double_t* fGenRespPrDeltaPrimePr; //! Generated responses for a single track
334   /*
335   Double_t* fGenRespElDeltaEl; //! Generated responses for a single track
336   Double_t* fGenRespElDeltaKa; //! Generated responses for a single track
337   Double_t* fGenRespElDeltaPi; //! Generated responses for a single track
338   Double_t* fGenRespElDeltaPr; //! Generated responses for a single track
339   Double_t* fGenRespKaDeltaEl; //! Generated responses for a single track
340   Double_t* fGenRespKaDeltaKa; //! Generated responses for a single track
341   Double_t* fGenRespKaDeltaPi; //! Generated responses for a single track
342   Double_t* fGenRespKaDeltaPr; //! Generated responses for a single track
343   Double_t* fGenRespPiDeltaEl; //! Generated responses for a single track
344   Double_t* fGenRespPiDeltaKa; //! Generated responses for a single track
345   Double_t* fGenRespPiDeltaPi; //! Generated responses for a single track
346   Double_t* fGenRespPiDeltaPr; //! Generated responses for a single track
347   Double_t* fGenRespMuDeltaEl; //! Generated responses for a single track
348   Double_t* fGenRespMuDeltaKa; //! Generated responses for a single track
349   Double_t* fGenRespMuDeltaPi; //! Generated responses for a single track
350   Double_t* fGenRespMuDeltaPr; //! Generated responses for a single track
351   Double_t* fGenRespPrDeltaEl; //! Generated responses for a single track
352   Double_t* fGenRespPrDeltaKa; //! Generated responses for a single track
353   Double_t* fGenRespPrDeltaPi; //! Generated responses for a single track
354   Double_t* fGenRespPrDeltaPr; //! Generated responses for a single track
355   */
356   
357   TH1D* fhEventsProcessed; //! Histo holding the number of processed events (i.e. passing trigger selection, vtx and zvtx cuts
358   TH1D* fhEventsTriggerSel; //! Histo holding the number of events passing trigger selection
359   TH1D* fhEventsTriggerSelVtxCut; //! Histo holding the number of events passing trigger selection and vtx cut
360   
361   TH2D* fhSkippedTracksForSignalGeneration; //! Number of tracks that have been skipped for the signal generation
362   THnSparseD* fhMCgeneratedYieldsPrimaries; //! Histo holding the generated (no reco, no cuts) primary particle yields in considered eta range
363   
364   TH2D* fh2FFJetPtRec;            //! Number of reconstructed jets vs. jetPt and centrality
365   TH2D* fh2FFJetPtGen;            //! Number of generated jets vs. jetPt and centrality
366   
367   TProfile* fh1Xsec;              //! pythia cross section and trials
368   TH1D*     fh1Trials;            //! sum of trials
369   
370   AliCFContainer* fContainerEff; //! Container for efficiency determination
371   
372   THnSparseD* fPtResolution[AliPID::kSPECIES]; //! Pt Resolution for the individual species
373   
374   TObjArray* fOutputContainer;  //! output data container
375   
376   TObjArray* fPtResolutionContainer;  //! output data container for pt resolution
377   
378   AliAnalysisTaskPID(const AliAnalysisTaskPID&); // not implemented
379   AliAnalysisTaskPID& operator=(const AliAnalysisTaskPID&); // not implemented
380   
381   ClassDef(AliAnalysisTaskPID, 16);
382 };
383
384
385 //_____________________________________________________________________________
386 inline Bool_t AliAnalysisTaskPID::FillEfficiencyContainer(const Double_t* values, AliAnalysisTaskPID::EffSteps step,
387                                                           Double_t weight) 
388 {
389   // Fill efficiency container at step "step" with the values
390   
391   if (!fDoEfficiency)
392     return kFALSE;
393   
394   if (!fContainerEff) {
395     AliError("Efficiency container not initialised -> cannot be filled!");
396     return kFALSE;
397   }
398   
399   fContainerEff->Fill(values, step, weight);    
400   
401   return kTRUE;
402 }
403
404
405 //_____________________________________________________________________________
406 inline Bool_t AliAnalysisTaskPID::FillGeneratedYield(const Double_t* values, Double_t weight)
407 {
408   // Fill histos with generated primary yields with provided values
409   
410   if (!fDoPID)
411     return kFALSE;
412   
413   if (!fhMCgeneratedYieldsPrimaries) {
414     AliError("Histo for generated primary yield not initialised -> cannot be filled!");
415     return kFALSE;
416   }
417   
418   fhMCgeneratedYieldsPrimaries->Fill(values, weight);
419     
420   return kTRUE;
421 }
422
423
424 //_____________________________________________________________________________
425 inline Bool_t AliAnalysisTaskPID::FillGenJets(Double_t centralityPercentile, Double_t jetPt, Double_t norm)
426 {
427   if (!fDoPID && !fDoEfficiency)
428     return kFALSE;
429   
430   if (!fh2FFJetPtGen)
431     return kFALSE;
432   
433   if (norm > 0.)
434     fh2FFJetPtGen->Fill(centralityPercentile, jetPt, 1. / norm);
435   else
436     fh2FFJetPtGen->Fill(centralityPercentile, jetPt);
437   
438   return kTRUE;
439 }
440
441
442 //_____________________________________________________________________________
443 inline Bool_t AliAnalysisTaskPID::FillRecJets(Double_t centralityPercentile, Double_t jetPt, Double_t norm)
444 {
445   if (!fDoPID && !fDoEfficiency)
446     return kFALSE;
447   
448   if (!fh2FFJetPtRec)
449     return kFALSE;
450   
451   if (norm > 0.)
452     fh2FFJetPtRec->Fill(centralityPercentile, jetPt, 1. / norm);
453   else
454     fh2FFJetPtRec->Fill(centralityPercentile, jetPt);
455   
456   return kTRUE;
457 }
458
459
460 //_____________________________________________________________________________
461 inline Bool_t AliAnalysisTaskPID::FillPtResolution(Int_t mcID, const Double_t* values)
462 {
463   // Fill histos with pT resolution with provided values
464   
465   if (!fDoPtResolution || mcID < 0 || mcID >= AliPID::kSPECIES)
466     return kFALSE;
467   
468   if (!fPtResolution[mcID]) {
469     AliError(Form("Histo for pT resolution (species: %s) not initialised -> cannot be filled!", AliPID::ParticleName(mcID)));
470     return kFALSE;
471   }
472   
473   fPtResolution[mcID]->Fill(values);
474     
475   return kTRUE;
476 }
477  
478
479 //_____________________________________________________________________________
480 inline Bool_t AliAnalysisTaskPID::IncrementEventCounter(Double_t centralityPercentile, AliAnalysisTaskPID::EventCounterType type)
481 {
482   // Increment the number of events for the given centrality percentile and the corresponding counter type
483   
484   if (type == kTriggerSel) {
485     if (!fhEventsTriggerSel) {
486       AliError("Histogram for number of events (kTriggerSel) not initialised -> cannot be incremented!");
487       return kFALSE;
488     }
489     
490     fhEventsTriggerSel->Fill(centralityPercentile);
491   }
492   else if (type == kTriggerSelAndVtxCut) {
493     if (!fhEventsTriggerSelVtxCut) {
494       AliError("Histogram for number of events (kTriggerSelAndVtxCut) not initialised -> cannot be incremented!");
495       return kFALSE;
496     }
497     
498     fhEventsTriggerSelVtxCut->Fill(centralityPercentile);
499   }
500   else if (type == kTriggerSelAndVtxCutAndZvtxCut) {
501     if (!fhEventsProcessed) {
502       AliError("Histogram for number of events (kTriggerSelAndVtxCutAndZvtxCut) not initialised -> cannot be incremented!");
503       return kFALSE;
504     }
505     
506     fhEventsProcessed->Fill(centralityPercentile);
507   }
508   
509   
510   return kTRUE;
511 };
512
513
514 //_____________________________________________________________________________
515 inline Bool_t AliAnalysisTaskPID::SetEtaAbsCutRange(Double_t lowerLimit, Double_t upperLimit)
516 {
517   if (lowerLimit >= upperLimit) {
518     AliError(Form("Requested lower |eta| cut limit >= upper |eta| cut limit. Old eta cut range will be used (low %f, high %f).",
519                   fEtaAbsCutLow, fEtaAbsCutUp));
520     return kFALSE;
521   }
522   
523   fEtaAbsCutLow = lowerLimit;
524   fEtaAbsCutUp = upperLimit;
525   
526   return kTRUE;
527 };
528
529
530 //_____________________________________________________________________________
531 inline Double_t AliAnalysisTaskPID::GetConvolutedGaussTransitionPar(Int_t index) const
532 {
533   if (index < 0 || index >= 3) {
534     printf("Invalid index %d!\n", index);
535     return -1;
536   }
537   return fConvolutedGaussTransitionPars[index];
538 }
539
540
541 //_____________________________________________________________________________
542 inline Int_t AliAnalysisTaskPID::GetParticleFractionHistoNbinsTrackPt() const
543 {
544   if (!fFractionHists[AliPID::kPion])
545     return -1;
546   
547   return fFractionHists[AliPID::kPion]->GetNbinsX();
548 }
549
550
551 //_____________________________________________________________________________
552 inline Int_t AliAnalysisTaskPID::GetParticleFractionHistoNbinsJetPt() const
553 {
554   if (!fFractionHists[AliPID::kPion])
555     return -1;
556   
557   return fFractionHists[AliPID::kPion]->GetNbinsY();
558 }
559
560
561 //_____________________________________________________________________________
562 inline Int_t AliAnalysisTaskPID::GetParticleFractionHistoNbinsCentrality() const
563 {
564   if (!fFractionHists[AliPID::kPion])
565     return -1;
566   
567   return fFractionHists[AliPID::kPion]->GetNbinsZ();
568 }
569
570
571 //_____________________________________________________________________________
572 inline Double_t AliAnalysisTaskPID::GetCentralityPercentile(AliVEvent* evt) const
573 {
574   if (!evt)
575     return -1;
576   
577   Double_t centralityPercentile = -1.;
578   if (fStoreCentralityPercentile)
579     centralityPercentile = evt->GetCentrality()->GetCentralityPercentile(fCentralityEstimator.Data());
580   
581   return centralityPercentile;
582 }
583
584
585 //_____________________________________________________________________________
586 inline void AliAnalysisTaskPID::PostOutputData()
587 {
588   PostData(1, fOutputContainer);
589   
590   if (fDoEfficiency)
591     PostData(2, fContainerEff);
592   
593   if (fDoPtResolution)
594     PostData(3, fPtResolutionContainer);
595 }
596
597 #endif