CMS 3D CMS Logo

List of all members | Classes | Public Member Functions | Private Member Functions | Private Attributes
TtFullLepKinSolver Class Reference

#include <TtFullLepKinSolver.h>

Classes

struct  NeutrinoSolution
 

Public Member Functions

TtDilepEvtSolution addKinSolInfo (TtDilepEvtSolution *asol)
 
NeutrinoSolution getNuSolution (const TLorentzVector &LV_l, const TLorentzVector &LV_l_, const TLorentzVector &LV_b, const TLorentzVector &LV_b_)
 
void SetConstraints (const double xx=0, const double yy=0)
 
 TtFullLepKinSolver ()
 default constructor More...
 
 TtFullLepKinSolver (const double, const double, const double, const std::vector< double > &, const double=80.4, const double=4.8)
 constructor with parameters to configure the top-mass scan and the neutrino spectrum More...
 
void useWeightFromMC (bool useMC)
 
 ~TtFullLepKinSolver ()
 destructor More...
 

Private Member Functions

int cubic (const double *c_coeff, double *c_sol) const
 
void FindCoeff (const TLorentzVector &al, const TLorentzVector &l, const TLorentzVector &b_al, const TLorentzVector &b_l, const double mt, const double mat, const double pxboost, const double pyboost, double *q_coeff)
 
int quartic (double *q_coeff, double *q_sol) const
 
double sqr (const double x) const
 
void SWAP (double &realone, double &realtwo) const
 
void TopRec (const TLorentzVector &al, const TLorentzVector &l, const TLorentzVector &b_al, const TLorentzVector &b_l, const double sol)
 
double WeightSolfromMC () const
 
double WeightSolfromShape () const
 use the parametrized event shape to obtain the solution weight. More...
 

Private Attributes

double C
 
double D
 
TF2 * EventShape_
 Event shape. More...
 
double F
 
TLorentzVector genLV_n
 provisional More...
 
TLorentzVector genLV_n_
 
double k16
 
double k26
 
double k36
 
double k46
 
double k51
 
double k56
 
double k61
 
TLorentzVector LV_n
 
TLorentzVector LV_n_
 
TLorentzVector LV_t
 
TLorentzVector LV_t_
 
TLorentzVector LV_tt_t
 
TLorentzVector LV_tt_t_
 
double m1
 
double m2
 
double m3
 
const double mb
 
const double mw
 
double n1
 
double n2
 
double n3
 
double pom
 
double pxmiss_
 
double pymiss_
 
const double topmass_begin
 
const double topmass_end
 
const double topmass_step
 
bool useMCforBest_
 flag to swith from WeightSolfromMC() to WeightSolfromShape() More...
 

Detailed Description

Definition at line 25 of file TtFullLepKinSolver.h.

Constructor & Destructor Documentation

◆ TtFullLepKinSolver() [1/2]

TtFullLepKinSolver::TtFullLepKinSolver ( )

default constructor

Definition at line 4 of file TtFullLepKinSolver.cc.

References EventShape_.

5  : topmass_begin(0), topmass_end(0), topmass_step(0), mw(80.4), mb(4.8), pxmiss_(0), pymiss_(0) {
6  // That crude parametrisation has been obtained from a fit of O(1000) pythia events.
7  // It is normalized to 1.
8  EventShape_ = new TF2("landau2D", "[0]*TMath::Landau(x,[1],[2],0)*TMath::Landau(y,[3],[4],0)", 0, 500, 0, 500);
9  EventShape_->SetParameters(30.7137, 56.2880, 23.0744, 59.1015, 24.9145);
10 }
const double topmass_end
const double topmass_step
TF2 * EventShape_
Event shape.
const double topmass_begin

◆ TtFullLepKinSolver() [2/2]

TtFullLepKinSolver::TtFullLepKinSolver ( const double  b,
const double  e,
const double  s,
const std::vector< double > &  nupars,
const double  mW = 80.4,
const double  mB = 4.8 
)

constructor with parameters to configure the top-mass scan and the neutrino spectrum

Definition at line 12 of file TtFullLepKinSolver.cc.

References EventShape_.

14  : topmass_begin(b), topmass_end(e), topmass_step(s), mw(mW), mb(mB), pxmiss_(0), pymiss_(0) {
15  EventShape_ = new TF2("landau2D", "[0]*TMath::Landau(x,[1],[2],0)*TMath::Landau(y,[3],[4],0)", 0, 500, 0, 500);
16  EventShape_->SetParameters(nupars[0], nupars[1], nupars[2], nupars[3], nupars[4]);
17 }
const double topmass_end
const double topmass_step
TF2 * EventShape_
Event shape.
double b
Definition: hdecay.h:118
const double topmass_begin

◆ ~TtFullLepKinSolver()

TtFullLepKinSolver::~TtFullLepKinSolver ( )

destructor

Definition at line 22 of file TtFullLepKinSolver.cc.

References EventShape_.

22 { delete EventShape_; }
TF2 * EventShape_
Event shape.

Member Function Documentation

◆ addKinSolInfo()

TtDilepEvtSolution TtFullLepKinSolver::addKinSolInfo ( TtDilepEvtSolution asol)

Definition at line 24 of file TtFullLepKinSolver.cc.

