rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2013_I1230812

$Z$ + jets in $pp$ at 7 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1230812
Status: VALIDATED
Authors:
  • Katharina Bierwagen
  • Frank Siegert
References:
  • arXiv: 1304.7098
  • J. High Energy Phys. 07 (2013) 032
Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • Z+jets, electronic Z-decays (data are a weighted combination of electron/muon).

Measurements of the production of jets of particles in association with a $Z$ boson in $pp$ collisions at $\sqrt{s}$ = 7 TeV are presented, using data corresponding to an integrated luminosity of 4.6/fb collected by the ATLAS experiment at the Large Hadron Collider. Inclusive and differential jet cross sections in $Z$ events, with Z decaying into electron or muon pairs, are measured for jets with transverse momentum $p_T > 30$ GeV and rapidity $|y| < 4.4$. This Rivet module implements the event selection for the weighted combination of both decay channels to produce the average cross section for a single lepton flavour, and uses the data from that combination (as in the paper plots). In the default mode this is how it will run, assuming mixed electronic and muonic events. If LMODE is set to EL (MU), the plots for the individual decay channels (with a slightly different fiducial phase space) are made.

Source code: ATLAS_2013_I1230812.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/PromptFinalState.hh"
  7
  8namespace Rivet {
  9
 10
 11  /// Z + jets in pp at 7 TeV (combined channel / base class)
 12  /// @note This base class contains a "mode" variable for combined, e, and mu channel derived classes
 13  class ATLAS_2013_I1230812 : public Analysis {
 14  public:
 15
 16    /// @name Constructors etc.
 17    //@{
 18
 19    /// Constructor
 20    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2013_I1230812);
 21    //@}
 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 = 0;
 29      if ( getOption("LMODE") == "EL" ) _mode = 1;
 30      if ( getOption("LMODE") == "MU" ) _mode = 2;
 31
 32      // Determine the e/mu decay channels used (NB Prompt leptons only).
 33      /// @todo Note that Zs are accepted with any rapidity: the cuts are on the e/mu: is this correct?
 34      Cut pt20 = Cuts::pT >= 20*GeV;
 35      Cut eta_e = _mode? Cuts::abseta < 1.37 || Cuts::absetaIn(1.52, 2.47) : Cuts::abseta < 2.5;
 36      Cut eta_m = _mode? Cuts::abseta < 2.4 : Cuts::abseta < 2.5;
 37      ZFinder zfinder_el(FinalState(eta_e), pt20, PID::ELECTRON, 66*GeV, 116*GeV);
 38      ZFinder zfinder_mu(FinalState(eta_m), pt20, PID::MUON, 66*GeV, 116*GeV);
 39      declare(zfinder_el, "zfinder_el");
 40      declare(zfinder_mu, "zfinder_mu");
 41
 42      // Define veto FS in order to prevent Z-decay products entering the jet algorithm
 43      VetoedFinalState had_fs;
 44      had_fs.addVetoOnThisFinalState(getProjection<ZFinder>("zfinder_el"));
 45      had_fs.addVetoOnThisFinalState(getProjection<ZFinder>("zfinder_mu"));
 46      FastJets jets(had_fs, FastJets::ANTIKT, 0.4, JetAlg::Muons::ALL, JetAlg::Invisibles::ALL);
 47      declare(jets, "jets");
 48
 49      book(_h_njet_incl              ,  1, 1, _mode+1);
 50      book(_h_njet_incl_ratio        ,  2, 1, _mode+1, true);
 51      book(_h_njet_excl              ,  3, 1, _mode+1);
 52      book(_h_njet_excl_ratio        ,  4, 1, _mode+1, true);
 53      book(_h_njet_excl_pt150        ,  5, 1, _mode+1);
 54      book(_h_njet_excl_pt150_ratio  ,  6, 1, _mode+1, true);
 55      book(_h_njet_excl_vbf          ,  7, 1, _mode+1);
 56      book(_h_njet_excl_vbf_ratio    ,  8, 1, _mode+1, true);
 57      book(_h_ptlead                 ,  9, 1, _mode+1);
 58      book(_h_ptseclead              , 10, 1, _mode+1);
 59      book(_h_ptthirdlead            , 11, 1, _mode+1);
 60      book(_h_ptfourthlead           , 12, 1, _mode+1);
 61      book(_h_ptlead_excl            , 13, 1, _mode+1);
 62      book(_h_pt_ratio               , 14, 1, _mode+1);
 63      book(_h_pt_z                   , 15, 1, _mode+1);
 64      book(_h_pt_z_excl              , 16, 1, _mode+1);
 65      book(_h_ylead                  , 17, 1, _mode+1);
 66      book(_h_yseclead               , 18, 1, _mode+1);
 67      book(_h_ythirdlead             , 19, 1, _mode+1);
 68      book(_h_yfourthlead            , 20, 1, _mode+1);
 69      book(_h_deltay                 , 21, 1, _mode+1);
 70      book(_h_mass                   , 22, 1, _mode+1);
 71      book(_h_deltaphi               , 23, 1, _mode+1);
 72      book(_h_deltaR                 , 24, 1, _mode+1);
 73      book(_h_ptthirdlead_vbf        , 25, 1, _mode+1);
 74      book(_h_ythirdlead_vbf         , 26, 1, _mode+1);
 75      book(_h_ht                     , 27, 1, _mode+1);
 76      book(_h_st                     , 28, 1, _mode+1);
 77    }
 78
 79
 80    /// Perform the per-event analysis
 81    void analyze(const Event& event) {
 82
 83      FourMomentum z, lp, lm;
 84      const ZFinder& zfinder_el = apply<ZFinder>(event, "zfinder_el");
 85      const ZFinder& zfinder_mu = apply<ZFinder>(event, "zfinder_mu");
 86
 87      bool e_ok = zfinder_el.constituents().size() == 2 && zfinder_mu.constituents().size() ==0;
 88      bool m_ok = zfinder_el.constituents().size() == 0 && zfinder_mu.constituents().size() ==2;
 89
 90      if (_mode == 0 &&  !e_ok && !m_ok ) vetoEvent;
 91      if (_mode == 1 && !e_ok) vetoEvent;
 92      if (_mode == 2 && !m_ok) vetoEvent;
 93
 94      if (zfinder_el.constituents().size() == 2) {
 95        z = zfinder_el.boson().momentum();
 96        lp = zfinder_el.constituents()[0].momentum();
 97        lm = zfinder_el.constituents()[1].momentum();
 98      }
 99      else if (zfinder_mu.constituents().size() == 2) {
100        z = zfinder_mu.boson().momentum();
101        lp = zfinder_mu.constituents()[0].momentum();
102        lm = zfinder_mu.constituents()[1].momentum();
103      }
104      else  vetoEvent;
105
106      if (deltaR(lp, lm) < 0.2) vetoEvent;
107
108      Jets jets = apply<FastJets>(event, "jets").jetsByPt(Cuts::pT > 30*GeV && Cuts::absrap < 4.4);
109      ifilter_discard(jets, deltaRLess(lp, 0.5));
110      ifilter_discard(jets, deltaRLess(lm, 0.5));
111
112      // Fill jet multiplicities
113      for (size_t ijet = 0; ijet <= jets.size(); ++ijet) {
114        _h_njet_incl->fill(ijet);
115      }
116      _h_njet_excl->fill(jets.size());
117
118      // Require at least one jet
119      if (jets.size() >= 1) {
120        // Leading jet histos
121        const double ptlead   = jets[0].pT()/GeV;
122        const double yabslead = fabs(jets[0].rapidity());
123        const double ptz   = z.pT()/GeV;
124        _h_ptlead->fill(ptlead);
125        _h_ylead ->fill(yabslead);
126        _h_pt_z  ->fill(ptz);
127        // Fill jet multiplicities
128        if (ptlead > 150)  _h_njet_excl_pt150->fill(jets.size());
129
130        // Loop over selected jets, fill inclusive distributions
131        double st = 0;
132        double ht = lp.pT()/GeV + lm.pT()/GeV;
133        for (size_t ijet = 0; ijet < jets.size(); ++ijet) {
134          ht += jets[ijet].pT()/GeV;
135          st += jets[ijet].pT()/GeV;
136        }
137        _h_ht->fill(ht);
138        _h_st->fill(st);
139
140        // Require exactly one jet
141        if (jets.size() == 1) {
142          _h_ptlead_excl->fill(ptlead);
143          _h_pt_z_excl  ->fill(ptz);
144        }
145      }
146
147
148      // Require at least two jets
149      if (jets.size() >= 2) {
150        // Second jet histos
151        const double ptlead      = jets[0].pT()/GeV;
152        const double pt2ndlead   = jets[1].pT()/GeV;
153        const double ptratio     = pt2ndlead/ptlead;
154        const double yabs2ndlead = fabs(jets[1].rapidity());
155        _h_ptseclead->fill(pt2ndlead);
156        _h_yseclead->fill( yabs2ndlead);
157        _h_pt_ratio->fill( ptratio);
158
159        // Dijet histos
160        const double deltaphi = fabs(deltaPhi(jets[1], jets[0]));
161        const double deltarap = fabs(jets[0].rapidity() - jets[1].rapidity()) ;
162        const double deltar   = fabs(deltaR(jets[0], jets[1], RAPIDITY));
163        const double mass     = (jets[0].momentum() + jets[1].momentum()).mass()/GeV;
164        _h_mass->fill(    mass);
165        _h_deltay->fill(  deltarap);
166        _h_deltaphi->fill(deltaphi);
167        _h_deltaR->fill(  deltar);
168
169        if (mass > 350 && deltarap > 3)  _h_njet_excl_vbf->fill(jets.size());
170      }
171
172      // Require at least three jets
173      if (jets.size() >= 3) {
174        // Third jet histos
175        const double pt3rdlead   = jets[2].pT()/GeV;
176        const double yabs3rdlead = fabs(jets[2].rapidity());
177        _h_ptthirdlead->fill(pt3rdlead);
178        _h_ythirdlead->fill( yabs3rdlead);
179
180        //Histos after VBF preselection
181        const double deltarap = fabs(jets[0].rapidity() - jets[1].rapidity()) ;
182        const double mass     = (jets[0].momentum() + jets[1].momentum()).mass();
183        if (mass > 350 && deltarap > 3) {
184          _h_ptthirdlead_vbf->fill(pt3rdlead);
185          _h_ythirdlead_vbf->fill( yabs3rdlead);
186        }
187      }
188
189      // Require at least four jets
190      if (jets.size() >= 4) {
191        // Fourth jet histos
192        const double pt4thlead   = jets[3].pT()/GeV;
193        const double yabs4thlead = fabs(jets[3].rapidity());
194        _h_ptfourthlead->fill(pt4thlead);
195        _h_yfourthlead->fill( yabs4thlead);
196      }
197    }
198
199    /// @name Ratio calculator util functions
200    //@{
201
202    /// Calculate the efficiency error, being careful about div-by-zero
203    double err_incl(const HistoBin1D &M, const HistoBin1D &N, bool hasWeights) {
204      double r = safediv(M.sumW(), N.sumW());
205      if (hasWeights) { // use F. James's approximation for weighted events
206        return sqrt( safediv((1 - 2 * r) * M.sumW2() + r * r * N.sumW2(), N.sumW() * N.sumW()) );
207      }
208      return sqrt( safediv(r * (1 - r), N.sumW()) );
209    }
210
211    /// Calculate the ratio error, being careful about div-by-zero
212    double err_excl(const HistoBin1D &A, const HistoBin1D &B) {
213      double r = safediv(A.sumW(), B.sumW());
214      double dAsquared = safediv(A.sumW2(), A.sumW() * A.sumW()); // squared relative error of A
215      double dBsquared = safediv(B.sumW2(), B.sumW() * B.sumW()); // squared relative error of B
216      return r * sqrt(dAsquared + dBsquared);
217    }
218
219    //@}
220
221
222    void finalize() {
223      bool hasWeights = _h_njet_incl->effNumEntries() != _h_njet_incl->numEntries();
224      for (size_t i = 0; i < 6; ++i) {
225        _h_njet_incl_ratio->point(i).setY(safediv(_h_njet_incl->bin(i + 1).sumW(), _h_njet_incl->bin(i).sumW()),
226                                          err_incl(_h_njet_incl->bin(i + 1), _h_njet_incl->bin(i), hasWeights));
227        _h_njet_excl_ratio->point(i).setY(safediv(_h_njet_excl->bin(i + 1).sumW(), _h_njet_excl->bin(i).sumW()),
228                                          err_excl(_h_njet_excl->bin(i + 1), _h_njet_excl->bin(i)));
229        if (i >= 1) {
230          _h_njet_excl_pt150_ratio->point(i - 1).setY(safediv(_h_njet_excl_pt150->bin(i).sumW(), _h_njet_excl_pt150->bin(i - 1).sumW()),
231                                                      err_excl(_h_njet_excl_pt150->bin(i), _h_njet_excl_pt150->bin(i - 1)));
232          if (i >= 2) {
233            _h_njet_excl_vbf_ratio->point(i - 2).setY(safediv(_h_njet_excl_vbf->bin(i).sumW(), _h_njet_excl_vbf->bin(i - 1).sumW()),
234                                                      err_excl(_h_njet_excl_vbf->bin(i), _h_njet_excl_vbf->bin(i - 1)));
235          }
236        }
237      }
238
239      double sf = _mode? 1.0 : 0.5;
240      const double xs = sf * crossSectionPerEvent()/picobarn;
241
242      scale(_h_njet_incl, xs); scale(_h_njet_excl, xs); scale(_h_njet_excl_pt150, xs); 
243      scale(_h_njet_excl_vbf, xs); scale(_h_ptlead, xs); scale(_h_ptseclead, xs); 
244      scale(_h_ptthirdlead, xs); scale(_h_ptfourthlead, xs); scale(_h_ptlead_excl, xs);
245      scale(_h_pt_ratio, xs); scale(_h_pt_z, xs); scale(_h_pt_z_excl, xs);
246      scale(_h_ylead, xs); scale(_h_yseclead, xs); scale(_h_ythirdlead, xs); 
247      scale(_h_yfourthlead, xs); scale(_h_deltay, xs); scale(_h_mass, xs); 
248      scale(_h_deltaphi, xs); scale(_h_deltaR, xs); scale(_h_ptthirdlead_vbf, xs); 
249      scale(_h_ythirdlead_vbf, xs); scale(_h_ht, xs); scale(_h_st, xs);
250    }
251
252    //@}
253
254
255  protected:
256
257    size_t _mode;
258
259
260  private:
261
262    Scatter2DPtr _h_njet_incl_ratio;
263    Scatter2DPtr _h_njet_excl_ratio;
264    Scatter2DPtr _h_njet_excl_pt150_ratio;
265    Scatter2DPtr _h_njet_excl_vbf_ratio;
266    Histo1DPtr _h_njet_incl;
267    Histo1DPtr _h_njet_excl;
268    Histo1DPtr _h_njet_excl_pt150;
269    Histo1DPtr _h_njet_excl_vbf;
270    Histo1DPtr _h_ptlead;
271    Histo1DPtr _h_ptseclead;
272    Histo1DPtr _h_ptthirdlead;
273    Histo1DPtr _h_ptfourthlead;
274    Histo1DPtr _h_ptlead_excl;
275    Histo1DPtr _h_pt_ratio;
276    Histo1DPtr _h_pt_z;
277    Histo1DPtr _h_pt_z_excl;
278    Histo1DPtr _h_ylead;
279    Histo1DPtr _h_yseclead;
280    Histo1DPtr _h_ythirdlead;
281    Histo1DPtr _h_yfourthlead;
282    Histo1DPtr _h_deltay;
283    Histo1DPtr _h_mass;
284    Histo1DPtr _h_deltaphi;
285    Histo1DPtr _h_deltaR;
286    Histo1DPtr _h_ptthirdlead_vbf;
287    Histo1DPtr _h_ythirdlead_vbf;
288    Histo1DPtr _h_ht;
289    Histo1DPtr _h_st;
290  };
291
292
293  RIVET_DECLARE_PLUGIN(ATLAS_2013_I1230812);
294
295}