加入吉林一号影像

tools/viusal_tifs_gps.py

@@ -0,0 +1,39 @@
+import rasterio
+from rasterio.warp import transform
+import matplotlib.pyplot as plt
+import glob
+import os
+
+def get_bounds_lonlat(dataset):
+    # 直接用元数据获取边界
+    left, bottom, right, top = dataset.bounds
+    # 四角顺序：左下、左上、右上、右下
+    xs = [left, left, right, right]
+    ys = [bottom, top, top, bottom]
+    # 转为经纬度
+    lons, lats = transform(dataset.crs, 'EPSG:4326', xs, ys)
+    return list(zip(lons, lats))
+
+if __name__ == "__main__":
+    tif_dir = "test_geosat1_imgs"
+    pattern = os.path.join(tif_dir, "Jilin*.tif")
+    all_bounds = []
+    for tif in glob.glob(pattern):
+        with rasterio.open(tif) as ds:
+            corners_lonlat = get_bounds_lonlat(ds)
+            all_bounds.append((os.path.basename(tif), corners_lonlat))
+
+    plt.figure(figsize=(10, 10))
+    for name, corners in all_bounds:
+        # 闭合多边形
+        corners_closed = corners + [corners[0]]
+        lons, lats = zip(*corners_closed)
+        plt.plot(lons, lats, marker='o', label=name)
+        for lon, lat in corners:
+            plt.text(lon, lat, f'({lon:.4f},{lat:.4f})', color='red', fontsize=8)
+    plt.xlabel("Longitude")
+    plt.ylabel("Latitude")
+    plt.title("Jilin* 影像边界经纬度分布")
+    plt.legend()
+    plt.grid(True)
+    plt.show()