References reco::LeafCandidate::energy(), FindCoeff(), genLV_n, genLV_n_, TtDilepEvtSolution::getCalJetB(), TtDilepEvtSolution::getCalJetBbar(), TtDilepEvtSolution::getElectronm(), TtDilepEvtSolution::getElectronp(), TtDilepEvtSolution::getGenN(), TtDilepEvtSolution::getGenNbar(), TtDilepEvtSolution::getMuonm(), TtDilepEvtSolution::getMuonp(), TtDilepEvtSolution::getTaum(), TtDilepEvtSolution::getTaup(), TtDilepEvtSolution::getWmDecay(), TtDilepEvtSolution::getWpDecay(), TtSemiLepEvtBuilder_cfi::mt, reco::LeafCandidate::px(), pxmiss_, reco::LeafCandidate::py(), pymiss_, reco::LeafCandidate::pz(), quartic(), TtDilepEvtSolution::setRecTopMass(), TtDilepEvtSolution::setRecWeightMax(), topmass_begin, topmass_end, topmass_step, TopRec(), useMCforBest_, mps_merge::weight, WeightSolfromMC(), and WeightSolfromShape().

Referenced by TtDilepEvtSolutionMaker::produce().

24  {
25  TtDilepEvtSolution fitsol(*asol);
26 
27  //antilepton and lepton
28  TLorentzVector LV_e, LV_e_;
29  //b and bbar quark
30  TLorentzVector LV_b, LV_b_;
31 
32  bool hasMCinfo = true;
33  if (fitsol.getGenN()) { // protect against non-dilept genevents
34  genLV_n = TLorentzVector(
35  fitsol.getGenN()->px(), fitsol.getGenN()->py(), fitsol.getGenN()->pz(), fitsol.getGenN()->energy());
36  } else
37  hasMCinfo = false;
38 
39  if (fitsol.getGenNbar()) { // protect against non-dilept genevents
40  genLV_n_ = TLorentzVector(
41  fitsol.getGenNbar()->px(), fitsol.getGenNbar()->py(), fitsol.getGenNbar()->pz(), fitsol.getGenNbar()->energy());
42  } else
43  hasMCinfo = false;
44  // if MC is to be used to select the best top mass and is not available,
45  // then nothing can be done. Stop here.
46  if (useMCforBest_ && !hasMCinfo)
47  return fitsol;
48 
49  // first lepton
50  if (fitsol.getWpDecay() == "muon") {
51  LV_e = TLorentzVector(
52  fitsol.getMuonp().px(), fitsol.getMuonp().py(), fitsol.getMuonp().pz(), fitsol.getMuonp().energy());
53  } else if (fitsol.getWpDecay() == "electron") {
54  LV_e = TLorentzVector(fitsol.getElectronp().px(),
55  fitsol.getElectronp().py(),
56  fitsol.getElectronp().pz(),
57  fitsol.getElectronp().energy());
58  } else if (fitsol.getWpDecay() == "tau") {
59  LV_e =
60  TLorentzVector(fitsol.getTaup().px(), fitsol.getTaup().py(), fitsol.getTaup().pz(), fitsol.getTaup().energy());
61  }
62 
63  // second lepton
64  if (fitsol.getWmDecay() == "muon") {
65  LV_e_ = TLorentzVector(
66  fitsol.getMuonm().px(), fitsol.getMuonm().py(), fitsol.getMuonm().pz(), fitsol.getMuonm().energy());
67  } else if (fitsol.getWmDecay() == "electron") {
68  LV_e_ = TLorentzVector(fitsol.getElectronm().px(),
69  fitsol.getElectronm().py(),
70  fitsol.getElectronm().pz(),
71  fitsol.getElectronm().energy());
72  } else if (fitsol.getWmDecay() == "tau") {
73  LV_e_ =
74  TLorentzVector(fitsol.getTaum().px(), fitsol.getTaum().py(), fitsol.getTaum().pz(), fitsol.getTaum().energy());
75  }
76 
77  // first jet
78  LV_b = TLorentzVector(
79  fitsol.getCalJetB().px(), fitsol.getCalJetB().py(), fitsol.getCalJetB().pz(), fitsol.getCalJetB().energy());
80 
81  // second jet
82  LV_b_ = TLorentzVector(fitsol.getCalJetBbar().px(),
83  fitsol.getCalJetBbar().py(),
84  fitsol.getCalJetBbar().pz(),
85  fitsol.getCalJetBbar().energy());
86 
87  //loop on top mass parameter
88  double weightmax = -1e30;
89  double mtmax = 0;
90  for (double mt = topmass_begin; mt < topmass_end + 0.5 * topmass_step; mt += topmass_step) {
91  //cout << "mt = " << mt << endl;
92  double q_coeff[5], q_sol[4];
93  FindCoeff(LV_e, LV_e_, LV_b, LV_b_, mt, mt, pxmiss_, pymiss_, q_coeff);
94  int NSol = quartic(q_coeff, q_sol);
95 
96  //loop on all solutions
97  for (int isol = 0; isol < NSol; isol++) {
98  TopRec(LV_e, LV_e_, LV_b, LV_b_, q_sol[isol]);
100  if (weight > weightmax) {
101  weightmax = weight;
102  mtmax = mt;
103  }
104  }
105 
106  //for (int i=0;i<5;i++) cout << " q_coeff["<<i<< "]= " << q_coeff[i];
107  //cout << endl;
108 
109  //for (int i=0;i<4;i++) cout << " q_sol["<<i<< "]= " << q_sol[i];
110  //cout << endl;
111  //cout << "NSol_" << NSol << endl;
112  }
113 
114  fitsol.setRecTopMass(mtmax);
115  fitsol.setRecWeightMax(weightmax);
116 
117  return fitsol;
118 }
double WeightSolfromMC() const
const double topmass_end
void FindCoeff(const TLorentzVector &al, const TLorentzVector &l, const TLorentzVector &b_al, const TLorentzVector &b_l, const double mt, const double mat, const double pxboost, const double pyboost, double *q_coeff)
Definition: weight.py:1
bool useMCforBest_
flag to swith from WeightSolfromMC() to WeightSolfromShape()
double WeightSolfromShape() const
use the parametrized event shape to obtain the solution weight.
TLorentzVector genLV_n
provisional
const double topmass_step
int quartic(double *q_coeff, double *q_sol) const
TLorentzVector genLV_n_
void TopRec(const TLorentzVector &al, const TLorentzVector &l, const TLorentzVector &b_al, const TLorentzVector &b_l, const double sol)
const double topmass_begin

