Register Now
Member Count: 627,818 - April 24, 2014  [Get Time]
Login
forums   
Round Tables
News Discussions
Algorithm Matches
Marathon Matches
NASA Tournament Lab
Software Forums
TopCoder Cookbook
High School Matches
Sponsor Discussions
Search Watch Thread  |  My Post History  |  My Watches  |  User Settings
View: Flat (newest first)  | Threaded  | Tree
Previous Thread  |  Next Thread
Forums TopCoder Cookbook Algorithm Competitions - New Recipes Working with Polygons
Working with Polygons | Reply
Name
Working with Polygons

Problem
This recipe covers some common tasks encounter when working with polygons.

Solution
Here we consider only simple polygons, i.e. without self-crossings and self-contacts. However, non-simple polygons encounter relatively rare.
So, let’s look at concrete tasks.

How to calculate an area of a polygon?
There is a formula for calculating the area of a triangle given its three vertices based on vector multiplication. It’s this: SABC = 1/2 * (xA(yB – yC) + xB(yC – yA) + xC(yA – yB)). Strictly speaking, this formula is for calculating oriented area, i.e. it will be positive if vertices A, B, C are located counterclockwise and negative if clockwise. If you don’t know their orientation and you don’t need to know it, just take an absolute value of this area. Using this formula we can obtain the area of any polygon considering it as a union of triangles. So SP1…PN = SP1P2P3 + SP1P3P4 + … + SP1PN-1PN. Note that this formula works correctly even with non-convex polygons, subtracting areas of some triangles automatically from a total sum.

How to check is a polygon convex?
That’s performed using the same formula for area calculating. A polygon P1 P2 … PN is convex if all oriented areas SP1P2P3, SP2P3P4, …, SPN-1PNP1, SPNP1P2 are of the same sign, or, speaking more strictly, all nonnegative oriented areas among these are of the same sign. If we know that all vertices of the polygon are distinct, P0ero value of some area SP(i), P(i mod N + 1), P(i mod N + 2) tells us that an angle P(i mod N + 1) is equal to 180°.

How to check is a point inside of a polygon, outside of it, or on its boundary?
First of all, to check is the point lie on the boundary of the polygon that’s sufficient to check consequently is it lie on any of segments forming the boundary. For this we refer to recipe about dealing with lines and segments.
Second, if the polygon is convex, we can check if all triangles P0P1P2, …, P0PN-1PN, P0PNP1 have the same sign of oriented area. If it is true, this means the point P0 lies inside the polygon, if false – outside.
And, third, in the case when we don’t know is the polygon convex, we should use a somewhat more complicated algorithm. Let’s consider a line parallel to OX axis, on which point P0 being under consideration lies. For every segment of the boundary, find a point of intersection of this segment and that line (if any) and count the number of these points lying to the right from point P0. If this number is odd, the point lies inside the polygon. But there are some disgusting cases, with which we should treat carefully. First, what if a segment lies on this line and lies to the right of the point? OK, just skip this segment. Second, what if this line goes through a vertex of a polygon? Well, here we should consider, for example, only upper ends of segments, but not lower – and for doing this we need a special test routine, because this discrimination is not the same as discrimination of “starting” and “ending” ends of segments in the order of a roundabout way of the polygon’s perimeter. On the other hand, we may get rid of taking into account these cases, if we can draw not an axis-parallel line, but a line that is guaranteed not to be parallel to any one of polygon’s sides. For example, if we deal with points with integer coordinates, we can draw a line through P0 and a point with coordinates (xP0 + c; yP0 + c +1) where c is a constant bigger than maximum difference of any coordinates between any two points given in the problem.
Re: Working with Polygons (response to post by Ferlon) | Reply
1. Do we really need to place the formula for triangle area here? It would be better to write it as half of cross product explicitly like [(b-a) x (c-a)] / 2. The formula of course does less arithmetic operations, but it is not a significant difference for algorithm competition. On the other hand, formula involving cross product is easier to remember and understand.

2. The example with determining the place for the point and polygon is awful. There is a lot of "redundant code" there that freaks people away from it. You can delete/inline manually 60% of all those tiny functions.
Moreover, the algorithm you describe here is something I'm strongly against because it involves tricky special cases.
There is another way: given point (X, Y) and a polygon create a segment with end points (X, Y) --- (X+10000, Y+10001). The second point is considered to be infinitely far away. Since vector has components with gcd=1, the segment cannot intersect any of the polygon vertices. So you can trace this segment without handling any special cases.
Re: Working with Polygons (response to post by Ferlon) | Reply
Discussion
An example of finding area is ConvexPolygon problem. Just do it. Here is the code:
#include <vector>
#include <cmath>
using namespace std;
struct point {
	double x,y;
	point (double xx=0, double yy=0) {
		x=xx;
		y=yy;
	}
};
double sq(const point &a, const point &b, const point &c) {
	return (a.x*(b.y-c.y)+b.x*(c.y-a.y)+c.x*(a.y-b.y))/2;
}
vector<point> a;
class ConvexPolygon {
public:
	double findArea(vector<int> x, vector<int> y) {
		int i;
		int n=x.size();
		for (i=0; i<n; ++i) a.push_back(point(x[i],y[i]));
		double s=0;
		for (i=1; i<=n-2; ++i) s+=sq(a[0],a[i],a[i+1]);
		s=fabs(s);		
		return s;
	}
};

