]>
Commit | Line | Data |
---|---|---|
d89b8229 | 1 | /************************************************************************* |
2 | * * | |
3 | * Task for fast embedding * | |
4 | * read extra input from AOD * | |
5 | * * | |
6 | *************************************************************************/ | |
7 | ||
8 | ||
9 | /************************************************************************** | |
10 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
11 | * * | |
12 | * Author: The ALICE Off-line Project. * | |
13 | * Contributors are mentioned in the code where appropriate. * | |
14 | * * | |
15 | * Permission to use, copy, modify and distribute this software and its * | |
16 | * documentation strictly for non-commercial purposes is hereby granted * | |
17 | * without fee, provided that the above copyright notice appears in all * | |
18 | * copies and that both the copyright notice and this permission notice * | |
19 | * appear in the supporting documentation. The authors make no claims * | |
20 | * about the suitability of this software for any purpose. It is * | |
21 | * provided "as is" without express or implied warranty. * | |
22 | **************************************************************************/ | |
23 | ||
24 | /* $Id: */ | |
25 | ||
26 | #include <TFile.h> | |
27 | #include <TTree.h> | |
28 | #include <TChain.h> | |
29 | #include <TClonesArray.h> | |
30 | #include <TDirectory.h> | |
31 | #include <TSystem.h> | |
32 | #include <TRef.h> | |
33 | #include <TRandom3.h> | |
34 | #include <TH1F.h> | |
35 | #include <TH2F.h> | |
36 | #include <TProfile.h> | |
37 | #include <TKey.h> | |
38 | ||
39 | ||
40 | #include "AliAnalysisTaskFastEmbedding.h" | |
41 | #include "AliAnalysisManager.h" | |
42 | #include "AliESDEvent.h" | |
43 | #include "AliESDtrack.h" | |
44 | #include "AliAODEvent.h" | |
45 | #include "AliAODTrack.h" | |
46 | #include "AliAODJet.h" | |
47 | #include "AliAODMCParticle.h" | |
48 | #include "AliAODMCHeader.h" | |
49 | #include "AliInputEventHandler.h" | |
50 | ||
51 | #include "AliLog.h" | |
52 | ||
53 | ClassImp(AliAnalysisTaskFastEmbedding) | |
54 | ||
55 | //__________________________________________________________________________ | |
56 | AliAnalysisTaskFastEmbedding::AliAnalysisTaskFastEmbedding() | |
57 | : AliAnalysisTaskSE() | |
58 | ,fESD(0) | |
59 | ,fAODout(0) | |
60 | ,fAODevent(0) | |
61 | ,fAODtree(0) | |
62 | ,fAODfile(0) | |
63 | ,mcHeader(0) | |
64 | ,rndm(0) | |
65 | ,fInputEntries(0) | |
66 | ,fAODPathArray(0) | |
67 | ,fAODEntriesArray(0) | |
68 | ,fAODPath("AliAOD.root") | |
69 | ,fAODEntries(-1) | |
70 | ,fAODEntriesSum(0) | |
71 | ,fAODEntriesMax(0) | |
72 | ,fOfflineTrgMask(AliVEvent::kAny) | |
73 | ,fMinContribVtx(1) | |
74 | ,fVtxZMin(-8.) | |
75 | ,fVtxZMax(8.) | |
76 | ,fEvtClassMin(0) | |
77 | ,fEvtClassMax(4) | |
78 | ,fCentMin(0.) | |
79 | ,fCentMax(100.) | |
80 | ,fNInputTracksMin(0) | |
81 | ,fNInputTracksMax(-1) | |
82 | ,fTrackBranch("aodExtraTracks") | |
83 | ,fMCparticlesBranch("aodExtraMCparticles") | |
84 | ,fJetBranch("") | |
85 | ,fFileId(-1) | |
86 | ,fAODEntry(0) | |
87 | ,fCountEvents(-1) | |
88 | ,fEmbedMode(0) | |
89 | ,fEvtSelecMode(0) | |
90 | ,fEvtSelMinJetPt(-1) | |
91 | ,fEvtSelMaxJetPt(-1) | |
92 | ,fEvtSelMinJetEta(-999.) | |
93 | ,fEvtSelMaxJetEta( 999.) | |
94 | ,fEvtSelMinJetPhi(0.) | |
95 | ,fEvtSelMaxJetPhi(TMath::Pi()*2.) | |
96 | ,fToyMinNbOfTracks(1) | |
97 | ,fToyMaxNbOfTracks(1) | |
98 | ,fToyMinTrackPt(50.) | |
99 | ,fToyMaxTrackPt(50.) | |
100 | ,fToyDistributionTrackPt(0.) | |
101 | ,fToyMinTrackEta(-.5) | |
102 | ,fToyMaxTrackEta(.5) | |
103 | ,fToyMinTrackPhi(0.) | |
104 | ,fToyMaxTrackPhi(2*TMath::Pi()) | |
105 | ,fToyFilterMap(0) | |
106 | ,fTrackFilterMap(0) | |
107 | ,fNPtHard(10) | |
108 | ,fPtHard(0) | |
109 | ,fPtHardBin(-1) | |
110 | ,fAODJets(0x0) | |
111 | ,fNevents(0) | |
112 | ,fXsection(0) | |
113 | ,fAvgTrials(0) | |
114 | ,fHistList(0) | |
115 | ,fHistEvtSelection(0) | |
116 | ,fh1Xsec(0) | |
117 | ,fh1Trials(0) | |
118 | ,fh1TrialsEvtSel(0) | |
119 | ,fh2PtHard(0) | |
120 | ,fh2PtHardEvtSel(0) | |
121 | ,fh2PtHardTrials(0) | |
122 | ,fh1TrackPt(0) | |
123 | ,fh2TrackEtaPhi(0) | |
124 | ,fh1TrackN(0) | |
125 | ,fh1JetPt(0) | |
126 | ,fh2JetEtaPhi(0) | |
127 | ,fh1JetN(0) | |
128 | ,fh1MCTrackPt(0) | |
129 | ,fh2MCTrackEtaPhi(0) | |
130 | ,fh1MCTrackN(0) | |
131 | ,fh1AODfile(0) | |
132 | ,fh2AODevent(0) | |
133 | { | |
134 | // default constructor | |
135 | ||
136 | } | |
137 | ||
138 | ||
139 | //__________________________________________________________________________ | |
140 | AliAnalysisTaskFastEmbedding::AliAnalysisTaskFastEmbedding(const char *name) | |
141 | : AliAnalysisTaskSE(name) | |
142 | ,fESD(0) | |
143 | ,fAODout(0) | |
144 | ,fAODevent(0) | |
145 | ,fAODtree(0) | |
146 | ,fAODfile(0) | |
147 | ,mcHeader(0) | |
148 | ,rndm(0) | |
149 | ,fInputEntries(0) | |
150 | ,fAODPathArray(0) | |
151 | ,fAODEntriesArray(0) | |
152 | ,fAODPath("AliAOD.root") | |
153 | ,fAODEntries(-1) | |
154 | ,fAODEntriesSum(0) | |
155 | ,fAODEntriesMax(0) | |
156 | ,fOfflineTrgMask(AliVEvent::kAny) | |
157 | ,fMinContribVtx(1) | |
158 | ,fVtxZMin(-8.) | |
159 | ,fVtxZMax(8.) | |
160 | ,fEvtClassMin(0) | |
161 | ,fEvtClassMax(4) | |
162 | ,fCentMin(0.) | |
163 | ,fCentMax(100.) | |
164 | ,fNInputTracksMin(0) | |
165 | ,fNInputTracksMax(-1) | |
166 | ,fTrackBranch("aodExtraTracks") | |
167 | ,fMCparticlesBranch("aodExtraMCparticles") | |
168 | ,fJetBranch("") | |
169 | ,fFileId(-1) | |
170 | ,fAODEntry(0) | |
171 | ,fCountEvents(-1) | |
172 | ,fEmbedMode(0) | |
173 | ,fEvtSelecMode(0) | |
174 | ,fEvtSelMinJetPt(-1) | |
175 | ,fEvtSelMaxJetPt(-1) | |
176 | ,fEvtSelMinJetEta(-999.) | |
177 | ,fEvtSelMaxJetEta( 999.) | |
178 | ,fEvtSelMinJetPhi(0.) | |
179 | ,fEvtSelMaxJetPhi(TMath::Pi()*2.) | |
180 | ,fToyMinNbOfTracks(1) | |
181 | ,fToyMaxNbOfTracks(1) | |
182 | ,fToyMinTrackPt(50.) | |
183 | ,fToyMaxTrackPt(50.) | |
184 | ,fToyDistributionTrackPt(0.) | |
185 | ,fToyMinTrackEta(-.5) | |
186 | ,fToyMaxTrackEta(.5) | |
187 | ,fToyMinTrackPhi(0.) | |
188 | ,fToyMaxTrackPhi(2*TMath::Pi()) | |
189 | ,fToyFilterMap(0) | |
190 | ,fTrackFilterMap(0) | |
191 | ,fNPtHard(10) | |
192 | ,fPtHard(0) | |
193 | ,fPtHardBin(-1) | |
194 | ,fAODJets(0x0) | |
195 | ,fNevents(0) | |
196 | ,fXsection(0) | |
197 | ,fAvgTrials(0) | |
198 | ,fHistList(0) | |
199 | ,fHistEvtSelection(0) | |
200 | ,fh1Xsec(0) | |
201 | ,fh1Trials(0) | |
202 | ,fh1TrialsEvtSel(0) | |
203 | ,fh2PtHard(0) | |
204 | ,fh2PtHardEvtSel(0) | |
205 | ,fh2PtHardTrials(0) | |
206 | ,fh1TrackPt(0) | |
207 | ,fh2TrackEtaPhi(0) | |
208 | ,fh1TrackN(0) | |
209 | ,fh1JetPt(0) | |
210 | ,fh2JetEtaPhi(0) | |
211 | ,fh1JetN(0) | |
212 | ,fh1MCTrackPt(0) | |
213 | ,fh2MCTrackEtaPhi(0) | |
214 | ,fh1MCTrackN(0) | |
215 | ,fh1AODfile(0) | |
216 | ,fh2AODevent(0) | |
217 | { | |
218 | // constructor | |
219 | DefineOutput(1, TList::Class()); | |
220 | } | |
221 | ||
222 | ||
223 | //__________________________________________________________________________ | |
224 | AliAnalysisTaskFastEmbedding::AliAnalysisTaskFastEmbedding(const AliAnalysisTaskFastEmbedding ©) | |
225 | : AliAnalysisTaskSE() | |
226 | ,fESD(copy.fESD) | |
227 | ,fAODout(copy.fAODout) | |
228 | ,fAODevent(copy.fAODevent) | |
229 | ,fAODtree(copy.fAODtree) | |
230 | ,fAODfile(copy.fAODfile) | |
231 | ,mcHeader(copy.mcHeader) | |
232 | ,rndm(copy.rndm) | |
233 | ,fInputEntries(copy.fInputEntries) | |
234 | ,fAODPathArray(copy.fAODPathArray) | |
235 | ,fAODEntriesArray(copy.fAODEntriesArray) | |
236 | ,fAODPath(copy.fAODPath) | |
237 | ,fAODEntries(copy.fAODEntries) | |
238 | ,fAODEntriesSum(copy.fAODEntriesSum) | |
239 | ,fAODEntriesMax(copy.fAODEntriesMax) | |
240 | ,fOfflineTrgMask(copy.fOfflineTrgMask) | |
241 | ,fMinContribVtx(copy.fMinContribVtx) | |
242 | ,fVtxZMin(copy.fVtxZMin) | |
243 | ,fVtxZMax(copy.fVtxZMax) | |
244 | ,fEvtClassMin(copy.fEvtClassMin) | |
245 | ,fEvtClassMax(copy.fEvtClassMax) | |
246 | ,fCentMin(copy.fCentMin) | |
247 | ,fCentMax(copy.fCentMax) | |
248 | ,fNInputTracksMin(copy.fNInputTracksMin) | |
249 | ,fNInputTracksMax(copy.fNInputTracksMax) | |
250 | ,fTrackBranch(copy.fTrackBranch) | |
251 | ,fMCparticlesBranch(copy.fMCparticlesBranch) | |
252 | ,fJetBranch(copy.fJetBranch) | |
253 | ,fFileId(copy.fFileId) | |
254 | ,fAODEntry(copy.fAODEntry) | |
255 | ,fCountEvents(copy.fCountEvents) | |
256 | ,fEmbedMode(copy.fEmbedMode) | |
257 | ,fEvtSelecMode(copy.fEvtSelecMode) | |
258 | ,fEvtSelMinJetPt(copy.fEvtSelMinJetPt) | |
259 | ,fEvtSelMaxJetPt(copy.fEvtSelMaxJetPt) | |
260 | ,fEvtSelMinJetEta(copy.fEvtSelMinJetEta) | |
261 | ,fEvtSelMaxJetEta(copy.fEvtSelMaxJetEta) | |
262 | ,fEvtSelMinJetPhi(copy.fEvtSelMinJetPhi) | |
263 | ,fEvtSelMaxJetPhi(copy.fEvtSelMaxJetPhi) | |
264 | ,fToyMinNbOfTracks(copy.fToyMinNbOfTracks) | |
265 | ,fToyMaxNbOfTracks(copy.fToyMaxNbOfTracks) | |
266 | ,fToyMinTrackPt(copy.fToyMinTrackPt) | |
267 | ,fToyMaxTrackPt(copy.fToyMaxTrackPt) | |
268 | ,fToyDistributionTrackPt(copy.fToyDistributionTrackPt) | |
269 | ,fToyMinTrackEta(copy.fToyMinTrackEta) | |
270 | ,fToyMaxTrackEta(copy.fToyMaxTrackEta) | |
271 | ,fToyMinTrackPhi(copy.fToyMinTrackPhi) | |
272 | ,fToyMaxTrackPhi(copy.fToyMaxTrackPhi) | |
273 | ,fToyFilterMap(copy.fToyFilterMap) | |
274 | ,fTrackFilterMap(copy.fTrackFilterMap) | |
275 | ,fNPtHard(copy.fNPtHard) | |
276 | ,fPtHard(copy.fPtHard) | |
277 | ,fPtHardBin(copy.fPtHardBin) | |
278 | ,fAODJets(copy.fAODJets) | |
279 | ,fNevents(copy.fNevents) | |
280 | ,fXsection(copy.fXsection) | |
281 | ,fAvgTrials(copy.fAvgTrials) | |
282 | ,fHistList(copy.fHistList) | |
283 | ,fHistEvtSelection(copy.fHistEvtSelection) | |
284 | ,fh1Xsec(copy.fh1Xsec) | |
285 | ,fh1Trials(copy.fh1Trials) | |
286 | ,fh1TrialsEvtSel(copy.fh1TrialsEvtSel) | |
287 | ,fh2PtHard(copy.fh2PtHard) | |
288 | ,fh2PtHardEvtSel(copy.fh2PtHardEvtSel) | |
289 | ,fh2PtHardTrials(copy.fh2PtHardTrials) | |
290 | ,fh1TrackPt(copy.fh1TrackPt) | |
291 | ,fh2TrackEtaPhi(copy.fh2TrackEtaPhi) | |
292 | ,fh1TrackN(copy.fh1TrackN) | |
293 | ,fh1JetPt(copy.fh1JetPt) | |
294 | ,fh2JetEtaPhi(copy.fh2JetEtaPhi) | |
295 | ,fh1JetN(copy.fh1JetN) | |
296 | ,fh1MCTrackPt(copy.fh1MCTrackPt) | |
297 | ,fh2MCTrackEtaPhi(copy.fh2MCTrackEtaPhi) | |
298 | ,fh1MCTrackN(copy.fh1MCTrackN) | |
299 | ,fh1AODfile(copy.fh1AODfile) | |
300 | ,fh2AODevent(copy.fh2AODevent) | |
301 | { | |
302 | // copy constructor | |
303 | } | |
304 | ||
305 | ||
306 | //__________________________________________________________________________ | |
307 | AliAnalysisTaskFastEmbedding& AliAnalysisTaskFastEmbedding::operator=(const AliAnalysisTaskFastEmbedding& o) | |
308 | { | |
309 | // assignment | |
310 | ||
311 | if(this!=&o){ | |
312 | AliAnalysisTaskSE::operator=(o); | |
313 | fESD = o.fESD; | |
314 | fAODout = o.fAODout; | |
315 | fAODevent = o.fAODevent; | |
316 | fAODtree = o.fAODtree; | |
317 | fAODfile = o.fAODfile; | |
318 | mcHeader = o.mcHeader; | |
319 | rndm = o.rndm; | |
320 | fInputEntries = o.fInputEntries; | |
321 | fAODPathArray = o.fAODPathArray; | |
322 | fAODEntriesArray = o.fAODEntriesArray; | |
323 | fAODPath = o.fAODPath; | |
324 | fAODEntries = o.fAODEntries; | |
325 | fAODEntriesSum = o.fAODEntriesSum; | |
326 | fAODEntriesMax = o.fAODEntriesMax; | |
327 | fOfflineTrgMask = o.fOfflineTrgMask; | |
328 | fMinContribVtx = o.fMinContribVtx; | |
329 | fVtxZMin = o.fVtxZMin; | |
330 | fVtxZMax = o.fVtxZMax; | |
331 | fEvtClassMin = o.fEvtClassMin; | |
332 | fEvtClassMax = o.fEvtClassMax; | |
333 | fCentMin = o.fCentMin; | |
334 | fCentMax = o.fCentMax; | |
335 | fNInputTracksMin = o.fNInputTracksMin; | |
336 | fNInputTracksMax = o.fNInputTracksMax; | |
337 | fTrackBranch = o.fTrackBranch; | |
338 | fMCparticlesBranch = o.fMCparticlesBranch; | |
339 | fJetBranch = o.fJetBranch; | |
340 | fFileId = o.fFileId; | |
341 | fAODEntry = o.fAODEntry; | |
342 | fCountEvents = o.fCountEvents; | |
343 | fEmbedMode = o.fEmbedMode; | |
344 | fEvtSelecMode = o.fEvtSelecMode; | |
345 | fEvtSelMinJetPt = o.fEvtSelMinJetPt; | |
346 | fEvtSelMaxJetPt = o.fEvtSelMaxJetPt; | |
347 | fEvtSelMinJetEta = o.fEvtSelMinJetEta; | |
348 | fEvtSelMaxJetEta = o.fEvtSelMaxJetEta; | |
349 | fEvtSelMinJetPhi = o.fEvtSelMinJetPhi; | |
350 | fEvtSelMaxJetPhi = o.fEvtSelMaxJetPhi; | |
351 | fToyMinNbOfTracks = o.fToyMinNbOfTracks; | |
352 | fToyMaxNbOfTracks = o.fToyMaxNbOfTracks; | |
353 | fToyMinTrackPt = o.fToyMinTrackPt; | |
354 | fToyMaxTrackPt = o.fToyMaxTrackPt; | |
355 | fToyDistributionTrackPt = o.fToyDistributionTrackPt; | |
356 | fToyMinTrackEta = o.fToyMinTrackEta; | |
357 | fToyMaxTrackEta = o.fToyMaxTrackEta; | |
358 | fToyMinTrackPhi = o.fToyMinTrackPhi; | |
359 | fToyMaxTrackPhi = o.fToyMaxTrackPhi; | |
360 | fToyFilterMap = o.fToyFilterMap; | |
361 | fTrackFilterMap = o.fTrackFilterMap; | |
362 | fNPtHard = o.fNPtHard; | |
363 | fPtHard = o.fPtHard; | |
364 | fPtHardBin = o.fPtHardBin; | |
365 | fAODJets = o.fAODJets; | |
366 | fNevents = o.fNevents; | |
367 | fXsection = o.fXsection; | |
368 | fAvgTrials = o.fAvgTrials; | |
369 | fHistList = o.fHistList; | |
370 | fHistEvtSelection = o.fHistEvtSelection; | |
371 | fh1Xsec = o.fh1Xsec; | |
372 | fh1Trials = o.fh1Trials; | |
373 | fh1TrialsEvtSel = o.fh1TrialsEvtSel; | |
374 | fh2PtHard = o.fh2PtHard; | |
375 | fh2PtHardEvtSel = o.fh2PtHardEvtSel; | |
376 | fh2PtHardTrials = o.fh2PtHardTrials; | |
377 | fh1TrackPt = o.fh1TrackPt; | |
378 | fh2TrackEtaPhi = o.fh2TrackEtaPhi; | |
379 | fh1TrackN = o.fh1TrackN; | |
380 | fh1JetPt = o.fh1JetPt; | |
381 | fh2JetEtaPhi = o.fh2JetEtaPhi; | |
382 | fh1JetN = o.fh1JetN; | |
383 | fh1MCTrackPt = o.fh1MCTrackPt; | |
384 | fh2MCTrackEtaPhi = o.fh2MCTrackEtaPhi; | |
385 | fh1MCTrackN = o.fh1MCTrackN; | |
386 | fh1AODfile = o.fh1AODfile; | |
387 | fh2AODevent = o.fh2AODevent; | |
388 | } | |
389 | ||
390 | return *this; | |
391 | } | |
392 | ||
393 | ||
394 | //__________________________________________________________________________ | |
395 | AliAnalysisTaskFastEmbedding::~AliAnalysisTaskFastEmbedding() | |
396 | { | |
397 | // destructor | |
398 | delete rndm; | |
399 | } | |
400 | ||
401 | ||
402 | //__________________________________________________________________________ | |
403 | void AliAnalysisTaskFastEmbedding::UserCreateOutputObjects() | |
404 | { | |
405 | // create output objects | |
406 | if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::UserCreateOutputObjects()"); | |
407 | AliLog::SetClassDebugLevel("AliAnalysisTaskFastEmbedding", AliLog::kInfo); | |
408 | ||
409 | OpenFile(1); | |
410 | if(!fHistList) fHistList = new TList(); | |
411 | fHistList->SetOwner(kTRUE); | |
412 | ||
413 | ||
414 | // set seed | |
415 | rndm = new TRandom3(); | |
416 | Int_t id = GetJobID(); | |
417 | if(id>-1) rndm->SetSeed(id); | |
418 | else rndm->SetSeed(); // a TTUID is generated and used for seed | |
419 | AliInfo(Form("TRandom3 seed: %d", rndm->GetSeed())); | |
420 | ||
421 | ||
422 | ||
423 | // embed mode with AOD | |
424 | if(fEmbedMode==kAODFull || fEmbedMode==kAODJetTracks || fEmbedMode==kAODJet4Mom){ | |
425 | ||
426 | // open input AOD | |
427 | fFileId = OpenAODfile(); | |
428 | fNevents = 0; // force to open another aod in UserExec() | |
429 | if(fFileId<0){ | |
430 | AliError(""); | |
431 | PostData(1, fHistList); | |
432 | return; | |
433 | } | |
434 | } //end: embed mode with AOD | |
435 | ||
436 | ||
437 | // connect output aod | |
438 | // create a new branch for extra tracks | |
439 | fAODout = AODEvent(); | |
440 | if(!fAODout){ | |
441 | AliError("Output AOD not found."); | |
442 | PostData(1, fHistList); | |
443 | return; | |
444 | } | |
445 | if(!fAODout->FindListObject(fTrackBranch.Data()) && strlen(fTrackBranch.Data())){ | |
446 | AliInfo(Form("Add AOD branch %s", fTrackBranch.Data())); | |
447 | TClonesArray *tracks = new TClonesArray("AliAODTrack",0); | |
448 | tracks->SetName(fTrackBranch.Data()); | |
449 | AddAODBranch("TClonesArray", &tracks); | |
450 | } | |
451 | // create new branch for extra mcparticle if available as input | |
452 | if(fAODevent && fAODevent->FindListObject("mcparticles") && strlen(fMCparticlesBranch.Data())){ | |
453 | AliInfo(Form("Add AOD branch %s", fMCparticlesBranch.Data())); | |
454 | TClonesArray *mcparticles = new TClonesArray("AliAODMCParticle",0); | |
455 | mcparticles->SetName(fMCparticlesBranch.Data()); | |
456 | AddAODBranch("TClonesArray", &mcparticles); | |
457 | } | |
458 | ||
459 | ||
460 | ||
461 | //qa histograms | |
462 | Bool_t oldStatus = TH1::AddDirectoryStatus(); | |
463 | TH1::AddDirectory(kFALSE); | |
464 | ||
465 | fHistEvtSelection = new TH1I("fHistEvtSelection", "event selection", 6, -0.5, 5.5); | |
466 | fHistEvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED"); | |
467 | fHistEvtSelection->GetXaxis()->SetBinLabel(2,"events IN"); | |
468 | fHistEvtSelection->GetXaxis()->SetBinLabel(3,"event selection (rejected)"); | |
469 | fHistEvtSelection->GetXaxis()->SetBinLabel(4,"vertex cut (rejected)"); | |
470 | fHistEvtSelection->GetXaxis()->SetBinLabel(5,"centrality (rejected)"); | |
471 | fHistEvtSelection->GetXaxis()->SetBinLabel(6,"multiplicity (rejected)"); | |
472 | ||
473 | fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root;p_{T,hard} bin;<#sigma>",fNPtHard+1,-1.5,fNPtHard-0.5); | |
474 | fh1Trials = new TH1F("fh1Trials","trials (simulation) from pyxsec.root;p_{T,hard} bin;#sum{ntrials}",fNPtHard+1,-1.5,fNPtHard-0.5); | |
475 | fh1TrialsEvtSel = new TH1F("fh1TrialsEvtSel","trials (event selection) from pyxsec.root;p_{T,hard} bin;#sum{ntrials}",fNPtHard+1,-1.5,fNPtHard-0.5); | |
476 | fh2PtHard = new TH2F("fh2PtHard","PYTHIA Pt hard;p_{T,hard} bin;p_{T,hard}",fNPtHard+1,-1.5,fNPtHard-0.5,350,-.5,349.5); | |
477 | fh2PtHardEvtSel = new TH2F("fh2PtHardEvtSel","PYTHIA Pt hard (event selection);p_{T,hard} bin;p_{T,hard}",fNPtHard+1,-1.5,fNPtHard-0.5,350,-.5,349.5); | |
478 | fh2PtHardTrials = new TH2F("fh2PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard} bin;#sum{p_{T,hard}}",fNPtHard+1,-1.5,fNPtHard-0.5,350,-.5,349.5); | |
479 | ||
480 | fHistList->Add(fHistEvtSelection); | |
481 | fHistList->Add(fh1Xsec); | |
482 | fHistList->Add(fh1Trials); | |
483 | fHistList->Add(fh1TrialsEvtSel); | |
484 | fHistList->Add(fh2PtHard); | |
485 | fHistList->Add(fh2PtHardEvtSel); | |
486 | fHistList->Add(fh2PtHardTrials); | |
487 | ||
488 | fh1TrackPt = new TH1F("fh1TrackPt","pT of extra tracks;p_{T};entries", 250, 0., 250.); | |
489 | fh2TrackEtaPhi = new TH2F("fh2TrackEtaPhi","eta-phi distribution of extra tracks;#eta;#phi", 20, -1., 1., 60, 0., 2*TMath::Pi()); | |
490 | fh1TrackN = new TH1F("fh1TrackN", "nb. of extra tracks per event;nb. of tracks;entries",300, 0., 300.); | |
491 | ||
492 | fHistList->Add(fh1TrackPt); | |
493 | fHistList->Add(fh2TrackEtaPhi); | |
494 | fHistList->Add(fh1TrackN); | |
495 | ||
496 | if(fEmbedMode==kAODFull || fEmbedMode==kAODJetTracks || fEmbedMode==kAODJet4Mom){ | |
497 | ||
498 | fh1JetPt = new TH1F("fh1JetPt", "pT of extra jets;p_{T};entries", 120, 0., 120.); | |
499 | fh2JetEtaPhi = new TH2F("fh2JetEtaPhi", "eta-phi distribution of extra jets;#eta;#phi", 20, -1., 1., 60, 0., 2*TMath::Pi()); | |
500 | fh1JetN = new TH1F("fh1JetN", "nb. of extra jets per event;nb. of jets;entries",20,0.,20.); | |
501 | ||
502 | fHistList->Add(fh1JetPt); | |
503 | fHistList->Add(fh2JetEtaPhi); | |
504 | fHistList->Add(fh1JetN); | |
505 | } | |
506 | ||
507 | ||
508 | if(fAODevent && fAODevent->FindListObject("mcparticles") && strlen(fMCparticlesBranch.Data())){ | |
509 | ||
510 | fh1MCTrackPt = new TH1F("fh1MCTrackPt","pT of MC extra tracks;p_{T};entries", 250, 0., 250.); | |
511 | fh2MCTrackEtaPhi = new TH2F("fh2MCTrackEtaPhi","eta-phi distribution of MC extra tracks;#eta;#phi", 20, -1., 1., 60, 0., 2*TMath::Pi()); | |
512 | fh1MCTrackN = new TH1F("fh1MCTrackN", "nb. of MC extra tracks per event;nb. of tracks;entries",300, 0., 300.); | |
513 | ||
514 | fHistList->Add(fh1MCTrackPt); | |
515 | fHistList->Add(fh2MCTrackEtaPhi); | |
516 | fHistList->Add(fh1MCTrackN); | |
517 | ||
518 | } | |
519 | ||
520 | fh1AODfile = new TH1I("fh1AODfile", "overview of opened AOD files from the array", 2300, -0.5, 2299.5); | |
521 | fh2AODevent = new TH2I("fh1AODevent","selected events;file;event", 2500,-.5,2499.5,5000,-.5,4999.5); | |
522 | fHistList->Add(fh1AODfile); | |
523 | fHistList->Add(fh2AODevent); | |
524 | ||
525 | // =========== Switch on Sumw2 for all histos =========== | |
526 | for (Int_t i=0; i<fHistList->GetEntries(); ++i) { | |
527 | TH1 *h1 = dynamic_cast<TH1*>(fHistList->At(i)); | |
528 | if (h1){ | |
529 | h1->Sumw2(); | |
530 | continue; | |
531 | } | |
532 | } | |
533 | ||
534 | TH1::AddDirectory(oldStatus); | |
535 | ||
536 | PostData(1, fHistList); | |
537 | } | |
538 | ||
539 | ||
540 | //__________________________________________________________________________ | |
541 | void AliAnalysisTaskFastEmbedding::Init() | |
542 | { | |
543 | // Initialization | |
544 | if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::Init()"); | |
545 | ||
546 | } | |
547 | ||
548 | //__________________________________________________________________________ | |
549 | Bool_t AliAnalysisTaskFastEmbedding::UserNotify() | |
550 | { | |
b284a58f | 551 | // User defined Notify(), called once |
d89b8229 | 552 | if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::UserNotify()"); |
553 | ||
554 | // get total nb of events in tree (of this subjob) | |
555 | AliInputEventHandler* inputHandler = (AliInputEventHandler*)((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler()); | |
556 | fInputEntries = inputHandler->GetTree()->GetEntriesFast(); | |
557 | AliInfo(Form("Total nb. of events: %d", fInputEntries)); | |
558 | ||
559 | return kTRUE; | |
560 | ||
561 | } | |
562 | ||
563 | ||
564 | //__________________________________________________________________________ | |
565 | void AliAnalysisTaskFastEmbedding::UserExec(Option_t *) | |
566 | { | |
b284a58f | 567 | // Analysis, called once per event |
d89b8229 | 568 | if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::UserExec()"); |
569 | ||
570 | if(!fAODout){ | |
571 | AliError("Need output AOD, but is not connected."); | |
572 | PostData(1, fHistList); | |
573 | return; | |
574 | } | |
575 | ||
576 | fESD=dynamic_cast<AliESDEvent*>(InputEvent()); | |
577 | if (!fESD) { | |
578 | AliError("ESD not available"); | |
579 | PostData(1, fHistList); | |
580 | return; | |
581 | } | |
582 | ||
583 | // -- event selection -- | |
584 | fHistEvtSelection->Fill(1); // number of events before event selection | |
585 | ||
586 | // physics selection | |
587 | AliInputEventHandler* inputHandler = (AliInputEventHandler*) | |
588 | ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler()); | |
589 | if(!(inputHandler->IsEventSelected() & fOfflineTrgMask)){ | |
590 | if(fDebug) Printf(" Trigger Selection: event REJECTED ... "); | |
591 | fHistEvtSelection->Fill(2); | |
592 | PostData(1, fHistList); | |
593 | return; | |
594 | } | |
595 | ||
596 | // vertex selection | |
597 | AliAODVertex* primVtx = fAODout->GetPrimaryVertex(); | |
598 | Int_t nTracksPrim = primVtx->GetNContributors(); | |
599 | if ((nTracksPrim < fMinContribVtx) || | |
600 | (primVtx->GetZ() < fVtxZMin) || | |
601 | (primVtx->GetZ() > fVtxZMax) ){ | |
602 | if(fDebug) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ()); | |
603 | fHistEvtSelection->Fill(3); | |
604 | PostData(1, fHistList); | |
605 | return; | |
606 | } | |
607 | ||
608 | /** takes wrong value, since jet service tasks runs later | |
609 | // event class selection (from jet helper task) | |
610 | Int_t eventClass = AliAnalysisHelperJetTasks::EventClass(); | |
611 | if(fDebug) Printf("Event class %d", eventClass); | |
612 | if (eventClass < fEvtClassMin || eventClass > fEvtClassMax){ | |
613 | fHistEvtSelection->Fill(4); | |
614 | PostData(1, fHistList); | |
615 | return; | |
616 | }*/ | |
617 | ||
618 | // centrality selection | |
619 | AliCentrality *cent = 0x0; | |
620 | Float_t centValue = 0.; | |
621 | if(fESD) cent = fESD->GetCentrality(); | |
622 | if(cent) centValue = cent->GetCentralityPercentile("V0M"); | |
623 | if(fDebug) printf("centrality: %f\n", centValue); | |
624 | if (centValue < fCentMin || centValue > fCentMax){ | |
625 | fHistEvtSelection->Fill(4); | |
626 | PostData(1, fHistList); | |
627 | return; | |
628 | } | |
629 | ||
630 | ||
631 | /* ** not implemented ** | |
632 | // multiplicity due to input tracks | |
633 | Int_t nInputTracks = GetNInputTracks(); | |
634 | ||
635 | if (nInputTracks < fNInputTracksMin || (fNInputTracksMax > -1 && nInputTracks > fNInputTracksMax)){ | |
636 | fHistEvtSelection->Fill(5); | |
637 | PostData(1, fHistList); | |
638 | return; | |
639 | } | |
640 | */ | |
641 | ||
642 | fHistEvtSelection->Fill(0); // accepted events | |
643 | // -- end event selection -- | |
644 | ||
645 | // connect aod out | |
646 | TClonesArray *tracks = (TClonesArray*)(fAODout->FindListObject(fTrackBranch.Data())); | |
647 | if(!tracks){ | |
648 | AliError("Extra track branch not found in output."); | |
649 | PostData(1, fHistList); | |
650 | return; | |
651 | } | |
652 | tracks->Delete(); | |
653 | Int_t nAODtracks=0; | |
654 | ||
655 | TRef dummy; | |
656 | ||
657 | // === embed mode with AOD === | |
658 | if(fEmbedMode==kAODFull || fEmbedMode==kAODJetTracks || fEmbedMode==kAODJet4Mom){ | |
659 | if(!fAODevent){ | |
660 | AliError("Need input AOD, but is not connected."); | |
661 | PostData(1, fHistList); | |
662 | return; | |
663 | } | |
664 | ||
665 | ||
666 | ||
667 | Bool_t useEntry = kFALSE; | |
668 | while(!useEntry){ // protection need, if no event fulfills requierment | |
669 | ||
670 | fAODEntry++; // go to next event | |
671 | fCountEvents++; | |
672 | if(fAODEntry>=fNevents) fAODEntry=0; | |
673 | ||
674 | // open new aod file | |
675 | if(fCountEvents>=fNevents){ | |
676 | if(!fAODPathArray){ | |
677 | AliDebug(AliLog::kDebug, "No list of AODs, keep current AOD."); | |
678 | } | |
679 | else { | |
680 | AliDebug(AliLog::kDebug, "Select new AOD file ..."); | |
681 | ||
682 | fFileId = OpenAODfile(); | |
683 | if(fFileId<0) { | |
684 | PostData(1, fHistList); | |
685 | return; | |
686 | } | |
687 | fh1AODfile->Fill(fFileId); | |
688 | if(fAODEntries!=fNevents){ | |
689 | AliError("File with list of AODs and file with nb. of entries inconsistent"); | |
690 | PostData(1, fHistList); | |
691 | return; | |
692 | } | |
693 | } | |
694 | } | |
695 | ||
696 | ||
697 | ||
698 | // get next event | |
699 | fAODtree->GetEvent(fAODEntry); | |
700 | ||
701 | // get pt hard | |
702 | if(mcHeader){ | |
703 | fPtHard = mcHeader->GetPtHard(); | |
704 | GetPtHard(kTRUE,fPtHard); // make value available for other tasks (e.g. jet response task) | |
705 | AliDebug(AliLog::kDebug,Form("pt hard %.2f", fPtHard)); | |
706 | } | |
707 | else{ | |
708 | AliDebug(AliLog::kDebug,"No mcHeader"); | |
709 | } | |
710 | fPtHardBin = GetPtHardBin(fPtHard); | |
711 | ||
712 | //Printf("pt hard (bin) %.2f (%d), xSection %.2e, trials %.2f",fPtHard, fPtHardBin, fXsection, fAvgTrials); | |
713 | ||
714 | // fill event variables for each event | |
715 | fh1Xsec->Fill(fPtHardBin,fXsection); | |
716 | fh2PtHard->Fill(fPtHardBin,fPtHard); | |
717 | ||
718 | fh1Trials->Fill(fPtHardBin, fAvgTrials); | |
719 | fh1TrialsEvtSel->Fill(fPtHardBin); | |
720 | fh2PtHardTrials->Fill(fPtHardBin,fPtHard,fAvgTrials); | |
721 | ||
722 | // jet pt selection | |
723 | if(fEvtSelecMode==kEventsJetPt){ | |
724 | Int_t nJets = fAODJets->GetEntries(); | |
725 | //Printf("nb. jets: %d",nJets); | |
726 | for(Int_t ij=0; ij<nJets; ++ij){ | |
727 | AliAODJet *jet = dynamic_cast<AliAODJet*>(fAODJets->At(ij)); | |
728 | if(!jet) continue; | |
729 | ||
730 | if( (fEvtSelMinJetPt==-1. || jet->Pt()>=fEvtSelMinJetPt) | |
731 | && (fEvtSelMaxJetPt==-1. || jet->Pt()<=fEvtSelMaxJetPt) | |
732 | && (jet->Eta()>=fEvtSelMinJetEta && jet->Eta()<=fEvtSelMaxJetEta) | |
733 | && (jet->Phi()>=fEvtSelMinJetPhi && jet->Phi()<=fEvtSelMaxJetPhi)){ | |
734 | useEntry = kTRUE; | |
735 | break; | |
736 | } | |
737 | } | |
738 | } | |
739 | ||
740 | // no selection | |
741 | if(fEvtSelecMode==kEventsAll){ | |
742 | useEntry = kTRUE; | |
743 | } | |
744 | ||
745 | } | |
746 | AliDebug(AliLog::kDebug,Form("Use entry %d from extra AOD.", fAODEntry)); | |
747 | ||
748 | fh2PtHardEvtSel->Fill(fPtHardBin,fPtHard); | |
749 | fh2AODevent->Fill(fFileId,fAODEntry); | |
750 | ||
751 | TClonesArray *mcpartIN = (TClonesArray*)(fAODevent->FindListObject("mcparticles")); | |
752 | TClonesArray *mcpartOUT = 0x0; | |
753 | if(mcpartIN){ | |
754 | mcpartOUT = (TClonesArray*)(fAODout->FindListObject(fMCparticlesBranch.Data())); | |
755 | if(!mcpartOUT){ | |
756 | AliError("Extra MC particles branch not found in output."); | |
757 | PostData(1, fHistList); | |
758 | return; | |
759 | } | |
760 | mcpartOUT->Delete(); | |
761 | } else { | |
762 | AliInfo("No extra MC particles found."); | |
763 | } | |
764 | ||
765 | ||
766 | if(fEmbedMode==kAODFull || fEmbedMode==kAODJetTracks){ // take all tracks or jet tracks | |
767 | // loop over input tracks | |
768 | // add to output aod | |
769 | Int_t nTracks = 0; | |
770 | Int_t nJets = fAODJets->GetEntries(); | |
771 | Int_t nSelectedJets = 0; | |
772 | AliAODJet *leadJet = 0x0; // used for jet tracks only | |
773 | ||
774 | if(fEmbedMode==kAODFull){ | |
775 | nTracks = fAODevent->GetNumberOfTracks(); | |
776 | ||
777 | for(Int_t ij=0; ij<nJets; ++ij){ | |
778 | AliAODJet *jet = dynamic_cast<AliAODJet*>(fAODJets->At(ij)); | |
779 | if(!jet) continue; | |
780 | if( (fEvtSelMinJetPt==-1. || jet->Pt()>=fEvtSelMinJetPt) | |
781 | && (fEvtSelMaxJetPt==-1. || jet->Pt()<=fEvtSelMaxJetPt) | |
782 | && (jet->Eta()>=fEvtSelMinJetEta && jet->Eta()<=fEvtSelMaxJetEta) | |
783 | && (jet->Phi()>=fEvtSelMinJetPhi && jet->Eta()<=fEvtSelMaxJetPhi)){ | |
784 | ||
785 | fh1JetPt->Fill(jet->Pt()); | |
786 | fh2JetEtaPhi->Fill(jet->Eta(), jet->Phi()); | |
787 | nSelectedJets++; | |
788 | ||
789 | } | |
790 | } | |
791 | fh1JetN->Fill(nSelectedJets); | |
792 | } | |
793 | ||
794 | if(fEmbedMode==kAODJetTracks){ | |
795 | // look for leading jet within selection | |
796 | // get reference tracks for this jet | |
797 | Float_t leadJetPt = 0.; | |
798 | for(Int_t ij=0; ij<nJets; ++ij){ | |
799 | AliAODJet *jet = dynamic_cast<AliAODJet*>(fAODJets->At(ij)); | |
800 | if(!jet) continue; | |
801 | if( (fEvtSelMinJetPt==-1. || jet->Pt()>=fEvtSelMinJetPt) | |
802 | && (fEvtSelMaxJetPt==-1. || jet->Pt()<=fEvtSelMaxJetPt) | |
803 | && (jet->Eta()>=fEvtSelMinJetEta && jet->Eta()<=fEvtSelMaxJetEta) | |
804 | && (jet->Phi()>=fEvtSelMinJetPhi && jet->Eta()<=fEvtSelMaxJetPhi)){ | |
805 | if(jet->Pt()>leadJetPt){ | |
806 | leadJetPt = jet->Pt(); | |
807 | leadJet = jet; | |
808 | } | |
809 | } | |
810 | } | |
811 | if(leadJet){ | |
812 | nTracks = leadJet->GetRefTracks()->GetEntriesFast(); | |
813 | fh1JetPt->Fill(leadJet->Pt()); | |
814 | fh2JetEtaPhi->Fill(leadJet->Eta(), leadJet->Pt()); | |
815 | fh1JetN->Fill(1); | |
816 | } | |
817 | } | |
818 | ||
819 | fh1TrackN->Fill((Float_t)nTracks); | |
820 | ||
821 | for(Int_t it=0; it<nTracks; ++it){ | |
822 | AliAODTrack *tmpTr = 0x0; | |
823 | if(fEmbedMode==kAODFull) tmpTr = fAODevent->GetTrack(it); | |
824 | if(fEmbedMode==kAODJetTracks) tmpTr = dynamic_cast<AliAODTrack*>(leadJet->GetRefTracks()->At(it)); | |
825 | if(!tmpTr) continue; | |
826 | ||
827 | tmpTr->SetStatus(AliESDtrack::kEmbedded); | |
828 | ||
829 | new ((*tracks)[nAODtracks++]) AliAODTrack(*tmpTr); | |
830 | dummy = (*tracks)[nAODtracks-1]; | |
831 | ||
832 | if(fTrackFilterMap<=0 || tmpTr->TestFilterBit(fTrackFilterMap)){ | |
833 | fh1TrackPt->Fill(tmpTr->Pt()); | |
834 | fh2TrackEtaPhi->Fill(tmpTr->Eta(), tmpTr->Phi()); | |
835 | } | |
836 | } | |
837 | ||
838 | if(mcpartIN){ | |
839 | Int_t nMCpart = mcpartIN->GetEntriesFast(); | |
840 | ||
841 | Int_t nPhysicalPrimary=0; | |
842 | Int_t nAODmcpart=0; | |
843 | for(Int_t ip=0; ip<nMCpart; ++ip){ | |
844 | AliAODMCParticle *tmpPart = (AliAODMCParticle*) mcpartIN->At(ip); | |
845 | ||
846 | if(fEmbedMode==kAODJetTracks){ | |
847 | // jet track? else continue | |
848 | // not implemented yet | |
849 | continue; | |
850 | } | |
851 | ||
852 | new((*mcpartOUT)[nAODmcpart++]) AliAODMCParticle(*tmpPart); | |
853 | dummy = (*mcpartOUT)[nAODmcpart-1]; | |
854 | ||
855 | if(tmpPart->IsPhysicalPrimary() && tmpPart->Charge()!=0. && tmpPart->Charge()!=-99. ){ | |
856 | fh1MCTrackPt->Fill(tmpPart->Pt()); | |
857 | fh2MCTrackEtaPhi->Fill(tmpPart->Eta(), tmpPart->Phi()); | |
858 | nPhysicalPrimary++; | |
859 | } | |
860 | } | |
861 | fh1MCTrackN->Fill((Float_t)nPhysicalPrimary); | |
862 | ||
863 | } | |
864 | } // end: embed all tracks or jet tracks | |
865 | ||
866 | ||
867 | if(fEmbedMode==kAODJet4Mom){ | |
868 | ||
869 | // loop over jets | |
870 | Int_t nJets = fAODJets->GetEntries(); | |
871 | fh1TrackN->Fill((Float_t)nJets); | |
872 | for(Int_t ij=0; ij<nJets; ++ij){ | |
873 | AliAODJet *jet = dynamic_cast<AliAODJet*>(fAODJets->At(ij)); | |
874 | if(!jet) continue; | |
875 | AliAODTrack *tmpTr = (AliAODTrack*)jet; | |
876 | tmpTr->SetFlags(AliESDtrack::kEmbedded); | |
877 | ||
878 | new ((*tracks)[nAODtracks++]) AliAODTrack(*tmpTr); | |
879 | dummy = (*tracks)[nAODtracks-1]; | |
880 | ||
881 | fh1TrackPt->Fill(tmpTr->Pt()); | |
882 | fh2TrackEtaPhi->Fill(tmpTr->Eta(), tmpTr->Phi()); | |
883 | } | |
884 | ||
885 | } // end: embed jets as 4-momenta | |
886 | ||
887 | ||
888 | } //end: embed mode with AOD | |
889 | ||
890 | ||
891 | // === embed mode with toy === | |
892 | if(fEmbedMode==kToyTracks){ | |
893 | Int_t nT = (Int_t)(rndm->Uniform(fToyMinNbOfTracks, fToyMaxNbOfTracks)+0.5); | |
894 | ||
895 | fh1TrackN->Fill((Float_t)nT); | |
896 | ||
897 | for(Int_t i=0; i<nT; ++i){ | |
898 | ||
899 | Double_t pt = -1.; | |
900 | if(fToyMinTrackPt!=fToyMaxTrackPt){ | |
901 | if(fToyDistributionTrackPt==0){ | |
902 | pt = rndm->Uniform(fToyMinTrackPt, fToyMaxTrackPt); | |
903 | } else { | |
904 | while(pt<fToyMinTrackPt||pt>fToyMaxTrackPt){ | |
905 | pt = rndm->Exp(fToyDistributionTrackPt); // no power law yet!! | |
906 | pt += fToyMinTrackPt; | |
907 | } | |
908 | } | |
909 | } else { | |
910 | pt = fToyMinTrackPt; | |
911 | } | |
912 | Double_t eta = rndm->Uniform(fToyMinTrackEta,fToyMaxTrackEta); | |
913 | Double_t phi = rndm->Uniform(fToyMinTrackPhi,fToyMaxTrackPhi); | |
914 | phi = TVector2::Phi_0_2pi(phi); | |
915 | ||
916 | if(fDebug) Printf("Add track #%d: pt %.2f, eta %.2f, phi %.2f", i, pt, eta, phi); | |
917 | ||
918 | Double_t theta = 2*TMath::ATan(1/TMath::Exp(eta)); | |
919 | Float_t mom[3]; | |
920 | mom[0] = pt; | |
921 | mom[1] = phi; | |
922 | mom[2] = theta; | |
923 | ||
924 | Float_t xyz[3]; | |
925 | xyz[0] = -999.; | |
926 | xyz[1] = -999.; | |
927 | xyz[2] = -999.; | |
928 | ||
929 | AliAODTrack *tmpTr = new AliAODTrack( -999, // id | |
930 | -999, // label | |
931 | mom, // momentum[3] | |
932 | kFALSE, // cartesian | |
933 | xyz, // position | |
934 | kFALSE, // DCA | |
935 | NULL, // covMatrix, | |
936 | -99, // charge | |
937 | 0, // itsClusMap | |
938 | NULL, // pid | |
939 | NULL, // prodVertex | |
940 | kFALSE, // used for vertex fit | |
941 | kFALSE, // used for prim vtx fit | |
942 | AliAODTrack::kUndef, // type | |
943 | fToyFilterMap, // select info | |
944 | -999. // chi2 per NDF | |
945 | ); | |
946 | tmpTr->SetFlags(AliESDtrack::kEmbedded); | |
947 | ||
948 | new((*tracks)[nAODtracks++]) AliAODTrack(*tmpTr); | |
949 | dummy = (*tracks)[nAODtracks-1]; | |
950 | ||
951 | fh1TrackPt->Fill(pt); | |
952 | fh2TrackEtaPhi->Fill(eta,phi); | |
953 | ||
954 | delete tmpTr; | |
955 | } | |
956 | } //end: embed mode with toy | |
957 | ||
958 | ||
959 | PostData(1, fHistList); | |
960 | } | |
961 | ||
962 | ||
963 | //__________________________________________________________________________ | |
964 | void AliAnalysisTaskFastEmbedding::Terminate(Option_t *) | |
965 | { | |
b284a58f | 966 | // terminate |
d89b8229 | 967 | if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::Terminate()"); |
968 | ||
969 | if(fAODfile && fAODfile->IsOpen()) | |
970 | fAODfile->Close(); | |
971 | ||
972 | } | |
973 | ||
974 | //__________________________________________________________________________ | |
975 | Int_t AliAnalysisTaskFastEmbedding::GetJobID() | |
976 | { | |
b284a58f | 977 | // gives back the alien subjob id, if available, else -1 |
978 | ||
d89b8229 | 979 | Int_t id=-1; |
980 | ||
981 | const char* env = gSystem->Getenv("ALIEN_PROC_ID"); // GRID | |
982 | //if(!env || !strlen(env)) env = gSystem->Getenv("LSB_JOBINDEX"); // GSI | |
983 | ||
984 | if(env && strlen(env)){ | |
985 | id= atoi(env); | |
986 | AliInfo(Form("Job index %d", id)); | |
987 | } | |
988 | else{ | |
989 | AliInfo("Job index not found. Okay if running locally."); | |
990 | } | |
991 | ||
992 | return id; | |
993 | } | |
994 | ||
995 | //__________________________________________________________________________ | |
b284a58f | 996 | Int_t AliAnalysisTaskFastEmbedding::SelectAODfile() |
997 | { | |
998 | // select an AOD file from fAODPathArray | |
d89b8229 | 999 | |
1000 | Int_t nFiles = fAODPathArray->GetEntries(); | |
1001 | ||
1002 | // choose randomly file and event | |
1003 | Int_t n = -1; | |
1004 | Float_t tmpProp = -1.; | |
1005 | Int_t pendingEvents = fInputEntries-Entry(); | |
1006 | //Printf("input entries %d, entry %d, pending events %d", fInputEntries, (Int_t)Entry(), pendingEvents); | |
b284a58f | 1007 | if(fAODEntriesArray){ |
d89b8229 | 1008 | while(rndm->Rndm()>tmpProp){ |
1009 | n = rndm->Integer(nFiles); | |
1010 | fAODEntries = fAODEntriesArray->At(n); | |
1011 | Int_t tmpEntries = fAODEntries<pendingEvents ? pendingEvents : fAODEntries; | |
1012 | tmpProp = fAODEntriesMax ? (Float_t)tmpEntries/fAODEntriesMax : 1.; | |
1013 | } | |
1014 | fAODEntry = (Int_t)(rndm->Uniform(fAODEntries)); | |
b284a58f | 1015 | } else { |
1016 | AliWarning("Number of entries in extra AODs not set!"); | |
1017 | n = rndm->Integer(nFiles); | |
1018 | fAODEntry = 0; | |
1019 | } | |
d89b8229 | 1020 | |
1021 | AliInfo(Form("Select AOD file %d", n)); | |
1022 | TObjString *objStr = (TObjString*) fAODPathArray->At(n); | |
1023 | if(!objStr){ | |
1024 | AliError("Could not get path of aod file."); | |
1025 | return -1; | |
1026 | } | |
1027 | fAODPath = objStr->GetString(); | |
1028 | ||
1029 | return n; | |
1030 | } | |
1031 | ||
1032 | //__________________________________________________________________________ | |
b284a58f | 1033 | Int_t AliAnalysisTaskFastEmbedding::OpenAODfile(Int_t trial) |
1034 | { | |
1035 | // select and open an AOD file from fAODPathArray | |
d89b8229 | 1036 | |
1037 | if(fAODPathArray) fFileId = SelectAODfile(); | |
1038 | if(fFileId<0) return -1; | |
1039 | ||
1040 | TDirectory *owd = gDirectory; | |
b284a58f | 1041 | if (fAODfile && fAODfile->IsOpen()) fAODfile->Close(); |
d89b8229 | 1042 | fAODfile = TFile::Open(fAODPath.Data(),"TIMEOUT=180"); |
1043 | owd->cd(); | |
1044 | if(!fAODfile){ | |
1045 | ||
1046 | fFileId = -1; | |
1047 | if(fAODPathArray){ | |
1048 | if(trial<=3){ | |
1049 | AliError(Form("Could not open AOD file %s (trial %d), try again ...", fAODPath.Data(), trial)); | |
1050 | fFileId = OpenAODfile(trial+1); | |
1051 | } else { | |
1052 | AliError(Form("Could not open AOD file %s, give up ...",fAODPath.Data())); | |
1053 | return -1; | |
1054 | } | |
1055 | } else { | |
1056 | AliError(Form("Could not open AOD file %s.", fAODPath.Data())); | |
1057 | } | |
1058 | ||
1059 | return fFileId; | |
1060 | } | |
1061 | ||
1062 | fAODtree = (TTree*)fAODfile->Get("aodTree"); | |
1063 | ||
1064 | if(!fAODtree){ | |
1065 | AliError("AOD tree not found."); | |
1066 | return -1; | |
1067 | } | |
1068 | ||
1069 | /* | |
1070 | fAODtree->SetBranchStatus("*",0); | |
1071 | fAODtree->SetBranchStatus("header",1); | |
1072 | fAODtree->SetBranchStatus("tracks*",1); | |
1073 | fAODtree->SetBranchStatus("jets*",1); | |
1074 | fAODtree->SetBranchStatus("clusters*",1); | |
1075 | fAODtree->SetBranchStatus("mcparticles*",1); | |
1076 | fAODtree->SetBranchStatus("mcHeader*",1); | |
1077 | */ | |
1078 | ||
1079 | delete fAODevent; | |
1080 | fAODevent = new AliAODEvent(); | |
1081 | fAODevent->ReadFromTree(fAODtree); | |
1082 | ||
1083 | ||
1084 | // fetch header, jets, etc. from new file | |
1085 | fNevents = fAODtree->GetEntries(); | |
1086 | mcHeader = (AliAODMCHeader*)fAODevent->FindListObject(AliAODMCHeader::StdBranchName()); | |
1087 | ||
1088 | if(fJetBranch.Length()) fAODJets = dynamic_cast<TClonesArray*>(fAODevent->FindListObject(fJetBranch.Data())); | |
1089 | else fAODJets = fAODevent->GetJets(); | |
1090 | if(!fAODJets){ | |
1091 | AliError("Could not find jets in AOD. Check jet branch when indicated."); | |
1092 | return -1; | |
1093 | } | |
1094 | ||
1095 | TFile *curfile = fAODtree->GetCurrentFile(); | |
1096 | if (!curfile) { | |
1097 | AliError("No current file."); | |
1098 | return -1; | |
1099 | } | |
1100 | ||
1101 | Float_t trials = 1.; | |
1102 | Float_t xsec = 0.; | |
1103 | PythiaInfoFromFile(curfile->GetName(),xsec,trials); | |
1104 | fXsection = xsec; | |
1105 | ||
1106 | // construct a poor man average trials | |
1107 | Float_t nEntries = (Float_t)fAODtree->GetTree()->GetEntries(); | |
1108 | if(trials>=nEntries && nEntries>0.)fAvgTrials = trials/nEntries; | |
1109 | ||
1110 | if(fFileId>=0){ | |
1111 | AliInfo(Form("Read successfully AOD event from file %d",fFileId)); | |
1112 | } | |
1113 | ||
1114 | fCountEvents=0; // new file, reset counter | |
1115 | ||
1116 | return fFileId; // file position in AOD path array, if array available | |
1117 | } | |
1118 | ||
1119 | ||
1120 | //____________________________________________________________________________ | |
1121 | Float_t AliAnalysisTaskFastEmbedding::GetPtHard(Bool_t bSet, Float_t newValue){ | |
1122 | ||
b284a58f | 1123 | // static stored, available for other tasks in train |
1124 | ||
d89b8229 | 1125 | static Float_t ptHard = -1.; |
1126 | if(bSet) ptHard = newValue; | |
1127 | ||
1128 | return ptHard; | |
1129 | } | |
1130 | ||
1131 | //____________________________________________________________________________ | |
1132 | Int_t AliAnalysisTaskFastEmbedding::GetPtHardBin(Double_t ptHard){ | |
1133 | ||
b284a58f | 1134 | // returns pt hard bin (for LHC10e14, LHC11a1x, LHC11a2x simulations) |
1135 | ||
d89b8229 | 1136 | const Int_t nBins = 10; |
1137 | Double_t binLimits[nBins] = { 5., 11., 21., 36., 57., 84., 117., 156., 200., 249. }; // lower limits | |
1138 | ||
1139 | Int_t bin = -1; | |
1140 | while(bin<nBins-1 && binLimits[bin+1]<ptHard){ | |
1141 | bin++; | |
1142 | } | |
1143 | ||
1144 | return bin; | |
1145 | } | |
1146 | ||
1147 | //____________________________________________________________________________ | |
1148 | Bool_t AliAnalysisTaskFastEmbedding::PythiaInfoFromFile(const char* currFile,Float_t &fXsec,Float_t &fTrials){ | |
1149 | // | |
1150 | // get the cross section and the trails either from pyxsec.root or from pysec_hists.root | |
1151 | // This should provide the path to the AOD/ESD file | |
1152 | // (taken from PWG4/JetTasks/AliAnalysisHelperJetTasks) | |
1153 | ||
1154 | TString file(currFile); | |
1155 | fXsec = 0; | |
1156 | fTrials = 1; | |
1157 | ||
1158 | if(file.Contains("root_archive.zip#")){ | |
1159 | Ssiz_t pos1 = file.Index("root_archive",12,TString::kExact); | |
1160 | Ssiz_t pos = file.Index("#",1,pos1,TString::kExact); | |
1161 | file.Replace(pos+1,20,""); | |
1162 | } | |
1163 | else { | |
1164 | // not an archive take the basename.... | |
1165 | file.ReplaceAll(gSystem->BaseName(file.Data()),""); | |
1166 | } | |
1167 | ||
1168 | TFile *fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec.root")); // problem that we cannot really test the existance of a file in a archive so we have to lvie with open error message from root | |
1169 | if(!fxsec){ | |
1170 | // next trial fetch the histgram file | |
1171 | fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec_hists.root")); | |
1172 | if(!fxsec){ | |
1173 | // not a severe condition but inciate that we have no information | |
1174 | AliDebug(AliLog::kDebug,Form("Neither pyxsec.root nor pyxsec_hists.root found in %s",file.Data())); | |
1175 | return kFALSE; | |
1176 | } | |
1177 | else{ | |
1178 | // find the tlist we want to be independtent of the name so use the Tkey | |
1179 | TKey* key = (TKey*)fxsec->GetListOfKeys()->At(0); | |
1180 | if(!key){ | |
1181 | fxsec->Close(); | |
1182 | return kFALSE; | |
1183 | } | |
1184 | TList *list = dynamic_cast<TList*>(key->ReadObj()); | |
1185 | if(!list){ | |
1186 | fxsec->Close(); | |
1187 | return kFALSE; | |
1188 | } | |
1189 | fXsec = ((TProfile*)list->FindObject("h1Xsec"))->GetBinContent(1); | |
1190 | fTrials = ((TH1F*)list->FindObject("h1Trials"))->GetBinContent(1); | |
1191 | fxsec->Close(); | |
1192 | } | |
1193 | } // no tree pyxsec.root | |
1194 | else { | |
1195 | TTree *xtree = (TTree*)fxsec->Get("Xsection"); | |
1196 | if(!xtree){ | |
1197 | fxsec->Close(); | |
1198 | return kFALSE; | |
1199 | } | |
1200 | UInt_t ntrials = 0; | |
1201 | Double_t xsection = 0; | |
1202 | xtree->SetBranchAddress("xsection",&xsection); | |
1203 | xtree->SetBranchAddress("ntrials",&ntrials); | |
1204 | xtree->GetEntry(0); | |
1205 | fTrials = ntrials; | |
1206 | fXsec = xsection; | |
1207 | fxsec->Close(); | |
1208 | } | |
1209 | return kTRUE; | |
1210 | } | |
1211 | ||
1212 |