AliPhysics  de71be2 (de71be2)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
AliAnalysisTaskEmcal.h
Go to the documentation of this file.
1 #ifndef ALIANALYSISTASKEMCAL_H
2 #define ALIANALYSISTASKEMCAL_H
3 
4 class TClonesArray;
5 class TString;
6 class TList;
7 class AliEmcalParticle;
8 class AliMCParticle;
9 class AliVCluster;
10 class AliVTrack;
11 class AliVParticle;
12 class AliVCaloCells;
13 class TH1;
14 class TProfile;
15 class AliEMCALGeometry;
16 class AliGenPythiaEventHeader;
17 class AliVCaloTrigger;
18 class AliAnalysisUtils;
19 class AliEMCALTriggerPatchInfo;
20 class AliAODTrack;
21 class AliEmcalPythiaInfo;
22 
23 #include "Rtypes.h"
24 
25 #include "AliParticleContainer.h"
26 #include "AliMCParticleContainer.h"
27 #include "AliTrackContainer.h"
28 #include "AliClusterContainer.h"
29 
30 #include "AliAnalysisTaskSE.h"
31 
32 class AliAnalysisTaskEmcal : public AliAnalysisTaskSE {
33  public:
34 
35  enum BeamType {
36  kNA = -1,
37  kpp = 0,
38  kAA = 1,
39  kpA = 2
40  };
41 
42  enum TriggerType {
43  kND = -1,
44  kJ1 = 0,
45  kJ2 = 1,
46  kG1 = 2,
47  kG2 = 3,
48  kL0 = 4
49  };
50 
51  enum TriggerCategory { // Online trigger categories
55  kTriggerRecalcJet = 3, // Recalculated max trigger patch; does not need to be above trigger threshold
57  };
58 
60  kNoSpecialTreatment, // No special treatment for EMCal triggers
61  kOverlapWithLowThreshold // The overlap between low and high threshold trigger is assigned to the lower threshold only
62  };
63 
65  AliAnalysisTaskEmcal(const char *name, Bool_t histo=kFALSE);
66  virtual ~AliAnalysisTaskEmcal();
67 
69  AliTrackContainer *AddTrackContainer(const char *n);
75  AliParticleContainer *GetParticleContainer(const char* name) const;
76  AliClusterContainer *GetClusterContainer(Int_t i=0) const;
77  AliClusterContainer *GetClusterContainer(const char* name) const;
79  AliMCParticleContainer *GetMCParticleContainer(const char* name) const { return dynamic_cast<AliMCParticleContainer*>(GetParticleContainer(name)); }
80  AliTrackContainer *GetTrackContainer(Int_t i=0) const { return dynamic_cast<AliTrackContainer*>(GetParticleContainer(i)) ; }
81  AliTrackContainer *GetTrackContainer(const char* name) const { return dynamic_cast<AliTrackContainer*>(GetParticleContainer(name)) ; }
82  void RemoveParticleContainer(Int_t i=0) { fParticleCollArray.RemoveAt(i) ; }
83  void RemoveClusterContainer(Int_t i=0) { fClusterCollArray.RemoveAt(i) ; }
84  void SetCaloCellsName(const char *n) { fCaloCellsName = n ; }
86  void SetCaloTriggersName(const char *n) { fCaloTriggersName = n ; }
87  void SetCentRange(Double_t min, Double_t max) { fMinCent = min ; fMaxCent = max ; }
88  void SetCentralityEstimator(const char *c) { fCentEst = c ; }
89  void SetClusName(const char *n) { AddClusterContainer(n) ; }
90  void SetClusPtCut(Double_t cut, Int_t c=0);
91  void SetClusTimeCut(Double_t min, Double_t max, Int_t c=0);
92  void SetEventPlaneVsEmcal(Double_t ep) { fEventPlaneVsEmcal = ep ; }
94  void SetHistoBins(Int_t nbins, Double_t min, Double_t max) { fNbins = nbins; fMinBinPt = min; fMaxBinPt = max ; }
95  void SetIsEmbedded(Bool_t i) { fIsEmbedded = i ; }
96  void SetIsPythia(Bool_t i) { fIsPythia = i ; }
98  void SetMCLabelShift(Int_t s) { fMCLabelShift = s ; }
99  void SetMinMCLabel(Int_t s) { fMinMCLabel = s ; }
100  void SetMinNTrack(Int_t min) { fMinNTrack = min ; }
101  void SetMinPtTrackInEmcal(Double_t min) { fMinPtTrackInEmcal = min ; }
102  virtual void SetNCentBins(Int_t n) { fNcentBins = n ; }
103  void SetNeedEmcalGeom(Bool_t n) { fNeedEmcalGeom = n ; }
104  void SetOffTrigger(UInt_t t) { fOffTrigger = t ; }
105  void SetTrackEtaLimits(Double_t min, Double_t max, Int_t c=0);
106  void SetTrackPhiLimits(Double_t min, Double_t max, Int_t c=0);
107  void SetTrackPtCut(Double_t cut, Int_t c=0);
108  void SetTracksName(const char *n) { AddParticleContainer(n) ; }
109  void SetTrigClass(const char *n) { fTrigClass = n ; }
111  void SetUseAliAnaUtils(Bool_t b, Bool_t bRejPilup = kTRUE) { fUseAliAnaUtils = b ; fRejectPileup = bRejPilup ; }
112  void SetVzRange(Double_t min, Double_t max) { fMinVz = min ; fMaxVz = max ; }
117  void SetPythiaInfoName(const char *n) { fPythiaInfoName = n ; }
118  const TString& GetPythiaInfoName() const { return fPythiaInfoName ; }
119  const AliEmcalPythiaInfo *GetPythiaInfo() const { return fPythiaInfo ; }
120 
121  protected:
122  void LoadPythiaInfo(AliVEvent *event);
123  void SetRejectionReasonLabels(TAxis* axis);
124  Bool_t AcceptCluster(AliVCluster *clus, Int_t c = 0) const;
125  Bool_t AcceptTrack(AliVParticle *track, Int_t c = 0) const;
126  void AddObjectToEvent(TObject *obj, Bool_t attempt = kFALSE);
127  AliVParticle *GetAcceptParticleFromArray(Int_t p, Int_t c=0) const;
128  AliVCluster *GetAcceptClusterFromArray(Int_t cl, Int_t c=0) const;
129  TClonesArray *GetArrayFromEvent(const char *name, const char *clname=0);
131  TClonesArray *GetParticleArray(Int_t i=0) const;
132  TClonesArray *GetClusterArray(Int_t i=0) const;
133  Int_t GetNParticles(Int_t i=0) const;
134  Int_t GetNClusters(Int_t i=0) const;
135  AliEMCALTriggerPatchInfo *GetMainTriggerPatch(TriggerCategory triggersel = kTriggerLevel1Jet, Bool_t doOfflinSimple = kFALSE);
136  Bool_t HasTriggerType(TriggerType triggersel);
137  ULong_t GetTriggerList();
138  Bool_t PythiaInfoFromFile(const char* currFile, Float_t &fXsec, Float_t &fTrials, Int_t &pthard);
139  Bool_t IsTrackInEmcalAcceptance(AliVParticle* part, Double_t edges=0.9) const;
140 
141  void GeneratePythiaInfoObject(AliMCEvent* mcEvent);
142 
143  // Overloaded AliAnalysisTaskSE methods
145  void UserExec(Option_t *option);
146  Bool_t UserNotify();
147 
148  // Virtual functions, to be overloaded in derived classes
149  virtual void ExecOnce();
150  virtual Bool_t FillGeneralHistograms();
151  virtual Bool_t IsEventSelected();
152  virtual Bool_t RetrieveEventObjects();
153  virtual Bool_t FillHistograms() { return kTRUE ; }
154  virtual Bool_t Run() { return kTRUE ; }
155 
156  // Static utilities
157  static void GetEtaPhiDiff(const AliVTrack *t, const AliVCluster *v, Double_t &phidiff, Double_t &etadiff);
158  static Byte_t GetTrackType(const AliVTrack *t);
159  static Byte_t GetTrackType(const AliAODTrack *aodTrack, UInt_t filterBit1, UInt_t filterBit2);
160  static Double_t DeltaPhi(Double_t phia, Double_t phib, Double_t rMin = -TMath::Pi()/2, Double_t rMax = 3*TMath::Pi()/2);
161  static Double_t* GenerateFixedBinArray(Int_t n, Double_t min, Double_t max);
162  static void GenerateFixedBinArray(Int_t n, Double_t min, Double_t max, Double_t* array);
163  static Double_t GetParallelFraction(AliVParticle* part1, AliVParticle* part2);
164  static Double_t GetParallelFraction(const TVector3& vect1, AliVParticle* part2);
165 
166  static Double_t fgkEMCalDCalPhiDivide; // phi value used to distinguish between DCal and EMCal
167 
168  // Task configuration
169  TString fPythiaInfoName; // name of pythia info object
170  BeamType fForceBeamType; // forced beam type
171  Bool_t fGeneralHistograms; // whether or not it should fill some general histograms
172  Bool_t fInitialized; // whether or not the task has been already initialized
173  Bool_t fCreateHisto; // whether or not create histograms
174  TString fCaloCellsName; // name of calo cell collection
175  TString fCaloTriggersName; // name of calo triggers collection
176  TString fCaloTriggerPatchInfoName; // trigger patch info array name
177  Double_t fMinCent; // min centrality for event selection
178  Double_t fMaxCent; // max centrality for event selection
179  Double_t fMinVz; // min vertex for event selection
180  Double_t fMaxVz; // max vertex for event selection
181  Double_t fTrackPtCut; // cut on track pt in event selection
182  Int_t fMinNTrack; // minimum nr of tracks in event with pT>fTrackPtCut
183  Bool_t fUseAliAnaUtils; // used for LHC13* data: z-vtx, Ncontributors, z-vtx resolution cuts
184  Bool_t fRejectPileup; // Reject pilup using function AliAnalysisUtils::IsPileUpEvent()
185  Bool_t fTklVsClusSPDCut; // Apply tracklet-vs-cluster SPD cut to reject background events in pp
186  UInt_t fOffTrigger; // offline trigger for event selection
187  TString fTrigClass; // trigger class name for event selection
188  TriggerType fTriggerTypeSel; // trigger type to select based on trigger patches
189  Int_t fNbins; // no. of pt bins
190  Double_t fMinBinPt; // min pt in histograms
191  Double_t fMaxBinPt; // max pt in histograms
192  Double_t fMinPtTrackInEmcal; // min pt track in emcal
193  Double_t fEventPlaneVsEmcal; // select events which have a certain event plane wrt the emcal
194  Double_t fMinEventPlane; // minimum event plane value
195  Double_t fMaxEventPlane; // maximum event plane value
196  TString fCentEst; // name of V0 centrality estimator
197  Bool_t fIsEmbedded; // trigger, embedded signal
198  Bool_t fIsPythia; // trigger, if it is a PYTHIA production
199  Int_t fSelectPtHardBin; // select one pt hard bin for analysis
200  Int_t fMinMCLabel; // minimum MC label value for the tracks/clusters being considered MC particles
201  Int_t fMCLabelShift; // if MC label > fMCLabelShift, MC label -= fMCLabelShift
202  Int_t fNcentBins; // how many centrality bins
203  Bool_t fNeedEmcalGeom; // whether or not the task needs the emcal geometry
204  TObjArray fParticleCollArray; // particle/track collection array
205  TObjArray fClusterCollArray; // cluster collection array
206  ULong_t fTriggers; // list of fired triggers
207  EMCalTriggerMode_t fEMCalTriggerMode; // EMCal trigger selection mode
208  Bool_t fUseNewCentralityEstimation; // Use new centrality estimation (for 2015 data)
209  Bool_t fGeneratePythiaInfoObject; // Generate Pythia info object
210 
211  // Service fields
212  AliAnalysisUtils *fAliAnalysisUtils;
213  Bool_t fIsEsd;
214  AliEMCALGeometry *fGeom;
215  TClonesArray *fTracks;
216  TClonesArray *fCaloClusters;
217  AliVCaloCells *fCaloCells;
218  AliVCaloTrigger *fCaloTriggers;
219  TClonesArray *fTriggerPatchInfo;
220  Double_t fCent;
221  Int_t fCentBin;
222  Double_t fEPV0;
223  Double_t fEPV0A;
224  Double_t fEPV0C;
225  Double_t fVertex[3];
226  Int_t fNVertCont;
228  AliGenPythiaEventHeader *fPythiaHeader;
229  Double_t fPtHard;
230  Int_t fPtHardBin;
231  Int_t fNTrials;
232  Float_t fXsection;
234 
235  // Output
236  TList *fOutput;
241  TH1 *fHistTrials;
242  TH1 *fHistEvents;
243  TProfile *fHistXsection;
244  TH1 *fHistPtHard;
250 
251  private:
252  AliAnalysisTaskEmcal(const AliAnalysisTaskEmcal&); // not implemented
253  AliAnalysisTaskEmcal &operator=(const AliAnalysisTaskEmcal&); // not implemented
254 
255  ClassDef(AliAnalysisTaskEmcal, 14) // EMCAL base analysis task
256 };
257 
258 //________________________________________________________________________
259 inline Double_t AliAnalysisTaskEmcal::DeltaPhi(Double_t phia, Double_t phib, Double_t rangeMin, Double_t rangeMax)
260 {
261  // Calculate Delta Phi.
262 
263  Double_t dphi = -999;
264  const Double_t tpi = TMath::TwoPi();
265 
266  if (phia < 0) phia += tpi;
267  else if (phia > tpi) phia -= tpi;
268  if (phib < 0) phib += tpi;
269  else if (phib > tpi) phib -= tpi;
270  dphi = phib - phia;
271  if (dphi < rangeMin) dphi += tpi;
272  else if (dphi > rangeMax) dphi -= tpi;
273 
274  return dphi;
275 }
276 
277 //________________________________________________________________________
278 inline void AliAnalysisTaskEmcal::GenerateFixedBinArray(Int_t n, Double_t min, Double_t max, Double_t* array)
279 {
280  Double_t binWidth = (max-min)/n;
281  array[0] = min;
282  for (Int_t i = 1; i <= n; i++) {
283  array[i] = array[i-1]+binWidth;
284  }
285 }
286 
287 //________________________________________________________________________
288 inline Double_t* AliAnalysisTaskEmcal::GenerateFixedBinArray(Int_t n, Double_t min, Double_t max)
289 {
290  Double_t *array = new Double_t[n+1];
291  GenerateFixedBinArray(n, min, max, array);
292  return array;
293 }
294 
295 #endif
void SetCentRange(Double_t min, Double_t max)
void SetPythiaInfoName(const char *n)
AliMCParticleContainer * GetMCParticleContainer(const char *name) const
TH1 * fHistTrials
x section from pythia header
static Double_t DeltaPhi(Double_t phia, Double_t phib, Double_t rMin=-TMath::Pi()/2, Double_t rMax=3 *TMath::Pi()/2)
AliEmcalPythiaInfo * fPythiaInfo
x-section from pythia header
Bool_t AcceptTrack(AliVParticle *track, Int_t c=0) const
Bool_t HasTriggerType(TriggerType triggersel)
Int_t fNTrials
event pt hard bin
void SetMinNTrack(Int_t min)
void AdoptParticleContainer(AliParticleContainer *cont)
void SetUseAliAnaUtils(Bool_t b, Bool_t bRejPilup=kTRUE)
void SetUseSPDTrackletVsClusterBG(Bool_t b)
Double_t fPtHard
event Pythia header
void SetTrackPtCut(Double_t cut, Int_t c=0)
static void GetEtaPhiDiff(const AliVTrack *t, const AliVCluster *v, Double_t &phidiff, Double_t &etadiff)
Double_t fEPV0
event centrality bin
void SetUseNewCentralityEstimation(Bool_t b)
Bool_t AcceptCluster(AliVCluster *clus, Int_t c=0) const
AliAnalysisTaskEmcal & operator=(const AliAnalysisTaskEmcal &)
Int_t fCentBin
event centrality
TH1 * fHistEventsAfterSel
total number of trials per pt hard bin after selection
void SetTriggerTypeSel(TriggerType t)
void SetVzRange(Double_t min, Double_t max)
TH1 * fHistEventPlane
z vertex position
TList * fOutput
event parton info
void SetTracksName(const char *n)
AliTrackContainer * GetTrackContainer(const char *name) const
TH1 * fHistEvents
trials from pyxsec.root
void SetClusPtCut(Double_t cut, Int_t c=0)
AliClusterContainer * AddClusterContainer(const char *n)
Double_t fEPV0C
event plane V0A
TH1 * fHistCentrality
pt hard distribution
void SetCaloTriggerPatchInfoName(const char *n)
void AdoptClusterContainer(AliClusterContainer *cont)
void SetTrackEtaLimits(Double_t min, Double_t max, Int_t c=0)
AliParticleContainer * GetParticleContainer(Int_t i=0) const
void SetEMCalTriggerMode(EMCalTriggerMode_t m)
TProfile * fHistXsectionAfterSel
total number of events per pt hard bin after selection
EMCalTriggerMode_t fEMCalTriggerMode
virtual Bool_t FillHistograms()
void SetForceBeamType(BeamType f)
Int_t GetNParticles(Int_t i=0) const
TClonesArray * fCaloClusters
tracks
Bool_t IsTrackInEmcalAcceptance(AliVParticle *part, Double_t edges=0.9) const
TH1 * fHistTriggerClasses
book keep reasons for rejecting event
void GeneratePythiaInfoObject(AliMCEvent *mcEvent)
AliEMCALGeometry * fGeom
whether it's an ESD analysis
AliGenPythiaEventHeader * fPythiaHeader
event beam type
void SetTrackPhiLimits(Double_t min, Double_t max, Int_t c=0)
AliParticleContainer * AddParticleContainer(const char *n)
AliAnalysisUtils * fAliAnalysisUtils
AliClusterContainer * GetClusterContainer(Int_t i=0) const
virtual Bool_t FillGeneralHistograms()
TClonesArray * GetParticleArray(Int_t i=0) const
BeamType fBeamType
event vertex number of contributors
Double_t fCent
trigger patch info array
Int_t GetNClusters(Int_t i=0) const
Int_t fNVertCont
event vertex
void SetClusName(const char *n)
void RemoveClusterContainer(Int_t i=0)
AliMCParticleContainer * AddMCParticleContainer(const char *n)
static Double_t GetParallelFraction(AliVParticle *part1, AliVParticle *part2)
virtual Bool_t RetrieveEventObjects()
void RemoveParticleContainer(Int_t i=0)
TProfile * fHistXsection
total number of events per pt hard bin
Bool_t PythiaInfoFromFile(const char *currFile, Float_t &fXsec, Float_t &fTrials, Int_t &pthard)
void UserExec(Option_t *option)
virtual void SetNCentBins(Int_t n)
AliVCaloCells * fCaloCells
clusters
void SetTrigClass(const char *n)
void SetIsEmbedded(Bool_t i)
TClonesArray * GetArrayFromEvent(const char *name, const char *clname=0)
virtual Bool_t IsEventSelected()
TH1 * fHistPtHard
x section from pyxsec.root
static Double_t * GenerateFixedBinArray(Int_t n, Double_t min, Double_t max)
Int_t fPtHardBin
event pt hard
void SetHistoBins(Int_t nbins, Double_t min, Double_t max)
AliMCParticleContainer * GetMCParticleContainer(Int_t i=0) const
TClonesArray * fTracks
emcal geometry
TH1 * fHistTrialsAfterSel
incoming and selected events
void LoadPythiaInfo(AliVEvent *event)
Bool_t fIsEsd
vertex selection (optional)
Double_t fVertex[3]
event plane V0C
AliTrackContainer * AddTrackContainer(const char *n)
AliTrackContainer * GetTrackContainer(Int_t i=0) const
Store some informaion about a Pythia eventThis class is used to store some information about a Pythia...
TH1 * fHistEventRejection
event plane distribution
void SetMakeGeneralHistograms(Bool_t g)
TClonesArray * fTriggerPatchInfo
calo triggers
TClonesArray * GetClusterArray(Int_t i=0) const
void SetNeedEmcalGeom(Bool_t n)
Double_t fEPV0A
event plane V0
void SetEventPlaneVsEmcal(Double_t ep)
void SetCaloTriggersName(const char *n)
const AliEmcalPythiaInfo * GetPythiaInfo() const
void AddObjectToEvent(TObject *obj, Bool_t attempt=kFALSE)
AliVCaloTrigger * fCaloTriggers
cells
void SetRejectionReasonLabels(TAxis *axis)
TH1 * fHistZVertex
event centrality distribution
static Byte_t GetTrackType(const AliVTrack *t)
void SetGeneratePythiaInfoObject(Bool_t b)
const Int_t nbins
void SetClusTimeCut(Double_t min, Double_t max, Int_t c=0)
Float_t fXsection
event trials
void SetCaloCellsName(const char *n)
void SetCentralityEstimator(const char *c)
TH1 * fHistEventCount
output list
void SetMinPtTrackInEmcal(Double_t min)
AliVParticle * GetAcceptParticleFromArray(Int_t p, Int_t c=0) const
AliVCluster * GetAcceptClusterFromArray(Int_t cl, Int_t c=0) const
const TString & GetPythiaInfoName() const
AliEMCALTriggerPatchInfo * GetMainTriggerPatch(TriggerCategory triggersel=kTriggerLevel1Jet, Bool_t doOfflinSimple=kFALSE)
static Double_t fgkEMCalDCalPhiDivide