Commit efa83459 authored by Hoek, Steven's avatar Hoek, Steven
Browse files

Algorithm in method get_nearest_edge_points improved. Warnings were built in...

Algorithm in method get_nearest_edge_points improved. Warnings were built in because found edge points may not actually lie on the edge nearest to the point.
parent 0b72648f
"""VectorLib - a Python library for representing and handling vector based spatial data"""
from math import sqrt
from math import sqrt, copysign
from delaunay_triangulation.triangulate import delaunay
from delaunay_triangulation import typing
import warnings
def distance(xy1, xy2):
# Be prepared to handle sequenced coordinates as well as x & y attributes
......@@ -87,6 +88,12 @@ class SimplePoint(object):
def __init__(self, x, y):
self.x = x
self.y = y
def __str__(self):
return "(" + str(self.x) + ", " + str(self.y) + ")"
def __repr__(self):
return self.__str__()
class SimpleShape(object):
# TODO: things may not work in case of island polygons
......@@ -118,7 +125,7 @@ class SimpleShape(object):
_len = len(x_coords)
return [sum(x_coords)/_len, sum(y_coords)/_len]
def get_centroid(self):
def get_centroid(self) -> typing.Vertex:
result = (-1.0, -1.0)
try:
# Get all the points and then do the triangulation
......@@ -135,7 +142,7 @@ class SimpleShape(object):
total_area += area
total_xweighted += centroid[0] * area
total_yweighted += centroid[1] * area
result = (total_xweighted / total_area, total_yweighted / total_area)
result = typing.Vertex(total_xweighted / total_area, total_yweighted / total_area)
except Exception as e:
print(e)
......@@ -161,7 +168,9 @@ class SimpleShape(object):
# Check that the result lies on the edge of the polygon
A = (min(pt1[0], pt2[0])) <= result.x and (max(pt1[0], pt2[0]) >= result.x)
B = (min(pt1[1], pt2[1])) <= result.y and (max(pt1[1], pt2[1]) >= result.y)
if (not A) or (not B): result = typing.Vertex(pt1[0], pt1[1])
if (not A) and (not B):
# The intersection point is on the line but outside the x- and y-range
result = typing.Vertex(pt1[0], pt1[1]) # Still nearest point
return result
def get_nearest_edge_points(self, otherpoint) -> (typing.Vertex, typing.Vertex):
......@@ -169,11 +178,29 @@ class SimpleShape(object):
if hasattr(otherpoint, "__len__"): x0, y0 = otherpoint[0], otherpoint[1]
elif hasattr(otherpoint, "x") and hasattr(otherpoint, "y"): x0, y0 = otherpoint.x, otherpoint.y
# Get the 2 points of the polygon that are nearest to the given point
# Get 2 points of the polygon that are near to the given point, on both sides
points = sum(self.get_point_series(), [])
mydicts = [{"point":(p.x, p.y), "distance":distance((p.x, p.y), (x0, y0))} for p in points]
n = len(points)
mydicts = [{"point":(p.x, p.y), "distance":distance((p.x,p.y), (x0,y0)), "index":i} for i, p in zip(range(n), points)]
sorted_points = sorted(mydicts, key=lambda p: p['distance'])
pt1, pt2 = sorted_points[0]["point"], sorted_points[1]["point"]
# We assume that this is the nearest point, but there may even be a line through points further away that runs closer by:
spt1 = sorted_points[0]
pt1 = spt1["point"]
# The 2 nearest points are not necessarily the right ones. Which line starting / ending in pt1 runs closest by?
for j in [-1, 1]:
spt2 = list(filter(lambda p: p["index"] == spt1["index"]+j, mydicts))[0]
pt2 = spt2["point"]
line = typing.get_standard_line(typing.Vertex(pt1[0], pt1[1]), typing.Vertex(pt2[0], pt2[1]))
if line.A != 0:
plumb_line = self.get_plumb_line(line, otherpoint)
edgepoint = typing.get_line_intersection(line, plumb_line)
else:
edgepoint = typing.Vertex(x0, pt1[1])
A = (min(pt1[0], pt2[0])) <= edgepoint.x and (max(pt1[0], pt2[0]) >= edgepoint.x)
B = (min(pt1[1], pt2[1])) <= edgepoint.y and (max(pt1[1], pt2[1]) >= edgepoint.y)
if A and B: break
return typing.Vertex(pt1[0], pt1[1]), typing.Vertex(pt2[0], pt2[1])
def get_plumb_line(self, line, point) -> typing.StandardLine:
......@@ -240,23 +267,23 @@ class SimpleShape(object):
eps = 0.0000000001
if abs(pt1.x - pt2.x) < eps and abs(pt1.y - pt2.y) < eps: # TODO: find a good way to handle this
raise ValueError("Two consecutive points of the shape are located too close to each other!")
A = (min(pt1[0], pt2[0])) <= x and (max(pt1[0], pt2[0]) >= x)
B = (min(pt1[1], pt2[1])) <= y and (max(pt1[1], pt2[1]) >= y)
if A and B:
# At least the edgepoint is located inside the boundary box around pt1 and pt2
# Make sure a point inside the shape is found that is distance away from the edge
line1 = typing.get_standard_line(typing.Vertex(pt1[0], pt1[1]), typing.Vertex(pt2[0], pt2[1]))
if line1.A == 0:
if self.contains_point((x, y + distance)): result = typing.Vertex(x, y + distance)
else: result = typing.Vertex(x, y - distance)
else:
line2 = self.get_plumb_line(line1, edgepoint)
slope = line2.A / line2.B
dx = distance / sqrt(1 + slope**2)
dy = slope * dx
if self.contains_point((x + dx, y + dy)): result = typing.Vertex(x + dx, y + dy)
else: result = typing.Vertex(x - dx, y - dy)
else: raise ValueError("Given point is not located on the edge!")
A = (min(pt1[0], pt2[0]) <= x) and (max(pt1[0], pt2[0]) >= x)
B = (min(pt1[1], pt2[1]) <= y) and (max(pt1[1], pt2[1]) >= y)
if not (A and B): warnings.warn("It seems that the given point is not located on the nearest edge!")
# Make sure a point inside the shape is found that is distance away from the edge
line1 = typing.get_standard_line(typing.Vertex(pt1[0], pt1[1]), typing.Vertex(pt2[0], pt2[1]))
if line1.A == 0:
if self.contains_point((x, y + distance)): result = typing.Vertex(x, y + distance)
else: result = typing.Vertex(x, y - distance)
else:
line2 = self.get_plumb_line(line1, edgepoint)
slope = line2.A / line2.B
dx = distance / sqrt(1 + slope**2)
dy = slope * dx
if self.contains_point((x + dx, y + dy)): result = typing.Vertex(x + dx, y + dy)
else: result = typing.Vertex(x - dx, y - dy)
except Exception as e:
print(e)
finally:
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment