Rivet analyses referenceATLAS_2014_I1307243Measurements of jet vetoes and azimuthal decorrelations in dijet events produced in pp collisions at $\sqrt{s}$ = 7 TeV using the ATLAS detectorExperiment: ATLAS (LHC) Inspire ID: 1307243 Status: VALIDATED Authors:
Beam energies: (3500.0, 3500.0) GeV Run details:
Additional jet activity in dijet events is measured using $pp$ collisions at ATLAS at a centre-of-mass energy of 7 TeV, for jets reconstructed using the anti-kt algorithm with radius parameter R=0.6. This is done using variables such as the fraction of dijet events without an additional jet in the rapidity interval bounded by the dijet subsystem and correlations between the azimuthal angles of the dijets. They are presented, both with and without a veto on additional jet activity in the rapidity interval, as a function of the mean transverse momentum of the dijets and of the rapidity interval size. The double differential dijet cross section is also measured as a function of the interval size and the azimuthal angle between the dijets. These variables probe differences in the approach to resummation of large logarithms when performing QCD calculations. The data are compared to POWHEG, interfaced to the PYTHIA 8 and HERWIG parton shower generators, as well as to HEJ with and without interfacing it to the ARIADNE parton shower generator. None of the theoretical predictions agree with the data across the full phase-space considered; however, POWHEG+PYTHIA 8 and HEJ+ARIADNE are found to provide the best agreement with the data. These measurements use the full data sample collected with the ATLAS detector in 7 TeV $pp$ collisions at the LHC and correspond to integrated luminosities of 36.1 pb$^-1$ and 4.5 fb$^-1$ for data collected during 2010 and 2011 respectively. Source code: ATLAS_2014_I1307243.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Tools/Logging.hh"
4#include "Rivet/Projections/FastJets.hh"
5#include "Rivet/Tools/BinnedHistogram.hh"
6
7namespace Rivet {
8
9
10 /// @brief ATLAS azimuthal decorrelation with jet veto analysis
11 /// @author James Robinson <james.robinson@cern.ch>
12 class ATLAS_2014_I1307243 : public Analysis {
13 public:
14
15 /// Constructor
16 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2014_I1307243);
17
18
19 /// Book histograms and initialise projections before the run
20 void init() {
21
22 _dy_max = 8; _years = { 2010, 2011 };
23 _Qnoughts = { 20., 30., 40., 50., 60., 70., 80., 90., 100. };
24 /// Initialise and register projections here
25 FastJets fastJets(FinalState(), FastJets::ANTIKT, 0.6, JetAlg::Muons::ALL, JetAlg::Invisibles::ALL);
26 declare(fastJets, "AntiKt6JetsWithInvisibles");
27
28
29 /// Book histograms
30 for (const string cat : { "inclusive", "gap" }) {
31 const size_t offset = (cat == "gap") ? 1 : 0;
32
33 // Temporary inclusive and gap histograms
34 book(_aux_dy[cat], "_" + cat + "_dy", refData(1, 1, 1));
35 book(_aux_pTbar[cat], "_" + cat + "_pTbar", refData(2, 1, 1));
36
37 book(_h_C2C1_dy[cat], 7 + 4 * offset, 1, 1, false);
38 book(_h_C2C1_pTbar[cat], 8 + 4 * offset, 1, 1, false);
39
40 // Azimuthal moment histograms
41 book(_p_cosDeltaPhi_dy[cat], 5 + 4 * offset, 1, 1);
42 book(_p_cosDeltaPhi_pTbar[cat], 6 + 4 * offset, 1, 1);
43 book(_p_cosTwoDeltaPhi_dy[cat], 37 + 2 * offset, 1, 1);
44 book(_p_cosTwoDeltaPhi_pTbar[cat], 38 + 2 * offset, 1, 1);
45
46 // Gap fraction vs. Q0 and cross-section in dy slices
47 _s_gapFrac_Q0.resize(_dy_max);
48 for (size_t dyLow = 0; dyLow < _dy_max; ++dyLow ) {
49 const string hname("_" + cat + "_dySlice_" + toString(dyLow) + "_" + toString(dyLow+1) + "_Q0");
50 { Histo1DPtr tmp; _aux_Q0_dySlices[cat].add(dyLow, dyLow+1, book(tmp, hname, refData(29+dyLow, 1, 1))); }
51 { Histo1DPtr tmp; _h_dphi_dySlices[cat].add(dyLow, dyLow+1, book(tmp, 13+(_dy_max*offset)+dyLow, 1, 1)); }
52 if (!offset) book(_s_gapFrac_Q0[dyLow], 29 + dyLow, 1, 1); // only book once
53 }
54
55 }
56
57 // Number of jets in rapidity interval
58 book(_s_gapFrac_dy, 1, 1, 1);
59 book(_s_gapFrac_pTbar, 2, 1, 1);
60 book(_p_nGapJets_dy, 3, 1, 1);
61 book(_p_nGapJets_pTbar, 4, 1, 1);
62 }
63
64
65 /// Perform the per-event analysis
66 void analyze(const Event& event) {
67
68 for (size_t reg = 0; reg < 2; ++reg) {
69
70 // Retrieve all anti-kt R=0.6 jets
71 const double maxRap = reg? 2.4 : 4.4;
72 const Jets& akt6Jets = apply<JetAlg>(event, "AntiKt6JetsWithInvisibles").jetsByPt(Cuts::absrap < maxRap);
73 // If there are fewer than 2 jets then bail
74 if ( akt6Jets.size() < 2 ) { vetoEvent; }
75
76 // Require jets to be above {60, 50} GeV
77 if ( akt6Jets[0].pT() < 60*GeV || akt6Jets[1].pT() < 50*GeV ) { vetoEvent; }
78
79 // Identify gap boundaries
80 const double yMin = std::min(akt6Jets[0].rap(), akt6Jets[1].rap());
81 const double yMax = std::max(akt6Jets[0].rap(), akt6Jets[1].rap());
82
83 // Determine azimuthal decorrelation quantities
84 const double dy = yMax - yMin;
85 const double dphi = mapAngle0ToPi(deltaPhi(akt6Jets[0], akt6Jets[1]));
86 const double pTbar = 0.5*(akt6Jets[0].pT() + akt6Jets[1].pT())/GeV;
87
88 // Impose minimum dy for the 2011 phase space
89 if ( _years[reg] == 2011 && dy < 1.0 ) { vetoEvent; }
90
91 // Determine gap quantities
92 size_t nGapJets = 0;
93 double maxGapQ0 = 0.;
94 const double vetoScale = reg? 30*GeV : 20*GeV;
95 for (const Jet& jet : akt6Jets) {
96 if (!inRange(jet.rap(), yMin, yMax, OPEN, OPEN)) continue;
97 const double pT = jet.pT()/GeV;
98 if (pT > vetoScale) { ++nGapJets; }
99 if (pT > maxGapQ0) { maxGapQ0 = pT; }
100 }
101
102 // Fill histograms
103 fillHists(_years[reg], nGapJets, { dy, pTbar, dphi, maxGapQ0 });
104 }
105 return;
106 }
107
108 void fillHists(const size_t region, const size_t nGapJets, const vector<double>& vars) {
109 assert(vars.size() == 4);
110 const double dy = vars[0];
111 const double pTbar = vars[1];
112 const double dphi = vars[2];
113 const double maxGapQ0 = vars[3];
114 // Determine gap category
115 vector<string> categories = {"inclusive"};
116 if (!nGapJets) { categories += string("gap"); }
117
118 // Fill histograms relevant for comparison with 2010 data
119 if (region == _years[0]) {
120 // Fill inclusive and gap histograms
121 for (const string& cat : categories) {
122 _aux_dy[cat]->fill(dy);
123 _h_dphi_dySlices[cat].fill(dy, dphi/M_PI);
124 _p_cosDeltaPhi_dy[cat]->fill(dy, cos(M_PI - dphi));
125 _p_cosTwoDeltaPhi_dy[cat]->fill(dy, cos(2*dphi));
126 }
127 // Fill profiled nGapJets
128 _p_nGapJets_dy->fill(dy, nGapJets);
129 // Fill Q0 histograms - can fill multiple points per event
130 for (const double Q0 : _Qnoughts) {
131 _aux_Q0_dySlices["inclusive"].fill(dy, Q0);
132 if (maxGapQ0 <= Q0) { _aux_Q0_dySlices["gap"].fill(dy, Q0); }
133 }
134
135 // Fill histograms relevant for comparison with 2011 data
136 }
137 else if (region == _years[1]) {
138 // Fill inclusive and gap histograms
139 for (const string& cat : categories) {
140 _aux_pTbar[cat]->fill(pTbar);
141 _p_cosDeltaPhi_pTbar[cat]->fill(pTbar, cos(M_PI - dphi));
142 _p_cosTwoDeltaPhi_pTbar[cat]->fill(pTbar, cos(2*dphi));
143 }
144 // Fill profiled nGapJets
145 _p_nGapJets_pTbar->fill(pTbar, nGapJets);
146 }
147 return;
148 }
149
150 /// Normalise histograms etc., after the run
151 void finalize() {
152 // Normalise cross-section plots to correct cross-section
153 const double ySpan = 1.0; // all dy spans are 1
154 const double sf = crossSection() / picobarn / sumOfWeights();
155 for (const string cat : { "inclusive", "gap" }) {
156 _h_dphi_dySlices[cat].scale(sf/ySpan/M_PI, this);
157 // Create C2/C1 scatter from profiles
158 divide(_p_cosTwoDeltaPhi_dy[cat], _p_cosDeltaPhi_dy[cat], _h_C2C1_dy[cat]);
159 divide(_p_cosTwoDeltaPhi_pTbar[cat], _p_cosDeltaPhi_pTbar[cat], _h_C2C1_pTbar[cat]);
160 }
161
162 // Fill simple gap fractions
163 efficiency(_aux_dy["gap"], _aux_dy["inclusive"], _s_gapFrac_dy);
164 efficiency(_aux_pTbar["gap"], _aux_pTbar["inclusive"], _s_gapFrac_pTbar);
165
166 // Register and fill Q0 gap fractions
167 for (size_t dyLow = 0; dyLow < _dy_max; ++dyLow) {
168 efficiency(_aux_Q0_dySlices["gap"].histos()[dyLow], _aux_Q0_dySlices["inclusive"].histos()[dyLow], _s_gapFrac_Q0[dyLow]);
169 }
170 }
171
172
173 private:
174
175 /// Member variables
176 vector<size_t> _years;
177 vector<double> _Qnoughts;
178
179 size_t _dy_max;
180
181 /// auxiliary histograms for gap fractions
182 map<string, Histo1DPtr> _aux_dy;
183 map<string, Histo1DPtr> _aux_pTbar;
184 map<string, BinnedHistogram> _aux_Q0_dySlices;
185
186 // Gap fractions
187 Scatter2DPtr _s_gapFrac_dy, _s_gapFrac_pTbar;
188 vector<Scatter2DPtr> _s_gapFrac_Q0;
189
190 // Number of jets in rapidity interval
191 Profile1DPtr _p_nGapJets_dy;
192 Profile1DPtr _p_nGapJets_pTbar;
193
194 // Azimuthal moment histograms
195 map<string, Profile1DPtr> _p_cosDeltaPhi_dy;
196 map<string, Profile1DPtr> _p_cosDeltaPhi_pTbar;
197 map<string, Profile1DPtr> _p_cosTwoDeltaPhi_dy;
198 map<string, Profile1DPtr> _p_cosTwoDeltaPhi_pTbar;
199 map<string, Scatter2DPtr> _h_C2C1_dy;
200 map<string, Scatter2DPtr> _h_C2C1_pTbar;
201
202 // Cross-section vs. deltaPhi in deltaY slices
203 map<string, BinnedHistogram> _h_dphi_dySlices;
204 };
205
206 // The hook for the plugin system
207 RIVET_DECLARE_PLUGIN(ATLAS_2014_I1307243);
208
209}
|