I would like to determine a polygon and implement an algorithm which would check if a point is inside or outside the polygon.
Does anyone know if there is any example available of any similar algorithm?
I would like to determine a polygon and implement an algorithm which would check if a point is inside or outside the polygon.
Does anyone know if there is any example available of any similar algorithm?
C# code
bool IsPointInPolygon(List<Loc> poly, Loc point)
{
int i, j;
bool c = false;
for (i = 0, j = poly.Count - 1; i < poly.Count; j = i++)
{
if ((((poly[i].Lt <= point.Lt) && (point.Lt < poly[j].Lt))
|| ((poly[j].Lt <= point.Lt) && (point.Lt < poly[i].Lt)))
&& (point.Lg < (poly[j].Lg - poly[i].Lg) * (point.Lt - poly[i].Lt)
/ (poly[j].Lt - poly[i].Lt) + poly[i].Lg))
{
c = !c;
}
}
return c;
}
Location class
public class Loc
{
private double lt;
private double lg;
public double Lg
{
get { return lg; }
set { lg = value; }
}
public double Lt
{
get { return lt; }
set { lt = value; }
}
public Loc(double lt, double lg)
{
this.lt = lt;
this.lg = lg;
}
}
After searching the web and trying various implementations and porting them from C++ to C# I finally got my code straight:
public static bool PointInPolygon(LatLong p, List<LatLong> poly)
{
int n = poly.Count();
poly.Add(new LatLong { Lat = poly[0].Lat, Lon = poly[0].Lon });
LatLong[] v = poly.ToArray();
int wn = 0; // the winding number counter
// loop through all edges of the polygon
for (int i = 0; i < n; i++)
{ // edge from V[i] to V[i+1]
if (v[i].Lat <= p.Lat)
{ // start y <= P.y
if (v[i + 1].Lat > p.Lat) // an upward crossing
if (isLeft(v[i], v[i + 1], p) > 0) // P left of edge
++wn; // have a valid up intersect
}
else
{ // start y > P.y (no test needed)
if (v[i + 1].Lat <= p.Lat) // a downward crossing
if (isLeft(v[i], v[i + 1], p) < 0) // P right of edge
--wn; // have a valid down intersect
}
}
if (wn != 0)
return true;
else
return false;
}
private static int isLeft(LatLong P0, LatLong P1, LatLong P2)
{
double calc = ((P1.Lon - P0.Lon) * (P2.Lat - P0.Lat)
- (P2.Lon - P0.Lon) * (P1.Lat - P0.Lat));
if (calc > 0)
return 1;
else if (calc < 0)
return -1;
else
return 0;
}
The isLeft function was giving me rounding problems and I spent hours without realizing that I was doing the conversion wrong, so forgive me for the lame if block at the end of that function.
BTW, this is the original code and article: http://softsurfer.com/Archive/algorithm_0103/algorithm_0103.htm
I think there is a simpler and more efficient solution.
Here is the code in C++. I should be simple to convert it to C#.
int pnpoly(int npol, float *xp, float *yp, float x, float y)
{
int i, j, c = 0;
for (i = 0, j = npol-1; i < npol; j = i++) {
if ((((yp[i] <= y) && (y < yp[j])) ||
((yp[j] <= y) && (y < yp[i]))) &&
(x < (xp[j] - xp[i]) * (y - yp[i]) / (yp[j] - yp[i]) + xp[i]))
c = !c;
}
return c;
}
By far the best explanation and implementation can be found at Point In Polygon Winding Number Inclusion
There is even a C++ implementation at the end of the well explained article. This site also contains some great algorithms/solutions for other geometry based problems.
I have modified and used the C++ implementation and also created a C# implementation. You definitely want to use the Winding Number algorithm as it is more accurate than the edge crossing algorithm and it is very fast.
The complete solution in asp.Net C#, you can see the complete detail here, you can see how to find point(lat,lon) whether its inside or Outside the Polygon using the latitude and longitudes ? Article Reference Link
private static bool checkPointExistsInGeofencePolygon(string latlnglist, string lat, string lng) {
List<Loc> objList = new List<Loc>();
// sample string should be like this strlatlng = "39.11495,-76.873259|39.114588,-76.872808|39.112921,-76.870373|";
string[] arr = latlnglist.Split('|');
for (int i = 0; i <= arr.Length - 1; i++)
{
string latlng = arr[i];
string[] arrlatlng = latlng.Split(',');
Loc er = new Loc(Convert.ToDouble(arrlatlng[0]), Convert.ToDouble(arrlatlng[1]));
objList.Add(er);
}
Loc pt = new Loc(Convert.ToDouble(lat), Convert.ToDouble(lng));
if (IsPointInPolygon(objList, pt) == true)
{
return true;
}
else
{
return false;
}
}
private static bool IsPointInPolygon(List<Loc> poly, Loc point)
{
int i, j;
bool c = false;
for (i = 0, j = poly.Count - 1; i < poly.Count; j = i++)
{
if ((((poly[i].Lt <= point.Lt) && (point.Lt < poly[j].Lt)) |
((poly[j].Lt <= point.Lt) && (point.Lt < poly[i].Lt))) &&
(point.Lg < (poly[j].Lg - poly[i].Lg) * (point.Lt - poly[i].Lt) / (poly[j].Lt - poly[i].Lt) + poly[i].Lg))
c = !c;
}
return c;
}
Just a heads up (using answer as I can't comment), if you want to use point-in-polygon for geo fencing, then you need to change your algorithm to work with spherical coordinates. -180 longitude is the same as 180 longitude and point-in-polygon will break in such situation.
Relating to kobers answer I worked it out with more readable clean code and changed the longitudes that crosses the date border:
public bool IsPointInPolygon(List<PointPosition> polygon, double latitude, double longitude)
{
bool isInIntersection = false;
int actualPointIndex = 0;
int pointIndexBeforeActual = polygon.Count - 1;
var offset = calculateLonOffsetFromDateLine(polygon);
longitude = longitude < 0.0 ? longitude + offset : longitude;
foreach (var actualPointPosition in polygon)
{
var p1Lat = actualPointPosition.Latitude;
var p1Lon = actualPointPosition.Longitude;
var p0Lat = polygon[pointIndexBeforeActual].Latitude;
var p0Lon = polygon[pointIndexBeforeActual].Longitude;
if (p1Lon < 0.0) p1Lon += offset;
if (p0Lon < 0.0) p0Lon += offset;
// Jordan curve theorem - odd even rule algorithm
if (isPointLatitudeBetweenPolyLine(p0Lat, p1Lat, latitude)
&& isPointRightFromPolyLine(p0Lat, p0Lon, p1Lat, p1Lon, latitude, longitude))
{
isInIntersection = !isInIntersection;
}
pointIndexBeforeActual = actualPointIndex;
actualPointIndex++;
}
return isInIntersection;
}
private double calculateLonOffsetFromDateLine(List<PointPosition> polygon)
{
double offset = 0.0;
var maxLonPoly = polygon.Max(x => x.Longitude);
var minLonPoly = polygon.Min(x => x.Longitude);
if (Math.Abs(minLonPoly - maxLonPoly) > 180)
{
offset = 360.0;
}
return offset;
}
private bool isPointLatitudeBetweenPolyLine(double polyLinePoint1Lat, double polyLinePoint2Lat, double poiLat)
{
return polyLinePoint2Lat <= poiLat && poiLat < polyLinePoint1Lat || polyLinePoint1Lat <= poiLat && poiLat < polyLinePoint2Lat;
}
private bool isPointRightFromPolyLine(double polyLinePoint1Lat, double polyLinePoint1Lon, double polyLinePoint2Lat, double polyLinePoint2Lon, double poiLat, double poiLon)
{
// lon <(lon1-lon2)*(latp-lat2)/(lat1-lat2)+lon2
return poiLon < (polyLinePoint1Lon - polyLinePoint2Lon) * (poiLat - polyLinePoint2Lat) / (polyLinePoint1Lat - polyLinePoint2Lat) + polyLinePoint2Lon;
}
I add one detail to help people who live in the... south of earth!! If you're in Brazil (that's my case), our GPS coord are all negatives. And all these algo give wrong results.
The easiest way if to use the absolute values of the Lat and Long of all point. And in that case Jan Kobersky's algo is perfect.
the polygon is defined as a sequential list of point pairs A, B, C .... A. no side A-B, B-C ... crosses any other side
Determine box Xmin, Xmax, Ymin, Ymax
case 1 the test point P lies outside the box
case 2 the test point P lies inside the box:
Determine the 'diameter' D of the box {[Xmin,Ymin] - [Xmax, Ymax]} ( and add a little extra to avoid possible confusion with D being on a side)
Determine the gradients M of all sides
Find a gradient Mt most different from all gradients M
The test line runs from P at gradient Mt a distance D.
Set the count of intersections to zero
For each of the sides A-B, B-C test for the intersection of P-D with a side from its start up to but NOT INCLUDING its end. Increment the count of intersections if required. Note that a zero distance from P to the intersection indicates that P is ON a side
An odd count indicates P is inside the polygon
I translated c# method in Php and I added many comments to understand code.
Description of PolygonHelps:
Check if a point is inside or outside of a polygon. This procedure uses gps coordinates and it works when polygon has a little geographic area.
INPUT:
$poly: array of Point: polygon vertices list; [{Point}, {Point}, ...];
$point: point to check; Point: {"lat" => "x.xxx", "lng" => "y.yyy"}
When $c is false, the number of intersections with polygon is even, so the point is outside of polygon;
When $c is true, the number of intersections with polygon is odd, so the point is inside of polygon;
$n is the number of vertices in polygon;
For each vertex in polygon, method calculates line through current vertex and previous vertex and check if the two lines have an intersection point.
$c changes when intersection point exists.
So, method can return true if point is inside the polygon, else return false.
class PolygonHelps {
public static function isPointInPolygon(&$poly, $point){
$c = false;
$n = $j = count($poly);
for ($i = 0, $j = $n - 1; $i < $n; $j = $i++){
if ( ( ( ( $poly[$i]->lat <= $point->lat ) && ( $point->lat < $poly[$j]->lat ) )
|| ( ( $poly[$j]->lat <= $point->lat ) && ( $point->lat < $poly[$i]->lat ) ) )
&& ( $point->lng < ( $poly[$j]->lng - $poly[$i]->lng )
* ( $point->lat - $poly[$i]->lat )
/ ( $poly[$j]->lat - $poly[$i]->lat )
+ $poly[$i]->lng ) ){
$c = !$c;
}
}
return $c;
}
}
Jan's answer is great.
Here is the same code using the GeoCoordinate class instead.
using System.Device.Location;
...
public static bool IsPointInPolygon(List<GeoCoordinate> poly, GeoCoordinate point)
{
int i, j;
bool c = false;
for (i = 0, j = poly.Count - 1; i < poly.Count; j = i++)
{
if ((((poly[i].Latitude <= point.Latitude) && (point.Latitude < poly[j].Latitude))
|| ((poly[j].Latitude <= point.Latitude) && (point.Latitude < poly[i].Latitude)))
&& (point.Longitude < (poly[j].Longitude - poly[i].Longitude) * (point.Latitude - poly[i].Latitude)
/ (poly[j].Latitude - poly[i].Latitude) + poly[i].Longitude))
c = !c;
}
return c;
}
Check if a point is inside a polygon or not -
Consider the polygon which has vertices a1,a2,a3,a4,a5. The following set of steps should help in ascertaining whether point P lies inside the polygon or outside.
Compute the vector area of the triangle formed by edge a1->a2 and the vectors connecting a2 to P and P to a1. Similarly, compute the vector area of the each of the possible triangles having one side as the side of the polygon and the other two connecting P to that side.
For a point to be inside a polygon, each of the triangles need to have positive area. Even if one of the triangles have a negative area then the point P stands out of the polygon.
In order to compute the area of a triangle given vectors representing its 3 edges, refer to http://www.jtaylor1142001.net/calcjat/Solutions/VCrossProduct/VCPATriangle.htm
The problem is easier if your polygon is convex. If so, you can do a simple test for each line to see if the point is on the inside or outside of that line (extending to infinity in both directions). Otherwise, for concave polygons, draw an imaginary ray from your point out to infinity (in any direction). Count how many times it crosses a boundary line. Odd means the point is inside, even means the point is outside.
This last algorithm is trickier than it looks. You will have to be very careful about what happens when your imaginary ray exactly hits one of the polygon's vertices.
If your imaginary ray goes in the -x direction, you can choose only to count lines that include at least one point whose y coordinate is strictly less than the y coordinate of your point. This is how you get most of the weird edge cases to work correctly.
If you have a simple polygon (none of the lines cross) and you don't have holes you can also triangulate the polygon, which you are probably going to do anyway in a GIS app to draw a TIN, then test for points in each triangle. If you have a small number of edges to the polygon but a large number of points this is fast.
For an interesting point in triangle see link text
Otherwise definately use the winding rule rather than edge crossing, edge crossing has real problems with points on edges, which if your data is generated form a GPS with limited precision is very likely.
You can try this simple class https://github.com/xopbatgh/sb-polygon-pointer
It is easy to deal with it
$polygonBox = [
[55.761515, 37.600375],
[55.759428, 37.651156],
[55.737112, 37.649566],
[55.737649, 37.597301],
];
$sbPolygonEngine = new sbPolygonEngine($polygonBox);
$isCrosses = $sbPolygonEngine->isCrossesWith(55.746768, 37.625605);
// $isCrosses is boolean
(answer was returned from deleted by myself because initially it was formatted wrong)