Rivet analyses referenceATLAS_2011_I930220Inclusive and dijet cross-sections of $b$-jets in $pp$ collisions at $\sqrt{s} = 7$ TeV with ATLASExperiment: ATLAS (LHC) Inspire ID: 930220 Status: VALIDATED Authors:
Beam energies: (3500.0, 3500.0) GeV Run details:
The inclusive and dijet production cross-sections have been measured for jets containing $b$-hadrons ($b$-jets) in proton--proton collisions at a centre-of-mass energy of $\sqrt{s} = 7$ TeV, using the ATLAS detector at the LHC. The measurements use data corresponding to an integrated luminosity of 34 pb$^-1$. The $b$-jets are identified using either a lifetime-based method, where secondary decay vertices of $b$-hadrons in jets are reconstructed using information from the tracking detectors, or a muon-based method where the presence of a muon is used to identify semileptonic decays of $b$-hadrons inside jets. The inclusive $b$-jet cross-section is measured as a function of transverse momentum in the range $20 < pT < 400$ GeV and rapidity in the range $|y| < 2.1$. The $b\bar{b}$-dijet cross-section is measured as a function of the dijet invariant mass in the range $110 < $m_{jj}$ < 760$ GeV, the azimuthal angle difference between the two jets and the angular variable $\chi$ in two dijet mass regions. The results are compared with next-to-leading-order QCD predictions. Good agreement is observed between the measured cross-sections and the predictions obtained using POWHEG+Pythia6. MC@NLO+Herwig shows good agreement with the measured $b\bar{b}$-dijet cross-section. However, it does not reproduce the measured inclusive cross-section well, particularly for central $b$-jets with large transverse momenta. Source code: ATLAS_2011_I930220.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FastJets.hh"
4#include "Rivet/Projections/HeavyHadrons.hh"
5
6namespace Rivet {
7
8
9 /// @brief ATLAS inclusive b-jet pT spectrum, di-jet mass and di-jet chi
10 class ATLAS_2011_I930220: public Analysis {
11 public:
12
13 ATLAS_2011_I930220()
14 : Analysis("ATLAS_2011_I930220")
15 { }
16
17
18 void init() {
19 FinalState fs((Cuts::etaIn(-3.5, 3.5)));
20 declare(fs, "FinalState");
21 FastJets fj(fs, JetAlg::ANTIKT, 0.4);
22 fj.useInvisibles();
23 declare(fj, "Jets");
24 declare(HeavyHadrons(Cuts::abseta < 3.5 && Cuts::pT > 5*GeV), "BHadrons");
25
26 vector<double> ybins{ 0.0, 0.3, 0.8, 1.2, 2.1 };
27 book(_bjetpT_SV0, ybins, {"d01-x01-y01", "d02-x01-y01", "d03-x01-y01", "d04-x01-y01" });
28 book(_bjetpT_SV0_All , 5, 1, 1);
29 book(_bjetpT_pTRel , 6, 1, 1);
30 book(_dijet_mass , 7, 1, 1);
31 book(_dijet_phi , 8, 1, 1);
32 book(_dijet_chi_110_370, 9, 1, 1);
33 book(_dijet_chi_370_850, 10, 1, 1);
34
35 book(_chiCounter1, "_chiCounter1");
36 book(_chiCounter2, "_chiCounter2");
37 book(_phiCounter, "_phiCounter1");
38 }
39
40
41 void analyze(const Event& evt) {
42
43 const Particles& bHadrons = apply<HeavyHadrons>(evt, "BHadrons").bHadrons();
44 const Jets& jets = apply<JetFinder>(evt, "Jets").jetsByPt(Cuts::pT > 15*GeV);
45
46 FourMomentum leadingJet, subleadingJet;
47 int leadJet = 0, subJet = 0;
48 for (const Jet& j : jets) {
49 bool hasB = false;
50 for (const Particle& b : bHadrons)
51 if (deltaR(j, b) < 0.3) { hasB = true; break; }
52
53 // Identify and classify the leading and subleading jets
54 if (j.absrap() < 2.1) { ///< Move this into the jets defn
55 if (!leadJet) {
56 leadingJet = j.momentum();
57 leadJet = (hasB && j.pT() > 40*GeV) ? 2 : 1;
58 }
59 else if (leadJet && !subJet) {
60 subleadingJet = j.momentum();
61 subJet = (hasB && j.pT() > 40*GeV) ? 2 : 1;
62 }
63 if (hasB) {
64 _bjetpT_SV0->fill(j.absrap(), j.pT()/GeV);
65 _bjetpT_SV0_All->fill(j.pT()/GeV);
66 _bjetpT_pTRel->fill(j.pT()/GeV);
67 }
68 }
69 }
70
71 // Di-b-jet plots require both the leading and subleading jets to be b-tagged and have pT > 40 GeV
72 if (leadJet == 2 && subJet == 2) {
73 const double mass = FourMomentum( leadingJet + subleadingJet ).mass();
74 _dijet_mass->fill(mass/GeV);
75
76 // Plot dphi for high-mass di-b-jets
77 if (mass > 110*GeV) {
78 _phiCounter->fill();
79 const double d_phi = deltaPhi( leadingJet.phi(), subleadingJet.phi() );
80 _dijet_phi->fill(fabs(d_phi));
81 }
82
83 // Plot chi for low y_boost di-b-jets (in two high-mass bins)
84 const double y_boost = 0.5 * (leadingJet.rapidity() + subleadingJet.rapidity());
85 const double chi = exp( fabs( leadingJet.rapidity() - subleadingJet.rapidity() ) );
86 if ( fabs(y_boost) < 1.1 ) {
87 if (inRange(mass/GeV, 110, 370)) {
88 _chiCounter1->fill();
89 _dijet_chi_110_370->fill(chi);
90 } else if (inRange(mass/GeV, 370, 850)) {
91 _chiCounter2->fill();
92 _dijet_chi_370_850->fill(chi);
93 }
94 }
95 }
96 }
97
98
99 void finalize() {
100 // Normalizing to cross-section and mass
101 // Additional factors represent the division by rapidity
102 const double xsec = crossSectionPerEvent()/(picobarn);
103 const double chiScale1 = 1 / dbl(*_chiCounter1) / 260.0;
104 const double chiScale2 = 1 / dbl(*_chiCounter2) / 480.0;
105 const double phiScale = 1 / dbl(*_phiCounter);
106
107 scale(_bjetpT_SV0, 0.5*xsec);
108 scale(_bjetpT_SV0_All, xsec);
109 scale(_bjetpT_pTRel, xsec);
110 scale(_dijet_mass, xsec);
111 scale(_dijet_phi, phiScale );
112 scale(_dijet_chi_110_370, chiScale1);
113 scale(_dijet_chi_370_850, chiScale2);
114
115 divByGroupWidth(_bjetpT_SV0);
116 }
117
118
119 private:
120
121 Histo1DGroupPtr _bjetpT_SV0;
122
123 Histo1DPtr _bjetpT_SV0_All;
124 Histo1DPtr _bjetpT_pTRel;
125 Histo1DPtr _dijet_mass;
126 Histo1DPtr _dijet_phi;
127 Histo1DPtr _dijet_chi_110_370;
128 Histo1DPtr _dijet_chi_370_850;
129
130 CounterPtr _chiCounter1;
131 CounterPtr _chiCounter2;
132 CounterPtr _phiCounter;
133 };
134
135
136 RIVET_DECLARE_PLUGIN(ATLAS_2011_I930220);
137
138}
|