checkmate is hosted by Hepforge, IPPP Durham
AnalysisBase.cc
Go to the documentation of this file.
1 #include "AnalysisBase.h"
2 
3 AnalysisBase::AnalysisBase(std::string inFile, std::string outFol, std::string outPre, double xs, double xserr, std::map<std::string, std::string> branches, std::map<std::string, std::vector<int> > flags) {
4  // Setting internal parameters
5  outputFolder = outFol;
6  outputPrefix = outPre;
7  inputFile = inFile;
8  xsect = xs;
9  xsecterr = xserr;
10  luminosity = 0;
11  ignoreElectrons = ignoreElectrons_LooseIsolation = ignoreElectrons_MediumEfficiency = ignoreElectrons_TightEfficiency = ignoreMuons = ignoreMuons_LooseIsolation = ignoreMuons_CombinedPlusEfficiency = ignoreMuons_CombinedEfficiency = ignorePhotons = ignorePhotons_LooseIsolation = ignoreTracks = ignoreTowers = false;
12 
13  // Reading the Delphes Tree
14  chain = new TChain("Delphes");
15  chain->Add(inFile.c_str());
16  treeReader = new ExRootTreeReader(chain);
17  result = new ExRootResult();
18  branchEvent = treeReader->UseBranch(branches["Event"].c_str());
19  branchTrack = treeReader->UseBranch(branches["Track"].c_str());
20  branchTower = treeReader->UseBranch(branches["Tower"].c_str());
21  branchJet = treeReader->UseBranch(branches["Jet"].c_str());
22  branchJet2 = treeReader->UseBranch(branches["Jet2"].c_str());
23  branchElectron = treeReader->UseBranch(branches["Electron"].c_str());
24  branchPhoton = treeReader->UseBranch(branches["Photon"].c_str());
25  branchMuon = treeReader->UseBranch(branches["Muon"].c_str());
26  branchMissingET = treeReader->UseBranch(branches["MissingET"].c_str());
27 
28  // Saving which flagvalues correspond to this analysis (given by CheckMATE)
29  electronIsolationFlags = flags["electron_isolation"];
30  muonIsolationFlags = flags["muon_isolation"];
31  photonIsolationFlags = flags["photon_isolation"];
32  jetBTagFlags = flags["jet_btags"];
33  jet2BTagFlags = flags["jet2_btags"];
34 
35  // Setting the random number generator
36  if (!flags["randomseed"].empty())
37  srand(flags["randomseed"][0]);
38  else
39  srand(time(0));
40 
41  std::cout << rand() << std::endl;
42 }
43 
45  // Free all pointers
46  delete chain;
47  delete treeReader;
48  delete result;
49 
50  for (int i = 0; i < fStreams.size(); i++) {
51  fStreams[i]->close();
52  delete fStreams[i];
53  }
54 }
55 
56 void AnalysisBase::bookSignalRegions(std::string listOfRegions) {
57  std::string currKey = "";
58  // Sum letter by letter and define map-key as soon as ; is reached
59  for(int i = 0; i < strlen(listOfRegions.c_str()); i++) {
60  char c = listOfRegions[i];
61  if (c == ';') {
62  signalRegions[currKey] = signalRegions2[currKey] = 0.0;
63  currKey = "";
64  }
65  else
66  currKey += c;
67  }
68  // The last key might not be separated by ;
69  if(currKey != "")
70  signalRegions[currKey] = signalRegions2[currKey] = 0.0;
71 }
72 
73 void AnalysisBase::bookControlRegions(std::string listOfRegions) {
74  std::string currKey = "";
75  // Sum letter by letter and define map-key as soon as ; is reached
76  for(int i = 0; i < strlen(listOfRegions.c_str()); i++) {
77  char c = listOfRegions[i];
78  if (c == ';') {
79  controlRegions[currKey] = controlRegions2[currKey] = 0.0;
80  currKey = "";
81  }
82  else
83  currKey += c;
84  }
85  // The last key might not be separated by ;
86  if(currKey != "")
87  controlRegions[currKey] = controlRegions2[currKey] = 0.0;
88 }
89 void AnalysisBase::bookCutflowRegions(std::string listOfRegions) {
90  std::string currKey = "";
91  // Sum letter by letter and define map-key as soon as ; is reached
92  for(int i = 0; i < strlen(listOfRegions.c_str()); i++) {
93  char c = listOfRegions[i];
94  if (c == ';') {
95  cutflowRegions[currKey] = cutflowRegions2[currKey] = 0.0;
96  currKey = "";
97  }
98  else
99  currKey += c;
100  }
101  // The last key might not be separated by ;
102  if(currKey != "")
103  cutflowRegions[currKey] = cutflowRegions2[currKey] = 0.0;
104 }
105 
106 
108  // These const vectors are needed for efficiency/isolation filters
109  std::vector<int> looseIsolationFlags;
110  looseIsolationFlags.push_back(0);
111  std::vector<int> mediumEfficiencyFlags;
112  mediumEfficiencyFlags.push_back(0);
113  std::vector<int> tightEfficiencyFlags;
114  tightEfficiencyFlags.push_back(1);
115 
116  initialize();
117  sumOfWeights = 0;
118  nEvents = treeReader->GetEntries();
119  double missingET_ET = 0, missingET_Phi = 0, missingET_Ex = 0, missingET_Ey = 0;
120  for(Int_t entry = 0; entry < nEvents; ++entry) {
121  treeReader->ReadEntry(entry);
122  weight = ((LHEFEvent*)branchEvent->At(0))->Weight;
123  // If the events do not have any weights / the wrong weight branch is used, Delphes will save NAN in the corresponding branch
124  if (weight != weight) {
125  weight = ((HepMCEvent*)branchEvent->At(0))->Weight;
126  if (weight != weight)
127  weight=1.;
128  }
129 
130 
131  sumOfWeights += weight;
132  sumOfWeights2 += weight*weight;
133 
134  if(!ignoreElectrons) {
135  electrons.clear();
136  if(branchElectron) {
137  for(int i = 0; i < branchElectron->GetEntries(); i++)
138  electrons.push_back((Electron*)branchElectron->At(i));
139  }
140  // All electrons need to fulfill minimum isolation criteria
141  // (Note that this flag is not accessible by the user from within the analyses!)
142  if(!ignoreElectrons_LooseIsolation) {
143  electronsLoose = filterFlags(electrons, "isolation", looseIsolationFlags);
144  if(!ignoreElectrons_MediumEfficiency) {
145  electronsMedium = filterFlags(electronsLoose, "efficiency", mediumEfficiencyFlags);
146  if(!ignoreElectrons_TightEfficiency)
147  electronsTight = filterFlags(electronsMedium, "efficiency", tightEfficiencyFlags);
148  }
149  }
150  }
151 
152  if(!ignoreMuons) {
153  muons.clear();
154  if(branchMuon) {
155  for(int i = 0; i < branchMuon->GetEntries(); i++)
156  muons.push_back((Muon*)branchMuon->At(i));
157  }
158  // All muons need to fulfill minimum isolation criteria
159  if(!ignoreMuons_LooseIsolation) {
160  muons = filterFlags(muons, "isolation", looseIsolationFlags);
161 
162  if(!ignoreMuons_CombinedPlusEfficiency) {
163  muonsCombinedPlus = filterFlags(muons, "efficiency", mediumEfficiencyFlags);
164  if(!ignoreMuons_CombinedEfficiency)
165  muonsCombined = filterFlags(muonsCombinedPlus, "efficiency", tightEfficiencyFlags);
166  }
167  }
168  }
169 
170  jets.clear();
171  for(int i = 0; i < branchJet->GetEntries(); i++)
172  jets.push_back((Jet*)branchJet->At(i));
173 
174  jets2.clear();
175  if(branchJet2) {
176  for(int i = 0; i < branchJet2->GetEntries(); i++)
177  jets2.push_back((Jet*)branchJet2->At(i));
178  }
179 
180  if(!ignorePhotons) {
181  photons.clear();
182  if(branchPhoton) {
183  for(int i = 0; i < branchPhoton->GetEntries(); i++)
184  photons.push_back((Photon*)branchPhoton->At(i));
185  }
186  // All photons need to fulfill minimum isolation criteria
187  if(!ignorePhotons_LooseIsolation)
188  photons = filterFlags(photons, "isolation", looseIsolationFlags);
189  }
190 
191  if(!ignoreTracks) {
192  tracks.clear();
193  if(branchTrack) {
194  for(int i = 0; i < branchTrack->GetEntries(); i++)
195  tracks.push_back((Track*)branchTrack->At(i));
196  }
197  }
198 
199  if(!ignoreTowers) {
200  towers.clear();
201  if(branchTower) {
202  for(int i = 0; i < branchTower->GetEntries(); i++)
203  towers.push_back((Tower*)branchTower->At(i));
204  }
205  }
206 
207  missingET = new ETMiss((MissingET*)branchMissingET->At(0));
208 
209  analyze();
210  }
211  finalize();
212 
213  // Save results of signal-, control- and/or cutflow-regions in specific files
214  if(!cutflowRegions.empty()) {
215  int cutflowOutput = bookFile(analysis+"_cutflow.dat");
216  *fStreams[cutflowOutput] << "Cut Sum_W Sum_W2 Acc N_Norm\n";
217  for(std::map<std::string, double>::iterator cutflow_iter=cutflowRegions.begin(); cutflow_iter!=cutflowRegions.end(); ++cutflow_iter)
218  *fStreams[cutflowOutput] << cutflow_iter->first << " " << cutflow_iter->second << " " << cutflowRegions2[cutflow_iter->first] << " " << cutflow_iter->second/sumOfWeights << " " << normalize(cutflow_iter->second) << "\n";
219  }
220  if(!signalRegions.empty()) {
221  int signalOutput = bookFile(analysis+"_signal.dat");
222  *fStreams[signalOutput] << "SR Sum_W Sum_W2 Acc N_Norm\n";
223  for(std::map<std::string, double>::iterator signal_iter=signalRegions.begin(); signal_iter!=signalRegions.end(); ++signal_iter)
224  *fStreams[signalOutput] << signal_iter->first << " " << signal_iter->second << " " << signalRegions2[signal_iter->first] << " " << signal_iter->second/sumOfWeights << " " << normalize(signal_iter->second) << "\n";
225  }
226  if(!controlRegions.empty()) {
227  int controlOutput = bookFile(analysis+"_control.dat");
228  *fStreams[controlOutput] << "CR Sum_W Sum_W2 Acc N_Norm\n";
229  for(std::map<std::string, double>::iterator control_iter=controlRegions.begin(); control_iter!=controlRegions.end(); ++control_iter)
230  *fStreams[controlOutput] << control_iter->first << " " << control_iter->second << " " << controlRegions2[control_iter->first] << " " << control_iter->second/sumOfWeights << " " << normalize(control_iter->second) << "\n";
231  }
232 }
233 
234 
235 int AnalysisBase::bookFile(std::string name, bool noheader) {
236  // Assemble absolute filename
237  std::string filename = outputFolder+"/"+outputPrefix+"_"+name;
238  // Open stream
239  std::ofstream* file = new ofstream(filename.c_str());
240  // Write standard information
241  if (!noheader) {
242  *file << information << "\n";
243  *file << "@Inputfile: " << inputFile << "\n";
244  *file << "@XSect: " << xsect << " fb\n";
245  *file << "@ Error: " << xsecterr << " fb\n";
246  *file << "@MCEvents: " << nEvents << "\n";
247  *file << "@ SumOfWeights: " << sumOfWeights << "\n";
248  *file << "@ SumOfWeights2: " << sumOfWeights2 << "\n";
249  *file << "@ NormEvents: " << normalize(sumOfWeights) << "\n\n";
250  }
251  fStreams.push_back(file);
252  fNames.push_back(filename);
253  // The returned number denotes the corresponding index for the fStreams vector
254  return fStreams.size()-1;
255 }
256 
257 
258 bool AnalysisBase::checkTauTag(Jet* candidate, std::string efficiency) {
259  if ((efficiency == "loose") && (candidate->TauFlags >= 1))
260  return true;
261  else if ((efficiency == "medium") && (candidate->TauFlags >= 2))
262  return true;
263  else if ((efficiency == "tight") && (candidate->TauFlags == 3))
264  return true;
265  return false;
266 }
267 
268 
269 bool AnalysisBase::checkBTag(Jet* candidate, int relative_flag, std::string option) {
270  // First, decode the member's isoflag into a vector of valid flags
271  int absolute_flag = 0;
272  if(option == "secondJet") {
273  if (relative_flag >= jet2BTagFlags.size()) {
274  std::cerr << "Error: There is no second jet b tag " << relative_flag << std::endl;
275  std::cerr << "Exiting... "<< std::endl;
276  exit(1);
277  }
278  absolute_flag = jet2BTagFlags[relative_flag];
279  }
280  else {
281  if (relative_flag >= jetBTagFlags.size()) {
282  std::cerr << "Error: There is no b tag " << relative_flag << std::endl;
283  std::cerr << "Exiting... "<< std::endl;
284  exit(1);
285  }
286  absolute_flag = jetBTagFlags[relative_flag];
287  }
288 
289  int code = candidate->BFlags;
290  std::vector<int> candidate_flags = code_to_flags(code);
291  // Return bool whether the flag is in the candidate's flags
292  bool result = (std::find(candidate_flags.begin(), candidate_flags.end(), absolute_flag) != candidate_flags.end());
293  return result;
294 }
295 
296 void AnalysisBase::ignore(std::string ignore_what) {
297  if(ignore_what == "electrons")
298  ignoreElectrons = true;
299  else if(ignore_what == "electrons_looseIsolation")
300  ignoreElectrons_LooseIsolation = true;
301  else if(ignore_what == "electronsMedium")
302  ignoreElectrons_MediumEfficiency = true;
303  else if(ignore_what == "electronsTight")
304  ignoreElectrons_TightEfficiency = true;
305  else if(ignore_what == "muons")
306  ignoreMuons = true;
307  else if(ignore_what == "muons_looseIsolation")
308  ignoreMuons_LooseIsolation = true;
309  else if(ignore_what == "muonsCombinedPlus")
310  ignoreMuons_CombinedPlusEfficiency = true;
311  else if(ignore_what == "muonsCombined")
312  ignoreMuons_CombinedEfficiency = true;
313  else if(ignore_what == "photons")
314  ignorePhotons = true;
315  else if(ignore_what == "photons_looseIsolation")
316  ignorePhotons_LooseIsolation = true;
317  else if(ignore_what == "tracks")
318  ignoreTracks = true;
319  else if(ignore_what == "towers")
320  ignoreTowers = true;
321  else {
322  std::cerr << "Cannot ignore " << ignore_what << std::endl;
323  std::cerr << "Exiting." << std::endl;
324  exit(1);
325  }
326 }
327 
328 double AnalysisBase::mT(const TLorentzVector & vis, const TLorentzVector & invis) {
329 
330  return sqrt(2.*vis.Pt()*invis.Et()*(1.-cos(fabs(vis.DeltaPhi(invis)))));
331 }
332 
333 double AnalysisBase::mT2(const TLorentzVector & vis1, const TLorentzVector & vis2, double m_inv, const TLorentzVector & invis) {
334  // Setup mt2 evaluation object.
335  mt2_bisect::mt2 mt2_event;
336  TLorentzVector zeroVector = TLorentzVector(0. ,0. ,0. ,0.);
337  // If no invis is given, use missingET. Note that pmiss[0] is irrelvant, which is why we set it to -1.
338  double pmiss[] = {-1, missingET->P4().Px(), missingET->P4().Py()};
339  if (invis != zeroVector) {
340  pmiss[0] = -1;
341  pmiss[1] = invis.Px();
342  pmiss[2] = invis.Py();
343  }
344 
345  // Construct arrays that mt2_bisect needs as input and start evaluation
346  double p1[3] = {vis1.M(), vis1.Px(), vis1.Py()};
347  double p2[3] = {vis2.M(), vis2.Px(), vis2.Py()};
348  mt2_event.set_momenta( p1, p2, pmiss );
349  mt2_event.set_mn( m_inv );
350  return mt2_event.get_mt2();
351 }
352 
353 double AnalysisBase::mCT(const TLorentzVector & v1,const TLorentzVector & v2)
354 {
355  mctlib::mct mct_event;
356  double v1t[4] = {v1.E(),v1.Px(),v1.Py(),v1.Pz()};
357  double v2t[4] = {v2.E(),v2.Px(),v2.Py(),v2.Pz()};
358  return mct_event.mctnorm(v1t,v2t); //returns 'normal' mCT
359 }
360 
361 double AnalysisBase::mCTcorr(const TLorentzVector & v1, const TLorentzVector & v2,const TLorentzVector & vds,const TLorentzVector & invis,const double ecm,const double mxlo)
362 {
363  mctlib::mct mct_event;
364  double v1t[4] = {v1.E(),v1.Px(),v1.Py(),v1.Pz()};
365  double v2t[4] = {v2.E(),v2.Px(),v2.Py(),v2.Pz()};
366  double vdst[4] = {vds.E(),vds.Px(),vds.Py(),vds.Pz()};
367  double ptmt[2] = {invis.Px(),invis.Py()};
368  return mct_event.mctcorr(v1t,v2t,vdst,ptmt,ecm,mxlo);
369 }
370 
371 double AnalysisBase::mCTy(const TLorentzVector & v1, const TLorentzVector & v2,const TLorentzVector & vds,const TLorentzVector & invis)
372 {
373  mctlib::mct mct_event;
374  double v1t[4] = {v1.E(),v1.Px(),v1.Py(),v1.Pz()};
375  double v2t[4] = {v2.E(),v2.Px(),v2.Py(),v2.Pz()};
376  double vdst[4] = {vds.E(),vds.Px(),vds.Py(),vds.Pz()};
377  double ptmt[2] = {invis.Px(),invis.Py()};
378  return mct_event.mcy(v1t,v2t,vdst,ptmt);
379 }
380 
381 double AnalysisBase::mT2_bl(const TLorentzVector & pl_in, const TLorentzVector & pb1_in, const TLorentzVector & pb2_in, const TLorentzVector & invis) {
382  // Setup mt2_bl evaluation object.
383  mt2bl_bisect::mt2bl mt2bl_event;
384  TLorentzVector zeroVector = TLorentzVector(0. ,0. ,0. ,0.);
385 
386  double pl[4] = {pl_in.E(), pl_in.Px(), pl_in.Py(), pl_in.Pz()}; // El, plx, ply, plz, (visible lepton)
387  double pb1[4] = {pb1_in.E(), pb1_in.Px(), pb1_in.Py(), pb1_in.Pz()}; // Eb1, pb1x, pb1y, pb1z (bottom on the same side as the visible lepton)
388  double pb2[4] = {pb2_in.E(), pb2_in.Px(), pb2_in.Py(), pb2_in.Pz()}; // Eb2, pb2x, pb2y, pb2z (other bottom, paired with the invisible W)
389 
390  // If no invis is given, use missingET. Note that pmiss[0] is irrelvant, which is why we set it to -1.
391  double pmiss[] = {-1, missingET->P4().Px(), missingET->P4().Py()};
392  if (invis != zeroVector) {
393  pmiss[0] = -1;
394  pmiss[1] = invis.Px();
395  pmiss[2] = invis.Py();
396  }
397 
398  // Construct arrays that mt2_bisect needs as input and start evaluation
399  mt2bl_event.set_momenta(pl,pb1,pb2,pmiss);
400 
401  return mt2bl_event.get_mt2bl();
402 }
403 
404 double AnalysisBase::alphaT(const std::vector<Jet*> & jets, const double thresh_ET) {
405  // alphaT code supplied by CMS
406  // Need to pass jets and energy threshold of jets to be considered for variable (allows different ET to that of baseline jets)
407 
408  double HT = 0.;
409  double mHTNorm = 0.;
410  std::vector<Jet*> jetsThresh;
411  TLorentzVector vecHT;
412  for (int i = 0; i < jets.size(); i++) {
413  double ET = jets[i]->P4().Et();
414  if (ET > thresh_ET) {
415  jetsThresh.push_back(jets[i]);
416  HT += ET;
417  vecHT += jets[i]->P4();
418  }
419  }
420  double mHT = vecHT.Pt();
421 
422  std::vector<double> diff( 1<<(jetsThresh.size()-1) , 0. );
423  for (unsigned i=0; i < diff.size(); i++) {
424  for (unsigned j=0; j < jetsThresh.size(); j++) {
425  diff[i] += jetsThresh[j]->P4().Et() * ( 1 - 2 * (int(i>>j)&1) ) ;
426  }
427  }
428  double DHT = fabs( *min_element( diff.begin(), diff.end(), fabs_less() ) );
429  double alphaT = 0.5 * ( HT - DHT ) / sqrt( HT*HT - mHT*mHT );
430 
431  return alphaT;
432 }
433 
434 std::vector<double> AnalysisBase::razor(const std::vector<TLorentzVector> & finalStateObj, const TLorentzVector & invis) {
435  // Razor code supplied by CMS (Maurizio Pierini)
436  // Need to pass objects (4 vectors that could be jets or leptons) and missingET
437  // Result is a 2 dimensional vector, [0] = MR, [1] = R
438 
439  TLorentzVector j1, j2;
440  bool foundGood = false;
441  int N_comb = 1;
442  for(int i = 0; i < finalStateObj.size(); i++)
443  N_comb *= 2;
444 
445  double M_min = 9999999999.0;
446  int j_count;
447  for(int i = 1; i < N_comb-1; i++) {
448  TLorentzVector j_temp1, j_temp2;
449  int itemp = i;
450  j_count = N_comb/2;
451  int count = 0;
452  while(j_count > 0){
453  if(itemp/j_count == 1)
454  j_temp1 += finalStateObj[count];
455  else
456  j_temp2 += finalStateObj[count];
457  itemp -= j_count*(itemp/j_count);
458  j_count /= 2;
459  count++;
460  }
461  double M_temp = j_temp1.M2()+j_temp2.M2();
462  // smallest mass
463  if (M_temp < M_min) {
464  M_min = M_temp;
465  j1 = j_temp1;
466  j2 = j_temp2;
467  }
468  }
469  if (j2.Pt() > j1.Pt()) {
470  TLorentzVector temp = j1;
471  j1 = j2;
472  j2 = temp;
473  }
474 
475  double A = j1.P();
476  double B = j2.P();
477  double az = j1.Pz();
478  double bz = j2.Pz();
479  TVector3 j1T, j2T;
480  j1T.SetXYZ(j1.Px(),j1.Py(),0.0);
481  j2T.SetXYZ(j2.Px(),j2.Py(),0.0);
482  double ATBT = (j1T+j2T).Mag2();
483  double MR = sqrt((A+B)*(A+B)-(az+bz)*(az+bz)-(j2T.Dot(j2T)-j1T.Dot(j1T))*(j2T.Dot(j2T)-j1T.Dot(j1T))/(j1T+j2T).Mag2());
484  double mybeta = (j2T.Dot(j2T)-j1T.Dot(j1T))/sqrt(ATBT*((A+B)*(A+B)-(az+bz)*(az+bz)));
485  double mygamma = 1./sqrt(1.-mybeta*mybeta);
486  //gamma times MRstar
487  MR = MR*mygamma;
488 
489  double MRT = missingET->P4().Et()*(j1.Pt()+j2.Pt()) - missingET->P4().Vect().Dot(j1.Vect()+j2.Vect());
490  MRT /= 2.;
491  MRT = sqrt(MRT);
492 
493  double R = MRT/MR;
494 
495  std::vector<double> razor;
496  razor.push_back(MR);
497  razor.push_back(R);
498 
499  return razor;
500 }
501 
std::vector< std::ofstream * > fStreams
Container of output file streams booked via bookFile()
Definition: AnalysisBase.h:496
double mT2(const TLorentzVector &vis1, const TLorentzVector &vis2, double m_inv, const TLorentzVector &invis=TLorentzVector(0., 0., 0., 0.))
Evaluates .
void bookSignalRegions(std::string listOfRegions)
Function to book signal regions.
Definition: AnalysisBase.cc:56
double alphaT(const std::vector< Jet * > &jets, const double thresh_ET=0.)
Evaluates .
void set_momenta(double *pl0, double *pb10, double *pb20, double *pmiss0)
Definition: mt2bl_bisect.cc:78
TLorentzVector P4()
Returns a TLorentzVector object to the full 4 momentum.
Definition: ETMiss.h:73
std::vector< Muon * > muonsCombined
Container of 'muonsCombinedPlus' objects that pass 'CB' efficiency.
Definition: AnalysisBase.h:107
std::vector< Photon * > photons
Container of all truth photons after detector smearing and loose isolation condition.
Definition: AnalysisBase.h:110
double weight
Current event weight usable for e.g. histograms.
Definition: AnalysisBase.h:562
double normalize(double x)
Normalises number to luminosity and cross section.
Definition: AnalysisBase.h:547
virtual void initialize()
Function to prepare the analysis.
Definition: AnalysisBase.h:80
double mT(const TLorentzVector &vis, const TLorentzVector &invis)
Evaluates the transverse mass .
std::vector< Electron * > electronsLoose
Container of 'electrons' objects that pass loose isolation condition.
Definition: AnalysisBase.h:102
std::vector< Jet * > jets2
Container of all reconstructed 'second type' jets if defined in analysis settings.
Definition: AnalysisBase.h:109
double mctnorm(const double v1[4], const double v2[4])
Definition: mctlib.cc:87
double mctcorr(const double v1[4], const double v2[4], const double vds[4], const double ptm[2], const double ecm=14000.0, const double mxlo=0.0)
Definition: mctlib.cc:10
~AnalysisBase()
Destructor function to free pointers and close opened file streams.
Definition: AnalysisBase.cc:44
double mCT(const TLorentzVector &v1, const TLorentzVector &v2)
Evaluates normal .
void bookCutflowRegions(std::string listOfRegions)
Function to book cutflow regions.
Definition: AnalysisBase.cc:89
double mT2_bl(const TLorentzVector &pl_in, const TLorentzVector &pb1_in, const TLorentzVector &pb2_in, const TLorentzVector &invis=TLorentzVector(0., 0., 0., 0.))
Evaluates .
ETMiss * missingET
Reconstructed missingET vector excluding muons.
Definition: AnalysisBase.h:113
std::vector< Electron * > electrons
Container of all truth electrons after detector smearing in acceptance range.
Definition: AnalysisBase.h:101
ExRootResult * result
Object to ExRootAnalysis for internal studies.
Definition: AnalysisBase.h:549
Int_t TauFlags
std::vector< Electron * > electronsTight
Container of 'electronsLoose' objects that pass 'medium' efficiency cut.
Definition: AnalysisBase.h:104
void set_momenta(double *pa0, double *pb0, double *pmiss0)
Definition: mt2_bisect.cc:61
virtual void finalize()
Function to finish the analyis.
Definition: AnalysisBase.h:84
virtual void analyze()
Function containing the event wise analysis code.
Definition: AnalysisBase.h:82
const double B
Definition: Units.h:11
void loopOverEvents()
Internal function that starts the event loop in main().
std::vector< Muon * > muons
Container of all truth muons after detector smearing in acceptance range.
Definition: AnalysisBase.h:105
double get_mt2()
Definition: mt2_bisect.cc:49
A class to parametrise missingET.
Definition: ETMiss.h:19
std::vector< Jet * > jets
Container of all reconstructed jets.
Definition: AnalysisBase.h:108
std::vector< Tower * > towers
Container of all calorimeter towers.
Definition: AnalysisBase.h:112
AnalysisBase(std::string inFile, std::string outFol, std::string outPre, double xs, double xserr, std::map< std::string, std::string > branches, std::map< std::string, std::vector< int > > flags)
Constructor function to load the Delphes ROOT file.
Definition: AnalysisBase.cc:3
bool checkBTag(Jet *candidate, int relative_flag=0, std::string option="")
Checks if candidate jet fulfills given b-jet identification.
int bookFile(std::string name, bool noheader=false)
Function to book file streams accessible via fStreams and fNames.
void set_mn(double mn)
Definition: mt2_bisect.cc:129
Int_t BFlags
std::vector< Electron * > electronsMedium
Container of 'electronsLoose' objects that pass 'medium' efficiency cut.
Definition: AnalysisBase.h:103
std::vector< double > razor(const std::vector< TLorentzVector > &obj, const TLorentzVector &invis=TLorentzVector(0., 0., 0., 0.))
Evaluates razor.
bool checkTauTag(Jet *candidate, std::string efficiency)
Checks if candidate jet fulfills given tau identification cut.
double mCTcorr(const TLorentzVector &v1, const TLorentzVector &v2, const TLorentzVector &vds, const TLorentzVector &invis, const double ecm=8000.0, const double mxlo=0.0)
Evaluates boost corrected .
std::vector< Muon * > muonsCombinedPlus
Container of objects that pass loose isolation condition and 'CBplusST' efficiency.
Definition: AnalysisBase.h:106
void ignore(std::string ignore_what)
Does not read out unneeded ROOT information.
void bookControlRegions(std::string listOfRegions)
Function to book control regions.
Definition: AnalysisBase.cc:73
double mcy(const double v1[4], const double v2[4], const double vds[4], const double ptm[2])
Definition: mctlib.cc:196
double mCTy(const TLorentzVector &v1, const TLorentzVector &v2, const TLorentzVector &vds, const TLorentzVector &invis)
Evaluates transverse.
std::vector< Track * > tracks
Container of all reconstructed tracks.
Definition: AnalysisBase.h:111
std::vector< std::string > fNames
Container of output file names booked via bookFile()
Definition: AnalysisBase.h:497