四川Shp矢量数据处理指南:从坐标系、行政区划到GeoJSON导出
简介四川地图SHP文件是一套覆盖四川省的矢量地理数据集面向GIS从业者、规划研究人员及空间数据分析学习者可用于地图制图、空间查询与专题分析。压缩包共200个文件包含shp几何文件、dbf属性表、shx索引、prj投影定义以及xml元数据和mxd工程文档总大小118.18MB。其中shp/shx/dbf组合存储点线面要素及属性prj文件确保坐标系统正确mxd便于直接在ArcMap中加载查看。内容涵盖行政边界、道路、水系、乡镇村道以及各类兴趣点如政府机关、医疗、餐饮、购物等能够支撑缓冲区分析、叠图分析或应急响应等场景。目前已有3588人学习下载适合需要快速获取四川省基础地理数据并开展空间可视化与分析的读者。1. 四川地图shp一张图背后藏着省界、市州界、区县界三套数据搜“四川地图shp文件”结果里十个有八个是从全国矢量数据里裁出来的直接拿来用往往能画出图却接不上业务数据。问题通常集中在四个点用的是哪个行政层级、坐标系是不是统一基准、属性表编码对不对、区划有没有更新。四川本身横跨第一级和第二级地形阶梯市州之间边界复杂县级行政区数量多且经历过撤县设区、合并改名对这些问题格外敏感一份四川地图shp要在分析和开发里被反复依赖第一步就是把数据层级和坐标基准这两个开关掰对。这篇文章讲一套我反复验证过的流程选层级、选坐标系、用Python读取绘制、做统计关联、对照四个常见坑排查最后转出Web端能直接用的轻量GeoJSON。2. 选对四川Shp的行政层级和坐标系21个市州、183个区县怎么取舍2.1 省级、市州级、区县级Shp的差异字段、精度和应用场景拿到一份四川shp先别急着画打开属性表看记录数。要素数量直接暴露数据的行政级别省级shp通常只有一条记录是整个四川的一个大面适合做轮廓底图、地图裁剪边界或者和省外区域做对比市州级shp包含21个左右的市州要素适合做专题着色和区域统计区县级shp包含183个左右的县级行政区要素适合做县域经济分析、人口分布和点位落图。很多人下载的时候只看到“四川”两个字不确认级别结果做县级统计时拿的是市州面merge之后全是空值。行政层级要素个数四川典型用途省级1条地图轮廓、裁剪边界、邻省对比市州级21条左右专题图着色、区域汇总、大屏展示区县级183个左右县级统计、经济人口分析、点位匹配来源上我一般会优先找两类一类是测绘部门公开发布的行政区域界线矢量成果这类文件通常带元数据说明和更新时间另一类是开发者整理的全国省市县合集字段结构统一适合用province字段筛选出四川。后者风险在于发布人可能把已经撤销合并的老县留在文件里所以拿到后要核对行政区划代码。比较稳的做法是下全国合集按省级字段过滤出四川字段统一后续换省做图也不需要重写代码。筛选时要注意有的合集里“四川”字段写成“省份”“省名称”或拼音SICHUAN有的是多部件要素把甘孜、阿坝这种大面积自治州拆成几十条记录。不管哪种过滤后一定要reset_index并且检查记录数是否符合该级别预期否则后面绘图会莫名出现一堆零碎小面。2.2 坐标系选择与投影转换WGS84、CGCS2000还是Albers等积坐标系是shp使用里最容易被忽略的坑。多数开源数据集以WGS84经纬度存储而国内测绘成果常用CGCS2000两者在大比例尺地图上差别可以忽略但和在线底图叠加时任何形如“底图偏了几百米”的问题根源几乎都和坐标基准有关。我的选择逻辑很简单如果数据只用于Web可视化全程保留WGS84经纬度也就是EPSG:4326如果要做面积统计、距离计算、缓冲区分析就转投影坐标系。四川经度跨度约11度、纬度跨度约8度东西跨度明显超过一个UTM分带用UTM 47N或48N都会有较大变形所以面积计算我会优先选Albers等积投影中央经线取105°E标准纬线取25°N和47°N这是中国区域制图常用的配置。使用场景推荐坐标系原因Web底图叠加、GeoJSON输出WGS84经纬度EPSG:4326与在线底图对齐无需额外转换打印出图、成果提交CGCS2000EPSG:4490符合国内制图规范面积、距离、缓冲区计算Albers等积投影四川东西跨度大等积投影面积变形小拿到数据后先看有没有.prj文件。如果缺失GIS软件会默认按WGS84去猜这时坐标偏移问题就出现了。用文本编辑器打开.prj能看到坐标系完整定义如果没有.prj可以从坐标数值判断四川经纬度大约是东经97°到108°、北纬26°到34°如果看到的坐标是几百万的大数字那是投影坐标不是经纬度不能直接当经纬度用。3. 用Python读取并画四川Shp从过滤全国数据到导出地图3.1 先读数据再看边界GeoPandas读取与几何有效性检查Python侧处理shp我基本只用GeoPandas。读取的第一步不是plot而是先确认读进来的记录数、字段名和边界范围这三项能暴露大部分数据问题。下面这段是从全国市州合集里筛四川的常规操作import geopandas as gpd # 读取全国市州界数据encoding参数仅影响属性表不影响几何 gdf_nation gpd.read_file(nation_city.shp, encodingutf-8) print(gdf_nation.columns.tolist()) # 按省级字段筛选四川省字段名以实际数据为准 gdf gdf_nation[gdf_nation[province] 四川省].copy() gdf gdf.reset_index(dropTrue) print(记录数:, len(gdf)) print(字段:, gdf.columns.tolist()) # total_bounds返回(minx, miny, maxx, maxy) print(边界范围:, gdf.total_bounds)这段代码里encoding只决定.dbf属性表的读取方式和几何坐标无关如果属性表是GBK编码而这里写了utf-8记录数不会受影响但中文名会乱码。total_bounds非常关键四川的正常范围应该是西到约97、南到约26、东到约108、北到约34如果发现边界值完全对不上说明坐标基准判断错了后面所有叠加分析都会偏。几何有效性也建议写进流程。shp来源复杂很多数据在编辑过程中产生了自相交多边形这类错误面在图上肉眼很难发现但会在面积计算和叠加分析时给出离谱结果# 统计无效几何数量 invalid gdf[~gdf.geometry.is_valid] print(无效几何数量:, len(invalid)) # 若有无效几何用shapely的make_valid修复 from shapely.validation import make_valid if len(invalid) 0: gdf.loc[invalid.index, geometry] gdf.loc[invalid.index, geometry].apply(make_valid) print(修复后无效几何数量:, (~gdf.geometry.is_valid).sum())make_valid会把自相交的多边形拆成合法几何个别情况下会返回MultiPolygon这会改变要素类型但不会影响后续绘图和面积计算。修复后再跑一次is_valid是为了确认没有遗漏这条检查建议作为固定步骤因为我在实际项目里遇到四川区县级数据时无效几何的比例常在1%到5%之间。3.2 把省界和市州界画在一起叠加绘制与样式参数调整数据检查完画图就是水到渠成的事。绘制四川市州界最简单的代码是这样import matplotlib.pyplot as plt fig, ax plt.subplots(figsize(8, 8)) # facecolor是填充色edgecolor是边线色linewidth是边线宽度 gdf.plot(axax, facecolor#f0f0f0, edgecolor#444444, linewidth1.2) ax.set_title(四川省市州界地图) ax.set_axis_off() plt.tight_layout() plt.savefig(sichuan_city.png, dpi300) plt.show()如果你手里同时有省级轮廓和市州边界想叠在一张图上注意两条数据必须处于相同坐标系。我把省级面设为无填充、粗边线市州面只画细边线这样既能看到全省轮廓又能看到内部划分# 读取省级轮廓 gdf_prov gpd.read_file(sichuan_province.shp, encodingutf-8) # 先画市州面再叠省界线 ax gdf.plot(facecolor#f0f0f0, edgecolor#888888, linewidth0.6, axax) gdf_prov.plot(axax, facecolornone, edgecolor#000000, linewidth2.0)样式调整里比较实操的参数是linewidth和dpi。屏幕预览用linewidth0.6、dpi100足够但导出印刷或汇报用图建议市州线加粗到1.2以上、dpi设置300否则PDF里放大后边线会出现明显的锯齿。还有一点如果两个shp的坐标系不同Geopandas在同一个ax上叠加时不会自动重投影你得先手动统一否则画出来的图就会有一层是歪的。4. 基于四川Shp做统计和关联面积核算与区划代码合并4.1 面积怎么算才准投影转换后统计各县面积很多人直接在WGS84经纬度坐标下调用geometry.area算出来的是一个没有物理意义的数字因为经度方向的长度在不同纬度下被拉伸了。四川南部和北部的纬度差约8度直接用经纬度算面积误差会非常明显。正确的做法是先转Albers等积投影再计算面积。# Albers等积投影定义中央经线105E标准纬线25N和47N albers (projaea lat_125 lat_247 lat_00 lon_0105 x_00 y_00 ellpsWGS84 unitsm no_defs) gdf_proj gdf.to_crs(albers) # area结果单位是平方米除以1e6转成平方公里 gdf[area_km2] gdf_proj.geometry.area / 1_000_000 print(gdf[[name, area_km2]].sort_values(area_km2, ascendingFalse).head(5)) print(总面积:, gdf[area_km2].sum())这里有个很好的自检方法四川总面积约48.6万平方公里如果你对市州级shp求面积总和结果应该在48.6万左右浮动。如果偏得离谱一是检查坐标是否还是经纬度状态二是检查是否有重复要素或无效几何。我遇到过一份数据里甘孜、阿坝两个大面各出现了两次面积总和直接多出三分之一这种问题只有通过总和校验才能暴露。面积结果建议直接写回原始表后续做人均GDP、人口密度这类计算就不用反复重投影了。写回时保留原始几何不要用投影后的几何覆盖原字段。4.2 业务数据怎么挂到地图上用行政区划代码而不是名称把Excel或CSV里的统计数据和shp关联最稳的键是行政区划代码而不是县名。四川县级名称经历过不少调整有的县改成了区有的县级市升了级别名称一变旧数据按name关联就会全部落空。行政区划代码是6位数字省级两位、市州两位、区县两位变更时整体调整只要两边用的都是同一套标准匹配准确率接近百分之百。import pandas as pd # 业务表至少有adcode和value两列 biz pd.read_csv(county_data.csv) # 统一转成6位字符串避免数字前导0被读掉 biz[adcode] biz[adcode].astype(str).str.zfill(6) gdf[adcode] gdf[adcode].astype(str).str.zfill(6) merged gdf.merge(biz, onadcode, howleft) # 统计匹配不上的区县 unmatched merged[merged[value].isna()] print(无法匹配的区县数量:, len(unmatched)) print(unmatched[name].tolist())必须提醒一个细节有些shp里的adcode字段本身就已过时比如曾经存在的县已经合并代码在旧表里查无此项。这时要区分两种情况一种是代码变更查最新的行政区划代码表把旧码更新为新码另一种是实质合并需要根据合并规则手动映射。最忌直接删掉匹配不上的行因为缺失区县恰恰是分析中最容易出问题的那部分。实际执行时我还会加一步检查merge之后记录数是否和原shp一致。如果用了howinner匹配不上的区县会静默消失记录数减少而你毫无察觉所以这里用howleft让不匹配的区县保留名称和空值方便人工排查。5. 四川Shp四大疑难排查坐标偏移、乱码、区划过期、几何自相交5.1 和在线底图对不上坐标整体偏移了几百米现象四川shp在本地画着正常但叠加在线地图后省界整体朝某个方向平移了几百米甚至一公里市州边界也跟着偏完全无法对上道路和地名。原因分两类一类是数据发布时没有保留.prj文件GIS软件默认按WGS84读取实际原始坐标系并不是WGS84另一类是数据来自在线地图导出经纬度已经做过坐标偏移处理属性表里的坐标值并非真实经纬度。解决先用文本编辑器打开.prj查看坐标系定义如果没有.prj找一个四川境内标志明显的地物把shp坐标和在线底图坐标做点对点比较判断偏移方向和距离。最直接的办法是换用同一坐标系来源的底图或者找到该数据的原始发布说明。对业务项目来说与其在shp上做坐标纠偏不如统一用WGS84基准的在线底图省掉后续所有坐标系扯皮。5.2 属性表中文全是乱码的解决先探测dbf编码再读现象用GeoPandas读取四川区县shpname字段显示成“镣”“鍘”之类毫无意义的汉字组合。原因绝大多数国产shp属性表采用GBK或GB18030编码而GeoPandas在缺少.cpg文件时默认按UTF-8解码中文全部乱码。解决读取前用chardet探测一下编码再传给read_fileimport chardet from pathlib import Path # 只读取dbf文件头部字节来探测编码 raw Path(sichuan_county.dbf).read_bytes()[:10000] result chardet.detect(raw) print(result) # 将探测到的编码传给read_file gdf gpd.read_file(sichuan_county.shp, encodingresult.get(encoding, utf-8))如果用chardet探测出GB2312或GBK直接写encodinggbk或encodinggb18030都行。问题解决后建议顺手把数据重导出为UTF-8编码避免下一个接手的人再踩一遍坑gdf.to_file(sichuan_county_utf8.shp, encodingutf-8)这里要注意重导出会同时改写.shp、.dbf和.prj等文件务必另存为新文件不要覆盖原始数据。5.3 边界过期导致数据合并不上核对最新行政区划代码现象业务数据里有某县2020年以后的GDP但shp属性表里查不到这个县map展示时始终缺一块区域。原因近几年四川经历了撤县设区、设立县级市等多轮区划调整一些老县被撤销或合并。旧版shp发布于调整之前自然没有新单位。解决每次拿到shp先看属性表的更新年份和来源说明再对照当前行政区划代码表逐条核对。重点检查两类一是代码变了但名称相近二是老代码已消失、新代码需要查找合并关系。我一般保留两份文件一份旧版带old_name字段一份新版带当前名称这样历史分析和当前展示互不干扰。对时效性要求高的项目建议在读取流程里直接加一段代码校验# 假设code_latest.csv是最新行政区划代码表 latest pd.read_csv(code_latest.csv, dtype{adcode: str}) gdf[adcode] gdf[adcode].astype(str).str.zfill(6) diff set(gdf[adcode]) - set(latest[adcode]) print(需要核对的旧代码数量:, len(diff))5.4 多边形自相交与环方向异常几何修复的快速做法现象地图上某个区县填充色出现异常空洞或者打印时该区域有一块像被挖掉了一样还有的情况是面积计算明显小于常识值。原因多边形在数字化过程中产生了自相交环即边线自己穿过自己破坏了面的拓扑结构也有部分数据的内外环方向反转导致渲染时被当作空洞。解决shapely的make_valid是最直接的修复手段。前面3.1已经给过修复代码这里补充一点修复后要再跑一次面积校验如果一个县的面积和当地公开统计对不上多半还有未修复的拓扑错误。环方向问题在GeoPandas渲染中较少引发空洞但导出到CAD或打印出版时就会出现处理方式是用统一方向函数强制外环逆时针、内环顺时针不过绝大多数Web项目不需要这一步只做make_valid就够。6. 转存轻量GeoJSON并做一次自检把Shp用于Web端的一套顺手流程6.1 简化节点并导出GeoJSONGeoPandas把shp直接转成GeoJSON很简单但原始shp的坐标节点往往过于密集转出来的文件大、前端加载慢。我的做法是先简化再导出简化容差在经纬度坐标系下可以这样估算0.005度大约对应0.5公里展示市州轮廓足够如果细化到区县边界容差建议降到0.001。# 只保留必要字段再导出缩小文件体积 simple gdf[[adcode, name, geometry]].copy() # simplify用Douglas-Peucker算法preserve_topology避免相邻面出现缝隙 simple[geometry] simple.geometry.simplify(tolerance0.005, preserve_topologyTrue) simple.to_file(sichuan_city.geojson, driverGeoJSON)6.2 用随机采样点验证简化后的边界简化最常见的副作用是相邻多边形之间出现微小缝隙或重叠肉眼在小比例尺下看不出放大后却很显眼。转完GeoJSON后我用随机采样点做一次快速自检确保边界仍然闭合。import random from shapely.geometry import Point minx, miny, maxx, maxy simple.total_bounds for _ in range(10): pt Point(random.uniform(minx, maxx), random.uniform(miny, maxy)) inside simple.geometry.contains(pt).any() print(f({pt.x:.3f}, {pt.y:.3f}) inside{inside})如果外层省份没有出现大面积空洞各市州的随机点结果基本一致挑一两个落在市州交界线上的点尤其能考验简化质量。这个验证做完文件就可以放心交给前端了。我也曾经图方便直接用原始未简化的shp丢给前端结果5MB的GeoJSON让页面卡了好几秒后来才养成“先简化再验证再交付”的固定步骤。希望这套流程能帮到正在和四川shp打交道的你。本文还有配套的精品资源点击获取