#include <TLorentzVector.h>
#include <TClonesArray.h>
#include "TDatabasePDG.h"
+#include <TGrid.h>
#include "AliAnalysisTaskJetCluster.h"
#include "AliAnalysisManager.h"
fUseAODMCInput(kFALSE),
fUseBackgroundCalc(kFALSE),
fEventSelection(kFALSE),
+ fRequireVZEROAC(kFALSE),
+ fRequireTZEROvtx(kFALSE),
+ fUseHFcuts(kFALSE),
fFilterMask(0),
fFilterMaskBestPt(0),
fFilterType(0),
fAvgTrials(1),
fExternalWeight(1),
fTrackEtaWindow(0.9),
+ fRequireITSRefit(0),
+ fApplySharedClusterCut(0),
fRecEtaWindow(0.5),
fTrackPtCut(0.),
fJetOutputMinPt(0.150),
fVtxR2Cut(1),
fCentCutUp(0),
fCentCutLo(0),
+ fStoreRhoLeadingTrackCorr(kFALSE),
fNonStdBranch(""),
fBackgroundBranch(""),
fNonStdFile(""),
fhEffH3(0x0),
fUseTrPtResolutionSmearing(kFALSE),
fUseDiceEfficiency(kFALSE),
- fDiceEfficiencyMinPt(0.),
+ fDiceEfficiencyMinPt(-1.),
fUseTrPtResolutionFromOADB(kFALSE),
fUseTrEfficiencyFromOADB(kFALSE),
fPathTrPtResolution(""),
fPathTrEfficiency(""),
fChangeEfficiencyFraction(0.),
+ fEfficiencyFixed(1.),
fRparam(1.0),
fAlgorithm(fastjet::kt_algorithm),
fStrategy(fastjet::Best),
fh2PtNchNRan(0x0),
fh2TracksLeadingJetPhiPtRan(0x0),
fh2TracksLeadingJetPhiPtWRan(0x0),
+ fh3CentvsRhoLeadingTrackPt(0x0),
+ fh3CentvsSigmaLeadingTrackPt(0x0),
+ fh3MultvsRhoLeadingTrackPt(0x0),
+ fh3MultvsSigmaLeadingTrackPt(0x0),
+ fh3CentvsRhoLeadingTrackPtQ1(0x0),
+ fh3CentvsRhoLeadingTrackPtQ2(0x0),
+ fh3CentvsRhoLeadingTrackPtQ3(0x0),
+ fh3CentvsRhoLeadingTrackPtQ4(0x0),
+ fh3CentvsSigmaLeadingTrackPtQ1(0x0),
+ fh3CentvsSigmaLeadingTrackPtQ2(0x0),
+ fh3CentvsSigmaLeadingTrackPtQ3(0x0),
+ fh3CentvsSigmaLeadingTrackPtQ4(0x0),
+ fh3MultvsRhoLeadingTrackPtQ1(0x0),
+ fh3MultvsRhoLeadingTrackPtQ2(0x0),
+ fh3MultvsRhoLeadingTrackPtQ3(0x0),
+ fh3MultvsRhoLeadingTrackPtQ4(0x0),
+ fh3MultvsSigmaLeadingTrackPtQ1(0x0),
+ fh3MultvsSigmaLeadingTrackPtQ2(0x0),
+ fh3MultvsSigmaLeadingTrackPtQ3(0x0),
+ fh3MultvsSigmaLeadingTrackPtQ4(0x0),
+ fh3CentvsDeltaRhoLeadingTrackPtQ1(0x0),
+ fh3CentvsDeltaRhoLeadingTrackPtQ2(0x0),
+ fh3CentvsDeltaRhoLeadingTrackPtQ3(0x0),
+ fh3CentvsDeltaRhoLeadingTrackPtQ4(0x0),
fh2PtGenPtSmeared(0x0),
fp1Efficiency(0x0),
fp1PtResolution(0x0),
fUseAODTrackInput(kFALSE),
fUseAODMCInput(kFALSE),
fUseBackgroundCalc(kFALSE),
- fEventSelection(kFALSE), fFilterMask(0),
+ fEventSelection(kFALSE),
+ fRequireVZEROAC(kFALSE),
+ fRequireTZEROvtx(kFALSE),
+ fUseHFcuts(kFALSE),
+ fFilterMask(0),
fFilterMaskBestPt(0),
fFilterType(0),
fJetTypes(kJet),
fAvgTrials(1),
fExternalWeight(1),
fTrackEtaWindow(0.9),
+ fRequireITSRefit(0),
+ fApplySharedClusterCut(0),
fRecEtaWindow(0.5),
fTrackPtCut(0.),
fJetOutputMinPt(0.150),
fVtxR2Cut(1),
fCentCutUp(0),
fCentCutLo(0),
+ fStoreRhoLeadingTrackCorr(kFALSE),
fNonStdBranch(""),
fBackgroundBranch(""),
fNonStdFile(""),
fhEffH3(0x0),
fUseTrPtResolutionSmearing(kFALSE),
fUseDiceEfficiency(kFALSE),
- fDiceEfficiencyMinPt(0.),
+ fDiceEfficiencyMinPt(-1.),
fUseTrPtResolutionFromOADB(kFALSE),
fUseTrEfficiencyFromOADB(kFALSE),
fPathTrPtResolution(""),
fPathTrEfficiency(""),
fChangeEfficiencyFraction(0.),
+ fEfficiencyFixed(1.),
fRparam(1.0),
fAlgorithm(fastjet::kt_algorithm),
fStrategy(fastjet::Best),
fh2PtNchNRan(0x0),
fh2TracksLeadingJetPhiPtRan(0x0),
fh2TracksLeadingJetPhiPtWRan(0x0),
+ fh3CentvsRhoLeadingTrackPt(0x0),
+ fh3CentvsSigmaLeadingTrackPt(0x0),
+ fh3MultvsRhoLeadingTrackPt(0x0),
+ fh3MultvsSigmaLeadingTrackPt(0x0),
+ fh3CentvsRhoLeadingTrackPtQ1(0x0),
+ fh3CentvsRhoLeadingTrackPtQ2(0x0),
+ fh3CentvsRhoLeadingTrackPtQ3(0x0),
+ fh3CentvsRhoLeadingTrackPtQ4(0x0),
+ fh3CentvsSigmaLeadingTrackPtQ1(0x0),
+ fh3CentvsSigmaLeadingTrackPtQ2(0x0),
+ fh3CentvsSigmaLeadingTrackPtQ3(0x0),
+ fh3CentvsSigmaLeadingTrackPtQ4(0x0),
+ fh3MultvsRhoLeadingTrackPtQ1(0x0),
+ fh3MultvsRhoLeadingTrackPtQ2(0x0),
+ fh3MultvsRhoLeadingTrackPtQ3(0x0),
+ fh3MultvsRhoLeadingTrackPtQ4(0x0),
+ fh3MultvsSigmaLeadingTrackPtQ1(0x0),
+ fh3MultvsSigmaLeadingTrackPtQ2(0x0),
+ fh3MultvsSigmaLeadingTrackPtQ3(0x0),
+ fh3MultvsSigmaLeadingTrackPtQ4(0x0),
+ fh3CentvsDeltaRhoLeadingTrackPtQ1(0x0),
+ fh3CentvsDeltaRhoLeadingTrackPtQ2(0x0),
+ fh3CentvsDeltaRhoLeadingTrackPtQ3(0x0),
+ fh3CentvsDeltaRhoLeadingTrackPtQ4(0x0),
fh2PtGenPtSmeared(0x0),
fp1Efficiency(0x0),
fp1PtResolution(0x0),
fh2TracksLeadingJetPhiPtC[i] = 0;
fh2TracksLeadingJetPhiPtWC[i] = 0;
}
- DefineOutput(1, TList::Class());
+ DefineOutput(1, TList::Class());
}
if(fJetTypes&kJetRan){
fTCAJetsOutRan = new TClonesArray("AliAODJet", 0);
fTCAJetsOutRan->SetName(Form("%s_%s",fNonStdBranch.Data(),"random"));
- if(fUseDiceEfficiency || fUseTrPtResolutionSmearing)
- fTCAJetsOutRan->SetName(Form("%s_%sDetector%d%dFr%d",fNonStdBranch.Data(),"random",fUseTrPtResolutionSmearing,fUseDiceEfficiency,(int)(fChangeEfficiencyFraction*100.)));
+ if(fUseDiceEfficiency || fUseTrPtResolutionSmearing) {
+ if( fEfficiencyFixed < 1.)
+ fTCAJetsOutRan->SetName(Form("%s_%sDetector%d%dEffFixed%d",fNonStdBranch.Data(),"random",fUseTrPtResolutionSmearing,fUseDiceEfficiency,(int)(fEfficiencyFixed*100.)));
+ else
+ fTCAJetsOutRan->SetName(Form("%s_%sDetector%d%dFr%d",fNonStdBranch.Data(),"random",fUseTrPtResolutionSmearing,fUseDiceEfficiency,(int)(fChangeEfficiencyFraction*100.)));
+ }
AddAODBranch("TClonesArray",&fTCAJetsOutRan,fNonStdFile.Data());
}
if(!AODEvent()->FindListObject(Form("%s_%s",AliAODJetEventBackground::StdBranchName(),fNonStdBranch.Data()))){
fAODJetBackgroundOut = new AliAODJetEventBackground();
fAODJetBackgroundOut->SetName(Form("%s_%s",AliAODJetEventBackground::StdBranchName(),fNonStdBranch.Data()));
- if(fUseDiceEfficiency || fUseTrPtResolutionSmearing)
- fAODJetBackgroundOut->SetName(Form("%s_%sDetector%d%dFr%d",AliAODJetEventBackground::StdBranchName(),fNonStdBranch.Data(),fUseTrPtResolutionSmearing,fUseDiceEfficiency,(int)(fChangeEfficiencyFraction*100.)));
-
+ if(fUseDiceEfficiency || fUseTrPtResolutionSmearing) {
+ if( fEfficiencyFixed < 1.)
+ fAODJetBackgroundOut->SetName(Form("%s_%sDetector%d%dEffFixed%d",AliAODJetEventBackground::StdBranchName(),fNonStdBranch.Data(),fUseTrPtResolutionSmearing,fUseDiceEfficiency,(int)(fEfficiencyFixed*100.)));
+ else
+ fAODJetBackgroundOut->SetName(Form("%s_%sDetector%d%dFr%d",AliAODJetEventBackground::StdBranchName(),fNonStdBranch.Data(),fUseTrPtResolutionSmearing,fUseDiceEfficiency,(int)(fChangeEfficiencyFraction*100.)));
+ }
AddAODBranch("AliAODJetEventBackground",&fAODJetBackgroundOut,fNonStdFile.Data());
}
}
// create the branch for the random cones with the same R
TString cName = Form("%sRandomConeSkip%02d",fNonStdBranch.Data(),fNSkipLeadingCone);
- if(fUseDiceEfficiency || fUseTrPtResolutionSmearing)
- cName = Form("%sDetector%d%dFr%d_RandomConeSkip%02d",fNonStdBranch.Data(),fUseTrPtResolutionSmearing,fUseDiceEfficiency,(int)(fChangeEfficiencyFraction*100.),fNSkipLeadingCone);
-
+ if(fUseDiceEfficiency || fUseTrPtResolutionSmearing) {
+ if( fEfficiencyFixed < 1.)
+ cName = Form("%sDetector%d%dEffFixed%d_RandomConeSkip%02d",fNonStdBranch.Data(),fUseTrPtResolutionSmearing,fUseDiceEfficiency,(int)(fEfficiencyFixed*100.),fNSkipLeadingCone);
+ else
+ cName = Form("%sDetector%d%dFr%d_RandomConeSkip%02d",fNonStdBranch.Data(),fUseTrPtResolutionSmearing,fUseDiceEfficiency,(int)(fChangeEfficiencyFraction*100.),fNSkipLeadingCone);
+ }
if(fNRandomCones>0){
if(fJetTypes&kRC){
if(!AODEvent()->FindListObject(cName.Data())){
nBinEta,binLimitsEta,nBinPt,binLimitsPt);
fh2TrackEtaPt = new TH2F("fh2TrackEtaPt","pt vs eta all jets;#eta;p_{T}",
- nBinEta,binLimitsEta,nBinPt,binLimitsPt);
+ nBinEta,binLimitsEta,nBinPt,binLimitsPt);
fh2LeadingTrackEtaPt = new TH2F("fh2LeadingTrackEtaPt","pT vs eta leading jets;#eta;p_{T}",
- nBinEta,binLimitsEta,nBinPt,binLimitsPt);
+ nBinEta,binLimitsEta,nBinPt,binLimitsPt);
fh2JetsLeadingPhiEta = new TH2F("fh2JetsLeadingPhiEta","delta eta vs delta phi to leading jet;#Delta#phi;#Delta#eta",
- nBinPhi,binLimitsPhi,nBinEta,binLimitsEta);
+ nBinPhi,binLimitsPhi,nBinEta,binLimitsEta);
fh2JetsLeadingPhiPt = new TH2F("fh2JetsLeadingPhiPt","leading p_T vs delta phi to leading jet;#Delta#phi;p_{T} (GeV/c)",
- nBinPhi,binLimitsPhi,nBinPt,binLimitsPt);
+ nBinPhi,binLimitsPhi,nBinPt,binLimitsPt);
fh2TracksLeadingPhiEta = new TH2F("fh2TracksLeadingPhiEta","delta eta vs delta phi to leading track;#Delta#phi;#Delta#eta",
nBinPhi,binLimitsPhi,nBinEta,binLimitsEta);
fh2TracksLeadingPhiPt = new TH2F("fh2TracksLeadingPhiPt","leading p_T vs delta phi to leading jet;#Delta#phi;p_{T} (GeV/c)",
- nBinPhi,binLimitsPhi,nBinPt,binLimitsPt);
+ nBinPhi,binLimitsPhi,nBinPt,binLimitsPt);
fh2TracksLeadingJetPhiPt = new TH2F("fh2TracksLeadingJetPhiPt","leading p_T vs delta phi to leading jet;#Delta#phi;p_{T} (GeV/c)",
- nBinPhi,binLimitsPhi,nBinPt,binLimitsPt);
+ nBinPhi,binLimitsPhi,nBinPt,binLimitsPt);
fh2TracksLeadingJetPhiPtRan = new TH2F("fh2TracksLeadingJetPhiPtRan","leading p_T vs delta phi to leading jet;#Delta#phi;p_{T} (GeV/c)",
- nBinPhi,binLimitsPhi,nBinPt,binLimitsPt);
+ nBinPhi,binLimitsPhi,nBinPt,binLimitsPt);
fh2JetsLeadingPhiPtW = new TH2F("fh2JetsLeadingPhiPtW","leading p_T vs delta phi p_T weigted to leading jet;#Delta#phi;p_{T} (GeV/c)",
- nBinPhi,binLimitsPhi,nBinPt,binLimitsPt);
+ nBinPhi,binLimitsPhi,nBinPt,binLimitsPt);
fh2TracksLeadingPhiPtW = new TH2F("fh2TracksLeadingPhiPtW","leading p_T vs delta phi to leading jet (p_T weighted);#Delta#phi;p_{T} (GeV/c)",
- nBinPhi,binLimitsPhi,nBinPt,binLimitsPt);
+ nBinPhi,binLimitsPhi,nBinPt,binLimitsPt);
fh2TracksLeadingJetPhiPtW = new TH2F("fh2TracksLeadingJetPhiPtW","leading p_T vs delta phi to leading jet;#Delta#phi;p_{T} (GeV/c)",
nBinPhi,binLimitsPhi,nBinPt,binLimitsPt);
fh2TracksLeadingJetPhiPtWRan = new TH2F("fh2TracksLeadingJetPhiPtWRan","leading p_T vs delta phi to leading jet;#Delta#phi;p_{T} (GeV/c)",
- nBinPhi,binLimitsPhi,nBinPt,binLimitsPt);
+ nBinPhi,binLimitsPhi,nBinPt,binLimitsPt);
+
+ if(fStoreRhoLeadingTrackCorr) {
+ fh3CentvsRhoLeadingTrackPt = new TH3F("fh3CentvsRhoLeadingTrackPt","centrality vs background density full event; centrality; #rho", 50,0.,100.,500,0.,250.,100,0.,100.);
+ fh3CentvsSigmaLeadingTrackPt = new TH3F("fh3CentvsSigmaLeadingTrackPt","centrality vs sigma full event; centrality; #sigma(#rho)", 50,0.,100.,50,0.,50.,100,0.,100.);
+ fh3MultvsRhoLeadingTrackPt = new TH3F("fh3MultvsRhoLeadingTrackPt","multiplicity vs background density full event; multiplicity; #rho", 100,0.,5000.,500,0.,250.,100,0.,100.);
+ fh3MultvsSigmaLeadingTrackPt = new TH3F("fh3MultvsSigmaLeadingTrackPt","multiplicity vs sigma full event; multiplicity; #sigma(#rho)", 100,0.,5000.,50,0.,50.,100,0.,100.);
+
+
+ fh3CentvsRhoLeadingTrackPtQ1 = new TH3F("fh3CentvsRhoLeadingTrackPtQ1","centrality vs background density Q1; centrality; #rho; leading p_{t}^{track}", 50,0.,100.,500,0.,250.,100,0.,100.);
+ fh3CentvsRhoLeadingTrackPtQ2 = new TH3F("fh3CentvsRhoLeadingTrackPtQ2","centrality vs background density Q2; centrality; #rho; leading p_{t}^{track}", 50,0.,100.,500,0.,250.,100,0.,100.);
+ fh3CentvsRhoLeadingTrackPtQ3 = new TH3F("fh3CentvsRhoLeadingTrackPtQ3","centrality vs background density Q3; centrality; #rho; leading p_{t}^{track}", 50,0.,100.,500,0.,250.,100,0.,100.);
+ fh3CentvsRhoLeadingTrackPtQ4 = new TH3F("fh3CentvsRhoLeadingTrackPtQ4","centrality vs background density Q4; centrality; #rho; leading p_{t}^{track}", 50,0.,100.,500,0.,250.,100,0.,100.);
+
+ fh3CentvsSigmaLeadingTrackPtQ1 = new TH3F("fh3CentvsSigmaLeadingTrackPtQ1","centrality vs background density Q1; centrality; #sigma(#rho); leading p_{t}^{track}", 50,0.,100.,50,0.,50.,100,0.,100.);
+ fh3CentvsSigmaLeadingTrackPtQ2 = new TH3F("fh3CentvsSigmaLeadingTrackPtQ2","centrality vs background density Q2; centrality; #sigma(#rho); leading p_{t}^{track}", 50,0.,100.,50,0.,50.,100,0.,100.);
+ fh3CentvsSigmaLeadingTrackPtQ3 = new TH3F("fh3CentvsSigmaLeadingTrackPtQ3","centrality vs background density Q3; centrality; #sigma(#rho); leading p_{t}^{track}", 50,0.,100.,50,0.,50.,100,0.,100.);
+ fh3CentvsSigmaLeadingTrackPtQ4 = new TH3F("fh3CentvsSigmaLeadingTrackPtQ4","centrality vs background density Q4; centrality; #sigma(#rho); leading p_{t}^{track}", 50,0.,100.,50,0.,50.,100,0.,100.);
+
+ fh3MultvsRhoLeadingTrackPtQ1 = new TH3F("fh3MultvsRhoLeadingTrackPtQ1","multiplicity vs background density Q1; multiplicity; #rho; leading p_{t}^{track}", 100,0.,5000.,500,0.,250.,100,0.,100.);
+ fh3MultvsRhoLeadingTrackPtQ2 = new TH3F("fh3MultvsRhoLeadingTrackPtQ2","multiplicity vs background density Q2; multiplicity; #rho; leading p_{t}^{track}", 100,0.,5000.,500,0.,250.,100,0.,100.);
+ fh3MultvsRhoLeadingTrackPtQ3 = new TH3F("fh3MultvsRhoLeadingTrackPtQ3","multiplicity vs background density Q3; multiplicity; #rho; leading p_{t}^{track}", 100,0.,5000.,500,0.,250.,100,0.,100.);
+ fh3MultvsRhoLeadingTrackPtQ4 = new TH3F("fh3MultvsRhoLeadingTrackPtQ4","multiplicity vs background density Q4; multiplicity; #rho; leading p_{t}^{track}", 100,0.,5000.,500,0.,250.,100,0.,100.);
+
+ fh3MultvsSigmaLeadingTrackPtQ1 = new TH3F("fh3MultvsSigmaLeadingTrackPtQ1","multiplicity vs background density Q1; multiplicity; #sigma(#rho); leading p_{t}^{track}", 100,0.,5000.,50,0.,50.,100,0.,100.);
+ fh3MultvsSigmaLeadingTrackPtQ2 = new TH3F("fh3MultvsSigmaLeadingTrackPtQ2","multiplicity vs background density Q2; multiplicity; #sigma(#rho); leading p_{t}^{track}", 100,0.,5000.,50,0.,50.,100,0.,100.);
+ fh3MultvsSigmaLeadingTrackPtQ3 = new TH3F("fh3MultvsSigmaLeadingTrackPtQ3","multiplicity vs background density Q3; multiplicity; #sigma(#rho); leading p_{t}^{track}", 100,0.,5000.,50,0.,50.,100,0.,100.);
+ fh3MultvsSigmaLeadingTrackPtQ4 = new TH3F("fh3MultvsSigmaLeadingTrackPtQ4","multiplicity vs background density Q4; multiplicity; #sigma(#rho); leading p_{t}^{track}", 100,0.,5000.,50,0.,50.,100,0.,100.);
+
+
+ fh3CentvsDeltaRhoLeadingTrackPtQ1 = new TH3F("fh3CentvsDeltaRhoLeadingTrackPtQ1","centrality vs background density Q1; centrality; #rho_{quadrant}-#rho_{full event}; leading p_{t}^{track}", 50,0.,100.,200,-10.,10.,100,0.,100.);
+ fh3CentvsDeltaRhoLeadingTrackPtQ2 = new TH3F("fh3CentvsDeltaRhoLeadingTrackPtQ2","centrality vs background density Q2; centrality; #rho_{quadrant}-#rho_{full event}; leading p_{t}^{track}", 50,0.,100.,200,-10.,10.,100,0.,100.);
+ fh3CentvsDeltaRhoLeadingTrackPtQ3 = new TH3F("fh3CentvsDeltaRhoLeadingTrackPtQ3","centrality vs background density Q3; centrality; #rho_{quadrant}-#rho_{full event}; leading p_{t}^{track}", 50,0.,100.,200,-10.,10.,100,0.,100.);
+ fh3CentvsDeltaRhoLeadingTrackPtQ4 = new TH3F("fh3CentvsDeltaRhoLeadingTrackPtQ4","centrality vs background density Q4; centrality; #rho_{quadrant}-#rho_{full event}; leading p_{t}^{track}", 50,0.,100.,200,-10.,10.,100,0.,100.);
+
+ fHistList->Add(fh3CentvsRhoLeadingTrackPt);
+ fHistList->Add(fh3CentvsSigmaLeadingTrackPt);
+ fHistList->Add(fh3MultvsRhoLeadingTrackPt);
+ fHistList->Add(fh3MultvsSigmaLeadingTrackPt);
+
+ fHistList->Add(fh3CentvsRhoLeadingTrackPtQ1);
+ fHistList->Add(fh3CentvsRhoLeadingTrackPtQ2);
+ fHistList->Add(fh3CentvsRhoLeadingTrackPtQ3);
+ fHistList->Add(fh3CentvsRhoLeadingTrackPtQ4);
+
+ fHistList->Add(fh3CentvsSigmaLeadingTrackPtQ1);
+ fHistList->Add(fh3CentvsSigmaLeadingTrackPtQ2);
+ fHistList->Add(fh3CentvsSigmaLeadingTrackPtQ3);
+ fHistList->Add(fh3CentvsSigmaLeadingTrackPtQ4);
+
+ fHistList->Add(fh3MultvsRhoLeadingTrackPtQ1);
+ fHistList->Add(fh3MultvsRhoLeadingTrackPtQ2);
+ fHistList->Add(fh3MultvsRhoLeadingTrackPtQ3);
+ fHistList->Add(fh3MultvsRhoLeadingTrackPtQ4);
+
+ fHistList->Add(fh3MultvsSigmaLeadingTrackPtQ1);
+ fHistList->Add(fh3MultvsSigmaLeadingTrackPtQ2);
+ fHistList->Add(fh3MultvsSigmaLeadingTrackPtQ3);
+ fHistList->Add(fh3MultvsSigmaLeadingTrackPtQ4);
+
+ fHistList->Add(fh3CentvsDeltaRhoLeadingTrackPtQ1);
+ fHistList->Add(fh3CentvsDeltaRhoLeadingTrackPtQ2);
+ fHistList->Add(fh3CentvsDeltaRhoLeadingTrackPtQ3);
+ fHistList->Add(fh3CentvsDeltaRhoLeadingTrackPtQ4);
+ }
+
//Detector level effects histos
fh2PtGenPtSmeared = new TH2F("fh2PtGenPtSmeared","fh2PtGenPtSmeared",nBinPt,binLimitsPt,nBinPt,binLimitsPt);
fHistList->Add(fh1BiARandomConesRan[i]);
}
}
- for(int i = 0;i < kMaxCent;i++){
- fHistList->Add(fh2JetsLeadingPhiPtC[i]);
- fHistList->Add(fh2JetsLeadingPhiPtWC[i]);
- fHistList->Add(fh2TracksLeadingJetPhiPtC[i]);
- fHistList->Add(fh2TracksLeadingJetPhiPtWC[i]);
- }
+ for(int i = 0;i < kMaxCent;i++){
+ fHistList->Add(fh2JetsLeadingPhiPtC[i]);
+ fHistList->Add(fh2JetsLeadingPhiPtWC[i]);
+ fHistList->Add(fh2TracksLeadingJetPhiPtC[i]);
+ fHistList->Add(fh2TracksLeadingJetPhiPtWC[i]);
+ }
fHistList->Add(fh2NRecJetsPt);
fHistList->Add(fh2NRecTracksPt);
//Check if information is provided detector level effects
if(!fMomResH1 || !fMomResH2 || !fMomResH3) fUseTrPtResolutionSmearing = kFALSE;
- if(!fhEffH1 || !fhEffH2 || !fhEffH3) fUseDiceEfficiency = kFALSE;
-
+ if( fEfficiencyFixed < 1. ) {
+ if (!fUseDiceEfficiency)
+ fUseDiceEfficiency = 1; // 1 is the default; 2 can be set by user
+ }
+ else {
+ if(!fhEffH1 || !fhEffH2 || !fhEffH3 ) fUseDiceEfficiency = kFALSE;
+ }
+
Bool_t selectEvent = false;
Bool_t physicsSelection = true;// handled by the framework(fInputHandler->IsEventSelected()&AliVEvent::kMB)==AliVEvent::kMB;
selectEvent = true;
}
}
+
+
+ Bool_t T0 = false;
+ Bool_t V0 = false;
+ const AliAODVZERO *vzero = fAOD->GetVZEROData();
+ if(vzero){
+ if((vzero->GetTriggerChargeA()>0)&&(vzero->GetTriggerChargeC()>0)){
+ V0 = true;
+ }
+ }
+
+ const AliAODTZERO *tzero = fAOD->GetTZEROData();
+ if(tzero){
+ if(TMath::Abs(tzero->GetT0VertexRaw())<100){
+ T0 = true;
+ }
+ }
+ if(fRequireVZEROAC&&fRequireTZEROvtx)selectEvent = selectEvent&&V0&&T0;
+ else if(fRequireTZEROvtx)selectEvent = selectEvent&&T0;
+ else if(fRequireVZEROAC)selectEvent = selectEvent&&V0;
+
if(!selectEvent){
PostData(1, fHistList);
Int_t nT = GetListOfTracks(&recParticles,fTrackTypeRec);
Float_t nCh = recParticles.GetEntries();
+ Float_t nGen=genParticles.GetEntries();
fh1Nch->Fill(nCh);
if(fDebug>2)Printf("%s:%d Selected Rec tracks: %d %d",(char*)__FILE__,__LINE__,nT,recParticles.GetEntries());
nT = GetListOfTracks(&genParticles,fTrackTypeGen);
// Dice to decide if particle is kept or not - toy model for efficiency
//
- Double_t rnd = fRandom->Uniform(1.);
- Double_t pT = vp->Pt();
+ Double_t sumEff = 0.;
+ Double_t pT = 0.;
Double_t eff[3] = {0.};
- Double_t pTtmp = pT;
- if(pT>10.) pTtmp = 10.;
- Double_t eff1 = fhEffH1->GetBinContent(fhEffH1->FindBin(pTtmp));
- Double_t eff2 = fhEffH2->GetBinContent(fhEffH2->FindBin(pTtmp));
- Double_t eff3 = fhEffH3->GetBinContent(fhEffH3->FindBin(pTtmp));
Int_t cat[3] = {0};
- //Sort efficiencies from large to small
- if(eff1>eff2 && eff1>eff3) {
- eff[0] = eff1;
- cat[0] = 1;
- if(eff2>eff3) {
- eff[1] = eff2;
- eff[2] = eff3;
- cat[1] = 2;
- cat[2] = 3;
- } else {
- eff[1] = eff3;
- eff[2] = eff2;
- cat[1] = 3;
- cat[2] = 2;
+
+ Double_t rnd = fRandom->Uniform(1.);
+ if( fEfficiencyFixed<1. ) {
+ sumEff = fEfficiencyFixed;
+ if (fUseDiceEfficiency == 2) {
+ sumEff = (nCh - fEfficiencyFixed*nGen) / nCh; // rescale eff; fEfficiencyFixed is wrt to nGen, but dicing is fraction of nCh
+ }
+ } else {
+
+ pT = vp->Pt();
+ Double_t pTtmp = pT;
+ if(pT>10.) pTtmp = 10.;
+ Double_t eff1 = fhEffH1->GetBinContent(fhEffH1->FindBin(pTtmp));
+ Double_t eff2 = fhEffH2->GetBinContent(fhEffH2->FindBin(pTtmp));
+ Double_t eff3 = fhEffH3->GetBinContent(fhEffH3->FindBin(pTtmp));
+
+ //Sort efficiencies from large to small
+ if(eff1>eff2 && eff1>eff3) {
+ eff[0] = eff1;
+ cat[0] = 1;
+ if(eff2>eff3) {
+ eff[1] = eff2;
+ eff[2] = eff3;
+ cat[1] = 2;
+ cat[2] = 3;
+ } else {
+ eff[1] = eff3;
+ eff[2] = eff2;
+ cat[1] = 3;
+ cat[2] = 2;
+ }
}
- }
- else if(eff2>eff1 && eff2>eff3) {
- eff[0] = eff2;
- cat[0] = 2;
- if(eff1>eff3) {
- eff[1] = eff1;
- eff[2] = eff3;
- cat[1] = 1;
- cat[2] = 3;
- } else {
- eff[1] = eff3;
- eff[2] = eff1;
- cat[1] = 3;
- cat[2] = 1;
+ else if(eff2>eff1 && eff2>eff3) {
+ eff[0] = eff2;
+ cat[0] = 2;
+ if(eff1>eff3) {
+ eff[1] = eff1;
+ eff[2] = eff3;
+ cat[1] = 1;
+ cat[2] = 3;
+ } else {
+ eff[1] = eff3;
+ eff[2] = eff1;
+ cat[1] = 3;
+ cat[2] = 1;
+ }
}
- }
- else if(eff3>eff1 && eff3>eff2) {
- eff[0] = eff3;
- cat[0] = 3;
- if(eff1>eff2) {
- eff[1] = eff1;
- eff[2] = eff2;
- cat[1] = 1;
- cat[2] = 2;
- } else {
- eff[1] = eff2;
- eff[2] = eff1;
- cat[1] = 2;
- cat[2] = 1;
+ else if(eff3>eff1 && eff3>eff2) {
+ eff[0] = eff3;
+ cat[0] = 3;
+ if(eff1>eff2) {
+ eff[1] = eff1;
+ eff[2] = eff2;
+ cat[1] = 1;
+ cat[2] = 2;
+ } else {
+ eff[1] = eff2;
+ eff[2] = eff1;
+ cat[1] = 2;
+ cat[2] = 1;
+ }
}
+
+ sumEff = eff[0]+eff[1]+eff[2];
}
-
- Double_t sumEff = eff[0]+eff[1]+eff[2];
fp1Efficiency->Fill(vp->Pt(),sumEff);
if(rnd>sumEff && pT > fDiceEfficiencyMinPt) continue;
// now create the object that holds info about ghosts
/*
- if(!fUseBackgroundCalc&& fNonStdBranch.Length()==0){
+ if(!fUseBackgroundCalc&& fNonStdBranch.Length()==0){
// reduce CPU time...
fGhostArea = 0.5;
fActiveAreaRepeats = 0;
- }
+ }
*/
fastjet::GhostedAreaSpec ghostSpec(fGhostEtamax, fActiveAreaRepeats, fGhostArea);
fh1NJetsRec->Fill(sortedJets.size());
- // loop over all jets an fill information, first on is the leading jet
+ // loop over all jets an fill information, first on is the leading jet
Int_t nRecOver = inclusiveJets.size();
Int_t nRec = inclusiveJets.size();
if(externalBackground){
// carefull has to be filled in a task before
// todo, ReArrange to the botom
- pTback = externalBackground->GetBackground(2)*leadingJet.EffectiveAreaCharged();
+ pTback = externalBackground->GetBackground(2)*leadingJet.EffectiveAreaCharged();
}
pt = leadingJet.Pt() - pTback;
// correlation of leading jet with tracks
Float_t tmpPtBack = 0;
if(externalBackground){
// carefull has to be filled in a task before
- // todo, ReArrange to the botom
+ // todo, ReArrange to the botom
tmpPtBack = externalBackground->GetBackground(2)*tmpRec.EffectiveAreaCharged();
}
tmpPt = tmpPt - tmpPtBack;
}
AliAODTrack *aodT = 0;
- if(partLead&&aodOutJet){
- //set pT of leading constituent of jet
- aodOutJet->SetPtLeading(partLead->Pt());
- aodT = dynamic_cast<AliAODTrack*>(partLead);
- if(aodT){
- if(aodT->TestFilterBit(fFilterMaskBestPt)){
- aodOutJet->SetTrigger(AliAODJet::kHighTrackPtBest);
+ if(partLead){
+ if(aodOutJet){
+ //set pT of leading constituent of jet
+ aodOutJet->SetPtLeading(partLead->Pt());
+ aodT = dynamic_cast<AliAODTrack*>(partLead);
+ if(aodT){
+ if(aodT->TestFilterBit(fFilterMaskBestPt)){
+ aodOutJet->SetTrigger(AliAODJet::kHighTrackPtBest);
+ }
}
}
}
-
- // correlation
- Float_t tmpPhi = tmpRec.Phi();
- Float_t tmpEta = tmpRec.Eta();
- if(tmpPhi<0)tmpPhi+=TMath::Pi()*2.;
- if(j==0){
- fh1PtJetsLeadingRecIn->Fill(tmpPt);
- fh2LeadingJetPhiEta->Fill(tmpPhi,tmpEta);
- fh2LeadingJetEtaPt->Fill(tmpEta,tmpPt);
- fh1NConstLeadingRec->Fill(constituents.size());
- fh2NConstLeadingPt->Fill(tmpPt,constituents.size());
- continue;
- }
- fh2JetPhiEta->Fill(tmpRec.Phi(),tmpEta);
- fh2JetEtaPt->Fill(tmpEta,tmpPt);
- Float_t dPhi = phi - tmpPhi;
- if(dPhi>TMath::Pi())dPhi = dPhi - 2.*TMath::Pi();
- if(dPhi<(-1.*TMath::Pi()))dPhi = dPhi + 2.*TMath::Pi();
- Float_t dEta = eta - tmpRec.Eta();
- fh2JetsLeadingPhiEta->Fill(dPhi,dEta);
- fh2JetsLeadingPhiPt->Fill(dPhi,pt);
- if(cenClass>=0){
- fh2JetsLeadingPhiPtC[cenClass]->Fill(dPhi,pt);
- fh2JetsLeadingPhiPtWC[cenClass]->Fill(dPhi,pt,tmpPt);
- }
- fh2JetsLeadingPhiPtW->Fill(dPhi,pt,tmpPt);
+
+ // correlation
+ Float_t tmpPhi = tmpRec.Phi();
+ Float_t tmpEta = tmpRec.Eta();
+ if(tmpPhi<0)tmpPhi+=TMath::Pi()*2.;
+ if(j==0){
+ fh1PtJetsLeadingRecIn->Fill(tmpPt);
+ fh2LeadingJetPhiEta->Fill(tmpPhi,tmpEta);
+ fh2LeadingJetEtaPt->Fill(tmpEta,tmpPt);
+ fh1NConstLeadingRec->Fill(constituents.size());
+ fh2NConstLeadingPt->Fill(tmpPt,constituents.size());
+ continue;
+ }
+ fh2JetPhiEta->Fill(tmpRec.Phi(),tmpEta);
+ fh2JetEtaPt->Fill(tmpEta,tmpPt);
+ Float_t dPhi = phi - tmpPhi;
+ if(dPhi>TMath::Pi())dPhi = dPhi - 2.*TMath::Pi();
+ if(dPhi<(-1.*TMath::Pi()))dPhi = dPhi + 2.*TMath::Pi();
+ Float_t dEta = eta - tmpRec.Eta();
+ fh2JetsLeadingPhiEta->Fill(dPhi,dEta);
+ fh2JetsLeadingPhiPt->Fill(dPhi,pt);
+ if(cenClass>=0){
+ fh2JetsLeadingPhiPtC[cenClass]->Fill(dPhi,pt);
+ fh2JetsLeadingPhiPtWC[cenClass]->Fill(dPhi,pt,tmpPt);
+ }
+ fh2JetsLeadingPhiPtW->Fill(dPhi,pt,tmpPt);
}// loop over reconstructed jets
- delete recIter;
-
-
-
- // Add the random cones...
- if(fNRandomCones>0&&fTCARandomConesOut){
- // create a random jet within the acceptance
- Double_t etaMax = fTrackEtaWindow - fRparam;
- Int_t nCone = 0;
- Int_t nConeRan = 0;
- Double_t pTC = 1; // small number
- for(int ir = 0;ir < fNRandomCones;ir++){
- Double_t etaC = etaMax*2.*(fRandom->Rndm()-0.5); // +- etamax
- Double_t phiC = fRandom->Rndm()*2.*TMath::Pi(); // 0 - 2pi
- // massless jet
- Double_t thetaC = 2.*TMath::ATan(TMath::Exp(-etaC));
- Double_t pZC = pTC/TMath::Tan(thetaC);
- Double_t pXC = pTC * TMath::Cos(phiC);
- Double_t pYC = pTC * TMath::Sin(phiC);
- Double_t pC = TMath::Sqrt(pTC*pTC+pZC*pZC);
- AliAODJet tmpRecC (pXC,pYC,pZC, pC);
- bool skip = false;
- for(int jj = 0; jj < TMath::Min(nRec,fNSkipLeadingCone);jj++){// test for overlap with leading jets
- AliAODJet jet (sortedJets[jj].px(), sortedJets[jj].py(), sortedJets[jj].pz(), sortedJets[jj].E());
- if(jet.DeltaR(& tmpRecC)<2.*fRparam+0.2){
- skip = true;
- break;
- }
- }
- // test for overlap with previous cones to avoid double counting
- for(int iic = 0;iic<ir;iic++){
- AliAODJet *iicone = (AliAODJet*)fTCARandomConesOut->At(iic);
- if(iicone){
- if(iicone->DeltaR(&tmpRecC)<2.*fRparam){
- skip = true;
- break;
- }
- }
- }
- if(skip)continue;
- tmpRecC.SetBgEnergy(0,0); // this is use as temporary storage of the summed p_T below
- tmpRecC.SetPtLeading(-1.);
- if(fTCARandomConesOut)new ((*fTCARandomConesOut)[nCone++]) AliAODJet(tmpRecC);
- if(fTCARandomConesOutRan)new ((*fTCARandomConesOutRan)[nConeRan++]) AliAODJet(tmpRecC);
- }// loop over random cones creation
+ delete recIter;
+
+
+
+ // Add the random cones...
+ if(fNRandomCones>0&&fTCARandomConesOut){
+ // create a random jet within the acceptance
+ Double_t etaMax = fTrackEtaWindow - fRparam;
+ Int_t nCone = 0;
+ Int_t nConeRan = 0;
+ Double_t pTC = 1; // small number
+ for(int ir = 0;ir < fNRandomCones;ir++){
+ Double_t etaC = etaMax*2.*(fRandom->Rndm()-0.5); // +- etamax
+ Double_t phiC = fRandom->Rndm()*2.*TMath::Pi(); // 0 - 2pi
+ // massless jet
+ Double_t thetaC = 2.*TMath::ATan(TMath::Exp(-etaC));
+ Double_t pZC = pTC/TMath::Tan(thetaC);
+ Double_t pXC = pTC * TMath::Cos(phiC);
+ Double_t pYC = pTC * TMath::Sin(phiC);
+ Double_t pC = TMath::Sqrt(pTC*pTC+pZC*pZC);
+ AliAODJet tmpRecC (pXC,pYC,pZC, pC);
+ bool skip = false;
+ for(int jj = 0; jj < TMath::Min(nRec,fNSkipLeadingCone);jj++){// test for overlap with leading jets
+ AliAODJet jet (sortedJets[jj].px(), sortedJets[jj].py(), sortedJets[jj].pz(), sortedJets[jj].E());
+ if(jet.DeltaR(& tmpRecC)<2.*fRparam+0.2){
+ skip = true;
+ break;
+ }
+ }
+ // test for overlap with previous cones to avoid double counting
+ for(int iic = 0;iic<ir;iic++){
+ AliAODJet *iicone = (AliAODJet*)fTCARandomConesOut->At(iic);
+ if(iicone){
+ if(iicone->DeltaR(&tmpRecC)<2.*fRparam){
+ skip = true;
+ break;
+ }
+ }
+ }
+ if(skip)continue;
+ tmpRecC.SetBgEnergy(0,0); // this is use as temporary storage of the summed p_T below
+ tmpRecC.SetPtLeading(-1.);
+ if(fTCARandomConesOut)new ((*fTCARandomConesOut)[nCone++]) AliAODJet(tmpRecC);
+ if(fTCARandomConesOutRan)new ((*fTCARandomConesOutRan)[nConeRan++]) AliAODJet(tmpRecC);
+ }// loop over random cones creation
- // loop over the reconstructed particles and add up the pT in the random cones
- // maybe better to loop over randomized particles not in the real jets...
- // but this by definition brings dow average energy in the whole event
- AliAODJet vTmpRanR(1,0,0,1);
- for(int i = 0; i < recParticles.GetEntries(); i++){
- AliVParticle *vp = (AliVParticle*)recParticles.At(i);
- if(fTCARandomConesOut){
- for(int ir = 0;ir < fNRandomCones;ir++){
- AliAODJet *jC = (AliAODJet*)fTCARandomConesOut->At(ir);
- if(jC&&jC->DeltaR(vp)<fRparam){
- if(vp->Pt()>fMaxTrackPtInJet)jC->SetTrigger(AliAODJet::kHighTrackPtTriggered);
- jC->SetBgEnergy(jC->ChargedBgEnergy()+vp->Pt(),0);
- if(vp->Pt() > jC->GetPtLeading()) jC->SetPtLeading(vp->Pt());
- }
- }
- }// add up energy in cone
-
- // the randomized input changes eta and phi, but keeps the p_T
- if(i>=fNSkipLeadingRan){// eventually skip the leading particles
- Double_t pTR = vp->Pt();
- Double_t etaR = 2.*fTrackEtaWindow* fRandom->Rndm() - fTrackEtaWindow;
- Double_t phiR = 2.* TMath::Pi() * fRandom->Rndm();
+ // loop over the reconstructed particles and add up the pT in the random cones
+ // maybe better to loop over randomized particles not in the real jets...
+ // but this by definition brings dow average energy in the whole event
+ AliAODJet vTmpRanR(1,0,0,1);
+ for(int i = 0; i < recParticles.GetEntries(); i++){
+ AliVParticle *vp = (AliVParticle*)recParticles.At(i);
+ if(fTCARandomConesOut){
+ for(int ir = 0;ir < fNRandomCones;ir++){
+ AliAODJet *jC = (AliAODJet*)fTCARandomConesOut->At(ir);
+ if(jC&&jC->DeltaR(vp)<fRparam){
+ if(vp->Pt()>fMaxTrackPtInJet)jC->SetTrigger(AliAODJet::kHighTrackPtTriggered);
+ jC->SetBgEnergy(jC->ChargedBgEnergy()+vp->Pt(),0);
+ if(vp->Pt() > jC->GetPtLeading()) jC->SetPtLeading(vp->Pt());
+ }
+ }
+ }// add up energy in cone
+
+ // the randomized input changes eta and phi, but keeps the p_T
+ if(i>=fNSkipLeadingRan){// eventually skip the leading particles
+ Double_t pTR = vp->Pt();
+ Double_t etaR = 2.*fTrackEtaWindow* fRandom->Rndm() - fTrackEtaWindow;
+ Double_t phiR = 2.* TMath::Pi() * fRandom->Rndm();
- Double_t thetaR = 2.*TMath::ATan(TMath::Exp(-etaR));
- Double_t pZR = pTR/TMath::Tan(thetaR);
+ Double_t thetaR = 2.*TMath::ATan(TMath::Exp(-etaR));
+ Double_t pZR = pTR/TMath::Tan(thetaR);
- Double_t pXR = pTR * TMath::Cos(phiR);
- Double_t pYR = pTR * TMath::Sin(phiR);
- Double_t pR = TMath::Sqrt(pTR*pTR+pZR*pZR);
- vTmpRanR.SetPxPyPzE(pXR,pYR,pZR,pR);
- if(fTCARandomConesOutRan){
- for(int ir = 0;ir < fTCARandomConesOutRan->GetEntriesFast();ir++){
- AliAODJet *jC = (AliAODJet*)fTCARandomConesOutRan->At(ir);
- if(jC&&jC->DeltaR(&vTmpRanR)<fRparam){
- if(vTmpRanR.Pt()>fMaxTrackPtInJet)jC->SetTrigger(AliAODJet::kHighTrackPtTriggered);
- jC->SetBgEnergy(jC->ChargedBgEnergy()+vTmpRanR.Pt(),0);
- if(vTmpRanR.Pt() > jC->GetPtLeading()) jC->SetPtLeading(vTmpRanR.Pt());
- }
- }
- }
- }
- }// loop over recparticles
+ Double_t pXR = pTR * TMath::Cos(phiR);
+ Double_t pYR = pTR * TMath::Sin(phiR);
+ Double_t pR = TMath::Sqrt(pTR*pTR+pZR*pZR);
+ vTmpRanR.SetPxPyPzE(pXR,pYR,pZR,pR);
+ if(fTCARandomConesOutRan){
+ for(int ir = 0;ir < fTCARandomConesOutRan->GetEntriesFast();ir++){
+ AliAODJet *jC = (AliAODJet*)fTCARandomConesOutRan->At(ir);
+ if(jC&&jC->DeltaR(&vTmpRanR)<fRparam){
+ if(vTmpRanR.Pt()>fMaxTrackPtInJet)jC->SetTrigger(AliAODJet::kHighTrackPtTriggered);
+ jC->SetBgEnergy(jC->ChargedBgEnergy()+vTmpRanR.Pt(),0);
+ if(vTmpRanR.Pt() > jC->GetPtLeading()) jC->SetPtLeading(vTmpRanR.Pt());
+ }
+ }
+ }
+ }
+ }// loop over recparticles
- Float_t jetArea = fRparam*fRparam*TMath::Pi();
- if(fTCARandomConesOut){
- for(int ir = 0;ir < fTCARandomConesOut->GetEntriesFast();ir++){
- // rescale the momentum vectors for the random cones
+ Float_t jetArea = fRparam*fRparam*TMath::Pi();
+ if(fTCARandomConesOut){
+ for(int ir = 0;ir < fTCARandomConesOut->GetEntriesFast();ir++){
+ // rescale the momentum vectors for the random cones
- AliAODJet *rC = (AliAODJet*)fTCARandomConesOut->At(ir);
- if(rC){
- Double_t etaC = rC->Eta();
- Double_t phiC = rC->Phi();
- // massless jet, unit vector
- pTC = rC->ChargedBgEnergy();
- if(pTC<=0)pTC = 0.001; // for almost empty events
- Double_t thetaC = 2.*TMath::ATan(TMath::Exp(-etaC));
- Double_t pZC = pTC/TMath::Tan(thetaC);
- Double_t pXC = pTC * TMath::Cos(phiC);
- Double_t pYC = pTC * TMath::Sin(phiC);
- Double_t pC = TMath::Sqrt(pTC*pTC+pZC*pZC);
- rC->SetPxPyPzE(pXC,pYC,pZC, pC);
- rC->SetBgEnergy(0,0);
- rC->SetEffArea(jetArea,0);
- }
- }
- }
- if(fTCARandomConesOutRan){
- for(int ir = 0;ir < fTCARandomConesOutRan->GetEntriesFast();ir++){
- AliAODJet* rC = (AliAODJet*)fTCARandomConesOutRan->At(ir);
- // same wit random
- if(rC){
- Double_t etaC = rC->Eta();
- Double_t phiC = rC->Phi();
- // massless jet, unit vector
- pTC = rC->ChargedBgEnergy();
- if(pTC<=0)pTC = 0.001;// for almost empty events
- Double_t thetaC = 2.*TMath::ATan(TMath::Exp(-etaC));
- Double_t pZC = pTC/TMath::Tan(thetaC);
- Double_t pXC = pTC * TMath::Cos(phiC);
- Double_t pYC = pTC * TMath::Sin(phiC);
- Double_t pC = TMath::Sqrt(pTC*pTC+pZC*pZC);
- rC->SetPxPyPzE(pXC,pYC,pZC, pC);
- rC->SetBgEnergy(0,0);
- rC->SetEffArea(jetArea,0);
- }
- }
- }
- }// if(fNRandomCones
+ AliAODJet *rC = (AliAODJet*)fTCARandomConesOut->At(ir);
+ if(rC){
+ Double_t etaC = rC->Eta();
+ Double_t phiC = rC->Phi();
+ // massless jet, unit vector
+ pTC = rC->ChargedBgEnergy();
+ if(pTC<=0)pTC = 0.001; // for almost empty events
+ Double_t thetaC = 2.*TMath::ATan(TMath::Exp(-etaC));
+ Double_t pZC = pTC/TMath::Tan(thetaC);
+ Double_t pXC = pTC * TMath::Cos(phiC);
+ Double_t pYC = pTC * TMath::Sin(phiC);
+ Double_t pC = TMath::Sqrt(pTC*pTC+pZC*pZC);
+ rC->SetPxPyPzE(pXC,pYC,pZC, pC);
+ rC->SetBgEnergy(0,0);
+ rC->SetEffArea(jetArea,0);
+ }
+ }
+ }
+ if(fTCARandomConesOutRan){
+ for(int ir = 0;ir < fTCARandomConesOutRan->GetEntriesFast();ir++){
+ AliAODJet* rC = (AliAODJet*)fTCARandomConesOutRan->At(ir);
+ // same wit random
+ if(rC){
+ Double_t etaC = rC->Eta();
+ Double_t phiC = rC->Phi();
+ // massless jet, unit vector
+ pTC = rC->ChargedBgEnergy();
+ if(pTC<=0)pTC = 0.001;// for almost empty events
+ Double_t thetaC = 2.*TMath::ATan(TMath::Exp(-etaC));
+ Double_t pZC = pTC/TMath::Tan(thetaC);
+ Double_t pXC = pTC * TMath::Cos(phiC);
+ Double_t pYC = pTC * TMath::Sin(phiC);
+ Double_t pC = TMath::Sqrt(pTC*pTC+pZC*pZC);
+ rC->SetPxPyPzE(pXC,pYC,pZC, pC);
+ rC->SetBgEnergy(0,0);
+ rC->SetEffArea(jetArea,0);
+ }
+ }
+ }
+ }// if(fNRandomCones
- //background estimates:all bckg jets(0) & wo the 2 hardest(1)
-
+ //background estimates:all bckg jets(0) & wo the 2 hardest(1)
+ if(fAODJetBackgroundOut){
+ vector<fastjet::PseudoJet> jets2=sortedJets;
+ if(jets2.size()>2) jets2.erase(jets2.begin(),jets2.begin()+2);
+ Double_t bkg1=0;
+ Double_t sigma1=0.;
+ Double_t meanarea1=0.;
+ Double_t bkg2=0;
+ Double_t sigma2=0.;
+ Double_t meanarea2=0.;
+ Double_t bkg4=0;
+ Double_t sigma4=0.;
+ Double_t meanarea4=0.;
- if(fAODJetBackgroundOut){
- vector<fastjet::PseudoJet> jets2=sortedJets;
- if(jets2.size()>2) jets2.erase(jets2.begin(),jets2.begin()+2);
- Double_t bkg1=0;
- Double_t sigma1=0.;
- Double_t meanarea1=0.;
- Double_t bkg2=0;
- Double_t sigma2=0.;
- Double_t meanarea2=0.;
+ clustSeq.get_median_rho_and_sigma(jets2, range, true, bkg1, sigma1, meanarea1, true);
+ fAODJetBackgroundOut->SetBackground(0,bkg1,sigma1,meanarea1);
- clustSeq.get_median_rho_and_sigma(jets2, range, true, bkg1, sigma1, meanarea1, true);
- fAODJetBackgroundOut->SetBackground(0,bkg1,sigma1,meanarea1);
+ // fh1BiARandomCones[0]->Fill(omCone-(bkg1*areaRandomCone));
+ // fh1BiARandomConesRan[0]->Fill(ptRandomConeRan-(bkg1*areaRandomCone));
+ clustSeq.get_median_rho_and_sigma(sortedJets, range, true, bkg4, sigma4, meanarea4, true);
+ fAODJetBackgroundOut->SetBackground(3,bkg4,sigma4,meanarea4);
- // fh1BiARandomCones[0]->Fill(omCone-(bkg1*areaRandomCone));
- // fh1BiARandomConesRan[0]->Fill(ptRandomConeRan-(bkg1*areaRandomCone));
+ //////////////////////
- clustSeq.get_median_rho_and_sigma(jets2, range, false, bkg2, sigma2, meanarea2, true);
- fAODJetBackgroundOut->SetBackground(1,bkg2,sigma2,meanarea2);
- // fh1BiARandomCones[1]->Fill(ptRandomCone-(bkg2*areaRandomCone));
- // fh1BiARandomConesRan[1]->Fill(ptRandomConeRan-(bkg2*areaRandomCone));
+ clustSeq.get_median_rho_and_sigma(jets2, range, false, bkg2, sigma2, meanarea2, true);
+ fAODJetBackgroundOut->SetBackground(1,bkg2,sigma2,meanarea2);
+ // fh1BiARandomCones[1]->Fill(ptRandomCone-(bkg2*areaRandomCone));
+ // fh1BiARandomConesRan[1]->Fill(ptRandomConeRan-(bkg2*areaRandomCone));
+
+ }
+ }
+
+ if(fStoreRhoLeadingTrackCorr) {
+ vector<fastjet::PseudoJet> jets3=sortedJets;
+ if(jets3.size()>2) jets3.erase(jets3.begin(),jets3.begin()+2);
+
+ Double_t bkg2=0;
+ Double_t sigma2=0.;
+ Double_t meanarea2=0.;
+
+ clustSeq.get_median_rho_and_sigma(jets3, range, false, bkg2, sigma2, meanarea2, true);
+ fAODJetBackgroundOut->SetBackground(1,bkg2,sigma2,meanarea2);
+
+ //Get leading track in event
+ AliVParticle *leading = (AliVParticle*)recParticles.At(0);
+
+ fh3CentvsRhoLeadingTrackPt->Fill(cent,bkg2,leading->Pt());
+ fh3CentvsSigmaLeadingTrackPt->Fill(cent,sigma2,leading->Pt());
+ fh3MultvsRhoLeadingTrackPt->Fill(nCh,bkg2,leading->Pt());
+ fh3MultvsSigmaLeadingTrackPt->Fill(nCh,sigma2,leading->Pt());
+
+
+ //Calculate rho with 4-vector area method for different quadrants with respect to the leading track in the event
+ //exclude 2 leading jets from event
+ //Quadrant 1: |phi_leadingTrack - phi_bkgCluster| < pi/2./2. - R (Near side to leading track)
+ //Quadrant 2: pi/2 - (pi/2./2. - R) < |phi_leadingTrack - phi_bkgCluster| < pi/2 + (pi/2./2. - R) (Orthogonal to leading track)
+ //Quadrant 3: pi - (pi/2/2 - R) < |phi_leadingTrack - phi_bkgCluster| < pi + (pi/2./2. - R) (Away side to leading track)
+ //Quadrant 4: 3/2*pi - (pi/2./2. - R) < |phi_leadingTrack - phi_bkgCluster| < 3/2*pi + (pi/2./2. - R) (Orthogonal to leading track)
+
+ //Assuming R=0.2 for background clusters
+
+ Double_t bkg2Q[4] = {0.};
+ Double_t sigma2Q[4] = {0.};
+ Double_t meanarea2Q[4] = {0.};
+
+ //Get phi of leading track in event
+ Float_t phiLeadingTrack = leading->Phi();
+ Float_t phiStep = (TMath::Pi()/2./2. - 0.2);
+
+ //Quadrant1 - near side
+ phiMin = phiLeadingTrack - phiStep;
+ phiMax = phiLeadingTrack + phiStep;
+ fastjet::RangeDefinition rangeQ1(rapMin,rapMax, phiMin, phiMax);
+ clustSeq.get_median_rho_and_sigma(jets3, rangeQ1, false, bkg2Q[0], sigma2Q[0], meanarea2Q[0], true);
+
+ //Quadrant2 - orthogonal
+ Double_t phiQ2 = phiLeadingTrack + TMath::Pi()/2.;
+ if(phiQ2>TMath::TwoPi()) phiQ2 = phiQ2 - TMath::TwoPi();
+ phiMin = phiQ2 - phiStep;
+ phiMax = phiQ2 + phiStep;
+ fastjet::RangeDefinition rangeQ2(rapMin,rapMax, phiMin, phiMax);
+ clustSeq.get_median_rho_and_sigma(jets3, rangeQ2, false, bkg2Q[1], sigma2Q[1], meanarea2Q[1], true);
+
+ //Quadrant3 - away side
+ Double_t phiQ3 = phiLeadingTrack + TMath::Pi();
+ if(phiQ3>TMath::TwoPi()) phiQ3 = phiQ3 - TMath::TwoPi();
+ phiMin = phiQ3 - phiStep;
+ phiMax = phiQ3 + phiStep;
+ fastjet::RangeDefinition rangeQ3(rapMin,rapMax, phiMin, phiMax);
+ clustSeq.get_median_rho_and_sigma(jets3, rangeQ3, false, bkg2Q[2], sigma2Q[2], meanarea2Q[2], true);
+
+ //Quadrant4 - othogonal
+ Double_t phiQ4 = phiLeadingTrack + 3./2.*TMath::Pi();
+ if(phiQ4>TMath::TwoPi()) phiQ4 = phiQ4 - TMath::TwoPi();
+ phiMin = phiQ4 - phiStep;
+ phiMax = phiQ4 + phiStep;
+ fastjet::RangeDefinition rangeQ4(rapMin,rapMax, phiMin, phiMax);
+ clustSeq.get_median_rho_and_sigma(jets3, rangeQ4, false, bkg2Q[3], sigma2Q[3], meanarea2Q[3], true);
+
+ fh3CentvsRhoLeadingTrackPtQ1->Fill(cent,bkg2Q[0],leading->Pt());
+ fh3CentvsRhoLeadingTrackPtQ2->Fill(cent,bkg2Q[1],leading->Pt());
+ fh3CentvsRhoLeadingTrackPtQ3->Fill(cent,bkg2Q[2],leading->Pt());
+ fh3CentvsRhoLeadingTrackPtQ4->Fill(cent,bkg2Q[3],leading->Pt());
+
+ fh3CentvsSigmaLeadingTrackPtQ1->Fill(cent,sigma2Q[0],leading->Pt());
+ fh3CentvsSigmaLeadingTrackPtQ2->Fill(cent,sigma2Q[1],leading->Pt());
+ fh3CentvsSigmaLeadingTrackPtQ3->Fill(cent,sigma2Q[2],leading->Pt());
+ fh3CentvsSigmaLeadingTrackPtQ4->Fill(cent,sigma2Q[3],leading->Pt());
+
+ fh3MultvsRhoLeadingTrackPtQ1->Fill(nCh,bkg2Q[0],leading->Pt());
+ fh3MultvsRhoLeadingTrackPtQ2->Fill(nCh,bkg2Q[1],leading->Pt());
+ fh3MultvsRhoLeadingTrackPtQ3->Fill(nCh,bkg2Q[2],leading->Pt());
+ fh3MultvsRhoLeadingTrackPtQ4->Fill(nCh,bkg2Q[3],leading->Pt());
+
+ fh3MultvsSigmaLeadingTrackPtQ1->Fill(nCh,sigma2Q[0],leading->Pt());
+ fh3MultvsSigmaLeadingTrackPtQ2->Fill(nCh,sigma2Q[1],leading->Pt());
+ fh3MultvsSigmaLeadingTrackPtQ3->Fill(nCh,sigma2Q[2],leading->Pt());
+ fh3MultvsSigmaLeadingTrackPtQ4->Fill(nCh,sigma2Q[3],leading->Pt());
+
+ fh3CentvsDeltaRhoLeadingTrackPtQ1->Fill(cent,bkg2Q[0]-bkg2,leading->Pt());
+ fh3CentvsDeltaRhoLeadingTrackPtQ2->Fill(cent,bkg2Q[1]-bkg2,leading->Pt());
+ fh3CentvsDeltaRhoLeadingTrackPtQ3->Fill(cent,bkg2Q[2]-bkg2,leading->Pt());
+ fh3CentvsDeltaRhoLeadingTrackPtQ4->Fill(cent,bkg2Q[3]-bkg2,leading->Pt());
- }
}
// do the same for tracks and jets
if(nTrackOver>0){
- TIterator *recIter = recParticles.MakeIterator();
- AliVParticle *tmpRec = (AliVParticle*)(recIter->Next());
- Float_t pt = tmpRec->Pt();
-
- // Printf("Leading track p_t %3.3E",pt);
- for(int i = 1;i <= fh2NRecTracksPt->GetNbinsX();i++){
- Float_t ptCut = fh2NRecTracksPt->GetXaxis()->GetBinCenter(i);
- while(pt<ptCut&&tmpRec){
- nTrackOver--;
- tmpRec = (AliVParticle*)(recIter->Next());
- if(tmpRec){
- pt = tmpRec->Pt();
- }
- }
- if(nTrackOver<=0)break;
- fh2NRecTracksPt->Fill(ptCut,nTrackOver);
- }
+ TIterator *recIter = recParticles.MakeIterator();
+ AliVParticle *tmpRec = (AliVParticle*)(recIter->Next());
+ Float_t pt = tmpRec->Pt();
+
+ // Printf("Leading track p_t %3.3E",pt);
+ for(int i = 1;i <= fh2NRecTracksPt->GetNbinsX();i++){
+ Float_t ptCut = fh2NRecTracksPt->GetXaxis()->GetBinCenter(i);
+
while(pt<ptCut&&tmpRec){
+ nTrackOver--;
+ tmpRec = (AliVParticle*)(recIter->Next());
+ if(tmpRec){
+ pt = tmpRec->Pt();
+ }
+ }
+ if(nTrackOver<=0)break;
+ fh2NRecTracksPt->Fill(ptCut,nTrackOver);
+ }
- recIter->Reset();
- AliVParticle *leading = (AliVParticle*)recParticles.At(0);
- Float_t phi = leading->Phi();
- if(phi<0)phi+=TMath::Pi()*2.;
- Float_t eta = leading->Eta();
- pt = leading->Pt();
- while((tmpRec = (AliVParticle*)(recIter->Next()))){
- Float_t tmpPt = tmpRec->Pt();
- Float_t tmpEta = tmpRec->Eta();
- fh1PtTracksRecIn->Fill(tmpPt);
- fh2TrackEtaPt->Fill(tmpEta,tmpPt);
- if(tmpRec==leading){
- fh1PtTracksLeadingRecIn->Fill(tmpPt);
- fh2LeadingTrackEtaPt->Fill(tmpEta,tmpPt);
- continue;
- }
+ recIter->Reset();
+ AliVParticle *leading = (AliVParticle*)recParticles.At(0);
+ Float_t phi = leading->Phi();
+ if(phi<0)phi+=TMath::Pi()*2.;
+ Float_t eta = leading->Eta();
+ pt = leading->Pt();
+ while((tmpRec = (AliVParticle*)(recIter->Next()))){
+ Float_t tmpPt = tmpRec->Pt();
+ Float_t tmpEta = tmpRec->Eta();
+ fh1PtTracksRecIn->Fill(tmpPt);
+ fh2TrackEtaPt->Fill(tmpEta,tmpPt);
+ if(tmpRec==leading){
+ fh1PtTracksLeadingRecIn->Fill(tmpPt);
+ fh2LeadingTrackEtaPt->Fill(tmpEta,tmpPt);
+ continue;
+ }
// correlation
- Float_t tmpPhi = tmpRec->Phi();
+ Float_t tmpPhi = tmpRec->Phi();
- if(tmpPhi<0)tmpPhi+=TMath::Pi()*2.;
- Float_t dPhi = phi - tmpPhi;
- if(dPhi>TMath::Pi())dPhi = dPhi - 2.*TMath::Pi();
- if(dPhi<(-1.*TMath::Pi()))dPhi = dPhi + 2.*TMath::Pi();
- Float_t dEta = eta - tmpRec->Eta();
- fh2TracksLeadingPhiEta->Fill(dPhi,dEta);
- fh2TracksLeadingPhiPt->Fill(dPhi,pt);
- fh2TracksLeadingPhiPtW->Fill(dPhi,pt,tmpPt);
- }
- delete recIter;
- }
+ if(tmpPhi<0)tmpPhi+=TMath::Pi()*2.;
+ Float_t dPhi = phi - tmpPhi;
+ if(dPhi>TMath::Pi())dPhi = dPhi - 2.*TMath::Pi();
+ if(dPhi<(-1.*TMath::Pi()))dPhi = dPhi + 2.*TMath::Pi();
+ Float_t dEta = eta - tmpRec->Eta();
+ fh2TracksLeadingPhiEta->Fill(dPhi,dEta);
+ fh2TracksLeadingPhiPt->Fill(dPhi,pt);
+ fh2TracksLeadingPhiPtW->Fill(dPhi,pt,tmpPt);
+ }
+ delete recIter;
+ }
- // find the random jets
+ // find the random jets
- fastjet::ClusterSequenceArea clustSeqRan(inputParticlesRecRan, jetDef, areaDef);
+ fastjet::ClusterSequenceArea clustSeqRan(inputParticlesRecRan, jetDef, areaDef);
- // fill the jet information from random track
- const vector <fastjet::PseudoJet> &inclusiveJetsRan = clustSeqRan.inclusive_jets();
- const vector <fastjet::PseudoJet> &sortedJetsRan = sorted_by_pt(inclusiveJetsRan);
+ // fill the jet information from random track
+ const vector <fastjet::PseudoJet> &inclusiveJetsRan = clustSeqRan.inclusive_jets();
+ const vector <fastjet::PseudoJet> &sortedJetsRan = sorted_by_pt(inclusiveJetsRan);
fh1NJetsRecRan->Fill(sortedJetsRan.size());
- // loop over all jets an fill information, first on is the leading jet
+ // loop over all jets an fill information, first on is the leading jet
- Int_t nRecOverRan = inclusiveJetsRan.size();
- Int_t nRecRan = inclusiveJetsRan.size();
+ Int_t nRecOverRan = inclusiveJetsRan.size();
+ Int_t nRecRan = inclusiveJetsRan.size();
- if(inclusiveJetsRan.size()>0){
- AliAODJet leadingJet (sortedJetsRan[0].px(), sortedJetsRan[0].py(), sortedJetsRan[0].pz(), sortedJetsRan[0].E());
- Float_t pt = leadingJet.Pt();
+ if(inclusiveJetsRan.size()>0){
+ AliAODJet leadingJet (sortedJetsRan[0].px(), sortedJetsRan[0].py(), sortedJetsRan[0].pz(), sortedJetsRan[0].E());
+ Float_t pt = leadingJet.Pt();
- Int_t iCount = 0;
- TLorentzVector vecarearanb;
+ Int_t iCount = 0;
+ TLorentzVector vecarearanb;
- for(int i = 1;i <= fh2NRecJetsPtRan->GetNbinsX();i++){
- Float_t ptCut = fh2NRecJetsPtRan->GetXaxis()->GetBinCenter(i);
+ for(int i = 1;i <= fh2NRecJetsPtRan->GetNbinsX();i++){
+ Float_t ptCut = fh2NRecJetsPtRan->GetXaxis()->GetBinCenter(i);
while(pt<ptCut&&iCount<nRecRan){
nRecOverRan--;
iCount++;
}
Int_t nAodOutJetsRan = 0;
- AliAODJet *aodOutJetRan = 0;
+ AliAODJet *aodOutJetRan = 0;
for(int j = 0; j < nRecRan;j++){
AliAODJet tmpRec (sortedJetsRan[j].px(), sortedJetsRan[j].py(), sortedJetsRan[j].pz(), sortedJetsRan[j].E());
Float_t tmpPt = tmpRec.Pt();
fh2PtNchNRan->Fill(nCh,tmpPt,constituents.size());
- if(tmpPt>fJetOutputMinPt&&fTCAJetsOutRan){
- aodOutJetRan = new ((*fTCAJetsOutRan)[nAodOutJetsRan++]) AliAODJet(tmpRec);
- Double_t arearan=clustSeqRan.area(sortedJetsRan[j]);
- aodOutJetRan->GetRefTracks()->Clear();
- aodOutJetRan->SetEffArea(arearan,0);
- fastjet::PseudoJet vecarearan=clustSeqRan.area_4vector(sortedJetsRan[j]);
- vecarearanb.SetPxPyPzE(vecarearan.px(),vecarearan.py(),vecarearan.pz(),vecarearan.e());
- aodOutJetRan->SetVectorAreaCharged(&vecarearanb);
+ if(tmpPt>fJetOutputMinPt&&fTCAJetsOutRan){
+ aodOutJetRan = new ((*fTCAJetsOutRan)[nAodOutJetsRan++]) AliAODJet(tmpRec);
+ Double_t arearan=clustSeqRan.area(sortedJetsRan[j]);
+ aodOutJetRan->GetRefTracks()->Clear();
+ aodOutJetRan->SetEffArea(arearan,0);
+ fastjet::PseudoJet vecarearan=clustSeqRan.area_4vector(sortedJetsRan[j]);
+ vecarearanb.SetPxPyPzE(vecarearan.px(),vecarearan.py(),vecarearan.pz(),vecarearan.e());
+ aodOutJetRan->SetVectorAreaCharged(&vecarearanb);
- }
+ }
// correlation
Float_t tmpPhi = tmpRec.Phi();
if(fAODJetBackgroundOut){
- Double_t bkg3=0.;
- Double_t sigma3=0.;
- Double_t meanarea3=0.;
- clustSeqRan.get_median_rho_and_sigma(sortedJetsRan ,range, false, bkg3, sigma3, meanarea3, true);
- fAODJetBackgroundOut->SetBackground(2,bkg3,sigma3,meanarea3);
- // float areaRandomCone = rRandomCone2 *TMath::Pi();
- /*
- fh1BiARandomCones[2]->Fill(ptRandomCone-(bkg3*areaRandomCone));
- fh1BiARandomConesRan[2]->Fill(ptRandomConeRan-(bkg3*areaRandomCone));
- */
+ Double_t bkg3=0.;
+ Double_t sigma3=0.;
+ Double_t meanarea3=0.;
+ clustSeqRan.get_median_rho_and_sigma(sortedJetsRan ,range, false, bkg3, sigma3, meanarea3, true);
+ fAODJetBackgroundOut->SetBackground(2,bkg3,sigma3,meanarea3);
+ // float areaRandomCone = rRandomCone2 *TMath::Pi();
+ /*
+ fh1BiARandomCones[2]->Fill(ptRandomCone-(bkg3*areaRandomCone));
+ fh1BiARandomConesRan[2]->Fill(ptRandomConeRan-(bkg3*areaRandomCone));
+ */
}
- }
+ }
- // do the event selection if activated
- if(fJetTriggerPtCut>0){
- bool select = false;
- Float_t minPt = fJetTriggerPtCut;
- /*
- // hard coded for now ...
- // 54.50 44.5 29.5 18.5 for anti-kt rejection 1E-3
- if(cent<10)minPt = 50;
- else if(cent<30)minPt = 42;
- else if(cent<50)minPt = 28;
- else if(cent<80)minPt = 18;
- */
- float rho = 0;
- if(externalBackground)rho = externalBackground->GetBackground(2);
- if(fTCAJetsOut){
- for(int i = 0;i < fTCAJetsOut->GetEntriesFast();i++){
- AliAODJet *jet = (AliAODJet*)fTCAJetsOut->At(i);
- Float_t ptSub = jet->Pt() - rho *jet->EffectiveAreaCharged();
- if(ptSub>=minPt){
- select = true;
- break;
- }
- }
- }
+ // do the event selection if activated
+ if(fJetTriggerPtCut>0){
+ bool select = false;
+ Float_t minPt = fJetTriggerPtCut;
+ /*
+ // hard coded for now ...
+ // 54.50 44.5 29.5 18.5 for anti-kt rejection 1E-3
+ if(cent<10)minPt = 50;
+ else if(cent<30)minPt = 42;
+ else if(cent<50)minPt = 28;
+ else if(cent<80)minPt = 18;
+ */
+ float rho = 0;
+ if(externalBackground)rho = externalBackground->GetBackground(2);
+ if(fTCAJetsOut){
+ for(int i = 0;i < fTCAJetsOut->GetEntriesFast();i++){
+ AliAODJet *jet = (AliAODJet*)fTCAJetsOut->At(i);
+ Float_t ptSub = jet->Pt() - rho *jet->EffectiveAreaCharged();
+ if(ptSub>=minPt){
+ select = true;
+ break;
+ }
+ }
+ }
- if(select){
- static AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
- fh1CentralitySelect->Fill(cent);
- fh1ZSelect->Fill(zVtx);
- aodH->SetFillAOD(kTRUE);
- }
- }
- if (fDebug > 2){
- if(fTCAJetsOut)Printf("%s:%d Rec Jets %d",(char*)__FILE__,__LINE__,fTCAJetsOut->GetEntriesFast());
- if(fTCAJetsOutRan)Printf("%s:%d Rec Jets Ran %d",(char*)__FILE__,__LINE__,fTCAJetsOutRan->GetEntriesFast());
- if(fTCARandomConesOut)Printf("%s:%d RC %d",(char*)__FILE__,__LINE__,fTCARandomConesOut->GetEntriesFast());
- if(fTCARandomConesOutRan)Printf("%s:%d RC Ran %d",(char*)__FILE__,__LINE__,fTCARandomConesOutRan->GetEntriesFast());
- }
- PostData(1, fHistList);
+ if(select){
+
static AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
+ fh1CentralitySelect->Fill(cent);
+ fh1ZSelect->Fill(zVtx);
+ aodH->SetFillAOD(kTRUE);
+ }
+ }
+ if (fDebug > 2){
+ if(fTCAJetsOut)Printf("%s:%d Rec Jets %d",(char*)__FILE__,__LINE__,fTCAJetsOut->GetEntriesFast());
+ if(fTCAJetsOutRan)Printf("%s:%d Rec Jets Ran %d",(char*)__FILE__,__LINE__,fTCAJetsOutRan->GetEntriesFast());
+ if(fTCARandomConesOut)Printf("%s:%d RC %d",(char*)__FILE__,__LINE__,fTCARandomConesOut->GetEntriesFast());
+ if(fTCARandomConesOutRan)Printf("%s:%d RC Ran %d",(char*)__FILE__,__LINE__,fTCARandomConesOutRan->GetEntriesFast());
+ }
+ PostData(1, fHistList);
}
void AliAnalysisTaskJetCluster::Terminate(Option_t */*option*/)
//
// Terminate analysis
//
- if (fDebug > 1) printf("AnalysisJetCluster: Terminate() \n");
+ if (fDebug > 1) printf("AnalysisJetCluster: Terminate() \n");
- if(fMomResH1Fit) delete fMomResH1Fit;
- if(fMomResH2Fit) delete fMomResH2Fit;
- if(fMomResH3Fit) delete fMomResH3Fit;
+ if(fMomResH1Fit) delete fMomResH1Fit;
+ if(fMomResH2Fit) delete fMomResH2Fit;
+ if(fMomResH3Fit) delete fMomResH3Fit;
}
// get list of tracks/particles for different types
//
- if(fDebug>2)Printf("%s:%d Selecting tracks with %d",(char*)__FILE__,__LINE__,type);
+ if(fDebug>2) Printf("%s:%d Selecting tracks with %d",(char*)__FILE__,__LINE__,type);
Int_t iCount = 0;
- if(type==kTrackAOD || type==kTrackAODextra || type==kTrackAODextraonly){
- if(type!=kTrackAODextraonly) {
+ if(type==kTrackAOD || type==kTrackAODextra || type==kTrackAODextraonly || type==kTrackAODMCextra || type==kTrackAODMCextraonly){
+
+ if(type!=kTrackAODextraonly && type!=kTrackAODMCextraonly) {
AliAODEvent *aod = 0;
if(fUseAODTrackInput)aod = dynamic_cast<AliAODEvent*>(InputEvent());
else aod = AODEvent();
if(fDebug>10)Printf("%s:%d Not matching filter %d/%d %d/%d",(char*)__FILE__,__LINE__,it,aod->GetNumberOfTracks(),fFilterMask,tr->GetFilterMap());
continue;
}
+
+ // heavy flavor jets
+ if(fFilterMask==528 && fUseHFcuts){
+ Double_t ntpcClus = tr->GetTPCNcls();
+ Double_t trPt=tr->Pt();
+ TFormula NTPCClsCut("f1NClustersTPCLinearPtDep","70.+30./20.*x");
+
+ if (trPt <= 20. && (ntpcClus < NTPCClsCut.Eval(trPt))) continue;
+ else if (trPt > 20. && ntpcClus < 100) continue;
+
+ if (AvoidDoubleCountingHF(aod,tr)) continue;
+ }
+ //
+
+ if(fRequireITSRefit){if((tr->GetStatus()&AliESDtrack::kITSrefit)==0)continue;}
+ if (fApplySharedClusterCut) {
+ Double_t frac = Double_t(tr->GetTPCnclsS()) /Double_t(tr->GetTPCncls());
+ if (frac > 0.4) continue;
+ }
if(TMath::Abs(tr->Eta())>fTrackEtaWindow){
if(fDebug>10)Printf("%s:%d Not matching eta %d/%d",(char*)__FILE__,__LINE__,it,aod->GetNumberOfTracks());
continue;
else if(fFilterType == 1)bGood = trackAOD->IsHybridTPCConstrainedGlobal();
else if(fFilterType == 2)bGood = trackAOD->IsHybridGlobalConstrainedGlobal();
if((fFilterMask>0)&&((!trackAOD->TestFilterBit(fFilterMask)||(!bGood))))continue;
- if(TMath::Abs(trackAOD->Eta())>fTrackEtaWindow) continue;
+ if(fRequireITSRefit){if((trackAOD->GetStatus()&AliESDtrack::kITSrefit)==0)continue;}
+ if (fApplySharedClusterCut) {
+ Double_t frac = Double_t(trackAOD->GetTPCnclsS()) /Double_t(trackAOD->GetTPCncls());
+ if (frac > 0.4) continue;
+ }
+
+
+ if(TMath::Abs(trackAOD->Eta())>fTrackEtaWindow) continue;
if(trackAOD->Pt()<fTrackPtCut) continue;
+ if(fDebug) printf("pt extra track %.2f \n", trackAOD->Pt());
list->Add(trackAOD);
iCount++;
}
}
+
+ if(type==kTrackAODMCextra || type==kTrackAODMCextraonly) { //embed generator level particles
+ AliAODEvent *aod = 0;
+ if(fUseAODTrackInput)aod = dynamic_cast<AliAODEvent*>(InputEvent());
+ else aod = AODEvent();
+ if(!aod){
+ return iCount;
+ }
+ TClonesArray *aodExtraTracks = dynamic_cast<TClonesArray*>(aod->FindListObject("aodExtraMCparticles"));
+ if(!aodExtraTracks)return iCount;
+ for(int it =0; it<aodExtraTracks->GetEntries(); it++) {
+ AliVParticle *track = dynamic_cast<AliVParticle*> ((*aodExtraTracks)[it]);
+ AliAODMCParticle *partmc = dynamic_cast<AliAODMCParticle*> ((*aodExtraTracks)[it]);
+ if (!track) {
+ if(fDebug>10) printf("track %d does not exist\n",it);
+ continue;
+ }
+
+ if(!partmc) continue;
+ if(!partmc->IsPhysicalPrimary())continue;
+
+ if (track->Pt()<fTrackPtCut) {
+ if(fDebug>10) printf("track %d has too low pt %.2f\n",it,track->Pt());
+ continue;
+ }
+
+ /*
+ AliAODTrack *trackAOD = dynamic_cast<AliAODTrack*>((*aodExtraTracks)[it]);//(track);
+
+ if(!trackAOD) {
+ if(fDebug>10) printf("trackAOD %d does not exist\n",it);
+ continue;
+ }
+
+ trackAOD->SetFlags(AliESDtrack::kEmbedded);
+ trackAOD->SetFilterMap(fFilterMask);
+ */
+ if(fDebug>10) printf("pt extra track %.2f \n", track->Pt());
+
+ if(TMath::Abs(track->Eta())>fTrackEtaWindow) continue;
+ if(track->Pt()<fTrackPtCut) continue;
+ list->Add(track);
+
+ iCount++;
+ }
+ }
+
}
else if (type == kTrackKineAll||type == kTrackKineCharged){
AliMCEvent* mcEvent = MCEvent();
}
}
}// AODMCparticle
+ else if (type == kTrackAODMCHF){
+
+ AliAODEvent *aod = 0;
+ if(fUseAODMCInput)aod = dynamic_cast<AliAODEvent*>(InputEvent());
+ else aod = AODEvent();
+ if(!aod)return iCount;
+ TClonesArray *tca = dynamic_cast<TClonesArray*>(aod->FindListObject(AliAODMCParticle::StdBranchName()));
+ if(!tca)return iCount;
+ for(int it = 0;it < tca->GetEntriesFast();++it){
+ AliAODMCParticle *part = (AliAODMCParticle*)(tca->At(it));
+ if(!part) continue;
+ int partpdg= part->PdgCode();
+ if(!part->IsPhysicalPrimary() && !IsBMeson(partpdg) && !IsDMeson(partpdg) )continue;
+
+ if (IsBMeson(partpdg) || IsDMeson(partpdg)) {
+ iCount+= AddDaughters( list , part,tca);
+ }
+ else {
+
+ if(part->Pt()<fTrackPtCut) continue;
+ if(TMath::Abs(part->Eta())>fTrackEtaWindow) continue;
+ if(part->Charge()==0) continue;
+
+ if((part->Pt()>=fTrackPtCut) && (TMath::Abs(part->Eta())<=fTrackEtaWindow) && (part->Charge()!=0))list->Add(part);
+ iCount++;
+ }
+ }
+ }
+
list->Sort();
return iCount;
}
+Int_t AliAnalysisTaskJetCluster::AddDaughters(TList * list, AliAODMCParticle *part, TClonesArray * tca){
+ Int_t count=0;
+ Int_t nDaugthers = part->GetNDaughters();
+ for(Int_t i=0;i< nDaugthers;++i){
+ AliAODMCParticle *partdaughter = (AliAODMCParticle*)(tca->At(i));
+ if(!partdaughter) continue;
+ if(partdaughter->Pt()<fTrackPtCut)continue;
+ if(TMath::Abs(partdaughter->Eta())>fTrackEtaWindow)continue;
+ if(partdaughter->Charge()==0)continue;
+
+ if(!IsDMeson(partdaughter->PdgCode()) && !IsBMeson(partdaughter->PdgCode()) ){
+ list->Add(partdaughter);
+ count++;
+ }
+ else count+=AddDaughters(list,part,tca);
+ }
+return count;
+}
+
+
+Bool_t AliAnalysisTaskJetCluster::AvoidDoubleCountingHF(AliAODEvent *aod, AliAODTrack *tr1){
+
+ if(!(tr1->TestFilterBit(BIT(9)))) return kFALSE;
+
+ Int_t idtr1 = tr1->GetID();
+
+ for(int jt = 0;jt < aod->GetNumberOfTracks();++jt){
+
+ const AliAODTrack *tr2 = aod->GetTrack(jt);
+ Int_t idtr2 = tr2->GetID();
+
+ if (!(tr2->TestFilterBit(BIT(4)))) continue;
+ if (idtr1==(idtr2+1)*-1.) return kTRUE;
+
+ }
+ return kFALSE;
+}
void AliAnalysisTaskJetCluster::LoadTrPtResolutionRootFileFromOADB() {
- TFile *f = new TFile(fPathTrPtResolution.Data());
+ if (!gGrid) {
+ AliInfo("Trying to connect to AliEn ...");
+ TGrid::Connect("alien://");
+ }
+
+ TFile *f = TFile::Open(fPathTrPtResolution.Data());
if(!f)return;
SetSmearResolution(kTRUE);
SetMomentumResolutionHybrid(fProfPtPtSigma1PtGlobSt,fProfPtPtSigma1PtGlobCnoITS,fProfPtPtSigma1PtGlobCnoSPD);
-
-
- if(f) delete f;
+
}
void AliAnalysisTaskJetCluster::LoadTrEfficiencyRootFileFromOADB() {
- TFile *f = new TFile(fPathTrEfficiency.Data());
+ if (!gGrid) {
+ AliInfo("Trying to connect to AliEn ...");
+ TGrid::Connect("alien://");
+ }
+
+ TFile *f = TFile::Open(fPathTrEfficiency.Data());
if(!f)return;
TH1D *hEffPosGlobSt = (TH1D*)f->Get("hEffPosGlobSt");
else
SetEfficiencyHybrid(hEffPosGlobSt,hEffPosGlobCnoITS,hEffPosGlobCnoSPD);
- if(f) delete f;
-
}
void AliAnalysisTaskJetCluster::SetMomentumResolutionHybrid(TProfile *p1, TProfile *p2, TProfile *p3) {
// set mom res profiles
//
- fMomResH1 = (TProfile*)p1->Clone("fMomResH1");
- fMomResH2 = (TProfile*)p2->Clone("fMomResH2");
- fMomResH3 = (TProfile*)p3->Clone("fMomResH3");
+ if(fMomResH1) delete fMomResH1;
+ if(fMomResH2) delete fMomResH2;
+ if(fMomResH3) delete fMomResH3;
+
+ fMomResH1 = new TProfile(*p1);//(TProfile*)p1->Clone("fMomResH1");
+ fMomResH2 = new TProfile(*p2);//(TProfile*)p2->Clone("fMomResH2");
+ fMomResH3 = new TProfile(*p3);//(TProfile*)p3->Clone("fMomResH3");
+
}
void AliAnalysisTaskJetCluster:: SetEfficiencyHybrid(TH1 *h1, TH1 *h2, TH1 *h3) {
}
/*
-Int_t AliAnalysisTaskJetCluster::AddParticlesFastJet(TList &particles,vector<fastjet::PseudoJet> &inputParticles){
+ Int_t AliAnalysisTaskJetCluster::AddParticlesFastJet(TList &particles,vector<fastjet::PseudoJet> &inputParticles){
for(int i = 0; i < particles.GetEntries(); i++){
- AliVParticle *vp = (AliVParticle*)particles.At(i);
- // Carefull energy is not well determined in real data, should not matter for p_T scheme?
- fastjet::PseudoJet jInp(vp->Px(),vp->Py(),vp->Pz(),vp->E());
- jInp.set_user_index(i);
- inputParticles.push_back(jInp);
+ AliVParticle *vp = (AliVParticle*)particles.At(i);
+ // Carefull energy is not well determined in real data, should not matter for p_T scheme?
+ fastjet::PseudoJet jInp(vp->Px(),vp->Py(),vp->Pz(),vp->E());
+ jInp.set_user_index(i);
+ inputParticles.push_back(jInp);
}
return 0;
-}
+ }
*/
+
+
+bool AliAnalysisTaskJetCluster::IsBMeson(int pc){
+ int bPdG[] = {511,521,10511,10521,513,523,10513,10523,20513,20523,20513,20523,515,525,531,
+ 10531,533,10533,20533,535,541,10541,543,10543,20543,545,551,10551,100551,
+ 110551,200551,210551,553,10553,20553,30553,100553,110553,120553,130553,200553,210553,220553,
+ 300553,9000533,9010553,555,10555,20555,100555,110555,120555,200555,557,100557};
+ for(int i=0;i< (int)(sizeof(bPdG)/sizeof(int));++i) if(abs(pc) == bPdG[i]) return true;
+return false;
+}
+bool AliAnalysisTaskJetCluster::IsDMeson(int pc){
+ int bPdG[] = {411,421,10411,10421,413,423,10413,10423,20431,20423,415,
+ 425,431,10431,433,10433,20433,435,441,10441,100441,443,10443,20443,
+ 100443,30443,9000443,9010443,9020443,445,100445};
+ for(int i=0;i< (int)(sizeof(bPdG)/sizeof(int));++i) if(abs(pc) == bPdG[i]) return true;
+return false;
+}
+
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff