Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 13:33:36

0001 #include <string>
0002 
0003 #include "RecoBTag/BTagTools/interface/SignedDecayLength3D.h"
0004 
0005 #include "DataFormats/GeometrySurface/interface/Line.h"
0006 
0007 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
0008 #include "TrackingTools/GeomPropagators/interface/AnalyticalTrajectoryExtrapolatorToLine.h"
0009 #include "DataFormats/GeometryCommonDetAlgo/interface/Measurement1D.h"
0010 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
0011 
0012 #include "CLHEP/Vector/ThreeVector.h"
0013 #include "CLHEP/Vector/LorentzVector.h"
0014 #include "CLHEP/Matrix/Vector.h"
0015 
0016 using namespace std;
0017 using namespace reco;
0018 
0019 pair<bool, Measurement1D> SignedDecayLength3D::apply(const TransientTrack& transientTrack,
0020                                                      const GlobalVector& direction,
0021                                                      const Vertex& vertex) {
0022   double theError = 0.;
0023   bool theIsValid;
0024 
0025   //TrajectoryStateOnSurface TSOS = (aRecTrack).impactPointStateOnSurface();
0026   TrajectoryStateOnSurface TSOS = transientTrack.impactPointState();
0027   const FreeTrajectoryState* FTS = TSOS.freeTrajectoryState();
0028 
0029   TrajectoryStateOnSurface theTSOS = closestApproachToJet(*FTS, vertex, direction, transientTrack.field());
0030   theIsValid = theTSOS.isValid();
0031 
0032   if (theIsValid) {
0033     GlobalVector J = direction.unit();
0034     GlobalPoint vertexPosition(vertex.x(), vertex.y(), vertex.z());
0035 
0036     double theValue = J.dot(theTSOS.globalPosition() - vertexPosition);
0037 
0038     //error calculation
0039 
0040     AlgebraicVector3 j;
0041     j[0] = J.x();
0042     j[1] = J.y();
0043     j[2] = J.z();
0044     AlgebraicVector6 jj;
0045     jj[0] = J.x();
0046     jj[1] = J.y();
0047     jj[2] = J.z();
0048     jj[3] = 0.;
0049     jj[4] = 0.;
0050     jj[5] = 0.;  ///chech it!!!!!!!!!!!!!!!!!!!!!!!
0051     double E1 = ROOT::Math::Similarity(jj, theTSOS.cartesianError().matrix());
0052     //  double E2 = (aJet.vertex().positionError().matrix()).similarity(j);
0053     double E2 = ROOT::Math::Similarity(j, vertex.covariance());
0054 
0055     theError = sqrt(E1 + E2);
0056 
0057     //cout<< "Error ="<< theError<<endl;
0058     Measurement1D A(theValue, theError);
0059     return pair<bool, Measurement1D>(theIsValid, A);
0060   } else {
0061     return pair<bool, Measurement1D>(theIsValid, Measurement1D(0., 0.));
0062   }  // endif (isValid)
0063 }  // end constructor declaration
0064 
0065 TrajectoryStateOnSurface SignedDecayLength3D::closestApproachToJet(const FreeTrajectoryState& aFTS,
0066                                                                    const Vertex& vertex,
0067                                                                    const GlobalVector& aJetDirection,
0068                                                                    const MagneticField* field) {
0069   GlobalVector J = aJetDirection.unit();
0070 
0071   Line::PositionType pos(GlobalPoint(vertex.x(), vertex.y(), vertex.z()));
0072   Line::DirectionType dir(J);
0073   Line Jet(pos, dir);
0074 
0075   AnalyticalTrajectoryExtrapolatorToLine TETL(field);
0076 
0077   return TETL.extrapolate(aFTS, Jet);
0078 }