gambit is hosted by Hepforge, IPPP Durham
GAMBIT  v1.5.0-2191-ga4742ac
a Global And Modular Bsm Inference Tool
Analysis_ATLAS_13TeV_1LEPStop_36invfb.cpp
Go to the documentation of this file.
1 #include "gambit/cmake/cmake_variables.hpp"
2 #ifndef EXCLUDE_ROOT
3 
4 #include <vector>
5 #include <cmath>
6 #include <memory>
7 #include <iomanip>
8 
13 #include "TLorentzVector.h"
14 #include "TVector3.h"
15 #include "HEPUtils/FastJet.h"
16 #include "TRandom3.h"
17 
18 using namespace std;
19 
20 /* The ATLAS 1 lepton direct stop analysis
21 
22  Based on: https://arxiv.org/abs/1711.11520
23 
24  Code by Martin White (based on ATLAS public code snippet on HepData)
25 
26  KNOWN ISSUES
27 
28  1) Have not added the BDT signal regions (despite having BDT code from ATLAS). They cover a specific kinematic region where the m_stop - m_chi1 mass difference is m_top, which we already know Pythia does badly with.
29 
30  2) We have no equivalent of the ATLAS fakeJER method. Am assuming a 3% JER on every jet for now.
31 
32  3) Have used TLorentzVectors for boosting. Could probably be done without ROOT?
33 
34 */
35 
36 namespace Gambit {
37  namespace ColliderBit {
38 
39  // Need two different functions here for use with std::sort
40  bool sortByPT_1l(const HEPUtils::Jet* jet1, const HEPUtils::Jet* jet2) { return (jet1->pT() > jet2->pT()); }
41  bool sortByPT_1l_sharedptr(std::shared_ptr<HEPUtils::Jet> jet1, std::shared_ptr<HEPUtils::Jet> jet2) { return sortByPT_1l(jet1.get(), jet2.get()); }
42 
43  // Need two different functions here for use with std::sort
44  bool sortByMass_1l(const HEPUtils::Jet* jet1, const HEPUtils::Jet* jet2) { return (jet1->mass() > jet2->mass()); }
45  bool sortByMass_1l_sharedptr(std::shared_ptr<HEPUtils::Jet> jet1, std::shared_ptr<HEPUtils::Jet> jet2) { return sortByMass_1l(jet1.get(), jet2.get()); }
46 
47  double calcMT_1l(HEPUtils::P4 jetMom,HEPUtils::P4 metMom)
48  {
49  //std::cout << "metMom.px() " << metMom.px() << " jetMom PT " << jetMom.pT() << std::endl;
50  double met=sqrt(metMom.px()*metMom.px()+metMom.py()*metMom.py());
51  double dphi = metMom.deltaPhi(jetMom);
52  double mt=sqrt(2*jetMom.pT()*met*(1-std::cos(dphi)));
53  return mt;
54  }
55 
56 
57  class Analysis_ATLAS_13TeV_1LEPStop_36invfb : public Analysis
58  {
59  private:
60 
61  // Numbers passing cuts
62  std::map<string, EventCounter> _counters = {
63  {"tN_med", EventCounter("tN_med")},
64  {"tN_high", EventCounter("tN_high")},
65  {"bWN", EventCounter("bWN")},
66  {"bC2x_diag", EventCounter("bC2x_diag")},
67  {"bC2x_med", EventCounter("bC2x_med")},
68  {"bCbv", EventCounter("bCbv")},
69  {"DM_low_loose", EventCounter("DM_low_loose")},
70  {"DM_low", EventCounter("DM_low")},
71  {"DM_high", EventCounter("DM_high")},
72  {"bffN", EventCounter("bffN")},
73  {"bCsoft_diag", EventCounter("bCsoft_diag")},
74  {"bCsoft_med", EventCounter("bCsoft_med")},
75  {"bCsoft_high", EventCounter("bCsoft_high")},
76  };
77 
78  vector<int> cutFlowVector;
79  vector<string> cutFlowVector_str;
80  int NCUTS; //=16;
81 
82 
83  void LeptonLeptonOverlapRemoval(vector<const HEPUtils::Particle*> &lep1vec, vector<const HEPUtils::Particle*> &lep2vec, double DeltaRMax)
84  {
85 
86  //Routine to do jet-lepton check
87  //Discards jets if they are within DeltaRMax of a lepton
88 
89  vector<const HEPUtils::Particle*> Survivors;
90 
91  for(unsigned int itlep1 = 0; itlep1 < lep1vec.size(); itlep1++)
92  {
93  bool overlap = false;
94  HEPUtils::P4 lep1mom=lep1vec.at(itlep1)->mom();
95  for(unsigned int itlep2 = 0; itlep2 < lep2vec.size(); itlep2++)
96  {
97  HEPUtils::P4 lep2mom=lep2vec.at(itlep2)->mom();
98  double dR;
99 
100  dR=lep1mom.deltaR_eta(lep2mom);
101 
102  if(fabs(dR) <= DeltaRMax) overlap=true;
103  }
104  if(overlap) continue;
105  Survivors.push_back(lep1vec.at(itlep1));
106  }
107  lep1vec=Survivors;
108 
109  return;
110  }
111 
112 
113  void JetLeptonOverlapRemoval(vector<const HEPUtils::Jet*> &jetvec, vector<const HEPUtils::Particle*> &lepvec, double DeltaRMax)
114  {
115  //Routine to do jet-lepton check
116  //Discards jets if they are within DeltaRMax of a lepton
117 
118  vector<const HEPUtils::Jet*> Survivors;
119 
120  for(unsigned int itjet = 0; itjet < jetvec.size(); itjet++)
121  {
122  bool overlap = false;
123  HEPUtils::P4 jetmom=jetvec.at(itjet)->mom();
124  for(unsigned int itlep = 0; itlep < lepvec.size(); itlep++)
125  {
126  HEPUtils::P4 lepmom=lepvec.at(itlep)->mom();
127  double dR;
128 
129  dR=jetmom.deltaR_eta(lepmom);
130 
131  if(fabs(dR) <= DeltaRMax) overlap=true;
132  }
133  if(overlap) continue;
134  Survivors.push_back(jetvec.at(itjet));
135  }
136  jetvec=Survivors;
137 
138  return;
139  }
140 
141 
142  void LeptonJetOverlapRemoval(vector<const HEPUtils::Particle*> &lepvec, vector<const HEPUtils::Jet*> &jetvec)
143  {
144  //Routine to do lepton-jet check
145  //Discards leptons if they are within dR of a jet as defined in analysis paper
146 
147  vector<const HEPUtils::Particle*> Survivors;
148 
149  for(unsigned int itlep = 0; itlep < lepvec.size(); itlep++)
150  {
151  bool overlap = false;
152  HEPUtils::P4 lepmom=lepvec.at(itlep)->mom();
153  for(unsigned int itjet= 0; itjet < jetvec.size(); itjet++)
154  {
155  HEPUtils::P4 jetmom=jetvec.at(itjet)->mom();
156  double dR;
157  double DeltaRMax = std::max(0.1,std::min(0.4, 0.04 + 10 / lepmom.pT()));
158  dR=jetmom.deltaR_eta(lepmom);
159 
160  if(fabs(dR) <= DeltaRMax) overlap=true;
161  }
162  if(overlap) continue;
163  Survivors.push_back(lepvec.at(itlep));
164  }
165  lepvec=Survivors;
166 
167  return;
168  }
169 
170 
171  public:
172 
173  // Required detector sim
174  static constexpr const char* detector = "ATLAS";
175 
176  Analysis_ATLAS_13TeV_1LEPStop_36invfb()
177  {
178 
179  set_analysis_name("ATLAS_13TeV_1LEPStop_36invfb");
180  set_luminosity(36.);
181 
182  NCUTS=150;
183 
184  for(int i=0;i<NCUTS;i++)
185  {
186  cutFlowVector.push_back(0);
187  cutFlowVector_str.push_back("");
188  }
189 
190  }
191 
192  struct ClusteringHistory : public FJNS::PseudoJet::UserInfoBase
193  {
194  enum Status
195  {
196  GOOD,
197  JET_TOO_SMALL,
198  JET_TOO_LARGE,
199  TOO_MANY_ITERATIONS,
200  NONE,
201  };
202 
203  struct Step
204  {
205  double pt;
206  double r;
207  size_t constit;
208  Status status;
209  };
210 
211  size_t id; // a per-event unique jet id that is needed for the event dump
212  std::vector<Step> steps;
213 
214  static ClusteringHistory* AddStep(ClusteringHistory& history, const Step& step)
215  {
216  auto newHistory = new ClusteringHistory(history);
217  newHistory->steps.push_back(step);
218  return newHistory;
219  }
220  };
221 
222  // Return the history of a PseudoJet object, handling all the ugly casting.
223  ClusteringHistory& GetHistory(const FJNS::PseudoJet& jet)
224  {
225  auto shared_ptr = jet.user_info_shared_ptr();
226  return *dynamic_cast<ClusteringHistory*>(shared_ptr.get());
227  }
228 
229  static std::vector<FJNS::PseudoJet> SortedByNConstit(std::vector<FJNS::PseudoJet> jets)
230  {
231  std::sort(jets.begin(), jets.end(), [](const FJNS::PseudoJet& a, const FJNS::PseudoJet& b) {
232  if (a.constituents().size() != b.constituents().size())
233  return a.constituents().size() > b.constituents().size();
234  return a.pt() > b.pt();
235  });
236 
237  return jets;
238  }
239 
240  inline double optimalRadius(const double pT, const double m) { return 2 * m / pT; }
241  inline double minRadius(const double pT, const double m) { return optimalRadius(pT, m) - 0.3; }
242  inline double maxRadius(const double pT, const double m) { return optimalRadius(pT, m) + 0.5; }
243 
244 
245  std::pair<bool, FJNS::PseudoJet> RecursiveRecluster(const FJNS::PseudoJet& candidate, double candRadius,
246  const double mass, size_t step)
247  {
248  if (minRadius(candidate.pt(), mass) > candRadius)
249  {
250  GetHistory(candidate).steps.back().status = ClusteringHistory::JET_TOO_SMALL;
251  return std::make_pair(false, candidate);
252  }
253  else if (maxRadius(candidate.pt(), mass) < candRadius)
254  {
255  const double newR = std::max(maxRadius(candidate.pt(), mass), candRadius / 2.);
256  GetHistory(candidate).steps.back().status = ClusteringHistory::JET_TOO_LARGE;
257 
258  if (step > 10)
259  {
260  GetHistory(candidate).steps.back().status = ClusteringHistory::TOO_MANY_ITERATIONS;
261  return std::make_pair(false, candidate);
262  }
263 
264  FJNS::JetDefinition jetDef(FJNS::antikt_algorithm, newR);
265  auto cs = new FJNS::ClusterSequence(candidate.constituents(), jetDef);
266 
267  std::vector<FJNS::PseudoJet> reclusteredJets;
268  reclusteredJets = SortedByNConstit(cs->inclusive_jets());
269 
270  if (reclusteredJets.size() == 0)
271  {
272  delete cs;
273  return std::make_pair(false, FJNS::PseudoJet());
274  }
275 
276  cs->delete_self_when_unused();
277  auto newCandidate = reclusteredJets[0];
278 
279  auto newHistory = ClusteringHistory::AddStep(
280  GetHistory(candidate),
281  {newCandidate.pt(), newR, newCandidate.constituents().size(), ClusteringHistory::NONE});
282  newCandidate.set_user_info(newHistory);
283 
284  return RecursiveRecluster(newCandidate, newR, mass, step + 1);
285  }
286  else
287  {
288  GetHistory(candidate).steps.back().status = ClusteringHistory::GOOD;
289  return std::make_pair(true, candidate);
290  }
291  }
292 
293 
294  HEPUtils::P4 reclusteredParticle(vector<const HEPUtils::Jet*> jets, vector<const HEPUtils::Jet*> bjets,
295  const double mass, const bool useBJets)
296  {
297 
298  //AnalysisObject p = AnalysisObject(0., 0., 0., 0., 0, 0, AnalysisObjectType::JET, 0, 0);
299  HEPUtils::P4 p;
300  double r0 = 3.0;
301 
302  vector<const HEPUtils::Jet*> usejets;
303  for(const HEPUtils::Jet* jet : jets)
304  {
305  usejets.push_back(jet);
306  }
307 
308  if (useBJets && bjets.size())
309  {
310  for(const HEPUtils::Jet* bjet : bjets)
311  {
312  usejets.push_back(bjet);
313  }
314  }
315 
316  std::vector<FJNS::PseudoJet> initialJets;
317 
318  for (const HEPUtils::Jet* jet : usejets)
319  {
320  FJNS::PseudoJet Pjet(jet->mom().px(), jet->mom().py(), jet->mom().pz(), jet->mom().E());
321  initialJets.push_back(Pjet);
322  }
323 
324  FJNS::JetDefinition jetDef(FJNS::antikt_algorithm, r0);
325  FJNS::ClusterSequence cs(initialJets, jetDef);
326 
327  auto candidates = FJNS::sorted_by_pt(cs.inclusive_jets());
328 
329  std::vector<FJNS::PseudoJet> selectedJets;
330  selectedJets.reserve(candidates.size());
331  std::vector<FJNS::PseudoJet> badJets;
332  badJets.reserve(candidates.size());
333 
334  size_t i = 0;
335  for (auto& cand : candidates)
336  {
337  auto history = new ClusteringHistory();
338  history->id = i;
339  history->steps.push_back({cand.pt(), r0, cand.constituents().size(), ClusteringHistory::NONE});
340  cand.set_user_info(history);
341  ++i;
342  }
343 
344  for (const auto& cand : candidates)
345  {
346  bool selected = false;
347  FJNS::PseudoJet jet;
348 
349  std::tie(selected, jet) = RecursiveRecluster(cand, r0, mass, 0);
350 
351  if (selected)
352  selectedJets.push_back(jet);
353  else
354  badJets.push_back(jet);
355  }
356 
357  if (selectedJets.size() < 1)
358  {
359  return p;
360  }
361 
362  vector<std::shared_ptr<HEPUtils::Jet>> aoSelectedJets;
363  for (const FJNS::PseudoJet& j : selectedJets) aoSelectedJets.push_back(std::make_shared<HEPUtils::Jet>(HEPUtils::mk_p4(j)));
364 
365  //for (const auto jet : selectedJets)
366  // aoSelectedJets.push_back(
367  // AnalysisObject(jet.px(), jet.py(), jet.pz(), jet.E(), 0, 0, AnalysisObjectType::COMBINED, 0, 0));
368 
369  std::sort(aoSelectedJets.begin(), aoSelectedJets.end(), sortByPT_1l_sharedptr);
370  p = aoSelectedJets[0]->mom();
371 
372  return p;
373  }
374 
375 
376  void run(const HEPUtils::Event* event)
377  {
378 
379  // Missing energy
380  HEPUtils::P4 metVec = event->missingmom();
381  double Met = event->met();
382 
383  // Construct baseline electron objects
384  vector<const HEPUtils::Particle*> baselineElectrons;
385  for (const HEPUtils::Particle* electron : event->electrons())
386  {
387  if (electron->pT() > 5. && electron->abseta() < 2.47)
388  {
389  baselineElectrons.push_back(electron);
390  }
391  }
392 
393  // Apply electron efficiency
394  ATLAS::applyElectronEff(baselineElectrons);
395 
396  // Construct baseline muon objects
397  vector<const HEPUtils::Particle*> baselineMuons;
398  for (const HEPUtils::Particle* muon : event->muons())
399  {
400  if (muon->pT() > 4. && muon->abseta() < 2.7)
401  {
402  baselineMuons.push_back(muon);
403  }
404  }
405 
406  // Apply muon efficiency
407  ATLAS::applyMuonEff(baselineMuons);
408 
409  // Construct set of all light baseline leptons
410  vector<const HEPUtils::Particle*> baselineLeptons = baselineElectrons;
411  baselineLeptons.insert(baselineLeptons.end(), baselineMuons.begin(), baselineMuons.end() );
412 
413  // Construct baseline tau objects
414  vector<const HEPUtils::Particle*> baselineTaus;
415  for (const HEPUtils::Particle* tau : event->taus())
416  {
417  if (tau->pT() > 20. && fabs(tau->eta()) < 2.5) baselineTaus.push_back(tau);
418  }
419  // Apply tau efficiency
420  ATLAS::applyTauEfficiencyR1(baselineTaus);
421 
422  // Photons
423  vector<const HEPUtils::Particle*> signalPhotons;
424  for (const HEPUtils::Particle* photon : event->photons())
425  {
426  signalPhotons.push_back(photon);
427  }
428 
429  // Jets
430  vector<const HEPUtils::Jet*> bJets;
431  vector<const HEPUtils::Jet*> nonBJets;
432  vector<const HEPUtils::Jet*> trueBJets; //for debugging
433 
434  // Get b jets
436  const std::vector<double> a = {0,10.};
437  const std::vector<double> b = {0,10000.};
438  const std::vector<double> c = {0.77}; // set b-tag efficiency to 77%
439  HEPUtils::BinnedFn2D<double> _eff2d(a,b,c);
440  for (const HEPUtils::Jet* jet : event->jets())
441  {
442  bool hasTag=has_tag(_eff2d, fabs(jet->eta()), jet->pT());
443  if (jet->pT() > 20. && fabs(jet->eta()) < 4.9)
444  {
445  if(jet->btag() && hasTag && fabs(jet->eta()) < 2.5 && jet->pT() > 20.)
446  {
447  bJets.push_back(jet);
448  }
449  else
450  {
451  nonBJets.push_back(jet);
452  }
453  }
454  }
455 
456  // Note: use paper description instead of code snippet
457  // This is not identical to the overlap removal in the paper
458  // Probably good enough though
459  LeptonLeptonOverlapRemoval(baselineMuons,baselineElectrons,0.01); // mimics shared track requirement
460  JetLeptonOverlapRemoval(nonBJets,baselineElectrons,0.2);
461  LeptonJetOverlapRemoval(baselineElectrons,nonBJets);
462  LeptonJetOverlapRemoval(baselineElectrons,bJets);
463  LeptonJetOverlapRemoval(baselineMuons,nonBJets);
464  LeptonJetOverlapRemoval(baselineMuons,bJets);
465  LeptonLeptonOverlapRemoval(baselineTaus,baselineElectrons,0.1);
466 
467  // Fill a jet-pointer-to-bool map to make it easy to check
468  // if a given jet is treated as a b-jet in this analysis
469  map<const HEPUtils::Jet*,bool> analysisBtags;
470  for (const HEPUtils::Jet* jet : bJets) {
471  analysisBtags[jet] = true;
472  }
473  for (const HEPUtils::Jet* jet : nonBJets) {
474  analysisBtags[jet] = false;
475  }
476 
477  // Signal object containers
478  vector<const HEPUtils::Particle*> signalElectrons;
479  vector<const HEPUtils::Particle*> signalSoftElectrons;
480  vector<const HEPUtils::Particle*> signalMuons;
481  vector<const HEPUtils::Particle*> signalSoftMuons;
482  vector<const HEPUtils::Particle*> signalLeptons;
483  vector<const HEPUtils::Particle*> signalSoftLeptons;
484  vector<const HEPUtils::Particle*> electronsForVeto;
485  vector<const HEPUtils::Particle*> muonsForVeto;
486 
487  vector<const HEPUtils::Jet*> signalJets;
488  vector<const HEPUtils::Jet*> signalBJets;
489  vector<const HEPUtils::Jet*> signalNonBJets;
490 
491  // Now apply signal jet cuts
492  for (const HEPUtils::Jet* jet : bJets)
493  {
494  if(jet->pT() > 25. && fabs(jet->eta())<2.5)
495  {
496  signalJets.push_back(jet);
497  signalBJets.push_back(jet);
498  }
499  }
500 
501  for (const HEPUtils::Jet* jet : nonBJets)
502  {
503  if(jet->pT() > 25. && fabs(jet->eta())<2.5)
504  {
505  signalJets.push_back(jet);
506  signalNonBJets.push_back(jet);
507  }
508  }
509 
510  // Note that the isolation requirements and tight selection are currently missing
511 
512  for (const HEPUtils::Particle* electron : baselineElectrons)
513  {
514  signalSoftElectrons.push_back(electron);
515  signalSoftLeptons.push_back(electron);
516  if(electron->pT() > 25.)
517  {
518  signalElectrons.push_back(electron);
519  signalLeptons.push_back(electron);
520  }
521  }
522 
523  for (const HEPUtils::Particle* muon : baselineMuons)
524  {
525  signalSoftMuons.push_back(muon);
526  signalSoftLeptons.push_back(muon);
527  if(muon->pT() > 25.)
528  {
529  signalMuons.push_back(muon);
530  signalLeptons.push_back(muon);
531  }
532  }
533 
534  // We now have the signal electrons, muons, jets and b jets- move on to the analysis
535 
536  int nJets=signalJets.size();
537  int nBJets = signalBJets.size();
538 
539  // Minimal event selection
540  bool cut_minSelection=false;
541  if((Met > 100. && (baselineElectrons.size()+baselineMuons.size()) == 1 &&
542  ((signalSoftLeptons.size() == 1 || signalLeptons.size() == 1)) && nJets > 1)) cut_minSelection=true;
543 
544  vector<const HEPUtils::Jet*> mostBjetLike;
545  vector<const HEPUtils::Jet*> signalNotBjetLike;
546  vector<const HEPUtils::Jet*> signalNotBjet;
547 
548  // create containers with exactly 2 jets being considered to be b-jets and the inverse
549  int bJet1 = -1, bJet2 = -1;
550 
551  for (unsigned int i = 0; i < signalJets.size(); ++i)
552  {
553  if (!analysisBtags.at(signalJets[i])) continue;
554  if (bJet1 == -1)
555  bJet1 = i;
556  else if (bJet2 == -1)
557  {
558  bJet2 = i;
559  break;
560  }
561  }
562  if (bJet2 == -1)
563  {
564  for (unsigned int i = 0; i < signalJets.size(); ++i)
565  {
566  if (analysisBtags.at(signalJets[i])) continue;
567  if (bJet1 == -1)
568  bJet1 = i;
569  else if (bJet2 == -1)
570  {
571  bJet2 = i;
572  break;
573  }
574  }
575  }
576 
577  if(signalJets.size()>1)
578  {
579  mostBjetLike.push_back(signalJets.at(bJet1));
580  mostBjetLike.push_back(signalJets.at(bJet2));
581  }
582 
583  for (int i = 0; i < (int)signalJets.size(); ++i)
584  {
585  if (!analysisBtags.at(signalJets[i])) signalNotBjet.push_back(signalJets.at(i));
586  if (i == bJet1 || i == bJet2) continue;
587  signalNotBjetLike.push_back(signalJets.at(i));
588  }
589 
590  /* ensure object collections to be pT sorted */
591  std::sort(signalJets.begin(), signalJets.end(), sortByPT_1l);
592  //if (baseTaus.size() > 0) sortObjectsByPt(baseTaus);
593 
594  // Now make a collection to hold the JER for each jet
595  // Have obtained the values from Matthias' BuckFast code
596  // https://atlas.web.cern.ch/Atlas/GROUPS/PHYSICS/CONFNOTES/ATLAS-CONF-2015-017/
597  // Parameterisation can be still improved, but eta dependence is minimal
598  const std::vector<double> binedges_eta = {0,10.};
599  const std::vector<double> binedges_pt = {0,50.,70.,100.,150.,200.,1000.,10000.};
600  const std::vector<double> JetsJER = {0.145,0.115,0.095,0.075,0.07,0.05,0.04};
601  static HEPUtils::BinnedFn2D<double> _resJets2D(binedges_eta,binedges_pt,JetsJER);
602  vector<double> signalJER;
603 
604  for(unsigned int i = 0; i < signalJets.size(); ++i)signalJER.push_back(_resJets2D.get_at(signalJets[i]->abseta(), signalJets[i]->pT()));
605 
606  float sigmaAbsHtMiss = 0;
607  float Ht = 0;
608  /* calculate vecHtMiss */
609  HEPUtils::P4 vecHtMiss;
610  HEPUtils::P4 leptonHtMiss;
611 
612  bool preselLowMet=false;
613  bool preselHighMet=false;
614  double MetPerp = 0.;
615  double HtSigMiss=0.;
616  double absDPhiJMet0 = 0.;
617  double absDPhiJMet1 = 0.;
618  double absDPhiJiMet = 0.;
619  double mT=0.;
620  double topReclM=0.;
621  double WReclM=0.;
622  double amT2=0.;
623  double dRbl=9999.;
624  double mT2Tau=0.;
625  double dPhiMetLep;
626  double pTLepOverMet=999.;
627 
628  if(cut_minSelection)
629  {
630 
631  for (unsigned int i = 0; i < baselineLeptons.size(); ++i)
632  {
633  vecHtMiss -= baselineLeptons[i]->mom();
634  leptonHtMiss -= baselineLeptons[i]->mom();
635  }
636 
637  for (unsigned int i = 0; i < signalJets.size(); ++i)
638  vecHtMiss -= signalJets[i]->mom();
639 
640  /* calculate Ht and HtSig */
641  for (const HEPUtils::Jet* jet : signalJets) Ht += jet->pT();
642 
643  TRandom3 myRandom;
644  myRandom.SetSeed(signalJets[0]->pT());
645 
646  int PEs = 100;
647  // double smear_factor;
648  double ETmissmean = 0, ETmissRMS = 0;
649  for (int j = 0; j < PEs; ++j)
650  {
651  double jetHtx = leptonHtMiss.px();
652  double jetHty = leptonHtMiss.py();
653 
654  for (unsigned int i = 0; i < signalJets.size(); ++i)
655  {
656  //std::normal_distribution<> dx(signalJets[i]->mom().px(), signalJets[i]->mom().px() * signalJER[i]);
657  //std::normal_distribution<> dy(signalJets[i]->mom().py(), signalJets[i]->mom().px() * signalJER[i]);
658  //jetHtx -= dx(Random::rng());
659  //jetHty -= dy(Random::rng());
660  jetHtx -= myRandom.Gaus(signalJets[i]->mom().px(), signalJets[i]->mom().px() * signalJER[i]);
661  jetHty -= myRandom.Gaus(signalJets[i]->mom().py(), signalJets[i]->mom().px() * signalJER[i]);
662  }
663  double ETtemp = sqrt(jetHtx * jetHtx + jetHty * jetHty);
664  ETmissmean += ETtemp;
665  ETmissRMS += ETtemp * ETtemp;
666  }
667 
668  ETmissmean = ETmissmean / PEs;
669  sigmaAbsHtMiss = sqrt((ETmissRMS / PEs) - ETmissmean * ETmissmean);
670 
671  HtSigMiss = (ETmissmean - 100.) / sigmaAbsHtMiss;
672 
673  double absDPhiJMet[4] = {fabs(signalJets[0]->mom().deltaPhi(metVec)), fabs(signalJets[1]->mom().deltaPhi(metVec)),
674  signalJets.size() > 2 ? fabs(signalJets[2]->mom().deltaPhi(metVec)) : NAN,
675  signalJets.size() > 3 ? fabs(signalJets[3]->mom().deltaPhi(metVec)) : NAN};
676 
677  if(nJets>0)absDPhiJMet0 = absDPhiJMet[0];
678  if(nJets>1)absDPhiJMet1 = absDPhiJMet[1];
679 
680  for (int i = 1; i < 4; i++)
681  if (absDPhiJMet[i] < absDPhiJiMet) absDPhiJiMet = absDPhiJMet[i];
682 
683  mT = calcMT_1l(baselineLeptons[0]->mom(), metVec);
684  dPhiMetLep = fabs(metVec.deltaPhi(baselineLeptons[0]->mom()));
685 
686  // Calculate MT2 tau using the leading tau in the event
687  mT2Tau = 120.;
688  if(baselineTaus.size() > 0)
689  {
690  double pa_tau[3] = { 0, baselineTaus[0]->mom().px(), baselineTaus[0]->mom().py() };
691  double pb_tau[3] = { 0, baselineLeptons[0]->mom().px(), baselineLeptons[0]->mom().py() };
692  double pmiss_tau[3] = { 0, metVec.px(), metVec.py() };
693  double mn_tau = 0.;
694  mt2_bisect::mt2 mt2_event_tau;
695  mt2_event_tau.set_momenta(pa_tau,pb_tau,pmiss_tau);
696  mt2_event_tau.set_mn(mn_tau);
697 
698  mT2Tau = mt2_event_tau.get_mt2();
699  }
700 
701  pTLepOverMet = baselineLeptons[0]->pT() / Met;
702  preselHighMet = Met > 230 && mT > 30;
703  preselLowMet = baselineLeptons[0]->pT() > 27 && signalBJets.size() > 0 && signalJets[0]->pT() > 50. && Met > 100 && mT > 90;
704 
705  // Apply tight selection if lepton is an electron
706  // Am using same selection as 8 TeV (probably needs updating)
707  // Note that we have already applied a 1 lepton cut
708  if (baselineElectrons.size()==1 && baselineMuons.size()==0)
709  {
710  vector<const HEPUtils::Particle*> tightElectrons;
711  tightElectrons.push_back(baselineElectrons[0]);
712  ATLAS::applyTightIDElectronSelection(tightElectrons);
713  preselLowMet = preselLowMet && (tightElectrons.size()==1);
714  }
715 
716  // Now calculate amT2 using two different assignments of the b jets and the leptons
717 
718  HEPUtils::P4 lepton_plus_bjet0;
719  HEPUtils::P4 lepton_plus_bjet1;
720 
721  lepton_plus_bjet0 = baselineLeptons[0]->mom()+mostBjetLike[0]->mom();
722  lepton_plus_bjet1 = baselineLeptons[0]->mom()+mostBjetLike[1]->mom();
723 
724  double pa_a[3] = { 0, lepton_plus_bjet0.px(), lepton_plus_bjet0.py() };
725  double pb_a[3] = { 80, mostBjetLike[1]->mom().px(), mostBjetLike[1]->mom().py() };
726  double pmiss_a[3] = { 0, metVec.px(), metVec.py() };
727  double mn_a = 0.;
728 
729  mt2_bisect::mt2 mt2_event_a;
730 
731  mt2_event_a.set_momenta(pa_a,pb_a,pmiss_a);
732  mt2_event_a.set_mn(mn_a);
733 
734  double mt2a = mt2_event_a.get_mt2();
735 
736  double pa_b[3] = { 0, lepton_plus_bjet1.px(), lepton_plus_bjet1.py() };
737  double pb_b[3] = { 80, mostBjetLike[0]->mom().px(), mostBjetLike[0]->mom().py() };
738  double pmiss_b[3] = { 0, metVec.px(), metVec.py() };
739  double mn_b = 0.;
740 
741  mt2_bisect::mt2 mt2_event_b;
742 
743  mt2_event_b.set_momenta(pa_b,pb_b,pmiss_b);
744  mt2_event_b.set_mn(mn_b);
745  double mt2b = mt2_event_b.get_mt2();
746 
747  amT2 = std::min(mt2a,mt2b);
748  dRbl = baselineLeptons[0]->mom().deltaR_eta(mostBjetLike[0]->mom());
749 
750  /* Reconstruct top by a chi2 based method */
751  float mW = 80.;
752  float mTop = 170.;
753  float chi2min = 9e99;
754  //AnalysisObject* topChi2 = new AnalysisObject(0., 0., 0., 0., 0, 0, AnalysisObjectType::JET, 0, 0);
755  HEPUtils::P4 topChi2;
756  int jetComb[3] = {0, 0, 0};
757  vector<double> signalBJER;
758 
759  for(unsigned int i = 0; i < mostBjetLike.size(); ++i)signalBJER.push_back(_resJets2D.get_at(mostBjetLike[i]->abseta(), mostBjetLike[i]->pT()));
760  float f;
761 
762  for (int i = 0; i < (int)signalJets.size(); ++i)
763  {
764  if (i == bJet1 || i == bJet2) continue;
765  for (int j = i + 1; j < (int)signalJets.size(); ++j)
766  {
767  if (j == bJet1 || j == bJet2) continue;
768  for (unsigned int k = 0; k < mostBjetLike.size() && k < 2; ++k)
769  {
770  f = pow((signalJets[i]->mom() + signalJets[j]->mom() + mostBjetLike[k]->mom()).m() - mTop, 2) /
771  (pow((signalJets[i]->mom() + signalJets[j]->mom() + mostBjetLike[k]->mom()).m(), 2) *
772  (pow(signalJER[i], 2) + pow(signalJER[j], 2) + pow(signalBJER[k], 2))) +
773  pow((signalJets[i]->mom() + signalJets[j]->mom()).m() - mW, 2) /
774  (pow((signalJets[i]->mom() + signalJets[j]->mom()).m(), 2) * (pow(signalJER[i], 2) + pow(signalJER[j], 2)));
775  if (f < chi2min)
776  {
777  chi2min = f;
778  jetComb[0] = i;
779  jetComb[1] = j;
780  jetComb[2] = k;
781  }
782  }
783  }
784  }
785  topChi2 = signalJets[jetComb[0]]->mom() + signalJets[jetComb[1]]->mom() + mostBjetLike[jetComb[2]]->mom();
786 
787 
788  HEPUtils::P4 top1;
789  top1 = baselineLeptons[0]->mom() + (jetComb[2] == 0 ? mostBjetLike[1]->mom() : mostBjetLike[0]->mom());
790 
791  /* calculate MetPerp */
792 
793  TLorentzVector ttbar;
794  ttbar.SetPxPyPzE((topChi2 + top1).px(),(topChi2 + top1).py(),(topChi2 + top1).pz(),(topChi2 + top1).E());
795  TLorentzVector top1Rest;
796  top1Rest.SetPxPyPzE(top1.px(),top1.py(),top1.pz(),top1.E());
797  TLorentzVector metRest;
798  metRest.SetPxPyPzE(metVec.px(),metVec.py(),metVec.pz(),metVec.E());
799 
800 
801  ttbar.Boost(-ttbar.Px() / ttbar.E(), -ttbar.Py() / ttbar.E(), -ttbar.Pz() / ttbar.E());
802 
803 
804  top1Rest.Boost(-ttbar.Px() / ttbar.E(), -ttbar.Py() / ttbar.E(), -ttbar.Pz() / ttbar.E());
805  metRest.Boost(-ttbar.Px() / ttbar.E(), -ttbar.Py() / ttbar.E(), -ttbar.Pz() / ttbar.E());
806  MetPerp = metRest.Vect().XYvector().Norm(top1Rest.Vect().XYvector()).Mod();
807 
808  // Now we have to do the fancy jet reclustering to get reconstructed W and top particles
809 
810  HEPUtils::P4 WRecl = reclusteredParticle(signalNotBjet, mostBjetLike, mW, false); //signalNotBjet+mostBjetLike is inconsistent but bjets are not used anyway
811  HEPUtils::P4 topRecl = reclusteredParticle(signalNotBjetLike, mostBjetLike, 175., true);
812 
813  topReclM=0;
814 
815  if (nBJets > 0 && nJets > 3 && preselHighMet)topReclM=topRecl.m();
816  WReclM = WRecl.m();
817 
818  }
819 
820  // Should now be ready to do signal selections
821 
822  bool is_tN_med=false;
823  bool is_tN_high=false;
824  bool is_bWN=false;
825  bool is_bC2x_diag=false;
826  bool is_bC2x_med=false;
827  bool is_bCbv=false;
828  // bool is_DM_low_loose=false; // <-- We currently don't use this
829  bool is_DM_low=false;
830  bool is_DM_high=false;
831 
832  bool is_bffN=false;
833  bool is_bCsoft_diag=false;
834  bool is_bCsoft_med=false;
835  bool is_bCsoft_high=false;
836 
837  // non-soft lepton selections
838  if (signalLeptons.size() == 1)
839  {
840  //tN_med
841  if (nJets > 3 && nBJets > 0 && preselLowMet && signalJets[0]->pT() > 60 && signalJets[1]->pT() > 50 &&
842  signalJets[3]->pT() > 40 && Met > 250 && MetPerp > 230 && HtSigMiss > 14 && mT > 160 && amT2 > 175 &&
843  topReclM > 150 && dRbl < 2.0 && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80)
844  is_tN_med=true;
845 
846  //tN_high
847  if (nJets > 3 && nBJets > 0 && preselHighMet && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 80 &&
848  signalJets[2]->pT() > 50 && signalJets[3]->pT() > 30 && Met > 550 && HtSigMiss > 27 && mT > 160 &&
849  amT2 > 175 && topReclM > 130 && dRbl < 2.0 && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 &&
850  mT2Tau > 80)
851  is_tN_high=true;
852 
853  //bWN
854  if (nJets > 3 && nBJets > 0 && preselHighMet && signalJets[0]->pT() > 50 && Met > 300 && mT > 130 &&
855  amT2 < 110 && dPhiMetLep < 2.5 && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80)
856  is_bWN=true;
857 
858  //bC2x_diag
859  if (nJets > 3 && nBJets > 1 && preselHighMet && signalJets[2]->pT() > 75 && signalJets[3]->pT() > 30 &&
860  signalBJets[1]->pT() > 30 && Met > 230 && HtSigMiss > 13 && mT > 180 && amT2 > 175 &&
861  absDPhiJMet0 > 0.7 && absDPhiJMet1 > 0.7 && WReclM > 50 && mT2Tau > 80)
862  is_bC2x_diag=true;
863 
864  //bC2x_med
865  if (nJets > 3 && nBJets > 1 && preselHighMet && signalJets[0]->pT() > 200 && signalJets[1]->pT() > 140 &&
866  signalBJets[1]->pT() > 140 && Met > 230 && HtSigMiss > 10 && mT > 120 && amT2 > 300 &&
867  absDPhiJMet0 > 0.9 && absDPhiJMet1 > 0.9 && WReclM > 50 && mT2Tau > 80)
868  is_bC2x_med=true;
869 
870  //bCbv
871  if (nJets > 1 && nBJets == 0 && preselHighMet && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 80 &&
872  Met > 360 && HtSigMiss > 16 && mT > 200 && absDPhiJMet0 > 2.0 && absDPhiJMet1 > 0.8 &&
873  WReclM >= 70 && WReclM <= 100 && dPhiMetLep > 1.2 && baselineLeptons[0]->pT() > 60)
874  is_bCbv=true;
875 
876  // We currently don't use this.
877  // //DM_low_loose
878  // if (nJets > 3 && nBJets > 0 && preselHighMet && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40 &&
879  // Met > 300 && mT > 120 && HtSigMiss > 14 && amT2 > 140 && dPhiMetLep > 0.8 && absDPhiJiMet > 1.4)
880  // is_DM_low_loose=true;
881 
882  //DM_low
883  if (nJets > 3 && nBJets > 0 && preselHighMet && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 85 &&
884  signalJets[2]->pT() > 65 && signalBJets[0]->pT() > 60 && Met > 320 && mT > 170 && HtSigMiss > 14 &&
885  amT2 > 160 && topReclM > 130 && dPhiMetLep > 1.2 && absDPhiJiMet > 1.0 && mT2Tau > 80)
886  is_DM_low=true;
887 
888  //DM_high
889  if (nJets > 3 && nBJets > 0 && preselHighMet && signalJets[0]->pT() > 125 && signalJets[1]->pT() > 75 &&
890  signalJets[2]->pT() > 65 && Met > 380 && mT > 225 && amT2 > 190 && topReclM > 130 &&
891  dPhiMetLep > 1.2 && absDPhiJiMet > 1.0)
892  is_DM_high=true;
893  }
894 
895  // Soft-lepton selections
896  bool preselSoftLep=false;
897  double Wpt = 0.;
898 
899  double minDPhiMetBJet = 99999.;
900  for(size_t i=0;i<signalBJets.size();i++)
901  {
902  double dPhi_tmp = fabs(signalBJets[i]->mom().deltaPhi(metVec));
903  if(dPhi_tmp<minDPhiMetBJet)minDPhiMetBJet=dPhi_tmp;
904  }
905 
906  double dRbb=0;
907  if(nBJets>1)dRbb=mostBjetLike[0]->mom().deltaR_eta(mostBjetLike[1]->mom());
908 
909 
910  if (signalSoftLeptons.size() == 1)
911  {
912  preselSoftLep = Met > 230;
913 
914  // Apply tight selection if lepton is an electron
915  // Am using same selection as 8 TeV (probably needs updating)
916  // Note that we have already applied a 1 lepton cut
917  if (signalSoftElectrons.size()==1 && signalSoftMuons.size()==0)
918  {
919  vector<const HEPUtils::Particle*> tightElectrons;
920  tightElectrons.push_back(signalSoftElectrons[0]);
921  ATLAS::applyTightIDElectronSelection(tightElectrons);
922  preselSoftLep = preselSoftLep && (tightElectrons.size()==1);
923  }
924 
925  Wpt = (signalSoftLeptons[0]->mom()+metVec).pT();
926  }
927 
928  //bffN
929  if (nJets > 1 && nBJets > 0 && preselSoftLep && signalJets[0]->pT() > 400 && Met > 300 && mT < 160 &&
930  pTLepOverMet < 0.02 && minDPhiMetBJet < 1.5 && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 &&
931  topReclM < 150 && !analysisBtags.at(signalJets[0]))is_bffN=true;
932 
933  //bCsoft_diag
934  if (nJets > 1 && nBJets > 0 && preselSoftLep && signalJets[0]->pT() > 400 && Met > 300 && mT < 50 &&
935  pTLepOverMet < 0.02 && minDPhiMetBJet < 1.5 && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 &&
936  topReclM < 150 && !analysisBtags.at(signalJets[0]))is_bCsoft_diag=true;
937 
938  //bCsoft_med
939  if (nJets > 2 && nBJets > 1 && preselSoftLep && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 60 &&
940  signalJets[2]->pT() > 40 && signalBJets[0]->pT() > 120 && signalBJets[1]->pT() > 60 && Met > 230 &&
941  mT < 160 && pTLepOverMet < 0.03 && amT2 > 200 && minDPhiMetBJet > 0.8 && absDPhiJMet0 > 0.4 &&
942  absDPhiJMet1 > 0.4 && Wpt > 400)is_bCsoft_med=true;
943 
944  //bCsoft_high
945  if (nJets > 1 && nBJets > 1 && preselSoftLep && signalJets[1]->pT() > 100 && signalBJets[1]->pT() > 100 &&
946  Met > 230 && mT < 160 && pTLepOverMet < 0.03 && amT2 > 300 && minDPhiMetBJet > 0.4 &&
947  absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && Wpt > 500 && dRbb > 0.8)is_bCsoft_high=true;
948 
949  //bool isSRD_high=false;
950 
951  cutFlowVector_str[0] = "No cuts ";
952  cutFlowVector_str[1] = "Derivation skim";
953  cutFlowVector_str[2] = ">=1 baseline lepton ";
954  cutFlowVector_str[3] = ">=1 signal lepton ";
955  cutFlowVector_str[4] = "==1 signal lepton ";
956  cutFlowVector_str[5] = "==1 baseline lepton ";
957  cutFlowVector_str[6] = "XE trigger, >=4 jets, met > 230 GeV ";
958  cutFlowVector_str[7] = "deltaPhi(j1,met) > 0.4 ";
959  cutFlowVector_str[8] = "deltaPhi(j2,met) > 0.4 ";
960  cutFlowVector_str[9] = "mT2tau > 80 GeV ";
961  cutFlowVector_str[10] = "tN_med: j0 pT > 60 GeV";
962  cutFlowVector_str[11] = "tN_med: j1 pT > 50 GeV";
963  cutFlowVector_str[12] = "tN_med: j2 pT > 40 GeV ";
964  cutFlowVector_str[13] = "tN_med: j3 pT > 40 GeV ";
965  cutFlowVector_str[14] = "tN_med: met > 250 GeV ";
966  cutFlowVector_str[15] = "tN_med: metPerp > 230 GeV";
967  cutFlowVector_str[16] = "tN_med: HTmissSig > 14 ";
968  cutFlowVector_str[17] = "tN_med: mT > 160 GeV";
969  cutFlowVector_str[18] = "tN_med: amt2 > 175 GeV";
970  cutFlowVector_str[19] = "tN_med: >=1 b jet ";
971  cutFlowVector_str[20] = "tN_med: deltaR(b,l) < 2.0";
972  cutFlowVector_str[21] = "tN_med: mtop_recl > 150 GeV";
973  cutFlowVector_str[22] = "tN_high: j0 pT > 100 GeV ";
974  cutFlowVector_str[23] = "tN_high: j1 pT > 80 GeV ";
975  cutFlowVector_str[24] = "tN_high: j2 pT > 50 GeV";
976  cutFlowVector_str[25] = "tN_high: j3 pT > 30 GeV";
977  cutFlowVector_str[26] = "tN_high: met > 550 GeV ";
978  cutFlowVector_str[27] = "tN_high: HTmissSig > 27";
979  cutFlowVector_str[28] = "tN_high: mT > 160 GeV ";
980  cutFlowVector_str[29] = "tN_high: amT2 > 175 GeV ";
981  cutFlowVector_str[30] = "tN_high: >= 1 b jet";
982  cutFlowVector_str[31] = "tN_high: deltaR(b,l) < 2.0";
983  cutFlowVector_str[32] = "tN_high: mtop_recl > 130 GeV";
984  cutFlowVector_str[33] = "bWN: jet0 pT > 50 GeV";
985  cutFlowVector_str[34] = "bWN: Met > 300 GeV ";
986  cutFlowVector_str[35] = "bWN: mT > 130 GeV";
987  cutFlowVector_str[36] = "bWN: amT2 < 110 GeV";
988  cutFlowVector_str[37] = "bWN: >=1 b jet";
989  cutFlowVector_str[38] = "bWN: deltaPhi(l,ptmiss) < 2.5";
990  cutFlowVector_str[39] = "bffN: Soft lepton preselection";
991  cutFlowVector_str[40] = "bffN: Met > 300 GeV";
992  cutFlowVector_str[41] = "bffN: jet0 pT > 400 GeV";
993  cutFlowVector_str[42] = "bffN: mT < 160 GeV";
994  cutFlowVector_str[43] = "bffN: leading jet not b-tagged";
995  cutFlowVector_str[44] = "bffN: min(DPhi(ptmiss, b-jet)) < 1.5";
996  cutFlowVector_str[45] = "bffN: pTl/Met < 0.05 ";
997  cutFlowVector_str[46] = "bffN: top veto (or mtop_recl < 150 GeV) ";
998  cutFlowVector_str[47] = "bffN: pTl/met < 0.02 ";
999  cutFlowVector_str[48] = "bC2x_diag: jet0 pT > 75 GeV";
1000  cutFlowVector_str[49] = "bC2x_diag: jet2 pT > 75 GeV ";
1001  cutFlowVector_str[50] = "bC2x_diag: jet3 pT > 75 GeV ";
1002  cutFlowVector_str[51] = "bC2x_diag: jet4 pT > 30 GeC ";
1003  cutFlowVector_str[52] = "bC2x_diag: >=2 b jet ";
1004  cutFlowVector_str[53] = "bC2x_diag: bjet1 pT > 30 GeV ";
1005  cutFlowVector_str[54] = "bC2x_diag: dPhi(j1,ptmiss) > 0.7 ";
1006  cutFlowVector_str[55] = "bC2x_diag: dPhi(j2,ptmiss) > 0.7 ";
1007  cutFlowVector_str[56] = "bC2x_diag: HTmissSig > 13 ";
1008  cutFlowVector_str[57] = "bC2x_diag: mT > 180 GeV ";
1009  cutFlowVector_str[58] = "bC2x_diag: amT2 > 175 GeV ";
1010  cutFlowVector_str[59] = "bC2x_diag: mWrecl > 50 GeV ";
1011  cutFlowVector_str[60] = "bC2x_med: jet0 pT > 200 GeV ";
1012  cutFlowVector_str[61] = "bC2x_med: jet1 pT > 140 GeV ";
1013  cutFlowVector_str[62] = "bC2x_med: >=2 b jet ";
1014  cutFlowVector_str[63] = "bC2x_med: bjet0 pT > 140 GeV ";
1015  cutFlowVector_str[64] = "bC2x_med: bjet1 pT > 140 GeV";
1016  cutFlowVector_str[65] = "bC2x_med: dPhi(j1,ptmiss) > 0.9 ";
1017  cutFlowVector_str[66] = "bC2x_med: dPhi(j2,ptmiss) > 0.9 " ;
1018  cutFlowVector_str[67] = "bC2x_med: HTmissSig > 10";
1019  cutFlowVector_str[68] = "bC2x_med: mT > 120 GeV";
1020  cutFlowVector_str[69] = "bC2x_med: amT2 > 300 GeV ";
1021  cutFlowVector_str[70] = "bC2x_med: mWrecl > 50 GeV ";
1022  cutFlowVector_str[71] = "bCbv: jet0 pT > 120 GeV";
1023  cutFlowVector_str[72] = "bCbv: jet1 pT > 80 GeV ";
1024  cutFlowVector_str[73] = "bCbv: ==0 b jets ";
1025  cutFlowVector_str[74] = "bCbv: lepton pt > 60 GeV ";
1026  cutFlowVector_str[75] = "bCbv: dPhi(j1,ptmiss) > 2.0 ";
1027  cutFlowVector_str[76] = "bCbv: dPhi(j2,ptmiss) > 0.8 ";
1028  cutFlowVector_str[77] = "bCbv: Met > 360 GeV ";
1029  cutFlowVector_str[78] = "bCbv: HtmissSig > 16";
1030  cutFlowVector_str[79] = "bCbv: mT > 200 GeV";
1031  cutFlowVector_str[80] = "bCbv: mWrecl in [70,100] GeV";
1032  cutFlowVector_str[81] = "bCbv: dPhi(l,ptmiss) > 1.2";
1033  cutFlowVector_str[82] = "bCsoft_diag: Soft lepton preselection";
1034  cutFlowVector_str[83] = "bCsoft_diag: Met > 300 GeV";
1035  cutFlowVector_str[84] = "bCsoft_diag: jet0 pT > 400 GeV ";
1036  cutFlowVector_str[85] = "bCsoft_diag: mT < 160 GeV ";
1037  cutFlowVector_str[86] = "bCsoft_diag: leading jet not b-tagged ";
1038  cutFlowVector_str[87] = "bCsoft_diag: mT < 50 GeV ";
1039  cutFlowVector_str[88] = "bCsoft_diag: min(dPhi(ptmiss, b-jet)) < 1.5 ";
1040  cutFlowVector_str[89] = "bCsoft_diag: pTl/Met < 0.05 ";
1041  cutFlowVector_str[90] = "bCsoft_diag: top veto (or mtop_recl < 150 GeV) ";
1042  cutFlowVector_str[91] = "bCsoft_diag: pTl/Met < 0.02";
1043  cutFlowVector_str[92] = "bCsoft_med: Soft lepton preselection ";
1044  cutFlowVector_str[93] = "bCsoft_med: >=3 jets";
1045  cutFlowVector_str[94] = "bCsoft_med: pTW > 400 GeV";
1046  cutFlowVector_str[95] = "bCsoft_med: jet0 pT > 120 GeV";
1047  cutFlowVector_str[96] = "bCsoft_med: jet1 pT > 60 GeV";
1048  cutFlowVector_str[97] = "bCsoft_med: jet2 pT > 40 GeV";
1049  cutFlowVector_str[98] = "bCsoft_med: mT < 160 GeV";
1050  cutFlowVector_str[99] = "bCsoft_med: amT2 > 200 GeV";
1051  cutFlowVector_str[100] = "bCsoft_med: >=2 b jet";
1052  cutFlowVector_str[101] = "bCsoft_med: bjet0 pT > 120 GeV ";
1053  cutFlowVector_str[102] = "bCsoft_med: bjet1 pT > 60 GeV";
1054  cutFlowVector_str[103] = "bCsoft_med: min(dPhi(ptmiss, b-jet)) > 0.8 ";
1055  cutFlowVector_str[104] = "bCsoft_med: pTl/Met < 0.1";
1056  cutFlowVector_str[105] = "bCsoft_med: pTl/Met < 0.03";
1057  cutFlowVector_str[106] = "bCsoft_high: XE trigger, >=2 jets, Met > 230 GeV";
1058  cutFlowVector_str[107] = "bCsoft_high: dPhi(j1,ptmiss) > 0.4";
1059  cutFlowVector_str[108] = "bCsoft_high: dPhi(j2,ptmiss) > 0.4";
1060  cutFlowVector_str[109] = "bCsoft_high: jet0 pt > 100 GeV";
1061  cutFlowVector_str[110] = "bCsoft_high: jet1 pt > 100 GeV";
1062  cutFlowVector_str[111] = "bCsoft_high: mT < 160 GeV";
1063  cutFlowVector_str[112] = "bCsoft_high: pTW > 500 GeV";
1064  cutFlowVector_str[113] = "bCsoft_high: dRbb > 0.8";
1065  cutFlowVector_str[114] = "bCsoft_high: min(dPhi(ptmiss, b-jet))";
1066  cutFlowVector_str[115] = "bCsoft_high: >=2 b jet";
1067  cutFlowVector_str[116] = "bCsoft_high: bjet0 pT > 100 GeV";
1068  cutFlowVector_str[117] = "bCsoft_high: bjet1 pT > 100 GeV";
1069  cutFlowVector_str[118] = "bCsoft_high: amT2 > 300 GeV";
1070 
1071  for(int j=0;j<NCUTS;j++)
1072  {
1073  if(
1074  (j==0) ||
1075 
1076  (j==1 ) ||
1077 
1078  (j==2 && baselineLeptons.size()>0) ||
1079 
1080  (j==3 && baselineLeptons.size()>0 && signalLeptons.size()>0) ||
1081 
1082  (j==4 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1) ||
1083 
1084  (j==5 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1) ||
1085 
1086  (j==6 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230.) ||
1087 
1088  (j==7 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4) ||
1089 
1090  (j==8 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4) ||
1091 
1092  (j==9 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80) ||
1093 
1094  //tN_med cutflow
1095 
1096  (j==10 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 60 ) ||
1097 
1098  (j==11 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 60 && signalJets[1]->pT() > 60) ||
1099 
1100  (j==12 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 60 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40) ||
1101 
1102  (j==13 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 60 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40 && signalJets[3]->pT() > 40) ||
1103 
1104  (j==14 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 60 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40 && signalJets[3]->pT() > 40 && Met > 250.) ||
1105 
1106  (j==15 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 60 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40 && signalJets[3]->pT() > 40 && Met > 250. && MetPerp > 230) ||
1107 
1108  (j==16 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 60 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40 && signalJets[3]->pT() > 40 && Met > 250. && MetPerp > 230 && HtSigMiss > 14) ||
1109 
1110  (j==17 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 60 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40 && signalJets[3]->pT() > 40 && Met > 250. && MetPerp > 230 && HtSigMiss > 14 && mT > 160) ||
1111 
1112  (j==18 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 60 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40 && signalJets[3]->pT() > 40 && Met > 250. && MetPerp > 230 && HtSigMiss > 14 && mT > 160 && amT2 > 175) ||
1113 
1114  (j==19 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 60 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40 && signalJets[3]->pT() > 40 && Met > 250. && MetPerp > 230 && HtSigMiss > 14 && mT > 160 && amT2 > 175 && nBJets >=1) ||
1115 
1116  (j==20 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 60 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40 && signalJets[3]->pT() > 40 && Met > 250. && MetPerp > 230 && HtSigMiss > 14 && mT > 160 && amT2 > 175 && nBJets >=1 && dRbl < 2.0) ||
1117 
1118  (j==21 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 60 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40 && signalJets[3]->pT() > 40 && Met > 250. && MetPerp > 230 && HtSigMiss > 14 && mT > 160 && amT2 > 175 && nBJets >=1 && dRbl < 2.0 && topReclM > 150) ||
1119 
1120  // tN_high cutflow
1121 
1122  (j==22 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 100) ||
1123 
1124  (j==23 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 80) ||
1125 
1126  (j==24 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 80 && signalJets[2]->pT() > 50) ||
1127 
1128  (j==25 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 80 && signalJets[2]->pT() > 50 && signalJets[3]->pT() > 30) ||
1129 
1130  (j==26 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 80 && signalJets[2]->pT() > 50 && signalJets[3]->pT() > 30 && Met > 550.) ||
1131 
1132  (j==27 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 80 && signalJets[2]->pT() > 50 && signalJets[3]->pT() > 30 && Met > 550. && HtSigMiss > 27) ||
1133 
1134  (j==28 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 80 && signalJets[2]->pT() > 50 && signalJets[3]->pT() > 30 && Met > 550. && HtSigMiss > 27 && mT > 160) ||
1135 
1136  (j==29 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 80 && signalJets[2]->pT() > 50 && signalJets[3]->pT() > 30 && Met > 550. && HtSigMiss > 27 && mT > 160 && amT2 > 175.) ||
1137 
1138  (j==30 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 80 && signalJets[2]->pT() > 50 && signalJets[3]->pT() > 30 && Met > 550. && HtSigMiss > 27 && mT > 160 && amT2 > 175. && nBJets >=1) ||
1139 
1140  (j==31 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 80 && signalJets[2]->pT() > 50 && signalJets[3]->pT() > 30 && Met > 550. && HtSigMiss > 27 && mT > 160 && amT2 > 175. && nBJets >=1 && dRbl < 2.0) ||
1141 
1142  (j==32 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 80 && signalJets[2]->pT() > 50 && signalJets[3]->pT() > 30 && Met > 550. && HtSigMiss > 27 && mT > 160 && amT2 > 175. && nBJets >=1 && dRbl < 2.0 && topReclM > 130.) ||
1143 
1144  // bWN cutflow
1145 
1146  (j==33 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 50.) ||
1147 
1148  (j==34 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 50. && Met > 300) ||
1149 
1150  (j==35 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 50. && Met > 300 && mT > 130) ||
1151 
1152  (j==36 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 50. && Met > 300 && mT > 130 && amT2 < 110) ||
1153 
1154  (j==37 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 50. && Met > 300 && mT > 130 && amT2 < 110 && nBJets>=1) ||
1155 
1156  (j==38 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 50. && Met > 300 && mT > 130 && amT2 < 110 && nBJets>=1 && dPhiMetLep < 2.5 ) ||
1157 
1158  // bffN cutflow
1159 
1160  (j==39 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep) ||
1161 
1162  (j==40 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && Met > 300.) ||
1163 
1164  (j==41 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && Met > 300. && nJets > 0 && signalJets[0]->pT() > 400.) ||
1165 
1166  (j==42 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && Met > 300. && nJets > 0 && signalJets[0]->pT() > 400. && mT < 160.) ||
1167 
1168  (j==43 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && Met > 300. && nJets > 0 && signalJets[0]->pT() > 400. && mT < 160. && !analysisBtags.at(signalJets[0])) ||
1169 
1170  (j==44 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && Met > 300. && nJets > 0 && signalJets[0]->pT() > 400. && mT < 160. && !analysisBtags.at(signalJets[0]) && minDPhiMetBJet < 1.5) ||
1171 
1172  (j==45 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && Met > 300. && nJets > 0 && signalJets[0]->pT() > 400. && mT < 160. && !analysisBtags.at(signalJets[0]) && minDPhiMetBJet < 1.5 && pTLepOverMet < 0.05) ||
1173 
1174  (j==46 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && Met > 300. && nJets > 0 && signalJets[0]->pT() > 400. && mT < 160. && !analysisBtags.at(signalJets[0]) && minDPhiMetBJet < 1.5 && pTLepOverMet < 0.05 && topReclM < 150) ||
1175 
1176  (j==47 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && Met > 300. && nJets > 0 && signalJets[0]->pT() > 400. && mT < 160. && !analysisBtags.at(signalJets[0]) && minDPhiMetBJet < 1.5 && pTLepOverMet < 0.05 && topReclM < 150 && pTLepOverMet < 0.02) ||
1177 
1178  // bC2x_diag cutflow
1179 
1180  (j==48 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 75) ||
1181 
1182  (j==49 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 75 && signalJets[1]->pT() > 75) ||
1183 
1184  (j==50 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 75 && signalJets[1]->pT() > 75 && signalJets[2]->pT() > 75) ||
1185 
1186  (j==51 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 75 && signalJets[1]->pT() > 75 && signalJets[2]->pT() > 75 && signalJets[3]->pT() > 30) ||
1187 
1188  (j==52 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 75 && signalJets[1]->pT() > 75 && signalJets[2]->pT() > 75 && signalJets[3]->pT() > 30 && nBJets>=2) ||
1189 
1190  (j==53 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 75 && signalJets[1]->pT() > 75 && signalJets[2]->pT() > 75 && signalJets[3]->pT() > 30 && nBJets>=2 && signalBJets[1]->pT() > 30) ||
1191 
1192  (j==54 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 75 && signalJets[1]->pT() > 75 && signalJets[2]->pT() > 75 && signalJets[3]->pT() > 30 && nBJets>=2 && signalBJets[1]->pT() > 30 && absDPhiJMet0 > 0.7) ||
1193 
1194  (j==55 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 75 && signalJets[1]->pT() > 75 && signalJets[2]->pT() > 75 && signalJets[3]->pT() > 30 && nBJets>=2 && signalBJets[1]->pT() > 30 && absDPhiJMet0 > 0.7 && absDPhiJMet1 > 0.7) ||
1195 
1196  (j==56 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 75 && signalJets[1]->pT() > 75 && signalJets[2]->pT() > 75 && signalJets[3]->pT() > 30 && nBJets>=2 && signalBJets[1]->pT() > 30 && absDPhiJMet0 > 0.7 && absDPhiJMet1 > 0.7 && HtSigMiss > 13) ||
1197 
1198  (j==57 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 75 && signalJets[1]->pT() > 75 && signalJets[2]->pT() > 75 && signalJets[3]->pT() > 30 && nBJets>=2 && signalBJets[1]->pT() > 30 && absDPhiJMet0 > 0.7 && absDPhiJMet1 > 0.7 && HtSigMiss > 13 && mT > 180) ||
1199 
1200  (j==58 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 75 && signalJets[1]->pT() > 75 && signalJets[2]->pT() > 75 && signalJets[3]->pT() > 30 && nBJets>=2 && signalBJets[1]->pT() > 30 && absDPhiJMet0 > 0.7 && absDPhiJMet1 > 0.7 && HtSigMiss > 13 && mT > 180 && amT2 > 175) ||
1201 
1202  (j==59 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 75 && signalJets[1]->pT() > 75 && signalJets[2]->pT() > 75 && signalJets[3]->pT() > 30 && nBJets>=2 && signalBJets[1]->pT() > 30 && absDPhiJMet0 > 0.7 && absDPhiJMet1 > 0.7 && HtSigMiss > 13 && mT > 180 && amT2 > 175 && WReclM > 50) ||
1203 
1204  // bC2x_med
1205 
1206  (j==60 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 200) ||
1207 
1208  (j==61 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 200 && signalJets[1]->pT() > 140) ||
1209 
1210  (j==62 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 200 && signalJets[1]->pT() > 140 && nBJets >=2) ||
1211 
1212  (j==63 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 200 && signalJets[1]->pT() > 140 && nBJets >=2 && signalBJets[0]->pT() > 140) ||
1213 
1214  (j==64 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 200 && signalJets[1]->pT() > 140 && nBJets >=2 && signalBJets[0]->pT() > 140 && signalBJets[1]->pT() > 140) ||
1215 
1216  (j==65 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 200 && signalJets[1]->pT() > 140 && nBJets >=2 && signalBJets[0]->pT() > 140 && signalBJets[1]->pT() > 140 && absDPhiJMet0 > 0.9) ||
1217 
1218  (j==66 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 200 && signalJets[1]->pT() > 140 && nBJets >=2 && signalBJets[0]->pT() > 140 && signalBJets[1]->pT() > 140 && absDPhiJMet0 > 0.9 && absDPhiJMet1 > 0.9) ||
1219 
1220  (j==67 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 200 && signalJets[1]->pT() > 140 && nBJets >=2 && signalBJets[0]->pT() > 140 && signalBJets[1]->pT() > 140 && absDPhiJMet0 > 0.9 && absDPhiJMet1 > 0.9 && HtSigMiss > 10) ||
1221 
1222  (j==68 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 200 && signalJets[1]->pT() > 140 && nBJets >=2 && signalBJets[0]->pT() > 140 && signalBJets[1]->pT() > 140 && absDPhiJMet0 > 0.9 && absDPhiJMet1 > 0.9 && HtSigMiss > 10 && mT > 120) ||
1223 
1224  (j==69 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 200 && signalJets[1]->pT() > 140 && nBJets >=2 && signalBJets[0]->pT() > 140 && signalBJets[1]->pT() > 140 && absDPhiJMet0 > 0.9 && absDPhiJMet1 > 0.9 && HtSigMiss > 10 && mT > 120 && amT2 > 300) ||
1225 
1226  (j==70 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=4 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && mT2Tau > 80 && signalJets[0]->pT() > 200 && signalJets[1]->pT() > 140 && nBJets >=2 && signalBJets[0]->pT() > 140 && signalBJets[1]->pT() > 140 && absDPhiJMet0 > 0.9 && absDPhiJMet1 > 0.9 && HtSigMiss > 10 && mT > 120 && amT2 > 300 && WReclM > 50) ||
1227 
1228  // bCbv cutflow
1229 
1230  (j==71 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && signalJets[0]->pT() > 120) ||
1231 
1232  (j==72 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 80) ||
1233 
1234  (j==73 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 80 && nBJets==0) ||
1235 
1236  (j==74 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 80 && nBJets==0 && baselineLeptons[0]->pT() > 60) ||
1237 
1238  (j==75 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 80 && nBJets==0 && baselineLeptons[0]->pT() > 60 && absDPhiJMet0 > 2.0 ) ||
1239 
1240  (j==76 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 80 && nBJets==0 && baselineLeptons[0]->pT() > 60 && absDPhiJMet0 > 2.0 && absDPhiJMet1 > 0.8) ||
1241 
1242  (j==77 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 80 && nBJets==0 && baselineLeptons[0]->pT() > 60 && absDPhiJMet0 > 2.0 && absDPhiJMet1 > 0.8 && Met > 360) ||
1243 
1244  (j==78 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 80 && nBJets==0 && baselineLeptons[0]->pT() > 60 && absDPhiJMet0 > 2.0 && absDPhiJMet1 > 0.8 && Met > 360 && HtSigMiss > 16) ||
1245 
1246  (j==79 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 80 && nBJets==0 && baselineLeptons[0]->pT() > 60 && absDPhiJMet0 > 2.0 && absDPhiJMet1 > 0.8 && Met > 360 && HtSigMiss > 16 && mT > 200) ||
1247 
1248  (j==80 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 80 && nBJets==0 && baselineLeptons[0]->pT() > 60 && absDPhiJMet0 > 2.0 && absDPhiJMet1 > 0.8 && Met > 360 && HtSigMiss > 16 && mT > 200 && WReclM >= 70 && WReclM <= 100) ||
1249 
1250  (j==81 && baselineLeptons.size()>0 && signalLeptons.size()>0 && signalLeptons.size()==1 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 80 && nBJets==0 && baselineLeptons[0]->pT() > 60 && absDPhiJMet0 > 2.0 && absDPhiJMet1 > 0.8 && Met > 360 && HtSigMiss > 16 && mT > 200 && WReclM >= 70 && WReclM <= 100 && dPhiMetLep > 1.2 ) ||
1251 
1252  // bCsoft_diag
1253 
1254  (j==82 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep) ||
1255 
1256  (j==83 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && Met > 300.) ||
1257 
1258  (j==84 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && Met > 300. && signalJets[0]->pT() > 400) ||
1259 
1260  (j==85 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && Met > 300. && signalJets[0]->pT() > 400 && mT < 160) ||
1261 
1262  (j==86 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && Met > 300. && signalJets[0]->pT() > 400 && mT < 160 && !analysisBtags.at(signalJets[0])) ||
1263 
1264  (j==87 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && Met > 300. && signalJets[0]->pT() > 400 && mT < 160 && !analysisBtags.at(signalJets[0]) && mT < 50) ||
1265 
1266  (j==88 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && Met > 300. && signalJets[0]->pT() > 400 && mT < 160 && !analysisBtags.at(signalJets[0]) && mT < 50 && minDPhiMetBJet < 1.5) ||
1267 
1268  (j==89 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && Met > 300. && signalJets[0]->pT() > 400 && mT < 160 && !analysisBtags.at(signalJets[0]) && mT < 50 && minDPhiMetBJet < 1.5 && pTLepOverMet < 0.05) ||
1269 
1270  (j==90 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && Met > 300. && signalJets[0]->pT() > 400 && mT < 160 && !analysisBtags.at(signalJets[0]) && mT < 50 && minDPhiMetBJet < 1.5 && pTLepOverMet < 0.05 && topReclM < 150 ) ||
1271 
1272  (j==91 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && Met > 300. && signalJets[0]->pT() > 400 && mT < 160 && !analysisBtags.at(signalJets[0]) && mT < 50 && minDPhiMetBJet < 1.5 && pTLepOverMet < 0.05 && topReclM < 150 && pTLepOverMet < 0.02) ||
1273 
1274  // bCsoft_med cutflow
1275 
1276  (j==92 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep) ||
1277 
1278  (j==93 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && nJets >=3) ||
1279 
1280  (j==94 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && nJets >=3 && Wpt > 400) ||
1281 
1282  (j==95 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && nJets >=3 && Wpt > 400 && signalJets[0]->pT() > 120) ||
1283 
1284  (j==96 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && nJets >=3 && Wpt > 400 && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 60) ||
1285 
1286  (j==97 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && nJets >=3 && Wpt > 400 && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40) ||
1287  (j==98 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && nJets >=3 && Wpt > 400 && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40 && mT < 160) ||
1288 
1289  (j==99 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && nJets >=3 && Wpt > 400 && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40 && mT < 160 && amT2 > 200) ||
1290 
1291  (j==100 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && nJets >=3 && Wpt > 400 && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40 && mT < 160 && amT2 > 200 && nBJets > 1) ||
1292 
1293  (j==101 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && nJets >=3 && Wpt > 400 && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40 && mT < 160 && amT2 > 200 && nBJets > 1 && signalBJets[0]->pT() > 120) ||
1294 
1295  (j==102 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && nJets >=3 && Wpt > 400 && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40 && mT < 160 && amT2 > 200 && nBJets > 1 && signalBJets[0]->pT() > 120 && signalBJets[1]->pT() > 60) ||
1296 
1297  (j==103 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && nJets >=3 && Wpt > 400 && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40 && mT < 160 && amT2 > 200 && nBJets > 1 && signalBJets[0]->pT() > 120 && signalBJets[1]->pT() > 60 && minDPhiMetBJet > 0.8) ||
1298 
1299  (j==104 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && nJets >=3 && Wpt > 400 && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40 && mT < 160 && amT2 > 200 && nBJets > 1 && signalBJets[0]->pT() > 120 && signalBJets[1]->pT() > 60 && minDPhiMetBJet > 0.8 && pTLepOverMet < 0.1) ||
1300 
1301  (j==105 && baselineLeptons.size()==1 && nJets > 1 && nBJets > 0 && preselSoftLep && nJets >=3 && Wpt > 400 && signalJets[0]->pT() > 120 && signalJets[1]->pT() > 60 && signalJets[2]->pT() > 40 && mT < 160 && amT2 > 200 && nBJets > 1 && signalBJets[0]->pT() > 120 && signalBJets[1]->pT() > 60 && minDPhiMetBJet > 0.8 && pTLepOverMet < 0.1 && pTLepOverMet < 0.03) ||
1302 
1303  // bCsoft_high cutflow
1304 
1305  (j==106 && baselineLeptons.size()>0 && baselineLeptons.size()==1 && nJets >=2 && Met > 230.) ||
1306 
1307  (j==107 && baselineLeptons.size()>0 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && absDPhiJMet0 > 0.4) ||
1308 
1309  (j==108 && baselineLeptons.size()>0 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4) ||
1310  (j==109 && baselineLeptons.size()>0 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && signalJets[0]->pT() > 100) ||
1311 
1312  (j==110 && baselineLeptons.size()>0 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 100) ||
1313 
1314  (j==111 && baselineLeptons.size()>0 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 100 && mT < 160) ||
1315 
1316  (j==112 && baselineLeptons.size()>0 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 100 && mT < 160 && Wpt > 500) ||
1317 
1318  (j==113 && baselineLeptons.size()>0 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 100 && mT < 160 && Wpt > 500 && dRbb > 0.8) ||
1319 
1320  (j==114 && baselineLeptons.size()>0 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 100 && mT < 160 && Wpt > 500 && dRbb > 0.8 && minDPhiMetBJet > 0.4) ||
1321 
1322  (j==115 && baselineLeptons.size()>0 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 100 && mT < 160 && Wpt > 500 && dRbb > 0.8 && minDPhiMetBJet > 0.4 && nBJets>=2 ) ||
1323 
1324  (j==116 && baselineLeptons.size()>0 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 100 && mT < 160 && Wpt > 500 && dRbb > 0.8 && minDPhiMetBJet > 0.4 && nBJets>=2 && signalBJets[0]->pT() > 100) ||
1325 
1326  (j==117 && baselineLeptons.size()>0 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 100 && mT < 160 && Wpt > 500 && dRbb > 0.8 && minDPhiMetBJet > 0.4 && nBJets>=2 && signalBJets[0]->pT() > 100 && signalBJets[1]->pT() > 100) ||
1327 
1328  (j==118 && baselineLeptons.size()>0 && baselineLeptons.size()==1 && nJets >=2 && Met > 230. && absDPhiJMet0 > 0.4 && absDPhiJMet1 > 0.4 && signalJets[0]->pT() > 100 && signalJets[1]->pT() > 100 && mT < 160 && Wpt > 500 && dRbb > 0.8 && minDPhiMetBJet > 0.4 && nBJets>=2 && signalBJets[0]->pT() > 100 && signalBJets[1]->pT() > 100 && amT2 > 300)
1329  )
1330  {
1331  cutFlowVector[j]++;
1332  }
1333 
1334  }
1335 
1336  if(is_tN_med) _counters.at("tN_med").add_event(event);
1337  if(is_tN_high) _counters.at("tN_high").add_event(event);
1338  if(is_bWN) _counters.at("bWN").add_event(event);
1339  if(is_bC2x_diag) _counters.at("bC2x_diag").add_event(event);
1340  if(is_bC2x_med) _counters.at("bC2x_med").add_event(event);
1341  if(is_bCbv) _counters.at("bCbv").add_event(event);
1342  if(is_DM_low) _counters.at("DM_low_loose").add_event(event);
1343  if(is_DM_low) _counters.at("DM_low").add_event(event);
1344  if(is_DM_high) _counters.at("DM_high").add_event(event);
1345 
1346  if(is_bffN) _counters.at("bffN").add_event(event);
1347  if(is_bCsoft_diag) _counters.at("bCsoft_diag").add_event(event);
1348  if(is_bCsoft_med) _counters.at("bCsoft_med").add_event(event);
1349  if(is_bCsoft_high) _counters.at("bCsoft_high").add_event(event);
1350 
1351  return;
1352 
1353  }
1354 
1356  void combine(const Analysis* other)
1357  {
1358  const Analysis_ATLAS_13TeV_1LEPStop_36invfb* specificOther
1359  = dynamic_cast<const Analysis_ATLAS_13TeV_1LEPStop_36invfb*>(other);
1360 
1361  for (auto& pair : _counters) { pair.second += specificOther->_counters.at(pair.first); }
1362 
1363  if (NCUTS != specificOther->NCUTS) NCUTS = specificOther->NCUTS;
1364  for (int j=0; j<NCUTS; j++)
1365  {
1366  cutFlowVector[j] += specificOther->cutFlowVector[j];
1367  cutFlowVector_str[j] = specificOther->cutFlowVector_str[j];
1368  }
1369  }
1370 
1371 
1372  void collect_results()
1373  {
1374 
1375  // // For debugging:
1376  // double scale_by=1.;
1377  // cout << "------------------------------------------------------------------------------------------------------------------------------ "<<endl;
1378  // cout << "CUT FLOW: ATLAS 13 TeV 1 lep stop paper "<<endl;
1379  // cout << "------------------------------------------------------------------------------------------------------------------------------"<<endl;
1380  // cout<< right << setw(40) << "CUT" << setw(20) << "RAW" << setw(20) << "SCALED"
1381  // << setw(20) << "%" << setw(20) << "clean adj RAW"<< setw(20) << "clean adj %" << endl;
1382  // for (int j=0; j<NCUTS; j++)
1383  // {
1384  // cout << right << setw(40) << cutFlowVector_str[j].c_str() << setw(20)
1385  // << cutFlowVector[j] << setw(20) << cutFlowVector[j]*scale_by << setw(20)
1386  // << 100.*cutFlowVector[j]/cutFlowVector[0] << "%" << setw(20)
1387  // << cutFlowVector[j]*scale_by << setw(20) << 100.*cutFlowVector[j]/cutFlowVector[0]<< "%" << endl;
1388  // }
1389  // cout << "------------------------------------------------------------------------------------------------------------------------------ "<<endl;
1390 
1392 
1393  add_result(SignalRegionData(_counters.at("tN_med"), 50., { 36.3, 6.6}));
1394  add_result(SignalRegionData(_counters.at("tN_high"), 8., { 3.8, 1.0}));
1395  add_result(SignalRegionData(_counters.at("bWN"), 68., { 71, 16}));
1396  add_result(SignalRegionData(_counters.at("bC2x_diag"), 22., { 21.3, 5.0}));
1397  add_result(SignalRegionData(_counters.at("bC2x_med"), 4., { 5.8, 1.6}));
1398  add_result(SignalRegionData(_counters.at("bCbv"), 25., { 25.1, 3.8}));
1399  add_result(SignalRegionData(_counters.at("DM_low_loose"), 65., { 48.3, 8.2}));
1400  add_result(SignalRegionData(_counters.at("DM_low"), 13., { 13.8, 3.6}));
1401  add_result(SignalRegionData(_counters.at("DM_high"), 5., { 7.4, 2.1}));
1402  add_result(SignalRegionData(_counters.at("bffN"), 70., { 60.5, 6.1}));
1403  add_result(SignalRegionData(_counters.at("bCsoft_diag"), 33., { 24.7, 3.1}));
1404  add_result(SignalRegionData(_counters.at("bCsoft_med"), 19., { 13.7, 2.1}));
1405  add_result(SignalRegionData(_counters.at("bCsoft_high"), 2., { 1.8, 0.3}));
1406 
1407  return;
1408  }
1409 
1410 
1411  protected:
1412  void analysis_specific_reset()
1413  {
1414  for (auto& pair : _counters) { pair.second.reset(); }
1415 
1416  std::fill(cutFlowVector.begin(), cutFlowVector.end(), 0);
1417  }
1418 
1419  };
1420 
1421  DEFINE_ANALYSIS_FACTORY(ATLAS_13TeV_1LEPStop_36invfb)
1422 
1423  }
1424 }
1425 
1426 #endif
bool sortByMass_1l_sharedptr(std::shared_ptr< HEPUtils::Jet > jet1, std::shared_ptr< HEPUtils::Jet > jet2)
bool has_tag(const HEPUtils::BinnedFn2D< double > &effmap, double eta, double pt)
Randomly get a tag result (can be anything) from a 2D |eta|-pT efficiency map.
Definition: Utils.hpp:207
#define DEFINE_ANALYSIS_FACTORY(ANAME)
For analysis factory function definition.
Definition: Analysis.hpp:124
bool sortByMass_1l(const HEPUtils::Jet *jet1, const HEPUtils::Jet *jet2)
void run(bool, bool, bool, int, double, double, int, int, int, int, int, double, char[], int, int[], bool, bool, bool, bool, double, int, double(*)(double *, int, int, void *), void(*)(int, int, int, double *, double *, double *, double, double, double, void *), void *)
const double mt
Definition: topness.h:39
STL namespace.
HEPUtils::P4 mk_p4(const Vec4T &p)
Definition: Py8Utils.hpp:45
START_MODEL b
Definition: demo.hpp:270
const double mW
Definition: topness.h:40
double calcMT_1l(HEPUtils::P4 jetMom, HEPUtils::P4 metMom)
void set_momenta(double *pa0, double *pb0, double *pmiss0)
Definition: mt2_bisect.cpp:62
A simple class for counting events of type HEPUtils::Event.
A class for collider analyses within ColliderBit.
Definition: Analysis.hpp:41
double get_mt2()
Definition: mt2_bisect.cpp:50
A threadsafe interface to the STL random number generators.
bool sortByPT_1l_sharedptr(std::shared_ptr< HEPUtils::Jet > jet1, std::shared_ptr< HEPUtils::Jet > jet2)
START_MODEL dNur_CMB r
A simple container for the result of one signal region from one analysis.
void applyTightIDElectronSelection(std::vector< const HEPUtils::Particle *> &electrons)
Efficiency function for Tight ID electrons.
void applyMuonEff(std::vector< const HEPUtils::Particle *> &muons)
Randomly filter the supplied particle list by parameterised muon efficiency.
DS5_MSPCTM DS_INTDOF int
void set_mn(double mn)
Definition: mt2_bisect.cpp:130
double pow(const double &a)
Outputs a^i.
def combine(region_file, pip_file)
Definition: colouring.py:169
void applyElectronEff(std::vector< const HEPUtils::Particle *> &electrons)
Randomly filter the supplied particle list by parameterised electron efficiency.
TODO: see if we can use this one:
Definition: Analysis.hpp:33
void applyTauEfficiencyR1(std::vector< const HEPUtils::Particle *> &taus)
Randomly filter the supplied particle list by parameterised Run 1 tau efficiency. ...
Class for ColliderBit analyses.
double E(const double x, const double y)
Functions that do super fast ATLAS detector simulation based on four-vector smearing.
bool sortByPT_1l(const HEPUtils::Jet *jet1, const HEPUtils::Jet *jet2)