rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2020_I1788444

Differential cross-sections for Z + b-jets
Experiment: ATLAS (LHC)
Inspire ID: 1788444
Status: VALIDATED
Authors:
  • Federico Sforza
  • Deepak Kar
References: Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • Z+jets at 13 TeV

This paper presents a measurement of the production cross-section of a $Z$ boson in association with $b$-jets, in proton-proton collisions at $\sqrt{s} = 13$ TeV with the ATLAS experiment at the Large Hadron Collider using data corresponding to an integrated luminosity of 35.6 fb$^{-1}$. Inclusive and differential cross-sections are measured for events containing a $Z$ boson decaying into electrons or muons and produced in association with at least one or at least two $b$-jets with transverse momentum $p_T> 20$ GeV and rapidity $|y|<2.5$. Predictions from several Monte Carlo generators based on leading-order (LO) or next-to-leading-order (NLO) matrix elements interfaced with a parton-shower simulation and testing different flavour schemes for the choice of initial-state partons are compared with measured cross-sections. The 5-flavour number scheme predictions at NLO accuracy agree better with data than 4-flavour number scheme ones. The 4-flavour number scheme predictions underestimate data in events with at least one $b$-jet. N.B.: Data correspond to average of electron and muon channel. Use LMODE=EL,MU to run on individual channels.

