252 lines
8.7 KiB
Python
252 lines
8.7 KiB
Python
import ctypes
|
|
import platform
|
|
import os
|
|
import math
|
|
from typing import Any
|
|
from .s0_base import get_node_links, get_link_nodes, is_pipe
|
|
from .s5_pipes import get_pipe
|
|
from .database import read, try_read, read_all, write
|
|
from .s24_coordinates import node_has_coord, get_node_coord
|
|
|
|
|
|
def from_postgis_polygon(polygon: str) -> list[tuple[float, float]]:
|
|
boundary = polygon.lower().removeprefix('polygon((').removesuffix('))').split(',')
|
|
xys = []
|
|
for pt in boundary:
|
|
xy = pt.split(' ')
|
|
xys.append((float(xy[0]), float(xy[1])))
|
|
return xys
|
|
|
|
|
|
def to_postgis_polygon(boundary: list[tuple[float, float]]) -> str:
|
|
polygon = ''
|
|
for pt in boundary:
|
|
polygon += f'{pt[0]} {pt[1]},'
|
|
return f'polygon(({polygon[:-1]}))'
|
|
|
|
|
|
def get_nodes_in_boundary(name: str, boundary: list[tuple[float, float]]) -> list[str]:
|
|
api = 'get_nodes_in_boundary'
|
|
write(name, f"delete from temp_region where id = '{api}'")
|
|
write(name, f"insert into temp_region (id, boundary) values ('{api}', '{to_postgis_polygon(boundary)}')")
|
|
|
|
nodes: list[str] = []
|
|
for row in read_all(name, f"select c.node from coordinates as c, temp_region as r where ST_Intersects(c.coord, r.boundary) and r.id = '{api}'"):
|
|
nodes.append(row['node'])
|
|
|
|
write(name, f"delete from temp_region where id = '{api}'")
|
|
|
|
return nodes
|
|
|
|
|
|
def get_nodes_in_region(name: str, region_id: str) -> list[str]:
|
|
nodes: list[str] = []
|
|
for row in read_all(name, f"select c.node from coordinates as c, region as r where ST_Intersects(c.coord, r.boundary) and r.id = '{region_id}'"):
|
|
nodes.append(row['node'])
|
|
return nodes
|
|
|
|
|
|
def calculate_convex_hull(name: str, nodes: list[str]) -> list[tuple[float, float]]:
|
|
write(name, f'delete from temp_node')
|
|
for node in nodes:
|
|
write(name, f"insert into temp_node values ('{node}')")
|
|
|
|
# TODO: check none
|
|
polygon = read(name, f'select st_astext(st_convexhull(st_collect(array(select coord from coordinates where node in (select * from temp_node))))) as boundary' )['boundary']
|
|
write(name, f'delete from temp_node')
|
|
|
|
return from_postgis_polygon(polygon)
|
|
|
|
|
|
def _verify_platform():
|
|
_platform = platform.system()
|
|
if _platform != "Windows":
|
|
raise Exception(f'Platform {_platform} unsupported (not yet)')
|
|
|
|
|
|
def _normal(v: tuple[float, float]) -> tuple[float, float]:
|
|
l = math.sqrt(v[0] * v[0] + v[1] * v[1])
|
|
return (v[0] / l, v[1] / l)
|
|
|
|
|
|
def _angle(v: tuple[float, float]) -> float:
|
|
if v[0] >= 0 and v[1] >= 0:
|
|
return math.asin(v[1])
|
|
elif v[0] <= 0 and v[1] >= 0:
|
|
return math.pi - math.asin(v[1])
|
|
elif v[0] <= 0 and v[1] <= 0:
|
|
return math.asin(-v[1]) + math.pi
|
|
elif v[0] >= 0 and v[1] <= 0:
|
|
return math.pi * 2 - math.asin(-v[1])
|
|
return 0
|
|
|
|
|
|
def _angle_of_node_link(node: str, link: str, nodes, links) -> float:
|
|
n1 = node
|
|
n2 = links[link]['node1'] if n1 == links[link]['node2'] else links[link]['node2']
|
|
x1, y1 = nodes[n1]['x'], nodes[n1]['y']
|
|
x2, y2 = nodes[n2]['x'], nodes[n2]['y']
|
|
v = _normal((x2 - x1, y2 - y1))
|
|
return _angle(v)
|
|
|
|
|
|
class Topology:
|
|
def __init__(self, db: str, nodes: list[str]) -> None:
|
|
self._nodes: dict[str, Any] = {}
|
|
self._max_x_node = ''
|
|
self._node_list: list[str] = []
|
|
for node in nodes:
|
|
if not node_has_coord(db, node):
|
|
continue
|
|
if get_node_links(db, node) == 0:
|
|
continue
|
|
self._nodes[node] = get_node_coord(db, node) | { 'links': [] }
|
|
self._node_list.append(node)
|
|
if self._max_x_node == '' or self._nodes[node]['x'] > self._nodes[self._max_x_node]['x']:
|
|
self._max_x_node = node
|
|
|
|
self._links = {}
|
|
self._link_list: list[str] = []
|
|
for node in self._nodes:
|
|
for link in get_node_links(db, node):
|
|
candidate = True
|
|
link_nodes = get_link_nodes(db, link)
|
|
for link_node in link_nodes:
|
|
if link_node not in self._nodes:
|
|
candidate = False
|
|
break
|
|
if candidate:
|
|
length = get_pipe(db, link)['length'] if is_pipe(db, link) else 0.0
|
|
self._links[link] = { 'node1' : link_nodes[0], 'node2' : link_nodes[1], 'length' : length }
|
|
self._link_list.append(link)
|
|
if link not in self._nodes[link_nodes[0]]['links']:
|
|
self._nodes[link_nodes[0]]['links'].append(link)
|
|
if link not in self._nodes[link_nodes[1]]['links']:
|
|
self._nodes[link_nodes[1]]['links'].append(link)
|
|
|
|
def nodes(self):
|
|
return self._nodes
|
|
|
|
def node_list(self):
|
|
return self._node_list
|
|
|
|
def max_x_node(self):
|
|
return self._max_x_node
|
|
|
|
def links(self):
|
|
return self._links
|
|
|
|
def link_list(self):
|
|
return self._link_list
|
|
|
|
|
|
def calculate_boundary(name: str, nodes: list[str]) -> list[tuple[float, float]]:
|
|
topology = Topology(name, nodes)
|
|
t_nodes = topology.nodes()
|
|
t_links = topology.links()
|
|
|
|
cursor = topology.max_x_node()
|
|
in_angle = 0
|
|
|
|
paths: list[str] = []
|
|
while True:
|
|
#print(cursor)
|
|
paths.append(cursor)
|
|
sorted_links = []
|
|
overlapped_link = ''
|
|
for link in t_nodes[cursor]['links']:
|
|
angle = _angle_of_node_link(cursor, link, t_nodes, t_links)
|
|
if angle == in_angle:
|
|
overlapped_link = link
|
|
continue
|
|
sorted_links.append((angle, link))
|
|
|
|
# work into a branch, return
|
|
if len(sorted_links) == 0:
|
|
cursor = paths[-2]
|
|
if cursor == topology.max_x_node():
|
|
paths.append(cursor)
|
|
break
|
|
in_angle = in_angle = _angle_of_node_link(cursor, overlapped_link, t_nodes, t_links)
|
|
continue
|
|
|
|
sorted_links = sorted(sorted_links, key=lambda s:s[0])
|
|
out_link = sorted_links[0][1]
|
|
for angle, link in sorted_links:
|
|
if angle > in_angle:
|
|
out_link = link
|
|
break
|
|
|
|
cursor = t_links[out_link]['node1'] if cursor == t_links[out_link]['node2'] else t_links[out_link]['node2']
|
|
# end up trip :)
|
|
if cursor == topology.max_x_node():
|
|
paths.append(cursor)
|
|
break
|
|
|
|
in_angle = _angle_of_node_link(cursor, out_link, t_nodes, t_links)
|
|
|
|
boundary: list[tuple[float, float]] = []
|
|
for node in paths:
|
|
boundary.append((t_nodes[node]['x'], t_nodes[node]['y']))
|
|
|
|
return boundary
|
|
|
|
|
|
'''
|
|
# CClipper2.dll
|
|
# int inflate_paths(double* path, size_t size, double delta, int jt, int et, double miter_limit, int precision, double arc_tolerance, double** out_path, size_t* out_size);
|
|
# int simplify_paths(double* path, size_t size, double epsilon, int is_closed_path, double** out_path, size_t* out_size);
|
|
# void free_paths(double** paths);
|
|
'''
|
|
def inflate_boundary(name: str, boundary: list[tuple[float, float]], delta: float = 0.5) -> list[tuple[float, float]]:
|
|
if boundary[0] == boundary[-1]:
|
|
del(boundary[-1])
|
|
|
|
lib = ctypes.CDLL(os.path.join(os.getcwd(), 'api', 'CClipper2.dll'))
|
|
|
|
c_size = ctypes.c_size_t(len(boundary) * 2)
|
|
c_path = (ctypes.c_double * c_size.value)()
|
|
i = 0
|
|
for xy in boundary:
|
|
c_path[i] = xy[0]
|
|
i += 1
|
|
c_path[i] = xy[1]
|
|
i += 1
|
|
c_delta = ctypes.c_double(delta)
|
|
JoinType_Square, JoinType_Round, JoinType_Miter = 0, 1, 2
|
|
c_jt = ctypes.c_int(JoinType_Square)
|
|
EndType_Polygon, EndType_Joined, EndType_Butt, EndType_Square, EndType_Round = 0, 1, 2, 3, 4
|
|
c_et = ctypes.c_int(EndType_Polygon)
|
|
c_miter_limit = ctypes.c_double(2.0)
|
|
c_precision = ctypes.c_int(2)
|
|
c_arc_tolerance = ctypes.c_double(0.0)
|
|
c_out_path = ctypes.POINTER(ctypes.c_double)()
|
|
c_out_size = ctypes.c_size_t(0)
|
|
|
|
lib.inflate_paths(c_path, c_size, c_delta, c_jt, c_et, c_miter_limit, c_precision, c_arc_tolerance, ctypes.byref(c_out_path), ctypes.byref(c_out_size))
|
|
if c_out_size.value == 0:
|
|
lib.free_paths(ctypes.byref(c_out_path))
|
|
return []
|
|
|
|
# TODO: simplify_paths :)
|
|
|
|
result: list[tuple[float, float]] = []
|
|
for i in range(0, c_out_size.value, 2):
|
|
result.append((c_out_path[i], c_out_path[i + 1]))
|
|
result.append(result[0])
|
|
|
|
lib.free_paths(ctypes.byref(c_out_path))
|
|
return result
|
|
|
|
|
|
def inflate_region(name: str, region_id: str, delta: float = 0.5) -> list[tuple[float, float]]:
|
|
r = try_read(name, f"select id, st_astext(boundary) as boundary_geom from region where id = '{region_id}'")
|
|
if r == None:
|
|
return []
|
|
boundary = from_postgis_polygon(str(r['boundary_geom']))
|
|
return inflate_boundary(name, boundary, delta)
|
|
|
|
|
|
if __name__ == '__main__':
|
|
_verify_platform()
|