Python应用指南：获取行政区最小外接矩形

ops/2024/11/15 6:55:59/

有的时候我们用行政区的面层进行裁剪的时候，特别是在基于行政区裁剪的渔网，因为行政边界本身就是不规则的原因，在边界的渔网匹配数据会出现匹配不上的问题，其实如果周边有数据完全可以采用最近匹配的原则，但是就是因为行政区裁太可丁可卯了，会丢失一部分数据，这个时候可以通过形成一个最小外接矩形的形式，可以让数据裁剪保持一个微妙的范围，不至于太大，又不至于没有二次处理的空间，第二个方面就是比如我们需要行政区里面的POI，以高德为例，因为行政区本身不是规则的，我们通过建立外接矩阵，这样获取其中的数据就不会出现丢失的情况了；

行政区边界获取部分，这里采用了阿里云旗下的数据可视化平台，其数据源本身也是来着高德地图，数据质量还是比较可靠的；

全国行政区划边界GeoJSON数据下载平台阿里云：DataV.GeoAtlas地理小工具系列 (aliyun.com)

这里把下载到本地的geojson或者json的文件路径改成自己的即可；

完整代码#运行环境Python 3.11

python">import geopandas as gp# 结果变量
resultMinLon = 99999.0  # 最小经度
resultMaxLon = -99999.0  # 最大经度
resultMinLat = 99999.0  # 最小纬度
resultMaxLat = -99999.0  # 最大纬度# 读取厦门市边界数据
# 注意：请确保文件路径正确
loadGeoData = gp.read_file("D:\\data\\厦门市.json")# 遍历每一行数据
for i in range(len(loadGeoData)):# 读取几何数据loadGeometry = loadGeoData.loc[i, "geometry"]# 检查几何数据类型，可能是单个多边形（Polygon）或多个多边形（MultiPolygon）if loadGeometry.geom_type == 'Polygon':# 单个多边形# 遍历多边形的所有顶点for z in range(len(loadGeometry.exterior.coords)):# 获取当前顶点的经度和纬度loadLon, loadLat = loadGeometry.exterior.coords[z]# 更新最小经度if loadLon < resultMinLon:resultMinLon = loadLon# 更新最大经度if loadLon > resultMaxLon:resultMaxLon = loadLon# 更新最小纬度if loadLat < resultMinLat:resultMinLat = loadLat# 更新最大纬度if loadLat > resultMaxLat:resultMaxLat = loadLatelif loadGeometry.geom_type == 'MultiPolygon':# 多个多边形# 遍历每个多边形for polygon in loadGeometry.geoms:# 遍历多边形的所有顶点for z in range(len(polygon.exterior.coords)):# 获取当前顶点的经度和纬度loadLon, loadLat = polygon.exterior.coords[z]# 更新最小经度if loadLon < resultMinLon:resultMinLon = loadLon# 更新最大经度if loadLon > resultMaxLon:resultMaxLon = loadLon# 更新最小纬度if loadLat < resultMinLat:resultMinLat = loadLat# 更新最大纬度if loadLat > resultMaxLat:resultMaxLat = loadLat# 结果输出
print("厦门市最小经度：", resultMinLon)
print("厦门市最大经度：", resultMaxLon)
print("厦门市最小纬度：", resultMinLat)
print("厦门市最大纬度：", resultMaxLat)

这里会直接打印出来两个对角线坐标，另外，因为数据源是高德的，用的是（GCJ-02）坐标系，最好先转成WGS84地理坐标系或者CGCS2000地理坐标系。如果需要可视化的话,可以跑一下下面的代码，与上面的区别就是通过调用了matplotlib包，使得结果可视化；

完整代码#运行环境Python 3.11

