Rivet analyses referenceATLAS_2014_I1327229Inclusive multilepton search at 8 TeVExperiment: ATLAS (LHC) Inspire ID: 1327229 Status: VALIDATED Authors:
Beam energies: (4000.0, 4000.0) GeV Run details:
A generic search for anomalous production of events with at least three charged leptons is presented. The data sample consists of $pp$ collisions at $\sqrt{s} = 8$\,TeV collected in 2012 by the ATLAS experiment at the CERN Large Hadron Collider, and corresponds to an integrated luminosity of 20.3\,$\text{fb}^{-1}$. Events are required to have at least three selected lepton candidates, at least two of which must be electrons or muons, while the third may be a hadronically decaying tau. Selected events are categorized based on their lepton flavour content and signal regions are constructed using several kinematic variables of interest. No significant deviations from Standard Model predictions are observed. Model-independent upper limits on contributions from beyond the Standard Model phenomena are provided for each signal region, along with prescription to re-interpret the limits for any model. Constraints are also placed on models predicting doubly charged Higgs bosons and excited leptons. For doubly charged Higgs bosons decaying to $e\tau$ or $\muon\tau$, lower limits on the mass are set at 400\,GeV at 95\,\% confidence level. For excited leptons, constraints are provided as functions of both the mass of the excited state and the compositeness scale $\Lambda$, with the strongest mass constraints arising in regions where the mass equals $\Lambda$. In such scenarios, lower mass limits are set at 3.0\,TeV for excited electrons and muons, 2.5\,TeV for excited taus, and 1.6\,TeV for every excited-neutrino flavour. Source code: ATLAS_2014_I1327229.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/ChargedFinalState.hh"
5#include "Rivet/Projections/VisibleFinalState.hh"
6#include "Rivet/Projections/VetoedFinalState.hh"
7#include "Rivet/Projections/IdentifiedFinalState.hh"
8#include "Rivet/Projections/UnstableParticles.hh"
9#include "Rivet/Projections/FastJets.hh"
10#include "Rivet/Tools/Random.hh"
11
12namespace Rivet {
13
14
15 /// Inclusive multilepton search at 8 TeV
16 class ATLAS_2014_I1327229 : public Analysis {
17 public:
18
19 /// Constructor
20 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2014_I1327229);
21
22
23 /// Book histograms and initialise projections before the run
24 void init() {
25
26 // To calculate the acceptance without having the fiducial lepton efficiencies included, this part can be turned off
27 _use_fiducial_lepton_efficiency = true;
28
29 // Read in all signal regions
30 _signal_regions = getSignalRegions();
31
32 // Set number of events per signal region to 0
33 for (size_t i = 0; i < _signal_regions.size(); i++)
34 book(_eventCountsPerSR[_signal_regions[i]], "_eventCountsPerSR_" + _signal_regions[i]);
35
36 // Final state including all charged and neutral particles
37 const FinalState fs((Cuts::etaIn(-5.0, 5.0) && Cuts::pT >= 1*GeV));
38 declare(fs, "FS");
39
40 // Final state including all charged particles
41 declare(ChargedFinalState((Cuts::etaIn(-2.5, 2.5) && Cuts::pT >= 1*GeV)), "CFS");
42
43 // Final state including all visible particles (to calculate MET, Jets etc.)
44 declare(VisibleFinalState((Cuts::etaIn(-5.0,5.0))),"VFS");
45
46 // Final state including all AntiKt 04 Jets
47 VetoedFinalState vfs;
48 vfs.addVetoPairId(PID::MUON);
49 declare(FastJets(vfs, JetAlg::ANTIKT, 0.4), "AntiKtJets04");
50
51 // Final state including all unstable particles (including taus)
52 declare(UnstableParticles(Cuts::abseta < 5.0 && Cuts::pT > 5*GeV),"UFS");
53
54 // Final state including all electrons
55 IdentifiedFinalState elecs(Cuts::abseta < 2.47 && Cuts::pT > 10*GeV);
56 elecs.acceptIdPair(PID::ELECTRON);
57 declare(elecs, "elecs");
58
59 // Final state including all muons
60 IdentifiedFinalState muons(Cuts::abseta < 2.5 && Cuts::pT > 10*GeV);
61 muons.acceptIdPair(PID::MUON);
62 declare(muons, "muons");
63
64
65
66 /// Book histograms:
67 book(_h_HTlep_all ,"HTlep_all", 30,0,3000);
68 book(_h_HTjets_all ,"HTjets_all", 30,0,3000);
69 book(_h_MET_all ,"MET_all", 30,0,1500);
70 book(_h_Meff_all ,"Meff_all", 50,0,5000);
71 book(_h_min_pT_all ,"min_pT_all", 50, 0, 2000);
72 book(_h_mT_all ,"mT_all", 50, 0, 2000);
73
74 book(_h_e_n ,"e_n", 10, -0.5, 9.5);
75 book(_h_mu_n ,"mu_n", 10, -0.5, 9.5);
76 book(_h_tau_n ,"tau_n", 10, -0.5, 9.5);
77
78 book(_h_pt_1_3l ,"pt_1_3l", 100, 0, 2000);
79 book(_h_pt_2_3l ,"pt_2_3l", 100, 0, 2000);
80 book(_h_pt_3_3l ,"pt_3_3l", 100, 0, 2000);
81 book(_h_pt_1_2ltau ,"pt_1_2ltau", 100, 0, 2000);
82 book(_h_pt_2_2ltau ,"pt_2_2ltau", 100, 0, 2000);
83 book(_h_pt_3_2ltau ,"pt_3_2ltau", 100, 0, 2000);
84
85 book(_h_excluded ,"excluded", 2, -0.5, 1.5);
86 }
87
88
89 /// Perform the per-event analysis
90 void analyze(const Event& event) {
91
92 // Muons
93 Particles muon_candidates;
94 const Particles charged_tracks = apply<ChargedFinalState>(event, "CFS").particles();
95 const Particles visible_particles = apply<VisibleFinalState>(event, "VFS").particles();
96 for (const Particle& mu : apply<IdentifiedFinalState>(event, "muons").particlesByPt() ) {
97
98 // Calculate pTCone30 variable (pT of all tracks within dR<0.3 - pT of muon itself)
99 double pTinCone = -mu.pT();
100 for (const Particle& track : charged_tracks ) {
101 if (deltaR(mu.momentum(),track.momentum()) < 0.3 )
102 pTinCone += track.pT();
103 }
104
105 // Calculate eTCone30 variable (pT of all visible particles within dR<0.3)
106 double eTinCone = 0.;
107 for (const Particle& visible_particle : visible_particles) {
108 if (visible_particle.abspid() != PID::MUON && inRange(deltaR(mu.momentum(),visible_particle.momentum()), 0.1, 0.3))
109 eTinCone += visible_particle.pT();
110 }
111
112 // Apply reconstruction efficiency and simulate reconstruction
113 int muon_id = 13;
114 if (mu.hasAncestorWith(Cuts::pid == PID::TAU) || mu.hasAncestorWith(Cuts::pid == -PID::TAU)) muon_id = 14;
115 const double eff = (_use_fiducial_lepton_efficiency) ? apply_reco_eff(muon_id,mu) : 1.0;
116 const bool keep_muon = rand01()<=eff;
117
118 // Keep muon if pTCone30/pT < 0.15 and eTCone30/pT < 0.2 and reconstructed
119 if (keep_muon && pTinCone/mu.pT() <= 0.1 && eTinCone/mu.pT() < 0.1)
120 muon_candidates.push_back(mu);
121 }
122
123 // Electrons
124 Particles electron_candidates;
125 for (const Particle& e : apply<IdentifiedFinalState>(event, "elecs").particlesByPt() ) {
126 // Neglect electrons in crack regions
127 if (inRange(e.abseta(), 1.37, 1.52)) continue;
128
129 // Calculate pTCone30 variable (pT of all tracks within dR<0.3 - pT of electron itself)
130 double pTinCone = -e.pT();
131 for (const Particle& track : charged_tracks) {
132 if (deltaR(e.momentum(), track.momentum()) < 0.3 ) pTinCone += track.pT();
133 }
134
135 // Calculate eTCone30 variable (pT of all visible particles (except muons) within dR<0.3)
136 double eTinCone = 0.;
137 for (const Particle& visible_particle : visible_particles) {
138 if (visible_particle.abspid() != PID::MUON && inRange(deltaR(e.momentum(),visible_particle.momentum()), 0.1, 0.3))
139 eTinCone += visible_particle.pT();
140 }
141
142 // Apply reconstruction efficiency and simulate reconstruction
143 int elec_id = 11;
144 if (e.hasAncestorWith(Cuts::pid == 15) || e.hasAncestorWith(Cuts::pid == -15)) elec_id = 12;
145 const double eff = (_use_fiducial_lepton_efficiency) ? apply_reco_eff(elec_id,e) : 1.0;
146 const bool keep_elec = rand01()<=eff;
147
148 // Keep electron if pTCone30/pT < 0.13 and eTCone30/pT < 0.2 and reconstructed
149 if (keep_elec && pTinCone/e.pT() <= 0.1 && eTinCone/e.pT() < 0.1)
150 electron_candidates.push_back(e);
151 }
152
153
154 // Taus
155 Particles tau_candidates;
156 for (const Particle& tau : apply<UnstableParticles>(event, "UFS").particles() ) {
157 // Only pick taus out of all unstable particles
158 if ( tau.abspid() != PID::TAU) continue;
159 // Check that tau has decayed into daughter particles
160 if (tau.genParticle()->end_vertex() == 0) continue;
161 // Calculate visible tau momentum using the tau neutrino momentum in the tau decay
162 FourMomentum daughter_tau_neutrino_momentum = get_tau_neutrino_momentum(tau);
163 Particle tau_vis = tau;
164 tau_vis.setMomentum(tau.momentum()-daughter_tau_neutrino_momentum);
165 // keep only taus in certain eta region and above 15 GeV of visible tau pT
166 if ( tau_vis.pT()/GeV <= 15.0 || tau_vis.abseta() > 2.5) continue;
167
168 // Get prong number (number of tracks) in tau decay and check if tau decays leptonically
169 unsigned int nprong = 0;
170 bool lep_decaying_tau = false;
171 get_prong_number(tau.genParticle(),nprong,lep_decaying_tau);
172
173 // Apply reconstruction efficiency and simulate reconstruction
174 int tau_id = 15;
175 if (nprong == 1) tau_id = 15;
176 else if (nprong == 3) tau_id = 16;
177
178
179 const double eff = (_use_fiducial_lepton_efficiency) ? apply_reco_eff(tau_id,tau_vis) : 1.0;
180 const bool keep_tau = rand01()<=eff;
181
182 // Keep tau if nprong = 1, it decays hadronically and it is reconstructed
183 if ( !lep_decaying_tau && nprong == 1 && keep_tau) tau_candidates.push_back(tau_vis);
184 }
185
186 // Jets (all anti-kt R=0.4 jets with pT > 30 GeV and eta < 4.9
187 Jets jet_candidates = apply<FastJets>(event, "AntiKtJets04").jetsByPt(Cuts::pT > 30.0*GeV && Cuts::abseta < 4.9);
188
189 // ETmiss
190 Particles vfs_particles = apply<VisibleFinalState>(event, "VFS").particles();
191 FourMomentum pTmiss;
192 for (const Particle& p : vfs_particles)
193 pTmiss -= p.momentum();
194 double eTmiss = pTmiss.pT()/GeV;
195
196
197 // -------------------------
198 // Overlap removal
199
200 // electron - electron
201 Particles electron_candidates_2;
202 for(size_t ie = 0; ie < electron_candidates.size(); ++ie) {
203 const Particle& e = electron_candidates[ie];
204 bool away = true;
205 // If electron pair within dR < 0.1: remove electron with lower pT
206 for(size_t ie2 = 0; ie2 < electron_candidates_2.size(); ++ie2) {
207 if (deltaR(e.momentum(),electron_candidates_2[ie2].momentum()) < 0.1 ) {
208 away = false;
209 break;
210 }
211 }
212 // If isolated keep it
213 if ( away )
214 electron_candidates_2.push_back( e );
215 }
216
217 // jet - electron
218 Jets recon_jets;
219 for (const Jet& jet : jet_candidates) {
220 bool away = true;
221 // If jet within dR < 0.2 of electron: remove jet
222 for (const Particle& e : electron_candidates_2) {
223 if (deltaR(e.momentum(), jet.momentum()) < 0.2 ) {
224 away = false;
225 break;
226 }
227 }
228
229 // jet - tau
230 if ( away ) {
231 // If jet within dR < 0.2 of tau: remove jet
232 for (const Particle& tau : tau_candidates) {
233 if (deltaR(tau.momentum(), jet.momentum()) < 0.2 ) {
234 away = false;
235 break;
236 }
237 }
238 }
239 // If isolated keep it
240 if ( away )
241 recon_jets.push_back( jet );
242 }
243
244 // electron - jet
245 Particles recon_leptons, recon_e;
246 for (size_t ie = 0; ie < electron_candidates_2.size(); ++ie) {
247 const Particle& e = electron_candidates_2[ie];
248 // If electron within 0.2 < dR < 0.4 from any jets: remove electron
249 bool away = true;
250 for (const Jet& jet : recon_jets) {
251 if (deltaR(e.momentum(), jet.momentum()) < 0.4 ) {
252 away = false;
253 break;
254 }
255 }
256 // electron - muon
257 // If electron within dR < 0.1 of a muon: remove electron
258 if (away) {
259 for (const Particle& mu : muon_candidates) {
260 if (deltaR(mu.momentum(),e.momentum()) < 0.1) {
261 away = false;
262 break;
263 }
264 }
265 }
266 // If isolated keep it
267 if ( away ) {
268 recon_e.push_back( e );
269 recon_leptons.push_back( e );
270 }
271 }
272
273 // tau - electron
274 Particles recon_tau;
275 for (const Particle& tau : tau_candidates) {
276 bool away = true;
277 // If tau within dR < 0.2 of an electron: remove tau
278 for (const Particle & e : recon_e) {
279 if (deltaR(tau.momentum(),e.momentum()) < 0.2 ) {
280 away = false;
281 break;
282 }
283 }
284 // tau - muon
285 // If tau within dR < 0.2 of a muon: remove tau
286 if (away) {
287 for (const Particle& mu : muon_candidates) {
288 if (deltaR(tau.momentum(), mu.momentum()) < 0.2 ) {
289 away = false;
290 break;
291 }
292 }
293 }
294 // If isolated keep it
295 if (away) recon_tau.push_back( tau );
296 }
297
298 // muon - jet
299 Particles recon_mu, trigger_mu;
300 // If muon within dR < 0.4 of a jet: remove muon
301 for (const Particle& mu : muon_candidates ) {
302 bool away = true;
303 for (const Jet& jet : recon_jets) {
304 if (deltaR(mu.momentum(), jet.momentum()) < 0.4 ) {
305 away = false;
306 break;
307 }
308 }
309 if (away) {
310 recon_mu.push_back( mu );
311 recon_leptons.push_back( mu );
312 if (mu.abseta() < 2.4) trigger_mu.push_back( mu );
313 }
314 }
315
316 // End overlap removal
317 // ---------------------
318
319 // Jet cleaning
320 if (rand01() <= 0.42) {
321 for (const Jet& jet : recon_jets ) {
322 const double eta = jet.rapidity();
323 const double phi = jet.azimuthalAngle(MINUSPI_PLUSPI);
324 if(jet.pT() > 25*GeV && inRange(eta,-0.1,1.5) && inRange(phi,-0.9,-0.5)) vetoEvent;
325 }
326 }
327
328 // Event selection
329 // Require at least 3 charged tracks in event
330 if (charged_tracks.size() < 3) vetoEvent;
331
332 // And at least one e/mu passing trigger
333 if( !( !recon_e.empty() && recon_e[0].pT()>26.*GeV) &&
334 !( !trigger_mu.empty() && trigger_mu[0].pT()>26.*GeV) ) {
335 MSG_DEBUG("Hardest lepton fails trigger");
336 vetoEvent;
337 }
338
339 // And only accept events with at least 2 electrons and muons and at least 3 leptons in total
340 if (recon_mu.size() + recon_e.size() + recon_tau.size() < 3 || recon_leptons.size() < 2) vetoEvent;
341
342 // Sort leptons by decreasing pT
343 isortByPt(recon_leptons);
344 isortByPt(recon_tau);
345
346 // Calculate HTlep, fill lepton pT histograms & store chosen combination of 3 leptons
347 double HTlep = 0.;
348 Particles chosen_leptons;
349 if (recon_leptons.size() > 2) {
350 _h_pt_1_3l->fill(recon_leptons[0].pT()/GeV);
351 _h_pt_2_3l->fill(recon_leptons[1].pT()/GeV);
352 _h_pt_3_3l->fill(recon_leptons[2].pT()/GeV);
353 HTlep = (recon_leptons[0].pT() + recon_leptons[1].pT() + recon_leptons[2].pT())/GeV;
354 chosen_leptons.push_back( recon_leptons[0] );
355 chosen_leptons.push_back( recon_leptons[1] );
356 chosen_leptons.push_back( recon_leptons[2] );
357 }
358 else {
359 _h_pt_1_2ltau->fill(recon_leptons[0].pT()/GeV);
360 _h_pt_2_2ltau->fill(recon_leptons[1].pT()/GeV);
361 _h_pt_3_2ltau->fill(recon_tau[0].pT()/GeV);
362 HTlep = recon_leptons[0].pT()/GeV + recon_leptons[1].pT()/GeV + recon_tau[0].pT()/GeV;
363 chosen_leptons.push_back( recon_leptons[0] );
364 chosen_leptons.push_back( recon_leptons[1] );
365 chosen_leptons.push_back( recon_tau[0] );
366 }
367
368 // Calculate mT and mTW variable
369 Particles mT_leptons;
370 Particles mTW_leptons;
371 for (size_t i1 = 0; i1 < 3; i1 ++) {
372 for (size_t i2 = i1+1; i2 < 3; i2 ++) {
373 double OSSF_inv_mass = isOSSF_mass(chosen_leptons[i1],chosen_leptons[i2]);
374 if (OSSF_inv_mass != 0.) {
375 for (size_t i3 = 0; i3 < 3 ; i3 ++) {
376 if (i3 != i2 && i3 != i1) {
377 mT_leptons.push_back(chosen_leptons[i3]);
378 if ( fabs(91.0 - OSSF_inv_mass) < 20. )
379 mTW_leptons.push_back(chosen_leptons[i3]);
380 }
381 }
382 }
383 else {
384 mT_leptons.push_back(chosen_leptons[0]);
385 mTW_leptons.push_back(chosen_leptons[0]);
386 }
387 }
388 }
389
390 isortByPt(mT_leptons);
391 isortByPt(mTW_leptons);
392
393 double mT = sqrt(2*pTmiss.pT()/GeV*mT_leptons[0].pT()/GeV*(1-cos(pTmiss.phi()-mT_leptons[0].phi())));
394 double mTW = sqrt(2*pTmiss.pT()/GeV*mTW_leptons[0].pT()/GeV*(1-cos(pTmiss.phi()-mTW_leptons[0].phi())));
395
396 // Calculate Min pT variable
397 double min_pT = chosen_leptons[2].pT()/GeV;
398
399 // Number of prompt e/mu and had taus
400 _h_e_n->fill(recon_e.size());
401 _h_mu_n->fill(recon_mu.size());
402 _h_tau_n->fill(recon_tau.size());
403
404 // Calculate HTjets variable
405 double HTjets = 0.;
406 for (const Jet& jet : recon_jets)
407 HTjets += jet.pT()/GeV;
408
409 // Calculate meff variable
410 double meff = eTmiss + HTjets;
411 Particles all_leptons;
412 for (const Particle& e : recon_e ) {
413 meff += e.pT()/GeV;
414 all_leptons.push_back( e );
415 }
416 for (const Particle& mu : recon_mu) {
417 meff += mu.pT()/GeV;
418 all_leptons.push_back( mu );
419 }
420 for (const Particle& tau : recon_tau) {
421 meff += tau.pT()/GeV;
422 all_leptons.push_back( tau );
423 }
424
425 // Fill histograms of kinematic variables
426 _h_HTlep_all->fill(HTlep);
427 _h_HTjets_all->fill(HTjets);
428 _h_MET_all->fill(eTmiss);
429 _h_Meff_all->fill(meff);
430 _h_min_pT_all->fill(min_pT);
431 _h_mT_all->fill(mT);
432
433 // Determine signal region (3l / 2ltau , onZ / offZ OSSF / offZ no-OSSF)
434 // 3l vs. 2ltau
435 string basic_signal_region;
436 if (recon_mu.size() + recon_e.size() > 2)
437 basic_signal_region += "3l_";
438 else if ( (recon_mu.size() + recon_e.size() == 2) && (recon_tau.size() > 0))
439 basic_signal_region += "2ltau_";
440
441 // Is there an OSSF pair or a three lepton combination with an invariant mass close to the Z mass
442 int onZ = isonZ(chosen_leptons);
443 if (onZ == 1) basic_signal_region += "onZ";
444 else if (onZ == 0) {
445 bool OSSF = isOSSF(chosen_leptons);
446 if (OSSF) basic_signal_region += "offZ_OSSF";
447 else basic_signal_region += "offZ_noOSSF";
448 }
449
450 // Check in which signal regions this event falls and adjust event counters
451 // INFO: The b-jet signal regions of the paper are not included in this Rivet implementation
452 fillEventCountsPerSR(basic_signal_region,onZ,HTlep,eTmiss,HTjets,meff,min_pT,mTW);
453 }
454
455
456 /// Normalise histograms etc., after the run
457 void finalize() {
458
459 // Normalize to an integrated luminosity of 1 fb-1
460 double norm = crossSection()/femtobarn/sumOfWeights();
461
462 string best_signal_region = "";
463 double ratio_best_SR = 0.;
464
465 // Loop over all signal regions and find signal region with best sensitivity (ratio signal events/visible cross-section)
466 for (size_t i = 0; i < _signal_regions.size(); i++) {
467 double signal_events = _eventCountsPerSR[_signal_regions[i]]->val() * norm;
468 // Use expected upper limits to find best signal region:
469 double UL95 = getUpperLimit(_signal_regions[i],false);
470 double ratio = signal_events / UL95;
471 if (ratio > ratio_best_SR) {
472 best_signal_region = _signal_regions.at(i);
473 ratio_best_SR = ratio;
474 }
475 }
476
477 double signal_events_best_SR = _eventCountsPerSR[best_signal_region]->val() * norm;
478 double exp_UL_best_SR = getUpperLimit(best_signal_region, false);
479 double obs_UL_best_SR = getUpperLimit(best_signal_region, true);
480
481
482 // Print out result
483 cout << "----------------------------------------------------------------------------------------" << '\n';
484 cout << "Number of total events: " << sumOfWeights() << '\n';
485 cout << "Best signal region: " << best_signal_region << '\n';
486 cout << "Normalized number of signal events in this best signal region (per fb-1): " << signal_events_best_SR << '\n';
487 cout << "Efficiency*Acceptance: " << _eventCountsPerSR[best_signal_region]->val()/sumOfWeights() << '\n';
488 cout << "Cross-section [fb]: " << crossSection()/femtobarn << '\n';
489 cout << "Expected visible cross-section (per fb-1): " << exp_UL_best_SR << '\n';
490 cout << "Ratio (signal events / expected visible cross-section): " << ratio_best_SR << '\n';
491 cout << "Observed visible cross-section (per fb-1): " << obs_UL_best_SR << '\n';
492 cout << "Ratio (signal events / observed visible cross-section): " << signal_events_best_SR/obs_UL_best_SR << '\n';
493 cout << "----------------------------------------------------------------------------------------" << '\n';
494
495 cout << "Using the EXPECTED limits (visible cross-section) of the analysis: " << '\n';
496 if (signal_events_best_SR > exp_UL_best_SR) {
497 cout << "Since the number of signal events > the visible cross-section, this model/grid point is EXCLUDED with 95% C.L." << '\n';
498 _h_excluded->fill(1);
499 }
500 else {
501 cout << "Since the number of signal events < the visible cross-section, this model/grid point is NOT EXCLUDED." << '\n';
502 _h_excluded->fill(0);
503 }
504 cout << "----------------------------------------------------------------------------------------" << '\n';
505
506 cout << "Using the OBSERVED limits (visible cross-section) of the analysis: " << '\n';
507 if (signal_events_best_SR > obs_UL_best_SR) {
508 cout << "Since the number of signal events > the visible cross-section, this model/grid point is EXCLUDED with 95% C.L." << '\n';
509 _h_excluded->fill(1);
510 }
511 else {
512 cout << "Since the number of signal events < the visible cross-section, this model/grid point is NOT EXCLUDED." << '\n';
513 _h_excluded->fill(0);
514 }
515 cout << "----------------------------------------------------------------------------------------" << '\n';
516 cout << "INFO: The b-jet signal regions of the paper are not included in this Rivet implementation." << '\n';
517 cout << "----------------------------------------------------------------------------------------" << '\n';
518
519
520 /// Normalize to cross section
521
522 if (norm != 0) {
523 scale(_h_HTlep_all, norm);
524 scale(_h_HTjets_all, norm);
525 scale(_h_MET_all, norm);
526 scale(_h_Meff_all, norm);
527 scale(_h_min_pT_all, norm);
528 scale(_h_mT_all, norm);
529
530 scale(_h_pt_1_3l, norm);
531 scale(_h_pt_2_3l, norm);
532 scale(_h_pt_3_3l, norm);
533 scale(_h_pt_1_2ltau, norm);
534 scale(_h_pt_2_2ltau, norm);
535 scale(_h_pt_3_2ltau, norm);
536
537 scale(_h_e_n, norm);
538 scale(_h_mu_n, norm);
539 scale(_h_tau_n, norm);
540
541 scale(_h_excluded, norm);
542 }
543
544 }
545
546
547 /// Helper functions
548 /// @{
549 /// Function giving a list of all signal regions
550 vector<string> getSignalRegions() {
551
552 // List of basic signal regions
553 vector<string> basic_signal_regions;
554 basic_signal_regions.push_back("3l_offZ_OSSF");
555 basic_signal_regions.push_back("3l_offZ_noOSSF");
556 basic_signal_regions.push_back("3l_onZ");
557 basic_signal_regions.push_back("2ltau_offZ_OSSF");
558 basic_signal_regions.push_back("2ltau_offZ_noOSSF");
559 basic_signal_regions.push_back("2ltau_onZ");
560
561 // List of kinematic variables
562 vector<string> kinematic_variables;
563 kinematic_variables.push_back("HTlep");
564 kinematic_variables.push_back("METStrong");
565 kinematic_variables.push_back("METWeak");
566 kinematic_variables.push_back("Meff");
567 kinematic_variables.push_back("MeffStrong");
568 kinematic_variables.push_back("MeffMt");
569 kinematic_variables.push_back("MinPt");
570
571 vector<string> signal_regions;
572 // Loop over all kinematic variables and basic signal regions
573 for (size_t i0 = 0; i0 < kinematic_variables.size(); i0++) {
574 for (size_t i1 = 0; i1 < basic_signal_regions.size(); i1++) {
575 // Is signal region onZ?
576 int onZ = (basic_signal_regions[i1].find("onZ") != string::npos) ? 1 : 0;
577 // Get cut values for this kinematic variable
578 vector<int> cut_values = getCutsPerSignalRegion(kinematic_variables[i0], onZ);
579 // Loop over all cut values
580 for (size_t i2 = 0; i2 < cut_values.size(); i2++) {
581 // Push signal region into vector
582 signal_regions.push_back( kinematic_variables[i0] + "_" + basic_signal_regions[i1] + "_cut_" + toString(cut_values[i2]) );
583 }
584 }
585 }
586 return signal_regions;
587 }
588
589
590
591 /// Function giving all cut values per kinematic variable
592 vector<int> getCutsPerSignalRegion(const string& signal_region, int onZ = 0) {
593 vector<int> cutValues;
594
595 // Cut values for HTlep
596 if (signal_region.compare("HTlep") == 0) {
597 cutValues.push_back(0);
598 cutValues.push_back(200);
599 cutValues.push_back(500);
600 cutValues.push_back(800);
601 }
602 // Cut values for MinPt
603 else if (signal_region.compare("MinPt") == 0) {
604 cutValues.push_back(0);
605 cutValues.push_back(50);
606 cutValues.push_back(100);
607 cutValues.push_back(150);
608 }
609 // Cut values for METStrong (HTjets > 150 GeV) and METWeak (HTjets < 150 GeV)
610 else if (signal_region.compare("METStrong") == 0 || signal_region.compare("METWeak") == 0) {
611 cutValues.push_back(0);
612 cutValues.push_back(100);
613 cutValues.push_back(200);
614 cutValues.push_back(300);
615 }
616 // Cut values for Meff
617 if (signal_region.compare("Meff") == 0) {
618 cutValues.push_back(0);
619 cutValues.push_back(600);
620 cutValues.push_back(1000);
621 cutValues.push_back(1500);
622 }
623 // Cut values for MeffStrong (MET > 100 GeV)
624 if ((signal_region.compare("MeffStrong") == 0 || signal_region.compare("MeffMt") == 0) && onZ ==1) {
625 cutValues.push_back(0);
626 cutValues.push_back(600);
627 cutValues.push_back(1200);
628 }
629
630 return cutValues;
631 }
632
633 /// function fills map _eventCountsPerSR by looping over all signal regions
634 /// and looking if the event falls into this signal region
635 void fillEventCountsPerSR(const string& basic_signal_region, int onZ,
636 double HTlep, double eTmiss, double HTjets,
637 double meff, double min_pT, double mTW) {
638
639 // Get cut values for HTlep, loop over them and add event if cut is passed
640 vector<int> cut_values = getCutsPerSignalRegion("HTlep", onZ);
641 for (size_t i = 0; i < cut_values.size(); i++) {
642 if (HTlep > cut_values[i])
643 _eventCountsPerSR[("HTlep_" + basic_signal_region + "_cut_" + toString(cut_values[i]))]->fill();
644 }
645
646 // Get cut values for MinPt, loop over them and add event if cut is passed
647 cut_values = getCutsPerSignalRegion("MinPt", onZ);
648 for (size_t i = 0; i < cut_values.size(); i++) {
649 if (min_pT > cut_values[i])
650 _eventCountsPerSR[("MinPt_" + basic_signal_region + "_cut_" + toString(cut_values[i]))]->fill();
651 }
652
653 // Get cut values for METStrong, loop over them and add event if cut is passed
654 cut_values = getCutsPerSignalRegion("METStrong", onZ);
655 for (size_t i = 0; i < cut_values.size(); i++) {
656 if (eTmiss > cut_values[i] && HTjets > 150.)
657 _eventCountsPerSR[("METStrong_" + basic_signal_region + "_cut_" + toString(cut_values[i]))]->fill();
658 }
659
660 // Get cut values for METWeak, loop over them and add event if cut is passed
661 cut_values = getCutsPerSignalRegion("METWeak", onZ);
662 for (size_t i = 0; i < cut_values.size(); i++) {
663 if (eTmiss > cut_values[i] && HTjets <= 150.)
664 _eventCountsPerSR[("METWeak_" + basic_signal_region + "_cut_" + toString(cut_values[i]))]->fill();
665 }
666
667 // Get cut values for Meff, loop over them and add event if cut is passed
668 cut_values = getCutsPerSignalRegion("Meff", onZ);
669 for (size_t i = 0; i < cut_values.size(); i++) {
670 if (meff > cut_values[i])
671 _eventCountsPerSR[("Meff_" + basic_signal_region + "_cut_" + toString(cut_values[i]))]->fill();
672 }
673
674 // Get cut values for MeffStrong, loop over them and add event if cut is passed
675 cut_values = getCutsPerSignalRegion("MeffStrong", onZ);
676 for (size_t i = 0; i < cut_values.size(); i++) {
677 if (meff > cut_values[i] && eTmiss > 100.)
678 _eventCountsPerSR[("MeffStrong_" + basic_signal_region + "_cut_" + toString(cut_values[i]))]->fill();
679 }
680
681 // Get cut values for MeffMt, loop over them and add event if cut is passed
682 cut_values = getCutsPerSignalRegion("MeffMt", onZ);
683 for (size_t i = 0; i < cut_values.size(); i++) {
684 if (meff > cut_values[i] && mTW > 100. && onZ == 1)
685 _eventCountsPerSR[("MeffMt_" + basic_signal_region + "_cut_" + toString(cut_values[i]))]->fill();
686 }
687
688 }
689
690 /// Function returning 4-momentum of daughter-particle if it is a tau neutrino
691 FourMomentum get_tau_neutrino_momentum(const Particle& p) {
692 assert(p.abspid() == PID::TAU);
693 ConstGenVertexPtr dv = p.genParticle()->end_vertex();
694 assert(dv != nullptr);
695 // Loop over all daughter particles
696 for(ConstGenParticlePtr pp: HepMCUtils::particles(dv, Relatives::CHILDREN)){
697 if (abs(pp->pdg_id()) == PID::NU_TAU) return FourMomentum(pp->momentum());
698 }
699 return FourMomentum();
700 }
701
702 /// Function calculating the prong number of taus
703 void get_prong_number(ConstGenParticlePtr p, unsigned int& nprong, bool& lep_decaying_tau) {
704 assert(p != nullptr);
705 ConstGenVertexPtr dv = p->end_vertex();
706 assert(dv != nullptr);
707 for(ConstGenParticlePtr pp: HepMCUtils::particles(dv, Relatives::CHILDREN)){
708 // If they have status 1 and are charged they will produce a track and the prong number is +1
709 if (pp->status() == 1 ) {
710 const int id = pp->pdg_id();
711 if (Rivet::PID::charge(id) != 0 ) ++nprong;
712 // Check if tau decays leptonically
713 if (( abs(id) == PID::ELECTRON || abs(id) == PID::MUON || abs(id) == PID::TAU ) && abs(p->pdg_id()) == PID::TAU) lep_decaying_tau = true;
714 }
715 // If the status of the daughter particle is 2 it is unstable and the further decays are checked
716 else if (pp->status() == 2 ) {
717 get_prong_number(pp,nprong,lep_decaying_tau);
718 }
719 }
720 }
721
722 /// Function giving fiducial lepton efficiency
723 double apply_reco_eff(int flavor, const Particle& p) {
724 double pt = p.pT()/GeV;
725 double eta = p.eta();
726
727 double eff = 0.;
728
729 if (flavor == 11) { // weight prompt electron -- now including data/MC ID SF in eff.
730 double avgrate = 0.685;
731 const static double wz_ele[] = {0.0256,0.522,0.607,0.654,0.708,0.737,0.761,0.784,0.815,0.835,0.851,0.841,0.898};
732 // double ewz_ele[] = {0.000257,0.00492,0.00524,0.00519,0.00396,0.00449,0.00538,0.00513,0.00773,0.00753,0.0209,0.0964,0.259};
733 int ibin = 0;
734 if(pt > 10 && pt < 15) ibin = 0;
735 if(pt > 15 && pt < 20) ibin = 1;
736 if(pt > 20 && pt < 25) ibin = 2;
737 if(pt > 25 && pt < 30) ibin = 3;
738 if(pt > 30 && pt < 40) ibin = 4;
739 if(pt > 40 && pt < 50) ibin = 5;
740 if(pt > 50 && pt < 60) ibin = 6;
741 if(pt > 60 && pt < 80) ibin = 7;
742 if(pt > 80 && pt < 100) ibin = 8;
743 if(pt > 100 && pt < 200) ibin = 9;
744 if(pt > 200 && pt < 400) ibin = 10;
745 if(pt > 400 && pt < 600) ibin = 11;
746 if(pt > 600) ibin = 12;
747 double eff_pt = 0.;
748 eff_pt = wz_ele[ibin];
749
750 eta = fabs(eta);
751
752 const static double wz_ele_eta[] = {0.65,0.714,0.722,0.689,0.635,0.615};
753 // double ewz_ele_eta[] = {0.00642,0.00355,0.00335,0.004,0.00368,0.00422};
754 ibin = 0;
755 if(eta > 0 && eta < 0.1) ibin = 0;
756 if(eta > 0.1 && eta < 0.5) ibin = 1;
757 if(eta > 0.5 && eta < 1.0) ibin = 2;
758 if(eta > 1.0 && eta < 1.5) ibin = 3;
759 if(eta > 1.5 && eta < 2.0) ibin = 4;
760 if(eta > 2.0 && eta < 2.5) ibin = 5;
761 double eff_eta = 0.;
762 eff_eta = wz_ele_eta[ibin];
763
764 eff = (eff_pt * eff_eta) / avgrate;
765 }
766
767 if (flavor == 12) { // weight electron from tau
768 double avgrate = 0.476;
769 const static double wz_ele[] = {0.00855,0.409,0.442,0.55,0.632,0.616,0.615,0.642,0.72,0.617};
770 // double ewz_ele[] = {0.000573,0.0291,0.0366,0.0352,0.0363,0.0474,0.0628,0.0709,0.125,0.109};
771 int ibin = 0;
772 if(pt > 10 && pt < 15) ibin = 0;
773 if(pt > 15 && pt < 20) ibin = 1;
774 if(pt > 20 && pt < 25) ibin = 2;
775 if(pt > 25 && pt < 30) ibin = 3;
776 if(pt > 30 && pt < 40) ibin = 4;
777 if(pt > 40 && pt < 50) ibin = 5;
778 if(pt > 50 && pt < 60) ibin = 6;
779 if(pt > 60 && pt < 80) ibin = 7;
780 if(pt > 80 && pt < 100) ibin = 8;
781 if(pt > 100) ibin = 9;
782 double eff_pt = 0.;
783 eff_pt = wz_ele[ibin];
784
785 eta = fabs(eta);
786
787 const static double wz_ele_eta[] = {0.546,0.5,0.513,0.421,0.47,0.433};
788 //double ewz_ele_eta[] = {0.0566,0.0257,0.0263,0.0263,0.0303,0.0321};
789 ibin = 0;
790 if(eta > 0 && eta < 0.1) ibin = 0;
791 if(eta > 0.1 && eta < 0.5) ibin = 1;
792 if(eta > 0.5 && eta < 1.0) ibin = 2;
793 if(eta > 1.0 && eta < 1.5) ibin = 3;
794 if(eta > 1.5 && eta < 2.0) ibin = 4;
795 if(eta > 2.0 && eta < 2.5) ibin = 5;
796 double eff_eta = 0.;
797 eff_eta = wz_ele_eta[ibin];
798
799 eff = (eff_pt * eff_eta) / avgrate;
800 }
801
802 if (flavor == 13) { // weight prompt muon
803 int ibin = 0;
804 if(pt > 10 && pt < 15) ibin = 0;
805 if(pt > 15 && pt < 20) ibin = 1;
806 if(pt > 20 && pt < 25) ibin = 2;
807 if(pt > 25 && pt < 30) ibin = 3;
808 if(pt > 30 && pt < 40) ibin = 4;
809 if(pt > 40 && pt < 50) ibin = 5;
810 if(pt > 50 && pt < 60) ibin = 6;
811 if(pt > 60 && pt < 80) ibin = 7;
812 if(pt > 80 && pt < 100) ibin = 8;
813 if(pt > 100 && pt < 200) ibin = 9;
814 if(pt > 200 && pt < 400) ibin = 10;
815 if(pt > 400) ibin = 11;
816 if(fabs(eta) < 0.1) {
817 const static double wz_mu[] = {0.00705,0.402,0.478,0.49,0.492,0.499,0.527,0.512,0.53,0.528,0.465,0.465};
818 //double ewz_mu[] = {0.000298,0.0154,0.017,0.0158,0.0114,0.0123,0.0155,0.0133,0.0196,0.0182,0.0414,0.0414};
819 double eff_pt = 0.;
820 eff_pt = wz_mu[ibin];
821 eff = eff_pt;
822 }
823 if(fabs(eta) > 0.1) {
824 const static double wz_mu[] = {0.0224,0.839,0.887,0.91,0.919,0.923,0.925,0.925,0.922,0.918,0.884,0.834};
825 //double ewz_mu[] = {0.000213,0.00753,0.0074,0.007,0.00496,0.00534,0.00632,0.00583,0.00849,0.00804,0.0224,0.0963};
826 double eff_pt = 0.;
827 eff_pt = wz_mu[ibin];
828 eff = eff_pt;
829 }
830 }
831
832 if (flavor == 14) { // weight muon from tau
833 int ibin = 0;
834 if(pt > 10 && pt < 15) ibin = 0;
835 if(pt > 15 && pt < 20) ibin = 1;
836 if(pt > 20 && pt < 25) ibin = 2;
837 if(pt > 25 && pt < 30) ibin = 3;
838 if(pt > 30 && pt < 40) ibin = 4;
839 if(pt > 40 && pt < 50) ibin = 5;
840 if(pt > 50 && pt < 60) ibin = 6;
841 if(pt > 60 && pt < 80) ibin = 7;
842 if(pt > 80 && pt < 100) ibin = 8;
843 if(pt > 100) ibin = 9;
844
845 if(fabs(eta) < 0.1) {
846 const static double wz_mu[] = {0.0,0.664,0.124,0.133,0.527,0.283,0.495,0.25,0.5,0.331};
847 //double ewz_mu[] = {0.0,0.192,0.0437,0.0343,0.128,0.107,0.202,0.125,0.25,0.191};
848 double eff_pt = 0.;
849 eff_pt = wz_mu[ibin];
850 eff = eff_pt;
851 }
852 if(fabs(eta) > 0.1) {
853 const static double wz_mu[] = {0.0,0.617,0.655,0.676,0.705,0.738,0.712,0.783,0.646,0.745};
854 //double ewz_mu[] = {0.0,0.043,0.0564,0.0448,0.0405,0.0576,0.065,0.0825,0.102,0.132};
855 double eff_pt = 0.;
856 eff_pt = wz_mu[ibin];
857 eff = eff_pt;
858 }
859 }
860
861 if (flavor == 15) { // weight hadronic tau 1p
862 double avgrate = 0.16;
863 const static double wz_tau1p[] = {0.0,0.0311,0.148,0.229,0.217,0.292,0.245,0.307,0.227,0.277};
864 //double ewz_tau1p[] = {0.0,0.00211,0.0117,0.0179,0.0134,0.0248,0.0264,0.0322,0.0331,0.0427};
865 int ibin = 0;
866 if(pt > 10 && pt < 15) ibin = 0;
867 if(pt > 15 && pt < 20) ibin = 1;
868 if(pt > 20 && pt < 25) ibin = 2;
869 if(pt > 25 && pt < 30) ibin = 3;
870 if(pt > 30 && pt < 40) ibin = 4;
871 if(pt > 40 && pt < 50) ibin = 5;
872 if(pt > 50 && pt < 60) ibin = 6;
873 if(pt > 60 && pt < 80) ibin = 7;
874 if(pt > 80 && pt < 100) ibin = 8;
875 if(pt > 100) ibin = 9;
876 double eff_pt = 0.;
877 eff_pt = wz_tau1p[ibin];
878
879 const static double wz_tau1p_eta[] = {0.166,0.15,0.188,0.175,0.142,0.109};
880 //double ewz_tau1p_eta[] ={0.0166,0.00853,0.0097,0.00985,0.00949,0.00842};
881 ibin = 0;
882 if(eta > 0.0 && eta < 0.1) ibin = 0;
883 if(eta > 0.1 && eta < 0.5) ibin = 1;
884 if(eta > 0.5 && eta < 1.0) ibin = 2;
885 if(eta > 1.0 && eta < 1.5) ibin = 3;
886 if(eta > 1.5 && eta < 2.0) ibin = 4;
887 if(eta > 2.0 && eta < 2.5) ibin = 5;
888 double eff_eta = 0.;
889 eff_eta = wz_tau1p_eta[ibin];
890
891 eff = (eff_pt * eff_eta) / avgrate;
892 }
893
894 return eff;
895 }
896
897
898 /// Function giving observed and expected upper limits (on the visible cross-section)
899 double getUpperLimit(const string& signal_region, bool observed) {
900
901 map<string,double> upperLimitsObserved;
902 map<string,double> upperLimitsExpected;
903
904 upperLimitsObserved["HTlep_3l_offZ_OSSF_cut_0"] = 2.435;
905 upperLimitsObserved["HTlep_3l_offZ_OSSF_cut_200"] = 0.704;
906 upperLimitsObserved["HTlep_3l_offZ_OSSF_cut_500"] = 0.182;
907 upperLimitsObserved["HTlep_3l_offZ_OSSF_cut_800"] = 0.147;
908 upperLimitsObserved["HTlep_2ltau_offZ_OSSF_cut_0"] = 13.901;
909 upperLimitsObserved["HTlep_2ltau_offZ_OSSF_cut_200"] = 1.677;
910 upperLimitsObserved["HTlep_2ltau_offZ_OSSF_cut_500"] = 0.141;
911 upperLimitsObserved["HTlep_2ltau_offZ_OSSF_cut_800"] = 0.155;
912 upperLimitsObserved["HTlep_3l_offZ_noOSSF_cut_0"] = 1.054;
913 upperLimitsObserved["HTlep_3l_offZ_noOSSF_cut_200"] = 0.341;
914 upperLimitsObserved["HTlep_3l_offZ_noOSSF_cut_500"] = 0.221;
915 upperLimitsObserved["HTlep_3l_offZ_noOSSF_cut_800"] = 0.140;
916 upperLimitsObserved["HTlep_2ltau_offZ_noOSSF_cut_0"] = 4.276;
917 upperLimitsObserved["HTlep_2ltau_offZ_noOSSF_cut_200"] = 0.413;
918 upperLimitsObserved["HTlep_2ltau_offZ_noOSSF_cut_500"] = 0.138;
919 upperLimitsObserved["HTlep_2ltau_offZ_noOSSF_cut_800"] = 0.150;
920 upperLimitsObserved["HTlep_3l_onZ_cut_0"] = 29.804;
921 upperLimitsObserved["HTlep_3l_onZ_cut_200"] = 3.579;
922 upperLimitsObserved["HTlep_3l_onZ_cut_500"] = 0.466;
923 upperLimitsObserved["HTlep_3l_onZ_cut_800"] = 0.298;
924 upperLimitsObserved["HTlep_2ltau_onZ_cut_0"] = 205.091;
925 upperLimitsObserved["HTlep_2ltau_onZ_cut_200"] = 3.141;
926 upperLimitsObserved["HTlep_2ltau_onZ_cut_500"] = 0.290;
927 upperLimitsObserved["HTlep_2ltau_onZ_cut_800"] = 0.157;
928 upperLimitsObserved["METStrong_3l_offZ_OSSF_cut_0"] = 1.111;
929 upperLimitsObserved["METStrong_3l_offZ_OSSF_cut_100"] = 0.354;
930 upperLimitsObserved["METStrong_3l_offZ_OSSF_cut_200"] = 0.236;
931 upperLimitsObserved["METStrong_3l_offZ_OSSF_cut_300"] = 0.150;
932 upperLimitsObserved["METStrong_2ltau_offZ_OSSF_cut_0"] = 1.881;
933 upperLimitsObserved["METStrong_2ltau_offZ_OSSF_cut_100"] = 0.406;
934 upperLimitsObserved["METStrong_2ltau_offZ_OSSF_cut_200"] = 0.194;
935 upperLimitsObserved["METStrong_2ltau_offZ_OSSF_cut_300"] = 0.134;
936 upperLimitsObserved["METStrong_3l_offZ_noOSSF_cut_0"] = 0.770;
937 upperLimitsObserved["METStrong_3l_offZ_noOSSF_cut_100"] = 0.295;
938 upperLimitsObserved["METStrong_3l_offZ_noOSSF_cut_200"] = 0.149;
939 upperLimitsObserved["METStrong_3l_offZ_noOSSF_cut_300"] = 0.140;
940 upperLimitsObserved["METStrong_2ltau_offZ_noOSSF_cut_0"] = 2.003;
941 upperLimitsObserved["METStrong_2ltau_offZ_noOSSF_cut_100"] = 0.806;
942 upperLimitsObserved["METStrong_2ltau_offZ_noOSSF_cut_200"] = 0.227;
943 upperLimitsObserved["METStrong_2ltau_offZ_noOSSF_cut_300"] = 0.138;
944 upperLimitsObserved["METStrong_3l_onZ_cut_0"] = 6.383;
945 upperLimitsObserved["METStrong_3l_onZ_cut_100"] = 0.959;
946 upperLimitsObserved["METStrong_3l_onZ_cut_200"] = 0.549;
947 upperLimitsObserved["METStrong_3l_onZ_cut_300"] = 0.182;
948 upperLimitsObserved["METStrong_2ltau_onZ_cut_0"] = 10.658;
949 upperLimitsObserved["METStrong_2ltau_onZ_cut_100"] = 0.637;
950 upperLimitsObserved["METStrong_2ltau_onZ_cut_200"] = 0.291;
951 upperLimitsObserved["METStrong_2ltau_onZ_cut_300"] = 0.227;
952 upperLimitsObserved["METWeak_3l_offZ_OSSF_cut_0"] = 1.802;
953 upperLimitsObserved["METWeak_3l_offZ_OSSF_cut_100"] = 0.344;
954 upperLimitsObserved["METWeak_3l_offZ_OSSF_cut_200"] = 0.189;
955 upperLimitsObserved["METWeak_3l_offZ_OSSF_cut_300"] = 0.148;
956 upperLimitsObserved["METWeak_2ltau_offZ_OSSF_cut_0"] = 12.321;
957 upperLimitsObserved["METWeak_2ltau_offZ_OSSF_cut_100"] = 0.430;
958 upperLimitsObserved["METWeak_2ltau_offZ_OSSF_cut_200"] = 0.137;
959 upperLimitsObserved["METWeak_2ltau_offZ_OSSF_cut_300"] = 0.134;
960 upperLimitsObserved["METWeak_3l_offZ_noOSSF_cut_0"] = 0.562;
961 upperLimitsObserved["METWeak_3l_offZ_noOSSF_cut_100"] = 0.153;
962 upperLimitsObserved["METWeak_3l_offZ_noOSSF_cut_200"] = 0.154;
963 upperLimitsObserved["METWeak_3l_offZ_noOSSF_cut_300"] = 0.141;
964 upperLimitsObserved["METWeak_2ltau_offZ_noOSSF_cut_0"] = 2.475;
965 upperLimitsObserved["METWeak_2ltau_offZ_noOSSF_cut_100"] = 0.244;
966 upperLimitsObserved["METWeak_2ltau_offZ_noOSSF_cut_200"] = 0.141;
967 upperLimitsObserved["METWeak_2ltau_offZ_noOSSF_cut_300"] = 0.142;
968 upperLimitsObserved["METWeak_3l_onZ_cut_0"] = 24.769;
969 upperLimitsObserved["METWeak_3l_onZ_cut_100"] = 0.690;
970 upperLimitsObserved["METWeak_3l_onZ_cut_200"] = 0.198;
971 upperLimitsObserved["METWeak_3l_onZ_cut_300"] = 0.138;
972 upperLimitsObserved["METWeak_2ltau_onZ_cut_0"] = 194.360;
973 upperLimitsObserved["METWeak_2ltau_onZ_cut_100"] = 0.287;
974 upperLimitsObserved["METWeak_2ltau_onZ_cut_200"] = 0.144;
975 upperLimitsObserved["METWeak_2ltau_onZ_cut_300"] = 0.130;
976 upperLimitsObserved["Meff_3l_offZ_OSSF_cut_0"] = 2.435;
977 upperLimitsObserved["Meff_3l_offZ_OSSF_cut_600"] = 0.487;
978 upperLimitsObserved["Meff_3l_offZ_OSSF_cut_1000"] = 0.156;
979 upperLimitsObserved["Meff_3l_offZ_OSSF_cut_1500"] = 0.140;
980 upperLimitsObserved["Meff_2ltau_offZ_OSSF_cut_0"] = 13.901;
981 upperLimitsObserved["Meff_2ltau_offZ_OSSF_cut_600"] = 0.687;
982 upperLimitsObserved["Meff_2ltau_offZ_OSSF_cut_1000"] = 0.224;
983 upperLimitsObserved["Meff_2ltau_offZ_OSSF_cut_1500"] = 0.155;
984 upperLimitsObserved["Meff_3l_offZ_noOSSF_cut_0"] = 1.054;
985 upperLimitsObserved["Meff_3l_offZ_noOSSF_cut_600"] = 0.249;
986 upperLimitsObserved["Meff_3l_offZ_noOSSF_cut_1000"] = 0.194;
987 upperLimitsObserved["Meff_3l_offZ_noOSSF_cut_1500"] = 0.145;
988 upperLimitsObserved["Meff_2ltau_offZ_noOSSF_cut_0"] = 4.276;
989 upperLimitsObserved["Meff_2ltau_offZ_noOSSF_cut_600"] = 0.772;
990 upperLimitsObserved["Meff_2ltau_offZ_noOSSF_cut_1000"] = 0.218;
991 upperLimitsObserved["Meff_2ltau_offZ_noOSSF_cut_1500"] = 0.204;
992 upperLimitsObserved["Meff_3l_onZ_cut_0"] = 29.804;
993 upperLimitsObserved["Meff_3l_onZ_cut_600"] = 2.933;
994 upperLimitsObserved["Meff_3l_onZ_cut_1000"] = 0.912;
995 upperLimitsObserved["Meff_3l_onZ_cut_1500"] = 0.225;
996 upperLimitsObserved["Meff_2ltau_onZ_cut_0"] = 205.091;
997 upperLimitsObserved["Meff_2ltau_onZ_cut_600"] = 1.486;
998 upperLimitsObserved["Meff_2ltau_onZ_cut_1000"] = 0.641;
999 upperLimitsObserved["Meff_2ltau_onZ_cut_1500"] = 0.204;
1000 upperLimitsObserved["MeffStrong_3l_offZ_OSSF_cut_0"] = 0.479;
1001 upperLimitsObserved["MeffStrong_3l_offZ_OSSF_cut_600"] = 0.353;
1002 upperLimitsObserved["MeffStrong_3l_offZ_OSSF_cut_1200"] = 0.187;
1003 upperLimitsObserved["MeffStrong_2ltau_offZ_OSSF_cut_0"] = 0.617;
1004 upperLimitsObserved["MeffStrong_2ltau_offZ_OSSF_cut_600"] = 0.320;
1005 upperLimitsObserved["MeffStrong_2ltau_offZ_OSSF_cut_1200"] = 0.281;
1006 upperLimitsObserved["MeffStrong_3l_offZ_noOSSF_cut_0"] = 0.408;
1007 upperLimitsObserved["MeffStrong_3l_offZ_noOSSF_cut_600"] = 0.240;
1008 upperLimitsObserved["MeffStrong_3l_offZ_noOSSF_cut_1200"] = 0.150;
1009 upperLimitsObserved["MeffStrong_2ltau_offZ_noOSSF_cut_0"] = 0.774;
1010 upperLimitsObserved["MeffStrong_2ltau_offZ_noOSSF_cut_600"] = 0.417;
1011 upperLimitsObserved["MeffStrong_2ltau_offZ_noOSSF_cut_1200"] = 0.266;
1012 upperLimitsObserved["MeffStrong_3l_onZ_cut_0"] = 1.208;
1013 upperLimitsObserved["MeffStrong_3l_onZ_cut_600"] = 0.837;
1014 upperLimitsObserved["MeffStrong_3l_onZ_cut_1200"] = 0.269;
1015 upperLimitsObserved["MeffStrong_2ltau_onZ_cut_0"] = 0.605;
1016 upperLimitsObserved["MeffStrong_2ltau_onZ_cut_600"] = 0.420;
1017 upperLimitsObserved["MeffStrong_2ltau_onZ_cut_1200"] = 0.141;
1018 upperLimitsObserved["MeffMt_3l_onZ_cut_0"] = 1.832;
1019 upperLimitsObserved["MeffMt_3l_onZ_cut_600"] = 0.862;
1020 upperLimitsObserved["MeffMt_3l_onZ_cut_1200"] = 0.222;
1021 upperLimitsObserved["MeffMt_2ltau_onZ_cut_0"] = 1.309;
1022 upperLimitsObserved["MeffMt_2ltau_onZ_cut_600"] = 0.481;
1023 upperLimitsObserved["MeffMt_2ltau_onZ_cut_1200"] = 0.146;
1024 upperLimitsObserved["MinPt_3l_offZ_OSSF_cut_0"] = 2.435;
1025 upperLimitsObserved["MinPt_3l_offZ_OSSF_cut_50"] = 0.500;
1026 upperLimitsObserved["MinPt_3l_offZ_OSSF_cut_100"] = 0.203;
1027 upperLimitsObserved["MinPt_3l_offZ_OSSF_cut_150"] = 0.128;
1028 upperLimitsObserved["MinPt_2ltau_offZ_OSSF_cut_0"] = 13.901;
1029 upperLimitsObserved["MinPt_2ltau_offZ_OSSF_cut_50"] = 0.859;
1030 upperLimitsObserved["MinPt_2ltau_offZ_OSSF_cut_100"] = 0.158;
1031 upperLimitsObserved["MinPt_2ltau_offZ_OSSF_cut_150"] = 0.155;
1032 upperLimitsObserved["MinPt_3l_offZ_noOSSF_cut_0"] = 1.054;
1033 upperLimitsObserved["MinPt_3l_offZ_noOSSF_cut_50"] = 0.295;
1034 upperLimitsObserved["MinPt_3l_offZ_noOSSF_cut_100"] = 0.148;
1035 upperLimitsObserved["MinPt_3l_offZ_noOSSF_cut_150"] = 0.137;
1036 upperLimitsObserved["MinPt_2ltau_offZ_noOSSF_cut_0"] = 4.276;
1037 upperLimitsObserved["MinPt_2ltau_offZ_noOSSF_cut_50"] = 0.314;
1038 upperLimitsObserved["MinPt_2ltau_offZ_noOSSF_cut_100"] = 0.134;
1039 upperLimitsObserved["MinPt_2ltau_offZ_noOSSF_cut_150"] = 0.140;
1040 upperLimitsObserved["MinPt_3l_onZ_cut_0"] = 29.804;
1041 upperLimitsObserved["MinPt_3l_onZ_cut_50"] = 1.767;
1042 upperLimitsObserved["MinPt_3l_onZ_cut_100"] = 0.690;
1043 upperLimitsObserved["MinPt_3l_onZ_cut_150"] = 0.301;
1044 upperLimitsObserved["MinPt_2ltau_onZ_cut_0"] = 205.091;
1045 upperLimitsObserved["MinPt_2ltau_onZ_cut_50"] = 1.050;
1046 upperLimitsObserved["MinPt_2ltau_onZ_cut_100"] = 0.155;
1047 upperLimitsObserved["MinPt_2ltau_onZ_cut_150"] = 0.146;
1048 upperLimitsObserved["nbtag_3l_offZ_OSSF_cut_0"] = 2.435;
1049 upperLimitsObserved["nbtag_3l_offZ_OSSF_cut_1"] = 0.865;
1050 upperLimitsObserved["nbtag_3l_offZ_OSSF_cut_2"] = 0.474;
1051 upperLimitsObserved["nbtag_2ltau_offZ_OSSF_cut_0"] = 13.901;
1052 upperLimitsObserved["nbtag_2ltau_offZ_OSSF_cut_1"] = 1.566;
1053 upperLimitsObserved["nbtag_2ltau_offZ_OSSF_cut_2"] = 0.426;
1054 upperLimitsObserved["nbtag_3l_offZ_noOSSF_cut_0"] = 1.054;
1055 upperLimitsObserved["nbtag_3l_offZ_noOSSF_cut_1"] = 0.643;
1056 upperLimitsObserved["nbtag_3l_offZ_noOSSF_cut_2"] = 0.321;
1057 upperLimitsObserved["nbtag_2ltau_offZ_noOSSF_cut_0"] = 4.276;
1058 upperLimitsObserved["nbtag_2ltau_offZ_noOSSF_cut_1"] = 2.435;
1059 upperLimitsObserved["nbtag_2ltau_offZ_noOSSF_cut_2"] = 1.073;
1060 upperLimitsObserved["nbtag_3l_onZ_cut_0"] = 29.804;
1061 upperLimitsObserved["nbtag_3l_onZ_cut_1"] = 3.908;
1062 upperLimitsObserved["nbtag_3l_onZ_cut_2"] = 0.704;
1063 upperLimitsObserved["nbtag_2ltau_onZ_cut_0"] = 205.091;
1064 upperLimitsObserved["nbtag_2ltau_onZ_cut_1"] = 9.377;
1065 upperLimitsObserved["nbtag_2ltau_onZ_cut_2"] = 0.657;
1066 upperLimitsExpected["HTlep_3l_offZ_OSSF_cut_0"] = 2.893;
1067 upperLimitsExpected["HTlep_3l_offZ_OSSF_cut_200"] = 1.175;
1068 upperLimitsExpected["HTlep_3l_offZ_OSSF_cut_500"] = 0.265;
1069 upperLimitsExpected["HTlep_3l_offZ_OSSF_cut_800"] = 0.155;
1070 upperLimitsExpected["HTlep_2ltau_offZ_OSSF_cut_0"] = 14.293;
1071 upperLimitsExpected["HTlep_2ltau_offZ_OSSF_cut_200"] = 1.803;
1072 upperLimitsExpected["HTlep_2ltau_offZ_OSSF_cut_500"] = 0.159;
1073 upperLimitsExpected["HTlep_2ltau_offZ_OSSF_cut_800"] = 0.155;
1074 upperLimitsExpected["HTlep_3l_offZ_noOSSF_cut_0"] = 0.836;
1075 upperLimitsExpected["HTlep_3l_offZ_noOSSF_cut_200"] = 0.340;
1076 upperLimitsExpected["HTlep_3l_offZ_noOSSF_cut_500"] = 0.218;
1077 upperLimitsExpected["HTlep_3l_offZ_noOSSF_cut_800"] = 0.140;
1078 upperLimitsExpected["HTlep_2ltau_offZ_noOSSF_cut_0"] = 4.132;
1079 upperLimitsExpected["HTlep_2ltau_offZ_noOSSF_cut_200"] = 0.599;
1080 upperLimitsExpected["HTlep_2ltau_offZ_noOSSF_cut_500"] = 0.146;
1081 upperLimitsExpected["HTlep_2ltau_offZ_noOSSF_cut_800"] = 0.148;
1082 upperLimitsExpected["HTlep_3l_onZ_cut_0"] = 32.181;
1083 upperLimitsExpected["HTlep_3l_onZ_cut_200"] = 4.879;
1084 upperLimitsExpected["HTlep_3l_onZ_cut_500"] = 0.473;
1085 upperLimitsExpected["HTlep_3l_onZ_cut_800"] = 0.266;
1086 upperLimitsExpected["HTlep_2ltau_onZ_cut_0"] = 217.801;
1087 upperLimitsExpected["HTlep_2ltau_onZ_cut_200"] = 3.676;
1088 upperLimitsExpected["HTlep_2ltau_onZ_cut_500"] = 0.235;
1089 upperLimitsExpected["HTlep_2ltau_onZ_cut_800"] = 0.150;
1090 upperLimitsExpected["METStrong_3l_offZ_OSSF_cut_0"] = 1.196;
1091 upperLimitsExpected["METStrong_3l_offZ_OSSF_cut_100"] = 0.423;
1092 upperLimitsExpected["METStrong_3l_offZ_OSSF_cut_200"] = 0.208;
1093 upperLimitsExpected["METStrong_3l_offZ_OSSF_cut_300"] = 0.158;
1094 upperLimitsExpected["METStrong_2ltau_offZ_OSSF_cut_0"] = 2.158;
1095 upperLimitsExpected["METStrong_2ltau_offZ_OSSF_cut_100"] = 0.461;
1096 upperLimitsExpected["METStrong_2ltau_offZ_OSSF_cut_200"] = 0.186;
1097 upperLimitsExpected["METStrong_2ltau_offZ_OSSF_cut_300"] = 0.138;
1098 upperLimitsExpected["METStrong_3l_offZ_noOSSF_cut_0"] = 0.495;
1099 upperLimitsExpected["METStrong_3l_offZ_noOSSF_cut_100"] = 0.284;
1100 upperLimitsExpected["METStrong_3l_offZ_noOSSF_cut_200"] = 0.150;
1101 upperLimitsExpected["METStrong_3l_offZ_noOSSF_cut_300"] = 0.146;
1102 upperLimitsExpected["METStrong_2ltau_offZ_noOSSF_cut_0"] = 1.967;
1103 upperLimitsExpected["METStrong_2ltau_offZ_noOSSF_cut_100"] = 0.732;
1104 upperLimitsExpected["METStrong_2ltau_offZ_noOSSF_cut_200"] = 0.225;
1105 upperLimitsExpected["METStrong_2ltau_offZ_noOSSF_cut_300"] = 0.147;
1106 upperLimitsExpected["METStrong_3l_onZ_cut_0"] = 7.157;
1107 upperLimitsExpected["METStrong_3l_onZ_cut_100"] = 1.342;
1108 upperLimitsExpected["METStrong_3l_onZ_cut_200"] = 0.508;
1109 upperLimitsExpected["METStrong_3l_onZ_cut_300"] = 0.228;
1110 upperLimitsExpected["METStrong_2ltau_onZ_cut_0"] = 12.441;
1111 upperLimitsExpected["METStrong_2ltau_onZ_cut_100"] = 0.534;
1112 upperLimitsExpected["METStrong_2ltau_onZ_cut_200"] = 0.243;
1113 upperLimitsExpected["METStrong_2ltau_onZ_cut_300"] = 0.218;
1114 upperLimitsExpected["METWeak_3l_offZ_OSSF_cut_0"] = 2.199;
1115 upperLimitsExpected["METWeak_3l_offZ_OSSF_cut_100"] = 0.391;
1116 upperLimitsExpected["METWeak_3l_offZ_OSSF_cut_200"] = 0.177;
1117 upperLimitsExpected["METWeak_3l_offZ_OSSF_cut_300"] = 0.144;
1118 upperLimitsExpected["METWeak_2ltau_offZ_OSSF_cut_0"] = 12.431;
1119 upperLimitsExpected["METWeak_2ltau_offZ_OSSF_cut_100"] = 0.358;
1120 upperLimitsExpected["METWeak_2ltau_offZ_OSSF_cut_200"] = 0.150;
1121 upperLimitsExpected["METWeak_2ltau_offZ_OSSF_cut_300"] = 0.135;
1122 upperLimitsExpected["METWeak_3l_offZ_noOSSF_cut_0"] = 0.577;
1123 upperLimitsExpected["METWeak_3l_offZ_noOSSF_cut_100"] = 0.214;
1124 upperLimitsExpected["METWeak_3l_offZ_noOSSF_cut_200"] = 0.155;
1125 upperLimitsExpected["METWeak_3l_offZ_noOSSF_cut_300"] = 0.140;
1126 upperLimitsExpected["METWeak_2ltau_offZ_noOSSF_cut_0"] = 2.474;
1127 upperLimitsExpected["METWeak_2ltau_offZ_noOSSF_cut_100"] = 0.382;
1128 upperLimitsExpected["METWeak_2ltau_offZ_noOSSF_cut_200"] = 0.144;
1129 upperLimitsExpected["METWeak_2ltau_offZ_noOSSF_cut_300"] = 0.146;
1130 upperLimitsExpected["METWeak_3l_onZ_cut_0"] = 26.305;
1131 upperLimitsExpected["METWeak_3l_onZ_cut_100"] = 1.227;
1132 upperLimitsExpected["METWeak_3l_onZ_cut_200"] = 0.311;
1133 upperLimitsExpected["METWeak_3l_onZ_cut_300"] = 0.188;
1134 upperLimitsExpected["METWeak_2ltau_onZ_cut_0"] = 205.198;
1135 upperLimitsExpected["METWeak_2ltau_onZ_cut_100"] = 0.399;
1136 upperLimitsExpected["METWeak_2ltau_onZ_cut_200"] = 0.166;
1137 upperLimitsExpected["METWeak_2ltau_onZ_cut_300"] = 0.140;
1138 upperLimitsExpected["Meff_3l_offZ_OSSF_cut_0"] = 2.893;
1139 upperLimitsExpected["Meff_3l_offZ_OSSF_cut_600"] = 0.649;
1140 upperLimitsExpected["Meff_3l_offZ_OSSF_cut_1000"] = 0.252;
1141 upperLimitsExpected["Meff_3l_offZ_OSSF_cut_1500"] = 0.150;
1142 upperLimitsExpected["Meff_2ltau_offZ_OSSF_cut_0"] = 14.293;
1143 upperLimitsExpected["Meff_2ltau_offZ_OSSF_cut_600"] = 0.657;
1144 upperLimitsExpected["Meff_2ltau_offZ_OSSF_cut_1000"] = 0.226;
1145 upperLimitsExpected["Meff_2ltau_offZ_OSSF_cut_1500"] = 0.154;
1146 upperLimitsExpected["Meff_3l_offZ_noOSSF_cut_0"] = 0.836;
1147 upperLimitsExpected["Meff_3l_offZ_noOSSF_cut_600"] = 0.265;
1148 upperLimitsExpected["Meff_3l_offZ_noOSSF_cut_1000"] = 0.176;
1149 upperLimitsExpected["Meff_3l_offZ_noOSSF_cut_1500"] = 0.146;
1150 upperLimitsExpected["Meff_2ltau_offZ_noOSSF_cut_0"] = 4.132;
1151 upperLimitsExpected["Meff_2ltau_offZ_noOSSF_cut_600"] = 0.678;
1152 upperLimitsExpected["Meff_2ltau_offZ_noOSSF_cut_1000"] = 0.243;
1153 upperLimitsExpected["Meff_2ltau_offZ_noOSSF_cut_1500"] = 0.184;
1154 upperLimitsExpected["Meff_3l_onZ_cut_0"] = 32.181;
1155 upperLimitsExpected["Meff_3l_onZ_cut_600"] = 3.219;
1156 upperLimitsExpected["Meff_3l_onZ_cut_1000"] = 0.905;
1157 upperLimitsExpected["Meff_3l_onZ_cut_1500"] = 0.261;
1158 upperLimitsExpected["Meff_2ltau_onZ_cut_0"] = 217.801;
1159 upperLimitsExpected["Meff_2ltau_onZ_cut_600"] = 1.680;
1160 upperLimitsExpected["Meff_2ltau_onZ_cut_1000"] = 0.375;
1161 upperLimitsExpected["Meff_2ltau_onZ_cut_1500"] = 0.178;
1162 upperLimitsExpected["MeffStrong_3l_offZ_OSSF_cut_0"] = 0.571;
1163 upperLimitsExpected["MeffStrong_3l_offZ_OSSF_cut_600"] = 0.386;
1164 upperLimitsExpected["MeffStrong_3l_offZ_OSSF_cut_1200"] = 0.177;
1165 upperLimitsExpected["MeffStrong_2ltau_offZ_OSSF_cut_0"] = 0.605;
1166 upperLimitsExpected["MeffStrong_2ltau_offZ_OSSF_cut_600"] = 0.335;
1167 upperLimitsExpected["MeffStrong_2ltau_offZ_OSSF_cut_1200"] = 0.249;
1168 upperLimitsExpected["MeffStrong_3l_offZ_noOSSF_cut_0"] = 0.373;
1169 upperLimitsExpected["MeffStrong_3l_offZ_noOSSF_cut_600"] = 0.223;
1170 upperLimitsExpected["MeffStrong_3l_offZ_noOSSF_cut_1200"] = 0.150;
1171 upperLimitsExpected["MeffStrong_2ltau_offZ_noOSSF_cut_0"] = 0.873;
1172 upperLimitsExpected["MeffStrong_2ltau_offZ_noOSSF_cut_600"] = 0.428;
1173 upperLimitsExpected["MeffStrong_2ltau_offZ_noOSSF_cut_1200"] = 0.210;
1174 upperLimitsExpected["MeffStrong_3l_onZ_cut_0"] = 2.034;
1175 upperLimitsExpected["MeffStrong_3l_onZ_cut_600"] = 1.093;
1176 upperLimitsExpected["MeffStrong_3l_onZ_cut_1200"] = 0.293;
1177 upperLimitsExpected["MeffStrong_2ltau_onZ_cut_0"] = 0.690;
1178 upperLimitsExpected["MeffStrong_2ltau_onZ_cut_600"] = 0.392;
1179 upperLimitsExpected["MeffStrong_2ltau_onZ_cut_1200"] = 0.156;
1180 upperLimitsExpected["MeffMt_3l_onZ_cut_0"] = 2.483;
1181 upperLimitsExpected["MeffMt_3l_onZ_cut_600"] = 0.845;
1182 upperLimitsExpected["MeffMt_3l_onZ_cut_1200"] = 0.255;
1183 upperLimitsExpected["MeffMt_2ltau_onZ_cut_0"] = 1.448;
1184 upperLimitsExpected["MeffMt_2ltau_onZ_cut_600"] = 0.391;
1185 upperLimitsExpected["MeffMt_2ltau_onZ_cut_1200"] = 0.146;
1186 upperLimitsExpected["MinPt_3l_offZ_OSSF_cut_0"] = 2.893;
1187 upperLimitsExpected["MinPt_3l_offZ_OSSF_cut_50"] = 0.703;
1188 upperLimitsExpected["MinPt_3l_offZ_OSSF_cut_100"] = 0.207;
1189 upperLimitsExpected["MinPt_3l_offZ_OSSF_cut_150"] = 0.143;
1190 upperLimitsExpected["MinPt_2ltau_offZ_OSSF_cut_0"] = 14.293;
1191 upperLimitsExpected["MinPt_2ltau_offZ_OSSF_cut_50"] = 0.705;
1192 upperLimitsExpected["MinPt_2ltau_offZ_OSSF_cut_100"] = 0.149;
1193 upperLimitsExpected["MinPt_2ltau_offZ_OSSF_cut_150"] = 0.155;
1194 upperLimitsExpected["MinPt_3l_offZ_noOSSF_cut_0"] = 0.836;
1195 upperLimitsExpected["MinPt_3l_offZ_noOSSF_cut_50"] = 0.249;
1196 upperLimitsExpected["MinPt_3l_offZ_noOSSF_cut_100"] = 0.135;
1197 upperLimitsExpected["MinPt_3l_offZ_noOSSF_cut_150"] = 0.136;
1198 upperLimitsExpected["MinPt_2ltau_offZ_noOSSF_cut_0"] = 4.132;
1199 upperLimitsExpected["MinPt_2ltau_offZ_noOSSF_cut_50"] = 0.339;
1200 upperLimitsExpected["MinPt_2ltau_offZ_noOSSF_cut_100"] = 0.149;
1201 upperLimitsExpected["MinPt_2ltau_offZ_noOSSF_cut_150"] = 0.145;
1202 upperLimitsExpected["MinPt_3l_onZ_cut_0"] = 32.181;
1203 upperLimitsExpected["MinPt_3l_onZ_cut_50"] = 2.260;
1204 upperLimitsExpected["MinPt_3l_onZ_cut_100"] = 0.438;
1205 upperLimitsExpected["MinPt_3l_onZ_cut_150"] = 0.305;
1206 upperLimitsExpected["MinPt_2ltau_onZ_cut_0"] = 217.801;
1207 upperLimitsExpected["MinPt_2ltau_onZ_cut_50"] = 1.335;
1208 upperLimitsExpected["MinPt_2ltau_onZ_cut_100"] = 0.162;
1209 upperLimitsExpected["MinPt_2ltau_onZ_cut_150"] = 0.149;
1210 upperLimitsExpected["nbtag_3l_offZ_OSSF_cut_0"] = 2.893;
1211 upperLimitsExpected["nbtag_3l_offZ_OSSF_cut_1"] = 0.923;
1212 upperLimitsExpected["nbtag_3l_offZ_OSSF_cut_2"] = 0.452;
1213 upperLimitsExpected["nbtag_2ltau_offZ_OSSF_cut_0"] = 14.293;
1214 upperLimitsExpected["nbtag_2ltau_offZ_OSSF_cut_1"] = 1.774;
1215 upperLimitsExpected["nbtag_2ltau_offZ_OSSF_cut_2"] = 0.549;
1216 upperLimitsExpected["nbtag_3l_offZ_noOSSF_cut_0"] = 0.836;
1217 upperLimitsExpected["nbtag_3l_offZ_noOSSF_cut_1"] = 0.594;
1218 upperLimitsExpected["nbtag_3l_offZ_noOSSF_cut_2"] = 0.298;
1219 upperLimitsExpected["nbtag_2ltau_offZ_noOSSF_cut_0"] = 4.132;
1220 upperLimitsExpected["nbtag_2ltau_offZ_noOSSF_cut_1"] = 2.358;
1221 upperLimitsExpected["nbtag_2ltau_offZ_noOSSF_cut_2"] = 0.958;
1222 upperLimitsExpected["nbtag_3l_onZ_cut_0"] = 32.181;
1223 upperLimitsExpected["nbtag_3l_onZ_cut_1"] = 3.868;
1224 upperLimitsExpected["nbtag_3l_onZ_cut_2"] = 0.887;
1225 upperLimitsExpected["nbtag_2ltau_onZ_cut_0"] = 217.801;
1226 upperLimitsExpected["nbtag_2ltau_onZ_cut_1"] = 9.397;
1227 upperLimitsExpected["nbtag_2ltau_onZ_cut_2"] = 0.787;
1228
1229
1230
1231 if (observed) return upperLimitsObserved[signal_region];
1232 else return upperLimitsExpected[signal_region];
1233 }
1234
1235
1236 /// Function checking if there is an OSSF lepton pair or a combination of 3 leptons with an invariant mass close to the Z mass
1237 int isonZ (const Particles& particles) {
1238 int onZ = 0;
1239 double best_mass_2 = 999.;
1240 double best_mass_3 = 999.;
1241
1242 // Loop over all 2 particle combinations to find invariant mass of OSSF pair closest to Z mass
1243 for (const Particle& p1 : particles) {
1244 for (const Particle& p2 : particles) {
1245 double mass_difference_2_old = fabs(91.0 - best_mass_2);
1246 double mass_difference_2_new = fabs(91.0 - (p1.momentum() + p2.momentum()).mass()/GeV);
1247
1248 // If particle combination is OSSF pair calculate mass difference to Z mass
1249 if ((p1.pid()*p2.pid() == -121 || p1.pid()*p2.pid() == -169)) {
1250
1251 // Get invariant mass closest to Z mass
1252 if (mass_difference_2_new < mass_difference_2_old)
1253 best_mass_2 = (p1.momentum() + p2.momentum()).mass()/GeV;
1254 // In case there is an OSSF pair take also 3rd lepton into account (e.g. from FSR and photon to electron conversion)
1255 for (const Particle& p3 : particles ) {
1256 double mass_difference_3_old = fabs(91.0 - best_mass_3);
1257 double mass_difference_3_new = fabs(91.0 - (p1.momentum() + p2.momentum() + p3.momentum()).mass()/GeV);
1258 if (mass_difference_3_new < mass_difference_3_old)
1259 best_mass_3 = (p1.momentum() + p2.momentum() + p3.momentum()).mass()/GeV;
1260 }
1261 }
1262 }
1263 }
1264
1265 // Pick the minimum invariant mass of the best OSSF pair combination and the best 3 lepton combination
1266 double best_mass = min(best_mass_2,best_mass_3);
1267 // if this mass is in a 20 GeV window around the Z mass, the event is classified as onZ
1268 if ( fabs(91.0 - best_mass) < 20. ) onZ = 1;
1269 return onZ;
1270 }
1271
1272 /// function checking if two leptons are an OSSF lepton pair and giving out the invariant mass (0 if no OSSF pair)
1273 double isOSSF_mass (const Particle& p1, const Particle& p2) {
1274 double inv_mass = 0.;
1275 // Is particle combination OSSF pair?
1276 if ((p1.pid()*p2.pid() == -121 || p1.pid()*p2.pid() == -169)) {
1277 // Get invariant mass
1278 inv_mass = (p1.momentum() + p2.momentum()).mass()/GeV;
1279 }
1280 return inv_mass;
1281 }
1282
1283 /// Function checking if there is an OSSF lepton pair
1284 bool isOSSF (const Particles& particles) {
1285 for (size_t i1=0 ; i1 < 3 ; i1 ++) {
1286 for (size_t i2 = i1+1 ; i2 < 3 ; i2 ++) {
1287 if ((particles[i1].pid()*particles[i2].pid() == -121 || particles[i1].pid()*particles[i2].pid() == -169)) {
1288 return true;
1289 }
1290 }
1291 }
1292 return false;
1293 }
1294
1295 /// @}
1296
1297 private:
1298
1299 /// Histograms
1300 /// @{
1301 Histo1DPtr _h_HTlep_all, _h_HTjets_all, _h_MET_all, _h_Meff_all, _h_min_pT_all, _h_mT_all;
1302 Histo1DPtr _h_pt_1_3l, _h_pt_2_3l, _h_pt_3_3l, _h_pt_1_2ltau, _h_pt_2_2ltau, _h_pt_3_2ltau;
1303 Histo1DPtr _h_e_n, _h_mu_n, _h_tau_n;
1304 Histo1DPtr _h_excluded;
1305 /// @}
1306
1307 /// Fiducial efficiencies to model the effects of the ATLAS detector
1308 bool _use_fiducial_lepton_efficiency;
1309
1310 /// List of signal regions and event counts per signal region
1311 vector<string> _signal_regions;
1312 map<string, CounterPtr> _eventCountsPerSR;
1313
1314 };
1315
1316
1317 RIVET_DECLARE_PLUGIN(ATLAS_2014_I1327229);
1318
1319}
|