rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CMS_2022_I2079374

Measurement of the mass dependence of the transverse momentum of lepton pairs in Drell--Yan production in proton-proton collisions at 13 TeV
Experiment: CMS (LHC)
Inspire ID: 2079374
Status: VALIDATED
Authors:
  • cms-pag-conveners-smp@cern.ch
  • Louis Moureaux
  • Buğra Bilin
  • Itana Bubanja
  • Philippe Gras
References:
  • CMS-SMP-20-003
  • arXiv: 2205.04897
  • Submitted to Eur. Phys. J. C
Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • Run MC generators with Z decaying leptonically into both electrons and muons at 13 TeV CoM energy. If only one of the two decay channels is included, set LMODE accordingly.

The double differential cross sections of the Drell-Yan lepton pair ($\ell^+\ell^-$, dielectron or dimuon) production are measured, as functions of the invariant mass $m_{\ell\ell}$, transverse momentum $p_T^{\ell\ell}$, and $\varphi^*_\eta$. The $\varphi^*_\eta$ observable, derived from angular measurements of the leptons and highly correlated with $p_T^{\ell\ell}$, is used to probe the low-$p_T^{\ell\ell}$ region in a complementary way. Dilepton masses up to 1 TeV are investigated. Additionally, a measurement is performed requiring at least one jet in the final state. To benefit from partial cancellation of the systematic uncertainty, the ratios of the differential cross sections for various $m_{\ell\ell}$ ranges to those in the Z mass peak interval are presented. The collected data correspond to an integrated luminosity of 36.3 fb$^{-1}$ of proton-proton collisions recorded with the CMS detector at the LHC at a centre-of-mass energy of 13 TeV. Measurements are compared with predictions based on perturbative quantum chromodynamics, including soft-gluon resummation.

