From 4092854b0bea9e1fc02f29222c6f7cd876e565a5 Mon Sep 17 00:00:00 2001
From: shuishen <1109946754@qq.com>
Date: Mon, 17 Aug 2026 09:43:20 +0800
Subject: [PATCH] feat:项目基础整理

---
 capabilities/15-trajectory-analysis/run_trajectory_analysis.py |   57 ++++++++++++++++++++++++++++++++++++++++++++++++++++-----
 1 files changed, 52 insertions(+), 5 deletions(-)

diff --git a/capabilities/15-trajectory-analysis/run_trajectory_analysis.py b/capabilities/15-trajectory-analysis/run_trajectory_analysis.py
index 947c0f8..1afb80d 100644
--- a/capabilities/15-trajectory-analysis/run_trajectory_analysis.py
+++ b/capabilities/15-trajectory-analysis/run_trajectory_analysis.py
@@ -25,6 +25,7 @@
 matplotlib.use("Agg")
 import matplotlib.pyplot as plt
 from matplotlib.lines import Line2D
+from matplotlib.patches import Patch
 
 
 DEFAULT_THRESHOLDS = {
@@ -105,7 +106,7 @@
     return CRS.from_epsg(epsg)
 
 
-def load_case(manifest_path: Path) -> tuple[dict[str, Any], pd.DataFrame, gpd.GeoDataFrame, gpd.GeoDataFrame, int, CRS]:
+def load_case(manifest_path: Path) -> tuple[dict[str, Any], pd.DataFrame, gpd.GeoDataFrame, gpd.GeoDataFrame, gpd.GeoDataFrame | None, int, CRS]:
     manifest = load_json(manifest_path)
     for field in ("case_id", "crs", "observations", "reference_routes", "zones"):
         if not manifest.get(field):
@@ -142,6 +143,12 @@
         raise ValueError("Every track must contain at least two distinct timestamps.")
     routes = gpd.read_file(routes_path)
     zones = gpd.read_file(zones_path)
+    flyable_zones: gpd.GeoDataFrame | None = None
+    if manifest.get("flyable_zones"):
+        flyable_path = resolve_input_path(manifest_path, str(manifest["flyable_zones"]))
+        flyable_zones = gpd.read_file(flyable_path)
+        if flyable_zones.crs is None:
+            raise ValueError("Flyable-zone GeoJSON must declare a CRS.")
     if "track_id" not in routes.columns:
         raise ValueError("Reference route GeoJSON must contain a track_id property.")
     if "zone_id" not in zones.columns or "zone_type" not in zones.columns:
@@ -153,7 +160,15 @@
     if missing_routes:
         raise ValueError(f"Missing reference routes for: {', '.join(sorted(missing_routes))}")
     metric_crs = choose_metric_crs(float(frame["longitude"].mean()), float(frame["latitude"].mean()))
-    return manifest, frame, routes.to_crs(metric_crs), zones.to_crs(metric_crs), duplicate_count, metric_crs
+    return (
+        manifest,
+        frame,
+        routes.to_crs(metric_crs),
+        zones.to_crs(metric_crs),
+        flyable_zones.to_crs(metric_crs) if flyable_zones is not None else None,
+        duplicate_count,
+        metric_crs,
+    )
 
 
 def true_runs(flags: list[bool], times: list[pd.Timestamp], max_gap: float) -> Iterable[tuple[int, int]]:
@@ -410,6 +425,7 @@
     tracks: dict[str, gpd.GeoDataFrame],
     routes_wgs84: gpd.GeoDataFrame,
     zones_wgs84: gpd.GeoDataFrame,
+    flyable_zones_wgs84: gpd.GeoDataFrame | None,
     metadata: dict[str, Any],
 ) -> None:
     output_dir.mkdir(parents=True, exist_ok=True)
@@ -435,7 +451,22 @@
         for event in events
     ]
     write_geojson(output_dir / "events.geojson", event_features)
-    render_map(output_dir / "analysis.png", tracks, routes_wgs84, zones_wgs84, events, manifest["case_id"])
+    route_features = json.loads(routes_wgs84.to_json())["features"]
+    write_geojson(output_dir / "reference_routes.geojson", route_features)
+    zone_features = json.loads(zones_wgs84.to_json())["features"]
+    write_geojson(output_dir / "zones.geojson", zone_features)
+    if flyable_zones_wgs84 is not None:
+        flyable_features = json.loads(flyable_zones_wgs84.to_json())["features"]
+        write_geojson(output_dir / "flyable_zones.geojson", flyable_features)
+    render_map(
+        output_dir / "analysis.png",
+        tracks,
+        routes_wgs84,
+        zones_wgs84,
+        flyable_zones_wgs84,
+        events,
+        manifest["case_id"],
+    )
     (output_dir / "run_metadata.json").write_text(
         json.dumps(metadata, ensure_ascii=False, indent=2), encoding="utf-8"
     )
