10000 Improved Code for Segments Intersect by TarasKuzyo · Pull Request #11553 · matplotlib/matplotlib · GitHub
[go: up one dir, main page]

Skip to content

Improved Code for Segments Intersect #11553

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 16 commits into from
Mar 2, 2019
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
77 changes: 77 additions & 0 deletions lib/matplotlib/tests/test_path.py
Original file line number Diff line number Diff line change
Expand Up @@ -273,6 +273,83 @@ def test_path_deepcopy():
copy.deepcopy(path2)


def test_path_intersect_path():
# test for the range of intersection angles
base_angles = np.array([0, 15, 30, 45, 60, 75, 90, 105, 120, 135])
angles = np.concatenate([base_angles, base_angles + 1, base_angles - 1])
eps_array = [1e-5, 1e-8, 1e-10, 1e-12]

for phi in angles:

transform = transforms.Affine2D().rotate(np.deg2rad(phi))

# a and b intersect at angle phi
a = Path([(-2, 0), (2, 0)])
b = transform.transform_path(a)
assert a.intersects_path(b) and b.intersects_path(a)

# a and b touch at angle phi at (0, 0)
a = Path([(0, 0), (2, 0)])
b = transform.transform_path(a)
assert a.intersects_path(b) and b.intersects_path(a)

# a and b are orthogonal and intersect at (0, 3)
a = transform.transform_path(Path([(0, 1), (0, 3)]))
b = transform.transform_path(Path([(1, 3), (0, 3)]))
assert a.intersects_path(b) and b.intersects_path(a)

# a and b are collinear and intersect at (0, 3)
a = transform.transform_path(Path([(0, 1), (0, 3)]))
b = transform.transform_path(Path([(0, 5), (0, 3)]))
assert a.intersects_path(b) and b.intersects_path(a)

# self-intersect
assert a.intersects_path(a)

# a contains b
a = transform.transform_path(Path([(0, 0), (5, 5)]))
b = transform.transform_path(Path([(1, 1), (3, 3)]))
assert a.intersects_path(b) and b.intersects_path(a)

# a and b are collinear but do not intersect
a = transform.transform_path(Path([(0, 1), (0, 5)]))
b = transform.transform_path(Path([(3, 0), (3, 3)]))
assert not a.intersects_path(b) and not b.intersects_path(a)

# a and b are on the same line but do not intersect
a = transform.transform_path(Path([(0, 1), (0, 5)]))
b = transform.transform_path(Path([(0, 6), (0, 7)]))
assert not a.intersects_path(b) and not b.intersects_path(a)

# Note: 1e-13 is the absolute tolerance error used for
# `isclose` function from src/_path.h

# a and b are parallel but do not touch
for eps in eps_array:
a = transform.transform_path(Path([(0, 1), (0, 5)]))
b = transform.transform_path(Path([(0 + eps, 1), (0 + eps, 5)]))
assert not a.intersects_path(b) and not b.intersects_path(a)

# a and b are on the same line but do not intersect (really close)
for eps in eps_array:
a = transform.transform_path(Path([(0, 1), (0, 5)]))
b = transform.transform_path(Path([(0, 5 + eps), (0, 7)]))
assert not a.intersects_path(b) and not b.intersects_path(a)

# a and b are on the same line and intersect (really close)
for eps in eps_array:
a = transform.transform_path(Path([(0, 1), (0, 5)]))
b = transform.transform_path(Path([(0, 5 - eps), (0, 7)]))
assert a.intersects_path(b) and b.intersects_path(a)

# b is the same as a but with an extra point
a = transform.transform_path(Path([(0, 1), (0, 5)]))
b = transform.transform_path(Path([(0, 1), (0, 2), (0, 5)]))
assert a.intersects_path(b) and b.intersects_path(a)

return


@pytest.mark.parametrize('offset', range(-720, 361, 45))
def test_full_arc(offset):
low = offset
Expand Down
34 changes: 32 additions & 2 deletions src/_path.h
8000
Original file line number Diff line number Diff line change
Expand Up @@ -813,6 +813,14 @@ int count_bboxes_overlapping_bbox(agg::rect_d &a, BBoxArray &bboxes)
return count;
}


inline bool isclose(double a, double b, double rtol, double atol)
{
// as per python's math.isclose
return fabs(a-b) <= fmax(rtol * fmax(fabs(a), fabs(b)), atol);
}


inline bool segments_intersect(const double &x1,
const double &y1,
const double &x2,
Expand All @@ -822,8 +830,27 @@ inline bool segments_intersect(const double &x1,
const double &x4,
const double &y4)
{
// relative and absolute tolerance values are chosen empirically
// it looks the atol value matters here bacause of round-off errors
const double rtol = 1e-10;
const double atol = 1e-13;
// determinant
double den = ((y4 - y3) * (x2 - x1)) - ((x4 - x3) * (y2 - y1));
if (den == 0.0) {

if (isclose(den, 0.0, rtol, atol)) { // collinear segments
if (x1 == x2 && x2 == x3) { // segments have infinite slope (vertical lines)
// and lie on the same line
return (fmin(y1, y2) <= fmin(y3, y4) && fmin(y3, y4) <= fmax(y1, y2)) ||
(fmin(y3, y4) <= fmin(y1, y2) && fmin(y1, y2) <= fmax(y3, y4));
}
else {
double intercept = (y1*x2 - y2*x1)*(x4 - x3) - (y3*x4 - y4*x3)*(x1 - x2);
if (isclose(intercept, 0.0, rtol, atol)) { // segments lie on the same line
return (fmin(x1, x2) <= fmin(x3, x4) && fmin(x3, x4) <= fmax(x1, x2)) ||
(fmin(x3, x4) <= fmin(x1, x2) && fmin(x1, x2) <= fmax(x3, x4));
}
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This looks better to me, though I didn't check your math. Note that you could just have the one return statement using y only if you wanted it all to be more elegant. Thanks for sticking with this!

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I don't think we can use y's only because of the same reasoning why we cannot use only x's.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If intercept == 0, it doesn't matter if you check y or x, everything is on the same line...

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In the intercept == 0, we excluded x = const lines, but if we replace x's with y's in the return statement the code will fail for y = const lines.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Ah, OK, fair enough... Thanks!

}

return false;
}

Expand All @@ -833,7 +860,10 @@ inline bool segments_intersect(const double &x1,
double u1 = n1 / den;
double u2 = n2 / den;

return (u1 >= 0.0 && u1 <= 1.0 && u2 >= 0.0 && u2 <= 1.0);
return ((u1 > 0.0 || isclose(u1, 0.0, rtol, atol)) &&
(u1 < 1.0 || isclose(u1, 1.0, rtol, atol)) &&
(u2 > 0.0 || isclose(u2, 0.0, rtol, atol)) &&
(u2 < 1.0 || isclose(u2, 1.0, rtol, atol)));
}

template <class PathIterator1, class PathIterator2>
Expand Down
0