directly call the vector sincos function

DetectorDescription/GeoModel/GeoSpecialShapes/GeoSpecialShapes/LArWheelCalculator.h

@@ -126,9 +126,6 @@ class LArWheelCalculator
     std::pair<int, int> GetPhiGapAndSide(const CLHEP::Hep3Vector &p) const;
     double AmplitudeOfSurface(const CLHEP::Hep3Vector& P, int side, int fan_number) const;
 
-    /// Set the sincos calculator
-    // void SetVectorizedSincos() {m_sincos_calculator = &m_vsincos_par.eval;};
-
     /// @}
 
   private:
@@ -220,10 +217,6 @@ class LArWheelCalculator
     vsincos_par m_vsincos_par{};
 #endif
 
-  private:
-    typedef void (LArWheelCalculator::*SincosCalculator)(const double, double &, double &) const;
-    SincosCalculator m_sincos_calculator;
-    void EvalSincos(const double P, double &x, double &y) const {(this->*m_sincos_calculator)(P, x, y);};
 };
 
 #ifndef XAOD_STANDALONE

DetectorDescription/GeoModel/GeoSpecialShapes/GeoSpecialShapes/vec_parametrized_sincos.h

@@ -3,9 +3,10 @@
 */
 
 /**
- * @brief vectorized version of parametrized sincos
- * @ author Miha Muskinja, Chistos
- */
+  * @brief vectorized version of parametrized sincos
+  * see ATLASSIM-4753 for details
+  * @author Miha Muskinja, Chistos Anastopoulos
+*/
 #ifndef VEC_PARAMETRIZED_SINCOS_H
 #define VEC_PARAMETRIZED_SINCOS_H
 
@@ -14,61 +15,61 @@
 struct vsincos_par
 {
 
-  CxxUtils::vec<double, 4> param_0 = {};
-  CxxUtils::vec<double, 4> param_1 = {};
-  CxxUtils::vec<double, 4> param_2 = {};
+    CxxUtils::vec<double, 4> param_0 = {};
+    CxxUtils::vec<double, 4> param_1 = {};
+    CxxUtils::vec<double, 4> param_2 = {};
 
 #if HAVE_FUNCTION_MULTIVERSIONING
 #if defined(__x86_64__)
 
-  __attribute__((target("avx2,fma")))
+    __attribute__((target("avx2,fma")))
+
+    void
+    evaluate(const double r,
+             double &ATH_RESTRICT sin_a,
+             double &ATH_RESTRICT cos_a) const ATH_RESTRICT
+    {
+        const double r2 = r * r;
+        CxxUtils::vec<double, 4> P = r2 * param_0 + param_1;
+        P = r2 * P + param_2;
+        CxxUtils::vec<double, 4> P2 = {P[1], P[0], P[3], P[2]};
+        CxxUtils::vec<double, 4> res = r * P2 + P;
+        sin_a = res[0];
+        cos_a = res[2];
+    }
 
-  void
-  evaluate(const double r,
-           double& ATH_RESTRICT sin_a,
-           double& ATH_RESTRICT cos_a) const ATH_RESTRICT
-  {
-    const double r2 = r * r;
-    CxxUtils::vec<double, 4> P = r2 * param_0 + param_1;
-    P = r2 * P + param_2;
-    CxxUtils::vec<double, 4> P2 = { P[1], P[0], P[3], P[2] };
-    CxxUtils::vec<double, 4> res = r * P2 + P;
-    sin_a = res[0];
-    cos_a = res[2];
-  }
+    __attribute__((target("avx2")))
 
-  __attribute__((target("avx2")))
+    void
+    eval(const double r,
+         double &ATH_RESTRICT sin_a,
+         double &ATH_RESTRICT cos_a) const ATH_RESTRICT
+    {
+        const double r2 = r * r;
+        CxxUtils::vec<double, 4> P = r2 * param_0 + param_1;
+        P = r2 * P + param_2;
+        CxxUtils::vec<double, 4> P2 = {P[1], P[0], P[3], P[2]};
+        CxxUtils::vec<double, 4> res = r * P2 + P;
+        sin_a = res[0];
+        cos_a = res[2];
+    }
 
-  void
-  eval(const double r,
-       double& ATH_RESTRICT sin_a,
-       double& ATH_RESTRICT cos_a) const ATH_RESTRICT
-  {
-    const double r2 = r * r;
-    CxxUtils::vec<double, 4> P = r2 * param_0 + param_1;
-    P = r2 * P + param_2;
-    CxxUtils::vec<double, 4> P2 = { P[1], P[0], P[3], P[2] };
-    CxxUtils::vec<double, 4> res = r * P2 + P;
-    sin_a = res[0];
-    cos_a = res[2];
-  }
+    __attribute__((target("default")))
 
-  __attribute__((target("default")))
 #endif // x86
 #endif // FMV
 
-  void
-  eval(const double r,
-       double& ATH_RESTRICT sin_a,
-       double& ATH_RESTRICT cos_a) const ATH_RESTRICT
-  {
-    const double r2 = r * r;
-    CxxUtils::vec<double, 4> P = r2 * param_0 + param_1;
-    P = r2 * P + param_2;
-    sin_a = r * P[1] + P[0];
-    cos_a = r * P[3] + P[2];
-  }
+    void
+    eval(const double r,
+         double &ATH_RESTRICT sin_a,
+         double &ATH_RESTRICT cos_a) const ATH_RESTRICT
+    {
+        const double r2 = r * r;
+        CxxUtils::vec<double, 4> P = r2 * param_0 + param_1;
+        P = r2 * P + param_2;
+        sin_a = r * P[1] + P[0];
+        cos_a = r * P[3] + P[2];
+    }
 };
 
 #endif
