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#include "Rivet/Tools/BinnedHistogram.hh"
6
7namespace Rivet {
8
9 // This analysis is a derived from the class Analysis:
10 class CMS_2013_I1208923 : public Analysis {
11
12 public:
13 // Constructor
14 CMS_2013_I1208923()
15 : Analysis("CMS_2013_I1208923") {
16
17 }
18
19 // Book histograms and initialize projections:
20 void init() {
21 const FinalState fs;
22 declare(fs, "FS");
23
24 // Initialize the projections
25 declare(FastJets(fs, FastJets::ANTIKT, 0.7), "Jets");
26
27 // Book histograms
28 {Histo1DPtr tmp; _h_sigma.add(0.0, 0.5, book(tmp, 1, 1, 1));}
29 {Histo1DPtr tmp; _h_sigma.add(0.5, 1.0, book(tmp, 1, 1, 2));}
30 {Histo1DPtr tmp; _h_sigma.add(1.0, 1.5, book(tmp, 1, 1, 3));}
31 {Histo1DPtr tmp; _h_sigma.add(1.5, 2.0, book(tmp, 1, 1, 4));}
32 {Histo1DPtr tmp; _h_sigma.add(2.0, 2.5, book(tmp, 1, 1, 5));}
33
34 {Histo1DPtr tmp; _h_invMass.add(0.0, 0.5, book(tmp, 2, 1, 1));}
35 {Histo1DPtr tmp; _h_invMass.add(0.5, 1.0, book(tmp, 2, 1, 2));}
36 {Histo1DPtr tmp; _h_invMass.add(1.0, 1.5, book(tmp, 2, 1, 3));}
37 {Histo1DPtr tmp; _h_invMass.add(1.5, 2.0, book(tmp, 2, 1, 4));}
38 {Histo1DPtr tmp; _h_invMass.add(2.0, 2.5, book(tmp, 2, 1, 5));}
39 }
40
41 // Analysis
42 void analyze(const Event &event) {
43 const double weight = 1.0;
44 const FastJets &fJets = apply<FastJets>(event, "Jets");
45
46 // Fill the jet pT spectra
47 const Jets& jets = fJets.jetsByPt(Cuts::pt>100.*GeV && Cuts::absrap <2.5);
48 for (const Jet &j : jets) {
49 _h_sigma.fill(fabs(j.momentum().rapidity()), j.momentum().pT() / GeV, weight);
50 }
51
52 // Require two jets
53 const Jets& dijets = fJets.jetsByPt(Cuts::pt>30.*GeV && Cuts::absrap < 2.5);
54 if (dijets.size() > 1) {
55 if (dijets[0].momentum().pT() / GeV > 60.) {
56 // Fill the invariant mass histogram
57 double ymax = max(dijets[0].momentum().absrapidity(), dijets[1].momentum().absrapidity());
58 double invMass = FourMomentum(dijets[0].momentum() + dijets[1].momentum()).mass();
59 _h_invMass.fill(fabs(ymax), invMass, weight);
60 }
61 }
62
63 }
64
65
66 // Scale histograms by the production cross section
67 void finalize() {
68 _h_sigma.scale( crossSection() / sumOfWeights() / 2.0, this);
69 _h_invMass.scale(crossSection() / sumOfWeights() / 2.0, this);
70 }
71
72 private:
73 BinnedHistogram _h_sigma;
74 BinnedHistogram _h_invMass;
75 };
76
77 // This global object acts as a hook for the plugin system.
78 RIVET_DECLARE_PLUGIN(CMS_2013_I1208923);
79}
|