]>
Commit | Line | Data |
---|---|---|
62567a80 | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | ||
16 | //--------------------------------------------------------------------- | |
17 | // JetAnalysis class | |
83a444b1 | 18 | // Analyse Jets (already found jets) |
19 | // Authors: andreas.morsch@cern.ch, jgcn@mail.cern.ch | |
20 | // mercedes.lopez.noriega@cern.ch | |
62567a80 | 21 | //--------------------------------------------------------------------- |
22 | ||
83a444b1 | 23 | #include <Riostream.h> |
62567a80 | 24 | #include "AliJetAnalysis.h" |
25 | ClassImp(AliJetAnalysis) | |
83a444b1 | 26 | |
27 | // root | |
28 | #include <Riostream.h> | |
62567a80 | 29 | #include <TH1F.h> |
83a444b1 | 30 | #include <TH2F.h> |
62567a80 | 31 | #include <TProfile.h> |
62567a80 | 32 | #include <TFile.h> |
33 | #include <TTree.h> | |
34 | #include <TStyle.h> | |
35 | #include <TSystem.h> | |
36 | #include <TLorentzVector.h> | |
83a444b1 | 37 | #include <TMath.h> |
38 | // aliroot | |
62567a80 | 39 | #include "AliJetProductionDataPDC2004.h" |
40 | #include "AliJet.h" | |
62567a80 | 41 | #include "AliUA1JetHeader.h" |
42 | #include "AliLeading.h" | |
7aec0427 | 43 | #include "AliJetReaderHeader.h" |
83a444b1 | 44 | |
1b7d5d7e | 45 | AliJetAnalysis::AliJetAnalysis(): |
46 | fReaderHeader(0x0), | |
47 | fDirectory(0x0), | |
48 | fBkgdDirectory(0x0), | |
a6e0ebfe | 49 | fFile("anaJets.root"), |
1b7d5d7e | 50 | fEventMin(0), |
51 | fEventMax(-1), | |
52 | fRunMin(0), | |
53 | fRunMax(11), | |
54 | fminMult(0), | |
55 | fPercentage(-1.0), | |
56 | fPartPtCut(0.0), | |
57 | fdrJt(1.0), | |
58 | fdrdNdxi(0.7), | |
59 | fdrdEdr(1.0), | |
60 | fEfactor(1.0), | |
61 | fp0(0.0), | |
62 | fPtJ(0.0), | |
63 | fEJ(0.0), | |
64 | fEtaJ(0.0), | |
65 | fPhiJ(0.0), | |
66 | fjv3X(0.0), | |
67 | fjv3Y(0.0), | |
68 | fjv3Z(0.0), | |
69 | fPythia(kFALSE), | |
70 | fDoPart(kTRUE), | |
71 | fDoGenJ(kTRUE), | |
72 | fDoRecJ(kTRUE), | |
73 | fDoKine(kTRUE), | |
74 | fDoCorr(kTRUE), | |
75 | fDoCorr50(kFALSE), | |
76 | fDoShap(kTRUE), | |
77 | fDoFrag(kTRUE), | |
78 | fDoTrig(kTRUE), | |
79 | fDoJt(kTRUE), | |
80 | fDodNdxi(kTRUE), | |
81 | fDoBkgd(kTRUE), | |
82 | fWeight(1.0), | |
83 | fWShapR(0.0), | |
84 | fWFragR(0.0), | |
85 | fWdEdr(0.0), | |
86 | fWJt(0.0), | |
87 | fWdNdxi(0.0), | |
88 | fPart(0), | |
89 | fGenJ(0), | |
90 | fRecJ(0), | |
91 | fRecB(0), | |
92 | fRKineEneH(0), | |
93 | fRKinePtH(0), | |
94 | fRKinePhiH(0), | |
95 | fRKineEtaH(0), | |
96 | fGKineEneH(0), | |
97 | fGKinePtH(0), | |
98 | fGKinePhiH(0), | |
99 | fGKineEtaH(0), | |
100 | fPKineEneH(0), | |
101 | fPKinePtH(0), | |
102 | fPKinePhiH(0), | |
103 | fPKineEtaH(0), | |
104 | fPGCorrEneH(0), | |
105 | fPGCorrPtH(0), | |
106 | fPGCorrEtaH(0), | |
107 | fPGCorrPhiH(0), | |
108 | fPRCorrEneH(0), | |
109 | fPRCorrPtH(0), | |
110 | fPRCorrEtaH(0), | |
111 | fPRCorrPhiH(0), | |
112 | fRGCorrEneH(0), | |
113 | fRGCorrPtH(0), | |
114 | fRGCorrEtaH(0), | |
115 | fRGCorrPhiH(0), | |
116 | fPRCorr50EneH(0), | |
117 | fPRCorr50PtH(0), | |
118 | fPRCorr50EtaH(0), | |
119 | fPRCorr50PhiH(0), | |
120 | fRGCorr50EneH(0), | |
121 | fRGCorr50PtH(0), | |
122 | fRGCorr50EtaH(0), | |
123 | fRGCorr50PhiH(0), | |
124 | fRFragSelH(0), | |
125 | fRFragRejH(0), | |
126 | fRFragAllH(0), | |
127 | fRShapSelH(0), | |
128 | fRShapRejH(0), | |
129 | fRShapAllH(0), | |
130 | fGTriggerEneH(0), | |
131 | fRTriggerEneH(0), | |
132 | fGPTriggerEneH(0), | |
133 | fPTriggerEneH(0), | |
134 | fdEdrH(0), | |
135 | fdEdrB(0), | |
136 | fPtEneH2(0), | |
137 | fdEdrW(0), | |
138 | fJtH(0), | |
139 | fJtB(0), | |
140 | fJtW(0), | |
141 | fdNdxiH(0), | |
142 | fdNdxiB(0), | |
143 | fdNdxiW(0), | |
144 | fPtEneH(0) | |
83a444b1 | 145 | { |
1b7d5d7e | 146 | // Default constructor |
83a444b1 | 147 | |
83a444b1 | 148 | // initialize weight for dE/dr histo |
149 | SetdEdrWeight(); | |
150 | } | |
151 | ||
152 | AliJetAnalysis::~AliJetAnalysis() | |
153 | { | |
154 | // Destructor | |
155 | } | |
156 | ||
157 | ||
158 | //////////////////////////////////////////////////////////////////////// | |
159 | // define histogrames | |
160 | ||
161 | void AliJetAnalysis::DefineHistograms() | |
162 | { | |
163 | // Define the histograms to be filled | |
164 | if (fDoKine) DefineKineH(); | |
165 | if (fDoCorr) DefineCorrH(); | |
166 | if (fDoCorr50) DefineCorr50H(); | |
167 | if (fDoShap) DefineShapH(); | |
168 | if (fDoFrag) DefineFragH(); | |
169 | if (fDoTrig) DefineTrigH(); | |
170 | if (fDoJt) DefineJtH(); | |
171 | if (fDodNdxi) DefinedNdxiH(); | |
172 | } | |
173 | ||
174 | void AliJetAnalysis::DefineKineH() | |
175 | { | |
176 | // leading particle | |
177 | if (fDoPart) { | |
178 | fPKineEneH = new TH1F("PKineEne","Energy of leading particle",50,0,200); | |
179 | SetProperties(fPKineEneH,"Energy (GeV)","Entries"); | |
180 | fPKinePtH = new TH1F("PKinePt","Pt of leading particle",50,0,200); | |
181 | SetProperties(fPKinePtH,"P_{T} (GeV)","Entries"); | |
182 | fPKinePhiH = new TH1F("PKinePhiH","Azimuthal angle of leading particle", | |
183 | 90,0.,2.0*TMath::Pi()); | |
184 | SetProperties(fPKinePhiH,"#phi","Entries"); | |
185 | fPKineEtaH = new TH1F("PKineEtaH","Pseudorapidity of leading particle", | |
186 | 40,-1.0,1.0); | |
187 | SetProperties(fPKineEtaH,"#eta","Entries"); | |
188 | } | |
189 | // leading generated jet | |
190 | if (fDoGenJ) { | |
191 | fGKineEneH = new TH1F("GKineEne","Energy of generated jet",50,0,200); | |
192 | SetProperties(fGKineEneH,"Energy (GeV)","Entries"); | |
193 | fGKinePtH = new TH1F("GKinePt","Pt of generated jet",50,0,200); | |
194 | SetProperties(fGKinePtH,"P_{T} (GeV)","Entries"); | |
195 | fGKinePhiH = new TH1F("GKinePhiH","Azimuthal angle of generated jet", | |
196 | 90,0.,2.0*TMath::Pi()); | |
197 | SetProperties(fGKinePhiH,"#phi","Entries"); | |
198 | fGKineEtaH = new TH1F("GKineEtaH","Pseudorapidity of generated jet", | |
199 | 40,-1.0,1.0); | |
200 | SetProperties(fGKineEtaH,"#eta","Entries"); | |
201 | } | |
202 | // leading reconstructed jet | |
203 | if (fDoRecJ) { | |
204 | fRKineEneH = new TH1F("RKineEne","Energy of reconstructed jet",50,0,200); | |
205 | SetProperties(fRKineEneH,"Energy (GeV)","Entries"); | |
206 | fRKinePtH = new TH1F("RKinePt","Pt of reconstructed jet",50,0,200); | |
207 | SetProperties(fRKinePtH,"P_{T} (GeV)","Entries"); | |
208 | fRKinePhiH = new TH1F("RKinePhiH","Azimuthal angle of reconstructed jet", | |
209 | 90,0.,2.0*TMath::Pi()); | |
210 | SetProperties(fRKinePhiH,"#phi","Entries"); | |
211 | fRKineEtaH = new TH1F("RKineEtaH","Pseudorapidity of reconstructed jet", | |
212 | 40,-1.0,1.0); | |
213 | SetProperties(fRKineEtaH,"#eta","Entries"); | |
214 | } | |
215 | } | |
62567a80 | 216 | |
83a444b1 | 217 | void AliJetAnalysis::DefineCorrH() |
62567a80 | 218 | { |
83a444b1 | 219 | // correlation |
220 | if (fDoPart && fDoGenJ) { | |
221 | fPGCorrEneH = new TH2F("PGCorrEne","Energy correlation part-gen jet", | |
222 | 40,0,200,40,0,200); | |
223 | SetProperties(fPGCorrEneH,"Part Energy (GeV)","Gen Jet Energy (GeV)"); | |
224 | fPGCorrPtH = new TH2F("PGCorrPt","Pt correlation part-gen jet", | |
225 | 40,0,200,40,0,200); | |
226 | SetProperties(fPGCorrPtH,"Part P_{T} (GeV)","Gen Jet P_{T} (GeV)"); | |
227 | fPGCorrEtaH = new TH2F("PGCorrEta","Pseudorapidity correlation part-gen jet", | |
228 | 40,-1.0,1.0,40,-1.0,1.0); | |
229 | SetProperties(fPGCorrEtaH,"Part #eta","Gen Jet #eta"); | |
230 | fPGCorrPhiH = new TH2F("PGCorrPhi","Azimuthal angle correlation part-gen jet", | |
231 | 90,0.,2.0*TMath::Pi(),90,0.,2.0*TMath::Pi()); | |
232 | SetProperties(fPGCorrPhiH,"Part #phi","Gen Jet #phi"); | |
233 | } | |
234 | if (fDoPart && fDoRecJ) { | |
235 | fPRCorrEneH = new TH2F("PRCorrEne","Energy correlation part-rec jet", | |
236 | 40,0,200,40,0,200); | |
237 | SetProperties(fPRCorrEneH,"Part Jet Energy (GeV)","Rec Jet Energy (GeV)"); | |
238 | fPRCorrPtH = new TH2F("PRCorrPt","Pt correlation part-rec jet", | |
239 | 40,0,200,40,0,200); | |
240 | SetProperties(fPRCorrPtH,"Part Jet P_{T} (GeV)","Rec Jet P_{T} (GeV)"); | |
241 | fPRCorrEtaH = new TH2F("PRCorrEta","Pseudorapidity correlation part-rec jet", | |
242 | 40,-1.0,1.0,40,-1.0,1.0); | |
243 | SetProperties(fPRCorrEtaH,"part #eta","Rec Jet #eta"); | |
244 | fPRCorrPhiH = new TH2F("PRCorrPhi","Azimuthal angle correlation part-rec jet", | |
245 | 90,0.,2.0*TMath::Pi(),90,0.,2.0*TMath::Pi()); | |
246 | SetProperties(fPRCorrPhiH,"Part #phi","Rec Jet #phi"); | |
247 | } | |
248 | if (fDoGenJ && fDoRecJ) { | |
249 | fRGCorrEneH = new TH2F("RGCorrEne","Energy correlation rec jet-gen jet", | |
250 | 40,0,200,40,0,200); | |
251 | SetProperties(fRGCorrEneH,"Rec Jet Energy (GeV)","Gen Jet Energy (GeV)"); | |
252 | fRGCorrPtH = new TH2F("RGCorrPt","Pt correlation rec jet-gen jet", | |
253 | 40,0,200,40,0,200); | |
254 | SetProperties(fRGCorrPtH,"Rec Jet P_{T} (GeV)","Gen Jet P_{T} (GeV)"); | |
255 | fRGCorrEtaH = new TH2F("RGCorrEta","Pseudorapidity correlation rec jet-gen jet", | |
256 | 40,-1.0,1.0,40,-1.0,1.0); | |
257 | SetProperties(fRGCorrEtaH,"Rec Jet #eta","Gen Jet #eta"); | |
258 | fRGCorrPhiH = new TH2F("RGCorrPhi","Azimuthal angle correlation rec jet-gen jet", | |
259 | 90,0.,2.0*TMath::Pi(),90,0.,2.0*TMath::Pi()); | |
260 | SetProperties(fRGCorrPhiH,"Rec Jet #phi","Gen Jet #phi"); | |
261 | } | |
62567a80 | 262 | } |
263 | ||
83a444b1 | 264 | void AliJetAnalysis::DefineCorr50H() |
62567a80 | 265 | { |
83a444b1 | 266 | // correlation |
267 | if (fDoPart && fDoRecJ) { | |
268 | fPRCorr50EneH = new TH2F("PRCorr50Ene","Energy correlation part-rec jet", | |
269 | 40,0,200,40,0,200); | |
270 | SetProperties(fPRCorr50EneH,"Part Jet Energy (GeV)","Rec Jet Energy (GeV)"); | |
271 | fPRCorr50PtH = new TH2F("PRCorr50Pt","Pt correlation part-rec jet", | |
272 | 40,0,200,40,0,200); | |
273 | SetProperties(fPRCorr50PtH,"Part Jet P_{T} (GeV)","Rec Jet P_{T} (GeV)"); | |
274 | fPRCorr50EtaH = new TH2F("PRCorr50Eta","Pseudorapidity correlation part-rec jet", | |
275 | 40,-1.0,1.0,40,-1.0,1.0); | |
276 | SetProperties(fPRCorr50EtaH,"part #eta","Rec Jet #eta"); | |
277 | fPRCorr50PhiH = new TH2F("PRCorr50Phi","Azimuthal angle correlation part-rec jet", | |
278 | 90,0.,2.0*TMath::Pi(),90,0.,2.0*TMath::Pi()); | |
279 | SetProperties(fPRCorr50PhiH,"Part #phi","Rec Jet #phi"); | |
280 | } | |
281 | ||
282 | if (fDoGenJ && fDoRecJ) { | |
283 | fRGCorr50EneH = new TH2F("RGCorr50Ene","Energy correlation rec jet-gen jet", | |
284 | 40,0,200,40,0,200); | |
285 | SetProperties(fRGCorr50EneH,"Rec Jet Energy (GeV)","Gen Jet Energy (GeV)"); | |
286 | fRGCorr50PtH = new TH2F("RGCorr50Pt","Pt correlation rec jet-gen jet", | |
287 | 40,0,200,40,0,200); | |
288 | SetProperties(fRGCorr50PtH,"Rec Jet P_{T} (GeV)","Gen Jet P_{T} (GeV)"); | |
289 | fRGCorr50EtaH = new TH2F("RGCorr50Eta","Pseudorapidity correlation rec jet-gen jet", | |
290 | 40,-1.0,1.0,40,-1.0,1.0); | |
291 | SetProperties(fRGCorr50EtaH,"Rec Jet #eta","Gen Jet #eta"); | |
292 | fRGCorr50PhiH = new TH2F("RGCorr50Phi","Azimuthal angle correlation rec jet-gen jet", | |
293 | 90,0.,2.0*TMath::Pi(),90,0.,2.0*TMath::Pi()); | |
294 | SetProperties(fRGCorr50PhiH,"Rec Jet #phi","Gen Jet #phi"); | |
295 | } | |
296 | } | |
297 | ||
298 | void AliJetAnalysis::DefineShapH() | |
299 | { | |
300 | // leading reconstructed jet | |
301 | if (fDoRecJ) { | |
302 | fdEdrH = new TH2F("fdEdrH","dE/dr histo",20,0,1,40,0,200); | |
303 | SetProperties(fdEdrH,"r","Rec Jet P_{T}"); | |
304 | fdEdrB = new TH2F("fdEdrB","dE/dr bkgdhisto",20,0,1,40,0,200); | |
305 | SetProperties(fdEdrB,"r","Rec P_{T}"); | |
306 | fPtEneH2 = new TH2F("fPtEneH2","fPtEneH2",40,0,200,40,0,200); | |
307 | SetProperties(fPtEneH2,"Rec Jet E","Rec Jet P_{T}"); | |
308 | fdEdrW = new TH1F("fdEdrW","weights for dE/dr",40,0,200); | |
309 | SetProperties(fdEdrW,"Rec Jet P_{T}","weights"); | |
62567a80 | 310 | |
83a444b1 | 311 | fRShapSelH = new TH1F("RShapSel","Shape of generated jets (sel part)",20,0.,1.); |
312 | SetProperties(fRShapSelH,"r","1/N_{JET}#Sigma P_{T}(0,r)/P_{T}(0,R_{JET})"); | |
313 | fRShapRejH = new TH1F("RShapRej","Shape of generated jets (rej part)",20,0.,1.); | |
314 | SetProperties(fRShapRejH,"r","1/N_{JET}#Sigma P_{T}(0,r)/P_{T}(0,R_{JET})"); | |
315 | fRShapAllH = new TH1F("RShapAll","Shape of generated jets (all part)",20,0.,1.); | |
316 | SetProperties(fRShapAllH,"r","1/N_{JET}#Sigma P_{T}(0,r)/P_{T}(0,R_{JET})"); | |
317 | } | |
318 | } | |
62567a80 | 319 | |
83a444b1 | 320 | void AliJetAnalysis::DefineJtH() |
321 | { | |
322 | // Define the histogram for J_T | |
323 | if (fDoRecJ) { | |
324 | fJtH = new TH2F("fJtH","J_{T} histogram",80,0.,4.,40,0.,200.); | |
325 | SetProperties(fJtH,"J_{T}","Rec Jet P_{T}"); | |
326 | fJtB = new TH2F("fJtB","J_{T} bkgd histogram",80,0.,4.,40,0.,200.); | |
327 | SetProperties(fJtB,"J_{T}","Rec P_{T}"); | |
328 | fJtW = new TH1F("fJtW","J_{T} weight",40,0,200); | |
329 | SetProperties(fJtW,"J_{T}W","weight"); | |
330 | } | |
331 | } | |
62567a80 | 332 | |
83a444b1 | 333 | void AliJetAnalysis::DefinedNdxiH() |
334 | { | |
335 | // Define the histogram for dN/dxi | |
336 | if (fDoRecJ) { | |
337 | fdNdxiH = new TH2F("fdNdxiH","dN/d#xi histo",200,0,10,40,0,200); | |
338 | SetProperties(fdNdxiH,"xi","Rec Jet P_{T}"); | |
339 | fdNdxiB = new TH2F("fdNdxiB","dN/d#xi bkgd histo",200,0,10,40,0,200); | |
340 | SetProperties(fdNdxiB,"xi","Rec P_{T}"); | |
341 | fdNdxiW = new TH1F("fdNdxiW","dN/d#xi histo",40,0,200); | |
342 | SetProperties(fdNdxiW,"Rec Jet P_{T}","weights"); | |
343 | fPtEneH = new TH2F("fPtEneH","fPtEneH",40,0,200,40,0,200); | |
344 | SetProperties(fPtEneH,"Rec Jet E","Rec Jet P_{T}"); | |
345 | } | |
346 | } | |
62567a80 | 347 | |
83a444b1 | 348 | void AliJetAnalysis::DefineFragH() |
349 | { | |
350 | // leading reconstructed jet | |
351 | if (fDoRecJ) { | |
352 | fRFragSelH = new TH1F("RFragSel","Frag Fun of reconstructed jets (sel part)",20,0.,10.); | |
353 | SetProperties(fRFragSelH,"#xi=ln(E_{JET}/E_{i})","1/N_{JET}dN_{ch}/d#xi"); | |
354 | fRFragRejH = new TH1F("RFragRej","Frag Fun of reconstructed jets (rej part)",20,0.,10.); | |
355 | SetProperties(fRFragRejH,"#xi=ln(E_{JET}/E_{i})","1/N_{JET}dN_{ch}/d#xi"); | |
356 | fRFragAllH = new TH1F("RFragAll","Frag Fun of reconstructed jets (all part)",20,0.,10.); | |
357 | SetProperties(fRFragAllH,"#xi=ln(E_{JET}/E_{i})","1/N_{JET}dN_{ch}/d#xi"); | |
358 | } | |
359 | } | |
62567a80 | 360 | |
83a444b1 | 361 | void AliJetAnalysis::DefineTrigH() |
362 | { | |
363 | // generated energy | |
364 | fGTriggerEneH = new TProfile("GTriggerEne","Generated energy (trigger bias)", | |
365 | 20, 0., 200., 0., 1., "S"); | |
366 | fGTriggerEneH->SetXTitle("E_{gen}"); | |
367 | fGTriggerEneH->SetYTitle("E_{rec}/E_{gen}"); | |
368 | // reconstructed energy | |
369 | fRTriggerEneH = new TProfile("RTriggerEne","Reconstructed energy (trigger bias)", | |
370 | 20, 0., 200., 0., 1., "S"); | |
371 | fRTriggerEneH->SetXTitle("E_{rec}"); | |
372 | fRTriggerEneH->SetYTitle("E_{rec}/E_{gen}"); | |
373 | // generated energy vs generated/leading | |
374 | fGPTriggerEneH = new TProfile("GPTriggerEne","Generated energy (trigger bias)", | |
375 | 20, 0., 200., 0., 1., "S"); | |
376 | fGPTriggerEneH->SetXTitle("E_{gen}"); | |
377 | fGPTriggerEneH->SetYTitle("E_{L}/E_{gen}"); | |
378 | // leading particle energy | |
379 | fPTriggerEneH = new TProfile("PTriggerEne","Leading particle energy (trigger bias)", | |
380 | 20, 0., 200., 0., 1., "S"); | |
381 | fPTriggerEneH->SetXTitle("E_{L}/0.2"); | |
382 | fPTriggerEneH->SetYTitle("E_{L}/E_{gen}"); | |
62567a80 | 383 | |
83a444b1 | 384 | } |
62567a80 | 385 | |
83a444b1 | 386 | void AliJetAnalysis::SetProperties(TH1* h,const char* x, const char* y) const |
387 | { | |
388 | //Set properties of histos (x and y title and error propagation) | |
389 | h->SetXTitle(x); | |
390 | h->SetYTitle(y); | |
391 | h->Sumw2(); | |
392 | } | |
393 | ||
394 | //////////////////////////////////////////////////////////////////////// | |
395 | // compute weights for dE/dr histogram | |
62567a80 | 396 | |
83a444b1 | 397 | void AliJetAnalysis::SetdEdrWeight() |
398 | { | |
399 | // Due to the limited acceptance, each bin in the dE/dr has a different | |
400 | // weight given by the ratio of the area of a disk dR to the area | |
401 | // within the eta acceptance. The weight depends on the eta position | |
402 | // of the jet. Here a look up table for the weights is computed. It | |
403 | // assumes |etaJet|<0.5 and |eta_lego|<0.9. It makes bins of 0.05 | |
404 | // units in eta. Note that this is fixed by the bin width chosen for | |
405 | // the histogram --> this weights are tailored for the specific | |
406 | // histogram definition used here! | |
407 | ||
408 | // two dimensional table: first index, bin in eta of jet, second | |
409 | // index bin of dE/dr histo | |
410 | ||
411 | Int_t nBins = 20; | |
412 | Float_t xMin = 0.0; | |
413 | Float_t xMax = 1.0; | |
414 | Float_t binSize = (xMax-xMin)/nBins; | |
62567a80 | 415 | |
83a444b1 | 416 | Float_t da,ds,r1,r2,h1,h2,a1,a2,theta1,theta2,area1,area2; |
417 | Int_t ji; | |
418 | for (Int_t i=0;i<(nBins/2);i++) { | |
419 | // index of first histo bin needing a scale factor | |
420 | ji=(nBins-2)-i; | |
421 | for(Int_t j=0;j<nBins;j++) { | |
422 | // area of ring. | |
423 | da = TMath::Pi()*(binSize*binSize)*(1.0+2.0*j); | |
424 | ds = 1.0; | |
425 | if (j>=ji) { // compute missing area using segments of circle | |
426 | r1=j*binSize; | |
427 | r2=(j+1)*binSize; | |
428 | h1=(j-ji)*binSize; | |
429 | h2=(j+1-ji)*binSize; | |
430 | a1=2.0*TMath::Sqrt(2.0*h1*r1-h1*h1); | |
431 | a2=2.0*TMath::Sqrt(2.0*h2*r2-h2*h2); | |
432 | theta1=2*TMath::ACos((r1-h1)/r1); | |
433 | theta2=2*TMath::ACos((r2-h2)/r2); | |
434 | area1=binSize*(r1*r1*theta1-a1*(r1-h1)); | |
435 | area2=binSize*(r2*r2*theta2-a2*(r2-h2)); | |
436 | ds = (da-(area2-area1))/da; | |
437 | } | |
438 | fWeightdEdr[i][j]=ds/da; | |
439 | } | |
440 | } | |
441 | } | |
442 | ||
443 | // get weight for dE/dr histogram | |
444 | Float_t AliJetAnalysis::GetdEdrWeight(Float_t etaJet, Float_t r) | |
445 | { | |
446 | // Return the correponding weight for the dE/dr plot | |
447 | Int_t nBins = 20; | |
448 | Float_t xMin = 0.0; | |
449 | Float_t xMax = 1.0; | |
450 | Float_t binSize = (xMax-xMin)/nBins; | |
451 | ||
452 | Float_t eta = TMath::Abs(etaJet); | |
453 | if ((eta > 0.5) || (r > fdrdEdr)) return 0.0; | |
454 | Int_t binJet = (Int_t) TMath::Floor(eta/binSize); | |
455 | Int_t binR = (Int_t) TMath::Floor(r/binSize); | |
456 | Float_t w = fWeightdEdr[binJet][binR]; | |
457 | return w; | |
458 | } | |
459 | ||
460 | ||
461 | //////////////////////////////////////////////////////////////////////// | |
462 | // fill histograms | |
463 | ||
464 | void AliJetAnalysis::FillHistograms() | |
465 | { | |
466 | // fill histograms | |
62567a80 | 467 | |
468 | // Run data | |
83a444b1 | 469 | AliJetProductionDataPDC2004* runData = new AliJetProductionDataPDC2004(); |
470 | ||
471 | // Loop over runs | |
472 | TFile* jFile = 0x0; | |
473 | ||
474 | for (Int_t iRun = fRunMin; iRun <= fRunMax; iRun++) { | |
475 | // Open file | |
476 | char fn[256]; | |
477 | sprintf(fn, "%s/%s.root", fDirectory, (runData->GetRunTitle(iRun)).Data()); | |
478 | jFile = new TFile(fn); | |
479 | printf(" Analyzing run: %d %s\n", iRun,fn); | |
480 | ||
481 | // Get reader header and events to be looped over | |
482 | AliJetReaderHeader *jReaderH = (AliJetReaderHeader*)(jFile->Get(fReaderHeader)); | |
483 | if (fEventMin == -1) fEventMin = jReaderH->GetFirstEvent(); | |
484 | if (fEventMax == -1) { | |
485 | fEventMax = jReaderH->GetLastEvent(); | |
486 | } else { | |
487 | fEventMax = TMath::Min(fEventMax, jReaderH->GetLastEvent()); | |
488 | } | |
489 | ||
490 | AliUA1JetHeader *jH = (AliUA1JetHeader *) (jFile->Get("AliUA1JetHeader")); | |
62567a80 | 491 | |
83a444b1 | 492 | // Calculate weight |
493 | fWeight = runData->GetWeight(iRun)/ ( (Float_t) (fEventMax - fEventMin + 1)); | |
62567a80 | 494 | |
83a444b1 | 495 | // Loop over events |
496 | for (Int_t i = fEventMin; i < fEventMax; i++) { | |
497 | if (i%100 == 0) printf(" Analyzing event %d / %d \n",i,fEventMax); | |
498 | ||
499 | // Get next tree with AliJet | |
500 | char nameT[100]; | |
501 | sprintf(nameT, "TreeJ%d",i); | |
502 | TTree *jetT =(TTree *)(jFile->Get(nameT)); | |
503 | if (fDoRecJ) jetT->SetBranchAddress("FoundJet", &fRecJ); | |
504 | if (fDoGenJ) jetT->SetBranchAddress("GenJet", &fGenJ); | |
505 | if (fDoPart) jetT->SetBranchAddress("LeadingPart", &fPart); | |
506 | jetT->GetEntry(0); | |
507 | ||
508 | int nJets = fRecJ->GetNJets(); | |
509 | ||
510 | TArrayI inJet = fRecJ->GetInJet(); | |
511 | if(inJet.GetSize()>fminMult){ // removes events with low multiplicity | |
512 | if (fDoKine) FillKineH(); | |
513 | if (fDoCorr) FillCorrH(); | |
514 | if (fDoCorr50) FillCorr50H(); | |
515 | if (fDoShap) FillShapH(jH->GetRadius()); | |
516 | if (fDoFrag) FillFragH(); | |
517 | if (fDoTrig) FillTrigH(); | |
518 | if (fDoJt) FillJtH(); | |
519 | if (fDodNdxi) FilldNdxiH(); | |
520 | } | |
521 | delete jetT; // jet should be deleted before creating a new one | |
522 | if(inJet.GetSize()>fminMult){ // removes events with low multiplicity | |
523 | if (fDoBkgd && nJets>0) FillBkgd(i,iRun); | |
524 | } | |
525 | } // end loop over events in one file | |
526 | if (jFile) jFile->Close(); | |
527 | delete jFile; | |
528 | } // end loop over files | |
529 | } | |
530 | ||
531 | void AliJetAnalysis::FillKineH() | |
532 | { | |
533 | // leading particle | |
534 | if (fDoPart && fPart->LeadingFound()){ | |
535 | fPKineEneH->Fill(fPart->GetE(),fWeight); | |
536 | fPKinePtH->Fill(fPart->GetPt(),fWeight); | |
537 | fPKinePhiH->Fill(fPart->GetPhi(),fWeight); | |
538 | fPKineEtaH->Fill(fPart->GetEta(),fWeight); | |
539 | } | |
540 | // leading generated jet | |
541 | if (fDoGenJ && fGenJ->GetNJets()> 0){ | |
542 | fGKineEneH->Fill(fGenJ->GetE(0),fWeight); | |
543 | fGKinePtH->Fill(fGenJ->GetPt(0),fWeight); | |
544 | fGKinePhiH->Fill(fGenJ->GetPhi(0),fWeight); | |
545 | fGKineEtaH->Fill(fGenJ->GetEta(0),fWeight); | |
546 | } | |
547 | // leading reconstructed jet | |
548 | if (fDoRecJ && fRecJ->GetNJets()> 0) { | |
549 | TArrayF p=fRecJ->GetPtFromSignal(); | |
550 | if (p[0]>fPercentage) { | |
551 | fRKineEneH->Fill(fRecJ->GetE(0)/fEfactor,fWeight); | |
552 | fRKinePtH->Fill(fRecJ->GetPt(0),fWeight); | |
553 | fRKinePhiH->Fill(fRecJ->GetPhi(0),fWeight); | |
554 | fRKineEtaH->Fill(fRecJ->GetEta(0),fWeight); | |
555 | } | |
556 | } | |
557 | } | |
558 | ||
559 | void AliJetAnalysis::FillCorrH() | |
560 | { | |
561 | // Fill correlation histograms | |
562 | if (fDoPart && fPart->LeadingFound() && fDoGenJ && fGenJ->GetNJets()> 0) | |
563 | Correlation(fPart->GetLeading(),fGenJ->GetJet(0), | |
564 | fPGCorrEneH,fPGCorrPtH,fPGCorrEtaH,fPGCorrPhiH); | |
565 | if (fDoPart && fPart->LeadingFound() && fDoRecJ && fRecJ->GetNJets()> 0) { | |
566 | TArrayF p=fRecJ->GetPtFromSignal(); | |
567 | if (p[0]>fPercentage) | |
568 | Correlation(fPart->GetLeading(),fRecJ->GetJet(0), | |
569 | fPRCorrEneH,fPRCorrPtH,fPRCorrEtaH,fPRCorrPhiH); | |
570 | } | |
571 | if (fDoGenJ && fGenJ->GetNJets()> 0 && fDoRecJ && fRecJ->GetNJets()> 0) { | |
572 | TArrayF p=fRecJ->GetPtFromSignal(); | |
573 | if (p[0]>fPercentage) | |
574 | Correlation(fRecJ->GetJet(0), fGenJ->GetJet(0), | |
575 | fRGCorrEneH,fRGCorrPtH,fRGCorrEtaH,fRGCorrPhiH); | |
576 | } | |
577 | } | |
578 | ||
579 | void AliJetAnalysis::Correlation(TLorentzVector *lv1,TLorentzVector *lv2, | |
580 | TH2F *h1, TH2F *h2, TH2F *h3, TH2F *h4) | |
581 | { | |
582 | // Correlation histograms | |
583 | h1->Fill(lv1->E(),lv2->E(),fWeight); | |
584 | h2->Fill(lv1->Pt(),lv2->Pt(),fWeight); | |
585 | h3->Fill(lv1->Eta(),lv2->Eta(),fWeight); | |
586 | Float_t p1, p2; | |
587 | p1 = ((lv1->Phi() < 0) ? (lv1->Phi()) + 2. * TMath::Pi() : lv1->Phi()); | |
588 | p2 = ((lv2->Phi() < 0) ? (lv2->Phi()) + 2. * TMath::Pi() : lv2->Phi()); | |
589 | h4->Fill(p1,p2,fWeight); | |
590 | } | |
591 | ||
592 | void AliJetAnalysis::FillCorr50H() | |
593 | { | |
594 | // Fill correlation histograms when one particle has > 50% of jet energy | |
595 | if (fDoRecJ && fRecJ->GetNJets()> 0) { | |
596 | TArrayF p = fRecJ->GetPtFromSignal(); | |
597 | if (p[0]>fPercentage) { | |
598 | if (fDoPart && fPart->LeadingFound()) | |
599 | Correlation50(fRecJ, fPart->GetLeading(),fRecJ->GetJet(0), | |
600 | fPRCorr50EneH,fPRCorr50PtH,fPRCorr50EtaH,fPRCorr50PhiH); | |
601 | if (fDoGenJ && fGenJ->GetNJets()> 0) | |
602 | Correlation50(fRecJ, fRecJ->GetJet(0), fGenJ->GetJet(0), | |
603 | fRGCorr50EneH,fRGCorr50PtH,fRGCorr50EtaH,fRGCorr50PhiH); | |
604 | } | |
605 | } | |
606 | } | |
607 | ||
608 | void AliJetAnalysis::Correlation50(AliJet *j,TLorentzVector *lv1,TLorentzVector *lv2, | |
609 | TH2F *h1, TH2F *h2, TH2F *h3, TH2F *h4) | |
610 | { | |
611 | // Correlation histograms when one particle has > 50% of jet energy | |
612 | TArrayF ptin = j->GetPtIn(); | |
613 | TArrayF etain = j->GetEtaIn(); | |
614 | TArrayI inJet = j->GetInJet(); | |
615 | ||
616 | Int_t flag = 0; | |
617 | for(Int_t i=0;i<(inJet.GetSize());i++) { | |
618 | if (inJet[i] == 1) { | |
619 | Float_t t1 = TMath::Tan(2.0*TMath::ATan(TMath::Exp(-etain[i]))); | |
620 | Float_t x = (ptin[i]*TMath::Sqrt(1.+1./(t1*t1)))/(j->GetE(0)); | |
621 | if (x>0.5) flag++; | |
622 | } | |
623 | } | |
624 | if (flag>1) cout << " Flag = " << flag << endl; | |
625 | if (flag == 1) { | |
626 | h1->Fill(lv1->E(),lv2->E(),fWeight); | |
627 | h2->Fill(lv1->Pt(),lv2->Pt(),fWeight); | |
628 | h3->Fill(lv1->Eta(),lv2->Eta(),fWeight); | |
629 | Float_t p1, p2; | |
630 | p1 = ((lv1->Phi() < 0) ? | |
631 | (lv1->Phi()) + 2. * TMath::Pi() : lv1->Phi()); | |
632 | p2 = ((lv2->Phi() < 0) ? | |
633 | (lv2->Phi()) + 2. * TMath::Pi() : lv2->Phi()); | |
634 | h4->Fill(p1,p2,fWeight); | |
635 | } | |
636 | } | |
637 | ||
638 | void AliJetAnalysis::FillJtH() | |
639 | { | |
640 | // Fill J_T histogram | |
641 | if (fRecJ->GetNJets()> 0) { | |
642 | fjv3X = 0.0; fjv3Y = 0.0; fjv3Z = 0.0; | |
643 | TArrayF p=fRecJ->GetPtFromSignal(); | |
644 | if (p[0]>fPercentage) { | |
645 | // initialize | |
646 | const TVector3 kJv3 = fRecJ->GetJet(0)->Vect(); | |
647 | fjv3X = kJv3.X(); fjv3Y = kJv3.Y(); fjv3Z = kJv3.Z(); | |
648 | TVector3 trk; | |
649 | TArrayI inJet = fRecJ->GetInJet(); | |
650 | TArrayF etain = fRecJ->GetEtaIn(); | |
651 | TArrayF ptin = fRecJ->GetPtIn(); | |
652 | TArrayF phiin = fRecJ->GetPhiIn(); | |
653 | Float_t deta, dphi,jt, dr; | |
654 | for(Int_t i=0;i<inJet.GetSize();i++) { | |
655 | deta = etain[i] - fRecJ->GetEta(0); | |
656 | dphi = phiin[i] - fRecJ->GetPhi(0); | |
657 | if (dphi < -TMath::Pi()) dphi= -dphi - 2.0 * TMath::Pi(); | |
658 | if (dphi > TMath::Pi()) dphi = 2.0 * TMath::Pi() - dphi; | |
659 | dr = TMath::Sqrt(deta * deta + dphi * dphi); | |
660 | if ((dr<fdrJt) && (ptin[i] > fPartPtCut)) { | |
661 | trk.SetPtEtaPhi(ptin[i],etain[i],phiin[i]); | |
662 | jt = TMath::Sin(trk.Angle(kJv3))*trk.Mag(); | |
663 | fJtH->Fill(jt,fRecJ->GetPt(0),fWeight); | |
62567a80 | 664 | } |
83a444b1 | 665 | } |
666 | fJtW->Fill(fRecJ->GetPt(0),fWeight); | |
667 | fWJt+=fWeight; | |
668 | } | |
669 | } | |
670 | } | |
671 | ||
672 | void AliJetAnalysis::FilldNdxiH() | |
673 | { | |
674 | // Fill dN/d#xi histograms | |
675 | if (fRecJ->GetNJets()> 0) { | |
676 | TArrayF p=fRecJ->GetPtFromSignal(); | |
677 | if (p[0]>fPercentage) { | |
678 | fp0 = p[0]; | |
679 | TArrayI inJet = fRecJ->GetInJet(); | |
680 | TArrayF etain = fRecJ->GetEtaIn(); | |
681 | TArrayF ptin = fRecJ->GetPtIn(); | |
682 | TArrayF phiin = fRecJ->GetPhiIn(); | |
683 | Float_t xi,t1,ene,dr,deta,dphi; | |
684 | for(Int_t i=0;i<inJet.GetSize();i++) { | |
685 | deta = etain[i] - fRecJ->GetEta(0); | |
686 | dphi = phiin[i] - fRecJ->GetPhi(0); | |
687 | if (dphi < -TMath::Pi()) dphi= -dphi - 2.0 * TMath::Pi(); | |
688 | if (dphi > TMath::Pi()) dphi = 2.0 * TMath::Pi() - dphi; | |
689 | dr = TMath::Sqrt(deta * deta + dphi * dphi); | |
690 | if ((dr<fdrdNdxi) && (ptin[i] > fPartPtCut)) { | |
691 | t1 = TMath::Tan(2.0*TMath::ATan(TMath::Exp(-etain[i]))); | |
692 | ene = ptin[i]*TMath::Sqrt(1.+1./(t1*t1)); | |
693 | xi = (Float_t) TMath::Log((fRecJ->GetE(0)/fEfactor)/ene); | |
694 | fdNdxiH->Fill(xi,fRecJ->GetPt(0),fWeight); | |
695 | } | |
696 | } | |
697 | fdNdxiW->Fill(fRecJ->GetPt(0),fWeight); | |
698 | fPtEneH->Fill(fRecJ->GetE(0)/fEfactor,fRecJ->GetPt(0),fWeight); | |
699 | fWdNdxi+=fWeight; | |
700 | } | |
701 | } | |
702 | } | |
703 | ||
704 | void AliJetAnalysis::FillBkgd(Int_t eventN, Int_t runN) | |
705 | { | |
706 | // Background calculated with hijing events (no pythia) | |
707 | if (fp0>fPercentage) { | |
708 | fPtJ=0.,fEJ=0.,fEtaJ=0.,fPhiJ=0.; | |
709 | // Background calculated with hijing events (no pythia) | |
710 | AliJetProductionDataPDC2004* runDataB = new AliJetProductionDataPDC2004(); | |
711 | TFile* jFileB =0x0;; | |
712 | char fnB[256]; | |
62567a80 | 713 | |
83a444b1 | 714 | sprintf(fnB, "%s/%s.root", fBkgdDirectory, (runDataB->GetRunTitle(runN)).Data()); |
715 | jFileB = new TFile(fnB); | |
716 | ||
717 | char nameB[100]; | |
718 | sprintf(nameB, "TreeJ%d",eventN); | |
719 | TTree *jetB =(TTree *)(jFileB->Get(nameB)); | |
720 | jetB->SetBranchAddress("FoundJet", &fRecB); | |
721 | jetB->GetEntry(0); | |
62567a80 | 722 | |
83a444b1 | 723 | TArrayI inJetB = fRecB->GetInJet(); |
724 | TArrayF etainB = fRecB->GetEtaIn(); | |
725 | TArrayF ptinB = fRecB->GetPtIn(); | |
726 | TArrayF phiinB = fRecB->GetPhiIn(); | |
727 | fPtJ = fRecJ->GetPt(0); | |
728 | fEJ = fRecJ->GetE(0); | |
729 | fEtaJ = fRecJ->GetEta(0); | |
730 | fPhiJ = fRecJ->GetPhi(0); | |
731 | Float_t t1,ene,xi,detaB,dphiB,drB,jt,wB; | |
732 | TVector3 trkB; | |
733 | TVector3 jv3B; | |
734 | jv3B.SetX(fjv3X); jv3B.SetY(fjv3Y); jv3B.SetZ(fjv3Z); | |
735 | ||
736 | for(Int_t k=0;k<inJetB.GetSize();k++){ | |
737 | if(ptinB[k] > fPartPtCut){ | |
738 | detaB = etainB[k] - fEtaJ; | |
739 | dphiB = phiinB[k] - fPhiJ; | |
740 | if (dphiB < -TMath::Pi()) dphiB= -dphiB - 2.0 * TMath::Pi(); | |
741 | if (dphiB > TMath::Pi()) dphiB = 2.0 * TMath::Pi() - dphiB; | |
742 | drB = TMath::Sqrt(detaB * detaB + dphiB * dphiB); | |
743 | t1 = TMath::Tan(2.0*TMath::ATan(TMath::Exp(-etainB[k]))); | |
744 | ene = ptinB[k]*TMath::Sqrt(1.+1./(t1*t1)); | |
745 | trkB.SetPtEtaPhi(ptinB[k],etainB[k],phiinB[k]); | |
746 | // --- dN/dxi | |
747 | if (drB<fdrdNdxi) { | |
748 | xi = (Float_t) TMath::Log(fEJ/ene); | |
749 | fdNdxiB->Fill(xi,fPtJ,fWeight); | |
750 | } | |
751 | // --- Jt | |
752 | if (drB<fdrJt) { | |
753 | jt = TMath::Sin(trkB.Angle(jv3B))*(trkB.Mag()); | |
754 | fJtB->Fill(jt,fPtJ,fWeight); | |
755 | } | |
756 | // --- dE/dr | |
757 | if (drB<fdrdEdr) { | |
758 | wB = GetdEdrWeight(fEtaJ,drB)*fWeight*ene; | |
759 | fdEdrB->Fill(drB,fPtJ,wB); | |
760 | } | |
761 | } | |
762 | } | |
763 | delete jetB; | |
764 | if (jFileB) jFileB->Close(); | |
765 | delete jFileB; | |
766 | } | |
767 | } | |
62567a80 | 768 | |
83a444b1 | 769 | void AliJetAnalysis::FillShapH(Float_t r) |
770 | { | |
771 | // Fill jet shape histograms | |
772 | if (fDoRecJ && fRecJ->GetNJets()> 0) { | |
773 | TArrayF p=fRecJ->GetPtFromSignal(); | |
774 | if (p[0]>fPercentage) { | |
775 | Shape(fRecJ,fRShapSelH,fRShapRejH,fRShapAllH,fdEdrH,fPtEneH2,fdEdrW,r); | |
776 | fWdEdr+=fWeight; | |
777 | fWShapR+=fWeight; | |
778 | } | |
779 | } | |
780 | } | |
62567a80 | 781 | |
83a444b1 | 782 | void AliJetAnalysis::Shape(AliJet *j,TH1F* hs, TH1F* hr, TH1F* ha, |
783 | TH2F *hdedr, TH2F *hptene, TH1F* wdedr, Float_t radius) | |
784 | { | |
785 | // initialize | |
786 | TArrayI inJet = j->GetInJet(); | |
787 | TArrayF etain = j->GetEtaIn(); | |
788 | TArrayF ptin = j->GetPtIn(); | |
789 | TArrayF phiin = j->GetPhiIn(); | |
62567a80 | 790 | |
83a444b1 | 791 | // first compute dE/dr histo |
792 | Float_t etaj = j->GetEta(0); | |
793 | Float_t ene,w,deta,dphi,dr,t1; | |
794 | for(Int_t i=0;i<inJet.GetSize();i++) { | |
795 | deta = etain[i] - j->GetEta(0); | |
796 | dphi = phiin[i] - j->GetPhi(0); | |
797 | if (dphi < -TMath::Pi()) dphi= -dphi - 2.0 * TMath::Pi(); | |
798 | if (dphi > TMath::Pi()) dphi = 2.0 * TMath::Pi() - dphi; | |
799 | dr = TMath::Sqrt(deta * deta + dphi * dphi); | |
800 | if ((dr<fdrdEdr) && (ptin[i] > fPartPtCut)) { | |
801 | t1 = TMath::Tan(2.0*TMath::ATan(TMath::Exp(-etain[i]))); | |
802 | ene = ptin[i]*TMath::Sqrt(1.+1./(t1*t1)); | |
803 | w = GetdEdrWeight(etaj,dr)*fWeight*ene; | |
804 | hdedr->Fill(dr,j->GetPt(0),w); | |
805 | } | |
806 | } | |
807 | hptene->Fill(fRecJ->GetE(0),fRecJ->GetPt(0),fWeight); | |
808 | wdedr->Fill(j->GetPt(0),fWeight); | |
62567a80 | 809 | |
83a444b1 | 810 | // now compute shape histos |
811 | Int_t nBins = ha->GetNbinsX(); | |
812 | Float_t xMin = ha->GetXaxis()->GetXmin(); | |
813 | Float_t xMax = ha->GetXaxis()->GetXmax(); | |
814 | Float_t binSize,halfBin; | |
815 | binSize = (xMax-xMin)/nBins; | |
816 | halfBin = binSize/2; | |
817 | Float_t rptS[20], rptR[20], rptA[20]; | |
818 | for(Int_t i=0;i<nBins;i++) rptS[i]=rptR[i]=rptA[i]=0.0; | |
819 | // fill bins in r for leading jet | |
820 | for(Int_t i=0;i<inJet.GetSize();i++) { | |
821 | if (inJet[i] == 1 || (inJet[i] == -1 && fPythia)) { | |
822 | deta = etain[i] - j->GetEta(0); | |
823 | dphi = phiin[i] - j->GetPhi(0); | |
824 | if (dphi < -TMath::Pi()) dphi= -dphi - 2.0 * TMath::Pi(); | |
825 | if (dphi > TMath::Pi()) dphi = 2.0 * TMath::Pi() - dphi; | |
826 | dr = TMath::Sqrt(deta * deta + dphi * dphi); | |
827 | Float_t rR = dr/radius; | |
828 | Int_t bin = (Int_t) TMath::Floor(rR/binSize); | |
829 | rptA[bin]+=ptin[i]/(j->GetPt(0)); | |
830 | if (inJet[i] == 1) rptS[bin]+=ptin[i]/(j->GetPt(0)); | |
831 | if (fPythia && inJet[i] == -1) | |
832 | rptR[bin]+=ptin[i]/(j->GetPt(0)); | |
833 | } | |
834 | } | |
62567a80 | 835 | |
83a444b1 | 836 | // compute shape and fill histogram |
837 | Float_t ptS,ptR,ptA,r; | |
838 | ptS=ptR=ptA=0.0; | |
839 | for (Int_t i=0;i<nBins;i++) { | |
840 | ptS+=rptS[i]; | |
841 | ptR+=rptR[i]; | |
842 | ptA+=rptA[i]; | |
843 | r=(i+1)*binSize-halfBin; | |
844 | hs->Fill(r,ptS*fWeight); | |
845 | if(fPythia) { | |
846 | hr->Fill(r,ptR*fWeight); | |
847 | ha->Fill(r,ptA*fWeight); | |
848 | } | |
849 | } | |
850 | } | |
851 | ||
852 | void AliJetAnalysis::FillFragH() | |
853 | { | |
854 | // Fill fragmentation histogram | |
855 | if (fDoRecJ && fRecJ->GetNJets()> 0) { | |
856 | TArrayF p=fRecJ->GetPtFromSignal(); | |
857 | if (p[0]>fPercentage) { | |
858 | FragFun(fRecJ,fRFragSelH,fRFragRejH,fRFragAllH); | |
859 | fWFragR+=fWeight; | |
860 | } | |
861 | } | |
862 | } | |
62567a80 | 863 | |
83a444b1 | 864 | void AliJetAnalysis::FragFun(AliJet *j,TH1F* hs, TH1F* hr, TH1F* ha) |
865 | { | |
866 | // Calculate de fragmentation function | |
867 | TArrayI inJet = j->GetInJet(); | |
868 | TArrayF etain = j->GetEtaIn(); | |
869 | TArrayF ptin = j->GetPtIn(); | |
62567a80 | 870 | |
83a444b1 | 871 | Float_t xi,t1,ene; |
62567a80 | 872 | |
83a444b1 | 873 | for(Int_t i=0;i<(inJet.GetSize());i++) { |
874 | if (inJet[i] == 1 || (inJet[i] == -1 && fPythia)) { | |
875 | t1 = TMath::Tan(2.0*TMath::ATan(TMath::Exp(-etain[i]))); | |
876 | ene = ptin[i]*TMath::Sqrt(1.+1./(t1*t1)); | |
877 | xi = (Float_t) TMath::Log((j->GetE(0))/ene); | |
878 | if (fPythia) ha->Fill(xi,fWeight); | |
879 | if (inJet[i] == 1) hs->Fill(xi,fWeight); | |
880 | if (fPythia && inJet[i] == -1) hr->Fill(xi,fWeight); | |
881 | } | |
882 | } | |
883 | } | |
62567a80 | 884 | |
83a444b1 | 885 | void AliJetAnalysis:: FillTrigH() |
886 | { | |
887 | // Fill trigger bias histograms | |
888 | if(fDoTrig && fGenJ->GetNJets()>0 && fRecJ->GetNJets()>0){ | |
889 | TArrayF p=fRecJ->GetPtFromSignal(); | |
890 | if (p[0]>fPercentage) { | |
891 | float genEne = fGenJ->GetE(0); | |
892 | float recEne = fRecJ->GetE(0); | |
893 | float eneRatio = (recEne/fEfactor)/genEne; | |
894 | ||
895 | fGTriggerEneH->Fill(genEne, eneRatio, fWeight); | |
896 | fRTriggerEneH->Fill(recEne/fEfactor, eneRatio, fWeight); | |
897 | ||
898 | if (fPart->LeadingFound()){ | |
899 | float leaEne = fPart->GetE(); | |
900 | float eneRatio2 = leaEne/genEne; | |
901 | fGPTriggerEneH->Fill(genEne, eneRatio2, fWeight); | |
902 | fPTriggerEneH->Fill(leaEne/0.2, eneRatio2, fWeight); | |
903 | } | |
904 | } | |
905 | } | |
906 | } | |
907 | ||
908 | ||
909 | //////////////////////////////////////////////////////////////////////// | |
910 | // Normalize histogrames | |
911 | ||
912 | void AliJetAnalysis::NormHistograms() | |
913 | { | |
914 | // normalize shape histograms | |
915 | if (fDoShap) { | |
916 | if (fDoRecJ && fWShapR>0.) { // leading reconstructed jet | |
917 | fRShapSelH->Scale(1.0/fWShapR); | |
918 | fRShapRejH->Scale(1.0/fWShapR); | |
919 | fRShapAllH->Scale(1.0/fWShapR); | |
920 | } | |
921 | } | |
62567a80 | 922 | |
83a444b1 | 923 | // normalize fragmentation function histograms |
924 | if (fDoFrag) { | |
925 | if (fDoRecJ && fWFragR>0.) { // leading reconstructed jet | |
926 | fRFragSelH->Scale(2.0/fWFragR); | |
927 | fRFragRejH->Scale(2.0/fWFragR); | |
928 | fRFragAllH->Scale(2.0/fWFragR); | |
929 | } | |
930 | } | |
931 | } | |
932 | ||
933 | //////////////////////////////////////////////////////////////////////// | |
934 | ||
935 | void AliJetAnalysis::PlotHistograms() | |
936 | { | |
937 | // Plot histogramas (to be done...) | |
938 | if (fDoKine) PlotKineH(); | |
939 | if (fDoCorr) PlotCorrH(); | |
940 | if (fDoCorr50) PlotCorr50H(); | |
941 | if (fDoShap) PlotShapH(); | |
942 | if (fDoFrag) PlotFragH(); | |
943 | if (fDoTrig) PlotTrigH(); | |
944 | } | |
945 | ||
946 | void AliJetAnalysis::PlotKineH() const | |
947 | { | |
948 | // missing | |
83a444b1 | 949 | } |
950 | ||
951 | void AliJetAnalysis::PlotCorrH() const | |
952 | { | |
953 | // missing | |
83a444b1 | 954 | } |
955 | void AliJetAnalysis::PlotCorr50H() const | |
956 | { | |
957 | // missing | |
83a444b1 | 958 | } |
959 | ||
960 | void AliJetAnalysis::PlotShapH() const | |
961 | { | |
962 | // missing | |
83a444b1 | 963 | } |
964 | ||
965 | void AliJetAnalysis::PlotFragH() const | |
966 | { | |
967 | // missing | |
83a444b1 | 968 | } |
969 | ||
970 | void AliJetAnalysis::PlotTrigH() | |
971 | { | |
972 | // missing | |
973 | ||
974 | } | |
975 | ||
976 | //////////////////////////////////////////////////////////////////////// | |
977 | ||
978 | void AliJetAnalysis::SaveHistograms() | |
979 | { | |
980 | // Save histograms | |
981 | TFile *fOut = new TFile(fFile,"recreate"); | |
982 | fOut->cd(); | |
983 | if (fDoKine) SaveKineH(); | |
984 | if (fDoCorr) SaveCorrH(); | |
985 | if (fDoCorr50) SaveCorr50H(); | |
986 | if (fDoShap) SaveShapH(); | |
987 | if (fDoFrag) SaveFragH(); | |
988 | if (fDoTrig) SaveTrigH(); | |
989 | if (fDoJt) SaveJtH(); | |
990 | if (fDodNdxi) SavedNdxiH(); | |
991 | fOut->Close(); | |
992 | } | |
993 | ||
994 | void AliJetAnalysis::SaveKineH() | |
995 | { | |
996 | // Save kinematic histograms | |
997 | if (fDoPart) { | |
998 | fPKineEneH->Write(); | |
999 | fPKinePtH->Write(); | |
1000 | fPKinePhiH->Write(); | |
1001 | fPKineEtaH->Write(); | |
1002 | } | |
62567a80 | 1003 | |
83a444b1 | 1004 | if (fDoGenJ) { |
1005 | fGKineEneH->Write(); | |
1006 | fGKinePtH->Write(); | |
1007 | fGKinePhiH->Write(); | |
1008 | fGKineEtaH->Write(); | |
1009 | } | |
62567a80 | 1010 | |
83a444b1 | 1011 | if (fDoRecJ) { |
1012 | fRKineEneH->Write(); | |
1013 | fRKinePtH->Write(); | |
1014 | fRKinePhiH->Write(); | |
1015 | fRKineEtaH->Write(); | |
1016 | } | |
1017 | } | |
1018 | ||
1019 | void AliJetAnalysis::SaveCorrH() | |
1020 | { | |
1021 | // Save correlation histograms | |
1022 | if (fDoPart && fDoGenJ) { | |
1023 | fPGCorrEneH->Write(); | |
1024 | fPGCorrPtH->Write(); | |
1025 | fPGCorrEtaH->Write(); | |
1026 | fPGCorrPhiH->Write(); | |
1027 | } | |
62567a80 | 1028 | |
83a444b1 | 1029 | if (fDoPart && fDoRecJ) { |
1030 | fPRCorrEneH->Write(); | |
1031 | fPRCorrPtH->Write(); | |
1032 | fPRCorrEtaH->Write(); | |
1033 | fPRCorrPhiH->Write(); | |
1034 | } | |
62567a80 | 1035 | |
83a444b1 | 1036 | if (fDoGenJ && fDoRecJ) { |
1037 | fRGCorrEneH->Write(); | |
1038 | fRGCorrPtH->Write(); | |
1039 | fRGCorrEtaH->Write(); | |
1040 | fRGCorrPhiH->Write(); | |
1041 | } | |
1042 | } | |
1043 | ||
1044 | void AliJetAnalysis::SaveCorr50H() | |
1045 | { | |
1046 | // Save correlation histograms (special case) | |
1047 | if (fDoPart && fDoRecJ) { | |
1048 | fPRCorr50EneH->Write(); | |
1049 | fPRCorr50PtH->Write(); | |
1050 | fPRCorr50EtaH->Write(); | |
1051 | fPRCorr50PhiH->Write(); | |
1052 | } | |
1053 | if (fDoGenJ && fDoRecJ) { | |
1054 | fRGCorr50EneH->Write(); | |
1055 | fRGCorr50PtH->Write(); | |
1056 | fRGCorr50EtaH->Write(); | |
1057 | fRGCorr50PhiH->Write(); | |
1058 | } | |
1059 | } | |
1060 | ||
1061 | void AliJetAnalysis::SaveShapH() | |
1062 | { | |
1063 | // Save jet shape histograms | |
1064 | if (fDoRecJ) { | |
1065 | fRShapSelH->Write(); | |
1066 | fdEdrH->Write(); | |
1067 | if(fDoBkgd) fdEdrB->Write(); | |
1068 | fPtEneH2->Write(); | |
1069 | fdEdrW->Write(); | |
1070 | if (fPythia){ | |
1071 | fRShapRejH->Write(); | |
1072 | fRShapAllH->Write(); | |
1073 | } | |
1074 | } | |
1075 | } | |
1076 | ||
1077 | void AliJetAnalysis::SaveJtH() | |
1078 | { | |
1079 | // Save J_T histograms | |
1080 | if (fDoRecJ) { | |
1081 | fJtH->Write(); | |
1082 | if(fDoBkgd) fJtB->Write(); | |
1083 | fJtW->Write(); | |
1084 | } | |
1085 | } | |
1086 | ||
1087 | void AliJetAnalysis::SavedNdxiH() | |
1088 | { | |
1089 | // Save dN/d#xi histograms | |
1090 | if (fDoRecJ) { | |
1091 | fdNdxiH->Write(); | |
1092 | if(fDoBkgd) fdNdxiB->Write(); | |
1093 | fPtEneH->Write(); | |
1094 | fdNdxiW->Write(); | |
1095 | } | |
1096 | } | |
1097 | ||
1098 | void AliJetAnalysis::SaveFragH() | |
1099 | { | |
1100 | // Save fragmentation histograms | |
1101 | if (fDoRecJ) { | |
1102 | fRFragSelH->Write(); | |
1103 | if(fPythia){ | |
1104 | fRFragRejH->Write(); | |
1105 | fRFragAllH->Write(); | |
1106 | } | |
1107 | } | |
1108 | } | |
1109 | ||
1110 | void AliJetAnalysis::SaveTrigH() | |
1111 | { | |
1112 | // Save trigger bias histograms | |
1113 | if(fDoTrig){ | |
1114 | fGTriggerEneH->Write(); | |
1115 | fRTriggerEneH->Write(); | |
1116 | fGPTriggerEneH->Write(); | |
1117 | fPTriggerEneH->Write(); | |
1118 | } | |
1119 | } | |
1120 | ||
1121 | //////////////////////////////////////////////////////////////////////// | |
1122 | // main Analysis function | |
1123 | ||
1124 | void AliJetAnalysis::Analyze() | |
1125 | ||
1126 | { | |
1127 | // Kinematics for | |
1128 | // leading particle | |
1129 | // leading generated jet | |
1130 | // leading reconstructed jet | |
1131 | ||
1132 | // Correlations amd resolutions | |
1133 | // a) correlations in energy, pt, phi, eta | |
1134 | // b) resolutions in energy, pt, phi, eta, r | |
1135 | // leading particle and leading generated jet | |
1136 | // leading particle and leading reconstructed jet | |
1137 | // leading generated jet and leading reconstructed jet | |
1138 | ||
1139 | // Fragmentation functions and Shapes | |
1140 | // a) integrated over all pt | |
1141 | // b) in 3 flavors: | |
1142 | // b.1) only for user selected particles in jet | |
1143 | // b.2) only for user rejected particles in jet | |
1144 | // b.3) for all particles in jet | |
1145 | ||
1146 | DefineHistograms(); | |
1147 | FillHistograms(); | |
1148 | NormHistograms(); | |
1149 | PlotHistograms(); | |
1150 | SaveHistograms(); | |
1151 | } | |
62567a80 | 1152 |