-

DetectorDescription/GeoModel/GeoSpecialShapes/src/LArWheelCalculator.cxx

@@ -443,9 +443,6 @@ LArWheelCalculator::LArWheelCalculator(LArG4::LArWheelCalculator_t a_wheelType,
      fclose(O);
      exit(0);
   */
-
-   // set the sincos calculator
-   m_sincos_calculator = &LArWheelCalculator::parameterized_sincos;
 }
 
 /* converts module gap number into wheel gap number */

DetectorDescription/GeoModel/GeoSpecialShapes/src/LArWheelCalculator_Impl/DistanceCalculatorSaggingOff.cxx

@@ -48,7 +48,7 @@ namespace LArWheelCalculator_Impl
     const double cos_a = scalpha.cs, sin_a = scalpha.sn;
 #else // parameterized sine
     double cos_a, sin_a;
-    lwc()->EvalSincos(P.y(), sin_a, cos_a);
+    lwc()->m_vsincos_par.eval(P.y(), sin_a, cos_a);
 #endif
     // determination of the nearest quarter-wave number
     int nqwave = (z < 0.) ? 0 : int(z / lwc()->m_QuarterWaveLength);
@@ -160,7 +160,7 @@ namespace LArWheelCalculator_Impl
     double cos_a = scalpha.cs, sin_a = scalpha.sn;
 #else // parameterized sine
     double cos_a, sin_a;
-    lwc()->EvalSincos(P.y(), sin_a, cos_a);
+    lwc()->m_vsincos_par.eval(P.y(), sin_a, cos_a);
 #endif
 
     bool sqw = false;
@@ -251,7 +251,7 @@ namespace LArWheelCalculator_Impl
     const double cos_a = scalpha.cs, sin_a = scalpha.sn;
 #else // parameterized sine
     double cos_a, sin_a;
-    lwc()->EvalSincos(P.y(), sin_a, cos_a);
+    lwc()->m_vsincos_par.eval(P.y(), sin_a, cos_a);
 #endif
 
     int nqwave;
@@ -339,7 +339,7 @@ namespace LArWheelCalculator_Impl
     double cos_a = scalpha.cs, sin_a = scalpha.sn;
 #else // parameterized sine
     double cos_a, sin_a;
-    lwc()->EvalSincos(P.y(), sin_a, cos_a);
+    lwc()->m_vsincos_par.eval(P.y(), sin_a, cos_a);
 #endif
 
     bool sqw = false;
@@ -434,7 +434,7 @@ namespace LArWheelCalculator_Impl
     // parameterized sine
 #else
     double cos_a, sin_a;
-    lwc()->EvalSincos(P.y(), sin_a, cos_a);
+    lwc()->m_vsincos_par.eval(P.y(), sin_a, cos_a);
 #endif
 
     // determination of the nearest quarter-wave number

DetectorDescription/GeoModel/GeoSpecialShapes/src/LArWheelCalculator_Impl/sincos_poly.cxx

@@ -99,6 +99,8 @@ void LArWheelCalculator::fill_sincos_parameterization()
       m_cos_parametrization[i] = cos_parametrization[S][i];
     }
 
+    // Parameterization for the vectorized sincos calculation
+    // see ATLASSIM-4753 for details
     // s4, s5, c4, c5
     // s2, s3, c2, c3
     // s0, s1, c0, c1
@@ -155,6 +157,8 @@ void LArWheelCalculator::fill_sincos_parameterization()
     cos_parametrization[S][i] = params_cos[i];
   }
 
+  // Parameterization for the vectorized sincos calculation
+  // see ATLASSIM-4753 for details
   // s4, s5, c4, c5
   // s2, s3, c2, c3
   // s0, s1, c0, c1
@@ -191,8 +195,12 @@ void LArWheelCalculator::fill_sincos_parameterization()
     parameterized_sincos(r, sin_a, cos_a);
     m_vsincos_par.eval(r, sin_a_v, cos_a_v);
 #if DEBUGPRINT
-    std::cout << "default: " << r << " " << sin_a << " " << cos_a << std::endl;
+    std::streamsize ss = std::cout.precision();
+    std::cout.precision(16);
+    std::cout << "def: " << r << " " << sin_a << " " << cos_a << std::endl;
     std::cout << "vec: " << r << " " << sin_a_v << " " << cos_a_v << std::endl;
+    std::cout << "dif: " << r << " " << (sin_a - sin_a_v) / sin_a << " " << (cos_a - cos_a_v) / cos_a << std::endl;
+    std::cout.precision(ss);
 #endif
     double ds = fabs(scalpha.sn - sin_a);
     if(ds > dsin){