|
Parallelogram: | ![]() |
![]() |
![]() |
Triangle: | ![]() |
![]() |
However, except in special situations, finding the height of a triangle at an arbitrary orientation usually requires also computing the perpendicular distance of the top vertex from the base.
For example, if one knows the lengths of two sides, a and b, of a triangle and also the angle q between them, then Euclid says this is enough to determine the triangle and its area. Using trigonometry, the height of the triangle over the base b is given by h = a sin q, and thus the area is:
![]() |
![]() |
![]() |
Another frequently used computation is derived from the fact that triangles with equal sides are congruent, and thus have the same area. This observation from Euclid (~300 BC) culminated in Heron's formula (~50 AD) for area as a function of the lengths of its three sides [Note: some historians attribute this result to Archimedes (~250 BC)]; namely:
![]() |
![]() |
![]() |
where a,b,c are the lengths of the sides, and s is the semiperimeter. There are interesting algebraic variations of this formula; such as:
which avoids calculating the 3 square roots to explicitly get the lengths a,b,c from the triangle's vertex coordinates. Other variations on Heron's formula can be found at Kevin Brown's Heron's Formula page and Eric Weisstein's Triangle page.
The remaining classical triangle congruence is when two angles and one side are known. Knowing two angles gives all three, so we can assume the angles q and j are both adjacent to the known base b. Then the formula for area is:
![]() |
![]() |
![]() |
More recently, starting in the 17-th century with Descartes and Fermat, linear algebra produced new simple formulas for area. In 3 dimensional space (3D), the area of a parallelogram and triangle can be expressed as the magnitude of the cross-product of two edge vectors, since |v×w| = |v||w||sin q| where q is the angle between the two vectors v and w. Thus for a 3D triangle with vertices V0V1V2 putting v=V1-V0 and w=V2-V0, one gets:
![]() |
![]() |
![]() |
In 2 dimensional space (2D), a vector can be viewed as embedded in 3D by adding a third component which is set = 0. This lets one take the cross-product of 2D vectors, and use it to compute area. Given a triangle with vertices Vi=(xi,yi)=(xi,yi,0) for i=0,2, we can compute that:
And the absolute value of the third z-component is twice the absolute area of the triangle. However, it is useful to not take the absolute value here, and instead let the area be a signed quantity.
![]() |
![]() |
![]() |
This formula for area is a very efficient computation with no roots or trigonometric functions involved - just 2 multiplications and 5 additions, and possibly 1 division by 2 (which can sometimes be avoided).
Note that the signed area will be positive if V0V1V2 are oriented counterclockwise around the triangle, and will be negative if the triangle is oriented clockwise; and so, this area computation can be used to test for a triangle's orientation. The signed area can also be used to test whether the point V2 is to the left (positive) or the right (negative) of the directed line segment V0V1 . So this is a very useful primitive, and it's great to have such an efficient formula for it.
The Greeks singled out certain quadrilaterals (also called quadrangles) for special treatment, including the square, the rectangle, the parallelogram, and the trapezium. Then, given an arbitrary quadrilateral, they showed how to construct a parallelogram [Euclid, Book I, Prop 45] or square [Euclid, Book II, Prop 14] with an equal area. And the area of a parallelogram was equal to its base times its height. But there was no general formula for the quadrilateral's area.
An extension of Heron's triangle area formula to quadrilaterals was discovered by the Hindu geometer Brahmagupta (~620 AD) [Coxeter,1967, Section 3.2]. However, it only works for cyclic quadrilaterals where all four vertices lie on the same circle. For a cyclic quadrilateral Q, let the lengths of the four sides be a, b, c, d, and the semiperimeter be s=(a+b+c+d)/2. Then, the area of Q is given by:
which is an amazing symmetric formula. If one side is zero length, say d=0, then we have a triangle (which is always cyclic) and this formula reduces to Heron's one.
In modern linear algebra, as already noted, the area of a planar parallelogram is the magnitude of the cross product of two adjacent edge vectors. So, for any 3D planar parallelogram V0V1V2V3, we have:
![]() |
![]() |
![]() |
In 2D, with vertices Vi=(xi,yi)=(xi,yi,0) for i=0,3, this becomes:
which is again a signed area just as we had for triangles.
Next, for an arbitrary quadrilateral, one can compute its area using a parallelogram discovered by Pierre Varignon (first published in 1731). It is amazing that the Greeks missed Varignon's simple result which was discovered 2000 years after Euclid! Given any quadrilateral, one can take the midpoints of its 4 edges to get 4 vertices which form a new quadrilateral. It is then easy to show that this midpoint quadrilateral is always a parallelogram, called the "Varignon parallelogram", and that its area is exactly one-half the area of the original quadrilateral [Coxeter,1967, Section 3.1]. So, for a quadrilateral Q=V0V1V2V3, let this parallelogram have midpoint vertices M0M1M2M3 as shown in the diagram:
From elementary geometry, we know that in triangle V0V1V2 the midpoint line M0M1 is parallel to the base V0V2. In triangle V0V1V2V3, the line M3M2 is parallel to that same base V0V2. Thus, M0M1 and M3M2 are parallel to each other. Similarly, M0M3 and M1M2 are parallel, which shows that M0M1M2M3 is a parallelogram. The area relation is also easy to demonstrate, and we can compute the quadrilateral's area as:
which equals one-half the magnitude of the cross-product of the two diagonals of the quadrilateral. This result was noted by [Van Gelder,1995] who used a different proof.. This formula holds for any 3D planar quadrilateral. When restricted to 2D with Vi=(xi,yi), it becomes a formula for a signed area:
This formula for an arbitrary quadrilateral is just as efficient as the one for an arbitrary triangle, using only 2 multiplications and 5 additions. For simple quadrilaterals, the area is positive when the vertices are oriented counterclockwise, and negative when they are clockwise. However, it also works for nonsimple quadrilaterals and is equal to the difference in area of the two regions the quadrilateral bounds. For example, in the following diagram where I is the self-intersection point of a nonsimple quadrilateral Q=V0V1V2V3, we have:
![]() |
![]() |
![]() |
A 2D polygon can be decomposed into triangles. For computing area, there is a very easy decomposition method for simple polygons (i.e. ones without self intersections). Let a polygon W be defined by its vertices Vi=(xi,yi) for i=0,n with Vn=V0. Also, let P be any point; and for each edge ViVi+1 of W, form the triangle Di=DPViVi+1. Then, the area of W is equal to the sum of the signed areas of all the triangles Di for i=0,n-1; and we have:
![]() |
![]() |
![]() |
Notice that, for a counterclockwise oriented polygon, when the point P is on the "inside" left side of an edge ViVi+1, then the area of Di is positive; whereas, when P is on the "outside" right side of an edge ViVi+1, then Di has a negative area. If instead the polygon is oriented clockwise, then the signs are reversed, and inside triangles become negative.
For example, in the above diagram, the triangles D2=DPV2V3 and Dn-1=DPVn-1Vn have positive area, and contribute positively to the total area of polygon W. However, as one easily observes, only part of D2 and Dn-1 are actually inside W and there is a part of each triangle that is also exterior. On the other hand, the triangles D0 and D1 have negative area and this cancels out the exterior excesses of positive area triangles. In the final analysis, exterior areas all get cancelled, and one is left with exactly the area of the polygon W.
One can make the formula more explicit by picking a specific point P and expanding the terms. By selecting P=(0,0), the area of each triangle reduces to 2A(Di) = (xiyi+1-xi+1yi). This yields:
![]() |
![]() |
![]() |
A little algebra shows that the second and
third summations are equal to the first. For a polygon with n vertices, the first summation uses 2n multiplications and (2n-1) additions; the second uses n multiplications and (3n-1) additions; and the third uses only n multiplications and (2n-1) additions. So, the third is
preferred for efficiency, but to avoid any overhead from computing the
index i(mod n), one must extend
the polygon array up to Vn+1=V1.
This computation gives a signed area for a polygon; and, similar to the signed area of a triangle, is positive when the vertices are oriented counterclockwise around the polygon, and negative when oriented clockwise. So, this computation can be used to test for a polygon's global orientation. However, there are other more efficient algorithms for determining polygon orientation. The easiest is to find the rightmost lowest vertex of the polygon, and then test the orientation of the entering and leaving edges at this vertex. This test can be made by checking if the end vertex of the leaving edge is to the left of the entering edge, which means that the orientation is counterclockwise.
An important generalization is for planar polygons embedded in 3D space [Goldman, 1994]. We have already shown that the area of a 3D triangle D=V0V1V2 is given by half the magnitude of the cross product of two of its edge vectors; namely, ½ |V0V1 × V0V2|.
There is a classic standard formula for the area of a 3D polygon [Goldman, 1994] that extends the cross-product formula for a triangle. It can be derived from Stokes Theorem. However, we show here how to derive it from a 3D triangular decomposition that is geometrically more intuitive.
A general 3D planar polygon W has vertices Vi=(xi,yi,zi)
for i=0,n, with Vn=V0. All
the vertices lie on the same 3D plane p which has a unit normal vector
n. Now, the same as in the 2D
case, let P be any 3D point (not
generally on the plane p);
and for each edge ei=ViVi+1
of W, form the 3D triangle
Di=DPViVi+1. We
would like to relate the sum of the areas of all these triangles to the
area of the polygon W in the
plane p. But what we
have is a pyramidal cone with P as an apex over
the polygon W as a
base. We need to project the triangular sides of this cone onto the
plane p of the base polygon,
and compute signed areas of the projected triangles. If we can do
this, then the sum of the projected areas will equal the total area of the
planar polygon.
To achieve this, start by associating to
each triangle Di an area vector ai=½(PVi ×
PVi+1), perpendicular to Di,
whose magnitude we know is equal to that triangle's area.
Next, drop a perpendicular from P to a point
P0 on p, and consider the projected triangle Ti=DP0ViVi+1.
Drop a perpendicular P0Bi from P0 to Bi on the edge ei=ViVi+1.
Since PP0 is also perpendicular to
ei, the three points
PP0Bi define a
plane normal to ei,
and PBi is a perpendicular
from P to the edge ei. Thus, |PBi| is the height of Di, and |P0Bi| is the height of
Ti. Further, the angle between
these two altitudes = q = the angle
between n and ai since a
90° rotation results in congruence. This gives:
![]() |
This signed area computation is positive if the vertices of Ti are oriented counterclockwise when we look at the plane p from the side pointed to by n. The same as in the 2D case, we can now add together the signed areas of all the triangles Ti to get the area of the polygon W. Writing this down, we have:
![]() |
Finally, by selecting P=(0,0,0), we have PVi=Vi and this yields the simple formula:
![]() |
![]() |
![]() |
which uses 6n+3 multiplications and 4n+2 additions
Similar to the 2D case, this is a signed area which is positive when the vertices are oriented counterclockwise around the polygon when viewed from the side of p pointed to by n.
[Van Gelder, 1995] has shown how to significantly speed up this computation by using a decomposition into quadrilaterals instead of triangles. He first noted that the area of a 3D planar quadrilateral Q=V0V1V2V3 can be computed in terms of the cross-product of its diagonals; that is:
which reduces four expensive cross-product computations to just one!
Then, any large polygon W with n>4 vertices can be decomposed into quadrilaterals formed by V0 and three other sequential vertices V2i-1, V2i, and V2i+1 for i=1,h where h= the greatest integer £(n-1)/2. For n odd, the decomposition ends with a triangle. This gives:
where k=0 for n odd, and k=n-1 for n even. This formula reduces the number of expensive cross-products by a factor of two (replacing them with vector subtractions). In total there are 3n+3 multiplications and 5n+1 additions making this formula roughly twice as fast as the classic standard one.
[Van Gelder, 1995] also states that this method can be applied to 2D polygons, but he does not write down the details. Working this out produces a formula that uses n multiplications and 3n-1 additions, which is not as fast as the prior formulas we have given. We simply note this here, and do not pursue it further.
One can speed up the computation of 3D planar polygon area by projecting the polygon onto a 2D plane [Snyder & Barr, 1987]. Then, the area can be computed in 2D using our fastest formula, and the 3D area is recovered using an area scaling factor. This method is implemented by projecting onto an axis-aligned plane by ignoring one of the three coordinates. To avoid degeneracy and optimize robustness, we look at the plane's normal vector n (see the April 2001 Algorithm About Planes), and choose the component with the greatest absolute as the one to ignore. Let Projc() be the projection that ignores the coordinate c = x, y, or z. Then, the ratio of areas for the projected polygon Projc(W) and original planar polygon W with normal n = (nx,ny,nz) is:
Thus, the 3D planar area can be recovered by a single extra multiplication, and in total this algorithm uses n+5 multiplications, 2n+1 additions, 1 square root (when n is not a unit normal), plus a small overhead choosing the coordinate to ignore. This is a very significant improvement of the classic standard formula, achieving almost a 6 times speed-up. We give an efficient implementation below in the routine area3D_Polygon().
Here are some sample "C++" implementations of these formulas as algorithms. We just give the 2D case with integer coordinates, and use the simplest structures for a point, a triangle, and a polygon which may differ in your application. We represent a polygon as an array of points, but it is often more convenient to have it as a linked list of vertices (to allow insertion or deletion during drawing operations), and the polygon routines can be easily modified to scan through the linked list (see [O'Rourke, 1998] for an example of this approach).
// Copyright 2000,
softSurfer (www.softsurfer.com)
// This code may be freely used and
modified for any purpose
// providing that this copyright notice is
included with it.
// SoftSurfer makes no warranty for this code, and
cannot be held
// liable for any real or imagined damage resulting from
its use.
// Users of this code must verify correctness for their
application.
// a Point (or
vector) is defined by its coordinates
typedef struct {int x, y, z;}
Point; // exclude z for 2D
// a Triangle is
given by three points: Point V0, V1, V2
// a Polygon is given
by:
// int n = number of
vertex points
// Point* V[] =
an array of points with V[n]=V[0], V[n+1]=V[1]
// Note: for efficiency low-level functions are declared to be inline.
//
isLeft(): tests if a point is Left|On|Right of an infinite
line.
// Input: three points P0, P1, and
P2
// Return: >0 for P2 left of the line through
P0 and
P1
//
=0 for P2 on the
line
//
<0 for P2 right of the line
inline int
isLeft( Point P0, Point
P1, Point P2 )
{
return ( (P1.x - P0.x) * (P2.y -
P0.y)
- (P2.x - P0.x) * (P1.y - P0.y)
);
}
//===================================================================
// orientation2D_Triangle(): test
the orientation of a triangle
// Input: three
vertex points V0, V1, V2
// Return: >0 for
counterclockwise
//
=0 for none
(degenerate)
//
<0 for clockwise
inline int
orientation2D_Triangle( Point V0,
Point V1, Point V2 )
{
return isLeft(V0, V1,
V2);
}
//===================================================================
// area2D_Triangle(): compute the area of
a triangle
// Input: three vertex points V0,
V1, V2
// Return: the (float) area of T
inline
float
area2D_Triangle( Point V0, Point V1, Point V2
)
{
return (float)isLeft(V0, V1, V2) /
2.0;
}
//===================================================================
// orientation2D_Polygon(): tests
the orientation of a simple polygon
// Input:
int n = the number of vertices in the
polygon
//
Point* V = an array of n+1 vertices with V[n]=V[0]
//
Return: >0 for counterclockwise
//
=0 for none
(degenerate)
//
<0 for clockwise
// Note: this algorithm is faster
than computing the signed area.
int
orientation2D_Polygon( int n,
Point* V )
{
// first find rightmost lowest
vertex of the polygon
int rmin =
0;
int xmin = V[0].x;
int ymin
= V[0].y;
for (int i=1; i<n; i++)
{
if (V[i].y >
ymin)
continue;
if (V[i].y ==
ymin) { // just as
low
if (V[i].x < xmin) // and to
left
continue;
}
rmin =
i; // a new
rightmost lowest vertex
xmin
= V[i].x;
ymin =
V[i].y;
}
// test
orientation at this rmin vertex
// ccw <=> the
edge leaving is left of the entering edge
if (rmin
== 0)
return isLeft( V[n-1],
V[0], V[1] );
else
return isLeft(
V[rmin-1], V[rmin], V[rmin+1] );
}
//===================================================================
// area2D_Polygon(): computes the area of a
2D polygon
// Input: int n = the number of
vertices in the
polygon
//
Point* V = an array of n+2 vertices
//
with V[n]=V[0] and V[n+1]=V[1]
// Return: the (float)
area of the polygon
float
area2D_Polygon( int n, Point* V
)
{
float area = 0;
int i, j, k; //
indices
for (i=1, j=2, k=0; i<=n; i++, j++,
k++) {
area += V[i].x *
(V[j].y - V[k].y);
}
return
area /
2.0;
}
//===================================================================
// area3D_Polygon(): computes the area of a
3D planar polygon
// Input: int n = the number
of vertices in the
polygon
//
Point* V = an array of n+2 vertices in a
plane
//
with V[n]=V[0] and
V[n+1]=V[1]
//
Point N = unit normal vector of the polygon's
plane
// Return: the (float) area of the
polygon
float
area3D_Polygon( int n, Point* V, Point N
)
{
float area = 0;
float
an, ax, ay, az; // abs value of normal and its
coords
int
coord; //
coord to ignore: 1=x, 2=y, 3=z
int i, j,
k; // loop
indices
// select largest abs coordinate to
ignore for projection
ax = (N.x>0 ? N.x :
-N.x); // abs x-coord
ay =
(N.y>0 ? N.y : -N.y); // abs
y-coord
az = (N.z>0 ? N.z :
-N.z); // abs z-coord
coord =
3;
// ignore z-coord
if (ax > ay)
{
if (ax > az) coord =
1; // ignore x-coord
}
else if (ay > az) coord = 2; //
ignore y-coord
// compute area of the 2D
projection
for (i=1, j=2, k=0; i<=n; i++, j++,
k++)
switch (coord)
{
case
1:
area += (V[i].y * (V[j].z -
V[k].z));
continue;
case
2:
area += (V[i].x * (V[j].z -
V[k].z));
continue;
case
3:
area += (V[i].x * (V[j].y -
V[k].y));
continue;
}
// scale to get area before
projection
an = sqrt( ax*ax + ay*ay + az*az);
// length of normal vector
switch (coord)
{
case
1:
area *= (an /
(2*ax));
break;
case
2:
area *= (an /
(2*ay));
break;
case
3:
area *= (an /
(2*az));
}
return
area;
}
//===================================================================
Kevin Brown, MathPages : Heron's Formula
Donald Coxeter & Samuel Greitzer, Geometry Revisited (1967)
Ronald Goldman, "Area of Planar Polygons and Volume of Polyhedra" in Graphics Gems II (1994)
Joseph O'Rourke, Computational Geometry in C (2nd Edition), Section 1.3 "Area of a Polygon" (1998)
J.P. Snyder and A.H. Barr, "Ray Tracing Complex Models Containing Surface Tessellations", ACM Comp Graphics 21, (1987)
Allen Van Gelder, "Efficient Computation of Polygon Area and Polyhedron Volume" in Graphics Gems V (1995)
Eric Weisstein, MathWorld Site: Triangle or in The CRC Concise Encyclopedia of Mathematics (1998)
Copyright ©
2001-2002 softSurfer. All rights reserved. |