CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Types | Public Member Functions | Static Public Member Functions
tauImpactParameter::MultiProngTauSolver Class Reference

#include <MultiProngTauSolver.h>

Inheritance diagram for tauImpactParameter::MultiProngTauSolver:
tauImpactParameter::TauA1NuConstrainedFitter

Public Types

enum  Ambiguity { zero, minus, plus, NAmbiguity }
 

Public Member Functions

 MultiProngTauSolver ()
 
virtual ~MultiProngTauSolver ()
 

Static Public Member Functions

static void analyticESolver (TLorentzVector &nu_plus, TLorentzVector &nu_minus, const TLorentzVector &A1, bool &isReal)
 
static LorentzVectorParticle estimateNu (const LorentzVectorParticle &a1, const TVector3 &pv, int ambiguity, TLorentzVector &tau)
 
static void numericalESolver (TLorentzVector &nu_plus, TLorentzVector &nu_minus, const TLorentzVector &A1, bool &isReal)
 
static void quadratic (double &x_plus, double &x_minus, double a, double b, double c, bool &isReal)
 
static TVectorT< double > rotateToTauFrame (const TVectorT< double > &inpar)
 
static bool setTauDirectionatThetaGJMax (const TLorentzVector &a1, double &theta, double &phi, double scale=1.0)
 
static void solveByRotation (const TVector3 &TauDir, const TLorentzVector &A1, TLorentzVector &Tau_plus, TLorentzVector &Tau_minus, TLorentzVector &nu_plus, TLorentzVector &nu_minus, bool &isReal, bool rotateback=true)
 
static double thetaGJMax (const TLorentzVector &a1)
 

Detailed Description

Definition at line 17 of file MultiProngTauSolver.h.

Member Enumeration Documentation

Constructor & Destructor Documentation

tauImpactParameter::MultiProngTauSolver::MultiProngTauSolver ( )
inline

Definition at line 22 of file MultiProngTauSolver.h.

22 {};
virtual tauImpactParameter::MultiProngTauSolver::~MultiProngTauSolver ( )
inlinevirtual

Definition at line 23 of file MultiProngTauSolver.h.

23 {};

Member Function Documentation

void MultiProngTauSolver::analyticESolver ( TLorentzVector &  nu_plus,
TLorentzVector &  nu_minus,
const TLorentzVector &  A1,
bool &  isReal 
)
static

Definition at line 16 of file MultiProngTauSolver.cc.

References a, b, EnergyCorrector::c, quadratic(), mathSSE::sqrt(), and tauImpactParameter::PDGInfo::tau_mass().

Referenced by solveByRotation().

16  {
17  double a=(A1.Pz()*A1.Pz())/(A1.E()*A1.E())-1.0;
18  double K=(PDGInfo::tau_mass()*PDGInfo::tau_mass()-A1.M2()-2.0*A1.Pt()*A1.Pt())/(2.0*A1.E());
19  double b=2.0*K*A1.Pz()/A1.E();
20  double c=K*K-A1.Pt()*A1.Pt();
21  double z_plus(0),z_minus(0);
22  quadratic(z_plus,z_minus,a,b,c,isReal);
23  nu_plus=TLorentzVector(-A1.Px(),-A1.Py(),z_plus,sqrt(z_plus*z_plus+A1.Pt()*A1.Pt()));
24  nu_minus=TLorentzVector(-A1.Px(),-A1.Py(),z_minus,sqrt(z_minus*z_minus+A1.Pt()*A1.Pt()));
25 }
T sqrt(T t)
Definition: SSEVec.h:18
double b
Definition: hdecay.h:120
double a
Definition: hdecay.h:121
static double tau_mass()
Definition: PDGInfo.h:14
static void quadratic(double &x_plus, double &x_minus, double a, double b, double c, bool &isReal)
LorentzVectorParticle MultiProngTauSolver::estimateNu ( const LorentzVectorParticle a1,
const TVector3 &  pv,
int  ambiguity,
TLorentzVector &  tau 
)
static

