#include <TRandom3.h>
#include <TH1F.h>
#include <TH2F.h>
+#include <TProfile.h>
+#include <TKey.h>
+#include <TGrid.h>
#include "AliAnalysisTaskFastEmbedding.h"
#include "AliAnalysisManager.h"
+#include "AliESDEvent.h"
+#include "AliESDtrack.h"
#include "AliAODEvent.h"
#include "AliAODTrack.h"
#include "AliAODJet.h"
#include "AliAODMCParticle.h"
+#include "AliAODMCHeader.h"
+#include "AliInputEventHandler.h"
#include "AliLog.h"
//__________________________________________________________________________
AliAnalysisTaskFastEmbedding::AliAnalysisTaskFastEmbedding()
- : AliAnalysisTaskSE()
- ,fAODout(0)
- ,fAODevent(0)
- ,fAODtree(0)
- ,fAODfile(0)
- ,rndm(0)
- ,fAODPathArray(0)
- ,fAODPath("AliAOD.root")
- ,fTrackBranch("aodExtraTracks")
- ,fMCparticlesBranch("aodExtraMCparticles")
- ,fJetBranch("")
- ,fEntry(0)
- ,fEmbedMode(0)
- ,fEvtSelecMode(0)
- ,fEvtSelMinJetPt(-1)
- ,fEvtSelMaxJetPt(-1)
- ,fToyMinNbOfTracks(1)
- ,fToyMaxNbOfTracks(1)
- ,fToyMinTrackPt(50.)
- ,fToyMaxTrackPt(50.)
- ,fToyDistributionTrackPt(0.)
- ,fToyMinTrackEta(-.5)
- ,fToyMaxTrackEta(.5)
- ,fToyMinTrackPhi(0.)
- ,fToyMaxTrackPhi(2*TMath::Pi())
- ,fToyFilterMap(0)
- ,fHistList(0)
- ,fh1TrackPt(0)
- ,fh2TrackEtaPhi(0)
- ,fh1TrackN(0)
- ,fh1MCTrackPt(0)
- ,fh2MCTrackEtaPhi(0)
- ,fh1MCTrackN(0)
-
-
+: AliAnalysisTaskSE()
+,fESD(0)
+,fAODout(0)
+,fAODevent(0)
+,fAODtree(0)
+,fAODfile(0)
+,mcHeader(0)
+,rndm(0)
+,fInputEntries(0)
+,fAODPathArray(0)
+,fAODEntriesArray(0)
+,fAODPath("AliAOD.root")
+,fAODEntries(-1)
+,fAODEntriesSum(0)
+,fAODEntriesMax(0)
+,fOfflineTrgMask(AliVEvent::kAny)
+,fMinContribVtx(1)
+,fVtxZMin(-8.)
+,fVtxZMax(8.)
+,fEvtClassMin(0)
+,fEvtClassMax(4)
+,fCentMin(0.)
+,fCentMax(100.)
+,fNInputTracksMin(0)
+,fNInputTracksMax(-1)
+,fTrackBranch("aodExtraTracks")
+,fMCparticlesBranch("aodExtraMCparticles")
+,fJetBranch("")
+,fFileId(-1)
+,fAODEntry(0)
+,fCountEvents(-1)
+,fEmbedMode(0)
+,fEvtSelecMode(0)
+,fEvtSelMinJetPt(-1)
+,fEvtSelMaxJetPt(-1)
+,fEvtSelMinJetEta(-999.)
+,fEvtSelMaxJetEta( 999.)
+,fEvtSelMinJetPhi(0.)
+,fEvtSelMaxJetPhi(TMath::Pi()*2.)
+,fExtraEffPb(1)
+,fToyMinNbOfTracks(1)
+,fToyMaxNbOfTracks(1)
+,fToyMinTrackPt(50.)
+,fToyMaxTrackPt(50.)
+,fToyDistributionTrackPt(0.)
+,fToyMinTrackEta(-.5)
+,fToyMaxTrackEta(.5)
+,fToyMinTrackPhi(0.)
+,fToyMaxTrackPhi(2*TMath::Pi())
+,fToyFilterMap(0)
+,fTrackFilterMap(0)
+,fNPtHard(10)
+,fPtHard(0)
+,fPtHardBin(-1)
+,fAODJets(0x0)
+,fNevents(0)
+,fXsection(0)
+,fAvgTrials(0)
+,fHistList(0)
+,fHistEvtSelection(0)
+,fh1Xsec(0)
+,fh1Trials(0)
+,fh1TrialsEvtSel(0)
+,fh2PtHard(0)
+,fh2PtHardEvtSel(0)
+,fh2PtHardTrials(0)
+,fh1TrackPt(0)
+,fh2TrackEtaPhi(0)
+,fh1TrackN(0)
+,fh1JetPt(0)
+,fh2JetEtaPhi(0)
+,fh1JetN(0)
+,fh1MCTrackPt(0)
+,fh2MCTrackEtaPhi(0)
+,fh1MCTrackN(0)
+,fh1AODfile(0)
+,fh2AODevent(0)
{
- // default constructor
+ // default constructor
}
//__________________________________________________________________________
AliAnalysisTaskFastEmbedding::AliAnalysisTaskFastEmbedding(const char *name)
- : AliAnalysisTaskSE(name)
- ,fAODout(0)
- ,fAODevent(0)
- ,fAODtree(0)
- ,fAODfile(0)
- ,rndm(0)
- ,fAODPathArray(0)
- ,fAODPath("AliAOD.root")
- ,fTrackBranch("aodExtraTracks")
- ,fMCparticlesBranch("aodExtraMCparticles")
- ,fJetBranch("")
- ,fEntry(0)
- ,fEmbedMode(0)
- ,fEvtSelecMode(0)
- ,fEvtSelMinJetPt(-1)
- ,fEvtSelMaxJetPt(-1)
- ,fToyMinNbOfTracks(1)
- ,fToyMaxNbOfTracks(1)
- ,fToyMinTrackPt(50.)
- ,fToyMaxTrackPt(50.)
- ,fToyDistributionTrackPt(0.)
- ,fToyMinTrackEta(-.5)
- ,fToyMaxTrackEta(.5)
- ,fToyMinTrackPhi(0.)
- ,fToyMaxTrackPhi(2*TMath::Pi())
- ,fToyFilterMap(0)
- ,fHistList(0)
- ,fh1TrackPt(0)
- ,fh2TrackEtaPhi(0)
- ,fh1TrackN(0)
- ,fh1MCTrackPt(0)
- ,fh2MCTrackEtaPhi(0)
- ,fh1MCTrackN(0)
+: AliAnalysisTaskSE(name)
+,fESD(0)
+,fAODout(0)
+,fAODevent(0)
+,fAODtree(0)
+,fAODfile(0)
+,mcHeader(0)
+,rndm(0)
+,fInputEntries(0)
+,fAODPathArray(0)
+,fAODEntriesArray(0)
+,fAODPath("AliAOD.root")
+,fAODEntries(-1)
+,fAODEntriesSum(0)
+,fAODEntriesMax(0)
+,fOfflineTrgMask(AliVEvent::kAny)
+,fMinContribVtx(1)
+,fVtxZMin(-8.)
+,fVtxZMax(8.)
+,fEvtClassMin(0)
+,fEvtClassMax(4)
+,fCentMin(0.)
+,fCentMax(100.)
+,fNInputTracksMin(0)
+,fNInputTracksMax(-1)
+,fTrackBranch("aodExtraTracks")
+,fMCparticlesBranch("aodExtraMCparticles")
+,fJetBranch("")
+,fFileId(-1)
+,fAODEntry(0)
+,fCountEvents(-1)
+,fEmbedMode(0)
+,fEvtSelecMode(0)
+,fEvtSelMinJetPt(-1)
+,fEvtSelMaxJetPt(-1)
+,fEvtSelMinJetEta(-999.)
+,fEvtSelMaxJetEta( 999.)
+,fEvtSelMinJetPhi(0.)
+,fEvtSelMaxJetPhi(TMath::Pi()*2.)
+,fExtraEffPb(1)
+,fToyMinNbOfTracks(1)
+,fToyMaxNbOfTracks(1)
+,fToyMinTrackPt(50.)
+,fToyMaxTrackPt(50.)
+,fToyDistributionTrackPt(0.)
+,fToyMinTrackEta(-.5)
+,fToyMaxTrackEta(.5)
+,fToyMinTrackPhi(0.)
+,fToyMaxTrackPhi(2*TMath::Pi())
+,fToyFilterMap(0)
+,fTrackFilterMap(0)
+,fNPtHard(10)
+,fPtHard(0)
+,fPtHardBin(-1)
+,fAODJets(0x0)
+,fNevents(0)
+,fXsection(0)
+,fAvgTrials(0)
+,fHistList(0)
+,fHistEvtSelection(0)
+,fh1Xsec(0)
+,fh1Trials(0)
+,fh1TrialsEvtSel(0)
+,fh2PtHard(0)
+,fh2PtHardEvtSel(0)
+,fh2PtHardTrials(0)
+,fh1TrackPt(0)
+,fh2TrackEtaPhi(0)
+,fh1TrackN(0)
+,fh1JetPt(0)
+,fh2JetEtaPhi(0)
+,fh1JetN(0)
+,fh1MCTrackPt(0)
+,fh2MCTrackEtaPhi(0)
+,fh1MCTrackN(0)
+,fh1AODfile(0)
+,fh2AODevent(0)
{
- // constructor
- DefineOutput(1, TList::Class());
+ // constructor
+ DefineOutput(1, TList::Class());
}
//__________________________________________________________________________
AliAnalysisTaskFastEmbedding::AliAnalysisTaskFastEmbedding(const AliAnalysisTaskFastEmbedding ©)
- : AliAnalysisTaskSE()
- ,fAODout(copy.fAODout)
- ,fAODevent(copy.fAODevent)
- ,fAODtree(copy.fAODtree)
- ,fAODfile(copy.fAODfile)
- ,rndm(copy.rndm)
- ,fAODPathArray(copy.fAODPathArray)
- ,fAODPath(copy.fAODPath)
- ,fTrackBranch(copy.fTrackBranch)
- ,fMCparticlesBranch(copy.fMCparticlesBranch)
- ,fJetBranch(copy.fJetBranch)
- ,fEntry(copy.fEntry)
- ,fEmbedMode(copy.fEmbedMode)
- ,fEvtSelecMode(copy.fEvtSelecMode)
- ,fEvtSelMinJetPt(copy.fEvtSelMinJetPt)
- ,fEvtSelMaxJetPt(copy.fEvtSelMaxJetPt)
- ,fToyMinNbOfTracks(copy.fToyMinNbOfTracks)
- ,fToyMaxNbOfTracks(copy.fToyMaxNbOfTracks)
- ,fToyMinTrackPt(copy.fToyMinTrackPt)
- ,fToyMaxTrackPt(copy.fToyMaxTrackPt)
- ,fToyDistributionTrackPt(copy.fToyDistributionTrackPt)
- ,fToyMinTrackEta(copy.fToyMinTrackEta)
- ,fToyMaxTrackEta(copy.fToyMaxTrackEta)
- ,fToyMinTrackPhi(copy.fToyMinTrackPhi)
- ,fToyMaxTrackPhi(copy.fToyMaxTrackPhi)
- ,fToyFilterMap(copy.fToyFilterMap)
- ,fHistList(copy.fHistList)
- ,fh1TrackPt(copy.fh1TrackPt)
- ,fh2TrackEtaPhi(copy.fh2TrackEtaPhi)
- ,fh1TrackN(copy.fh1TrackN)
- ,fh1MCTrackPt(copy.fh1MCTrackPt)
- ,fh2MCTrackEtaPhi(copy.fh2MCTrackEtaPhi)
- ,fh1MCTrackN(copy.fh1MCTrackN)
+: AliAnalysisTaskSE()
+,fESD(copy.fESD)
+,fAODout(copy.fAODout)
+,fAODevent(copy.fAODevent)
+,fAODtree(copy.fAODtree)
+,fAODfile(copy.fAODfile)
+,mcHeader(copy.mcHeader)
+,rndm(copy.rndm)
+,fInputEntries(copy.fInputEntries)
+,fAODPathArray(copy.fAODPathArray)
+,fAODEntriesArray(copy.fAODEntriesArray)
+,fAODPath(copy.fAODPath)
+,fAODEntries(copy.fAODEntries)
+,fAODEntriesSum(copy.fAODEntriesSum)
+,fAODEntriesMax(copy.fAODEntriesMax)
+,fOfflineTrgMask(copy.fOfflineTrgMask)
+,fMinContribVtx(copy.fMinContribVtx)
+,fVtxZMin(copy.fVtxZMin)
+,fVtxZMax(copy.fVtxZMax)
+,fEvtClassMin(copy.fEvtClassMin)
+,fEvtClassMax(copy.fEvtClassMax)
+,fCentMin(copy.fCentMin)
+,fCentMax(copy.fCentMax)
+,fNInputTracksMin(copy.fNInputTracksMin)
+,fNInputTracksMax(copy.fNInputTracksMax)
+,fTrackBranch(copy.fTrackBranch)
+,fMCparticlesBranch(copy.fMCparticlesBranch)
+,fJetBranch(copy.fJetBranch)
+,fFileId(copy.fFileId)
+,fAODEntry(copy.fAODEntry)
+,fCountEvents(copy.fCountEvents)
+,fEmbedMode(copy.fEmbedMode)
+,fEvtSelecMode(copy.fEvtSelecMode)
+,fEvtSelMinJetPt(copy.fEvtSelMinJetPt)
+,fEvtSelMaxJetPt(copy.fEvtSelMaxJetPt)
+,fEvtSelMinJetEta(copy.fEvtSelMinJetEta)
+,fEvtSelMaxJetEta(copy.fEvtSelMaxJetEta)
+,fEvtSelMinJetPhi(copy.fEvtSelMinJetPhi)
+,fEvtSelMaxJetPhi(copy.fEvtSelMaxJetPhi)
+,fExtraEffPb(copy.fExtraEffPb)
+,fToyMinNbOfTracks(copy.fToyMinNbOfTracks)
+,fToyMaxNbOfTracks(copy.fToyMaxNbOfTracks)
+,fToyMinTrackPt(copy.fToyMinTrackPt)
+,fToyMaxTrackPt(copy.fToyMaxTrackPt)
+,fToyDistributionTrackPt(copy.fToyDistributionTrackPt)
+,fToyMinTrackEta(copy.fToyMinTrackEta)
+,fToyMaxTrackEta(copy.fToyMaxTrackEta)
+,fToyMinTrackPhi(copy.fToyMinTrackPhi)
+,fToyMaxTrackPhi(copy.fToyMaxTrackPhi)
+,fToyFilterMap(copy.fToyFilterMap)
+,fTrackFilterMap(copy.fTrackFilterMap)
+,fNPtHard(copy.fNPtHard)
+,fPtHard(copy.fPtHard)
+,fPtHardBin(copy.fPtHardBin)
+,fAODJets(copy.fAODJets)
+,fNevents(copy.fNevents)
+,fXsection(copy.fXsection)
+,fAvgTrials(copy.fAvgTrials)
+,fHistList(copy.fHistList)
+,fHistEvtSelection(copy.fHistEvtSelection)
+,fh1Xsec(copy.fh1Xsec)
+,fh1Trials(copy.fh1Trials)
+,fh1TrialsEvtSel(copy.fh1TrialsEvtSel)
+,fh2PtHard(copy.fh2PtHard)
+,fh2PtHardEvtSel(copy.fh2PtHardEvtSel)
+,fh2PtHardTrials(copy.fh2PtHardTrials)
+,fh1TrackPt(copy.fh1TrackPt)
+,fh2TrackEtaPhi(copy.fh2TrackEtaPhi)
+,fh1TrackN(copy.fh1TrackN)
+,fh1JetPt(copy.fh1JetPt)
+,fh2JetEtaPhi(copy.fh2JetEtaPhi)
+,fh1JetN(copy.fh1JetN)
+,fh1MCTrackPt(copy.fh1MCTrackPt)
+,fh2MCTrackEtaPhi(copy.fh2MCTrackEtaPhi)
+,fh1MCTrackN(copy.fh1MCTrackN)
+,fh1AODfile(copy.fh1AODfile)
+,fh2AODevent(copy.fh2AODevent)
{
- // copy constructor
+ // copy constructor
}
//__________________________________________________________________________
AliAnalysisTaskFastEmbedding& AliAnalysisTaskFastEmbedding::operator=(const AliAnalysisTaskFastEmbedding& o)
{
- // assignment
-
- if(this!=&o){
- AliAnalysisTaskSE::operator=(o);
- fAODout = o.fAODout;
- fAODevent = o.fAODevent;
- fAODtree = o.fAODtree;
- fAODfile = o.fAODfile;
- rndm = o.rndm;
- fAODPathArray = o.fAODPathArray;
- fAODPath = o.fAODPath;
- fTrackBranch = o.fTrackBranch;
- fMCparticlesBranch = o.fMCparticlesBranch;
- fJetBranch = o.fJetBranch;
- fEntry = o.fEntry;
- fEmbedMode = o.fEmbedMode;
- fEvtSelecMode = o.fEvtSelecMode;
- fEvtSelMinJetPt = o.fEvtSelMinJetPt;
- fEvtSelMaxJetPt = o.fEvtSelMaxJetPt;
- fToyMinNbOfTracks = o.fToyMinNbOfTracks;
- fToyMaxNbOfTracks = o.fToyMaxNbOfTracks;
- fToyMinTrackPt = o.fToyMinTrackPt;
- fToyMaxTrackPt = o.fToyMaxTrackPt;
- fToyDistributionTrackPt = o.fToyDistributionTrackPt;
- fToyMinTrackEta = o.fToyMinTrackEta;
- fToyMaxTrackEta = o.fToyMaxTrackEta;
- fToyMinTrackPhi = o.fToyMinTrackPhi;
- fToyMaxTrackPhi = o.fToyMaxTrackPhi;
- fToyFilterMap = o.fToyFilterMap;
- fHistList = o.fHistList;
- fh1TrackPt = o.fh1TrackPt;
- fh2TrackEtaPhi = o.fh2TrackEtaPhi;
- fh1TrackN = o.fh1TrackN;
- fh1MCTrackPt = o.fh1MCTrackPt;
- fh2MCTrackEtaPhi = o.fh2MCTrackEtaPhi;
- fh1MCTrackN = o.fh1MCTrackN;
- }
-
- return *this;
+ // assignment
+
+ if(this!=&o){
+ AliAnalysisTaskSE::operator=(o);
+ fESD = o.fESD;
+ fAODout = o.fAODout;
+ fAODevent = o.fAODevent;
+ fAODtree = o.fAODtree;
+ fAODfile = o.fAODfile;
+ mcHeader = o.mcHeader;
+ rndm = o.rndm;
+ fInputEntries = o.fInputEntries;
+ fAODPathArray = o.fAODPathArray;
+ fAODEntriesArray = o.fAODEntriesArray;
+ fAODPath = o.fAODPath;
+ fAODEntries = o.fAODEntries;
+ fAODEntriesSum = o.fAODEntriesSum;
+ fAODEntriesMax = o.fAODEntriesMax;
+ fOfflineTrgMask = o.fOfflineTrgMask;
+ fMinContribVtx = o.fMinContribVtx;
+ fVtxZMin = o.fVtxZMin;
+ fVtxZMax = o.fVtxZMax;
+ fEvtClassMin = o.fEvtClassMin;
+ fEvtClassMax = o.fEvtClassMax;
+ fCentMin = o.fCentMin;
+ fCentMax = o.fCentMax;
+ fNInputTracksMin = o.fNInputTracksMin;
+ fNInputTracksMax = o.fNInputTracksMax;
+ fTrackBranch = o.fTrackBranch;
+ fMCparticlesBranch = o.fMCparticlesBranch;
+ fJetBranch = o.fJetBranch;
+ fFileId = o.fFileId;
+ fAODEntry = o.fAODEntry;
+ fCountEvents = o.fCountEvents;
+ fEmbedMode = o.fEmbedMode;
+ fEvtSelecMode = o.fEvtSelecMode;
+ fEvtSelMinJetPt = o.fEvtSelMinJetPt;
+ fEvtSelMaxJetPt = o.fEvtSelMaxJetPt;
+ fEvtSelMinJetEta = o.fEvtSelMinJetEta;
+ fEvtSelMaxJetEta = o.fEvtSelMaxJetEta;
+ fEvtSelMinJetPhi = o.fEvtSelMinJetPhi;
+ fEvtSelMaxJetPhi = o.fEvtSelMaxJetPhi;
+ fExtraEffPb = o.fExtraEffPb;
+ fToyMinNbOfTracks = o.fToyMinNbOfTracks;
+ fToyMaxNbOfTracks = o.fToyMaxNbOfTracks;
+ fToyMinTrackPt = o.fToyMinTrackPt;
+ fToyMaxTrackPt = o.fToyMaxTrackPt;
+ fToyDistributionTrackPt = o.fToyDistributionTrackPt;
+ fToyMinTrackEta = o.fToyMinTrackEta;
+ fToyMaxTrackEta = o.fToyMaxTrackEta;
+ fToyMinTrackPhi = o.fToyMinTrackPhi;
+ fToyMaxTrackPhi = o.fToyMaxTrackPhi;
+ fToyFilterMap = o.fToyFilterMap;
+ fTrackFilterMap = o.fTrackFilterMap;
+ fNPtHard = o.fNPtHard;
+ fPtHard = o.fPtHard;
+ fPtHardBin = o.fPtHardBin;
+ fAODJets = o.fAODJets;
+ fNevents = o.fNevents;
+ fXsection = o.fXsection;
+ fAvgTrials = o.fAvgTrials;
+ fHistList = o.fHistList;
+ fHistEvtSelection = o.fHistEvtSelection;
+ fh1Xsec = o.fh1Xsec;
+ fh1Trials = o.fh1Trials;
+ fh1TrialsEvtSel = o.fh1TrialsEvtSel;
+ fh2PtHard = o.fh2PtHard;
+ fh2PtHardEvtSel = o.fh2PtHardEvtSel;
+ fh2PtHardTrials = o.fh2PtHardTrials;
+ fh1TrackPt = o.fh1TrackPt;
+ fh2TrackEtaPhi = o.fh2TrackEtaPhi;
+ fh1TrackN = o.fh1TrackN;
+ fh1JetPt = o.fh1JetPt;
+ fh2JetEtaPhi = o.fh2JetEtaPhi;
+ fh1JetN = o.fh1JetN;
+ fh1MCTrackPt = o.fh1MCTrackPt;
+ fh2MCTrackEtaPhi = o.fh2MCTrackEtaPhi;
+ fh1MCTrackN = o.fh1MCTrackN;
+ fh1AODfile = o.fh1AODfile;
+ fh2AODevent = o.fh2AODevent;
+ }
+
+ return *this;
}
//__________________________________________________________________________
AliAnalysisTaskFastEmbedding::~AliAnalysisTaskFastEmbedding()
{
- // destructor
- delete rndm;
+ // destructor
+ delete rndm;
}
//__________________________________________________________________________
void AliAnalysisTaskFastEmbedding::UserCreateOutputObjects()
{
- // create output objects
- if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::UserCreateOutputObjects()");
+ // create output objects
+ if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::UserCreateOutputObjects()");
+ AliLog::SetClassDebugLevel("AliAnalysisTaskFastEmbedding", AliLog::kInfo);
+
+ OpenFile(1);
+ if(!fHistList) fHistList = new TList();
+ fHistList->SetOwner(kTRUE);
+
+
+ // set seed
+ rndm = new TRandom3();
+ Int_t id = GetJobID();
+ if(id>-1) rndm->SetSeed(id);
+ else rndm->SetSeed(); // a TTUID is generated and used for seed
+ AliInfo(Form("TRandom3 seed: %d", rndm->GetSeed()));
+
+
+
+ // embed mode with AOD
+ if(fEmbedMode==kAODFull || fEmbedMode==kAODJetTracks || fEmbedMode==kAODJet4Mom){
+
+ // open input AOD
+ fFileId = OpenAODfile();
+ fNevents = 0; // force to open another aod in UserExec()
+ if(fFileId<0){
+ AliError("");
+ PostData(1, fHistList);
+ return;
+ }
+ } //end: embed mode with AOD
+
+
+ // connect output aod
+ // create a new branch for extra tracks
+ fAODout = AODEvent();
+ if(!fAODout){
+ AliError("Output AOD not found.");
+ PostData(1, fHistList);
+ return;
+ }
+ if(!fAODout->FindListObject(fTrackBranch.Data()) && strlen(fTrackBranch.Data())){
+ AliInfo(Form("Add AOD branch %s", fTrackBranch.Data()));
+ TClonesArray *tracks = new TClonesArray("AliAODTrack",0);
+ tracks->SetName(fTrackBranch.Data());
+ AddAODBranch("TClonesArray", &tracks);
+ }
+ // create new branch for extra mcparticle if available as input
+ if(fAODevent && fAODevent->FindListObject("mcparticles") && strlen(fMCparticlesBranch.Data())){
+ AliInfo(Form("Add AOD branch %s", fMCparticlesBranch.Data()));
+ TClonesArray *mcparticles = new TClonesArray("AliAODMCParticle",0);
+ mcparticles->SetName(fMCparticlesBranch.Data());
+ AddAODBranch("TClonesArray", &mcparticles);
+ }
- rndm = new TRandom3();
- Int_t id = GetJobID();
- if(id>-1) rndm->SetSeed(id);
- else rndm->SetSeed(); // a TTUID is generated and used for seed
- AliInfo(Form("TRandom3 seed: %d", rndm->GetSeed()));
- // embed mode with AOD
- if(fEmbedMode==kAODFull || fEmbedMode==kAODJetTracks || fEmbedMode==kAODJet4Mom){
+ //qa histograms
+ Bool_t oldStatus = TH1::AddDirectoryStatus();
+ TH1::AddDirectory(kFALSE);
+
+ fHistEvtSelection = new TH1I("fHistEvtSelection", "event selection", 6, -0.5, 5.5);
+ fHistEvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
+ fHistEvtSelection->GetXaxis()->SetBinLabel(2,"events IN");
+ fHistEvtSelection->GetXaxis()->SetBinLabel(3,"event selection (rejected)");
+ fHistEvtSelection->GetXaxis()->SetBinLabel(4,"vertex cut (rejected)");
+ fHistEvtSelection->GetXaxis()->SetBinLabel(5,"centrality (rejected)");
+ fHistEvtSelection->GetXaxis()->SetBinLabel(6,"multiplicity (rejected)");
+
+ fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root;p_{T,hard} bin;<#sigma>",fNPtHard+1,-1.5,fNPtHard-0.5);
+ fh1Trials = new TH1F("fh1Trials","trials (simulation) from pyxsec.root;p_{T,hard} bin;#sum{ntrials}",fNPtHard+1,-1.5,fNPtHard-0.5);
+ fh1TrialsEvtSel = new TH1F("fh1TrialsEvtSel","trials (event selection) from pyxsec.root;p_{T,hard} bin;#sum{ntrials}",fNPtHard+1,-1.5,fNPtHard-0.5);
+ 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);
+ 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);
+ 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);
+
+ fHistList->Add(fHistEvtSelection);
+ fHistList->Add(fh1Xsec);
+ fHistList->Add(fh1Trials);
+ fHistList->Add(fh1TrialsEvtSel);
+ fHistList->Add(fh2PtHard);
+ fHistList->Add(fh2PtHardEvtSel);
+ fHistList->Add(fh2PtHardTrials);
+
+ fh1TrackPt = new TH1F("fh1TrackPt","pT of extra tracks;p_{T};entries", 250, 0., 250.);
+ fh2TrackEtaPhi = new TH2F("fh2TrackEtaPhi","eta-phi distribution of extra tracks;#eta;#phi", 20, -1., 1., 60, 0., 2*TMath::Pi());
+ fh1TrackN = new TH1F("fh1TrackN", "nb. of extra tracks per event;nb. of tracks;entries",300, 0., 300.);
+
+ fHistList->Add(fh1TrackPt);
+ fHistList->Add(fh2TrackEtaPhi);
+ fHistList->Add(fh1TrackN);
+
+ if(fEmbedMode==kAODFull || fEmbedMode==kAODJetTracks || fEmbedMode==kAODJet4Mom){
+
+ fh1JetPt = new TH1F("fh1JetPt", "pT of extra jets;p_{T};entries", 120, 0., 120.);
+ fh2JetEtaPhi = new TH2F("fh2JetEtaPhi", "eta-phi distribution of extra jets;#eta;#phi", 20, -1., 1., 60, 0., 2*TMath::Pi());
+ fh1JetN = new TH1F("fh1JetN", "nb. of extra jets per event;nb. of jets;entries",20,0.,20.);
+
+ fHistList->Add(fh1JetPt);
+ fHistList->Add(fh2JetEtaPhi);
+ fHistList->Add(fh1JetN);
+ }
- // open input AOD
- if(fAODPathArray){
- Int_t rc = SelectAODfile();
- if(rc<0) return;
- }
- Int_t rc = OpenAODfile();
- if(rc<0) return;
+ if(fAODevent && fAODevent->FindListObject("mcparticles") && strlen(fMCparticlesBranch.Data())){
- } //end: embed mode with AOD
+ fh1MCTrackPt = new TH1F("fh1MCTrackPt","pT of MC extra tracks;p_{T};entries", 250, 0., 250.);
+ fh2MCTrackEtaPhi = new TH2F("fh2MCTrackEtaPhi","eta-phi distribution of MC extra tracks;#eta;#phi", 20, -1., 1., 60, 0., 2*TMath::Pi());
+ fh1MCTrackN = new TH1F("fh1MCTrackN", "nb. of MC extra tracks per event;nb. of tracks;entries",300, 0., 300.);
+
+ fHistList->Add(fh1MCTrackPt);
+ fHistList->Add(fh2MCTrackEtaPhi);
+ fHistList->Add(fh1MCTrackN);
+
+ }
+
+ fh1AODfile = new TH1I("fh1AODfile", "overview of opened AOD files from the array", 23, -0.5, 2299.5);
+ fh2AODevent = new TH2I("fh1AODevent","selected events;file;event", 25,-.5,2499.5,50,-.5,4999.5);
+ fHistList->Add(fh1AODfile);
+ fHistList->Add(fh2AODevent);
+
+ // =========== Switch on Sumw2 for all histos ===========
+ for (Int_t i=0; i<fHistList->GetEntries(); ++i) {
+ TH1 *h1 = dynamic_cast<TH1*>(fHistList->At(i));
+ if (h1){
+ h1->Sumw2();
+ continue;
+ }
+ }
+ TH1::AddDirectory(oldStatus);
- // connect output aod
- // create a new branch for extra tracks
- fAODout = AODEvent();
- if(!fAODout){
- AliError("Output AOD not found.");
- return;
- }
- if(!fAODout->FindListObject(fTrackBranch.Data()) && strlen(fTrackBranch.Data())){
- AliInfo(Form("Add AOD branch %s", fTrackBranch.Data()));
- TClonesArray *tracks = new TClonesArray("AliAODTrack",0);
- tracks->SetName(fTrackBranch.Data());
- AddAODBranch("TClonesArray", &tracks);
- }
- // create new branch for extra mcparticle if available as input
- if(fAODevent && fAODevent->FindListObject("mcparticles") && strlen(fMCparticlesBranch.Data())){
- AliInfo(Form("Add AOD branch %s", fMCparticlesBranch.Data()));
- TClonesArray *mcparticles = new TClonesArray("AliAODMCParticle",0);
- mcparticles->SetName(fMCparticlesBranch.Data());
- AddAODBranch("TClonesArray", &mcparticles);
- }
+ PostData(1, fHistList);
+}
+//__________________________________________________________________________
+void AliAnalysisTaskFastEmbedding::Init()
+{
+ // Initialization
+ if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::Init()");
+}
- //qa histograms
+//__________________________________________________________________________
+Bool_t AliAnalysisTaskFastEmbedding::UserNotify()
+{
+ // User defined Notify(), called once
+ if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::UserNotify()");
- OpenFile(1);
- fHistList = new TList();
+ // get total nb of events in tree (of this subjob)
+ AliInputEventHandler* inputHandler = (AliInputEventHandler*)((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
+ fInputEntries = inputHandler->GetTree()->GetEntriesFast();
+ AliInfo(Form("Total nb. of events: %d", fInputEntries));
- Bool_t oldStatus = TH1::AddDirectoryStatus();
- TH1::AddDirectory(kFALSE);
+ return kTRUE;
- fh1TrackPt = new TH1F("fh1TrackPt","pT of extra tracks;p_{T};entries", 300, 0., 300.);
- fh2TrackEtaPhi = new TH2F("fh2TrackEtaPhi","eta-phi distribution of extra tracks;#eta;#phi", 20, -1., 1., 60, 0., 2*TMath::Pi());
- fh1TrackN = new TH1F("fh1TrackN", "nb. of extra tracks per event;nb. of tracks;entries",300, 0., 3000.);
+}
- fHistList->Add(fh1TrackPt);
- fHistList->Add(fh2TrackEtaPhi);
- fHistList->Add(fh1TrackN);
-
- if(fAODevent && fAODevent->FindListObject("mcparticles") && strlen(fMCparticlesBranch.Data())){
+//__________________________________________________________________________
+void AliAnalysisTaskFastEmbedding::UserExec(Option_t *)
+{
+ // Analysis, called once per event
+ if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::UserExec()");
+
+ if(!fAODout){
+ AliError("Need output AOD, but is not connected.");
+ PostData(1, fHistList);
+ return;
+ }
- fh1MCTrackPt = new TH1F("fh1MCTrackPt","pT of MC extra tracks;p_{T};entries", 300, 0., 300.);
- fh2MCTrackEtaPhi = new TH2F("fh2MCTrackEtaPhi","eta-phi distribution of MC extra tracks;#eta;#phi", 20, -1., 1., 60, 0., 2*TMath::Pi());
- fh1MCTrackN = new TH1F("fh1MCTrackN", "nb. of MC extra tracks per event;nb. of tracks;entries",300, 0., 3000.);
+ fESD=dynamic_cast<AliESDEvent*>(InputEvent());
+ if (!fESD) {
+ AliError("ESD not available");
+ PostData(1, fHistList);
+ return;
+ }
- fHistList->Add(fh1MCTrackPt);
- fHistList->Add(fh2MCTrackEtaPhi);
- fHistList->Add(fh1MCTrackN);
+ // -- event selection --
+ fHistEvtSelection->Fill(1); // number of events before event selection
+
+ // physics selection
+ AliInputEventHandler* inputHandler = (AliInputEventHandler*)
+ ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
+ if(!(inputHandler->IsEventSelected() & fOfflineTrgMask)){
+ if(fDebug) Printf(" Trigger Selection: event REJECTED ... ");
+ fHistEvtSelection->Fill(2);
+ PostData(1, fHistList);
+ return;
+ }
+
+ // vertex selection
+ AliAODVertex* primVtx = fAODout->GetPrimaryVertex();
+ Int_t nTracksPrim = primVtx->GetNContributors();
+ if ((nTracksPrim < fMinContribVtx) ||
+ (primVtx->GetZ() < fVtxZMin) ||
+ (primVtx->GetZ() > fVtxZMax) ){
+ if(fDebug) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
+ fHistEvtSelection->Fill(3);
+ PostData(1, fHistList);
+ return;
+ }
- }
+ /** takes wrong value, since jet service tasks runs later
+ // event class selection (from jet helper task)
+ Int_t eventClass = AliAnalysisHelperJetTasks::EventClass();
+ if(fDebug) Printf("Event class %d", eventClass);
+ if (eventClass < fEvtClassMin || eventClass > fEvtClassMax){
+ fHistEvtSelection->Fill(4);
+ PostData(1, fHistList);
+ return;
+ }*/
+
+ // centrality selection
+ AliCentrality *cent = 0x0;
+ Float_t centValue = 0.;
+ if(fESD) cent = fESD->GetCentrality();
+ if(cent) centValue = cent->GetCentralityPercentile("V0M");
+ if(fDebug) printf("centrality: %f\n", centValue);
+ if (centValue < fCentMin || centValue > fCentMax){
+ fHistEvtSelection->Fill(4);
+ PostData(1, fHistList);
+ return;
+ }
- TH1::AddDirectory(oldStatus);
-}
+ /* ** not implemented **
+ // multiplicity due to input tracks
+ Int_t nInputTracks = GetNInputTracks();
+ if (nInputTracks < fNInputTracksMin || (fNInputTracksMax > -1 && nInputTracks > fNInputTracksMax)){
+ fHistEvtSelection->Fill(5);
+ PostData(1, fHistList);
+ return;
+ }
+ */
-//__________________________________________________________________________
-void AliAnalysisTaskFastEmbedding::Init()
-{
- // Initialization
- if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::Init()");
+ fHistEvtSelection->Fill(0); // accepted events
+ // -- end event selection --
-}
+ // connect aod out
+ TClonesArray *tracks = (TClonesArray*)(fAODout->FindListObject(fTrackBranch.Data()));
+ if(!tracks){
+ AliError("Extra track branch not found in output.");
+ PostData(1, fHistList);
+ return;
+ }
+ tracks->Delete();
+ Int_t nAODtracks=0;
+
+ TRef dummy;
+
+ // === embed mode with AOD ===
+ if(fEmbedMode==kAODFull || fEmbedMode==kAODJetTracks || fEmbedMode==kAODJet4Mom){
+ if(!fAODevent){
+ AliError("Need input AOD, but is not connected.");
+ PostData(1, fHistList);
+ return;
+ }
+
+
+
+ Bool_t useEntry = kFALSE;
+ while(!useEntry){ // protection need, if no event fulfills requirement
+
+ fAODEntry++; // go to next event
+ fCountEvents++;
+ if(fAODEntry>=fNevents) fAODEntry=0;
+
+ // open new aod file
+ if(fCountEvents>=fNevents){
+ if(!fAODPathArray){
+ AliDebug(AliLog::kDebug, "No list of AODs, keep current AOD.");
+ }
+ else {
+ AliDebug(AliLog::kDebug, "Select new AOD file ...");
+
+ fFileId = OpenAODfile();
+ if(fFileId<0) {
+ PostData(1, fHistList);
+ return;
+ }
+ fh1AODfile->Fill(fFileId);
+ if(fAODEntries!=fNevents){
+ AliError("File with list of AODs and file with nb. of entries inconsistent");
+ PostData(1, fHistList);
+ return;
+ }
+ }
+ }
+
+
+
+ // get next event
+ fAODtree->GetEvent(fAODEntry);
+
+ // get pt hard
+ if(mcHeader){
+ fPtHard = mcHeader->GetPtHard();
+ GetPtHard(kTRUE,fPtHard); // make value available for other tasks (e.g. jet response task)
+ AliDebug(AliLog::kDebug,Form("pt hard %.2f", fPtHard));
+ }
+ else{
+ AliDebug(AliLog::kDebug,"No mcHeader");
+ }
+ fPtHardBin = GetPtHardBin(fPtHard);
+
+ //Printf("pt hard (bin) %.2f (%d), xSection %.2e, trials %.2f",fPtHard, fPtHardBin, fXsection, fAvgTrials);
+
+ // fill event variables for each event
+ fh1Xsec->Fill(fPtHardBin,fXsection);
+ fh2PtHard->Fill(fPtHardBin,fPtHard);
+
+ fh1Trials->Fill(fPtHardBin, fAvgTrials);
+ fh1TrialsEvtSel->Fill(fPtHardBin);
+ fh2PtHardTrials->Fill(fPtHardBin,fPtHard,fAvgTrials);
+
+ // jet pt selection
+ if(fEvtSelecMode==kEventsJetPt){
+ Int_t nJets = fAODJets->GetEntries();
+ //Printf("nb. jets: %d",nJets);
+ for(Int_t ij=0; ij<nJets; ++ij){
+ AliAODJet *jet = dynamic_cast<AliAODJet*>(fAODJets->At(ij));
+ if(!jet) continue;
+ if( (fEvtSelMinJetPt==-1. || jet->Pt()>=fEvtSelMinJetPt)
+ && (fEvtSelMaxJetPt==-1. || jet->Pt()<=fEvtSelMaxJetPt)
+ && (jet->Eta()>=fEvtSelMinJetEta && jet->Eta()<=fEvtSelMaxJetEta)
+ && (jet->Phi()>=fEvtSelMinJetPhi && jet->Phi()<=fEvtSelMaxJetPhi)){
+ useEntry = kTRUE;
+ break;
+ }
+ }
+ }
+
+ // no selection
+ if(fEvtSelecMode==kEventsAll){
+ useEntry = kTRUE;
+ }
+
+ }
+ AliDebug(AliLog::kDebug,Form("Use entry %d from extra AOD.", fAODEntry));
+
+ fh2PtHardEvtSel->Fill(fPtHardBin,fPtHard);
+ fh2AODevent->Fill(fFileId,fAODEntry);
+
+ TClonesArray *mcpartIN = (TClonesArray*)(fAODevent->FindListObject("mcparticles"));
+ TClonesArray *mcpartOUT = 0x0;
+ if(mcpartIN){
+ mcpartOUT = (TClonesArray*)(fAODout->FindListObject(fMCparticlesBranch.Data()));
+ if(!mcpartOUT){
+ AliError("Extra MC particles branch not found in output.");
+ PostData(1, fHistList);
+ return;
+ }
+ mcpartOUT->Delete();
+ } else {
+ AliInfo("No extra MC particles found.");
+ }
+
+
+ if(fEmbedMode==kAODFull || fEmbedMode==kAODJetTracks){ // take all tracks or jet tracks
+ // loop over input tracks
+ // add to output aod
+ Int_t nTracks = 0;
+ Int_t nJets = fAODJets->GetEntries();
+ Int_t nSelectedJets = 0;
+ AliAODJet *leadJet = 0x0; // used for jet tracks only
+
+ if(fEmbedMode==kAODFull){
+ nTracks = fAODevent->GetNumberOfTracks();
+
+ for(Int_t ij=0; ij<nJets; ++ij){
+ AliAODJet *jet = dynamic_cast<AliAODJet*>(fAODJets->At(ij));
+ if(!jet) continue;
+ if( (fEvtSelMinJetPt==-1. || jet->Pt()>=fEvtSelMinJetPt)
+ && (fEvtSelMaxJetPt==-1. || jet->Pt()<=fEvtSelMaxJetPt)
+ && (jet->Eta()>=fEvtSelMinJetEta && jet->Eta()<=fEvtSelMaxJetEta)
+ && (jet->Phi()>=fEvtSelMinJetPhi && jet->Eta()<=fEvtSelMaxJetPhi)){
+
+ fh1JetPt->Fill(jet->Pt());
+ fh2JetEtaPhi->Fill(jet->Eta(), jet->Phi());
+ nSelectedJets++;
+
+ }
+ }
+ fh1JetN->Fill(nSelectedJets);
+ }
+
+ if(fEmbedMode==kAODJetTracks){
+ // look for leading jet within selection
+ // get reference tracks for this jet
+ Float_t leadJetPt = 0.;
+ for(Int_t ij=0; ij<nJets; ++ij){
+ AliAODJet *jet = dynamic_cast<AliAODJet*>(fAODJets->At(ij));
+ if(!jet) continue;
+ if( (fEvtSelMinJetPt==-1. || jet->Pt()>=fEvtSelMinJetPt)
+ && (fEvtSelMaxJetPt==-1. || jet->Pt()<=fEvtSelMaxJetPt)
+ && (jet->Eta()>=fEvtSelMinJetEta && jet->Eta()<=fEvtSelMaxJetEta)
+ && (jet->Phi()>=fEvtSelMinJetPhi && jet->Eta()<=fEvtSelMaxJetPhi)){
+ if(jet->Pt()>leadJetPt){
+ leadJetPt = jet->Pt();
+ leadJet = jet;
+ }
+ }
+ }
+ if(leadJet){
+ nTracks = leadJet->GetRefTracks()->GetEntriesFast();
+ fh1JetPt->Fill(leadJet->Pt());
+ fh2JetEtaPhi->Fill(leadJet->Eta(), leadJet->Pt());
+ fh1JetN->Fill(1);
+ }
+ }
+
+ fh1TrackN->Fill((Float_t)nTracks);
+
+ for(Int_t it=0; it<nTracks; ++it){
+ AliAODTrack *tmpTr = 0x0;
+ if(fEmbedMode==kAODFull) tmpTr = fAODevent->GetTrack(it);
+ if(fEmbedMode==kAODJetTracks) tmpTr = dynamic_cast<AliAODTrack*>(leadJet->GetRefTracks()->At(it));
+ if(!tmpTr) continue;
+ Double_t rd=rndm->Uniform(0.,1.);
+ if(rd>fExtraEffPb) continue;
+ tmpTr->SetStatus(AliESDtrack::kEmbedded);
+
+ new ((*tracks)[nAODtracks++]) AliAODTrack(*tmpTr);
+ dummy = (*tracks)[nAODtracks-1];
+
+ if(fTrackFilterMap<=0 || tmpTr->TestFilterBit(fTrackFilterMap)){
+ fh1TrackPt->Fill(tmpTr->Pt());
+ fh2TrackEtaPhi->Fill(tmpTr->Eta(), tmpTr->Phi());
+ }
+ }
-//__________________________________________________________________________
-void AliAnalysisTaskFastEmbedding::UserExec(Option_t *)
-{
- if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::UserExec()");
-
- if(!fAODout){
- AliError("Need output AOD, but is not connected.");
- return;
- }
-
- // connect aod out
- TClonesArray *tracks = (TClonesArray*)(fAODout->FindListObject(fTrackBranch.Data()));
- if(!tracks){
- AliError("Extra track branch not found in output.");
- return;
- }
- tracks->Delete();
- Int_t nAODtracks=0;
-
-
- TRef dummy;
-
- // === embed mode with AOD ===
- if(fEmbedMode==kAODFull || fEmbedMode==kAODJetTracks || fEmbedMode==kAODJet4Mom){
- if(!fAODevent){
- AliError("Need input AOD, but is not connected.");
- return;
- }
-
- // fetch jets
- TClonesArray *aodJets = 0;
- if(fJetBranch.Length()) aodJets = dynamic_cast<TClonesArray*>(fAODevent->FindListObject(fJetBranch.Data()));
- else aodJets = fAODevent->GetJets();
- if(!aodJets){
- AliError("Could not find jets in AOD. Check jet branch when indicated.");
- return;
- }
-
-
- Int_t nEvents = fAODtree->GetEntries();
-
- Bool_t useEntry = kFALSE;
- while(!useEntry){ // protection need, if no event fulfills requierment
- if(fEntry>nEvents){
- fEntry=0;
- if(!fAODPathArray){
- AliDebug(AliLog::kDebug, "Last event in AOD reached, start from entry 0 again.");
- }
- else {
- AliDebug(AliLog::kDebug, "Last event in AOD reached, select new AOD file ...");
- Int_t rc = SelectAODfile();
- if(rc<0) return;
-
- rc = OpenAODfile();
- if(rc<0) return;
- }
- }
-
- fAODtree->GetEvent(fEntry);
-
- // jet pt selection
- if(fEvtSelecMode==kEventsJetPt){
- Int_t nJets = aodJets->GetEntries();
- for(Int_t ij=0; ij<nJets; ++ij){
- AliAODJet *jet = dynamic_cast<AliAODJet*>(aodJets->At(ij));
- if(!jet) continue;
-
- if( (fEvtSelMinJetPt==-1. || jet->Pt()>=fEvtSelMinJetPt)
- && (fEvtSelMaxJetPt==-1. || jet->Pt()<=fEvtSelMaxJetPt)){
- useEntry = kTRUE;
- break;
- }
- }
- }
-
- // no selection
- if(fEvtSelecMode==kEventsAll){
- useEntry = kTRUE;
- }
-
- fEntry++;
- }
- AliInfo(Form("Use entry %d from extra AOD.", fEntry-1));
-
-
- TClonesArray *mcpartIN = (TClonesArray*)(fAODevent->FindListObject("mcparticles"));
- TClonesArray *mcpartOUT = 0x0;
- if(mcpartIN){
- mcpartOUT = (TClonesArray*)(fAODout->FindListObject(fMCparticlesBranch.Data()));
- mcpartOUT->Delete();
- } else {
- AliInfo("No extra MC particles found.");
- }
-
-
- if(fEmbedMode==kAODFull || fEmbedMode==kAODJetTracks){ // take all tracks or jet tracks
- // loop over input tracks
- // add to output aod
- Int_t nTracks = fAODevent->GetNumberOfTracks();
- fh1TrackN->Fill((Float_t)nTracks);
-
- for(Int_t it=0; it<nTracks; ++it){
- AliAODTrack *tmpTr = fAODevent->GetTrack(it);
- // ?? test filter bit, or something else??
+ if(mcpartIN){
+ Int_t nMCpart = mcpartIN->GetEntriesFast();
+ Int_t nPhysicalPrimary=0;
+ Int_t nAODmcpart=0;
+ for(Int_t ip=0; ip<nMCpart; ++ip){
+ AliAODMCParticle *tmpPart = (AliAODMCParticle*) mcpartIN->At(ip);
+ if(!tmpPart) continue;
if(fEmbedMode==kAODJetTracks){
// jet track? else continue
// not implemented yet
continue;
}
- new ((*tracks)[nAODtracks++]) AliAODTrack(*tmpTr);
- dummy = (*tracks)[nAODtracks-1];
-
- fh1TrackPt->Fill(tmpTr->Pt());
- fh2TrackEtaPhi->Fill(tmpTr->Eta(), tmpTr->Phi());
- }
-
- if(mcpartIN){
- Int_t nMCpart = mcpartIN->GetEntriesFast();
-
- Int_t nPhysicalPrimary=0;
- Int_t nAODmcpart=0;
- for(Int_t ip=0; ip<nMCpart; ++ip){
- AliAODMCParticle *tmpPart = (AliAODMCParticle*) mcpartIN->At(ip);
-
- if(fEmbedMode==kAODJetTracks){
- // jet track? else continue
- // not implemented yet
- continue;
- }
-
- new((*mcpartOUT)[nAODmcpart++]) AliAODMCParticle(*tmpPart);
- dummy = (*mcpartOUT)[nAODmcpart-1];
-
- if(tmpPart->IsPhysicalPrimary() && tmpPart->Charge()!=0. && tmpPart->Charge()!=-99. ){
- fh1MCTrackPt->Fill(tmpPart->Pt());
- fh2MCTrackEtaPhi->Fill(tmpPart->Eta(), tmpPart->Phi());
- nPhysicalPrimary++;
- }
- }
- fh1MCTrackN->Fill((Float_t)nPhysicalPrimary);
-
- }
- } // end: embed all tracks or jet tracks
-
-
- if(fEmbedMode==kAODJet4Mom){
-
- // loop over jets
- Int_t nJets = aodJets->GetEntries();
- fh1TrackN->Fill((Float_t)nJets);
- for(Int_t ij=0; ij<nJets; ++ij){
- AliAODJet *jet = dynamic_cast<AliAODJet*>(aodJets->At(ij));
- if(!jet) continue;
- AliAODTrack *tmpTr = (AliAODTrack*)jet;
-
- new ((*tracks)[nAODtracks++]) AliAODTrack(*tmpTr);
- dummy = (*tracks)[nAODtracks-1];
-
- fh1TrackPt->Fill(tmpTr->Pt());
- fh2TrackEtaPhi->Fill(tmpTr->Eta(), tmpTr->Phi());
- }
-
- } // end: embed jets as 4-momenta
-
-
- } //end: embed mode with AOD
-
-
- // === embed mode with toy ===
- if(fEmbedMode==kToyTracks){
- Int_t nT = (Int_t)(rndm->Uniform(fToyMinNbOfTracks, fToyMaxNbOfTracks)+0.5);
-
- fh1TrackN->Fill((Float_t)nT);
-
- for(Int_t i=0; i<nT; ++i){
-
- Double_t pt = -1.;
- if(fToyMinTrackPt!=fToyMaxTrackPt){
- if(fToyDistributionTrackPt==0){
- pt = rndm->Uniform(fToyMinTrackPt, fToyMaxTrackPt);
- } else {
- while(pt<fToyMinTrackPt||pt>fToyMaxTrackPt){
- pt = rndm->Exp(fToyDistributionTrackPt); // no power law yet!!
- pt += fToyMinTrackPt;
- }
- }
- } else {
- pt = fToyMinTrackPt;
- }
- Double_t eta = rndm->Uniform(fToyMinTrackEta,fToyMaxTrackEta);
- Double_t phi = rndm->Uniform(fToyMinTrackPhi,fToyMaxTrackPhi);
- phi = TVector2::Phi_0_2pi(phi);
-
- Printf("Add track #%d: pt %.2f, eta %.2f, phi %.2f", i, pt, eta, phi);
-
- Double_t theta = 2*TMath::ATan(1/TMath::Exp(eta));
- Float_t mom[3];
- mom[0] = pt;
- mom[1] = phi;
- mom[2] = theta;
-
- Float_t xyz[3];
- xyz[0] = -999.;
- xyz[1] = -999.;
- xyz[2] = -999.;
-
- AliAODTrack *tmpTr = new AliAODTrack( -999, // id
- -999, // label
- mom, // momentum[3]
- kFALSE, // cartesian
- xyz, // position
- kFALSE, // DCA
- NULL, // covMatrix,
- -99, // charge
- 0, // itsClusMap
- NULL, // pid
- NULL, // prodVertex
- kFALSE, // used for vertex fit
- kFALSE, // used for prim vtx fit
- AliAODTrack::kUndef, // type
- fToyFilterMap, // select info
- -999. // chi2 per NDF
- );
- tmpTr->SetFlags(1<<27);
-
- new((*tracks)[nAODtracks++]) AliAODTrack(*tmpTr);
- dummy = (*tracks)[nAODtracks-1];
-
- fh1TrackPt->Fill(pt);
- fh2TrackEtaPhi->Fill(eta,phi);
-
- delete tmpTr;
- }
- } //end: embed mode with toy
-
-
- PostData(1, fHistList);
+ if(tmpPart->IsPhysicalPrimary() && tmpPart->Charge()!=0. && tmpPart->Charge()!=-99. && tmpPart->Pt()>0.){
+ new((*mcpartOUT)[nAODmcpart++]) AliAODMCParticle(*tmpPart);
+ dummy = (*mcpartOUT)[nAODmcpart-1];
+
+ if(fDebug>10) printf("added track %d with pT=%.2f to extra branch\n",nAODmcpart,tmpPart->Pt());
+
+ fh1MCTrackPt->Fill(tmpPart->Pt());
+ fh2MCTrackEtaPhi->Fill(tmpPart->Eta(), tmpPart->Phi());
+ nPhysicalPrimary++;
+
+ }
+ }
+ fh1MCTrackN->Fill((Float_t)nPhysicalPrimary);
+
+ }
+ } // end: embed all tracks or jet tracks
+
+
+ if(fEmbedMode==kAODJet4Mom){
+
+ // loop over jets
+ Int_t nJets = fAODJets->GetEntries();
+ fh1TrackN->Fill((Float_t)nJets);
+ for(Int_t ij=0; ij<nJets; ++ij){
+ AliAODJet *jet = dynamic_cast<AliAODJet*>(fAODJets->At(ij));
+ if(!jet) continue;
+ AliAODTrack *tmpTr = (AliAODTrack*)jet;
+ tmpTr->SetFlags(AliESDtrack::kEmbedded);
+
+ new ((*tracks)[nAODtracks++]) AliAODTrack(*tmpTr);
+ dummy = (*tracks)[nAODtracks-1];
+
+ fh1TrackPt->Fill(tmpTr->Pt());
+ fh2TrackEtaPhi->Fill(tmpTr->Eta(), tmpTr->Phi());
+ }
+
+ } // end: embed jets as 4-momenta
+
+
+ } //end: embed mode with AOD
+
+
+ // === embed mode with toy ===
+ if(fEmbedMode==kToyTracks){
+ Int_t nT = (Int_t)(rndm->Uniform(fToyMinNbOfTracks, fToyMaxNbOfTracks)+0.5);
+
+ fh1TrackN->Fill((Float_t)nT);
+
+ for(Int_t i=0; i<nT; ++i){
+
+ Double_t pt = -1.;
+ if(fToyMinTrackPt!=fToyMaxTrackPt){
+ if(fToyDistributionTrackPt==0){
+ pt = rndm->Uniform(fToyMinTrackPt, fToyMaxTrackPt);
+ } else {
+ while(pt<fToyMinTrackPt||pt>fToyMaxTrackPt){
+ pt = rndm->Exp(fToyDistributionTrackPt); // no power law yet!!
+ pt += fToyMinTrackPt;
+ }
+ }
+ } else {
+ pt = fToyMinTrackPt;
+ }
+ Double_t eta = rndm->Uniform(fToyMinTrackEta,fToyMaxTrackEta);
+ Double_t phi = rndm->Uniform(fToyMinTrackPhi,fToyMaxTrackPhi);
+ phi = TVector2::Phi_0_2pi(phi);
+
+ if(fDebug) Printf("Add track #%d: pt %.2f, eta %.2f, phi %.2f", i, pt, eta, phi);
+
+ Double_t theta = 2*TMath::ATan(1/TMath::Exp(eta));
+ Float_t mom[3];
+ mom[0] = pt;
+ mom[1] = phi;
+ mom[2] = theta;
+
+ Float_t xyz[3];
+ xyz[0] = -999.;
+ xyz[1] = -999.;
+ xyz[2] = -999.;
+
+ AliAODTrack *tmpTr = new AliAODTrack( -999, // id
+ -999, // label
+ mom, // momentum[3]
+ kFALSE, // cartesian
+ xyz, // position
+ kFALSE, // DCA
+ NULL, // covMatrix,
+ -99, // charge
+ 0, // itsClusMap
+ // NULL, // pid
+ NULL, // prodVertex
+ kFALSE, // used for vertex fit
+ kFALSE, // used for prim vtx fit
+ AliAODTrack::kUndef, // type
+ fToyFilterMap, // select info
+ -999. // chi2 per NDF
+ );
+ tmpTr->SetFlags(AliESDtrack::kEmbedded);
+
+ new((*tracks)[nAODtracks++]) AliAODTrack(*tmpTr);
+ dummy = (*tracks)[nAODtracks-1];
+
+ fh1TrackPt->Fill(pt);
+ fh2TrackEtaPhi->Fill(eta,phi);
+
+ delete tmpTr;
+ }
+ } //end: embed mode with toy
+
+
+ PostData(1, fHistList);
}
//__________________________________________________________________________
void AliAnalysisTaskFastEmbedding::Terminate(Option_t *)
{
- if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::Terminate()");
+ // terminate
+ if(fDebug > 1) Printf("AliAnalysisTaskFastEmbedding::Terminate()");
- if(fAODfile && fAODfile->IsOpen())
- fAODfile->Close();
+ if(fAODfile && fAODfile->IsOpen())
+ fAODfile->Close();
}
//__________________________________________________________________________
Int_t AliAnalysisTaskFastEmbedding::GetJobID()
{
+ // gives back the alien subjob id, if available, else -1
+
Int_t id=-1;
const char* env = gSystem->Getenv("ALIEN_PROC_ID"); // GRID
//if(!env || !strlen(env)) env = gSystem->Getenv("LSB_JOBINDEX"); // GSI
if(env && strlen(env)){
- id= atoi(env);
- AliInfo(Form("Job index %d", id));
+ id= atoi(env);
+ AliInfo(Form("Job index %d", id));
}
else{
- AliInfo("Job index not found. Okay if running locally.");
+ AliInfo("Job index not found. Okay if running locally.");
}
return id;
}
//__________________________________________________________________________
+Int_t AliAnalysisTaskFastEmbedding::SelectAODfile()
+{
+ // select an AOD file from fAODPathArray
+
+ Int_t nFiles = fAODPathArray->GetEntries();
+
+ // choose randomly file and event
+ Int_t n = -1;
+ Float_t tmpProp = -1.;
+ Int_t pendingEvents = fInputEntries-Entry();
+ //Printf("input entries %d, entry %d, pending events %d", fInputEntries, (Int_t)Entry(), pendingEvents);
+ if(fAODEntriesArray){
+ while(rndm->Rndm()>tmpProp){
+ n = rndm->Integer(nFiles);
+ fAODEntries = fAODEntriesArray->At(n);
+ Int_t tmpEntries = fAODEntries<pendingEvents ? pendingEvents : fAODEntries;
+ tmpProp = fAODEntriesMax ? (Float_t)tmpEntries/fAODEntriesMax : 1.;
+ }
+ fAODEntry = (Int_t)(rndm->Uniform(fAODEntries));
+ }
+ else {
+ AliWarning("Number of entries in extra AODs not set!");
+ n = rndm->Integer(nFiles);
+ fAODEntry = 0;
+ }
-Int_t AliAnalysisTaskFastEmbedding::SelectAODfile(){
+ AliInfo(Form("Select AOD file %d", n));
+ TObjString *objStr = (TObjString*) fAODPathArray->At(n);
+ if(!objStr){
+ AliError("Could not get path of aod file.");
+ return -1;
+ }
+ fAODPath = objStr->GetString();
- Int_t nFiles = fAODPathArray->GetEntries();
- Int_t n = rndm->Integer(nFiles);
- AliInfo(Form("Select AOD file %d", n));
- TObjString *objStr = (TObjString*) fAODPathArray->At(n);
- if(!objStr){
- AliError("Could not get path of aod file.");
- return -1;
- }
- fAODPath = objStr->GetString();
+ return n;
+}
+
+//__________________________________________________________________________
+Int_t AliAnalysisTaskFastEmbedding::OpenAODfile(Int_t trial)
+{
+ // select and open an AOD file from fAODPathArray
+
+ if(fAODPathArray) fFileId = SelectAODfile();
+ if(fFileId<0) return -1;
+
+ if (!gGrid) {
+ AliInfo("Trying to connect to AliEn ...");
+ TGrid::Connect("alien://");
+ }
+
+ TDirectory *owd = gDirectory;
+ if (fAODfile && fAODfile->IsOpen()) fAODfile->Close();
+ fAODfile = TFile::Open(fAODPath.Data(),"TIMEOUT=180");
+ owd->cd();
+ if(!fAODfile){
+
+ fFileId = -1;
+ if(fAODPathArray){
+ if(trial<=3){
+ AliError(Form("Could not open AOD file %s (trial %d), try again ...", fAODPath.Data(), trial));
+ fFileId = OpenAODfile(trial+1);
+ } else {
+ AliError(Form("Could not open AOD file %s, give up ...",fAODPath.Data()));
+ return -1;
+ }
+ } else {
+ AliError(Form("Could not open AOD file %s.", fAODPath.Data()));
+ }
+
+ return fFileId;
+ }
+
+ fAODtree = (TTree*)fAODfile->Get("aodTree");
+
+ if(!fAODtree){
+ AliError("AOD tree not found.");
+ return -1;
+ }
- return n;
+ /*
+ fAODtree->SetBranchStatus("*",0);
+ fAODtree->SetBranchStatus("header",1);
+ fAODtree->SetBranchStatus("tracks*",1);
+ fAODtree->SetBranchStatus("jets*",1);
+ fAODtree->SetBranchStatus("clusters*",1);
+ fAODtree->SetBranchStatus("mcparticles*",1);
+ fAODtree->SetBranchStatus("mcHeader*",1);
+ */
+
+ delete fAODevent;
+ fAODevent = new AliAODEvent();
+ fAODevent->ReadFromTree(fAODtree);
+
+
+ // fetch header, jets, etc. from new file
+ fNevents = fAODtree->GetEntries();
+ mcHeader = (AliAODMCHeader*)fAODevent->FindListObject(AliAODMCHeader::StdBranchName());
+
+ if(fJetBranch.Length()) fAODJets = dynamic_cast<TClonesArray*>(fAODevent->FindListObject(fJetBranch.Data()));
+ else fAODJets = fAODevent->GetJets();
+ if(!fAODJets){
+ AliError("Could not find jets in AOD. Check jet branch when indicated.");
+ return -1;
+ }
+
+ TFile *curfile = fAODtree->GetCurrentFile();
+ if (!curfile) {
+ AliError("No current file.");
+ return -1;
+ }
+
+ Float_t trials = 1.;
+ Float_t xsec = 0.;
+ PythiaInfoFromFile(curfile->GetName(),xsec,trials);
+ fXsection = xsec;
+
+ // construct a poor man average trials
+ Float_t nEntries = (Float_t)fAODtree->GetTree()->GetEntries();
+ if(trials>=nEntries && nEntries>0.)fAvgTrials = trials/nEntries;
+
+ if(fFileId>=0){
+ AliInfo(Form("Read successfully AOD event from file %d",fFileId));
+ }
+
+ fCountEvents=0; // new file, reset counter
+
+ return fFileId; // file position in AOD path array, if array available
+}
+
+
+//____________________________________________________________________________
+Float_t AliAnalysisTaskFastEmbedding::GetPtHard(Bool_t bSet, Float_t newValue){
+
+ // static stored, available for other tasks in train
+
+ static Float_t ptHard = -1.;
+ if(bSet) ptHard = newValue;
+
+ return ptHard;
}
-Int_t AliAnalysisTaskFastEmbedding::OpenAODfile(){
+//____________________________________________________________________________
+Int_t AliAnalysisTaskFastEmbedding::GetPtHardBin(Double_t ptHard){
- TDirectory *owd = gDirectory;
- fAODfile = TFile::Open(fAODPath.Data());
- owd->cd();
- if(!fAODfile){
- AliError("Could not open AOD file.");
- return -1;
- }
+ // returns pt hard bin (for LHC10e14, LHC11a1x, LHC11a2x simulations)
- fAODtree = (TTree*)fAODfile->Get("aodTree");
+ const Int_t nBins = 10;
+ Double_t binLimits[nBins] = { 5., 11., 21., 36., 57., 84., 117., 156., 200., 249. }; // lower limits
- if(!fAODtree){
- AliError("AOD tree not found.");
- return -1;
- }
+ Int_t bin = -1;
+ while(bin<nBins-1 && binLimits[bin+1]<ptHard){
+ bin++;
+ }
- fAODevent = new AliAODEvent();
- fAODevent->ReadFromTree(fAODtree);
- if(!fAODevent){
- AliError("AOD event not found.");
- return -1;
- }
+ return bin;
+}
- return 1;
+//____________________________________________________________________________
+Bool_t AliAnalysisTaskFastEmbedding::PythiaInfoFromFile(const char* currFile,Float_t &fXsec,Float_t &fTrials){
+ //
+ // get the cross section and the trails either from pyxsec.root or from pysec_hists.root
+ // This should provide the path to the AOD/ESD file
+ // (taken from PWG4/JetTasks/AliAnalysisHelperJetTasks)
+
+ TString file(currFile);
+ fXsec = 0;
+ fTrials = 1;
+
+ if(file.Contains("root_archive.zip#")){
+ Ssiz_t pos1 = file.Index("root_archive",12,0,TString::kExact);
+ Ssiz_t pos = file.Index("#",1,pos1,TString::kExact);
+ Ssiz_t pos2 = file.Index(".root",5,TString::kExact);
+ file.Replace(pos+1,pos2-pos1,"");
+ // file.Replace(pos+1,20,"");
+ }
+ else {
+ // not an archive take the basename....
+ file.ReplaceAll(gSystem->BaseName(file.Data()),"");
+ }
+
+ 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
+ if(!fxsec){
+ // next trial fetch the histgram file
+ fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec_hists.root"));
+ if(!fxsec){
+ // not a severe condition but inciate that we have no information
+ AliDebug(AliLog::kDebug,Form("Neither pyxsec.root nor pyxsec_hists.root found in %s",file.Data()));
+ return kFALSE;
+ }
+ else{
+ // find the tlist we want to be independtent of the name so use the Tkey
+ TKey* key = (TKey*)fxsec->GetListOfKeys()->At(0);
+ if(!key){
+ fxsec->Close();
+ return kFALSE;
+ }
+ TList *list = dynamic_cast<TList*>(key->ReadObj());
+ if(!list){
+ fxsec->Close();
+ return kFALSE;
+ }
+ fXsec = ((TProfile*)list->FindObject("h1Xsec"))->GetBinContent(1);
+ fTrials = ((TH1F*)list->FindObject("h1Trials"))->GetBinContent(1);
+ fxsec->Close();
+ }
+ } // no tree pyxsec.root
+ else {
+ TTree *xtree = (TTree*)fxsec->Get("Xsection");
+ if(!xtree){
+ fxsec->Close();
+ return kFALSE;
+ }
+ UInt_t ntrials = 0;
+ Double_t xsection = 0;
+ xtree->SetBranchAddress("xsection",&xsection);
+ xtree->SetBranchAddress("ntrials",&ntrials);
+ xtree->GetEntry(0);
+ fTrials = ntrials;
+ fXsec = xsection;
+ fxsec->Close();
+ }
+ return kTRUE;
}
+
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff