Rivet analyses referenceCMS_2017_I1631985Measurement of differential cross sections in the kinematic angular variable $\phi$* for inclusive Z boson production in pp collisions at 8 TeVExperiment: CMS (LHC) Inspire ID: 1631985 Status: VALIDATED Authors:
Beam energies: (4000.0, 4000.0) GeV Run details:
Measurements of differential cross sections d$\sigma$/d$\phi^*$ and double-differential cross sections d$\sigma^{2}$/d$\phi^*$d|y| for inclusive Z boson production are presented using the dielectron and dimuon final states at the born level. The kinematic observable $\phi$* correlates with the dilepton transverse momentum but has better resolution, and y is the dilepton rapidity. The analysis is based on data collected with the CMS experiment at a centre-of-mass energy of 8 TeV corresponding to an integrated luminosity of 19.7 fb$^{-1}$. The normalised cross section (1/$\sigma$) d$\sigma$/d$\phi^*$, within the fiducial kinematic region, is measured with a precision of better than 0.5% for $\phi$* < 1. The measurements are compared to theoretical predictions and they agree, typically, within few percent. Source code: CMS_2017_I1631985.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/FastJets.hh"
5#include "Rivet/Projections/LeptonFinder.hh"
6#include "Rivet/Projections/DileptonFinder.hh"
7
8namespace Rivet {
9
10
11 /// @brief Differential Z cross section measurement in phi* at 8 TeV
12 class CMS_2017_I1631985 : public Analysis {
13 public:
14
15 /// Constructor
16 RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2017_I1631985);
17
18
19 /// @name Analysis methods
20 /// @{
21
22 /// Book histograms and initialise projections before the run
23 void init() {
24 // default to combined ee+mumu
25 _mode = 2;
26 if ( getOption("LMODE") == "EL" ) _mode = 0;
27 if ( getOption("LMODE") == "MU" ) _mode = 1;
28 if ( getOption("LMODE") == "EMU" ) _mode = 2;
29
30 _twodim = false;
31 if ( getOption("TWODIM") == "YES" ) _twodim = true;
32
33 // The experimental result is unfolded to the 'born' level (pre-FSR), which is discouraged
34 // The implementation in Rivet here is using dressed leptons (or disable generator QED FSR)
35 Cut cut = Cuts::abseta < 2.4 && Cuts::pT > 20*GeV;
36 DileptonFinder zeeFind(91.2*GeV, 0.1, cut && Cuts::abspid == PID::ELECTRON, Cuts::massIn(60*GeV, 120*GeV));
37 declare(zeeFind, "ZeeFind");
38 DileptonFinder zmmFind(91.2*GeV, 0.1, cut && Cuts::abspid == PID::MUON , Cuts::massIn(60*GeV, 120*GeV));
39 declare(zmmFind, "ZmmFind");
40
41 // Book histograms
42 // take binning from reference data using HEPData ID (digits in "d01-x01-y01" etc.)
43 book(_h_Zll_phiStar, 1, 1, 1);
44 book(_h_Zll_phiStar_norm, 2, 1, 1);
45 if (_twodim) {
46 book(_h2D_Zll_phiStar_y, 3, 1, 1);
47 book(_h2D_Zll_phiStar_y_norm, 4, 1, 1);
48 }
49
50 vector<double> edges = {0.0, 0.4, 0.8, 1.2, 1.6, 2.0, 2.4};
51 book(_b_Zll_phiStar, edges);
52 book(_b_Zll_phiStar_norm, edges);
53 for (size_t i = 0; i < _b_Zll_phiStar->numBins(); ++i) {
54 book(_b_Zll_phiStar->bin(i+1), 9, 1, i+1);
55 book(_b_Zll_phiStar_norm->bin(i+1), 10, 1, i+1);
56 }
57 }
58
59
60 /// Perform the per-event analysis
61 void analyze(const Event& event) {
62
63 const DileptonFinder& zeeFS = apply<DileptonFinder>(event, "ZeeFind");
64 const DileptonFinder& zmumuFS = apply<DileptonFinder>(event, "ZmmFind");
65
66 const Particles& zees = zeeFS.bosons();
67 const Particles& zmumus = zmumuFS.bosons();
68
69 if (zees.size() + zmumus.size() != 1) {
70 MSG_DEBUG("Did not find exactly one good Z candidate");
71 vetoEvent;
72 }
73 if (zees.size() == 1 && _mode == 1)
74 vetoEvent;
75 if (zmumus.size() == 1 && _mode == 0)
76 vetoEvent;
77
78
79 bool ee_event=false;
80 if( zees.size() == 1 ) ee_event = true;
81
82 const Particles& theLeptons = ee_event ? zeeFS.constituents() : zmumuFS.constituents();
83 const Particle& leadingLepton = theLeptons[0].pt() > theLeptons[1].pt() ? theLeptons[0] : theLeptons[1];
84
85 // asymmetric cut
86 if( leadingLepton.pt() < 30.0 ) vetoEvent;
87 if( leadingLepton.abseta() > 2.1 ) vetoEvent;
88
89 // calculate phi*
90 const Particle& lminus = theLeptons[0].charge() < 0 ? theLeptons[0] : theLeptons[1];
91 const Particle& lplus = theLeptons[0].charge() < 0 ? theLeptons[1] : theLeptons[0];
92 const double thetaStar = acos(tanh( 0.5 * (lminus.eta() - lplus.eta()) ));
93 const double dPhi = M_PI - deltaPhi(lminus, lplus);
94 const double phiStar = tan(0.5 * dPhi) * sin(thetaStar);
95
96 const Particle& zcand = ee_event ? zees[0] : zmumus[0];
97
98 _h_Zll_phiStar->fill(phiStar);
99 _h_Zll_phiStar_norm->fill(phiStar);
100
101 double absRap = zcand.absrap();
102 if (_twodim) {
103 _h2D_Zll_phiStar_y->fill(phiStar, absRap);
104 _h2D_Zll_phiStar_y_norm->fill(phiStar, absRap);
105 }
106
107 _b_Zll_phiStar->fill(absRap, phiStar);
108 _b_Zll_phiStar_norm->fill(absRap, phiStar);
109
110 }
111
112 /// Normalise histograms etc., after the run
113 void finalize() {
114 double norm = (sumOfWeights() != 0) ? crossSection()/picobarn/sumOfWeights() : 1.0;
115
116 // when running for both ee and mm channel, need to average to get lepton xsec
117 if (_mode == 2) norm /= 2.;
118
119 scale(_h_Zll_phiStar, norm);
120 scale(_b_Zll_phiStar, norm);
121 divByGroupWidth(_b_Zll_phiStar);
122 normalize(_h_Zll_phiStar_norm);
123
124 if (_twodim) {
125 scale(_h2D_Zll_phiStar_y, norm);
126 normalize(_h2D_Zll_phiStar_y_norm);
127 }
128
129 // normalized using the sum of 2D
130 normalizeGroup(_b_Zll_phiStar_norm);
131 divByGroupWidth(_b_Zll_phiStar_norm);
132 }
133
134 /// @}
135
136 Histo1DPtr _h_Zll_phiStar, _h_Zll_phiStar_norm;
137 Histo1DGroupPtr _b_Zll_phiStar, _b_Zll_phiStar_norm;
138 Histo2DPtr _h2D_Zll_phiStar_y, _h2D_Zll_phiStar_y_norm;
139
140 size_t _mode;
141 bool _twodim;
142 };
143
144
145 RIVET_DECLARE_PLUGIN(CMS_2017_I1631985);
146
147}
|