Fix the order of principal curvature (#1767)

Co-authored-by: lchaooo <[email protected]>

include/igl/principal_curvature.cpp

@@ -40,8 +40,8 @@ class CurvatureCalculator
 {
 public:
   /* Row number i represents the i-th vertex, whose columns are:
-   curv[i][0] : K2
-   curv[i][1] : K1
+   curv[i][0] : K1
+   curv[i][1] : K2
    curvDir[i][0] : PD1
    curvDir[i][1] : PD2
    */
@@ -425,20 +425,20 @@ IGL_INLINE void CurvatureCalculator::finalEigenStuff(int i, const std::vector<Ei
   if (c_val[0] > c_val[1])
   {
     curv[i]=std::vector<double>(2);
-    curv[i][0]=c_val(1);
-    curv[i][1]=c_val(0);
+    curv[i][0]=c_val(0);
+    curv[i][1]=c_val(1);
     curvDir[i]=std::vector<Eigen::Vector3d>(2);
-    curvDir[i][0]=v2global;
-    curvDir[i][1]=v1global;
+    curvDir[i][0]=v1global;
+    curvDir[i][1]=v2global;
   }
   else
   {
     curv[i]=std::vector<double>(2);
-    curv[i][0]=c_val(0);
-    curv[i][1]=c_val(1);
+    curv[i][0]=c_val(1);
+    curv[i][1]=c_val(0);
     curvDir[i]=std::vector<Eigen::Vector3d>(2);
-    curvDir[i][0]=v1global;
-    curvDir[i][1]=v2global;
+    curvDir[i][0]=v2global;
+    curvDir[i][1]=v1global;
   }
   // ---- end Eigen stuff
 }

tests/include/igl/principal_curvature.cpp

@@ -0,0 +1,29 @@
+#include <test_common.h>
+#include <igl/principal_curvature.h>
+#include <igl/cylinder.h>
+
+TEST_CASE("principal_curvature: cylinder", "[igl]")
+{
+  using namespace igl;
+  const int axis_devisions = 20;
+  const int height_devisions = 20;
+  Eigen::MatrixXd V;
+  Eigen::MatrixXi F;
+  cylinder(axis_devisions,height_devisions,V,F);
+  Eigen::MatrixXd PD1,PD2;
+  Eigen::VectorXd PV1,PV2;
+  //PV1: maximal curvature value for each vertex.
+  //PV2: minimal curvature value for each vertex.
+  igl::principal_curvature(V,F,PD1,PD2,PV1,PV2);
+  REQUIRE (PD1.rows() == V.rows());
+  REQUIRE (PD2.rows() == V.rows()); 
+  REQUIRE (PD1.cols() == 3);
+  REQUIRE (PD2.cols() == 3);
+  REQUIRE (PV1.size() == V.rows());
+  REQUIRE (PV2.size() == V.rows());
+  for(int i = 0; i<PV1.size(); ++i)
+  {
+    //max curvature is greater than or equal to min curvature
+    REQUIRE (PV1[i]>=PV2[i]);
+  }
+}

tutorial/203_CurvatureDirections/main.cpp

@@ -51,12 +51,12 @@ int main(int argc, char *argv[])
   // Average edge length for sizing
   const double avg = igl::avg_edge_length(V,F);
 
-  // Draw a blue segment parallel to the minimal curvature direction
+  // Draw a red segment parallel to the maximal curvature direction
   const RowVector3d red(0.8,0.2,0.2),blue(0.2,0.2,0.8);
-  viewer.data().add_edges(V + PD1*avg, V - PD1*avg, blue);
+  viewer.data().add_edges(V + PD1*avg, V - PD1*avg, red);
 
-  // Draw a red segment parallel to the maximal curvature direction
-  viewer.data().add_edges(V + PD2*avg, V - PD2*avg, red);
+  // Draw a blue segment parallel to the minimal curvature direction
+  viewer.data().add_edges(V + PD2*avg, V - PD2*avg, blue);
 
   // Hide wireframe
   viewer.data().show_lines = false;