]> 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
Fix for true centrality estimator with hijing npart (Alberica)
[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);
846
847 if(fEmbedMode==kAODJetTracks){
848 // jet track? else continue
849 // not implemented yet
850 continue;
851 }
852
853 new((*mcpartOUT)[nAODmcpart++]) AliAODMCParticle(*tmpPart);
854 dummy = (*mcpartOUT)[nAODmcpart-1];
855
856 if(tmpPart->IsPhysicalPrimary() && tmpPart->Charge()!=0. && tmpPart->Charge()!=-99. ){
857 fh1MCTrackPt->Fill(tmpPart->Pt());
858 fh2MCTrackEtaPhi->Fill(tmpPart->Eta(), tmpPart->Phi());
859 nPhysicalPrimary++;
860 }
861 }
862 fh1MCTrackN->Fill((Float_t)nPhysicalPrimary);
863
864 }
865 } // end: embed all tracks or jet tracks
866
867
868 if(fEmbedMode==kAODJet4Mom){
869
870 // loop over jets
871 Int_t nJets = fAODJets->GetEntries();
872 fh1TrackN->Fill((Float_t)nJets);
873 for(Int_t ij=0; ij<nJets; ++ij){
874 AliAODJet *jet = dynamic_cast<AliAODJet*>(fAODJets->At(ij));
875 if(!jet) continue;
876 AliAODTrack *tmpTr = (AliAODTrack*)jet;
877 tmpTr->SetFlags(AliESDtrack::kEmbedded);
878
879 new ((*tracks)[nAODtracks++]) AliAODTrack(*tmpTr);
880 dummy = (*tracks)[nAODtracks-1];
881
882 fh1TrackPt->Fill(tmpTr->Pt());
883 fh2TrackEtaPhi->Fill(tmpTr->Eta(), tmpTr->Phi());
884 }
885
886 } // end: embed jets as 4-momenta
887
888
889 } //end: embed mode with AOD
890
891
892 // === embed mode with toy ===
893 if(fEmbedMode==kToyTracks){
894 Int_t nT = (Int_t)(rndm->Uniform(fToyMinNbOfTracks, fToyMaxNbOfTracks)+0.5);
895
896 fh1TrackN->Fill((Float_t)nT);
897
898 for(Int_t i=0; i<nT; ++i){
899
900 Double_t pt = -1.;
901 if(fToyMinTrackPt!=fToyMaxTrackPt){
902 if(fToyDistributionTrackPt==0){
903 pt = rndm->Uniform(fToyMinTrackPt, fToyMaxTrackPt);
904 } else {
905 while(pt<fToyMinTrackPt||pt>fToyMaxTrackPt){
906 pt = rndm->Exp(fToyDistributionTrackPt); // no power law yet!!
907 pt += fToyMinTrackPt;
908 }
909 }
910 } else {
911 pt = fToyMinTrackPt;
912 }
913 Double_t eta = rndm->Uniform(fToyMinTrackEta,fToyMaxTrackEta);
914 Double_t phi = rndm->Uniform(fToyMinTrackPhi,fToyMaxTrackPhi);
915 phi = TVector2::Phi_0_2pi(phi);
916
917 if(fDebug) Printf("Add track #%d: pt %.2f, eta %.2f, phi %.2f", i, pt, eta, phi);
918
919 Double_t theta = 2*TMath::ATan(1/TMath::Exp(eta));
920 Float_t mom[3];
921 mom[0] = pt;
922 mom[1] = phi;
923 mom[2] = theta;
924
925 Float_t xyz[3];
926 xyz[0] = -999.;
927 xyz[1] = -999.;
928 xyz[2] = -999.;
929
930 AliAODTrack *tmpTr = new AliAODTrack( -999, // id
931 -999, // label
932 mom, // momentum[3]
933 kFALSE, // cartesian
934 xyz, // position
935 kFALSE, // DCA
936 NULL, // covMatrix,
937 -99, // charge
938 0, // itsClusMap
939 NULL, // pid
940 NULL, // prodVertex
941 kFALSE, // used for vertex fit
942 kFALSE, // used for prim vtx fit
943 AliAODTrack::kUndef, // type
944 fToyFilterMap, // select info
945 -999. // chi2 per NDF
946 );
947 tmpTr->SetFlags(AliESDtrack::kEmbedded);
948
949 new((*tracks)[nAODtracks++]) AliAODTrack(*tmpTr);
950 dummy = (*tracks)[nAODtracks-1];
951
952 fh1TrackPt->Fill(pt);
953 fh2TrackEtaPhi->Fill(eta,phi);
954
955 delete tmpTr;
956 }
957 } //end: embed mode with toy
958
959
960 PostData(1, fHistList);
b725dccf 961}
962
963
964//__________________________________________________________________________
965void AliAnalysisTaskFastEmbedding::Terminate(Option_t *)
966{
7a88ef30 967 // terminate
31b9d515 968 if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::Terminate()");
b725dccf 969
31b9d515 970 if(fAODfile && fAODfile->IsOpen())
971 fAODfile->Close();
b725dccf 972
973}
974
975//__________________________________________________________________________
b725dccf 976Int_t AliAnalysisTaskFastEmbedding::GetJobID()
977{
7a88ef30 978 // gives back the alien subjob id, if available, else -1
979
0ed51e0d 980 Int_t id=-1;
b725dccf 981
0ed51e0d 982 const char* env = gSystem->Getenv("ALIEN_PROC_ID"); // GRID
983 //if(!env || !strlen(env)) env = gSystem->Getenv("LSB_JOBINDEX"); // GSI
b725dccf 984
985 if(env && strlen(env)){
31b9d515 986 id= atoi(env);
987 AliInfo(Form("Job index %d", id));
b725dccf 988 }
989 else{
31b9d515 990 AliInfo("Job index not found. Okay if running locally.");
b725dccf 991 }
992
993 return id;
0ed51e0d 994}
b725dccf 995
996//__________________________________________________________________________
7a88ef30 997Int_t AliAnalysisTaskFastEmbedding::SelectAODfile()
998{
999 // select an AOD file from fAODPathArray
b725dccf 1000
31b9d515 1001 Int_t nFiles = fAODPathArray->GetEntries();
1002
1003 // choose randomly file and event
1004 Int_t n = -1;
1005 Float_t tmpProp = -1.;
1006 Int_t pendingEvents = fInputEntries-Entry();
1007 //Printf("input entries %d, entry %d, pending events %d", fInputEntries, (Int_t)Entry(), pendingEvents);
7a88ef30 1008 if(fAODEntriesArray){
b5fee0f9 1009 while(rndm->Rndm()>tmpProp){
1010 n = rndm->Integer(nFiles);
1011 fAODEntries = fAODEntriesArray->At(n);
1012 Int_t tmpEntries = fAODEntries<pendingEvents ? pendingEvents : fAODEntries;
1013 tmpProp = fAODEntriesMax ? (Float_t)tmpEntries/fAODEntriesMax : 1.;
1014 }
1015 fAODEntry = (Int_t)(rndm->Uniform(fAODEntries));
7a88ef30 1016 } else {
1017 AliWarning("Number of entries in extra AODs not set!");
1018 n = rndm->Integer(nFiles);
1019 fAODEntry = 0;
1020 }
31b9d515 1021
1022 AliInfo(Form("Select AOD file %d", n));
1023 TObjString *objStr = (TObjString*) fAODPathArray->At(n);
1024 if(!objStr){
1025 AliError("Could not get path of aod file.");
1026 return -1;
1027 }
1028 fAODPath = objStr->GetString();
1029
1030 return n;
1a2bb6d5 1031}
1032
46465e39 1033//__________________________________________________________________________
7a88ef30 1034Int_t AliAnalysisTaskFastEmbedding::OpenAODfile(Int_t trial)
1035{
1036 // select and open an AOD file from fAODPathArray
46465e39 1037
31b9d515 1038 if(fAODPathArray) fFileId = SelectAODfile();
1039 if(fFileId<0) return -1;
1040
b5fee0f9 1041 if (!gGrid) {
1042 AliInfo("Trying to connect to AliEn ...");
1043 TGrid::Connect("alien://");
1044 }
1045
31b9d515 1046 TDirectory *owd = gDirectory;
7a88ef30 1047 if (fAODfile && fAODfile->IsOpen()) fAODfile->Close();
31b9d515 1048 fAODfile = TFile::Open(fAODPath.Data(),"TIMEOUT=180");
1049 owd->cd();
1050 if(!fAODfile){
1051
1052 fFileId = -1;
1053 if(fAODPathArray){
1054 if(trial<=3){
1055 AliError(Form("Could not open AOD file %s (trial %d), try again ...", fAODPath.Data(), trial));
1056 fFileId = OpenAODfile(trial+1);
1057 } else {
1058 AliError(Form("Could not open AOD file %s, give up ...",fAODPath.Data()));
1059 return -1;
1060 }
1061 } else {
1062 AliError(Form("Could not open AOD file %s.", fAODPath.Data()));
1063 }
1064
1065 return fFileId;
1066 }
1067
1068 fAODtree = (TTree*)fAODfile->Get("aodTree");
1069
1070 if(!fAODtree){
1071 AliError("AOD tree not found.");
1072 return -1;
1073 }
1074
1075 /*
1076 fAODtree->SetBranchStatus("*",0);
1077 fAODtree->SetBranchStatus("header",1);
1078 fAODtree->SetBranchStatus("tracks*",1);
1079 fAODtree->SetBranchStatus("jets*",1);
1080 fAODtree->SetBranchStatus("clusters*",1);
1081 fAODtree->SetBranchStatus("mcparticles*",1);
1082 fAODtree->SetBranchStatus("mcHeader*",1);
1083 */
1084
1085 delete fAODevent;
1086 fAODevent = new AliAODEvent();
1087 fAODevent->ReadFromTree(fAODtree);
1088
1089
1090 // fetch header, jets, etc. from new file
1091 fNevents = fAODtree->GetEntries();
1092 mcHeader = (AliAODMCHeader*)fAODevent->FindListObject(AliAODMCHeader::StdBranchName());
1093
1094 if(fJetBranch.Length()) fAODJets = dynamic_cast<TClonesArray*>(fAODevent->FindListObject(fJetBranch.Data()));
1095 else fAODJets = fAODevent->GetJets();
1096 if(!fAODJets){
1097 AliError("Could not find jets in AOD. Check jet branch when indicated.");
1098 return -1;
1099 }
1100
1101 TFile *curfile = fAODtree->GetCurrentFile();
1102 if (!curfile) {
1103 AliError("No current file.");
1104 return -1;
1105 }
1106
1107 Float_t trials = 1.;
1108 Float_t xsec = 0.;
1109 PythiaInfoFromFile(curfile->GetName(),xsec,trials);
1110 fXsection = xsec;
1111
1112 // construct a poor man average trials
1113 Float_t nEntries = (Float_t)fAODtree->GetTree()->GetEntries();
1114 if(trials>=nEntries && nEntries>0.)fAvgTrials = trials/nEntries;
1115
1116 if(fFileId>=0){
1117 AliInfo(Form("Read successfully AOD event from file %d",fFileId));
1118 }
1119
1120 fCountEvents=0; // new file, reset counter
1121
1122 return fFileId; // file position in AOD path array, if array available
1123}
1124
1125
1126//____________________________________________________________________________
1127Float_t AliAnalysisTaskFastEmbedding::GetPtHard(Bool_t bSet, Float_t newValue){
1128
7a88ef30 1129 // static stored, available for other tasks in train
1130
31b9d515 1131 static Float_t ptHard = -1.;
1132 if(bSet) ptHard = newValue;
1133
1134 return ptHard;
1a2bb6d5 1135}
b725dccf 1136
31b9d515 1137//____________________________________________________________________________
1138Int_t AliAnalysisTaskFastEmbedding::GetPtHardBin(Double_t ptHard){
1139
7a88ef30 1140 // returns pt hard bin (for LHC10e14, LHC11a1x, LHC11a2x simulations)
1141
31b9d515 1142 const Int_t nBins = 10;
1143 Double_t binLimits[nBins] = { 5., 11., 21., 36., 57., 84., 117., 156., 200., 249. }; // lower limits
1144
1145 Int_t bin = -1;
1146 while(bin<nBins-1 && binLimits[bin+1]<ptHard){
1147 bin++;
1148 }
1149
1150 return bin;
1151}
1152
1153//____________________________________________________________________________
1154Bool_t AliAnalysisTaskFastEmbedding::PythiaInfoFromFile(const char* currFile,Float_t &fXsec,Float_t &fTrials){
1155 //
1156 // get the cross section and the trails either from pyxsec.root or from pysec_hists.root
1157 // This should provide the path to the AOD/ESD file
1158 // (taken from PWG4/JetTasks/AliAnalysisHelperJetTasks)
1159
1160 TString file(currFile);
1161 fXsec = 0;
1162 fTrials = 1;
1163
1164 if(file.Contains("root_archive.zip#")){
b5fee0f9 1165 Ssiz_t pos1 = file.Index("root_archive",12,0,TString::kExact);
1166 Ssiz_t pos = file.Index("#",1,pos1,TString::kExact);
1167 Ssiz_t pos2 = file.Index(".root",5,TString::kExact);
1168 file.Replace(pos+1,pos2-pos1,"");
1169 // file.Replace(pos+1,20,"");
31b9d515 1170 }
1171 else {
1172 // not an archive take the basename....
1173 file.ReplaceAll(gSystem->BaseName(file.Data()),"");
1174 }
1175
1176 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
1177 if(!fxsec){
1178 // next trial fetch the histgram file
1179 fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec_hists.root"));
1180 if(!fxsec){
1181 // not a severe condition but inciate that we have no information
1182 AliDebug(AliLog::kDebug,Form("Neither pyxsec.root nor pyxsec_hists.root found in %s",file.Data()));
1183 return kFALSE;
1184 }
1185 else{
1186 // find the tlist we want to be independtent of the name so use the Tkey
1187 TKey* key = (TKey*)fxsec->GetListOfKeys()->At(0);
1188 if(!key){
1189 fxsec->Close();
1190 return kFALSE;
1191 }
1192 TList *list = dynamic_cast<TList*>(key->ReadObj());
1193 if(!list){
1194 fxsec->Close();
1195 return kFALSE;
1196 }
1197 fXsec = ((TProfile*)list->FindObject("h1Xsec"))->GetBinContent(1);
1198 fTrials = ((TH1F*)list->FindObject("h1Trials"))->GetBinContent(1);
1199 fxsec->Close();
1200 }
1201 } // no tree pyxsec.root
1202 else {
1203 TTree *xtree = (TTree*)fxsec->Get("Xsection");
1204 if(!xtree){
1205 fxsec->Close();
1206 return kFALSE;
1207 }
1208 UInt_t ntrials = 0;
1209 Double_t xsection = 0;
1210 xtree->SetBranchAddress("xsection",&xsection);
1211 xtree->SetBranchAddress("ntrials",&ntrials);
1212 xtree->GetEntry(0);
1213 fTrials = ntrials;
1214 fXsec = xsection;
1215 fxsec->Close();
1216 }
1217 return kTRUE;
1218}
1219
1220