Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fix the ear clipping #2743

Open
wants to merge 2 commits into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
9 changes: 5 additions & 4 deletions surface/include/pcl/surface/ear_clipping.h
Original file line number Diff line number Diff line change
Expand Up @@ -83,11 +83,12 @@ namespace pcl
void
triangulate (const Vertices& vertices, PolygonMesh& output);

/** \brief Compute the signed area of a polygon.
* \param[in] vertices the vertices representing the polygon
/** \brief Triangulate one polygon, assume the vertices are clockwise.
* \param[in] vertices the set of vertices
* \param[out] output the resultant polygonal mesh
*/
float
area (const std::vector<uint32_t>& vertices);
size_t
triangulate (std::vector<uint32_t>& vertices, PolygonMesh& output);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@taketwo what's our policy on API changes of non-virtual protected methods? Purist side: still breaks the API since our consumers might have derived classes from this and use this method. Pragmatic side: unlikely.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think generallry we can approach "protected" API stability in a relaxed way, i.e. accept breakage without deprecation cycle. That said, isn't it super easy to keep area() method around and thus definitely not break anyone's code?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I haven't checked yet. I decided to read a little on ear clipping triangulation before having a look with some proper eyes to the PR. I was failing to get the jargon of the OP. ^^

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The code is not used at all, so there's no reason to keep it around. I would remove it.


