简介：这份中国喀斯特岩溶空间分布矢量数据集面向GIS从业者、地质地理研究者及环境规划人员，用于获取全国岩溶地块的边界与岩性属性信息，支撑地貌分析、农业规划、水利工程与旅游开发等场景。资源包共8个文件，约1.2MB，以SHP矢量数据为核心，配套SHX、SBX、SBN空间索引文件，DBF属性表，PRJ空间参考，CPG编码说明及XML元数据，构成一套可直接在GIS软件中加载检索的完整数据。属性字段包含rock_type岩性分类（连续与不连续碳酸盐岩）、Shape_Area面积、Shape_Len周长及RTypeLabel文本标签，便于量化评估岩溶地貌空间范围与地表地下水系统规模。目前已有394人学习下载，适合需要中国岩溶分布基础底图的中高级GIS用户参考使用。
1. 喀斯特岩溶空间分布矢量数据集：一份能直接进 GIS 的 SHP 底图
做西南地区水文、地质灾害或生态评价的同行大概都遇到过同一个尴尬：手头有 DEM、有降雨栅格、有行政区划，唯独缺一层靠谱的岩溶分布边界。自己从地质图数字化，光是配准和判读就能耗掉一周，还未必对得上。这份中国喀斯特岩溶空间分布矢量数据集，就是把这个环节直接省掉——它以 SHP 格式提供岩溶空间分布范围，能直接拖进 ArcGIS、QGIS 或 PostGIS 参与叠加分析。适合做区域岩溶发育评价、地下水脆弱性制图、碳酸盐岩区工程选址的人。它解决的不是"有没有数据"，而是"这层数据能不能和你的其他图层对齐、能不能直接算面积、能不能导出成别的格式继续用"。
2. 先搞懂 SHP 这套文件到底装了什么：从坐标系到属性表
2.1 为什么岩溶研究绕不开矢量边界
岩溶区和非岩溶区的处理逻辑完全不同。碳酸盐岩分布区要考虑溶蚀、塌陷、地下管道流，非岩溶区按常规坡面流处理就行。所以第一步永远是"把研究区裁成岩溶和非岩溶两块"。栅格分类当然也能做，但岩溶边界本质是一条地质界线，用矢量表达更准，也方便后续按行政区、流域做统计。SHP 作为 Esri 定义的开放格式，几乎被所有 GIS 软件原生支持，这也是它比 GeoJSON、KML 更适合做底图的原因——你不用担心某个软件读不了。
这份数据集的核心价值在于"空间分布"四个字。它不是一张图片，而是带几何和属性的矢量图层，你可以查询、裁剪、合并、算面积，也能导出成 WKT、TXT 或转成 3dtiles 做三维展示。对做县域尺度评价的人来说，直接按行政区划边界做相交，就能得到每个县的岩溶面积占比。
2.2 SHP 文件组的构成与字段含义
一个完整的 SHP 从来不是单个文件。拿到压缩包解压后，常见的是这样一组：
| 文件后缀 | 作用 | 缺失后果 |
|---|---|---|
| .shp | 存储几何形状 | 没有图形，只剩空表 |
| .shx | 几何索引 | 部分软件无法定位要素 |
| .dbf | 属性表 | 丢失字段，无法分类统计 |
| .prj | 坐标系定义 | 坐标变成"无名"，无法叠加 |
| .cpg | 字符编码 | 中文属性可能乱码 |
我一般拿到数据先看有没有 .prj。没有的话，得问清楚是地理坐标（经纬度）还是投影坐标（米），否则一叠加就偏到姥姥家。属性表里通常会有岩溶类型或分区编码字段，具体字段名以实际数据为准，常见的是用整型编码区分裸露型、覆盖型、埋藏型。做统计前先打开属性表看一眼字段类型，字符串和数值混用是后面分组统计翻车的重灾区。
2.3 坐标系判断：一个能省半天返工的检查动作
很多人栽在坐标系上。判断方法很简单：在 QGIS 里加载后看右下角坐标值。如果是 73~135 这种量级，基本是经纬度（CGCS2000 或 WGS84）；如果是六七位数，那是投影坐标，单位是米。两种情况下做面积计算的姿势完全不同——经纬度直接算面积会得到"平方度"，毫无意义。
# 用 geopandas 快速检查坐标系和字段 import geopandas as gpd gdf = gpd.read_file("karst_distribution.shp") print("坐标系:", gdf.crs) # 看是 EPSG:4326 还是投影带 print("要素数量:", len(gdf)) print("字段列表:", list(gdf.columns)) print("几何类型:", gdf.geom_type.unique())这段代码干三件事：确认 CRS 是否已定义、看数据规模、看几何是面还是多面。如果gdf.crs返回 None，说明 .prj 丢了，得手动指定。geom_type如果是 MultiPolygon，说明有飞地或破碎图斑，后面做融合时要留意。参数上没什么可调的，纯粹是体检。跑完这一步，你对这份数据能不能直接用就有底了。
3. 把 SHP 接进实际工作流：裁剪、统计与格式转换
3.1 按研究区裁剪并统计岩溶面积占比
拿到全国或大区域数据，第一步几乎都是裁到自己的研究区。假设你手头有县域行政区划边界 shp，想算每个县的岩溶面积占比，思路是"先投影到等面积坐标系，再相交，再分组求和"。
import geopandas as gpd # 1. 读取岩溶数据和行政区划 karst = gpd.read_file("karst_distribution.shp") county = gpd.read_file("county_boundary.shp") # 2. 统一到等面积投影，避免面积失真（CGCS2000 高斯投影或 Albers） karst = karst.to_crs("EPSG:4527") # 按实际区域选带号 county = county.to_crs("EPSG:4527") # 3. 相交，得到每个县内的岩溶图斑 inter = gpd.overlay(karst, county, how="intersection") # 4. 计算岩溶面积并汇总到县 inter["karst_area"] = inter.geometry.area result = inter.groupby("county_name")["karst_area"].sum().reset_index() # 5. 除以县总面积得到占比 county["total_area"] = county.geometry.area result = result.merge( county[["county_name", "total_area"]], on="county_name" ) result["ratio"] = result["karst_area"] / result["total_area"] print(result.head())逻辑说明：to_crs是关键，不投影直接算面积会得到错误结果。overlay的how="intersection"只保留两者重叠部分。分组字段county_name要换成你数据里真实的县名字段。参数上，投影带号必须按研究区经度选，选错带号面积会有千分之几的系统偏差。这一步跑通，你就有了每个县的岩溶面积占比表，可以直接进后续评价模型。
3.2 导出为 WKT、TXT 与转 3dtiles 的取舍
有些场景不需要 SHP 本身，而是要它的坐标文本。比如把岩溶边界喂给某个只认 WKT 的数据库，或者导出成 TXT 做批量处理。geopandas 一行就能转：
# 导出为 WKT 文本 karst["wkt"] = karst.geometry.apply(lambda g: g.wkt) karst[["wkt"]].to_csv("karst_wkt.txt", index=False, header=False) # 导出为带经纬度的 TXT（取质心） cent = karst.geometry.centroid karst["lon"] = cent.x karst["lat"] = cent.y karst[["lon", "lat"]].to_csv("karst_centroid.txt", index=False)WKT 适合入库和程序解析，TXT 质心适合做点位标注。要注意的是，转 WKT 前坐标系必须是经纬度，否则 WKT 里的数值是投影米，别人拿到对不上。至于 shp 转 3dtiles，那是三维可视化的事，需要先把面拉伸成体，再用专门工具切片，属于进阶操作，普通二维评价用不上，别为了炫技绕远路。
3.3 和 DEM 提取的 shp 叠加时注意什么
热词里有人问"arcgis 从 dem 提取 shp"，这通常是提取等高线或流域边界。把岩溶分布和这类派生 shp 叠加时，最大的坑是精度不匹配。DEM 派生的边界往往有锯齿，岩溶边界是平滑地质线，两者相交会产生大量碎小图斑。我的习惯是叠加前先对 DEM 派生 shp 做一次平滑或简化，容差按 DEM 分辨率定，比如 30 米 DEM 用 15~30 米容差。否则后面统计会被碎图斑拖慢，结果也不干净。
4. 避坑与排查：SHP 用起来最容易翻车的五个地方
4.1 中文属性乱码
现象：打开属性表，原本的岩溶类型变成一堆问号或方块。原因：.dbf 的字符编码和软件默认编码不一致，老数据常见 GBK，新软件默认 UTF-8。解决：确认有没有 .cpg 文件，没有就手动建一个，内容写 UTF-8 或 GBK；或者在 QGIS 里用"图层编码"选项强制指定。我一般先用文本编辑器打开 .cpg 看一眼，比在软件里反复试快。
4.2 叠加后面积对不上
现象：同一块岩溶区，单独算面积和相交后汇总的面积差很多。原因：多半是坐标系不统一，一个经纬度一个投影，软件做了隐式转换但没提示。解决：叠加前强制把所有图层to_crs到同一个投影坐标系，别信软件的自动对齐。养成"叠加前先打印所有图层 CRS"的习惯，能省掉大量返工。
4.3 图斑破碎导致统计缓慢
现象：相交后要素数量暴涨，groupby 跑几分钟不出结果。原因：原始数据图斑细碎，加上边界不吻合，产生海量小多边形。解决：统计前先做一次融合（dissolve），按岩溶类型合并，或者用buffer(0)修复几何后简化。参数上，简化容差别超过最小关注图斑的尺度，否则会丢信息。
4.4 导出 WKT 后坐标是米不是度
现象：导出的 WKT 喂给别的系统，位置偏到国外。原因：导出时数据还在投影坐标系，WKT 里是米。解决：导出前to_crs("EPSG:4326")转回经纬度。这个坑我踩过不止一次，后来固定成流程：凡是导出文本坐标，先转 4326。
4.5 渔网分割后属性丢失
现象：用渔网分割 shp 做格网统计，结果每个格子没有岩溶类型字段。原因：渔网是新建图层，分割时没把原属性带过去。解决：用overlay而不是简单裁剪，或者分割后做空间连接把属性挂回来。渔网分割 shp 本身没问题，问题出在属性传递这一步，很多人以为几何切了属性自然跟着，其实不会。
5. 进阶：把岩溶分布做成可复用的空间权重与批量出图
走到这一步，数据本身已经能用了，但真正拉开效率差距的是"复用"。我一般会把岩溶分布转成一个空间权重矩阵，喂给后续的空间自相关或地理加权回归。做法是：先按研究单元（县、格网）算岩溶面积占比，再基于邻接关系构建权重。
import libpysal from libpysal.weights import Queen # 基于行政区划构建邻接权重 w = Queen.from_dataframe(county) w.transform = "r" # 行标准化 # 把岩溶占比作为变量，做空间自相关 from esda.moran import Moran moran = Moran(result["ratio"].values, w) print("Moran's I:", moran.I, "p值:", moran.p_sim)Queen是共边或共点都算邻接，transform="r"做行标准化让权重可比。Moran输出的是岩溶分布有没有空间聚集性。参数上，如果研究单元大小差异大，考虑用距离权重替代邻接权重。这一步把静态的 SHP 变成了能参与统计推断的变量，价值比单纯出图高一个量级。
批量出图也是同理。把裁剪、统计、出图写成函数，换一个县名就跑一遍，几十个县一晚上出完。我习惯把投影带号、字段名、输出路径都做成参数，避免硬编码。最后提醒一句：这份数据是空间分布底图，不是实时地质调查成果，用它做评价可以，但涉及具体工程选址，还是得叠加上实地勘察资料交叉验证。从那以后我每次拿到新的 SHP，都强制先跑一遍坐标系和字段体检，再动手分析。希望帮到你。
本文还有配套的精品资源，点击获取
