diff --git a/del_items.py b/tools/del_items.py similarity index 100% rename from del_items.py rename to tools/del_items.py diff --git a/ingest_stac_items.py b/tools/ingest_stac_items.py similarity index 100% rename from ingest_stac_items.py rename to tools/ingest_stac_items.py diff --git a/tools/viusal_tifs_gps.py b/tools/viusal_tifs_gps.py new file mode 100644 index 0000000..560ddcf --- /dev/null +++ b/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() \ No newline at end of file