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