]>
Commit | Line | Data |
---|---|---|
1c5acb87 | 1 | #ifndef ALIANAPARTICLEJETLEADINGCONECORRELATION_H |
2 | #define ALIANAPARTICLEJETLEADINGCONECORRELATION_H | |
3 | /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
4 | * See cxx source for full Copyright notice */ | |
1c5acb87 | 5 | |
6 | //_________________________________________________________________________ | |
7 | // Class that contains the algorithm for the reconstruction of jet, cone around leading particle | |
8 | // The seed is a backward particle (direct photon) | |
78219bac | 9 | // 1) Take the a trigger particle found stored in AliAODPWG4ParticleCorrelation, |
1c5acb87 | 10 | // 2) Search for the highest pt leading particle opposite to the trigger within a phi, pt window |
11 | // 3) Take all particles around leading in a cone R with pt larger than threshold and construct the jet | |
12 | // | |
13 | // Class created from old AliPHOSGammaJet | |
14 | // (see AliRoot versions previous Release 4-09) | |
15 | // | |
16 | //-- Author: Gustavo Conesa (INFN-LNF) | |
17 | ||
18 | // --- ROOT system --- | |
19 | class TH2F; | |
20 | ||
21 | //---- Analysis system ---- | |
22 | class AliAODTrack; | |
0ae57829 | 23 | class AliVCluster; |
1c5acb87 | 24 | class AliCaloTrackReader; |
25 | class AliNeutralMesonSelection; | |
1c5acb87 | 26 | |
745913ae | 27 | #include "AliAnaCaloTrackCorrBaseClass.h" |
1c5acb87 | 28 | |
745913ae | 29 | class AliAnaParticleJetLeadingConeCorrelation : public AliAnaCaloTrackCorrBaseClass { |
1c5acb87 | 30 | |
31 | public: | |
1c5acb87 | 32 | AliAnaParticleJetLeadingConeCorrelation() ; // default ctor |
78219bac | 33 | virtual ~AliAnaParticleJetLeadingConeCorrelation() ; //virtual dtor |
78219bac | 34 | |
1c5acb87 | 35 | TList * GetCreateOutputObjects(); |
36 | ||
37 | void InitParameters(); | |
38 | ||
39 | void Print(const Option_t * opt) const; | |
40 | ||
41 | Bool_t AreJetsRecalculated() const {return fReMakeJet ; } | |
42 | void SwitchOnJetsRecalculation(){fReMakeJet = kTRUE; } | |
43 | void SwitchOffJetsRecalculation(){fReMakeJet = kFALSE; } | |
44 | ||
45 | Bool_t AreJetsOnlyInCTS() const {return fJetsOnlyInCTS ; } | |
46 | void SwitchOnJetsOnlyInCTS(){fJetsOnlyInCTS = kTRUE; } | |
47 | void SwitchOffJetsOnlyInCTS(){fJetsOnlyInCTS = kFALSE; } | |
48 | ||
49 | Bool_t AreSeveralConeAndPtCuts() const {return fSeveralConeAndPtCuts ; } | |
50 | void SwitchOnSeveralConeAndPtCuts(){fSeveralConeAndPtCuts = kTRUE ;} | |
51 | void SwitchOffSeveralConeAndPtCuts(){fSeveralConeAndPtCuts = kFALSE ;} | |
52 | ||
53 | Bool_t IsPbPb() const {return fPbPb ; } | |
54 | void SetppCollisions(){fPbPb = kFALSE; } | |
55 | void SetPbPbCollisions(){fPbPb = kTRUE; } | |
56 | ||
57 | Double_t GetDeltaPhiMaxCut() const {return fDeltaPhiMaxCut ; } | |
58 | Double_t GetDeltaPhiMinCut() const {return fDeltaPhiMinCut ; } | |
59 | Double_t GetLeadingRatioMaxCut() const {return fLeadingRatioMaxCut ; } | |
60 | Double_t GetLeadingRatioMinCut() const {return fLeadingRatioMinCut ; } | |
61 | ||
62 | Double_t GetPtTriggerSelectionCut() const {return fPtTriggerSelectionCut ; } | |
63 | Double_t GetJetRatioMaxCut() const {return fJetRatioMaxCut ; } | |
64 | Double_t GetJetRatioMinCut() const {return fJetRatioMinCut ; } | |
65 | ||
66 | void SetPtTriggerSelectionCut(Double_t cut){fPtTriggerSelectionCut = cut; } | |
67 | void SetJetSelectionMode(UInt_t select){ fSelect= select ; } | |
68 | ||
69 | Int_t GetJetNCones() const {return fJetNCone ; } | |
70 | Int_t GetJetNPtThres() const {return fJetNPt ; } | |
71 | Float_t GetJetCone() const {return fJetCone ; } | |
72 | Float_t GetJetPtThreshold() const {return fJetPtThreshold ; } | |
73 | Float_t GetJetPtThresPbPb() const {return fJetPtThresPbPb ; } | |
74 | Float_t GetJetCones(Int_t i) const {return fJetCones[i] ; } | |
75 | Float_t GetJetPtThreshold(Int_t i) const {return fJetPtThres[i] ; } | |
76 | TString GetJetConeName(Int_t i) const {return fJetNameCones[i] ; } | |
77 | TString GetJetPtThresName(Int_t i) const {return fJetNamePtThres[i] ; } | |
78 | ||
79 | ||
80 | void SetDeltaPhiCutRange(Double_t phimin, Double_t phimax) | |
81 | {fDeltaPhiMaxCut =phimax; fDeltaPhiMinCut =phimin;} | |
82 | void SetLeadingRatioCutRange(Double_t ratiomin, Double_t ratiomax) | |
83 | {fLeadingRatioMaxCut =ratiomax; fLeadingRatioMinCut = ratiomin ; } | |
84 | ||
85 | void SetJetNCones(Int_t n){fJetNCone = n ; } | |
86 | void SetJetNPtThresholds(Int_t n){fJetNPt = n ; } | |
87 | void SetJetCones(Int_t i, Float_t cone, TString sc) {fJetCones[i] = cone ; fJetNameCones[i] = sc; }; | |
88 | void SetCone(Float_t cone) {fJetCone = cone; } | |
89 | void SetJetPtThreshold(Float_t pt){fJetPtThreshold = pt; }; | |
90 | void SetJetPtThresPbPb(Float_t pt){fJetPtThresPbPb = pt; }; | |
91 | void SetJetPtThresholds(Int_t i,Float_t pt, TString spt){fJetPtThres[i] = pt ; fJetNamePtThres[i] = spt; }; | |
92 | ||
93 | void SetJetRatioCutRange(Double_t ratiomin, Double_t ratiomax) | |
94 | {fJetRatioMaxCut =ratiomax; fJetRatioMinCut = ratiomin ; } | |
95 | void SetJetCTSRatioCutRange(Double_t ratiomin, Double_t ratiomax) | |
96 | {fJetCTSRatioMaxCut =ratiomax; fJetCTSRatioMinCut = ratiomin ; } | |
97 | ||
477d6cee | 98 | Bool_t OnlyIsolated() const {return fSelectIsolated ; } |
99 | void SelectIsolated(Bool_t select) {fSelectIsolated = select ; } | |
100 | ||
1c5acb87 | 101 | private: |
102 | ||
103 | Double_t CalculateJetRatioLimit(const Double_t ptTrig, const Double_t *param, const Double_t *x) const ; | |
104 | ||
105 | void FillJetHistos(AliAODPWG4ParticleCorrelation * particle, const TLorentzVector leading, const TLorentzVector jet, const TString type, const TString lastname); | |
106 | ||
107 | TList * GetOutputContainer() const {return fOutCont; } | |
108 | ||
109 | Bool_t IsJetSelected(const Double_t ptTrig, const Double_t ptjet) const ; | |
110 | Bool_t IsParticleInJetCone(const Double_t eta, Double_t phi, const Double_t etal, Double_t phil) const ; | |
111 | ||
233e0df8 | 112 | void GetLeadingCharge(AliAODPWG4ParticleCorrelation* const particle, TLorentzVector & pLeading) const ; |
decca433 | 113 | void GetLeadingPi0 (AliAODPWG4ParticleCorrelation* const particle, TLorentzVector & pLeading) ; |
114 | Bool_t GetLeadingParticle(AliAODPWG4ParticleCorrelation *particle, TLorentzVector & pLeading) ; | |
1c5acb87 | 115 | |
116 | void MakeAnalysisFillAOD(); | |
117 | void MakeAnalysisFillHistograms(); | |
67a9899b | 118 | void MakeAODJet(AliAODPWG4ParticleCorrelation * particle, const TLorentzVector pLeading) ; |
1c5acb87 | 119 | void MakeJetFromAOD(AliAODPWG4ParticleCorrelation * particle, const TLorentzVector pLeading, |
120 | TLorentzVector & jet, TLorentzVector & bkg) const ; | |
121 | ||
0ae57829 | 122 | Bool_t SelectCluster(AliVCluster * calo, Double_t *vertex, TLorentzVector & mom, Int_t & pdg) ; |
1c5acb87 | 123 | |
124 | private: | |
125 | ||
126 | Bool_t fJetsOnlyInCTS ; // Jets measured only in TPC+ITS. | |
127 | Bool_t fPbPb; // PbPb event | |
128 | Bool_t fSeveralConeAndPtCuts; // To play with the jet cone size and pt th. | |
129 | Bool_t fReMakeJet ; //Re make the jet reconstruction from AODParticleCorrelation input | |
130 | ||
131 | //Leading particle selection parameters | |
132 | Double_t fDeltaPhiMaxCut ; // Minimum Delta Phi Gamma-Leading | |
133 | Double_t fDeltaPhiMinCut ; // Maximum Delta Phi Gamma-Leading | |
134 | Double_t fLeadingRatioMaxCut ; // Leading /gamma Ratio cut maximum | |
135 | Double_t fLeadingRatioMinCut ; // Leading/gamma Ratio cut minimum | |
136 | ||
137 | //Jet selection parameters | |
138 | //Fixed cuts (old) | |
139 | Double_t fJetCTSRatioMaxCut ; // Jet(CTS) /gamma Ratio cut maximum | |
140 | Double_t fJetCTSRatioMinCut ; // Jet(CTS) /gamma Ratio cut maximum | |
141 | Double_t fJetRatioMaxCut ; // Jet(EMCAL+CTS)/gamma Ratio cut maximum | |
142 | Double_t fJetRatioMinCut ; // Jet(EMCAL+CTS)/gamma Ratio cut minimum | |
143 | ||
144 | //Cuts depending on jet pt | |
145 | Double_t fJetE1[2]; //Rec. jet energy parameters | |
146 | Double_t fJetE2[2]; //Rec. jet energy parameters | |
147 | Double_t fJetSigma1[2];//Rec. sigma of jet energy parameters | |
148 | Double_t fJetSigma2[2];//Rec. sigma of jet energy parameters | |
149 | Double_t fBkgMean[6]; //Background mean energy | |
150 | Double_t fBkgRMS[6]; //Background RMS | |
151 | Double_t fJetXMin1[6]; //X Factor to set jet min limit for pp | |
152 | Double_t fJetXMin2[6]; //X Factor to set jet min limit for PbPb | |
153 | Double_t fJetXMax1[6]; //X Factor to set jet max limit for pp | |
154 | Double_t fJetXMax2[6]; //X Factor to set jet max limit for PbPb | |
155 | ||
156 | Int_t fJetNCone ; // Number of jet cones sizes, maximum 5 | |
157 | Int_t fJetNPt ; // Number of jet particle pT threshold, maximum 5 | |
158 | Double_t fJetCone ; // Jet cone sizes under study (!fSeveralConeAndPtCuts) | |
159 | Double_t fJetCones[5]; // Jet cone sizes under study (fSeveralConeAndPtCuts) | |
160 | TString fJetNameCones[5]; // String name of cone to append to histos | |
161 | Double_t fJetPtThreshold; // Jet pT threshold under study(!fSeveralConeAndPtCuts) | |
162 | Double_t fJetPtThresPbPb; // Jet pT threshold under study(!fSeveralConeAndPtCuts) | |
163 | Double_t fJetPtThres[5]; // Jet pT threshold under study(fSeveralConeAndPtCuts) | |
164 | TString fJetNamePtThres[5]; // String name of pt th to append to histos | |
165 | Double_t fPtTriggerSelectionCut; // Jet pt to change to low pt jets analysis | |
166 | UInt_t fSelect ; //kTRUE: Selects all jets, no limits. | |
477d6cee | 167 | Bool_t fSelectIsolated ; // Select only trigger particles isolated |
168 | ||
1c5acb87 | 169 | //Histograms |
170 | //Leading particle distributions | |
171 | TList * fOutCont ; //! Container for histograms | |
172 | ||
173 | TH2F * fhChargedLeadingPt ; //! Pt(Pt trigger) distribution of charged hadrons | |
174 | TH2F * fhChargedLeadingPhi ; //! Phi(Pt trigger) distribution of charged hadrons | |
175 | TH2F * fhChargedLeadingEta ; //! Eta(Pt trigger) distribution of charged hadrons | |
176 | TH2F * fhChargedLeadingDeltaPt ; //! Difference of charged hadron and trigger pT as function of trigger p | |
177 | TH2F * fhChargedLeadingDeltaPhi ; //! Difference of charged hadron and trigger phi as function of trigger pT | |
178 | TH2F * fhChargedLeadingDeltaEta ; //! Difference of charged particle and trigger eta as function of trigger pT | |
179 | TH2F * fhChargedLeadingRatioPt ; //! Ratio of Pt leading charge and trigger | |
180 | ||
181 | TH2F * fhNeutralLeadingPt ; //! Pt(Pt trigger) distribution of neutral hadrons | |
182 | TH2F * fhNeutralLeadingPhi ; //! Phi(Pt trigger) distribution of neutral hadrons | |
183 | TH2F * fhNeutralLeadingEta ; //! Eta(Pt trigger) distribution of neutral hadrons | |
184 | TH2F * fhNeutralLeadingDeltaPt ; //! Difference of neutral hadron and trigger pT as function of trigger pT | |
185 | TH2F * fhNeutralLeadingDeltaPhi ; //! Difference of neutral hadron and trigger phi as function of trigger pT | |
186 | TH2F * fhNeutralLeadingDeltaEta ; //! Difference of charged particle and trigger eta as function of trigger pT | |
187 | TH2F * fhNeutralLeadingRatioPt ; //! Ratio of Pt leading neutral and trigger | |
9415d854 | 188 | |
189 | TH2F * fhChargedLeadingXi ; //! Ln (pt leading charge / pt trigger) | |
190 | TH2F * fhNeutralLeadingXi ; //! Ln (pt leading neutral / pt trigger) | |
dde5a268 | 191 | |
192 | TH2F * fhChargedLeadingDeltaPhiRatioPt30 ; //! Difference of charged hadron and trigger phi as function of pT leading / trigger pT, pT Trigger > 30 GeV | |
193 | TH2F * fhNeutralLeadingDeltaPhiRatioPt30 ; //! Difference of neutral hadron and trigger phi as function of pT leading / trigger pT, pT Trigger > 30 GeV | |
194 | TH2F * fhChargedLeadingDeltaPhiRatioPt50 ; //! Difference of charged hadron and trigger phi as function of pT leading / trigger pT, pT Trigger > 50 GeV | |
195 | TH2F * fhNeutralLeadingDeltaPhiRatioPt50 ; //! Difference of neutral hadron and trigger phi as function of pT leading / trigger pT, pT Trigger > 50 GeV | |
9415d854 | 196 | |
1c5acb87 | 197 | // Jet distributions |
198 | // Fixed cone and pt threshold | |
199 | TH2F * fhJetPt ; //! leading pt jet vs pt trigger | |
200 | TH2F * fhJetRatioPt ; //! Ratio of pt jet and pt trigger | |
201 | TH2F * fhJetDeltaPhi ; //! Delta phi jet-trigger | |
202 | TH2F * fhJetDeltaEta ; //! Delta eta jet-trigger | |
203 | TH2F * fhJetLeadingRatioPt ; //! Ratio of pt leading and pt jet | |
204 | TH2F * fhJetLeadingDeltaPhi ; //! Delta phi jet-leading | |
205 | TH2F * fhJetLeadingDeltaEta ; //! Delta eta jet-leading | |
206 | TH2F * fhJetFFz; //! Accepted reconstructed jet fragmentation function, z=ptjet/pttrig | |
207 | TH2F * fhJetFFxi; //! Accepted reconstructed jet fragmentation function, xsi = ln(pttrig/ptjet) | |
208 | TH2F * fhJetFFpt; //! Jet particle pt distribution in cone | |
209 | TH2F * fhJetNTracksInCone ; //! jet multiplicity in cone | |
210 | ||
211 | TH2F * fhBkgPt ; //! leading pt bakground vs pt trigger | |
212 | TH2F * fhBkgRatioPt ; //! Ratio of pt background and pt trigger | |
213 | TH2F * fhBkgDeltaPhi ; //! Delta phi background-trigger | |
214 | TH2F * fhBkgDeltaEta ; //! Delta eta background-trigger | |
215 | TH2F * fhBkgLeadingRatioPt ; //! Ratio of pt leading and pt background | |
216 | TH2F * fhBkgLeadingDeltaPhi ; //! Delta phi background-leading | |
217 | TH2F * fhBkgLeadingDeltaEta ; //! Delta eta background-leading | |
218 | TH2F * fhBkgFFz; //! Accepted reconstructed background fragmentation function, z=ptjet/pttrig | |
219 | TH2F * fhBkgFFxi; //! Accepted reconstructed background fragmentation function, xsi = ln(pttrig/ptjet) | |
220 | TH2F * fhBkgFFpt; //! Background particle pt distribution in cone | |
221 | TH2F * fhBkgNTracksInCone ; //! Background multiplicity in cone | |
222 | ||
223 | // Variable cone and pt threshold | |
224 | ||
225 | TH2F * fhJetPts[5][5]; //! leading pt jet vs pt trigger | |
226 | TH2F * fhJetRatioPts[5][5]; //! Ratio of pt jet and pt trigger | |
227 | TH2F * fhJetDeltaPhis[5][5]; //! Delta phi jet-trigger | |
228 | TH2F * fhJetDeltaEtas[5][5]; //! Delta eta jet-trigger | |
229 | TH2F * fhJetLeadingRatioPts[5][5]; //! Ratio of pt leading and pt jet | |
230 | TH2F * fhJetLeadingDeltaPhis[5][5]; //! Delta phi jet-leading | |
231 | TH2F * fhJetLeadingDeltaEtas[5][5]; //! Delta eta jet-leading | |
232 | TH2F * fhJetFFzs[5][5]; //! Accepted reconstructed jet fragmentation function, z=ptjet/pttrig | |
233 | TH2F * fhJetFFxis[5][5]; //! Accepted reconstructed jet fragmentation function, xsi = ln(pttrig/ptjet) | |
234 | TH2F * fhJetFFpts[5][5]; //! Jet particle pt distribution in cone | |
235 | TH2F * fhJetNTracksInCones[5][5]; //! jet multiplicity in cone | |
236 | ||
237 | TH2F * fhBkgPts[5][5]; //! leading pt bakground vs pt trigger | |
238 | TH2F * fhBkgRatioPts[5][5]; //! Ratio of pt background and pt trigger | |
239 | TH2F * fhBkgDeltaPhis[5][5]; //! Delta phi background-trigger | |
240 | TH2F * fhBkgDeltaEtas[5][5]; //! Delta eta background-trigger | |
241 | TH2F * fhBkgLeadingRatioPts[5][5]; //! Ratio of pt leading and pt background | |
242 | TH2F * fhBkgLeadingDeltaPhis[5][5]; //! Delta phi background-leading | |
243 | TH2F * fhBkgLeadingDeltaEtas[5][5]; //! Delta eta background-leading | |
244 | TH2F * fhBkgFFzs[5][5]; //! Accepted reconstructed background fragmentation function, z=ptjet/pttrig | |
245 | TH2F * fhBkgFFxis[5][5]; //! Accepted reconstructed background fragmentation function, xsi = ln(pttrig/ptjet) | |
246 | TH2F * fhBkgFFpts[5][5]; //! Background particle pt distribution in cone | |
247 | TH2F * fhBkgNTracksInCones[5][5]; //! Background multiplicity in cone | |
248 | ||
c5693f62 | 249 | AliAnaParticleJetLeadingConeCorrelation(const AliAnaParticleJetLeadingConeCorrelation & g) ; // cpy ctor |
250 | AliAnaParticleJetLeadingConeCorrelation & operator = (const AliAnaParticleJetLeadingConeCorrelation & g) ;//cpy assignment | |
1c5acb87 | 251 | |
252 | ClassDef(AliAnaParticleJetLeadingConeCorrelation,1) | |
253 | } ; | |
254 | ||
255 | ||
256 | #endif //ALIANAPARTICLEJETLEADINGCONECORRELATION_H | |
257 | ||
258 | ||
259 |