Rivet analyses referenceATLAS_2016_I1469071Measurement of the $WZ$ production cross section at 13 TeVExperiment: ATLAS (LHC) Inspire ID: 1469071 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
The production of $W^\pm Z$ events in proton--proton collisions at a centre-of-mass energy of 13 TeV is measured with the ATLAS detector at the LHC. The collected data correspond to an integrated luminosity of 3.2 fb${}^{-1}$ . The $W^\pm Z$ candidates are reconstructed using leptonic decays of the gauge bosons into electrons or muons. The measured inclusive cross section in the detector fiducial region for leptonic decay modes is $\sigma(W^\pm Z\to\ell^\prime\nu\ell\ell \text{fid.}) = 63.2\pm 3.2$(stat.)$\pm 2.6$(sys.)$\pm 1.5$(lumi.) fb. In comparison, the next-to-leading-order Standard Model prediction is 53.4-2.8+3.6 fb. The extrapolation of the measurement from the fiducial to the total phase space yields $\sigma(W^\pm Z\text{tot.})=50.6\pm 2.6$(stat.)$\pm 2.0$(sys.)$\pm 0.9$(th.)$\pm 1.2$(lumi.) pb, in agreement with a recent next-to-next-to-leading-order calculation of 48.2-1.0+1.1 pb. The cross section as a function of jet multiplicity is also measured, together with the charge-dependent $W^+Z$ and $W^-Z$ cross sections and their ratio. Uses SM neutrino-lepton flavour matching and a resonant shape algorithm assuming the Standard Model, to match the MC-based correction to the fiducial region applied in the paper. This routine is therefore only valid under the assumption of the Standard Model and cannot be used for BSM reinterpretation. Source code: ATLAS_2016_I1469071.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/FastJets.hh"
5#include "Rivet/Projections/IdentifiedFinalState.hh"
6#include "Rivet/Projections/PromptFinalState.hh"
7#include "Rivet/Projections/LeptonFinder.hh"
8#include "Rivet/Projections/VetoedFinalState.hh"
9
10namespace Rivet {
11
12
13 /// @brief Measurement of the WZ production cross section at 13 TeV
14 class ATLAS_2016_I1469071 : public Analysis {
15 public:
16
17 /// Constructor
18 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2016_I1469071);
19
20 /// @name Analysis methods
21 /// @{
22
23 /// Book histograms and initialise projections before the run
24 void init() {
25
26 // Lepton cuts
27 Cut FS_Zlept = Cuts::abseta < 2.5 && Cuts::pT > 15*GeV;
28
29 FinalState fs;
30 Cut fs_z = Cuts::abseta < 2.5 && Cuts::pT > 15*GeV;
31 Cut fs_j = Cuts::abseta < 4.5 && Cuts::pT > 25*GeV;
32
33 // Get photons to dress leptons
34 PromptFinalState photons(Cuts::abspid == PID::PHOTON);
35
36 // Electrons and muons in Fiducial PS
37 PromptFinalState leptons(FinalState(fs_z && (Cuts::abspid == PID::ELECTRON || Cuts::abspid == PID::MUON)));
38 leptons.acceptTauDecays(false);
39 LeptonFinder dressedleptons(leptons, photons, 0.1, FS_Zlept);
40 declare(dressedleptons, "LeptonFinder");
41
42 // Electrons and muons in Total PS
43 PromptFinalState leptons_total(Cuts::abspid == PID::ELECTRON || Cuts::abspid == PID::MUON);
44 leptons_total.acceptTauDecays(false);
45 LeptonFinder dressedleptonsTotal(leptons_total, photons, 0.1);
46 declare(dressedleptonsTotal, "LeptonFinderTotal");
47
48 // Promot neutrinos (yikes!)
49 IdentifiedFinalState nu_id;
50 nu_id.acceptNeutrinos();
51 PromptFinalState neutrinos(nu_id);
52 neutrinos.acceptTauDecays(false);
53 declare(neutrinos, "Neutrinos");
54 MSG_WARNING("\033[91;1mLIMITED VALIDITY - check info file for details!\033[m");
55
56 // Jets
57 VetoedFinalState veto;
58 veto.addVetoOnThisFinalState(dressedleptons);
59 FastJets jets(veto, JetAlg::ANTIKT, 0.4);
60 declare(jets, "Jets");
61
62 // Book histograms
63 book(_h["eee"] , 1, 1, 1);
64 book(_h["mee"] , 1, 1, 2);
65 book(_h["emm"] , 1, 1, 3);
66 book(_h["mmm"] , 1, 1, 4);
67 book(_h["fid"] , 1, 1, 5);
68 book(_h["eee_Plus"] , 2, 1, 1);
69 book(_h["mee_Plus"] , 2, 1, 2);
70 book(_h["emm_Plus"] , 2, 1, 3);
71 book(_h["mmm_Plus"] , 2, 1, 4);
72 book(_h["fid_Plus"] , 2, 1, 5);
73 book(_h["eee_Minus"], 3, 1, 1);
74 book(_h["mee_Minus"], 3, 1, 2);
75 book(_h["emm_Minus"], 3, 1, 3);
76 book(_h["mmm_Minus"], 3, 1, 4);
77 book(_h["fid_Minus"], 3, 1, 5);
78 book(_htot, 6, 1, 1);
79 book(_hnjets, 8, 1, 1);
80
81 }
82
83
84 void analyze(const Event& event) {
85 const DressedLeptons& dressedleptons = apply<LeptonFinder>(event, "LeptonFinder").dressedLeptons();
86 const DressedLeptons& dressedleptonsTotal = apply<LeptonFinder>(event, "LeptonFinderTotal").dressedLeptons();
87 const Particles& neutrinos = apply<PromptFinalState>(event, "Neutrinos").particlesByPt();
88 Jets jets = apply<JetFinder>(event, "Jets").jetsByPt( (Cuts::abseta < 4.5) && (Cuts::pT > 25*GeV) );
89
90 if (dressedleptonsTotal.size() < 3 || neutrinos.size() < 1) vetoEvent;
91
92 //---Total PS: assign leptons to W and Z bosons using Resonant shape algorithm
93 // NB: This resonant shape algorithm assumes the Standard Model and can therefore
94 // NOT be used for any kind of reinterpretation in terms of new-physics models..
95
96 // Try Z pair of leptons 01
97 double MassZ01 = 0, MassW2 = 0;
98 if ( (dressedleptonsTotal[0].pid() ==-(dressedleptonsTotal[1].pid())) && (dressedleptonsTotal[2].abspid()==neutrinos[0].abspid()-1)){
99 MassZ01 = (dressedleptonsTotal[0].momentum()+dressedleptonsTotal[1].momentum()).mass();
100 MassW2 = (dressedleptonsTotal[2].momentum()+neutrinos[0].momentum()).mass();
101 }
102 // Try Z pair of leptons 02
103 double MassZ02 = 0, MassW1 = 0;
104 if ( (dressedleptonsTotal[0].pid()==-(dressedleptonsTotal[2].pid())) && (dressedleptonsTotal[1].abspid()==neutrinos[0].abspid()-1)){
105 MassZ02 = (dressedleptonsTotal[0].momentum()+dressedleptonsTotal[2].momentum()).mass();
106 MassW1 = (dressedleptonsTotal[1].momentum()+neutrinos[0].momentum()).mass();
107 }
108 // Try Z pair of leptons 12
109 double MassZ12 = 0, MassW0 = 0;
110 if ( (dressedleptonsTotal[1].pid()==-(dressedleptonsTotal[2].pid())) && (dressedleptonsTotal[0].abspid()==neutrinos[0].abspid()-1)){
111 MassZ12 = (dressedleptonsTotal[1].momentum()+dressedleptonsTotal[2].momentum()).mass();
112 MassW0 = (dressedleptonsTotal[0].momentum()+neutrinos[0].momentum()).mass();
113 }
114 double WeightZ1 = 1/(pow(MassZ01*MassZ01 - MZ_PDG*MZ_PDG,2) + pow(MZ_PDG*GammaZ_PDG,2));
115 double WeightW1 = 1/(pow(MassW2*MassW2 - MW_PDG*MW_PDG,2) + pow(MW_PDG*GammaW_PDG,2));
116 double WeightTotal1 = WeightZ1*WeightW1;
117 double M1 = -1*WeightTotal1;
118
119 double WeightZ2 = 1/(pow(MassZ02*MassZ02- MZ_PDG*MZ_PDG,2) + pow(MZ_PDG*GammaZ_PDG,2));
120 double WeightW2 = 1/(pow(MassW1*MassW1- MW_PDG*MW_PDG,2) + pow(MW_PDG*GammaW_PDG,2));
121 double WeightTotal2 = WeightZ2*WeightW2;
122 double M2 = -1*WeightTotal2;
123
124 double WeightZ3 = 1/(pow(MassZ12*MassZ12 - MZ_PDG*MZ_PDG,2) + pow(MZ_PDG*GammaZ_PDG,2));
125 double WeightW3 = 1/(pow(MassW0*MassW0 - MW_PDG*MW_PDG,2) + pow(MW_PDG*GammaW_PDG,2));
126 double WeightTotal3 = WeightZ3*WeightW3;
127 double M3 = -1*WeightTotal3;
128
129 int i = -1, j = -1, k = -1;
130 if ((M1 < M2 && M1 < M3) || (MassZ01 != 0 && MassW2 != 0 && MassZ02 == 0 && MassZ12 == 0)) {
131 i = 0; j = 1; k = 2;
132 }
133 if ((M2 < M1 && M2 < M3) || (MassZ02 != 0 && MassW1 != 0 && MassZ01 == 0 && MassZ12 == 0)) {
134 i = 0; j = 2; k = 1;
135 }
136 if ((M3 < M1 && M3 < M2) || (MassZ12 != 0 && MassW0 != 0 && MassZ01 == 0 && MassZ02 == 0)) {
137 i = 1; j = 2; k = 0;
138 }
139 ///
140 if (i < 0 || j < 0 || k < 0) vetoEvent;
141
142 FourMomentum ZbosonTotal = dressedleptonsTotal[i].momentum()+dressedleptonsTotal[j].momentum();
143 if ( ZbosonTotal.mass() >= 66*GeV && ZbosonTotal.mass() <= 116*GeV ) _htot->fill(13000);
144
145 //---end Total PS
146
147
148 //---Fiducial PS: assign leptons to W and Z bosons using Resonant shape algorithm
149 if (dressedleptons.size() < 3) vetoEvent;
150
151 int EventType = -1;
152 int Nel = 0, Nmu = 0;
153
154 for (const DressedLepton& l : dressedleptons) {
155 if (l.abspid() == 11) ++Nel;
156 if (l.abspid() == 13) ++Nmu;
157 }
158
159 if ( (Nel == 3) && (Nmu==0) ) EventType = 3;
160 if ( (Nel == 2) && (Nmu==1) ) EventType = 2;
161 if ( (Nel == 1) && (Nmu==2) ) EventType = 1;
162 if ( (Nel == 0) && (Nmu==3) ) EventType = 0;
163
164 int EventCharge = -dressedleptons[0].charge() * dressedleptons[1].charge() * dressedleptons[2].charge();
165
166 MassZ01 = 0; MassZ02 = 0; MassZ12 = 0;
167 MassW0 = 0; MassW1 = 0; MassW2 = 0;
168
169 // try Z pair of leptons 01
170 if (dressedleptons[0].pid() == -dressedleptons[1].pid()) {
171 MassZ01 = (dressedleptons[0].momentum() + dressedleptons[1].momentum()).mass();
172 MassW2 = (dressedleptons[2].momentum() + neutrinos[0].momentum()).mass();
173 }
174 // try Z pair of leptons 02
175 if (dressedleptons[0].pid() == -dressedleptons[2].pid()) {
176 MassZ02 = (dressedleptons[0].momentum() + dressedleptons[2].momentum()).mass();
177 MassW1 = (dressedleptons[1].momentum() + neutrinos[0].momentum()).mass();
178 }
179 // try Z pair of leptons 12
180 if (dressedleptons[1].pid() == -dressedleptons[2].pid()) {
181 MassZ12 = (dressedleptons[1].momentum() + dressedleptons[2].momentum()).mass();
182 MassW0 = (dressedleptons[0].momentum() + neutrinos[0].momentum()).mass();
183 }
184 WeightZ1 = 1/(pow(MassZ01*MassZ01 - MZ_PDG*MZ_PDG,2) + pow(MZ_PDG*GammaZ_PDG,2));
185 WeightW1 = 1/(pow(MassW2*MassW2 - MW_PDG*MW_PDG,2) + pow(MW_PDG*GammaW_PDG,2));
186 WeightTotal1 = WeightZ1*WeightW1;
187 M1 = -1*WeightTotal1;
188
189 WeightZ2 = 1/(pow(MassZ02*MassZ02- MZ_PDG*MZ_PDG,2) + pow(MZ_PDG*GammaZ_PDG,2));
190 WeightW2 = 1/(pow(MassW1*MassW1- MW_PDG*MW_PDG,2) + pow(MW_PDG*GammaW_PDG,2));
191 WeightTotal2 = WeightZ2*WeightW2;
192 M2 = -1*WeightTotal2;
193
194 WeightZ3 = 1/(pow(MassZ12*MassZ12 - MZ_PDG*MZ_PDG,2) + pow(MZ_PDG*GammaZ_PDG,2));
195 WeightW3 = 1/(pow(MassW0*MassW0 - MW_PDG*MW_PDG,2) + pow(MW_PDG*GammaW_PDG,2));
196 WeightTotal3 = WeightZ3*WeightW3;
197 M3 = -1*WeightTotal3;
198
199 if( (M1 < M2 && M1 < M3) || (MassZ01 != 0 && MassW2 != 0 && MassZ02 == 0 && MassZ12 == 0) ) {
200 i = 0; j = 1; k = 2;
201 }
202 if((M2 < M1 && M2 < M3) || (MassZ02 != 0 && MassW1 != 0 && MassZ01 == 0 && MassZ12 == 0) ) {
203 i = 0; j = 2; k = 1;
204 }
205 if((M3 < M1 && M3 < M2) || (MassZ12 != 0 && MassW0 != 0 && MassZ01 == 0 && MassZ02 == 0) ) {
206 i = 1; j = 2; k = 0;
207 }
208
209 FourMomentum Zlepton1 = dressedleptons[i].momentum();
210 FourMomentum Zlepton2 = dressedleptons[j].momentum();
211 FourMomentum Wlepton = dressedleptons[k].momentum();
212 FourMomentum Zboson = dressedleptons[i].momentum()+dressedleptons[j].momentum();
213
214 double Wboson_mT = sqrt( 2 * Wlepton.pT() * neutrinos[0].pt() * (1 - cos(deltaPhi(Wlepton, neutrinos[0]))) );
215
216 if (fabs(Zboson.mass()/GeV - MZ_PDG) >= 10.) vetoEvent;
217 if (Wboson_mT <= 30*GeV) vetoEvent;
218 if (Wlepton.pT() <= 20*GeV) vetoEvent;
219 if (deltaR(Zlepton1, Zlepton2) < 0.2) vetoEvent;
220 if (deltaR(Zlepton1, Wlepton) < 0.3) vetoEvent;
221 if (deltaR(Zlepton2, Wlepton) < 0.3) vetoEvent;
222
223 if (EventType == 3) _h["eee"]->fill(13000);
224 if (EventType == 2) _h["mee"]->fill(13000);
225 if (EventType == 1) _h["emm"]->fill(13000);
226 if (EventType == 0) _h["mmm"]->fill(13000);
227 _h["fid"]->fill(13000);
228
229 if (EventCharge == 1) {
230 if (EventType == 3) _h["eee_Plus"]->fill(13000);
231 if (EventType == 2) _h["mee_Plus"]->fill(13000);
232 if (EventType == 1) _h["emm_Plus"]->fill(13000);
233 if (EventType == 0) _h["mmm_Plus"]->fill(13000);
234 _h["fid_Plus"]->fill(13000);
235 } else {
236 if (EventType == 3) _h["eee_Minus"]->fill(13000);
237 if (EventType == 2) _h["mee_Minus"]->fill(13000);
238 if (EventType == 1) _h["emm_Minus"]->fill(13000);
239 if (EventType == 0) _h["mmm_Minus"]->fill(13000);
240 _h["fid_Minus"]->fill(13000);
241 }
242
243 if (jets.size() == 0) _hnjets->fill("0.0"s);
244 else if (jets.size() == 1) _hnjets->fill("1.0"s);
245 else if (jets.size() == 2) _hnjets->fill("2.0"s);
246 else if (jets.size() == 3) _hnjets->fill("3.0"s);
247 else _hnjets->fill("$\\geq 4$"s);
248
249 }
250
251
252 void finalize() {
253
254 // Print summary info
255 const double sf_pb = crossSection() / picobarn / sumOfWeights();
256 const double sf_fb = crossSection() / femtobarn / sumOfWeights();
257 const float totalBR= 4*0.1086*0.033658; // W and Z leptonic branching fractions
258
259 scale(_h, sf_fb);
260 scale({ _h["fid"], _h["fid_Plus"], _h["fid_Minus"]}, 0.25);
261 scale(_hnjets, 0.25*sf_fb);
262 scale(_htot, sf_pb/totalBR);
263
264 }
265
266 /// @}
267
268
269 private:
270
271
272 /// @name Histograms
273 /// @{
274 map<string,BinnedHistoPtr<int>> _h;
275 BinnedHistoPtr<int> _htot;
276 BinnedHistoPtr<string> _hnjets;
277
278 /// @}
279
280 double MZ_PDG = 91.1876;
281 double MW_PDG = 83.385;
282 double GammaZ_PDG = 2.4952;
283 double GammaW_PDG = 2.085;
284
285 };
286
287 RIVET_DECLARE_PLUGIN(ATLAS_2016_I1469071);
288}
|
