how to find out if a point is in a polygon, both being in geographical coordinates

Question

I have simulated this algorithm http://forcefront.github.io/point-in-polygon/point-in-polygon.html.

But I think latitude and longitude are spheroidal coordinates.

Is that algorithm suitable for testing whether points are in a polygon on a map?

Answer 1

Here's a relatively simple system you can use if you have a map projection library, and if your polygon's points are quite close together (to minimise errors caused by projection of edges of your polygon, presumed to be great circles, to lines that are not great circles).

1. Project your polygon and your point using the Azimuthal Equidistant projection. Use the North Pole as the centre of the projection. If your point is the South Pole, use the South Pole as the centre of the projection. Alternatively, always use the point itself as the centre of projection.

2. Use a standard point-in-polygon routine to find out whether the projected point is in the projected polygon.

3. That's it!

This method correctly handles the discontinuity between 180 and -180 degrees of longitude, and the cases of a polygon that contains one or both poles.

If you don't have a map projection library, the formulae for the projection (using a spherical approximation to the earth's shape) are quite simple when the centre is a pole, and are explained in the Wikipedia article referred to.

Answer 2

It works in all cases, except when you

1) overlap the datum limit (limit where longitude jumps from 180 to -180) or when
2) overlap the poles.

So for 99.99% of all applications, just check that both conditions at geo data import and then use the linked algorithm.

Provide the info that polygons which do not meet these condition are ignored, and the datat provider has to split them before.

Answer 3

Mostly when you just have angles between -179 and 180 degree it should work. In case your angles is smaller or bigger then you need to convert it. You can use a conversion to world coordinate of every vertex and project it back to -179,180 degree. Here's how to compute a bounding box for example from New York in the U.S. to Bejing in China: Need to calculate latitude longitude from postal code database when location has multiple codes. Here is how you can normalize an angle Easy way to keeping angles between -179 and 180 degrees. You can also use leaflet JavaScript library with the wrap function.

Answer 4

A general-purpose point-in-polygon implementation would not account for the fact that you can wrap from a value like 179.999 degrees to -180 degrees, and we haven't been given sufficient information on the problem at hand. If your input is properly normalized (as in, it's placed in [-180, 180) for my case) and your polygons don't touch your map bounds, then you should have no problem.

When you look at an algorithm like this, the x and y coordinate systems are arbitrary; you could apply them to y and z, x and t, etc (though the x and t case doesn't necessary have any meaning that I can quickly think of).

If I haven't been clear in my description: If I gave you a polygon circling the north pole, it would be ambiguous as to whether or not that polygon included the north pole, or whether it included everything but the north pole; that hasn't been clearly relayed to me.

Likewise, if you look at a latitude/longitude projection of a sphere, a polygon intersect with lat or long = 180, and a general purpose point in polygon would not handle that.

Take my picture of the world to make things more clear. I'm not sure how you would define a polygon to be the green region on the left vs the region between the two half-circles (crossing through us and China).

Off the top of my head, you could probably solve these problems by splitting the two circles into two different polygons and testing them independently. You could represent the top half-circle as is, or the negation of it through a polygon with its points AND the 4 corners of the entire map.

Alternatively, you could shift the coordinate system whenever you test against your polygon (so that it doesn't wrap and its points are all normalized) though that would fail to account for polygons that span your globe numerous times (normalization would break your polygon).