Points in Polygons

[Shaggy Hiker]

I don't really need an answer to this question, since I already have a solution that will work for this year, but I'll put it out here in case anybody has some thoughts on the subject:

Suppose you have a polygon that is defined by a series of points. The exact number of points could be quite large, but it is ok to assume that the shape is a convex polygon. Given another point, how would you determine whether or not it was in the polygon?

Where am I going with this? I'm not really sure. This is robot code, as most questions I post here seem to be these days. The robot will be seeing with sonar, but the sonar will return nothing but echoes with a certain amount of error (about 2cm, I believe). Eventually, the robot will find a point and need to decide whether it is looking at object A, or the adjacent object B. If the two are close enough together, then it makes not a whit of difference, which is how I have solved the problem (as long as the new point is close enough that the echo could have come from either object, then I assign it to one or the other based on a half-assed algorithm, but that's fine, because if the two objects are really close together, this particular module doesn't need to be able to tell them apart).

Eventually, there will be another module that will be deciding whether or not to investigate suspicious points, and how much urgency it should attach to the investigation if it decides to investigate. At that point, it may be of value to come up with a definitive answer as to whether a point is within a polygon defined by a series of points. Might actually be working on that in a week or two.

[jemidiah]

If you could get the order of the points in the polygon, then you could cut your polygon into a series of triangles with two adjacent points on the polygon connected to the test point. After that you could sum the angle next to the test point for each polygon. Inside the polygon, you would get a sum of 2*pi. Outside would be greater (I think).

Ordering the points would be the tough part. Summing the angles could be done (as I seem to love to do) using dot products and n=number of sides inverse cosines. If you want speed that might be a limiting factor of this method. Otherwise I think it's pretty simple.

Edit: I take that back. I think since the polygon is convex, you can assume that the two adjacent points for any edge vertex are those that are closest in distance on the left and right of the line connecting the point with the middle (center of mass) of the polygon.

[Shaggy Hiker]

I considered splitting the polygon into triangles, and felt that the assumption of concavity would be reasonable, which makes that task easier. The problem that may arise is that the number of points on the polygon may well be absurd. Unfortunately, I won't even be able to test that for at least another month or two, so at this point, it's all just theoretical.

Normally, a polygon is made up of some reasonable number of points, but in my case, it is entirely possible that any object will be defined by a HUGE number of points....or not. In any case, the points can't be farther than 6cm apart, and are likely to be much closer.

Now that I think about it, this may be a problem more akin to character recognition rather than geometry. If there really are lots of points along each edge (and there really shouldn't be any in the interior, then I might be better off "walking the edge" to come up with a series of line segments that outline the object. I'll have to think about that a bit more.

[krtxmrtz]

...Given another point, how would you determine whether or not it was in the polygon?...

Here's an answer to that question.

You have to consider a horizontal straight line starting at the test point and going in the positive x direction and count how many times it intercepts the polygon. If it's an odd number of times the point is inside, else it's outside. If the polygon is convex these numbers will be 1 (inside) and 0 or 2 (outside at the right or left side of the polygon respectively), but the code below works for any polygon including concave ones.

Function InOrOut(p, q, k, p0, q0)
 
'Function to determine the position of point (p0,q0) relative
'to the closed contour defined by a number of points with
'coordinates p(i), q(i), i=1, 2, ..., k
'Point 1 must be the same as point k
'Function values():
'       1: point is inside the contour
'      -1: point is outside the contour
'       0: point lies on the contour border
'*********************************************************************
    Dim kross As Integer, i As Integer
    Dim pp As Single
    'Function initialization
    InOrOut = 0
    'Initialization of kross, a variable which keeps track
    'of how many times the horizontal semi-infinite straight
    'line starting at (p0,q0) and in the positive (right hand side)
    'x axis direction intercepts the contour
    '(The contour may have vertices with angles larger than 180 degrees)
    kross = 0
    'Loop over all contour sides
    For i = 1 To k - 1
        'If the side between the i and i+1 vertices lies entirely
        'above or below the point (p0,q0), then continue on
        'with the next side (there is no intercept)
        If (q(i) > q0 And q(i + 1) > q0) Or (q(i) < q0 And q(i + 1) < q0) Then GoTo NextItem
        'If the side is horizontal avoid the calculation of the intercept by interpolation
        'as there would be a division by zero
        If q(i) = q(i + 1) Then
            'It has to be determined if the point (p0,q0) lies on this horizontal segment
            If (p(i) > p0 And p(i + 1) > p0) Or (p(i) < p0 And p(i + 1) < p0) Then GoTo NextItem
            'If it doesn't, we're done!
            Exit Function
        End If
        'Calculation of pp, the coordinate of the point where the segment connecting the
        'sides i and i+1 and the horizontal straight line that goes through (p0,q0) intercept
        pp = p(i) + (q0 - q(i)) * ((p(i + 1) - p(i)) / (q(i + 1) - q(i)))
        'The sign of pp-p0 determines the position of the intercept relative to (p0,q0)
        If pp - p0 > 0 Then
            'Intercept to the right: increment counter (but only if the intercept does not lie
            'on the first vertex, to avoid counting the same intercept twice!)
            If q0 <> q(i) Then kross = kross + 1
        ElseIf pp - p0 = 0 Then
            Exit Function
        End If
        'If the intercept lies to the left, take no action and continue on
NextItem:
    Next
    'If the number of intercepts is even, the point (p0,q0) lies out of the contour
    'If it's odd then it lies inside
    If kross Mod 2 = 0 Then
        InOrOut = -1
    Else
        InOrOut = 1
    End If
End Function

[Shaggy Hiker]

What language was that written in? It looks like an early flavor of VB. I can readily convert it to my needs, and the point you made prior to the code is a good one. However, in this case, I am beginning to think that I will need to do something slightly different prior to any such attempt.

Thanks for posting that, it may well come in handy.

[krtxmrtz]

What language was that written in? It looks like an early flavor of VB...

Good old vb6!

[jemidiah]

That's a really clever idea. I like it, much simpler than mine with the difficulty of ordering points.

[Shaggy Hiker]

I actually got working on designing a module that dealt with gaps between objects (the reason behind the intersecting lines thread), and while working on that design, I realized that this question was irrelevant. For those points that are clearly part of the object, then I don't need to know if they are inside or outside, whereas for points that are NOT clearly part of one object or another close object, then the point gets ignored for certain calculations, rendering the whole question moot at this point.

I also decided that I can greatly improve speed for some calculations by figuring a rough bounding box for each object. I'm not yet sure whether or not I need to know where a point lies with relation to a true polygon, or if knowing how it lies relative to the bounding box will be enough, but for now, it doesn't matter either way.

I'm sure I will change my mind as I begin working on the next module after this one, but for now, all is well.