简介面向遥感与GIS学习者的高分一号影像水体反演完整流程资源以重庆部分区域为案例重点演示黑臭水体与水华识别方法全程基于ENVI和ArcGIS两大软件操作。包体内含1个docx文档压缩包约13.95MB文档按数据下载、预处理、反演及出图三大阶段组织。数据下载部分介绍陆地观测卫星数据服务网站的参数检索要点预处理涵盖辐射定标、FLAASH大气校正、RPC正射校正及水体掩膜裁剪反演阶段给出NDWI与MNDWI波段计算公式并借助ArcGIS完成符号化制图与空间统计分析。适合需要系统掌握遥感水质监测流程或准备相关课程设计、毕业设计的学生与从业者。该资源已有557人学习下载文档内附各工具的具体参数设置与操作截图可帮助读者避开常见坑点快速复现从原始影像到最终水体分布图的完整链路。1. 高分一号水体反演全流程从原始影像到黑臭水与水华分级图重庆部分区域的黑臭水体和水华反演不是一个“跑通就行”的活。很多人在 ENVI 里一步一步操作完到了 ArcGIS 出图阶段才发现水体指数算出来的结果要么噪点爆炸要么大片是黑色无效值要么水域边界把江心洲都吞了。问题十有八九不是反演公式写错而是前面的数据下载、辐射定标、大气校正、正射校正和掩膜裁剪某一个环节埋了雷。这篇笔记把从数据获取到最终成图的完整流程拆开讲透覆盖陆地观测卫星数据服务检索、高分一号 WFV 数据下载、ENVI 预处理参数设置、FLAASH 大气校正的细节、RPC 正射校正与去黑边操作最后落到 ENVI 波段计算器跑 NDWI/MNDWI 以及 ArcGIS 出图分级。适合正在做遥感水质课程设计、本科毕业论文或刚开始接触环境遥感反演的从业者照着复现。2. 数据下载与产品选型检索参数、云覆盖阈值与 GF-1 WFV 传感器选择2.1 为什么选高分一号 WFV 传感器高分一号搭载了两类载荷一类是 2 米分辨率全色加 8 米多光谱的 PMS另一类是 16 米分辨率四波段多光谱的 WFV 宽幅相机。水体反演用的是 WFV不是 PMS。原因很简单水体面积大但光谱差异细微需要的是覆盖范围广、重访周期短、能捕捉水色变化的宽幅数据16 米分辨率对城市河道和大面积湖库来说已经够用而且单景覆盖范围可以达到 800 公里幅宽做区域水体监测一次成像能把大片水域收进来。PMS 的 8 米多光谱虽然在空间细节上更好但覆盖范围小、单景面积有限用来做整个区域的连续水体监测需要镶嵌很多景工作量翻倍而且拼接缝容易干扰水质反演。GF-1 的 WFV 传感器有 WFV1、WFV2、WFV3、WFV4 四个相机各相机光谱范围基本一致都是四个波段蓝光约 450-520nm、绿光约 520-590nm、红光约 630-690nm、近红外约 770-890nm。这次重庆部分区域的案例用的是 WFV3 相机光谱响应曲线文件在 ENVI 里可以直接调取后面的 FLAASH 大气校正需要指定 gf1-wfv3 这个滤波函数文件。验证一下数据头文件里 Sensor Type 是否写对WFV1 和 WFV3 的光谱响应略有差异选错光谱函数会让大气校正后的反射率出现系统性偏移表面看不出问题NDWI 阈值一变就全乱了。2.2 陆地观测卫星数据服务检索参数数据下载步骤在陆地观测卫星数据服务网站完成选“通用卫星载荷检索”入口。检索时有几个参数直接决定拿到的数据能不能用检索参数推荐设置说明时间范围反演目标时段前后 15-30 天时间太长水体状态变化大时间太短可能一景都没有空间范围重庆目标区域行政区或经纬度框用经纬度矩形框更精确避免下载无关区域浪费硬盘云覆盖量≤ 10%最多放宽到 20%重庆多雾20% 以上的云量基本没法做大气校正传感器GF-1 WFV注明 WFV1/WFV3 等相机号不同相机光谱响应不同FLAASH 参数设置时要用对应文件产品级别L1A 级预处理产品这是几何粗校正产品适合接下来自己做辐射定标和大气校正检索出来的结果里同一个时间段往往有多景候选影像不要只看文件名点开查看每景的云覆盖百分比和缩略图。云覆盖量是整个文件的统计值容易出现“整体云量 5%目标河道正好被一片云盖住”的情况这时候老老实实换一景重复操作不然大气校正时被云覆盖的水体反射率会被强行插值反演结果会出现一条条“云影假水体”。2.3 下载后的第一步先做数据体检下载下来的高分一号数据一般是压缩包解压后里面包含一个 XML 元数据文件、一个 TIFF 或 GeoTIFF 影像文件以及 RPC 系数文件。我一般会先用 ENVI 打开影像检查三件事影像波段数量是否为 4 个、空间分辨率是否为 16 米、影像最上方是否有明显条带或坏线。检查元数据里的成像时间也很重要FLAASH 大气校正需要精确到小时的成像时刻如果元数据时间和下载列表里的时间对不上输出结果在校正时用的太阳高度角会不准水体反射率会整体偏高或偏低。这个体检过程只有 5 分钟但能避开后面两小时的白费功夫。3. 辐射定标与 FLAASH 大气校正参数表、光谱响应函数与常见错误3.1 辐射定标把 DN 值变成亮度值高分一号下载的 L1A 级产品存的是原始 DN 值也就是传感器记录的数字量化值。DN 值本身没有物理意义不同时间、不同光照条件下同一个地物的 DN 值可能差很多所以必须先做辐射定标把 DN 值转换成辐亮度或表观反射率。在 ENVI 的工具栏搜索 rad打开 Radiometric Calibration 工具选择打开的高分影像。这里有几个关键参数定标参数推荐设置原因Calibration TypeRadiance输出辐亮度供 FLAASH 使用Output Data TypeFloat浮点型避免整型截断导致波段差值失真Scale Factor1.0保持原始量级FLAASH 能正常读取Output StoredInterleaveBIL 或 BSQ 均可以FLAASH 推荐 BSQ点击 OK 后输出辐亮度文件文件名我用 rad.dat 这种明显标识方便后面找文件。定标结果可以快速拉一个波段查看亮度和对比度重庆地区植被覆盖率高正常影像应显示绿色植被区域明显偏亮、水体区域较暗如果整体亮度异常那大概率是数据类型设置不对。3.2 FLAASH 大气校正参数设置辐射定标完成后在 ENVI 工具栏搜索 fla点击 FLAASH Atmospheric Correction。这里的大气校正是整个流程里参数最多、最容易出问题的一步也是决定反演精度上限的一步。首先是 Radiance Image 选择刚才的辐射定标结果Output Reflectance File 设置输出反射率文件名这里我习惯命名为 ref.dat。接着需要填写中心点经纬度、传感器类型、飞行时间和大气模型。中心点的经纬度直接读 FLAASH 界面提供的 RFL 信息或者在影像几何信息里查看注意坐标要填对填到别的城市去了大气模型和太阳几何角度全部错位。大气模型的选择以影像获取季节和纬度带为基准重庆地处中低纬度气候湿润夏季选 Tropical冬季选 Mid-Latitude Summer 或 Mid-Latitude Winter。重庆的最大问题是水汽含量高大气水汽反演经常报错所以在 Multispectral Settings 里打开 Filter Function File 时一定要确认选择的 gf1-wfv3 光谱响应函数和图像传感器一致否则水汽反演这一列数据会出现大量 0 值。气溶胶模型城区及工业区选 Urban郊区植被密布选 Rural。重庆主城区水体周边建筑密集气溶胶成分复杂选 Urban 更稳。初始能见度 Initial Visibility 通常设置为 20-40 km重庆多雾天气能见度很可能只有 10-20 km如果设置太高大气校正输出的反射率在蓝绿波段容易出现负值。3.3 光谱响应函数加载方法在 FLAASH 界面点击 Multispectral Settings打开后按 Open然后进入 Spectral Library 找 gf1-wfv3 光谱响应曲线文件。这个文件是 ENVI 内置的高分一号 WFV3 传感器光谱响应函数误差不要自己拍脑袋去修改。选择完成后 FLAASH 会读取每个波段的响应函数来计算大气透过率和散射贡献。这里要清楚一点如果找不到 gf1-wfv3直接选了 gf1-wfv1大气校正虽然能跑完但输出的反射率在绿光波段和近红外波段会有一个系统性偏差NDWI 的值整体被抬高或压低后面阈值分割时用固定值就不准了。设置好后点击 ApplyFLAASH 开始大气校正过程。处理时间取决于影像大小和机器性能16 米分辨率的 WFV 数据一般需要 15-30 分钟。如果中途报错或者输出的反射率影像出现大面积负值优先检查两个地方一是输入影像是否真的是辐亮度而非反射率二是能见度和气溶胶模型是否过于激进。负值偶尔出现在个别像元是正常现象出现在水体区域则说明初始能见度设置偏小或气溶胶模型偏滞后。3.4 大气校正结果的验证方法大气校正做没做准不用等到最后反演阶段。快速验证方法是对比水体在近红外波段的反射率深水体的近红外反射率通常很低一般在 3% 到 5% 以下。如果近红外反射率在 8% 以上要么是水体本身含沙量很高要么是大气校正没做干净。另一个方法是看水体在蓝光波段的反射率是否比绿光低正常清水体符合这个规律如果蓝光高于绿光且数值异常偏大说明大气校正的瑞利散射部分没有消除干净。4. 正射校正、去黑边与水域裁剪RPC 模型参数与掩膜处理细节4.1 RPC 正射校正参数选择大气校正后的反射率影像还存在地形起伏引起的几何形变。重庆是典型山地城市地形高差大河道两侧有大量梯坎和山体如果不做正射校正山体阴影和水体会互相错位在 ArcGIS 里叠加水域边界时会出现边界错开的地方。ENVI 的正射校正工具是 RPC Orthorectification Using Reference Image在工具栏搜索 RPC 找到它。参数设置上输入文件是大气校正后的反射率影像参考影像可以用已有的经过精确校正的高分辨率影像或具有准确地理定位的底图高程数据用覆盖重庆区域的 DEM我用的是 ASTER GDEM 30 米分辨率高程数据对 16 米分辨率的 WFV 数据足够。重采样方法选三次卷积输出像元大小设置为 16 米保持原始分辨率。需要注意 RPC 正射校正时需要读取影像自带 RPC 系数文件下载的数据包里的 RPC 文件是关键没有它正射校正是跑不动的。如果 ENVI 提示找不到 RPC检查影像文件是否和 .rpc 同名且放同一目录某些下载数据解压后会丢失 RPC 文件导致无法处理这时候只能重新下载。4.2 去除黑色背景值Data Ignore Value 设置正射校正输出的结果影像边缘通常有一圈黑色背景这些像元值为 0但和真实水体里的 0 反射率完全不同。如果不处理后续波段计算和水域裁剪会把黑色背景当成“超低反射率纯黑水体”NDWI 计算时这些区域的值会被拉低水域提取结果会出现从影像边缘延伸进内部的假水道。去除黑边要双击输出结果图层选择 Edit Metadata然后点击 Add Metadata items找到 Data Ignore Value 并添加将值设置为 0点击 OK 确认。这一步是把黑色背景像元标记为忽略值后续统计和波段运算时 ENVI 自动跳过这些像元。做完这步再打开图层黑边部分不会参与后面的数值计算。4.3 Build Mask Apply Mask 裁剪水体反演黑臭水体和水华只关心水体区域陆地和植被的信号会干扰阈值提取。水域裁剪分两步第一步用 Build Mask 构建掩膜第二步用 Apply Mask 应用掩膜。Build Mask 里输入文件选正射校正后的影像然后载入重庆水域 shp 矢量数据。这里有个容易踩的坑Build Mask 默认认为加载进来的矢量区域是“感兴趣区”但需要确认掩膜输出的逻辑值方向。在面板里把 Shp 区域的值设为 1其他区域设为 0保存成 mask 文件。如果你的 shp 是线状水系或者包含很多岛礁注意检查 shp 拓扑是否闭合未闭合的 shp 生成掩膜时会出现锯齿状空洞。Apply Mask 操作里输入文件选择正射校正结果掩膜文件选择刚才生成的 mask输出文件命名为 water.dat。通过这个操作陆地区域全部被掩膜掉只留下水体像元。打开 water.dat 检查一下重庆的长江、嘉陵江以及周边水库和池塘应该清晰可见而山地和城区全部变成忽略值。如果水体区域出现了断裂说明输入的 bool mask 值方向反了把陆地保留了下来而水体被裁掉重新回到 Build Mask 检查。5. 避坑排查定标过暗、黑边残留、掩膜反向等五个实测问题5.1 辐射定标后影像整体偏暗水体细节看不见现象做完辐射定标打开影像发现整体亮度很低水体区域几乎全黑拉直方图也救不回来。 原因定标时输出数据类型选成了 Byte 整型辐亮度数值范围大整型存储会把大量中间灰度截断导致暗部被压缩成一团。 解决重新做辐射定标把 Output Data Type 设为 Float输出浮点型辐亮度。操作不麻烦但很多人第一步顺手默认就选了 Byte后续 FLAASH 读取这样的“假辐亮度”数据也会出错。记住一条ENVI 做定标和大气校正输出类型永远是 Float宁可后面存储空间多占一点也不能为了省空间把物理量级截断。5.2 FLAASH 大气校正后输出的反射率出现大量负值现象大气校正后打开影像水体区域大量像元反射率为负整幅图像是打了马赛克。 原因两种可能一是初始能见度设置过大把大气浑浊度估计得过于干净导致校正算法把大气散射分量过度扣除二是气溶胶模型选错把乡村模型用到城市污染严重的区域或者反过来把城市模型用到干净山区。 解决重新跑 FLAASH把能见度从 40 km 调到 20 km气溶胶模型选 Urban然后查看水体区域在近红外波段的取值。通常把能见度调低一档后负值大幅减少。注意重庆这种多雾低能见度城市不要照抄平原地区的能见度参数。5.3 正射校正后黑边没去掉掩膜裁剪把黑边也包进来了现象Apply Mask 后输出的水域结果不是按 shp 边界裁剪而是连同黑色背景一起形成了一个矩形区域水体边缘全是锯齿状黑块。 原因跳过 Edit Metadata 设置 Data Ignore Value 的步骤黑色背景像元值为 0在 Build Mask 的时候软件会把这个区域判定为有效数据导致掩膜和黑边一起作用。 解决回到正射校正结果重新设置 Data Ignore Value0保存后再做 Build Mask 和 Apply Mask。如果已经做了掩膜直接重新建掩膜重来不要懒黑边一旦参与掩膜计算后期每步都会受影响。5.4 波段计算时提示“Division by zero”或输出大量异常大值现象用 Band Math 计算 NDWI 时输出结果里出现很多极端大正数或大负数水体提取结果像雪花屏。 原因分母为 0 或接近 0。正射校正结果里黑色背景像元值虽然是 0但如果没设置 Data Ignore Value波段计算器会正常参与运算(b2b4) 等于 0 的区域在浮点运算里就变成无穷大。 解决两种方式一是在波段计算前确保设置了 Data Ignore Value0ENVI 在波段运算时会自动跳过忽略值二是把表达式写成 float 类型保护例如 (float(b2)-float(b4))/(float(b2)float(b4))。另外注意计算前先检查影像的坐标系和范围是否一致不一致的影像做波段运算结果也是错位相加。5.5 水域裁剪后长江主河道断裂成很多小段现象水域 shp 本身连续但裁剪后的河道出现多处断裂尤其是江心洲两侧和桥梁附近。 原因两类可能。一类是 shp 数据本身经过简化处理在江心洲和桥梁位置存在多边形的自相交或空洞另一类是掩膜生成时投影变换导致矢量边界和栅格像元不对齐像元重采样后细窄河道被吞掉。 解决对 shp 做拓扑修复在 ArcGIS 里用 Repair Geometry 工具处理后再放进 ENVI同时确认 shp 和影像的坐标系保持一致影像通常是 UTM 投影shp 也需要是同一个投影。投影不一致时先在 ArcGIS 里做投影转换不要指望 ENVI 自动处理好。6. 波段计算反演与 ArcGIS 制图从水体指数到黑臭水/水华分级出图6.1 NDWI 与 MNDWI 的波段计算器实现在 ENVI 工具栏搜索 band打开 Band Math 波段计算器输入 NDWI 表达式并以 Map Variable to Input File 映射到 water.dat 裁剪后的水体影像。高分一号 WFV3 的波段顺序是b1 蓝光、b2 绿光、b3 红光、b4 近红外。NDWI 公式用绿光和近红外波段# NDWI - 归一化差异水体指数GF-1 WFV3 对应波段2绿光、波段4近红外 (float(b2)-float(b4))/(float(b2)float(b4))表达式里的 float() 运算符做浮点转换避免整型相除时接近 0 的分母导致极值。NDWI 的物理含义是绿光反射率减去近红外反射率再除以二者之和水体在绿光反射率高于近红外所以水体区域的 NDWI 通常为正值而植被和土壤为负值。长江和嘉陵江主河道应该呈现明显正值区而支流和水库因为水深不同会有梯度差异。MNDWI 的表达式则把绿光波段与中红外波段之间的差异换成了绿光与近红外需要再确认一下 GF-1 WFV 没有中红外波段所以严格意义上 MNDWI 在 WFV 数据上是没法直接算的。但很多教材项目里说的“MNDWI”其实是用绿光和红光波段的差值替代表达式可以写成# MNDWI 替代表达式用绿光和红光波段的差异增强水体信号 (float(b2)-float(b3))/(float(b2)float(b3))注意这是一种经验替代不要把它当成严格意义上的 MNDWI 去发表论文。实际处理时我一般优先用 NDWI 做水体提取主判据用绿光和红光差值来辅助区分高浊度水体因为黑臭水体的典型特征是绿光和红光反射率关系发生改变。6.2 黑臭水体和水华的阈值分割NDWI 计算结果是一幅连续变化的灰度图需要设定阈值把水体信息分类。水体区域的 NDWI 阈值定位在 0 附近但重庆的长江水体含沙量较高近红外反射率偏大实际 NDWI 在 0 附近的像元很多。查看 NDWI 直方图在 0 到 0.2 之间会看到两个波峰第一个波峰对应含沙量较高的浑浊水体第二个波峰对应较清澈的水体阈值选在两个波峰之间的谷底处。黑臭水体的判别在光谱上表现为蓝紫光反射率显著升高绿光反射率相对降低红光反射率升高水体整体呈现黑暗甚至深褐色。实际操作中我习惯用两个波段组合来判别一是用蓝光波段和绿光波段的差值 (b1-b2) 提取发黑发暗的水体二是用红光波段和绿光波段的差值与 NDWI 值组合分级。在 ENVI Band Math 里输入分类表达式# 黑臭水体判别条件波段值均为反射率条件是经验阈值需要按区域微调 ((float(b1)-float(b2)) 0.02) AND ((float(b4)) 0.05) AND (NDWI -0.1)水华区域则表现为叶绿素浓度高近红外反射率异常升高所以用 NDVI 更好# 水华区域判别近红外高于红光叶绿素富集 (float(b4)-float(b3))/(float(b4)float(b3))水华区域 NDVI 明显高于周围水体用直方图分割找异常高值区一般阈值在 0.2 以上。使用 ArcGIS 的栅格计算器配合条件判断把正常水体、黑臭水体和疑似水华区域合并成一个三值栅格方便出图和统计。建议保留原 NDWI 连续值图层用于后续定量分析而分类结果只作为可视化展示。6.3 ArcGIS 制图与面积统计技巧反演结果最后要落到 ArcGIS 出图。先把 ENVI 输出的分类栅格导出成 GeoTIFF然后用 ArcMap 加载设置不同分类值的符号正常水体用浅蓝色黑臭水体用深褐色疑似水华区域用黄绿色。图例、比例尺、指北针和文字注记逐一加上重庆长江和嘉陵江交汇处的点位标注清楚图名写成“重庆部分区域黑臭水体与水华分布反演图”再导出成 JPEG 或 PDF 作为成果输出。面积统计在 ArcGIS 里用栅格转面工具把分类栅格转成多边形要素类然后新建字段计算 Area按分类值汇总面积即可。注意利用栅格转面前先把分类栅格的 Data Ignore Value 或 NoData 值清理好不然统计会把无效像元算进去。如果是论文或报告建议同时在 ArcGIS 里生成采样点的 NDWI 值提取用来和实地水质监测数据做相关性分析这样反演成果的可靠性更高。我自己的习惯是无论项目要求多紧从数据下载到反演出图每一步都保留过程文件。数据下载的 XML、定标后的 rad.dat、大气校正后的 ref.dat、正射校正后的 ortho.dat、裁剪后的 water.dat、NDWI 中间结果和最终分类栅格全部按阶段存好每一层的元数据里都备注当时用的参数。尤其是原始数据下载时那个云覆盖量和影像的规范性检查从那以后我每次做重庆区域的水体反演都会强制走一遍“先体检再开工”的流程绝不在没确认云量分布和传感器型号的情况下盲目开跑。做遥感反演翻车往往不在最后一步而在最开始那些你没耐心检查的参数上。希望这篇笔记能帮你少走点弯路顺利把整套流程复现出来。本文还有配套的精品资源点击获取