Rivet analyses referenceCMS_2019_I1753680Measurements of differential Z boson production cross sections in proton-proton collisions at 13 TeVExperiment: CMS (LHC) Inspire ID: 1753680 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
Measurements are presented of the differential cross sections for Z bosons produced in proton-proton collisions at $\sqrt{s} = 13 \text{TeV}$ and decaying to muons and electrons. The data analyzed were collected in 2016 with the CMS detector at the LHC and correspond to an integrated luminosity of $35.9/fb$. The measured fiducial inclusive product of cross section and branching fraction agrees with next-to-next-to-leading order quantum chromodynamics calculations. Differential cross sections of the transverse momentum \pt, the optimized angular variable $\phi^{*}_{\eta}$, and the rapidity of lepton pairs are measured. The data are corrected for detector effects and compared to theoretical predictions using fixed order, resummed, and parton shower calculations. The uncertainties of the measured normalized cross sections are smaller than 0.5% for $\phi^{*}_{\eta} < 0.5$ and for $p_{T}^{Z} < 50 \text{GeV}$. Source code: CMS_2019_I1753680.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/FastJets.hh"
5#include "Rivet/Projections/DressedLeptons.hh"
6#include "Rivet/Projections/ZFinder.hh"
7
8namespace Rivet {
9
10
11 /// @brief Measurements of differential Z boson production cross sections in proton-proton collisions at 13 TeV
12 class CMS_2019_I1753680 : public Analysis {
13 public:
14
15 /// Constructor
16 RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2019_I1753680);
17
18
19 /// @name Analysis methods
20 //@{
21
22 /// Book histograms and initialise projections before the run
23 void init() {
24
25 // Get options from the new option system
26 // default to combined.
27 _mode = 2;
28 if ( getOption("LMODE") == "EL" ) _mode = 0;
29 if ( getOption("LMODE") == "MU" ) _mode = 1;
30 if ( getOption("LMODE") == "EMU" ) _mode = 2;
31
32 // Initialise and register projections
33 FinalState fs;
34 Cut cut = Cuts::abseta < 2.4 && Cuts::pT > 25*GeV;
35
36 ZFinder zeeFind(fs, cut, PID::ELECTRON, 76.1876*GeV, 106.1876*GeV, 0.1, ZFinder::ChargedLeptons::PROMPT, ZFinder::ClusterPhotons::NODECAY, ZFinder::AddPhotons::YES );
37 declare(zeeFind, "ZeeFind");
38 ZFinder zmmFind(fs, cut, PID::MUON , 76.1876*GeV, 106.1876*GeV, 0.1, ZFinder::ChargedLeptons::PROMPT, ZFinder::ClusterPhotons::NODECAY, ZFinder::AddPhotons::YES );
39 declare(zmmFind, "ZmmFind");
40
41 // Book histograms
42 book(_h_Zmm_absY , 26, 1, 1);
43 book(_h_Zee_absY , 26, 1, 2);
44 book(_h_Zll_absY , 26, 1, 3);
45 book(_h_Zmm_pt , 27, 1, 1);
46 book(_h_Zee_pt , 27, 1, 2);
47 book(_h_Zll_pt , 27, 1, 3);
48 book(_h_Zmm_phiStar , 28, 1, 1);
49 book(_h_Zee_phiStar , 28, 1, 2);
50 book(_h_Zll_phiStar , 28, 1, 3);
51 book(_h_Zll_pt_Y0 , 29, 1, 1);
52 book(_h_Zll_pt_Y1 , 29, 1, 2);
53 book(_h_Zll_pt_Y2 , 29, 1, 3);
54 book(_h_Zll_pt_Y3 , 29, 1, 4);
55 book(_h_Zll_pt_Y4 , 29, 1, 5);
56
57 book(_h_Zll_pt_norm , 30, 1, 1);
58 book(_h_Zll_phiStar_norm , 31, 1, 1);
59 book(_h_Zll_absY_norm , 32, 1, 1);
60 book(_h_Zll_pt_Y0_norm , 33, 1, 1);
61 book(_h_Zll_pt_Y1_norm , 33, 1, 2);
62 book(_h_Zll_pt_Y2_norm , 33, 1, 3);
63 book(_h_Zll_pt_Y3_norm , 33, 1, 4);
64 book(_h_Zll_pt_Y4_norm , 33, 1, 5);
65
66 }
67
68
69 /// Perform the per-event analysis
70 void analyze(const Event& event) {
71
72 const ZFinder& zeeFS = apply<ZFinder>(event, "ZeeFind");
73 const ZFinder& zmumuFS = apply<ZFinder>(event, "ZmmFind");
74
75 const Particles& zees = zeeFS.bosons();
76 const Particles& zmumus = zmumuFS.bosons();
77
78 if (zees.size() + zmumus.size() != 1) {
79 MSG_DEBUG("Did not find exactly one good Z candidate");
80 vetoEvent;
81 }
82
83 //event identification depending on mass window
84 bool ee_event=false;
85 bool mm_event=false;
86
87 if (zees.size() == 1) {
88 ee_event = true;
89 }
90 if (zmumus.size() == 1) {
91 mm_event = true;
92 }
93
94 if (ee_event && _mode == 1)
95 vetoEvent;
96 if (mm_event && _mode == 0)
97 vetoEvent;
98
99 const Particles& theLeptons = ee_event ? zeeFS.constituents() : zmumuFS.constituents();
100 const Particle& lminus = theLeptons[0].charge() < 0 ? theLeptons[0] : theLeptons[1];
101 const Particle& lplus = theLeptons[0].charge() < 0 ? theLeptons[1] : theLeptons[0];
102
103 //calculate phi*
104 const double thetaStar = acos(tanh( 0.5 * (lminus.eta() - lplus.eta()) ));
105 const double dPhi = M_PI - deltaPhi(lminus, lplus);
106 const double phiStar = tan(0.5 * dPhi) * sin(thetaStar);
107
108 const Particle& zcand = ee_event ? zees[0] : zmumus[0];
109
110 if (ee_event) {
111 _h_Zee_absY->fill(zcand.absrap());
112 _h_Zee_pt->fill(zcand.pt());
113 _h_Zee_phiStar->fill(phiStar);
114 }
115 else if (mm_event) {
116 _h_Zmm_absY->fill(zcand.absrap());
117 _h_Zmm_pt->fill(zcand.pt());
118 _h_Zmm_phiStar->fill(phiStar);
119 }
120
121 _h_Zll_pt->fill(zcand.pt());
122 _h_Zll_pt_norm->fill(zcand.pt());
123 _h_Zll_phiStar->fill(phiStar);
124 _h_Zll_phiStar_norm->fill(phiStar);
125 _h_Zll_absY->fill(zcand.absrap());
126 _h_Zll_absY_norm->fill(zcand.absrap());
127
128 if (zcand.absrap()<0.4) {
129 _h_Zll_pt_Y0->fill(zcand.pt());
130 _h_Zll_pt_Y0_norm->fill(zcand.pt());
131 }
132 else if (zcand.absrap()<0.8) {
133 _h_Zll_pt_Y1->fill(zcand.pt());
134 _h_Zll_pt_Y1_norm->fill(zcand.pt());
135 }
136 else if (zcand.absrap()<1.2) {
137 _h_Zll_pt_Y2->fill(zcand.pt());
138 _h_Zll_pt_Y2_norm->fill(zcand.pt());
139 }
140 else if (zcand.absrap()<1.6) {
141 _h_Zll_pt_Y3->fill(zcand.pt());
142 _h_Zll_pt_Y3_norm->fill(zcand.pt());
143 }
144 else if (zcand.absrap()<2.4) {
145 _h_Zll_pt_Y4->fill(zcand.pt());
146 _h_Zll_pt_Y4_norm->fill(zcand.pt());
147 }
148
149 }
150
151 void normalizeToSum(Histo1DPtr hist) {
152 double sum = 0.;
153 for (size_t i = 0; i < hist->numBins(); ++i) {
154 sum += hist->bin(i).height();
155 }
156 scale(hist, 1./sum);
157 }
158
159 /// Normalise histograms etc., after the run
160 void finalize() {
161
162 double norm = (sumOfWeights() != 0) ? crossSection()/picobarn/sumOfWeights() : 1.0;
163
164 scale(_h_Zmm_pt, norm);
165 scale(_h_Zmm_absY, norm);
166 scale(_h_Zmm_phiStar, norm);
167
168 scale(_h_Zee_pt, norm);
169 scale(_h_Zee_absY, norm);
170 scale(_h_Zee_phiStar, norm);
171
172 // when running in combined mode, need to average to get lepton xsec
173 if (_mode == 2) norm /= 2.;
174
175 scale(_h_Zll_pt, norm);
176 scale(_h_Zll_absY, norm);
177 scale(_h_Zll_phiStar, norm);
178 scale(_h_Zll_pt_Y0, norm);
179 scale(_h_Zll_pt_Y1, norm);
180 scale(_h_Zll_pt_Y2, norm);
181 scale(_h_Zll_pt_Y3, norm);
182 scale(_h_Zll_pt_Y4, norm);
183
184 normalizeToSum(_h_Zll_pt_norm);
185 normalizeToSum(_h_Zll_absY_norm);
186 normalizeToSum(_h_Zll_phiStar_norm);
187 normalizeToSum(_h_Zll_pt_Y0_norm);
188 normalizeToSum(_h_Zll_pt_Y1_norm);
189 normalizeToSum(_h_Zll_pt_Y2_norm);
190 normalizeToSum(_h_Zll_pt_Y3_norm);
191 normalizeToSum(_h_Zll_pt_Y4_norm);
192
193 }
194
195 //@}
196
197 protected:
198
199 size_t _mode;
200
201 /// @name Histograms
202
203 private:
204
205 Histo1DPtr _h_Zmm_pt, _h_Zmm_phiStar, _h_Zmm_absY;
206 Histo1DPtr _h_Zee_pt, _h_Zee_phiStar, _h_Zee_absY;
207
208 Histo1DPtr _h_Zll_pt, _h_Zll_phiStar, _h_Zll_absY;
209 Histo1DPtr _h_Zll_pt_Y0, _h_Zll_pt_Y1, _h_Zll_pt_Y2, _h_Zll_pt_Y3, _h_Zll_pt_Y4;
210
211 Histo1DPtr _h_Zll_pt_norm, _h_Zll_phiStar_norm, _h_Zll_absY_norm;
212 Histo1DPtr _h_Zll_pt_Y0_norm, _h_Zll_pt_Y1_norm, _h_Zll_pt_Y2_norm, _h_Zll_pt_Y3_norm, _h_Zll_pt_Y4_norm;
213
214 };
215
216
217 // The hook for the plugin system
218 RIVET_DECLARE_PLUGIN(CMS_2019_I1753680);
219
220
221}
|