GhostTrackFitter/src/VertexGhostTrackState.cc

0001 #include <cmath>
0002
0003 #include <Math/SMatrix.h>
0004 #include <Math/MatrixFunctions.h>
0005
0006 #include "DataFormats/GeometryVector/interface/GlobalPoint.h"
0007 #include "DataFormats/GeometryVector/interface/GlobalVector.h"
0008
0009 #include "RecoVertex/GhostTrackFitter/interface/GhostTrackPrediction.h"
0010
0011 #include "RecoVertex/GhostTrackFitter/interface/VertexGhostTrackState.h"
0012
0013 using namespace reco;
0014
0015 namespace {
0016 #ifndef __clang__
0017   static inline double sqr(double arg) { return arg * arg; }
0018 #endif
0019   using namespace ROOT::Math;
0020
0021   typedef SVector<double, 3> Vector3;
0022
0023   typedef SMatrix<double, 3, 3, MatRepSym<double, 3> > Matrix3S;
0024   typedef SMatrix<double, 6, 6, MatRepSym<double, 6> > Matrix6S;
0025   typedef SMatrix<double, 3, 3> Matrix33;
0026   typedef SMatrix<double, 3, 6> Matrix36;
0027
0028   inline Vector3 conv(const GlobalVector &vec) {
0029     Vector3 result;
0030     result[0] = vec.x();
0031     result[1] = vec.y();
0032     result[2] = vec.z();
0033     return result;
0034   }
0035 }  // namespace
0036
0037 BasicGhostTrackState::Vertex VertexGhostTrackState::vertexStateOnGhostTrack(const GhostTrackPrediction &pred,
0038                                                                             bool withMeasurementError) const {
0039   using namespace ROOT::Math;
0040
0041   GlobalPoint origin = pred.origin();
0042   GlobalVector direction = pred.direction();
0043
0044   double rho2 = pred.rho2();
0045   double rho = std::sqrt(rho2);
0046   double lambda = (position_ - origin) * direction / rho2;
0047   GlobalPoint pos = origin + lambda * direction;
0048
0049   Vector3 b = conv(direction) / rho;
0050   Vector3 ca = conv(position_ - pos);
0051
0052   Matrix33 pA = TensorProd(b, b);
0053   Matrix33 pB = TensorProd(b, ca);
0054
0055   Matrix36 jacobian;
0056   jacobian.Place_at(-pA + Matrix33(SMatrixIdentity()), 0, 0);
0057   jacobian.Place_at(pB / rho, 0, 3);
0058   Matrix3S error = Similarity(jacobian, pred.cartesianError(lambda));
0059
0060   if (withMeasurementError)
0061     error += Similarity(pA, covariance_);
0062
0063   return Vertex(pos, error);
0064 }
0065
0066 BasicGhostTrackState::Vertex VertexGhostTrackState::vertexStateOnMeasurement(const GhostTrackPrediction &pred,
0067                                                                              bool withGhostTrackError) const {
0068   return Vertex(position_, covariance_);
0069 }