Rivet analyses referenceCMS_2018_I1703993Measurements of ttbar differential cross sections in proton-proton collisions at 13 TeV using events containing two leptonsExperiment: CMS (LHC) Inspire ID: 1703993 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
Absolute and normalised differential top quark pair production cross sections are measured in the dilepton decay channels. The differential cross sections have been obtained by the CMS experiment at the CERN LHC in 2016 in proton-proton collisions at a centre-of-mass energy of 13 TeV. The measurements are performed with data corresponding to an integrated luminosity of 35.9/fb. The cross sections are measured differentially as a function of the kinematic properties of the top quarks and the tt system at the particle and parton levels, and the top quark decay products at the particle level. The results are compared to Monte Carlo simulations from POWHEG interfaced with the parton shower generators PYTHIA 8 and HERWIG 7 up to NNLO (next-to-next leading order) accuracy. Source code: CMS_2018_I1703993.cc 1#include "Rivet/Analysis.hh"
2#include "Rivet/Projections/FinalState.hh"
3#include "Rivet/Projections/FastJets.hh"
4#include "Rivet/Projections/PartonicTops.hh"
5#include "Rivet/Projections/LeptonFinder.hh"
6#include "Rivet/Projections/IdentifiedFinalState.hh"
7#include "Rivet/Projections/InvisibleFinalState.hh"
8#include "Rivet/Projections/PromptFinalState.hh"
9#include "Rivet/Projections/VetoedFinalState.hh"
10#include "Rivet/Projections/ChargedLeptons.hh"
11#include "Rivet/Math/LorentzTrans.hh"
12
13namespace Rivet {
14
15
16 /// ttbar dilepton differential cross-sections in pp collisions at 13 TeV
17 class CMS_2018_I1703993 : public Analysis { //
18 public:
19
20 RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2018_I1703993); //
21
22 void init() {
23 // Parton level top quarks dilepton e/mu channels only
24 declare(PartonicTops(TopDecay::E_MU, PromptEMuFromTau::NO), "PartonTops"); // Partonic top decaying to e or mu
25
26 // Build particle level tops starting from FinalState
27 const FinalState fs(Cuts::pT > 0. && Cuts::abseta < 6.);
28
29 // Neutrinos
30 InvisibleFinalState neutrinos(OnlyPrompt::YES, TauDecaysAs::PROMPT, MuDecaysAs::PROMPT);
31 declare(neutrinos, "Neutrinos");
32
33 // Projection for electrons and muons
34 Cut lepton_cut = Cuts::pt > 20 * GeV && Cuts::abseta < 2.4;
35
36 // Dressed leptons
37 ChargedLeptons charged_leptons(fs);
38 IdentifiedFinalState photons(fs);
39 photons.acceptIdPair(PID::PHOTON);
40
41 PromptFinalState prompt_leptons(charged_leptons, TauDecaysAs::NONPROMPT);
42 PromptFinalState prompt_photons(photons, TauDecaysAs::PROMPT);
43
44 LeptonFinder dressed_leptons(prompt_leptons, prompt_photons, 0.1, lepton_cut);
45 declare(dressed_leptons, "LeptonFinder");
46
47 // Projection for jets
48 VetoedFinalState fs_jets(fs);
49 fs_jets.addVetoOnThisFinalState(dressed_leptons);
50 fs_jets.vetoNeutrinos();
51 declare(FastJets(fs_jets, JetAlg::ANTIKT, 0.4), "ak4jets");
52
53 // Book hists for particle level, normalized values
54 book(_h["top_pt_norm"], "d07-x01-y01"); // table 15, normalized values
55 book(_h["atop_pt_norm"], "d15-x01-y01"); // table 16, normalized values
56 book(_h["top1_pt_norm"], "d23-x01-y01"); // table 17, normalized values
57 book(_h["top2_pt_norm"], "d31-x01-y01"); // table 18, normalized values
58 book(_h["top_pt_ttrf_norm"], "d39-x01-y01"); // table 19, normalized values
59 book(_h["top_y_norm"], "d47-x01-y01"); // table 20, normalized values
60 book(_h["atop_y_norm"], "d55-x01-y01"); // table 21, normalized values
61 book(_h["top1_y_norm"], "d63-x01-y01"); // table 22, normalized values
62 book(_h["top2_y_norm"], "d71-x01-y01"); // table 23, normalized values
63 book(_h["tt_pt_norm"], "d79-x01-y01"); // table 24, normalized values
64 book(_h["tt_y_norm"], "d87-x01-y01"); // table 25, normalized values
65 book(_h["tt_m_norm"], "d95-x01-y01"); // table 26, normalized values
66 book(_h["tt_dy_norm"], "d103-x01-y01"); // table 27, normalized values
67 book(_h["tt_dphi_norm"], "d111-x01-y01"); // table 28, normalized values
68 book(_h["lep_pt_norm"], "d115-x01-y01"); // table 29, normalized values,lepton
69 book(_h["alep_pt_norm"], "d119-x01-y01"); // table 30, normalized values,antilepton
70 book(_h["lep1_pt_norm"], "d123-x01-y01"); // table 31, normalized values,leading
71 book(_h["lep2_pt_norm"], "d127-x01-y01"); // table 32, normalized values,trailing
72 book(_h["lep_eta_norm"], "d131-x01-y01"); // table 33, normalized values, eta lepton
73 book(_h["alep_eta_norm"], "d135-x01-y01"); // table 34, normalized values, eta antilepton
74 book(_h["lep1_eta_norm"], "d139-x01-y01"); // table 35, normalized values, eta l leading
75 book(_h["lep2_eta_norm"], "d143-x01-y01"); // table 36, normalized values, eta l trailing
76 book(_h["ll_pt_norm"], "d147-x01-y01"); // table 37, normalized values
77 book(_h["ll_m_norm"], "d151-x01-y01"); // table 38, normalized values
78 book(_h["ll_dphi_norm"], "d155-x01-y01"); // table 39, normalized values
79 book(_h["ll_deta_norm"], "d159-x01-y01"); // table 40, normalized values
80 book(_h["jet_sz_norm"], "d187-x01-y01"); // table 41(47), normalized values
81 book(_h["jet1_pt_norm"], "d163-x01-y01"); // table 42, normalized values, leading
82 book(_h["jet2_pt_norm"], "d167-x01-y01"); // table 43, normalized values
83 book(_h["b1_eta_norm"], "d171-x01-y01"); // table 44(43), normalized values, leading
84 book(_h["b2_eta_norm"], "d175-x01-y01"); // table 45(44), normalized values, trailing
85 book(_h["bb_pt_norm"], "d179-x01-y01"); // table 46(45), normalized values
86 book(_h["bb_m_norm"], "d183-x01-y01"); // table 47(46), normalized values
87
88 // book hists for particle level, absolute values
89 book(_h_abs["top_pt"], "d05-x01-y01"); // table 15
90 book(_h_abs["atop_pt"], "d13-x01-y01"); // table 16
91 book(_h_abs["top1_pt"], "d21-x01-y01"); // table 17
92 book(_h_abs["top2_pt"], "d29-x01-y01"); // table 18
93 book(_h_abs["top_pt_ttrf"], "d37-x01-y01"); // table 19
94 book(_h_abs["top_y"], "d45-x01-y01"); // table 20
95 book(_h_abs["atop_y"], "d53-x01-y01"); // table 21
96 book(_h_abs["top1_y"], "d61-x01-y01"); // table 22
97 book(_h_abs["top2_y"], "d69-x01-y01"); // table 23
98 book(_h_abs["tt_pt"], "d77-x01-y01"); // table 24
99 book(_h_abs["tt_y"], "d85-x01-y01"); // table 25
100 book(_h_abs["tt_m"], "d93-x01-y01"); // table 26
101 book(_h_abs["tt_dy"], "d101-x01-y01"); // table 27
102 book(_h_abs["tt_dphi"], "d109-x01-y01"); // table 28
103 book(_h_abs["lep_pt"], "d113-x01-y01"); // table 29, lepton
104 book(_h_abs["alep_pt"], "d117-x01-y01"); // table 30, antilepton
105 book(_h_abs["lep1_pt"], "d121-x01-y01"); // table 31, leading
106 book(_h_abs["lep2_pt"], "d125-x01-y01"); // table 32, trailing
107 book(_h_abs["lep_eta"], "d129-x01-y01"); // table 33, eta lepton
108 book(_h_abs["alep_eta"], "d133-x01-y01"); // table 34, eta antilepton
109 book(_h_abs["lep1_eta"], "d137-x01-y01"); // table 35, eta l leading
110 book(_h_abs["lep2_eta"], "d141-x01-y01"); // table 36, eta l trailing
111 book(_h_abs["ll_pt"], "d145-x01-y01"); // table 37
112 book(_h_abs["ll_m"], "d149-x01-y01"); // table 38
113 book(_h_abs["ll_dphi"], "d153-x01-y01"); // table 39
114 book(_h_abs["ll_deta"], "d157-x01-y01"); // table 40
115 book(_h_abs["jet_sz"], "d185-x01-y01"); // table 41(47)
116 book(_h_abs["jet1_pt"], "d161-x01-y01"); // table 42(41), leading
117 book(_h_abs["jet2_pt"], "d165-x01-y01"); // table 43(42)
118 book(_h_abs["b1_eta"], "d169-x01-y01"); // table 44(43), leading
119 book(_h_abs["b2_eta"], "d173-x01-y01"); // table 45(44), trailing
120 book(_h_abs["bb_pt"], "d177-x01-y01"); // table 46(45)
121 book(_h_abs["bb_m"], "d181-x01-y01"); // table 47(46)
122
123 //book hists for parton level, normalized values
124 book(_h_part["top_pt_norm"], "d03-x01-y01"); // table 1
125 book(_h_part["atop_pt_norm"], "d11-x01-y01"); // table 2
126 book(_h_part["top1_pt_norm"], "d19-x01-y01"); // table 3
127 book(_h_part["top2_pt_norm"], "d27-x01-y01"); // table 4
128 book(_h_part["top_pt_ttrf_norm"], "d35-x01-y01"); // table 5
129 book(_h_part["top_y_norm"], "d43-x01-y01"); // table 6
130 book(_h_part["atop_y_norm"], "d51-x01-y01"); // table 7
131 book(_h_part["top1_y_norm"], "d59-x01-y01"); // table 8
132 book(_h_part["top2_y_norm"], "d67-x01-y01"); // table 9
133 book(_h_part["tt_pt_norm"], "d75-x01-y01"); // table 10
134 book(_h_part["tt_y_norm"], "d83-x01-y01"); // table 11
135 book(_h_part["tt_m_norm"], "d91-x01-y01"); // table 12
136 book(_h_part["tt_dy_norm"], "d99-x01-y01"); // table 13
137 book(_h_part["tt_dphi_norm"], "d107-x01-y01"); // table 14
138
139 //book hists for parton level, absolute values
140 book(_h_part_abs["top_pt"], "d01-x01-y01"); // table 1
141 book(_h_part_abs["atop_pt"], "d09-x01-y01"); // table 2
142 book(_h_part_abs["top1_pt"], "d17-x01-y01"); // table 3
143 book(_h_part_abs["top2_pt"], "d25-x01-y01"); // table 4
144 book(_h_part_abs["top_pt_ttrf"], "d33-x01-y01"); // table 5
145 book(_h_part_abs["top_y"], "d41-x01-y01"); // table 6
146 book(_h_part_abs["atop_y"], "d49-x01-y01"); // table 7
147 book(_h_part_abs["top1_y"], "d57-x01-y01"); // table 8
148 book(_h_part_abs["top2_y"], "d65-x01-y01"); // table 9
149 book(_h_part_abs["tt_pt"], "d73-x01-y01"); // table 10
150 book(_h_part_abs["tt_y"], "d81-x01-y01"); // table 11
151 book(_h_part_abs["tt_m"], "d89-x01-y01"); // table 12
152 book(_h_part_abs["tt_dy"], "d97-x01-y01"); // table 13
153 book(_h_part_abs["tt_dphi"], "d105-x01-y01"); // table 14
154 }
155
156 void analyze(const Event& event) {
157 // Parton-level analysis
158 const Particles partonTops = apply<ParticleFinder>(event, "PartonTops").particlesByPt();
159 if (partonTops.size() == 2) {
160 Particle top1 = partonTops[0];
161 Particle top2 = partonTops[1];
162
163 Particle top = top1;
164 Particle atop = top2;
165 if (top.pid() < 0) {
166 std::swap(top, atop);
167 }
168
169 const FourMomentum tt = top1.momentum() + top2.momentum();
170 const LorentzTransform rf = LorentzTransform::mkFrameTransform(tt);
171 const FourMomentum top_rf = rf.transform(top);
172
173 //Parton level normalized values
174 _h_part["top_pt_norm"]->fill(top.pt()/GeV);
175 _h_part["atop_pt_norm"]->fill(atop.pt()/GeV);
176 _h_part["top1_pt_norm"]->fill(top1.pt()/GeV);
177 _h_part["top2_pt_norm"]->fill(top2.pt()/GeV);
178 _h_part["top_pt_ttrf_norm"]->fill(top_rf.pt()/GeV);
179 _h_part["top_y_norm"]->fill(top.rapidity());
180 _h_part["atop_y_norm"]->fill(atop.rapidity());
181 _h_part["top1_y_norm"]->fill(top1.rapidity());
182 _h_part["top2_y_norm"]->fill(top2.rapidity());
183 _h_part["tt_pt_norm"]->fill(tt.pt()/GeV);
184 _h_part["tt_y_norm"]->fill(tt.rapidity());
185 _h_part["tt_m_norm"]->fill(tt.mass()/GeV);
186 _h_part["tt_dy_norm"]->fill(deltaRap(top.absrap(), atop.absrap(), true));
187 _h_part["tt_dphi_norm"]->fill(deltaPhi(top.phi(), atop.phi()));
188
189 //Parton level absolute values
190 _h_part_abs["top_pt"]->fill(top.pt()/GeV);
191 _h_part_abs["atop_pt"]->fill(atop.pt()/GeV);
192 _h_part_abs["top1_pt"]->fill(top1.pt()/GeV);
193 _h_part_abs["top2_pt"]->fill(top2.pt()/GeV);
194 _h_part_abs["top_pt_ttrf"]->fill(top_rf.pt()/GeV);
195 _h_part_abs["top_y"]->fill(top.rapidity());
196 _h_part_abs["atop_y"]->fill(atop.rapidity());
197 _h_part_abs["top1_y"]->fill(top1.rapidity());
198 _h_part_abs["top2_y"]->fill(top2.rapidity());
199 _h_part_abs["tt_pt"]->fill(tt.pt()/GeV);
200 _h_part_abs["tt_y"]->fill(tt.rapidity());
201 _h_part_abs["tt_m"]->fill(tt.mass()/GeV);
202 _h_part_abs["tt_dy"]->fill(deltaRap(top.absrap(), atop.absrap(), true));
203 _h_part_abs["tt_dphi"]->fill(deltaPhi(top.phi(), atop.phi()));
204 }
205
206 //Particle-level analysis
207 // Select leptons
208 const DressedLeptons& dressedLeptons = apply<LeptonFinder>(event, "LeptonFinder").dressedLeptons();
209 if (dressedLeptons.size() != 2) vetoEvent;
210 sortByPt(dressedLeptons);
211
212 const FourMomentum& lepton1 = dressedLeptons[0].momentum();
213 const FourMomentum& lepton2 = dressedLeptons[1].momentum();
214 if ((lepton1 + lepton2).mass() < 20*GeV)
215 vetoEvent;
216
217 // Select neutrinos
218 const Particles neutrinos = apply<InvisibleFinalState>(event, "Neutrinos").particlesByPt();
219 if (neutrinos.size() < 2)
220 vetoEvent;
221
222 // Select bjets
223 const FastJets& fjJets = apply<FastJets>(event, "ak4jets");
224 const Jets jets = fjJets.jetsByPt(Cuts::abseta < 2.4 && Cuts::pT > 30 * GeV);
225 const Jets bJets = select(jets, hasBTag());
226
227 // There should at least two b jets.
228 if (bJets.size() < 2)
229 vetoEvent;
230
231 // Construct particle level top
232 FourMomentum nu1 = neutrinos[0].momentum();
233 FourMomentum nu2 = neutrinos[1].momentum();
234 if (std::abs((lepton1 + nu1).mass() - 80.4*GeV) + std::abs((lepton2 + nu2).mass() - 80.4*GeV) >
235 std::abs((lepton1 + nu2).mass() - 80.4*GeV) + std::abs((lepton2 + nu1).mass() - 80.4*GeV)) {
236 std::swap(nu1, nu2);
237 }
238 const FourMomentum w1 = lepton1 + nu1;
239 const FourMomentum w2 = lepton2 + nu2;
240
241 double dm = 1e9; // Reset once again for top combination.
242 int selB1 = -1, selB2 = -1;
243 for (unsigned int i = 0; i < bJets.size(); ++i) {
244 const auto& bjet1 = bJets.at(i);
245 if (deltaR(bjet1, lepton1) < 0.4)
246 continue;
247 if (deltaR(bjet1, lepton2) < 0.4)
248 continue;
249
250 const double dm1 = std::abs((w1 + bjet1).mass() - 172.5*GeV);
251 for (unsigned int j = 0; j < bJets.size(); ++j) {
252 if (i == j)
253 continue;
254 const auto& bjet2 = bJets.at(j);
255 if (deltaR(bjet2, lepton1) < 0.4)
256 continue;
257 if (deltaR(bjet2, lepton2) < 0.4)
258 continue;
259 const double dm2 = std::abs((w2 + bjet2).mass() - 172.5*GeV);
260 const double newDm = dm1 + dm2;
261
262 if (newDm < dm) {
263 dm = newDm;
264 selB1 = i;
265 selB2 = j;
266 }
267 }
268 }
269
270 if (dm >= 1e9)
271 vetoEvent;
272
273 FourMomentum bjet1 = bJets[selB1].momentum();
274 FourMomentum bjet2 = bJets[selB2].momentum();
275
276 const FourMomentum t1 = w1 + bjet1;
277 const FourMomentum t2 = w2 + bjet2;
278 const FourMomentum tt = t1 + t2;
279
280 int q1 = dressedLeptons[0].charge();
281 int q2 = dressedLeptons[1].charge();
282 if (q1 * q2 > 0)
283 vetoEvent;
284
285 int idx_lep = (q1 == -1) ? 0 : 1;
286 int idx_alep = (q2 == 1) ? 1 : 0;
287
288 FourMomentum top = t1;
289 FourMomentum atop = t2;
290 if (q1 == -1) {
291 std::swap(top, atop);
292 }
293
294 FourMomentum top1 = t1;
295 FourMomentum top2 = t2;
296 if (top1.pt() < top2.pt()) {
297 std::swap(top1, top2);
298 }
299
300 if (bjet1.pt() < bjet2.pt()) {
301 std::swap(bjet1, bjet2);
302 }
303
304 const LorentzTransform rf = LorentzTransform::mkFrameTransform(tt);
305 const FourMomentum top_rf = rf.transform(top);
306
307 const FourMomentum ll = lepton1 + lepton2;
308
309 const FourMomentum bb = bjet1 + bjet2;
310
311 //particle level normalized values
312 _h["top_pt_norm"]->fill(top.pt()/GeV);
313 _h["atop_pt_norm"]->fill(atop.pt()/GeV);
314 _h["top1_pt_norm"]->fill(top1.pt()/GeV);
315 _h["top2_pt_norm"]->fill(top2.pt()/GeV);
316 _h["top_pt_ttrf_norm"]->fill(top_rf.pt()/GeV);
317 _h["top_y_norm"]->fill(top.rapidity());
318 _h["atop_y_norm"]->fill(atop.rapidity());
319 _h["top1_y_norm"]->fill(top1.rapidity());
320 _h["top2_y_norm"]->fill(top2.rapidity());
321 _h["tt_pt_norm"]->fill(tt.pt()/GeV);
322 _h["tt_dy_norm"]->fill(deltaRap(top.absrap(), atop.absrap(), true));
323 _h["tt_y_norm"]->fill(tt.rapidity());
324 _h["tt_m_norm"]->fill(tt.mass()/GeV);
325 _h["tt_dphi_norm"]->fill(deltaPhi(top.phi(), atop.phi()));
326 _h["lep_pt_norm"]->fill(dressedLeptons[idx_lep].pt()/GeV);
327 _h["alep_pt_norm"]->fill(dressedLeptons[idx_alep].pt()/GeV);
328 _h["lep1_pt_norm"]->fill(lepton1.pt()/GeV);
329 _h["lep2_pt_norm"]->fill(lepton2.pt()/GeV);
330 _h["lep_eta_norm"]->fill(dressedLeptons[idx_lep].eta());
331 _h["alep_eta_norm"]->fill(dressedLeptons[idx_alep].eta());
332 _h["lep1_eta_norm"]->fill(lepton1.eta());
333 _h["lep2_eta_norm"]->fill(lepton2.eta());
334 _h["ll_pt_norm"]->fill(ll.pt()/GeV);
335 _h["ll_m_norm"]->fill(ll.mass()/GeV);
336 _h["ll_dphi_norm"]->fill(deltaPhi(dressedLeptons[idx_lep], dressedLeptons[idx_alep]));
337 _h["ll_deta_norm"]->fill(deltaEta(dressedLeptons[idx_lep].abseta(), dressedLeptons[idx_alep].abseta(), true));
338 _h["jet_sz_norm"]->fill(jets.size());
339 _h["jet1_pt_norm"]->fill(bjet1.pt()/GeV);
340 _h["jet2_pt_norm"]->fill(bjet2.pt()/GeV);
341 _h["b1_eta_norm"]->fill(bjet1.eta());
342 _h["b2_eta_norm"]->fill(bjet2.eta());
343 _h["bb_pt_norm"]->fill(bb.pt()/GeV);
344 _h["bb_m_norm"]->fill(bb.mass()/GeV);
345
346 // absolute values
347 _h_abs["top_pt"]->fill(top.pt()/GeV);
348 _h_abs["atop_pt"]->fill(atop.pt()/GeV);
349 _h_abs["top1_pt"]->fill(top1.pt()/GeV);
350 _h_abs["top2_pt"]->fill(top2.pt()/GeV);
351 _h_abs["top_pt_ttrf"]->fill(top_rf.pt()/GeV);
352 _h_abs["top_y"]->fill(top.rapidity());
353 _h_abs["atop_y"]->fill(atop.rapidity());
354 _h_abs["top1_y"]->fill(top1.rapidity());
355 _h_abs["top2_y"]->fill(top2.rapidity());
356 _h_abs["tt_pt"]->fill(tt.pt()/GeV);
357 _h_abs["tt_y"]->fill(tt.rapidity());
358 _h_abs["tt_m"]->fill(tt.mass()/GeV);
359 _h_abs["tt_dy"]->fill(deltaRap(top.absrap(), atop.absrap(), true));
360 _h_abs["tt_dphi"]->fill(deltaPhi(top.phi(), atop.phi()));
361 _h_abs["lep_pt"]->fill(dressedLeptons[idx_lep].pt()/GeV);
362 _h_abs["alep_pt"]->fill(dressedLeptons[idx_alep].pt()/GeV);
363 _h_abs["lep1_pt"]->fill(lepton1.pt()/GeV);
364 _h_abs["lep2_pt"]->fill(lepton2.pt()/GeV);
365 _h_abs["lep_eta"]->fill(dressedLeptons[idx_lep].eta());
366 _h_abs["alep_eta"]->fill(dressedLeptons[idx_alep].eta());
367 _h_abs["lep1_eta"]->fill(lepton1.eta());
368 _h_abs["lep2_eta"]->fill(lepton2.eta());
369 _h_abs["ll_pt"]->fill(ll.pt()/GeV);
370 _h_abs["ll_m"]->fill(ll.mass()/GeV);
371 _h_abs["ll_dphi"]->fill(deltaPhi(dressedLeptons[idx_lep], dressedLeptons[idx_alep]));
372 _h_abs["ll_deta"]->fill(deltaEta(dressedLeptons[idx_lep].abseta(), dressedLeptons[idx_alep].abseta(), true));
373 _h_abs["jet_sz"]->fill(jets.size());
374 _h_abs["jet1_pt"]->fill(bjet1.pt()/GeV);
375 _h_abs["jet2_pt"]->fill(bjet2.pt()/GeV);
376 _h_abs["b1_eta"]->fill(bjet1.eta());
377 _h_abs["b2_eta"]->fill(bjet2.eta());
378 _h_abs["bb_pt"]->fill(bb.pt()/GeV);
379 _h_abs["bb_m"]->fill(bb.mass()/GeV);
380 }
381
382 /// Normalise histograms etc., after the run
383 void finalize() {
384 normalize(_h);
385 scale(_h_abs, crossSection() / picobarn / sumOfWeights());
386
387 normalize(_h_part);
388 scale(_h_part_abs,
389 crossSection() / picobarn / sumOfWeights() /
390 0.04553956); // BR correction for E_MU: 0.04553956 //BR correction for MUON: 0.01129969
391 }
392
393 private:
394 map<string, Histo1DPtr> _h;
395 map<string, Histo1DPtr> _h_abs;
396 map<string, Histo1DPtr> _h_part;
397 map<string, Histo1DPtr> _h_part_abs;
398 };
399
400 RIVET_DECLARE_PLUGIN(CMS_2018_I1703993);
401
402} // namespace Rivet
|