ARTICLE DETAIL

资讯详情

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

鄱阳湖水系shp数据整理实战:从文件补全到水文分析

鄱阳湖水系shp数据整理实战:从文件补全到水文分析 简介面向ArcGIS初学者与流域研究人员提供长江流域鄱阳湖水系完整地形数据包。包内含湖区、河流等矢量图层及dem90m数字高程、山体阴影等栅格地形文件并配好了鄱阳湖流域.mxd工程文档双击打开后重新链接各图层到解压目录即可一键出图即使完全不熟悉GIS出图文件夹中也已备好jpg、PDF、EPS三种格式成品图可直接用于论文插图、汇报展示只是无法再调整制图区域。全包共71个文件以shp、dbf、prj、shx等矢量组件、adf栅格数据及xml元数据为主压缩包大小约62.91MB矢量栅格分层存放、目录结构直观。已有2854人学习下载该资源。数据来自网络搜集并进一步加工仅供学习科研使用准确性需自审如有侵权可联系删除。1. 拿到鄱阳湖水系shp却打不开地形图矢量数据缺了什么做水利、规划或环评的人多半都经历过这个场景从某个数据共享平台下载到“鄱阳湖水系流域地形图.shp”拖进ArcGIS一看要么只显示一条湖岸线要么图层顶着“未知坐标系”在地图角落飘着要么湖泊边界和影像底图错开了几十公里。这份shp不是单一文件它通常是一组由等高线、高程点、水系、湖面和流域边界组成的矢量包真正能在ArcGIS里做流域分析、裁剪出图和面积统计至少还差文件补全、坐标统一和几何整理三步。本文按一线做法把这套流程走一遍适合刚接触shp的新手也适合被坐标问题反复绊倒的老用户。2. 拆解鄱阳湖水系地形图shp图层构成、文件结构与坐标系落差2.1 一套地形shp实际包含的图层等高线、高程点、水系与边界标题里的“地形图”三个字在ArcGIS语境里通常不是单图层栅格而是一套矢量成果。常见做法是等高线polyline、高程点point、河流与湖泊polyline/polygon、流域边界polygon有的还带堤防、蓄滞洪区和水利工程点。不同单位出的数据图层命名五花八门“contour”“等高线”“DLTB”“DEM_Line”都有属性字段更是各写各的。图层几何类型常见字段用途等高线线ELEV / 高程 / Contour生成DEM、地形分析高程点点Height / GC插值补充、制图标注河流水系线Name / 名称河网提取、缓冲分析湖泊水库面Name / Area水域统计、淹没分析流域边界面Code / 流域名裁剪、分区统计新手拿到数据包不要急着双击加载先看文件总大小和文件个数。一个正经的shp数据包至少应该包含主文件名相同的多个后缀文件如果解压出来只有一个“.shp”后面十有八九要翻车。属性字段的差异比图层名更让人头疼。有的数据用“DLDM”表示地类代码有的用“GB”表示国标码还有直接用拼音缩写。打开属性表看一眼字段前缀规律基本能猜出数据是按什么标准生产的猜不出来就把字段英文名当作黑匣子先别删保留原始字段是最稳妥的做法。2.2 shp文件不是一个文件.shp、.shx、.dbf、.prj的身份关系shp格式由ESRI在几十年前定义设计上就是一个“文件集合”而不是单个文件。ArcGIS读取shp时会同时查找几何、索引、属性和坐标信息缺了任何一个轻则属性为空重则整个图层打不开。文件后缀作用丢失后果.shp空间几何信息线条和面的坐标文件不存在则完全打不开.shx几何索引加速读取部分软件打不开可能只显示外框.dbf属性表字段和值都在这图层能加载属性全空.prj坐标系描述文本格式显示为未知坐标系容易错位.cpg属性表字符编码声明中文属性大概率乱码快速体检一个shp包我一般做两步。第一步是在文件夹里把扩展名显示打开确认上述文件是否齐备第二步是用记事本打开.prj文件看看里面是GCS还是Projected。prj是纯文本内容如果包含“GCS_WGS_1984”“GCS_CN_2000”那这是地理坐标系如果包含“Gauss_Kruger”“3_Degree_GK”“CM_117E”这类字样那已经是投影坐标系后面处理方向会完全不同。文件名最好不要动。shp的辅助文件是按主文件名匹配的有人喜欢把“鄱阳湖水系地形图.shp”改成“poyang.shp”结果只改了主文件忘了同步改.shx和.dbf打开时ArcGIS直接报错。要改名就在文件管理器里全选这组文件一起改或者拷进地理数据库后再改。2.3 坐标系落差CGCS2000、西安80、WGS84与投影坐标系鄱阳湖流域的shp数据坐标系混乱几乎是常态。同一个数据包里有的图层是CGCS2000地理坐标有的是WGS84还有老资料用西安80更麻烦的是有人把经纬度数据直接当平面坐标用加载到ArcGIS里所有要素都挤在一个像素点上。坐标系名称类型常见出现场景GCS_WGS_1984地理坐标全球公开数据、在线底图默认GCS_CN_2000地理坐标新版基础测绘成果Xian_1980_3_Degree_GK_CM_117E投影坐标老版本地形图成果CGCS2000_3_Degree_GK_CM_117E投影坐标近年自然资源数据标准成果查看某个shp的坐标系在ArcMap里右键图层属性切到“源”选项卡就能看到在ArcGIS Pro里则是在图层属性里找“源信息”。如果显示“未知”不要急着用“投影”工具去修——投影工具要求输入数据本身已经带坐标系未知坐标系的正确工具是“定义投影”。这里有个高频误用用户拿着一个已经漂移的shp用了“定义投影”随手选一个坐标系结果数据没修正反而被套上错误坐标。定义投影只是给数据“补身份证”不是“校正位置”。不知道原始坐标系时先看数据来源说明再看同数据包其他图层的.prj最后才考虑用同名影像或底图辅助判断。3. 用ArcGIS整理鄱阳湖水系shp投影统一、裁剪脚本与字段清洗3.1 先定底图再做投影中央经线选117°E还是114°E开始处理前先决定整批数据的“最终参照系”。鄱阳湖主体在东经115°到117°之间跨了两个3度带做全流域制图时常见做法是统一到CGCS2000_3_Degree_GK_CM_117E部分靠近修水上游的区域离中央经线远变形会稍大但不影响一般制图与拓扑分析。如果目标是做面积统计建议改用等积投影。自然资源数据里常用Albers等积圆锥投影ArcGIS的投影工具里可以直接选“China_Albers_Equal_Area_Conic”。平面面积的统计结果在Web墨卡托底下算出来会偏大在等积投影下才接近真实值。使用场景推荐坐标系理由出图、叠加、目视检查CGCS2000 3 Degree GK CM 117E与常见测绘成果一致面积精确统计Albers等积投影面积变形小发布到Web端WGS84地理坐标在线底图同框架统一步骤在ArcToolbox里是“数据管理工具→投影和变换→要素→投影”一次只能处理一个图层。图层多的时候用批处理脚本更稳见3.2。3.2 用arcpy脚本批量投影与裁剪几十个图层一次跑完手动对几十个图层逐个做投影和裁剪容易漏也容易在中间环节选错坐标系。我一般用一个临时目录放原始数据再用脚本统一处理输出到独立的地理数据库。地理数据库的好处是文件结构完整不再有shp辅助文件缺失的问题。# 批量投影将工作空间内所有shp统一到CGCS2000 3度带117E import arcpy arcpy.env.workspace rD:\poyang\raw_shp # 原始shp所在文件夹 arcpy.env.overwriteOutput True out_gdb rD:\poyang\poyang.gdb # 输出地理数据库 target_sr arcpy.SpatialReference(CGCS2000 3 Degree GK CM 117E) for fc in arcpy.ListFeatureClasses(): out_fc out_gdb \\proj_ fc[:-4] # 去掉.shp后缀作为输出名 print(开始处理, fc) try: arcpy.Project_management(in_datasetfc, out_datasetout_fc, out_coor_systemtarget_sr) print(完成, fc) except Exception as e: print(失败, fc, e)这段脚本先通过arcpy.ListFeatureClasses()枚举工作空间内所有要素类再逐个调用Project_management。arcpy.SpatialReference里写坐标系名称要和ArcGIS安装的坐标系定义一致如果写成“CGCS2000 3 Degree GK CM 117E”找不到可以换成EPSG编号。输出时用fc[:-4]把.shp后缀去掉避免输出要素类名称带点号导致后续工具报错。投影完成后再按流域边界裁剪。这里要求边界必须是面要素如果拿到的是流域边界线先用“要素转面”工具处理。# 按鄱阳湖流域边界批量裁剪 import arcpy arcpy.env.workspace rD:\poyang\poyang.gdb arcpy.env.overwriteOutput True boundary rD:\poyang\poyang_boundary.shp # 流域边界面 for fc in arcpy.ListFeatureClasses(proj_*): desc arcpy.Describe(fc) out_fc clip_ desc.baseName try: arcpy.Clip_analysis(in_featuresfc, clip_featuresboundary, out_feature_classout_fc) count arcpy.GetCount_management(out_fc) print(desc.baseName, 裁剪完成要素数:, count[0]) except Exception as e: print(desc.baseName, 裁剪失败:, e)Clip_analysis会把输入要素严格限制在边界范围内边界以外的要素直接丢弃。裁剪后每个图层都用一个独立要素类后续做水文分析或出图时不必反复屏蔽范围外数据。GetCount_management返回的是结果对象用count[0]取出计数字符串方便在日志里确认是否出现“裁剪后要素数为0”的异常。3.3 字段清洗中文乱码、字段类型和空值处理坐标系和几何问题解决后属性表往往还有一堆坑。最常见的是中文乱码。老数据的.dbf文件用GBK或GB2312编码而ArcGIS Pro默认按UTF-8读取打开属性表就是“锟斤拷”三件套。这时候看数据包里的.cpg文件里面写着“UTF-8”就改成“GBK”写着“GBK”就改成“UTF-8”改完重启ArcGIS再看属性表。处理字段用ListFields先看全貌再决定是改名还是改类型import arcpy arcpy.env.workspace rD:\poyang\poyang.gdb fields arcpy.ListFields(clip_rivers) for field in fields: print(字段名:, field.name, | 类型:, field.type, | 长度:, field.length)输出结果里如果出现字段名为空或长度异常的项多半是dbf文件损坏或字段定义不规范。字段改名的常用工具是AlterField改类型也可以用同一工具但字段类型改动会丢原有内容尤其是浮点转整型时小数部分直接截断。更稳妥的做法是新建一个正确类型的字段用字段计算器把旧值拷过去核对完再删旧字段。这套做法麻烦但不会把原始数据改坏。字段值也不全是可靠的。水系名称字段里常见前后空格、全角括号、编码残留批量清理时用字段计算器或strip()函数跑一遍。属性清洗不是技术门槛但它决定后续出图的图例和注记是否还要返工值得在裁剪后立刻做掉。4. 从地形图shp到流域分析成果DEM生成、河网提取与分水岭划分4.1 等高线转DEM还是直接用已有DEM为什么推荐TopoToRaster很多网络下载的“流域地形图shp”里只有等高线没有现成DEM。这时候就需要用等高线生成栅格高程模型。ArcGIS里常用两种方式TIN转栅格或者TopoToRaster。TIN方式生成的DEM在山坡陡峭处容易出现三角面棱角水文分析时流向计算会沿着三角面边界跑出“假河道”。TopoToRaster是专为等高线和水系设计的插值方法它同时考虑等高线、河流和边界约束生成的表面在水文分析里表现更平滑也更贴合真实地形。做流域分析时我优先选TopoToRaster。# 用等高线生成30米DEM import arcpy arcpy.env.workspace rD:\poyang\poyang.gdb arcpy.env.extent MAXOF # 用等高线范围定边界 arcpy.env.cellSize 30 # 像元大小30米 contours clip_contours # 等高线要素类 dem arcpy.ddd.TopoToRaster( [[contours, ELEV]], # [要素类, 高程字段] cell_size30, extentarcpy.env.extent, boundaryclip_watershed # 流域边界约束必填 ) dem.save(poyang_dem_30m)TopoToRaster的输入是一个列表的列表每个内层列表写“要素类 高程字段”字段必须是数值型且单位为米。boundary传入流域边界能防止插值在边界外飞出去。像元大小不是越小越好1:5万地形图等高距通常是10到20米30米像元足够匹配原图精度强行设置成5米生成的DEM只是一个更光滑的曲面不会增加真实信息计算时间却翻几倍。4.2 水文分析流程填洼、流向、流量、河网提取参数怎么设有了DEM就可以做水文分析。标准流程是填洼、流向、流量累积、提取河网。鄱阳湖周边平原区地形平坦DEM里的伪洼地比山区更多填洼这一步不能跳过否则流向计算会产生大量不连通的断流。# 鄱阳湖DEM水文分析填洼到河网提取 import arcpy from arcpy.sa import * arcpy.env.workspace rD:\poyang\poyang.gdb arcpy.env.snapRaster poyang_dem_30m # 对齐到DEM像元 arcpy.env.overwriteOutput True dem Raster(poyang_dem_30m) # 1. 填洼消除DEM中的伪凹陷 fill Fill(dem) fill.save(fill_30m) # 2. 流向D8算法每个像元指向最陡下坡方向 fd FlowDirection(fill) fd.save(fd_30m) # 3. 流量累积统计每个像元上游汇水面积 fac FlowAccumulation(fd) fac.save(fac_30m) # 4. 提取河网流量大于阈值的像元记为河道 threshold 3000 stream Con(fac, 1, None, VALUE {}.format(threshold)) stream.save(stream_30m) # 5. 河网栅格转矢量shp stream_shp StreamToFeature(stream, fd, stream_net_30m)Fill默认填平所有洼地但鄱阳湖区真实的湖泊洼地也会被填掉。如果需要保留湖盆可以用Fill(dem, z_limit)加填洼深度限制比如只填小于5米的洼坑避免把主湖盆填成平地。FlowDirection采用D8算法把每个像元的流向分到8个相邻像元中最陡的方向。FlowAccumulation统计上游汇水像元数量数量越大代表越可能是河道。阈值3000不是定死的。30米分辨率下3000个像元大约对应2.7平方公里汇水面积平原区取小一点比如1000到2000能捕捉更多小水沟但也会产生大量平行短流山区取5000以上更干净。第一次跑建议用3000然后用符号系统按流量分级看河网形态再调阈值。“阈值这东西一半是公式一半是玄学”多跑两次就有手感。4.3 子流域划分从出水口到分水岭面shp河网提取完成后如果还想划分鄱阳湖各支流的子流域就要用出水口点和流向栅格做分水岭分析。出水口的选取直接影响子流域形状和面积一般选水文站、桥涵位置或干流汇合点不要随手在图上点一个地方。# 从出水口到子流域边界shp import arcpy from arcpy.sa import * arcpy.env.workspace rD:\poyang\poyang.gdb arcpy.env.snapRaster fill_30m # 与填洼栅格对齐 fd Raster(fd_30m) fac Raster(fac_30m) outlet_pt outlet_stations # 出水口点要素类 # 将出水口吸附到流量最大的像元 snap SnapPourPoint(outlet_pt, fac, 100, OUTLET) snap.save(snap_outlet) # 基于流向栅格划分子流域 ws Watershed(fd, snap) ws.save(sub_watershed) # 栅格转面shp arcpy.RasterToPolygon_conversion( sub_watershed, sub_basins_shp, NO_SIMPLIFY )SnapPourPoint的作用是把出水口位置吸附到流量累积值最大的相邻像元上第三个参数100是搜索半径单位是像元。如果出水口点离河道只有几十米但没落在高流量像元上不吸附的话分水岭会从错误位置开始边界偏得离谱。Watershed工具用流向栅格从出水口向上游找分水岭输出栅格的每个像元都带一个子流域编号最后转成面shp做后续统计。这里有一个细节NO_SIMPLIFY参数保留栅格边界的原始锯齿分析用途保留锯齿没问题如果出图要美观再用“简化面”工具做一次平滑。先保真再求美这是分析数据的底线。5. 鄱阳湖水系shp使用避坑指南文件打不开、坐标漂移、边界对不上的排查5.1 只有shp文件打开报错或只显示外框现象双击.shpArcCatalog里能预览但拖进ArcMap后提示“无法打开要素类”或者图层只显示一个矩形框属性表为空。原因shp是文件集合缺少.shx或.dbf会导致几何或属性读取失败。下载数据时只选了“shp”这一项或者传输过程中辅助文件被误删。解决回到源数据包补齐.shx、.dbf、.prj注意主文件名必须完全一致。如果源数据已经找不到在ArcCatalog里对该shp右键选“导出”→“要素类到地理数据库”ArcGIS有时能从残缺文件中抢救出几何信息但属性字段不保证完整。以后传输shp请打压缩包别只发单个文件。5.2 图层挤在地图左下角或放大后什么都看不见现象添加shp后全图范围跑到负坐标区域要素缩成一个点或一条短线。原因数据本身是经纬度坐标但缺少.prjArcGIS不知道这是经纬度把它当作平面米制坐标直接显示。经纬度的单位是度一个度在屏幕上被当成一个米画面自然缩没了。解决先确认数据真的是经纬度值。查看属性表里坐标字段如果X在115到117左右、Y在28到30左右这是经度纬度。这时用“定义投影”工具选择GCS_CN_2000或GCS_WGS_1984再按第3章的流程做投影。注意这里只能选地理坐标系不能直接选投影坐标系。5.3 湖泊边界与影像底图差出几百米甚至几公里现象叠加在线影像时shp里的湖岸线和影像上的湖岸线平行错开有的地方差几百米有的地方差几公里。原因数据集内部坐标系不一致或者有人把WGS84坐标的数据“定义”成了CGCS2000两个椭球差异导致整体偏移。另一个常见原因是3度带选错鄱阳湖区域用了114E的成果去对接117E的底图。解决把所有图层统一到同一个坐标系再叠加。检查每个图层的.prj凡是地理坐标系不同的先“投影”到统一地理坐标凡是投影带不同的先“投影”到同一中央经线。WGS84和CGCS2000在实际使用中差得很小肉眼基本看不出来如果差到公里级大概率是投影带错而不是椭球错。5.4 裁剪后“要素数为0”或湖泊面积明显偏小现象用流域边界裁剪湖泊shp结果图层里什么也没有或者裁剪后湖泊面积比原数据小了一大块。原因裁剪边界是线而不是面Clip_analysis不接受线要素还有一种情况是边界与目标图层坐标系不一致空间关系判断失败所有要素都被判定为“不相交”。再一种常见原因是shp几何本身有自相交ArcGIS在裁剪时把坏几何一并清掉了。解决线边界先用“要素转面”转成面统一坐标系后重跑裁剪。裁剪前对原数据执行“检查几何”和“修复几何”这两个工具在数据管理工具里跑一下能排除大部分几何损坏问题。修复后如果还是缺就用“按位置选择”把边界内的要素选中后导出换个思路绕开Clip的严格判定。5.5 ArcGIS启动时提示许可没反应现象双击ArcMap或ArcGIS Pro图标出现后没反应一会弹许可证错误或者直接闪退。原因新版许可服务没有启动机器码变动导致许可失效或者杀毒软件拦截了许可进程。这和shp数据本身无关但经常在“好不容易弄到数据”的节骨眼上冒出来让人火大。解决WinR输入services.msc找到ArcGIS License Manager相关的服务确认状态是“正在运行”不是就右键启动。启动后仍失败打开License Server Administrator重新读取许可文件。大多数“点启动没反应”都是服务停了不用急着重装软件。血泪经验是重装前先试服务重启能省半天时间。6. 数据验证与成果交付让鄱阳湖shp在ArcGIS里真正“能用”6.1 用在线底图叠加验证shp精度所有投影、裁剪做完最后一步是验证。把shp叠加到在线影像底图上重点看三处鄱阳湖主湖岸线是否贴合影像水边线赣江、抚河、信江、饶河、修水的主河道和影像走向是否一致等高线在山区部分是否沿山脊线走。放大到1:5万以内再下结论别在缩略图里觉得“差不多”量级不同。如果偏移出现在局部而不是整体问题多半出在原始数据本身局部是旧版测量成果拼接造成的。遇到这种情况不要硬调整个图层的坐标用“空间校正”工具做局部配准校正后再做一次几何检查。6.2 shp与dwg、GeoJSON、3DTiles互转的取舍交付对象不一定都用ArcGIS。设计院常用CAD那就用ArcToolbox里的“CAD转地理数据库”把shp转成dwgWeb端要展示转GeoJSON最通用生成3D Tiles则没有“一键转换”的说法ArcGIS原生也不直接导出3DTiles需要先给shp配高程信息再借助专门的转换工具处理。转GeoJSON可以用GDAL命令行一行就能完成ogr2ogr -f GeoJSON poyang_rivers.geojson poyang_rivers.shp -t_srs EPSG:4326 -lco ENCODINGUTF-8-t_srs EPSG:4326把结果转到WGS84经纬度Web地图直接用-lco ENCODINGUTF-8让属性表里的中文不乱码。转换前确认原始shp的坐标系已经有prj定义否则ogr2ogr不知道从哪里转起。6.3 面积计算的精度问题用椭球面积还是平面面积最后处理面积字段时很多人直接在属性表里右键计算几何选了平面面积结果和真实面积差出几个百分点。鄱阳湖这种大范围水域平面面积受投影变形影响明显应该用椭球面积。import arcpy arcpy.env.workspace rD:\poyang\poyang.gdb arcpy.management.CalculateGeometryAttributes( clip_lake_polygon, [[AREA_KM2, AREA_GEODESIC]], area_unitSQUARE_KILOMETERS, coordinate_systemGCS_CN_2000 )AREA_GEODESIC按椭球面计算不依赖投影带coordinate_system指定用CGCS2000地理坐标参与计算出来的面积更接近真实值。这个参数对省级尺度的水域统计尤其重要别嫌麻烦。我现在的习惯是任何shp交付前都做三件事补全.prj和.cpg另存一个地理坐标系副本在成果说明里写清楚投影带。这套习惯救过我很多次不然换个电脑、换个软件前面所有分析都得重来。希望帮到你。本文还有配套的精品资源点击获取
返回列表