CODEWITHSUNDEEP
graphicsgeometrycomputational-geometry
Andy
Andy

Reputation: 3789

How do I QUICKLY find the closest intersection in 2D between a ray and m polylines?

How do I find the closest intersection in 2D between a ray:

x = x0 + t*cos(a), y = y0 + t*sin(a)

and m polylines:

{(x1,y1), (x2,y2), ..., (xn,yn)}

QUICKLY?

I started by looping trough all linesegments and for each linesegment; {(x1,y1),(x2,y2)} solving:

x1 + u*(x2-x1) = x0 + t*cos(a)

y1 + u*(y2-y1) = y0 + t*sin(a)

by Cramer's rule, and afterward sorting the intersections on distance, but that was slow :-(

BTW: the polylines happens to be monotonically increasing in x.

Upvotes: 2

Views: 1469

Answers (2)

MvG
MvG

Reputation: 60858

Coordinate system transformation

I suggest you first transform your setup to something with easier coordinates:

  1. Take your point p = (x, y).
  2. Move it by (-x0, -y0) so that the ray now starts at the center.
  3. Rotate it by -a so that the ray now lies on the x axis.

So far the above operations have cost you four additions and four multiplications per point:

ca = cos(a) # computed only once
sa = sin(a) # likewise

x' = x - x0
y' = y - y0
x'' = x'*ca + y'*sa 
y'' = y'*ca - x'*sa

Checking for intersections

Now you know that a segment of the polyline will only intersect the ray if the sign of its y'' value changes, i.e. y1'' * y2'' < 0. You could even postpone the computation of the x'' values until after this check. Furthermore, the segment will only intersect the ray if the intersection of the segment with the x axis occurs for x > 0, which can only happen if either value is greater than zero, i.e. x1'' > 0 or x2'' > 0. If both x'' are greater than zero, then you know there is an intersection.

The following paragraph is kind of optional, don't worry if you don't understand it, there is an alternative noted later on.
If one x'' is positive but the other is negative, then you have to check further. Suppose that the sign of y'' changed from negative to positive, i.e. y1'' < 0 < y2''. The line from p1'' to p2'' will intersect the x axis at x > 0 if and only if the triangle formed by p1'', p2'' and the origin is oriented counter-clockwise. You can determine the orientation of that triangle by examining the sign of the determinant x1''*y2'' - x2''*y1'', it will be positive for a counter-clockwise triangle. If the direction of the sign change is different, the orientation has to be different as well. So to take this together, you can check whether 

(x1'' * y2'' - x2'' * y1'') * y2'' > 0

If that is the case, then you have an intersection. Notice that there were no costly divisions involved so far.

Computing intersections

As you want to not only decide whether an intersection exists, but actually find a specific one, you now have to compute that intersection. Let's call it p3. It must satisfy the equations

(x2'' - x3'')/(y2'' - y3'') = (x1'' - x3'')/(y1'' - y3'') and 
                       y3'' = 0

which results in

x3'' = (x1'' * y1'' - x2'' * y2'')/(y1'' - y2'')

Instead of the triangle orientation check from the previous paragraph, you could always compute this x3'' value and discard any results where it turns out to be negative. Less code, but more divisions. Benchmark if in doubt about performance.

To find the point closest to the origin of the ray, you take the result with minimal x3'' value, which you can then transform back into its original position:

x3 = x3''*ca + x0
y3 = x3''*sa + y0

There you are.

Note that all of the above assumed that all numbers were either positive or negative. If you have zeros, it depends on the exact interpretation of what you actually want to compute, how you want to handle these border cases.

Upvotes: 2

Ante
Ante

Reputation: 5448

To avoid checking intersection with all segments, some space partition is needed, like Quadtree, BSP tree. With space partition it is needed to check ray intersection with space partitions.

In this case, since points are sorted by x-coordinate, it is possible to make space partition with boxes (min x, min y)-(max x, max y) for parts of polyline. Root box is min-max of all points, and it is split in 2 boxes for first and second part of a polyline. Number of segments in parts is same or one box has one more segment. This box splitting is done recursively until only one segment is in a box.

To check ray intersection start with root box and check is it intersected with a ray, if it is than check 2 sub-boxes for an intersection and first test closer sub-box then farther sub-box.

Checking ray-box intersection is checking if ray is crossing axis aligned line between 2 positions. That is done for 4 box boundaries.

Upvotes: 1

Related Questions