ATLAS Offline Software
MissingMassCalculator.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 // vim: ts=8 sw=2
6 /*
7  Missing Mass Calculator
8 */
9 //
10 // to be done : tau 4-vect and type should be data member of MMC.
11 
12 // if histogram smoothing
13 //#define SMOOTH
14 
15 #include "DiTauMassTools/MissingMassCalculator.h" // this is for RootCore package
16 #include <fstream>
17 #include <iomanip>
18 #include <iostream>
19 #include <sstream>
20 // #include "MissingMassCalculator.h" // this is for standalone
21 // package
22 
23 #include <TObject.h>
24 // SpeedUp committed from revision 163876
25 #include <TF1.h>
26 #include <TFitResult.h>
27 #include <TFitResultPtr.h>
28 #include <TMatrixDSym.h>
29 #include <TObject.h>
30 #include <TVectorD.h>
31 #include "Math/VectorUtil.h"
32 
33 using namespace DiTauMassTools;
34 using ROOT::Math::PtEtaPhiMVector;
35 using ROOT::Math::PxPyPzMVector;
36 using ROOT::Math::XYVector;
39 
40 //______________________________constructor________________________________
42  MMCCalibrationSet::e aset, std::string m_paramFilePath)
43  : m_randomGen(), Prob(new MissingMassProb(aset, m_paramFilePath)) {
44  m_mmcCalibrationSet = aset;
46  preparedInput.m_beamEnergy = 6500.0; // for now LHC default is sqrt(S)=7 TeV
47  m_niter_fit1 = 20;
48  m_niter_fit2 = 30;
49  m_niter_fit3 = 10;
50  m_NsucStop = -1;
51  m_NiterRandom = -1; // if the user does not set it to positive value,will be set
52  // in SpaceWalkerInit
53  m_niterRandomLocal = -1; // niterandom which is really used
54  // to be used with RMSSTOP NiterRandom=10000000; // number of random
55  // iterations for lh. Multiplied by 10 for ll, divided by 10 for hh (to be
56  // optimised)
57  // RMSStop=200;// Stop criteria depending of rms of histogram
58  m_reRunWithBestMET = false;
59  m_RMSStop = -1; // disable
60 
61  m_RndmSeedAltering = 0; // can be changed to re-compute with different random seed
62  m_dRmax_tau = 0.4; // changed from 0.2
63  m_nsigma_METscan = -1; // number of sigmas for MET scan
64  m_nsigma_METscan_ll = 3.0; // number of sigmas for MET scan
65  m_nsigma_METscan_lh = 3.0; // number of sigmas for MET scan
66  m_nsigma_METscan_hh = 4.0; // number of sigmas for MET scan (4 for hh 2013)
67  m_nsigma_METscan_lfv_ll = 5.0; // number of sigmas for MET scan (LFV leplep)
68  m_nsigma_METscan_lfv_lh = 5.0; // number of sigmas for MET scan (LFV lephad)
69 
70  m_meanbinStop = -1; // meanbin stopping criterion (-1 if not used)
71  m_proposalTryMEt = -1; // loop on METproposal disable // FIXME should be cleaner
72  m_ProposalTryPhi = -1; // loop on Phiproposal disable
73  m_ProposalTryMnu = -1; // loop on MNuProposal disable
74  m_ProposalTryEtau = -1; // loop on ETauProposal disable
75 
76  Prob->SetUseTauProbability(true); // TauProbability is ON by default DRMERGE comment out for now
77  Prob->SetUseMnuProbability(false); // MnuProbability is OFF by default
78  Prob->SetUseDphiLL(false); // added by Tomas Davidek for lep-lep
79  m_dTheta3d_binMin = 0.0025;
80  m_dTheta3d_binMax = 0.02;
81  preparedInput.m_METresSyst = 0; // no MET resolution systematics by default (+/-1: up/down 1 sigma)
82  preparedInput.m_dataType = 1; // set to "data" by default
83  preparedInput.m_fUseTailCleanup = 1; // cleanup by default for lep-had Moriond 2012 analysis
84  preparedInput.m_fUseDefaults = 0; // use pre-set defaults for various configurations; if set it to 0
85  // if need to study various options
86  m_fUseEfficiencyRecovery = 0; // no re-fit by default
88 
89  preparedInput.m_METScanScheme = 1; // MET-scan scheme: 0- use JER; 1- use simple sumEt & missingHt
90  // for Njet=0 events in (lep-had winter 2012)
91  // MnuScanRange=1.777; // range of M(nunu) scan
92  m_MnuScanRange = 1.5; // better value (sacha)
93  preparedInput.m_LFVmode = -1; // by default consider case of H->mu+tau(->ele)
95 
97  m_iterTheta3d = 0;
98  m_debugThisIteration = false;
99  m_lfvLeplepRefit = true;
100  m_SaveLlhHisto = false;
101 
102  m_nsolmax = 4;
104 
105  m_nuvecsol1.resize(m_nsolmax);
106  m_nuvecsol2.resize(m_nsolmax);
107  m_tauvecsol1.resize(m_nsolmax);
108  m_tauvecsol2.resize(m_nsolmax);
109  m_tauvecprob1.resize(m_nsolmax);
110  m_tauvecprob2.resize(m_nsolmax);
111 
112  m_nsol = 0;
117 
118  m_nsolOld = 0;
123 
124  float hEmax = 3000.0; // maximum energy (GeV)
125  // number of bins
126  int hNbins = 1500; // original 2500 for mass, 10000 for P
127  // choice of hNbins also related to size of window for fitting (see
128  // maxFromHist)
129 
130  //--- define histograms for histogram method
131  //--- upper limits need to be revisied in the future!!! It may be not enough
132  // for some analyses
133 
134  m_fMfit_all = std::make_shared<TH1F>("MMC_h1", "M", hNbins, 0.0,
135  hEmax); // all solutions
136  m_fMfit_all->Sumw2(); // allow proper error bin calculation. Slightly slower but
137  // completely negligible
138 
139  // histogram without weight. useful for debugging. negligibly slow until now
141  std::make_shared<TH1F>("MMC_h1NoW", "M no weight", hNbins, 0.0, hEmax); // all solutions
142 
143  m_fPXfit1 = std::make_shared<TH1F>("MMC_h2", "Px1", 4 * hNbins, -hEmax,
144  hEmax); // Px for tau1
145  m_fPYfit1 = std::make_shared<TH1F>("MMC_h3", "Py1", 4 * hNbins, -hEmax,
146  hEmax); // Py for tau1
147  m_fPZfit1 = std::make_shared<TH1F>("MMC_h4", "Pz1", 4 * hNbins, -hEmax,
148  hEmax); // Pz for tau1
149  m_fPXfit2 = std::make_shared<TH1F>("MMC_h5", "Px2", 4 * hNbins, -hEmax,
150  hEmax); // Px for tau2
151  m_fPYfit2 = std::make_shared<TH1F>("MMC_h6", "Py2", 4 * hNbins, -hEmax,
152  hEmax); // Py for tau2
153  m_fPZfit2 = std::make_shared<TH1F>("MMC_h7", "Pz2", 4 * hNbins, -hEmax,
154  hEmax); // Pz for tau2
155 
156  m_fMfit_all->SetDirectory(0);
157 
158  m_fMfit_allNoWeight->SetDirectory(0);
159  m_fPXfit1->SetDirectory(0);
160  m_fPYfit1->SetDirectory(0);
161  m_fPZfit1->SetDirectory(0);
162  m_fPXfit2->SetDirectory(0);
163  m_fPYfit2->SetDirectory(0);
164  m_fPZfit2->SetDirectory(0);
165 
166  // max hist fitting function
167  m_fFitting =
168  new TF1("MMC_maxFitting", this, &MissingMassCalculator::maxFitting, 0., hEmax, 3);
169  // Sets initial parameter names
170  m_fFitting->SetParNames("Max", "Mean", "InvWidth2");
171 
172  if (preparedInput.m_fUseVerbose == 1) {
173  gDirectory->pwd();
174  gDirectory->ls();
175  }
176 
177  if (preparedInput.m_fUseVerbose == 1) {
178  gDirectory->pwd();
179  gDirectory->ls();
180  }
181 }
182 
184 
185 //_____________________________________________________________________________
186 // Main Method to run MissingMassCalculator
188  const xAOD::IParticle *part2,
189  const xAOD::MissingET *met,
190  const int &njets) {
191  m_reRunWithBestMET = false;
192 
194  if (preparedInput.m_fUseVerbose == 1) {
195  Info("DiTauMassTools", "------------- Raw Input for MissingMassCalculator --------------");
196  }
197  FinalizeSettings(part1, part2, met, njets); // rawInput, preparedInput );
199  if (preparedInput.m_fUseVerbose == 1) {
200  Info("DiTauMassTools", "------------- Prepared Input for MissingMassCalculator--------------");
202  }
203 
204  if (preparedInput.m_LFVmode < 0) {
205  // remove argument DiTauMassCalculatorV9Walk work directly on preparedInput
207 
208  // re-running MMC for on failed events
210  // most events where MMC failed happened to have dPhi>2.9. Run re-fit only
211  // on these events
212  if (preparedInput.m_DelPhiTT > 2.9) {
213  // preparedInput.MetVec.Set(-(preparedInput.vistau1+preparedInput.vistau2).Px(),-(preparedInput.vistau1+preparedInput.vistau2).Py());
214  // // replace MET by MPT
215 
216  XYVector dummy_met(-(preparedInput.m_vistau1 + preparedInput.m_vistau2).Px(),
218  preparedInput.m_METcovphi = dummy_met.Phi();
219  double dummy_METres =
222  dummy_METres * std::abs(cos(dummy_met.Phi() - preparedInput.m_MetVec.Phi()));
224  dummy_METres * std::abs(sin(dummy_met.Phi() - preparedInput.m_MetVec.Phi()));
225  if (preparedInput.m_METsigmaP < 5.0)
227  m_nsigma_METscan_lh = 6.0; // increase range of MET scan
228  m_nsigma_METscan_hh = 6.0; // increase range of MET scan
229 
230  OutputInfo.ClearOutput(preparedInput.m_fUseVerbose); // clear output stuff before re-running
231  OutputInfo.m_FitStatus = DitauMassCalculatorV9walk(); // run MMC again
232  }
233  }
234 
235  }
236 
237  // running MMC in LFV mode for reconstructing mass of X->lep+tau
238  else {
239  if (preparedInput.m_fUseVerbose == 1) {
240  Info("DiTauMassTools", "Calling DitauMassCalculatorV9lfv");
241  }
243  }
244 
245  if(m_SaveLlhHisto){
246  TFile *outFile = TFile::Open("MMC_likelihoods.root", "UPDATE");
247  outFile->cd();
249  if (!outFile->GetDirectory(path.c_str()))
250  outFile->mkdir(path.c_str());
251  outFile->cd(path.c_str());
252  m_fMfit_all->Write(m_fMfit_all->GetName(), TObject::kOverwrite);
253  m_fMEtP_all->Write(m_fMEtP_all->GetName(), TObject::kOverwrite);
254  m_fMEtL_all->Write(m_fMEtL_all->GetName(), TObject::kOverwrite);
255  m_fMnu1_all->Write(m_fMnu1_all->GetName(), TObject::kOverwrite);
256  m_fMnu2_all->Write(m_fMnu2_all->GetName(), TObject::kOverwrite);
257  m_fPhi1_all->Write(m_fPhi1_all->GetName(), TObject::kOverwrite);
258  m_fPhi2_all->Write(m_fPhi2_all->GetName(), TObject::kOverwrite);
259  m_fMfit_allNoWeight->Write(m_fMfit_allNoWeight->GetName(), TObject::kOverwrite);
260  m_fMfit_allGraph->Write("Graph", TObject::kOverwrite);
261  TH1D *nosol = new TH1D("nosol", "nosol", 7, 0, 7);
262  nosol->SetBinContent(1, m_testptn1);
263  nosol->SetBinContent(2, m_testptn2);
264  nosol->SetBinContent(3, m_testdiscri1);
265  nosol->SetBinContent(4, m_testdiscri2);
266  nosol->SetBinContent(5, m_nosol1);
267  nosol->SetBinContent(6, m_nosol1);
268  nosol->SetBinContent(7, m_iterNuPV3);
269  nosol->Write(nosol->GetName(), TObject::kOverwrite);
270  outFile->Write();
271  outFile->Close();
272  }
273 
274  DoOutputInfo();
275  PrintResults();
277  return 1;
278 }
279 
280 //-------- clearing ditau container
282  fStuff.Mditau_best = 0.0;
283  fStuff.Sign_best = 1.0E6;
284  fStuff.nutau1 = PtEtaPhiMVector(0., 0., 0., 0.);
285  fStuff.nutau2 = PtEtaPhiMVector(0., 0., 0., 0.);
286  fStuff.vistau1 = PtEtaPhiMVector(0., 0., 0., 0.);
287  fStuff.vistau2 = PtEtaPhiMVector(0., 0., 0., 0.);
288  fStuff.RMSoverMPV = 0.0;
289 
290  return;
291 }
292 
293 //---------------------------- Accessors to output parameters
294 //------------------------
295 // finalizes output information
297  if (OutputInfo.m_FitStatus > 0) {
298  if (preparedInput.m_fUseVerbose == 1) {
299  Info("DiTauMassTools", "Retrieving output from fDitauStuffFit");
300  }
301  // MAXW method : get from fDittauStuffFit
304  double q1 = (1. - 0.68) / 2.;
305  double q2 = 1. - q1;
306  double xq[2], yq[2];
307  xq[0] = q1;
308  xq[1] = q2;
309  m_fMfit_all->GetQuantiles(2, yq, xq);
321  XYVector metmaxw(OutputInfo.m_nuvec1[MMCFitMethod::MAXW].Px() +
326 
330 
331  PtEtaPhiMVector tlvdummy(0., 0., 0., 0.);
332  XYVector metdummy(0., 0.);
340 
341  // MLNU3P method : get from fDittauStuffHisto 4 momentum
355 
356  XYVector metmlnu3p(OutputInfo.m_nuvec1[MMCFitMethod::MLNU3P].Px() +
361 
363  }
364 
367  OutputInfo.m_NSolutions = m_fMfit_all->GetEntries();
368  OutputInfo.m_SumW = m_fMfit_all->GetSumOfWeights();
369 
370  //----------------- Check if input was re-ordered in FinalizeInputStuff() and
371  // restore the original order if needed
372  if (preparedInput.m_InputReorder == 1) {
373  PtEtaPhiMVector dummy_vec1(0.0, 0.0, 0.0, 0.0);
374  PtEtaPhiMVector dummy_vec2(0.0, 0.0, 0.0, 0.0);
375  for (int i = 0; i < 3; i++) {
376  // re-ordering neutrinos
377  dummy_vec1 = OutputInfo.m_nuvec1[i];
378  dummy_vec2 = OutputInfo.m_nuvec2[i];
379  OutputInfo.m_nuvec1[i] = dummy_vec2;
380  OutputInfo.m_nuvec2[i] = dummy_vec1;
381  // re-ordering tau's
382  dummy_vec1 = OutputInfo.m_objvec1[i];
383  dummy_vec2 = OutputInfo.m_objvec2[i];
384  OutputInfo.m_objvec1[i] = dummy_vec2;
385  OutputInfo.m_objvec2[i] = dummy_vec1;
386  }
387  }
388 
389  return;
390 }
391 
392 // Printout of final results
394  if (preparedInput.m_fUseVerbose != 1)
395  return;
396 
397  Info("DiTauMassTools",
398  ".........................Other input.....................................");
399  Info("DiTauMassTools", "%s",
400  ("Beam energy =" + std::to_string(preparedInput.m_beamEnergy) +
401  " sqrt(S) for collisions =" + std::to_string(2.0 * preparedInput.m_beamEnergy))
402  .c_str());
403  Info("DiTauMassTools", "%s",
404  ("CalibrationSet " + MMCCalibrationSet::name[m_mmcCalibrationSet])
405  .c_str());
406  Info("DiTauMassTools", "%s",
407  ("LFV mode " + std::to_string(preparedInput.m_LFVmode) + " seed=" + std::to_string(m_seed))
408  .c_str());
409  Info("DiTauMassTools", "%s", ("usetauProbability=" + std::to_string(Prob->GetUseTauProbability()) +
410  " useTailCleanup=" + std::to_string(preparedInput.m_fUseTailCleanup))
411  .c_str());
412 
413  if (preparedInput.m_InputReorder != 0) {
414  Info("DiTauMassTools",
415  "tau1 and tau2 were internally swapped (visible on prepared input printout)");
416  } else {
417  Info("DiTauMassTools", "tau1 and tau2 were NOT internally swapped");
418  }
419 
420  Info("DiTauMassTools", "%s",
421  (" MEtLMin=" + std::to_string(m_MEtLMin) + " MEtLMax=" + std::to_string(m_MEtLMax)).c_str());
422  Info("DiTauMassTools", "%s",
423  (" MEtPMin=" + std::to_string(m_MEtPMin) + " MEtPMax=" + std::to_string(m_MEtPMax)).c_str());
424  Info("DiTauMassTools", "%s",
425  (" Phi1Min=" + std::to_string(m_Phi1Min) + " Phi1Max=" + std::to_string(m_Phi1Max)).c_str());
426  Info("DiTauMassTools", "%s",
427  (" Phi2Min=" + std::to_string(m_Phi2Min) + " Phi2Max=" + std::to_string(m_Phi2Max)).c_str());
428  Info("DiTauMassTools", "%s",
429  (" Mnu1Min=" + std::to_string(m_Mnu1Min) + " Mnu1Max=" + std::to_string(m_Mnu1Max)).c_str());
430  Info("DiTauMassTools", "%s",
431  (" Mnu2Min=" + std::to_string(m_Mnu2Min) + " Mnu2Max=" + std::to_string(m_Mnu2Max)).c_str());
432 }
433 
434 // Printout of final results
436 
437  if (preparedInput.m_fUseVerbose != 1)
438  return;
439 
440  const PtEtaPhiMVector *origVisTau1 = 0;
441  const PtEtaPhiMVector *origVisTau2 = 0;
442 
443  if (preparedInput.m_InputReorder == 0) {
444  origVisTau1 = &preparedInput.m_vistau1;
445  origVisTau2 = &preparedInput.m_vistau2;
446  } else // input order was flipped
447  {
448  origVisTau1 = &preparedInput.m_vistau2;
449  origVisTau2 = &preparedInput.m_vistau1;
450  }
451 
452  PrintOtherInput();
453 
454  Info("DiTauMassTools",
455  "------------- Printing Final Results for MissingMassCalculator --------------");
456  Info("DiTauMassTools",
457  ".............................................................................");
458  Info("DiTauMassTools", "%s", ("Fit status=" + std::to_string(OutputInfo.m_FitStatus)).c_str());
459 
460  for (int imeth = 0; imeth < MMCFitMethod::MAX; ++imeth) {
461  Info("DiTauMassTools", "%s",
462  ("___ Results for " + MMCFitMethod::name[imeth] + "Method ___")
463  .c_str());
464  Info("DiTauMassTools", "%s",
465  (" signif=" + std::to_string(OutputInfo.m_FitSignificance[imeth])).c_str());
466  Info("DiTauMassTools", "%s", (" mass=" + std::to_string(OutputInfo.m_FittedMass[imeth])).c_str());
467  Info("DiTauMassTools", "%s", (" rms/mpv=" + std::to_string(OutputInfo.m_RMS2MPV)).c_str());
468 
469  if (imeth == MMCFitMethod::MLM) {
470  Info("DiTauMassTools", " no 4-momentum or MET from this method ");
471  continue;
472  }
473 
474  if (OutputInfo.m_FitStatus <= 0) {
475  Info("DiTauMassTools", " fit failed ");
476  }
477 
478  const PtEtaPhiMVector &tlvnu1 = OutputInfo.m_nuvec1[imeth];
479  const PtEtaPhiMVector &tlvnu2 = OutputInfo.m_nuvec2[imeth];
480  const PtEtaPhiMVector &tlvo1 = OutputInfo.m_objvec1[imeth];
481  const PtEtaPhiMVector &tlvo2 = OutputInfo.m_objvec2[imeth];
482  const XYVector &tvmet = OutputInfo.m_FittedMetVec[imeth];
483 
484  Info("DiTauMassTools", "%s",
485  (" Neutrino-1: P=" + std::to_string(tlvnu1.P()) + " Pt=" + std::to_string(tlvnu1.Pt()) +
486  " Eta=" + std::to_string(tlvnu1.Eta()) + " Phi=" + std::to_string(tlvnu1.Phi()) +
487  " M=" + std::to_string(tlvnu1.M()) + " Px=" + std::to_string(tlvnu1.Px()) +
488  " Py=" + std::to_string(tlvnu1.Py()) + " Pz=" + std::to_string(tlvnu1.Pz()))
489  .c_str());
490  Info("DiTauMassTools", "%s",
491  (" Neutrino-2: P=" + std::to_string(tlvnu2.P()) + " Pt=" + std::to_string(tlvnu2.Pt()) +
492  " Eta=" + std::to_string(tlvnu2.Eta()) + " Phi=" + std::to_string(tlvnu2.Phi()) +
493  " M=" + std::to_string(tlvnu2.M()) + " Px=" + std::to_string(tlvnu2.Px()) +
494  " Py=" + std::to_string(tlvnu2.Py()) + " Pz=" + std::to_string(tlvnu2.Pz()))
495  .c_str());
496  Info("DiTauMassTools", "%s",
497  (" Tau-1: P=" + std::to_string(tlvo1.P()) + " Pt=" + std::to_string(tlvo1.Pt()) +
498  " Eta=" + std::to_string(tlvo1.Eta()) + " Phi=" + std::to_string(tlvo1.Phi()) +
499  " M=" + std::to_string(tlvo1.M()) + " Px=" + std::to_string(tlvo1.Px()) +
500  " Py=" + std::to_string(tlvo1.Py()) + " Pz=" + std::to_string(tlvo1.Pz()))
501  .c_str());
502  Info("DiTauMassTools", "%s",
503  (" Tau-2: P=" + std::to_string(tlvo2.P()) + " Pt=" + std::to_string(tlvo2.Pt()) +
504  " Eta=" + std::to_string(tlvo2.Eta()) + " Phi=" + std::to_string(tlvo2.Phi()) +
505  " M=" + std::to_string(tlvo2.M()) + " Px=" + std::to_string(tlvo2.Px()) +
506  " Py=" + std::to_string(tlvo2.Py()) + " Pz=" + std::to_string(tlvo2.Pz()))
507  .c_str());
508 
509  Info("DiTauMassTools", "%s",
510  (" dR(nu1-visTau1)=" + std::to_string(DeltaR(tlvnu1,*origVisTau1))).c_str());
511  Info("DiTauMassTools", "%s",
512  (" dR(nu2-visTau2)=" + std::to_string(DeltaR(tlvnu2,*origVisTau2))).c_str());
513 
514  Info("DiTauMassTools", "%s",
515  (" Fitted MET =" + std::to_string(tvmet.R()) + " Phi=" + std::to_string(tlvnu1.Phi()) +
516  " Px=" + std::to_string(tvmet.X()) + " Py=" + std::to_string(tvmet.Y()))
517  .c_str());
518 
519  Info("DiTauMassTools", "%s", (" Resonance: P=" + std::to_string(OutputInfo.m_totalvec[imeth].P()) +
520  " Pt=" + std::to_string(OutputInfo.m_totalvec[imeth].Pt()) +
521  " Eta=" + std::to_string(OutputInfo.m_totalvec[imeth].Eta()) +
522  " Phi=" + std::to_string(OutputInfo.m_totalvec[imeth].Phi()) +
523  " M=" + std::to_string(OutputInfo.m_totalvec[imeth].M()) +
524  " Px=" + std::to_string(OutputInfo.m_totalvec[imeth].Px()) +
525  " Py=" + std::to_string(OutputInfo.m_totalvec[imeth].Py()) +
526  " Pz=" + std::to_string(OutputInfo.m_totalvec[imeth].Pz()))
527  .c_str());
528  }
529 
530  return;
531 }
532 
533 // returns P1, P2, and theta1 & theta2 solutions
534 // This compute the nu1 nu2 solution in the most efficient way. Wrt to
535 // NuPsolutionV2, the output nu1 nu2 4-vector have non zero mass (if relevant).
536 // It is not optimised for grid running so much less caching is done (which
537 // makes it more readable). Only quantities fixed within an event are cached.
538 // relies on a number of these variables to be initialised before the loop.
539 
540 int MissingMassCalculator::NuPsolutionV3(const double &mNu1, const double &mNu2,
541  const double &phi1, const double &phi2,
542  int &nsol1, int &nsol2) {
543 
544  // Pv1, Pv2 : visible tau decay product momentum
545  // Pn1 Pn2 : neutrino momentum
546  // phi1, phi2 : neutrino azymutal angles
547  // PTmiss2=PTmissy Cos[phi2] - PTmissx Sin[phi2]
548  // PTmiss2cscdphi=PTmiss2/Sin[phi1-phi2]
549  // Pv1proj=Pv1x Cos[phi1] + Pv1y Sin[phi1]
550  // M2noma1=Mtau^2-Mv1^2-Mn1^2
551  // ETv1^2=Ev1^2-Pv1z^2
552 
553  // discriminant : 16 Ev1^2 (M2noma1^2 + 4 M2noma1 PTmiss2cscdphi Pv1proj - 4
554  // (ETv1^2 (Mn1^2 + PTmiss2cscdphi^2) - PTmiss2cscdphi^2 Pv1proj^2))
555  // two solutions for epsilon = +/-1
556  // Pn1z=(1/(2 ETv1^2))(epsilon Ev1 Sqrt[ M2noma1^2 + 4 M2noma1 PTmiss2cscdphi
557  // Pv1proj - 4 (ETv1^2 (Mn1^2 + qPTmiss2cscdphi^2) - PTmiss2cscdphi^2
558  // Pv1proj^2)] + M2noma1 Pv1z + 2 PTmiss2cscdphi Pv1proj Pv1z)
559  // with conditions: M2noma1 + 2 PTmiss2cscdphi Pv1proj + 2 Pn1z Pv1z > 0
560  // PTn1 -> PTmiss2 Csc[phi1 - phi2]
561 
562  // if initialisation precompute some quantities
563  int solution_code = 0; // 0 with no solution, 1 with solution
564  nsol1 = 0;
565  nsol2 = 0;
566 
567  // Variables used to test PTn1 and PTn2 > 0
568 
569  const double &pTmissx = preparedInput.m_MEtX;
570  const double &pTmissy = preparedInput.m_MEtY;
571 
573  double pTmiss2 = pTmissy * m_cosPhi2 - pTmissx * m_sinPhi2;
574 
575  int dPhiSign = 0;
576  dPhiSign = fixPhiRange(phi1 - phi2) > 0 ? +1 : -1;
577 
578  // Test if PTn1 and PTn2 > 0. Then MET vector is between the two neutrino
579  // vector
580 
581  if (pTmiss2 * dPhiSign < 0) {
582  ++m_testptn1;
583  return solution_code;
584  }
585 
587  double pTmiss1 = pTmissy * m_cosPhi1 - pTmissx * m_sinPhi1;
588 
589  if (pTmiss1 * (-dPhiSign) < 0) {
590  ++m_testptn2;
591  return solution_code;
592  }
593 
594  // Variables used to calculate discri1
595 
596  double m2Vis1 = m_tauVec1M * m_tauVec1M;
598  m_m2Nu1 = mNu1 * mNu1;
599  double m2noma1 = m_mTau2 - m_m2Nu1 - m2Vis1;
600  double m4noma1 = m2noma1 * m2noma1;
601  double pv1proj = m_tauVec1Px * m_cosPhi1 + m_tauVec1Py * m_sinPhi1;
602  double p2v1proj = std::pow(pv1proj, 2);
603  double sinDPhi2 = m_cosPhi2 * m_sinPhi1 - m_sinPhi2 * m_cosPhi1; // sin(Phi1-Phi2)
604  double pTmiss2CscDPhi = pTmiss2 / sinDPhi2;
605  double &pTn1 = pTmiss2CscDPhi;
606  double pT2miss2CscDPhi = pTmiss2CscDPhi * pTmiss2CscDPhi;
607 
608  // Test on discri1
609  const double discri1 = m4noma1 + 4 * m2noma1 * pTmiss2CscDPhi * pv1proj -
610  4 * (m_ET2v1 * (m_m2Nu1 + pT2miss2CscDPhi) - (pT2miss2CscDPhi * p2v1proj));
611 
612  if (discri1 < 0) // discriminant negative -> no solution
613  {
614  ++m_testdiscri1;
615  return solution_code;
616  }
617 
618  // Variables used to calculate discri2
619  double m2Vis2 = m_tauVec2M * m_tauVec2M;
621  m_m2Nu2 = mNu2 * mNu2;
622  double m2noma2 = m_mTau2 - m_m2Nu2 - m2Vis2;
623  double m4noma2 = m2noma2 * m2noma2;
624  double pv2proj = m_tauVec2Px * m_cosPhi2 + m_tauVec2Py * m_sinPhi2;
625  double p2v2proj = std::pow(pv2proj, 2);
626  double sinDPhi1 = -sinDPhi2;
627  double pTmiss1CscDPhi = pTmiss1 / sinDPhi1;
628  double &pTn2 = pTmiss1CscDPhi;
629  double pT2miss1CscDPhi = pTmiss1CscDPhi * pTmiss1CscDPhi;
630 
631  const double discri2 = m4noma2 + 4 * m2noma2 * pTmiss1CscDPhi * pv2proj -
632  4 * (m_ET2v2 * (m_m2Nu2 + pT2miss1CscDPhi) - (pT2miss1CscDPhi * p2v2proj));
633 
634  if (discri2 < 0) // discriminant negative -> no solution
635  {
636  ++m_testdiscri2;
637  return solution_code;
638  }
639 
640  // this should be done only once we know there are solutions for nu2
642  m_Ev1 = sqrt(m_E2v1);
643  double sqdiscri1 = sqrt(discri1);
644  double first1 =
645  (m2noma1 * m_tauVec1Pz + 2 * pTmiss2CscDPhi * pv1proj * m_tauVec1Pz) / (2 * m_ET2v1);
646  double second1 = sqdiscri1 * m_Ev1 / (2 * m_ET2v1);
647 
648  // first solution
649  double pn1Z = first1 + second1;
650 
651  if (m2noma1 + 2 * pTmiss2CscDPhi * pv1proj + 2 * pn1Z * m_tauVec1Pz >
652  0) // Condition for solution to exist
653  {
654  m_nuvecsol1[nsol1].SetPxPyPzE(pTn1 * m_cosPhi1, pTn1 * m_sinPhi1, pn1Z,
655  sqrt(std::pow(pTn1, 2) + std::pow(pn1Z, 2) + m_m2Nu1));
656 
657  ++nsol1;
658  }
659 
660  pn1Z = first1 - second1;
661 
662  if (m2noma1 + 2 * pTmiss2CscDPhi * pv1proj + 2 * pn1Z * m_tauVec1Pz >
663  0) // Condition for solution to exist
664  {
665 
666  m_nuvecsol1[nsol1].SetPxPyPzE(pTn1 * m_cosPhi1, pTn1 * m_sinPhi1, pn1Z,
667  sqrt(std::pow(pTn1, 2) + std::pow(pn1Z, 2) + m_m2Nu1));
668 
669  ++nsol1;
670  }
671 
672  if (nsol1 == 0) {
673  ++m_nosol1;
674  return solution_code;
675  }
676 
678  m_Ev2 = sqrt(m_E2v2);
679  double sqdiscri2 = sqrt(discri2);
680  double first2 =
681  (m2noma2 * m_tauVec2Pz + 2 * pTmiss1CscDPhi * pv2proj * m_tauVec2Pz) / (2 * m_ET2v2);
682  double second2 = sqdiscri2 * m_Ev2 / (2 * m_ET2v2);
683 
684  // second solution
685  double pn2Z = first2 + second2;
686 
687  if (m2noma2 + 2 * pTmiss1CscDPhi * pv2proj + 2 * pn2Z * m_tauVec2Pz >
688  0) // Condition for solution to exist
689  {
690  m_nuvecsol2[nsol2].SetPxPyPzE(pTn2 * m_cosPhi2, pTn2 * m_sinPhi2, pn2Z,
691  sqrt(std::pow(pTn2, 2) + std::pow(pn2Z, 2) + m_m2Nu2));
692 
693  ++nsol2;
694  }
695 
696  pn2Z = first2 - second2;
697  ;
698 
699  if (m2noma2 + 2 * pTmiss1CscDPhi * pv2proj + 2 * pn2Z * m_tauVec2Pz >
700  0) // Condition for solution to exist
701  {
702  m_nuvecsol2[nsol2].SetPxPyPzE(pTn2 * m_cosPhi2, pTn2 * m_sinPhi2, pn2Z,
703  sqrt(std::pow(pTn2, 2) + std::pow(pn2Z, 2) + m_m2Nu2));
704 
705  ++nsol2;
706  }
707 
708  if (nsol2 == 0) {
709  ++m_nosol2;
710  return solution_code;
711  }
712 
713  // Verification if solution exist
714 
715  solution_code = 1;
716  ++m_iterNuPV3;
717 
718  // double check solutions from time to time
719  if (m_iterNuPV3 % 1000 == 1) {
720  double pnux = m_nuvecsol1[0].Px() + m_nuvecsol2[0].Px();
721  double pnuy = m_nuvecsol1[0].Py() + m_nuvecsol2[0].Py();
722  double mtau1plus = (m_nuvecsol1[0] + m_tauVec1).M();
723  double mtau1moins = (m_nuvecsol1[1] + m_tauVec1).M();
724  double mtau2plus = (m_nuvecsol2[0] + m_tauVec2).M();
725  double mtau2moins = (m_nuvecsol2[1] + m_tauVec2).M();
726  if (std::abs(pnux - pTmissx) > 0.001 || std::abs(pnuy - pTmissy) > 0.001) {
727  Info("DiTauMassTools", "%s", ("NuPsolutionV3 ERROR Pnux-Met.X or Pnuy-Met.Y > 0.001 : " +
728  std::to_string(pnux - pTmissx) + " and " +
729  std::to_string(pnuy - pTmissx) + " " + "Invalid solutions")
730  .c_str());
731  }
732  if (std::abs(mtau1plus - m_mTau) > 0.001 || std::abs(mtau1moins - m_mTau) > 0.001 ||
733  std::abs(mtau2plus - m_mTau) > 0.001 || std::abs(mtau2moins - m_mTau) > 0.001) {
734  Info("DiTauMassTools", "%s", ("NuPsolutionV3 ERROR tau mass not recovered : " +
735  std::to_string(mtau1plus) + " " + std::to_string(mtau1moins) + " " +
736  std::to_string(mtau2plus) + " " + std::to_string(mtau2moins))
737  .c_str());
738  }
739  }
740 
741  return solution_code;
742 }
743 
744 // returns solution for Lepton Flavor Violating X->lep+tau studies
745 int MissingMassCalculator::NuPsolutionLFV(const XYVector &met_vec,
746  const PtEtaPhiMVector &tau, const double &l_nu,
747  std::vector<PtEtaPhiMVector> &nu_vec) {
748  int solution_code = 0; // 0 with no solution, 1 with solution
749 
750  nu_vec.clear();
751  PxPyPzMVector nu(met_vec.X(), met_vec.Y(), 0.0, l_nu);
752  PxPyPzMVector nu2(met_vec.X(), met_vec.Y(), 0.0, l_nu);
753 
754  const double Mtau = 1.777;
755  // double msq = (Mtau*Mtau-tau.M()*tau.M())/2;
756  double msq = (Mtau * Mtau - tau.M() * tau.M() - l_nu * l_nu) /
757  2; // to take into account the fact that 2-nu systema has mass
758  double gamma = nu.Px() * nu.Px() + nu.Py() * nu.Py();
759  double beta = tau.Px() * nu.Px() + tau.Py() * nu.Py() + msq;
760  double a = tau.E() * tau.E() - tau.Pz() * tau.Pz();
761  double b = -2 * tau.Pz() * beta;
762  double c = tau.E() * tau.E() * gamma - beta * beta;
763  if ((b * b - 4 * a * c) < 0)
764  return solution_code; // no solution found
765  else
766  solution_code = 2;
767  double pvz1 = (-b + sqrt(b * b - 4 * a * c)) / (2 * a);
768  double pvz2 = (-b - sqrt(b * b - 4 * a * c)) / (2 * a);
769 
770  nu.SetCoordinates(met_vec.X(), met_vec.Y(), pvz1, l_nu);
771  nu2.SetCoordinates(met_vec.X(), met_vec.Y(), pvz2, l_nu);
772 
773  PtEtaPhiMVector return_nu(nu.Pt(), nu.Eta(), nu.Phi(), nu.M());
774  PtEtaPhiMVector return_nu2(nu2.Pt(), nu2.Eta(), nu2.Phi(), nu2.M());
775  nu_vec.push_back(return_nu);
776  nu_vec.push_back(return_nu2);
777  return solution_code;
778 }
779 
780 // like v9fast, but the parameter space scanning is now factorised out, to allow
781 // flexibility
783 
784  int nsuccesses = 0;
785 
786  int fit_code = 0; // 0==bad, 1==good
789  OutputInfo.m_AveSolRMS = 0.;
790 
791  m_fMfit_all->Reset();
792 
793  if(m_SaveLlhHisto){
794  m_fMEtP_all->Reset();
795  m_fMEtL_all->Reset();
796  m_fMnu1_all->Reset();
797  m_fMnu2_all->Reset();
798  m_fPhi1_all->Reset();
799  m_fPhi2_all->Reset();
800  }
801 
802  m_fMfit_allNoWeight->Reset();
803  m_fPXfit1->Reset();
804  m_fPYfit1->Reset();
805  m_fPZfit1->Reset();
806  m_fPXfit2->Reset();
807  m_fPYfit2->Reset();
808  m_fPZfit2->Reset();
809 
810  // these histograms are used for the floating stopping criterion
811  if (m_fUseFloatStopping) {
812  m_fMmass_split1->Reset();
813  m_fMEtP_split1->Reset();
814  m_fMEtL_split1->Reset();
815  m_fMnu1_split1->Reset();
816  m_fMnu2_split1->Reset();
817  m_fPhi1_split1->Reset();
818  m_fPhi2_split1->Reset();
819  m_fMmass_split2->Reset();
820  m_fMEtP_split2->Reset();
821  m_fMEtL_split2->Reset();
822  m_fMnu1_split2->Reset();
823  m_fMnu2_split2->Reset();
824  m_fPhi1_split2->Reset();
825  m_fPhi2_split2->Reset();
826  }
827 
828  m_prob_tmp = 0.0;
829 
830  m_iter1 = 0;
831 
832  m_totalProbSum = 0;
833  m_mtautauSum = 0;
834 
835  XYVector deltamet_vec;
836 
837  // initialize a spacewalker, which walks the parameter space according to some
838  // algorithm
839  SpaceWalkerInit();
840 
841  while (SpaceWalkerWalk()) {
842  bool paramInsideRange = false;
843  m_nsol = 0;
844 
845  paramInsideRange = checkAllParamInRange();
846 
847  // FIXME if no tau scanning, or symmetric matrices, rotatin is made twice
848  // which is inefficient
849  const double deltaMetx = m_MEtL * m_metCovPhiCos - m_MEtP * m_metCovPhiSin;
850  const double deltaMety = m_MEtL * m_metCovPhiSin + m_MEtP * m_metCovPhiCos;
851 
852  // deltaMetVec.Set(met_smear_x,met_smear_y);
854  preparedInput.m_inputMEtY + deltaMety);
855 
856  // save in global variable for speed sake
860 
861  if (paramInsideRange)
863 
864  // DR for markov chain need to enter handleSolution also when zero solutions
865  handleSolutions();
866  // be careful that with markov, current solution is from now on stored in
867  // XYZOldSolVec
868 
869  if (m_nsol <= 0)
870  continue;
871 
872  // for markov, nsuccess more difficult to define. Decide this is the number
873  // of independent point accepted (hence without weight)
874  nsuccesses = m_markovNAccept;
876 
877  m_iter1 += m_nsol;
878  fit_code = 1;
879 
880  } // while loop
881 
883  OutputInfo.m_NSuccesses = nsuccesses;
884 
885  if (nsuccesses > 0) {
886  OutputInfo.m_AveSolRMS /= nsuccesses;
887  } else {
888  OutputInfo.m_AveSolRMS = -1.;
889  }
890 
891  double Px1, Py1, Pz1;
892  double Px2, Py2, Pz2;
893  if (nsuccesses > 0) {
894 
895  // note that smoothing can slightly change the integral of the histogram
896 
897 #ifdef SMOOTH
898  m_fMfit_all->Smooth();
899  m_fMfit_allNoWeight->Smooth();
900  m_fPXfit1->Smooth();
901  m_fPYfit1->Smooth();
902  m_fPZfit1->Smooth();
903  m_fPXfit2->Smooth();
904  m_fPYfit2->Smooth();
905  m_fPZfit2->Smooth();
906 #endif
907 
908  // default max finding method defined in MissingMassCalculator.h
909  // note that window defined in terms of number of bin, so depend on binning
910  std::vector<double> histInfo(HistInfo::MAXHISTINFO);
912  double prob_hist = histInfo.at(HistInfo::PROB);
913 
914  if (prob_hist != 0.0)
915  m_fDitauStuffHisto.Sign_best = -log10(std::abs(prob_hist));
916  else {
917  // this mean the histogram is empty.
918  // possible but very rare if all entries outside histogram range
919  // fall back to maximum
922  }
923 
926  std::vector<double> histInfoOther(HistInfo::MAXHISTINFO);
927  //---- getting full tau1 momentum
928  Px1 = maxFromHist(m_fPXfit1, histInfoOther);
929  Py1 = maxFromHist(m_fPYfit1, histInfoOther);
930  Pz1 = maxFromHist(m_fPZfit1, histInfoOther);
931 
932  //---- getting full tau2 momentum
933  Px2 = maxFromHist(m_fPXfit2, histInfoOther);
934  Py2 = maxFromHist(m_fPYfit2, histInfoOther);
935  Pz2 = maxFromHist(m_fPZfit2, histInfoOther);
936 
937  //---- setting 4-vecs
938  PxPyPzMVector fulltau1, fulltau2;
939  fulltau1.SetCoordinates(Px1, Py1, Pz1, 1.777);
940  fulltau2.SetCoordinates(Px2, Py2, Pz2, 1.777);
941  // PtEtaPhiMVector fulltau1(_fulltau1.Pt(), _fulltau1.Eta(), _fulltau1.Phi(), _fulltau1.M());
942  //PtEtaPhiMVector fulltau2(_fulltau2.Pt(), _fulltau2.Eta(), _fulltau2.Phi(), _fulltau2.M());
943 
944  if (fulltau1.P() < preparedInput.m_vistau1.P())
945  fulltau1 = 1.01 * preparedInput.m_vistau1; // protection against cases when fitted tau
946  // momentum is smaller than visible tau momentum
947  if (fulltau2.P() < preparedInput.m_vistau2.P())
948  fulltau2 = 1.01 * preparedInput.m_vistau2; // protection against cases when fitted tau
949  // momentum is smaller than visible tau momentum
950  m_fDitauStuffHisto.vistau1 = preparedInput.m_vistau1; // FIXME should also be fitted if tau scan
952  m_fDitauStuffHisto.nutau1 = fulltau1 - preparedInput.m_vistau1; // these are the original tau vis
954  fulltau2 - preparedInput.m_vistau2; // FIXME neutrino mass not necessarily zero
955  }
956 
957  // Note that for v9walk, points outside the METx MEty disk are counted, while
958  // this was not the case for v9
959  if (preparedInput.m_fUseVerbose == 1) {
960  Info("DiTauMassTools", "Scanning ");
961  Info("DiTauMassTools", " Markov ");
962  Info("DiTauMassTools", "%s",
963  (" V9W niters=" + std::to_string(m_iter0) + " " + std::to_string(m_iter1)).c_str());
964  Info("DiTauMassTools", "%s", (" nFullScan " + std::to_string(m_markovNFullScan)).c_str());
965  Info("DiTauMassTools", "%s", (" nRejectNoSol " + std::to_string(m_markovNRejectNoSol)).c_str());
966  Info("DiTauMassTools", "%s", (" nRejectMetro " + std::to_string(m_markovNRejectMetropolis)).c_str());
967  Info("DiTauMassTools", "%s", (" nAccept " + std::to_string(m_markovNAccept)).c_str());
968  Info("DiTauMassTools", "%s",
969  (" probsum " + std::to_string(m_totalProbSum) + " msum " + std::to_string(m_mtautauSum))
970  .c_str());
971  }
972 
973  if (preparedInput.m_fUseVerbose == 1) {
974  if (fit_code == 0) {
975  Info("DiTauMassTools", "%s", ("!!!----> Warning-3 in "
976  "MissingMassCalculator::DitauMassCalculatorV9Walk() : fit status=" +
977  std::to_string(fit_code))
978  .c_str());
979  Info("DiTauMassTools", "%s", "....... No solution is found. Printing input info .......");
980 
981  Info("DiTauMassTools", "%s", (" vis Tau-1: Pt=" + std::to_string(preparedInput.m_vistau1.Pt()) +
982  " M=" + std::to_string(preparedInput.m_vistau1.M()) +
983  " eta=" + std::to_string(preparedInput.m_vistau1.Eta()) +
984  " phi=" + std::to_string(preparedInput.m_vistau1.Phi()) +
986  .c_str());
987  Info("DiTauMassTools", "%s", (" vis Tau-2: Pt=" + std::to_string(preparedInput.m_vistau2.Pt()) +
988  " M=" + std::to_string(preparedInput.m_vistau2.M()) +
989  " eta=" + std::to_string(preparedInput.m_vistau2.Eta()) +
990  " phi=" + std::to_string(preparedInput.m_vistau2.Phi()) +
992  .c_str());
993  Info("DiTauMassTools", "%s", (" MET=" + std::to_string(preparedInput.m_MetVec.R()) +
994  " Met_X=" + std::to_string(preparedInput.m_MetVec.X()) +
995  " Met_Y=" + std::to_string(preparedInput.m_MetVec.Y()))
996  .c_str());
997  Info("DiTauMassTools", " ---------------------------------------------------------- ");
998  }
999  }
1000 
1001  return fit_code;
1002 }
1003 
1005 
1006  // debugThisIteration=false;
1007  m_debugThisIteration = true;
1008 
1009  int fit_code = 0; // 0==bad, 1==good
1012  OutputInfo.m_NTrials = 0;
1014  OutputInfo.m_AveSolRMS = 0.;
1015 
1016  //------- Settings -------------------------------
1017  int NiterMET = m_niter_fit2; // number of iterations for each MET scan loop
1018  int NiterMnu = m_niter_fit3; // number of iterations for Mnu loop
1019  const double Mtau = 1.777;
1020  double Mnu_binSize = m_MnuScanRange / NiterMnu;
1021 
1022  double METresX = preparedInput.m_METsigmaL; // MET resolution in direction parallel to
1023  // leading jet, for MET scan
1024  double METresY = preparedInput.m_METsigmaP; // MET resolution in direction perpendicular to
1025  // leading jet, for MET scan
1026 
1027  //-------- end of Settings
1028 
1029  // if m_nsigma_METscan was not set by user, set to default values
1030  if(m_nsigma_METscan == -1){
1031  if (preparedInput.m_tauTypes == TauTypes::ll) { // both tau's are leptonic
1033  } else if (preparedInput.m_tauTypes == TauTypes::lh) { // lep had
1035  }
1036  }
1037 
1038  double N_METsigma = m_nsigma_METscan; // number of sigmas for MET scan
1039  double METresX_binSize = 2 * N_METsigma * METresX / NiterMET;
1040  double METresY_binSize = 2 * N_METsigma * METresY / NiterMET;
1041 
1042  int solution = 0;
1043 
1044  std::vector<PtEtaPhiMVector> nu_vec;
1045 
1046  m_totalProbSum = 0;
1047  m_mtautauSum = 0;
1048 
1049  double metprob = 1.0;
1050  double sign_tmp = 0.0;
1051  double tauprob = 1.0;
1052  double totalProb = 0.0;
1053 
1054  m_prob_tmp = 0.0;
1055 
1056  double met_smear_x = 0.0;
1057  double met_smear_y = 0.0;
1058  double met_smearL = 0.0;
1059  double met_smearP = 0.0;
1060 
1061  double angle1 = 0.0;
1062 
1063  if (m_fMfit_all) {
1064  m_fMfit_all->Reset();
1065  }
1066  if (m_fMfit_allNoWeight) {
1067  m_fMfit_allNoWeight->Reset();
1068  }
1069  if (m_fPXfit1) {
1070  m_fPXfit1->Reset();
1071  }
1072  if (m_fPYfit1) {
1073  m_fPYfit1->Reset();
1074  }
1075  if (m_fPZfit1) {
1076  m_fPZfit1->Reset();
1077  }
1078 
1079  int iter0 = 0;
1080  m_iter1 = 0;
1081  m_iter2 = 0;
1082  m_iter3 = 0;
1083  m_iter4 = 0;
1084 
1085  const double met_coscovphi = cos(preparedInput.m_METcovphi);
1086  const double met_sincovphi = sin(preparedInput.m_METcovphi);
1087 
1088  m_iang1low = 0;
1089  m_iang1high = 0;
1090 
1091  // double Mvis=(tau_vec1+tau_vec2).M();
1092  // PtEtaPhiMVector met4vec(0.0,0.0,0.0,0.0);
1093  // met4vec.SetPxPyPzE(met_vec.X(),met_vec.Y(),0.0,met_vec.R());
1094  // double Meff=(tau_vec1+tau_vec2+met4vec).M();
1095  // double met_det=met_vec.R();
1096 
1097  //---------------------------------------------
1098  if (preparedInput.m_tauTypes == TauTypes::ll) // dilepton case
1099  {
1100  if (preparedInput.m_fUseVerbose == 1) {
1101  Info("DiTauMassTools", "Running in dilepton mode");
1102  }
1103  double input_metX = preparedInput.m_MetVec.X();
1104  double input_metY = preparedInput.m_MetVec.Y();
1105 
1106  PtEtaPhiMVector tau_tmp(0.0, 0.0, 0.0, 0.0);
1107  PtEtaPhiMVector lep_tmp(0.0, 0.0, 0.0, 0.0);
1108  int tau_type_tmp;
1109  int tau_ind = 0;
1110 
1111  if (preparedInput.m_LFVmode == 1) // muon case: H->mu+tau(->ele) decays
1112  {
1113  if ((preparedInput.m_vistau1.M() > 0.05 &&
1114  preparedInput.m_vistau2.M() < 0.05) != refit) // choosing lepton from Higgs decay
1115  //When the mass calculator is rerun with refit==true the alternative lepton ordering is used
1116  {
1117  tau_tmp = preparedInput.m_vistau2;
1118  lep_tmp = preparedInput.m_vistau1;
1119  tau_type_tmp = preparedInput.m_type_visTau2;
1120  tau_ind = 2;
1121  } else {
1122  tau_tmp = preparedInput.m_vistau1;
1123  lep_tmp = preparedInput.m_vistau2;
1124  tau_type_tmp = preparedInput.m_type_visTau1;
1125  tau_ind = 1;
1126  }
1127  }
1128  if (preparedInput.m_LFVmode == 0) // electron case: H->ele+tau(->mu) decays
1129  {
1130  if ((preparedInput.m_vistau1.M() < 0.05 &&
1131  preparedInput.m_vistau2.M() > 0.05) != refit) // choosing lepton from Higgs decay
1132  //When the mass calculator is rerun with refit=true the alternative lepton ordering is used
1133  {
1134  tau_tmp = preparedInput.m_vistau2;
1135  lep_tmp = preparedInput.m_vistau1;
1136  tau_type_tmp = preparedInput.m_type_visTau2;
1137  tau_ind = 2;
1138  } else {
1139  tau_tmp = preparedInput.m_vistau1;
1140  lep_tmp = preparedInput.m_vistau2;
1141  tau_type_tmp = preparedInput.m_type_visTau1;
1142  tau_ind = 1;
1143  }
1144  }
1145 
1146  //------- Settings -------------------------------
1147  double Mlep = tau_tmp.M();
1148  // double dMnu_max=m_MnuScanRange-Mlep;
1149  // double Mnu_binSize=dMnu_max/NiterMnu;
1150  //-------- end of Settings
1151 
1152  // double M=Mtau;
1153  double M_nu = 0.0;
1154  double MnuProb = 1.0;
1155  //---------------------------------------------
1156  for (int i3 = 0; i3 < NiterMnu; i3++) //---- loop-3: virtual neutrino mass
1157  {
1158  M_nu = Mnu_binSize * i3;
1159  if (M_nu >= (Mtau - Mlep))
1160  continue;
1161  // M=sqrt(Mtau*Mtau-M_nu*M_nu);
1162  MnuProb = Prob->MnuProbability(preparedInput, M_nu,
1163  Mnu_binSize); // Mnu probability
1164  //---------------------------------------------
1165  for (int i4 = 0; i4 < NiterMET + 1; i4++) // MET_X scan
1166  {
1167  met_smearL = METresX_binSize * i4 - N_METsigma * METresX;
1168  for (int i5 = 0; i5 < NiterMET + 1; i5++) // MET_Y scan
1169  {
1170  met_smearP = METresY_binSize * i5 - N_METsigma * METresY;
1171  if (pow(met_smearL / METresX, 2) + pow(met_smearP / METresY, 2) > pow(N_METsigma, 2))
1172  continue; // use ellipse instead of square
1173  met_smear_x = met_smearL * met_coscovphi - met_smearP * met_sincovphi;
1174  met_smear_y = met_smearL * met_sincovphi + met_smearP * met_coscovphi;
1175  metvec_tmp.SetXY(input_metX + met_smear_x, input_metY + met_smear_y);
1176 
1177  solution = NuPsolutionLFV(metvec_tmp, tau_tmp, M_nu, nu_vec);
1178 
1179  ++iter0;
1180 
1181  if (solution < 1)
1182  continue;
1183  ++m_iter1;
1184 
1185  // if fast sin cos, result to not match exactly nupsolutionv2, so skip
1186  // test
1187  // SpeedUp no nested loop to compute individual probability
1188  int ngoodsol1 = 0;
1189 
1190  metprob = Prob->MetProbability(preparedInput, met_smearL, met_smearP, METresX, METresY);
1191  if (metprob <= 0)
1192  continue;
1193  for (unsigned int j1 = 0; j1 < nu_vec.size(); j1++) {
1194  if (tau_tmp.E() + nu_vec[j1].E() >= preparedInput.m_beamEnergy)
1195  continue;
1196  const double tau1_tmpp = (tau_tmp + nu_vec[j1]).P();
1197  angle1 = Angle(nu_vec[j1], tau_tmp);
1198 
1199  if (angle1 < dTheta3DLimit(tau_type_tmp, 0, tau1_tmpp)) {
1200  ++m_iang1low;
1201  continue;
1202  } // lower 99% bound
1203  if (angle1 > dTheta3DLimit(tau_type_tmp, 1, tau1_tmpp)) {
1204  ++m_iang1high;
1205  continue;
1206  } // upper 99% bound
1207  double tauvecprob1j =
1208  Prob->dTheta3d_probabilityFast(preparedInput, tau_type_tmp, angle1, tau1_tmpp);
1209  if (tauvecprob1j == 0.)
1210  continue;
1211  tauprob = Prob->TauProbabilityLFV(preparedInput, tau_type_tmp, tau_tmp, nu_vec[j1]);
1212  totalProb = tauvecprob1j * metprob * MnuProb * tauprob;
1213 
1214  m_tautau_tmp.SetPxPyPzE(0.0, 0.0, 0.0, 0.0);
1215  m_tautau_tmp += tau_tmp;
1216  m_tautau_tmp += lep_tmp;
1217  m_tautau_tmp += nu_vec[j1];
1218 
1219  const double mtautau = m_tautau_tmp.M();
1220 
1221  m_totalProbSum += totalProb;
1222  m_mtautauSum += mtautau;
1223 
1224  fit_code = 1; // at least one solution is found
1225 
1226  m_fMfit_all->Fill(mtautau, totalProb);
1227  m_fMfit_allNoWeight->Fill(mtautau, 1.);
1228  //----------------- using P*fit to fill Px,y,z_tau
1229  m_fPXfit1->Fill((tau_tmp + nu_vec[j1]).Px(), totalProb);
1230  m_fPYfit1->Fill((tau_tmp + nu_vec[j1]).Py(), totalProb);
1231  m_fPZfit1->Fill((tau_tmp + nu_vec[j1]).Pz(), totalProb);
1232 
1233  if (totalProb > m_prob_tmp) // fill solution with highest probability
1234  {
1235  sign_tmp = -log10(totalProb);
1236  m_prob_tmp = totalProb;
1237  m_fDitauStuffFit.Mditau_best = mtautau;
1238  m_fDitauStuffFit.Sign_best = sign_tmp;
1239  if (tau_ind == 1)
1240  m_fDitauStuffFit.nutau1 = nu_vec[j1];
1241  if (tau_ind == 2)
1242  m_fDitauStuffFit.nutau2 = nu_vec[j1];
1243  }
1244 
1245  ++ngoodsol1;
1246  }
1247 
1248  if (ngoodsol1 == 0)
1249  continue;
1250  m_iter2 += 1;
1251 
1252  m_iter3 += 1;
1253  }
1254  }
1255  }
1256  } else if (preparedInput.m_tauTypes == TauTypes::lh) // lepton+tau case
1257  {
1258  if (preparedInput.m_fUseVerbose == 1) {
1259  Info("DiTauMassTools", "Running in lepton+tau mode");
1260  }
1261  //------- Settings -------------------------------
1262 
1263  //----- Stuff below are for Winter 2012 lep-had analysis only; it has to be
1264  // replaced by a more common scheme once other channels are optimized
1265  // XYVector
1266  // mht_vec((tau_vec1+tau_vec2).Px(),(tau_vec1+tau_vec2).Py()); //
1267  // missing Ht vector for Njet25=0 events const double
1268  // mht=mht_vec.R();
1269  double input_metX = preparedInput.m_MetVec.X();
1270  double input_metY = preparedInput.m_MetVec.Y();
1271 
1272  // double mht_offset=0.0;
1273  // if(InputInfo.UseHT) // use missing Ht (for 0-jet events only for
1274  // now)
1275  // {
1276  // input_metX=-mht_vec.X();
1277  // input_metY=-mht_vec.Y();
1278  // }
1279  // else // use MET (for 0-jet and 1-jet events)
1280  // {
1281  // input_metX=met_vec.X();
1282  // input_metY=met_vec.Y();
1283  // }
1284 
1285  PtEtaPhiMVector tau_tmp(0.0, 0.0, 0.0, 0.0);
1286  PtEtaPhiMVector lep_tmp(0.0, 0.0, 0.0, 0.0);
1287  int tau_type_tmp;
1288  if (preparedInput.m_type_visTau1 == 8) {
1289  tau_tmp = preparedInput.m_vistau2;
1290  lep_tmp = preparedInput.m_vistau1;
1291  tau_type_tmp = preparedInput.m_type_visTau2;
1292  }
1293  if (preparedInput.m_type_visTau2 == 8) {
1294  tau_tmp = preparedInput.m_vistau1;
1295  lep_tmp = preparedInput.m_vistau2;
1296  tau_type_tmp = preparedInput.m_type_visTau1;
1297  }
1298 
1299  //---------------------------------------------
1300  for (int i4 = 0; i4 < NiterMET + 1; i4++) // MET_X scan
1301  {
1302  met_smearL = METresX_binSize * i4 - N_METsigma * METresX;
1303  for (int i5 = 0; i5 < NiterMET + 1; i5++) // MET_Y scan
1304  {
1305  met_smearP = METresY_binSize * i5 - N_METsigma * METresY;
1306  if (pow(met_smearL / METresX, 2) + pow(met_smearP / METresY, 2) > pow(N_METsigma, 2))
1307  continue; // use ellipse instead of square
1308  met_smear_x = met_smearL * m_metCovPhiCos - met_smearP * m_metCovPhiSin;
1309  met_smear_y = met_smearL * m_metCovPhiSin + met_smearP * m_metCovPhiCos;
1310  metvec_tmp.SetXY(input_metX + met_smear_x, input_metY + met_smear_y);
1311 
1312  solution = NuPsolutionLFV(metvec_tmp, tau_tmp, 0.0, nu_vec);
1313 
1314  ++iter0;
1315 
1316  if (solution < 1)
1317  continue;
1318  ++m_iter1;
1319 
1320  // if fast sin cos, result to not match exactly nupsolutionv2, so skip
1321  // test
1322  // SpeedUp no nested loop to compute individual probability
1323  int ngoodsol1 = 0;
1324 
1325  metprob = Prob->MetProbability(preparedInput, met_smearL, met_smearP, METresX, METresY);
1326  if (metprob <= 0)
1327  continue;
1328  for (unsigned int j1 = 0; j1 < nu_vec.size(); j1++) {
1329  if (tau_tmp.E() + nu_vec[j1].E() >= preparedInput.m_beamEnergy)
1330  continue;
1331  const double tau1_tmpp = (tau_tmp + nu_vec[j1]).P();
1332  angle1 = Angle(nu_vec[j1], tau_tmp);
1333 
1334  if (angle1 < dTheta3DLimit(tau_type_tmp, 0, tau1_tmpp)) {
1335  ++m_iang1low;
1336  continue;
1337  } // lower 99% bound
1338  if (angle1 > dTheta3DLimit(tau_type_tmp, 1, tau1_tmpp)) {
1339  ++m_iang1high;
1340  continue;
1341  } // upper 99% bound
1342  double tauvecprob1j =
1343  Prob->dTheta3d_probabilityFast(preparedInput, tau_type_tmp, angle1, tau1_tmpp);
1344  if (tauvecprob1j == 0.)
1345  continue;
1346  tauprob = Prob->TauProbabilityLFV(preparedInput, tau_type_tmp, tau_tmp, nu_vec[j1]);
1347  totalProb = tauvecprob1j * metprob * tauprob;
1348 
1349  m_tautau_tmp.SetPxPyPzE(0.0, 0.0, 0.0, 0.0);
1350  m_tautau_tmp += tau_tmp;
1351  m_tautau_tmp += lep_tmp;
1352  m_tautau_tmp += nu_vec[j1];
1353 
1354  const double mtautau = m_tautau_tmp.M();
1355 
1356  m_totalProbSum += totalProb;
1357  m_mtautauSum += mtautau;
1358 
1359  fit_code = 1; // at least one solution is found
1360 
1361  m_fMfit_all->Fill(mtautau, totalProb);
1362  m_fMfit_allNoWeight->Fill(mtautau, 1.);
1364  // m_fPXfit1->Fill((tau_tmp+nu_vec[j1]).Px(),totalProb);
1365  // m_fPYfit1->Fill((tau_tmp+nu_vec[j1]).Py(),totalProb);
1366  // m_fPZfit1->Fill((tau_tmp+nu_vec[j1]).Pz(),totalProb);
1367 
1368  if (totalProb > m_prob_tmp) // fill solution with highest probability
1369  {
1370  sign_tmp = -log10(totalProb);
1371  m_prob_tmp = totalProb;
1372  m_fDitauStuffFit.Mditau_best = mtautau;
1373  m_fDitauStuffFit.Sign_best = sign_tmp;
1374  if (preparedInput.m_type_visTau1 == 8) {
1375  m_fDitauStuffFit.vistau1 = lep_tmp;
1376  m_fDitauStuffFit.vistau2 = tau_tmp;
1377  m_fDitauStuffFit.nutau2 = nu_vec[j1];
1378  } else if (preparedInput.m_type_visTau2 == 8) {
1379  m_fDitauStuffFit.vistau2 = lep_tmp;
1380  m_fDitauStuffFit.vistau1 = tau_tmp;
1381  m_fDitauStuffFit.nutau1 = nu_vec[j1];
1382  }
1383  }
1384 
1385  ++ngoodsol1;
1386  }
1387 
1388  if (ngoodsol1 == 0)
1389  continue;
1390  m_iter2 += 1;
1391 
1392  m_iter3 += 1;
1393  }
1394  }
1395  } else {
1396  Info("DiTauMassTools", "Running in an unknown mode?!?!");
1397  }
1398 
1399  OutputInfo.m_NTrials = iter0;
1401 
1402  if (preparedInput.m_fUseVerbose == 1) {
1403  Info("DiTauMassTools", "%s",
1404  ("SpeedUp niters=" + std::to_string(iter0) + " " + std::to_string(m_iter1) + " " +
1406  " " + std::to_string(m_iang1high))
1407  .c_str());
1408  }
1409 
1410  if (m_fMfit_all->GetEntries() > 0 && m_iter3 > 0) {
1411 #ifdef SMOOTH
1412  m_fMfit_all->Smooth();
1413  m_fMfit_allNoWeight->Smooth();
1414  m_fPXfit1->Smooth();
1415  m_fPYfit1->Smooth();
1416  m_fPZfit1->Smooth();
1417 #endif
1418 
1419  // default max finding method defined in MissingMassCalculator.h
1420  // note that window defined in terms of number of bin, so depend on binning
1421  std::vector<double> histInfo(HistInfo::MAXHISTINFO);
1423  double prob_hist = histInfo.at(HistInfo::PROB);
1424 
1425  if (prob_hist != 0.0)
1426  m_fDitauStuffHisto.Sign_best = -log10(std::abs(prob_hist));
1427  else {
1428  // this mean the histogram is empty.
1429  // possible but very rare if all entries outside histogram range
1430  // fall back to maximum
1431  m_fDitauStuffHisto.Sign_best = -999.;
1433  }
1434 
1435  if (m_fDitauStuffHisto.Mditau_best > 0.0)
1437  std::vector<double> histInfoOther(HistInfo::MAXHISTINFO);
1438  //---- getting Nu1
1439  double Px1 = maxFromHist(m_fPXfit1, histInfoOther);
1440  double Py1 = maxFromHist(m_fPYfit1, histInfoOther);
1441  double Pz1 = maxFromHist(m_fPZfit1, histInfoOther);
1442  //---- setting 4-vecs
1443  PxPyPzMVector nu1_tmp(0.0, 0.0, 0.0, 0.0);
1444  PxPyPzMVector nu2_tmp(0.0, 0.0, 0.0, 0.0);
1445  if (preparedInput.m_type_visTau1 == 8) {
1446  nu1_tmp = preparedInput.m_vistau1;
1447  nu2_tmp.SetCoordinates(Px1, Py1, Pz1, 1.777);
1448  }
1449  if (preparedInput.m_type_visTau2 == 8) {
1450  nu2_tmp = preparedInput.m_vistau2;
1451  nu1_tmp.SetCoordinates(Px1, Py1, Pz1, 1.777);
1452  }
1455  }
1456  if (m_lfvLeplepRefit && fit_code==0 && !refit) {
1457  fit_code = DitauMassCalculatorV9lfv(true);
1458  return fit_code;
1459  }
1460 
1461 
1462 
1463  if (preparedInput.m_fUseVerbose == 1) {
1464  if (fit_code == 0) {
1465  Info(
1466  "DiTauMassTools", "%s",
1467  ("!!!----> Warning-3 in MissingMassCalculator::DitauMassCalculatorV9lfv() : fit status=" +
1468  std::to_string(fit_code))
1469  .c_str());
1470  Info("DiTauMassTools", "....... No solution is found. Printing input info .......");
1471 
1472  Info("DiTauMassTools", "%s", (" vis Tau-1: Pt="+std::to_string(preparedInput.m_vistau1.Pt())
1474  +" phi="+std::to_string(preparedInput.m_vistau1.Phi())
1475  +" type="+std::to_string(preparedInput.m_type_visTau1)).c_str());
1476  Info("DiTauMassTools", "%s", (" vis Tau-2: Pt="+std::to_string(preparedInput.m_vistau2.Pt())
1478  +" phi="+std::to_string(preparedInput.m_vistau2.Phi())
1479  +" type="+std::to_string(preparedInput.m_type_visTau2)).c_str());
1480  Info("DiTauMassTools", "%s", (" MET="+std::to_string(preparedInput.m_MetVec.R())+" Met_X="+std::to_string(preparedInput.m_MetVec.X())
1481  +" Met_Y="+std::to_string(preparedInput.m_MetVec.Y())).c_str());
1482  Info("DiTauMassTools", " ---------------------------------------------------------- ");
1483  }
1484  }
1485  return fit_code;
1486 }
1487 
1488 // function to fit maximum
1489 Double_t MissingMassCalculator::maxFitting(Double_t *x, Double_t *par)
1490 // Double_t maxFitting(Double_t *x, Double_t *par)
1491 {
1492  // parabola with parameters max, mean and invwidth
1493  const double mM = x[0];
1494  const double mMax = par[0];
1495  const double mMean = par[1];
1496  const double mInvWidth2 = par[2]; // if param positif distance between intersection of the
1497  // parabola with x axis: 1/Sqrt(mInvWidth2)
1498  const double fitval = mMax * (1 - 4 * mInvWidth2 * std::pow(mM - mMean, 2));
1499  return fitval;
1500 }
1501 
1502 // determine the maximum from the histogram
1503 // if input prob not default , compute also some probability
1504 // MaxHistStrategy : different method to find maximum
1505 // TODO should get the array on work on it
1506 // should also find the effective range of the hist
1507 
1508 double
1509 MissingMassCalculator::maxFromHist(TH1F *theHist, std::vector<double> &histInfo,
1510  const MaxHistStrategy::e maxHistStrategy,
1511  const int winHalfWidth, bool debug) {
1512  // namespace HistInfo
1513  // enum e {
1514  // PROB=0,INTEGRAL,CHI2,DISCRI,TANTHETA,TANTHETAW,FITLENGTH,RMS,RMSVSDISCRI,MAXHISTINFO
1515  // };
1516  double maxPos = 0.;
1517  double prob = 0.;
1518 
1519  for (std::vector<double>::iterator itr = histInfo.begin(); itr != histInfo.end(); ++itr) {
1520  *itr = -1;
1521  }
1522 
1523  histInfo[HistInfo::INTEGRAL] = theHist->Integral();
1524 
1525  if (maxHistStrategy == MaxHistStrategy::MAXBIN ||
1526  ((maxHistStrategy == MaxHistStrategy::MAXBINWINDOW ||
1527  maxHistStrategy == MaxHistStrategy::SLIDINGWINDOW) &&
1528  winHalfWidth == 0)) {
1529 
1530  // simple max search
1531  // original version, simple bin maximum
1532  int max_bin = theHist->GetMaximumBin();
1533  maxPos = theHist->GetBinCenter(max_bin);
1534 
1535  // FIXME GetEntries is unweighted
1536  prob = theHist->GetBinContent(max_bin) / double(theHist->GetEntries());
1537  if (prob > 1.)
1538  prob = 1.;
1539  histInfo[HistInfo::PROB] = prob;
1540  return maxPos;
1541  }
1542 
1543  int hNbins = theHist->GetNbinsX();
1544 
1545  if (maxHistStrategy == MaxHistStrategy::MAXBINWINDOW) {
1546  // average around maximum bin (nearly useless in fact)
1547  // could be faster
1548  int max_bin = theHist->GetMaximumBin();
1549  int iBinMin = max_bin - winHalfWidth;
1550  if (iBinMin < 0)
1551  iBinMin = 0;
1552  int iBinMax = max_bin + winHalfWidth;
1553  if (iBinMax > hNbins)
1554  iBinMax = hNbins - 1;
1555  double sumw = 0;
1556  double sumx = 0;
1557  for (int iBin = iBinMin; iBin <= iBinMax; ++iBin) {
1558  const double weight = theHist->GetBinContent(iBin);
1559  sumw += weight;
1560  sumx += weight * theHist->GetBinCenter(iBin);
1561  }
1562  maxPos = sumx / sumw;
1563 
1564  // FIXME GetEntries is unweighted
1565  prob = sumw / theHist->GetEntries();
1566  if (prob > 1.)
1567  prob = 1.;
1568 
1569  return maxPos;
1570  }
1571 
1572  // now compute sliding window anyway
1573  if (maxHistStrategy != MaxHistStrategy::SLIDINGWINDOW &&
1574  maxHistStrategy != MaxHistStrategy::FIT) {
1575  Error("DiTauMassTools", "%s",
1576  ("ERROR undefined maxHistStrategy:" + std::to_string(maxHistStrategy)).c_str());
1577  return -10.;
1578  }
1579 
1580  // first iteration to find the first and last non zero bin, and the histogram
1581  // integral (not same as Entries because of weights)
1582  int lastNonZeroBin = -1;
1583  int firstNonZeroBin = -1;
1584  double totalSumw = 0.;
1585  bool firstNullPart = true;
1586  for (int iBin = 0; iBin < hNbins; ++iBin) {
1587  const double weight = theHist->GetBinContent(iBin);
1588  if (weight > 0) {
1589  totalSumw += weight;
1590  lastNonZeroBin = iBin;
1591  if (firstNullPart) {
1592  firstNullPart = false;
1593  firstNonZeroBin = iBin;
1594  }
1595  }
1596  }
1597 
1598  // enlarge first and last non zero bin with window width to avoid side effect
1599  // (maximum close to the edge)
1600  firstNonZeroBin = std::max(0, firstNonZeroBin - winHalfWidth - 1);
1601  lastNonZeroBin = std::min(hNbins - 1, lastNonZeroBin + winHalfWidth + 1);
1602 
1603  // if null histogram quit
1604  if (firstNullPart)
1605  return maxPos;
1606 
1607  // determine the size of the sliding window in the fit case
1608 
1609  // sliding window
1610  const int nwidth = 2 * winHalfWidth + 1;
1611  double winsum = 0.;
1612 
1613  for (int ibin = 0; ibin < nwidth; ++ibin) {
1614  winsum += theHist->GetBinContent(ibin);
1615  }
1616  double winmax = winsum;
1617 
1618  int max_bin = 0.;
1619  int iBinL = firstNonZeroBin;
1620  int iBinR = iBinL + 2 * winHalfWidth;
1621  bool goingUp = true;
1622 
1623  do {
1624  ++iBinL;
1625  ++iBinR;
1626  const double deltawin = theHist->GetBinContent(iBinR) - theHist->GetBinContent(iBinL - 1);
1627 
1628  if (deltawin < 0) {
1629  if (goingUp) {
1630  // if were climbing and now loose more on the left
1631  // than win on the right. This was a local maxima
1632  if (winsum > winmax) {
1633  // global maximum one so far
1634  winmax = winsum;
1635  max_bin = (iBinR + iBinL) / 2 - 1;
1636  }
1637  goingUp = false; // now going down
1638  }
1639  } else {
1640  // do not care about minima, simply indicate we are going down
1641  goingUp = true;
1642  }
1643 
1644  winsum += deltawin;
1645 
1646  } while (iBinR < lastNonZeroBin);
1647 
1648  // now compute average
1649  int iBinMin = max_bin - winHalfWidth;
1650  if (iBinMin < 0)
1651  iBinMin = 0;
1652  int iBinMax = max_bin + winHalfWidth;
1653  if (iBinMax >= hNbins)
1654  iBinMax = hNbins - 1;
1655  double sumw = 0;
1656  double sumx = 0;
1657  for (int iBin = iBinMin; iBin <= iBinMax; ++iBin) {
1658  const double weight = theHist->GetBinContent(iBin);
1659  sumw += weight;
1660  sumx += weight * theHist->GetBinCenter(iBin);
1661  }
1662 
1663  double maxPosWin = -1.;
1664 
1665  if (sumw > 0.) {
1666  maxPosWin = sumx / sumw;
1667  }
1668  // prob if the fraction of events in the window
1669  prob = sumw / totalSumw;
1670 
1671  // Definitions of some useful parameters
1672 
1673  const double h_rms = theHist->GetRMS(1);
1674  histInfo[HistInfo::RMS] = h_rms;
1675 
1676  double num = 0;
1677  double numerator = 0;
1678  double denominator = 0;
1679  bool nullBin = false;
1680 
1681  for (int i = iBinMin; i < iBinMax; ++i) {
1682  double binError = theHist->GetBinError(i);
1683  if (binError < 1e-10) {
1684  nullBin = true;
1685  }
1686  double binErrorSquare = std::pow(binError, 2);
1687  num = theHist->GetBinContent(i) / (binErrorSquare);
1688  numerator = numerator + num;
1689  denominator = denominator + (1 / (binErrorSquare));
1690  }
1691  if (numerator < 1e-10 || denominator < 1e-10 || nullBin == true) {
1692  histInfo[HistInfo::MEANBIN] = -1;
1693  } else {
1694  histInfo[HistInfo::MEANBIN] = sqrt(1 / denominator) / (numerator / denominator);
1695  }
1696 
1697  // stop here if only looking for sliding window
1698  if (maxHistStrategy == MaxHistStrategy::SLIDINGWINDOW) {
1699  return maxPosWin;
1700  }
1701 
1702  maxPos = maxPosWin;
1703  // now FIT maxHistStrategy==MaxHistStrategy::FIT
1704 
1705  // now mass fit in range defined by sliding window
1706  // window will be around maxPos
1707  const double binWidth = theHist->GetBinCenter(2) - theHist->GetBinCenter(1);
1708  double fitWidth = (winHalfWidth + 0.5) * binWidth;
1709  // fit range 2 larger than original window range, 3 if less than 20% of the
1710  // histogram in slinding window
1711 
1712  if (prob > 0.2) {
1713  fitWidth *= 2.;
1714  } else {
1715  fitWidth *= 3.;
1716  }
1717  // fit option : Q == Quiet, no printout S result of the fit returned in
1718  // TFitResultPtr N do not draw the resulting function
1719 
1720  // if debug plot the fitted function
1721  TString fitOption = debug ? "QS" : "QNS";
1722  // root fit
1723  // Sets initial values
1724  m_fFitting->SetParameters(sumw / winHalfWidth, maxPos, 0.0025);
1725  // TFitResultPtr
1726  // fitRes=theHist->Fit("pol2",fitOption,"",maxPos-fitWidth,maxPos+fitWidth);
1727  TFitResultPtr fitRes =
1728  theHist->Fit(m_fFitting, fitOption, "", maxPos - fitWidth, maxPos + fitWidth);
1729 
1730  double maxPosFit = -1.;
1731 
1732  if (int(fitRes) == 0) {
1733  // root fit
1734  histInfo[HistInfo::CHI2] = fitRes->Chi2();
1735  const double mMax = fitRes->Parameter(0);
1736  const double mMean = fitRes->Parameter(1);
1737  const double mInvWidth2 = fitRes->Parameter(2);
1738  double mMaxError = fitRes->ParError(0);
1739  m_PrintmMaxError = mMaxError;
1740  double mMeanError = fitRes->ParError(1);
1741  m_PrintmMeanError = mMeanError;
1742  double mInvWidth2Error = fitRes->ParError(2);
1743  m_PrintmInvWidth2Error = mInvWidth2Error;
1744  mMeanError = 0.; // avoid warning
1745  mInvWidth2Error = 0.; // avoid warning
1746  const double c = mMax * (1 - 4 * mMean * mMean * mInvWidth2);
1747  const double b = 8 * mMax * mMean * mInvWidth2;
1748  const double a = -4 * mMax * mInvWidth2;
1749  // when built in polynomial fit
1750  // const double c=fitRes->Parameter(0);
1751  // const double b=fitRes->Parameter(1);
1752  // const double a=fitRes->Parameter(2);
1753 
1754  const double h_discri = b * b - 4 * a * c;
1755  histInfo[HistInfo::DISCRI] = h_discri;
1756  const double sqrth_discri = sqrt(h_discri);
1757  const double h_fitLength = sqrth_discri / a;
1758  histInfo[HistInfo::FITLENGTH] = h_fitLength;
1759  histInfo[HistInfo::TANTHETA] = 2 * a / sqrth_discri;
1760  histInfo[HistInfo::TANTHETAW] = 2 * a * sumw / sqrth_discri;
1761  histInfo[HistInfo::RMSVSDISCRI] = h_rms / h_fitLength;
1762  // compute maximum position (only if inverted parabola)
1763  if (a < 0)
1764  maxPosFit = -b / (2 * a);
1765  }
1766 
1767  // keep fit result only if within 80% of fit window, and fit succeeded
1768  if (maxPosFit >= 0. and std::abs(maxPosFit - maxPosWin) < 0.8 * fitWidth) {
1769  histInfo[HistInfo::PROB] = prob;
1770  return maxPosFit;
1771  } else {
1772  // otherwise keep the weighted average
1773  // negate prob just to flag such event
1774  prob = -prob;
1775  histInfo[HistInfo::PROB] = prob;
1776  return maxPosWin;
1777  }
1778 }
1779 
1780 // compute probability for any input value,can be called from a pure parameter
1781 // scan
1782 // deltametvec is along phijet
1783 // returns number of solution if positive, return code if negative, vector of
1784 // probability and mass
1785 int MissingMassCalculator::probCalculatorV9fast(const double &phi1, const double &phi2,
1786  const double &M_nu1,
1787  const double &M_nu2) {
1788  // bool debug=true;
1789 
1790  int nsol1;
1791  int nsol2;
1792 
1793  const int solution = NuPsolutionV3(M_nu1, M_nu2, phi1, phi2, nsol1, nsol2);
1794 
1795  if (solution != 1)
1796  return -4;
1797  // refineSolutions ( M_nu1,M_nu2,
1798  // met_smearL,met_smearP,metvec_tmp.R(),
1799  // nsol1, nsol2,m_Mvis,m_Meff);
1800  refineSolutions(M_nu1, M_nu2, nsol1, nsol2, m_Mvis, m_Meff);
1801 
1802  if (m_nsol <= 0)
1803  return 0;
1804 
1805  // success
1806 
1807  return m_nsol; // for backward compatibility
1808 }
1809 
1810 // nuvecsol1 and nuvecsol2 passed by MMC
1811 int MissingMassCalculator::refineSolutions(const double &M_nu1, const double &M_nu2,
1812  const int nsol1, const int nsol2,
1813  const double &Mvis, const double &Meff)
1814 
1815 {
1816  m_nsol = 0;
1817 
1818  if (int(m_probFinalSolVec.size()) < m_nsolfinalmax)
1819  Error("DiTauMassTools", "%s",
1820  ("refineSolutions ERROR probFinalSolVec.size() should be " + std::to_string(m_nsolfinalmax))
1821  .c_str());
1822  if (int(m_mtautauFinalSolVec.size()) < m_nsolfinalmax)
1823  Error("DiTauMassTools", "%s",
1824  ("refineSolutions ERROR mtautauSolVec.size() should be " + std::to_string(m_nsolfinalmax))
1825  .c_str());
1826  if (int(m_nu1FinalSolVec.size()) < m_nsolfinalmax)
1827  Error("DiTauMassTools", "%s",
1828  ("refineSolutions ERROR nu1FinalSolVec.size() should be " + std::to_string(m_nsolfinalmax))
1829  .c_str());
1830  if (int(m_nu2FinalSolVec.size()) < m_nsolfinalmax)
1831  Error("DiTauMassTools", "%s",
1832  ("refineSolutions ERROR nu2FinalSolVec.size() should be " + std::to_string(m_nsolfinalmax))
1833  .c_str());
1834  if (nsol1 > int(m_nsolmax))
1835  Error("DiTauMassTools", "%s", ("refineSolutions ERROR nsol1 " + std::to_string(nsol1) +
1836  " > nsolmax !" + std::to_string(m_nsolmax))
1837  .c_str());
1838  if (nsol2 > int(m_nsolmax))
1839  Error("DiTauMassTools", "%s", ("refineSolutions ERROR nsol1 " + std::to_string(nsol2) +
1840  " > nsolmax !" + std::to_string(m_nsolmax))
1841  .c_str());
1842 
1843  int ngoodsol1 = 0;
1844  int ngoodsol2 = 0;
1845  double constProb =
1846  Prob->apply(preparedInput, -99, -99, PtEtaPhiMVector(0, 0, 0, 0), PtEtaPhiMVector(0, 0, 0, 0),
1847  PtEtaPhiMVector(0, 0, 0, 0), PtEtaPhiMVector(0, 0, 0, 0), true, false, false);
1848 
1849  for (int j1 = 0; j1 < nsol1; ++j1) {
1850  PtEtaPhiMVector &nuvec1_tmpj = m_nuvecsol1[j1];
1851  PtEtaPhiMVector &tauvecsol1j = m_tauvecsol1[j1];
1852  double &tauvecprob1j = m_tauvecprob1[j1];
1853  tauvecprob1j = 0.;
1854  // take first or second solution
1855  // no time to call rndm, switch more or less randomely, according to an
1856  // oscillating switch perturbed by m_phi1
1857  if (nsol1 > 1) {
1858  if (j1 == 0) { // decide at the first solution which one we will take
1859  const int pickInt = std::abs(10000 * m_Phi1);
1860  const int pickDigit = pickInt - 10 * (pickInt / 10);
1861  if (pickDigit < 5)
1862  m_switch1 = !m_switch1;
1863  }
1864  m_switch1 = !m_switch1;
1865  }
1866 
1867  if (!m_switch1) {
1868  nuvec1_tmpj.SetCoordinates(nuvec1_tmpj.Pt(), nuvec1_tmpj.Eta(), nuvec1_tmpj.Phi(), M_nu1);
1869  tauvecsol1j.SetPxPyPzE(0., 0., 0., 0.);
1870  tauvecsol1j += nuvec1_tmpj;
1871  tauvecsol1j += m_tauVec1;
1872  if (tauvecsol1j.E() >= preparedInput.m_beamEnergy)
1873  continue;
1875  PtEtaPhiMVector(0, 0, 0, 0), nuvec1_tmpj,
1876  PtEtaPhiMVector(0, 0, 0, 0), false, true, false);
1877  ++ngoodsol1;
1878  }
1879 
1880  for (int j2 = 0; j2 < nsol2; ++j2) {
1881  PtEtaPhiMVector &nuvec2_tmpj = m_nuvecsol2[j2];
1882  PtEtaPhiMVector &tauvecsol2j = m_tauvecsol2[j2];
1883  double &tauvecprob2j = m_tauvecprob2[j2];
1884  if (j1 == 0) {
1885  tauvecprob2j = 0.;
1886  // take first or second solution
1887  // no time to call rndm, switch more or less randomely, according to an
1888  // oscillating switch perturbed by m_phi2
1889  if (nsol2 > 1) {
1890  if (j2 == 0) { // decide at the first solution which one we will take
1891  const int pickInt = std::abs(10000 * m_Phi2);
1892  const int pickDigit = pickInt - 10 * int(pickInt / 10);
1893  if (pickDigit < 5)
1894  m_switch2 = !m_switch2;
1895  }
1896  m_switch2 = !m_switch2;
1897  }
1898 
1899  if (!m_switch2) {
1900  nuvec2_tmpj.SetCoordinates(nuvec2_tmpj.Pt(), nuvec2_tmpj.Eta(), nuvec2_tmpj.Phi(), M_nu2);
1901  tauvecsol2j.SetPxPyPzE(0., 0., 0., 0.);
1902  tauvecsol2j += nuvec2_tmpj;
1903  tauvecsol2j += m_tauVec2;
1904  if (tauvecsol2j.E() >= preparedInput.m_beamEnergy)
1905  continue;
1906  tauvecprob2j = Prob->apply(preparedInput, -99, preparedInput.m_type_visTau2,
1907  PtEtaPhiMVector(0, 0, 0, 0), m_tauVec2,
1908  PtEtaPhiMVector(0, 0, 0, 0), nuvec2_tmpj, false, true, false);
1909  ++ngoodsol2;
1910  }
1911  }
1912  if (tauvecprob1j == 0.)
1913  continue;
1914  if (tauvecprob2j == 0.)
1915  continue;
1916 
1917  double totalProb = 1.;
1918 
1919  m_tautau_tmp.SetPxPyPzE(0., 0., 0., 0.);
1920  m_tautau_tmp += tauvecsol1j;
1921  m_tautau_tmp += tauvecsol2j;
1922  const double mtautau = m_tautau_tmp.M();
1923 
1924  if (TailCleanUp(m_tauVec1, nuvec1_tmpj, m_tauVec2, nuvec2_tmpj, mtautau, Mvis, Meff,
1925  preparedInput.m_DelPhiTT) == 0) {
1926  continue;
1927  }
1928 
1929  totalProb *=
1930  (constProb * tauvecprob1j * tauvecprob2j *
1932  m_tauVec1, m_tauVec2, nuvec1_tmpj, nuvec2_tmpj, false, false, true));
1933 
1934  if (totalProb <= 0) {
1936  Warning("DiTauMassTools", "%s",
1937  ("null proba solution, rejected "+std::to_string(totalProb)).c_str());
1938  } else {
1939  // only count solution with non zero probability
1940  m_totalProbSum += totalProb;
1941  m_mtautauSum += mtautau;
1942 
1943  if (m_nsol >= int(m_nsolfinalmax)) {
1944  Error("DiTauMassTools", "%s",
1945  ("refineSolutions ERROR nsol getting larger than nsolfinalmax!!! " +
1947  .c_str());
1948  Error("DiTauMassTools", "%s",
1949  (" j1 " + std::to_string(j1) + " j2 " + std::to_string(j2) + " nsol1 " +
1950  std::to_string(nsol1) + " nsol2 " + std::to_string(nsol2))
1951  .c_str());
1952  --m_nsol; // overwrite last solution. However this should really never
1953  // happen
1954  }
1955 
1956  // good solution found, copy in vector
1957  m_mtautauFinalSolVec[m_nsol] = mtautau;
1958  m_probFinalSolVec[m_nsol] = totalProb;
1959 
1960  PtEtaPhiMVector &nu1Final = m_nu1FinalSolVec[m_nsol];
1961  PtEtaPhiMVector &nu2Final = m_nu2FinalSolVec[m_nsol];
1962  // for (int iv=0;iv<4;++iv){
1963 
1964  nu1Final.SetPxPyPzE(nuvec1_tmpj.Px(), nuvec1_tmpj.Py(), nuvec1_tmpj.Pz(), nuvec1_tmpj.E());
1965  nu2Final.SetPxPyPzE(nuvec2_tmpj.Px(), nuvec2_tmpj.Py(), nuvec2_tmpj.Pz(), nuvec2_tmpj.E());
1966  // }
1967 
1968  ++m_nsol;
1969  } // else totalProb<=0
1970 
1971  } // loop j2
1972  } // loop j1
1973  if (ngoodsol1 == 0) {
1974  return -1;
1975  }
1976  if (ngoodsol2 == 0) {
1977  return -2;
1978  }
1979  return m_nsol;
1980 }
1981 
1982 int MissingMassCalculator::TailCleanUp(const PtEtaPhiMVector &vis1,
1983  const PtEtaPhiMVector &nu1,
1984  const PtEtaPhiMVector &vis2,
1985  const PtEtaPhiMVector &nu2, const double &mmc_mass,
1986  const double &vis_mass, const double &eff_mass,
1987  const double &dphiTT) {
1988 
1989  int pass_code = 1;
1991  return pass_code;
1992 
1993  // the Clean-up cuts are specifically for rel16 analyses.
1994  // the will change in rel17 analyses and after the MMC is updated
1995 
1996  if (preparedInput.m_tauTypes == TauTypes::ll) // lepton-lepton channel
1997  {
1998  const double MrecoMvis = mmc_mass / vis_mass;
1999  if (MrecoMvis > 2.6)
2000  return 0;
2001  const double MrecoMeff = mmc_mass / eff_mass;
2002  if (MrecoMeff > 1.9)
2003  return 0;
2004  const double e1p1 = nu1.E() / vis1.P();
2005  const double e2p2 = nu2.E() / vis2.P();
2006  if ((e1p1 + e2p2) > 4.5)
2007  return 0;
2008  if (e2p2 > 4.0)
2009  return 0;
2010  if (e1p1 > 3.0)
2011  return 0;
2012  }
2013 
2014  //-------- these are new cuts for lep-had analysis for Moriond
2015  if (preparedInput.m_tauTypes == TauTypes::lh) // lepton-hadron channel
2016  {
2017 
2023  return pass_code; // don't use TailCleanup for 8 & 13 TeV data
2024 
2025  //--------- leave code uncommented to avoid Compilation warnings
2026  if (Prob->GetUseHT()) {
2027  const double MrecoMvis = mmc_mass / vis_mass;
2028  const double MrecoMeff = mmc_mass / eff_mass;
2029  const double x = dphiTT > 1.5 ? dphiTT : 1.5;
2030  if ((MrecoMeff + MrecoMvis) > 5.908 - 1.881 * x + 0.2995 * x * x)
2031  return 0;
2032  }
2033  }
2034  return pass_code;
2035 }
2036 
2037 // note that if MarkovChain the input solutions can be modified
2039 
2040 {
2041 
2042  bool reject = true;
2043  double totalProbSumSol = 0.;
2044  double totalProbSumSolOld = 0.;
2045  bool firstPointWithSol = false;
2046 
2047  for (int isol = 0; isol < m_nsol; ++isol) {
2048  totalProbSumSol += m_probFinalSolVec[isol];
2049  }
2050 
2051  double uMC = -1.;
2052  bool notSureToKeep = true;
2053  // note : if no solution, the point is treated as having a zero probability
2054  if (m_fullParamSpaceScan) {
2055  reject = false; // accept anyway in this mode
2056  notSureToKeep = false; // do not need to test on prob
2057  if (m_nsol <= 0) {
2058  // if initial full scaning and no sol : continue
2059  m_markovNFullScan += 1;
2060  } else {
2061  // if we were in in full scan mode and we have a solution, switch it off
2062  m_fullParamSpaceScan = false;
2063  firstPointWithSol = true; // as this is the first point without a solution
2064  // there is no old sol
2065  m_iter0 = 0; // reset the counter so that separately the full scan pphase
2066  // and the markov phase use m_niterRandomLocal points
2067  // hack for hh : allow 10 times less iteration for markov than for the
2068  // fullscan phase
2070  m_niterRandomLocal /= 10;
2071  }
2072  }
2073  }
2074 
2075  if (notSureToKeep) {
2076  // apply Metropolis algorithm to decide to keep this point.
2077  // compute the probability of the previous point and the current one
2078  for (int isol = 0; isol < m_nsolOld; ++isol) {
2079  totalProbSumSolOld += m_probFinalSolOldVec[isol];
2080  }
2081 
2082  // accept anyway if null old probability (should only happen for the very
2083  // first point with a solution)
2084  if (!firstPointWithSol && totalProbSumSolOld <= 0.) {
2085  Error("DiTauMassTools", "%s",
2086  (" ERROR null old probability !!! " + std::to_string(totalProbSumSolOld) + " nsolOld " +
2088  .c_str());
2089  reject = false;
2090  } else if (totalProbSumSol > totalProbSumSolOld) {
2091  // if going up, accept anyway
2092  reject = false;
2093  // else if (totalProbSumSol < 1E-16) { // if null target probability,
2094  // reject anyway
2095  } else if (totalProbSumSol < totalProbSumSolOld * 1E-6) { // if ratio of probability <1e6, point
2096  // will be accepted only every 1E6
2097  // iteration, so can reject anyway
2098  reject = true;
2099  } else if (m_nsol <= 0) { // new parametrisation give prob too small to
2100  // trigger above condition if no solution is found
2101  reject = true;
2102  } else {
2103  // if going down, reject with a probability
2104  // 1-totalProbSum/totalProbSumOld)
2105  uMC = m_randomGen.Rndm();
2106  reject = (uMC > totalProbSumSol / totalProbSumSolOld);
2107  }
2108  } // if reject
2109 
2110  // proceed with the handling of the solutions wether the old or the new ones
2111 
2112  // optionally fill the vectors with the complete list of points (for all
2113  // walkstrategy)
2114 
2115  if (reject) {
2116  // current point reset to the previous one
2117  // Note : only place where m_MEtP etc... are modified outside spacewalkerXYZ
2118  m_MEtP = m_MEtP0;
2119  m_MEtL = m_MEtL0;
2120  m_Phi1 = m_Phi10;
2121  m_Phi2 = m_Phi20;
2122  m_eTau1 = m_eTau10;
2123  m_eTau2 = m_eTau20;
2124  if (m_scanMnu1)
2125  m_Mnu1 = m_Mnu10;
2126  if (m_scanMnu2)
2127  m_Mnu2 = m_Mnu20;
2128  }
2129 
2130  // default case : fill the histogram with solution, using current point
2131  bool fillSolution = true;
2132  bool oldToBeUsed = false;
2133 
2134  // now handle the reject or accept cases
2135  // the tricky thing is that for markov, we accept the old point as soon as a
2136  // new accepted point is found with a weight equal to one plus the number of
2137  // rejected point inbetween
2138 
2139  if (reject) {
2140  fillSolution = false; // do not fill solution, just count number of replication
2142  if (m_nsol <= 0) {
2143  m_markovNRejectNoSol += 1;
2144  } else {
2146  }
2147 
2148  } else {
2149  // if accept, will fill solution (except for very first point) but taking
2150  // the values from the previous point
2151  if (!m_fullParamSpaceScan) {
2152  m_markovNAccept += 1;
2153  }
2154  if (!firstPointWithSol) {
2155  fillSolution = true;
2156  oldToBeUsed = true;
2157  } else {
2158  fillSolution = false;
2159  }
2160  } // else reject
2161 
2162  // if do not fill solution exit now
2163  // for the first point with solution we need to copy the new sol into the old
2164  // one before leaving
2165  if (!fillSolution) {
2166  if (firstPointWithSol) {
2167  // current point is the future previous one
2168  m_nsolOld = m_nsol;
2169  for (int isol = 0; isol < m_nsol; ++isol) {
2172  m_nu1FinalSolOldVec[isol] = m_nu1FinalSolVec[isol];
2173  m_nu2FinalSolOldVec[isol] = m_nu2FinalSolVec[isol];
2174  }
2175  }
2176  return;
2177  }
2178 
2179  // compute RMS of the different solutions
2180  double solSum = 0.;
2181  double solSum2 = 0.;
2182 
2183  for (int isol = 0; isol < m_nsol; ++isol) {
2184  ++m_iter5;
2185  double totalProb;
2186  double mtautau;
2187  const PtEtaPhiMVector *pnuvec1_tmpj;
2188  const PtEtaPhiMVector *pnuvec2_tmpj;
2189 
2190  if (oldToBeUsed) {
2191  totalProb = m_probFinalSolOldVec[isol];
2192  mtautau = m_mtautauFinalSolOldVec[isol];
2193  pnuvec1_tmpj = &m_nu1FinalSolOldVec[isol];
2194  pnuvec2_tmpj = &m_nu2FinalSolOldVec[isol];
2195  } else {
2196  totalProb = m_probFinalSolVec[isol];
2197  mtautau = m_mtautauFinalSolVec[isol];
2198  pnuvec1_tmpj = &m_nu1FinalSolVec[isol];
2199  pnuvec2_tmpj = &m_nu2FinalSolVec[isol];
2200  }
2201  const PtEtaPhiMVector &nuvec1_tmpj = *pnuvec1_tmpj;
2202  const PtEtaPhiMVector &nuvec2_tmpj = *pnuvec2_tmpj;
2203 
2204  solSum += mtautau;
2205  solSum2 += mtautau * mtautau;
2206 
2207  double weight;
2208  // MarkovChain : accepted events already distributed according to
2209  // probability distribution, so weight is 1. acutally to have a proper
2210  // estimate of per bin error, instead of putting several time the same point
2211  // when metropolis alg reject one (or no solution), rather put it with the
2212  // multiplicity weight. Should only change the error bars might change if
2213  // weighted markov chain are used there is also an issue with the 4 very
2214  // close nearly identical solution
2216  1; // incremented only when a point is rejected, hence need to add 1
2217 
2218  m_fMfit_all->Fill(mtautau, weight);
2219 
2220  if(m_SaveLlhHisto){
2221  m_fMEtP_all->Fill(m_MEtP, weight);
2222  m_fMEtL_all->Fill(m_MEtL, weight);
2223  m_fMnu1_all->Fill(m_Mnu1, weight);
2224  m_fMnu2_all->Fill(m_Mnu2, weight);
2225  m_fPhi1_all->Fill(m_Phi1, weight);
2226  m_fPhi2_all->Fill(m_Phi2, weight);
2227  if (mtautau != 0. && weight != 0.)
2228  m_fMfit_allGraph->SetPoint(m_iter0, mtautau, -TMath::Log(weight));
2229  }
2230 
2231  m_fMfit_allNoWeight->Fill(mtautau, 1.);
2232 
2233  // m_fPXfit1->Fill(nuvec1_tmpj.Px(),weight);
2234  // m_fPYfit1->Fill(nuvec1_tmpj.Py(),weight);
2235  // m_fPZfit1->Fill(nuvec1_tmpj.Pz(),weight);
2236  // m_fPXfit2->Fill(nuvec2_tmpj.Px(),weight);
2237  // m_fPYfit2->Fill(nuvec2_tmpj.Py(),weight);
2238  // m_fPZfit2->Fill(nuvec2_tmpj.Pz(),weight);
2239 
2240  //----------------- using P*fit to fill Px,y,z_tau
2241  // Note that the original vistau are used there deliberately,
2242  // since they will be subtracted after histogram fitting
2243  // DR, kudos Antony Lesage : do not create temporary TLV within each Fill,
2244  // saves 10% CPU
2245  m_fPXfit1->Fill(preparedInput.m_vistau1.Px() + nuvec1_tmpj.Px(), totalProb);
2246  m_fPYfit1->Fill(preparedInput.m_vistau1.Py() + nuvec1_tmpj.Py(), totalProb);
2247  m_fPZfit1->Fill(preparedInput.m_vistau1.Pz() + nuvec1_tmpj.Pz(), totalProb);
2248  m_fPXfit2->Fill(preparedInput.m_vistau2.Px() + nuvec2_tmpj.Px(), totalProb);
2249  m_fPYfit2->Fill(preparedInput.m_vistau2.Py() + nuvec2_tmpj.Py(), totalProb);
2250  m_fPZfit2->Fill(preparedInput.m_vistau2.Pz() + nuvec2_tmpj.Pz(), totalProb);
2251 
2252  // fill histograms for floating stopping criterion, split randomly
2253  if (m_fUseFloatStopping) {
2254  if (m_randomGen.Rndm() <= 0.5) {
2255  m_fMmass_split1->Fill(mtautau, weight);
2256  m_fMEtP_split1->Fill(m_MEtP, weight);
2257  m_fMEtL_split1->Fill(m_MEtL, weight);
2258  m_fMnu1_split1->Fill(m_Mnu1, weight);
2259  m_fMnu2_split1->Fill(m_Mnu2, weight);
2260  m_fPhi1_split1->Fill(m_Phi1, weight);
2261  m_fPhi2_split1->Fill(m_Phi2, weight);
2262  } else {
2263  m_fMmass_split2->Fill(mtautau, weight);
2264  m_fMEtP_split2->Fill(m_MEtP, weight);
2265  m_fMEtL_split2->Fill(m_MEtL, weight);
2266  m_fMnu1_split2->Fill(m_Mnu1, weight);
2267  m_fMnu2_split2->Fill(m_Mnu2, weight);
2268  m_fPhi1_split2->Fill(m_Phi1, weight);
2269  m_fPhi2_split2->Fill(m_Phi2, weight);
2270  }
2271  }
2272 
2273  if (totalProb > m_prob_tmp) // fill solution with highest probability
2274  {
2275  m_prob_tmp = totalProb;
2276  m_fDitauStuffFit.Mditau_best = mtautau;
2277  m_fDitauStuffFit.Sign_best = -log10(totalProb);
2278  ;
2279  m_fDitauStuffFit.nutau1 = nuvec1_tmpj;
2280  m_fDitauStuffFit.nutau2 = nuvec2_tmpj;
2283  }
2284  } // loop on solutions
2285 
2286  m_markovCountDuplicate = 0; // now can reset the duplicate count
2287 
2288  if (oldToBeUsed) {
2289  // current point is the future previous one
2290  // TLV copy not super efficient but not dramatic
2291  m_nsolOld = m_nsol;
2292  for (int isol = 0; isol < m_nsol; ++isol) {
2295  m_nu1FinalSolOldVec[isol] = m_nu1FinalSolVec[isol];
2296  m_nu2FinalSolOldVec[isol] = m_nu2FinalSolVec[isol];
2297  }
2298  }
2299 
2300  // compute rms of solutions
2301  const double solRMS = sqrt(solSum2 / m_nsol - std::pow(solSum / m_nsol, 2));
2302  OutputInfo.m_AveSolRMS += solRMS;
2303 
2304  return;
2305 }
2306 
2308  // FIXME could use function pointer to switch between functions
2309  m_nsolOld = 0;
2310 
2311  double METresX = preparedInput.m_METsigmaL; // MET resolution in direction parallel to MET
2312  // resolution major axis, for MET scan
2313  double METresY = preparedInput.m_METsigmaP; // MET resolution in direction perpendicular to
2314  // to MET resolution major axis, for MET scan
2315 
2316  // precompute some quantities and store in m_ data members
2317  precomputeCache();
2324  }
2325 
2326  // if m_nsigma_METscan was not set by user, set to default values
2327  if(m_nsigma_METscan == -1){
2328  if (preparedInput.m_tauTypes == TauTypes::ll) // both tau's are leptonic
2329  {
2331  } else if (preparedInput.m_tauTypes == TauTypes::lh) // lep had
2332  {
2334  } else // hh
2335  {
2337  }
2338  }
2339 
2341 
2344 
2345  m_walkWeight = 1.;
2346 
2347  // dummy initial value to avoid printout with random values
2348  m_Phi10 = 0.;
2349  m_Phi20 = 0.;
2350  m_MEtL0 = 0.;
2351  m_MEtP0 = 0.;
2352  m_Mnu10 = 0.;
2353  m_Mnu20 = 0.;
2354 
2355  m_mTau = 1.777;
2356 
2357  // seeds the random generator in a reproducible way from the phi of both tau;
2358  double aux = std::abs(m_tauVec1Phi + double(m_tauVec2Phi) / 100. / TMath::Pi()) * 100;
2359  m_seed = (aux - floor(aux)) * 1E6 * (1 + m_RndmSeedAltering) + 13;
2360 
2361  m_randomGen.SetSeed(m_seed);
2362  // int Niter=Niter_fit1; // number of points for each dR loop
2363  // int NiterMET=Niter_fit2; // number of iterations for each MET scan loop
2364  // int NiterMnu=Niter_fit3; // number of iterations for Mnu loop
2365 
2366  // approximately compute the number of points from the grid scanning
2367  // divide by abritry number to recover timing with still better results
2368  // m_NiterRandom=(NiterMET+1)*(NiterMET+1)*4*Niter*Niter/10;
2369 
2373 
2377 
2378  m_Mnu1Min = 0.;
2379  m_scanMnu1 = false;
2380  m_Mnu1 = m_Mnu1Min;
2381 
2382  // for markov chain use factor 2
2384 
2385  // NiterRandom set by user (default is -1). If negative, defines the default
2386  // here. no more automatic scaling for ll hl hh
2387  if (m_NiterRandom <= 0) {
2388  m_niterRandomLocal = 100000; // number of iterations for Markov for lh
2390  m_niterRandomLocal *= 2; // multiplied for ll , unchecked
2392  m_niterRandomLocal *= 5; // divided for hh ,checked
2393  } else {
2395  }
2396 
2397  if (preparedInput.m_type_visTau1 == 8) {
2398  // m_Mnu1Max=m_mTau-m_tauVec1M;
2401  m_scanMnu1 = true;
2402  }
2403 
2404  m_Mnu2Min = 0.;
2405  m_scanMnu2 = false;
2406  m_Mnu2 = m_Mnu2Min;
2407  if (preparedInput.m_type_visTau2 == 8) {
2408  // m_Mnu2Max=m_mTau-m_tauVec2M;
2411  m_scanMnu2 = true;
2412  }
2413 
2414  m_MEtLMin = -m_nsigma_METscan * METresX;
2415  m_MEtLMax = +m_nsigma_METscan * METresX;
2417 
2418  m_MEtPMin = -m_nsigma_METscan * METresY;
2419  m_MEtPMax = +m_nsigma_METscan * METresY;
2421 
2422  m_eTau1Min = -1;
2423  m_eTau1Max = -1;
2424  m_eTau2Min = -1;
2425  m_eTau2Max = -1;
2426 
2427  m_switch1 = true;
2428  m_switch2 = true;
2429 
2431  m_rmsStop = m_RMSStop;
2432 
2433  m_iter0 = -1;
2434  m_iterNuPV3 = 0;
2435  m_testptn1 = 0;
2436  m_testptn2 = 0;
2437  m_testdiscri1 = 0;
2438  m_testdiscri2 = 0;
2439  m_nosol1 = 0;
2440  m_nosol2 = 0;
2441  m_iterNsuc = 0;
2442  if (m_meanbinStop > 0) {
2443  m_meanbinToBeEvaluated = true;
2444  } else {
2445  m_meanbinToBeEvaluated = false;
2446  }
2447 
2451  m_markovNAccept = 0;
2452  m_markovNFullScan = 0;
2453  // set full parameter space scannning for the first steps, until a solution is
2454  // found
2455  m_fullParamSpaceScan = true;
2456  // size of step. Needs to be tune. Start with simple heuristic.
2457  if (m_proposalTryMEt < 0) {
2458  m_MEtProposal = m_MEtPRange / 30.;
2459  } else {
2461  }
2462  if (m_ProposalTryPhi < 0) {
2463  m_PhiProposal = 0.04;
2464  } else {
2466  }
2467  // FIXME if m_Mnu1Range !ne m_Mnu2Range same proposal will be done
2468  if (m_scanMnu1) {
2469  if (m_ProposalTryMnu < 0) {
2470  m_MnuProposal = m_Mnu1Range / 10.;
2471  } else {
2473  }
2474  }
2475  if (m_scanMnu2) {
2476  if (m_ProposalTryMnu < 0) {
2477  m_MnuProposal = m_Mnu2Range / 10.;
2478  } else {
2480  }
2481  }
2482 }
2483 
2484 // iterator. walk has internal counters, should only be used in a while loop
2485 // so far only implement grid strategy
2486 // act on MMC data member to be fast
2488  preparedInput.m_MEtX = -999.;
2489  preparedInput.m_MEtY = -999.;
2490 
2491  ++m_iter0;
2492 
2493  if (m_meanbinToBeEvaluated && m_iterNsuc == 500) {
2494  Info("DiTauMassTools", " in m_meanbinToBeEvaluated && m_iterNsuc==500 ");
2495  // for markov chain m_iterNsuc is the number of *accepted* points, so there
2496  // can be several iterations without any increment of m_iterNsuc. Hence need
2497  // to make sure meanbin is evaluated only once
2498  m_meanbinToBeEvaluated = false;
2499 
2500  // Meanbin stopping criterion
2501  std::vector<double> histInfo(HistInfo::MAXHISTINFO);
2502  // SLIDINGWINDOW strategy to avoid doing the parabola fit now given it will
2503  // not be use
2505  double meanbin = histInfo.at(HistInfo::MEANBIN);
2506  if (meanbin < 0) {
2507  m_nsucStop = -1; // no meaningful meanbin switch back to niter criterion
2508  } else {
2509  double stopdouble = 500 * std::pow((meanbin / m_meanbinStop), 2);
2510  int stopint = stopdouble;
2511  m_nsucStop = stopint;
2512  }
2513  if (m_nsucStop < 500)
2514  return false;
2515  }
2516  // should be outside m_meanbinStop test
2517  if (m_iterNsuc == m_nsucStop)
2518  return false; // Critere d'arret pour nombre de succes
2519 
2520  if (m_iter0 == m_niterRandomLocal)
2521  return false; // for now simple stopping criterion on number of iteration
2522 
2523  // floating stopping criterion, reduces run-time for lh, hh by a factor ~2 and ll by roughly
2524  // factor ~3 check if every scanned variable and resulting mass thermalised after 10k iterations
2525  // and then every 1k iterations do this by checking that the means of the split distributions is
2526  // comparable within 5% of their sigma
2527  if (m_iter0 >= 10000 && (m_iter0 % 1000) == 0 && m_fUseFloatStopping) {
2528  if (std::abs(m_fMEtP_split1->GetMean() - m_fMEtP_split2->GetMean()) <= 0.05 * m_fMEtP_split1->GetRMS()) {
2529  if (std::abs(m_fMEtL_split1->GetMean() - m_fMEtL_split2->GetMean()) <=
2530  0.05 * m_fMEtL_split1->GetRMS()) {
2531  if (std::abs(m_fMnu1_split1->GetMean() - m_fMnu1_split2->GetMean()) <=
2532  0.05 * m_fMnu1_split1->GetRMS()) {
2533  if (std::abs(m_fMnu2_split1->GetMean() - m_fMnu2_split2->GetMean()) <=
2534  0.05 * m_fMnu2_split1->GetRMS()) {
2535  if (std::abs(m_fPhi1_split1->GetMean() - m_fPhi1_split2->GetMean()) <=
2536  0.05 * m_fPhi1_split1->GetRMS()) {
2537  if (std::abs(m_fPhi2_split1->GetMean() - m_fPhi2_split2->GetMean()) <=
2538  0.05 * m_fPhi2_split1->GetRMS()) {
2539  if (std::abs(m_fMmass_split1->GetMean() - m_fMmass_split2->GetMean()) <=
2540  0.05 * m_fMmass_split1->GetRMS()) {
2541  return false;
2542  }
2543  }
2544  }
2545  }
2546  }
2547  }
2548  }
2549  }
2550 
2551  if (m_fullParamSpaceScan) {
2552  // as long as no solution found need to randomise on the full parameter
2553  // space
2554 
2555  // cut the corners in MissingET (not optimised at all)
2556  // not needed if distribution is already gaussian
2557  do {
2558  m_MEtP = m_MEtPMin + m_MEtPRange * m_randomGen.Rndm();
2559  m_MEtL = m_MEtLMin + m_MEtLRange * m_randomGen.Rndm();
2560  } while (!checkMEtInRange());
2561 
2562  if (m_scanMnu1) {
2563  m_Mnu1 = m_Mnu1Min + m_Mnu1Range * m_randomGen.Rndm();
2564  }
2565 
2566  if (m_scanMnu2) {
2567  m_Mnu2 = m_Mnu2Min + m_Mnu2Range * m_randomGen.Rndm();
2568  }
2569 
2570  m_Phi1 = m_Phi1Min + m_Phi1Range * m_randomGen.Rndm();
2571  m_Phi2 = m_Phi2Min + m_Phi2Range * m_randomGen.Rndm();
2572 
2573  return true;
2574  }
2575 
2576  // here the real markov chain takes place : "propose" the new point
2577  // note that if one parameter goes outside range, this should not be fixed
2578  // here but later in handleSolution, otherwise would cause a bias
2579 
2580  // m_MEtP0 etc... also store the position of the previous Markov Chain step,
2581  // which is needed by the algorithm
2582  m_MEtP0 = m_MEtP;
2583  m_MEtL0 = m_MEtL;
2584 
2586 
2588 
2589  if (m_scanMnu1) {
2590  m_Mnu10 = m_Mnu1;
2592  }
2593 
2594  if (m_scanMnu2) {
2595  m_Mnu20 = m_Mnu2;
2597  }
2598 
2599  m_Phi10 = m_Phi1;
2601 
2602  m_Phi20 = m_Phi2;
2603 
2605 
2606  return true;
2607 }
2608 
2609 // compute cached values (this value do not change within one call of MMC,
2610 // except for tau e scanning) return true if cache was already uptodatexs
2612 
2613  // copy tau 4 vect. If tau E scanning, these vectors will be modified
2616 
2617  const XYVector &metVec = preparedInput.m_MetVec;
2618 
2619  bool same = true;
2634 
2635  same = updateDouble(1.777, m_mTau) && same;
2640 
2641  PtEtaPhiMVector Met4vec;
2642  Met4vec.SetPxPyPzE(preparedInput.m_MetVec.X(), preparedInput.m_MetVec.Y(), 0.0,
2643  preparedInput.m_MetVec.R());
2644  same = updateDouble((m_tauVec1 + m_tauVec2 + Met4vec).M(), m_Meff) && same;
2645 
2647  // note that if useHT met_vec is actually -HT
2648  same = updateDouble(metVec.X(), preparedInput.m_inputMEtX) && same;
2649  same = updateDouble(metVec.Y(), preparedInput.m_inputMEtY) && same;
2650  same = updateDouble(metVec.R(), preparedInput.m_inputMEtT) && same;
2651 
2652  return same;
2653 }
2654 
2655 // return true if all parameters are within their domain
2657 
2658  if (m_scanMnu1) {
2659  if (m_Mnu1 < m_Mnu1Min)
2660  return false;
2661  if (m_Mnu1 > m_Mnu1Max)
2662  return false;
2663  if (m_Mnu1 > m_mTau - m_tauVec1M)
2664  return false;
2665  }
2666 
2667  if (m_scanMnu2) {
2668  if (m_Mnu2 < m_Mnu2Min)
2669  return false;
2670  if (m_Mnu2 > m_Mnu2Max)
2671  return false;
2672  if (m_Mnu2 > m_mTau - m_tauVec2M)
2673  return false;
2674  }
2675 
2676  // FIXME note that since there is a coupling between Met and tau, should
2677  // rigorously test both together however since the 3 sigma range is just a
2678  // hack, it is probably OK
2679 
2680  if (m_Phi1 < m_Phi1Min)
2681  return false;
2682  if (m_Phi1 > m_Phi1Max)
2683  return false;
2684 
2685  if (m_Phi2 < m_Phi2Min)
2686  return false;
2687  if (m_Phi2 > m_Phi2Max)
2688  return false;
2689 
2690  if (!checkMEtInRange())
2691  return false;
2692 
2693  return true;
2694 }
2695 
2696 // return true if Met is within disk instead of withing square (cut the corners)
2698  // check MEt is in allowed range
2699  // range is 3sigma disk ("cutting the corners")
2703  return false;
2704  } else {
2705  return true;
2706  }
2707 }
2708 
2709 // ----- returns dTheta3D lower and upper boundaries:
2710 // limit_code=0: 99% lower limit
2711 // limit_code=1; 99% upper limit
2712 // limit_code=2; 95% upper limit
2713 double MissingMassCalculator::dTheta3DLimit(const int &tau_type, const int &limit_code,
2714  const double &P_tau) {
2715 
2716 #ifndef WITHDTHETA3DLIM
2717  // make the test ineffective if desired
2718  if (limit_code == 0)
2719  return 0.;
2720  if (limit_code == 1)
2721  return 10.;
2722  if (limit_code == 2)
2723  return 10.;
2724 #endif
2725 
2726  double limit = 1.0;
2727  if (limit_code == 0)
2728  limit = 0.0;
2729  double par[3] = {0.0, 0.0, 0.0};
2730  // ---- leptonic tau's
2731  if (tau_type == 8) {
2732  if (limit_code == 0) // lower 99% limit
2733  {
2734  par[0] = 0.3342;
2735  par[1] = -0.3376;
2736  par[2] = -0.001377;
2737  }
2738  if (limit_code == 1) // upper 99% limit
2739  {
2740  par[0] = 3.243;
2741  par[1] = -12.87;
2742  par[2] = 0.009656;
2743  }
2744  if (limit_code == 2) // upper 95% limit
2745  {
2746  par[0] = 2.927;
2747  par[1] = -7.911;
2748  par[2] = 0.007783;
2749  }
2750  }
2751  // ---- 1-prong tau's
2752  if (tau_type >= 0 && tau_type <= 2) {
2753  if (limit_code == 0) // lower 99% limit
2754  {
2755  par[0] = 0.2673;
2756  par[1] = -14.8;
2757  par[2] = -0.0004859;
2758  }
2759  if (limit_code == 1) // upper 99% limit
2760  {
2761  par[0] = 9.341;
2762  par[1] = -15.88;
2763  par[2] = 0.0333;
2764  }
2765  if (limit_code == 2) // upper 95% limit
2766  {
2767  par[0] = 6.535;
2768  par[1] = -8.649;
2769  par[2] = 0.00277;
2770  }
2771  }
2772  // ---- 3-prong tau's
2773  if (tau_type >= 3 && tau_type <= 5) {
2774  if (limit_code == 0) // lower 99% limit
2775  {
2776  par[0] = 0.2308;
2777  par[1] = -15.24;
2778  par[2] = -0.0009458;
2779  }
2780  if (limit_code == 1) // upper 99% limit
2781  {
2782  par[0] = 14.58;
2783  par[1] = -6.043;
2784  par[2] = -0.00928;
2785  }
2786  if (limit_code == 2) // upper 95% limit
2787  {
2788  par[0] = 8.233;
2789  par[1] = -0.3018;
2790  par[2] = -0.009399;
2791  }
2792  }
2793 
2794  if (std::abs(P_tau + par[1]) > 0.0)
2795  limit = par[0] / (P_tau + par[1]) + par[2];
2796  if (limit_code == 0) {
2797  if (limit < 0.0) {
2798  limit = 0.0;
2799  } else if (limit > 0.03) {
2800  limit = 0.03;
2801  }
2802  } else {
2803  if (limit < 0.0 || limit > 0.5 * TMath::Pi()) {
2804  limit = 0.5 * TMath::Pi();
2805  } else if (limit < 0.05 && limit > 0.0) {
2806  limit = 0.05; // parameterization only runs up to P~220 GeV in this regime
2807  // will set an upper bound of 0.05
2808  }
2809  }
2810 
2811  return limit;
2812 }
2813 
2814 // checks units of input variables, converts into [GeV] if needed, make all
2815 // possible corrections DR new : now a second structure preparedInput is derived
2816 // from the input one which only has direct user input
2818  const xAOD::IParticle *part2,
2819  const xAOD::MissingET *met,
2820  const int &njets) {
2821  const double GEV = 1000.;
2822  int mmcType1 = mmcType(part1);
2823  if (mmcType1 < 0)
2824  return; // return CP::CorrectionCode::Error;
2825 
2826  int mmcType2 = mmcType(part2);
2827  if (mmcType2 < 0)
2828  return; // return CP::CorrectionCode::Error;
2829 
2830  preparedInput.SetLFVmode(-2); // initialise LFV mode value for this event with being *not* LFV
2831  // if(getLFVMode(part1, part2, mmcType1, mmcType2) ==
2832  // CP::CorrectionCode::Error) {
2834  int LFVMode = getLFVMode(part1, part2, mmcType1, mmcType2);
2835  if (LFVMode == -1) {
2836  return; // return CP::CorrectionCode::Error;
2837  } else if (LFVMode != -2) {
2838  preparedInput.SetLFVmode(LFVMode);
2839  }
2840  }
2841 
2842  // this will be in MeV but MMC allows MeV
2843  // assume the mass is correct as well
2844  PtEtaPhiMVector tlvTau1(part1->pt(), part1->eta(), part1->phi(), part1->m());
2845  PtEtaPhiMVector tlvTau2(part2->pt(), part2->eta(), part2->phi(), part2->m());
2846 
2847  // Convert to GeV. In principle, MMC should cope with MeV but should check
2848  // thoroughly
2849  PtEtaPhiMVector fixedtau1;
2850  fixedtau1.SetCoordinates(tlvTau1.Pt() / GEV, tlvTau1.Eta(), tlvTau1.Phi(), tlvTau1.M() / GEV);
2851  PtEtaPhiMVector fixedtau2;
2852  fixedtau2.SetCoordinates(tlvTau2.Pt() / GEV, tlvTau2.Eta(), tlvTau2.Phi(), tlvTau2.M() / GEV);
2853 
2854  preparedInput.SetVisTauType(0, mmcType1);
2855  preparedInput.SetVisTauType(1, mmcType2);
2856  preparedInput.SetVisTauVec(0, fixedtau1);
2857  preparedInput.SetVisTauVec(1, fixedtau2);
2858 
2859  if (mmcType1 == 8 && mmcType2 == 8) {
2861  } else if (mmcType1 >= 0 && mmcType1 <= 5 && mmcType2 >= 0 && mmcType2 <= 5) {
2863  } else {
2865  }
2867  Info("DiTauMassTools", "%s", ("running for tau types "+std::to_string(preparedInput.m_type_visTau1)+" "+std::to_string(preparedInput.m_type_visTau2)).c_str());
2868  XYVector met_vec(met->mpx() / GEV, met->mpy() / GEV);
2869  preparedInput.SetMetVec(met_vec);
2871  Info("DiTauMassTools", "%s", ("passing SumEt="+std::to_string(met->sumet() / GEV)).c_str());
2872  preparedInput.SetSumEt(met->sumet() / GEV);
2873  preparedInput.SetNjet25(njets);
2874 
2875  // check that the calibration set has been chosen explicitly, otherwise abort
2877  Error("DiTauMassTools", "MMCCalibrationSet has not been set !. Please use "
2878  "fMMC.SetCalibrationSet(MMCCalibrationSet::MMC2019) or fMMC.SetCalibrationSet(MMCCalibrationSet::MMC2024)"
2879  ". Abort now. ");
2880  std::abort();
2881  }
2882  //----------- Re-ordering input info, to make sure there is no dependence of
2883  // results on input order
2884  // this might be needed because a random scan is used
2885  // highest pT tau is always first
2886  preparedInput.m_InputReorder = 0; // set flag to 0 by default, i.e. no re-ordering
2888  preparedInput.m_type_visTau2 == 8) // if hadron-lepton, reorder to have lepton first
2889  {
2891  1; // re-order to be done, this flag is to be checked in DoOutputInfo()
2892  } else if (!((preparedInput.m_type_visTau2 >= 0 && preparedInput.m_type_visTau2 <= 5) &&
2893  preparedInput.m_type_visTau1 == 8)) // if not lep-had nor had lep, reorder if tau1 is
2894  // after tau2 clockwise
2895  {
2896  if (fixPhiRange(preparedInput.m_vistau1.Phi() - preparedInput.m_vistau2.Phi()) > 0) {
2897  preparedInput.m_InputReorder = 1; // re-order to be done, this flag is to be
2898  // checked in DoOutputInfo()
2899  }
2900  }
2901 
2902  if (preparedInput.m_InputReorder == 1) // copy and re-order
2903  {
2907  }
2908  //--------- re-ordering is done ---------------------------------------
2909 
2911  std::abs(Phi_mpi_pi(preparedInput.m_vistau1.Phi() - preparedInput.m_vistau2.Phi()));
2912 
2913  for (unsigned int i = 0; i < preparedInput.m_jet4vecs.size(); i++) {
2914  // correcting sumEt, give priority to SetMetScanParamsUE()
2915  if (preparedInput.m_METScanScheme == 0) {
2916  if ((preparedInput.m_METsigmaP < 0.1 || preparedInput.m_METsigmaL < 0.1) &&
2918  preparedInput.m_jet4vecs[i].Pt() > 20.0) {
2919  if (preparedInput.m_fUseVerbose == 1) {
2920  Info("DiTauMassTools", "correcting sumET");
2921  }
2923  }
2924  }
2925  }
2926 
2927  // give priority to SetVisTauType, only do this if type_visTau1 and
2928  // type_visTau2 are not set
2929  /*if(type_visTau1<0 && type_visTau2<0 && Nprong_tau1>-1 && Nprong_tau2>-1)
2930  {
2931  if(Nprong_tau1==0) type_visTau1 = 8; // leptonic tau
2932  else if( Nprong_tau1==1) type_visTau1 = 0; // set to 1p0n for now, may use
2933 different solution later like explicit integer for this case that pantau info is
2934 not available? else if( Nprong_tau1==3) type_visTau1 = 3; // set to 3p0n for
2935 now, see above if(Nprong_tau2==0) type_visTau2 = 8; // leptonic tau else if(
2936 Nprong_tau2==1) type_visTau2 = 0; // set to 1p0n for now, see above else if(
2937 Nprong_tau2==3) type_visTau2=3; // set to 3p0n for now, see above
2938  }
2939  */
2940  // checking input mass of hadronic tau-1
2941  // one prong
2942  // // checking input mass of hadronic tau-1
2943  // DRMERGE LFV addition
2946  preparedInput.m_vistau1.M() != 1.1) {
2948  preparedInput.m_vistau1.Phi(), 1.1);
2949  }
2951  preparedInput.m_vistau1.M() != 1.35) {
2953  preparedInput.m_vistau1.Phi(), 1.35);
2954  }
2955  // checking input mass of hadronic tau-2
2957  preparedInput.m_vistau2.M() != 1.1) {
2959  preparedInput.m_vistau2.Phi(), 1.1);
2960  }
2962  preparedInput.m_vistau2.M() != 1.35) {
2964  preparedInput.m_vistau2.Phi(), 1.35);
2965  }
2966  } else {
2967  // DRMERGE end LFV addition
2969  preparedInput.m_vistau1.M() != 0.8) {
2971  preparedInput.m_vistau1.Phi(), 0.8);
2972  }
2973  // 3 prong
2975  preparedInput.m_vistau1.M() != 1.2) {
2977  preparedInput.m_vistau1.Phi(), 1.2);
2978  }
2979  // checking input mass of hadronic tau-2
2980  // one prong
2982  preparedInput.m_vistau2.M() != 0.8) {
2984  preparedInput.m_vistau2.Phi(), 0.8);
2985  }
2986  // 3 prong
2988  preparedInput.m_vistau2.M() != 1.2) {
2990  preparedInput.m_vistau2.Phi(), 1.2);
2991  }
2992  } // DRDRMERGE LFV else closing
2993 
2994  // correcting sumEt for electron pt, give priority to SetMetScanParamsUE()
2995  // DR20150615 in tag 00-00-11 and before. The following was done before the
2996  // mass of the hadronic tau was set which mean that sumEt was wrongly
2997  // corrected for the hadronic tau pt if the hadronic tau mass was set to zero
2998  // Sasha 08/12/15: don't do electron Pt subtraction for high mass studies; in
2999  // the future, need to check if lepton Pt needs to be subtracted for both ele
3000  // and muon
3001  if (preparedInput.m_METsigmaP < 0.1 || preparedInput.m_METsigmaL < 0.1) {
3002 
3003  // T. Davidek: hack for lep-lep -- subtract lepton pT both for muon and
3004  // electron
3008  preparedInput.m_vistau1.M() < 0.12 && preparedInput.m_vistau2.M() < 0.12) { // lep-lep channel
3013  } else {
3014  // continue with the original code
3017  if (preparedInput.m_fUseVerbose == 1) {
3018  Info("DiTauMassTools", "Substracting pt1 from sumEt");
3019  }
3021  }
3024  if (preparedInput.m_fUseVerbose == 1) {
3025  Info("DiTauMassTools", "Substracting pt2 from sumEt");
3026  }
3028  }
3029  }
3030  }
3031 
3032  // controling TauProbability settings for UPGRADE studies
3035  if ((preparedInput.m_vistau1.M() < 0.12 && preparedInput.m_vistau2.M() > 0.12) ||
3036  (preparedInput.m_vistau2.M() < 0.12 && preparedInput.m_vistau1.M() > 0.12)) {
3037  Prob->SetUseTauProbability(true); // lep-had case
3038  }
3039  if (preparedInput.m_vistau1.M() > 0.12 && preparedInput.m_vistau2.M() > 0.12) {
3040  Prob->SetUseTauProbability(false); // had-had case
3041  }
3042  }
3043 
3044  // change Beam Energy for different running conditions
3046 
3047  //--------------------- pre-set defaults for Run-2. To disable pre-set
3048  // defaults set fUseDefaults=0
3049  if (preparedInput.m_fUseDefaults == 1) {
3051  SetNsigmaMETscan_ll(4.0);
3052  SetNsigmaMETscan_lh(4.0);
3053  SetNsigmaMETscan_hh(4.0);
3055  if ((preparedInput.m_vistau1.M() < 0.12 && preparedInput.m_vistau2.M() > 0.12) ||
3056  (preparedInput.m_vistau2.M() < 0.12 && preparedInput.m_vistau1.M() > 0.12))
3057  Prob->SetUseTauProbability(false); // lep-had
3059  Prob->SetUseTauProbability(true); // had-had
3060  Prob->SetUseMnuProbability(false);
3061  }
3062  }
3063 
3064  // compute HTOffset if relevant
3065  if (Prob->GetUseHT()) // use missing Ht for Njet25=0 events
3066  {
3067  // dPhi(l-t) dependence of misHt-trueMET
3068  double HtOffset = 0.;
3069  // proper for hh
3071  // hh
3072  double x = preparedInput.m_DelPhiTT;
3073  HtOffset = 87.5 - 27.0 * x;
3074  }
3075 
3076  preparedInput.m_HtOffset = HtOffset;
3077 
3078  // if use HT, replace MET with HT
3080  preparedInput.m_MHtSigma2; // sigma of 2nd Gaussian for missing Ht resolution
3082 
3083  PtEtaPhiMVector tauSum = preparedInput.m_vistau1 + preparedInput.m_vistau2;
3084  preparedInput.m_MetVec.SetXY(-tauSum.Px(), -tauSum.Py()); // WARNING this replace metvec by -mht
3085  }
3086 }
3087 
3090  if(!m_SaveLlhHisto) return;
3091 
3092  float hEmax = 3000.0; // maximum energy (GeV)
3093  int hNbins = 1500;
3094  m_fMEtP_all = std::make_shared<TH1F>("MEtP_h1", "M", hNbins, -100.0,
3095  100.); // all solutions
3096  m_fMEtL_all = std::make_shared<TH1F>("MEtL_h1", "M", hNbins, -100.0,
3097  100.); // all solutions
3098  m_fMnu1_all = std::make_shared<TH1F>("Mnu1_h1", "M", hNbins, 0.0,
3099  hEmax); // all solutions
3100  m_fMnu2_all = std::make_shared<TH1F>("Mnu2_h1", "M", hNbins, 0.0,
3101  hEmax); // all solutions
3102  m_fPhi1_all = std::make_shared<TH1F>("Phi1_h1", "M", hNbins, -10.0,
3103  10.); // all solutions
3104  m_fPhi2_all = std::make_shared<TH1F>("Phi2_h1", "M", hNbins, -10.0,
3105  10.); // all solutions
3106  m_fMfit_allGraph = std::make_shared<TGraph>(); // all solutions
3107 
3108  m_fMEtP_all->Sumw2();
3109  m_fMEtL_all->Sumw2();
3110  m_fMnu1_all->Sumw2();
3111  m_fMnu2_all->Sumw2();
3112  m_fPhi1_all->Sumw2();
3113  m_fPhi2_all->Sumw2();
3114 
3115  m_fMEtP_all->SetDirectory(0);
3116  m_fMEtL_all->SetDirectory(0);
3117  m_fMnu1_all->SetDirectory(0);
3118  m_fMnu2_all->SetDirectory(0);
3119  m_fPhi1_all->SetDirectory(0);
3120  m_fPhi2_all->SetDirectory(0);
3121 }
3122 
3125  if(!m_fUseFloatStopping) return;
3126 
3127  float hEmax = 3000.0; // maximum energy (GeV)
3128  int hNbins = 1500;
3129  m_fMmass_split1 = std::make_shared<TH1F>("mass_h1_1", "M", hNbins, 0.0, hEmax);
3130  m_fMEtP_split1 = std::make_shared<TH1F>("MEtP_h1_1", "M", hNbins, -100.0, 100.0);
3131  m_fMEtL_split1 = std::make_shared<TH1F>("MEtL_h1_1", "M", hNbins, -100.0, 100.0);
3132  m_fMnu1_split1 = std::make_shared<TH1F>("Mnu1_h1_1", "M", hNbins, 0.0, hEmax);
3133  m_fMnu2_split1 = std::make_shared<TH1F>("Mnu2_h1_1", "M", hNbins, 0.0, hEmax);
3134  m_fPhi1_split1 = std::make_shared<TH1F>("Phi1_h1_1", "M", hNbins, -10.0, 10.0);
3135  m_fPhi2_split1 = std::make_shared<TH1F>("Phi2_h1_1", "M", hNbins, -10.0, 10.0);
3136  m_fMmass_split2 = std::make_shared<TH1F>("mass_h1_2", "M", hNbins, 0.0, hEmax);
3137  m_fMEtP_split2 = std::make_shared<TH1F>("MEtP_h1_2", "M", hNbins, -100.0, 100.0);
3138  m_fMEtL_split2 = std::make_shared<TH1F>("MEtL_h1_2", "M", hNbins, -100.0, 100.0);
3139  m_fMnu1_split2 = std::make_shared<TH1F>("Mnu1_h1_2", "M", hNbins, 0.0, hEmax);
3140  m_fMnu2_split2 = std::make_shared<TH1F>("Mnu2_h1_2", "M", hNbins, 0.0, hEmax);
3141  m_fPhi1_split2 = std::make_shared<TH1F>("Phi1_h1_2", "M", hNbins, -10.0, 10.0);
3142  m_fPhi2_split2 = std::make_shared<TH1F>("Phi2_h1_2", "M", hNbins, -10.0, 10.0);
3143 
3144  m_fMmass_split1->Sumw2();
3145  m_fMEtP_split1->Sumw2();
3146  m_fMEtL_split1->Sumw2();
3147  m_fMnu1_split1->Sumw2();
3148  m_fMnu2_split1->Sumw2();
3149  m_fPhi1_split1->Sumw2();
3150  m_fPhi2_split1->Sumw2();
3151  m_fMmass_split2->Sumw2();
3152  m_fMEtP_split2->Sumw2();
3153  m_fMEtL_split2->Sumw2();
3154  m_fMnu1_split2->Sumw2();
3155  m_fMnu2_split2->Sumw2();
3156  m_fPhi1_split2->Sumw2();
3157  m_fPhi2_split2->Sumw2();
3158 
3159  m_fMmass_split1->SetDirectory(0);
3160  m_fMEtP_split1->SetDirectory(0);
3161  m_fMEtL_split1->SetDirectory(0);
3162  m_fMnu1_split1->SetDirectory(0);
3163  m_fMnu2_split1->SetDirectory(0);
3164  m_fPhi1_split1->SetDirectory(0);
3165  m_fPhi2_split1->SetDirectory(0);
3166  m_fMmass_split2->SetDirectory(0);
3167  m_fMEtP_split2->SetDirectory(0);
3168  m_fMEtL_split2->SetDirectory(0);
3169  m_fMnu1_split2->SetDirectory(0);
3170  m_fMnu2_split2->SetDirectory(0);
3171  m_fPhi1_split2->SetDirectory(0);
3172  m_fPhi2_split2->SetDirectory(0);
3173 }
DiTauMassTools::MMCCalibrationSet::MMC2024
@ MMC2024
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:40
xAOD::iterator
JetConstituentVector::iterator iterator
Definition: JetConstituentVector.cxx:68
DiTauMassTools::MissingMassCalculator::DitauStuff::Mditau_best
double Mditau_best
Definition: MissingMassCalculator.h:54
DiTauMassTools::MissingMassCalculator::m_Phi1
double m_Phi1
Definition: MissingMassCalculator.h:141
DiTauMassTools::MissingMassProb::SetUseDphiLL
void SetUseDphiLL(bool val)
Definition: MissingMassProb.h:57
AllowedVariables::e
e
Definition: AsgElectronSelectorTool.cxx:37
DiTauMassTools::MissingMassCalculator::m_fFitting
TF1 * m_fFitting
Definition: MissingMassCalculator.h:231
DiTauMassTools::MissingMassCalculator::NuPsolutionLFV
int NuPsolutionLFV(const XYVector &met_vec, const PtEtaPhiMVector &tau, const double &m_nu, std::vector< PtEtaPhiMVector > &nu_vec)
Definition: MissingMassCalculator.cxx:745
DiTauMassTools::MissingMassCalculator::m_fMnu1_split2
std::shared_ptr< TH1F > m_fMnu1_split2
Definition: MissingMassCalculator.h:226
DiTauMassTools::MissingMassCalculator::m_testdiscri2
int m_testdiscri2
Definition: MissingMassCalculator.h:115
DiTauMassTools::MissingMassCalculator::m_eTau2
double m_eTau2
Definition: MissingMassCalculator.h:142
DiTauMassTools::MissingMassProb::setParamAngle
void setParamAngle(const PtEtaPhiMVector &tauvec, int tau, int tautype)
Definition: MissingMassProb.cxx:155
DiTauMassTools::TauTypes::lh
@ lh
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:53
DiTauMassTools::MissingMassCalculator::m_tautau_tmp
PtEtaPhiMVector m_tautau_tmp
Definition: MissingMassCalculator.h:88
DiTauMassTools::MissingMassCalculator::m_switch2
bool m_switch2
Definition: MissingMassCalculator.h:120
DiTauMassTools::MissingMassCalculator::SaveLlhHisto
void SaveLlhHisto(const bool val)
Definition: MissingMassCalculator.cxx:3088
DiTauMassTools::MissingMassCalculator::m_nuvecsol1
std::vector< PtEtaPhiMVector > m_nuvecsol1
Definition: MissingMassCalculator.h:77
DiTauMassTools::MaxHistStrategy::e
e
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:30
DiTauMassTools::MissingMassCalculator::SpaceWalkerInit
void SpaceWalkerInit()
Definition: MissingMassCalculator.cxx:2307
DiTauMassTools::MissingMassCalculator::m_MEtPMax
double m_MEtPMax
Definition: MissingMassCalculator.h:146
DiTauMassTools::MissingMassCalculator::m_fMnu2_split2
std::shared_ptr< TH1F > m_fMnu2_split2
Definition: MissingMassCalculator.h:227
athena.path
path
python interpreter configuration --------------------------------------—
Definition: athena.py:128
DiTauMassTools::MissingMassCalculator::precomputeCache
bool precomputeCache()
Definition: MissingMassCalculator.cxx:2611
DiTauMassTools::MissingMassCalculator::m_MEtL0
double m_MEtL0
Definition: MissingMassCalculator.h:144
DiTauMassTools::MissingMassCalculator::NuPsolutionV3
int NuPsolutionV3(const double &mNu1, const double &mNu2, const double &phi1, const double &phi2, int &nsol1, int &nsol2)
Definition: MissingMassCalculator.cxx:540
DiTauMassTools::MissingMassCalculator::DitauStuff::nutau1
PtEtaPhiMVector nutau1
Definition: MissingMassCalculator.h:56
DiTauMassTools::MMCFitMethod::MAX
@ MAX
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:46
DiTauMassTools::MMCCalibrationSet::MMC2015HIGHMASS
@ MMC2015HIGHMASS
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:40
DiTauMassTools::MissingMassCalculator::m_ProposalTryPhi
double m_ProposalTryPhi
Definition: MissingMassCalculator.h:135
DiTauMassTools::MissingMassProb::GetUseTauProbability
bool GetUseTauProbability()
Definition: MissingMassProb.h:52
DiTauMassTools::MissingMassCalculator::m_nosol2
int m_nosol2
Definition: MissingMassCalculator.h:117
DiTauMassTools::TauTypes::hh
@ hh
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:53
DiTauMassTools::MissingMassCalculator::m_E2v1
double m_E2v1
Definition: MissingMassCalculator.h:189
DiTauMassTools::MissingMassCalculator::m_fMnu2_all
std::shared_ptr< TH1F > m_fMnu2_all
Definition: MissingMassCalculator.h:202
DiTauMassTools::MissingMassCalculator::m_fDitauStuffHisto
DitauStuff m_fDitauStuffHisto
Definition: MissingMassCalculator.h:247
DiTauMassTools::MissingMassCalculator::m_iter5
int m_iter5
Definition: MissingMassCalculator.h:98
DiTauMassTools::MissingMassCalculator::m_tauvecprob1
std::vector< double > m_tauvecprob1
Definition: MissingMassCalculator.h:82
DiTauMassTools::MissingMassCalculator::m_PrintmMaxError
double m_PrintmMaxError
Definition: MissingMassCalculator.h:130
DiTauMassTools::MissingMassCalculator::m_m2Nu1
double m_m2Nu1
Definition: MissingMassCalculator.h:185
DiTauMassTools::MMCCalibrationSet::name
const std::string name[MAXMMCCALIBRATIONSET]
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:41
DiTauMassTools::MissingMassInput::m_METsigmaP
double m_METsigmaP
Definition: MissingMassInput.h:62
DiTauMassTools::MissingMassCalculator::m_Ev2
double m_Ev2
Definition: MissingMassCalculator.h:191
DiTauMassTools::MissingMassCalculator::m_Phi2Max
double m_Phi2Max
Definition: MissingMassCalculator.h:146
DiTauMassTools::MissingMassCalculator::m_tauVec2Phi
double m_tauVec2Phi
Definition: MissingMassCalculator.h:178
DiTauMassTools::MissingMassCalculator::m_fMEtL_split1
std::shared_ptr< TH1F > m_fMEtL_split1
Definition: MissingMassCalculator.h:218
DiTauMassTools::MissingMassInput::m_METScanScheme
int m_METScanScheme
Definition: MissingMassInput.h:74
DiTauMassTools::MissingMassCalculator::m_fPhi2_split1
std::shared_ptr< TH1F > m_fPhi2_split1
Definition: MissingMassCalculator.h:222
DiTauMassTools::MissingMassCalculator::m_lfvLeplepRefit
bool m_lfvLeplepRefit
Definition: MissingMassCalculator.h:90
DiTauMassTools::MissingMassCalculator::refineSolutions
int refineSolutions(const double &M_nu1, const double &M_nu2, const int nsol1, const int nsol2, const double &Mvis, const double &Meff)
Definition: MissingMassCalculator.cxx:1811
max
constexpr double max()
Definition: ap_fixedTest.cxx:33
DiTauMassTools::MMCFitMethod::MLM
@ MLM
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:46
DMTest::P
P_v1 P
Definition: P.h:23
DiTauMassTools::MissingMassCalculator::m_metCovPhiCos
double m_metCovPhiCos
Definition: MissingMassCalculator.h:150
DiTauMassTools::MissingMassCalculator::m_eventNumber
int m_eventNumber
Definition: MissingMassCalculator.h:106
min
constexpr double min()
Definition: ap_fixedTest.cxx:26
met::DeltaR
@ DeltaR
Definition: METRecoCommon.h:11
DiTauMassTools::MissingMassInput::m_metVec
XYVector m_metVec
Definition: MissingMassInput.h:79
DiTauMassTools::MissingMassCalculator::DitauMassCalculatorV9lfv
int DitauMassCalculatorV9lfv(bool refit)
Definition: MissingMassCalculator.cxx:1004
DiTauMassTools::MissingMassProb::MET
void MET(MissingMassInput &preparedInput)
Definition: MissingMassProb.cxx:1160
DiTauMassTools::MissingMassCalculator::m_niter_fit2
int m_niter_fit2
Definition: MissingMassCalculator.h:250
DiTauMassTools::MissingMassCalculator::m_tauVec1
PtEtaPhiMVector m_tauVec1
Definition: MissingMassCalculator.h:177
DiTauMassTools::MissingMassCalculator::DitauStuff::vistau2
PtEtaPhiMVector vistau2
Definition: MissingMassCalculator.h:59
DiTauMassTools::MissingMassInput::m_inputMEtT
double m_inputMEtT
Definition: MissingMassInput.h:81
CSV_InDetExporter.new
new
Definition: CSV_InDetExporter.py:145
DiTauMassTools::MissingMassCalculator::m_nsolOld
int m_nsolOld
Definition: MissingMassCalculator.h:157
DiTauMassTools::MissingMassCalculator::m_Mvis
double m_Mvis
Definition: MissingMassCalculator.h:193
DiTauMassTools::MissingMassProb::SetUseMnuProbability
void SetUseMnuProbability(bool val)
Definition: MissingMassProb.h:54
DiTauMassTools::MissingMassCalculator::PrintOtherInput
void PrintOtherInput()
Definition: MissingMassCalculator.cxx:393
DiTauMassTools::MissingMassInput::m_METsigmaL
double m_METsigmaL
Definition: MissingMassInput.h:63
DiTauMassTools::MissingMassCalculator::m_PhiProposal
double m_PhiProposal
Definition: MissingMassCalculator.h:149
DiTauMassTools::MissingMassCalculator::m_fPYfit2
std::shared_ptr< TH1F > m_fPYfit2
Definition: MissingMassCalculator.h:212
DiTauMassTools::MissingMassCalculator::m_fPhi1_all
std::shared_ptr< TH1F > m_fPhi1_all
Definition: MissingMassCalculator.h:203
DiTauMassTools::MissingMassCalculator::DitauStuff::vistau1
PtEtaPhiMVector vistau1
Definition: MissingMassCalculator.h:58
DiTauMassTools::MissingMassCalculator::m_Phi2
double m_Phi2
Definition: MissingMassCalculator.h:141
DiTauMassTools::MissingMassOutput::m_FittedMass
double m_FittedMass[MMCFitMethod::MAX]
Definition: MissingMassOutput.h:55
DiTauMassTools::MissingMassCalculator::m_meanbinToBeEvaluated
bool m_meanbinToBeEvaluated
Definition: MissingMassCalculator.h:122
DiTauMassTools::MissingMassCalculator::m_tauVec2
PtEtaPhiMVector m_tauVec2
Definition: MissingMassCalculator.h:177
DiTauMassTools::MissingMassCalculator::m_fPXfit2
std::shared_ptr< TH1F > m_fPXfit2
Definition: MissingMassCalculator.h:211
DiTauMassTools::MissingMassCalculator::m_nsigma_METscan_lfv_ll
double m_nsigma_METscan_lfv_ll
Definition: MissingMassCalculator.h:95
DiTauMassTools::MissingMassInput::m_vistau1
PtEtaPhiMVector m_vistau1
Definition: MissingMassInput.h:54
DiTauMassTools::MissingMassCalculator::m_ET2v2
double m_ET2v2
Definition: MissingMassCalculator.h:188
DiTauMassTools::MissingMassOutput::m_FitSignificance
double m_FitSignificance[MMCFitMethod::MAX]
Definition: MissingMassOutput.h:54
DiTauMassTools::MissingMassCalculator::m_MEtLMin
double m_MEtLMin
Definition: MissingMassCalculator.h:145
DiTauMassTools::MissingMassCalculator::m_mmcCalibrationSet
MMCCalibrationSet::e m_mmcCalibrationSet
Definition: MissingMassCalculator.h:65
DiTauMassTools::MissingMassCalculator::m_niterRandomLocal
int m_niterRandomLocal
Definition: MissingMassCalculator.h:71
DiTauMassTools::MissingMassCalculator::m_nsigma_METscan_hh
double m_nsigma_METscan_hh
Definition: MissingMassCalculator.h:95
DiTauMassTools::MissingMassCalculator::m_iang1low
int m_iang1low
Definition: MissingMassCalculator.h:98
DiTauMassTools::MissingMassCalculator::m_MEtP0
double m_MEtP0
Definition: MissingMassCalculator.h:144
DiTauMassTools::MissingMassCalculator::MissingMassCalculator
MissingMassCalculator(MMCCalibrationSet::e aset, std::string m_paramFilePath)
Definition: MissingMassCalculator.cxx:41
DiTauMassTools::MissingMassInput::SetLFVmode
void SetLFVmode(int val)
Definition: MissingMassInput.h:47
DiTauMassTools::MissingMassCalculator::m_fullParamSpaceScan
bool m_fullParamSpaceScan
Definition: MissingMassCalculator.h:155
DiTauMassTools::MissingMassOutput::m_FittedMassLowerError
double m_FittedMassLowerError[MMCFitMethod::MAX]
Definition: MissingMassOutput.h:57
DiTauMassTools::MissingMassCalculator::m_Mnu2Max
double m_Mnu2Max
Definition: MissingMassCalculator.h:146
DiTauMassTools::MissingMassCalculator::m_Phi20
double m_Phi20
Definition: MissingMassCalculator.h:144
drawFromPickle.cos
cos
Definition: drawFromPickle.py:36
DiTauMassTools::MissingMassCalculator::m_tauVec1M
double m_tauVec1M
Definition: MissingMassCalculator.h:179
DiTauMassTools::HistInfo::TANTHETA
@ TANTHETA
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:35
DiTauMassTools::MissingMassCalculator::m_cosPhi2
double m_cosPhi2
Definition: MissingMassCalculator.h:173
plotBeamSpotVxVal.sumw
int sumw
Definition: plotBeamSpotVxVal.py:236
DiTauMassTools::MissingMassCalculator::m_nsigma_METscan2
double m_nsigma_METscan2
Definition: MissingMassCalculator.h:94
DiTauMassTools::MissingMassOutput::m_hMfit_all
std::shared_ptr< TH1F > m_hMfit_all
Definition: MissingMassOutput.h:65
DiTauMassTools::MissingMassCalculator::m_MEtP
double m_MEtP
Definition: MissingMassCalculator.h:141
DiTauMassTools::MaxHistStrategy::MAXBINWINDOW
@ MAXBINWINDOW
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:30
DiTauMassTools::MissingMassCalculator::m_tauVec2Py
double m_tauVec2Py
Definition: MissingMassCalculator.h:181
DiTauMassTools::HistInfo::MEANBIN
@ MEANBIN
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:35
DiTauMassTools::MMCFitMethod::name
const std::string name[MAX]
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:47
xAOD::IParticle
Class providing the definition of the 4-vector interface.
Definition: Event/xAOD/xAODBase/xAODBase/IParticle.h:41
DiTauMassTools::MissingMassCalculator::m_iterNuPV3
int m_iterNuPV3
Definition: MissingMassCalculator.h:111
covarianceTool.prob
prob
Definition: covarianceTool.py:678
DiTauMassTools::MissingMassCalculator::m_nsucStop
int m_nsucStop
Definition: MissingMassCalculator.h:72
Phi_mpi_pi
__HOSTDEV__ double Phi_mpi_pi(double)
Definition: GeoRegion.cxx:7
python.TrigEgammaFastCaloHypoTool.same
def same(val, tool)
Definition: TrigEgammaFastCaloHypoTool.py:12
x
#define x
DiTauMassTools::MissingMassCalculator::m_fMEtP_split1
std::shared_ptr< TH1F > m_fMEtP_split1
Definition: MissingMassCalculator.h:217
DiTauMassTools::MissingMassCalculator::m_mTau2
double m_mTau2
Definition: MissingMassCalculator.h:140
DiTauMassTools::MissingMassCalculator::m_fMnu2_split1
std::shared_ptr< TH1F > m_fMnu2_split1
Definition: MissingMassCalculator.h:220
DiTauMassTools::MMCCalibrationSet::LFVMMC2012
@ LFVMMC2012
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:40
DiTauMassTools::MissingMassCalculator::m_fMEtP_all
std::shared_ptr< TH1F > m_fMEtP_all
Definition: MissingMassCalculator.h:199
DiTauMassTools::MissingMassCalculator::maxFitting
Double_t maxFitting(Double_t *x, Double_t *par)
Definition: MissingMassCalculator.cxx:1489
DiTauMassTools::MissingMassCalculator::SetNsigmaMETscan_ll
void SetNsigmaMETscan_ll(const double val)
Definition: MissingMassCalculator.h:387
DiTauMassTools::MissingMassCalculator::m_fUseFloatStopping
bool m_fUseFloatStopping
Definition: MissingMassCalculator.h:68
binWidth
void binWidth(TH1 *h)
Definition: listroot.cxx:80
DiTauMassTools::MissingMassCalculator::m_Mnu2
double m_Mnu2
Definition: MissingMassCalculator.h:141
DiTauMassTools::MissingMassCalculator::m_nsigma_METscan_lh
double m_nsigma_METscan_lh
Definition: MissingMassCalculator.h:95
DiTauMassTools::MissingMassCalculator::m_walkWeight
double m_walkWeight
Definition: MissingMassCalculator.h:172
DiTauMassTools::MissingMassCalculator::m_MEtPMin
double m_MEtPMin
Definition: MissingMassCalculator.h:145
DiTauMassTools::MissingMassCalculator::m_ET2v1
double m_ET2v1
Definition: MissingMassCalculator.h:187
DiTauMassTools::MissingMassCalculator::m_Mnu20
double m_Mnu20
Definition: MissingMassCalculator.h:144
DiTauMassTools::MissingMassCalculator::m_debugThisIteration
bool m_debugThisIteration
Definition: MissingMassCalculator.h:90
DiTauMassTools::fastSinCos
void fastSinCos(const double &phi, double &sinPhi, double &cosPhi)
Definition: PhysicsAnalysis/TauID/DiTauMassTools/Root/HelperFunctions.cxx:42
DiTauMassTools::MissingMassOutput::m_hMfit_allNoWeight
std::shared_ptr< TH1F > m_hMfit_allNoWeight
Definition: MissingMassOutput.h:66
DiTauMassTools::MissingMassCalculator::m_fPZfit1
std::shared_ptr< TH1F > m_fPZfit1
Definition: MissingMassCalculator.h:210
DiTauMassTools::updateDouble
bool updateDouble(const double in, double &out)
Definition: PhysicsAnalysis/TauID/DiTauMassTools/Root/HelperFunctions.cxx:75
DiTauMassTools::MissingMassInput::m_METcovphi
double m_METcovphi
Definition: MissingMassInput.h:61
dqt_zlumi_pandas.weight
int weight
Definition: dqt_zlumi_pandas.py:189
xAOD::EgammaParameters::deltaPhi1
@ deltaPhi1
difference between the cluster eta (1st sampling) and the eta of the track extrapolated to the 1st sa...
Definition: EgammaEnums.h:196
DiTauMassTools::MissingMassCalculator::m_markovNRejectNoSol
int m_markovNRejectNoSol
Definition: MissingMassCalculator.h:126
DiTauMassTools::MissingMassCalculator::m_fMEtL_split2
std::shared_ptr< TH1F > m_fMEtL_split2
Definition: MissingMassCalculator.h:225
DiTauMassTools::MissingMassInput::SetVisTauVec
void SetVisTauVec(int i, const PtEtaPhiMVector &vec)
Definition: MissingMassInput.cxx:123
DiTauMassTools::MissingMassInput::m_tauTypes
TauTypes::e m_tauTypes
Definition: MissingMassInput.h:85
DiTauMassTools::MissingMassCalculator::m_Mnu2Min
double m_Mnu2Min
Definition: MissingMassCalculator.h:145
DiTauMassTools::MissingMassCalculator::SpaceWalkerWalk
bool SpaceWalkerWalk()
Definition: MissingMassCalculator.cxx:2487
DiTauMassTools::MissingMassCalculator::m_Mnu1Range
double m_Mnu1Range
Definition: MissingMassCalculator.h:148
DiTauMassTools::HistInfo::MAXHISTINFO
@ MAXHISTINFO
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:35
DiTauMassTools::HistInfo::PROB
@ PROB
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:35
DiTauMassTools::MissingMassProb::apply
double apply(MissingMassInput &preparedInput, const int &tau_type1, const int &tau_type2, const PtEtaPhiMVector &tauvec1, const PtEtaPhiMVector &tauvec2, const PtEtaPhiMVector nuvec1, const PtEtaPhiMVector &nuvec2, bool constant=false, bool oneTau=false, bool twoTau=false)
Definition: MissingMassProb.cxx:534
DiTauMassTools::MissingMassCalculator::m_ProposalTryMnu
double m_ProposalTryMnu
Definition: MissingMassCalculator.h:136
part1
Definition: part1.py:1
DiTauMassTools
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:24
DiTauMassTools::MissingMassCalculator::m_nsigma_METscan_ll
double m_nsigma_METscan_ll
Definition: MissingMassCalculator.h:94
DiTauMassTools::MissingMassCalculator::m_fPXfit1
std::shared_ptr< TH1F > m_fPXfit1
Definition: MissingMassCalculator.h:208
DiTauMassTools::MissingMassCalculator::m_fMnu1_all
std::shared_ptr< TH1F > m_fMnu1_all
Definition: MissingMassCalculator.h:201
DiTauMassTools::MissingMassCalculator::m_fPhi1_split2
std::shared_ptr< TH1F > m_fPhi1_split2
Definition: MissingMassCalculator.h:228
DiTauMassTools::MissingMassCalculator::m_Mnu10
double m_Mnu10
Definition: MissingMassCalculator.h:144
DiTauMassTools::MaxDelPhi
double MaxDelPhi(int tau_type, double Pvis, double dRmax_tau)
Definition: PhysicsAnalysis/TauID/DiTauMassTools/Root/HelperFunctions.cxx:15
DiTauMassTools::MissingMassCalculator::m_switch1
bool m_switch1
Definition: MissingMassCalculator.h:119
met
Definition: IMETSignificance.h:24
DiTauMassTools::MissingMassCalculator::m_fMEtP_split2
std::shared_ptr< TH1F > m_fMEtP_split2
Definition: MissingMassCalculator.h:224
DiTauMassTools::MissingMassCalculator::m_Phi10
double m_Phi10
Definition: MissingMassCalculator.h:144
DiTauMassTools::MissingMassCalculator::m_nuvecsol2
std::vector< PtEtaPhiMVector > m_nuvecsol2
Definition: MissingMassCalculator.h:78
DiTauMassTools::MissingMassCalculator::m_nCallprobCalculatorV9fast
int m_nCallprobCalculatorV9fast
Definition: MissingMassCalculator.h:92
DiTauMassTools::MissingMassCalculator::m_Mnu1
double m_Mnu1
Definition: MissingMassCalculator.h:141
DiTauMassTools::MissingMassCalculator::preparedInput
MissingMassInput preparedInput
Definition: MissingMassCalculator.h:337
DiTauMassTools::MissingMassInput::m_htOffset
double m_htOffset
Definition: MissingMassInput.h:83
DiTauMassTools::MissingMassCalculator::m_Phi2Min
double m_Phi2Min
Definition: MissingMassCalculator.h:145
DiTauMassTools::MissingMassInput::m_SumEt
double m_SumEt
Definition: MissingMassInput.h:64
doubleTestComp.j1
j1
Definition: doubleTestComp.py:21
DiTauMassTools::MissingMassProb::MnuProbability
double MnuProbability(MissingMassInput &preparedInput, double mnu, double binsize)
Definition: MissingMassProb.cxx:939
DiTauMassTools::MissingMassCalculator::FinalizeSettings
void FinalizeSettings(const xAOD::IParticle *part1, const xAOD::IParticle *part2, const xAOD::MissingET *met, const int &njets)
Definition: MissingMassCalculator.cxx:2817
DiTauMassTools::MissingMassCalculator::m_fMnu1_split1
std::shared_ptr< TH1F > m_fMnu1_split1
Definition: MissingMassCalculator.h:219
TrigVtx::gamma
@ gamma
Definition: TrigParticleTable.h:26
DiTauMassTools::MissingMassInput::m_Nprong_tau2
int m_Nprong_tau2
Definition: MissingMassInput.h:59
DiTauMassTools::MissingMassCalculator::m_nu1FinalSolOldVec
std::vector< PtEtaPhiMVector > m_nu1FinalSolOldVec
Definition: MissingMassCalculator.h:160
DiTauMassTools::MissingMassCalculator::m_Ev1
double m_Ev1
Definition: MissingMassCalculator.h:192
ParseInputs.gDirectory
gDirectory
Definition: Final2012/ParseInputs.py:133
lumiFormat.i
int i
Definition: lumiFormat.py:85
DiTauMassTools::MaxHistStrategy::MAXBIN
@ MAXBIN
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:30
DiTauMassTools::MissingMassInput::m_MHtSigma2
double m_MHtSigma2
Definition: MissingMassInput.h:69
DiTauMassTools::MissingMassCalculator::m_MnuScanRange
double m_MnuScanRange
Definition: MissingMassCalculator.h:261
DiTauMassTools::MissingMassCalculator::m_tauVec1Pz
double m_tauVec1Pz
Definition: MissingMassCalculator.h:180
DiTauMassTools::MissingMassCalculator::m_ProposalTryEtau
double m_ProposalTryEtau
Definition: MissingMassCalculator.h:137
DiTauMassTools::MissingMassOutput::ClearOutput
void ClearOutput(bool fUseVerbose)
Definition: MissingMassOutput.cxx:24
DiTauMassTools::MissingMassCalculator::DitauStuff::nutau2
PtEtaPhiMVector nutau2
Definition: MissingMassCalculator.h:57
DiTauMassTools::MissingMassCalculator::~MissingMassCalculator
~MissingMassCalculator()
Definition: MissingMassCalculator.cxx:183
python.LArMinBiasAlgConfig.int
int
Definition: LArMinBiasAlgConfig.py:59
DiTauMassTools::MissingMassCalculator::m_PrintmInvWidth2Error
double m_PrintmInvWidth2Error
Definition: MissingMassCalculator.h:132
DiTauMassTools::MissingMassCalculator::m_niter_fit3
int m_niter_fit3
Definition: MissingMassCalculator.h:251
DiTauMassTools::MissingMassInput::m_MEtT
double m_MEtT
Definition: MissingMassInput.h:82
DiTauMassTools::MissingMassCalculator::m_fMfit_all
std::shared_ptr< TH1F > m_fMfit_all
Definition: MissingMassCalculator.h:198
DiTauMassTools::getLFVMode
int getLFVMode(const xAOD::IParticle *p1, const xAOD::IParticle *p2, int mmcType1, int mmcType2)
Definition: PhysicsAnalysis/TauID/DiTauMassTools/Root/HelperFunctions.cxx:83
DiTauMassTools::MissingMassCalculator::m_seed
int m_seed
Definition: MissingMassCalculator.h:107
DiTauMassTools::MissingMassInput::m_Nprong_tau1
int m_Nprong_tau1
Definition: MissingMassInput.h:58
DiTauMassTools::MissingMassCalculator::m_tauVec2P
double m_tauVec2P
Definition: MissingMassCalculator.h:182
DiTauMassTools::MissingMassCalculator::probCalculatorV9fast
int probCalculatorV9fast(const double &phi1, const double &phi2, const double &M_nu1, const double &M_nu2)
Definition: MissingMassCalculator.cxx:1785
DiTauMassTools::MissingMassCalculator::m_Phi1Max
double m_Phi1Max
Definition: MissingMassCalculator.h:146
DiTauMassTools::MMCCalibrationSet::MMC2019
@ MMC2019
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:40
DiTauMassTools::MissingMassCalculator::m_tauVec1Phi
double m_tauVec1Phi
Definition: MissingMassCalculator.h:178
xAOD::MissingET_v1
Principal data object for Missing ET.
Definition: MissingET_v1.h:25
DiTauMassTools::MissingMassProb::dTheta3d_probabilityFast
double dTheta3d_probabilityFast(MissingMassInput &preparedInput, const int &tau_type, const double &dTheta3d, const double &P_tau)
Definition: MissingMassProb.cxx:1036
DiTauMassTools::MissingMassCalculator::m_tauVec1E
double m_tauVec1E
Definition: MissingMassCalculator.h:183
DiTauMassTools::MissingMassProb::GetUseHT
bool GetUseHT()
Definition: MissingMassProb.h:49
DiTauMassTools::MaxHistStrategy::FIT
@ FIT
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:30
DiTauMassTools::MissingMassCalculator::m_randomGen
TRandom2 m_randomGen
Definition: MissingMassCalculator.h:63
DiTauMassTools::MissingMassCalculator::m_scanMnu2
bool m_scanMnu2
Definition: MissingMassCalculator.h:175
GEV
#define GEV
Definition: PrintPhotonSF.cxx:24
MissingMassCalculator.h
DiTauMassTools::MissingMassInput::SetNjet25
void SetNjet25(int val)
Definition: MissingMassInput.cxx:107
DiTauMassTools::MissingMassCalculator::SetNsigmaMETscan_hh
void SetNsigmaMETscan_hh(const double val)
Definition: MissingMassCalculator.h:389
DiTauMassTools::MissingMassCalculator::m_nsigma_METscan
double m_nsigma_METscan
Definition: MissingMassCalculator.h:94
DiTauMassTools::MMCCalibrationSet::e
e
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:40
DiTauMassTools::MissingMassProb::GetUseMnuProbability
bool GetUseMnuProbability()
Definition: MissingMassProb.h:55
DiTauMassTools::MissingMassCalculator::DitauMassCalculatorV9walk
int DitauMassCalculatorV9walk()
Definition: MissingMassCalculator.cxx:782
DiTauMassTools::MissingMassOutput::m_NSolutions
int m_NSolutions
Definition: MissingMassOutput.h:69
DiTauMassTools::MaxHistStrategy::SLIDINGWINDOW
@ SLIDINGWINDOW
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:30
DiTauMassTools::MissingMassCalculator::m_tauVec2M
double m_tauVec2M
Definition: MissingMassCalculator.h:179
DiTauMassTools::MissingMassInput::m_dataType
int m_dataType
Definition: MissingMassInput.h:60
xAOD::double
double
Definition: CompositeParticle_v1.cxx:159
DiTauMassTools::mmcType
int mmcType(const xAOD::IParticle *part)
Definition: PhysicsAnalysis/TauID/DiTauMassTools/Root/HelperFunctions.cxx:136
DiTauMassTools::MissingMassCalculator::m_nsol
int m_nsol
Definition: MissingMassCalculator.h:163
DiTauMassTools::MissingMassInput::m_fUseVerbose
bool m_fUseVerbose
Definition: MissingMassInput.h:78
DiTauMassTools::MissingMassCalculator::m_reRunWithBestMET
bool m_reRunWithBestMET
Definition: MissingMassCalculator.h:194
DiTauMassTools::MissingMassCalculator::m_sinPhi2
double m_sinPhi2
Definition: MissingMassCalculator.h:173
DiTauMassTools::MissingMassCalculator::m_iter4
int m_iter4
Definition: MissingMassCalculator.h:98
DiTauMassTools::MissingMassCalculator::m_scanMnu1
bool m_scanMnu1
Definition: MissingMassCalculator.h:175
ReadTripsProbsFromCool.denominator
denominator
Definition: ReadTripsProbsFromCool.py:96
DiTauMassTools::MissingMassCalculator::m_beamEnergy
double m_beamEnergy
Definition: MissingMassCalculator.h:96
WriteCalibToCool.swap
swap
Definition: WriteCalibToCool.py:94
DiTauMassTools::MissingMassCalculator::checkAllParamInRange
bool checkAllParamInRange()
Definition: MissingMassCalculator.cxx:2656
DiTauMassTools::MissingMassOutput::m_objvec2
PtEtaPhiMVector m_objvec2[MMCFitMethod::MAX]
Definition: MissingMassOutput.h:61
DiTauMassTools::MissingMassCalculator::m_fMfit_allNoWeight
std::shared_ptr< TH1F > m_fMfit_allNoWeight
Definition: MissingMassCalculator.h:206
DiTauMassTools::MissingMassCalculator::m_mtautauSum
double m_mtautauSum
Definition: MissingMassCalculator.h:104
DiTauMassTools::MissingMassOutput::m_RMS2MPV
double m_RMS2MPV
Definition: MissingMassOutput.h:64
DiTauMassTools::MissingMassCalculator::m_markovNRejectMetropolis
int m_markovNRejectMetropolis
Definition: MissingMassCalculator.h:127
DiTauMassTools::MissingMassCalculator::m_eTau1Max
double m_eTau1Max
Definition: MissingMassCalculator.h:153
DiTauMassTools::MissingMassCalculator::checkMEtInRange
bool checkMEtInRange()
Definition: MissingMassCalculator.cxx:2697
DiTauMassTools::MissingMassCalculator::m_fPhi2_all
std::shared_ptr< TH1F > m_fPhi2_all
Definition: MissingMassCalculator.h:204
DiTauMassTools::MissingMassOutput::m_objvec1
PtEtaPhiMVector m_objvec1[MMCFitMethod::MAX]
Definition: MissingMassOutput.h:59
DiTauMassTools::MissingMassCalculator::m_fMmass_split1
std::shared_ptr< TH1F > m_fMmass_split1
Definition: MissingMassCalculator.h:216
DiTauMassTools::MissingMassCalculator::m_iterTheta3d
int m_iterTheta3d
Definition: MissingMassCalculator.h:99
DiTauMassTools::MissingMassInput::m_HtOffset
double m_HtOffset
Definition: MissingMassInput.h:71
part2
Definition: part2.py:1
DiTauMassTools::MissingMassCalculator::m_dRmax_tau
double m_dRmax_tau
Definition: MissingMassCalculator.h:259
DiTauMassTools::MissingMassCalculator::m_nsigma_METscan_lfv_lh
double m_nsigma_METscan_lfv_lh
Definition: MissingMassCalculator.h:96
DiTauMassTools::MissingMassProb::setParamRatio
void setParamRatio(int tau, int tautype)
Definition: MissingMassProb.cxx:184
DiTauMassTools::MissingMassCalculator::m_nsolmax
int m_nsolmax
Definition: MissingMassCalculator.h:70
DiTauMassTools::MissingMassCalculator::m_markovNFullScan
int m_markovNFullScan
Definition: MissingMassCalculator.h:125
DiTauMassTools::MissingMassCalculator::m_markovNAccept
int m_markovNAccept
Definition: MissingMassCalculator.h:128
DQPostProcessTest.outFile
outFile
Comment Out Those You do not wish to run.
Definition: DQPostProcessTest.py:37
DiTauMassTools::MissingMassCalculator::m_testdiscri1
int m_testdiscri1
Definition: MissingMassCalculator.h:114
DiTauMassTools::MissingMassCalculator::m_iter3
int m_iter3
Definition: MissingMassCalculator.h:98
DiTauMassTools::MissingMassCalculator::DitauStuff::RMSoverMPV
double RMSoverMPV
Definition: MissingMassCalculator.h:60
DiTauMassTools::MissingMassCalculator::m_Phi2Range
double m_Phi2Range
Definition: MissingMassCalculator.h:148
trigbs_pickEvents.num
num
Definition: trigbs_pickEvents.py:76
DiTauMassTools::MissingMassOutput::m_NSuccesses
int m_NSuccesses
Definition: MissingMassOutput.h:68
DiTauMassTools::MissingMassCalculator::m_Phi1Range
double m_Phi1Range
Definition: MissingMassCalculator.h:148
DiTauMassTools::MissingMassCalculator::RunMissingMassCalculator
int RunMissingMassCalculator(const xAOD::IParticle *part1, const xAOD::IParticle *part2, const xAOD::MissingET *met, const int &njets)
Definition: MissingMassCalculator.cxx:187
DiTauMassTools::MissingMassCalculator::m_tauVec2Px
double m_tauVec2Px
Definition: MissingMassCalculator.h:181
DiTauMassTools::HistInfo::RMS
@ RMS
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:35
DiTauMassTools::MissingMassInput::m_MetVec
XYVector m_MetVec
Definition: MissingMassInput.h:53
DiTauMassTools::MissingMassInput::m_type_visTau2
int m_type_visTau2
Definition: MissingMassInput.h:57
DiTauMassTools::MissingMassCalculator::m_nu1FinalSolVec
std::vector< PtEtaPhiMVector > m_nu1FinalSolVec
Definition: MissingMassCalculator.h:166
DiTauMassTools::MissingMassCalculator::m_fPZfit2
std::shared_ptr< TH1F > m_fPZfit2
Definition: MissingMassCalculator.h:213
DiTauMassTools::MissingMassCalculator::m_niter_fit1
int m_niter_fit1
Definition: MissingMassCalculator.h:249
DiTauMassTools::MissingMassCalculator::DitauStuff
Definition: MissingMassCalculator.h:53
DiTauMassTools::MissingMassProb::MetProbability
double MetProbability(MissingMassInput &preparedInput, const double &met1, const double &met2, const double &MetSigma1, const double &MetSigma2)
Definition: MissingMassProb.cxx:553
DiTauMassTools::MissingMassOutput::m_FittedMassUpperError
double m_FittedMassUpperError[MMCFitMethod::MAX]
Definition: MissingMassOutput.h:56
debug
const bool debug
Definition: MakeUncertaintyPlots.cxx:53
DiTauMassTools::MissingMassInput::m_DelPhiTT
double m_DelPhiTT
Definition: MissingMassInput.h:67
createCoolChannelIdFile.par
par
Definition: createCoolChannelIdFile.py:29
VP1PartSpect::E
@ E
Definition: VP1PartSpectFlags.h:21
ActsTrk::to_string
std::string to_string(const DetectorType &type)
Definition: GeometryDefs.h:34
plotBeamSpotMon.b
b
Definition: plotBeamSpotMon.py:77
DiTauMassTools::MissingMassCalculator::m_eTau2Max
double m_eTau2Max
Definition: MissingMassCalculator.h:154
DiTauMassTools::MissingMassCalculator::m_proposalTryMEt
double m_proposalTryMEt
Definition: MissingMassCalculator.h:134
DiTauMassTools::MissingMassCalculator::m_probFinalSolOldVec
std::vector< double > m_probFinalSolOldVec
Definition: MissingMassCalculator.h:158
DiTauMassTools::MissingMassCalculator::m_tauvecsol1
std::vector< PtEtaPhiMVector > m_tauvecsol1
Definition: MissingMassCalculator.h:80
DiTauMassTools::MissingMassInput::PrintInputInfo
void PrintInputInfo()
Definition: MissingMassInput.cxx:58
DiTauMassTools::MissingMassCalculator::m_mtautauFinalSolOldVec
std::vector< double > m_mtautauFinalSolOldVec
Definition: MissingMassCalculator.h:159
DiTauMassTools::MissingMassProb::SetUseTauProbability
void SetUseTauProbability(bool val)
Definition: MissingMassProb.h:51
DiTauMassTools::MissingMassOutput::m_AveSolRMS
double m_AveSolRMS
Definition: MissingMassOutput.h:71
DiTauMassTools::MissingMassCalculator::m_tauVec1P
double m_tauVec1P
Definition: MissingMassCalculator.h:182
DiTauMassTools::MissingMassCalculator::m_eTau2Min
double m_eTau2Min
Definition: MissingMassCalculator.h:154
DiTauMassTools::MissingMassCalculator::m_fMfit_allGraph
std::shared_ptr< TGraph > m_fMfit_allGraph
Definition: MissingMassCalculator.h:205
xAOD::EgammaParameters::deltaPhi2
@ deltaPhi2
difference between the cluster phi (second sampling) and the phi of the track extrapolated to the sec...
Definition: EgammaEnums.h:204
DiTauMassTools::MissingMassCalculator::m_iter0
int m_iter0
Definition: MissingMassCalculator.h:110
DiTauMassTools::MissingMassCalculator::m_fPhi1_split1
std::shared_ptr< TH1F > m_fPhi1_split1
Definition: MissingMassCalculator.h:221
DiTauMassTools::MissingMassCalculator::TailCleanUp
int TailCleanUp(const PtEtaPhiMVector &vis1, const PtEtaPhiMVector &nu1, const PtEtaPhiMVector &vis2, const PtEtaPhiMVector &nu2, const double &mmc_mass, const double &vis_mass, const double &eff_mass, const double &dphiTT)
Definition: MissingMassCalculator.cxx:1982
DiTauMassTools::MissingMassCalculator::m_iang1high
int m_iang1high
Definition: MissingMassCalculator.h:98
DiTauMassTools::fixPhiRange
double fixPhiRange(const double &phi)
Definition: PhysicsAnalysis/TauID/DiTauMassTools/Root/HelperFunctions.cxx:21
DiTauMassTools::MissingMassCalculator::m_tauVec2Pz
double m_tauVec2Pz
Definition: MissingMassCalculator.h:181
DiTauMassTools::MissingMassCalculator::DitauStuff::Sign_best
double Sign_best
Definition: MissingMassCalculator.h:55
DiTauMassTools::MissingMassInput::m_fUseTailCleanup
bool m_fUseTailCleanup
Definition: MissingMassInput.h:77
DiTauMassTools::HistInfo::INTEGRAL
@ INTEGRAL
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:35
DiTauMassTools::MissingMassCalculator::m_MnuProposal
double m_MnuProposal
Definition: MissingMassCalculator.h:149
DiTauMassTools::MissingMassCalculator::m_fPhi2_split2
std::shared_ptr< TH1F > m_fPhi2_split2
Definition: MissingMassCalculator.h:229
DiTauMassTools::MissingMassOutput::m_nuvec1
PtEtaPhiMVector m_nuvec1[MMCFitMethod::MAX]
Definition: MissingMassOutput.h:58
DiTauMassTools::MissingMassCalculator::maxFromHist
double maxFromHist(TH1F *theHist, std::vector< double > &histInfo, const MaxHistStrategy::e maxHistStrategy=MaxHistStrategy::FIT, const int winHalfWidth=2, bool debug=false)
Definition: MissingMassCalculator.cxx:1509
DiTauMassTools::MissingMassCalculator::metvec_tmp
XYVector metvec_tmp
Definition: MissingMassCalculator.h:415
DiTauMassTools::MissingMassOutput::m_FittedMetVec
XYVector m_FittedMetVec[MMCFitMethod::MAX]
Definition: MissingMassOutput.h:63
DiTauMassTools::MissingMassCalculator::DoOutputInfo
void DoOutputInfo()
Definition: MissingMassCalculator.cxx:296
DiTauMassTools::MissingMassCalculator::m_RMSStop
int m_RMSStop
Definition: MissingMassCalculator.h:254
DiTauMassTools::MissingMassCalculator::m_NsucStop
int m_NsucStop
Definition: MissingMassCalculator.h:253
DiTauMassTools::MissingMassOutput::m_nuvec2
PtEtaPhiMVector m_nuvec2[MMCFitMethod::MAX]
Definition: MissingMassOutput.h:60
DiTauMassTools::MissingMassCalculator::m_eTau1Min
double m_eTau1Min
Definition: MissingMassCalculator.h:153
DiTauMassTools::MissingMassCalculator::m_testptn2
int m_testptn2
Definition: MissingMassCalculator.h:113
DiTauMassTools::HistInfo::DISCRI
@ DISCRI
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:35
DiTauMassTools::MissingMassCalculator::SetNsigmaMETscan_lh
void SetNsigmaMETscan_lh(const double val)
Definition: MissingMassCalculator.h:388
DiTauMassTools::MissingMassInput::m_MEtX
double m_MEtX
Definition: MissingMassInput.h:82
DiTauMassTools::MissingMassInput::m_vistau2
PtEtaPhiMVector m_vistau2
Definition: MissingMassInput.h:55
DiTauMassTools::MissingMassInput::SetVisTauType
void SetVisTauType(int i, int tautype)
Definition: MissingMassInput.cxx:117
DiTauMassTools::MissingMassCalculator::m_iter1
int m_iter1
Definition: MissingMassCalculator.h:98
DiTauMassTools::MissingMassCalculator::m_fUseEfficiencyRecovery
bool m_fUseEfficiencyRecovery
Definition: MissingMassCalculator.h:67
DiTauMassTools::MissingMassCalculator::m_MEtL
double m_MEtL
Definition: MissingMassCalculator.h:141
DiTauMassTools::MissingMassCalculator::m_NiterRandom
int m_NiterRandom
Definition: MissingMassCalculator.h:252
DiTauMassTools::MissingMassProb
Definition: MissingMassProb.h:28
DiTauMassTools::MissingMassCalculator::m_MEtPRange
double m_MEtPRange
Definition: MissingMassCalculator.h:148
DiTauMassTools::MissingMassOutput::m_NTrials
int m_NTrials
Definition: MissingMassOutput.h:67
DiTauMassTools::MissingMassCalculator::m_eTau1
double m_eTau1
Definition: MissingMassCalculator.h:142
DiTauMassTools::MissingMassInput::m_MEtY
double m_MEtY
Definition: MissingMassInput.h:82
DiTauMassTools::MissingMassCalculator::m_dTheta3d_binMin
double m_dTheta3d_binMin
Definition: MissingMassCalculator.h:257
DiTauMassTools::MissingMassCalculator::m_MEtProposal
double m_MEtProposal
Definition: MissingMassCalculator.h:149
DiTauMassTools::MissingMassCalculator::m_MEtLRange
double m_MEtLRange
Definition: MissingMassCalculator.h:148
DiTauMassTools::MissingMassCalculator::m_iterNsuc
int m_iterNsuc
Definition: MissingMassCalculator.h:118
a
TList * a
Definition: liststreamerinfos.cxx:10
DiTauMassTools::MissingMassInput::ClearInput
void ClearInput()
Definition: MissingMassInput.cxx:24
DiTauMassTools::MissingMassProb::setParamNuMass
void setParamNuMass()
Definition: MissingMassProb.cxx:147
DiTauMassTools::MMCFitMethod::MAXW
@ MAXW
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:46
Pythia8_RapidityOrderMPI.val
val
Definition: Pythia8_RapidityOrderMPI.py:14
DiTauMassTools::MissingMassProb::TauProbabilityLFV
double TauProbabilityLFV(MissingMassInput &preparedInput, const int &type1, const PtEtaPhiMVector &vis1, const PtEtaPhiMVector &nu1)
Definition: MissingMassProb.cxx:614
DiTauMassTools::MissingMassInput::m_inputMEtX
double m_inputMEtX
Definition: MissingMassInput.h:81
DiTauMassTools::MissingMassOutput::m_FitStatus
int m_FitStatus
Definition: MissingMassOutput.h:53
DiTauMassTools::MissingMassCalculator::m_eTau20
double m_eTau20
Definition: MissingMassCalculator.h:143
DiTauMassTools::MissingMassCalculator::m_fPYfit1
std::shared_ptr< TH1F > m_fPYfit1
Definition: MissingMassCalculator.h:209
DiTauMassTools::MissingMassCalculator::m_markovCountDuplicate
int m_markovCountDuplicate
Definition: MissingMassCalculator.h:124
DiTauMassTools::MMCFitMethod::MLNU3P
@ MLNU3P
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:46
DiTauMassTools::MMCCalibrationSet::MAXMMCCALIBRATIONSET
@ MAXMMCCALIBRATIONSET
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:40
DiTauMassTools::HistInfo::FITLENGTH
@ FITLENGTH
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:35
DiTauMassTools::HistInfo::RMSVSDISCRI
@ RMSVSDISCRI
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:35
DiTauMassTools::MissingMassCalculator::m_Mnu1Min
double m_Mnu1Min
Definition: MissingMassCalculator.h:145
DiTauMassTools::MissingMassInput::m_jet4vecs
std::vector< PtEtaPhiMVector > m_jet4vecs
Definition: MissingMassInput.h:65
DiTauMassTools::MissingMassInput::m_type_visTau1
int m_type_visTau1
Definition: MissingMassInput.h:56
DiTauMassTools::MMCCalibrationSet::MMC2016MC15C
@ MMC2016MC15C
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:40
DiTauMassTools::MissingMassCalculator::m_rmsStop
int m_rmsStop
Definition: MissingMassCalculator.h:73
DiTauMassTools::MissingMassInput::m_METresSyst
int m_METresSyst
Definition: MissingMassInput.h:75
DiTauMassTools::MissingMassCalculator::m_tauVec1Px
double m_tauVec1Px
Definition: MissingMassCalculator.h:180
DiTauMassTools::MissingMassCalculator::handleSolutions
void handleSolutions()
Definition: MissingMassCalculator.cxx:2038
DiTauMassTools::MissingMassCalculator::m_SaveLlhHisto
bool m_SaveLlhHisto
Definition: MissingMassCalculator.h:90
DiTauMassTools::MissingMassCalculator::m_m2Nu2
double m_m2Nu2
Definition: MissingMassCalculator.h:186
DiTauMassTools::MissingMassCalculator::SetUseFloatStopping
void SetUseFloatStopping(const bool val)
Definition: MissingMassCalculator.cxx:3123
DiTauMassTools::MissingMassCalculator::m_RndmSeedAltering
int m_RndmSeedAltering
Definition: MissingMassCalculator.h:255
DiTauMassTools::MissingMassInput::m_InputReorder
int m_InputReorder
Definition: MissingMassInput.h:73
DiTauMassTools::MissingMassOutput::m_totalvec
PtEtaPhiMVector m_totalvec[MMCFitMethod::MAX]
Definition: MissingMassOutput.h:62
DiTauMassTools::MissingMassCalculator::m_iter2
int m_iter2
Definition: MissingMassCalculator.h:98
DiTauMassTools::HistInfo::TANTHETAW
@ TANTHETAW
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:35
DiTauMassTools::MissingMassOutput::m_SumW
double m_SumW
Definition: MissingMassOutput.h:70
DiTauMassTools::MissingMassCalculator::m_nu2FinalSolVec
std::vector< PtEtaPhiMVector > m_nu2FinalSolVec
Definition: MissingMassCalculator.h:167
python.TrigEgammaMonitorHelper.TH1F
def TH1F(name, title, nxbins, bins_par2, bins_par3=None, path='', **kwargs)
Definition: TrigEgammaMonitorHelper.py:24
DiTauMassTools::MissingMassCalculator::m_probFinalSolVec
std::vector< double > m_probFinalSolVec
Definition: MissingMassCalculator.h:164
DiTauMassTools::MissingMassCalculator::m_Phi1Min
double m_Phi1Min
Definition: MissingMassCalculator.h:145
DiTauMassTools::MissingMassInput::m_beamEnergy
double m_beamEnergy
Definition: MissingMassInput.h:72
DiTauMassTools::MissingMassCalculator::m_meanbinStop
double m_meanbinStop
Definition: MissingMassCalculator.h:74
DiTauMassTools::MissingMassCalculator::m_mtautauFinalSolVec
std::vector< double > m_mtautauFinalSolVec
Definition: MissingMassCalculator.h:165
DiTauMassTools::MMCCalibrationSet::UPGRADE
@ UPGRADE
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:40
DiTauMassTools::MissingMassCalculator::m_fMmass_split2
std::shared_ptr< TH1F > m_fMmass_split2
Definition: MissingMassCalculator.h:223
DiTauMassTools::MissingMassCalculator::m_PrintmMeanError
double m_PrintmMeanError
Definition: MissingMassCalculator.h:131
DiTauMassTools::MissingMassCalculator::m_totalProbSum
double m_totalProbSum
Definition: MissingMassCalculator.h:103
drawFromPickle.sin
sin
Definition: drawFromPickle.py:36
DiTauMassTools::MissingMassCalculator::m_MEtLMax
double m_MEtLMax
Definition: MissingMassCalculator.h:146
DiTauMassTools::MissingMassCalculator::m_dTheta3d_binMax
double m_dTheta3d_binMax
Definition: MissingMassCalculator.h:258
DiTauMassTools::MissingMassInput::SetMetVec
void SetMetVec(const XYVector &vec)
Definition: MissingMassInput.cxx:130
MuonParameters::beta
@ beta
Definition: MuonParamDefs.h:144
DiTauMassTools::MissingMassInput::m_fUseDefaults
bool m_fUseDefaults
Definition: MissingMassInput.h:76
DiTauMassTools::MissingMassCalculator::m_tauvecprob2
std::vector< double > m_tauvecprob2
Definition: MissingMassCalculator.h:83
DiTauMassTools::MissingMassInput::SetSumEt
void SetSumEt(double sumEt)
Definition: MissingMassInput.cxx:112
DiTauMassTools::MissingMassCalculator::m_testptn1
int m_testptn1
Definition: MissingMassCalculator.h:112
DiTauMassTools::Angle
double Angle(const VectorType1 &vec1, const VectorType2 &vec2)
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:61
pow
constexpr int pow(int base, int exp) noexcept
Definition: ap_fixedTest.cxx:15
DiTauMassTools::HistInfo::CHI2
@ CHI2
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:35
DiTauMassTools::MissingMassCalculator::m_cosPhi1
double m_cosPhi1
Definition: MissingMassCalculator.h:173
DiTauMassTools::MissingMassCalculator::m_nsolfinalmax
int m_nsolfinalmax
Definition: MissingMassCalculator.h:70
DiTauMassTools::MissingMassCalculator::PrintResults
void PrintResults()
Definition: MissingMassCalculator.cxx:435
python.compressB64.c
def c
Definition: compressB64.py:93
updateCoolNtuple.limit
int limit
Definition: updateCoolNtuple.py:45
DiTauMassTools::MissingMassCalculator::m_prob_tmp
double m_prob_tmp
Definition: MissingMassCalculator.h:101
DiTauMassTools::MissingMassInput::m_inputMEtY
double m_inputMEtY
Definition: MissingMassInput.h:81
DiTauMassTools::MissingMassCalculator::m_sinPhi1
double m_sinPhi1
Definition: MissingMassCalculator.h:173
DiTauMassTools::MissingMassCalculator::m_Mnu2Range
double m_Mnu2Range
Definition: MissingMassCalculator.h:148
DiTauMassTools::MissingMassCalculator::m_tauvecsol2
std::vector< PtEtaPhiMVector > m_tauvecsol2
Definition: MissingMassCalculator.h:81
DiTauMassTools::MissingMassCalculator::dTheta3DLimit
double dTheta3DLimit(const int &tau_type, const int &limit_code, const double &P_tau)
Definition: MissingMassCalculator.cxx:2713
DiTauMassTools::MissingMassCalculator::m_Mnu1Max
double m_Mnu1Max
Definition: MissingMassCalculator.h:146
DiTauMassTools::MissingMassCalculator::m_metCovPhiSin
double m_metCovPhiSin
Definition: MissingMassCalculator.h:150
doubleTestComp.j2
j2
Definition: doubleTestComp.py:22
DiTauMassTools::MissingMassInput::m_LFVmode
int m_LFVmode
Definition: MissingMassInput.h:84
DiTauMassTools::MissingMassCalculator::OutputInfo
MissingMassOutput OutputInfo
Definition: MissingMassCalculator.h:338
DiTauMassTools::MissingMassCalculator::m_Meff
double m_Meff
Definition: MissingMassCalculator.h:193
DiTauMassTools::MissingMassCalculator::m_nu2FinalSolOldVec
std::vector< PtEtaPhiMVector > m_nu2FinalSolOldVec
Definition: MissingMassCalculator.h:161
DiTauMassTools::MissingMassCalculator::Prob
MissingMassProb * Prob
Definition: MissingMassCalculator.h:339
DiTauMassTools::MissingMassCalculator::m_tauVec1Py
double m_tauVec1Py
Definition: MissingMassCalculator.h:180
DiTauMassTools::MissingMassCalculator::m_mTau
double m_mTau
Definition: MissingMassCalculator.h:140
DiTauMassTools::MissingMassCalculator::m_eTau10
double m_eTau10
Definition: MissingMassCalculator.h:143
DiTauMassTools::MissingMassCalculator::m_fDitauStuffFit
DitauStuff m_fDitauStuffFit
Definition: MissingMassCalculator.h:246
DiTauMassTools::MissingMassCalculator::m_tauVec2E
double m_tauVec2E
Definition: MissingMassCalculator.h:184
DiTauMassTools::TauTypes::ll
@ ll
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:53
DiTauMassTools::MissingMassCalculator::m_E2v2
double m_E2v2
Definition: MissingMassCalculator.h:190
DiTauMassTools::MissingMassCalculator::ClearDitauStuff
void ClearDitauStuff(DitauStuff &fStuff)
Definition: MissingMassCalculator.cxx:281
DiTauMassTools::MissingMassCalculator::m_nosol1
int m_nosol1
Definition: MissingMassCalculator.h:116
DiTauMassTools::MissingMassCalculator::m_fMEtL_all
std::shared_ptr< TH1F > m_fMEtL_all
Definition: MissingMassCalculator.h:200