◆ cubic()

int TtFullLepKinSolver::cubic ( const double *  c_coeff,
double *  c_sol 
) const
private

Definition at line 372 of file TtFullLepKinSolver.cc.

References funct::cos(), h, mps_fire::i, AlCaHLTBitMon_ParallelJobs::p, phi, Pi, submitPVResolutionJobs::q, sqr(), mathSSE::sqrt(), SWAP(), and w().

Referenced by quartic().

372  {
373  unsigned nreal;
374  double w, p, q, dis, h, phi;
375 
376  if (coeffs[3] != 0.0) {
377  /* cubic problem? */
378  w = coeffs[2] / (3 * coeffs[3]);
379  p = sqr(coeffs[1] / (3 * coeffs[3]) - sqr(w)) * (coeffs[1] / (3 * coeffs[3]) - sqr(w));
380  q = -0.5 * (2 * sqr(w) * w - (coeffs[1] * w - coeffs[0]) / coeffs[3]);
381  dis = sqr(q) + p;
382  /* discriminant */
383  if (dis < 0.0) {
384  /* 3 real solutions */
385  h = q / sqrt(-p);
386  if (h > 1.0)
387  h = 1.0;
388  /* confine the argument of */
389  if (h < -1.0)
390  h = -1.0;
391  /* acos to [-1;+1] */
392  phi = acos(h);
393  p = 2 * TMath::Power(-p, 1.0 / 6.0);
394  for (unsigned i = 0; i < 3; i++)
395  koreny[i] = p * cos((phi + 2 * i * TMath::Pi()) / 3.0) - w;
396  if (koreny[1] < koreny[0])
397  SWAP(koreny[0], koreny[1]);
398  /* sort results */
399  if (koreny[2] < koreny[1])
400  SWAP(koreny[1], koreny[2]);
401  if (koreny[1] < koreny[0])
402  SWAP(koreny[0], koreny[1]);
403  nreal = 3;
404  } else {
405  /* only one real solution */
406  dis = sqrt(dis);
407  h = TMath::Power(fabs(q + dis), 1.0 / 3.0);
408  p = TMath::Power(fabs(q - dis), 1.0 / 3.0);
409  koreny[0] = ((q + dis > 0.0) ? h : -h) + ((q - dis > 0.0) ? p : -p) - w;
410  nreal = 1;
411  }
412 
413  /* Perform one step of a Newton iteration in order to minimize
414  round-off errors */
415  for (unsigned i = 0; i < nreal; i++) {
416  h = coeffs[1] + koreny[i] * (2 * coeffs[2] + 3 * koreny[i] * coeffs[3]);
417  if (h != 0.0)
418  koreny[i] -= (coeffs[0] + koreny[i] * (coeffs[1] + koreny[i] * (coeffs[2] + koreny[i] * coeffs[3]))) / h;
419  }
420  }
421 
422  else if (coeffs[2] != 0.0) {
423  /* quadratic problem? */
424  p = 0.5 * coeffs[1] / coeffs[2];
425  dis = sqr(p) - coeffs[0] / coeffs[2];
426  if (dis >= 0.0) {
427  /* two real solutions */
428  dis = sqrt(dis);
429  koreny[0] = -p - dis;
430  koreny[1] = -p + dis;
431  nreal = 2;
432  } else
433  /* no real solution */
434  nreal = 0;
435  }
436 
437  else if (coeffs[1] != 0.0) {
438  /* linear problem? */
439  koreny[0] = -coeffs[0] / coeffs[1];
440  nreal = 1;
441  }
442 
443  else
444  /* no equation */
445  nreal = 0;
446 
447  return nreal;
448 }
const double Pi
T w() const
void SWAP(double &realone, double &realtwo) const
T sqrt(T t)
Definition: SSEVec.h:19
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
double sqr(const double x) const
The Signals That Services Can Subscribe To This is based on ActivityRegistry h
Helper function to determine trigger accepts.
Definition: Activities.doc:4

◆ FindCoeff()

void TtFullLepKinSolver::FindCoeff ( const TLorentzVector &  al,
const TLorentzVector &  l,
const TLorentzVector &  b_al,
const TLorentzVector &  b_l,
const double  mt,
const double  mat,
const double  pxboost,
const double  pyboost,
double *  q_coeff 
)
private

Definition at line 157 of file TtFullLepKinSolver.cc.

