rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2016_I1502620

W and Z inclusive cross sections at 7 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1502620
Status: VALIDATED
Authors:
  • Jan Kretzschmar
  • Christian Gutschow
References: Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • p + p -> Z + X ( Z -> e e )

High-precision measurements by the ATLAS Collaboration are presented of inclusive $W^+\to\ell^+$, $W^-\to\ell^-$, $Z/\gamma^\ast\to\ell\ell$ ($\ell=e,\mu$) Drell-Yan production cross sections at the LHC. The data were collected in proton-proton collisions at $\sqrt{s}=7$ TeV with an integrated luminosity of 4.6 fb${}^{-1}$. Differential $W^+$ and $W^-$ cross sections are measured in a lepton pseudorapidity range $|\eta_\ell| < 2.5$. Differential $Z/\gamma^\ast$ cross sections are measured as a function of the absolute dilepton rapidity, for $|y_{\ell\ell}| < 3.6$, for three intervals of dilepton mass, $m_{\ell\ell}$, extending from 46 to 150 GeV. The integrated and differential electron- and muon-channel cross sections are combined and compared to theoretical predictions using recent sets of parton distribution functions. The data, together with the final inclusive $e^\pm p$ scattering cross-section data from H1 and ZEUS, are interpreted in a next-to-next-to-leading-order QCD analysis, and a new set of parton distribution functions, ATLAS-epWZ16, is obtained. The ratio of strange-to-light sea-quark densities in the proton is determined more accurately than in previous determinations based on collider data only, and is established to be close to unity in the sensitivity range of the data. A new measurement of the CKM matrix element $|V_{cs}|$ is also provided. N.B.: Use :LMODE to choose specify the signature decay channel directly. Default, run everything, assuming all decay modes are generated Contains Z (W): fill only Z (W) signature histograms Contains EL (MU): assume only electron (muon) decay modes are being generated.