Another example of a very similar problem is Surveyor.

An example of checking is a point inside, outside or on a boundary of a polygon is PointInPolygon problem. Because all edges of the polygon here are parallel to coordinate axes, this problem can be solved in simpler way, but we forget about this favorable circumstance to demonstrate the general method. Also we demonstrate using of some geometrical routines. Take your attention that we don’t use a secure comparison of doubles in cases when we know in advance they differ more than our “epsilon”, our allowable absolute error. Here is the code:
#include <vector>
#include <string>
#include <sstream>
#include <cmath>
using namespace std;
const double eps=1e-6;
struct point {
	double x,y;
	point (double xx=0, double yy=0) {
		x=xx;
		y=yy;
	}
	point (const string &s) {
		stringstream ss(s);
		double x,y;
		ss >> x >> y;
		*this=point(x,y);
	}
};
struct line {
	double a,b,c;
};
bool eq(double a, double b) {
	return fabs(a-b)<=eps;
}
bool atline(const point &a, const line &l) {
	return eq(l.a*a.x+l.b*a.y+l.c,0);
}
void point2toline(const point &p1, const point &p2, line &l) {
	l.a=p1.y-p2.y;
	l.b=p2.x-p1.x;
	l.c=p1.x*p2.y-p2.x*p1.y;
}
bool atline(const point &a, const point &p1, const point &p2) {
	line l;
	point2toline(p1,p2,l);
	return atline(a,l);
}
bool lseq(double a, double b) {
	return b-a>=-eps;
}
bool ls(double a, double b) {
	return b-a>eps;
}
bool lseq(const point &a, const point &b) {
	if (!eq(a.x,b.x)) return a.x<b.x;
	return lseq(a.y,b.y);
}
bool ls(const point &a, const point &b) {
	if (!eq(a.x,b.x)) return a.x<b.x;
	return ls(a.y,b.y);
}
bool inrectangle(const point &a, const point &lf, const point &rg) {
	return lseq(lf,a) && lseq(a,rg) || lseq(rg,a) && lseq(a,lf);
}
bool atsegm(const point &a, const point &lf, const point &rg) {
	return atline(a,lf,rg) && inrectangle(a,lf,rg);
}
bool eq(const point &a, const point &b) {
	return eq(a.x,b.x) && eq(a.y,b.y);
}
bool crossline(const line &l1, const line &l2) {
	return !eq(l1.a*l2.b-l2.a*l1.b,0);
}
void line2topoint(const line &l1, const line &l2, point &p) {
	double d=l1.a*l2.b-l1.b*l2.a;
	p.x=(l1.b*l2.c-l2.b*l1.c)/d;
	p.y=(l1.c*l2.a-l2.c*l1.a)/d;
}
point p;
vector<point> a;
int n;
class PointInPolygon {
public:
	string testPoint (vector<string> vs, int tx, int ty) {
		n=vs.size();
		p=point(tx,ty);
		int i;
		for (i=0; i<n; ++i) a.push_back(point(vs[i]));
		for (i=0; i<n; ++i) if (atsegm(p,a[i],a[(i+1)%n])) return "BOUNDARY";
		line l;
		l.a=0;
		l.b=1;
		l.c=-p.y;
		point pt;
		int cnt=0;
		line ln;
		for (i=0; i<n; ++i) {
			point2toline(a[i],a[(i+1)%n],ln);
			if (crossline(l,ln)) {
				line2topoint(l,ln,pt);
				if (p.x<pt.x && inrectangle(pt,a[i],a[(i+1)%n])) {
					if (!eq(pt,a[i]) && !eq(pt,a[(i+1)%n])) ++cnt; else {
						if (eq(pt,a[i]) && ls(a[(i+1)%n].y,a[i].y) ||
							eq(pt,a[(i+1)%n]) && ls(a[i].y,a[(i+1)%n].y)) ++cnt;
					}
				}
			}
		}
		if (cnt&1) return "INTERIOR"; else return "EXTERIOR";
	}
};

Last to add. As with other geometrical problems, here a big part of implementation work is handling corner cases and, sometimes, caring about precision.
Re: Working with Polygons (response to post by syg96) | Reply
1. People should learn an explicit vector multiplication formula from a recipe Working with Lines and Segments, so it doesn't matter here.

2. Yes, I've added a description of this alternative algorithm (which made me split the recipe).
But you should remember that we don't teach people to solve exactly this problem, we teach them to solve geometrical problems in general, and we haven't enough recipes to show all these tiny functions directly. And I suppose that's better to carry these standard routines into separate functions than write these operations manually every time they encounter in the solution, especially speaking about really implementation-complicated problems.
Forums TopCoder Cookbook Algorithm Competitions - New Recipes Working with Polygons
Previous Thread  |  Next Thread
RSS