References C, D, F, createfilelist::int, k16, k26, k36, k46, k51, k56, k61, cmsLHEtoEOSManager::l, m1, m2, m3, mb, TtSemiLepEvtBuilder_cfi::mt, mw, n1, n2, n3, pom, funct::pow(), and sqr().

Referenced by addKinSolInfo(), and getNuSolution().

165  {
166  double E, apom1, apom2, apom3;
167  double k11, k21, k31, k41, cpom1, cpom2, cpom3, l11, l21, l31, l41, l51, l61, k1, k2, k3, k4, k5, k6;
168  double l1, l2, l3, l4, l5, l6, k15, k25, k35, k45;
169 
170  C = -al.Px() - b_al.Px() - l.Px() - b_l.Px() + px_miss;
171  D = -al.Py() - b_al.Py() - l.Py() - b_l.Py() + py_miss;
172 
173  // right side of first two linear equations - missing pT
174 
175  E = (sqr(mt) - sqr(mw) - sqr(mb)) / (2 * b_al.E()) - sqr(mw) / (2 * al.E()) - al.E() +
176  al.Px() * b_al.Px() / b_al.E() + al.Py() * b_al.Py() / b_al.E() + al.Pz() * b_al.Pz() / b_al.E();
177  F = (sqr(mat) - sqr(mw) - sqr(mb)) / (2 * b_l.E()) - sqr(mw) / (2 * l.E()) - l.E() + l.Px() * b_l.Px() / b_l.E() +
178  l.Py() * b_l.Py() / b_l.E() + l.Pz() * b_l.Pz() / b_l.E();
179 
180  m1 = al.Px() / al.E() - b_al.Px() / b_al.E();
181  m2 = al.Py() / al.E() - b_al.Py() / b_al.E();
182  m3 = al.Pz() / al.E() - b_al.Pz() / b_al.E();
183 
184  n1 = l.Px() / l.E() - b_l.Px() / b_l.E();
185  n2 = l.Py() / l.E() - b_l.Py() / b_l.E();
186  n3 = l.Pz() / l.E() - b_l.Pz() / b_l.E();
187 
188  pom = E - m1 * C - m2 * D;
189  apom1 = sqr(al.Px()) - sqr(al.E());
190  apom2 = sqr(al.Py()) - sqr(al.E());
191  apom3 = sqr(al.Pz()) - sqr(al.E());
192 
193  k11 = 1 / sqr(al.E()) *
194  (pow(mw, 4) / 4 + sqr(C) * apom1 + sqr(D) * apom2 + apom3 * sqr(pom) / sqr(m3) +
195  sqr(mw) * (al.Px() * C + al.Py() * D + al.Pz() * pom / m3) + 2 * al.Px() * al.Py() * C * D +
196  2 * al.Px() * al.Pz() * C * pom / m3 + 2 * al.Py() * al.Pz() * D * pom / m3);
197  k21 = 1 / sqr(al.E()) *
198  (-2 * C * m3 * n3 * apom1 + 2 * apom3 * n3 * m1 * pom / m3 - sqr(mw) * m3 * n3 * al.Px() +
199  sqr(mw) * m1 * n3 * al.Pz() - 2 * al.Px() * al.Py() * D * m3 * n3 + 2 * al.Px() * al.Pz() * C * m1 * n3 -
200  2 * al.Px() * al.Pz() * n3 * pom + 2 * al.Py() * al.Pz() * D * m1 * n3);
201  k31 = 1 / sqr(al.E()) *
202  (-2 * D * m3 * n3 * apom2 + 2 * apom3 * n3 * m2 * pom / m3 - sqr(mw) * m3 * n3 * al.Py() +
203  sqr(mw) * m2 * n3 * al.Pz() - 2 * al.Px() * al.Py() * C * m3 * n3 + 2 * al.Px() * al.Pz() * C * m2 * n3 -
204  2 * al.Py() * al.Pz() * n3 * pom + 2 * al.Py() * al.Pz() * D * m2 * n3);
205  k41 = 1 / sqr(al.E()) *
206  (2 * apom3 * m1 * m2 * sqr(n3) + 2 * al.Px() * al.Py() * sqr(m3) * sqr(n3) -
207  2 * al.Px() * al.Pz() * m2 * m3 * sqr(n3) - 2 * al.Py() * al.Pz() * m1 * m3 * sqr(n3));
208  k51 = 1 / sqr(al.E()) *
209  (apom1 * sqr(m3) * sqr(n3) + apom3 * sqr(m1) * sqr(n3) - 2 * al.Px() * al.Pz() * m1 * m3 * sqr(n3));
210  k61 = 1 / sqr(al.E()) *
211  (apom2 * sqr(m3) * sqr(n3) + apom3 * sqr(m2) * sqr(n3) - 2 * al.Py() * al.Pz() * m2 * m3 * sqr(n3));
212 
213  cpom1 = sqr(l.Px()) - sqr(l.E());
214  cpom2 = sqr(l.Py()) - sqr(l.E());
215  cpom3 = sqr(l.Pz()) - sqr(l.E());
216 
217  l11 = 1 / sqr(l.E()) * (pow(mw, 4) / 4 + cpom3 * sqr(F) / sqr(n3) + sqr(mw) * l.Pz() * F / n3);
218  l21 =
219  1 / sqr(l.E()) *
220  (-2 * cpom3 * F * m3 * n1 / n3 + sqr(mw) * (l.Px() * m3 * n3 - l.Pz() * n1 * m3) + 2 * l.Px() * l.Pz() * F * m3);
221  l31 =
222  1 / sqr(l.E()) *
223  (-2 * cpom3 * F * m3 * n2 / n3 + sqr(mw) * (l.Py() * m3 * n3 - l.Pz() * n2 * m3) + 2 * l.Py() * l.Pz() * F * m3);
224  l41 = 1 / sqr(l.E()) *
225  (2 * cpom3 * n1 * n2 * sqr(m3) + 2 * l.Px() * l.Py() * sqr(m3) * sqr(n3) -
226  2 * l.Px() * l.Pz() * n2 * n3 * sqr(m3) - 2 * l.Py() * l.Pz() * n1 * n3 * sqr(m3));
227  l51 = 1 / sqr(l.E()) *
228  (cpom1 * sqr(m3) * sqr(n3) + cpom3 * sqr(n1) * sqr(m3) - 2 * l.Px() * l.Pz() * n1 * n3 * sqr(m3));
229  l61 = 1 / sqr(l.E()) *
230  (cpom2 * sqr(m3) * sqr(n3) + cpom3 * sqr(n2) * sqr(m3) - 2 * l.Py() * l.Pz() * n2 * n3 * sqr(m3));
231 
232  k1 = k11 * k61;
233  k2 = k61 * k21 / k51;
234  k3 = k31;
235  k4 = k41 / k51;
236  k5 = k61 / k51;
237  k6 = 1;
238 
239  l1 = l11 * k61;
240  l2 = l21 * k61 / k51;
241  l3 = l31;
242  l4 = l41 / k51;
243  l5 = l51 * k61 / (sqr(k51));
244  l6 = l61 / k61;
245 
246  k15 = k1 * l5 - l1 * k5;
247  k25 = k2 * l5 - l2 * k5;
248  k35 = k3 * l5 - l3 * k5;
249  k45 = k4 * l5 - l4 * k5;
250 
251  k16 = k1 * l6 - l1 * k6;
252  k26 = k2 * l6 - l2 * k6;
253  k36 = k3 * l6 - l3 * k6;
254  k46 = k4 * l6 - l4 * k6;
255  k56 = k5 * l6 - l5 * k6;
256 
257  koeficienty[0] = k15 * sqr(k36) - k35 * k36 * k16 - k56 * sqr(k16);
258  koeficienty[1] =
259  2 * k15 * k36 * k46 + k25 * sqr(k36) + k35 * (-k46 * k16 - k36 * k26) - k45 * k36 * k16 - 2 * k56 * k26 * k16;
260  koeficienty[2] = k15 * sqr(k46) + 2 * k25 * k36 * k46 + k35 * (-k46 * k26 - k36 * k56) -
261  k56 * (sqr(k26) + 2 * k56 * k16) - k45 * (k46 * k16 + k36 * k26);
262  koeficienty[3] = k25 * sqr(k46) - k35 * k46 * k56 - k45 * (k46 * k26 + k36 * k56) - 2 * sqr(k56) * k26;
263  koeficienty[4] = -k45 * k46 * k56 - pow(k56, 3);
264 
265  // normalization of coefficients
266  int moc = (int(log10(fabs(koeficienty[0]))) + int(log10(fabs(koeficienty[4])))) / 2;
267 
268  koeficienty[0] = koeficienty[0] / TMath::Power(10, moc);
269  koeficienty[1] = koeficienty[1] / TMath::Power(10, moc);
270  koeficienty[2] = koeficienty[2] / TMath::Power(10, moc);
271  koeficienty[3] = koeficienty[3] / TMath::Power(10, moc);
272  koeficienty[4] = koeficienty[4] / TMath::Power(10, moc);
273 }
double sqr(const double x) const
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:29

