Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Trouble with boolean_point_in_polygon() on irregular polygons. #100

Open
mcgeochd opened this issue Jun 12, 2021 · 1 comment
Open

Trouble with boolean_point_in_polygon() on irregular polygons. #100

mcgeochd opened this issue Jun 12, 2021 · 1 comment

Comments

@mcgeochd
Copy link

mcgeochd commented Jun 12, 2021

I have a program that uses scipy spatial SpericalVoronoi to generate irregular non-overlapping polygons covering the surface of a sphere and I'm trying to use turfy to rasterize the polygons with boolean_point_in_polygon(). Some example outputs are:

https://imgur.com/a/nIQTSNF
https://imgur.com/a/8DpitYe

The top image is a method of rasterizing that calculates the distance from a raster square to every edge of every polygon (an edge in this case refers to a single great circle segment between vertices) by projecting the point onto the line segment, finds the edge the raster point is closest to, and determines which side of this edge its on using the sign of the cross track distance (using a custom implementation, not the turpy implementation), and therefore which polygon it belongs to. This has some issues, as can be seen with the inlcusions of the wrong coloured polygon. This occurs when two edges are equally close to a raster point, i.e. the vertex shared by two edges is the closest point, and the wrong edge is selected, resulting in an incorrect rasterization.

I'm trying to get around this by using turfpy to rasterize the polygons, the results of which are shown in the lower image. The trouble is that the geojson polygon needs the polygon vertices to be supplied in clockwise order. While I can create a directed edge for each polygon, I'm having trouble detecting whether a given directed edge is defined clockwise or not, if it isn't I can just reverse the list order. I've tried using Python implementations of two methods described here without much luck: https://stackoverflow.com/questions/1165647/how-to-determine-if-a-list-of-polygon-points-are-in-clockwise-order/1180256#1180256

Signed Area

area = 0
for i in range(len(verts)):
    v1 = verts[i]
    v2 = verts[(i+1)%len(verts)]
    area += (v2[0] - v1[0]) * (v2[1] + v1[1])
if area < 0: verts = verts[::-1]
self.poly = Polygon([verts])

Convex Hull Point

minLat = 1e18
A = None
for v in verts:
    if v[1] < minLat:
        minLat = v[1]
        A = v
    elif v[1] == minLat and v[0] > A[0]:
        A = v
Aindex = indexOf(verts, A)
A = np.array(A)
B = np.array(verts[(Aindex-1)%len(verts)])
C = np.array(verts[(Aindex+1)%len(verts)])
AB = B - A
AC = C - A
cross = np.cross(AB, AC)
if cross < 0: verts = verts[::-1]
self.poly = Polygon([verts])

If anybody knows what is going wrong, or what algorithm turfpy uses for boolean_point_in_polygon so I can write an implementation that is direction agnostic, that would be super helpful.

@mcgeochd
Copy link
Author

I fixed a problem where the GeoJson polygons have to have the same first and last point, which I had missed previously, but the output is still full of errors: https://imgur.com/a/lV2ZqsY

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

1 participant