Rivet analyses referenceCMS_2013_I1208923Jet-pT and dijet mass at sqrt(s) = 7 TeVExperiment: CMS (LHC) Inspire ID: 1208923 Status: VALIDATED Authors:
Beams: p+ p+ Beam energies: (3500.0, 3500.0) GeV Run details:
The single jet differential cross section has been measured as a function of the jet momentum and the dijet differential cross section as a function of the dijet mass. The data has been recorded by the CMS detector at the center of mass energy of 7 TeV. To reconstruct the jets, the anti-$k_T$ algorithm with a cone radius of $R = 0.7$ (AK7) has been used. The results are split into five rapidity bins, ranging from 0.0 to 2.5 with a width of 0.5 each. Source code: CMS_2013_I1208923.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 // This analysis is a derived from the class Analysis:
9 class CMS_2013_I1208923 : public Analysis {
10
11 public:
12
13 // Constructor
14 RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2013_I1208923);
15
16 // Book histograms and initialize projections:
17 void init() {
18 const FinalState fs;
19 declare(fs, "FS");
20
21 // Initialize the projections
22 declare(FastJets(fs, JetAlg::ANTIKT, 0.7), "Jets");
23
24 // Book histograms
25 book(_h_sigma, {0., 0.5, 1., 1.5, 2., 2.5},
26 {"d01-x01-y01", "d01-x01-y02", "d01-x01-y03", "d01-x01-y04", "d01-x01-y05"});
27 book(_h_invMass, {0., 0.5, 1., 1.5, 2., 2.5},
28 {"d02-x01-y01", "d02-x01-y02", "d02-x01-y03", "d02-x01-y04", "d02-x01-y05"});
29
30 }
31
32 // Analysis
33 void analyze(const Event &event) {
34
35 const FastJets &fJets = apply<FastJets>(event, "Jets");
36
37 // Fill the jet pT spectra
38 const Jets& jets = fJets.jetsByPt(Cuts::pt>100.*GeV && Cuts::absrap <2.5);
39 for (const Jet& j : jets) {
40 _h_sigma->fill(fabs(j.rapidity()), j.pT() / GeV);
41 }
42
43 // Require two jets
44 const Jets& dijets = fJets.jetsByPt(Cuts::pt>30.*GeV && Cuts::absrap < 2.5);
45 if (dijets.size() > 1) {
46 if (dijets[0].momentum().pT() / GeV > 60.) {
47 // Fill the invariant mass histogram
48 double ymax = max(dijets[0].momentum().absrapidity(), dijets[1].momentum().absrapidity());
49 double invMass = FourMomentum(dijets[0].momentum() + dijets[1].momentum()).mass();
50 _h_invMass->fill(fabs(ymax), invMass);
51 }
52 }
53
54 }
55
56
57 // Scale histograms by the production cross section
58 void finalize() {
59 scale(_h_sigma, crossSection()/picobarn / sumOfWeights() / 2.0);
60 scale(_h_invMass, crossSection()/picobarn / sumOfWeights() / 2.0);
61 divByGroupWidth({_h_sigma, _h_invMass});
62 }
63
64 private:
65
66 Histo1DGroupPtr _h_sigma;
67 Histo1DGroupPtr _h_invMass;
68
69 };
70
71 RIVET_DECLARE_PLUGIN(CMS_2013_I1208923);
72}
|