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