◆ getNuSolution()

TtFullLepKinSolver::NeutrinoSolution TtFullLepKinSolver::getNuSolution ( const TLorentzVector &  LV_l,
const TLorentzVector &  LV_l_,
const TLorentzVector &  LV_b,
const TLorentzVector &  LV_b_ 
)

Definition at line 125 of file TtFullLepKinSolver.cc.

References FindCoeff(), LV_n, LV_n_, TtSemiLepEvtBuilder_cfi::mt, TtFullLepKinSolver::NeutrinoSolution::neutrino, TtFullLepKinSolver::NeutrinoSolution::neutrinoBar, pxmiss_, pymiss_, quartic(), topmass_begin, topmass_end, topmass_step, TopRec(), TtFullLepKinSolver::NeutrinoSolution::weight, mps_merge::weight, and WeightSolfromShape().

Referenced by TtFullLepKinSolutionProducer::produce().

128  {
129  math::XYZTLorentzVector maxLV_n = math::XYZTLorentzVector(0, 0, 0, 0);
130  math::XYZTLorentzVector maxLV_n_ = math::XYZTLorentzVector(0, 0, 0, 0);
131 
132  //loop on top mass parameter
133  double weightmax = -1;
134  for (double mt = topmass_begin; mt < topmass_end + 0.5 * topmass_step; mt += topmass_step) {
135  double q_coeff[5], q_sol[4];
136  FindCoeff(LV_l, LV_l_, LV_b, LV_b_, mt, mt, pxmiss_, pymiss_, q_coeff);
137  int NSol = quartic(q_coeff, q_sol);
138 
139  //loop on all solutions
140  for (int isol = 0; isol < NSol; isol++) {
141  TopRec(LV_l, LV_l_, LV_b, LV_b_, q_sol[isol]);
142  double weight = WeightSolfromShape();
143  if (weight > weightmax) {
144  weightmax = weight;
145  maxLV_n.SetPxPyPzE(LV_n.Px(), LV_n.Py(), LV_n.Pz(), LV_n.E());
146  maxLV_n_.SetPxPyPzE(LV_n_.Px(), LV_n_.Py(), LV_n_.Pz(), LV_n_.E());
147  }
148  }
149  }
151  nuSol.neutrino = reco::LeafCandidate(0, maxLV_n);
152  nuSol.neutrinoBar = reco::LeafCandidate(0, maxLV_n_);
153  nuSol.weight = weightmax;
154  return nuSol;
155 }
const double topmass_end
void FindCoeff(const TLorentzVector &al, const TLorentzVector &l, const TLorentzVector &b_al, const TLorentzVector &b_l, const double mt, const double mat, const double pxboost, const double pyboost, double *q_coeff)
Definition: weight.py:1
double WeightSolfromShape() const
use the parametrized event shape to obtain the solution weight.
const double topmass_step
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
int quartic(double *q_coeff, double *q_sol) const
void TopRec(const TLorentzVector &al, const TLorentzVector &l, const TLorentzVector &b_al, const TLorentzVector &b_l, const double sol)
const double topmass_begin

