Crossposted at

http://stackoverflow.com/questions/16531963/algorithm-for-points-inside-a-polyline-inexplicably-fails-autocad No answers yet.

I have an algorithm (which I did not write) for determining if a point lies within a polygon (straight lines only).

I have been testing my command on the same 51 polygons repeatedly. 99% it will work perfectly. Every once in a while it will fail on 1 or more of the polygons, returning false for over 2000 points I am creating inside the bounding box of the polyline. I have seen it fail when the polyline is a simple rectangle and all of the points lie distributed in a grid within the polyline. It should have returned true over 2000 times in that case. It will never fail for just 1 of the points, it will fail all of them. I have confirmed that the points are being correctly created where I expect them to be and that the vertices of the polygon are where I expect them to be. When it fails, the angle variable for the last point is at exactly double PI.

I am not doing any multi-threading. The only possibly 'funny' thing I am doing is COM Interop with Excel. This is happening after the transaction has been committed for the part with this algorithm, and I am sure I am cleaning up all my COM objects. I have not been able to reproduce the failure without the COM Interop part but I don't think I've tested it enough yet to have enough absence of evidence.

Any ideas what could be wrong?

bool IsInsidePolygon(Polyline polygon, Point3d pt)

{

int n = polygon.NumberOfVertices;

double angle = 0;

Point pt1, pt2;

for (int i = 0; i < n; i++)

{

pt1.X = polygon.GetPoint2dAt(i).X - pt.X;

pt1.Y = polygon.GetPoint2dAt(i).Y - pt.Y;

pt2.X = polygon.GetPoint2dAt((i + 1) % n).X - pt.X;

pt2.Y = polygon.GetPoint2dAt((i + 1) % n).Y - pt.Y;

angle += Angle2D(pt1.X, pt1.Y, pt2.X, pt2.Y);

}

if (Math.Abs(angle) < Math.PI)

return false;

else

return true;

}

public struct Point

{

public double X, Y;

};

public static double Angle2D(double x1, double y1, double x2, double y2)

{

double dtheta, theta1, theta2;

theta1 = Math.Atan2(y1, x1);

theta2 = Math.Atan2(y2, x2);

dtheta = theta2 - theta1;

while (dtheta > Math.PI)

dtheta -= (Math.PI * 2);

while (dtheta < -Math.PI)

dtheta += (Math.PI * 2);

return (dtheta);

}