Definition at line 101 of file MultiProngTauSolver.cc.

References tauImpactParameter::Particle::bField(), funct::cos(), i, j, tauImpactParameter::LorentzVectorParticle::m, minus, tauImpactParameter::LorentzVectorParticle::NLorentzandVertexPar, tauImpactParameter::PDGInfo::nu_tau, tauImpactParameter::LorentzVectorParticle::NVertex, tauImpactParameter::LorentzVectorParticle::p4(), tauImpactParameter::LorentzVectorParticle::parameter(), phi(), plus, MetAnalyzer::pv(), tauImpactParameter::LorentzVectorParticle::px, tauImpactParameter::LorentzVectorParticle::py, tauImpactParameter::LorentzVectorParticle::pz, setTauDirectionatThetaGJMax(), funct::sin(), solveByRotation(), theta(), findQualityFiles::v, tauImpactParameter::LorentzVectorParticle::vertex(), tauImpactParameter::LorentzVectorParticle::vertexCov(), tauImpactParameter::LorentzVectorParticle::vx, tauImpactParameter::LorentzVectorParticle::vy, tauImpactParameter::LorentzVectorParticle::vz, and zero.

101  {
102  TLorentzVector lorentzA1=a1.p4();
103  TVector3 sv=a1.vertex();
104  TVector3 tauFlghtDir=sv-pv;
105  TLorentzVector nuGuess;
106 
107  TVector3 startingtauFlghtDir=tauFlghtDir.Unit();
108  if(ambiguity==zero){
109  double theta=tauFlghtDir.Theta();
110  double phi=tauFlghtDir.Phi();
111  setTauDirectionatThetaGJMax(lorentzA1,theta,phi);
112  startingtauFlghtDir=TVector3(sin(theta)*cos(phi),sin(theta)*sin(phi),cos(theta));
113  }
114  TLorentzVector tau1,tau2,nu1,nu2;
115  bool isReal;
116  solveByRotation(startingtauFlghtDir,lorentzA1,tau1,tau2,nu1,nu2,isReal);
117  if(ambiguity==plus){ nuGuess=nu1; tau=tau1; }
118  if(ambiguity==minus){ nuGuess=nu2; tau=tau1; }
119  if(ambiguity==zero){ nuGuess=nu1; tau=tau1; }
120  TVectorT<double> par(LorentzVectorParticle::NLorentzandVertexPar);
124  par(LorentzVectorParticle::px)=nuGuess.Px();
125  par(LorentzVectorParticle::py)=nuGuess.Py();
126  par(LorentzVectorParticle::pz)=nuGuess.Pz();
127  par(LorentzVectorParticle::m) =nuGuess.M();
128  TMatrixTSym<double> Cov(LorentzVectorParticle::NLorentzandVertexPar);
129  TMatrixTSym<double> pvCov=a1.vertexCov();
131  for(int j=0; j<=i; j++){
132  if(i<LorentzVectorParticle::NVertex) Cov(i,j)=pvCov(i,j);
133  else Cov(i,j)=0;
134  }
135  double v=0;
137  if(v<1000.0) v=1000.0; // try lowing to test impact
138  Cov(i,i)+=v;
139  }
140  return LorentzVectorParticle(par,Cov,PDGInfo::nu_tau,0,a1.bField());
141 }
TMatrixTSym< double > vertexCov() const
int i
Definition: DBlmapReader.cc:9
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
Geom::Theta< T > theta() const
static bool setTauDirectionatThetaGJMax(const TLorentzVector &a1, double &theta, double &phi, double scale=1.0)
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
int j
Definition: DBlmapReader.cc:9
virtual double bField() const
Definition: Particle.h:26
static void solveByRotation(const TVector3 &TauDir, const TLorentzVector &A1, TLorentzVector &Tau_plus, TLorentzVector &Tau_minus, TLorentzVector &nu_plus, TLorentzVector &nu_minus, bool &isReal, bool rotateback=true)
void MultiProngTauSolver::numericalESolver ( TLorentzVector &  nu_plus,
TLorentzVector &  nu_minus,
const TLorentzVector &  A1,
bool &  isReal 
)
static