◆ quartic()

int TtFullLepKinSolver::quartic ( double *  q_coeff,
double *  q_sol 
) const
private

Definition at line 314 of file TtFullLepKinSolver.cc.

References b0, b1, b2, c, cubic(), d0, d1, h, mps_fire::i, multPhiCorr_741_25nsDY_cfi::px, sqr(), mathSSE::sqrt(), submitPVValidationJobs::t, w(), and z.

Referenced by addKinSolInfo(), and getNuSolution().

314  {
315  double w, b0, b1, b2;
316  double c[4];
317  double d0, d1, h, t, z;
318  double* px;
319 
320  if (koeficienty[4] == 0.0)
321  return cubic(koeficienty, koreny);
322  /* quartic problem? */
323  w = koeficienty[3] / (4 * koeficienty[4]);
324  /* offset */
325  b2 = -6 * sqr(w) + koeficienty[2] / koeficienty[4];
326  /* koeficienty. of shifted polynomial */
327  b1 = (8 * sqr(w) - 2 * koeficienty[2] / koeficienty[4]) * w + koeficienty[1] / koeficienty[4];
328  b0 = ((-3 * sqr(w) + koeficienty[2] / koeficienty[4]) * w - koeficienty[1] / koeficienty[4]) * w +
329  koeficienty[0] / koeficienty[4];
330 
331  c[3] = 1.0;
332  /* cubic resolvent */
333  c[2] = b2;
334  c[1] = -4 * b0;
335  c[0] = sqr(b1) - 4 * b0 * b2;
336 
337  cubic(c, koreny);
338  z = koreny[0];
339  //double z1=1.0,z2=2.0,z3=3.0;
340  //TMath::RootsCubic(c,z1,z2,z3);
341  //if (z2 !=0) z = z2;
342  //if (z1 !=0) z = z1;
343  /* only lowermost root needed */
344 
345  int nreal = 0;
346  px = koreny;
347  t = sqrt(0.25 * sqr(z) - b0);
348  for (int i = -1; i <= 1; i += 2) {
349  d0 = -0.5 * z + i * t;
350  /* coeffs. of quadratic factor */
351  d1 = (t != 0.0) ? -i * 0.5 * b1 / t : i * sqrt(-z - b2);
352  h = 0.25 * sqr(d1) - d0;
353  if (h >= 0.0) {
354  h = sqrt(h);
355  nreal += 2;
356  *px++ = -0.5 * d1 - h - w;
357  *px++ = -0.5 * d1 + h - w;
358  }
359  }
360 
361  // if (nreal==4) {
362  /* sort results */
363  // if (koreny[2]<koreny[0]) SWAP(koreny[0], koreny[2]);
364  // if (koreny[3]<koreny[1]) SWAP(koreny[1], koreny[3]);
365  // if (koreny[1]<koreny[0]) SWAP(koreny[0], koreny[1]);
366  // if (koreny[3]<koreny[2]) SWAP(koreny[2], koreny[3]);
367  // if (koreny[2]<koreny[1]) SWAP(koreny[1], koreny[2]);
368  // }
369  return nreal;
370 }
T w() const
int cubic(const double *c_coeff, double *c_sol) const
T sqrt(T t)
Definition: SSEVec.h:19
static constexpr float d0
double sqr(const double x) const
static constexpr float b2
static constexpr float b0
static constexpr float d1
The Signals That Services Can Subscribe To This is based on ActivityRegistry h
Helper function to determine trigger accepts.
Definition: Activities.doc:4
static constexpr float b1

◆ SetConstraints()

void TtFullLepKinSolver::SetConstraints ( const double  xx = 0,
const double  yy = 0 
)

◆ sqr()

double TtFullLepKinSolver::sqr ( const double  x) const
inlineprivate

Definition at line 80 of file TtFullLepKinSolver.h.

References x.

Referenced by cubic(), FindCoeff(), and quartic().

