rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2016_I1478355

Measurement of the bbar dijet dijet cross section at 7 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1478355
Status: VALIDATED
Authors:
  • Lucia Keszeghova
  • Deepak Kar
References:
  • Eur.Phys.J. C76 (2016) 670
  • DOI:10.1140/epjc/s10052-016-4521-y
  • arXiv: 1607.08430
Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • Dijet production, high pThatmin

The dijet production cross section for jets containing a $b$-hadron ($b$-jets) has been measured in proton-proton collisions with a centre-of-mass energy of $\sqrt{s} = 7$ TeV, using the ATLAS detector at the LHC. The data used correspond to an integrated luminosity of 4.2 fb$^{-1}$. The cross section is measured for events with two identified $b$-jets with a transverse momentum $p_\mathrm{T} > 20$GeV and a minimum separation in the $\eta$-$\phi$ plane of $\Delta R$ = 0.4. At least one of the jets in the event is required to have $p_\mathrm{T} > 270$GeV. The cross section is measured differentially as a function of dijet invariant mass, dijet transverse momentum, boost of the dijet system, and the rapidity difference, azimuthal angle and angular distance between the $b$-jets. The results are compared to different predictions of leading order and next-to-leading order perturbative quantum chromodynamics matrix elements supplemented with models for parton-showers and hadronization.

Source code: ATLAS_2016_I1478355.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 Measurement of the bbar dijet dijet cross section at 7 TeV
 11  class ATLAS_2016_I1478355 : public Analysis {
 12  public:
 13
 14    /// Constructor
 15    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2016_I1478355);
 16
 17    /// @name Analysis methods
 18    ///@{
 19
 20    /// Book histograms and initialise projections before the run
 21    void init() {
 22
 23      // Initialise and register projections
 24      FinalState fs(Cuts::abseta < 3.2);
 25      FastJets fj(fs, FastJets::ANTIKT, 0.4);
 26      fj.useInvisibles();
 27      declare(fj, "Jets");
 28      declare(HeavyHadrons(Cuts::abseta < 3.2 && Cuts::pT > 5*GeV), "BHadrons");
 29
 30      book(_h["m_bb"], 1, 1, 1);
 31      book(_h["Delta_phi"], 2, 1, 1);
 32      book(_h["y_diff"], 3, 1, 1);
 33      book(_h["Delta_R"], 4, 1, 1);
 34      book(_h["pT_bb"], 5, 1, 1);
 35      book(_h["y_B"], 6, 1, 1);
 36    }
 37
 38
 39    /// Perform the per-event analysis
 40    void analyze(const Event& event) {
 41
 42      // Retrieve clustered jets, sorted by pT, with a minimum pT cut
 43      const Jets& jets = apply<FastJets>(event, "Jets").jetsByPt(20*GeV);
 44      const Particles& bHadrons = apply<HeavyHadrons>(event, "BHadrons").bHadrons();
 45
 46      if (jets.size() < 1)  vetoEvent;
 47      if (jets[0].pT() > 270*GeV){
 48
 49        Jets foundBjets;
 50        for (const Jet& j : jets){
 51          for (const Particle& b : bHadrons){
 52            if ((deltaR(j, b) < 0.3) && (j.pT() > 20*GeV) && (fabs(j.eta()) < 2.5)) foundBjets.push_back(j);
 53          }
 54        }
 55
 56        Jet firstBjet;
 57        Jet secondBjet;
 58        bool bJetPair = false;
 59        for (const Jet& b1 : foundBjets){
 60          for (const Jet& b2 : foundBjets){
 61            if(deltaR(b1, b2) > 0.4){
 62              bJetPair = true;
 63              firstBjet = b1;
 64              secondBjet = b2;
 65              break;
 66            }
 67          }
 68          if (bJetPair) break;
 69        }
 70
 71        if (bJetPair) {
 72          const double mass = (firstBjet.momentum() + secondBjet.momentum()).mass();
 73          _h["m_bb"]->fill(mass/GeV);
 74          const double Delta_phi = deltaPhi(firstBjet.phi(), secondBjet.phi());
 75          _h["Delta_phi"]->fill(Delta_phi);
 76          const double y_diff = fabs(firstBjet.rapidity() - secondBjet.rapidity())/2;
 77          _h["y_diff"]->fill(y_diff);
 78          const double Delta_R = deltaR(firstBjet.momentum(), secondBjet.momentum());
 79          _h["Delta_R"]->fill(Delta_R);
 80          const double pT = (firstBjet.momentum() + secondBjet.momentum()).pT();
 81          _h["pT_bb"]->fill(pT/GeV);
 82          const double y_B = (firstBjet.rapidity() + secondBjet.rapidity())/2;
 83          _h["y_B"]->fill(y_B);
 84        }
 85      }
 86    }
 87
 88    /// Normalise histograms etc., after the run
 89    void finalize() {
 90
 91      const double xsec = crossSectionPerEvent()/picobarn; 
 92      scale(_h, xsec);
 93    }
 94
 95    ///@}
 96
 97    /// @name Histograms
 98    ///@{
 99     map<string, Histo1DPtr> _h;
100    ///@}
101
102  };
103
104  DECLARE_RIVET_PLUGIN(ATLAS_2016_I1478355);
105}