Go to the documentation of this file.00001 #include "TrackingTools/PatternTools/interface/TwoTrackMinimumDistanceLineLine.h"
00002 #include "TrackingTools/TrajectoryParametrization/interface/GlobalTrajectoryParameters.h"
00003 #include "DataFormats/GeometryVector/interface/GlobalVector.h"
00004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00005 #include "MagneticField/Engine/interface/MagneticField.h"
00006
00007
00008 bool TwoTrackMinimumDistanceLineLine::calculate(
00009 const GlobalTrajectoryParameters & theG,
00010 const GlobalTrajectoryParameters & theH)
00011 {
00012 GlobalPoint gOrig = theG.position();
00013 GlobalPoint hOrig = theH.position();
00014 bool isLineG = theG.magneticField().inTesla(gOrig).z() == 0 || theG.charge() == 0;
00015 bool isLineH = theH.magneticField().inTesla(hOrig).z() == 0 || theH.charge() == 0;
00016 if ( ! (isLineG && isLineH) )
00017 {
00018 edm::LogWarning ("TwoTrackMinimumDistanceLineLine")
00019 << "charge of input track is not zero or field non zero"
00020 << "\n positions: "<<gOrig<<" , "<<hOrig
00021 << "\n Bz fields: "<<theG.magneticField().inTesla(gOrig).z()<<" , "<<theH.magneticField().inTesla(hOrig).z()
00022 << "\n charges: "<<theG.charge()<<" , "<<theH.charge();
00023 return true;
00024 };
00025
00026 GlobalVector gVec = theG.momentum();
00027 const double gMag= theG.momentum().mag(); double gMag2 = gMag*gMag;
00028 GlobalVector hVec = theH.momentum();
00029 const double hMag= theH.momentum().mag(); double hMag2 = hMag*hMag;
00030
00031 if ( gMag == 0. || hMag == 0. )
00032 {
00033 edm::LogWarning ("TwoTrackMinimumDistanceLineLine")
00034 << "momentum of input trajectory is zero.";
00035 return true;
00036 };
00037
00038 phiG = theG.momentum().phi();
00039 phiH = theH.momentum().phi();
00040
00041 double gVec_Dot_hVec = gVec.dot(hVec);
00042 double norm = gVec_Dot_hVec*gVec_Dot_hVec - gMag2*hMag2;
00043
00044 if ( norm == 0 )
00045 {
00046 edm::LogWarning ("TwoTrackMinimumDistanceLineLine")
00047 << "Tracks are parallel.";
00048 return true;
00049 }
00050
00051 GlobalVector posDiff = gOrig - hOrig;
00052
00053 double tG = (posDiff.dot(gVec) * hMag2 - gVec_Dot_hVec * posDiff.dot(hVec))/ norm;
00054 double tH = (gVec_Dot_hVec * posDiff.dot(gVec) - posDiff.dot(hVec) * gMag2)/ norm;
00055
00056 gPos = GlobalPoint ( gOrig + tG * gVec);
00057 hPos = GlobalPoint ( hOrig + tH * hVec);
00058
00059
00060 pathG = tG * gMag;
00061 pathH = tH * hMag;
00062 return false;
00063 }
00064
00065 std::pair<GlobalPoint, GlobalPoint> TwoTrackMinimumDistanceLineLine::points() const
00066 {
00067 return std::pair<GlobalPoint, GlobalPoint > ( gPos, hPos );
00068 }
00069
00070 std::pair<double, double> TwoTrackMinimumDistanceLineLine::pathLength() const
00071 {
00072 return std::pair<double, double> ( pathG, pathH);
00073 }