/** \brief Check if the triangle (u,v,w) is an ear.
* \param[in] u the first triangle vertex
Expand Down
106 changes: 43 additions & 63 deletions surface/src/ear_clipping.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -66,87 +66,66 @@ pcl::EarClipping::performProcessing (PolygonMesh& output)
void
pcl::EarClipping::triangulate (const Vertices& vertices, PolygonMesh& output)
{
const int n_vertices = static_cast<const int> (vertices.vertices.size ());
const size_t n_vertices = vertices.vertices.size ();

if (n_vertices < 3)
return;
else if (n_vertices == 3)
{
output.polygons.push_back( vertices );
output.polygons.emplace_back( vertices );
return;
}

std::vector<uint32_t> remaining_vertices (n_vertices);
if (area (vertices.vertices) > 0) // clockwise?
remaining_vertices = vertices.vertices;
else
for (int v = 0; v < n_vertices; v++)
std::vector<uint32_t> remaining_vertices = vertices.vertices;
size_t count = triangulate(remaining_vertices, output);

// if the input vertices order is anti-clockwise, it always left a
// convex polygon and start infinite loops, which means will left more
// than 3 points.
if(remaining_vertices.size() < 3)return;

output.polygons.erase(output.polygons.cend(), output.polygons.cend() + count);
remaining_vertices.resize(n_vertices);
for (size_t v = 0; v < n_vertices; v++)
remaining_vertices[v] = vertices.vertices[n_vertices - 1 - v];
triangulate(remaining_vertices, output);
}

/////////////////////////////////////////////////////////////////////////////////////////////
size_t
pcl::EarClipping::triangulate (std::vector<uint32_t>& vertices, PolygonMesh& output)
{
// triangles count
size_t count = 0;

// Avoid closed loops.
if (remaining_vertices.front () == remaining_vertices.back ())
remaining_vertices.erase (remaining_vertices.end () - 1);
if (vertices.front () == vertices.back ())
vertices.erase (vertices.end () - 1);

// null_iterations avoids infinite loops if the polygon is not simple.
for (int u = static_cast<int> (remaining_vertices.size ()) - 1, null_iterations = 0;
remaining_vertices.size () > 2 && null_iterations < static_cast<int >(remaining_vertices.size () * 2);
++null_iterations, u = (u+1) % static_cast<int> (remaining_vertices.size ()))
for (int u = static_cast<int> (vertices.size ()) - 1, null_iterations = 0;
vertices.size () > 2 && null_iterations < static_cast<int >(vertices.size () * 2);
++null_iterations, u = (u+1) % static_cast<int> (vertices.size ()))
{
int v = (u + 1) % static_cast<int> (remaining_vertices.size ());
int w = (u + 2) % static_cast<int> (remaining_vertices.size ());
int v = (u + 1) % static_cast<int> (vertices.size ());
int w = (u + 2) % static_cast<int> (vertices.size ());

if (isEar (u, v, w, remaining_vertices))
if (vertices.size() == 3 || isEar (u, v, w, vertices))
{
Vertices triangle;
triangle.vertices.resize (3);
triangle.vertices[0] = remaining_vertices[u];
triangle.vertices[1] = remaining_vertices[v];
triangle.vertices[2] = remaining_vertices[w];
output.polygons.push_back (triangle);
remaining_vertices.erase (remaining_vertices.begin () + v);
triangle.vertices[0] = vertices[u];
triangle.vertices[1] = vertices[v];
triangle.vertices[2] = vertices[w];
output.polygons.emplace_back (triangle);
vertices.erase (vertices.begin () + v);
null_iterations = 0;
count++;
}
}
return count;
}


/////////////////////////////////////////////////////////////////////////////////////////////
float
pcl::EarClipping::area (const std::vector<uint32_t>& vertices)
{
//if the polygon is projected onto the xy-plane, the area of the polygon is determined
//by the trapeze formula of Gauss. However this fails, if the projection is one 'line'.
//Therefore the following implementation determines the area of the flat polygon in 3D-space
//using Stoke's law: http://code.activestate.com/recipes/578276-3d-polygon-area/

int n = static_cast<int> (vertices.size ());
float area = 0.0f;
Eigen::Vector3f prev_p, cur_p;
Eigen::Vector3f total (0,0,0);
Eigen::Vector3f unit_normal;

if (n > 3)
{
for (int prev = n - 1, cur = 0; cur < n; prev = cur++)
{
prev_p = points_->points[vertices[prev]].getVector3fMap();
cur_p = points_->points[vertices[cur]].getVector3fMap();

total += prev_p.cross( cur_p );
}

//unit_normal is unit normal vector of plane defined by the first three points
prev_p = points_->points[vertices[1]].getVector3fMap() - points_->points[vertices[0]].getVector3fMap();
cur_p = points_->points[vertices[2]].getVector3fMap() - points_->points[vertices[0]].getVector3fMap();
unit_normal = (prev_p.cross(cur_p)).normalized();

area = total.dot( unit_normal );
}

return area * 0.5f;
}


/////////////////////////////////////////////////////////////////////////////////////////////
bool
pcl::EarClipping::isEar (int u, int v, int w, const std::vector<uint32_t>& vertices)
Expand All @@ -157,12 +136,13 @@ pcl::EarClipping::isEar (int u, int v, int w, const std::vector<uint32_t>& verti
p_w = points_->points[vertices[w]].getVector3fMap();

const float eps = 1e-15f;
Eigen::Vector3f p_uv, p_uw;
p_uv = p_v - p_u;
p_uw = p_w - p_u;
Eigen::Vector3f p_vu, p_vw;
p_vu = p_u - p_v;
p_vw = p_w - p_v;

// Avoid flat triangles.
if ((p_uv.cross(p_uw)).norm() < eps)
// Avoid flat triangles and concave vertex
Eigen::Vector3f cross = p_vu.cross(p_vw);
if ((cross[2] > 0) || (cross.norm() < eps))
return (false);

Eigen::Vector3f p;
Expand Down
22 changes: 11 additions & 11 deletions test/surface/test_ear_clipping.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -126,23 +126,23 @@ TEST (PCL, EarClippingCubeTest)
vertices.vertices.resize(4);

const int squares[][4] = { {1, 5, 6, 2},
{2, 6, 7, 3},
{3, 7, 4, 0},
{0, 4, 5, 1},
{4, 7, 6, 5},
{0, 1, 2, 3} };
{2, 6, 7, 3},
{3, 7, 4, 0},
{0, 4, 5, 1},
{4, 7, 6, 5},
{0, 1, 2, 3} };

const int truth[][3] = { {2, 1, 5},
{6, 2, 5},
{3, 2, 6},
{3, 2, 6},
{7, 3, 6},
{0, 3, 7},
{0, 3, 7},
{4, 0, 7},
{1, 0, 4},
{1, 0, 4},
{5, 1, 4},
{5, 4, 7},
{6, 5, 7},
{3, 0, 1},
{4, 5, 6},
{7, 4, 6},
{3, 0, 1},
{2, 3, 1} };

PolygonMesh::Ptr mesh (new PolygonMesh);
Expand Down