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