80 { return (x * x); }

◆ SWAP()

void TtFullLepKinSolver::SWAP ( double &  realone,
double &  realtwo 
) const
private

Definition at line 450 of file TtFullLepKinSolver.cc.

References printConversionInfo::aux.

Referenced by cubic().

450  {
451  if (realtwo < realone) {
452  double aux = realtwo;
453  realtwo = realone;
454  realone = aux;
455  }
456 }

◆ TopRec()

void TtFullLepKinSolver::TopRec ( const TLorentzVector &  al,
const TLorentzVector &  l,
const TLorentzVector &  b_al,
const TLorentzVector &  b_l,
const double  sol 
)
private

Definition at line 275 of file TtFullLepKinSolver.cc.

References printConversionInfo::aux, C, D, F, k16, k26, k36, k46, k51, k56, k61, cmsLHEtoEOSManager::l, LV_n, LV_n_, LV_t, LV_t_, LV_tt_t, LV_tt_t_, m1, m2, m3, n1, n2, n3, pom, funct::pow(), and mkfit::Const::sol.

Referenced by addKinSolInfo(), and getNuSolution().

279  {
280  TVector3 t_ttboost;
281  TLorentzVector aux;
282  double pxp, pyp, pzp, pup, pvp, pwp;
283 
284  pxp = sol * (m3 * n3 / k51);
285  pyp = -(m3 * n3 / k61) * (k56 * pow(sol, 2) + k26 * sol + k16) / (k36 + k46 * sol);
286  pzp = -1 / n3 * (n1 * pxp + n2 * pyp - F);
287  pwp = 1 / m3 * (m1 * pxp + m2 * pyp + pom);
288  pup = C - pxp;
289  pvp = D - pyp;
290 
291  LV_n_.SetXYZM(pxp, pyp, pzp, 0.0);
292  LV_n.SetXYZM(pup, pvp, pwp, 0.0);
293 
294  LV_t_ = b_l + l + LV_n_;
295  LV_t = b_al + al + LV_n;
296 
297  aux = (LV_t_ + LV_t);
298  t_ttboost = -aux.BoostVector();
299  LV_tt_t_ = LV_t_;
300  LV_tt_t = LV_t;
301  LV_tt_t_.Boost(t_ttboost);
302  LV_tt_t.Boost(t_ttboost);
303 }
constexpr float sol
Definition: Config.h:48
TLorentzVector LV_tt_t
TLorentzVector LV_tt_t_
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:29

◆ useWeightFromMC()

void TtFullLepKinSolver::useWeightFromMC ( bool  useMC)
inline

Definition at line 43 of file TtFullLepKinSolver.h.

References PatBasicFWLiteJetAnalyzer_Selector_cfg::useMC, and useMCforBest_.

Referenced by TtDilepEvtSolutionMaker::produce().

43 { useMCforBest_ = useMC; }
bool useMCforBest_
flag to swith from WeightSolfromMC() to WeightSolfromShape()

◆ WeightSolfromMC()

double TtFullLepKinSolver::WeightSolfromMC ( ) const
private

Definition at line 305 of file TtFullLepKinSolver.cc.

References genLV_n, genLV_n_, LV_n, LV_n_, and mps_merge::weight.

Referenced by addKinSolInfo().

305  {
306  double weight = 1;
307  weight = ((LV_n.E() > genLV_n.E()) ? genLV_n.E() / LV_n.E() : LV_n.E() / genLV_n.E()) *
308  ((LV_n_.E() > genLV_n_.E()) ? genLV_n_.E() / LV_n_.E() : LV_n_.E() / genLV_n_.E());
309  return weight;
310 }
Definition: weight.py:1
TLorentzVector genLV_n
provisional
TLorentzVector genLV_n_

◆ WeightSolfromShape()

double TtFullLepKinSolver::WeightSolfromShape ( ) const
private

use the parametrized event shape to obtain the solution weight.

Definition at line 312 of file TtFullLepKinSolver.cc.

References EventShape_, LV_n, and LV_n_.

Referenced by addKinSolInfo(), and getNuSolution().

312 { return EventShape_->Eval(LV_n.E(), LV_n_.E()); }
TF2 * EventShape_
Event shape.

Member Data Documentation

◆ C

double TtFullLepKinSolver::C
private

Definition at line 98 of file TtFullLepKinSolver.h.

Referenced by FindCoeff(), and TopRec().

◆ D

double TtFullLepKinSolver::D
private

Definition at line 99 of file TtFullLepKinSolver.h.

Referenced by FindCoeff(), and TopRec().

◆ EventShape_

TF2* TtFullLepKinSolver::EventShape_
private

Event shape.

Definition at line 124 of file TtFullLepKinSolver.h.

Referenced by TtFullLepKinSolver(), WeightSolfromShape(), and ~TtFullLepKinSolver().

◆ F

double TtFullLepKinSolver::F
private

Definition at line 100 of file TtFullLepKinSolver.h.

Referenced by FindCoeff(), and TopRec().

◆ genLV_n

TLorentzVector TtFullLepKinSolver::genLV_n
private

provisional

Definition at line 119 of file TtFullLepKinSolver.h.

Referenced by addKinSolInfo(), and WeightSolfromMC().

◆ genLV_n_

TLorentzVector TtFullLepKinSolver::genLV_n_
private

Definition at line 119 of file TtFullLepKinSolver.h.

