]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - JETAN/AliAnalysisTaskFastEmbedding.cxx
Coding rule violations corrected
[u/mrichter/AliRoot.git] / JETAN / AliAnalysisTaskFastEmbedding.cxx
... / ...
CommitLineData
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
37
38#include "AliAnalysisTaskFastEmbedding.h"
39#include "AliAnalysisManager.h"
40#include "AliESDtrack.h"
41#include "AliAODEvent.h"
42#include "AliAODTrack.h"
43#include "AliAODJet.h"
44#include "AliAODMCParticle.h"
45
46#include "AliLog.h"
47
48ClassImp(AliAnalysisTaskFastEmbedding)
49
50//__________________________________________________________________________
51AliAnalysisTaskFastEmbedding::AliAnalysisTaskFastEmbedding()
52 : AliAnalysisTaskSE()
53 ,fAODout(0)
54 ,fAODevent(0)
55 ,fAODtree(0)
56 ,fAODfile(0)
57 ,rndm(0)
58 ,fAODPathArray(0)
59 ,fAODPath("AliAOD.root")
60 ,fTrackBranch("aodExtraTracks")
61 ,fMCparticlesBranch("aodExtraMCparticles")
62 ,fJetBranch("")
63 ,fEntry(0)
64 ,fEmbedMode(0)
65 ,fEvtSelecMode(0)
66 ,fEvtSelMinJetPt(-1)
67 ,fEvtSelMaxJetPt(-1)
68 ,fEvtSelMinJetEta(-999.)
69 ,fEvtSelMaxJetEta( 999.)
70 ,fEvtSelMinJetPhi(0.)
71 ,fEvtSelMaxJetPhi(TMath::Pi()*2.)
72 ,fToyMinNbOfTracks(1)
73 ,fToyMaxNbOfTracks(1)
74 ,fToyMinTrackPt(50.)
75 ,fToyMaxTrackPt(50.)
76 ,fToyDistributionTrackPt(0.)
77 ,fToyMinTrackEta(-.5)
78 ,fToyMaxTrackEta(.5)
79 ,fToyMinTrackPhi(0.)
80 ,fToyMaxTrackPhi(2*TMath::Pi())
81 ,fToyFilterMap(0)
82 ,fHistList(0)
83 ,fh1TrackPt(0)
84 ,fh2TrackEtaPhi(0)
85 ,fh1TrackN(0)
86 ,fh1JetPt(0)
87 ,fh2JetEtaPhi(0)
88 ,fh1JetN(0)
89 ,fh1MCTrackPt(0)
90 ,fh2MCTrackEtaPhi(0)
91 ,fh1MCTrackN(0)
92 ,fh1AODfile(0)
93
94
95{
96 // default constructor
97
98}
99
100
101//__________________________________________________________________________
102AliAnalysisTaskFastEmbedding::AliAnalysisTaskFastEmbedding(const char *name)
103 : AliAnalysisTaskSE(name)
104 ,fAODout(0)
105 ,fAODevent(0)
106 ,fAODtree(0)
107 ,fAODfile(0)
108 ,rndm(0)
109 ,fAODPathArray(0)
110 ,fAODPath("AliAOD.root")
111 ,fTrackBranch("aodExtraTracks")
112 ,fMCparticlesBranch("aodExtraMCparticles")
113 ,fJetBranch("")
114 ,fEntry(0)
115 ,fEmbedMode(0)
116 ,fEvtSelecMode(0)
117 ,fEvtSelMinJetPt(-1)
118 ,fEvtSelMaxJetPt(-1)
119 ,fEvtSelMinJetEta(-999.)
120 ,fEvtSelMaxJetEta( 999.)
121 ,fEvtSelMinJetPhi(0.)
122 ,fEvtSelMaxJetPhi(TMath::Pi()*2.)
123 ,fToyMinNbOfTracks(1)
124 ,fToyMaxNbOfTracks(1)
125 ,fToyMinTrackPt(50.)
126 ,fToyMaxTrackPt(50.)
127 ,fToyDistributionTrackPt(0.)
128 ,fToyMinTrackEta(-.5)
129 ,fToyMaxTrackEta(.5)
130 ,fToyMinTrackPhi(0.)
131 ,fToyMaxTrackPhi(2*TMath::Pi())
132 ,fToyFilterMap(0)
133 ,fHistList(0)
134 ,fh1TrackPt(0)
135 ,fh2TrackEtaPhi(0)
136 ,fh1TrackN(0)
137 ,fh1JetPt(0)
138 ,fh2JetEtaPhi(0)
139 ,fh1JetN(0)
140 ,fh1MCTrackPt(0)
141 ,fh2MCTrackEtaPhi(0)
142 ,fh1MCTrackN(0)
143 ,fh1AODfile(0)
144{
145 // constructor
146 DefineOutput(1, TList::Class());
147}
148
149
150//__________________________________________________________________________
151AliAnalysisTaskFastEmbedding::AliAnalysisTaskFastEmbedding(const AliAnalysisTaskFastEmbedding &copy)
152 : AliAnalysisTaskSE()
153 ,fAODout(copy.fAODout)
154 ,fAODevent(copy.fAODevent)
155 ,fAODtree(copy.fAODtree)
156 ,fAODfile(copy.fAODfile)
157 ,rndm(copy.rndm)
158 ,fAODPathArray(copy.fAODPathArray)
159 ,fAODPath(copy.fAODPath)
160 ,fTrackBranch(copy.fTrackBranch)
161 ,fMCparticlesBranch(copy.fMCparticlesBranch)
162 ,fJetBranch(copy.fJetBranch)
163 ,fEntry(copy.fEntry)
164 ,fEmbedMode(copy.fEmbedMode)
165 ,fEvtSelecMode(copy.fEvtSelecMode)
166 ,fEvtSelMinJetPt(copy.fEvtSelMinJetPt)
167 ,fEvtSelMaxJetPt(copy.fEvtSelMaxJetPt)
168 ,fEvtSelMinJetEta(copy.fEvtSelMinJetEta)
169 ,fEvtSelMaxJetEta(copy.fEvtSelMaxJetEta)
170 ,fEvtSelMinJetPhi(copy.fEvtSelMinJetPhi)
171 ,fEvtSelMaxJetPhi(copy.fEvtSelMaxJetPhi)
172 ,fToyMinNbOfTracks(copy.fToyMinNbOfTracks)
173 ,fToyMaxNbOfTracks(copy.fToyMaxNbOfTracks)
174 ,fToyMinTrackPt(copy.fToyMinTrackPt)
175 ,fToyMaxTrackPt(copy.fToyMaxTrackPt)
176 ,fToyDistributionTrackPt(copy.fToyDistributionTrackPt)
177 ,fToyMinTrackEta(copy.fToyMinTrackEta)
178 ,fToyMaxTrackEta(copy.fToyMaxTrackEta)
179 ,fToyMinTrackPhi(copy.fToyMinTrackPhi)
180 ,fToyMaxTrackPhi(copy.fToyMaxTrackPhi)
181 ,fToyFilterMap(copy.fToyFilterMap)
182 ,fHistList(copy.fHistList)
183 ,fh1TrackPt(copy.fh1TrackPt)
184 ,fh2TrackEtaPhi(copy.fh2TrackEtaPhi)
185 ,fh1TrackN(copy.fh1TrackN)
186 ,fh1JetPt(copy.fh1JetPt)
187 ,fh2JetEtaPhi(copy.fh2JetEtaPhi)
188 ,fh1JetN(copy.fh1JetN)
189 ,fh1MCTrackPt(copy.fh1MCTrackPt)
190 ,fh2MCTrackEtaPhi(copy.fh2MCTrackEtaPhi)
191 ,fh1MCTrackN(copy.fh1MCTrackN)
192 ,fh1AODfile(copy.fh1AODfile)
193{
194 // copy constructor
195}
196
197
198//__________________________________________________________________________
199AliAnalysisTaskFastEmbedding& AliAnalysisTaskFastEmbedding::operator=(const AliAnalysisTaskFastEmbedding& o)
200{
201 // assignment
202
203 if(this!=&o){
204 AliAnalysisTaskSE::operator=(o);
205 fAODout = o.fAODout;
206 fAODevent = o.fAODevent;
207 fAODtree = o.fAODtree;
208 fAODfile = o.fAODfile;
209 rndm = o.rndm;
210 fAODPathArray = o.fAODPathArray;
211 fAODPath = o.fAODPath;
212 fTrackBranch = o.fTrackBranch;
213 fMCparticlesBranch = o.fMCparticlesBranch;
214 fJetBranch = o.fJetBranch;
215 fEntry = o.fEntry;
216 fEmbedMode = o.fEmbedMode;
217 fEvtSelecMode = o.fEvtSelecMode;
218 fEvtSelMinJetPt = o.fEvtSelMinJetPt;
219 fEvtSelMaxJetPt = o.fEvtSelMaxJetPt;
220 fEvtSelMinJetEta = o.fEvtSelMinJetEta;
221 fEvtSelMaxJetEta = o.fEvtSelMaxJetEta;
222 fEvtSelMinJetPhi = o.fEvtSelMinJetPhi;
223 fEvtSelMaxJetPhi = o.fEvtSelMaxJetPhi;
224 fToyMinNbOfTracks = o.fToyMinNbOfTracks;
225 fToyMaxNbOfTracks = o.fToyMaxNbOfTracks;
226 fToyMinTrackPt = o.fToyMinTrackPt;
227 fToyMaxTrackPt = o.fToyMaxTrackPt;
228 fToyDistributionTrackPt = o.fToyDistributionTrackPt;
229 fToyMinTrackEta = o.fToyMinTrackEta;
230 fToyMaxTrackEta = o.fToyMaxTrackEta;
231 fToyMinTrackPhi = o.fToyMinTrackPhi;
232 fToyMaxTrackPhi = o.fToyMaxTrackPhi;
233 fToyFilterMap = o.fToyFilterMap;
234 fHistList = o.fHistList;
235 fh1TrackPt = o.fh1TrackPt;
236 fh2TrackEtaPhi = o.fh2TrackEtaPhi;
237 fh1TrackN = o.fh1TrackN;
238 fh1JetPt = o.fh1JetPt;
239 fh2JetEtaPhi = o.fh2JetEtaPhi;
240 fh1JetN = o.fh1JetN;
241 fh1MCTrackPt = o.fh1MCTrackPt;
242 fh2MCTrackEtaPhi = o.fh2MCTrackEtaPhi;
243 fh1MCTrackN = o.fh1MCTrackN;
244 fh1AODfile = o.fh1AODfile;
245 }
246
247 return *this;
248}
249
250
251//__________________________________________________________________________
252AliAnalysisTaskFastEmbedding::~AliAnalysisTaskFastEmbedding()
253{
254 // destructor
255 delete rndm;
256}
257
258
259//__________________________________________________________________________
260void AliAnalysisTaskFastEmbedding::UserCreateOutputObjects()
261{
262 // create output objects
263 if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::UserCreateOutputObjects()");
264
265 // connect output aod
266 // create a new branch for extra tracks
267 fAODout = AODEvent();
268 if(!fAODout){
269 AliError("Output AOD not found.");
270 return;
271 }
272 if(!fAODout->FindListObject(fTrackBranch.Data()) && strlen(fTrackBranch.Data())){
273 AliInfo(Form("Add AOD branch %s", fTrackBranch.Data()));
274 TClonesArray *tracks = new TClonesArray("AliAODTrack",0);
275 tracks->SetName(fTrackBranch.Data());
276 AddAODBranch("TClonesArray", &tracks);
277 }
278 // create new branch for extra mcparticle if available as input
279 if(fAODevent && fAODevent->FindListObject("mcparticles") && strlen(fMCparticlesBranch.Data())){
280 AliInfo(Form("Add AOD branch %s", fMCparticlesBranch.Data()));
281 TClonesArray *mcparticles = new TClonesArray("AliAODMCParticle",0);
282 mcparticles->SetName(fMCparticlesBranch.Data());
283 AddAODBranch("TClonesArray", &mcparticles);
284 }
285
286
287
288
289 //qa histograms
290
291 OpenFile(1);
292 if(!fHistList) fHistList = new TList();
293 fHistList->SetOwner(kTRUE);
294
295 Bool_t oldStatus = TH1::AddDirectoryStatus();
296 TH1::AddDirectory(kFALSE);
297
298 fh1TrackPt = new TH1F("fh1TrackPt","pT of extra tracks;p_{T};entries", 250, 0., 250.);
299 fh2TrackEtaPhi = new TH2F("fh2TrackEtaPhi","eta-phi distribution of extra tracks;#eta;#phi", 20, -1., 1., 60, 0., 2*TMath::Pi());
300 fh1TrackN = new TH1F("fh1TrackN", "nb. of extra tracks per event;nb. of tracks;entries",300, 0., 300.);
301
302 fHistList->Add(fh1TrackPt);
303 fHistList->Add(fh2TrackEtaPhi);
304 fHistList->Add(fh1TrackN);
305
306 if(fEmbedMode==kAODFull || fEmbedMode==kAODJetTracks || fEmbedMode==kAODJet4Mom){
307
308 fh1JetPt = new TH1F("fh1JetPt", "pT of extra jets;p_{T};entries", 250, 0., 250.);
309 fh2JetEtaPhi = new TH2F("fh2JetEtaPhi", "eta-phi distribution of extra jets;#eta;#phi", 20, -1., 1., 60, 0., 2*TMath::Pi());
310 fh1JetN = new TH1F("fh1JetN", "nb. of extra jets per event;nb. of jets;entries",20,0.,20.);
311
312 fHistList->Add(fh1JetPt);
313 fHistList->Add(fh2JetEtaPhi);
314 fHistList->Add(fh1JetN);
315 }
316
317
318 if(fAODevent && fAODevent->FindListObject("mcparticles") && strlen(fMCparticlesBranch.Data())){
319
320 fh1MCTrackPt = new TH1F("fh1MCTrackPt","pT of MC extra tracks;p_{T};entries", 250, 0., 250.);
321 fh2MCTrackEtaPhi = new TH2F("fh2MCTrackEtaPhi","eta-phi distribution of MC extra tracks;#eta;#phi", 20, -1., 1., 60, 0., 2*TMath::Pi());
322 fh1MCTrackN = new TH1F("fh1MCTrackN", "nb. of MC extra tracks per event;nb. of tracks;entries",300, 0., 300.);
323
324 fHistList->Add(fh1MCTrackPt);
325 fHistList->Add(fh2MCTrackEtaPhi);
326 fHistList->Add(fh1MCTrackN);
327
328 }
329
330 fh1AODfile = new TH1I("fh1AODfile", "overview of opened AOD files from the array", 100, 0, 100);
331 fHistList->Add(fh1AODfile);
332
333 // =========== Switch on Sumw2 for all histos ===========
334 for (Int_t i=0; i<fHistList->GetEntries(); ++i) {
335 TH1 *h1 = dynamic_cast<TH1*>(fHistList->At(i));
336 if (h1){
337 h1->Sumw2();
338 continue;
339 }
340 }
341
342 TH1::AddDirectory(oldStatus);
343
344
345 // set seed
346 rndm = new TRandom3();
347 Int_t id = GetJobID();
348 if(id>-1) rndm->SetSeed(id);
349 else rndm->SetSeed(); // a TTUID is generated and used for seed
350 AliInfo(Form("TRandom3 seed: %d", rndm->GetSeed()));
351
352 // embed mode with AOD
353 if(fEmbedMode==kAODFull || fEmbedMode==kAODJetTracks || fEmbedMode==kAODJet4Mom){
354
355 // open input AOD
356 Int_t rc = OpenAODfile();
357 if(rc<0) return;
358 fh1AODfile->Fill(rc);
359
360 } //end: embed mode with AOD
361
362
363
364 PostData(1, fHistList);
365}
366
367
368//__________________________________________________________________________
369void AliAnalysisTaskFastEmbedding::Init()
370{
371 // Initialization
372 if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::Init()");
373
374}
375
376
377//__________________________________________________________________________
378void AliAnalysisTaskFastEmbedding::UserExec(Option_t *)
379{
380 if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::UserExec()");
381
382 if(!fAODout){
383 AliError("Need output AOD, but is not connected.");
384 PostData(1, fHistList);
385 return;
386 }
387
388 // connect aod out
389 TClonesArray *tracks = (TClonesArray*)(fAODout->FindListObject(fTrackBranch.Data()));
390 if(!tracks){
391 AliError("Extra track branch not found in output.");
392 PostData(1, fHistList);
393 return;
394 }
395 tracks->Delete();
396 Int_t nAODtracks=0;
397
398
399 TRef dummy;
400
401 // === embed mode with AOD ===
402 if(fEmbedMode==kAODFull || fEmbedMode==kAODJetTracks || fEmbedMode==kAODJet4Mom){
403 if(!fAODevent){
404 AliError("Need input AOD, but is not connected.");
405 PostData(1, fHistList);
406 return;
407 }
408
409 // fetch jets
410 TClonesArray *aodJets = 0;
411 if(fJetBranch.Length()) aodJets = dynamic_cast<TClonesArray*>(fAODevent->FindListObject(fJetBranch.Data()));
412 else aodJets = fAODevent->GetJets();
413 if(!aodJets){
414 AliError("Could not find jets in AOD. Check jet branch when indicated.");
415 PostData(1, fHistList);
416 return;
417 }
418
419
420 Int_t nEvents = fAODtree->GetEntries();
421
422 Bool_t useEntry = kFALSE;
423 while(!useEntry){ // protection need, if no event fulfills requierment
424 if(fEntry>=nEvents){
425 fEntry=0;
426 if(!fAODPathArray){
427 AliDebug(AliLog::kDebug, "Last event in AOD reached, start from entry 0 again.");
428 }
429 else {
430 AliDebug(AliLog::kDebug, "Last event in AOD reached, select new AOD file ...");
431
432 Int_t rc = OpenAODfile();
433 if(rc<0) {
434 PostData(1, fHistList);
435 return;
436 }
437 fh1AODfile->Fill(rc);
438
439 // new file => we must use the new jet array
440 if(fJetBranch.Length()) aodJets = dynamic_cast<TClonesArray*>(fAODevent->FindListObject(fJetBranch.Data()));
441 else aodJets = fAODevent->GetJets();
442 if(!aodJets){
443 AliError("Could not find jets in AOD. Check jet branch when indicated.");
444 PostData(1, fHistList);
445 return;
446 }
447 }
448 }
449
450 fAODtree->GetEvent(fEntry);
451
452 // jet pt selection
453 if(fEvtSelecMode==kEventsJetPt){
454 Int_t nJets = aodJets->GetEntries();
455 for(Int_t ij=0; ij<nJets; ++ij){
456 AliAODJet *jet = dynamic_cast<AliAODJet*>(aodJets->At(ij));
457 if(!jet) continue;
458
459 if( (fEvtSelMinJetPt==-1. || jet->Pt()>=fEvtSelMinJetPt)
460 && (fEvtSelMaxJetPt==-1. || jet->Pt()<=fEvtSelMaxJetPt)
461 && (jet->Eta()>=fEvtSelMinJetEta && jet->Eta()<=fEvtSelMaxJetEta)
462 && (jet->Phi()>=fEvtSelMinJetPhi && jet->Phi()<=fEvtSelMaxJetPhi)){
463 useEntry = kTRUE;
464 break;
465 }
466 }
467 }
468
469 // no selection
470 if(fEvtSelecMode==kEventsAll){
471 useEntry = kTRUE;
472 }
473
474 fEntry++;
475 }
476 AliInfo(Form("Use entry %d from extra AOD.", fEntry-1));
477
478
479 TClonesArray *mcpartIN = (TClonesArray*)(fAODevent->FindListObject("mcparticles"));
480 TClonesArray *mcpartOUT = 0x0;
481 if(mcpartIN){
482 mcpartOUT = (TClonesArray*)(fAODout->FindListObject(fMCparticlesBranch.Data()));
483 mcpartOUT->Delete();
484 } else {
485 AliInfo("No extra MC particles found.");
486 }
487
488
489 if(fEmbedMode==kAODFull || fEmbedMode==kAODJetTracks){ // take all tracks or jet tracks
490 // loop over input tracks
491 // add to output aod
492 Int_t nTracks = 0;
493 Int_t nJets = aodJets->GetEntries();
494 Int_t nSelectedJets = 0;
495 AliAODJet *leadJet = 0x0; // used for jet tracks only
496
497 if(fEmbedMode==kAODFull){
498 nTracks = fAODevent->GetNumberOfTracks();
499
500 for(Int_t ij=0; ij<nJets; ++ij){
501 AliAODJet *jet = dynamic_cast<AliAODJet*>(aodJets->At(ij));
502 if(!jet) continue;
503 if( (fEvtSelMinJetPt==-1. || jet->Pt()>=fEvtSelMinJetPt)
504 && (fEvtSelMaxJetPt==-1. || jet->Pt()<=fEvtSelMaxJetPt)
505 && (jet->Eta()>=fEvtSelMinJetEta && jet->Eta()<=fEvtSelMaxJetEta)
506 && (jet->Phi()>=fEvtSelMinJetPhi && jet->Eta()<=fEvtSelMaxJetPhi)){
507
508 fh1JetPt->Fill(jet->Pt());
509 fh2JetEtaPhi->Fill(jet->Eta(), jet->Phi());
510 nSelectedJets++;
511
512 }
513 }
514 fh1JetN->Fill(nSelectedJets);
515 }
516
517 if(fEmbedMode==kAODJetTracks){
518 // look for leading jet within selection
519 // get reference tracks for this jet
520 Float_t leadJetPt = 0.;
521 for(Int_t ij=0; ij<nJets; ++ij){
522 AliAODJet *jet = dynamic_cast<AliAODJet*>(aodJets->At(ij));
523 if(!jet) continue;
524 if( (fEvtSelMinJetPt==-1. || jet->Pt()>=fEvtSelMinJetPt)
525 && (fEvtSelMaxJetPt==-1. || jet->Pt()<=fEvtSelMaxJetPt)
526 && (jet->Eta()>=fEvtSelMinJetEta && jet->Eta()<=fEvtSelMaxJetEta)
527 && (jet->Phi()>=fEvtSelMinJetPhi && jet->Eta()<=fEvtSelMaxJetPhi)){
528 if(jet->Pt()>leadJetPt){
529 leadJetPt = jet->Pt();
530 leadJet = jet;
531 }
532 }
533 }
534 if(leadJet){
535 nTracks = leadJet->GetRefTracks()->GetEntriesFast();
536 fh1JetPt->Fill(leadJet->Pt());
537 fh2JetEtaPhi->Fill(leadJet->Eta(), leadJet->Pt());
538 fh1JetN->Fill(1);
539 }
540 }
541
542 fh1TrackN->Fill((Float_t)nTracks);
543
544 for(Int_t it=0; it<nTracks; ++it){
545 AliAODTrack *tmpTr = 0x0;
546 if(fEmbedMode==kAODFull) tmpTr = fAODevent->GetTrack(it);
547 if(fEmbedMode==kAODJetTracks) tmpTr = dynamic_cast<AliAODTrack*>(leadJet->GetRefTracks()->At(it));
548 if(!tmpTr) continue;
549
550 tmpTr->SetStatus(AliESDtrack::kEmbedded);
551
552 new ((*tracks)[nAODtracks++]) AliAODTrack(*tmpTr);
553 dummy = (*tracks)[nAODtracks-1];
554
555 fh1TrackPt->Fill(tmpTr->Pt());
556 fh2TrackEtaPhi->Fill(tmpTr->Eta(), tmpTr->Phi());
557 }
558
559 if(mcpartIN){
560 Int_t nMCpart = mcpartIN->GetEntriesFast();
561
562 Int_t nPhysicalPrimary=0;
563 Int_t nAODmcpart=0;
564 for(Int_t ip=0; ip<nMCpart; ++ip){
565 AliAODMCParticle *tmpPart = (AliAODMCParticle*) mcpartIN->At(ip);
566
567 if(fEmbedMode==kAODJetTracks){
568 // jet track? else continue
569 // not implemented yet
570 continue;
571 }
572
573 new((*mcpartOUT)[nAODmcpart++]) AliAODMCParticle(*tmpPart);
574 dummy = (*mcpartOUT)[nAODmcpart-1];
575
576 if(tmpPart->IsPhysicalPrimary() && tmpPart->Charge()!=0. && tmpPart->Charge()!=-99. ){
577 fh1MCTrackPt->Fill(tmpPart->Pt());
578 fh2MCTrackEtaPhi->Fill(tmpPart->Eta(), tmpPart->Phi());
579 nPhysicalPrimary++;
580 }
581 }
582 fh1MCTrackN->Fill((Float_t)nPhysicalPrimary);
583
584 }
585 } // end: embed all tracks or jet tracks
586
587
588 if(fEmbedMode==kAODJet4Mom){
589
590 // loop over jets
591 Int_t nJets = aodJets->GetEntries();
592 fh1TrackN->Fill((Float_t)nJets);
593 for(Int_t ij=0; ij<nJets; ++ij){
594 AliAODJet *jet = dynamic_cast<AliAODJet*>(aodJets->At(ij));
595 if(!jet) continue;
596 AliAODTrack *tmpTr = (AliAODTrack*)jet;
597 tmpTr->SetFlags(AliESDtrack::kEmbedded);
598
599 new ((*tracks)[nAODtracks++]) AliAODTrack(*tmpTr);
600 dummy = (*tracks)[nAODtracks-1];
601
602 fh1TrackPt->Fill(tmpTr->Pt());
603 fh2TrackEtaPhi->Fill(tmpTr->Eta(), tmpTr->Phi());
604 }
605
606 } // end: embed jets as 4-momenta
607
608
609 } //end: embed mode with AOD
610
611
612 // === embed mode with toy ===
613 if(fEmbedMode==kToyTracks){
614 Int_t nT = (Int_t)(rndm->Uniform(fToyMinNbOfTracks, fToyMaxNbOfTracks)+0.5);
615
616 fh1TrackN->Fill((Float_t)nT);
617
618 for(Int_t i=0; i<nT; ++i){
619
620 Double_t pt = -1.;
621 if(fToyMinTrackPt!=fToyMaxTrackPt){
622 if(fToyDistributionTrackPt==0){
623 pt = rndm->Uniform(fToyMinTrackPt, fToyMaxTrackPt);
624 } else {
625 while(pt<fToyMinTrackPt||pt>fToyMaxTrackPt){
626 pt = rndm->Exp(fToyDistributionTrackPt); // no power law yet!!
627 pt += fToyMinTrackPt;
628 }
629 }
630 } else {
631 pt = fToyMinTrackPt;
632 }
633 Double_t eta = rndm->Uniform(fToyMinTrackEta,fToyMaxTrackEta);
634 Double_t phi = rndm->Uniform(fToyMinTrackPhi,fToyMaxTrackPhi);
635 phi = TVector2::Phi_0_2pi(phi);
636
637 if(fDebug) Printf("Add track #%d: pt %.2f, eta %.2f, phi %.2f", i, pt, eta, phi);
638
639 Double_t theta = 2*TMath::ATan(1/TMath::Exp(eta));
640 Float_t mom[3];
641 mom[0] = pt;
642 mom[1] = phi;
643 mom[2] = theta;
644
645 Float_t xyz[3];
646 xyz[0] = -999.;
647 xyz[1] = -999.;
648 xyz[2] = -999.;
649
650 AliAODTrack *tmpTr = new AliAODTrack( -999, // id
651 -999, // label
652 mom, // momentum[3]
653 kFALSE, // cartesian
654 xyz, // position
655 kFALSE, // DCA
656 NULL, // covMatrix,
657 -99, // charge
658 0, // itsClusMap
659 NULL, // pid
660 NULL, // prodVertex
661 kFALSE, // used for vertex fit
662 kFALSE, // used for prim vtx fit
663 AliAODTrack::kUndef, // type
664 fToyFilterMap, // select info
665 -999. // chi2 per NDF
666 );
667 tmpTr->SetFlags(AliESDtrack::kEmbedded);
668
669 new((*tracks)[nAODtracks++]) AliAODTrack(*tmpTr);
670 dummy = (*tracks)[nAODtracks-1];
671
672 fh1TrackPt->Fill(pt);
673 fh2TrackEtaPhi->Fill(eta,phi);
674
675 delete tmpTr;
676 }
677 } //end: embed mode with toy
678
679
680 PostData(1, fHistList);
681}
682
683
684//__________________________________________________________________________
685void AliAnalysisTaskFastEmbedding::Terminate(Option_t *)
686{
687 if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::Terminate()");
688
689 if(fAODfile && fAODfile->IsOpen())
690 fAODfile->Close();
691
692}
693
694//__________________________________________________________________________
695Int_t AliAnalysisTaskFastEmbedding::GetJobID()
696{
697 Int_t id=-1;
698
699 const char* env = gSystem->Getenv("ALIEN_PROC_ID"); // GRID
700 //if(!env || !strlen(env)) env = gSystem->Getenv("LSB_JOBINDEX"); // GSI
701
702 if(env && strlen(env)){
703 id= atoi(env);
704 AliInfo(Form("Job index %d", id));
705 }
706 else{
707 AliInfo("Job index not found. Okay if running locally.");
708 }
709
710 return id;
711}
712
713//__________________________________________________________________________
714
715Int_t AliAnalysisTaskFastEmbedding::SelectAODfile(){
716
717 Int_t nFiles = fAODPathArray->GetEntries();
718 Int_t n = rndm->Integer(nFiles);
719 AliInfo(Form("Select AOD file %d", n));
720 TObjString *objStr = (TObjString*) fAODPathArray->At(n);
721 if(!objStr){
722 AliError("Could not get path of aod file.");
723 return -1;
724 }
725 fAODPath = objStr->GetString();
726
727 return n;
728}
729
730//__________________________________________________________________________
731
732Int_t AliAnalysisTaskFastEmbedding::OpenAODfile(Int_t trial){
733
734 if(trial>3){
735 AliError("Could not open AOD files, give up ...");
736 return -1;
737 }
738
739 Int_t rc = 0;
740 if(fAODPathArray) rc = SelectAODfile();
741 if(rc<0) return -1;
742
743 TDirectory *owd = gDirectory;
744 if (fAODfile)
745 fAODfile->Close();
746 fAODfile = TFile::Open(fAODPath.Data());
747 owd->cd();
748 if(!fAODfile){
749
750 rc = -1;
751 if(fAODPathArray){
752 AliError(Form("Could not open AOD file %s, try another from the list ...", fAODPath.Data()));
753 rc = OpenAODfile(trial+1);
754 } else {
755 AliError(Form("Could not open AOD file %s.", fAODPath.Data()));
756 }
757
758 return rc;
759 }
760
761 fAODtree = (TTree*)fAODfile->Get("aodTree");
762
763 if(!fAODtree){
764 AliError("AOD tree not found.");
765 return -1;
766 }
767
768 delete fAODevent;
769 fAODevent = new AliAODEvent();
770 fAODevent->ReadFromTree(fAODtree);
771 return rc; // file position in AOD path array, if array available
772}
773