使用多边形形状文件从大陆海岸线上删除岛屿 - GeoPandas

Question

我想从澳大利亚海岸线上删除所有断开的岛屿。我已经将该国作为一个Shapefile导入，可以使用硬编码边界删除区域。但是，在整个海岸线上追踪以去除所有断开的岛屿的效率是低效的。使用地理杂志有更快的方法吗？

import geoplot as gplt
import geopandas as gpd

Oz = gpd.read_file('OZ.shp')

#polygon = Polygon([(140, -40), (150, -40), (150, -45), (140,-45)])
#Oz_clip = gpd.clip(Oz, polygon)

Answer 1

• 使用过这里的几何图形：https://www.abs.gov.au/statistics/standards/australian-statistical-geography-standard-asgs-edition-3/jul2021-jun2026/access-and-downloads/digital-boundary-files
• 这包含一个由 6624 个多边形组成的多边形，
• 已将多边形扩展为多边形，然后创建一个包含区域的地理数据框 现在，您可以对每个多边形
• 的过滤

import geopandas as gpd

gdf = gpd.read_file(
    "https://www.abs.gov.au/statistics/standards/australian-statistical-geography-standard-asgs-edition-3/jul2021-jun2026/access-and-downloads/digital-boundary-files/AUS_2021_AUST_SHP_GDA2020.zip"
)

# expand polygons and move to UTM geometry
polys = (
    gpd.GeoSeries(
        gdf.dropna()["geometry"].apply(lambda g: g.geoms).explode(), crs=gdf.crs
    )
    .to_crs(gdf.estimate_utm_crs())
    .simplify(100)
    .reset_index(drop=True)
)

# now build geodataframe limited by biggest area polygons
gdfm = (
    gpd.GeoDataFrame(data={"area": polys.area / 10**6}, geometry=polys, crs=polys.crs)
    .sort_values("area", ascending=0)
    # .head(3) # only one tasmania will be excluded
    .loc[lambda d: d["area"].gt(5000)]  # filter by size of polygon
    .reset_index()
)

# visualise
gdfm.explore(height=300, width=500)

    index   area    geometry
0   6411    7.787225e+06    POLYGON ((1023611.393 5771135.050, 1021232.002...
1   6373    6.490288e+04    POLYGON ((1497821.632 5153803.332, 1498006.176...
2   554 5.826181e+03    POLYGON ((26191.908 8735865.816, 26640.856 873..

进行最大面积多边形head()或面积超过阈值的多边形