Rivet analyses referenceATLAS_2019_I1768911Z pT and Z phi* at 13 TeVExperiment: ATLAS (LHC) Inspire ID: 1768911 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
This paper describes precision measurements of the transverse momentum $p_\mathrm{T}^{\ell\ell}$ ($\ell=e,\mu$) and of the angular variable $\phi^{*}_{\eta}$ distributions of Drell-Yan lepton pairs in a mass range of 66-116 GeV. The analysis uses data from 36.1 fb$^{-1}$ of proton-proton collisions at a centre-of-mass energy of $\sqrt{s}=13$ TeV collected by the ATLAS experiment at the LHC in 2015 and 2016. Measurements in electron-pair and muon-pair final states are performed in the same fiducial volumes, corrected for detector effects, and combined. Compared to previous measurements in proton-proton collisions at $\sqrt{s}=$7 and 8 TeV, this new measurement probes perturbative QCD at a higher centre-of-mass energy with a different composition of initial states. It reaches a precision of 0.2% for the normalized spectra at low values of $p_\mathrm{T}^{\ell\ell}$. The data are compared with different QCD predictions, where it is found that predictions based on resummation approaches can describe the full spectrum within uncertainties. Decay channel (e or mu) is selected by LMODE=EL,MU. The comparison is always to the combined measurement. Source code: ATLAS_2019_I1768911.cc 1#include "Rivet/Analysis.hh"
2#include "Rivet/Projections/FinalState.hh"
3#include "Rivet/Projections/DileptonFinder.hh"
4
5namespace Rivet {
6
7
8 /// @brief Z pT and Z phi* at 13 TeV
9 class ATLAS_2019_I1768911 : public Analysis {
10 public:
11
12 /// Constructor
13 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2019_I1768911);
14
15
16 /// @name Analysis methods
17 /// @{
18
19 /// Book histograms and initialise projections before the run
20 void init() {
21
22 // Get options
23 _mode = 0;
24 if ( getOption("LMODE") == "EL" ) _mode = 1;
25 if ( getOption("LMODE") == "MU" ) _mode = 2;
26
27 // Configure projections
28 Cut cuts = Cuts::abseta < 2.5 && Cuts::pT > 27*GeV;
29 DileptonFinder zmmFinder(91.2*GeV, 0.1, cuts && Cuts::abspid == PID::MUON, Cuts::massIn(66*GeV, 116*GeV));
30 declare(zmmFinder, "DileptonFinder_mu");
31 DileptonFinder zeeFinder(91.2*GeV, 0.1, cuts && Cuts::abspid == PID::ELECTRON, Cuts::massIn(66*GeV, 116*GeV));
32 declare(zeeFinder, "DileptonFinder_el");
33
34 // Book histograms
35 book(_h["zpt_combined_dressed_normalised"], 27, 1, 1);
36 book(_h["zphistar_combined_dressed_normalised"], 28, 1, 1);
37 }
38
39
40 /// Perform the per-event analysis
41 void analyze(const Event& event) {
42
43 // Get leptonic Z boson
44 const DileptonFinder& zmmFinder = apply<DileptonFinder>(event, "DileptonFinder_mu");
45 const DileptonFinder& zeeFinder = apply<DileptonFinder>(event, "DileptonFinder_el");
46 if (_mode == 2 && zmmFinder.bosons().size() != 1 && zeeFinder.bosons().size()) vetoEvent;
47 if (_mode == 1 && zeeFinder.bosons().size() != 1 && zmmFinder.bosons().size()) vetoEvent;
48 if (_mode == 0 && (zeeFinder.bosons().size() + zmmFinder.bosons().size()) != 1) vetoEvent;
49 const Particle& Zboson = zeeFinder.bosons().size()? zeeFinder.boson() : zmmFinder.boson();
50
51 // cut on Z boson leptons and calculate pTll and phistar
52 const Particles& leptons = zeeFinder.bosons().size()? zeeFinder.constituents() : zmmFinder.constituents();
53 if (leptons.size() != 2 || leptons[0].charge3() * leptons[1].charge3() > 0) vetoEvent;
54
55 const double zpt = Zboson.pT()/GeV;
56 const Particle& lminus = leptons[0].charge() < 0 ? leptons[0] : leptons[1];
57 const Particle& lplus = leptons[0].charge() < 0 ? leptons[1] : leptons[0];
58 const double phi_acop = M_PI - deltaPhi(lminus, lplus);
59 const double costhetastar = tanh( 0.5 * (lminus.eta() - lplus.eta()) );
60 const double sin2thetastar = (costhetastar > 1) ? 0.0 : (1.0 - sqr(costhetastar));
61 const double phistar = tan(0.5 * phi_acop) * sqrt(sin2thetastar);
62
63 _h["zpt_combined_dressed_normalised"]->fill(zpt);
64 _h["zphistar_combined_dressed_normalised"]->fill(phistar);
65 }
66
67
68 /// Normalise histograms etc., after the run
69 void finalize() { normalize(_h); }
70 /// @}
71
72
73 protected:
74
75 size_t _mode;
76
77
78 private:
79
80 /// @name Histograms
81 /// @{
82 map<string, Histo1DPtr> _h;
83 /// @}
84
85 };
86
87
88 RIVET_DECLARE_PLUGIN(ATLAS_2019_I1768911);
89}
|