Point in polygon algorithm for SQL Server

AlexB picture AlexB · Oct 27, 2013 · Viewed 8.2k times · Source

I'm trying to write a SQL query which determine if a given point is into a polygon. (I'm using SQL Server 2008 R2).

I was following this tutorial (just copy / paste it and change some table name) and it approximatively works BUT it is not precise at all. For example, let's considerate a given point which coordinates are :
P = 45.7664, 4.87383.

If you draw a little polygon (an approximate square) around this point with 4 vertices coordinates :
S = 45.97215 4.693909, 45.687 4.674683, 45.73302 5.460205, 46.05227 5.366821, 45.97215 4.693909

the procedure given in the link below answers the point is NOT in polygon, whereas it is... This is the output (with my own formatting text) :
false output
(Polygon 20 is the polygon above)

But if you enlarged the square (10 times bigger in my test), the procedure answers my point is in the square.

So, I'm seeking for another algorithm, more precise.
Here is my VERTICE table, containing all vertices coordinates of each polygon of my DB Vertice table

I need to check, for ALL polygons (there are a few thousands) if a given point passed in parameter is in a polygon (and if yes, which one(s)). I can do the loop by myself, but I miss a correct Point in polygon algorithm.

Could someone help me ? Thank you very much.

SUMMARY

The corresponding SQL fiddle : http://sqlfiddle.com/#!3/0caa4/1

Answer

Ben Thul picture Ben Thul · Oct 27, 2013

Your problem is that you've defined the polygon backwards. Let's explore this with the native SQL geography types

declare @s geography, @t geography, @p geography;
select @s = geography::STPolyFromText('POLYGON((45.97215 4.693909, 45.687 4.674683, 45.73302 5.460205, 46.05227 5.366821, 45.97215 4.693909))', 4326);
select @t = geography::STPolyFromText('POLYGON((46.05227 5.366821, 45.73302 5.460205, 45.687 4.674683, 45.97215 4.693909, 46.05227 5.366821))', 4326);
select @p = geography::STPointFromText('POINT(45.7664 4.87383)', 4326);

select @p.STIntersects(@s), @p.STIntersects(@t);
select @p.STBuffer(10), @s, @t;

As you can see, @s holds pretty much the whole world. In the result set viewer, if you zoom into where your "square" should be, you'll find a hole. Contrast that with @t which is the area that you expect. Also note that I ran this in SQL 2012 which improved on a limitation that existed pre-2012 that says that a geospatial instance can't cross a hemisphere boundary. If you run the above on a 2008 instance, you'll likely get an error to that effect for @s. Comment out the line defining @s and it'll run.

There's a "right-hand" rule when defining geo(graphy/metry) polygons. Imagine that you were in a car visiting the points in the order you've specified. The area that you're defining is what you'd see if you were looking out the right side of the car.