]> git.uio.no Git - u/mrichter/AliRoot.git/blame - JETAN/AliAnalysisTaskFastEmbedding.cxx
Add print on the found particle tag
[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>
b725dccf 38
39
40#include "AliAnalysisTaskFastEmbedding.h"
41#include "AliAnalysisManager.h"
31b9d515 42#include "AliESDEvent.h"
46465e39 43#include "AliESDtrack.h"
b725dccf 44#include "AliAODEvent.h"
45#include "AliAODTrack.h"
46#include "AliAODJet.h"
47#include "AliAODMCParticle.h"
31b9d515 48#include "AliAODMCHeader.h"
49#include "AliInputEventHandler.h"
b725dccf 50
51#include "AliLog.h"
52
53ClassImp(AliAnalysisTaskFastEmbedding)
54
55//__________________________________________________________________________
56AliAnalysisTaskFastEmbedding::AliAnalysisTaskFastEmbedding()
31b9d515 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)
b725dccf 133{
31b9d515 134 // default constructor
b725dccf 135
136}
137
138
139//__________________________________________________________________________
140AliAnalysisTaskFastEmbedding::AliAnalysisTaskFastEmbedding(const char *name)
31b9d515 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)
b725dccf 217{
31b9d515 218 // constructor
219 DefineOutput(1, TList::Class());
b725dccf 220}
221
222
223//__________________________________________________________________________
224AliAnalysisTaskFastEmbedding::AliAnalysisTaskFastEmbedding(const AliAnalysisTaskFastEmbedding &copy)
31b9d515 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)
b725dccf 301{
31b9d515 302 // copy constructor
b725dccf 303}
304
305
306//__________________________________________________________________________
307AliAnalysisTaskFastEmbedding& AliAnalysisTaskFastEmbedding::operator=(const AliAnalysisTaskFastEmbedding& o)
308{
31b9d515 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;
b725dccf 391}
392
393
394//__________________________________________________________________________
395AliAnalysisTaskFastEmbedding::~AliAnalysisTaskFastEmbedding()
396{
31b9d515 397 // destructor
398 delete rndm;
b725dccf 399}
400
401
402//__________________________________________________________________________
403void AliAnalysisTaskFastEmbedding::UserCreateOutputObjects()
404{
31b9d515 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);
b725dccf 537}
538
539
540//__________________________________________________________________________
541void AliAnalysisTaskFastEmbedding::Init()
542{
31b9d515 543 // Initialization
544 if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::Init()");
545
546}
547
548//__________________________________________________________________________
549Bool_t AliAnalysisTaskFastEmbedding::UserNotify()
550{
551
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;
b725dccf 560
b725dccf 561}
562
563
564//__________________________________________________________________________
565void AliAnalysisTaskFastEmbedding::UserExec(Option_t *)
566{
31b9d515 567 if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::UserExec()");
568
569 if(!fAODout){
570 AliError("Need output AOD, but is not connected.");
571 PostData(1, fHistList);
572 return;
573 }
574
575 fESD=dynamic_cast<AliESDEvent*>(InputEvent());
576 if (!fESD) {
577 AliError("ESD not available");
578 PostData(1, fHistList);
579 return;
580 }
581
582 // -- event selection --
583 fHistEvtSelection->Fill(1); // number of events before event selection
584
585 // physics selection
586 AliInputEventHandler* inputHandler = (AliInputEventHandler*)
587 ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
588 if(!(inputHandler->IsEventSelected() & fOfflineTrgMask)){
589 if(fDebug) Printf(" Trigger Selection: event REJECTED ... ");
590 fHistEvtSelection->Fill(2);
591 PostData(1, fHistList);
592 return;
593 }
594
595 // vertex selection
596 AliAODVertex* primVtx = fAODout->GetPrimaryVertex();
597 Int_t nTracksPrim = primVtx->GetNContributors();
598 if ((nTracksPrim < fMinContribVtx) ||
599 (primVtx->GetZ() < fVtxZMin) ||
600 (primVtx->GetZ() > fVtxZMax) ){
601 if(fDebug) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
602 fHistEvtSelection->Fill(3);
603 PostData(1, fHistList);
604 return;
605 }
606
607 /** takes wrong value, since jet service tasks runs later
608 // event class selection (from jet helper task)
609 Int_t eventClass = AliAnalysisHelperJetTasks::EventClass();
610 if(fDebug) Printf("Event class %d", eventClass);
611 if (eventClass < fEvtClassMin || eventClass > fEvtClassMax){
612 fHistEvtSelection->Fill(4);
613 PostData(1, fHistList);
614 return;
615 }*/
616
617 // centrality selection
618 AliCentrality *cent = 0x0;
619 Float_t centValue = 0.;
620 if(fESD) cent = fESD->GetCentrality();
621 if(cent) centValue = cent->GetCentralityPercentile("V0M");
622 if(fDebug) printf("centrality: %f\n", centValue);
623 if (centValue < fCentMin || centValue > fCentMax){
624 fHistEvtSelection->Fill(4);
625 PostData(1, fHistList);
626 return;
627 }
628
629
630 /* ** not implemented **
631 // multiplicity due to input tracks
632 Int_t nInputTracks = GetNInputTracks();
633
634 if (nInputTracks < fNInputTracksMin || (fNInputTracksMax > -1 && nInputTracks > fNInputTracksMax)){
635 fHistEvtSelection->Fill(5);
636 PostData(1, fHistList);
637 return;
638 }
639 */
640
641 fHistEvtSelection->Fill(0); // accepted events
642 // -- end event selection --
643
644 // connect aod out
645 TClonesArray *tracks = (TClonesArray*)(fAODout->FindListObject(fTrackBranch.Data()));
646 if(!tracks){
647 AliError("Extra track branch not found in output.");
648 PostData(1, fHistList);
649 return;
650 }
651 tracks->Delete();
652 Int_t nAODtracks=0;
653
654 TRef dummy;
655
656 // === embed mode with AOD ===
657 if(fEmbedMode==kAODFull || fEmbedMode==kAODJetTracks || fEmbedMode==kAODJet4Mom){
658 if(!fAODevent){
659 AliError("Need input AOD, but is not connected.");
aa117124 660 PostData(1, fHistList);
661 return;
31b9d515 662 }
663
664
665
666 Bool_t useEntry = kFALSE;
667 while(!useEntry){ // protection need, if no event fulfills requierment
668
669 fAODEntry++; // go to next event
670 fCountEvents++;
671 if(fAODEntry>=fNevents) fAODEntry=0;
672
673 // open new aod file
674 if(fCountEvents>=fNevents){
675 if(!fAODPathArray){
676 AliDebug(AliLog::kDebug, "No list of AODs, keep current AOD.");
677 }
678 else {
679 AliDebug(AliLog::kDebug, "Select new AOD file ...");
680
681 fFileId = OpenAODfile();
682 if(fFileId<0) {
683 PostData(1, fHistList);
684 return;
685 }
686 fh1AODfile->Fill(fFileId);
687 if(fAODEntries!=fNevents){
688 AliError("File with list of AODs and file with nb. of entries inconsistent");
689 PostData(1, fHistList);
690 return;
691 }
692 }
693 }
694
695
696
697 // get next event
698 fAODtree->GetEvent(fAODEntry);
699
700 // get pt hard
701 if(mcHeader){
702 fPtHard = mcHeader->GetPtHard();
703 GetPtHard(kTRUE,fPtHard); // make value available for other tasks (e.g. jet response task)
704 AliDebug(AliLog::kDebug,Form("pt hard %.2f", fPtHard));
705 }
706 else{
707 AliDebug(AliLog::kDebug,"No mcHeader");
708 }
709 fPtHardBin = GetPtHardBin(fPtHard);
710
711 //Printf("pt hard (bin) %.2f (%d), xSection %.2e, trials %.2f",fPtHard, fPtHardBin, fXsection, fAvgTrials);
712
713 // fill event variables for each event
714 fh1Xsec->Fill(fPtHardBin,fXsection);
715 fh2PtHard->Fill(fPtHardBin,fPtHard);
716
717 fh1Trials->Fill(fPtHardBin, fAvgTrials);
718 fh1TrialsEvtSel->Fill(fPtHardBin);
719 fh2PtHardTrials->Fill(fPtHardBin,fPtHard,fAvgTrials);
720
721 // jet pt selection
722 if(fEvtSelecMode==kEventsJetPt){
723 Int_t nJets = fAODJets->GetEntries();
724 //Printf("nb. jets: %d",nJets);
725 for(Int_t ij=0; ij<nJets; ++ij){
726 AliAODJet *jet = dynamic_cast<AliAODJet*>(fAODJets->At(ij));
727 if(!jet) continue;
728
729 if( (fEvtSelMinJetPt==-1. || jet->Pt()>=fEvtSelMinJetPt)
730 && (fEvtSelMaxJetPt==-1. || jet->Pt()<=fEvtSelMaxJetPt)
731 && (jet->Eta()>=fEvtSelMinJetEta && jet->Eta()<=fEvtSelMaxJetEta)
732 && (jet->Phi()>=fEvtSelMinJetPhi && jet->Phi()<=fEvtSelMaxJetPhi)){
733 useEntry = kTRUE;
734 break;
735 }
736 }
737 }
738
739 // no selection
740 if(fEvtSelecMode==kEventsAll){
741 useEntry = kTRUE;
742 }
743
744 }
745 AliDebug(AliLog::kDebug,Form("Use entry %d from extra AOD.", fAODEntry));
746
747 fh2PtHardEvtSel->Fill(fPtHardBin,fPtHard);
748 fh2AODevent->Fill(fFileId,fAODEntry);
749
750 TClonesArray *mcpartIN = (TClonesArray*)(fAODevent->FindListObject("mcparticles"));
751 TClonesArray *mcpartOUT = 0x0;
752 if(mcpartIN){
753 mcpartOUT = (TClonesArray*)(fAODout->FindListObject(fMCparticlesBranch.Data()));
754 if(!mcpartOUT){
755 AliError("Extra MC particles branch not found in output.");
756 PostData(1, fHistList);
757 return;
758 }
759 mcpartOUT->Delete();
760 } else {
761 AliInfo("No extra MC particles found.");
762 }
763
764
765 if(fEmbedMode==kAODFull || fEmbedMode==kAODJetTracks){ // take all tracks or jet tracks
766 // loop over input tracks
767 // add to output aod
768 Int_t nTracks = 0;
769 Int_t nJets = fAODJets->GetEntries();
770 Int_t nSelectedJets = 0;
771 AliAODJet *leadJet = 0x0; // used for jet tracks only
772
773 if(fEmbedMode==kAODFull){
774 nTracks = fAODevent->GetNumberOfTracks();
775
776 for(Int_t ij=0; ij<nJets; ++ij){
777 AliAODJet *jet = dynamic_cast<AliAODJet*>(fAODJets->At(ij));
778 if(!jet) continue;
779 if( (fEvtSelMinJetPt==-1. || jet->Pt()>=fEvtSelMinJetPt)
46465e39 780 && (fEvtSelMaxJetPt==-1. || jet->Pt()<=fEvtSelMaxJetPt)
781 && (jet->Eta()>=fEvtSelMinJetEta && jet->Eta()<=fEvtSelMaxJetEta)
782 && (jet->Phi()>=fEvtSelMinJetPhi && jet->Eta()<=fEvtSelMaxJetPhi)){
b725dccf 783
31b9d515 784 fh1JetPt->Fill(jet->Pt());
785 fh2JetEtaPhi->Fill(jet->Eta(), jet->Phi());
786 nSelectedJets++;
787
788 }
789 }
790 fh1JetN->Fill(nSelectedJets);
791 }
792
793 if(fEmbedMode==kAODJetTracks){
794 // look for leading jet within selection
795 // get reference tracks for this jet
796 Float_t leadJetPt = 0.;
797 for(Int_t ij=0; ij<nJets; ++ij){
798 AliAODJet *jet = dynamic_cast<AliAODJet*>(fAODJets->At(ij));
1a2bb6d5 799 if(!jet) continue;
31b9d515 800 if( (fEvtSelMinJetPt==-1. || jet->Pt()>=fEvtSelMinJetPt)
801 && (fEvtSelMaxJetPt==-1. || jet->Pt()<=fEvtSelMaxJetPt)
802 && (jet->Eta()>=fEvtSelMinJetEta && jet->Eta()<=fEvtSelMaxJetEta)
803 && (jet->Phi()>=fEvtSelMinJetPhi && jet->Eta()<=fEvtSelMaxJetPhi)){
804 if(jet->Pt()>leadJetPt){
805 leadJetPt = jet->Pt();
806 leadJet = jet;
807 }
808 }
809 }
810 if(leadJet){
811 nTracks = leadJet->GetRefTracks()->GetEntriesFast();
812 fh1JetPt->Fill(leadJet->Pt());
813 fh2JetEtaPhi->Fill(leadJet->Eta(), leadJet->Pt());
814 fh1JetN->Fill(1);
815 }
816 }
b725dccf 817
31b9d515 818 fh1TrackN->Fill((Float_t)nTracks);
b725dccf 819
31b9d515 820 for(Int_t it=0; it<nTracks; ++it){
821 AliAODTrack *tmpTr = 0x0;
822 if(fEmbedMode==kAODFull) tmpTr = fAODevent->GetTrack(it);
823 if(fEmbedMode==kAODJetTracks) tmpTr = dynamic_cast<AliAODTrack*>(leadJet->GetRefTracks()->At(it));
824 if(!tmpTr) continue;
825
826 tmpTr->SetStatus(AliESDtrack::kEmbedded);
827
828 new ((*tracks)[nAODtracks++]) AliAODTrack(*tmpTr);
829 dummy = (*tracks)[nAODtracks-1];
830
831 if(fTrackFilterMap<=0 || tmpTr->TestFilterBit(fTrackFilterMap)){
b725dccf 832 fh1TrackPt->Fill(tmpTr->Pt());
833 fh2TrackEtaPhi->Fill(tmpTr->Eta(), tmpTr->Phi());
31b9d515 834 }
835 }
836
837 if(mcpartIN){
838 Int_t nMCpart = mcpartIN->GetEntriesFast();
839
840 Int_t nPhysicalPrimary=0;
841 Int_t nAODmcpart=0;
842 for(Int_t ip=0; ip<nMCpart; ++ip){
843 AliAODMCParticle *tmpPart = (AliAODMCParticle*) mcpartIN->At(ip);
844
845 if(fEmbedMode==kAODJetTracks){
846 // jet track? else continue
847 // not implemented yet
848 continue;
849 }
850
851 new((*mcpartOUT)[nAODmcpart++]) AliAODMCParticle(*tmpPart);
852 dummy = (*mcpartOUT)[nAODmcpart-1];
853
854 if(tmpPart->IsPhysicalPrimary() && tmpPart->Charge()!=0. && tmpPart->Charge()!=-99. ){
855 fh1MCTrackPt->Fill(tmpPart->Pt());
856 fh2MCTrackEtaPhi->Fill(tmpPart->Eta(), tmpPart->Phi());
857 nPhysicalPrimary++;
858 }
859 }
860 fh1MCTrackN->Fill((Float_t)nPhysicalPrimary);
861
862 }
863 } // end: embed all tracks or jet tracks
864
865
866 if(fEmbedMode==kAODJet4Mom){
867
868 // loop over jets
869 Int_t nJets = fAODJets->GetEntries();
870 fh1TrackN->Fill((Float_t)nJets);
871 for(Int_t ij=0; ij<nJets; ++ij){
872 AliAODJet *jet = dynamic_cast<AliAODJet*>(fAODJets->At(ij));
873 if(!jet) continue;
874 AliAODTrack *tmpTr = (AliAODTrack*)jet;
875 tmpTr->SetFlags(AliESDtrack::kEmbedded);
876
877 new ((*tracks)[nAODtracks++]) AliAODTrack(*tmpTr);
878 dummy = (*tracks)[nAODtracks-1];
879
880 fh1TrackPt->Fill(tmpTr->Pt());
881 fh2TrackEtaPhi->Fill(tmpTr->Eta(), tmpTr->Phi());
882 }
883
884 } // end: embed jets as 4-momenta
885
886
887 } //end: embed mode with AOD
888
889
890 // === embed mode with toy ===
891 if(fEmbedMode==kToyTracks){
892 Int_t nT = (Int_t)(rndm->Uniform(fToyMinNbOfTracks, fToyMaxNbOfTracks)+0.5);
893
894 fh1TrackN->Fill((Float_t)nT);
895
896 for(Int_t i=0; i<nT; ++i){
897
898 Double_t pt = -1.;
899 if(fToyMinTrackPt!=fToyMaxTrackPt){
900 if(fToyDistributionTrackPt==0){
901 pt = rndm->Uniform(fToyMinTrackPt, fToyMaxTrackPt);
902 } else {
903 while(pt<fToyMinTrackPt||pt>fToyMaxTrackPt){
904 pt = rndm->Exp(fToyDistributionTrackPt); // no power law yet!!
905 pt += fToyMinTrackPt;
906 }
907 }
908 } else {
909 pt = fToyMinTrackPt;
910 }
911 Double_t eta = rndm->Uniform(fToyMinTrackEta,fToyMaxTrackEta);
912 Double_t phi = rndm->Uniform(fToyMinTrackPhi,fToyMaxTrackPhi);
913 phi = TVector2::Phi_0_2pi(phi);
914
915 if(fDebug) Printf("Add track #%d: pt %.2f, eta %.2f, phi %.2f", i, pt, eta, phi);
916
917 Double_t theta = 2*TMath::ATan(1/TMath::Exp(eta));
918 Float_t mom[3];
919 mom[0] = pt;
920 mom[1] = phi;
921 mom[2] = theta;
922
923 Float_t xyz[3];
924 xyz[0] = -999.;
925 xyz[1] = -999.;
926 xyz[2] = -999.;
927
928 AliAODTrack *tmpTr = new AliAODTrack( -999, // id
929 -999, // label
930 mom, // momentum[3]
931 kFALSE, // cartesian
932 xyz, // position
933 kFALSE, // DCA
934 NULL, // covMatrix,
935 -99, // charge
936 0, // itsClusMap
937 NULL, // pid
938 NULL, // prodVertex
939 kFALSE, // used for vertex fit
940 kFALSE, // used for prim vtx fit
941 AliAODTrack::kUndef, // type
942 fToyFilterMap, // select info
943 -999. // chi2 per NDF
944 );
945 tmpTr->SetFlags(AliESDtrack::kEmbedded);
946
947 new((*tracks)[nAODtracks++]) AliAODTrack(*tmpTr);
948 dummy = (*tracks)[nAODtracks-1];
949
950 fh1TrackPt->Fill(pt);
951 fh2TrackEtaPhi->Fill(eta,phi);
952
953 delete tmpTr;
954 }
955 } //end: embed mode with toy
956
957
958 PostData(1, fHistList);
b725dccf 959}
960
961
962//__________________________________________________________________________
963void AliAnalysisTaskFastEmbedding::Terminate(Option_t *)
964{
31b9d515 965 if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::Terminate()");
b725dccf 966
31b9d515 967 if(fAODfile && fAODfile->IsOpen())
968 fAODfile->Close();
b725dccf 969
970}
971
972//__________________________________________________________________________
b725dccf 973Int_t AliAnalysisTaskFastEmbedding::GetJobID()
974{
0ed51e0d 975 Int_t id=-1;
b725dccf 976
0ed51e0d 977 const char* env = gSystem->Getenv("ALIEN_PROC_ID"); // GRID
978 //if(!env || !strlen(env)) env = gSystem->Getenv("LSB_JOBINDEX"); // GSI
b725dccf 979
980 if(env && strlen(env)){
31b9d515 981 id= atoi(env);
982 AliInfo(Form("Job index %d", id));
b725dccf 983 }
984 else{
31b9d515 985 AliInfo("Job index not found. Okay if running locally.");
b725dccf 986 }
987
988 return id;
0ed51e0d 989}
b725dccf 990
991//__________________________________________________________________________
992
1a2bb6d5 993Int_t AliAnalysisTaskFastEmbedding::SelectAODfile(){
b725dccf 994
31b9d515 995 Int_t nFiles = fAODPathArray->GetEntries();
996
997 // choose randomly file and event
998 Int_t n = -1;
999 Float_t tmpProp = -1.;
1000 Int_t pendingEvents = fInputEntries-Entry();
1001 //Printf("input entries %d, entry %d, pending events %d", fInputEntries, (Int_t)Entry(), pendingEvents);
1002 while(rndm->Rndm()>tmpProp){
1003 n = rndm->Integer(nFiles);
1004 fAODEntries = fAODEntriesArray->At(n);
1005 Int_t tmpEntries = fAODEntries<pendingEvents ? pendingEvents : fAODEntries;
1006 tmpProp = fAODEntriesMax ? (Float_t)tmpEntries/fAODEntriesMax : 1.;
1007 }
1008 fAODEntry = (Int_t)(rndm->Uniform(fAODEntries));
1009
1010 AliInfo(Form("Select AOD file %d", n));
1011 TObjString *objStr = (TObjString*) fAODPathArray->At(n);
1012 if(!objStr){
1013 AliError("Could not get path of aod file.");
1014 return -1;
1015 }
1016 fAODPath = objStr->GetString();
1017
1018 return n;
1a2bb6d5 1019}
1020
46465e39 1021//__________________________________________________________________________
1022
1023Int_t AliAnalysisTaskFastEmbedding::OpenAODfile(Int_t trial){
1024
31b9d515 1025 if(fAODPathArray) fFileId = SelectAODfile();
1026 if(fFileId<0) return -1;
1027
1028 TDirectory *owd = gDirectory;
1029 if (fAODfile) fAODfile->Close();
1030 fAODfile = TFile::Open(fAODPath.Data(),"TIMEOUT=180");
1031 owd->cd();
1032 if(!fAODfile){
1033
1034 fFileId = -1;
1035 if(fAODPathArray){
1036 if(trial<=3){
1037 AliError(Form("Could not open AOD file %s (trial %d), try again ...", fAODPath.Data(), trial));
1038 fFileId = OpenAODfile(trial+1);
1039 } else {
1040 AliError(Form("Could not open AOD file %s, give up ...",fAODPath.Data()));
1041 return -1;
1042 }
1043 } else {
1044 AliError(Form("Could not open AOD file %s.", fAODPath.Data()));
1045 }
1046
1047 return fFileId;
1048 }
1049
1050 fAODtree = (TTree*)fAODfile->Get("aodTree");
1051
1052 if(!fAODtree){
1053 AliError("AOD tree not found.");
1054 return -1;
1055 }
1056
1057 /*
1058 fAODtree->SetBranchStatus("*",0);
1059 fAODtree->SetBranchStatus("header",1);
1060 fAODtree->SetBranchStatus("tracks*",1);
1061 fAODtree->SetBranchStatus("jets*",1);
1062 fAODtree->SetBranchStatus("clusters*",1);
1063 fAODtree->SetBranchStatus("mcparticles*",1);
1064 fAODtree->SetBranchStatus("mcHeader*",1);
1065 */
1066
1067 delete fAODevent;
1068 fAODevent = new AliAODEvent();
1069 fAODevent->ReadFromTree(fAODtree);
1070
1071
1072 // fetch header, jets, etc. from new file
1073 fNevents = fAODtree->GetEntries();
1074 mcHeader = (AliAODMCHeader*)fAODevent->FindListObject(AliAODMCHeader::StdBranchName());
1075
1076 if(fJetBranch.Length()) fAODJets = dynamic_cast<TClonesArray*>(fAODevent->FindListObject(fJetBranch.Data()));
1077 else fAODJets = fAODevent->GetJets();
1078 if(!fAODJets){
1079 AliError("Could not find jets in AOD. Check jet branch when indicated.");
1080 return -1;
1081 }
1082
1083 TFile *curfile = fAODtree->GetCurrentFile();
1084 if (!curfile) {
1085 AliError("No current file.");
1086 return -1;
1087 }
1088
1089 Float_t trials = 1.;
1090 Float_t xsec = 0.;
1091 PythiaInfoFromFile(curfile->GetName(),xsec,trials);
1092 fXsection = xsec;
1093
1094 // construct a poor man average trials
1095 Float_t nEntries = (Float_t)fAODtree->GetTree()->GetEntries();
1096 if(trials>=nEntries && nEntries>0.)fAvgTrials = trials/nEntries;
1097
1098 if(fFileId>=0){
1099 AliInfo(Form("Read successfully AOD event from file %d",fFileId));
1100 }
1101
1102 fCountEvents=0; // new file, reset counter
1103
1104 return fFileId; // file position in AOD path array, if array available
1105}
1106
1107
1108//____________________________________________________________________________
1109Float_t AliAnalysisTaskFastEmbedding::GetPtHard(Bool_t bSet, Float_t newValue){
1110
1111 static Float_t ptHard = -1.;
1112 if(bSet) ptHard = newValue;
1113
1114 return ptHard;
1a2bb6d5 1115}
b725dccf 1116
31b9d515 1117//____________________________________________________________________________
1118Int_t AliAnalysisTaskFastEmbedding::GetPtHardBin(Double_t ptHard){
1119
1120 const Int_t nBins = 10;
1121 Double_t binLimits[nBins] = { 5., 11., 21., 36., 57., 84., 117., 156., 200., 249. }; // lower limits
1122
1123 Int_t bin = -1;
1124 while(bin<nBins-1 && binLimits[bin+1]<ptHard){
1125 bin++;
1126 }
1127
1128 return bin;
1129}
1130
1131//____________________________________________________________________________
1132Bool_t AliAnalysisTaskFastEmbedding::PythiaInfoFromFile(const char* currFile,Float_t &fXsec,Float_t &fTrials){
1133 //
1134 // get the cross section and the trails either from pyxsec.root or from pysec_hists.root
1135 // This should provide the path to the AOD/ESD file
1136 // (taken from PWG4/JetTasks/AliAnalysisHelperJetTasks)
1137
1138 TString file(currFile);
1139 fXsec = 0;
1140 fTrials = 1;
1141
1142 if(file.Contains("root_archive.zip#")){
1143 Ssiz_t pos1 = file.Index("root_archive",12,TString::kExact);
1144 Ssiz_t pos = file.Index("#",1,pos1,TString::kExact);
1145 file.Replace(pos+1,20,"");
1146 }
1147 else {
1148 // not an archive take the basename....
1149 file.ReplaceAll(gSystem->BaseName(file.Data()),"");
1150 }
1151
1152 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
1153 if(!fxsec){
1154 // next trial fetch the histgram file
1155 fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec_hists.root"));
1156 if(!fxsec){
1157 // not a severe condition but inciate that we have no information
1158 AliDebug(AliLog::kDebug,Form("Neither pyxsec.root nor pyxsec_hists.root found in %s",file.Data()));
1159 return kFALSE;
1160 }
1161 else{
1162 // find the tlist we want to be independtent of the name so use the Tkey
1163 TKey* key = (TKey*)fxsec->GetListOfKeys()->At(0);
1164 if(!key){
1165 fxsec->Close();
1166 return kFALSE;
1167 }
1168 TList *list = dynamic_cast<TList*>(key->ReadObj());
1169 if(!list){
1170 fxsec->Close();
1171 return kFALSE;
1172 }
1173 fXsec = ((TProfile*)list->FindObject("h1Xsec"))->GetBinContent(1);
1174 fTrials = ((TH1F*)list->FindObject("h1Trials"))->GetBinContent(1);
1175 fxsec->Close();
1176 }
1177 } // no tree pyxsec.root
1178 else {
1179 TTree *xtree = (TTree*)fxsec->Get("Xsection");
1180 if(!xtree){
1181 fxsec->Close();
1182 return kFALSE;
1183 }
1184 UInt_t ntrials = 0;
1185 Double_t xsection = 0;
1186 xtree->SetBranchAddress("xsection",&xsection);
1187 xtree->SetBranchAddress("ntrials",&ntrials);
1188 xtree->GetEntry(0);
1189 fTrials = ntrials;
1190 fXsec = xsection;
1191 fxsec->Close();
1192 }
1193 return kTRUE;
1194}
1195
1196