Hit-testing / finding intersection in a distorted rectangular mesh

Question

I have a rectangular mesh that is represented as a bidimensional array of bidimensional points. Something like this (C#-inspired pseudo-code):

Point2D[,] mesh = CreateSampleMesh();
mesh.Plot();

Important characteristics are:

1. The mesh is topologically rectangular, but geometrically "irregular" as the image shows;
2. The X coordinate of each row and Y coordinate of each column is monotonical (that is, it does not "reverse directions");
3. Every rectangle is convex (a necessary consequence of item 2, I guess);

My question is: given a candidate Point2D, how do I find which " rectangle" it's inside, if any? Here is the method signature I would suggest:

GridPosition position = GetIndexOfBottomLeftRectangleContaining(Point2D candidatePoint);
if (position != null)
{
    int pi = Position.I;
    int pj = Positino.J;
}

Answer 1

You may find useful the algorithm for finding whether a point is inside a convex polygon or not. Note that the second answer there suggests a method which costs four dot products in 2D, which assuming that multiplications take constant time, means time complexity of 4 * O(D) = 4 * O(2) = O(1).

Now you could go insane and check in parallel all the polygons you have there, which of course wouldn't be O(N²), but still an overkill probably.

For that reason, you need to partition space, with any data structure you like, where every partition would contain one of your polygons.

A quadtree:

is a tree data structure in which each internal node has exactly four children. Quadtrees are the two-dimensional analog of octrees and are most often used to partition a two-dimensional space by recursively subdividing it into four quadrants or regions.

seems to come in handy in that case, where every leaf would contain one polygon (or few of them, you choose).

Once a query lands on a leaf, then you would go and check if that point is inside the polygon(s) of this leaf and give the answer.

In particular, I would try a Region Quadtree, which "represents a partition of space in two dimensions by decomposing the region into four equal quadrants, subquadrants, and so on with each leaf node containing data corresponding to a specific subregion."

The query time is logarithmic, or linear in worst case scenario, as you can read in these slides.

Answer 2

You are lucky to have monotonic chains (polylines), which are in a way of sorting the plane.

First, design a procedure that will tell if you are left or right of a vertical chain. To achieve this, you find the line segment facing the point vertically, and draw an horizontal line. The intersection point tells you if you are left or right. This takes a dichotomic search among the node ordinates.

Use this procedure to find the two closest chains on either sides, that define a stripe where the point lies. Now a final dichotomic search in this stripe, by comparing to the horizontal edges will tell you the exact tile.

For an array of N² nodes, this will take time proportional to Log²N and require no extra data structure !

Answer 3

You could use some kind of 2D binary search.

Example: If you have a point p(px,py) and your mesh is stored in a 5x5 array 'mesh[i,j]', you would start with i=5/2=3 and j=5/2=3. Here is some basic example code for illustration.:

int i = 3;
int j = 3;
int iOld = i
int jOld = j;
do {
    iOld = i
    jOld = j;

    Point m = mesh[i,j]; 
    //Compare p with m: 
    if (p.x > m.x) {
        i = (i+5)/2;
    } else {
        i = (i/2)/2;
    }
    if (p.y > m.y) {
        j = (j+5)/2;
    } else {
        j = (j/2)/2;
    }
} while (i != iOld && j != jOld);

This code is only for illustration! You may have to verify that i and j are calculated correctly when you get close to the node, and that they don't oscillate, otherwise the exit condition may never become true. Otherwise this should find your point in the mesh pretty fast.

Edit: This approach assumes that distortion is 'benign', meaning that nodes are not contorted such that they move into a different square (which requires that the grid lines would intersect). If grid lines can cross, then this approach won't work.