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#include "Rivet/Tools/BinnedHistogram.hh"
6
7namespace Rivet {
8
9
10 /// @brief ATLAS inclusive b-jet pT spectrum, di-jet mass and di-jet chi
11 class ATLAS_2011_I930220: public Analysis {
12 public:
13
14 ATLAS_2011_I930220()
15 : Analysis("ATLAS_2011_I930220")
16 { }
17
18
19 void init() {
20 FinalState fs((Cuts::etaIn(-3.5, 3.5)));
21 declare(fs, "FinalState");
22 FastJets fj(fs, FastJets::ANTIKT, 0.4);
23 fj.useInvisibles();
24 declare(fj, "Jets");
25 declare(HeavyHadrons(Cuts::abseta < 3.5 && Cuts::pT > 5*GeV), "BHadrons");
26
27 double ybins[] = { 0.0, 0.3, 0.8, 1.2, 2.1 };
28 for (size_t i = 0; i < 4; ++i) {
29 Histo1DPtr tmp;
30 _bjetpT_SV0.add(ybins[i], ybins[i+1], book(tmp, i+1, 1, 1));
31 }
32 book(_bjetpT_SV0_All ,5, 1, 1);
33 book(_bjetpT_pTRel ,6, 1, 1);
34 book(_dijet_mass ,7, 1, 1);
35 book(_dijet_phi ,8, 1, 1);
36 book(_dijet_chi_110_370 ,9, 1, 1);
37 book(_dijet_chi_370_850 ,10, 1, 1);
38
39 book(_chiCounter1, "_chiCounter1");
40 book(_chiCounter2, "_chiCounter2");
41 book(_phiCounter, "_phiCounter1");
42 }
43
44
45 void analyze(const Event& evt) {
46
47 const Particles& bHadrons = apply<HeavyHadrons>(evt, "BHadrons").bHadrons();
48 const Jets& jets = apply<JetAlg>(evt, "Jets").jetsByPt(15*GeV);
49
50 FourMomentum leadingJet, subleadingJet;
51 int leadJet = 0, subJet = 0;
52 for (const Jet& j : jets) {
53 bool hasB = false;
54 for (const Particle& b : bHadrons)
55 if (deltaR(j, b) < 0.3) { hasB = true; break; }
56
57 // Identify and classify the leading and subleading jets
58 if (j.absrap() < 2.1) { ///< Move this into the jets defn
59 if (!leadJet) {
60 leadingJet = j.momentum();
61 leadJet = (hasB && j.pT() > 40*GeV) ? 2 : 1;
62 }
63 else if (leadJet && !subJet) {
64 subleadingJet = j.momentum();
65 subJet = (hasB && j.pT() > 40*GeV) ? 2 : 1;
66 }
67 if (hasB) {
68 _bjetpT_SV0.fill(j.absrap(), j.pT()/GeV);
69 _bjetpT_SV0_All->fill(j.pT()/GeV);
70 _bjetpT_pTRel->fill(j.pT()/GeV);
71 }
72 }
73 }
74
75 // Di-b-jet plots require both the leading and subleading jets to be b-tagged and have pT > 40 GeV
76 if (leadJet == 2 && subJet == 2) {
77 const double mass = FourMomentum( leadingJet + subleadingJet ).mass();
78 _dijet_mass->fill(mass/GeV);
79
80 // Plot dphi for high-mass di-b-jets
81 if (mass > 110*GeV) {
82 _phiCounter->fill();
83 const double d_phi = deltaPhi( leadingJet.phi(), subleadingJet.phi() );
84 _dijet_phi->fill(fabs(d_phi));
85 }
86
87 // Plot chi for low y_boost di-b-jets (in two high-mass bins)
88 const double y_boost = 0.5 * (leadingJet.rapidity() + subleadingJet.rapidity());
89 const double chi = exp( fabs( leadingJet.rapidity() - subleadingJet.rapidity() ) );
90 if ( fabs(y_boost) < 1.1 ) {
91 if (inRange(mass/GeV, 110, 370)) {
92 _chiCounter1->fill();
93 _dijet_chi_110_370->fill(chi);
94 } else if (inRange(mass/GeV, 370, 850)) {
95 _chiCounter2->fill();
96 _dijet_chi_370_850->fill(chi);
97 }
98 }
99 }
100 }
101
102
103 void finalize() {
104 // Normalizing to cross-section and mass
105 // Additional factors represent the division by rapidity
106 const double xsec = crossSectionPerEvent()/(picobarn);
107 const double chiScale1 = 1 / dbl(*_chiCounter1) / 260.0;
108 const double chiScale2 = 1 / dbl(*_chiCounter2) / 480.0;
109 const double phiScale = 1 / dbl(*_phiCounter);
110
111 _bjetpT_SV0.scale(xsec/2, this);
112 scale(_bjetpT_SV0_All, xsec);
113 scale(_bjetpT_pTRel, xsec);
114 scale(_dijet_mass, xsec);
115 scale(_dijet_phi, phiScale );
116 scale(_dijet_chi_110_370, chiScale1);
117 scale(_dijet_chi_370_850, chiScale2);
118 }
119
120
121 private:
122
123 BinnedHistogram _bjetpT_SV0;
124
125 Histo1DPtr _bjetpT_SV0_All;
126 Histo1DPtr _bjetpT_pTRel;
127 Histo1DPtr _dijet_mass;
128 Histo1DPtr _dijet_phi;
129 Histo1DPtr _dijet_chi_110_370;
130 Histo1DPtr _dijet_chi_370_850;
131
132 CounterPtr _chiCounter1;
133 CounterPtr _chiCounter2;
134 CounterPtr _phiCounter;
135 };
136
137
138 // The hook for the plugin system
139 RIVET_DECLARE_PLUGIN(ATLAS_2011_I930220);
140
141}
|