python">import geopandas as gp
import matplotlib.pyplot as plt
from matplotlib.patches import Rectangle# 配置 matplotlib 支持中文显示
plt.rcParams['font.family'] = ['SimHei']  # 使用黑体字体
plt.rcParams['axes.unicode_minus'] = False  # 正常显示负号# 结果变量
resultMinLon = 99999.0  # 最小经度
resultMaxLon = -99999.0  # 最大经度
resultMinLat = 99999.0  # 最小纬度
resultMaxLat = -99999.0  # 最大纬度# 读取厦门市边界数据
# 注意：请确保文件路径正确
loadGeoData = gp.read_file("D:\\data\\厦门市.json")# 遍历每一行数据
for i in range(len(loadGeoData)):# 读取几何数据loadGeometry = loadGeoData.loc[i, "geometry"]# 检查几何数据类型，可能是单个多边形（Polygon）或多个多边形（MultiPolygon）if loadGeometry.geom_type == 'Polygon':# 单个多边形# 遍历多边形的所有顶点for z in range(len(loadGeometry.exterior.coords)):# 获取当前顶点的经度和纬度loadLon, loadLat = loadGeometry.exterior.coords[z]# 更新最小经度if loadLon < resultMinLon:resultMinLon = loadLon# 更新最大经度if loadLon > resultMaxLon:resultMaxLon = loadLon# 更新最小纬度if loadLat < resultMinLat:resultMinLat = loadLat# 更新最大纬度if loadLat > resultMaxLat:resultMaxLat = loadLatelif loadGeometry.geom_type == 'MultiPolygon':# 多个多边形# 遍历每个多边形for polygon in loadGeometry.geoms:# 遍历多边形的所有顶点for z in range(len(polygon.exterior.coords)):# 获取当前顶点的经度和纬度loadLon, loadLat = polygon.exterior.coords[z]# 更新最小经度if loadLon < resultMinLon:resultMinLon = loadLon# 更新最大经度if loadLon > resultMaxLon:resultMaxLon = loadLon# 更新最小纬度if loadLat < resultMinLat:resultMinLat = loadLat# 更新最大纬度if loadLat > resultMaxLat:resultMaxLat = loadLat# 结果输出
print("厦门最小经度：", resultMinLon)
print("厦门最大经度：", resultMaxLon)
print("厦门最小纬度：", resultMinLat)
print("厦门最大纬度：", resultMaxLat)# 可视化
fig, ax = plt.subplots(figsize=(10, 10))# 绘制厦门市的边界
loadGeoData.plot(ax=ax, edgecolor='black', facecolor='lightblue')# 绘制最小外接矩形
min_rect = Rectangle((resultMinLon, resultMinLat), resultMaxLon - resultMinLon, resultMaxLat - resultMinLat,linewidth=2, edgecolor='red', facecolor='none', linestyle='--')
ax.add_patch(min_rect)# 标注最小和最大经度与纬度点
ax.scatter(resultMinLon, resultMinLat, color='red', label='最小边界')
ax.scatter(resultMaxLon, resultMaxLat, color='green', label='最大边界')
ax.annotate(f'({resultMinLon:.6f}, {resultMinLat:.6f})', (resultMinLon, resultMinLat), textcoords="offset points",xytext=(0, 10), ha='center')
ax.annotate(f'({resultMaxLon:.6f}, {resultMaxLat:.6f})', (resultMaxLon, resultMaxLat), textcoords="offset points",xytext=(0, 10), ha='center')# 添加图例
ax.legend()# 设置坐标轴范围
ax.set_xlim(resultMinLon - 0.01, resultMaxLon + 0.01)
ax.set_ylim(resultMinLat - 0.01, resultMaxLat + 0.01)# 显示图形
plt.title('厦门市边界及最小外接矩形')
plt.xlabel('经度')
plt.ylabel('纬度')
plt.grid(True)
plt.show()

下面这部分代码增加了一个功能：就是把这个矩形导出来，这里是把这个矩形导出来形成shp文件，修改文件路径的基础上再改一下保存位置路径即可；

完整代码#运行环境Python 3.11