@@ -457,17 +488,20 @@
     tracks: dict[str, gpd.GeoDataFrame],
     routes: gpd.GeoDataFrame,
     zones: gpd.GeoDataFrame,
+    flyable_zones: gpd.GeoDataFrame | None,
     events: list[dict[str, Any]],
     case_id: str,
 ) -> None:
     figure, axis = plt.subplots(figsize=(10, 7), dpi=140)
     zones.plot(ax=axis, facecolor="#f0e442", edgecolor="#8c6d1f", alpha=0.28, linewidth=1.5)
-    routes.plot(ax=axis, color="#666666", linestyle="--", linewidth=1.2)
+    if flyable_zones is not None and not flyable_zones.empty:
+        flyable_zones.plot(ax=axis, facecolor="#36a269", edgecolor="#167344", alpha=0.22, linewidth=1.3)
     palette = ["#0072b2", "#009e73", "#d55e00", "#cc79a7", "#56b4e9"]
     for index, (track_id, track) in enumerate(sorted(tracks.items())):
         color = palette[index % len(palette)]
         axis.plot(track.geometry.x, track.geometry.y, color=color, linewidth=2.2, marker="o", markersize=2.8, label=track_id)
         axis.annotate(track_id, (track.geometry.x.iloc[-1], track.geometry.y.iloc[-1]), fontsize=7, color=color)
+    routes.plot(ax=axis, color="#444444", linestyle="--", linewidth=1.5, zorder=4, label="reference route")
     for event in events:
         event_type = event["event_type"]
         axis.scatter(
@@ -481,6 +515,11 @@
             zorder=5,
         )
     handles, labels = axis.get_legend_handles_labels()
+    handles.append(Patch(facecolor="#f0e442", edgecolor="#8c6d1f", alpha=0.45, label="restricted area"))
+    labels.append("restricted area")
+    if flyable_zones is not None and not flyable_zones.empty:
+        handles.append(Patch(facecolor="#36a269", edgecolor="#167344", alpha=0.45, label="flyable area"))
+        labels.append("flyable area")
     present_event_types = {event["event_type"] for event in events}
     handles.extend(
         Line2D(
@@ -503,6 +542,13 @@
     axis.set_ylabel("Latitude (WGS84)")
     axis.grid(alpha=0.18)
     axis.ticklabel_format(useOffset=False)
+    visible_geometries = [geometry for track in tracks.values() for geometry in track.geometry]
+    visible_geometries.extend(geometry for geometry in routes.geometry if geometry is not None)
+    bounds = gpd.GeoSeries(visible_geometries, crs="EPSG:4326").total_bounds
+    width = max(bounds[2] - bounds[0], 0.0005)
+    height = max(bounds[3] - bounds[1], 0.0005)
+    axis.set_xlim(bounds[0] - width * 0.08, bounds[2] + width * 0.08)
+    axis.set_ylim(bounds[1] - height * 0.08, bounds[3] + height * 0.08)
     figure.tight_layout()
     figure.savefig(path)
     plt.close(figure)
@@ -510,7 +556,7 @@
 
 def process_manifest(manifest_path: Path, output_dir: Path) -> dict[str, Any]:
     started = time.perf_counter()
-    manifest, frame, routes_metric, zones_metric, duplicates, metric_crs = load_case(manifest_path)
+    manifest, frame, routes_metric, zones_metric, flyable_zones_metric, duplicates, metric_crs = load_case(manifest_path)
     thresholds = dict(DEFAULT_THRESHOLDS)
     for key, value in manifest.get("thresholds", {}).items():
         if key not in thresholds:
@@ -566,6 +612,7 @@
         tracks,
         routes_metric.to_crs("EPSG:4326"),
         zones_metric.to_crs("EPSG:4326"),
+        flyable_zones_metric.to_crs("EPSG:4326") if flyable_zones_metric is not None else None,
         metadata,
     )
     return metadata

--
Gitblit v1.9.3