Definition at line 27 of file MultiProngTauSolver.cc.

References i, min(), mathSSE::sqrt(), metsig::tau, tauImpactParameter::PDGInfo::tau_mass(), and SiStripMonitorClusterAlca_cfi::zmin.

27  {
28  double rmin(-100), rmax(100), step(0.01), mtau2(PDGInfo::tau_mass()*PDGInfo::tau_mass()), z1(-999), z2(-999), zmin(-999), min(9999), prev(9999);
29  double z=rmin;
30  TLorentzVector nu,tau;
31  for(int i=0;i<=(int)(rmax-rmin)/step;i++){
32  nu.SetPxPyPzE(-A1.Px(),-A1.Py(),z,sqrt(z*z+A1.Pt()*A1.Pt()));
33  tau=A1+nu;
34  double m2=tau.M2();
35  if(m2-mtau2<0 && prev-mtau2>=0) z1=z;
36  if(m2-mtau2>0 && prev-mtau2<=0) z2=z;
37  if(min>m2){ zmin=z; min=m2;}
38  prev=m2;
39  z+=step;
40  }
41  if(z1!=-999 && z2!=-999){
42  nu_plus=TLorentzVector(-A1.Px(),-A1.Py(),z1,sqrt(z1*z1+A1.Pt()*A1.Pt()));
43  nu_minus=TLorentzVector(-A1.Px(),-A1.Py(),z2,sqrt(z2*z2+A1.Pt()*A1.Pt()));
44  }
45  else{
46  nu_plus=TLorentzVector(-A1.Px(),-A1.Py(),zmin,sqrt(zmin*zmin+A1.Pt()*A1.Pt()));
47  nu_minus=TLorentzVector(-A1.Px(),-A1.Py(),zmin,sqrt(zmin*zmin+A1.Pt()*A1.Pt()));
48  }
49 }
int i
Definition: DBlmapReader.cc:9
T sqrt(T t)
Definition: SSEVec.h:18
T min(T a, T b)
Definition: MathUtil.h:58
step
static double tau_mass()
Definition: PDGInfo.h:14
void MultiProngTauSolver::quadratic ( double &  x_plus,
double &  x_minus,
double  a,
double  b,
double  c,
bool &  isReal 
)
static

Definition at line 8 of file MultiProngTauSolver.cc.

References EnergyCorrector::c, dttmaxenums::R, and mathSSE::sqrt().

Referenced by analyticESolver().

8  {
9  double R=b*b-4*a*c;
10  isReal=true;
11  if(R<0){isReal=false;}// flag cases when R<0 but compute quadratic equation with |R|
12  x_minus=(-b+sqrt(fabs(R)))/(2.0*a); // opposite sign is smaller
13  x_plus=(-b-sqrt(fabs(R)))/(2.0*a);
14 }
T sqrt(T t)
Definition: SSEVec.h:18
double b
Definition: hdecay.h:120
double a
Definition: hdecay.h:121
TVectorT< double > MultiProngTauSolver::rotateToTauFrame ( const TVectorT< double > &  inpar)
static

Definition at line 143 of file MultiProngTauSolver.cc.

References funct::cos(), and funct::sin().

143  {
144  TVectorT<double> outpar(3);
145  TVector3 res(inpar(0),inpar(1),inpar(2));
146  TVector3 Uz(sin(inpar(4))*cos(inpar(3)),sin(inpar(4))*sin(inpar(3)),cos(inpar(4)));
147  res.RotateUz(Uz);
148  /* double phi=inpar(3,0);
149  double theta=inpar(4,0);
150  res.RotateZ(-phi);
151  TVector3 Y(0,1,0);
152  TVector3 thetadir=res.Cross(Y);
153  thetadir.RotateY(-theta);
154  res.RotateY(-theta);
155  res.RotateZ(thetadir.Phi());*/
156  outpar(0)=res.X();
157  outpar(1)=res.Y();
158  outpar(2)=res.Z();
159  return outpar;
160 }
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
bool MultiProngTauSolver::setTauDirectionatThetaGJMax ( const TLorentzVector &  a1,
double &  theta,
double &  phi,
double  scale = 1.0 
)
static