Source code: ATLAS_2016_I1502620.cc
  1#include "Rivet/Analysis.hh"
  2#include "Rivet/Projections/FinalState.hh"
  3#include "Rivet/Projections/MissingMomentum.hh"
  4#include "Rivet/Projections/LeptonFinder.hh"
  5#include "Rivet/Projections/DileptonFinder.hh"
  6
  7namespace Rivet {
  8
  9
 10  /// @brief inclusive W/Z cross sections at 7 TeV
 11  class ATLAS_2016_I1502620 : public Analysis {
 12  public:
 13
 14    /// Constructor
 15    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2016_I1502620);
 16
 17
 18    /// @name Analysis methods
 19    /// @{
 20
 21    /// Book histograms and initialise projections before the run
 22    void init() {
 23
 24      // Get options from the option system
 25      _mode = 0;
 26      _runZ = true;
 27      _runW = true;
 28      if ( getOption("LMODE") == "EL" ||
 29           getOption("LMODE") == "ZEL" ||
 30           getOption("LMODE") == "WEL" )
 31        _mode = 1;
 32      if ( getOption("LMODE") == "MU" ||
 33           getOption("LMODE") == "ZMU" ||
 34           getOption("LMODE") == "WMU" )
 35        _mode = 2;
 36      if ( getOption("LMODE") == "Z" ||
 37           getOption("LMODE") == "ZEL" ||
 38           getOption("LMODE") == "ZMU" )
 39        _runW = false;
 40      if ( getOption("LMODE") == "W" ||
 41           getOption("LMODE") == "WEL" ||
 42           getOption("LMODE") == "WMU" )
 43        _runZ = false;
 44
 45
 46      // Initialise and register projections
 47      declare("MET", MissingMomentum());
 48      LeptonFinder ef(0.1, Cuts::pT > 25*GeV && Cuts::abspid == PID::ELECTRON);
 49      declare(ef, "Elecs");
 50      LeptonFinder mf(0.1, Cuts::pT > 25*GeV && Cuts::abspid == PID::MUON);
 51      declare(mf, "Muons");
 52
 53      DileptonFinder zfe(91.2*GeV, 0.1, Cuts::pT > 20*GeV && Cuts::abspid == PID::ELECTRON, Cuts::massIn(46.0*GeV, 150*GeV));
 54      declare(zfe, "Zee");
 55      DileptonFinder zfm(91.2*GeV, 0.1, Cuts::pT > 20*GeV && Cuts::abspid == PID::MUON, Cuts::massIn(46.0*GeV, 150*GeV));
 56      declare(zfm, "Zmm");
 57
 58
 59      // Book histograms
 60      if (_runW) {
 61        book(_h_Wp_eta,  9, 1, 1);
 62        book(_h_Wm_eta, 10, 1, 1);
 63        book(_h_W_asym, 35, 1, 1);
 64      }
 65      if (_runZ) {
 66        book(_h_Zcenlow_y_dressed,  11, 1, 1);
 67        book(_h_Zcenpeak_y_dressed, 12, 1, 1);
 68        book(_h_Zcenhigh_y_dressed, 13, 1, 1);
 69        book(_h_Zfwdpeak_y_dressed, 14, 1, 1);
 70        book(_h_Zfwdhigh_y_dressed, 15, 1, 1);
 71      }
 72    }
 73
 74
 75    /// Perform the per-event analysis
 76    void analyze(const Event& event) {
 77
 78      // W reco, starting with MET
 79      const P4& pmiss = apply<MissingMom>(event, "MET").missingMom();
 80
 81      // Identify the closest-matching l+MET to m == mW
 82      const Particles& es = apply<LeptonFinder>(event, "Elecs").particles();
 83      const Particles es_mtfilt = select(es, [&](const Particle& e){ return mT(e, pmiss) > 40*GeV; });
 84      const int iefound = closestMatchIndex(es_mtfilt, pmiss, Kin::mass, 80.4*GeV);
 85      const Particles& mus = apply<LeptonFinder>(event, "Muons").particles();
 86      const Particles mus_mtfilt = select(mus, [&](const Particle& m){ return mT(m, pmiss) > 40*GeV; });
 87      const int imfound = closestMatchIndex(mus_mtfilt, pmiss, Kin::mass, 80.4*GeV);
 88
 89      // Require only one W candidate
 90      if (pmiss.pT() < 25*GeV && (int(iefound >= 0) + int(imfound >= 0)) == 1 && _runW) {
 91        Particle lep;
 92        if (_mode != 2 && iefound >= 0) lep = es_mtfilt[iefound];
 93        else if (_mode != 1 && imfound >= 0) lep = mus_mtfilt[imfound];
 94
 95        if (lep.charge3() == 3) _h_Wp_eta->fill(lep.abseta());
 96        else if (lep.charge3() == -3) _h_Wm_eta->fill(lep.abseta());
 97      }
 98
 99
100      // Now the Z stuff.
101      const DileptonFinder& zee = apply<DileptonFinder>(event, "Zee");
102      const DileptonFinder& zmm = apply<DileptonFinder>(event, "Zmm");
103      if ((zee.bosons().size() + zmm.bosons().size()) == 1 && _runZ) {
104        Particle zboson;
105        if (_mode != 2 && zee.bosons().size() == 1) zboson = zee.boson();
106        else if (_mode !=1 && zmm.bosons().size() == 1 ) zboson = zmm.boson();
107
108        const Particles& leptons = zboson.constituents();
109        if (leptons.size() > 1) {
110          const double zrap  = zboson.absrap();
111          const double zmass = zboson.mass();
112          const double eta1 = leptons[0].abseta();
113          const double eta2 = leptons[1].abseta();
114
115          // Separation into central/forward and three mass bins
116          if (eta1 < 2.5 && eta2 < 2.5) {
117            if (zmass < 66.0*GeV)        _h_Zcenlow_y_dressed->fill(zrap);
118            else if (zmass < 116.0*GeV)  _h_Zcenpeak_y_dressed->fill(zrap);
119            else                         _h_Zcenhigh_y_dressed->fill(zrap);
120          }
121          else if ((eta1 < 2.5 && 2.5 < eta2 && eta2 < 4.9) || (eta2 < 2.5 && 2.5 < eta1 && eta1 < 4.9)) {
122            if (zmass > 66.0*GeV) {
123              if (zmass < 116.0*GeV)  _h_Zfwdpeak_y_dressed->fill(zrap);
124              else                    _h_Zfwdhigh_y_dressed->fill(zrap);
125            }
126          }
127        }
128      }
129    }
130
131
132    /// Normalise histograms etc., after the run
133    void finalize() {
134
135      // Construct asymmetry: (dsig+/deta - dsig-/deta) / (dsig+/deta + dsig-/deta)
136      if (_runW) asymm(_h_Wp_eta, _h_Wm_eta, _h_W_asym);
137
138      ///  Normalise, scale and otherwise manipulate histograms here
139      double lfac = 1.0;
140      // If we running on both electrons and muons, divide by two -> xsec for one flavour
141      if (_mode == 0) lfac = 0.5;
142      const double sf = lfac * 0.5 * crossSection() /picobarn / sumOfWeights(); //< 0.5 accounts for rapidity bin width
143
144      if (_runW){
145        scale(_h_Wp_eta, sf);
146        scale(_h_Wm_eta, sf);
147      }
148
149      if (_runZ){
150        scale(_h_Zcenlow_y_dressed, sf);
151        scale(_h_Zcenpeak_y_dressed, sf);
152        scale(_h_Zcenhigh_y_dressed, sf);
153        scale(_h_Zfwdpeak_y_dressed, sf);
154        scale(_h_Zfwdhigh_y_dressed, sf);
155      }
156    }
157
158    /// @}
159
160
161  protected:
162
163    size_t _mode;
164    bool _runZ, _runW;
165
166
167  private:
168
169    /// @name Histograms
170    /// @{
171    Histo1DPtr _h_Wp_eta, _h_Wm_eta;
172    Estimate1DPtr _h_W_asym;
173
174    Histo1DPtr _h_Zcenlow_y_dressed;
175    Histo1DPtr _h_Zcenpeak_y_dressed;
176    Histo1DPtr _h_Zcenhigh_y_dressed;
177    Histo1DPtr _h_Zfwdpeak_y_dressed;
178    Histo1DPtr _h_Zfwdhigh_y_dressed;
179    /// @}
180
181  };
182
183
184  RIVET_DECLARE_PLUGIN(ATLAS_2016_I1502620);
185
186}