ARTICLE DETAIL

资讯详情

深耕网站视觉设计与运营推广的一线实战洞察。

洱海SHP文件从数据体检到出图:坐标系转换、KML与天地图叠加全流程

洱海SHP文件从数据体检到出图:坐标系转换、KML与天地图叠加全流程 简介一份面向GIS操作与空间分析的洱海地区矢量底图可直接作为ArcGIS、QGIS等软件的基础图层适用于环境规划、城乡规划、测绘及科研场景。压缩包共31个文件、约8.86MB涵盖shp、dbf、shx、prj等Shapefile标准组件并含adf、xml等栅格与元数据文件shp负责几何信息dbf存储属性数据prj定义坐标系统adf对应DEM高程数据能支撑制图、查询与空间分析。资源包含洱海湖泊边界、岸线、水系线、周边城镇等要素可直观展示湖区及邻近区域的地理格局。已有235人浏览学习。导入GIS后可进行湖泊面积与岸线长度统计、缓冲区分析、生态保护区划定、旅游及建设用地布局等操作为洱海保护与区域可持续发展提供可靠的数据底座。1. 洱海SHP文件到底能做什么先用一张底图搞定画图、测量和出图做洱海流域治理、环湖截污或水质监测布点第一步往往不是跑分析而是先有一张能用的底图。网上下到的洱海SHP文件版本很多有的带流域边界有的只有湖面轮廓有的坐标系齐全有的连 .prj 投影文件都没有。把一份只有几兆的矢量数据真正变成 GIS 操作底图要解决数据体检、投影转换、裁切合并、叠加在线地图、输出分享这一整条链路。这篇笔记就把链路上每个环节拆开讲透。有人把洱海SHP当成黑匣子打开就往图上叠结果叠到在线地图上偏移几百米还以为是数据不行其实多数时候问题出在坐标系文件缺失或投影选错。文章所有操作都兼顾 ArcGIS Pro 和 QGIS关键步骤给 Python 代码可以直接改路径跑通。适合刚入门 GIS、拿到 SHP 不知从哪下手的新手也适合做规划环评、想省掉重复修数据的熟手。2. 拿到洱海SHP以后先做体检四个附属文件、坐标系统和属性表缺一不可2.1 SHP 文件不只是一个文件.shp/.shx/.dbf/.prj 缺哪个都会翻车SHP 的正式名字叫 shapefile它不是单文件而是一组文件的集合。一个完整的洱海 SHP 至少要有四个配套文件.shp存几何图形.shx存图形索引.dbf存属性表.prj存坐标系信息。很多网盘里流传的“洱海SHP”压缩包点开以后只有 .shp 和 .dbf没有 .prj这种数据一旦拖进 ArcGIS Pro 或者 QGIS软件会猜一个坐标系猜对了算运气猜错了直接导致整个底图叠不到天地图上。拿到数据第一步不是急着画图是体检。用 Python 的 geopandas 读一次把最关键的几个信息打出来import geopandas as gpd shp_path rD:\gis_data\erhai\erhai.shp gdf gpd.read_file(shp_path) print(要素数量:, len(gdf)) print(几何类型:, gdf.geom_type.unique()) print(坐标系:, gdf.crs) print(字段列表:, gdf.columns.tolist()) print(空间范围:, gdf.total_bounds)这段代码做了四件事确认要素个数是否正常确认是面数据还是线数据确认坐标系有没有被正确读取最后把空间范围打出来。gdf.crs如果显示None说明数据缺 .prj 或 .prj 内容无法识别后面所有叠加、测距、转 KML 都会埋雷。如果确认 crs 是空的先别慌用手动指定坐标系的方式抢救。多数从公开渠道下载的洱海SHP是 WGS84 经纬度坐标也就是 EPSG:4326可以这样写from pyproj import CRS if gdf.crs is None: gdf.crs CRS.from_epsg(4326) gdf.to_file(rD:\gis_data\erhai\erhai_fixed.shp, encodingutf-8)注意to_file保存时指定encodingutf-8否则属性表里的中文地名洱海、挖色、双廊在 ArcGIS Pro 里打开可能显示成乱码。这个“乱码”问题看起来是显示问题实际上会连累字段筛选后面做标注时才发现就晚了。我一般会在体检脚本里把这一句也加上宁可多写一步也不在出图前浪费时间。2.2 先看坐标值再动数据地理坐标和投影坐标的判别打开洱海SHP之后先看图层的坐标范围再判断坐标系类型。洱海大概在东经 100 度、北纬 25 度附近。如果total_bounds打印出来是(100.0, 25.5, 100.4, 25.9)这种小数数据是经纬度的地理坐标系如果打印出来是(620000, 2800000, 640000, 2900000)这种六位七位的大数是投影坐标系单位是米。地理坐标系和投影坐标系不能混着用。直接用经纬度数据测面积算出来是“平方度”那是一个没有任何实际意义的数字。直接拿经纬度数据叠加天地图虽然能对上位置但做不了缓冲区、做不了渔网分割、算不了真实面积。所以体检完之后要立刻决定这个项目里数据在哪个坐标系下工作。洱海区域我一般推荐用 WGS84 / UTM zone 47NEPSG 编号是 32647。原因很简单UTM 47N 的中央经线是东经 99 度正好覆盖洱海它用米做单位做面积计算、生成渔网网格、做缓冲区分析都顺手和 WGS84 互相转换误差极小转回经纬度给在线底图用也不会产生额外漂移。可以用这样的代码快速确认和转换gdf_wgs84 gdf.to_crs(epsg4326) print(转成经纬度后范围:, gdf_wgs84.total_bounds) gdf_utm gdf.to_crs(epsg32647) print(转成UTM后范围:, gdf_utm.total_bounds) print(估算面积(平方公里):, round(gdf_utm.geometry.area.sum() / 1e6, 2))这里to_crs是投影转化epsg32647指定目标坐标系geometry.area返回的是平方度还是平方米取决于当前坐标系。只有转成 UTM 之后算的面积才是平方米除以 1e6 得到平方公里。这一步做完底图才具备“可测量”的基础条件。2.3 属性表藏着关键信息筛选、字段命名与面积单位的检查洱海SHP的属性表字段质量参差不齐。有的文件里只有一个NAME字段有的带TYPE、AREA、PERIMETER这类早期 ArcInfo 遗留字段。我见过一份洱海SHP里面同时有湖面、岛屿、环湖路三类要素全堆在一个图层里不筛选根本没法用。读取属性表内容前几行先看里面到底有什么print(gdf.head(3).T) if TYPE in gdf.columns: print(gdf[TYPE].value_counts())如果TYPE字段里有“湖面”“岛屿”等类别做底图时通常只保留湖面主体gdf_water gdf[gdf[TYPE].astype(str).str.contains(湖)] print(筛选后要素数:, len(gdf_water))筛选时用str.contains而不是因为不同版本的洱海SHP里“湖面”可能叫“湖泊水面”叫“水面”甚至叫“水体”用包含匹配不容易漏。属性表检查还有一个细节看AREA字段的数字量级。如果字段值只有零点几、几说明导出时用的是度面积字段没有参考价值后面出图时必须基于 UTM 重新计算别直接用旧字段。体检做完洱海SHP的基本原则就立住了坐标系明确、几何类型正常、属性能筛选。这时候才可以进入下一环节把它加工成真正能画能测的操作底图。3. 把洱海SHP变成能用的底图投影转换、局部裁切与KML互转3.1 GIS文件投影转化把经纬度转成米制坐标系的两种做法投影转化是 GIS 操作里最高频的动作没有之一。很多人拿到洱海SHP直接开始画图画完要量面积才发现单位不对这就是没在项目开始前把坐标系理清楚。常见做法是两种在 ArcGIS Pro 里用工具在代码里用to_crs。ArcGIS Pro 的路径是分析工具 - 数据管理工具 - 投影和变换 - 投影输入洱海图层输出坐标系选WGS 1984 UTM Zone 47N点击运行。这个工具界面上能实时看到处理范围适合单次操作。代码方式适合批处理和工程化流程前面已经写了基本用法这里补一个完整保存版本import geopandas as gpd gdf gpd.read_file(rD:\gis_data\erhai\erhai.shp) if gdf.crs is None: gdf gdf.set_crs(epsg4326) gdf_utm gdf.to_crs(epsg32647) gdf_utm.to_file(rD:\gis_data\erhai\erhai_utm.shp, encodingutf-8)set_crs和to_crs是两回事。set_crs是“我不知道坐标系我告诉你是哪个”它只改写元数据不改几何坐标to_crs是“我知道坐标系帮我转换到另一个坐标系”它会真正重新计算几何坐标。新手经常把两者混用导致投影转化后位置飞到了海里这是最典型的一个翻车点。投影转化时注意一个小习惯保存成新文件不要覆盖原文件。洱海SHP原始版本可能是某个课题组辛苦整理过的保留一份原版就是后悔药。我自己的项目目录里永远分成raw和processed两个子文件夹raw 绝不动。3.2 如何在SHP图中去掉一部分矢量裁剪边界与过滤要素“如何在shp图中去掉一部分矢量”这个需求在洱海项目里很常见。比如只需要洱海南岸做截污规划或者要把湖心小岛从面数据里抠掉。拆开看有两类操作按空间范围裁剪和按属性过滤。按空间裁剪用 geopandas 的clip。这里有个前提裁剪框和待裁剪数据必须在同一个坐标系下否则边界对不齐。先建一个比洱海外包围盒大 500 米的矩形框再裁from shapely.geometry import box # 在UTM坐标系下操作 minx, miny, maxx, maxy gdf_utm.total_bounds crop_box box(minx - 500, miny - 500, maxx 500, maxy 500) clipped gpd.clip(gdf_utm, crop_box) print(裁剪后要素数:, len(clipped))clip是按图形的空间位置做几何叠加输出保留的是原始要素落在裁剪框内的部分。如果裁剪框画小了没有被框住的部分会被整个丢掉所以留 500 米缓冲是常规操作避免洱海边线因为拓扑误差被切出毛边。按属性过滤则简单得多适合“只要湖面不要岛屿”这种需求gdf_land gdf_utm[~gdf_utm[TYPE].astype(str).str.contains(岛)]这里~是取反筛选出不包含“岛”字的要素。两个操作都完成后建议做一次几何检查用gdf.is_valid.sum()看看有没有拓扑错误数量不对先修复再导出。3.3 ArcGIS SHP转KML的三个步骤以及KML转SHP的回来路ArcGIS SHP转KML是热搜词里出现频率很高的操作因为 KML 可以直接丢进 Google Earth、奥维地图也能作为天地图在线标注的交换格式。转出 KML 有三个步骤少一步都会出问题。第一步确保数据在 WGS84 经纬度坐标系。KML 标准强制使用 WGS84如果洱海SHP还在 UTM 坐标系必须先转回 4326gdf_4326 gdf_utm.to_crs(epsg4326)第二步属性表里的中文要干净。KML 里的字段名和值会直接出现在地图软件的信息气泡里字段名最好改成英文或拼音避免一些国产地图 App 读取乱码。第三步指定 driver 导出gdf_4326.to_file(rD:\gis_data\erhai\erhai.kml, driverKML)导出完成后用 Google Earth 打开看一次位置。我曾经导出一份洱海SHP在 ArcGIS Pro 里看完全正常导到 KML 后在 Google Earth 里整体偏移了约 60 米。排查结果是原始数据虽然显示为 WGS84实际是 CGCS2000 坐标系两个坐标系在大地基准上有几十厘米到几十米的差异投影转化时没有处理基准转换偏移就顺着链路带到了 KML。KML转SHP是反向操作。在 ArcGIS Pro 里用转换工具 - KML转图层或者用 Python 直接读再写kml_gdf gpd.read_file(rD:\gis_data\erhai\erhai.kml) kml_gdf.to_file(rD:\gis_data\erhai\erhai_from_kml.shp, encodingutf-8)这里要提醒一个坑gpd.read_file读 KML 时Name字段会被读进来但 KML 里的样式颜色、线宽不会保留。做线划数据转换时样式信息丢失是常态别在转换后去找原来的配色。KML 更适合做位置交换不适合做数据野外工作的最终存储格式。4. 做一张能交付的操作底图合并、面积字段与天地图叠加出图4.1 ArcGIS Pro 里做SHP合并Merge 工具和字段统一多个洱海SHP文件需要合成一个图层时比如手上有洱海北岸和南岸两个分块文件要合并成完整的湖面用 ArcGIS Pro 的合并Merge工具分析工具 - 数据管理工具 - 合并添加所有输入图层输出要素类选一个目标位置。这个工具不挑字段结构但字段结构越一致合并后越省事。如果两个文件的字段一个是NAME一个是name合并后的属性表会出现两列看起来像重复数据。所以合并之前先统一字段名我用 Python 处理import pandas as pd gdf1 gpd.read_file(rD:\gis_data\erhai\north.shp) gdf2 gpd.read_file(rD:\gis_data\erhai\south.shp) for g in (gdf1, gdf2): g.columns [c.lower() for c in g.columns] merged gpd.GeoDataFrame( pd.concat([gdf1, gdf2], ignore_indexTrue), crsgdf1.crs ) print(len(merged))pd.concat把两个 GeoDataFrame 直接拼起来ignore_indexTrue避免索引打架。合并前一定要确认两边的坐标系一致一个 UTM 一个经纬度直接合并产生的图形会变成一个在洱海、一个在洱海西边几百公里外。crsgdf1.crs是给合并结果指定坐标系如果两个输入 crs 不同程序不会报错但这个隐患会在你量面积时炸开。4.2 高精度面积计算从度到米的两次转换保留两位小数面积计算是洱海底图使用频率最高的功能。上报材料里动辄“洱海水域面积 256.xx 平方公里”这个数字精度依赖坐标系的正确性。记住一个铁律先在投影坐标系下计算再保留两位小数。gdf_utm[area_km2] (gdf_utm.geometry.area / 1e6).round(2) print(gdf_utm[[NAME, area_km2]].head(10))这段代码给属性表新增了一个area_km2字段单位是平方公里保留两位小数。用round(2)而不是在 Excel 里设置单元格格式是因为属性表里存的是实数Excel 的显示格式不影响实际数值精度导出后该是多少还是多少。这类面积字段在填环保统计表、永久基本农田面积统计等场景里同样适用两位数小数的规范在这里统一处理不要后面手工再凑。ArcGIS Pro 里的等价操作是右键图层 -属性表- 新建字段 - 字段计算器选择面积几何单位选平方千米。注意字段计算器算面积时同样要求当前数据是投影坐标系否则工具会灰掉不让你选“面积几何”。4.3 GIS导入天地图底图解决在线地图加载不了的落地步骤洱海SHP做操作底图最理想的搭档是一张在线天地图影像或路网本地矢量负责边界和属性影像负责空间参照。ArcGIS Pro 里导入天地图有几个版本差异加载不了是热搜词里反复出现的问题十有八九是下面三个原因。第一个原因是没有配置天地图 key。天地图的在线服务需要申请 token拿到后要写进服务地址地址格式形如http://...tk你的key。在 ArcGIS Pro 里选择地图 - 底图 - 添加底图 - 从路径添加数据把带 token 的服务地址粘进去。第二个原因是服务地址用了 http 而项目设置强制 httpsPro 默认会拦截混合内容报错信息常常是“无法连接”。把地址栏里的 http 改成 https或者调整项目的安全选项能解决一大半问题。第三个原因在 QGIS 上反而少见QGIS 的XYZ Tiles加载天地图比 ArcGIS Pro 宽容对服务地址格式错误有更明确的报错提示。如果在线底图始终加载不出来还有个离线替代方案用之前导出的洱海SHP配上一张下载好的影像瓦片在 ArcGIS Pro 里用地理配准工具把影像对齐到洱海SHP边界上。这个操作不需要连外网适合内网环境的项目组。4.4 图例标签换行与出图版式让底图自解释操作底图最终要出成图片发给别人图例标签的排版直接影响可信度。ArcGIS Pro 默认图例标签是单行长字段值会顶出图框。图例标签如何换行这个需求Pro 里的做法是用字段表达式改写标注文本。# 示例标注显示为两行 # 第一行: 洱海 # 第二行: 256.32 km²在标注表达式中写成NAME vbNewLine 面积: Round([area_km2], 2) km²QGIS 的表达式写法稍微不同用的是字符串拼接NAME || \n || 面积: || round(area_km2, 2) || km²vbNewLine在 ArcGIS 里表示换行符QGIS 里用\n。两种软件对空格的容忍度不一样ArcGIS 的拼接会自动忽略表达式里的多余空格QGIS 则保持原样。标签换行之后再设置字体大小和背景色确保叠加在天地图上时文字可读。到这里一张自解释的操作底图才算真正做完。5. 洱海SHP的常见翻车现场与排查思路5.1 发给别人说打不开SHP文件怎么保存发送才对现象把“洱海.shp”从文件夹里拖出来只用邮件发给合作方对方说打开图层全是空的甚至直接报错无法打开。原因shapefile 是多文件格式只发送 .shp 文件等于发了半个图层。几何数据在 .shp 里索引在 .shx 里属性在 .dbf 里坐标系在 .prj 里缺了任何一个接收方都可能打不开或打开后缺属性、错坐标。解决在 ArcGIS Pro 里右键图层 -共享-打包图层生成一个.lpkx图层包直接发或者把整个文件夹压缩成 zip 再发。我用 Python 也写过一个小函数自动检查并打包zip -r erhai_shp.zip erhai.shp erhai.shx erhai.dbf erhai.prj最简单的原则SHP 的六个配套文件含 .shp、.shx、.dbf、.prj、.cpg、.sbn能一起发就连同 .cpg 字符集定义一起发中文属性能不能正常显示往往就差这个 .cpg 文件。5.2 叠在线底图整片偏移坐标系基准的“黑匣子”现象洱海SHP叠加天地图影像后边界整体向东或者向南偏移偏移量从几十米到几百米不等放大后看得清清楚楚。原因很多公开SHP的坐标基准并不是 WGS84而是 CGCS2000 或者早期北京54、西安80。这三个基准之间虽然有固定转换关系但转换参数因地区而异不做严格的七参数转换就会偏移。解决先用“添加 XY 坐标”或者读取 .prj 文本确认原始基准如果是 CGCS2000在 ArcGIS Pro 里用投影工具把“地理坐标变换”参数从默认改成CGCS2000_To_WGS_1984。这是坐标系问题的黑匣子阶段我不会在没有查清原数据的基础上硬套参数查清基准再动手花 10 分钟省下后面 2 个小时的对图时间。5.3 面积算出来是天文数字先投影再计算现象用 ArcGIS Pro 的字段计算器选了面积洱海面积算出 100 多亿单位看起来像万平方公里。原因图层还是经纬度坐标系几何面积的计算单位是平方度不是平方米。在经纬度下算面积本来就是无效操作工具没有拦截结果就失控。解决回到第 3 章先把洱海SHP转成 UTM 47N再算面积。我每次出面积数据前都会自查一步gdf_utm.geometry.area.sum()的结果如果小于几个亿基本不对洱海湖面实际大约 250 平方公里量级算出结果明显不在这个范围就是坐标系没转。5.4 转KML后位置不对导出前必须做的一步现象SHP 转 KML 后在 Google Earth 里打开湖面位置跟卫星影像对不上偏移方向还随缩放级别变化。原因KML 标准要求 WGS84 坐标系但原始SHP如果是投影坐标系加 CGCS2000 基准ArcGIS 的“KML转出”工具按默认 WGS84 处理基准不管底层基准差异直接重投影输出位置上就会带固定偏移。解决导出 KML 前显式执行一次to_crs(epsg4326)并且把转换参数里的基准变换设置好。代码里写清楚转换目标比让工具自己去猜少踩一半坑。转完先别急着交付用 Google Earth 闪烁检查一次确认边界压线了再发。5.5 在线地图瓦片刷不出来天地图token和服务地址检查清单现象洱海SHP底图上叠天地图灰屏、转圈、白底图层列表里天地图显示已加载但界面空白。原因天地图服务地址需要 token服务地址的协议是 http且部分网络环境下访问天地图域名本身就慢这三个因素叠加导致瓦片全部加载失败。解决按顺序排查先在浏览器里直接访问带 token 的天地图服务地址确认能出图再到 ArcGIS Pro 里检查数据源看服务类型选的是WMTS还是XYZ Tiles最后把项目设置里的网络代理关掉再试一次。这个排查思路同样适用于其他在线底图不只是天地图。在项目最忙的时候我会直接用 QGIS 临时顶替QGIS 的瓦片缓存机制对在线服务更宽容先把图出了再说。6. 进阶玩法用渔网分割洱海SHP做成监测网格与3D Tiles6.1 用渔网工具生成规则的监测网格单元洱海水质监测有个常见需求把湖面切成规则网格每个网格代表一个监测单元。ArcGIS Pro 的创建渔网工具在分析工具 - 要素类 - 创建渔网下输入洱海SHP的范围指定像元宽度和高度比如 1 千米乘 1 千米生成后与洱海面做相交就得到湖边界的网格集合。用 Python 做同样的事好处是可以把网格边长做成参数反复调整时不用点鼠标from shapely.geometry import box cell_size 1000 # 米 minx, miny, maxx, maxy gdf_utm.total_bounds grid_cells [] x minx - (minx % cell_size) while x maxx: y miny - (miny % cell_size) while y maxy: grid_cells.append(box(x, y, x cell_size, y cell_size)) y cell_size x cell_size grid gpd.GeoDataFrame(grid_cells, columns[grid_geometry], crsEPSG:32647) intersected gpd.overlay(grid, gdf_utm, howintersection) print(有效监测网格数:, len(intersected))gpd.overlay的howintersection返回网格和湖面相交的部分只有落在水里的网格会被保留。后面给每个网格加编号再关联监测指标数据就能画出一张水质空间分布图。网格大小按监测断面密度来定1 公里网格适合小范围精细监测5 公里网格适合整个洱海流域的大尺度评估这个参数没有绝对标准但网格越小统计噪声越大生成的文件也越大。6.2 从SHP到3D Tiles一场从桌面GIS到WebGIS的接力洱海SHP还可以继续往前走转成 3D Tiles 后在浏览器里加载做成流域管理平台的数字底座。标准路线是SHP 转 GeoJSONGeoJSON 再转 3D Tiles。SHP 转 GeoJSON 用gdf.to_file(..., driverGeoJSON)一行搞定GeoJSON 转 3D Tiles 可以用常见的切片工具完成不同工具的输出参数差异较大这里不做展开只强调两个关键点。第一个转 3D Tiles 前必须把坐标系转成 EPSG:4978地心坐标系否则切片后的模型位置会脱离地球表面。这个转换在部分转换工具里是隐式的但结果飘在天上时多半是这一步没做对。第二个SHP 只有边界线和面没有高程转 3D Tiles 后是贴在表面上的“压膜数据”如果平台需要表现湖底地形得先叠加 DEM 做拉伸这一步要在 GIS 里完成而不是在切片阶段。从我自己的项目习惯来说现在拿到任何一份洱海SHP都会先按第 2 章跑一遍体检脚本确认坐标系、属性和几何没问题再谈投影、转KML还是做网格。这个习惯帮我挡掉了绝大多数“数据是坏的”判断很多所谓坏数据其实只是缺一个 .prj 文件。希望这份流程对你也有用。本文还有配套的精品资源点击获取
返回列表