shapely / shapely

Manipulation and analysis of geometric objects
https://shapely.readthedocs.io/en/stable/
BSD 3-Clause "New" or "Revised" License
3.88k stars 566 forks source link

Possible way to fix issue when using intersection points to split linestring #630

Open Andreachen0707 opened 6 years ago

Andreachen0707 commented 6 years ago

Scenario: I first use function intersection to find the intersection points between a line and the boundary of a polygon. Then, using these intersection points to split boundary and line. However, due to the floating point calculation issues, function split_line_with_point cannot work properly, i.e the intersection points are not within the line.

A possible way to fix this issue is to change the criteria:

Change line 377 in ops.py if not line.relate_pattern(splitter, '0********'): to if (line.distance(splitter)>1e-8):

In this way, intersection points can be used to split line or boundary of a polygon.

danvk commented 6 years ago

This tripped me up today. Here's a concrete example:

> a.wkt
'LINESTRING (1766.36403451032 307.2453045763735, 1807.216199655729 331.480621841313)'
> b.wkt
'LINESTRING (1501.061748537531 586.9106713524088, 1743.485916080718 383.0477914251387, 1751.524430191241 369.4980262969621, 1815.0568436799 262.4061723989435, 2101.901143337514 -221.1336479363963)'
> a.intersects(b)
True
> pt = a.intersection(b)
> [*pt.coords]
[(1782.7049005687165, 316.9394314824876)]

So far so good. But then…

> a.intersects(pt), b.intersects(pt)
(False, False)
> len(shapely.ops.split(b, pt))
1

i.e. shapely does not consider the intersection point to be on the line and cannot use it for splitting.

Andreachen0707 commented 6 years ago

This tripped me up today. Here's a concrete example:

> a.wkt
'LINESTRING (1766.36403451032 307.2453045763735, 1807.216199655729 331.480621841313)'
> b.wkt
'LINESTRING (1501.061748537531 586.9106713524088, 1743.485916080718 383.0477914251387, 1751.524430191241 369.4980262969621, 1815.0568436799 262.4061723989435, 2101.901143337514 -221.1336479363963)'
> a.intersects(b)
True
> pt = a.intersection(b)
> [*pt.coords]
[(1782.7049005687165, 316.9394314824876)]

So far so good. But then…

> a.intersects(pt), b.intersects(pt)
(False, False)
> len(shapely.ops.split(b, pt))
1

i.e. shapely does not consider the intersection point to be on the line and cannot use it for splitting.

Thank you so much for providing this example! Hope they can solve it in the future. Does my way help you solve the issue?

danvk commented 6 years ago

I copy/pasted this code into a tolerant_split method (skipping the check you mentioned). It seems to be working well!

jujaro commented 3 years ago

Got the same issue with a simple example:

> from shapely.geometry import LineString
> from shapely.ops import split
> total_geom = LineString([(0, 0), (1, 1), (2, 1), (3, 0)])
> point = total_geom.interpolate(2 / 3, normalized=True)
> print(split(total_geom, point))

GEOMETRYCOLLECTION (LINESTRING (0 0, 1 1, 2 1, 3 0))

I guess adding some tolerance would help? sort of instead of:

if not line.relate_pattern(splitter, '0********'):

something like:

if (line.distance(splitter)>tolerance):