geopandas公交站点覆盖率分析

目标

分析研究区域内公交站点的覆盖率

使用的文件：

• 研究区域边界.shp（polygon类型数据）
• 站点.shp（point类型数据）

代码：

# 导入所需的库
import geopandas as gpd
import pandas as pd
import matplotlib.pyplot as plt
import numpy as np

# 导入文件
shp_border = '研究区域边界.shp'
shp_pts = '站点.shp'

# 读取文件数据
data_border = read_file(shp_border, encoding='utf-8')
data_pts = read_file(shp_pts, encoding='utf-8')

# 对数据进行坐标转换
reprojected_data_border = data_border.to_crs(crs="+proj =utm +zone=48 +ellps=WGS84 +datum=WGS84 +units=m +no_defs")
reprojected_data_pts = data_pts.to_crs(crs="+proj =utm +zone=48 +ellps=WGS84 +datum=WGS84 +units=m +no_defs")

proj =utm：声明投影方法为墨卡托
zone=48：声明对应背48区
ellps=WGS84：声明椭球面，使用WGS84
datum=WGS84：声明基准面为WGS84
units=m：声明坐标系单位设置为米

# 在站点周边生成半径为 500m 的缓冲区
coverage = reprojected_data_pts.buffer(500)
# 将所有站点的缓冲区合并为一个图形
cov_union = coverage.unary_union
# 取得缓冲区与区域边界的交集（即在边界内的缓冲区）
holes = reprojected_data_border['geometry'].intersection(cov_union)
# 绘图
ax = holes.plot(color='#ed2f6a')
ax = reprojected_data_border.boundary.plot(color='#74787a', ax=ax)
ax.axis('off')

第10行plot()中的参数 ax=ax 表示将这张图与上一张图绘制在一个画布上

绘制效果图：

# 统计分析
s1 = holes.area.sum()     # 计算区域内缓冲区的总面积
s2 = reprojected_data_border.area.sum()     # 计算研究区域的总面积
ratio = s1 / s2     # 计算站点 500m 覆盖率

若直接输出 holes.area，输出结果为各行政区中的缓冲区面积

geopandas官方说明文档：

https://geopandas.org/index.html