Fast point in polygon test using binary space partitioning tree

Question

I am trying to code the Ray casting algorithm (http://en.wikipedia.org/wiki/Point_in_polygon) in a manner that will be effective when:

• the algorithm will be called many times with the same polygon but a different point
• the polygon has a large number of vertices

The polygon is simple but not necessarily convex.

I will use a horizontal ray from the point (xp,yp) going right: xr = xp + a, a>=0, yr = yp.

Instead of having to loop trough all edges of the polygon and check if the point is within the y range of the edge I want to use some Binary Space Partitioning tree techique(http://en.wikipedia.org/wiki/Binary_space_partitioning) to quickly finds such edges.

I imagine a binary tree of boxes:

struct Box {
    float               min, max;
    Box               *leftChild; // [min,mid]
    Box              *rightChild; // [mid,max]
    std::vector<int> edgeIndices; // edges that are completely within the box
};

1. All edges are placed in the root box which has min = min y(polygon) and max = max y(polygon). This becomes the current box.
2. The current box (min,max) is split into to: left: (min,mid) and right: (mid,max). Edges that are completely within either left or right are moved into that child box.

Step 2 is repeated for left and right child until each box contains less than a predefined number of edges or the depth of that child exceeds a predefined max depth.

Now in order to find which edges intersect a value yp: start at root box and traverse downward adding together all edges along the path.

How do I find the optimal mid values? I want the tree to be as flat and balanced as possible. That is the number of edges should be approximately the same in the left and the right child.

I could e.g. sort the edges on min y or max y and use the average y of the median edge as split point (mid).

Answer 1

If you have a non-convex polygon where most edges cross a certain y-value, then when you query a point with that y-value you will get essentially linear time in the number of edges to solve your point in polygon query. If you have a lot of preprocessing time and preprocessing space, then you can sort vertices according to y-coordinate and define a horizontal cutting line in between each pair of adjacent y-coordinate vertices, and order edges according to the x-values where they cross the horizontal cutting line. This will give a data structure of size at worst O(n^2) but potentially much smaller. Then, given a query point, you can binary search to find the two vertices whose y-coordinates are the closest above and below the y-coordinate of the query point, and then binary search on these edges to find where the query point's x-coordinate lies in the list of edges, and then report inside or outside depending on whether the x-position of the point is even or odd.in the sorted list of edges. Query time is O(log n)

Answer 2

You could use BSP tree with nodes being 2D planes dividing space into two areas. In this approach leaves are always convex polygons that are either outside or inside of the main polygon. This is a similar way of how Quake and all the FPS clones were made (look for PVS - Possible Visible Set). Note that creating this structure will require a lot of polygon cutting as each node divides the main polygon into to subpolygons. Using this structure is pretty simple: starting from the root you check your point with node's plane and choose one of the children that contains it. Finally you reach leaf and are either inside or outside. The essential issue however, is how to choose the dividing planes. Unfortunately I have no answer for this. Obviously one can chose the polygon edges instead of arbitrary lines and take into account several other factors. Most common are: minimizing future splits number or dividing the space most evenly. But frankly I am not really sure for this structure as it is meant for games and other interactive walk-throughs. I think it would be suitable with polys with really large number of edges.

Answer 3

You can use a treemap instead of the bsp tree. A treemap is split along the x-and y axis and you can sort the edges in descending order.

Answer 4

The structure you are building is in fact a bounding volume hierarchy, which can be used for ray-tracing. The idea is to have the edges in the leafs of the bounding volume hierarchy and to descend the tree when the ray interect a box.

http://en.wikipedia.org/wiki/Bounding_volume_hierarchy

In order to have an efficient descent into the tree you could use the surface area heuristic (SAH). It can be used for a bounding volume hierarchy too. The idea is to estimate the cost of traversing the hierarchy using the probability of either descending into the left or the rigth side and then to minimize this cost. The SAH can also implement a stop condition when the cost of descending further into the hierarchy would be higher than just looping througth all the remaining edges.

This question have link to a thesis describing the algorithm: KDTree Splitting