Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-04-02 03:10:29

0001 #ifndef CUDADataFormats_Track_TrajectoryStateSOAT_H
0002 #define CUDADataFormats_Track_TrajectoryStateSOAT_H
0003 
0004 #include <Eigen/Dense>
0005 #include "HeterogeneousCore/CUDAUtilities/interface/eigenSoA.h"
0006 
0007 template <int32_t S>
0008 struct TrajectoryStateSoAT {
0009   using Vector5f = Eigen::Matrix<float, 5, 1>;
0010   using Vector15f = Eigen::Matrix<float, 15, 1>;
0011 
0012   using Vector5d = Eigen::Matrix<double, 5, 1>;
0013   using Matrix5d = Eigen::Matrix<double, 5, 5>;
0014 
0015   static constexpr int32_t stride() { return S; }
0016 
0017   eigenSoA::MatrixSoA<Vector5f, S> state;
0018   eigenSoA::MatrixSoA<Vector15f, S> covariance;
0019 
0020   template <typename V3, typename M3, typename V2, typename M2>
0021   __host__ __device__ inline void copyFromCircle(
0022       V3 const& cp, M3 const& ccov, V2 const& lp, M2 const& lcov, float b, int32_t i) {
0023     state(i) << cp.template cast<float>(), lp.template cast<float>();
0024     state(i)(2) *= b;
0025     auto cov = covariance(i);
0026     cov(0) = ccov(0, 0);
0027     cov(1) = ccov(0, 1);
0028     cov(2) = b * float(ccov(0, 2));
0029     cov(4) = cov(3) = 0;
0030     cov(5) = ccov(1, 1);
0031     cov(6) = b * float(ccov(1, 2));
0032     cov(8) = cov(7) = 0;
0033     cov(9) = b * b * float(ccov(2, 2));
0034     cov(11) = cov(10) = 0;
0035     cov(12) = lcov(0, 0);
0036     cov(13) = lcov(0, 1);
0037     cov(14) = lcov(1, 1);
0038   }
0039 
0040   template <typename V5, typename M5>
0041   __host__ __device__ inline void copyFromDense(V5 const& v, M5 const& cov, int32_t i) {
0042     state(i) = v.template cast<float>();
0043     for (int j = 0, ind = 0; j < 5; ++j)
0044       for (auto k = j; k < 5; ++k)
0045         covariance(i)(ind++) = cov(j, k);
0046   }
0047 
0048   template <typename V5, typename M5>
0049   __host__ __device__ inline void copyToDense(V5& v, M5& cov, int32_t i) const {
0050     v = state(i).template cast<typename V5::Scalar>();
0051     for (int j = 0, ind = 0; j < 5; ++j) {
0052       cov(j, j) = covariance(i)(ind++);
0053       for (auto k = j + 1; k < 5; ++k)
0054         cov(k, j) = cov(j, k) = covariance(i)(ind++);
0055     }
0056   }
0057 };
0058 
0059 #endif  // CUDADataFormats_Track_TrajectoryStateSOAT_H