gambit is hosted by Hepforge, IPPP Durham
GAMBIT  v1.5.0-2191-ga4742ac
a Global And Modular Bsm Inference Tool
Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb.cpp
Go to the documentation of this file.
1 
2 #include "gambit/cmake/cmake_variables.hpp"
3 #ifndef EXCLUDE_ROOT
4 #ifndef EXCLUDE_RESTFRAMES
5 
6 #include <vector>
7 #include <cmath>
8 #include <memory>
9 #include <algorithm>
10 #include <iomanip>
11 
15 
16 #include "RestFrames/RestFrames.hh"
17 #include "TLorentzVector.h"
18 
19 using namespace std;
20 
21 /* The ATLAS 13 TeV 3 lepton low mass recursive jigsaw search
22 
23  Based on code kindly supplied by Abhishek Sharma
24 
25  Note that use of ROOT is compulsory for the RestFrames package
26 
27  Based on: https://arxiv.org/pdf/1806.02293.pdf
28 
29  Code adapted by Martin White
30 
31  KNOWN ISSUES
32 
33  1) Need to check overlap removal step when the paper comes out. For now, have assumed it is the same as the stop analysis.
34 
35 */
36 
37 namespace Gambit {
38  namespace ColliderBit {
39 
40  bool sortByPT_RJ3L(const HEPUtils::Jet* jet1, const HEPUtils::Jet* jet2) { return (jet1->pT() > jet2->pT()); }
41  bool sortLepByPT_RJ3L(const HEPUtils::Particle* lep1, const HEPUtils::Particle* lep2) { return (lep1->pT() > lep2->pT());}
42  //bool sortByMass(const HEPUtils::Jet* jet1, const HEPUtils::Jet* jet2) { return (jet1->mass() > jet2->mass()); }
43 
44  bool SortLeptons(const pair<TLorentzVector,int> lv1, const pair<TLorentzVector,int> lv2)
45  //bool VariableConstruction::SortLeptons(const lep lv1, const lep lv2)
46  {
47  return lv1.first.Pt() > lv2.first.Pt();
48  }
49 
50  bool SortJets(const TLorentzVector jv1, const TLorentzVector jv2)
51  {
52  return jv1.Pt() > jv2.Pt();
53  }
54 
55 
56  class Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb : public Analysis {
57 
58  protected:
59  // Numbers passing cuts
60  std::map<string, EventCounter> _counters = {
61  {"2L2JHIGH", EventCounter("2L2JHIGH")},
62  {"2L2JINT", EventCounter("2L2JINT")},
63  {"2L2JLOW", EventCounter("2L2JLOW")},
64  {"2L2JCOMP", EventCounter("2L2JCOMP")},
65  {"3LHIGH", EventCounter("3LHIGH")},
66  {"3LINT", EventCounter("3LINT")},
67  {"3LLOW", EventCounter("3LLOW")},
68  {"3LCOMP", EventCounter("3LCOMP")},
69  };
70 
71  private:
72 
73  vector<int> cutFlowVector;
74  vector<string> cutFlowVector_str;
75  int NCUTS; //=16;
76 
77  // Recursive jigsaw objects (using RestFrames)
78 
82 
83  unique_ptr<RestFrames::LabRecoFrame> LAB_2L2J;
84  unique_ptr<RestFrames::DecayRecoFrame> C1N2_2L2J;
85  unique_ptr<RestFrames::DecayRecoFrame> C1a_2L2J;
86  unique_ptr<RestFrames::DecayRecoFrame> N2b_2L2J;
87 
88  unique_ptr<RestFrames::DecayRecoFrame> Wa_2L2J;
89  unique_ptr<RestFrames::DecayRecoFrame> Zb_2L2J;
90 
91  unique_ptr<RestFrames::VisibleRecoFrame> J1_2L2J;
92  unique_ptr<RestFrames::VisibleRecoFrame> J2_2L2J;
93  unique_ptr<RestFrames::VisibleRecoFrame> L1_2L2J;
94  unique_ptr<RestFrames::VisibleRecoFrame> L2_2L2J;
95 
96  unique_ptr<RestFrames::InvisibleRecoFrame> X1a_2L2J;
97  unique_ptr<RestFrames::InvisibleRecoFrame> X1b_2L2J;
98 
99  unique_ptr<RestFrames::InvisibleGroup> INV_2L2J;
100 
101  unique_ptr<RestFrames::SetMassInvJigsaw> X1_mass_2L2J;
102  unique_ptr<RestFrames::SetRapidityInvJigsaw> X1_eta_2L2J;
103 
104  unique_ptr<RestFrames::ContraBoostInvJigsaw> X1X1_contra_2L2J;
105 
109 
110  unique_ptr<RestFrames::LabRecoFrame> LAB_3L;
111  unique_ptr<RestFrames::DecayRecoFrame> C1N2_3L;
112  unique_ptr<RestFrames::DecayRecoFrame> C1a_3L;
113  unique_ptr<RestFrames::DecayRecoFrame> N2b_3L;
114 
115  // unique_ptr<RestFrames::DecayRecoFrame> Wa_3L;
116  // unique_ptr<RestFrames::DecayRecoFrame> Zb_3L;
117 
118  unique_ptr<RestFrames::VisibleRecoFrame> L1a_3L;
119  unique_ptr<RestFrames::VisibleRecoFrame> L1b_3L;
120  unique_ptr<RestFrames::VisibleRecoFrame> L2b_3L;
121 
122  unique_ptr<RestFrames::InvisibleRecoFrame> X1a_3L;
123  unique_ptr<RestFrames::InvisibleRecoFrame> X1b_3L;
124 
125  unique_ptr<RestFrames::InvisibleGroup> INV_3L;
126 
127  unique_ptr<RestFrames::SetMassInvJigsaw> X1_mass_3L;
128  unique_ptr<RestFrames::SetRapidityInvJigsaw> X1_eta_3L;
129 
130  unique_ptr<RestFrames::ContraBoostInvJigsaw> X1X1_contra_3L;
131 
132  // combinatoric (transverse) tree
133  // for jet assignment
134  unique_ptr<RestFrames::LabRecoFrame> LAB_comb;
135  unique_ptr<RestFrames::DecayRecoFrame> CM_comb;
136  unique_ptr<RestFrames::DecayRecoFrame> S_comb;
137  unique_ptr<RestFrames::VisibleRecoFrame> ISR_comb;
138  unique_ptr<RestFrames::VisibleRecoFrame> J_comb;
139  unique_ptr<RestFrames::VisibleRecoFrame> L_comb;
140  unique_ptr<RestFrames::InvisibleRecoFrame> I_comb;
141  unique_ptr<RestFrames::InvisibleGroup> INV_comb;
142  unique_ptr<RestFrames::SetMassInvJigsaw> InvMass_comb;
143  unique_ptr<RestFrames::CombinatoricGroup> JETS_comb;
144  unique_ptr<RestFrames::MinMassesCombJigsaw> SplitJETS_comb;
145 
146  // 2L+NJ tree (Z->ll + W/Z->qq)
147  unique_ptr<RestFrames::LabRecoFrame> LAB_2LNJ;
148  unique_ptr<RestFrames::DecayRecoFrame> CM_2LNJ;
149  unique_ptr<RestFrames::DecayRecoFrame> S_2LNJ;
150  unique_ptr<RestFrames::VisibleRecoFrame> ISR_2LNJ;
151 
152  unique_ptr<RestFrames::DecayRecoFrame> Ca_2LNJ;
153  unique_ptr<RestFrames::DecayRecoFrame> Z_2LNJ;
154  unique_ptr<RestFrames::VisibleRecoFrame> L1_2LNJ;
155  unique_ptr<RestFrames::VisibleRecoFrame> L2_2LNJ;
156 
157  unique_ptr<RestFrames::DecayRecoFrame> Cb_2LNJ;
158  unique_ptr<RestFrames::SelfAssemblingRecoFrame> JSA_2LNJ;
159  unique_ptr<RestFrames::VisibleRecoFrame> J_2LNJ;
160 
161  unique_ptr<RestFrames::InvisibleRecoFrame> Ia_2LNJ;
162  unique_ptr<RestFrames::InvisibleRecoFrame> Ib_2LNJ;
163 
164  unique_ptr<RestFrames::InvisibleGroup> INV_2LNJ;
165  unique_ptr<RestFrames::SetMassInvJigsaw> InvMass_2LNJ;
166  unique_ptr<RestFrames::SetRapidityInvJigsaw> InvRapidity_2LNJ;
167  unique_ptr<RestFrames::ContraBoostInvJigsaw> SplitINV_2LNJ;
168  unique_ptr<RestFrames::CombinatoricGroup> JETS_2LNJ;
169 
170  // 2L+1L tree (Z->ll + Z/W->l)
171  unique_ptr<RestFrames::LabRecoFrame> LAB_2L1L;
172  unique_ptr<RestFrames::DecayRecoFrame> CM_2L1L;
173  unique_ptr<RestFrames::DecayRecoFrame> S_2L1L;
174  unique_ptr<RestFrames::VisibleRecoFrame> ISR_2L1L;
175 
176  unique_ptr<RestFrames::DecayRecoFrame> Ca_2L1L;
177  unique_ptr<RestFrames::DecayRecoFrame> Z_2L1L;
178  unique_ptr<RestFrames::VisibleRecoFrame> L1_2L1L;
179  unique_ptr<RestFrames::VisibleRecoFrame> L2_2L1L;
180 
181  unique_ptr<RestFrames::DecayRecoFrame> Cb_2L1L;
182  unique_ptr<RestFrames::VisibleRecoFrame> Lb_2L1L;
183 
184  unique_ptr<RestFrames::InvisibleRecoFrame> Ia_2L1L;
185  unique_ptr<RestFrames::InvisibleRecoFrame> Ib_2L1L;
186 
187  unique_ptr<RestFrames::InvisibleGroup> INV_2L1L;
188  unique_ptr<RestFrames::SetMassInvJigsaw> InvMass_2L1L;
189  unique_ptr<RestFrames::SetRapidityInvJigsaw> InvRapidity_2L1L;
190  unique_ptr<RestFrames::ContraBoostInvJigsaw> SplitINV_2L1L;
191 
192  // Debug histos
193 
194  void JetLeptonOverlapRemoval(vector<const HEPUtils::Jet*> &jetvec, vector<const HEPUtils::Particle*> &lepvec, double DeltaRMax) {
195  //Routine to do jet-lepton check
196  //Discards jets if they are within DeltaRMax of a lepton
197 
198  vector<const HEPUtils::Jet*> Survivors;
199 
200  for(unsigned int itjet = 0; itjet < jetvec.size(); itjet++) {
201  bool overlap = false;
202  HEPUtils::P4 jetmom=jetvec.at(itjet)->mom();
203  for(unsigned int itlep = 0; itlep < lepvec.size(); itlep++) {
204  HEPUtils::P4 lepmom=lepvec.at(itlep)->mom();
205  double dR;
206 
207  dR=jetmom.deltaR_eta(lepmom);
208 
209  if(fabs(dR) <= DeltaRMax) overlap=true;
210  }
211  if(overlap) continue;
212  Survivors.push_back(jetvec.at(itjet));
213  }
214  jetvec=Survivors;
215 
216  return;
217  }
218 
219  void LeptonJetOverlapRemoval(vector<const HEPUtils::Particle*> &lepvec, vector<const HEPUtils::Jet*> &jetvec, double DeltaRMax) {
220  //Routine to do lepton-jet check
221  //Discards leptons if they are within DeltaRMax of a jet
222 
223  vector<const HEPUtils::Particle*> Survivors;
224 
225  for(unsigned int itlep = 0; itlep < lepvec.size(); itlep++) {
226  bool overlap = false;
227  HEPUtils::P4 lepmom=lepvec.at(itlep)->mom();
228  for(unsigned int itjet= 0; itjet < jetvec.size(); itjet++) {
229  HEPUtils::P4 jetmom=jetvec.at(itjet)->mom();
230  double dR;
231 
232  dR=jetmom.deltaR_eta(lepmom);
233 
234  if(fabs(dR) <= DeltaRMax) overlap=true;
235  }
236  if(overlap) continue;
237  Survivors.push_back(lepvec.at(itlep));
238  }
239  lepvec=Survivors;
240 
241  return;
242  }
243 
244 
245  public:
246 
247  // Required detector sim
248  static constexpr const char* detector = "ATLAS";
249 
250  Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb() {
251 
252  set_analysis_name("ATLAS_13TeV_RJ3L_lowmass_36invfb");
253  set_luminosity(36.);
254 
255  NCUTS=70;
256 
257  for(int i=0;i<NCUTS;i++){
258  cutFlowVector.push_back(0);
259  cutFlowVector_str.push_back("");
260  }
261 
262 
263  // Recursive jigsaw stuff
264  #pragma omp critical (init_ATLAS_13TeV_RJ3L_lowmass_36invfb)
265  {
266 
267  // // DEBUG:
268  // RestFrames::SetLogPrint(RestFrames::LogDebug, true);
269  // RestFrames::SetLogPrint(RestFrames::LogVerbose, true);
270 
271  LAB_2L2J.reset(new RestFrames::LabRecoFrame("LAB_2L2J","lab2L2J"));
272  C1N2_2L2J.reset(new RestFrames::DecayRecoFrame("C1N2_2L2J","#tilde{#chi}^{ #pm}_{1} #tilde{#chi}^{ 0}_{2}"));
273  C1a_2L2J.reset(new RestFrames::DecayRecoFrame("C1a_2L2J","#tilde{#chi}^{ #pm}_{1}"));
274  N2b_2L2J.reset(new RestFrames::DecayRecoFrame("N2b_2L2J","#tilde{#chi}^{ 0}_{2}"));
275 
276  Wa_2L2J.reset(new RestFrames::DecayRecoFrame("Wa_2L2J","W_{a}"));
277  Zb_2L2J.reset(new RestFrames::DecayRecoFrame("Zb_2L2J","Z_{b}"));
278 
279  J1_2L2J.reset(new RestFrames::VisibleRecoFrame("J1_2L2J","#it{j}_{1}"));
280  J2_2L2J.reset(new RestFrames::VisibleRecoFrame("J2_2L2J","#it{j}_{2}"));
281  L1_2L2J.reset(new RestFrames::VisibleRecoFrame("L1_2L2J","#it{l}_{1}"));
282  L2_2L2J.reset(new RestFrames::VisibleRecoFrame("L2_2L2J","#it{l}_{2}"));
283 
284  X1a_2L2J.reset(new RestFrames::InvisibleRecoFrame("X1a_2L2J","#tilde{#chi}^{ 0}_{1 a}"));
285  X1b_2L2J.reset(new RestFrames::InvisibleRecoFrame("X1b_2L2J","#tilde{#chi}^{ 0}_{1 b}"));
286 
287  LAB_2L2J->SetChildFrame(*C1N2_2L2J);
288 
289  C1N2_2L2J->AddChildFrame(*C1a_2L2J);
290  C1N2_2L2J->AddChildFrame(*N2b_2L2J);
291 
292  C1a_2L2J->AddChildFrame(*Wa_2L2J);
293  C1a_2L2J->AddChildFrame(*X1a_2L2J);
294 
295  N2b_2L2J->AddChildFrame(*Zb_2L2J);
296  N2b_2L2J->AddChildFrame(*X1b_2L2J);
297 
298  Wa_2L2J->AddChildFrame(*J1_2L2J);
299  Wa_2L2J->AddChildFrame(*J2_2L2J);
300 
301  Zb_2L2J->AddChildFrame(*L1_2L2J);
302  Zb_2L2J->AddChildFrame(*L2_2L2J);
303 
304 
305  if(!LAB_2L2J->InitializeTree())
306  {
307  str errmsg;
308  errmsg = "Some problem occurred when calling LAB_2L2J->InitializeTree() from the Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb analysis class.\n";
310  }
311 
312 
314  //Setting the invisible
316  INV_2L2J.reset(new RestFrames::InvisibleGroup("INV_2L2J","#tilde{#chi}_{1}^{ 0} Jigsaws"));
317  INV_2L2J->AddFrame(*X1a_2L2J);
318  INV_2L2J->AddFrame(*X1b_2L2J);
319 
320  // Set di-LSP mass to minimum Lorentz-invariant expression
321  X1_mass_2L2J.reset(new RestFrames::SetMassInvJigsaw("X1_mass_2L2J", "Set M_{#tilde{#chi}_{1}^{ 0} #tilde{#chi}_{1}^{ 0}} to minimum"));
322  INV_2L2J->AddJigsaw(*X1_mass_2L2J);
323 
324  // Set di-LSP rapidity to that of visible particles
325  X1_eta_2L2J.reset(new RestFrames::SetRapidityInvJigsaw("X1_eta_2L2J", "#eta_{#tilde{#chi}_{1}^{ 0} #tilde{#chi}_{1}^{ 0}} = #eta_{2jet+2#it{l}}"));
326  INV_2L2J->AddJigsaw(*X1_eta_2L2J);
327  X1_eta_2L2J->AddVisibleFrames(C1N2_2L2J->GetListVisibleFrames());
328 
329 
330  X1X1_contra_2L2J.reset(new RestFrames::ContraBoostInvJigsaw("X1X1_contra_2L2J","Contraboost invariant Jigsaw"));
331  INV_2L2J->AddJigsaw(*X1X1_contra_2L2J);
332  X1X1_contra_2L2J->AddVisibleFrames(C1a_2L2J->GetListVisibleFrames(), 0);
333  X1X1_contra_2L2J->AddVisibleFrames(N2b_2L2J->GetListVisibleFrames(), 1);
334  X1X1_contra_2L2J->AddInvisibleFrame(*X1a_2L2J, 0);
335  X1X1_contra_2L2J->AddInvisibleFrame(*X1b_2L2J, 1);
336 
337  if(!LAB_2L2J->InitializeAnalysis())
338  {
339  str errmsg;
340  errmsg = "Some problem occurred when calling LAB_2L2J->InitializeAnalysis() from the Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb analysis class.\n";
342  }
343 
344 
345  LAB_3L.reset(new RestFrames::LabRecoFrame("LAB_3L","lab"));
346  C1N2_3L.reset(new RestFrames::DecayRecoFrame("C1N2_3L","#tilde{#chi}^{ #pm}_{1} #tilde{#chi}^{ 0}_{2}"));
347  C1a_3L.reset(new RestFrames::DecayRecoFrame("C1a_3L","#tilde{#chi}^{ #pm}_{1}"));
348  N2b_3L.reset(new RestFrames::DecayRecoFrame("N2b_3L","#tilde{#chi}^{ 0}_{2}"));
349 
350  L1a_3L.reset(new RestFrames::VisibleRecoFrame("L1a_3L","#it{l}_{1a}"));
351  L1b_3L.reset(new RestFrames::VisibleRecoFrame("L1b_3L","#it{l}_{1b}"));
352  L2b_3L.reset(new RestFrames::VisibleRecoFrame("L2b_3L","#it{l}_{2b}"));
353 
354  X1a_3L.reset(new RestFrames::InvisibleRecoFrame("X1a_3L","#tilde{#chi}^{ 0}_{1 a} + #nu_{a}"));
355  X1b_3L.reset(new RestFrames::InvisibleRecoFrame("X1b_3L","#tilde{#chi}^{ 0}_{1 b}"));
356 
357 
358  LAB_3L->SetChildFrame(*C1N2_3L);
359 
360  C1N2_3L->AddChildFrame(*C1a_3L);
361  C1N2_3L->AddChildFrame(*N2b_3L);
362 
363  C1a_3L->AddChildFrame(*L1a_3L);
364  C1a_3L->AddChildFrame(*X1a_3L);
365 
366  N2b_3L->AddChildFrame(*L1b_3L);
367  N2b_3L->AddChildFrame(*L2b_3L);
368  N2b_3L->AddChildFrame(*X1b_3L);
369 
370 
371  if(!LAB_3L->InitializeTree())
372  {
373  str errmsg;
374  errmsg = "Some problem occurred when calling LAB_3L->InitializeTree() from the Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb analysis class.\n";
376  }
377 
378 
379  //setting the invisible components
380  INV_3L.reset(new RestFrames::InvisibleGroup("INV_3L","Invisible system LSP mass Jigsaw"));
381  INV_3L->AddFrame(*X1a_3L);
382  INV_3L->AddFrame(*X1b_3L);
383 
384 
385  // Set di-LSP mass to minimum Lorentz-invariant expression
386  X1_mass_3L.reset(new RestFrames::SetMassInvJigsaw("X1_mass_3L", "Set M_{#tilde{#chi}_{1}^{ 0} #tilde{#chi}_{1}^{ 0}} to minimum"));
387  INV_3L->AddJigsaw(*X1_mass_3L);
388 
389  // Set di-LSP rapidity to that of visible particles and neutrino
390  X1_eta_3L.reset(new RestFrames::SetRapidityInvJigsaw("X1_eta_3L", "#eta_{#tilde{#chi}_{1}^{ 0} #tilde{#chi}_{1}^{ 0}} = #eta_{3#it{l}}"));
391  INV_3L->AddJigsaw(*X1_eta_3L);
392  X1_eta_3L->AddVisibleFrames(C1N2_3L->GetListVisibleFrames());
393 
394 
395  X1X1_contra_3L.reset(new RestFrames::ContraBoostInvJigsaw("X1X1_contra_3L","Contraboost invariant Jigsaw"));
396  INV_3L->AddJigsaw(*X1X1_contra_3L);
397  X1X1_contra_3L->AddVisibleFrames(C1a_3L->GetListVisibleFrames(),0);
398  X1X1_contra_3L->AddVisibleFrames(N2b_3L->GetListVisibleFrames(),1);
399  X1X1_contra_3L->AddInvisibleFrames(C1a_3L->GetListInvisibleFrames(),0);
400  X1X1_contra_3L->AddInvisibleFrames(N2b_3L->GetListInvisibleFrames(),1);
401 
402  if(!LAB_3L->InitializeAnalysis())
403  {
404  str errmsg;
405  errmsg = "Some problem occurred when calling LAB_3L->InitializeAnalysis() from the Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb analysis class.\n";
407  }
408 
409 
411  // RestFrames stuff
412 
413  // combinatoric (transverse) tree
414  // for jet assignment
415  LAB_comb.reset(new RestFrames::LabRecoFrame("LAB_comb","LAB"));
416  CM_comb.reset(new RestFrames::DecayRecoFrame("CM_comb","CM"));
417  S_comb.reset(new RestFrames::DecayRecoFrame("S_comb","S"));
418  ISR_comb.reset(new RestFrames::VisibleRecoFrame("ISR_comb","ISR"));
419  J_comb.reset(new RestFrames::VisibleRecoFrame("J_comb","Jets"));
420  L_comb.reset(new RestFrames::VisibleRecoFrame("L_comb","#it{l}'s"));
421  I_comb.reset(new RestFrames::InvisibleRecoFrame("I_comb","Inv"));
422 
423  LAB_comb->SetChildFrame(*CM_comb);
424  CM_comb->AddChildFrame(*ISR_comb);
425  CM_comb->AddChildFrame(*S_comb);
426  S_comb->AddChildFrame(*L_comb);
427  S_comb->AddChildFrame(*J_comb);
428  S_comb->AddChildFrame(*I_comb);
429 
430  if(!LAB_comb->InitializeTree())
431  {
432  str errmsg;
433  errmsg = "Some problem occurred when calling LAB_comb->InitializeTree() from the Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb analysis class.\n";
435  }
436 
437  // 2L+NJ tree (Z->ll + W/Z->qq)
438  LAB_2LNJ.reset(new RestFrames::LabRecoFrame("LAB_2LNJ","LAB"));
439  CM_2LNJ.reset(new RestFrames::DecayRecoFrame("CM_2LNJ","CM"));
440  S_2LNJ.reset(new RestFrames::DecayRecoFrame("S_2LNJ","S"));
441  ISR_2LNJ.reset(new RestFrames::VisibleRecoFrame("ISR_2LNJ","ISR"));
442  Ca_2LNJ.reset(new RestFrames::DecayRecoFrame("Ca_2LNJ","C_{a}"));
443  Z_2LNJ.reset(new RestFrames::DecayRecoFrame("Z_2LNJ","Z"));
444  L1_2LNJ.reset(new RestFrames::VisibleRecoFrame("L1_2LNJ","#it{l}_{1}"));
445  L2_2LNJ.reset(new RestFrames::VisibleRecoFrame("L2_2LNJ","#it{l}_{2}"));
446  Cb_2LNJ.reset(new RestFrames::DecayRecoFrame("Cb_2LNJ","C_{b}"));
447  JSA_2LNJ.reset(new RestFrames::SelfAssemblingRecoFrame("JSA_2LNJ", "J"));
448  J_2LNJ.reset(new RestFrames::VisibleRecoFrame("J_2LNJ","Jets"));
449  Ia_2LNJ.reset(new RestFrames::InvisibleRecoFrame("Ia_2LNJ","I_{a}"));
450  Ib_2LNJ.reset(new RestFrames::InvisibleRecoFrame("Ib_2LNJ","I_{b}"));
451 
452  LAB_2LNJ->SetChildFrame(*CM_2LNJ);
453  CM_2LNJ->AddChildFrame(*ISR_2LNJ);
454  CM_2LNJ->AddChildFrame(*S_2LNJ);
455  S_2LNJ->AddChildFrame(*Ca_2LNJ);
456  S_2LNJ->AddChildFrame(*Cb_2LNJ);
457  Ca_2LNJ->AddChildFrame(*Z_2LNJ);
458  Ca_2LNJ->AddChildFrame(*Ia_2LNJ);
459  Cb_2LNJ->AddChildFrame(*JSA_2LNJ);
460  Cb_2LNJ->AddChildFrame(*Ib_2LNJ);
461  Z_2LNJ->AddChildFrame(*L1_2LNJ);
462  Z_2LNJ->AddChildFrame(*L2_2LNJ);
463  JSA_2LNJ->AddChildFrame(*J_2LNJ);
464 
465  if(!LAB_2LNJ->InitializeTree())
466  {
467  str errmsg;
468  errmsg = "Some problem occurred when calling LAB_2LNJ->InitializeTree() from the Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb analysis class.\n";
470  }
471 
472 
473  // 2L+1L tree (Z->ll + Z/W->l)
474  LAB_2L1L.reset(new RestFrames::LabRecoFrame("LAB_2L1L","LAB"));
475  CM_2L1L.reset(new RestFrames::DecayRecoFrame("CM_2L1L","CM"));
476  S_2L1L.reset(new RestFrames::DecayRecoFrame("S_2L1L","S"));
477  ISR_2L1L.reset(new RestFrames::VisibleRecoFrame("ISR_2L1L","ISR"));
478  Ca_2L1L.reset(new RestFrames::DecayRecoFrame("Ca_2L1L","C_{a}"));
479  Z_2L1L.reset(new RestFrames::DecayRecoFrame("Z_2L1L","Z"));
480  L1_2L1L.reset(new RestFrames::VisibleRecoFrame("L1_2L1L","#it{l}_{1}"));
481  L2_2L1L.reset(new RestFrames::VisibleRecoFrame("L2_2L1L","#it{l}_{2}"));
482  Cb_2L1L.reset(new RestFrames::DecayRecoFrame("Cb_2L1L","C_{b}"));
483  Lb_2L1L.reset(new RestFrames::VisibleRecoFrame("Lb_2L1L","#it{l}_{b}"));
484  Ia_2L1L.reset(new RestFrames::InvisibleRecoFrame("Ia_2L1L","I_{a}"));
485  Ib_2L1L.reset(new RestFrames::InvisibleRecoFrame("Ia_2L1L","I_{b}"));
486 
487  LAB_2L1L->SetChildFrame(*CM_2L1L);
488  CM_2L1L->AddChildFrame(*ISR_2L1L);
489  CM_2L1L->AddChildFrame(*S_2L1L);
490  S_2L1L->AddChildFrame(*Ca_2L1L);
491  S_2L1L->AddChildFrame(*Cb_2L1L);
492  Ca_2L1L->AddChildFrame(*Z_2L1L);
493  Ca_2L1L->AddChildFrame(*Ia_2L1L);
494  Z_2L1L->AddChildFrame(*L1_2L1L);
495  Z_2L1L->AddChildFrame(*L2_2L1L);
496  Cb_2L1L->AddChildFrame(*Lb_2L1L);
497  Cb_2L1L->AddChildFrame(*Ib_2L1L);
498 
499  if(!LAB_2L1L->InitializeTree())
500  {
501  str errmsg;
502  errmsg = "Some problem occurred when calling LAB_2L1L->InitializeTree() from the Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb analysis class.\n";
504  }
505 
506 
508 
509  // combinatoric (transverse) tree
510  // for jet assignment
511  INV_comb.reset(new RestFrames::InvisibleGroup("INV_comb","Invisible System"));
512  INV_comb->AddFrame(*I_comb);
513 
514  InvMass_comb.reset(new RestFrames::SetMassInvJigsaw("InvMass_comb", "Invisible system mass Jigsaw"));
515  INV_comb->AddJigsaw(*InvMass_comb);
516 
517  JETS_comb.reset(new RestFrames::CombinatoricGroup("JETS_comb","Jets System"));
518  JETS_comb->AddFrame(*ISR_comb);
519  JETS_comb->SetNElementsForFrame(*ISR_comb, 1);
520  JETS_comb->AddFrame(*J_comb);
521  JETS_comb->SetNElementsForFrame(*J_comb, 0);
522 
523  SplitJETS_comb.reset(new RestFrames::MinMassesCombJigsaw("SplitJETS_comb", "Minimize M_{ISR} and M_{S} Jigsaw"));
524  JETS_comb->AddJigsaw(*SplitJETS_comb);
525  SplitJETS_comb->AddCombFrame(*ISR_comb, 0);
526  SplitJETS_comb->AddCombFrame(*J_comb, 1);
527  SplitJETS_comb->AddObjectFrame(*ISR_comb,0);
528  SplitJETS_comb->AddObjectFrame(*S_comb,1);
529 
530  if(!LAB_comb->InitializeAnalysis())
531  {
532  str errmsg;
533  errmsg = "Some problem occurred when calling LAB_comb->InitializeAnalysis() from the Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb analysis class.\n";
535  }
536 
537  // 2L+NJ tree (Z->ll + W/Z->qq)
538  INV_2LNJ.reset(new RestFrames::InvisibleGroup("INV_2LNJ","Invisible System"));
539  INV_2LNJ->AddFrame(*Ia_2LNJ);
540  INV_2LNJ->AddFrame(*Ib_2LNJ);
541 
542  InvMass_2LNJ.reset(new RestFrames::SetMassInvJigsaw("InvMass_2LNJ", "Invisible system mass Jigsaw"));
543  INV_2LNJ->AddJigsaw(*InvMass_2LNJ);
544  InvRapidity_2LNJ.reset(new RestFrames::SetRapidityInvJigsaw("InvRapidity_2LNJ", "Set inv. system rapidity"));
545  INV_2LNJ->AddJigsaw(*InvRapidity_2LNJ);
546  InvRapidity_2LNJ->AddVisibleFrames(S_2LNJ->GetListVisibleFrames());
547  SplitINV_2LNJ.reset(new RestFrames::ContraBoostInvJigsaw("SplitINV_2LNJ", "INV -> I_{a}+ I_{b} jigsaw"));
548  INV_2LNJ->AddJigsaw(*SplitINV_2LNJ);
549  SplitINV_2LNJ->AddVisibleFrames(Ca_2LNJ->GetListVisibleFrames(), 0);
550  SplitINV_2LNJ->AddVisibleFrames(Cb_2LNJ->GetListVisibleFrames(), 1);
551  SplitINV_2LNJ->AddInvisibleFrame(*Ia_2LNJ, 0);
552  SplitINV_2LNJ->AddInvisibleFrame(*Ib_2LNJ, 1);
553 
554  JETS_2LNJ.reset(new RestFrames::CombinatoricGroup("JETS_comb","Jets System"));
555  JETS_2LNJ->AddFrame(*J_2LNJ);
556  JETS_2LNJ->SetNElementsForFrame(*J_2LNJ, 0);
557 
558  if(!LAB_2LNJ->InitializeAnalysis())
559  {
560  str errmsg;
561  errmsg = "Some problem occurred when calling LAB_2LNJ->InitializeAnalysis() from the Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb analysis class.\n";
563  }
564 
565  // 2L+1L tree (Z->ll + Z/W->l)
566  INV_2L1L.reset(new RestFrames::InvisibleGroup("INV_2L1L","Invisible System"));
567  INV_2L1L->AddFrame(*Ia_2L1L);
568  INV_2L1L->AddFrame(*Ib_2L1L);
569 
570  InvMass_2L1L.reset(new RestFrames::SetMassInvJigsaw("InvMass_2L1L", "Invisible system mass Jigsaw"));
571  INV_2L1L->AddJigsaw(*InvMass_2L1L);
572  InvRapidity_2L1L.reset(new RestFrames::SetRapidityInvJigsaw("InvRapidity_2L1L", "Set inv. system rapidity"));
573  INV_2L1L->AddJigsaw(*InvRapidity_2L1L);
574  InvRapidity_2L1L->AddVisibleFrames(S_2L1L->GetListVisibleFrames());
575  SplitINV_2L1L.reset(new RestFrames::ContraBoostInvJigsaw("SplitINV_2L1L", "INV -> I_{a}+ I_{b} jigsaw"));
576  INV_2L1L->AddJigsaw(*SplitINV_2L1L);
577  SplitINV_2L1L->AddVisibleFrames(Ca_2L1L->GetListVisibleFrames(), 0);
578  SplitINV_2L1L->AddVisibleFrames(Cb_2L1L->GetListVisibleFrames(), 1);
579  SplitINV_2L1L->AddInvisibleFrame(*Ia_2L1L, 0);
580  SplitINV_2L1L->AddInvisibleFrame(*Ib_2L1L, 1);
581 
582  if(!LAB_2L1L->InitializeAnalysis())
583  {
584  str errmsg;
585  errmsg = "Some problem occurred when calling LAB_2L1L->InitializeAnalysis() from the Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb analysis class.\n";
587  }
588 
589  }
590  }
591 
592 
593  void run(const HEPUtils::Event* event) {
594 
595  // Clear
596  LAB_2L2J->ClearEvent();
597  LAB_comb->ClearEvent();
598  LAB_2LNJ->ClearEvent();
599  LAB_2L1L->ClearEvent();
600  LAB_3L->ClearEvent();
601 
602  // Missing energy
603  HEPUtils::P4 ptot = event->missingmom();
604  TVector3 ETMiss;
605  ETMiss.SetXYZ(ptot.px(),ptot.py(),0.0);
606 
607  // Baseline lepton objects
608  vector<const HEPUtils::Particle*> baselineElectrons, baselineMuons, baselineTaus;
609 
610  for (const HEPUtils::Particle* electron : event->electrons()) {
611  if (electron->pT() > 10. && electron->abseta() < 2.47) baselineElectrons.push_back(electron);
612  }
613 
614  // Apply electron efficiency
615  ATLAS::applyElectronEff(baselineElectrons);
616 
617  for (const HEPUtils::Particle* muon : event->muons()) {
618  if (muon->pT() > 10. && muon->abseta() < 2.4) baselineMuons.push_back(muon);
619  }
620 
621  // Apply muon efficiency
622  ATLAS::applyMuonEff(baselineMuons);
623 
624  // Photons
625  vector<const HEPUtils::Particle*> signalPhotons;
626  for (const HEPUtils::Particle* photon : event->photons()) {
627  signalPhotons.push_back(photon);
628  }
629 
630 
631  // No taus used in 13 TeV analysis?
632  //for (const HEPUtils::Particle* tau : event->taus()) {
633  //if (tau->pT() > 10. && tau->abseta() < 2.47) baselineTaus.push_back(tau);
634  //}
635  //ATLAS::applyTauEfficiencyR1(baselineTaus);
636 
637 
638  // Jets
639  vector<const HEPUtils::Jet*> bJets;
640  vector<const HEPUtils::Jet*> nonBJets;
641 
642  // Get b jets
644  const std::vector<double> a = {0,10.};
645  const std::vector<double> b = {0,10000.};
646  const std::vector<double> c = {0.77}; // set b-tag efficiency to 77%
647  HEPUtils::BinnedFn2D<double> _eff2d(a,b,c);
648  for (const HEPUtils::Jet* jet : event->jets()) {
649  bool hasTag=has_tag(_eff2d, jet->abseta(), jet->pT());
650  if (jet->pT() > 20. && fabs(jet->eta()) < 2.4) {
651  if(jet->btag() && hasTag){
652  bJets.push_back(jet);
653  } else {
654  nonBJets.push_back(jet);
655  }
656  }
657  }
658 
659  // Overlap removal is the same as the 8 TeV analysis
660  JetLeptonOverlapRemoval(nonBJets,baselineElectrons,0.2);
661  LeptonJetOverlapRemoval(baselineElectrons,nonBJets,0.4);
662  LeptonJetOverlapRemoval(baselineElectrons,bJets,0.4);
663  LeptonJetOverlapRemoval(baselineMuons,nonBJets,0.4);
664  LeptonJetOverlapRemoval(baselineMuons,bJets,0.4);
665 
666  // Fill a jet-pointer-to-bool map to make it easy to check
667  // if a given jet is treated as a b-jet in this analysis
668  map<const HEPUtils::Jet*,bool> analysisBtags;
669  for (const HEPUtils::Jet* jet : bJets) {
670  analysisBtags[jet] = true;
671  }
672  for (const HEPUtils::Jet* jet : nonBJets) {
673  analysisBtags[jet] = false;
674  }
675 
676  // Signal object containers
677  vector<const HEPUtils::Particle*> signalElectrons;
678  vector<const HEPUtils::Particle*> signalMuons;
679  vector<const HEPUtils::Particle*> signalLeptons;
680  vector<const HEPUtils::Particle*> electronsForVeto;
681  vector<const HEPUtils::Particle*> muonsForVeto;
682 
683  vector<const HEPUtils::Jet*> signalJets;
684  vector<const HEPUtils::Jet*> signalBJets;
685  vector<const HEPUtils::Jet*> signalNonBJets;
686 
687  for (const HEPUtils::Jet* jet : bJets) {
688  // pT > 20 and |eta| < 2.4 already required for jets in bJets
689  signalJets.push_back(jet);
690  signalBJets.push_back(jet);
691  }
692 
693  for (const HEPUtils::Jet* jet : nonBJets) {
694  // pT > 20 and |eta| < 2.4 already required for jets in nonBJets
695  signalJets.push_back(jet);
696  signalNonBJets.push_back(jet);
697  }
698 
699  //Put signal jets in pT order
700  std::sort(signalJets.begin(), signalJets.end(), sortByPT_RJ3L);
701  std::sort(signalBJets.begin(), signalBJets.end(), sortByPT_RJ3L);
702  std::sort(signalNonBJets.begin(), signalNonBJets.end(), sortByPT_RJ3L);
703 
704  for (const HEPUtils::Particle* electron : baselineElectrons) {
705  signalElectrons.push_back(electron);
706  signalLeptons.push_back(electron);
707  }
708 
709  for (const HEPUtils::Particle* muon : baselineMuons) {
710  signalMuons.push_back(muon);
711  signalLeptons.push_back(muon);
712  }
713 
714  std::sort(signalLeptons.begin(), signalLeptons.end(), sortLepByPT_RJ3L);
715 
716  // We now have the signal electrons, muons, jets and b jets- move on to the analysis
717  bool m_is2Lep=false;
718  bool m_is2Lep2Jet=false;
719  bool m_is2L2JInt=false;
720 
721  bool m_is3Lep=false;
722  bool m_is3LInt=false;
723  // bool m_is3Lep2Jet=false;
724  // bool m_is3Lep3Jet=false;
725 
726  // bool m_is4Lep=false;
727  // bool m_is4Lep2Jet=false;
728  // bool m_is4Lep3Jet=false;
729 
730  bool m_foundSFOS=false;
731 
732  // double m_H2PP_visible = -999.;
733  // double m_H2PP_invisible = -999.;
734  // double m_IaPP = -999.;
735  // double m_IbPP = -999.;
736  // double m_IaPa = -999.;
737  // double m_IbPb = -999.;
738  // double m_IaLAB = -999;
739  // double m_IbLAB = -999;
740  // double m_H4PP_Lept1A = -999.;
741  // double m_H4PP_Lept1B = -999.;
742  // double m_H4PP_Lept2B = -999.;
743  // double m_mu = -999;
744  // double m_pileUp_weight = -999;
745 
746 
748  // int m_nBaselineLeptons = -999;
749  // int m_nSignalLeptons = -999;
750 
751  double m_lept1Pt = -999;
752  // double m_lept1Eta = -999;
753  // double m_lept1Phi =-999;
754  double m_lept1sign=-999;
755  // double m_lept1origin = -999;
756  // double m_lept1type = -999;
757 
758  double m_lept2Pt =-999;
759  // double m_lept2Eta=-999;
760  // double m_lept2Phi =-999;
761  double m_lept2sign =-999;
762  // double m_lept2origin = -999;
763  // double m_lept2type = -999;
764 
765  double m_lept3Pt =-999;
766  // double m_lept3Eta =-999;
767  // double m_lept3Phi =-999;
768  double m_lept3sign =-999;
769  // double m_lept3origin = -999;
770  // double m_lept3type = -999;
771 
772  // double m_lept4Pt =-999;
773  // double m_lept4Eta =-999;
774  // double m_lept4Phi =-999;
775  // double m_lept4sign =-999;
776  // double m_lept4origin = -999;
777  // double m_lept4type = -999;
778  // double m_Zlep1Pt = -999;
779  // double m_Zlep1Phi = -999;
780  // double m_Zlep1Eta = -999;
781  // double m_Zlep1No = -999;
782  // double m_Zlep1sign = -999;
783 
784  // double m_Zlep2Pt = -999;
785  // double m_Zlep2sign = -999;
786  // double m_Zlep2Phi = -999;
787  // double m_Zlep2Eta = -999;
788  // double m_Zlep2No = -999;
789 
790  // double m_WlepPt = -999;
791  // double m_WlepPhi = -999;
792  // double m_WlepEta = -999;
793  // double m_WlepNo = -999;
794  // double m_Wlepsign = -999;
795 
796  // VR setup
797  // double m_lept1Pt_VR = -999;
798  // double m_lept1Eta_VR = -999;
799  // double m_lept1Phi_VR = -999;
800  // double m_lept1sign_VR = -999;
801 
802  // double m_lept2Pt_VR = -999;
803  // double m_lept2Eta_VR = -999;
804  // double m_lept2Phi_VR = -999;
805  // double m_lept2sign_VR = -999;
806 
807  //Jet Variables
808  int m_nJets=0;
809  // int m_nBtagJets=0;
810 
811  double m_jet1Pt =-999;
812  // double m_jet1Eta =-999;
813  // double m_jet1Phi =-999;
814  // double m_jet1M=-999;
815  // double m_jet1origin=-999;
816  // double m_jet1type=-999;
817 
818  double m_jet2Pt=-999;
819  // double m_jet2Eta=-999;
820  // double m_jet2Phi=-999;
821  // double m_jet2M=-999;
822  // double m_jet2origin=-999;
823  // double m_jet2type=-999;
824 
825  // double m_jet3Pt=-999;
826  // double m_jet3Eta=-999;
827  // double m_jet3Phi=-999;
828  // double m_jet3M=-999;
829  // double m_jet3origin=-999;
830  // double m_jet3type=-999;
831 
832  // double m_jet4Pt=-999;
833  // double m_jet4Eta=-999;
834  // double m_jet4Phi=-999;
835  // double m_jet4M=-999;
836  // double m_jet4origin=-999;
837  // double m_jet4type=-999;
838 
839  //Di-Lepton System: Calculated for OS Pairs
840  double m_mll=-999;
841  // double m_mt2=-999;
842  // double m_dRll=-999;
843  // double m_ptll=-999;
844  // double m_Zeta=-999;
845 
846  //Tri-Lepton System:
847  // double m_mlll=-999;
848  // double m_ptlll=-999;
849  double m_mTW=-999;
850  // double m_mTW_alt = -999;
851  // double m_mll_alt = -999;
852  //Di-Jet system: Calculated for the Two Leading Jets
853  double m_mjj=-999;
854  // double m_dRjj=-999;
855  // double m_ptjj=-999;
856  // double m_mj2j3 = -999;
857  //calculation of overall jet mass
858  // double m_mJ=-999;
859  // double m_mjjW=-999;//closest to the W-boson mass
860 
861  //Cleaning Variable: If MET is in the same direction as the Jet
862  double m_minDphi=-999;
863  // Some lab frame angles and stuff
864  // double m_dphill = -999;
865  // double m_dphilep1MET = -999;
866  // double m_dphilep2MET = -999;
867  // double m_dphilep3MET = -999;
868  // double m_dphiJMET = -999;
869  // double m_dphilll = -999;
870  // double m_dphilllMET = -999;
871  // double m_dphillMET = -999;
872  // double m_dphijj = -999;
873  // double m_dphijet1MET = -999;
874  // double m_dphijet2MET = -999;
875  // double m_dphijjMET = -999;
876  // double m_dphil3MET = -999;
877  // double m_MET=-999;
878  // double m_MET_phi = -999;
879  // double m_METTST = -999;
880  // double m_METTST_phi = -999;
881  // double m_Meff=-999;
882  // double m_LT=-999;
883 
884  // double m_MDR=-999;
885  // double m_PP_VisShape=-999;
886  // double m_gaminvPP=-999;
887  // double m_MP=-999;
888 
889  // double m_mC1=-999;
890  // double m_mN2=-999;
891 
892  // double m_mTW_Pa=-999;
893  // double m_mTW_PP=-999;
894 
895  // double m_mTZ_Pb=-999;
896  // double m_mTZ_PP=-999;
897 
898  // 3L CA
899  // double m_min_mt = -999;
900  // double m_pt_lll = -999;
901  // double m_mTl3 = -999;
902  //##############################//
903  //# Recursive Jigsaw Variables #//
904  //##############################//
905 
906  //Scale Variables
907  double m_H2PP=-999;
908  // double m_HT2PP=-999;
909  // double m_H3PP=-999;
910  // double m_HT3PP=-999;
911  double m_H4PP=-999;
912  double m_HT4PP=-999;
913  double m_H5PP=-999;
914  double m_HT5PP=-999;
915  // double m_H6PP=-999;
916  // double m_HT6PP=-999;
917 
918  double m_H2Pa=-999;
919  double m_H2Pb=-999;
920  double m_minH2P=-999;
921  // double m_R_H2Pa_H2Pb=-999;
922  double m_H3Pa=-999;
923  double m_H3Pb=-999;
924  double m_minH3P=-999;
925  // double m_R_H3Pa_H3Pb=-999;
926  double m_R_minH2P_minH3P=-999;
927  // double m_minR_pT2i_HT3Pi=-999;
928  // double m_maxR_H1PPi_H2PPi=-999;
929 
930  //Anglular Variables
931  // double m_cosPP=-999;
932  // double m_cosPa=-999;
933  // double m_cosPb=-999;
934  double m_dphiVP=-999;
935  // double m_dphiPPV=-999;
936  // double m_dphiPC1=-999;
937  // double m_dphiPN2=-999;
938 
939  // double m_sangle=-999;
940  // double m_dangle=-999;
941 
942  //Ratio Variables
943  // double m_RPZ_HT4PP=-999;
944  double m_RPT_HT4PP=-999;
945  // double m_R_HT4PP_H4PP=-999;
946 
947  // double m_RPZ_HT5PP=-999;
948  double m_RPT_HT5PP=-999;
949  // double m_R_HT5PP_H5PP=-999;
950  // double m_W_PP = -999;
951  // double m_WZ_PP = -999;
952 
954  double m_PTCM=-999;
955  double m_PTISR=-999;
956  double m_PTI=-999;
957  double m_RISR=-999;
958  // double m_cosCM=-999;
959  // double m_cosS=-999;
960  // double m_MISR=-999;
961  // double m_dphiCMI=-999;
962  // double m_dphiSI=-999;
963  double m_dphiISRI=-999;
964  // double m_HN2S=-999;
965  // double m_R_Ib_Ia=-999;
966  // double m_H11S = -999.;
967  // double m_HN1Ca = -999.;
968  // double m_HN1Cb = -999.;
969  // double m_H11Ca = -999.;
970  // double m_H11Cb = -999.;
971  // double m_cosC = -999.;
972  // double m_Is_Z = -999.;
973  // double m_Is_OS = -999;
974  double m_MZ = -999.;
975  double m_MJ = -999.;
976  // double m_mTWComp =-999.;
977  // double m_cosZ = -999.;
978  // double m_cosJ = -999.;
979  int m_NjS = 0;
980  int m_NjISR = 0;
981  int m_NbS = 0;
982  int m_NbISR = 0;
983 
984  // double m_MZ_VR = -999;
985  // double m_MJ_VR = -999;
986  // double m_PTCM_VR = -999;
987  // double m_PTISR_VR = -999;
988  // double m_PTI_VR = -999;
989  // double m_RISR_VR = -999;
990  // double m_dphiISRI_VR = -999;
991  // int m_NjS_VR = 0;
992  // int m_NjISR_VR = 0;
993 
994 
995  // double m_H2PP_VR = -999;
996  // double m_H5PP_VR = -999;
997  // double m_HT5PP_VR = -999;
998  // double m_RPT_HT5PP_VR = -999;
999  // double m_dphiVP_VR = -999;
1000  // double m_R_minH2P_minH3P_VR=-999;
1001 
1002  // double m_DPhi_METW = -999;
1003  //compressed
1004  // double m_WmassOnZ = -999;
1005  // double m_WptOnZ = -999;
1006  // double m_DPhi_METZ = -999;
1007  // double m_NonWJet_pT = -999;
1008  // double m_DPhi_METJetLeading = -999;
1009  // double m_DR_WOnZ2Jet = -999;
1010  // double m_DPhi_METNonWJet = -999;
1011  // double m_DPhi_METWonZ = -999;
1012 
1013  // Testing for low mass 3L
1014  // double m_M_I = -999;
1015  // double m_p_z_I = -999;
1016  // double m_p_z_Ia = -999;
1017  // double m_p_z_Ib = -999;
1018  // double m_boostx = -999;
1019  // double m_boosty = -999;
1020  // double m_boostz = -999;
1021 
1022  // Classify events
1023 
1024  m_nJets = signalJets.size();
1025 
1026  //if(signalLeptons.size()==2)std::cout << "m_nJets " << m_nJets << " signalLeptons.size() " << signalLeptons.size() << " pt1 " << signalLeptons[0]->pT() << " pt2 " << signalLeptons[1]->pT() << std::endl;
1027 
1028  if (signalLeptons.size()==2) m_is2Lep = true;
1029  else if (signalLeptons.size()==3) {m_is3Lep = true; //cout << "3L here" << endl;
1030  }
1031  // else if (signalLeptons.size()==4) m_is4Lep = true;
1032  //else return;
1033 
1034  if(m_is2Lep && m_nJets>1 ) m_is2Lep2Jet = true;
1035  if(m_is2Lep && m_nJets>2 && m_nJets<5) m_is2L2JInt = true;
1036  if(m_is3Lep && m_nJets>0 ) m_is3LInt = true;
1037  // if(m_is3Lep && m_nJets>1) m_is3Lep2Jet = true; //
1038  // if(m_is3Lep && m_nJets>2) m_is3Lep3Jet = true; //
1039  // if(m_is4Lep && m_nJets>1) m_is4Lep2Jet = true; //
1040  // if(m_is4Lep && m_nJets>2) m_is4Lep3Jet = true; //
1041 
1042  if(signalLeptons.size()==3)m_is3Lep=true;
1043 
1044  TLorentzVector metLV;
1045  //TLorentzVector bigFatJet;
1046  metLV.SetPxPyPzE(ptot.px(),ptot.py(),0.,sqrt(ptot.px()*ptot.px()+ptot.py()*ptot.py()));
1047 
1048  //Put the Jets in a more useful form
1049  vector<TLorentzVector> myJets;
1050  for(unsigned int ijet=0; ijet<signalJets.size();ijet++)
1051  {
1052  TLorentzVector tmp;
1053  tmp.SetPtEtaPhiM(signalJets[ijet]->pT(),signalJets[ijet]->eta(),signalJets[ijet]->phi(),signalJets[ijet]->mass());
1054  myJets.push_back(tmp);
1055  }
1056 
1057 
1058  //Put the Leptons in a more useful form
1059  vector<pair<TLorentzVector,int> > myLeptons;
1060  //vector<lep> myLeptons;
1061  for(unsigned int ilep=0; ilep<signalLeptons.size(); ilep++)
1062  {
1063  pair<TLorentzVector,int> temp;
1064  TLorentzVector tlv_temp;
1065 
1066  tlv_temp.SetPtEtaPhiM(signalLeptons[ilep]->pT(),signalLeptons[ilep]->eta(),signalLeptons[ilep]->phi(),0.0);
1067  temp.first = tlv_temp;
1068  int lepton_charge=0;
1069  if(signalLeptons[ilep]->pid()<0)lepton_charge=-1;
1070  if(signalLeptons[ilep]->pid()>0)lepton_charge=1;
1071  temp.second = lepton_charge;
1072  //temp.third = lepton_origin->at(lep_signal_index[ilep]);
1073  //temp.fourth = lepton_type->at(lep_signal_index[ilep]);
1074  //temp = make_tuple(tlv_temp,lepton_charge->at(lep_signal_index[ilep]),lepton_origin->at(lep_signal_index_[ilep]),lepton_type->at(lepton_sig_index[ilep]));
1075  myLeptons.push_back(temp);
1076  }
1077 
1078  sort(myJets.begin(), myJets.end(), SortJets);
1079  sort(myLeptons.begin(), myLeptons.end(), SortLeptons);
1080 
1081  if(m_is2Lep2Jet)
1082  {
1083 
1084  // if(myLeptons[0].first.Pt()<25.0 || myLeptons[1].first.Pt()<25.0) return;
1085 
1086  //Setting the Standard Variables
1087  //Di-Lepton System:
1088  m_lept1Pt = myLeptons[0].first.Pt();
1089  // m_lept1Eta = myLeptons[0].first.Eta();
1090  // m_lept1Phi = myLeptons[0].first.Phi();
1091  m_lept1sign = myLeptons[0].second;
1092 
1093  m_lept2Pt = myLeptons[1].first.Pt();
1094  // m_lept2Eta = myLeptons[1].first.Eta();
1095  // m_lept2Phi = myLeptons[1].first.Phi();
1096  m_lept2sign = myLeptons[1].second;
1097 
1098  m_mll = (myLeptons[0].first+myLeptons[1].first).M();
1099  // m_ptll = (myLeptons[0].first+myLeptons[1].first).Pt();
1100  // m_dRll = myLeptons[0].first.DeltaR(myLeptons[1].first);
1101  // m_Zeta = fabs(myLeptons[0].first.Eta() - myLeptons[1].first.Eta());
1102 
1103  vector<TLorentzVector> vleptons;
1104  vleptons.push_back(myLeptons[0].first);
1105  vleptons.push_back(myLeptons[1].first);
1106  //m_mt2 = myTool.GetMt2(vleptons,metLV);
1107 
1108  //min{d#phi}
1109  double mindphi=100000;
1110  double dphi=0;
1111  TLorentzVector tempjet;
1112  for(unsigned int ijet=0; ijet<signalJets.size();ijet++)
1113  {
1114  tempjet.SetPtEtaPhiM(signalJets[ijet]->pT(),signalJets[ijet]->eta(),signalJets[ijet]->phi(),signalJets[ijet]->mass());
1115 
1116  dphi = fabs(metLV.DeltaPhi(tempjet));
1117 
1118  if(dphi<mindphi) mindphi=dphi;
1119  }
1120 
1121  m_minDphi = mindphi;//cleaning variable for missmeasured jets;
1122 
1123  //just use the two leading jets
1124  int indexJ1=0;
1125  int indexJ2=1;
1126 
1127  //Di-Jet System: Here we decide which jets to use as output. The leading and sub-leading jet pair, or the jet pair with invariant mass closest to the W-Mass
1128  //jet closest to the W-boson mass
1129  m_jet1Pt = myJets[indexJ1].Pt();
1130  // m_jet1Eta = myJets[indexJ1].Eta();
1131  // m_jet1Phi = myJets[indexJ1].Phi();
1132  // m_jet1M = myJets[indexJ1].M();
1133 
1134  m_jet2Pt = myJets[indexJ2].Pt();
1135  // m_jet2Eta = myJets[indexJ2].Eta();
1136  // m_jet2Phi = myJets[indexJ2].Phi();
1137  // m_jet2M = myJets[indexJ2].M();
1138 
1139  // if(m_nJets>2) {
1140  // // m_jet3Pt = myJets[2].Pt();
1141  // // m_jet3Eta = myJets[2].Eta();
1142  // // m_jet3Phi = myJets[2].Phi();
1143  // // m_jet3M = myJets[2].M();
1144  // // m_mj2j3 = (myJets[1] + myJets[2]).M();
1145  // if(m_nJets>3) {
1146  // // m_jet4Pt = myJets[3].Pt();
1147  // // m_jet4Eta = myJets[3].Eta();
1148  // // m_jet4Phi = myJets[3].Phi();
1149  // // m_jet4M = myJets[3].M();
1150  // }
1151  // }
1152 
1153  m_mjj = (myJets[indexJ1]+myJets[indexJ2]).M();
1154  // m_ptjj = (myJets[indexJ1]+myJets[indexJ2]).Pt();
1155  // m_dRjj = myJets[indexJ1].DeltaR(myJets[indexJ2]);
1156 
1157 
1159  //Variables for the conventional approach
1160  // m_DPhi_METW = fabs((myJets[indexJ1]+myJets[indexJ2]).DeltaPhi(metLV));
1161  //for the comrpessed tree
1162  // double min_dPhi = 1000;
1163  // int WindexJ1 = -999;
1164  // int WindexJ2 = 999;
1165  // for (int j0=0;j0<myJets.size();j0++) {
1166 
1167  // double my_min_dphi= fabs(myJets[j0].DeltaPhi(metLV+myLeptons[0].first+myLeptons[1].first));
1168  // if (min_dPhi>my_min_dphi) {
1169  // min_dPhi = my_min_dphi;
1170  // WindexJ1 = j0;
1171  // }
1172  // }
1173 
1174  // min_dPhi = 1000;
1175  // for (int j1=0;j1<myJets.size();j1++) {
1176  // double my_min_dphi= fabs(myJets[j1].DeltaPhi(metLV+myLeptons[0].first+myLeptons[1].first));
1177 
1178  // if (min_dPhi>my_min_dphi) {
1179  // if (j1!=WindexJ1) {
1180  // min_dPhi = my_min_dphi;
1181  // WindexJ2 = j1;
1182  // }
1183  // }
1184  // }
1185 
1186  // m_WmassOnZ=(myJets[WindexJ1]+myJets[WindexJ2]).M();
1187  // m_WptOnZ=(myJets[WindexJ1]+myJets[WindexJ2]).Pt();
1188  // m_DPhi_METZ=fabs((myLeptons[0].first+myLeptons[1].first).DeltaPhi(metLV));
1189  // TLorentzVector nonWjetsLV;
1190  // for(int kjet=0;kjet<myJets.size();kjet++) {
1191  // if(kjet!=WindexJ1 && kjet!=WindexJ2) {
1192  // nonWjetsLV+=myJets[kjet];
1193  // }
1194  // }
1195  // if(m_nJets>2) {
1196  // // m_NonWJet_pT=nonWjetsLV.Pt();
1197  // // m_DPhi_METJetLeading = fabs(myJets[0].DeltaPhi(metLV));
1198  // // m_DR_WOnZ2Jet = myJets[WindexJ1].DeltaR(myJets[WindexJ2]);
1199  // // m_DPhi_METNonWJet = fabs(nonWjetsLV.DeltaPhi(metLV));
1200  // // m_DPhi_METWonZ = fabs((myJets[WindexJ1]+myJets[WindexJ2]).DeltaPhi(metLV));
1201  // }
1203 
1204 
1205  L1_2L2J->SetLabFrameFourVector(myLeptons[0].first); // Set lepton 4-vectors
1206  L2_2L2J->SetLabFrameFourVector(myLeptons[1].first);
1207  J1_2L2J->SetLabFrameFourVector(myJets[indexJ1]); // Set jets 4-vectors
1208  J2_2L2J->SetLabFrameFourVector(myJets[indexJ2]);
1209  TVector3 MET = ETMiss; // Get the MET
1210  MET.SetZ(0.);
1211  INV_2L2J->SetLabFrameThreeVector(MET); // Set the MET in reco tree
1212  TLorentzVector lep1;
1213  TLorentzVector lep2;
1215  //Lotentz vectors have been set, now do the boosts
1217 
1218  // Analyze the event
1219  if(!LAB_2L2J->AnalyzeEvent())
1220  {
1221  str errmsg;
1222  errmsg = "Some problem occurred when calling LAB_2L2J->AnalyzeEvent() from the Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb analysis class.\n";
1224  return;
1225  }
1226 
1227 
1228  //cout << L1_2L2J->GetFourVector(*LAB_2L2J).Pt() << endl;
1229  if (L1_2L2J->GetFourVector(*LAB_2L2J).Pt() > L2_2L2J->GetFourVector(*LAB_2L2J).Pt()){
1230  // m_Zlep1Pt = L1_2L2J->GetFourVector(*LAB_2L2J).Pt();
1231  // m_Zlep1sign = myLeptons[0].second;
1232  // m_Zlep1No = 0;
1233  // m_Zlep2Pt = L2_2L2J->GetFourVector(*LAB_2L2J).Pt();
1234  // m_Zlep2sign = myLeptons[1].second;
1235  // m_Zlep2No = 1;
1236  lep1 = L1_2L2J->GetFourVector(*LAB_2L2J);
1237  lep2 = L2_2L2J->GetFourVector(*LAB_2L2J);
1238  }
1239  else {
1240  // m_Zlep1Pt = L2_2L2J->GetFourVector(*LAB_2L2J).Pt();
1241  // m_Zlep1sign = myLeptons[1].second;
1242  // m_Zlep1No = 1;
1243  // m_Zlep2Pt = L1_2L2J->GetFourVector(*LAB_2L2J).Pt();
1244  // m_Zlep2sign = myLeptons[0].second;
1245  // m_Zlep2No = 0;
1246 
1247  lep1 = L2_2L2J->GetFourVector(*LAB_2L2J);
1248  lep2 = L1_2L2J->GetFourVector(*LAB_2L2J);
1249  }
1250  // set the jet lab frame 4-vector
1251  TLorentzVector jet1 = J1_2L2J->GetFourVector(*LAB_2L2J);
1252  TLorentzVector jet2 = J2_2L2J->GetFourVector(*LAB_2L2J);
1253 
1254  // Some lab frame stuff
1255  // m_dphill = lep1.DeltaPhi(lep2);
1256  // m_dphilep1MET = fabs(lep1.DeltaPhi(metLV));
1257  // m_dphilep2MET = fabs(lep2.DeltaPhi(metLV));
1258  // m_dphillMET = fabs((lep1 + lep2).DeltaPhi(metLV));
1259  // m_dphijet1MET = fabs(jet1.DeltaPhi(metLV));
1260  // m_dphijet2MET = fabs(jet2.DeltaPhi(metLV));
1261  // m_dphijj = fabs(jet1.DeltaPhi(jet2));
1262  // m_dphijjMET = fabs((jet1 + jet2).DeltaPhi(metLV));
1263  //... then by setting the Variables
1264  TLorentzVector vP_V1aPP = J1_2L2J->GetFourVector(*C1N2_2L2J);
1265  TLorentzVector vP_V2aPP = J2_2L2J->GetFourVector(*C1N2_2L2J);
1266  TLorentzVector vP_V1bPP = L1_2L2J->GetFourVector(*C1N2_2L2J);
1267  TLorentzVector vP_V2bPP = L2_2L2J->GetFourVector(*C1N2_2L2J);
1268  TLorentzVector vP_IaPP = X1a_2L2J->GetFourVector(*C1N2_2L2J);
1269  TLorentzVector vP_IbPP = X1b_2L2J->GetFourVector(*C1N2_2L2J);
1270 
1271  TLorentzVector vP_V1aPa = J1_2L2J->GetFourVector(*C1a_2L2J);
1272  TLorentzVector vP_V2aPa = J2_2L2J->GetFourVector(*C1a_2L2J);
1273  TLorentzVector vP_IaPa = X1a_2L2J->GetFourVector(*C1a_2L2J);
1274  TLorentzVector vP_V1bPb = L1_2L2J->GetFourVector(*N2b_2L2J);
1275  TLorentzVector vP_V2bPb = L2_2L2J->GetFourVector(*N2b_2L2J);
1276  TLorentzVector vP_IbPb = X1b_2L2J->GetFourVector(*N2b_2L2J);
1277 
1278 
1279  //Variables w/ 4 objects
1280  //Four vector sum of all visible objets + four vector sum of inv objects
1281  m_H2PP = (vP_V1aPP + vP_V2aPP + vP_V1bPP + vP_V2bPP).P() + (vP_IaPP+vP_IbPP).P();//H(1,1)PP
1282  // m_HT2PP = (vP_V1aPP + vP_V2aPP + vP_V1bPP + vP_V2bPP).Pt() + (vP_IaPP+vP_IbPP).Pt();//HT(1,1)PP
1283  //Scalar sum of all visible objects + vector sum of invisible momenta
1284  m_H5PP = vP_V1aPP.P() + vP_V2aPP.P() + vP_V1bPP.P() + vP_V2bPP.P() + (vP_IaPP + vP_IbPP).P();//H(4,1)PP
1285  m_HT5PP = vP_V1aPP.Pt() + vP_V2aPP.Pt() + vP_V1bPP.Pt() + vP_V2bPP.Pt() + (vP_IaPP + vP_IbPP).Pt();//HT(4,1)PP
1286  //scalar sum of all objects
1287  // m_H6PP = vP_V1aPP.P() + vP_V2aPP.P() + vP_V1bPP.P() + vP_V2bPP.P() + vP_IaPP.P() + vP_IbPP.P();//H(4,2)PP
1288  // m_HT6PP = vP_V1aPP.Pt() + vP_V2aPP.Pt() + vP_V1bPP.Pt() + vP_V2bPP.Pt() + vP_IaPP.Pt() + vP_IbPP.Pt();
1289 
1290  m_H2Pa = (vP_V1aPa + vP_V2aPa).P() + vP_IaPa.P();
1291  m_H2Pb = (vP_V1bPb + vP_V2bPb).P() + vP_IbPb.P();
1292  m_H3Pa = vP_V1aPa.P() + vP_V2aPa.P() + vP_IaPa.P();
1293  m_H3Pb = vP_V1bPb.P() + vP_V2bPb.P() + vP_IbPb.P();
1294  m_minH2P = std::min(m_H2Pa,m_H2Pb);
1295  m_minH3P = std::min(m_H3Pa,m_H3Pb);
1296  // m_R_H2Pa_H2Pb = m_H2Pa/m_H2Pb;
1297  // m_R_H3Pa_H3Pb = m_H3Pa/m_H3Pb;
1298  m_R_minH2P_minH3P = m_minH2P/m_minH3P;
1299  // std::cout << " m_R_minH2P_minH3P " << m_R_minH2P_minH3P << " " << m_minH2P << " " << m_minH3P << std::endl;
1300  // double H3PTa = vP_V1aPa.Pt() + vP_V2aPa.Pt() + vP_IaPa.Pt();
1301 
1302  // m_minR_pT2i_HT3Pi = std::min(vP_V1aPa.Pt()/H3PTa,vP_V2aPa.Pt()/H3PTa);
1303 
1304 
1305  // m_R_HT5PP_H5PP = m_HT5PP/m_H5PP;
1306 
1307  // Invisible in the PP frame, Px frames and lab frame
1308  TLorentzVector vP_IaLAB = X1a_2L2J->GetFourVector(*LAB_2L2J);
1309  TLorentzVector vP_IbLAB = X1b_2L2J->GetFourVector(*LAB_2L2J);
1310  // m_IaLAB = vP_IaLAB.P();
1311  // m_IbLAB = vP_IbLAB.P();
1312  // m_IaPP = vP_IaPP.P();
1313  // m_IbPP = vP_IbPP.P();
1314  // m_IaPa = vP_IaPa.P();
1315  // m_IbPb = vP_IbPb.P();
1316 
1317  // double jetMetphiP = (vP_V1aPa+vP_V2aPa).DeltaPhi(vP_IaPa);
1318  // m_mTW_Pa = sqrt(2*(vP_V1aPa+vP_V2aPa).Pt()*vP_IaPa.Pt()*(1-cos(jetMetphiP)));
1319 
1320  // double jetMetphiPP = (vP_V1aPP+vP_V2aPP).DeltaPhi(vP_IaPP+vP_IbPP);
1321  // m_mTW_PP = sqrt(2*(vP_V1aPP+vP_V2aPP).Pt()*(vP_IaPP+vP_IbPP).Pt()*(1-cos(jetMetphiPP)));
1322 
1323  // double dilepMetphiP = (vP_V1bPb+vP_V2bPb).DeltaPhi(vP_IbPb);
1324  // m_mTZ_Pb = sqrt(2*(vP_V1bPb+vP_V2bPb).Pt()*vP_IbPb.Pt()*(1-cos(dilepMetphiP)));
1325 
1326  // double dilepMetphiPP = (vP_V1bPP+vP_V2bPP).DeltaPhi(vP_IaPP+vP_IbPP);
1327  // m_mTZ_PP = sqrt(2*(vP_V1bPP+vP_V2bPP).Pt()*(vP_IaPP+vP_IbPP).Pt()*(1-cos(dilepMetphiPP)));
1328 
1329 
1330  // double H1PPa = (vP_V1aPP + vP_V2aPP).P();
1331  // double H1PPb = (vP_V1bPP + vP_V2bPP).P();
1332  // double H2PPa = vP_V1aPP.P() + vP_V2aPP.P();
1333  // double H2PPb = vP_V1bPP.P() + vP_V2bPP.P();
1334  // m_maxR_H1PPi_H2PPi = std::max(H1PPa/H2PPa,H1PPb/H2PPb);
1335 
1336  // signal variables
1337  TLorentzVector vP_Va = C1a_2L2J->GetVisibleFourVector(*C1a_2L2J);
1338  TLorentzVector vP_Vb = N2b_2L2J->GetVisibleFourVector(*N2b_2L2J);
1339  // m_MP = (vP_Va.M2()-vP_Vb.M2())/(2.*(vP_Va.E()-vP_Vb.E()));
1340 
1341  // double P_P = C1a_2L2J->GetMomentum(*C1N2_2L2J);
1342 
1343  // double MPP = 2.*sqrt(P_P*P_P + m_MP*m_MP);
1344  TVector3 vP_PP = C1N2_2L2J->GetFourVector(*LAB_2L2J).Vect();
1345  double Pt_PP = vP_PP.Pt();
1346  // double Pz_PP = fabs(vP_PP.Pz());
1347  m_RPT_HT5PP = Pt_PP / (Pt_PP + m_HT5PP);
1348  // m_RPZ_HT5PP = Pz_PP / (Pz_PP + m_HT5PP);
1349 
1350  // m_PP_VisShape = C1N2_2L2J->GetVisibleShape();
1351 
1352  // m_gaminvPP = 2.*m_MP/MPP;
1353  // m_MDR = m_PP_VisShape*C1N2_2L2J->GetMass();
1354 
1355  // m_mC1 = C1a_2L2J->GetMass();
1356  // m_mN2 = N2b_2L2J->GetMass();
1357 
1358 
1359  //Angular properties of the sparticles system
1360  // m_cosPP = C1N2_2L2J->GetCosDecayAngle(); //decay angle of the PP system
1361  // m_cosPa = C1a_2L2J->GetCosDecayAngle(*X1a_2L2J);//decay angle of the C1a system
1362  // m_cosPb = N2b_2L2J->GetCosDecayAngle(*X1b_2L2J);//decay angle of the N2b system
1363 
1364  //difference in azimuthal angle between the total sum of visible ojects in the C1N2 frame
1365  // m_dphiPPV = C1N2_2L2J->GetDeltaPhiBoostVisible();
1366  m_dphiVP = C1N2_2L2J->GetDeltaPhiDecayVisible();
1367 
1368  //hemisphere variables
1369  // m_dphiPC1 = C1a_2L2J->GetDeltaPhiDecayPlanes(*Wa_2L2J);
1370  // m_dphiPN2 = N2b_2L2J->GetDeltaPhiDecayPlanes(*Zb_2L2J);
1371 
1372  // m_sangle =(m_cosPa+(m_dphiVP-acos(-1.)/2.)/(acos(-1.)/2.))/2.;
1373  // m_dangle =(m_cosPa-(m_dphiVP-acos(-1.)/2.)/(acos(-1.)/2.))/2.;
1374  }//end is 2L2J event
1375 
1376  if(m_is3Lep){
1377 
1378  // bool m_pass3L_presel;
1379 
1380  // if(myLeptons[0].first.Pt()<25.0 || myLeptons[1].first.Pt()<25.0 || myLeptons[2].first.Pt()<20.0) {
1381  // m_pass3L_presel=false;
1382  // }
1383  // else {
1384  // m_pass3L_presel=true;
1385  // }
1386 
1387  //if(!m_pass3L_presel)return;
1388 
1389  //Tri-Lepton System
1390  //Here we choose leptons based on where they "come from"
1391  //lept1 and lept2 are the lepton pair with invariant mass closest to the Z-Mass
1392  //lept3 is the remaining lepton
1393  //This is meant to emulate lept1 and lept2 being produced by the Z, while lept3 is produced by the W
1394 
1395  double diff = 10000000000.0;
1396  int Zlep1 = -99;
1397  int Zlep2 = -99;
1398  double Zmass = -999.0;
1399  bool foundSFOS = false;
1400 
1401  for(unsigned int i=0; i<myLeptons.size(); i++)
1402  {
1403  for(unsigned int j=i+1; j<myLeptons.size(); j++)
1404  {
1405  //Opposite-Sign
1406  if(myLeptons[i].second*myLeptons[j].second<0)
1407  {
1408  //Same-Flavor
1409  if(abs(myLeptons[i].second)==abs(myLeptons[j].second))
1410  {
1411  double mass = (myLeptons[i].first+myLeptons[j].first).M();
1412  double massdiff = fabs(mass-91.1876);
1413  if(massdiff<diff)
1414  {
1415  diff=massdiff;
1416  Zmass=mass;
1417  Zlep1 = i;
1418  Zlep2 = j;
1419  foundSFOS = true;
1420  }
1421  }
1422  }
1423  }
1424  }
1425 
1426  if(!foundSFOS) {
1427  m_foundSFOS=false;
1428  }
1429  else {
1430  m_foundSFOS=true;
1431  }
1432 
1433  if(m_foundSFOS){
1434 
1435  int Wlep1 = -999;
1436  if( (Zlep1==0 && Zlep2==1) || (Zlep1==1 && Zlep2==0) ) Wlep1=2;
1437  else if( (Zlep1==0 && Zlep2==2) || (Zlep1==2 && Zlep2==0) ) Wlep1=1;
1438  else if((Zlep1==1 && Zlep2==2) || (Zlep1==2 && Zlep2==1) ) Wlep1=0;
1439 
1440  //Knowing the indices, we perform assignments
1441  m_lept1Pt = myLeptons[0].first.Pt();
1442  m_lept1sign = myLeptons[0].second;
1443 
1444  m_lept2Pt = myLeptons[1].first.Pt();
1445  m_lept2sign = myLeptons[1].second;
1446 
1447  m_lept3Pt = myLeptons[2].first.Pt();
1448  m_lept3sign = myLeptons[2].second;
1449 
1450 
1451  m_mll = Zmass; //based on mass minimization
1452 
1453  vector<TLorentzVector> Leptons;
1454  Leptons.push_back(myLeptons[Wlep1].first);
1455  Leptons.push_back(myLeptons[Zlep1].first);
1456  Leptons.push_back(myLeptons[Zlep2].first);
1457 
1458 
1459  double wlepMetphi = myLeptons[Wlep1].first.DeltaPhi(metLV);
1460 
1461  m_mTW = sqrt(2*myLeptons[Wlep1].first.Pt()*metLV.Pt()*(1-cos(wlepMetphi)));
1462 
1463  INV_3L->SetLabFrameThreeVector(ETMiss); //set the MET in the event
1464 
1465 
1466  L1a_3L->SetLabFrameFourVector(Leptons[0]); // Set lepton from W
1467  L1b_3L->SetLabFrameFourVector(Leptons[1]); // Set lepton1 from Z
1468  L2b_3L->SetLabFrameFourVector(Leptons[2]); // Set lepton2 from Z
1469 
1470  if(!LAB_3L->AnalyzeEvent())
1471  {
1472  str errmsg;
1473  errmsg = "Some problem occurred when calling LAB_3L->AnalyzeEvent() from the Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb analysis class.\n";
1475  return;
1476  }
1477 
1478 
1479  TLorentzVector l1;
1480  TLorentzVector l2;
1481  TLorentzVector l3 = L1a_3L->GetFourVector(*LAB_3L);
1482 
1483  //if(DEBUG) cout << "WlepPt: " << m_WlepPt << " Wlepsign: " << m_Wlepsign << endl;
1484  if (L1b_3L->GetFourVector(*LAB_3L).Pt() > L2b_3L->GetFourVector(*LAB_3L).Pt()){
1485  l1 = L1b_3L->GetFourVector(*LAB_3L);
1486  l2 = L2b_3L->GetFourVector(*LAB_3L);
1487  }
1488  else {
1489  l2 = L1b_3L->GetFourVector(*LAB_3L);
1490  l1 = L2b_3L->GetFourVector(*LAB_3L);
1491  }
1492 
1493  // More lab frame stuff
1494 
1495  //if(DEBUG) cout << "Zlep1: " << m_Zlep1Pt << " " << m_Zlep1sign << " Zlep2Pt: " << m_Zlep2Pt << " " << m_Zlep2sign << endl;
1496  TLorentzVector vP_V1aPP = L1a_3L->GetFourVector(*C1N2_3L);
1497  TLorentzVector vP_V1bPP = L1b_3L->GetFourVector(*C1N2_3L);
1498  TLorentzVector vP_V2bPP = L2b_3L->GetFourVector(*C1N2_3L);
1499  TLorentzVector vP_I1aPP = X1a_3L->GetFourVector(*C1N2_3L);
1500  TLorentzVector vP_I1bPP = X1b_3L->GetFourVector(*C1N2_3L);
1501 
1502  TLorentzVector vP_V1aPa = L1a_3L->GetFourVector(*C1a_3L);
1503  TLorentzVector vP_I1aPa = X1a_3L->GetFourVector(*C1a_3L);
1504 
1505  TLorentzVector vP_V1bPb = L1b_3L->GetFourVector(*N2b_3L);
1506  TLorentzVector vP_V2bPb = L2b_3L->GetFourVector(*N2b_3L);
1507  TLorentzVector vP_I1bPb = X1b_3L->GetFourVector(*N2b_3L);
1508 
1509 
1510 
1511  //Variables w/ 4 objects
1512 
1514  //Four vector sum of all visible objets + four vector sum of inv objects
1515 
1516  //Scalar sum of all visible objects + vector sum of invisible momenta
1517  m_H4PP = vP_V1aPP.P() + vP_V1bPP.P() + vP_V2bPP.P() + (vP_I1aPP + vP_I1bPP).P();//H(3,1)PP
1518  m_HT4PP = vP_V1aPP.Pt() + vP_V1bPP.Pt() + vP_V2bPP.Pt() + (vP_I1aPP + vP_I1bPP).Pt();//HT(3,1)PP
1519 
1520  // Invisible components again
1521  TLorentzVector vP_IaLAB = X1a_3L->GetFourVector(*LAB_3L);
1522  TLorentzVector vP_IbLAB = X1b_3L->GetFourVector(*LAB_3L);
1523 
1524 
1525  // Testing for low mass 3L
1526  TLorentzVector p_Ia_Lab = X1a_3L->GetFourVector(*LAB_3L);
1527  TLorentzVector p_Ib_Lab = X1b_3L->GetFourVector(*LAB_3L);
1528  TVector3 lab_to_pp = C1N2_3L->GetBoostInParentFrame();
1529 
1531  m_H2Pa = (vP_V1aPa).P() + (vP_I1aPa).P(); //H(1,1)P
1532  m_H2Pb = (vP_V1bPb + vP_V2bPb).P() + vP_I1bPb.P();//H(1,1)P
1533 
1534  m_H3Pa = vP_V1aPa.P() + vP_I1aPa.P();//H(1,1)P
1535  m_H3Pb = vP_V1bPb.P() + vP_V2bPb.P() + vP_I1bPb.P();//H(2,1)P
1536 
1537  m_minH2P = std::min(m_H2Pa,m_H2Pb);
1538  m_minH3P = std::min(m_H3Pa,m_H3Pb);
1539  // m_R_H2Pa_H2Pb = m_H2Pa/m_H2Pb;
1540  // m_R_H3Pa_H3Pb = m_H3Pa/m_H3Pb;
1541  m_R_minH2P_minH3P = m_H2Pb/m_H3Pb;
1542 
1543  // double H1PPa = (vP_V1aPP).P();
1544  // double H1PPb = (vP_V1bPP + vP_V2bPP).P();
1545  // double H2PPa = vP_V1aPP.P() + vP_I1aPP.P();
1546  // double H2PPb = (vP_V1bPP+vP_V2bPP).P() + vP_I1bPP.P();
1547  // m_maxR_H1PPi_H2PPi = std::max(H1PPa/H2PPa,H1PPb/H2PPb);
1548 
1550  //m_dRll_I_PP = (vP_V1bPP+vP_V1bPP).DeltaR(vP_I1bPP);
1551  //m_R_Ib_Ia = (vP_V1bPP + vP_V2bPP + vP_I1bPP).P()/(vP_V1aPP+vP_I1aPP).P();
1552 
1553  // signal variables
1554  TLorentzVector vP_Va = C1a_3L->GetVisibleFourVector(*C1a_3L);
1555  TLorentzVector vP_Vb = N2b_3L->GetVisibleFourVector(*N2b_3L);
1556 
1557  TVector3 vP_PP = C1N2_3L->GetFourVector(*LAB_3L).Vect();
1558  double Pt_PP = vP_PP.Pt();
1559  m_RPT_HT4PP = Pt_PP / (Pt_PP + m_HT4PP);
1560 
1561 
1562  // mt_min here
1563 
1564  /*double min0 = -999;
1565  double min1 = -999;
1566  double min2 = -999;
1567  double lepmetphi0 = myLeptons[0].first.DeltaPhi(metLV);
1568  double lepmetphi1 = myLeptons[1].first.DeltaPhi(metLV);
1569  double lepmetphi2 = myLeptons[2].first.DeltaPhi(metLV);
1570 
1571  if (myLeptons[0].second == -myLeptons[1].second) min0 = sqrt(2*myLeptons[2].first.Pt()*metLV.Pt()*(1-cos(lepmetphi2)));
1572  if (myLeptons[0].second == -myLeptons[2].second) min1 = sqrt(2*myLeptons[1].first.Pt()*metLV.Pt()*(1-cos(lepmetphi1)));
1573  if (myLeptons[1].second == -myLeptons[2].second) min2 = sqrt(2*myLeptons[0].first.Pt()*metLV.Pt()*(1-cos(lepmetphi0)));
1574 
1575  if (min0 > 0 && min1 > 0) m_min_mt = min(min0,min1);
1576  else if (min0 > 0 && min2 > 0) m_min_mt = min(min0,min2);
1577  else if (min1 > 0 && min2 > 0) m_min_mt = min(min1,min2);
1578  else if (min0 > 0 && min1 < 0 && min2 < 0) m_min_mt = min0;
1579  else if (min1 > 0 && min0 < 0 && min2 < 0) m_min_mt = min1;
1580  else if (min2 > 0 && min0 < 0 && min1 <0) m_min_mt = min2;*/
1581 
1582  } // end of if(m_foundSFOS)
1583  } // end of m_is3Lep
1584 
1585  if(m_is3LInt || m_is2L2JInt) {
1586 
1587  //min{d#phi}
1588  double mindphi=100000;
1589  double dphi=0;
1590  TLorentzVector tempjet;
1591  for(unsigned int ijet=0; ijet<signalJets.size();ijet++)
1592  {
1593  tempjet.SetPtEtaPhiM(signalJets[ijet]->pT(),signalJets[ijet]->eta(),signalJets[ijet]->phi(),signalJets[ijet]->mass());
1594 
1595  dphi = fabs(metLV.DeltaPhi(tempjet));
1596 
1597  if(dphi<mindphi) mindphi=dphi;
1598  }
1599 
1600  m_minDphi = mindphi;//cleaning variable for missmeasured jets;
1601  //if( fabs(mindphi)<0.4) return;
1602 
1603 
1604 
1605  vector<RestFrames::RFKey> jetID;
1606  for(int i = 0; i < int(myJets.size()); i++){
1607 
1608  TLorentzVector jet = myJets[i];
1609 
1610  // transverse view of jet 4-vectors
1611  jet.SetPtEtaPhiM(jet.Pt(),0.0,jet.Phi(),jet.M());
1612  jetID.push_back(JETS_comb->AddLabFrameFourVector(jet));
1613  }
1614 
1615  TLorentzVector lepSys(0.,0.,0.,0.);
1616  for(int i = 0; i < int(myLeptons.size()); i++){
1617  TLorentzVector lep1;
1618  // transverse view of mu 4-vectors
1619  lep1.SetPtEtaPhiM(myLeptons[i].first.Pt(),0.0,myLeptons[i].first.Phi(),myLeptons[i].first.M());
1620  lepSys = lepSys + lep1;
1621  }
1622  L_comb->SetLabFrameFourVector(lepSys);
1623 
1624  INV_comb->SetLabFrameThreeVector(ETMiss);
1625 
1626  if(!LAB_comb->AnalyzeEvent())
1627  {
1628  str errmsg;
1629  errmsg = "Some problem occurred when calling LAB_comb->AnalyzeEvent() from the Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb analysis class.\n";
1631  return;
1632  }
1633 
1634  for(int i = 0; i < int(signalJets.size()); i++){
1635  if(JETS_comb->GetFrame(jetID[i]) == *J_comb){
1636  m_NjS++;
1637  if( analysisBtags.at(signalJets[i]) ) m_NbS++;
1638  } else {
1639  m_NjISR++;
1640  if( analysisBtags.at(signalJets[i]) ) m_NbISR++;
1641  }
1642  }
1643 
1644  // 2LNJ analysis
1645  if(m_is2L2JInt){
1646 
1647  // put jets in their place
1648  int NJ = jetID.size();
1649  TLorentzVector vISR(0.,0.,0.,0.);
1650  for(int i = 0; i < NJ; i++){
1651  if(JETS_comb->GetFrame(jetID[i]) == *J_comb){
1652  JETS_2LNJ->AddLabFrameFourVector(myJets[i]);
1653  } else {
1654  vISR += myJets[i];
1655  }
1656  }
1657 
1658  ISR_2LNJ->SetLabFrameFourVector(vISR);
1659 
1660  // put leptons in their place
1661  L1_2LNJ->SetLabFrameFourVector(myLeptons[0].first);
1662  L2_2LNJ->SetLabFrameFourVector(myLeptons[1].first);
1663 
1664  INV_2LNJ->SetLabFrameThreeVector(ETMiss);
1665 
1666  if(!LAB_2LNJ->AnalyzeEvent())
1667  {
1668  str errmsg;
1669  errmsg = "Some problem occurred when calling LAB_2LNJ->AnalyzeEvent() from the Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb analysis class.\n";
1671  return;
1672  }
1673 
1674  }
1675 
1676  // 2L1L analysis
1677  if(m_is3LInt){
1678 
1679  // put jets in their place
1680  int NJ = jetID.size();
1681  TLorentzVector vISR(0.,0.,0.,0.);
1682  for(int i = 0; i < NJ; i++){
1683  if(JETS_comb->GetFrame(jetID[i]) != *J_comb) vISR += myJets[i];
1684  }
1685 
1686  ISR_2L1L->SetLabFrameFourVector(vISR);
1687 
1688  // put leptons in their place
1689  // find min mass OS pair
1690  pair<int,int> iSFOS;
1691  double mSFOS = -1.;
1692  for(int i = 0; i < 2; i++){
1693  for(int j = i+1; j < 3; j++){
1694  if((signbit(myLeptons[i].second) && !signbit(myLeptons[j].second)) || (!signbit(myLeptons[i].second) && signbit(myLeptons[j].second))){
1695  if(mSFOS < 0. ||
1696  (myLeptons[i].first+myLeptons[j].first).M() < mSFOS){
1697  mSFOS = (myLeptons[i].first+myLeptons[j].first).M();
1698  iSFOS.first = i;
1699  iSFOS.second = j;
1700  }
1701  }
1702  }
1703  }
1704 
1705  for(int i = 0; i < 3; i++){
1706  if(i == iSFOS.first)
1707  L1_2L1L->SetLabFrameFourVector(myLeptons[i].first);
1708  if(i == iSFOS.second)
1709  L2_2L1L->SetLabFrameFourVector(myLeptons[i].first);
1710  if(i != iSFOS.first && i != iSFOS.second) {
1711  Lb_2L1L->SetLabFrameFourVector(myLeptons[i].first);
1712  //calculate the mTWComp with the remaining lepton
1713  TLorentzVector themetLV;
1714  themetLV.SetPxPyPzE(ETMiss.X(),ETMiss.Y(),0.,sqrt(ETMiss.X()*ETMiss.X()+ETMiss.Y()*ETMiss.Y()));
1715  // double wlepMetphi = myLeptons[i].first.DeltaPhi(themetLV);
1716  // m_mTWComp = sqrt(2*myLeptons[i].first.Pt()*themetLV.Pt()*(1-cos(wlepMetphi)));
1717  }
1718  }
1719 
1720  INV_2L1L->SetLabFrameThreeVector(ETMiss);
1721 
1722  if(!LAB_2L1L->AnalyzeEvent())
1723  {
1724  str errmsg;
1725  errmsg = "Some problem occurred when calling LAB_2L1L->AnalyzeEvent() from the Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb analysis class.\n";
1727  return;
1728  }
1729 
1730  }
1731 
1732  TLorentzVector vP_CM;
1733  TLorentzVector vP_ISR;
1734  TLorentzVector vP_I;
1735 
1736  if(m_is2L2JInt){
1737 
1738  vP_CM = CM_2LNJ->GetFourVector();
1739  vP_ISR = ISR_2LNJ->GetFourVector();
1740  vP_I = (*Ia_2LNJ+*Ib_2LNJ).GetFourVector();
1741 
1742  // m_cosCM = CM_2LNJ->GetCosDecayAngle();
1743  // m_cosS = S_2LNJ->GetCosDecayAngle();
1744  // m_MISR = ISR_2LNJ->GetMass();
1745  // m_dphiCMI = acos(-1.)-fabs(CM_2LNJ->GetDeltaPhiBoostVisible());
1746  // m_dphiSI = acos(-1.)-fabs(S_2LNJ->GetDeltaPhiBoostVisible());
1747 
1748  // m_HN2S = //Z_2LNJ->GetFourVector(*S_2LNJ).E() +
1749  // L1_2LNJ->GetFourVector(*S_2LNJ).E()+
1750  // L2_2LNJ->GetFourVector(*S_2LNJ).E()+
1751  // J_2LNJ->GetFourVector(*S_2LNJ).E() +
1752  // Ia_2LNJ->GetFourVector(*S_2LNJ).P() +
1753  // Ib_2LNJ->GetFourVector(*S_2LNJ).P();
1754  // m_H11S = 2.*(*Ia_2LNJ+*Ib_2LNJ).GetFourVector(*S_2LNJ).P();
1755  // m_HN1Ca = Z_2LNJ->GetFourVector(*Ca_2LNJ).E()+
1756  Ia_2LNJ->GetFourVector(*Ca_2LNJ).P();
1757  // m_HN1Cb = J_2LNJ->GetFourVector(*Cb_2LNJ).E()+
1758  // Ib_2LNJ->GetFourVector(*Cb_2LNJ).P();
1759  // m_H11Ca = 2.*Ia_2LNJ->GetFourVector(*Ca_2LNJ).P();
1760  // m_H11Cb = 2.*Ib_2LNJ->GetFourVector(*Cb_2LNJ).P();
1761  // m_cosC = Ca_2LNJ->GetCosDecayAngle();
1762 
1763  // if((signbit(myLeptons[0].second) && !signbit(myLeptons[1].second)) || (!signbit(myLeptons[0].second) && signbit(myLeptons[1].second))) m_Is_OS = 1;
1764  // if(myLeptons[0].second+myLeptons[1].second == 0) m_Is_Z = 1;
1765  m_MZ = Z_2LNJ->GetMass();
1766  m_MJ = J_2LNJ->GetMass();
1767 
1768  // m_cosZ = Z_2LNJ->GetCosDecayAngle();
1769  //if(m_NjS > 1)
1770  // m_cosJ = JSA_2LNJ->GetCosDecayAngle();
1771  // m_dphiJMET = fabs(J_2LNJ->GetFourVector(*LAB_2LNJ).DeltaPhi(metLV));
1772  }
1773 
1774  if(m_is3LInt){
1775 
1776  vP_CM = CM_2L1L->GetFourVector();
1777  vP_ISR = ISR_2L1L->GetFourVector();
1778  vP_I = (*Ia_2L1L+*Ib_2L1L).GetFourVector();
1779 
1780  // m_cosCM = CM_2L1L->GetCosDecayAngle();
1781  // m_cosS = S_2L1L->GetCosDecayAngle();
1782  // m_MISR = ISR_2L1L->GetMass();
1783  // m_dphiCMI = acos(-1.)-fabs(CM_2L1L->GetDeltaPhiBoostVisible());
1784  // m_dphiSI = acos(-1.)-fabs(S_2L1L->GetDeltaPhiBoostVisible());
1785 
1786  // m_HN2S = //Z_2L1L->GetFourVector(*S_2L1L).E() +
1787  // L1_2L1L->GetFourVector(*S_2L1L).E() +
1788  // L2_2L1L->GetFourVector(*S_2L1L).E() +
1789  // Lb_2L1L->GetFourVector(*S_2L1L).E() +
1790  // Ia_2L1L->GetFourVector(*S_2L1L).P() +
1791  // Ib_2L1L->GetFourVector(*S_2L1L).P();
1792  // m_H11S = 2.*(*Ia_2L1L+*Ib_2L1L).GetFourVector(*S_2L1L).P();
1793  // m_HN1Ca = Z_2L1L->GetFourVector(*Ca_2L1L).E()+
1794  Ia_2L1L->GetFourVector(*Ca_2L1L).P();
1795  // m_HN1Cb = Lb_2L1L->GetFourVector(*Cb_2L1L).E()+
1796  // Ib_2L1L->GetFourVector(*Cb_2L1L).P();
1797  // m_H11Ca = 2.*Ia_2L1L->GetFourVector(*Ca_2L1L).P();
1798  // m_H11Cb = 2.*Ib_2L1L->GetFourVector(*Cb_2L1L).P();
1799  // m_cosC = Ca_2L1L->GetCosDecayAngle();
1800  // m_Is_OS = 1;
1801  // if(myLeptons[0].second+myLeptons[1].second == 0 ||
1802  // myLeptons[0].second+myLeptons[2].second == 0 ||
1803  // myLeptons[1].second+myLeptons[2].second == 0) m_Is_Z=1;
1804  m_MZ = Z_2L1L->GetMass();
1805  // m_cosZ = Z_2L1L->GetCosDecayAngle();
1806  }
1807 
1808  m_PTCM = vP_CM.Pt();
1809 
1810  TVector3 boostZ = vP_CM.BoostVector();
1811  boostZ.SetX(0.);
1812  boostZ.SetY(0.);
1813 
1814  vP_CM.Boost(-boostZ);
1815  vP_ISR.Boost(-boostZ);
1816  vP_I.Boost(-boostZ);
1817 
1818  TVector3 boostT = vP_CM.BoostVector();
1819  vP_ISR.Boost(-boostT);
1820  vP_I.Boost(-boostT);
1821 
1822  TVector3 vPt_ISR = vP_ISR.Vect();
1823  TVector3 vPt_I = vP_I.Vect();
1824  vPt_ISR -= vPt_ISR.Dot(boostZ.Unit())*boostZ.Unit();
1825  vPt_I -= vPt_I.Dot(boostZ.Unit())*boostZ.Unit();
1826 
1827  m_PTISR = vPt_ISR.Mag();
1828  m_RISR = -vPt_I.Dot(vPt_ISR.Unit()) / m_PTISR;
1829  m_PTI = vPt_I.Mag();
1830  m_dphiISRI = fabs(vPt_ISR.Angle(vPt_I));
1831 
1832  }//end INTERMEDIATE
1833 
1834  // Cutflow check
1835 
1836  cutFlowVector_str[0] = "No cuts ";
1837  cutFlowVector_str[1] = "3LLOW: Preselection ";
1838  cutFlowVector_str[2] = "3LLOW: 75 GeV < mll < 105 GeV ";
1839  cutFlowVector_str[3] = "3LLOW: mTW > 100 GeV ";
1840  cutFlowVector_str[4] = "3LLOW: m_HT4PP/m_H4PP > 0.9 ";
1841  cutFlowVector_str[5] = "3LLOW: m_H4PP > 250 GeV ";
1842  cutFlowVector_str[6] = "3LLOW: pT_PP/(pT_PP + HT_PP(3,1)) ";
1843  cutFlowVector_str[7] = "2L2JLOW: Preselection ";
1844  cutFlowVector_str[8] = "2L2JLOW: mll ";
1845  cutFlowVector_str[9] = "2L2JLOW: mjj ";
1846  cutFlowVector_str[10] = "2L2JLOW: HT_PP(1,1)/HT_PP(4,1) ";
1847  cutFlowVector_str[11] = "2L2JLOW: pT_PP/(pT_PP + HT_PP(4,1)) ";
1848  cutFlowVector_str[12] = "2L2JLOW: minDPhi ";
1849  cutFlowVector_str[13] = "2L2JLOW: HPP(4,1) ";
1850  cutFlowVector_str[14] = "2L2JINT: Preselection ";
1851  cutFlowVector_str[15] = "2L2JINT: mll ";
1852  cutFlowVector_str[16] = "2L2JINT: mjj ";
1853  cutFlowVector_str[17] = "2L2JINT: HT_PP(1,1)/HT_PP(4,1) ";
1854  cutFlowVector_str[18] = "2L2JINT: pT_PP/(pT_PP + HT_PP(4,1)) ";
1855  cutFlowVector_str[19] = "2L2JINT: minDPhi ";
1856  cutFlowVector_str[20] = "2L2JINT: HPP(4,1) ";
1857  cutFlowVector_str[21] = "2L2JHIGH: Preselection ";
1858  cutFlowVector_str[22] = "2L2JHIGH: mll ";
1859  cutFlowVector_str[23] = "2L2JHIGH: mjj ";
1860  cutFlowVector_str[24] = "2L2JHIGH: m_R_minH2P_minH3P>0.8";
1861  cutFlowVector_str[25] = "2L2JHIGH: m_RPT_HT5PP < 0.05 ";
1862  cutFlowVector_str[26] = "2L2JHIGH: 0.3 < minDPhiVP > 2.8 ";
1863  cutFlowVector_str[27] = "2L2JHIGH: m_H5PP>800. ";
1864  cutFlowVector_str[28] = "2L2JCOMP: Preselection ";
1865  cutFlowVector_str[29] = "2L2JCOMP: mZ ";
1866  cutFlowVector_str[30] = "2L2JCOMP: mJ ";
1867  cutFlowVector_str[31] = "2L2JCOMP: dPhi_ISR_I ";
1868  cutFlowVector_str[32] = "2L2JCOMP: R_ISR ";
1869  cutFlowVector_str[33] = "2L2JCOMP: p_ISRT ";
1870  cutFlowVector_str[34] = "2L2JCOMP: p_IT ";
1871  cutFlowVector_str[35] = "2L2JCOMP: pT_CM ";
1872  cutFlowVector_str[36] = "3LHIGH: Preselection ";
1873  cutFlowVector_str[37] = "3LHIGH: mll ";
1874  cutFlowVector_str[38] = "3LHIGH: mTW ";
1875  cutFlowVector_str[39] = "3LHIGH: m_HT4PP/m_H4PP ";
1876  cutFlowVector_str[40] = "3LHIGH: HPb(1,1)/HPb(2,1) ";
1877  cutFlowVector_str[41] = "3LHIGH: m_H4PP ";
1878  cutFlowVector_str[42] = "3LHIGH: pT_PP/(pT_PP + HT_PP(3,1)) ";
1879  cutFlowVector_str[43] = "3LCOMP: Preselection ";
1880  cutFlowVector_str[44] = "3LCOMP: mll ";
1881  cutFlowVector_str[45] = "3LCOMP: mTW ";
1882  cutFlowVector_str[46] = "3LCOMP: dPhi_ISRI ";
1883  cutFlowVector_str[47] = "3LCOMP: R_ISR ";
1884  cutFlowVector_str[48] = "3LCOMP: p_ISRT ";
1885  cutFlowVector_str[49] = "3LCOMP: p_IT ";
1886  cutFlowVector_str[50] = "3LCOMP: pT_CM ";
1887  cutFlowVector_str[51] = "3LINT: Preselection ";
1888  cutFlowVector_str[52] = "3LINT: mll ";
1889  cutFlowVector_str[53] = "3LINT: mTW ";
1890  cutFlowVector_str[54] = "3LINT: m_HT4PP/m_H4PP ";
1891  cutFlowVector_str[55] = "3LINT: HPb(1,1)/HPb(2,1) ";
1892  cutFlowVector_str[56] = "3LINT: m_H4PP ";
1893  cutFlowVector_str[57] = "3LINT: pT_PP/(pT_PP + HT_PP(3,1)) ";
1894 
1895  //std::cout << " m_is3Lep " << m_is3Lep << " m_is2Lep2Jet " << m_is2Lep2Jet << " m_is2L2JInt " << m_is2L2JInt << " m_is3LInt " << m_is3LInt << " m_is3Lep2Jet " << m_is3Lep2Jet << " m_is3Lep3Jet " << m_is3Lep3Jet << " m_is4Lep2Jet " << m_is4Lep2Jet << " m_is4Lep3Jet " << m_is4Lep3Jet << std::endl;
1896 
1897  //if(m_is3Lep)std::cout << " m_is3Lep " << m_is3Lep << " m_lept1sign " << m_lept1sign << " m_lept2sign " << m_lept2sign << " m_lept1Pt " << m_lept1Pt << " m_lept2Pt " << m_lept2Pt << " m_lept3Pt " << m_lept3Pt << " signalBJets.size() " << signalBJets.size() << " signalJets.size() " << signalJets.size() << std::endl;
1898 
1899  for(int j=0;j<NCUTS;j++){
1900 
1901  if( (j==0) ||
1902 
1903  (j==1 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>40. && m_lept3Pt>30. && signalBJets.size()==0 && signalJets.size()==0) ||
1904 
1905  (j==2 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>40. && m_lept3Pt>30. && signalBJets.size()==0 && signalJets.size()==0 && m_mll>75. && m_mll<105.) ||
1906 
1907  (j==3 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>40. && m_lept3Pt>30. && signalBJets.size()==0 && signalJets.size()==0 && m_mll>75. && m_mll<105. && m_mTW>100.) ||
1908 
1909  (j==4 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>40. && m_lept3Pt>30. && signalBJets.size()==0 && signalJets.size()==0 && m_mll>75. && m_mll<105. && m_mTW>100. && m_HT4PP/m_H4PP > 0.9) ||
1910 
1911  (j==5 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>40. && m_lept3Pt>30. && signalBJets.size()==0 && signalJets.size()==0 && m_mll>75. && m_mll<105. && m_mTW>100. && m_HT4PP/m_H4PP > 0.9 && m_H4PP > 250.) ||
1912 
1913  (j==6 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>40. && m_lept3Pt>30. && signalBJets.size()==0 && signalJets.size()==0 && m_mll>75. && m_mll<105. && m_mTW>100. && m_HT4PP/m_H4PP > 0.9 && m_H4PP > 250. && m_RPT_HT4PP < 0.05) ||
1914 
1915  /*cutFlowVector_str[7] = "2L2JLOW: Preselection ";
1916  cutFlowVector_str[8] = "2L2JLOW: 80 GeV < mll < 100 GeV";
1917  cutFlowVector_str[9] = "2L2JLOW: 70 GeV < mjj < 90 GeV ";
1918  cutFlowVector_str[10] = "2L2JLOW: HT_PP(1,1)/HT_PP(4,1) ";
1919  cutFlowVector_str[11] = "2L2JLOW: pT_PP/(pT_PP + HT_PP(4,1)) ";
1920  cutFlowVector_str[12] = "2L2JLOW: minDPhi ";
1921  cutFlowVector_str[13] = "2L2JLOW: HPP(4,1) ";*/
1922 
1923  (j==7 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()==2) ||
1924 
1925  (j==8 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()==2 && m_mll>80. && m_mll<100.) ||
1926 
1927  (j==9 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()==2 && m_mll>80. && m_mll<100. && m_mjj>70. && m_mjj<90.) ||
1928 
1929  (j==10 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()==2 && m_mll>80. && m_mll<100. && m_mjj>70. && m_mjj<90. && m_H2PP/m_H5PP>0.35 && m_H2PP/m_H5PP<0.6) ||
1930 
1931  (j==11 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()==2 && m_mll>80. && m_mll<100. && m_mjj>70. && m_mjj<90. && m_H2PP/m_H5PP>0.35 && m_H2PP/m_H5PP<0.6 && m_RPT_HT5PP<0.05) ||
1932 
1933  (j==12 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()==2 && m_mll>80. && m_mll<100. && m_mjj>70. && m_mjj<90. && m_H2PP/m_H5PP>0.35 && m_H2PP/m_H5PP<0.6 && m_RPT_HT5PP<0.05 && m_minDphi>2.4) ||
1934 
1935  (j==13 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()==2 && m_mll>80. && m_mll<100. && m_mjj>70. && m_mjj<90. && m_H2PP/m_H5PP>0.35 && m_H2PP/m_H5PP<0.6 && m_RPT_HT5PP<0.05 && m_minDphi>2.4 && m_H5PP>400.) ||
1936 
1937  /*cutFlowVector_str[14] = "2L2JINT: Preselection ";
1938  cutFlowVector_str[15] = "2L2JINT: mll ";
1939  cutFlowVector_str[16] = "2L2JINT: mjj ";
1940  cutFlowVector_str[17] = "2L2JINT: HT_PP(1,1)/HT_PP(4,1) ";
1941  cutFlowVector_str[18] = "2L2JINT: pT_PP/(pT_PP + HT_PP(4,1)) ";
1942  cutFlowVector_str[19] = "2L2JINT: minDPhi ";
1943  cutFlowVector_str[20] = "2L2JINT: HPP(4,1) ";*/
1944 
1945  (j==14 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()>=2) ||
1946 
1947  (j==15 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()>=2 && m_mll>80. && m_mll<100.) ||
1948 
1949  (j==16 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()>=2 && m_mll>80. && m_mll<100. && m_mjj>60. && m_mjj<100.) ||
1950 
1951  (j==17 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()>=2 && m_mll>80. && m_mll<100. && m_mjj>60. && m_mjj<100. && m_R_minH2P_minH3P>0.8 ) ||
1952 
1953  (j==18 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()>=2 && m_mll>80. && m_mll<100. && m_mjj>60. && m_mjj<100. && m_R_minH2P_minH3P>0.8 && m_RPT_HT5PP<0.05) ||
1954 
1955  (j==19 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()>=2 && m_mll>80. && m_mll<100. && m_mjj>60. && m_mjj<100. && m_R_minH2P_minH3P>0.8 && m_RPT_HT5PP<0.05 && m_dphiVP>0.6) ||
1956 
1957  (j==20 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()>=2 && m_mll>80. && m_mll<100. && m_mjj>60. && m_mjj<100. && m_R_minH2P_minH3P>0.8 && m_RPT_HT5PP<0.05 && m_dphiVP>0.6 && m_H5PP>600.) ||
1958 
1959  /* cutFlowVector_str[21] = "2L2JHIGH: Preselection ";
1960  cutFlowVector_str[22] = "2L2JHIGH: mll ";
1961  cutFlowVector_str[23] = "2L2JHIGH: mjj ";
1962  cutFlowVector_str[24] = "2L2JHIGH: m_R_minH2P_minH3P>0.8";
1963  cutFlowVector_str[25] = "2L2JHIGH: m_RPT_HT5PP < 0.05 ";
1964  cutFlowVector_str[26] = "2L2JHIGH: 0.3 < minDPhiVP > 2.8 ";
1965  cutFlowVector_str[27] = "2L2JHIGH: m_H5PP>800. ";*/
1966 
1967  (j==21 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()>=2) ||
1968 
1969  (j==22 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()>=2 && m_mll>80. && m_mll<100.) ||
1970 
1971  (j==23 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()>=2 && m_mll>80. && m_mll<100. && m_mjj>60. && m_mjj<100.) ||
1972 
1973  (j==24 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()>=2 && m_mll>80. && m_mll<100. && m_mjj>60. && m_mjj<100. && m_R_minH2P_minH3P>0.8) ||
1974 
1975  (j==25 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()>=2 && m_mll>80. && m_mll<100. && m_mjj>60. && m_mjj<100. && m_R_minH2P_minH3P>0.8 && m_RPT_HT5PP< 0.05) ||
1976 
1977  (j==26 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()>=2 && m_mll>80. && m_mll<100. && m_mjj>60. && m_mjj<100. && m_R_minH2P_minH3P>0.8 && m_RPT_HT5PP< 0.05 && m_dphiVP>0.3 && m_dphiVP<2.8 ) ||
1978 
1979  (j==27 && m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()>=2 && m_mll>80. && m_mll<100. && m_mjj>60. && m_mjj<100. && m_R_minH2P_minH3P>0.8 && m_RPT_HT5PP< 0.05 && m_dphiVP>0.3 && m_dphiVP<2.8 && m_H5PP>800.) ||
1980 
1981  /*cutFlowVector_str[28] = "2L2JCOMP: Preselection ";
1982  cutFlowVector_str[29] = "2L2JCOMP: mZ ";
1983  cutFlowVector_str[30] = "2L2JCOMP: mJ ";
1984  cutFlowVector_str[31] = "2L2JCOMP: dPhi_ISR_I ";
1985  cutFlowVector_str[32] = "2L2JCOMP: R_ISR ";
1986  cutFlowVector_str[33] = "2L2JCOMP: p_ISRT ";
1987  cutFlowVector_str[34] = "2L2JCOMP: p_IT ";
1988  cutFlowVector_str[35] = "2L2JCOMP: pT_CM ";*/
1989 
1990  (j==28 && m_is2L2JInt && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && m_NjS==2 && m_NjISR>0) ||
1991 
1992  (j==29 && m_is2L2JInt && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && m_NjS==2 && m_NjISR>0 && m_MZ>80. && m_MZ<100.) ||
1993 
1994  (j==30 && m_is2L2JInt && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && m_NjS==2 && m_NjISR>0 && m_MZ>80. && m_MZ<100. && m_MJ>50. && m_MJ<110.) ||
1995 
1996  (j==31 && m_is2L2JInt && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && m_NjS==2 && m_NjISR>0 && m_MZ>80. && m_MZ<100. && m_MJ>50. && m_MJ<110. && m_dphiISRI>2.8) ||
1997 
1998  (j==32 && m_is2L2JInt && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && m_NjS==2 && m_NjISR>0 && m_MZ>80. && m_MZ<100. && m_MJ>50. && m_MJ<110. && m_dphiISRI>2.8 && m_RISR > 0.40 && m_RISR < 0.75 ) ||
1999 
2000  (j==33 && m_is2L2JInt && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && m_NjS==2 && m_NjISR>0 && m_MZ>80. && m_MZ<100. && m_MJ>50. && m_MJ<110. && m_dphiISRI>2.8 && m_RISR > 0.40 && m_RISR < 0.75 && m_PTISR > 180.) ||
2001 
2002  (j==34 && m_is2L2JInt && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && m_NjS==2 && m_NjISR>0 && m_MZ>80. && m_MZ<100. && m_MJ>50. && m_MJ<110. && m_dphiISRI>2.8 && m_RISR > 0.40 && m_RISR < 0.75 && m_PTISR > 180. && m_PTI > 100.) ||
2003 
2004  (j==35 && m_is2L2JInt && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && m_NjS==2 && m_NjISR>0 && m_MZ>80. && m_MZ<100. && m_MJ>50. && m_MJ<110. && m_dphiISRI>2.8 && m_RISR > 0.40 && m_RISR < 0.75 && m_PTISR > 180. && m_PTI > 100. && m_PTCM < 20.) ||
2005 
2006 
2007  /* cutFlowVector_str[36] = "3LHIGH: Preselection ";
2008  cutFlowVector_str[37] = "3LHIGH: 75 GeV < mll < 105 GeV ";
2009  cutFlowVector_str[38] = "3LHIGH: mTW ";
2010  cutFlowVector_str[39] = "3LHIGH: m_HT4PP/m_H4PP ";
2011  cutFlowVector_str[40] = "3LHIGH: HPb(1,1)/HPb(2,1) ";
2012  cutFlowVector_str[41] = "3LHIGH: m_H4PP ";
2013  cutFlowVector_str[42] = "3LHIGH: pT_PP/(pT_PP + HT_PP(3,1)) ";*/
2014 
2015 
2016  (j==36 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>60. && m_lept3Pt>40. && signalBJets.size()==0 && signalJets.size()<3) ||
2017 
2018  (j==37 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>60. && m_lept3Pt>40. && signalBJets.size()==0 && signalJets.size()<3 && m_mll>75. && m_mll<105.) ||
2019 
2020  (j==38 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>60. && m_lept3Pt>40. && signalBJets.size()==0 && signalJets.size()<3 && m_mll>75. && m_mll<105. && m_mTW>150.) ||
2021 
2022  (j==39 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>60. && m_lept3Pt>40. && signalBJets.size()==0 && signalJets.size()<3 && m_mll>75. && m_mll<105. && m_mTW>150. && m_HT4PP/m_H4PP > 0.75) ||
2023 
2024  (j==40 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>60. && m_lept3Pt>40. && signalBJets.size()==0 && signalJets.size()<3 && m_mll>75. && m_mll<105. && m_mTW>150. && m_HT4PP/m_H4PP > 0.75 && m_R_minH2P_minH3P>0.8) ||
2025 
2026  (j==41 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>60. && m_lept3Pt>40. && signalBJets.size()==0 && signalJets.size()<3 && m_mll>75. && m_mll<105. && m_mTW>150. && m_R_minH2P_minH3P>0.8 && m_HT4PP/m_H4PP > 0.75 && m_H4PP > 550. ) ||
2027 
2028  (j==42 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>60. && m_lept3Pt>40. && signalBJets.size()==0 && signalJets.size()<3 && m_mll>75. && m_mll<105. && m_mTW>150. && m_R_minH2P_minH3P>0.8 && m_HT4PP/m_H4PP > 0.75 && m_H4PP > 550. && m_RPT_HT4PP < 0.2) ||
2029 
2030  /*cutFlowVector_str[43] = "3LCOMP: Preselection ";
2031  cutFlowVector_str[44] = "3LCOMP: mll ";
2032  cutFlowVector_str[45] = "3LCOMP: mTW ";
2033  cutFlowVector_str[46] = "3LCOMP: dPhi_ISRI ";
2034  cutFlowVector_str[47] = "3LCOMP: R_ISR ";
2035  cutFlowVector_str[48] = "3LCOMP: p_ISRT ";
2036  cutFlowVector_str[49] = "3LCOMP: p_IT ";
2037  cutFlowVector_str[50] = "3LCOMP: pT_CM ";*/
2038 
2039  (j==43 && m_is3LInt && ((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign))) && m_lept1Pt>25. && m_lept2Pt>25. && m_lept3Pt>20. && signalBJets.size()==0 && signalJets.size()<4) ||
2040 
2041  (j==44 && m_is3LInt && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>25. && m_lept2Pt>25. && m_lept3Pt>20. && signalBJets.size()==0 && signalJets.size()<4 && m_mll>75. && m_mll<105.) ||
2042 
2043  (j==45 && m_is3LInt && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>25. && m_lept2Pt>25. && m_lept3Pt>20. && signalBJets.size()==0 && signalJets.size()<4 && m_mll>75. && m_mll<105. && m_mTW>100.) ||
2044 
2045  (j==46 && m_is3LInt && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>25. && m_lept2Pt>25. && m_lept3Pt>20. && signalBJets.size()==0 && signalJets.size()<4 && m_mll>75. && m_mll<105. && m_mTW>100. && m_dphiISRI>2.0) ||
2046 
2047  (j==47 && m_is3LInt && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>25. && m_lept2Pt>25. && m_lept3Pt>20. && signalBJets.size()==0 && signalJets.size()<4 && m_mll>75. && m_mll<105. && m_mTW>100. && m_dphiISRI>2.0 && m_RISR>0.55 && m_RISR<1.0) ||
2048 
2049  (j==48 && m_is3LInt && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>25. && m_lept2Pt>25. && m_lept3Pt>20. && signalBJets.size()==0 && signalJets.size()<4 && m_mll>75. && m_mll<105. && m_mTW>100. && m_dphiISRI>2.0 && m_RISR>0.55 && m_RISR<1.0 && m_PTISR>100.) ||
2050 
2051  (j==49 && m_is3LInt && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>25. && m_lept2Pt>25. && m_lept3Pt>20. && signalBJets.size()==0 && signalJets.size()<4 && m_mll>75. && m_mll<105. && m_mTW>100. && m_dphiISRI>2.0 && m_RISR>0.55 && m_RISR<1.0 && m_PTISR>100. && m_PTI>80.) ||
2052 
2053  (j==50 && m_is3LInt && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>25. && m_lept2Pt>25. && m_lept3Pt>20. && signalBJets.size()==0 && signalJets.size()<4 && m_mll>75. && m_mll<105. && m_mTW>100. && m_dphiISRI>2.0 && m_RISR>0.55 && m_RISR<1.0 && m_PTISR>100. && m_PTI>80. && m_PTCM<25.) ||
2054 
2055  /* cutFlowVector_str[51] = "3LINT: Preselection ";
2056  cutFlowVector_str[52] = "3LINT: mll ";
2057  cutFlowVector_str[53] = "3LINT: mTW ";
2058  cutFlowVector_str[54] = "3LINT: m_HT4PP/m_H4PP ";
2059  cutFlowVector_str[55] = "3LINT: HPb(1,1)/HPb(2,1) ";
2060  cutFlowVector_str[56] = "3LINT: m_H4PP ";
2061  cutFlowVector_str[57] = "3LINT: pT_PP/(pT_PP + HT_PP(3,1)) ";*/
2062 
2063  (j==51 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>50. && m_lept3Pt>30. && signalBJets.size()==0 && signalJets.size()<3) ||
2064 
2065  (j==52 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>50. && m_lept3Pt>30. && signalBJets.size()==0 && signalJets.size()<3 && m_mll>75. && m_mll<105.) ||
2066 
2067  (j==53 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>50. && m_lept3Pt>30. && signalBJets.size()==0 && signalJets.size()<3 && m_mll>75. && m_mll<105. && m_mTW>130.) ||
2068 
2069  (j==54 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>50. && m_lept3Pt>30. && signalBJets.size()==0 && signalJets.size()<3 && m_mll>75. && m_mll<105. && m_mTW>130. && m_HT4PP/m_H4PP > 0.8 ) ||
2070 
2071  (j==55 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>50. && m_lept3Pt>30. && signalBJets.size()==0 && signalJets.size()<3 && m_mll>75. && m_mll<105. && m_mTW>130. && m_HT4PP/m_H4PP > 0.8 && m_R_minH2P_minH3P>0.75) ||
2072 
2073  (j==56 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>50. && m_lept3Pt>30. && signalBJets.size()==0 && signalJets.size()<3 && m_mll>75. && m_mll<105. && m_mTW>130. && m_HT4PP/m_H4PP > 0.8 && m_R_minH2P_minH3P>0.75 && m_H4PP > 450. ) ||
2074 
2075  (j==57 && m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>50. && m_lept3Pt>30. && signalBJets.size()==0 && signalJets.size()<3 && m_mll>75. && m_mll<105. && m_mTW>130. && m_HT4PP/m_H4PP > 0.8 && m_R_minH2P_minH3P>0.75 && m_H4PP > 450. && m_RPT_HT4PP < 0.15)
2076 
2077  )cutFlowVector[j]++;
2078  }
2079 
2080  // Now apply the signal region cuts
2081 
2082  if(m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()>=2 && m_mll>80. && m_mll<100. && m_mjj>60. && m_mjj<100. && m_R_minH2P_minH3P>0.8 && m_RPT_HT5PP< 0.05 && m_dphiVP>0.3 && m_dphiVP<2.8 && m_H5PP>800.) _counters.at("2L2JHIGH").add_event(event);
2083 
2084  if(m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()>=2 && m_mll>80. && m_mll<100. && m_mjj>60. && m_mjj<100. && m_R_minH2P_minH3P>0.8 && m_RPT_HT5PP<0.05 && m_dphiVP>0.6 && m_dphiVP<2.6 && m_H5PP>600.) _counters.at("2L2JINT").add_event(event);
2085 
2086 
2087  if(m_is2Lep2Jet && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && signalJets.size()==2 && m_mll>80. && m_mll<100. && m_mjj>70. && m_mjj<90. && m_H2PP/m_H5PP>0.35 && m_H2PP/m_H5PP<0.6 && m_RPT_HT5PP<0.05 && m_minDphi>2.4 && m_H5PP>400.) _counters.at("2L2JLOW").add_event(event);
2088 
2089  if(m_is2L2JInt && m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign) && m_lept1Pt>25. && m_lept2Pt>25. && m_jet1Pt>30. && m_jet2Pt>30. && signalBJets.size()==0 && m_NjS==2 && m_NjISR>0 && m_MZ>80. && m_MZ<100. && m_MJ>50. && m_MJ<110. && m_dphiISRI>2.8 && m_RISR > 0.40 && m_RISR < 0.75 && m_PTISR > 180. && m_PTI > 100. && m_PTCM < 20.) _counters.at("2L2JCOMP").add_event(event);
2090 
2091  if(m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>60. && m_lept3Pt>40. && signalBJets.size()==0 && signalJets.size()<3 && m_mll>75. && m_mll<105. && m_mTW>150. && m_HT4PP/m_H4PP > 0.75 && m_R_minH2P_minH3P>0.8 && m_H4PP > 550. && m_RPT_HT4PP < 0.2) _counters.at("3LHIGH").add_event(event);
2092 
2093  if(m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>50. && m_lept3Pt>30. && signalBJets.size()==0 && signalJets.size()<3 && m_mll>75. && m_mll<105. && m_mTW>130. && m_HT4PP/m_H4PP > 0.8 && m_R_minH2P_minH3P>0.75 && m_H4PP > 450. && m_RPT_HT4PP < 0.15) _counters.at("3LINT").add_event(event);
2094 
2095  if(m_is3LInt && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>25. && m_lept2Pt>25. && m_lept3Pt>20. && signalBJets.size()==0 && signalJets.size()<4 && m_mll>75. && m_mll<105. && m_mTW>100. && m_dphiISRI>2.0 && m_RISR>0.55 && m_RISR<1.0 && m_PTISR>100. && m_PTI>80. && m_PTCM<25.) _counters.at("3LCOMP").add_event(event);
2096 
2097  if(m_is3Lep && (((m_lept1sign*m_lept2sign<0 && abs(m_lept1sign)==abs(m_lept2sign)) || (m_lept1sign*m_lept3sign<0 && abs(m_lept1sign)==abs(m_lept3sign)) || (m_lept2sign*m_lept3sign<0 && abs(m_lept2sign)==abs(m_lept3sign)))) && m_lept1Pt>60. && m_lept2Pt>40. && m_lept3Pt>30. && signalBJets.size()==0 && signalJets.size()==0 && m_mll>75. && m_mll<105. && m_mTW>100. && m_HT4PP/m_H4PP > 0.9 && m_H4PP > 250. && m_RPT_HT4PP < 0.05) _counters.at("3LLOW").add_event(event);
2098 
2099  return;
2100 
2101  } // end analyze method
2102 
2104  void combine(const Analysis* other)
2105  {
2106  const Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb* specificOther
2107  = dynamic_cast<const Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb*>(other);
2108 
2109  for (auto& pair : _counters) { pair.second += specificOther->_counters.at(pair.first); }
2110 
2111  if (NCUTS != specificOther->NCUTS) NCUTS = specificOther->NCUTS;
2112 
2113  for (int j=0; j<NCUTS; j++)
2114  {
2115  cutFlowVector[j] += specificOther->cutFlowVector[j];
2116  cutFlowVector_str[j] = specificOther->cutFlowVector_str[j];
2117  }
2118 
2119  }
2120 
2121 
2122  virtual void collect_results() {
2123 /*
2124  double scale_by=1.;
2125  cout << "------------------------------------------------------------------------------------------------------------------------------ "<<endl;
2126  cout << "CUT FLOW: ATLAS 13 TeV 3 lep low mass RJ signal region "<<endl;
2127  cout << "------------------------------------------------------------------------------------------------------------------------------"<<endl;
2128  cout << right << setw(40) << "CUT" << setw(20) << "RAW" << setw(20) << "SCALED"
2129  << setw(20) << "%" << setw(20) << "clean adj RAW"<< setw(20) << "clean adj %" << endl;
2130  for (int j=0; j<NCUTS; j++) {
2131  cout << right << setw(40) << cutFlowVector_str[j].c_str() << setw(20)
2132  << cutFlowVector[j] << setw(20) << cutFlowVector[j]*scale_by << setw(20)
2133  << 100.*cutFlowVector[j]/cutFlowVector[0] << "%" << setw(20)
2134  << cutFlowVector[j]*scale_by << setw(20) << 100.*cutFlowVector[j]/cutFlowVector[0]<< "%" << endl;
2135  }
2136  cout << "------------------------------------------------------------------------------------------------------------------------------ "<<endl;
2137 */
2138 
2139  add_result(SignalRegionData(_counters.at("2L2JHIGH"), 0, {1.9, 0.8}));
2140  add_result(SignalRegionData(_counters.at("2L2JINT"), 1, {2.4, 0.9}));
2141  add_result(SignalRegionData(_counters.at("2L2JLOW"), 19, {8.4, 5.8}));
2142  add_result(SignalRegionData(_counters.at("2L2JCOMP"), 11, {2.7, 2.7}));
2143  add_result(SignalRegionData(_counters.at("3LHIGH"), 2, {1.1, 0.5}));
2144  add_result(SignalRegionData(_counters.at("3LINT"), 1, {2.3, 0.5}));
2145  add_result(SignalRegionData(_counters.at("3LLOW"), 20, {10., 2.0}));
2146  add_result(SignalRegionData(_counters.at("3LCOMP"), 12, {3.9, 1.0}));
2147 
2148  return;
2149  }
2150 
2151 
2152  protected:
2153  void analysis_specific_reset() {
2154 
2155  for (auto& pair : _counters) { pair.second.reset(); }
2156 
2157  std::fill(cutFlowVector.begin(), cutFlowVector.end(), 0);
2158  }
2159 
2160  }; // end class Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb
2161 
2162  DEFINE_ANALYSIS_FACTORY(ATLAS_13TeV_RJ3L_lowmass_36invfb)
2163 
2164 
2165  //
2166  // Derived analysis class for the RJ3L_lowmass SRs
2167  //
2168  class Analysis_ATLAS_13TeV_RJ3L_2Lep2Jets_36invfb : public Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb {
2169  public:
2170  Analysis_ATLAS_13TeV_RJ3L_2Lep2Jets_36invfb() {
2171  set_analysis_name("ATLAS_13TeV_RJ3L_2Lep2Jets_36invfb");
2172  }
2173  virtual void collect_results() {
2174  add_result(SignalRegionData(_counters.at("2L2JHIGH"), 0, {1.9, 0.8}));
2175  add_result(SignalRegionData(_counters.at("2L2JINT"), 1, {2.4, 0.9}));
2176  add_result(SignalRegionData(_counters.at("2L2JLOW"), 19, {8.4, 5.8}));
2177  add_result(SignalRegionData(_counters.at("2L2JCOMP"), 11, {2.7, 2.7}));
2178  }
2179 
2180  };
2181  // Factory fn
2182  DEFINE_ANALYSIS_FACTORY(ATLAS_13TeV_RJ3L_2Lep2Jets_36invfb)
2183 
2184  class Analysis_ATLAS_13TeV_RJ3L_3Lep_36invfb : public Analysis_ATLAS_13TeV_RJ3L_lowmass_36invfb {
2185  public:
2186  Analysis_ATLAS_13TeV_RJ3L_3Lep_36invfb() {
2187  set_analysis_name("ATLAS_13TeV_RJ3L_3Lep_36invfb");
2188  }
2189  virtual void collect_results() {
2190  add_result(SignalRegionData(_counters.at("3LHIGH"), 2, {1.1, 0.5}));
2191  add_result(SignalRegionData(_counters.at("3LINT"), 1, {2.3, 0.5}));
2192  add_result(SignalRegionData(_counters.at("3LLOW"), 20, {10., 2.0}));
2193  add_result(SignalRegionData(_counters.at("3LCOMP"), 12, {3.9, 1.0}));
2194  }
2195 
2196  };
2197  // Factory fn
2198  DEFINE_ANALYSIS_FACTORY(ATLAS_13TeV_RJ3L_3Lep_36invfb)
2199 
2200 
2201  } // end namespace ColliderBit
2202 } // end namespace Gambit
2203 
2204 #endif
2205 #endif
Piped_exceptions piped_errors
Global instance of Piped_exceptions class for errors.
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 sortLepByPT_RJ3L(const HEPUtils::Particle *lep1, const HEPUtils::Particle *lep2)
void request(std::string origin, std::string message)
Request an exception.
Definition: exceptions.cpp:547
#define LOCAL_INFO
Definition: local_info.hpp:34
Piped_exceptions piped_warnings
Global instance of Piped_exceptions class for warnings.
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 *)
STL namespace.
START_MODEL b
Definition: demo.hpp:270
bool SortLeptons(const pair< TLorentzVector, int > lv1, const pair< TLorentzVector, int > lv2)
A simple class for counting events of type HEPUtils::Event.
A class for collider analyses within ColliderBit.
Definition: Analysis.hpp:41
std::string str
Shorthand for a standard string.
Definition: Analysis.hpp:35
A simple container for the result of one signal region from one analysis.
void applyMuonEff(std::vector< const HEPUtils::Particle *> &muons)
Randomly filter the supplied particle list by parameterised muon efficiency.
DS5_MSPCTM DS_INTDOF int
START_MODEL M
bool sortByPT_RJ3L(const HEPUtils::Jet *jet1, const HEPUtils::Jet *jet2)
def combine(region_file, pip_file)
Definition: colouring.py:169
bool SortJets(const TLorentzVector jv1, const TLorentzVector jv2)
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
Class for ColliderBit analyses.
Functions that do super fast ATLAS detector simulation based on four-vector smearing.