#include <RecoVertex/KinematicFit/interface/SimplePointingConstraint.h>
Public Member Functions | |
virtual SimplePointingConstraint * | clone () const |
Clone method. | |
virtual pair< AlgebraicMatrix, AlgebraicVector > | derivative (const vector< RefCountedKinematicParticle > par) const |
Vector of values and matrix of derivatives calculated using current state parameters as expansion point. | |
virtual pair< AlgebraicMatrix, AlgebraicVector > | derivative (const AlgebraicVector &exPoint) const |
virtual AlgebraicVector | deviations (int nStates) const |
Returns vector of sigma squared associated to the KinematicParameters of refitted particles Initial deviations are given by user for the constraining parameters (mass, momentum components etc). | |
virtual int | numberOfEquations () const |
Returns number of constraint equations used for fitting. | |
SimplePointingConstraint (const GlobalPoint &ref) | |
virtual pair< AlgebraicVector, AlgebraicVector > | value (const vector< RefCountedKinematicParticle > par) const |
Methods making value and derivative matrix using current state parameters as expansion 7-point. | |
virtual pair< AlgebraicVector, AlgebraicVector > | value (const AlgebraicVector &exPoint) const |
Vector of values and matrix of derivatives calculated at given expansion 7xNumberOfStates point. | |
Private Member Functions | |
pair< AlgebraicMatrix, AlgebraicVector > | makeDerivative (const AlgebraicVector &exPoint) const |
pair< AlgebraicVector, AlgebraicVector > | makeValue (const AlgebraicVector &exPoint) const |
Private Attributes | |
GlobalPoint | refPoint |
Example: if b-meson momentum is reconstructed at b-meson decay position (secondary vertex), making reconstructed momentum pointing the the primary vertex
Multiple track refit is not supported in current version
Kirill Prokofiev, March 2004 MultiState version: July 2004
Definition at line 21 of file SimplePointingConstraint.h.
SimplePointingConstraint::SimplePointingConstraint | ( | const GlobalPoint & | ref | ) | [inline] |
Definition at line 25 of file SimplePointingConstraint.h.
Referenced by clone().
00025 :refPoint(ref) 00026 {}
virtual SimplePointingConstraint* SimplePointingConstraint::clone | ( | ) | const [inline, virtual] |
Clone method.
Implements KinematicConstraint.
Definition at line 53 of file SimplePointingConstraint.h.
References SimplePointingConstraint().
00054 {return new SimplePointingConstraint(*this);}
pair< AlgebraicMatrix, AlgebraicVector > SimplePointingConstraint::derivative | ( | const vector< RefCountedKinematicParticle > | par | ) | const [virtual] |
Vector of values and matrix of derivatives calculated using current state parameters as expansion point.
Implements KinematicConstraint.
Definition at line 43 of file SimplePointingConstraint.cc.
References makeDerivative().
00044 { 00045 int nStates = par.size(); 00046 if(nStates == 0) throw VertexException("PointingKinematicConstraint::Empty vector of particles passed"); 00047 if(nStates != 1) throw VertexException("PointingKinematicConstraint::Current version does not support the multistate refit"); 00048 00049 AlgebraicMatrix dr(2,7,0); 00050 AlgebraicVector lPoint = asHepVector<7>(par.front()->currentState().kinematicParameters().vector()); 00051 00052 //2x7 derivative matrix for given state 00053 AlgebraicMatrix lDeriv = makeDerivative(lPoint).first; 00054 dr.sub(1,1,lDeriv); 00055 // cout<<"Derivative returned: "<<dr<<endl; 00056 // cout<<"For the value: "<<lPoint<<endl; 00057 return pair<AlgebraicMatrix,AlgebraicVector>(dr,lPoint); 00058 }
pair< AlgebraicMatrix, AlgebraicVector > SimplePointingConstraint::derivative | ( | const AlgebraicVector & | exPoint | ) | const [virtual] |
Implements KinematicConstraint.
Definition at line 25 of file SimplePointingConstraint.cc.
References makeDerivative().
00026 { 00027 if(exPoint.num_row() ==0 ) throw VertexException("PointingKinematicConstraint::value requested for zero Linearization point"); 00028 00029 //security check for extended cartesian parametrization 00030 int inSize = exPoint.num_row(); 00031 if((inSize%7) !=0) throw VertexException("PointingKinematicConstraint::linearization point has a wrong dimension"); 00032 int nStates = inSize/7; 00033 if(nStates != 1) throw VertexException("PointingKinematicConstraint::Current version does not support the multistate refit"); 00034 AlgebraicVector lPar = exPoint; 00035 00036 //2x7 derivative matrix for given particle 00037 AlgebraicMatrix lDeriv = makeDerivative(lPar).first; 00038 AlgebraicMatrix dr(2,7,0); 00039 dr.sub(1,1,lDeriv); 00040 return pair<AlgebraicMatrix,AlgebraicVector>(dr,lPar); 00041 }
AlgebraicVector SimplePointingConstraint::deviations | ( | int | nStates | ) | const [virtual] |
Returns vector of sigma squared associated to the KinematicParameters of refitted particles Initial deviations are given by user for the constraining parameters (mass, momentum components etc).
In case of multiple states exactly the same values are added to every particle parameters
Implements KinematicConstraint.
Definition at line 75 of file SimplePointingConstraint.cc.
00076 {return AlgebraicVector(7*nStates,0);}
pair< AlgebraicMatrix, AlgebraicVector > SimplePointingConstraint::makeDerivative | ( | const AlgebraicVector & | exPoint | ) | const [private] |
Definition at line 120 of file SimplePointingConstraint.cc.
References funct::pow(), refPoint, funct::sqrt(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().
Referenced by derivative().
00121 { 00122 AlgebraicMatrix dr(2,7,0); 00123 AlgebraicVector point = exPoint; 00124 double dx = point(1) - refPoint.x(); 00125 double dy = point(2) - refPoint.y(); 00126 double dz = point(3) - refPoint.z(); 00127 double px = point(4); 00128 double py = point(5); 00129 double pz = point(6); 00130 00131 00132 //half angle solution 00133 //d/dx_i 00134 dr(1,1) = (sqrt((1 + dx/sqrt(pow(dx,2) + pow(dy,2)))*(1 - px/sqrt(pow(px,2) + pow(py,2)))) - 00135 sqrt((1 - dx/sqrt(pow(dx,2) + pow(dy,2)))*(1 + px/sqrt(pow(px,2) + pow(py,2)))))/2.; 00136 00137 dr(1,2) = (((-(pow(dx,2)/pow(pow(dx,2) + pow(dy,2),1.5)) + 1/sqrt(pow(dx,2) + pow(dy,2)))* 00138 (1 - px/sqrt(pow(px,2) + pow(py,2))))/ 00139 (2.*sqrt((1 + dx/sqrt(pow(dx,2) + pow(dy,2)))*(1 - px/sqrt(pow(px,2) + pow(py,2))))) - 00140 ((pow(dx,2)/pow(pow(dx,2) + pow(dy,2),1.5) - 1/sqrt(pow(dx,2) + pow(dy,2)))* 00141 (1 + px/sqrt(pow(px,2) + pow(py,2))))/ 00142 (2.*sqrt((1 - dx/sqrt(pow(dx,2) + pow(dy,2)))*(1 + px/sqrt(pow(px,2) + pow(py,2))))))/2.; 00143 00144 00145 dr(1,3) = 0; 00146 00147 //d/dp_i 00148 //debug: x->p index xhange in denominator 00149 dr(1,4) = (-(dx*dy*(1 - px/sqrt(pow(px,2) + pow(py,2))))/ 00150 (2.*pow(pow(dx,2) + pow(dy,2),1.5)* 00151 sqrt((1 + dx/sqrt(pow(dx,2) + pow(dy,2)))*(1 - px/sqrt(pow(px,2) + pow(py,2))))) - 00152 (dx*dy*(1 + px/sqrt(pow(px,2) + pow(py,2))))/ 00153 (2.*pow(pow(dx,2) + pow(dy,2),1.5)* 00154 sqrt((1 - dx/sqrt(pow(dx,2) + pow(dy,2)))*(1 + px/sqrt(pow(px,2) + pow(py,2))))))/2.; 00155 00156 00157 dr(1,5) = (((1 + dx/sqrt(pow(dx,2) + pow(dy,2)))*px*py)/ 00158 (2.*pow(pow(px,2) + pow(py,2),1.5)* 00159 sqrt((1 + dx/sqrt(pow(dx,2) + pow(dy,2)))*(1 - px/sqrt(pow(px,2) + pow(py,2))))) + 00160 ((1 - dx/sqrt(pow(dx,2) + pow(dy,2)))*px*py)/ 00161 (2.*pow(pow(px,2) + pow(py,2),1.5)* 00162 sqrt((1 - dx/sqrt(pow(dx,2) + pow(dy,2)))*(1 + px/sqrt(pow(px,2) + pow(py,2))))))/2.; 00163 00164 00165 00166 dr(1,6) = 0; 00167 dr(1,7) = 0; 00168 00169 //2nd equation 00170 //d/dx_i 00171 00172 dr(2,1) =(((-((dx*sqrt(pow(dx,2) + pow(dy,2)))/pow(pow(dx,2) + pow(dy,2) + pow(dz,2),1.5)) + 00173 dx/(sqrt(pow(dx,2) + pow(dy,2))*sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2))))* 00174 (1 - sqrt(pow(px,2) + pow(py,2))/sqrt(pow(px,2) + pow(py,2) + pow(pz,2))))/ 00175 (2.*sqrt((1 + sqrt(pow(dx,2) + pow(dy,2))/sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2)))* 00176 (1 - sqrt(pow(px,2) + pow(py,2))/sqrt(pow(px,2) + pow(py,2) + pow(pz,2))))) - 00177 (((dx*sqrt(pow(dx,2) + pow(dy,2)))/pow(pow(dx,2) + pow(dy,2) + pow(dz,2),1.5) - 00178 dx/(sqrt(pow(dx,2) + pow(dy,2))*sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2))))* 00179 (1 + sqrt(pow(px,2) + pow(py,2))/sqrt(pow(px,2) + pow(py,2) + pow(pz,2))))/ 00180 (2.*sqrt((1 - sqrt(pow(dx,2) + pow(dy,2))/sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2)))* 00181 (1 + sqrt(pow(px,2) + pow(py,2))/sqrt(pow(px,2) + pow(py,2) + pow(pz,2))))))/2.; 00182 00183 00184 dr(2,2) = (((-((dy*sqrt(pow(dx,2) + pow(dy,2)))/pow(pow(dx,2) + pow(dy,2) + pow(dz,2),1.5)) + 00185 dy/(sqrt(pow(dx,2) + pow(dy,2))*sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2))))* 00186 (1 - sqrt(pow(px,2) + pow(py,2))/sqrt(pow(px,2) + pow(py,2) + pow(pz,2))))/ 00187 (2.*sqrt((1 + sqrt(pow(dx,2) + pow(dy,2))/sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2)))* 00188 (1 - sqrt(pow(px,2) + pow(py,2))/sqrt(pow(px,2) + pow(py,2) + pow(pz,2))))) - 00189 (((dy*sqrt(pow(dx,2) + pow(dy,2)))/pow(pow(dx,2) + pow(dy,2) + pow(dz,2),1.5) - 00190 dy/(sqrt(pow(dx,2) + pow(dy,2))*sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2))))* 00191 (1 + sqrt(pow(px,2) + pow(py,2))/sqrt(pow(px,2) + pow(py,2) + pow(pz,2))))/ 00192 (2.*sqrt((1 - sqrt(pow(dx,2) + pow(dy,2))/sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2)))* 00193 (1 + sqrt(pow(px,2) + pow(py,2))/sqrt(pow(px,2) + pow(py,2) + pow(pz,2))))))/2.; 00194 00195 00196 dr(2,3) = (-(sqrt(pow(dx,2) + pow(dy,2))*dz*(1 - sqrt(pow(px,2) + pow(py,2))/ 00197 sqrt(pow(px,2) + pow(py,2) + pow(pz,2))))/ 00198 (2.*pow(pow(dx,2) + pow(dy,2) + pow(dz,2),1.5)* 00199 sqrt((1 + sqrt(pow(dx,2) + pow(dy,2))/sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2)))* 00200 (1 - sqrt(pow(px,2) + pow(py,2))/sqrt(pow(px,2) + pow(py,2) + pow(pz,2))))) - 00201 (sqrt(pow(dx,2) + pow(dy,2))*dz*(1 + 00202 sqrt(pow(px,2) + pow(py,2))/sqrt(pow(px,2) + pow(py,2) + pow(pz,2))))/ 00203 (2.*pow(pow(dx,2) + pow(dy,2) + pow(dz,2),1.5)* 00204 sqrt((1 - sqrt(pow(dx,2) + pow(dy,2))/sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2)))* 00205 (1 + sqrt(pow(px,2) + pow(py,2))/sqrt(pow(px,2) + pow(py,2) + pow(pz,2))))))/2.; 00206 00207 00208 00209 //d/dp_i 00210 //debug: x->p index xhange in denominator 00211 00212 dr(2,4) = (((1 + sqrt(pow(dx,2) + pow(dy,2))/sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2)))* 00213 ((px*sqrt(pow(px,2) + pow(py,2)))/pow(pow(px,2) + pow(py,2) + pow(pz,2),1.5) - 00214 px/(sqrt(pow(px,2) + pow(py,2))*sqrt(pow(px,2) + pow(py,2) + pow(pz,2)))))/ 00215 (2.*sqrt((1 + sqrt(pow(dx,2) + pow(dy,2))/sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2)))* 00216 (1 - sqrt(pow(px,2) + pow(py,2))/sqrt(pow(px,2) + pow(py,2) + pow(pz,2))))) - 00217 ((1 - sqrt(pow(dx,2) + pow(dy,2))/sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2)))* 00218 (-((px*sqrt(pow(px,2) + pow(py,2)))/pow(pow(px,2) + pow(py,2) + pow(pz,2),1.5)) + 00219 px/(sqrt(pow(px,2) + pow(py,2))*sqrt(pow(px,2) + pow(py,2) + pow(pz,2)))))/ 00220 (2.*sqrt((1 - sqrt(pow(dx,2) + pow(dy,2))/sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2)))* 00221 (1 + sqrt(pow(px,2) + pow(py,2))/sqrt(pow(px,2) + pow(py,2) + pow(pz,2))))))/2.; 00222 00223 00224 dr(2,5) = (((1 + sqrt(pow(dx,2) + pow(dy,2))/sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2)))* 00225 ((py*sqrt(pow(px,2) + pow(py,2)))/pow(pow(px,2) + pow(py,2) + pow(pz,2),1.5) - 00226 py/(sqrt(pow(px,2) + pow(py,2))*sqrt(pow(px,2) + pow(py,2) + pow(pz,2)))))/ 00227 (2.*sqrt((1 + sqrt(pow(dx,2) + pow(dy,2))/sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2)))* 00228 (1 - sqrt(pow(px,2) + pow(py,2))/sqrt(pow(px,2) + pow(py,2) + pow(pz,2))))) - 00229 ((1 - sqrt(pow(dx,2) + pow(dy,2))/sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2)))* 00230 (-((py*sqrt(pow(px,2) + pow(py,2)))/pow(pow(px,2) + pow(py,2) + pow(pz,2),1.5)) + 00231 py/(sqrt(pow(px,2) + pow(py,2))*sqrt(pow(px,2) + pow(py,2) + pow(pz,2)))))/ 00232 (2.*sqrt((1 - sqrt(pow(dx,2) + pow(dy,2))/sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2)))* 00233 (1 + sqrt(pow(px,2) + pow(py,2))/sqrt(pow(px,2) + pow(py,2) + pow(pz,2))))))/2.; 00234 00235 00236 dr(2,6) = (((1 + sqrt(pow(dx,2) + pow(dy,2))/sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2)))* 00237 sqrt(pow(px,2) + pow(py,2))*pz)/ 00238 (2.*pow(pow(px,2) + pow(py,2) + pow(pz,2),1.5)* 00239 sqrt((1 + sqrt(pow(dx,2) + pow(dy,2))/sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2)))* 00240 (1 - sqrt(pow(px,2) + pow(py,2))/sqrt(pow(px,2) + pow(py,2) + pow(pz,2))))) + 00241 ((1 - sqrt(pow(dx,2) + pow(dy,2))/sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2)))* 00242 sqrt(pow(px,2) + pow(py,2))*pz)/ 00243 (2.*pow(pow(px,2) + pow(py,2) + pow(pz,2),1.5)* 00244 sqrt((1 - sqrt(pow(dx,2) + pow(dy,2))/sqrt(pow(dx,2) + pow(dy,2) + pow(dz,2)))* 00245 (1 + sqrt(pow(px,2) + pow(py,2))/sqrt(pow(px,2) + pow(py,2) + pow(pz,2))))))/2.; 00246 00247 00248 dr(2,7) = 0; 00249 00250 // cout<<"derivative matrix "<<dr<<endl; 00251 return pair<AlgebraicMatrix,AlgebraicVector>(dr,point); 00252 }
pair< AlgebraicVector, AlgebraicVector > SimplePointingConstraint::makeValue | ( | const AlgebraicVector & | exPoint | ) | const [private] |
Definition at line 81 of file SimplePointingConstraint.cc.
References refPoint, funct::sqrt(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().
Referenced by value().
00082 { 00083 // cout<<"Make value called"<<endl; 00084 AlgebraicVector vl(2,0); 00085 AlgebraicVector point = exPoint; 00086 double dx = point(1) - refPoint.x(); 00087 double dy = point(2) - refPoint.y(); 00088 double dz = point(3) - refPoint.z(); 00089 double px = point(4); 00090 double py = point(5); 00091 double pz = point(6); 00092 00093 00094 //half angle solution: sin((alpha - betha)/2) 00095 double cos_phi_p = px/sqrt(px*px + py*py); 00096 double cos_phi_x = dx/sqrt(dx*dx + dy*dy); 00097 // cout<<"mom cos phi"<<cos_phi_p<<endl; 00098 // cout<<"x cos phi"<<cos_phi_x<<endl; 00099 00100 00101 double cos_theta_p = sqrt(px*px + py*py)/sqrt(px*px + py*py + pz*pz); 00102 double cos_theta_x = sqrt(dx*dx + dy*dy)/sqrt(dx*dx + dy*dy + dz*dz); 00103 00104 float feq = sqrt((1-cos_phi_p)*(1+cos_phi_x)) - sqrt((1+cos_phi_p)*(1-cos_phi_x)); 00105 float seq = sqrt((1-cos_theta_p)*(1+cos_theta_x)) - sqrt((1+cos_theta_p)*(1-cos_theta_x)); 00106 00107 // cout<<"First factor: "<<feq/2<<endl; 00108 // cout<<"Second factor: "<<seq/2<<endl; 00109 00110 vl(1) = feq/2; 00111 vl(2) = seq/2; 00112 00113 // cout<<"Value "<<vl<<endl; 00114 //half angle corrected 00115 // vl(1) = (sin_x/(1+cos_x)) - (sin_p/(1+cos_p)); 00116 // vl(2) = (sin_xt/(1+cos_xt)) - (sin_pt/(1+cos_pt)); 00117 return pair<AlgebraicVector,AlgebraicVector>(vl,point); 00118 }
int SimplePointingConstraint::numberOfEquations | ( | ) | const [virtual] |
Returns number of constraint equations used for fitting.
Method is relevant for proper NDF calculations.
Implements KinematicConstraint.
Definition at line 78 of file SimplePointingConstraint.cc.
pair< AlgebraicVector, AlgebraicVector > SimplePointingConstraint::value | ( | const vector< RefCountedKinematicParticle > | par | ) | const [virtual] |
Methods making value and derivative matrix using current state parameters as expansion 7-point.
Constraint can be made equaly for single and multiple states
Implements KinematicConstraint.
Definition at line 60 of file SimplePointingConstraint.cc.
References makeValue().
00061 { 00062 int nStates = par.size(); 00063 if(nStates == 0) throw VertexException("PointingKinematicConstraint::Empty vector of particles passed"); 00064 if(nStates != 1) throw VertexException("PointingKinematicConstraint::Current version does not support the multistate refit"); 00065 AlgebraicVector vl(2,0); 00066 AlgebraicVector lPoint = asHepVector<7>(par.front()->currentState().kinematicParameters().vector()); 00067 vl(1) = makeValue(lPoint).first(1); 00068 vl(2) = makeValue(lPoint).first(2); 00069 // cout<<"Value returned: "<<vl<<endl; 00070 // cout<<"For the point: "<<lPoint<<endl; 00071 00072 return pair<AlgebraicVector,AlgebraicVector>(vl,lPoint); 00073 }
pair< AlgebraicVector, AlgebraicVector > SimplePointingConstraint::value | ( | const AlgebraicVector & | exPoint | ) | const [virtual] |
Vector of values and matrix of derivatives calculated at given expansion 7xNumberOfStates point.
Implements KinematicConstraint.
Definition at line 5 of file SimplePointingConstraint.cc.
References makeValue().
00006 { 00007 if(exPoint.num_row() ==0 ) throw VertexException("PointingKinematicConstraint::value requested for zero Linearization point"); 00008 00009 //security check for extended cartesian parametrization 00010 int inSize = exPoint.num_row(); 00011 if((inSize%7) !=0) throw VertexException("PointingKinematicConstraint::linearization point has a wrong dimension"); 00012 int nStates = inSize/7; 00013 if(nStates != 1) throw VertexException("PointingKinematicConstraint::Current version does not support the multistate refit"); 00014 00015 AlgebraicVector lPar = exPoint; 00016 AlgebraicVector vl(2,0); 00017 00018 //vector of values 1x2 for given particle 00019 AlgebraicVector lValue = makeValue(lPar).first; 00020 vl(1) =lValue(1); 00021 vl(2) =lValue(2); 00022 return pair<AlgebraicVector,AlgebraicVector>(vl,lPar); 00023 }
Definition at line 61 of file SimplePointingConstraint.h.
Referenced by makeDerivative(), and makeValue().