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