rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2014_I1306294

Measurement of Z boson in association with b-jets at 7 TeV in ATLAS (electron channel)
Experiment: ATLAS (LHC)
Inspire ID: 1306294
Status: VALIDATED
Authors:
  • Gavin Hesketh
References: Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • Z+b(b) production in $pp$ collisions at $\sqrt{s} = 7$ TeV, electronic Z-decays

Measurements of differential production cross-sections of a $Z$ boson in association with $b$-jets in $pp$ collisions at $\sqrt{s}=7$ TeV are reported. The data analysed correspond to an integrated luminosity of 4.6 fb$^{-1}$ recorded with the ATLAS detector at the Large Hadron Collider. Particle-level cross-sections are determined for events with a $Z$ boson decaying into an electron or muon pair, and containing $b$-jets. For events with at least one $b$-jet, the cross-section is presented as a function of the $Z$ boson transverse momentum and rapidity, together with the inclusive $b$-jet cross-section as a function of $b$-jet transverse momentum, rapidity and angular separations between the $b$-jet and the $Z$ boson. For events with at least two $b$-jets, the cross-section is determined as a function of the invariant mass and angular separation of the two highest transverse momentum $b$-jets, and as a function of the $Z$ boson transverse momentum and rapidity. Results are compared to leading-order and next-to-leading-order perturbative QCD calculations. This Rivet module implements the event selection for Z decaying into electrons. If you want to use muonic events, please refer to ATLAS\_2014\_I1306294\_MU