Source code: CMS_2022_I2079374.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FastJets.hh"
  4#include "Rivet/Projections/VisibleFinalState.hh"
  5#include "Rivet/Projections/FinalState.hh"
  6#include "Rivet/Projections/DressedLeptons.hh"
  7#include "Rivet/Projections/PromptFinalState.hh"
  8
  9namespace Rivet {
 10
 11  /// @brief Z pT over a wide mass range
 12  class CMS_2022_I2079374 : public Analysis {
 13    public:
 14
 15      /// Constructor
 16      RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2022_I2079374);
 17
 18      /// Book histograms and initialise projections before the run
 19      void init() {
 20        // Get options from the new option system
 21        // default to combined.
 22        _mode = 2;
 23        if (getOption("LMODE") == "EL") {
 24          _mode = 0;
 25        } else if (getOption("LMODE") == "MU") {
 26          _mode = 1;
 27        } else if (getOption("LMODE") == "EMU") {
 28          _mode = 2;
 29        } else {
 30          MSG_INFO("Assuming a mixed electron/muon sample. Set LMODE=EL/MU if your sample contains a single flavor.");
 31        }
 32
 33        FinalState fs;
 34        PromptFinalState pfs(fs);
 35
 36        PromptFinalState bareMuons(Cuts::abspid == PID::MUON);
 37        declare(DressedLeptons(pfs, bareMuons, 0.1, Cuts::abseta < 2.4 && Cuts::pT > 20*GeV, true),
 38                "muons");
 39
 40        PromptFinalState bareElectrons(Cuts::abspid == PID::ELECTRON);
 41        declare(DressedLeptons(pfs, bareElectrons, 0.1, Cuts::abseta < 2.4 && Cuts::pT > 20*GeV, true),
 42                "electrons");
 43
 44        FastJets jets(fs, FastJets::ANTIKT, 0.4);
 45        declare(jets, "jets");
 46
 47        // Histograms
 48        book(_pt_0jet_mass50_76, 1, 1, 1);
 49        book(_pt_0jet_mass76_106, 3, 1, 1);
 50        book(_pt_0jet_mass106_170, 5, 1, 1);
 51        book(_pt_0jet_mass170_350, 7, 1, 1);
 52        book(_pt_0jet_mass350_1000, 9, 1, 1);
 53
 54        book(_pt_1jet_mass50_76, 11, 1, 1);
 55        book(_pt_1jet_mass76_106, 13, 1, 1);
 56        book(_pt_1jet_mass106_170, 15, 1, 1);
 57        book(_pt_1jet_mass170_350, 17, 1, 1);
 58
 59        book(_phistar_mass50_76, 19, 1, 1);
 60        book(_phistar_mass76_106, 21, 1, 1);
 61        book(_phistar_mass106_170, 23, 1, 1);
 62        book(_phistar_mass170_350, 25, 1, 1);
 63        book(_phistar_mass350_1000, 27, 1, 1);
 64      }
 65
 66      /// Z boson finder.
 67      /// Note: we don't use the standard ZFinder class in order to stick to
 68      /// the definition of the publication that is simpler than the ZFinder
 69      /// algorithm
 70      /// @param leptons pt-ordered of electron or muon collection to use to build
 71      /// the Z boson
 72      std::unique_ptr<Particle> zfinder(const std::vector<DressedLepton> &leptons)
 73      {
 74        if (leptons.size() < 2) {
 75          return nullptr;
 76        }
 77        // Leading lepton pT cut
 78        if (leptons[0].pT() < 25.*GeV) {
 79          return nullptr;
 80        }
 81        if (leptons[0].charge() * leptons[1].charge() > 0) {
 82          return nullptr;
 83        }
 84
 85        std::unique_ptr<Particle> candidate(new Particle(PID::ZBOSON, leptons[0].mom() + leptons[1].mom()));
 86        if (candidate->mass() < 50*GeV || candidate->mass() > 1000*GeV) {
 87          return nullptr;
 88        }
 89        return candidate;
 90      }
 91
 92      /// Perform the per-event analysis
 93      void analyze(const Event& event)
 94      {
 95        // Fetch dressed leptons
 96        auto muons = apply<DressedLeptons>(event, "muons").dressedLeptons();
 97        auto electrons = apply<DressedLeptons>(event, "electrons").dressedLeptons();
 98
 99        const std::vector<DressedLepton> *dressedLeptons = nullptr;
100
101        //Look for Z->ee
102        std::unique_ptr<Particle> z = zfinder(electrons);
103        if (z != nullptr) {
104          dressedLeptons = &electrons;
105          if (_mode == 1) {
106            // Z->ee found, but running with LMODE=MU
107            vetoEvent;
108          }
109        } else if (_mode == 0) {
110          // LMODE=EL, but no Z->ee found
111          vetoEvent;
112        } else { //look for Z->mumu
113          z = zfinder(muons);
114          if (z.get() != nullptr) {
115            dressedLeptons = &muons;
116          } else { //no Z boson found
117            vetoEvent;
118          }
119        }
120
121        // Cluster jets
122        Jets jets = apply<FastJets>(event, "jets").jetsByPt(Cuts::absrap < 2.4 && Cuts::pT > 30*GeV);
123        idiscardIfAnyDeltaRLess(jets, *dressedLeptons, 0.4);
124
125        //------------------------------------------------
126
127        // Compute phi*
128        const Particle& lminus = dressedLeptons->at(0).charge() < 0 ? dressedLeptons->at(0) : dressedLeptons->at(1);
129        const Particle& lplus  = dressedLeptons->at(0).charge() < 0 ? dressedLeptons->at(1) : dressedLeptons->at(0);
130        double phi_acop = M_PI - deltaPhi(lminus, lplus);
131        double costhetastar = tanh( 0.5 * (lminus.eta() - lplus.eta()) );
132        double sin2thetastar = (costhetastar > 1) ? 0.0 : (1.0 - sqr(costhetastar));
133        double phistar = tan(0.5 * phi_acop) * sqrt(sin2thetastar);
134
135        if (z->mass() > 50. && z->mass() < 76.) {
136          _pt_0jet_mass50_76->fill(z->pT());
137          _phistar_mass50_76->fill(phistar);
138        } else if (z->mass() > 76. && z->mass() < 106.) {
139          _pt_0jet_mass76_106->fill(z->pT());
140          _phistar_mass76_106->fill(phistar);
141        } else if (z->mass() > 106. && z->mass() < 170.) {
142          _pt_0jet_mass106_170->fill(z->pT());
143          _phistar_mass106_170->fill(phistar);
144        } else if (z->mass() > 170. && z->mass() < 350.) {
145          _pt_0jet_mass170_350->fill(z->pT());
146          _phistar_mass170_350->fill(phistar);
147        } else if (z->mass() > 350. && z->mass() < 1000.) {
148          _pt_0jet_mass350_1000->fill(z->pT());
149          _phistar_mass350_1000->fill(phistar);
150        }
151
152        //---------------------- at least one jet --------------------------------
153        if (jets.size() < 1) {
154          return;
155        }
156
157        if (z->mass() > 50. && z->mass() < 76.) {
158          _pt_1jet_mass50_76->fill(z->pT());
159        } else if (z->mass() > 76. && z->mass() < 106.) {
160          _pt_1jet_mass76_106->fill(z->pT());
161        } else if (z->mass() > 106. && z->mass() < 170.) {
162          _pt_1jet_mass106_170->fill(z->pT());
163        } else if (z->mass() > 170. && z->mass() < 350.) {
164          _pt_1jet_mass170_350->fill(z->pT());
165        }
166      }
167
168      /// Calculates the ratio between two histograms
169      void calculateRatio(int d, const Histo1DPtr &numerator, const Histo1DPtr &denominator) {
170        Scatter2DPtr ratio;
171        book(ratio, d, 1, 1, true);
172
173        // Should be the case, but better check
174        assert(YODA::Histo1D(*ratio).sameBinning(*numerator));
175
176        // The denominator has a finer binning, so rebin it
177        auto rebinned_den = denominator->clone();
178        rebinned_den.rebin(numerator->xEdges());
179
180        divide(*numerator, rebinned_den, ratio);
181      }
182
183      /// Normalise histograms etc., after the run
184      void finalize() {
185        double norm = (sumOfWeights() != 0) ? crossSection()/sumOfWeights() : 1.0;
186
187        if (_mode == 2)  {
188          norm /= 2.;
189        }
190
191        scale(_pt_0jet_mass50_76, norm);
192        scale(_pt_0jet_mass76_106, norm);
193        scale(_pt_0jet_mass106_170, norm);
194        scale(_pt_0jet_mass170_350, norm);
195        scale(_pt_0jet_mass350_1000, norm);
196
197        scale(_pt_1jet_mass50_76, norm);
198        scale(_pt_1jet_mass76_106, norm);
199        scale(_pt_1jet_mass106_170, norm);
200        scale(_pt_1jet_mass170_350, norm);
201
202        scale(_phistar_mass50_76, norm);
203        scale(_phistar_mass76_106, norm);
204        scale(_phistar_mass106_170, norm);
205        scale(_phistar_mass170_350, norm);
206        scale(_phistar_mass350_1000, norm);
207
208        // Calculate the ratios
209        calculateRatio(29, _pt_0jet_mass50_76, _pt_0jet_mass76_106);
210        calculateRatio(31, _pt_0jet_mass106_170, _pt_0jet_mass76_106);
211        calculateRatio(33, _pt_0jet_mass170_350, _pt_0jet_mass76_106);
212        calculateRatio(35, _pt_0jet_mass350_1000, _pt_0jet_mass76_106);
213
214        calculateRatio(37, _pt_1jet_mass50_76, _pt_1jet_mass76_106);
215        calculateRatio(39, _pt_1jet_mass106_170, _pt_1jet_mass76_106);
216        calculateRatio(41, _pt_1jet_mass170_350, _pt_1jet_mass76_106);
217
218        calculateRatio(43, _phistar_mass50_76, _phistar_mass76_106);
219        calculateRatio(45, _phistar_mass106_170, _phistar_mass76_106);
220        calculateRatio(47, _phistar_mass170_350, _phistar_mass76_106);
221        calculateRatio(49, _phistar_mass350_1000, _phistar_mass76_106);
222      }
223
224    protected:
225      size_t _mode;
226
227    private:
228      /// @name Histograms
229      Histo1DPtr _pt_0jet_mass50_76,
230                 _pt_0jet_mass76_106,
231                 _pt_0jet_mass106_170,
232                 _pt_0jet_mass170_350,
233                 _pt_0jet_mass350_1000;
234      Histo1DPtr _pt_1jet_mass50_76,
235                 _pt_1jet_mass76_106,
236                 _pt_1jet_mass106_170,
237                 _pt_1jet_mass170_350;
238      Histo1DPtr _phistar_mass50_76,
239                 _phistar_mass76_106,
240                 _phistar_mass106_170,
241                 _phistar_mass170_350,
242                 _phistar_mass350_1000;
243    };
244
245  RIVET_DECLARE_PLUGIN(CMS_2022_I2079374);
246}