From 2b424e47100b43ced8026ae5d38ffaccd05ad020 Mon Sep 17 00:00:00 2001 From: "WQY\\qiong" Date: Fri, 9 Jun 2023 09:36:11 +0800 Subject: [PATCH] Improve boundary algorithm --- api/s32_region_util.py | 183 ++++++++++++++++++++++++++++++++++++++--- test_tjnetwork.py | 34 ++++---- 2 files changed, 188 insertions(+), 29 deletions(-) diff --git a/api/s32_region_util.py b/api/s32_region_util.py index 6a4f6cc..d5ec1c0 100644 --- a/api/s32_region_util.py +++ b/api/s32_region_util.py @@ -25,6 +25,13 @@ def to_postgis_polygon(boundary: list[tuple[float, float]]) -> str: return f'polygon(({polygon[:-1]}))' +def to_postgis_linestring(boundary: list[tuple[float, float]]) -> str: + line = '' + for pt in boundary: + line += f'{pt[0]} {pt[1]},' + return f'linestring({line[:-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}'") @@ -105,7 +112,7 @@ class Topology: 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._links: dict[str, Any] = {} self._link_list: list[str] = [] for node in self._nodes: for link in get_node_links(db, node): @@ -140,25 +147,21 @@ class Topology: 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() +def _calculate_boundary(name: str, cursor: str, t_nodes: dict[str, Any], t_links: dict[str, Any]) -> tuple[list[str], dict[str, list[str]], list[tuple[float, float]]]: in_angle = 0 - paths: list[str] = [] + vertices: list[str] = [] + path: dict[str, list[str]] = {} while True: # prevent duplicated node - if len(paths) > 0 and cursor == paths[-1]: + if len(vertices) > 0 and cursor == vertices[-1]: break # prevent duplicated path - if len(paths) >= 3 and paths[0] == paths[-1] and paths[1] == cursor: + if len(vertices) >= 3 and vertices[0] == vertices[-1] and vertices[1] == cursor: break - paths.append(cursor) + vertices.append(cursor) sorted_links = [] overlapped_link = '' @@ -171,7 +174,8 @@ def calculate_boundary(name: str, nodes: list[str]) -> list[tuple[float, float]] # work into a branch, return if len(sorted_links) == 0: - cursor = paths[-2] + path[overlapped_link] = [] + cursor = vertices[-2] in_angle = _angle_of_node_link(cursor, overlapped_link, t_nodes, t_links) continue @@ -182,13 +186,166 @@ def calculate_boundary(name: str, nodes: list[str]) -> list[tuple[float, float]] out_link = link break + path[out_link] = [] cursor = t_links[out_link]['node1'] if cursor == t_links[out_link]['node2'] else t_links[out_link]['node2'] in_angle = _angle_of_node_link(cursor, out_link, t_nodes, t_links) boundary: list[tuple[float, float]] = [] - for node in paths: + for node in vertices: boundary.append((t_nodes[node]['x'], t_nodes[node]['y'])) + return (vertices, path, boundary) + + +def _collect_new_links(in_links: dict[str, list[str]], t_nodes: dict[str, Any], t_links: dict[str, Any], new_nodes: dict[str, Any], new_links: dict[str, Any]) -> tuple[dict[str, Any], dict[str, Any]]: + for link, pts in in_links.items(): + node1 = t_links[link]['node1'] + node2 = t_links[link]['node2'] + x1, x2 = t_nodes[node1]['x'], t_nodes[node2]['x'] + y1, y2 = t_nodes[node1]['y'], t_nodes[node2]['y'] + + if node1 not in new_nodes: + new_nodes[node1] = { 'x': x1, 'y': y1, 'links': [] } + if node2 not in new_nodes: + new_nodes[node2] = { 'x': x2, 'y': y2, 'links': [] } + + x_delta = x2 - x1 + y_delta = y2 - y1 + use_x = abs(x_delta) > abs(y_delta) + + if len(pts) == 0: + new_links[link] = t_links[link] + else: + sorted_nodes: list[tuple[float, str]] = [] + sorted_nodes.append((0.0, node1)) + sorted_nodes.append((1.0, node2)) + i = 0 + for pt in pts: + x, y = new_nodes[pt]['x'], new_nodes[pt]['y'] + percent = ((x - x1) / x_delta) if use_x else ((y - y1) / y_delta) + sorted_nodes.append((percent, pt)) + i += 1 + sorted_nodes = sorted(sorted_nodes, key=lambda s:s[0]) + + for i in range(1, len(sorted_nodes)): + l = sorted_nodes[i - 1][1] + r = sorted_nodes[i][1] + new_link = f'LINK_[{l}]_[{r}]' + new_links[new_link] = { 'node1': l, 'node2': r } + + return (new_nodes, new_links) + + +def calculate_boundary(name: str, nodes: list[str]) -> list[tuple[float, float]]: + topology = Topology(name, nodes) + t_nodes = topology.nodes() + t_links = topology.links() + + vertices, path, boundary = _calculate_boundary(name, topology.max_x_node(), t_nodes, t_links) + + #return boundary + + api = 'calculate_boundary' + write(name, f"delete from temp_region where id = '{api}'") + # use linestring instead of polygon to reduce strict limitation + write(name, f"insert into temp_region (id, boundary) values ('{api}', '{to_postgis_linestring(boundary)}')") + + write(name, f'delete from temp_node') + for node in nodes: + write(name, f"insert into temp_node values ('{node}')") + + for row in read_all(name, f"select n.node from coordinates as c, temp_node as n, temp_region as r where c.node = n.node and ST_Intersects(c.coord, r.boundary) and r.id = '{api}'"): + node = row['node'] + write(name, f"delete from temp_node where node = '{node}'") + + outside_nodes: list[str] = [] + for row in read_all(name, "select node from temp_node"): + outside_nodes.append(row['node']) + + # no outside nodes, return + if len(outside_nodes) == 0: + write(name, f'delete from temp_node') + write(name, f"delete from temp_region where id = '{api}'") + return boundary + + new_nodes: dict[str, Any] = {} + new_links: dict[str, Any] = {} + + boundary_links: dict[str, list[str]] = {} + write(name, "delete from temp_link_2") + for node in outside_nodes: + for link in t_nodes[node]['links']: + node1 = t_links[link]['node1'] + node2 = t_links[link]['node2'] + if node1 in outside_nodes and node2 not in outside_nodes and node2 not in vertices and link: + if link not in boundary: + boundary_links[link] = [] + line = f"LINESTRING({t_nodes[node1]['x']} {t_nodes[node1]['y']}, {t_nodes[node2]['x']} {t_nodes[node2]['y']})" + write(name, f"insert into temp_link_2 values ('{link}', '{line}')") + if node2 in outside_nodes and node1 not in outside_nodes and node1 not in vertices: + if link not in boundary: + boundary_links[link] = [] + line = f"LINESTRING({t_nodes[node1]['x']} {t_nodes[node1]['y']}, {t_nodes[node2]['x']} {t_nodes[node2]['y']})" + write(name, f"insert into temp_link_2 values ('{link}', '{line}')") + if node1 in outside_nodes and node2 in outside_nodes: + x1, x2 = t_nodes[node1]['x'], t_nodes[node2]['x'] + y1, y2 = t_nodes[node1]['y'], t_nodes[node2]['y'] + if node1 not in new_nodes: + new_nodes[node1] = { 'x': x1, 'y': y1, 'links': [] } + if node2 not in new_nodes: + new_nodes[node2] = { 'x': x2, 'y': y2, 'links': [] } + if link not in new_links: + new_links[link] = t_links[link] + + # no boundary links, return + if len(boundary_links) == 0: + write(name, "delete from temp_link_2") + write(name, f'delete from temp_node') + write(name, f"delete from temp_region where id = '{api}'") + return boundary + + write(name, "delete from temp_link_1") + for link, _ in path.items(): + node1 = t_links[link]['node1'] + node2 = t_links[link]['node2'] + line = f"LINESTRING({t_nodes[node1]['x']} {t_nodes[node1]['y']}, {t_nodes[node2]['x']} {t_nodes[node2]['y']})" + write(name, f"insert into temp_link_1 (link, geom) values ('{link}', '{line}')") + + has_intersection = False + for row in read_all(name, f"select l1.link as l, l2.link as r, st_astext(st_intersection(l1.geom, l2.geom)) as p from temp_link_1 as l1, temp_link_2 as l2 where st_intersects(l1.geom, l2.geom)"): + has_intersection = True + + link1, link2, pt = str(row['l']), str(row['r']), str(row['p']) + pts = pt.lower().removeprefix('point(').removesuffix(')').split(' ') + xy = (float(pts[0]), float(pts[1])) + + new_node = f'NODE_[{link1}]_[{link2}]' + new_nodes[new_node] = { 'x': xy[0], 'y': xy[1], 'links': [] } + + path[link1].append(new_node) + boundary_links[link2].append(new_node) + + if not has_intersection: + write(name, "delete from temp_link_1") + write(name, "delete from temp_link_2") + write(name, 'delete from temp_node') + write(name, f"delete from temp_region where id = '{api}'") + return boundary + + new_nodes, new_links = _collect_new_links(path, t_nodes, t_links, new_nodes, new_links) + new_nodes, new_links = _collect_new_links(boundary_links, t_nodes, t_links, new_nodes, new_links) + + for link, values in new_links.items(): + new_nodes[values['node1']]['links'].append(link) + new_nodes[values['node2']]['links'].append(link) + + _, _, boundary = _calculate_boundary(name, topology.max_x_node(), new_nodes, new_links) + + write(name, "delete from temp_link_1") + write(name, "delete from temp_link_2") + write(name, 'delete from temp_node') + write(name, f"delete from temp_region where id = '{api}'") + return boundary diff --git a/test_tjnetwork.py b/test_tjnetwork.py index 108e778..1f7e717 100644 --- a/test_tjnetwork.py +++ b/test_tjnetwork.py @@ -5832,9 +5832,9 @@ class TestApi: dmas = calculate_district_metering_area_for_nodes(p, get_nodes(p), 3) assert len(dmas) == 3 - assert dmas[0] == ['173', '184', '185', '199', '2', '201', '203', '205', '206', '207', '208', '209', '211', '213', '215', '217', '219', '225', '229', '231', '237', '239', '241', '243', '247', '249', '251', '253', '255', '273', '275', '50'] - assert dmas[1] == ['1', '10', '101', '103', '109', '111', '113', '161', '163', '164', '166', '167', '169', '171', '179', '181', '183', '187', '189', '191', '193', '195', '197', '204', '265', '267', '269', '271', '35', '40', 'Lake'] - assert dmas[2] == ['105', '107', '115', '117', '119', '120', '121', '123', '125', '127', '129', '131', '139', '141', '143', '145', '147', '149', '15', '151', '153', '157', '159', '20', '257', '259', '261', '263', '3', '60', '601', '61', 'River'] + assert dmas[0] == ['117', '119', '120', '121', '123', '125', '127', '129', '131', '139', '141', '143', '145', '147', '149', '15', '151', '153', '157', '159', '161', '195', '20', '257', '259', '261', '263', '3', '60', '601', '61', 'River'] + assert dmas[1] == ['1', '163', '164', '166', '167', '169', '171', '173', '177', '179', '181', '183', '184', '185', '187', '189', '191', '193', '199', '201', '203', '204', '205', '207', '265', '267', '269', '271', '273', '275', '35', '40'] + assert dmas[2] == ['10', '101', '103', '105', '107', '109', '111', '113', '115', '197', '2', '206', '208', '209', '211', '213', '215', '217', '219', '225', '229', '231', '237', '239', '241', '243', '247', '249', '251', '253', '255', '50', 'Lake'] self.leave(p) @@ -5844,15 +5844,15 @@ class TestApi: read_inp(p, f'./inp/net3.inp', '3') open_project(p) - assert get_node_coord(p, '177') == { 'x' : 0.0, 'y' : 0.0 } - add_region(p, ChangeSet({'id': 'r', 'boundary': [(-10000.0, -10000.0), (10000.0, -10000.0), (10000.0, 10000.0), (-10000.0, 10000.0), (-10000.0, -10000.0)]})) nodes = get_nodes_in_region(p, 'r') assert len(nodes) == 97 dmas = calculate_district_metering_area_for_region(p, 'r', 3) - assert dmas == [['10', '60', '601', '61', '101', '103', '105', '107', '109', '111', '113', '115', '117', '119', '120', '121', '123', '125', '145', '147', '149', '151', '153', '157', '159', '197', '257', '259', '261', '263', 'River', 'Lake'], ['35', '40', '161', '163', '164', '166', '167', '169', '171', '173', '179', '181', '183', '184', '185', '187', '189', '191', '193', '195', '199', '201', '203', '204', '205', '207', '265', '267', '269', '271', '273', '275', '1', '177'], ['15', '20', '50', '127', '129', '131', '139', '141', '143', '206', '208', '209', '211', '213', '215', '217', '219', '225', '229', '231', '237', '239', '241', '243', '247', '249', '251', '253', '255', '2', '3']] + assert dmas[0] == ['15', '20', '60', '601', '61', '117', '119', '120', '121', '123', '125', '127', '129', '131', '139', '141', '143', '145', '147', '149', '151', '153', '157', '159', '161', '195', '257', '259', '261', '263', 'River', '3'] + assert dmas[1] == ['50', '171', '173', '184', '199', '201', '203', '205', '206', '207', '208', '209', '211', '213', '215', '217', '219', '225', '229', '231', '237', '239', '241', '243', '247', '249', '251', '253', '255', '273', '275', '2'] + assert dmas[2] == ['10', '35', '40', '101', '103', '105', '107', '109', '111', '113', '115', '163', '164', '166', '167', '169', '177', '179', '181', '183', '185', '187', '189', '191', '193', '197', '204', '265', '267', '269', '271', 'Lake', '1'] self.leave(p) @@ -5864,9 +5864,9 @@ class TestApi: dmas = calculate_district_metering_area_for_network(p, 3) assert len(dmas) == 3 - assert dmas[0] == ['173', '184', '185', '199', '2', '201', '203', '205', '206', '207', '208', '209', '211', '213', '215', '217', '219', '225', '229', '231', '237', '239', '241', '243', '247', '249', '251', '253', '255', '273', '275', '50'] - assert dmas[1] == ['1', '10', '101', '103', '109', '111', '113', '161', '163', '164', '166', '167', '169', '171', '179', '181', '183', '187', '189', '191', '193', '195', '197', '204', '265', '267', '269', '271', '35', '40', 'Lake'] - assert dmas[2] == ['105', '107', '115', '117', '119', '120', '121', '123', '125', '127', '129', '131', '139', '141', '143', '145', '147', '149', '15', '151', '153', '157', '159', '20', '257', '259', '261', '263', '3', '60', '601', '61', 'River'] + assert dmas[0] == ['117', '119', '120', '121', '123', '125', '127', '129', '131', '139', '141', '143', '145', '147', '149', '15', '151', '153', '157', '159', '161', '195', '20', '257', '259', '261', '263', '3', '60', '601', '61', 'River'] + assert dmas[1] == ['1', '163', '164', '166', '167', '169', '171', '173', '177', '179', '181', '183', '184', '185', '187', '189', '191', '193', '199', '201', '203', '204', '205', '207', '265', '267', '269', '271', '273', '275', '35', '40'] + assert dmas[2] == ['10', '101', '103', '105', '107', '109', '111', '113', '115', '197', '2', '206', '208', '209', '211', '213', '215', '217', '219', '225', '229', '231', '237', '239', '241', '243', '247', '249', '251', '253', '255', '50', 'Lake'] self.leave(p) @@ -6150,13 +6150,16 @@ class TestApi: assert cs[1]['id'] == 'DMA_[DMA_1_1]_2_2' cs = generate_sub_district_metering_area(p, 'DMA_1_2', 3).operations - assert len(cs) == 2 + assert len(cs) == 3 assert cs[0]['operation'] == API_ADD assert cs[0]['type'] == 'district_metering_area' - assert cs[0]['id'] == 'DMA_[DMA_1_2]_2_2' + assert cs[0]['id'] == 'DMA_[DMA_1_2]_2_1' assert cs[1]['operation'] == API_ADD assert cs[1]['type'] == 'district_metering_area' - assert cs[1]['id'] == 'DMA_[DMA_1_2]_2_3' + assert cs[1]['id'] == 'DMA_[DMA_1_2]_2_2' + assert cs[2]['operation'] == API_ADD + assert cs[2]['type'] == 'district_metering_area' + assert cs[2]['id'] == 'DMA_[DMA_1_2]_2_3' cs = generate_sub_district_metering_area(p, 'DMA_1_3', 2).operations assert len(cs) == 2 @@ -6168,10 +6171,10 @@ class TestApi: assert cs[1]['id'] == 'DMA_[DMA_1_3]_2_2' dmas = get_all_district_metering_area_ids(p) - assert len(dmas) == 9 + assert len(dmas) == 10 cs = generate_district_metering_area(p, 3).operations - assert len(cs) == 12 + assert len(cs) == 13 dmas = get_all_district_metering_area_ids(p) assert len(dmas) == 3 @@ -6585,8 +6588,7 @@ class TestApi: open_project(p) result = calculate_demand_to_network(p, 100.0) - assert result == {'10': 3.2923444181778216, '101': 4.194261304565972, '103': 1.226514223391597, '105': 1.561545046227298, '107': 0.7929449232512782, '109': 1.3772201298574833, '111': 1.7690554866687873, '113': 1.2381069854274345, '115': 1.6902247048250931, '117': 1.088560355165132, '119': 1.9174428407275061, '120': 1.2276734995951808, '121': 1.502421959844527, '123': 10.897196313687157, '125': 2.0240962514572103, '127': 1.141887060529984, '129': 2.3486935884606575, '131': 1.5024219598445272, '139': 1.1129051554403906, '141': 1.6137124753885663, '143': 0.7071584841860814, '145': 1.3238934244926313, '147': 0.7141141414075839, '149': 0.44052495736182123, '15': 0.38256114718263423, '151': 1.3099821100496263, '153': 1.328530529306966, '157': 1.1569576511765725, '159': 1.1384092319192327, '161': 1.0966752885902182, '163': 0.41270232847581145, '164': 0.14838735405871872, '166': 0.11360906795120652, '167': 0.01391131444300488, '169': 0.5949405476791755, '171': 0.48225890069083593, '173': 0.9390137249028294, '179': 0.32459733700344723, '181': 0.06723801980785693, '183': 0.4799403482836684, '184': 1.0734665789944717, '185': 0.4486167052628357, '187': 0.732639375140852, '189': 0.8638926269106031, '191': 0.9807476682318441, '193': 0.8601829430591352, '195': 0.6167349403065497, '197': 0.9135096484239872, '199': 1.2415848140381855, '20': 0.20496003279360525, '201': 0.2434480027525854, '203': 0.02782262888600976, '204': 0.33037053249729426, '205': 1.4780771595692686, '206': 0.22258103108807809, '207': 0.7141141414075839, '208': 0.3234380607998635, '209': 0.48573672930158707, '211': 0.9923404302676815, '213': 1.7331179243576913, '215': 1.378379406061067, '217': 1.2218771185772621, '219': 0.47530324346933345, '225': 0.3616941755181269, '229': 1.1476834415479027, '231': 0.4544362718048261, '237': 0.7836707136226083, '239': 0.22605885969882933, '241': 0.6213720451208846, '243': 0.5100815295768456, '247': 0.42777291912240006, '249': 0.4382064049546538, '251': 0.5912308638277075, '253': 0.2550407647884228, '255': 1.0468264118361172, '257': 0.5552933015166115, '259': 0.405746671254309, '261': 0.5552933015166115, '263': 1.0375522022074473, '265': 0.7813521612154408, '267': 0.9227838580526571, '269': 0.5986502315306435, '271': 0.535585606055688, '273': 0.8346788665802929, '275': 0.9181467532383222, '35': 0.00695565722150244, '40': 0.2988614052838882, '50': 0.23741976649394997, '60': 0.28564565656303353, '601': 0.00046371048143349603, '61': 10.549645307852751} - + assert result == {'10': 3.2914286561977604, '101': 4.1930946753955975, '103': 1.226173069808884, '105': 1.5611107041895715, '107': 0.7927243664927, '109': 1.3768370575925843, '111': 1.7685634258302052, '113': 1.237762607330707, '115': 1.689754570681808, '117': 1.0882575732991893, '119': 1.9169095061095407, '120': 1.2273320235610663, '121': 1.5020040628282736, '123': 10.894165270513714, '125': 2.0235332513103135, '127': 1.1415694458995753, '129': 2.34804030192136, '131': 1.5020040628282738, '139': 1.1125956020950176, '141': 1.6132636230377755, '143': 0.7069617888312091, '145': 1.323525184992198, '147': 0.713915511344303, '149': 0.4404024258292778, '15': 0.3824547382201623, '151': 1.3096177399660103, '153': 1.3281610000009272, '157': 1.1566358446779454, '159': 1.1380925846430283, '161': 1.0963702495644652, '163': 0.4125875357769024, '164': 0.14834608027933568, '166': 0.11357746771386638, '167': 0.01390744502618772, '169': 0.5947750656199615, '171': 0.482124760907841, '173': 0.9387525392676711, '177': 0.01390744502618772, '179': 0.3314607731241407, '181': 0.07417304013966784, '183': 0.4798068534034764, '184': 1.0731679954457753, '185': 0.44849192301951035, '187': 0.7324355923041763, '189': 0.8636523361262574, '191': 0.9804748743462343, '193': 0.859943684119274, '195': 0.6165633961609889, '197': 0.9132555567196603, '199': 1.241239468587254, '20': 0.2049030233858324, '201': 0.24338028795828512, '203': 0.02781489005237544, '204': 0.3302786402969147, '205': 1.4776660340324455, '206': 0.2225191204190035, '207': 0.713915511344303, '208': 0.32334809685886445, '209': 0.48560162216438785, '211': 0.9920644118680574, '213': 1.7326358595125533, '215': 1.3779960113447665, '217': 1.2215372548001548, '219': 0.4751710383947471, '225': 0.3615935706808807, '229': 1.1473642146604868, '231': 0.45430987085546554, '237': 0.7834527364752415, '239': 0.22599598167555046, '241': 0.6211992111697181, '243': 0.5099396509602164, '247': 0.4276539345552724, '249': 0.4380845183249132, '251': 0.5910664136129782, '253': 0.2549698254801082, '255': 1.0465352382206259, '257': 0.5551388472953265, '259': 0.4056338132638085, '261': 0.5551388472953265, '263': 1.0372636082031674, '265': 0.7811348289708769, '267': 0.9225271867371188, '269': 0.5984837176269449, '271': 0.5354366335082272, '273': 0.8344467015712631, '275': 0.9178913717283894, '35': 0.00695372251309386, '40': 0.29877827731259954, '50': 0.23735372844693708, '60': 0.28556620453772114, '601': 0.000463581500872924, '61': 10.546710935609457} self.leave(p)