Rivet analyses referenceATLAS_2018_I1634970ATLAS Inclusive jet and dijet cross section measurement at sqrt(s)=13TeVExperiment: ATLAS (LHC) Inspire ID: 1634970 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
Inclusive jet and dijet cross-sections are measured in proton-proton collisions at a centre-of-mass energy of 13 TeV. The measurement uses a dataset with an integrated luminosity of 3.2 fb$^{-1}$ recorded in 2015 with the ATLAS detector at the Large Hadron Collider. Jets are identified using the anti-$k_\text{t}$ algorithm with a radius parameter value of $R = 0.4$. The inclusive jet cross-sections are measured double-differentially as a function of the jet transverse momentum, covering the range from 100 GeV to 3.5 TeV, and the absolute jet rapidity up to $|y| = 3$. The double-differential dijet production cross-sections are presented as a function of the dijet mass, covering the range from 300 GeV to 9 TeV, and the half absolute rapidity separation between the two leading jets within $|y| < 3$, $y^\ast$, up to $y^\ast = 3$. Next-to-leading-order, and next-to-next-to-leading-order for the inclusive jet measurement, perturbative QCD calculations corrected for non-perturbative and electroweak effects are compared to the measured cross-sections. Source code: ATLAS_2018_I1634970.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/FastJets.hh"
5
6namespace Rivet {
7
8
9 /// @brief inclusive jet and dijet at 13 TeV
10 class ATLAS_2018_I1634970 : public Analysis {
11 public:
12
13 /// Constructor
14 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2018_I1634970);
15
16 /// Book histograms and initialise projections before the run
17 void init() {
18
19 const FinalState fs;
20 declare(fs,"FinalState");
21 FastJets fj04(fs, JetAlg::ANTIKT, 0.4, JetMuons::ALL, JetInvisibles::DECAY);
22 declare(fj04, "AntiKT04");
23
24 // |y| and ystar bins
25 vector<double> ybins{ 0.0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0};
26
27 // Book histograms
28 book(_pThistograms, ybins);
29 book(_mjjhistograms, ybins);
30 for (size_t i=1; i < _pThistograms->numBins()+1; ++i) {
31 book(_pThistograms->bin(i), i, 1, 1);
32 book(_mjjhistograms->bin(i), 6+i, 1, 1);
33 }
34 }
35
36
37 /// Perform the per-event analysis
38 void analyze(const Event& event) {
39
40 const Jets& kt4Jets = apply<FastJets>(event, "AntiKT04").jetsByPt(Cuts::pT > 75*GeV && Cuts::absrap < 3.0);
41
42 int nJets = kt4Jets.size();
43
44 // Inclusive jet selection
45 for(int ijet=0;ijet<nJets;++ijet){ // loop over jets
46 FourMomentum jet = kt4Jets[ijet].momentum();
47 // pT selection
48 if (jet.pt()>100.0*GeV){
49 // Fill distribution
50 const double absy = jet.absrap();
51 _pThistograms->fill(absy,jet.pt()/GeV);
52 }
53 }
54
55 // Dijet selection
56 if(nJets > 1){ // skip events with less than 2 jets passing pT>75GeV and |y|<3.0 cuts
57 FourMomentum jet0 = kt4Jets[0].momentum();
58 FourMomentum jet1 = kt4Jets[1].momentum();
59 const double rap0 = jet0.rapidity();
60 const double rap1 = jet1.rapidity();
61 const double ystar = fabs(rap0-rap1)/2;
62 const double mass = (jet0 + jet1).mass();
63 const double HT2 = jet0.pt()+jet1.pt();
64 if (HT2>200*GeV && ystar<3.0){
65 // Fill distribution
66 _mjjhistograms->fill(ystar, mass/GeV);
67 }
68 }
69
70 }
71
72
73 /// Normalise histograms etc., after the run
74 void finalize() {
75
76 const double xs_norm_factor = 0.5*crossSection() / picobarn / sumOfWeights();
77
78 scale(_pThistograms, xs_norm_factor);
79 scale(_mjjhistograms, crossSectionPerEvent()/picobarn);
80 divByGroupWidth({_pThistograms, _mjjhistograms});
81
82 }
83
84 private:
85
86 // The inclusive pT spectrum for akt4 jets
87 Histo1DGroupPtr _pThistograms;
88 // The dijet mass spectrum for akt4 jets
89 Histo1DGroupPtr _mjjhistograms;
90
91 };
92
93 RIVET_DECLARE_PLUGIN(ATLAS_2018_I1634970);
94
95}
|