File indexing completed on 2024-04-06 12:31:41
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 }