python">import geopandas as gp
import matplotlib.pyplot as plt
from matplotlib.patches import Rectangle
from shapely.geometry import box# 配置 matplotlib 支持中文显示
plt.rcParams['font.family'] = ['SimHei']  # 使用黑体字体
plt.rcParams['axes.unicode_minus'] = False  # 正常显示负号# 结果变量
resultMinLon = 99999.0  # 最小经度
resultMaxLon = -99999.0  # 最大经度
resultMinLat = 99999.0  # 最小纬度
resultMaxLat = -99999.0  # 最大纬度# 读取厦门市边界数据
# 注意：请确保文件路径正确
loadGeoData = gp.read_file("D:\\data\\厦门市.json")# 遍历每一行数据
for i in range(len(loadGeoData)):# 读取几何数据loadGeometry = loadGeoData.loc[i, "geometry"]# 检查几何数据类型，可能是单个多边形（Polygon）或多个多边形（MultiPolygon）if loadGeometry.geom_type == 'Polygon':# 单个多边形# 遍历多边形的所有顶点for z in range(len(loadGeometry.exterior.coords)):# 获取当前顶点的经度和纬度loadLon, loadLat = loadGeometry.exterior.coords[z]# 更新最小经度if loadLon < resultMinLon:resultMinLon = loadLon# 更新最大经度if loadLon > resultMaxLon:resultMaxLon = loadLon# 更新最小纬度if loadLat < resultMinLat:resultMinLat = loadLat# 更新最大纬度if loadLat > resultMaxLat:resultMaxLat = loadLatelif loadGeometry.geom_type == 'MultiPolygon':# 多个多边形# 遍历每个多边形for polygon in loadGeometry.geoms:# 遍历多边形的所有顶点for z in range(len(polygon.exterior.coords)):# 获取当前顶点的经度和纬度loadLon, loadLat = polygon.exterior.coords[z]# 更新最小经度if loadLon < resultMinLon:resultMinLon = loadLon# 更新最大经度if loadLon > resultMaxLon:resultMaxLon = loadLon# 更新最小纬度if loadLat < resultMinLat:resultMinLat = loadLat# 更新最大纬度if loadLat > resultMaxLat:resultMaxLat = loadLat# 结果输出
print("厦门最小经度：", resultMinLon)
print("厦门最大经度：", resultMaxLon)
print("厦门最小纬度：", resultMinLat)
print("厦门最大纬度：", resultMaxLat)# 创建最小外接矩形的几何对象
bounding_box = box(resultMinLon, resultMinLat, resultMaxLon, resultMaxLat)# 创建 GeoDataFrame
bbox_gdf = gp.GeoDataFrame(geometry=[bounding_box], crs=loadGeoData.crs)# 保存为 Shapefile
output_path = "D:\\data\\厦门最小外接矩形.shp"
bbox_gdf.to_file(output_path)# 可视化
fig, ax = plt.subplots(figsize=(10, 10))# 绘制厦门市的边界
loadGeoData.plot(ax=ax, edgecolor='black', facecolor='lightblue')# 绘制最小外接矩形
min_rect = Rectangle((resultMinLon, resultMinLat), resultMaxLon - resultMinLon, resultMaxLat - resultMinLat,linewidth=2, edgecolor='red', facecolor='none', linestyle='--')
ax.add_patch(min_rect)# 标注最小和最大经度与纬度点
ax.scatter(resultMinLon, resultMinLat, color='red', label='最小边界')
ax.scatter(resultMaxLon, resultMaxLat, color='green', label='最大边界')
ax.annotate(f'({resultMinLon:.6f}, {resultMinLat:.6f})', (resultMinLon, resultMinLat), textcoords="offset points",xytext=(0, 10), ha='center')
ax.annotate(f'({resultMaxLon:.6f}, {resultMaxLat:.6f})', (resultMaxLon, resultMaxLat), textcoords="offset points",xytext=(0, 10), ha='center')# 添加图例
ax.legend()# 设置坐标轴范围
ax.set_xlim(resultMinLon - 0.01, resultMaxLon + 0.01)
ax.set_ylim(resultMinLat - 0.01, resultMaxLat + 0.01)# 显示图形
plt.title('厦门市边界及最小外接矩形')
plt.xlabel('经度')
plt.ylabel('纬度')
plt.grid(True)
plt.show()

