TrajectoryState/test/testTSOS.cpp

0001 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
0002 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateAccessor.h"
0003 #include "DataFormats/GeometrySurface/interface/Surface.h"
0004 #include "DataFormats/GeometrySurface/interface/BoundPlane.h"
0005 #include "MagneticField/Engine/interface/MagneticField.h"
0006
0007 #include "TrackingTools/TrajectoryState/interface/BasicSingleTrajectoryState.h"
0008 #include "TrackingTools/TrajectoryState/interface/FreeTrajectoryState.h"
0009 #include "TrackingTools/AnalyticalJacobians/interface/JacobianCurvilinearToCartesian.h"
0010 #include "TrackingTools/AnalyticalJacobians/interface/JacobianCartesianToCurvilinear.h"
0011 #include "TrackingTools/TrajectoryState/interface/PerigeeConversions.h"
0012
0013 #include <iostream>
0014
0015 class ConstMagneticField final : public MagneticField {
0016 public:
0017   ConstMagneticField() { setNominalValue(); }
0018   GlobalVector inTesla(const GlobalPoint&) const override { return GlobalVector(0, 0, 4); }
0019 };
0020
0021 int main() {
0022   std::cout << "sizes tsos, bsts, fts" << std::endl;
0023   std::cout << sizeof(TrajectoryStateOnSurface) << std::endl;
0024   std::cout << sizeof(BasicSingleTrajectoryState) << std::endl;
0025   std::cout << sizeof(FreeTrajectoryState) << std::endl;
0026
0027   using namespace std;
0028
0029   MagneticField* field = new ConstMagneticField;
0030   GlobalPoint gp(0, 0, 0);
0031   GlobalVector gv(1, 1, 1);
0032   GlobalTrajectoryParameters gtp(gp, gv, 1, field);
0033   double v[15] = {0.01, -0.01, 0., 0., 0., 0.01, 0., 0., 0., 0.01, 0., 0., 1., 0., 1.};
0034   AlgebraicSymMatrix55 gerr(v, 15);
0035   BoundPlane* plane = new BoundPlane(gp, Surface::RotationType());
0036
0037   TrajectoryStateOnSurface ts(gtp, gerr, *plane);
0038
0039   cout << "ts.globalMomentum() " << ts.globalMomentum() << endl;
0040   cout << "ts.localMomentum()  " << ts.localMomentum() << endl;
0041   cout << "ts.transverseCurvature()  " << ts.transverseCurvature() << endl;
0042   cout << "ts inversePtErr " << TrajectoryStateAccessor(*ts.freeState()).inversePtError() << std::endl;
0043   cout << "ts curv err\n" << ts.curvilinearError().matrix() << std::endl;
0044   cout << "ts cart err\n" << ts.cartesianError().matrix() << std::endl;
0045   {
0046     JacobianCartesianToCurvilinear cart2Curv(ts.globalParameters());
0047     const AlgebraicMatrix56& jac = cart2Curv.jacobian();
0048
0049     CurvilinearTrajectoryError theCurvilinearError = ROOT::Math::Similarity(jac, ts.cartesianError().matrix());
0050     cout << "curv from cart \n" << theCurvilinearError.matrix() << std::endl;
0051
0052     auto a55 = PerigeeConversions::jacobianCurvilinear2Perigee(*ts.freeState());
0053     std::cout << " curv 2 per " << a55 << std::endl;
0054     a55 = PerigeeConversions::jacobianPerigee2Curvilinear(gtp);
0055     std::cout << " per 2 cuv " << a55 << std::endl;
0056     auto a66 = PerigeeConversions::jacobianParameters2Cartesian({1., 1., 1.}, gp, 1, field);
0057     std::cout << " per 2 cart " << a66 << std::endl;
0058   }
0059
0060   LocalPoint lp(0, 0, 0);
0061   LocalVector lv(1, 1, 1);
0062   LocalTrajectoryParameters ltp(lp, lv, 1);
0063   LocalTrajectoryError lerr(1., 1., 0.1, 0.1, 0.1);
0064   TrajectoryStateOnSurface ts2(ltp, lerr, *plane, field);
0065   cout << "ts2.globalMomentum() " << ts2.globalMomentum() << endl;
0066   cout << "ts2.localMomentum()  " << ts2.localMomentum() << endl;
0067   cout << "ts2.transverseCurvature()  " << ts2.transverseCurvature() << endl;
0068   cout << "ts2 inversePtErr " << TrajectoryStateAccessor(*ts2.freeState()).inversePtError() << std::endl;
0069   cout << "ts2 curv err\n" << ts2.curvilinearError().matrix() << std::endl;
0070   cout << "ts2 cart err\n" << ts2.cartesianError().matrix() << std::endl;
0071   {
0072     JacobianCartesianToCurvilinear cart2Curv(ts2.globalParameters());
0073     const AlgebraicMatrix56& jac = cart2Curv.jacobian();
0074
0075     CurvilinearTrajectoryError theCurvilinearError = ROOT::Math::Similarity(jac, ts2.cartesianError().matrix());
0076     cout << "curv from cart \n" << theCurvilinearError.matrix() << std::endl;
0077
0078     auto a55 = PerigeeConversions::jacobianCurvilinear2Perigee(*ts2.freeState());
0079     std::cout << " curv 2 per " << a55 << std::endl;
0080   }
0081 }