Polygon Triangulation
Polygon triangulation is the process of dividing a polygon into simple triangles. It is one of the fundamental algorithms in computational geometry. Polygons generally are complicated shapes and not always convex. This alone makes it difficult to draw a polygon in say OpenGL. Triangles on the other hand are always convex, very easy to draw and can be used in many other algorithms in computational geometry.
In this tutorial we will develop a simple algorithm to divide a polygon into a series of triangles. Since this is not a mathematics course we will not deal with the proof of the theorems used. At the end we will develop the code that performs the triangulation.
Introduction
Let's start with some definitions.
Polygon is an ordered set of points called vertices. The lines connecting the points are called edges.
Triangle is actually a polygon consisting only of three points.
Convex polygon is the polygon whose any edge does not intersect with the rest of the outline. In the contrary in concave polygon there are edges that partition the polygon in two or more parts as we can see in the image.
Over the years it was proved that all polygons can be triangulated. Many algorithms were invented to solve the problem. One of the simplest algorithms is presented in this article.
One of the first problems that people tried to solve by triangulating a polygon was the so called art gallery problem. The intention was to have as few cameras as possible that can cover all the area of an art gallery. We will not try to present the solution to this problem here. The main problem we are trying to solve is the rendering of concave polygons in OpenGL. They are very common in everyday images and we need an algorithm to decompose them into simpler shapes that the renderer can handle.
In this first article we will deal with simple polygons like the ones shown in the image. In a future article we will address the problem of polygons with holes.
The algorithm
The algorithm developed here is quite simple. The triangles generated might be bad for numerical analysis such as finite elements. Our goal though is served fine. We want something simple that anyone can understand that will decompose any polygon so that OpenGL can render it afterwards.
How it works
The algorithm is based on the fact that every polygon can be triangulated. So there must be at least one valid diagonal that lies inside it and generates a triangle. We construct a new polygon which can be rotated, the first point can go to the end and the second can be the new first, so that we always check the first three points. This way it is less complicated and the code is easier to follow. If these three points form a valid triangle then we store it and remove the SECOND point from the polygon. This creates a polygon with one less point and we start all over again until the original polygon becomes a triangle. If the triangle is invalid we rotate the polygon sending the first point to the end of the point chain and start again.
Mathematicians many years before us proved that this approach is bound to succeed by proving that every polygon can be triangulated using only diagonals. A simple search in the network can reveal many solutions to the triangulation problem and of course all the mathematical evidence you might seek.
Here is the code for the triangulation algorithm. It relies on the previous tutorial about point in polygon. We create a new version of intersection testing which ignores intersection if it occurs on the endpoints of the line segments because it is always the case in the case of a diagonal of a polygon.
The algorithm developed here does not cover polygns with holes. In a future article we will address this issue as well.
class clf_triangle2D {
public:
clf_triangle2D(){}
clf_triangle2D(const clf_triangle2D& t){
verts[0] = t.verts[0];
verts[1] = t.verts[1];
verts[2] = t.verts[2];
}
~clf_triangle2D(){}
clf_point2D verts[3];
};
// check for intersection, but exclude the vertices
bool intersect_no_vertex(const clf_point2D& p1, const clf_point2D& p2, const clf_point2D& p3, const clf_point2D& p4, clf_point2D* p0)
{
// first check horizontal and vertical gap
if (min(p1.x, p2.x)>max(p3.x,p4.x))
return false;
if (max(p1.x, p2.x)<min(p3.x,p4.x))
return false;
if (min(p1.y, p2.y)>max(p3.y,p4.y))
return false;
if (max(p1.y, p2.y)<min(p3.y,p4.y))
return false;
// first calculate the common denominator
float d = (p4.y-p3.y)*(p2.x-p1.x)-(p4.x-p3.x)*(p2.y-p1.y);
if (fabs(d) < 0.0000001) // almost 0 (compensate for computer accuracy)
return false;
// calculate fractional part for first segment
float na = ((p4.x-p3.x)*(p1.y-p3.y)-(p4.y-p3.y)*(p1.x-p3.x))/d;
if (na < 0.0000001f || na > 0.9999999f) // beyond limits of first segment, search no more
return false;
// calculate fractional part for second segment
float nb = ((p2.x-p1.x)*(p1.y-p3.y)-(p2.y-p1.y)*(p1.x-p3.x))/d;
if (nb < 0.0000001f || nb > 0.9999999f) // beyond limits of second segment, search no more
return false;
// we have intersection, calculate actual coordinates (if required)
if (p0)
{
p0->x = p1.x+na*(p2.x-p1.x);
p0->y = p1.y+na*(p2.y-p1.y);
// or using the second segment endpoints
// p0->x = p3.x+nb*(p4.x-p3.x);
// p0->y = p3.y+nb*(p4.y-p3.y);
}
return true;
}
// check if the line segment (p1,p2) is inside a polygon
bool line_in_polygon(clf_polygon2D& poly, clf_point2D& p1, clf_point2D& p2)
{
// the line that connects the first and the third point is entirely inside thepolygon
clf_point2D pt((p1.x+p2.x)/2, (p1.y+p2.y)/2); // center point
if (!point_in_polygon(poly, pt)) // if it is outside the polygon
return false; // triangle is invalid
// iterate through the edges of the polygon
clf_polygon2D::iterator next = poly.begin();
clf_polygon2D::iterator it;
int count = 0;
for (it=poly.begin() ; it!=poly.end(); ++it)
{
// -next- points to the 'next' point while -it- is the current
++next;
if (next == poly.end())
next = poly.begin();
// we know that the vertices are on the polygon!!!
if (intersect_no_vertex(*it, *next, p1, p2, NULL))
{
return false;
}
}
return true;
}
bool triangulate(clf_polygon2D& poly, std::vector<clf_triangle2D>& result)
{
// clear the output buffer
result.erase(result.begin(), result.end());
// generate index buffer, instead of the points we keep their respective indices
std::list<int> tp; // valid indices
std::vector<clf_point2D>::iterator it;
int i=0;
for (it=poly.begin(); it!=poly.end(); ++it)
{
tp.push_back(i);
i++;
}
// and here we go
clf_triangle2D t;
std::list<int>::iterator ti;
while (tp.size() > 3) // while we have a polygon
{
ti = tp.begin(); // go to its start
// and get the first three points
int i1 = *ti; ++ti;
int i2 = *ti; ++ti;
int i3 = *ti; ++ti;
// check if they form a valid triangle
// the line that connects the first and the third point is entirely inside the polygon
bool valid = false; // assume invalid triangle
// if so, does the line intersect with the polygon?
if (line_in_polygon(poly, poly[i1], poly[i3]))
{
// if it is inside we have a valid triangle
t.verts[0] = poly[i1];
t.verts[1] = poly[i2];
t.verts[2] = poly[i3];
result.push_back(t); // store it
valid = true; // and prepare for the next
}
ti = tp.begin(); // back to the start of the index buffer
if (valid) // if the triangle was valid move to the second point
++ti;
int tpi = *ti; // get the index of the current point (first or second)
tp.erase(ti); // remove the point from the buffer
if (!valid) // if the triangle was invalid
tp.push_back(tpi); // store the first point at the end of the point buffer
}
// now store the remaining triangle
{
ti = tp.begin();
int i1 = *ti; ++ti;
int i2 = *ti; ++ti;
int i3 = *ti; ++ti;
t.verts[0] = poly[i1];
t.verts[1] = poly[i2];
t.verts[2] = poly[i3];
result.push_back(t);
}
return true;
}