Rivet analyses referenceCMS_2012_I1089835Inclusive b-jet production in pp collisions at 7 TeVExperiment: CMS (LHC) Inspire ID: 1089835 Status: VALIDATED Authors:
Beam energies: (3500.0, 3500.0) GeV Run details:
The inclusive b-jet production cross section in pp collisions at a center-of-mass energy of 7 TeV is measured using data collected by the CMS experiment at the LHC. The cross section is presented as a function of the jet transverse momentum in the range 18 $<$ pT $<$ 200 GeV for several rapidity intervals. The results are also given as the ratio of the b-jet production cross section to the inclusive jet production cross section. The measurement is performed with two different analyses, which differ in their trigger selection and b-jet identification: a jet analysis that selects events with a b jet using a sample corresponding to an integrated luminosity of 34 inverse picobarns, and a muon analysis requiring a b jet with a muon based on an integrated luminosity of 3 inverse picobarns. In both approaches the b jets are identified by requiring a secondary vertex. The results from the two methods are in agreement with each other and with next-to-leading order calculations, as well as with predictions based on the PYTHIA event generator. Source code: CMS_2012_I1089835.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Tools/Logging.hh"
4#include "Rivet/Projections/FinalState.hh"
5#include "Rivet/Projections/FastJets.hh"
6
7namespace Rivet {
8
9 /// @brief Inclusive b-jet production in pp collisions at 7 TeV
10 class CMS_2012_I1089835 : public Analysis {
11 public:
12
13 /// @name Constructors etc.
14 //@{
15
16 /// Constructor
17 CMS_2012_I1089835()
18 : Analysis("CMS_2012_I1089835")
19 {
20
21 }
22
23 //@}
24
25
26 public:
27
28 /// @name Analysis methods
29 //@{
30
31 /// Book histograms and initialise projections before the run
32 void init() {
33
34 const FinalState cnfs;
35 declare(cnfs, "FS");
36 declare(FastJets(cnfs, FastJets::ANTIKT, 0.5), "Jets");
37
38 book(_h_dsigdpty05, 4, 1, 1);
39 book(_h_dsigdpty10, 5, 1, 1);
40 book(_h_dsigdpty15, 6, 1, 1);
41 book(_h_dsigdpty20, 7, 1, 1);
42 book(_h_dsigdpty22, 8, 1, 1);
43
44 }
45
46
47 /// Perform the per-event analysis
48 void analyze(const Event& event) {
49
50 const FastJets& fastjets = apply<FastJets>(event, "Jets");
51 const Jets jets = fastjets.jetsByPt(10.);
52
53 for (const Jet& j : jets) {
54
55 const double ptB= j.pT();
56 const double yB= j.rapidity();
57
58 if (j.bTagged()) {
59
60 if( fabs(yB) < 0.5) { _h_dsigdpty05->fill( ptB, 1.0 );}
61 else if( fabs(yB) >= 0.5 && fabs(yB) < 1.0) { _h_dsigdpty10->fill( ptB, 1.0 );}
62 else if( fabs(yB) >= 1.0 && fabs(yB) < 1.5) { _h_dsigdpty15->fill( ptB, 1.0 );}
63 else if( fabs(yB) >= 1.5 && fabs(yB) < 2.0) { _h_dsigdpty20->fill( ptB, 1.0 );}
64 else if( fabs(yB) >= 2.0 && fabs(yB) < 2.2) { _h_dsigdpty22->fill( ptB, 1.0 );}
65 }
66 }
67 }
68
69 /// Normalise histograms etc., after the run
70 void finalize() {
71
72 double invlumi = crossSection()/picobarn/sumOfWeights();
73
74 scale(_h_dsigdpty05, invlumi);
75 scale(_h_dsigdpty10, invlumi);
76 scale(_h_dsigdpty15, invlumi);
77 scale(_h_dsigdpty20, invlumi);
78 scale(_h_dsigdpty22, invlumi/0.4);
79
80 }
81
82 //@}
83
84
85 private:
86
87 /// @name Histograms
88 //@{
89 Histo1DPtr _h_dsigdpty05;
90 Histo1DPtr _h_dsigdpty10;
91 Histo1DPtr _h_dsigdpty15;
92 Histo1DPtr _h_dsigdpty20;
93 Histo1DPtr _h_dsigdpty22;
94 //@}
95
96
97 };
98
99
100
101 // The hook for the plugin system
102 DECLARE_RIVET_PLUGIN(CMS_2012_I1089835);
103
104}
|