Referenced by addKinSolInfo(), and WeightSolfromMC().

◆ k16

double TtFullLepKinSolver::k16
private

Definition at line 102 of file TtFullLepKinSolver.h.

Referenced by FindCoeff(), and TopRec().

◆ k26

double TtFullLepKinSolver::k26
private

Definition at line 103 of file TtFullLepKinSolver.h.

Referenced by FindCoeff(), and TopRec().

◆ k36

double TtFullLepKinSolver::k36
private

Definition at line 104 of file TtFullLepKinSolver.h.

Referenced by FindCoeff(), and TopRec().

◆ k46

double TtFullLepKinSolver::k46
private

Definition at line 105 of file TtFullLepKinSolver.h.

Referenced by FindCoeff(), and TopRec().

◆ k51

double TtFullLepKinSolver::k51
private

Definition at line 107 of file TtFullLepKinSolver.h.

Referenced by FindCoeff(), and TopRec().

◆ k56

double TtFullLepKinSolver::k56
private

Definition at line 106 of file TtFullLepKinSolver.h.

Referenced by FindCoeff(), and TopRec().

◆ k61

double TtFullLepKinSolver::k61
private

Definition at line 108 of file TtFullLepKinSolver.h.

Referenced by FindCoeff(), and TopRec().

◆ LV_n

TLorentzVector TtFullLepKinSolver::LV_n
private

Definition at line 117 of file TtFullLepKinSolver.h.

Referenced by getNuSolution(), TopRec(), WeightSolfromMC(), and WeightSolfromShape().

◆ LV_n_

TLorentzVector TtFullLepKinSolver::LV_n_
private

Definition at line 117 of file TtFullLepKinSolver.h.

Referenced by getNuSolution(), TopRec(), WeightSolfromMC(), and WeightSolfromShape().

◆ LV_t

TLorentzVector TtFullLepKinSolver::LV_t
private

Definition at line 117 of file TtFullLepKinSolver.h.

Referenced by TopRec().

◆ LV_t_

TLorentzVector TtFullLepKinSolver::LV_t_
private

Definition at line 117 of file TtFullLepKinSolver.h.

Referenced by TopRec().

◆ LV_tt_t

TLorentzVector TtFullLepKinSolver::LV_tt_t
private

Definition at line 117 of file TtFullLepKinSolver.h.

Referenced by TopRec().

◆ LV_tt_t_

TLorentzVector TtFullLepKinSolver::LV_tt_t_
private

Definition at line 117 of file TtFullLepKinSolver.h.

Referenced by TopRec().

◆ m1

double TtFullLepKinSolver::m1
private

Definition at line 109 of file TtFullLepKinSolver.h.

Referenced by FindCoeff(), and TopRec().

◆ m2

double TtFullLepKinSolver::m2
private

Definition at line 110 of file TtFullLepKinSolver.h.

Referenced by FindCoeff(), and TopRec().

◆ m3

double TtFullLepKinSolver::m3
private

Definition at line 111 of file TtFullLepKinSolver.h.

Referenced by FindCoeff(), and TopRec().

◆ mb

const double TtFullLepKinSolver::mb
private

Definition at line 94 of file TtFullLepKinSolver.h.

Referenced by FindCoeff().

◆ mw

const double TtFullLepKinSolver::mw
private

Definition at line 92 of file TtFullLepKinSolver.h.

Referenced by FindCoeff().

◆ n1

double TtFullLepKinSolver::n1
private

Definition at line 112 of file TtFullLepKinSolver.h.

Referenced by FindCoeff(), and TopRec().

◆ n2

double TtFullLepKinSolver::n2
private

Definition at line 113 of file TtFullLepKinSolver.h.

Referenced by FindCoeff(), and TopRec().

◆ n3

double TtFullLepKinSolver::n3
private

Definition at line 114 of file TtFullLepKinSolver.h.

Referenced by FindCoeff(), and TopRec().

◆ pom

double TtFullLepKinSolver::pom
private

Definition at line 101 of file TtFullLepKinSolver.h.

Referenced by FindCoeff(), and TopRec().

◆ pxmiss_

double TtFullLepKinSolver::pxmiss_
private

Definition at line 96 of file TtFullLepKinSolver.h.

Referenced by addKinSolInfo(), getNuSolution(), and SetConstraints().

◆ pymiss_

double TtFullLepKinSolver::pymiss_
private

Definition at line 96 of file TtFullLepKinSolver.h.

Referenced by addKinSolInfo(), getNuSolution(), and SetConstraints().

◆ topmass_begin

const double TtFullLepKinSolver::topmass_begin
private

Definition at line 86 of file TtFullLepKinSolver.h.

Referenced by addKinSolInfo(), and getNuSolution().

◆ topmass_end

const double TtFullLepKinSolver::topmass_end
private

Definition at line 88 of file TtFullLepKinSolver.h.

Referenced by addKinSolInfo(), and getNuSolution().

◆ topmass_step

const double TtFullLepKinSolver::topmass_step
private

Definition at line 90 of file TtFullLepKinSolver.h.

Referenced by addKinSolInfo(), and getNuSolution().

◆ useMCforBest_

bool TtFullLepKinSolver::useMCforBest_
private

flag to swith from WeightSolfromMC() to WeightSolfromShape()

Definition at line 122 of file TtFullLepKinSolver.h.

Referenced by addKinSolInfo(), and useWeightFromMC().