Source code: ATLAS_2020_I1788444.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/ZFinder.hh"
  4#include "Rivet/Projections/FastJets.hh"
  5#include "Rivet/Projections/VetoedFinalState.hh"
  6#include "Rivet/Projections/HeavyHadrons.hh"
  7#include "Rivet/Projections/PromptFinalState.hh"
  8
  9namespace Rivet {
 10
 11  /// @brief Z + b(b) in pp at 13 TeV
 12  class ATLAS_2020_I1788444 : public Analysis {
 13  public:
 14    /// @name Constructors etc.
 15    //@{
 16    /// Constructor
 17    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2020_I1788444);
 18    //@}
 19
 20    /// Book histograms and initialise projections before the run
 21    void init() {
 22
 23      _mode = 0;
 24      if ( getOption("LMODE") == "EL" ) _mode = 1;
 25      if ( getOption("LMODE") == "MU" ) _mode = 2;
 26
 27	    const FinalState fs;
 28	    // Define fiducial cuts for the leptons in the ZFinder
 29	    Cut lepcuts = (Cuts::pT > 27*GeV) & (Cuts::abseta < 2.5);
 30	    ZFinder zfinderE(fs, lepcuts, PID::ELECTRON, 76*GeV, 106*GeV);
 31	    ZFinder zfinderM(fs, lepcuts, PID::MUON, 76*GeV, 106*GeV);
 32	    declare(zfinderE, "zfinderE");
 33      declare(zfinderM, "zfinderM");
 34
 35	    declare(HeavyHadrons(), "HFHadrons");
 36
 37	    // // Photons
 38	    FinalState photons(Cuts::abspid == PID::PHOTON);
 39	    // Muons
 40	    PromptFinalState bare_mu(Cuts::abspid == PID::MUON, true);
 41	    DressedLeptons all_dressed_mu(photons, bare_mu, 0.1, Cuts::abseta < 2.5, true);
 42	    // Electrons
 43	    PromptFinalState bare_el(Cuts::abspid == PID::ELECTRON, true);
 44	    DressedLeptons all_dressed_el(photons, bare_el, 0.1, Cuts::abseta < 2.5, true);
 45
 46	    //Jet forming
 47	    VetoedFinalState vfs(FinalState(Cuts::abseta < 4.5));
 48	    vfs.addVetoOnThisFinalState(all_dressed_el);
 49	    vfs.addVetoOnThisFinalState(all_dressed_mu);
 50
 51	    FastJets jets(vfs, FastJets::ANTIKT, 0.4, JetAlg::Muons::ALL, JetAlg::Invisibles::DECAY);
 52	    declare(jets, "jets");
 53
 54	    // book histos - binning taken from data.yoda
 55      book(_h["i1b_ZpT"],2,1,1);
 56      book(_h["i1b_ZY"],4,1,1);
 57      book(_h["i1b_dPhiZb"],6,1,1);
 58      book(_h["i1b_dRZb"],8,1,1);
 59      book(_h["i1b_dYZb"],7,1,1);
 60      book(_h["i1b_bpT"],3,1,1);
 61      book(_h["i1b_bY"],5,1,1);
 62
 63      book(_h["i2b_ZpT"],13,1,1);
 64      book(_h["i2b_dPhibb"],9,1,1);
 65      book(_h["i2b_dRbb"],11,1,1);
 66      book(_h["i2b_dYbb"],10,1,1);
 67      book(_h["i2b_Mbb"],12,1,1);
 68      book(_h["i2b_pTbb"],14,1,1);
 69      book(_h["i2b_pTOnMbb"],15,1,1);
 70
 71      book(_h["ib_nBJets"],1,1,1);
 72    }
 73
 74    /////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
 75    /// Perform the per-event analysis
 76    void analyze(const Event& event) {
 77
 78      const ZFinder& zfinderE = apply<ZFinder>(event, "zfinderE");
 79      const Particles& els = zfinderE.constituents();
 80	    const ZFinder& zfinderM = apply<ZFinder>(event, "zfinderM");
 81	    const Particles& mus = zfinderM.constituents();
 82
 83      // default is to run average of Z->ee and Z->mm
 84      // use LMODE option to pick one channel
 85	    if ( (els.size() + mus.size()) != 2 )  vetoEvent;
 86
 87      if (      _mode == 0 && !(els.size()==2 || mus.size()==2) )  vetoEvent;
 88	    else if ( _mode == 1 && !(els.size() == 2 && mus.empty()) )  vetoEvent;
 89	    else if ( _mode == 2 && !(els.empty() && mus.size() == 2) )  vetoEvent;
 90
 91	    double Vpt = 0, Vy = 0, Veta = 0, Vphi = 0;
 92
 93	    if ( els.size()==2 ) {
 94        Vpt = zfinderE.boson().pt()/GeV;
 95        Vphi = zfinderE.boson().phi();
 96        Vy = zfinderE.boson().rapidity();
 97        Veta = zfinderE.boson().eta();
 98	    } else {
 99        Vpt = zfinderM.boson().pt()/GeV;
100        Vphi = zfinderM.boson().phi();
101        Vy = zfinderM.boson().rapidity();
102        Veta = zfinderM.boson().eta();
103	    }
104
105	    Jets jets = apply<JetAlg>(event, "jets").jetsByPt(Cuts::pT>20*GeV && Cuts::absrap < 2.5);
106      idiscardIfAnyDeltaRLess(jets, els, 0.4);
107      idiscardIfAnyDeltaRLess(jets, mus, 0.4);
108
109	    Jets btagged;
110	    const Particles allBs = apply<HeavyHadrons>(event, "HFHadrons").bHadrons(5.0*GeV);
111	    Particles matchedBs;
112
113	    for (const Jet& j : jets) {
114        Jet closest_j;
115        Particle closest_b;
116        double minDR_j_b = 10;
117
118        for (const Particle& bHad : allBs) {
119          bool alreadyMatched = false;
120          for (const Particle& bMatched : matchedBs) {
121            alreadyMatched |= bMatched.isSame(bHad);
122          }
123          if(alreadyMatched)  continue;
124
125          double DR_j_b = deltaR(j, bHad);
126          if ( DR_j_b <= 0.3 && DR_j_b < minDR_j_b) {
127            minDR_j_b = DR_j_b;
128            closest_j = j;
129            closest_b = bHad;
130          }
131        }
132
133        if(minDR_j_b < 0.3) {
134          btagged += closest_j;
135          matchedBs += closest_b;
136        }
137	    }
138	    //size_t njets = jets.size();
139	    size_t ntags = btagged.size();
140	    if (ntags < 1) vetoEvent;
141
142	    _h["ib_nBJets"]->fill(1); //inclusive 1-b
143
144	    double dYVb = fabs(Vy - btagged[0].rap());
145	    double dEtaVb = fabs(Veta - btagged[0].eta());
146	    double dPhiVb = deltaPhi(Vphi, btagged[0]);
147	    double dRVb = sqrt(dEtaVb*dEtaVb + dPhiVb*dPhiVb);
148
149       _h["i1b_ZpT"]   ->fill(Vpt/GeV);
150       _h["i1b_ZY"]    ->fill(fabs(Vy));
151       _h["i1b_dPhiZb"]->fill(dPhiVb);
152       _h["i1b_dRZb"]->fill(dRVb);
153       _h["i1b_dYZb"]->fill(dYVb);
154       _h["i1b_bpT"]->fill(btagged[0].pt()/GeV);
155       _h["i1b_bY"]->fill(btagged[0].absrap());
156
157	    if ( ntags>1 ) {
158        _h["ib_nBJets"]->fill(2); //inclusive 2-b
159
160        double dYbb   = fabs(btagged[0].rap() - btagged[1].rap());
161        double dPhibb = deltaPhi(btagged[0], btagged[1]);
162        double dRbb   = deltaR(btagged[0], btagged[1]);
163        double Mbb    = (btagged[0].mom() + btagged[1].mom()).mass()/GeV;
164        double Ptbb   = (btagged[0].mom() + btagged[1].mom()).pt()/GeV;
165
166        _h["i2b_ZpT"]->fill(Vpt);
167        _h["i2b_dPhibb"]->fill(dPhibb);
168        _h["i2b_dRbb"]->fill(dRbb);
169        _h["i2b_dYbb"]->fill(dYbb);
170        _h["i2b_Mbb"]->fill(Mbb);
171        _h["i2b_pTbb"]->fill(Ptbb);
172        _h["i2b_pTOnMbb"]->fill(Ptbb/Mbb);
173
174      }
175    }
176
177
178    void finalize() {
179      // routine accepts both Z->ee and Z->mm
180      // data corresponds to average
181      const double sf = _mode? 1.0 : 0.5;
182      scale(_h, sf * crossSectionPerEvent());
183    }
184
185    protected:
186
187      size_t _mode;
188
189    private:
190
191      map<string, Histo1DPtr> _h;
192
193  };
194
195  RIVET_DECLARE_PLUGIN(ATLAS_2020_I1788444);
196}