忽然发现arcgis/arcgispro也可以生成这个最小外接矩形，直接检索【最小边界几何】，几何类型选择【按面积矩形】或者【按宽度矩形】都可，看需求，组选型选【ALL】这样会基于全部结果形成最外层的矩形；

结果如下；

文章仅用于分享个人学习成果与个人存档之用，分享知识，如有侵权，请联系作者进行删除。所有信息均基于作者的个人理解和经验，不代表任何官方立场或权威解读。


http://www.ppmy.cn/ops/112036.html

相关文章

车载软件架构 --- SOA设计与应用（上）

车载软件架构 --- SOA设计与应用（上）

我是穿拖鞋的汉子，魔都中坚持长期主义的汽车电子工程师。 老规矩，分享一段喜欢的文字，避免自己成为高知识低文化的工程师： 屏蔽力是信息过载时代一个人的特殊竞争力，任何消耗你的人和事，多看一眼都是你的不对。非必要不费力证明自己，无利益不试图说服别人，是精神上的节…
阅读更多...
【重学 MySQL】二十七、七种 join 连接

【重学 MySQL】二十七、七种 join 连接

【重学 MySQL】二十七、七种 join 连接 union 的使用UNION 的基本用法示例UNION ALL 的用法 七种 join 连接代码实现语法格式小结 union 的使用 UNION 在 SQL 中用于合并两个或多个 SELECT 语句的结果集&#xff0c;并默认去除重复的行。如果希望包含重复行&#xff0c;可以使…
阅读更多...
基于Python实现一个庆祝国庆节的小程序

基于Python实现一个庆祝国庆节的小程序

功能&#xff1a; 添加互动功能&#xff1a;允许用户选择不同的祝福语或者查询不同的国庆节信息。动态背景音乐&#xff1a;播放国庆节相关的背景音乐。增加节日小测验&#xff1a;提供一些关于国庆节的趣味小测验&#xff0c;让用户参与。增强图形用户界面 (GUI)&#xff1a;…
阅读更多...
2024世界技能大赛某省选拔赛“网络安全项目”B模块--数字取证解析②（超详细~）

2024世界技能大赛某省选拔赛“网络安全项目”B模块--数字取证解析②（超详细~）

2024世界技能大赛某省选拔赛“网络安全项目”B模块--数字取证解析 PS: 关注鱼影安全第一部分 网络安全事件响应第二部分 数字取证调查任务 3： 网络数据包分析取证解析：第三部分 应用程序安全：需要环境可以私信博主~PS: 关注鱼影安全 模块 B 竞赛项目试题 本文件为：2024世界…
阅读更多...
android 发一个可以下载的的android studio历史版本

android 发一个可以下载的的android studio历史版本

1、Android Studio 下载文件归档 | Android Developers 2、上个图&#xff1a;
阅读更多...
Linux基础---10进程管理

Linux基础---10进程管理

一.查看和关闭进程 1.查看进程 基础指令: ps -efPID 进程编号&#xff0c;PPID 父进程编号&#xff0c; CMD命令名称 进阶指令–查看进程的树形结构&#xff1a; yum install psmisc -y #首先安装psmisc后可直接使用pstreepstree2.关闭进程 要想关闭某个或多个进程需要知道…
阅读更多...
Playwright 与 Selenium对比

Playwright 与 Selenium对比

通过这篇关于 Playwright 与 Selenium 的文章&#xff0c;我们将更容易理解 Playwright 和 Selenium 之间的关键区别&#xff0c;并找出哪个工具可能更适合您的需求。 在自动化测试工具方面&#xff0c;Playwright 和 Selenium 都是软件测试人员使用的强大的 Web 自动化工具。它…
阅读更多...
代码随想录打卡Day32

代码随想录打卡Day32

今天有点事&#xff0c;先做一题&#xff0c;剩下的明天补。 509. 斐波那契数 这道题目太简单了&#xff0c;递归几行代码就结束了&#xff0c;用动态规划做也可以&#xff0c;主要是学习一下动态规划五部曲。 这是递归的代码 class Solution { public:int fib(int n) {//确…
阅读更多...
最新文章

Copyright @ 2022~2023