Rivet analyses referenceATLAS_2016_I1467454High-mass Drell-Yan at 8 TeVExperiment: ATLAS (LHC) Inspire ID: 1467454 Status: VALIDATED Authors:
Beam energies: (4000.0, 4000.0) GeV Run details:
This paper presents a measurement of the double-differential cross section for the Drell-Yan $Z/\gamma^\ast \rightarrow \ell^+\ell^-$ and photon-induced $\gamma\gamma \rightarrow \ell^+\ell^-$ processes where $\ell$ is an electron or muon. The measurement is performed for invariant masses of the lepton pairs, $m_{\ell\ell}$, between 116 GeV and 1500 GeV using a sample of 20.3 fb$^{-1}$ of pp collisions data at centre-of-mass energy of $\sqrt{s} = 8$ TeV collected by the ATLAS detector at the LHC in 2012. The data are presented double differentially in invariant mass and absolute dilepton rapidity as well as in invariant mass and absolute pseudorapidity separation of the lepton pair. The single-differential cross section as a function of $m_{\ell\ell}$ is also reported. The electron and muon channel measurements are combined and a total experimental precision of better than 1\% is achieved at low $m_{\ell\ell}$. A comparison to next-to-next-to-leading order perturbative QCD predictions using several recent parton distribution functions and including next-to-leading order electroweak effects indicates the potential of the data to constrain parton distribution functions. In particular, a large impact of the data on the photon PDF is demonstrated. Decay channel (e or mu) is selected by LMODE=EL,MU Source code: ATLAS_2016_I1467454.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/DileptonFinder.hh"
5
6namespace Rivet {
7
8
9 /// @brief High-mass Drell-Yan at 8 TeV
10 class ATLAS_2016_I1467454 : public Analysis {
11 public:
12
13 /// Constructor
14 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2016_I1467454);
15 /// @}
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 new option system
25 _mode = 0;
26 if ( getOption("LMODE") == "EL" ) _mode = 0;
27 if ( getOption("LMODE") == "MU" ) _mode = 1;
28
29 Cut cuts = Cuts::abseta < 2.5 && Cuts::pT > 30*GeV;
30 DileptonFinder zfinder(91.2*GeV, 0.1, cuts && Cuts::abspid == (_mode ? PID::MUON : PID::ELECTRON), Cuts::massIn(116*GeV, 1500*GeV));
31 declare(zfinder, "DileptonFinder");
32
33 size_t ch = _mode? 11 : 0; // offset
34 book(_hist_mll, 18 + ch, 1, 1);
35
36 const vector<double> mll_bins{ 116., 150., 200., 300., 500., 1500. };
37 book(_hist_rap, mll_bins);
38 book(_hist_deta, mll_bins);
39 for (size_t i=0; i < _hist_rap->numBins(); ++i) {
40 book(_hist_rap->bin(i+1), 19 + ch + i, 1, 1);
41 book(_hist_deta->bin(i+1), 24 + ch + i, 1, 1);
42 }
43
44 }
45
46
47 /// Perform the per-event analysis
48 void analyze(const Event& event) {
49
50 const DileptonFinder& zfinder = apply<DileptonFinder>(event, "DileptonFinder");
51 if (zfinder.bosons().size() != 1) vetoEvent;
52
53 const Particle z0 = zfinder.bosons()[0];
54 /// @todo Could use z0.constituents()
55 const Particle el1 = zfinder.leptons()[0];
56 const Particle el2 = zfinder.leptons()[1];
57
58 if (el1.pT() > 40*GeV || el2.pT() > 40*GeV) {
59 const double mass = z0.mass();
60 _hist_mll->fill(mass/GeV);
61 _hist_rap->fill(mass/GeV, z0.absrap());
62 _hist_deta->fill(mass/GeV, deltaEta(el1,el2));
63 }
64 }
65
66
67 /// Normalise histograms etc., after the run
68 void finalize() {
69
70 const double sf = crossSection()/picobarn/sumOfWeights();
71 scale(_hist_mll, sf);
72 scale(_hist_rap, sf*0.5);
73 divByGroupWidth(_hist_rap);
74 scale(_hist_deta, sf*0.5);
75 divByGroupWidth(_hist_deta);
76
77 }
78
79 /// @}
80
81
82 /// Choose to work in electron or muon mode
83 size_t _mode;
84
85
86 /// @name Histograms
87 /// @{
88 Histo1DPtr _hist_mll;
89 Histo1DGroupPtr _hist_rap, _hist_deta;
90 /// @}
91
92 };
93
94 RIVET_DECLARE_PLUGIN(ATLAS_2016_I1467454);
95
96}
|