Rivet analyses

title: CMS_2016_I1486238

Studies of 2 $b$-jet + 2 jet production in proton-proton collisions at 7 TeV

Experiment: CMS (LHC)

Inspire ID: 1486238

Status: VALIDATED

Authors: - P. Gunnellini

References: - Expt page: CMS-FSQ-13-010 - CERN-EP-2016-191 - arXiv: 1609.03489 - Phys.Rev. D94 (2016) no.11, 112005

Beams: p+ p+

Beam energies: (3500.0, 3500.0)GeV

Run details: - Hard QCD events with \pT cut at generator level of 15 GeV. Needs 2 $b$-jets with $\pT > 20$ GeV in $|\eta|$ $<$ 2.4; 2 jets with $\pT > 20$ GeV in $|\eta|$ $<$ 4.7

Measurements are presented of the cross section for the production of at least four jets, of which at least two originate from $b$ quarks, in proton-proton collisions. Data collected with the CMS detector at the LHC at a center-of-mass energy of 7 TeV are used, corresponding to an integrated luminosity of 3 inverse picobarns. The cross section is measured as a function of the jet transverse momentum for \pT 20 GeV, and of the jet pseudorapidity for $|\eta|<2.4$ ($b$ jets) and 4.7 (untagged jets). The correlations in azimuthal angle and \pT between the jets are also studied. The $\eta$ and \pT distributions of the four jets and the correlations between them are well reproduced by event generators that combine perturbative QCD calculations at next-to-leading-order accuracy with contributions from parton showers and multiparton interactions.

Source code:CMS_2016_I1486238.cc

