Rivet analyses referenceMC_JET_IN_HIMonte Carlo validation, simple jet analysis in heavy ion collisions, using centrality framework.Experiment: () Status: VALIDATED Authors:
Beams: * * Beam energies: ANY Run details:
1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/FastJets.hh"
5#include "Rivet/Projections/WFinder.hh"
6#include "Rivet/Projections/ZFinder.hh"
7#include "fastjet/tools/Filter.hh"
8#include "fastjet/tools/Pruner.hh"
9#include "Rivet/Tools/AtlasCommon.hh"
10
11namespace Rivet {
12
13 class MC_JET_IN_HI : public Analysis {
14 public:
15
16 /// @name Constructors etc.
17 //@{
18
19 /// Constructor
20 RIVET_DEFAULT_ANALYSIS_CTOR(MC_JET_IN_HI);
21
22 //@}
23
24
25 public:
26
27 string ts(int in) {
28 std::stringstream ss;
29 ss << in;
30 return ss.str();
31 }
32 /// @name Analysis methods
33 //@{
34
35 /// Book histograms and initialise projections before the run
36 void init() {
37 // Declare centrality projection - we use the ATLAS PbPb definition
38 // to be able to compare to data.
39 declareCentrality(ATLAS::SumET_PBPB_Centrality(),"ATLAS_PBPB_CENTRALITY",
40 "sumETFwd","sumETFwd");
41 // The final state where jets are found
42 FinalState fs(Cuts::abseta < 2.5);
43 declare(fs, "FS");
44
45 ZFinder zfinder(fs, Cuts::abseta < 2.5 && Cuts::pT > 30*GeV, PID::MUON, 80*GeV, 100*GeV, 0.2,
46 ZFinder::ChargedLeptons::PROMPT, ZFinder::ClusterPhotons::NODECAY, ZFinder::AddPhotons::YES);
47 declare(zfinder, "ZFinder");
48
49 // Z+jet jet collections
50 declare(FastJets(zfinder.remainingFinalState(), FastJets::ANTIKT, 0.3), "JetsAK3");
51 declare(FastJets(zfinder.remainingFinalState(), FastJets::ANTIKT, 0.5), "JetsAK5");
52 declare(FastJets(zfinder.remainingFinalState(), FastJets::ANTIKT, 0.7), "JetsAK7");
53 declare(FastJets(zfinder.remainingFinalState(), FastJets::ANTIKT, 0.9), "JetsAK9");
54
55 jetFinders = {"JetsAK3", "JetsAK5", "JetsAK7", "JetsAK9"};
56
57 h_zpT.resize(jetFinders.size());
58 h_jetpT.resize(jetFinders.size());
59 for (size_t i = 0; i < jetFinders.size(); ++i) {
60 string s = jetFinders[i];
61 book(h_zpT[i], s + "zpT",logspace(50, 1.0,1000));
62 book(h_jetpT[i], s + "jetpT",logspace(50, 1.0,1000));
63 }
64 book(incSow, "incSow");
65
66 centData = {0., 0.2, 0.4, 0.6, 0.8,};
67 for (size_t i = 0; i < centData.size(); ++i) {
68 book(c_jetpT[centData[i]], "cjetpT" + ts(i),logspace(100, 10.0,1000));
69 book(c_zpT[centData[i]], "czpt" + ts(i),logspace(100, 10.0,1000));
70 book(sow[centData[i]], "sow_" + ts(i));
71 }
72 }
73
74 bool isBackToBack_zj(const ZFinder& zf, const fastjet::PseudoJet& psjet) {
75 const FourMomentum& z = zf.bosons()[0].momentum();
76 const FourMomentum jmom(psjet.e(), psjet.px(), psjet.py(), psjet.pz());
77 return (deltaPhi(z, jmom) > 7.*M_PI/8. );
78 }
79
80
81 /// Perform the per-event analysis
82 void analyze(const Event& event) {
83
84 // Get the Z
85 const ZFinder& zfinder = apply<ZFinder>(event, "ZFinder");
86 if (zfinder.bosons().size() != 1) vetoEvent;
87 Particle z = zfinder.bosons()[0];
88 Particle l1 = zfinder.constituents()[0];
89 Particle l2 = zfinder.constituents()[1];
90
91 // Require a high-pT Z (and constituents)
92 if (l1.pT() < 10*GeV || l2.pT() < 10*GeV || z.pT() < 60*GeV) vetoEvent;
93
94 // Get the centrality
95 const double c = apply<CentralityProjection>(event,"sumETFwd")();
96 auto jetpTItr = c_jetpT.upper_bound(c);
97 auto zpTItr = c_zpT.upper_bound(c);
98 auto sowItr = sow.upper_bound(c);
99 if (jetpTItr == c_jetpT.end() || zpTItr == c_zpT.end() ||
100 sowItr == sow.end()) vetoEvent;
101 sowItr->second->fill();
102 incSow->fill();
103 // Get the jets
104 for (size_t i = 0; i < jetFinders.size(); ++i ) {
105 const PseudoJets& psjets = apply<FastJets>(event,
106 jetFinders[i]).pseudoJetsByPt(30.0*GeV);
107 if (!psjets.empty()) {
108 // Get the leading jet and make sure it's back-to-back with the Z
109 const fastjet::PseudoJet& j0 = psjets[0];
110 if (isBackToBack_zj(zfinder, j0)) {
111 // Fill the centrality inclusive histograms
112 h_zpT[i]->fill(z.pT());
113 h_jetpT[i]->fill(j0.perp());
114 // Fill centrality dept histograms only for R = 0.3
115 if (i == 0) {
116 jetpTItr->second->fill(j0.perp());
117 zpTItr->second->fill(z.pT());
118 }
119 }
120 }
121 }
122 }
123
124
125 /// Normalise histograms etc., after the run
126 void finalize() {
127 for(size_t i = 0; i < jetFinders.size(); ++i) {
128 h_zpT[i]->scaleW(1./incSow->sumW());
129 h_jetpT[i]->scaleW(1./incSow->sumW());
130 }
131 for (size_t i = 0; i < centData.size(); ++i) {
132 c_jetpT[centData[i]]->scaleW(1./sow[centData[i]]->sumW());
133 c_zpT[centData[i]]->scaleW(1./sow[centData[i]]->sumW());
134 }
135 }
136
137
138 //@}
139
140
141 private:
142 vector<string> jetFinders;
143 // Centrality inclusive histograms
144 vector<Histo1DPtr> h_zpT;
145 vector<Histo1DPtr> h_jetpT;
146 CounterPtr incSow;
147 // Centrality intervals
148 vector<double> centData;
149 // Centrality binned histograms
150 map<double, Histo1DPtr> c_jetpT;
151 map<double, Histo1DPtr> c_zpT;
152 map<double, CounterPtr> sow;
153
154 };
155
156
157 // The hook for the plugin system
158 RIVET_DECLARE_PLUGIN(MC_JET_IN_HI);
159
160}
|