Source code: ATLAS_2014_I1306294.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/ZFinder.hh"
  5#include "Rivet/Projections/FastJets.hh"
  6#include "Rivet/Projections/HeavyHadrons.hh"
  7#include "Rivet/Projections/VetoedFinalState.hh"
  8
  9namespace Rivet {
 10
 11
 12
 13  /// Electroweak Wjj production at 8 TeV
 14  class ATLAS_2014_I1306294 : public Analysis {
 15  public:
 16
 17    /// Constructor
 18    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2014_I1306294);
 19
 20
 21    /// @name Analysis methods
 22    //@{
 23
 24    /// Book histograms and initialise projections before the run
 25    void init() {
 26
 27      // Get options from the new option system
 28      _mode = 1;
 29      if ( getOption("LMODE") == "EL" ) _mode = 1;
 30      if ( getOption("LMODE") == "MU" ) _mode = 2;
 31
 32      FinalState fs;
 33      Cut cuts = Cuts::abseta < 2.5 && Cuts::pT > 20*GeV;
 34
 35      ZFinder zfinder(fs, cuts, _mode==1? PID::ELECTRON : PID::MUON, 76.0*GeV, 106.0*GeV, 0.1, 
 36                      ZFinder::ChargedLeptons::ALL, ZFinder::ClusterPhotons::NODECAY, ZFinder::AddPhotons::NO);
 37      declare(zfinder, "ZFinder");
 38
 39      VetoedFinalState jet_fs(fs);
 40      jet_fs.addVetoOnThisFinalState(getProjection<ZFinder>("ZFinder"));
 41      FastJets jetpro1(jet_fs, FastJets::ANTIKT, 0.4, JetAlg::Muons::ALL, JetAlg::Invisibles::ALL);
 42      declare(jetpro1, "AntiKtJets04");
 43      declare(HeavyHadrons(), "BHadrons");
 44
 45      // Histograms with data binning
 46      book(_h_bjet_Pt      , 3, 1, 1);
 47      book(_h_bjet_Y       , 5, 1, 1);
 48      book(_h_bjet_Yboost  , 7, 1, 1);
 49      book(_h_bjet_DY20    , 9, 1, 1);
 50      book(_h_bjet_ZdPhi20 ,11, 1, 1);
 51      book(_h_bjet_ZdR20   ,13, 1, 1);
 52      book(_h_bjet_ZPt     ,15, 1, 1);
 53      book(_h_bjet_ZY      ,17, 1, 1);
 54      book(_h_2bjet_dR     ,21, 1, 1);
 55      book(_h_2bjet_Mbb    ,23, 1, 1);
 56      book(_h_2bjet_ZPt    ,25, 1, 1);
 57      book(_h_2bjet_ZY     ,27, 1, 1);
 58    }
 59
 60
 61    /// Perform the per-event analysis
 62    void analyze(const Event& e) {
 63
 64      // Check we have a Z:
 65      const ZFinder& zfinder = apply<ZFinder>(e, "ZFinder");
 66      if (zfinder.bosons().size() != 1) vetoEvent;
 67      
 68      const Particles boson_s =  zfinder.bosons();
 69      const Particle boson_f =  boson_s[0];
 70      const Particles zleps   =  zfinder.constituents();
 71
 72      // Stop processing the event if no true b-partons or hadrons are found
 73      const Particles allBs = apply<HeavyHadrons>(e, "BHadrons").bHadrons(5.0*GeV);
 74      Particles stableBs = filter_select(allBs, Cuts::abseta < 2.5);
 75      if (stableBs.empty()) vetoEvent;
 76
 77      // Get the b-jets
 78      const Jets& jets = apply<JetAlg>(e, "AntiKtJets04").jetsByPt(Cuts::pT >20.0*GeV && Cuts::abseta <2.4);
 79      Jets b_jets;
 80      for (const Jet& jet : jets) {
 81        //veto overlaps with Z leptons:
 82        bool veto = false;
 83        for (const Particle& zlep : zleps) {
 84          if (deltaR(jet, zlep) < 0.5) veto = true;
 85        }
 86        if (veto) continue;
 87
 88        for (const Particle& bhadron : stableBs) {
 89          if (deltaR(jet, bhadron) <= 0.3) {
 90            b_jets.push_back(jet);
 91            break; // match
 92          }
 93	}
 94      }
 95
 96      // Make sure we have at least 1
 97      if (b_jets.empty()) vetoEvent;
 98
 99      // Fill the plots
100      const double ZpT = boson_f.pT()/GeV;
101      const double ZY  = boson_f.absrap();
102
103      _h_bjet_ZPt->fill(ZpT);
104      _h_bjet_ZY ->fill(ZY);
105
106      for (const Jet& jet : b_jets) {
107        _h_bjet_Pt->fill(jet.pT()/GeV);
108        _h_bjet_Y ->fill(jet.absrap());
109
110        const double Yboost = 0.5 * fabs(boson_f.rapidity() + jet.rapidity());
111
112        _h_bjet_Yboost->fill(Yboost);
113
114        if(ZpT > 20.) {
115
116          const double ZBDY   = fabs( boson_f.rapidity() - jet.rapidity() );
117          const double ZBDPHI = fabs( deltaPhi(jet.phi(), boson_f.phi()) );
118          const double ZBDR   = deltaR(jet, boson_f, RAPIDITY);
119          _h_bjet_DY20->fill(   ZBDY);
120          _h_bjet_ZdPhi20->fill(ZBDPHI);
121          _h_bjet_ZdR20->fill(  ZBDR);
122        }
123
124      } //loop over b-jets
125
126      if (b_jets.size() < 2) return;
127
128      _h_2bjet_ZPt->fill(ZpT);
129      _h_2bjet_ZY ->fill(ZY);
130
131      const double BBDR = deltaR(b_jets[0], b_jets[1], RAPIDITY);
132      const double Mbb  = (b_jets[0].momentum() + b_jets[1].momentum()).mass();
133
134      _h_2bjet_dR ->fill(BBDR);
135      _h_2bjet_Mbb->fill(Mbb);
136
137    } // end of analysis loop
138
139
140    /// Normalise histograms etc., after the run
141    void finalize() {
142
143      const double normfac = crossSection() / sumOfWeights();
144
145      scale( _h_bjet_Pt,      normfac);
146      scale( _h_bjet_Y,       normfac);
147      scale( _h_bjet_Yboost,  normfac);
148      scale( _h_bjet_DY20,    normfac);
149      scale( _h_bjet_ZdPhi20, normfac);
150      scale( _h_bjet_ZdR20,   normfac);
151      scale( _h_bjet_ZPt,     normfac);
152      scale( _h_bjet_ZY,      normfac);
153      scale( _h_2bjet_dR,     normfac);
154      scale( _h_2bjet_Mbb,    normfac);
155      scale( _h_2bjet_ZPt,    normfac);
156      scale( _h_2bjet_ZY,     normfac);
157    }
158
159    //@}
160
161
162  protected:
163
164    // Data members like post-cuts event weight counters go here
165    size_t _mode;
166
167
168  private:
169
170    Histo1DPtr _h_bjet_Pt;
171    Histo1DPtr _h_bjet_Y;
172    Histo1DPtr _h_bjet_Yboost;
173    Histo1DPtr _h_bjet_DY20;
174    Histo1DPtr _h_bjet_ZdPhi20;
175    Histo1DPtr _h_bjet_ZdR20;
176    Histo1DPtr _h_bjet_ZPt;
177    Histo1DPtr _h_bjet_ZY;
178    Histo1DPtr _h_2bjet_dR;
179    Histo1DPtr _h_2bjet_Mbb;
180    Histo1DPtr _h_2bjet_ZPt;
181    Histo1DPtr _h_2bjet_ZY;
182
183  };
184
185
186  RIVET_DECLARE_PLUGIN(ATLAS_2014_I1306294);
187
188}