Rivet analyses referenceCMS_2017_I1635889Measurement of the underlying event using inclusive Z-boson production in proton-proton collisions at √s=13 TeVExperiment: CMS (LHC) Inspire ID: 1635889 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
Measurement of the underlying event activity, in proton-proton collisions at a centre-of-mass energy of 13TeV, performed using inclusive Z-boson events collected by the CMS experiment at the Large Hadron Collider. The analyzed data corresponds to an integrated luminosity of 2.1fb−1. The underlying event activity is quantified in terms of charged particle multiplicity and their scalar sum of their transverse momentum in different topological regions defined with respect to the resultatnt azimuthal direction of the two muons coming from the decay of a Z-boson. The distributions are unfolded to stable charged-particle level and are compared with predictions from various simulations, as well as with the measurements at a centre-of-mass energy of 1.96TeV and 7TeV by the CDF and CMS experiments. Source code: CMS_2017_I1635889.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/DileptonFinder.hh"
5#include "Rivet/Projections/FinalState.hh"
6#include "Rivet/Projections/ChargedFinalState.hh"
7#include "Rivet/Projections/VetoedFinalState.hh"
8
9namespace Rivet {
10
11
12 /// Underlying event activity in the Drell-Yan process at 13 TeV
13 class CMS_2017_I1635889 : public Analysis {
14 public:
15
16 /// Constructor
17 RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2017_I1635889);
18
19
20 /// Initialization
21 void init() {
22
23 /// @note Using a bare muon Z (but with a clustering radius!?)
24 Cut cut = Cuts::abseta < 2.4 && Cuts::pT > 10*GeV;
25 DileptonFinder zfinder(91.2*GeV, 0.2, cut && Cuts::abspid == PID::MUON, Cuts::massIn(81*GeV, 101*GeV));
26 declare(zfinder, "DileptonFinder");
27
28 ChargedFinalState cfs(zfinder.remainingFinalState());
29 declare(cfs, "cfs");
30
31 book(_h_Nchg_towards_pTmumu , 1, 1, 1);
32 book(_h_Nchg_transverse_pTmumu , 2, 1, 1);
33 book(_h_Nchg_away_pTmumu , 3, 1, 1);
34 book(_h_pTsum_towards_pTmumu , 4, 1, 1);
35 book(_h_pTsum_transverse_pTmumu , 5, 1, 1);
36 book(_h_pTsum_away_pTmumu , 6, 1, 1);
37 }
38
39
40 /// Perform the per-event analysis
41 void analyze(const Event& event) {
42 const DileptonFinder& zfinder = apply<DileptonFinder>(event, "DileptonFinder");
43
44 if (zfinder.bosons().size() != 1) vetoEvent;
45 if (zfinder.constituents()[0].pT()<20 && zfinder.constituents()[1].pT()<20)vetoEvent;
46 //std::cout<<"pt[0] = "<<zfinder.constituents()[0].pT()<<"pt[1] = "<<zfinder.constituents()[1].pT()<<std::endl;
47 double Zpt = zfinder.bosons()[0].pT()/GeV;
48 double Zphi = zfinder.bosons()[0].phi();
49 //double Zmass = zfinder.bosons()[0].mass()/GeV;
50
51 Particles particles = apply<ChargedFinalState>(event, "cfs").particlesByPt(Cuts::pT > 0.5*GeV && Cuts::abseta <2.0);
52
53 int nTowards = 0;
54 int nTransverse = 0;
55 int nAway = 0;
56 double ptSumTowards = 0;
57 double ptSumTransverse = 0;
58 double ptSumAway = 0;
59
60 for (const Particle& p : particles) {
61 double dphi = fabs(deltaPhi(Zphi, p.phi()));
62 double pT = p.pT();
63
64 if ( dphi < M_PI/3 ) {
65 nTowards++;
66 ptSumTowards += pT;
67 } else if ( dphi < 2.*M_PI/3 ) {
68 nTransverse++;
69 ptSumTransverse += pT;
70 } else {
71 nAway++;
72 ptSumAway += pT;
73 }
74
75 } // Loop over particles
76
77
78 const double area = 8./3.*M_PI;
79 _h_Nchg_towards_pTmumu-> fill(Zpt, 1./area * nTowards);
80 _h_Nchg_transverse_pTmumu-> fill(Zpt, 1./area * nTransverse);
81 _h_Nchg_away_pTmumu-> fill(Zpt, 1./area * nAway);
82 _h_pTsum_towards_pTmumu-> fill(Zpt, 1./area * ptSumTowards);
83 _h_pTsum_transverse_pTmumu-> fill(Zpt, 1./area * ptSumTransverse);
84 _h_pTsum_away_pTmumu-> fill(Zpt, 1./area * ptSumAway);
85
86
87 }
88
89
90 /// Normalise histograms etc., after the run
91 void finalize() {
92 }
93
94
95 private:
96
97
98 /// @name Histogram objects
99 /// @{
100 Profile1DPtr _h_Nchg_towards_pTmumu;
101 Profile1DPtr _h_Nchg_transverse_pTmumu;
102 Profile1DPtr _h_Nchg_away_pTmumu;
103 Profile1DPtr _h_pTsum_towards_pTmumu;
104 Profile1DPtr _h_pTsum_transverse_pTmumu;
105 Profile1DPtr _h_pTsum_away_pTmumu;
106 /// @}
107
108 };
109
110
111 RIVET_DECLARE_PLUGIN(CMS_2017_I1635889);
112
113}
|