better templates for marching tets (#2279) [ci skip]

* better templates

* fabs -> abs

* missing cast

include/igl/marching_tets.cpp

@@ -15,6 +15,7 @@
 #include <iostream>
 #include <algorithm>
 #include <utility>
+#include <cmath>
 
 template <typename DerivedTV,
           typename DerivedTT,
@@ -27,7 +28,7 @@ void igl::marching_tets(
     const Eigen::MatrixBase<DerivedTV>& TV,
     const Eigen::MatrixBase<DerivedTT>& TT,
     const Eigen::MatrixBase<DerivedS>& isovals,
-    double isovalue,
+    const typename DerivedS::Scalar isovalue,
     Eigen::PlainObjectBase<DerivedSV>& outV,
     Eigen::PlainObjectBase<DerivedSF>& outF,
     Eigen::PlainObjectBase<DerivedJ>& J,
@@ -163,11 +164,12 @@ void igl::marching_tets(
       {
         // Typedef to make sure we handle floats properly
         typedef Eigen::Matrix<typename DerivedTV::Scalar, 1, 3, Eigen::RowMajor, 1, 3> RowVector;
-        const RowVector v1 = TV.row(edge.first);
-        const RowVector v2 = TV.row(edge.second);
-        const double a = fabs(isovals(edge.first, 0) - isovalue);
-        const double b = fabs(isovals(edge.second, 0) - isovalue);
-        const double w = a / (a+b);
+        using Scalar = typename DerivedS::Scalar;
+        const RowVector v1 =  TV.row(edge.first).template cast<Scalar>();
+        const RowVector v2 = TV.row(edge.second).template cast<Scalar>();
+        const Scalar a = abs(isovals(edge.first, 0) - isovalue);
+        const Scalar b = abs(isovals(edge.second, 0) - isovalue);
+        const Scalar w = a / (a+b);
 
         // Create a casted copy in case BCType is a float and we need to downcast
         const BCType bc_w = static_cast<BCType>(w);
@@ -193,5 +195,5 @@ void igl::marching_tets(
 
 
 #ifdef IGL_STATIC_LIBRARY
-template void igl::marching_tets<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, 1, 0, -1, 1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, double>(Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> > const&, double, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&, Eigen::SparseMatrix<double, 0, int>&);
+template void igl::marching_tets<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, 1, 0, -1, 1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, double>(Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> > const&, const double, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&, Eigen::SparseMatrix<double, 0, int>&);
 #endif // IGL_STATIC_LIBRARY

include/igl/marching_tets.h

@@ -38,7 +38,7 @@ namespace igl {
       const Eigen::MatrixBase<DerivedTV>& TV,
       const Eigen::MatrixBase<DerivedTT>& TT,
       const Eigen::MatrixBase<DerivedS>& S,
-      double isovalue,
+      const typename DerivedS::Scalar isovalue,
       Eigen::PlainObjectBase<DerivedSV>& SV,
       Eigen::PlainObjectBase<DerivedSF>& SF,
       Eigen::PlainObjectBase<DerivedJ>& J,
@@ -73,11 +73,11 @@ namespace igl {
       const Eigen::MatrixBase<DerivedTV>& TV,
       const Eigen::MatrixBase<DerivedTT>& TT,
       const Eigen::MatrixBase<DerivedS>& S,
-      double isovalue,
+      const typename DerivedS::Scalar isovalue,
       Eigen::PlainObjectBase<DerivedSV>& SV,
       Eigen::PlainObjectBase<DerivedSF>& SF,
       Eigen::PlainObjectBase<DerivedJ>& J) {
-    Eigen::SparseMatrix<double> _BC;
+    Eigen::SparseMatrix<typename DerivedSV::Scalar> _BC;
     return igl::marching_tets(TV, TT, S, isovalue, SV, SF, J, _BC);
   }
   /// \overload
@@ -91,7 +91,7 @@ namespace igl {
       const Eigen::MatrixBase<DerivedTV>& TV,
       const Eigen::MatrixBase<DerivedTT>& TT,
       const Eigen::MatrixBase<DerivedS>& S,
-      double isovalue,
+      const typename DerivedS::Scalar isovalue,
       Eigen::PlainObjectBase<DerivedSV>& SV,
       Eigen::PlainObjectBase<DerivedSF>& SF,
       Eigen::SparseMatrix<BCType>& BC) {
@@ -108,11 +108,11 @@ namespace igl {
       const Eigen::MatrixBase<DerivedTV>& TV,
       const Eigen::MatrixBase<DerivedTT>& TT,
       const Eigen::MatrixBase<DerivedS>& S,
-      double isovalue,
+      const typename DerivedS::Scalar isovalue,
       Eigen::PlainObjectBase<DerivedSV>& SV,
       Eigen::PlainObjectBase<DerivedSF>& SF) {
     Eigen::VectorXi _J;
-    Eigen::SparseMatrix<double> _BC;
+    Eigen::SparseMatrix<typename DerivedSV::Scalar> _BC;
     return igl::marching_tets(TV, TT, S, isovalue, SV, SF, _J, _BC);
   }