简介面向城市规划、地理信息分析与区域研究者的合肥市功能区区划矢量数据集在常规区县行政边界之外补充了湖滨新区、经济开发区、经济开发区北区、高新区与新站区等特定功能区的空间划分可用于城市空间结构分析、功能分布制图及规划政策研究。压缩包共8个文件约276KB以SHP矢量主体配合DBF属性表、PRJ坐标系统、CPG编码说明及SBN/SBX索引等辅助文件构成兼顾几何数据、区域属性与投影定位导入ArcGIS、QGIS等软件即可查看与叠加分析。已有298人学习适合需要快速获取合肥功能区边界底图、开展空间统计或制作专题图的初级至中级GIS使用者参考。1. 先别急着双击功能区划 shp 数据的真实使用场景拿到一份功能区划 shp 数据第一反应不是双击打开而是先看 .prj、.dbf 字段宽度和 .cpg 编码。很多翻车现场不是分析模型错而是坐标系被静默改写面积统计从平方公里变成“玄学”。这份 GIS 数据集适合做底图核对、功能区叠加统计、专题图出图和空间服务入库不适合直接拿经纬度坐标算面积。新手可以按后面的加载、校验、导出流程走熟手重点看字段截断、拓扑自相交和投影选择这几个边界。它解决的是“手里有一份区划边界怎样安全地进管线”的问题而不是替代业务分类标准。2. 拆开 shp 文件组坐标系统、字段与 .cpg 编码先对齐2.1 文件组不是单个 .shp先认 .prj、.dbf、.cpg 的分工shapefile 在磁盘上不是一个孤立文件而是成套出现。只看 .shp 能打开不代表属性、投影和中文编码都完整。我一般会先列目录把同名的 .shp、.shx、.dbf、.prj、.cpg 成组核对再谈加载。缺少 .shx 时部分驱动还能读但空间索引会慢缺少 .dbf 时几何在、属性没了缺少 .prj 时软件可能默认按经纬度处理缺少 .cpg 时中文字段名或中文属性容易乱码。功能区划数据经常要交付给不同软件文件组齐不齐比单个文件能不能打开更重要。扩展名作用缺失后的典型表现.shp存几何坐标主文件缺失图层无法加载.shx几何索引能读但查询慢部分工具直接报错.dbf存属性表只剩图形字段全部丢失.prj存坐标参考CRS 为空面积和距离不可信.cpg存属性编码中文乱码字段名显示异常# 在数据目录里先看文件组是否齐全 ls -lh # 只统计扩展名快速确认 .shp/.shx/.dbf/.prj/.cpg 是否成组出现 find . -maxdepth 1 -type f | sed s/.*\.// | sort | uniq -c代码逻辑很直接第一条命令看文件是否同名成组第二条命令按扩展名计数避免只盯着 .shp。参数上-maxdepth 1表示只看当前目录防止把子目录里的旧版本混进来sed s/.*\.//取扩展名sort | uniq -c做计数。如果 .prj 和 .cpg 各只有一份而 .shp 有多份说明目录里可能混了多个版本先分目录再处理。注意shapefile 字段名最大 10 个字符这不是导出时才出现的问题读入后就要检查字段名是否已经被截断。2.2 坐标系统选型交换用地理坐标算面积用投影坐标功能区划数据的第一层选型是坐标系统。地理坐标适合交换、展示和 Web 地图底图对齐单位是度投影坐标适合面积、长度和叠加统计单位是米。用经纬度直接算面积结果会随纬度变化看起来像“玄学”本质是单位不对。常见做法是读入后先确认gdf.crs不为空展示时保留数据自带 CRS做面积统计前再转到数据生产方指定的投影坐标系或者用本地等面积投影。不要在不确认 .prj 的情况下强行set_crs那只是给错误坐标贴标签。import geopandas as gpd from pathlib import Path shp Path(data/function_zone.shp) gdf gpd.read_file(shp, encodingutf-8) # 中文乱码时改用 gbk 再试 print(CRS:, gdf.crs) print(几何类型:, gdf.geom_type.value_counts().to_dict()) print(字段:, gdf.columns.tolist()) print(前 3 行:) print(gdf.head(3)) if gdf.crs is None: raise ValueError(缺少 .prj 或 CRS 为空先确认坐标系不要直接算面积) # 仅演示面积统计路径正式项目用数据生产方指定的投影坐标系 gdf_metric gdf.to_crs(gdf.estimate_utm_crs()) gdf_metric[area_km2] gdf_metric.geometry.area / 1e6 print(gdf_metric[[area_km2]].describe())代码逻辑先读文件组再打印 CRS、几何类型、字段和样例行gdf.crs is None是硬性闸门防止把经纬度当平面坐标。encodingutf-8是属性表编码参数读出来乱码时换gbk不要靠改系统编码碰运气。estimate_utm_crs()只是演示一个可用的米制投影正式统计要按数据说明或本地投影参数走/1e6把平方米换成平方公里。参数上area单位跟随 CRS投影坐标下才是平方米地理坐标下是平方度不能混用。提示展示图层可以保留地理坐标统计图层另存一份投影坐标两个图层不要来回to_crs十几次既慢又容易把字段类型带乱。2.3 字段结构核对10 字符限制与数据字典功能区划数据的字段通常包含区划标识、名称、类型、父级标识、面积、更新时间这类信息但具体字段名要以实际数据为准。读入后先做字段体检字段名长度、空值数量、唯一值数量、最大字符串长度、数值字段是否混入文本。shapefile 的 .dbf 字段名最多 10 个字符超过就会被截断导出后很容易出现zone_name变成zone_nam这类对不上数的情况。字段类型也会漂移比如面积字段在源数据里是文本读进来后dtype是 object后面做汇总就会报错。import pandas as pd field_report [] for col in gdf.columns: if col gdf.geometry.name: continue series gdf[col] field_report.append({ field: col, dtype: str(series.dtype), null_count: int(series.isna().sum()), n_unique: int(series.nunique(dropnaTrue)), max_len: int(series.astype(str).str.len().max()), }) report pd.DataFrame(field_report).sort_values(field) print(report.to_string(indexFalse)) # 字段名超过 10 个字符时先做映射再导出 shapefile rename_map { zone_identifier: zone_id, zone_category_name: zone_type, parent_zone_identifier: parent_id, } safe_gdf gdf.rename(columnsrename_map)代码逻辑逐字段统计类型、空值、唯一值和最大长度先把“字段长什么样”说清楚再决定要不要重命名。rename_map是字段名映射表左边是源字段右边是导出用短名右边尽量不超过 10 个字符。参数上dropnaTrue让唯一值统计不被空值干扰astype(str)用于看最大字符串长度但不要把数值字段长期转成字符串。字段核对做完后最好把映射表和数据字典一起交付否则下一位使用者还得重新猜字段含义。常见字段角色示例短名核对重点区划唯一标识zone_id是否唯一、是否为空、是否被截断区划名称zone_name中文编码、最大长度、重名情况区划类型zone_type分类值是否在数据字典内父级标识parent_id是否能回连、是否存在孤儿记录面积字段area_km2单位是否明确、是否为文本3. 把区划数据塞进 GIS 管线CRS、几何有效性与拓扑预检3.1 加载后的三道预检CRS、几何类型、空几何数据能读进来只是第一步。加载后我固定跑三道预检CRS 是否为空几何类型是否统一空几何有多少。功能区划数据常见的是 Polygon 和 MultiPolygon 混存MultiPolygon 本身不是错误但如果后续脚本只处理 Polygon就会漏掉多部件要素。空几何更隐蔽面积统计时它贡献 0出图时看不见入库后却占记录。预检的目标不是把数据改漂亮而是先知道问题规模再决定清洗策略。import geopandas as gpd gdf gpd.read_file(data/function_zone.shp, encodingutf-8) print(CRS 是否为空:, gdf.crs is None) print(几何类型分布:, gdf.geom_type.value_counts().to_dict()) print(空几何数量:, int(gdf.geometry.is_empty.sum())) print(空几何记录索引:, gdf.index[gdf.geometry.is_empty].tolist()[:20]) # 统一几何类型前先保留原始类型便于追溯 gdf[geom_type_raw] gdf.geom_type gdf[geometry] gdf.geometry.make_valid() print(修复后无效几何:, int((~gdf.geometry.is_valid).sum()))代码逻辑先统计 CRS、几何类型、空几何再用make_valid()修复无效几何。geom_type_raw保留原始类型方便回查哪些要素原本是多部件。make_valid()可能把自相交多边形拆成多部件或几何集合修复后还要再看一次几何类型不能假设它一定还是原来的 Polygon。参数上is_empty判断几何是否为空is_valid判断拓扑是否有效两者不是一回事空几何不一定无效无效几何也不一定为空。注意buffer(0)是常见的修复偏方但它会改变几何边界面积和相邻关系可能跟着变能追溯的修复优先用make_valid()。3.2 拓扑检查自相交、重叠与悬挂节点功能区划数据最怕“看起来是相邻实际有重叠”。自相交会让叠加分析报拓扑异常重叠会让同一块面积被统计两次悬挂节点会让相邻边界对不齐。检查重叠时不能只看intersects因为相邻多边形共享边界也会返回 True但交集是线面积为 0。我的做法是先转到投影坐标再用空间索引筛候选对计算交集面积超过阈值才记为重叠。阈值不要设成 0浮点误差会让共享边界产生极小面积。import geopandas as gpd gdf gpd.read_file(data/function_zone.shp, encodingutf-8) gdf_m gdf.to_crs(gdf.estimate_utm_crs()) issues [] sindex gdf_m.sindex for i, geom in enumerate(gdf_m.geometry): for j in sindex.query(geom): if j i: continue inter geom.intersection(gdf_m.geometry.iloc[j]) if not inter.is_empty and inter.area 1e-6: # 1e-6 平方米按数据精度调整 issues.append({ left: i, right: j, overlap_m2: inter.area, }) print(重叠候选对数:, len(issues)) print(issues[:10])代码逻辑用sindex.query缩小候选范围避免全量两两求交j i防止重复对intersection后看面积面积大于阈值才算面重叠。参数上1e-6是平方米阈值数据精度高时可以调小数据精度低时调大不要直接写 0。相邻共享边界产生的是线交集面积为 0不会进入结果。悬挂节点和缝隙检查更依赖拓扑规则常见做法是在 QGIS 拓扑检查器里设“不能有悬挂点”“不能有缝隙”或者用 GDAL 的拓扑校验工具跑一遍再回到属性表定位。3.3 出图与字段标注先验证再谈样式出图不是最后一步而是验证手段。把区划按类型填色、把边界加白线、把异常要素单独标红比盯属性表更容易发现空洞、飞地和重叠。出图时先不追求配色先看三件事边界是否闭合、相邻面是否连续、图例分类是否和字段一致。字段标注要控制长度中文名称太长会叠在一起可以先截断或者只标关键要素。出图脚本里不要写死字段名用变量接换一份数据只改配置。import matplotlib.pyplot as plt # 按实际字段替换 zone_type字段不存在时先回看 2.3 的字段核对结果 fig, ax plt.subplots(figsize(10, 8)) gdf.plot( axax, columnzone_type, legendTrue, edgecolorwhite, linewidth0.5, missing_kwds{color: lightgray, label: 缺失分类}, ) ax.set_axis_off() plt.savefig(output/zone_preview.png, dpi200, bbox_inchestight) print(已输出预览图)代码逻辑用column做分类填色edgecolor描边missing_kwds把空分类单独标灰避免空值被默认色掩盖。dpi200保证预览图能看清边界bbox_inchestight去掉多余白边。参数上linewidth不要太大否则相邻边界会糊成一条legendTrue适合分类少的情况分类多时图例会挤可以改成单独出图例。出图后如果看到某块颜色断开先回查几何有效性和重叠检查不要急着调色板。4. 叠加统计与导出sjoin、overlay 和 PostGIS 入库的边界4.1 叠加分析前的对齐CRS、字段类型与几何有效性叠加分析翻车十有八九是 CRS 不一致或几何无效。点数据和区划面做空间连接先把点转到面同一 CRS面数据和网格做 overlay先确保两边都是有效几何。字段类型也要对齐区划标识在一张表里是字符串在另一张表里是整数连接后会出现001和1对不上。我的习惯是先做一个小样本抽 100 条点、10 个面跑通再上全量。全量跑之前把几何有效性、CRS、连接键类型三项打成一张检查表省得跑半小时后报错。import geopandas as gpd zones gpd.read_file(data/function_zone.shp, encodingutf-8) points gpd.read_file(data/business_points.geojson) # CRS 对齐以区划图层为准 points points.to_crs(zones.crs) # 连接键类型对齐两边都转成字符串避免 001 和 1 对不上 zones[zone_id] zones[zone_id].astype(str).str.zfill(6) points[zone_id] points[zone_id].astype(str).str.zfill(6) # 几何有效性闸门 zones[geometry] zones.geometry.make_valid() points[geometry] points.geometry.make_valid() print(zones.crs, points.crs) print(zones[zone_id].head().tolist(), points[zone_id].head().tolist())代码逻辑先对齐 CRS再对齐连接键类型最后修复几何。zfill(6)是示例表示按实际编码位数补零不要照抄位数。参数上astype(str)把数值和字符串统一成文本避免前导零丢失make_valid()处理自相交和重复点。连接键对齐后再跑空间连接结果才可信。4.2 用 sjoin 与 overlay 做功能区统计点落在功能区里用sjoin做空间连接predicatewithin表示点完全落在面内。边界上的点会被within排除如果业务要求边界点也算就改用intersects然后按点和区划的唯一键去重。面与网格叠加用overlayhowintersection保留相交部分keep_geom_typeTrue过滤掉线、点结果。叠加后算面积一定要用投影坐标不要拿经纬度直接算。import geopandas as gpd zones gpd.read_file(data/function_zone.shp, encodingutf-8) points gpd.read_file(data/business_points.geojson).to_crs(zones.crs) # 点落在面内统计边界点严格来说不算 within按业务决定是否改 intersects joined gpd.sjoin(points, zones, howinner, predicatewithin) count_by_zone joined.groupby(zone_id).size().reset_index(namepoint_count) print(count_by_zone.head()) # 面与网格叠加先统一到投影坐标再算面积 zones_m zones.to_crs(zones.estimate_utm_crs()) grid gpd.read_file(data/grid.shp).to_crs(zones_m.crs) inter gpd.overlay(zones_m, grid, howintersection, keep_geom_typeTrue) inter[inter_area_km2] inter.geometry.area / 1e6 summary inter.groupby([zone_id, grid_id])[inter_area_km2].sum().reset_index() print(summary.head())代码逻辑sjoin把点属性挂到面上groupby统计每个区划的点数overlay把区和网格切成相交面再按区划和网格汇总面积。参数上predicatewithin和intersects的差别在边界点keep_geom_typeTrue防止输出里混入线几何。zones.estimate_utm_crs()只是演示正式项目用本地投影。inter_area_km2必须在投影坐标下计算否则单位不对。提示sjoin后如果发现一个点挂了多个区划先查区划面是否重叠不要直接drop_duplicates把问题盖住。4.3 导出 GeoJSON、PostGIS 与字段名映射导出格式决定字段名和编码策略。GeoJSON 字段名没有 10 字符限制中文属性用 UTF-8适合交付和前端加载shapefile 字段名受限导出前要做短名映射PostGIS 适合入库和空间查询几何列、SRID、索引要一起建。导出前先make_valid()再检查空几何最后按目标格式选驱动。入库时不要用默认表名和默认几何列名显式指定后面查起来省事。import geopandas as gpd zones gpd.read_file(data/function_zone.shp, encodingutf-8) zones[geometry] zones.geometry.make_valid() # GeoJSON 导出字段名不受 10 字符限制中文用 UTF-8 zones.to_file(output/function_zone.geojson, driverGeoJSON, encodingutf-8) # shapefile 导出字段名先映射到短名避免截断 rename_map { zone_identifier: zone_id, zone_category_name: zone_type, parent_zone_identifier: parent_id, } zones.rename(columnsrename_map).to_file( output/function_zone_final.shp, driverESRI Shapefile, encodingutf-8, ) # PostGIS 入库连接串按实际环境替换 from sqlalchemy import create_engine engine create_engine(postgresql://user:passwordlocalhost:5432/gisdb) zones.to_postgis(function_zone, engine, if_existsreplace, indexFalse, chunksize500)代码逻辑先修几何再按格式导出。GeoJSON 用driverGeoJSONshapefile 用driverESRI Shapefile并配合字段映射PostGIS 用to_postgis写库。参数上encodingutf-8控制属性编码if_existsreplace会覆盖同名表正式环境先改成append或手工建表chunksize500控制批量写入大小。连接串不要写进代码仓库用环境变量接。入库后记得给几何列建 GiST 索引否则空间查询会慢。5. 避坑与排查五类 shp 异常的现象、原因与解决5.1 坐标系缺失被当成经纬度面积统计玄学现象图层能打开面积字段算出来忽大忽小和地图上看到的范围对不上。原因.prj 缺失或 CRS 为空软件把坐标当成经纬度面积单位变成平方度或者有人用set_crs强行贴了一个 EPSG坐标数值没变标签错了。解决先读gdf.crs为空就找数据生产方要 .prj 或坐标说明不要猜 EPSG。确认 CRS 后面积统计转到投影坐标单位统一到平方米或平方公里。已经入库的表查geometry_columns或ST_SRID别只看图层属性面板。if gdf.crs is None: raise ValueError(CRS 为空停止面积统计) # 确认后再投影不要用 set_crs 硬贴标签 gdf_m gdf.to_crs(gdf.estimate_utm_crs())5.2 中文属性乱码.cpg 缺失与编码判断现象字段名显示成问号或乱码中文区划名称变成一串不可读字符。原因.dbf 的属性编码和读取端不一致.cpg 缺失时驱动只能猜UTF-8 和 GBK 混用也会出现部分字段正常、部分字段乱码。解决先试encodingutf-8再试gbk读出来后看中文名称是否可读导出时显式写encodingutf-8并保留 .cpg。不要用文本编辑器直接改 .dbf容易把字段长度和编码一起改坏。交付前抽 10 条中文名称核对比看编码声明可靠。5.3 字段名截断与类型漂移导出后对不上数现象源数据里字段叫zone_category_name导出 shapefile 后变成zone_cate下游按原名取字段直接报 KeyError。原因shapefile 字段名最多 10 个字符驱动会截断数值字段里混了暂无整列变成 object。解决导出前做字段映射表短名不超过 10 个字符数值字段先做to_numeric(errorscoerce)把脏值转成空值并记录数量。映射表和数据字典一起交付别让下游猜。字段类型在入库前再核对一次PostGIS 的 numeric 和 text 不会自动帮你兜底。5.4 几何自相交与多部件叠加失败、面积偏大现象overlay报拓扑异常或者叠加后总面积大于原面积。原因多边形自相交、相邻面重叠、多部件要素重复计算。解决先make_valid()再用 3.2 的空间索引方法查重叠对按业务决定是合并、裁切还是回退源数据。多部件要素不要直接explode后不去重 explode 会把一个区划拆成多条记录统计口径就变了。面积汇总前按区划唯一标识分组确认一个标识只保留一条汇总记录。5.5 空几何、Z/M 值与图层锁定发布前最后一道坎现象导出 GeoJSON 时报几何类型不支持或者入库后地图上少了几块。原因空几何、Z/M 值、几何集合混在图层里Windows 下文件被 QGIS 或资源管理器占用导出覆盖失败。解决导出前统计is_empty删掉或补上空几何有 Z/M 值时先force_2d()或按目标格式转换导出前关闭占用文件的程序输出到新目录再替换。发布前跑一遍“记录数、有效几何数、空几何数、面积合计”四项对账少一项都别急着交付。6. 交付前跑一套校验脚本字段、CRS、几何、面积四项对账最后一章落到一个具体习惯交付前跑校验脚本不靠肉眼。脚本不复杂四个数必须对得上记录数、CRS 是否明确、有效几何数、面积合计。记录数看有没有丢要素CRS 看有没有丢 .prj有效几何数看make_valid()后是否还有无效项面积合计看投影后总面积是否和源数据说明在同一个量级。四项里任何一项对不上先停下来查不要先做样式。import geopandas as gpd gdf gpd.read_file(output/function_zone_final.shp, encodingutf-8) gdf_m gdf.to_crs(gdf.estimate_utm_crs()) check { record_count: len(gdf), crs: str(gdf.crs), valid_geom_count: int(gdf.geometry.is_valid.sum()), empty_geom_count: int(gdf.geometry.is_empty.sum()), area_km2_sum: round(float(gdf_m.geometry.area.sum() / 1e6), 3), } print(check)代码逻辑重新读导出后的成果不读中间过程文件避免“源数据没问题、导出坏了”被漏掉。record_count对记录数crs对坐标声明valid_geom_count和empty_geom_count对几何area_km2_sum对面积量级。参数上estimate_utm_crs()只用于快速对账正式面积以项目指定投影为准round(..., 3)只是让输出好读不要拿四舍五入后的数做精确比对。6.1 交付清单里我固定放的三个附件除了 shp 文件组我会放三个附件字段映射表、坐标说明、校验结果截图或文本。字段映射表写清源字段、导出短名、类型、是否必填坐标说明写清地理坐标和投影坐标各自用途校验结果写记录数、有效几何数、面积合计。这样下游拿到数据不用先猜字段也不用重新问坐标系。附件不追求长一页能看完最好。6.2 入库前再看一眼空间索引如果成果要进 PostGIS建完表先看几何列 SRID 和索引。没有 GiST 索引空间查询会从秒级掉到分钟级。常见做法是建表后执行CREATE INDEX ON function_zone USING GIST (geom);再用ST_SRID(geom)确认 SRID 和业务约定一致。索引名按表名起别用默认名后面迁移好找。从那以后我每次拿到 shp 都先跑字段、CRS、几何三连再谈样式和出图这份功能区划 shp 数据集也一样先按第 2 章核对文件组和字段再按第 3 章做拓扑预检最后用第 6 章的四项对账收口。需要这份数据集的话按资源目录取用先跑校验再入库能少走不少弯路。希望帮到你。本文还有配套的精品资源点击获取
