2025-summer/team-3
12
Fork 0
Code Issues Pull requests Projects Releases Activity

simplify cut_polygon_with_line()

Browse source
This commit is contained in:
Laura Klünder 2017-11-13 17:09:26 +01:00
parent 47ff4fe788
commit bcdc77d2b2
1 changed files with 73 additions and 66 deletions
Download patch file Download diff file Expand all files Collapse all files

139
src/c3nav/mapdata/utils/geometry.py
View file

@ -154,77 +154,84 @@ def cut_polygon_with_line(polygon: Polygon, line: LineString):
last = points.popleft() last = points.popleft()
while points: while points:
current = points.popleft() current = points.popleft()
if current.polygon == last.polygon:
polygon = polygons[current.polygon]
segment = cut_line_with_point(cut_line_with_point(line, last.point)[-1], current.point)[0]
if current.ring != last.ring:
# connect rings
ring1 = cut_line_with_point(polygon[last.ring], last.point)
ring2 = cut_line_with_point(polygon[current.ring], current.point)
new_ring = LinearRing(ring1[1].coords[:-1] + ring1[0].coords[:-1] + segment.coords[:-1] +
ring2[1].coords[:-1] + ring2[0].coords[:-1] + segment.coords[::-1])
if current.ring == 0 or last.ring == 0:
new_i = 0
polygon[0] = new_ring
interior = current.ring if last.ring == 0 else last.ring
polygon[interior] = None
mapping = {interior: new_i}
else:
new_i = len(polygon)
mapping = {last.ring: new_i, current.ring: new_i}
polygon.append(new_ring)
polygon[last.ring] = None
polygon[current.ring] = None
points = deque((cutpoint(item.point, item.polygon, mapping[item.ring]) # don't to anything between different polygons
if (item.polygon == current.polygon and item.ring in mapping) else item) if current.polygon != last.polygon:
for item in points) last = current
current = cutpoint(current.point, current.polygon, new_i) continue
elif current.ring == 0: polygon = polygons[current.polygon]
# cut polygon segment = cut_line_with_point(cut_line_with_point(line, last.point)[-1], current.point)[0]
new_i = len(polygons) if current.ring != last.ring:
exterior = cut_line_with_point(polygon[0], current.point) # connect rings
exterior = cut_line_with_point(LinearRing(exterior[1].coords[:-1] + ring1 = cut_line_with_point(polygon[last.ring], last.point)
exterior[0].coords[0:]), last.point) ring2 = cut_line_with_point(polygon[current.ring], current.point)
new_ring = LinearRing(ring1[1].coords[:-1] + ring1[0].coords[:-1] + segment.coords[:-1] +
ring2[1].coords[:-1] + ring2[0].coords[:-1] + segment.coords[::-1])
if current.ring == 0 or last.ring == 0:
new_i = 0
polygon[0] = new_ring
interior = current.ring if last.ring == 0 else last.ring
polygon[interior] = None
mapping = {interior: new_i}
else:
new_i = len(polygon)
mapping = {last.ring: new_i, current.ring: new_i}
polygon.append(new_ring)
polygon[last.ring] = None
polygon[current.ring] = None
point_forwards = exterior[1].coords[1] points = deque((cutpoint(item.point, item.polygon, mapping[item.ring])
point_backwards = exterior[0].coords[-2] if (item.polygon == current.polygon and item.ring in mapping) else item)
angle_forwards = math.atan2(point_forwards[0] - last.point.x, point_forwards[1] - last.point.y) for item in points)
angle_backwards = math.atan2(point_backwards[0] - last.point.x, point_backwards[1] - last.point.y) last = cutpoint(current.point, current.polygon, new_i)
angle_segment = math.atan2(current.point.x - last.point.x, current.point.y - last.point.y) continue
while angle_forwards <= angle_backwards: # cut polygon?
angle_forwards += 2*math.pi new_i = len(polygons)
if angle_segment < angle_backwards: exterior = cut_line_with_point(polygon[0], current.point)
while angle_segment < angle_backwards: exterior = cut_line_with_point(LinearRing(exterior[1].coords[:-1] +
angle_segment += 2*math.pi exterior[0].coords[0:]), last.point)
else:
while angle_segment > angle_forwards:
angle_segment -= 2*math.pi
# cut polygon only if segment is inside polygon point_forwards = exterior[1].coords[1]
if angle_backwards < angle_segment < angle_forwards: point_backwards = exterior[0].coords[-2]
exterior1 = LinearRing(exterior[0].coords[:-1] + segment.coords[0:]) angle_forwards = math.atan2(point_forwards[0] - last.point.x, point_forwards[1] - last.point.y)
exterior2 = LinearRing(exterior[1].coords[:-1] + segment.coords[::-1]) angle_backwards = math.atan2(point_backwards[0] - last.point.x, point_backwards[1] - last.point.y)
geom = Polygon(exterior1) angle_segment = math.atan2(current.point.x - last.point.x, current.point.y - last.point.y)
polygon[0] = exterior1
new_polygon = [exterior2]
polygons.append(new_polygon)
mapping = {}
for i, interior in enumerate(polygon[1:]):
if interior is not None and not geom.contains(interior):
mapping[i] = len(new_polygon)
new_polygon.append(interior)
points = deque((cutpoint(item.point, new_i, mapping[item.ring]) while angle_forwards <= angle_backwards:
if (item.polygon == current.polygon and item.ring in mapping) else item) angle_forwards += 2*math.pi
for item in points) if angle_segment < angle_backwards:
points = deque((cutpoint(item.point, new_i, 0) while angle_segment < angle_backwards:
if (item.polygon == current.polygon and item.ring == 0 and angle_segment += 2*math.pi
not exterior1.contains(item.point)) else item) else:
for item in points) while angle_segment > angle_forwards:
current = cutpoint(current.point, new_i, 0) angle_segment -= 2*math.pi
last = current
# don't cut polygon if segment is not inside polygon (= we don't cut through emptyness)
if not (angle_backwards < angle_segment < angle_forwards):
last = current
continue
exterior1 = LinearRing(exterior[0].coords[:-1] + segment.coords[0:])
exterior2 = LinearRing(exterior[1].coords[:-1] + segment.coords[::-1])
geom = Polygon(exterior1)
polygon[0] = exterior1
new_polygon = [exterior2]
polygons.append(new_polygon)
mapping = {}
for i, interior in enumerate(polygon[1:]):
if interior is not None and not geom.contains(interior):
mapping[i] = len(new_polygon)
new_polygon.append(interior)
points = deque((cutpoint(item.point, new_i, mapping[item.ring])
if (item.polygon == current.polygon and item.ring in mapping) else item)
for item in points)
points = deque((cutpoint(item.point, new_i, 0)
if (item.polygon == current.polygon and item.ring == 0 and
not exterior1.contains(item.point)) else item)
for item in points)
last = cutpoint(current.point, new_i, 0)
return tuple(Polygon(polygon[0], tuple(ring for ring in polygon[1:] if ring is not None)) return tuple(Polygon(polygon[0], tuple(ring for ring in polygon[1:] if ring is not None))
for polygon in polygons) for polygon in polygons)