------------------------------------------------------------------------------ Determine whether a point is inside a polygon. The function uses a winding number calculation generalized to the 3D plane one which the polygon resides. Returns 0 if point is not in the polygon; 1 if it is inside. Can also return -1 to indicate a degenerate polygon. This implementation is inspired by Dan Sunday's algor
| 727 | // inspired by Dan Sunday's algorithm found in the book Practical Geometry |
| 728 | // Algorithms. |
| 729 | int vtkPolygon::PointInPolygon(double x[3], int numPts, double* pts, double bounds[6], double* n) |
| 730 | { |
| 731 | // Do a quick bounds check to throw out trivial cases. |
| 732 | // winding plane. |
| 733 | if (x[0] < bounds[0] || x[0] > bounds[1] || x[1] < bounds[2] || x[1] > bounds[3] || |
| 734 | x[2] < bounds[4] || x[2] > bounds[5]) |
| 735 | { |
| 736 | return VTK_POLYGON_OUTSIDE; |
| 737 | } |
| 738 | |
| 739 | // Check that the normal is non-zero. |
| 740 | if (vtkMath::Norm(n) <= FLT_EPSILON) |
| 741 | { |
| 742 | return VTK_POLYGON_FAILURE; |
| 743 | } |
| 744 | |
| 745 | // Assess whether the point lies on the boundary of the polygon. Points on |
| 746 | // the boundary are considered inside the polygon. Need to define a small |
| 747 | // tolerance relative to the bounding box diagonal length of the polygon. |
| 748 | double tol2 = VTK_POLYGON_TOL * |
| 749 | ((bounds[1] - bounds[0]) * (bounds[1] - bounds[0]) + |
| 750 | (bounds[3] - bounds[2]) * (bounds[3] - bounds[2]) + |
| 751 | (bounds[5] - bounds[4]) * (bounds[5] - bounds[4])); |
| 752 | tol2 *= tol2; |
| 753 | tol2 = (tol2 == 0.0 ? FLT_EPSILON : tol2); |
| 754 | |
| 755 | for (int i = 0; i < numPts; i++) |
| 756 | { |
| 757 | // Check coincidence to polygon vertices |
| 758 | double* p0 = pts + 3 * i; |
| 759 | if (vtkMath::Distance2BetweenPoints(x, p0) <= tol2) |
| 760 | { |
| 761 | return VTK_POLYGON_INSIDE; |
| 762 | } |
| 763 | |
| 764 | // Check coincidence to polygon edges |
| 765 | double* p1 = pts + 3 * ((i + 1) % numPts); |
| 766 | double t; |
| 767 | double d2 = vtkLine::DistanceToLine(x, p0, p1, t); |
| 768 | if (d2 <= tol2 && 0.0 < t && t < 1.0) |
| 769 | { |
| 770 | return VTK_POLYGON_INSIDE; |
| 771 | } |
| 772 | } |
| 773 | |
| 774 | // If here, begin computation of the winding number. This method works for |
| 775 | // points/polygons arbitrarily oriented in 3D space. Hence a projection |
| 776 | // onto one of the x-y-z coordinate planes using the maximum normal |
| 777 | // component. The computation will be performed in the (axis0,axis1) plane. |
| 778 | int axis0, axis1; |
| 779 | if (fabs(n[0]) > fabs(n[1])) |
| 780 | { |
| 781 | if (fabs(n[0]) > fabs(n[2])) |
| 782 | { |
| 783 | axis0 = 1; |
| 784 | axis1 = 2; |
| 785 | } |
| 786 | else |
no test coverage detected