Definition at line 76 of file MultiProngTauSolver.cc.

References a, b, funct::cos(), funct::sin(), metsig::tau, theta(), and thetaGJMax().

Referenced by estimateNu().

76  {
77  double thetaGJMaxvar =thetaGJMax(a1);
78  TVector3 a1v(a1.Vect()); if(a1v.Mag()!=0) a1v*=1/a1v.Mag();
79  TVector3 tau(cos(phi)*sin(theta),sin(phi)*sin(theta),cos(theta));
80  double dphitheta=acos(a1v.Dot(tau)/(a1v.Mag()*tau.Mag()));
81  if(thetaGJMaxvar<dphitheta || scale<0){
82  if(scale<0) scale=1.0;
83  double a=(thetaGJMaxvar/dphitheta)-(1-scale);
84  double b=1-(thetaGJMaxvar/dphitheta)+(1-scale);
85  edm::LogInfo("RecoTauTag/ImpactParameter") << "SetTauDirectionatThetaGJMax before GF " << thetaGJMaxvar << " dot " << acos(a1v.Dot(tau)/(a1v.Mag()*tau.Mag())) << " a1 phi " << a1v.Phi() << " tau phi " << tau.Phi() << " a1 theta " <<a1v.Theta() << " tau theta " << tau.Theta() ;
86  tau*=a;
87  a1v*=b;
88  tau+=a1v;
89  theta=tau.Theta();
90  phi=tau.Phi();
91  edm::LogInfo("RecoTauTag/ImpactParameter") << "SetTauDirectionatThetaGJMax GF " << thetaGJMaxvar << " dot " << acos(a1v.Dot(tau)/(a1v.Mag()*tau.Mag())) << " phi " << phi << " theta " << theta ;
92  return true;
93  }
94  return false;
95 }
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
Geom::Theta< T > theta() const
static double thetaGJMax(const TLorentzVector &a1)
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
double b
Definition: hdecay.h:120
double a
Definition: hdecay.h:121
void MultiProngTauSolver::solveByRotation ( const TVector3 &  TauDir,
const TLorentzVector &  A1,
TLorentzVector &  Tau_plus,
TLorentzVector &  Tau_minus,
TLorentzVector &  nu_plus,
TLorentzVector &  nu_minus,
bool &  isReal,
bool  rotateback = true 
)
static

Definition at line 51 of file MultiProngTauSolver.cc.

References analyticESolver(), phi(), and theta().

Referenced by estimateNu(), tauImpactParameter::TauA1NuConstrainedFitter::fit(), and tauImpactParameter::TauA1NuConstrainedFitter::SolveAmbiguityAnalytically().

52  {
53  TLorentzVector A1rot=A1;
54  double phi(TauDir.Phi()),theta(TauDir.Theta());
55  A1rot.RotateZ(-phi);
56  A1rot.RotateY(-theta);
58  // NumericalESolver(nu_plus,nu_minus,A1rot); // for debugging analyticESolver (slow)
59  analyticESolver(nu_plus,nu_minus,A1rot,isReal);
61  if(rotateback){
62  nu_plus.RotateY(theta);
63  nu_plus.RotateZ(phi);
64  Tau_plus=A1+nu_plus;
65  //
66  nu_minus.RotateY(theta);
67  nu_minus.RotateZ(phi);
68  Tau_minus=A1+nu_minus;
69  }
70  else{
71  Tau_plus=A1rot+nu_plus;
72  Tau_minus=A1rot+nu_minus;
73  }
74 }
static void analyticESolver(TLorentzVector &nu_plus, TLorentzVector &nu_minus, const TLorentzVector &A1, bool &isReal)
Geom::Theta< T > theta() const
double MultiProngTauSolver::thetaGJMax ( const TLorentzVector &  a1)
static