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