support incremental AOI updates for map building

.gitignore

@@ -1,6 +1,8 @@
 ignore/
 cache/
 *.html
+*.json
+*.pkl
 *.ipynb
 # Byte-compiled / optimized / DLL files
 __pycache__/

mosstool/map/_map_util/add_aoi_pop.py

@@ -7,21 +7,12 @@
 import rasterio
 from shapely.geometry import MultiPolygon, Point, Polygon
 
+from .._map_util.aoiutils import geo_coords
 from .const import *
 
 __all__ = ["add_aoi_pop"]
 
 
-def _coords(geo):
-    if isinstance(geo, Polygon):
-        return list(geo.exterior.coords)
-    elif isinstance(geo, MultiPolygon):
-        all_coords = []
-        for p_geo in geo.geoms:
-            all_coords.extend(list(p_geo.exterior.coords))
-        return all_coords
-    else:
-        return list(geo.coords)
 
 
 def _gps_distance(
@@ -211,7 +202,7 @@ def _get_aoi_poly_pop_unit(aoi):
         return aoi, HAS_INSIDE_PIXEL
 
     # If aoi is too small so that no pixel falls within it, take the population of the pixel where it is located.
-    x, y = _coords(poly.centroid)[0][:2]
+    x, y = geo_coords(poly.centroid)[0][:2]
     try:
         t = (
             pixel_idx2point_pop[

mosstool/map/_map_util/aoi_matcher.py

@@ -19,6 +19,7 @@
 from .._util.angle import abs_delta_angle, delta_angle
 from .._util.line import (connect_line_string, get_line_angle,
                           get_start_vector, line_extend, offset_lane)
+from .aoiutils import geo_coords
 
 # ATTENTION: In order to achieve longer distance POI merging, the maximum recursion depth needs to be modified.
 sys.setrecursionlimit(50000)
@@ -47,7 +48,7 @@ def _split_merged_poi_unit(points):
         else:
             # Post-processing uses road network segmentation to form polygons
             lines = []  # dividing lines
-            (x, y) = convex.centroid.coords[:][0][:2]
+            (x, y) = geo_coords(convex.centroid)[0][:2]
             length = convex.length
             # {'id', 'geo', 'point', 'length'}
             for lane in road_lane_matcher:
@@ -323,7 +324,7 @@ def _process_matched_result(aoi, d_matched, w_matched):
 #     huge_candidate = []
 #     (x, y) = geo.centroid.coords[:][0][:2]
 #     length = geo.length
-#     bound_poss = [c[:2] for c in geo.exterior.coords[:-1]]
+#     bound_poss = [c[:2] for c in geo_coords(geo)[:-1]]
 #     for lane in matcher:  # {'id', 'geo', 'point', 'length'}
 #         mid_x, mid_y = lane["point"][:2]
 #         dis_upper_bound = (
@@ -464,19 +465,29 @@ def _str_tree_matcher_unit(
     matcher_lane_tree,
     dis_gate,
     huge_gate,
+    direction_geos: Optional[Dict[int, List[LineString]]] = None,
 ):
     global LENGTH_PER_DOOR, MAX_DOOR_NUM, AOI_GATE_OFFSET
     matched = []
-    bound_poss = [c[:2] for c in geo.exterior.coords[:-1]]
+    bound_poss = [c[:2] for c in geo_coords(geo)[:-1]]
     small_tree_ids = matcher_lane_tree.query(geo.buffer(dis_gate))
+    stop_lane_angles = (
+        [[get_line_angle(l) for l in lanes] for lanes in direction_geos.values()]
+        if direction_geos is not None
+        else None
+    )
     for tid in small_tree_ids:
         lane = matcher[tid]  # {'id', 'geo', 'point', 'length'}
         p_aoi, p_lane = ops.nearest_points(
             geo, lane["geo"]
         )  # Returns the calculated nearest points in the input geometries
         distance = p_aoi.distance(p_lane)
-        # if distance < dis_gate:
-        if True:
+        lane_angle = get_line_angle(lane["geo"])
+        if stop_lane_angles is None or any(
+            np.mean([abs_delta_angle(lane_angle, angle) for angle in angles])
+            < np.pi / 4
+            for angles in stop_lane_angles
+        ):
             # Project the point on the lane closest to the poly to the lane and return s
             s = lane["geo"].project(p_lane)
             if (
@@ -596,7 +607,7 @@ def _add_point_aoi_unit(arg):
     global W_DIS_GATE, W_HUGE_GATE
     aoi, aoi_type = arg
     geo = aoi["geo"]
-    x, y = geo.coords[:][0][:2]
+    x, y = geo_coords(geo)[0][:2]
     d_matched, w_matched = [], []
     for (
         matcher,
@@ -717,7 +728,7 @@ def _add_poly_aoi_unit(arg):
         base_aoi = {
             "id": 0,  # It is difficult to deal with the problem of uid += 1 during parallelization, so assign the id after parallelization is completed.
             "type": aoi_type,
-            "positions": [{"x": c[0], "y": c[1]} for c in geo.exterior.coords[:]],
+            "positions": [{"x": c[0], "y": c[1]} for c in geo_coords(geo)],
             "area": geo.area,
             "external": {
                 "osm_tencent_ids": [
@@ -749,7 +760,7 @@ def _add_aoi_stop_unit(arg):
     global LENGTH_PER_DOOR, MAX_DOOR_NUM, AOI_GATE_OFFSET
     aoi, aoi_type = arg
     geo = aoi["geo"]
-    bound_poss = [c[:2] for c in geo.exterior.coords[:-1]]
+    bound_poss = [c[:2] for c in geo_coords(geo)[:-1]]
     station_type = aoi["external"]["station_type"]
     d_matched, w_matched = [], []
     if station_type == "SUBWAY":
@@ -823,7 +834,7 @@ def _add_aoi_stop_unit(arg):
         base_aoi = {
             "id": 0,  # It is difficult to deal with the problem of uid += 1 during parallelization, so assign the id after parallelization is completed.
             "type": aoi_type,
-            "positions": [{"x": c[0], "y": c[1]} for c in geo.exterior.coords[:]],
+            "positions": [{"x": c[0], "y": c[1]} for c in geo_coords(geo)],
             "area": geo.area,
             "external": {
                 "osm_tencent_ids": [
@@ -925,10 +936,10 @@ def map_land_use(euluc: int) -> int:
                     polygon = candidate_poly
                 else:
                     polygon = MultiPoint(
-                        [pt for g in geo.geoms for pt in g.exterior.coords]
+                        [pt for g in geo.geoms for pt in geo_coords(g)]
                     ).convex_hull
             else:
-                polygon = MultiPoint([pt for pt in polygon.exterior.coords]).convex_hull
+                polygon = MultiPoint([pt for pt in geo_coords(polygon)]).convex_hull
             df.geometry[i] = polygon
 
     # Processing AOIs
@@ -959,7 +970,7 @@ def map_land_use(euluc: int) -> int:
                     polygon = candidate_poly
                 else:
                     polygon = MultiPoint(
-                        [pt for g in geo.geoms for pt in g.exterior.coords]
+                        [pt for g in geo.geoms for pt in geo_coords(g)]
                     ).convex_hull
             else:
                 aoi["land_use"] = map_land_use(-1)  # 默认
@@ -1630,7 +1641,7 @@ def _merge_covered_aoi(aois, workers):
     # Pre-compute geometric properties
     for aoi in aois:
         geo = aoi["geo"]
-        aoi["point"] = list(geo.centroid.coords)[0]  # Geometric center
+        aoi["point"] = geo_coords(geo.centroid)[0]  # Geometric center
         aoi["length"] = geo.length  # Perimeter
         aoi["area"] = geo.area  # area
         aoi["valid"] = geo.is_valid

mosstool/map/_map_util/aoiutils.py

@@ -11,14 +11,27 @@
 __all__ = [
     "generate_aoi_poi",
     "generate_sumo_aoi_poi",
+    "geo_coords",
 ]
 
 
-def _fix_polygon(polygon: Polygon):
-    if polygon.is_valid:
-        return polygon
+def geo_coords(geo):
+    if isinstance(geo, Polygon):
+        return list(geo.exterior.coords)
+    elif isinstance(geo, MultiPolygon):
+        all_coords = []
+        for p_geo in geo.geoms:
+            all_coords.extend(geo_coords(p_geo))
+        return all_coords
+    else:
+        return list(geo.coords)
+
+
+def _fix_polygon(input_poly: Polygon):
+    if input_poly.is_valid:
+        return input_poly
     else:
-        geo = polygon.buffer(0)
+        geo = input_poly.buffer(0)
         if isinstance(geo, Polygon) and geo.is_valid:
             polygon = geo
         elif isinstance(geo, MultiPolygon):
@@ -31,10 +44,10 @@ def _fix_polygon(polygon: Polygon):
                 polygon = candidate_poly
             else:
                 polygon = MultiPoint(
-                    [pt for g in geo.geoms for pt in g.exterior.coords]
+                    [pt for g in geo.geoms for pt in geo_coords(g)]
                 ).convex_hull
         else:
-            polygon = MultiPoint([pt for pt in polygon.exterior.coords]).convex_hull
+            polygon = MultiPoint([pt for pt in geo_coords(input_poly)]).convex_hull
         return polygon
 
 
@@ -44,7 +57,7 @@ def _fix_aois_poly(input_aois: list) -> list:
         coords = aoi["coords"]
         geo = _fix_polygon(Polygon(coords))
         if geo.is_valid and geo:
-            aoi["coords"] = [c for c in geo.exterior.coords]
+            aoi["coords"] = [c for c in geo_coords(geo)]
             aois.append(aoi)
     return aois
 
@@ -54,7 +67,7 @@ def _connect_aoi_unit1(poly):
     Find the children aoi of each poly in the merged geometry
     """
     global aois_small
-    x, y = list(poly.centroid.coords)[0]
+    x, y = geo_coords(poly.centroid)[0]
     length = poly.length
     children = []
     added_aoi_small = []
@@ -81,7 +94,7 @@ def _connect_aoi_unit2(arg):
         ]  # There is only 1 child, indicating that it is not connected to other aoi and restores the original shape
 
     poly_inner = MultiPoint(
-        [pt for aoi in children for pt in aoi["geo"].exterior.coords]
+        [pt for aoi in children for pt in geo_coords(aoi["geo"])]
     ).convex_hull
     poly = poly.intersection(poly_inner)
     if not isinstance(poly, Polygon):
@@ -113,7 +126,7 @@ def _connect_aoi_unit2(arg):
                 "inner_poi": inner_poi,
                 "inner_poi_catg": inner_poi_catg,
             },
-            "point": list(poly.centroid.coords)[
+            "point": geo_coords(poly.centroid)[
                 0
             ],  # For subsequent processing needs, calculate the merged geometric center
             "length": poly.length,  # perimeter
@@ -171,7 +184,7 @@ def _merge_aoi(input_aois: list, merge_aoi: bool = False, workers: int = 32):
             logging.warning(f"Invalid polygon {aoi['id']}")
             continue
         aoi["geo"] = geo
-        aoi["point"] = list(geo.centroid.coords)[0]  # Geometric center
+        aoi["point"] = geo_coords(geo.centroid)[0]  # Geometric center
         aoi["length"] = geo.length  # Perimeter
         aoi["area"] = geo.area  # area
         aois.append(aoi)
@@ -221,7 +234,7 @@ def _connect_aoi(input_aois: list, merge_aoi: bool = False, workers: int = 32):
             logging.warning(f"Invalid polygon {aoi['id']}")
             continue
         aoi["geo"] = geo
-        aoi["point"] = list(geo.centroid.coords)[0]  # Geometric center
+        aoi["point"] = geo_coords(geo.centroid)[0]  # Geometric center
         aoi["length"] = geo.length  # Perimeter
         aoi["area"] = geo.area  # area
         aois.append(aoi)
@@ -239,7 +252,7 @@ def _connect_aoi(input_aois: list, merge_aoi: bool = False, workers: int = 32):
     polys = [aoi["geo"] for aoi in aois_small]
     polys_scale = [scale(p, xfact=SCALE, yfact=SCALE, origin="centroid") for p in polys]
     geo_scale_connect = ops.unary_union(polys_scale)
-    args = list(geo_scale_connect.geoms)
+    args = list(geo_scale_connect.geoms)  # type: ignore
     with Pool(processes=workers) as pool:
         results = pool.map(
             _connect_aoi_unit1,
@@ -326,7 +339,7 @@ def point_extend(center: Point, length: float):
             geo = Polygon(coords)
         geo = _fix_polygon(geo)
         stop["geo"] = geo
-        stop["point"] = list(geo.centroid.coords)[0]  # Geometric center
+        stop["point"] = geo_coords(geo.centroid)[0]  # Geometric center
         stop["length"] = geo.length  # Perimeter
         stop["area"] = geo.area  # area
     return stops
@@ -381,7 +394,7 @@ def _match_poi_to_aoi(aois, pois, workers):
 def _post_compute_aoi_poi(aois, pois_isolate):
     # Update coordinates
     for a in aois:
-        coords = a["geo"].exterior.coords[:]
+        coords = geo_coords(a["geo"])
         a["coords"] = coords
     # poi becomes Aoi independently
     for p in pois_isolate:

mosstool/map/_map_util/map_aois_matchers.py

@@ -10,7 +10,7 @@
 
 from ...type import Map
 from ...util.format_converter import pb2dict
-from .aoi_matcher import _process_matched_result,_str_tree_matcher_unit
+from .aoi_matcher import _process_matched_result, _str_tree_matcher_unit
 from .const import *
 
 
@@ -23,12 +23,12 @@ def _add_aoi_unit(aoi):
     global d_matcher, w_matcher
     global D_DIS_GATE, D_HUGE_GATE
     global W_DIS_GATE, W_HUGE_GATE
-    global d_tree,w_tree
+    global d_tree, w_tree
     geo = _map_aoi2geo(aoi)
     # d_matched = _matcher_unit(geo, d_matcher, D_DIS_GATE, D_HUGE_GATE)
     # w_matched = _matcher_unit(geo, w_matcher, W_DIS_GATE, W_HUGE_GATE)
-    d_matched = _str_tree_matcher_unit(geo, d_matcher, d_tree,D_DIS_GATE, D_HUGE_GATE)
-    w_matched = _str_tree_matcher_unit(geo, w_matcher, w_tree,W_DIS_GATE, W_HUGE_GATE)
+    d_matched = _str_tree_matcher_unit(geo, d_matcher, d_tree, D_DIS_GATE, D_HUGE_GATE)
+    w_matched = _str_tree_matcher_unit(geo, w_matcher, w_tree, W_DIS_GATE, W_HUGE_GATE)
     aoi["external"] = {}
     if d_matched or w_matched:
         return _process_matched_result(
@@ -38,7 +38,7 @@ def _add_aoi_unit(aoi):
 
 def match_map_aois(net: Map, matchers: dict, workers: int):
     global d_matcher, w_matcher
-    global d_tree,w_tree
+    global d_tree, w_tree
     net_dict = pb2dict(net)
     orig_aois = net_dict["aois"]
     orig_pois = net_dict["pois"]
@@ -48,27 +48,38 @@ def match_map_aois(net: Map, matchers: dict, workers: int):
     )
     d_tree = STRtree([l["geo"] for l in d_matcher])
     w_tree = STRtree([l["geo"] for l in w_matcher])
-    results_aoi = []
+    results_aois = []
     for i in range(0, len(orig_aois), MAX_BATCH_SIZE):
         args_batch = orig_aois[i : i + MAX_BATCH_SIZE]
         with Pool(processes=workers) as pool:
-            results_aoi += pool.map(
+            results_aois += pool.map(
                 _add_aoi_unit,
                 args_batch,
                 chunksize=min(ceil(len(args_batch) / workers), 500),
             )
-    aois = [r for r in results_aoi if r]
+    aois = [r for r in results_aois if r]
     # filter pois
     all_poi_ids = set(pid for a in aois for pid in a["poi_ids"])
     pois = [p for p in orig_pois if p["id"] in all_poi_ids]
     # rearrange aoi id
     aois_dict = {}
+    _poi_uid_dict = {
+        poi["id"]: poi_id for poi_id, poi in enumerate(pois, start=POI_START_ID)
+    }
     for aoi_id, aoi in enumerate(aois, start=AOI_START_ID):
         aoi["id"] = aoi_id
         aois_dict[aoi_id] = aoi
+        # update poi_ids
+        aoi["poi_ids"] = [_poi_uid_dict[pid] for pid in aoi["poi_ids"]]
     # rearrange poi id
     pois_dict = {}
+    _valid_aois = [(a["id"], aoi) for aoi, a in zip(results_aois, orig_aois) if aoi]
+    _aoi_uid_dict = {
+        aid: AOI_START_ID + idx for idx, (aid, _) in enumerate(_valid_aois)
+    }
     for poi_id, poi in enumerate(pois, start=POI_START_ID):
         poi["id"] = poi_id
         pois_dict[poi_id] = poi
+        # update aoi_id
+        poi["aoi_id"] = _aoi_uid_dict[poi["aoi_id"]]
     return (aois_dict, pois_dict)