#include "TH2F.h"
#include "TString.h"
#include "THnSparse.h"
+#include "TProfile.h"
+#include "TFile.h"
+#include "TKey.h"
#include "AliAODInputHandler.h"
#include "AliAODHandler.h"
,fJetTypeGen(0)
,fJetTypeRecEff(0)
,fFilterMask(0)
+ ,fUsePhysicsSelection(kTRUE)
,fTrackPtCut(0)
,fTrackEtaMin(0)
,fTrackEtaMax(0)
,fh1VertexNContributors(0)
,fh1VertexZ(0)
,fh1EvtMult(0)
+ ,fh1Xsec(0)
+ ,fh1Trials(0)
+ ,fh1PtHard(0)
+ ,fh1PtHardTrials(0)
,fh1nRecJetsCuts(0)
,fh1nGenJets(0)
,fh1nRecEffJets(0)
,fJetTypeGen(0)
,fJetTypeRecEff(0)
,fFilterMask(0)
+ ,fUsePhysicsSelection(kTRUE)
,fTrackPtCut(0)
,fTrackEtaMin(0)
,fTrackEtaMax(0)
,fh1VertexNContributors(0)
,fh1VertexZ(0)
,fh1EvtMult(0)
+ ,fh1Xsec(0)
+ ,fh1Trials(0)
+ ,fh1PtHard(0)
+ ,fh1PtHardTrials(0)
,fh1nRecJetsCuts(0)
,fh1nGenJets(0)
,fh1nRecEffJets(0)
,fJetTypeGen(copy.fJetTypeGen)
,fJetTypeRecEff(copy.fJetTypeRecEff)
,fFilterMask(copy.fFilterMask)
+ ,fUsePhysicsSelection(copy.fUsePhysicsSelection)
,fTrackPtCut(copy.fTrackPtCut)
,fTrackEtaMin(copy.fTrackEtaMin)
,fTrackEtaMax(copy.fTrackEtaMax)
,fh1VertexNContributors(copy.fh1VertexNContributors)
,fh1VertexZ(copy.fh1VertexZ)
,fh1EvtMult(copy.fh1EvtMult)
+ ,fh1Xsec(copy.fh1Xsec)
+ ,fh1Trials(copy.fh1Trials)
+ ,fh1PtHard(copy.fh1PtHard)
+ ,fh1PtHardTrials(copy.fh1PtHardTrials)
,fh1nRecJetsCuts(copy.fh1nRecJetsCuts)
,fh1nGenJets(copy.fh1nGenJets)
,fh1nRecEffJets(copy.fh1nRecEffJets)
fJetTypeGen = o.fJetTypeGen;
fJetTypeRecEff = o.fJetTypeRecEff;
fFilterMask = o.fFilterMask;
+ fUsePhysicsSelection = o.fUsePhysicsSelection;
fTrackPtCut = o.fTrackPtCut;
fTrackEtaMin = o.fTrackEtaMin;
fTrackEtaMax = o.fTrackEtaMax;
fh1VertexNContributors = o.fh1VertexNContributors;
fh1VertexZ = o.fh1VertexZ;
fh1EvtMult = o.fh1EvtMult;
+ fh1Xsec = o.fh1Xsec;
+ fh1Trials = o.fh1Trials;
+ fh1PtHard = o.fh1PtHard;
+ fh1PtHardTrials = o.fh1PtHardTrials;
fh1nRecJetsCuts = o.fh1nRecJetsCuts;
fh1nGenJets = o.fh1nGenJets;
fh1nRecEffJets = o.fh1nRecEffJets;
}
+
+Bool_t AliAnalysisTaskFragmentationFunction::Notify()
+{
+ //
+ // Implemented Notify() to read the cross sections
+ // and number of trials from pyxsec.root
+ // (taken from AliAnalysisTaskJetSpectrum)
+ //
+ TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
+ Double_t xsection = 0;
+ UInt_t ntrials = 0;
+ Float_t ftrials = 0;
+ if(tree){
+ TFile *curfile = tree->GetCurrentFile();
+ if (!curfile) {
+ Error("Notify","No current file");
+ return kFALSE;
+ }
+ if(!fh1Xsec||!fh1Trials){
+ Printf("%s%d No Histogram fh1Xsec",(char*)__FILE__,__LINE__);
+ return kFALSE;
+ }
+
+ TString fileName(curfile->GetName());
+ if(fileName.Contains("AliESDs.root")){
+ fileName.ReplaceAll("AliESDs.root", "");
+ }
+ else if(fileName.Contains("AliAOD.root")){
+ fileName.ReplaceAll("AliAOD.root", "");
+ }
+ else if(fileName.Contains("AliAODs.root")){
+ fileName.ReplaceAll("AliAODs.root", "");
+ }
+ else if(fileName.Contains("galice.root")){
+ // for running with galice and kinematics alone...
+ fileName.ReplaceAll("galice.root", "");
+ }
+ TFile *fxsec = TFile::Open(Form("%s%s",fileName.Data(),"pyxsec.root"));
+ if(!fxsec){
+ if(fDebug>0)Printf("%s:%d %s not found in the Input",(char*)__FILE__,__LINE__,Form("%s%s",fileName.Data(),"pyxsec.root"));
+ // next trial fetch the histgram file
+ fxsec = TFile::Open(Form("%s%s",fileName.Data(),"pyxsec_hists.root"));
+ if(!fxsec){
+ // not a severe condition
+ if(fDebug>0)Printf("%s:%d %s not found in the Input",(char*)__FILE__,__LINE__,Form("%s%s",fileName.Data(),"pyxsec_hists.root"));
+ return kTRUE;
+ }
+ 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){
+ if(fDebug>0)Printf("%s:%d key not found in the file",(char*)__FILE__,__LINE__);
+ return kTRUE;
+ }
+ TList *list = dynamic_cast<TList*>(key->ReadObj());
+ if(!list){
+ if(fDebug>0)Printf("%s:%d key is not a tlist",(char*)__FILE__,__LINE__);
+ return kTRUE;
+ }
+ xsection = ((TProfile*)list->FindObject("h1Xsec"))->GetBinContent(1);
+ ftrials = ((TH1F*)list->FindObject("h1Trials"))->GetBinContent(1);
+ }
+ }
+ else{
+ TTree *xtree = (TTree*)fxsec->Get("Xsection");
+ if(!xtree){
+ Printf("%s:%d tree not found in the pyxsec.root",(char*)__FILE__,__LINE__);
+ }
+ xtree->SetBranchAddress("xsection",&xsection);
+ xtree->SetBranchAddress("ntrials",&ntrials);
+ ftrials = ntrials;
+ xtree->GetEntry(0);
+ }
+ fh1Xsec->Fill("<#sigma>",xsection);
+ fh1Trials->Fill("#sum{ntrials}",ftrials);
+ }
+
+ return kTRUE;
+}
+
+
+
//__________________________________________________________________
void AliAnalysisTaskFragmentationFunction::UserCreateOutputObjects()
{
fh1EvtSelection = new TH1F("fh1EvtSelection", "Event Selection", 6, -0.5, 5.5);
fh1VertexNContributors = new TH1F("fh1VertexNContributors", "Vertex N contributors", 11,-.5, 10.5);
fh1VertexZ = new TH1F("fh1VertexZ", "Vertex z distribution", 30, -15., 15.);
- fh1EvtMult = new TH1F("fh1EvtMult","Event multiplicity, track pT cut > 150 MeV/c, |#eta| < 0.9",100,0.,100.);
+ fh1EvtMult = new TH1F("fh1EvtMult","Event multiplicity, track pT cut > 150 MeV/c, |#eta| < 0.9",120,0.,120.);
+
+ fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
+ fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
+ fh1Trials = new TH1F("fh1Trials","trials from pyxsec.root",1,0,1);
+ fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
+ fh1PtHard = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",350,-.5,349.5);
+ fh1PtHardTrials = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",350,-.5,349.5);
+
fh1nRecJetsCuts = new TH1F("fh1nRecJetsCuts","reconstructed jets per event",10,-0.5,9.5);
fh1nGenJets = new TH1F("fh1nGenJets","generated jets per event",10,-0.5,9.5);
fh1nRecEffJets = new TH1F("fh1nRecEffJets","reconstruction effiency: jets per event",10,-0.5,9.5);
fFFHistosGen->AddToOutput(fCommonHistList);
fFFHistosGenLeading->AddToOutput(fCommonHistList);
fFFHistosGenLeadingTrack->AddToOutput(fCommonHistList);
+
+ fCommonHistList->Add(fh1Xsec);
+ fCommonHistList->Add(fh1Trials);
+ fCommonHistList->Add(fh1PtHard);
+ fCommonHistList->Add(fh1PtHardTrials);
}
}
if(saveLevel>1){
fh1EvtSelection->Fill(1.);
} else {
fh1EvtSelection->Fill(0.);
- if(inputHandler->InheritsFrom("AliESDInputHandler")){ // PhysicsSelection only with ESD input
+ if(inputHandler->InheritsFrom("AliESDInputHandler") && fUsePhysicsSelection){ // PhysicsSelection only with ESD input
if (fDebug > 1 ) Printf(" Trigger Selection: event REJECTED ... ");
PostData(1, fCommonHistList);
return;
}
if (fDebug > 1) Printf("%s:%d primary vertex selection: event ACCEPTED ...",(char*)__FILE__,__LINE__);
fh1EvtSelection->Fill(5.);
+
+
+ //___ get MC information __________________________________________________________________
+
+ Double_t ptHard = 0.;
+ Double_t nTrials = 1; // trials for MC trigger weight for real data
+
+ AliGenPythiaEventHeader* pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(fMCEvent);
+ if(!pythiaGenHeader){
+ if(fJetTypeGen != kJetsUndef && fTrackTypeGen != kTrackUndef){
+ Printf("%s:%d no pythiaGenHeader found", (char*)__FILE__,__LINE__);
+ return;
+ }
+ } else {
+ nTrials = pythiaGenHeader->Trials();
+ ptHard = pythiaGenHeader->GetPtHard();
+
+ fh1PtHard->Fill(ptHard);
+ fh1PtHardTrials->Fill(ptHard,nTrials);
+ }
//___ fetch jets __________________________________________________________________________
class TList;
class TH1F;
class TH2F;
+class TProfile;
#include "THnSparse.h"
#include "AliAnalysisTaskSE.h"
virtual void LocalInit() {Init();}
virtual void UserExec(Option_t *option);
virtual void Terminate(Option_t* );
+ virtual Bool_t Notify();
virtual void SetTrackTypeGen(Int_t i){fTrackTypeGen = i;}
virtual void SetJetTypeGen(Int_t i){fJetTypeGen = i;}
{fTrackPtCut = trackPt; fTrackEtaMin = trackEtaMin; fTrackEtaMax = trackEtaMax;
fTrackPhiMin = trackPhiMin; fTrackPhiMax = trackPhiMax;}
virtual void SetFilterMask(UInt_t i) {fFilterMask = i;}
+ virtual void UsePhysicsSelection(Bool_t b) {fUsePhysicsSelection = b;}
virtual void SetJetCuts(Float_t jetPt = 5., Float_t jetEtaMin = -0.5, Float_t jetEtaMax = 0.5,
Float_t jetPhiMin = 0., Float_t jetPhiMax = 2*TMath::Pi())
{fJetPtCut = jetPt; fJetEtaMin = jetEtaMin; fJetEtaMax = jetEtaMax;
void SetHighPtThreshold(Float_t pt = 5.) { fQATrackHighPtThreshold = pt; }
- void SetFFHistoBins(Int_t nJetPt = 55, Float_t jetPtMin = 5, Float_t jetPtMax = 60,
- Int_t nPt = 70, Float_t ptMin = 0., Float_t ptMax = 70.,
+ void SetFFHistoBins(Int_t nJetPt = 245, Float_t jetPtMin = 5, Float_t jetPtMax = 250,
+ Int_t nPt = 200, Float_t ptMin = 0., Float_t ptMax = 200.,
Int_t nXi = 70, Float_t xiMin = 0., Float_t xiMax = 7.,
Int_t nZ = 22, Float_t zMin = 0., Float_t zMax = 1.1)
{ fFFNBinsJetPt = nJetPt; fFFJetPtMin = jetPtMin; fFFJetPtMax = jetPtMax;
fFFNBinsXi = nXi; fFFXiMin = xiMin; fFFXiMax = xiMax;
fFFNBinsZ = nZ; fFFZMin = zMin; fFFZMax = zMax; }
- void SetQAJetHistoBins(Int_t nPt = 200, Float_t ptMin = 0., Float_t ptMax = 200.,
+ void SetQAJetHistoBins(Int_t nPt = 300, Float_t ptMin = 0., Float_t ptMax = 300.,
Int_t nEta = 20, Float_t etaMin = -1.0, Float_t etaMax = 1.0,
Int_t nPhi = 60, Float_t phiMin = 0., Float_t phiMax = 2*TMath::Pi())
{ fQAJetNBinsPt = nPt; fQAJetPtMin = ptMin; fQAJetPtMax = ptMax;
fQAJetNBinsEta = nEta; fQAJetEtaMin = etaMin; fQAJetEtaMax = etaMax;
fQAJetNBinsPhi = nPhi; fQAJetPhiMin = phiMin; fQAJetPhiMax = phiMax; }
- void SetQATrackHistoBins(Int_t nPt = 70, Float_t ptMin = 0., Float_t ptMax = 70.,
+ void SetQATrackHistoBins(Int_t nPt = 200, Float_t ptMin = 0., Float_t ptMax = 200.,
Int_t nEta = 20, Float_t etaMin = -1.0, Float_t etaMax = 1.0,
Int_t nPhi = 60, Float_t phiMin = 0., Float_t phiMax = 2*TMath::Pi())
{ fQATrackNBinsPt = nPt; fQATrackPtMin = ptMin; fQATrackPtMax = ptMax;
fQATrackNBinsEta = nEta; fQATrackEtaMin = etaMin; fQATrackEtaMax = etaMax;
fQATrackNBinsPhi = nPhi; fQATrackPhiMin = phiMin; fQATrackPhiMax = phiMax; }
- void SetIJHistoBins(Int_t nJetPt = 55, Float_t jetPtMin = 5, Float_t jetPtMax = 60,
- Int_t nPt = 70, Float_t ptMin = 0., Float_t ptMax = 70.,
+ void SetIJHistoBins(Int_t nJetPt = 245, Float_t jetPtMin = 5, Float_t jetPtMax = 250,
+ Int_t nPt = 200, Float_t ptMin = 0., Float_t ptMax = 200.,
Int_t nZ = 22, Float_t zMin = 0., Float_t zMax = 1.1,
Int_t nCosTheta = 100, Float_t costhetaMin = 0., Float_t costhetaMax = 1.,
Int_t nTheta = 200, Float_t thetaMin = -0.5, Float_t thetaMax = 1.5,
fIJNBinsTheta = nTheta; fIJThetaMin = thetaMin; fIJThetaMax = thetaMax;
fIJNBinsJt = nJt; fIJJtMin = jtMin; fIJJtMax = jtMax; }
- void SetDiJetHistoBins(Int_t nJetInvMass = 55, Float_t jetInvMassMin = 5, Float_t jetInvMassMax = 60,
- Int_t nJetPt = 55, Float_t jetPtMin = 5, Float_t jetPtMax = 60,
- Int_t nPt = 70, Float_t ptMin = 0., Float_t ptMax = 70.,
+ void SetDiJetHistoBins(Int_t nJetInvMass = 245, Float_t jetInvMassMin = 5, Float_t jetInvMassMax = 250,
+ Int_t nJetPt = 245, Float_t jetPtMin = 5, Float_t jetPtMax = 250,
+ Int_t nPt = 200, Float_t ptMin = 0., Float_t ptMax = 200.,
Int_t nXi = 70, Float_t xiMin = 0., Float_t xiMax = 7.,
Int_t nZ = 22, Float_t zMin = 0., Float_t zMax = 1.1)
{
fDiJetNBinsZ = nZ; fDiJetZMin = zMin; fDiJetZMax = zMax;
}
- void SetQADiJetHistoBins(Int_t nInvMass = 55, Float_t invMassMin = 5., Float_t invMassMax = 60.,
- Int_t nJetPt = 55, Float_t jetPtMin = 5, Float_t jetPtMax = 60,
+ void SetQADiJetHistoBins(Int_t nInvMass = 245, Float_t invMassMin = 5., Float_t invMassMax = 250.,
+ Int_t nJetPt = 245, Float_t jetPtMin = 5, Float_t jetPtMax = 250,
Int_t nDeltaPhi = 100, Float_t deltaPhiMin = 0., Float_t deltaPhiMax = TMath::Pi(),
Int_t nDeltaEta = 22, Float_t deltaEtaMin = 0., Float_t deltaEtaMax = 1.1,
Int_t nDeltaPt = 100, Float_t deltaPtMin = 0., Float_t deltaPtMax = 100.)
void FillJetTrackRecEffHisto(THnSparse* histo,Double_t jetPhi,Double_t jetEta,Double_t jetPt,TList* jetTrackList, TList* tracksGen,
TArrayI& indexAODTr,TArrayS& isGenPrim);
-
-
- private:
// Consts
enum {kTrackUndef=0, kTrackAOD, kTrackAODQualityCuts, kTrackAODCuts, kTrackKineAll, kTrackKineCharged, kTrackKineChargedAcceptance,
kTrackAODMCAll, kTrackAODMCCharged, kTrackAODMCChargedAcceptance};
enum {kJetsUndef=0, kJetsRec, kJetsRecAcceptance, kJetsGen, kJetsGenAcceptance, kJetsKine, kJetsKineAcceptance};
-
+
+
+ private:
Int_t GetListOfTracks(TList* list, Int_t type);
Int_t GetListOfJets(TList* list, Int_t type);
Int_t fJetTypeRecEff; // type of jets used for filling reconstruction efficiency histos
UInt_t fFilterMask; // filter bit for selected tracks
+ Bool_t fUsePhysicsSelection; // switch for event selection
// track cuts
Float_t fTrackPtCut; // track transverse momentum cut
TH1F *fh1VertexNContributors; //! NContributors to prim vertex
TH1F *fh1VertexZ; //! prim vertex z distribution
TH1F *fh1EvtMult; //! number of reconstructed tracks after cuts
+
+ TProfile* fh1Xsec; //! pythia cross section and trials
+ TH1F* fh1Trials; //! sum of trials
+ TH1F* fh1PtHard; //! pt hard of the event
+ TH1F* fh1PtHardTrials; //! pt hard of the event
+
TH1F *fh1nRecJetsCuts; //! number of jets from reconstructed tracks per event
TH1F *fh1nGenJets; //! number of jets from generated tracks per event
TH1F *fh1nRecEffJets; //! number of jets for reconstruction eff per event
THnSparseF *fhnSingleTrackRecEffHisto; //! track reconstruction efficiency
THnSparseF *fhnJetTrackRecEffHisto; //! reconstruction efficiency jet tracks
- ClassDef(AliAnalysisTaskFragmentationFunction, 4);
+ ClassDef(AliAnalysisTaskFragmentationFunction, 5);
};
#endif
// only reconstructed (default)
if(iFlag&(1<<0)) ff = AddTaskFragmentationFunction("jets", "", "", "", filterMask);
- // MC tracks in acceptance, MC jets in acceptance
- if(iFlag&(1<<1)) ff = AddTaskFragmentationFunction("jets", "AODMC2b", "AODMCb", "AODMCb", filterMask);
+ // charged MC tracks and jets
+ if(iFlag&(1<<1)) ff = AddTaskFragmentationFunction("jets", "jetsAODMC2_UA04", "AODMC", "AODMC2", filterMask);
+ // charged MC tracks and jets with acceptance cuts
+ if(iFlag&(1<<2)) ff = AddTaskFragmentationFunction("jets", "jetsAODMC2_UA04", "AODMCb", "AODMC2b", filterMask);
// kine tracks in acceptance, pythia jets in acceptance
- if(iFlag&(1<<2)) ff = AddTaskFragmentationFunction("jets", "", "KINEb", "KINEb", filterMask);
+ if(iFlag&(1<<3)) ff = AddTaskFragmentationFunction("jets", "", "KINEb", "KINEb", filterMask);
// reconstructed charged tracks after cuts, MC jets in acceptance
- if(iFlag&(1<<3)) ff = AddTaskFragmentationFunction("jets", "AODMC2b", "AODMCb", "AOD2b", filterMask);
+ if(iFlag&(1<<4)) ff = AddTaskFragmentationFunction("jets", "jetsMC2b", "AODMCb", "AOD2b", filterMask);
return ff;
}