```c++ // -- C++ --

include "Rivet/Analysis.hh"

include "Rivet/Projections/FastJets.hh"

include "Rivet/Projections/FinalState.hh"

define I_KNOW_THE_INITIAL_QUARKS_PROJECTION_IS_DODGY_BUT_NEED_TO_USE_IT

include "Rivet/Projections/InitialQuarks.hh"

namespace Rivet {

/// Studies of 2 b-jet + 2 jet production in proton-proton collisions at 7 TeV class CMS_2016_I1486238 : public Analysis { public:

/// Constructor
RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2016_I1486238);


/// @name Analysis methods
/// @{

/// Book histograms and initialise projections before the run
void init() {

  FastJets akt(FinalState(), JetAlg::ANTIKT, 0.5);
  declare(akt, "antikT");

  book(_h_Deltaphi_newway, 1, 1, 1);
  book(_h_deltaphiafterlight, 9, 1, 1);
  book(_h_SumPLight, 5, 1, 1);

  book(_h_LeadingBJetpt, 11, 1, 1);
  book(_h_SubleadingBJetpt, 15, 1, 1);
  book(_h_LeadingLightJetpt, 13, 1, 1);
  book(_h_SubleadingLightJetpt, 17, 1, 1);

  book(_h_LeadingBJeteta, 10, 1, 1);
  book(_h_SubleadingBJeteta, 14, 1, 1);
  book(_h_LeadingLightJeteta, 12, 1, 1);
  book(_h_SubleadingLightJeteta, 16, 1, 1);
}


/// Perform the per-event analysis
void analyze(const Event& event) {
  const Jets& jets = apply<JetFinder>(event, "antikT")
                         .jetsByPt(Cuts::absrap < 4.7 && Cuts::pT > 20 * GeV);
  if (jets.size() < 4) vetoEvent;

  // Initial quarks
  /// @note Quark-level tagging...
  Particles bquarks;
  for (ConstGenParticlePtr p : HepMCUtils::particles(event.genEvent())) {
    if (abs(p->pdg_id()) == PID::BQUARK) bquarks += Particle(p);
  }
  Jets bjets, ljets;
  for (const Jet& j : jets) {
    const bool btag = any(bquarks, deltaRLess(j, 0.3));
    // for (const Particle& b : bquarks) if (deltaR(j, b) < 0.3) btag = true;
    (btag && j.abseta() < 2.4 ? bjets : ljets).push_back(j);
  }

  // Fill histograms
  if (bjets.size() >= 2 && ljets.size() >= 2) {
    _h_LeadingBJetpt->fill(bjets[0].pT() / GeV);
    _h_SubleadingBJetpt->fill(bjets[1].pT() / GeV);
    _h_LeadingLightJetpt->fill(ljets[0].pT() / GeV);
    _h_SubleadingLightJetpt->fill(ljets[1].pT() / GeV);
    //
    _h_LeadingBJeteta->fill(bjets[0].eta());
    _h_SubleadingBJeteta->fill(bjets[1].eta());
    _h_LeadingLightJeteta->fill(ljets[0].eta());
    _h_SubleadingLightJeteta->fill(ljets[1].eta());

    const double lightdphi = deltaPhi(ljets[0], ljets[1]);
    _h_deltaphiafterlight->fill(lightdphi);

    const double vecsumlightjets = sqrt(
        sqr(ljets[0].px() + ljets[1].px())
        + sqr(ljets[0].py() + ljets[1].py())); //< @todo Just (lj0+lj1).pT()? Or use add_quad
    const double term2 = vecsumlightjets
        / (sqrt(sqr(ljets[0].px()) + sqr(ljets[0].py()))
           + sqrt(sqr(ljets[1].px()) + sqr(ljets[1].py()))); //< @todo lj0.pT() + lj1.pT()? Or add_quad
    _h_SumPLight->fill(term2);

    const double pxBsyst2 = bjets[0].px() + bjets[1].px();                 // @todo (bj0+bj1).px()
    const double pyBsyst2 = bjets[0].py() + bjets[1].py();                 // @todo (bj0+bj1).py()
    const double pxJetssyst2 = ljets[0].px() + ljets[1].px();              // @todo (lj0+lj1).px()
    const double pyJetssyst2 = ljets[0].py() + ljets[1].py();              // @todo (lj0+lj1).py()
    const double modulusB2 = sqrt(sqr(pxBsyst2) + sqr(pyBsyst2));          //< @todo add_quad
    const double modulusJets2 = sqrt(sqr(pxJetssyst2) + sqr(pyJetssyst2)); //< @todo add_quad
    const double cosphiBsyst2 = pxBsyst2 / modulusB2;
    const double cosphiJetssyst2 = pxJetssyst2 / modulusJets2;
    const double phiBsyst2 = ((pyBsyst2 > 0) ? 1 : -1) * acos(cosphiBsyst2); //< @todo sign(pyBsyst2)
    const double phiJetssyst2 = sign(pyJetssyst2) * acos(cosphiJetssyst2);
    const double Dphi2 = deltaPhi(phiBsyst2, phiJetssyst2);
    _h_Deltaphi_newway->fill(Dphi2);
  }
}


/// Normalise histograms etc., after the run
void finalize() {
  const double invlumi = crossSection() / picobarn / sumOfWeights();
  normalize(_h_SumPLight);
  normalize(_h_deltaphiafterlight);
  normalize(_h_Deltaphi_newway);
  scale(_h_LeadingLightJetpt, invlumi);
  scale(_h_SubleadingLightJetpt, invlumi);
  scale(_h_LeadingBJetpt, invlumi);
  scale(_h_SubleadingBJetpt, invlumi);
  scale(_h_LeadingLightJeteta, invlumi);
  scale(_h_SubleadingLightJeteta, invlumi);
  scale(_h_LeadingBJeteta, invlumi);
  scale(_h_SubleadingBJeteta, invlumi);
}

/// @}

private:

/// @name Histograms
/// @{

Histo1DPtr _h_deltaphiafterlight, _h_Deltaphi_newway, _h_SumPLight;
Histo1DPtr _h_LeadingBJetpt, _h_SubleadingBJetpt, _h_LeadingLightJetpt, _h_SubleadingLightJetpt;
Histo1DPtr _h_LeadingBJeteta, _h_SubleadingBJeteta, _h_LeadingLightJeteta, _h_SubleadingLightJeteta;

};

RIVET_DECLARE_PLUGIN(CMS_2016_I1486238);

} ```