简介面向遥感与GIS应用学习者及环境监测从业者的实操文档围绕高分一号影像在重庆部分区域的水体监测展开目标是从数据获取到成果出图完整走通黑臭水体与水华反演流程兼顾ArcGIS与ENVI协同使用。压缩包仅含1个docx文件大小约13.95MB以图文步骤组织操作说明便于对照软件界面复现实验。目前已有560人学习适合具备基础遥感概念、希望补齐工程化处理链条的读者。内容覆盖数据检索与下载、辐射定标、基于FLAASH的大气校正、RPC正射校正与黑色背景值处理、水域掩膜裁剪以及用Band Math计算NDWI/MNDWI提取水体并在ArcGIS中完成符号化、制图与面积统计。读者可理解各环节参数设置依据掌握从原始影像到专题图的完整方法并积累可迁移到其他高分影像水体识别任务的实践思路。1. 一景 WFV 影像到一张水情分级图这套链路能解决什么一景 16 米分辨率的高分一号 WFV 影像摊在屏幕上水面只是灰绿一片哪里在发黑发臭、哪里是藻类堆积肉眼能感觉到不对却说不清位置、边界和面积。这套流程就是把这种“感觉不对”变成能交付的成果从检索下载开始经辐射定标、FLAASH 大气校正、RPC 正射校正用水域矢量把陆地抠掉再用 ENVI 的波段运算反演出黑臭水体和水华的空间分布最后进 ArcGIS 分级出图并统计面积。它适合做水环境遥感监测的从业者、要交毕业设计的学生也适合已经会用 ENVI 却每次卡在大气校正参数上的人。这不是一键工具而是一条有明确参数、也有明确边界的作业链路参数对了一景影像从下载到成图能压在一个工作日内参数错了后面每一步都在放大误差。跑通一次换区域、换时相都能复用。2. 数据获取与辐射定标云量、时间窗和定标参数怎么定死这一章把两件事放在一起讲——下数据和把 DN 值变成有物理意义的辐射亮度。原因很直接检索时选错时相或者定标系数根本没生效错误会一路传导到最后的反演图上而且在前面的环节几乎看不出来。2.1 检索下单时间窗、云量和传感器编号要一起定死在某陆地观测卫星数据服务平台上按通用卫星载荷检索界面看起来选项很多真正决定成败的就四个时间窗、云量、传感器、覆盖范围。时间窗的选法取决于你要看什么——黑臭水体优先选枯水期或非汛期这时水体流动性差、外源输入相对稳定光谱特征更容易和正常水体拉开差距水华则要选夏秋季藻类高发期最好是连续晴好天气之后三五天。如果一景影像想同时看两种现象春末或初秋是折中方案但要有心理准备两类目标的阈值都不太好标定。云量筛选是这里最容易想当然的地方。平台给的云量是整景统计值一景 WFV 幅宽很大整景云量 5% 完全可能刚好压在你那条江上。我的习惯是云量卡在 10% 以内下载后先目视扫一遍水面区域有薄云或云影就直接换时相别指望后面用算法修补。传感器选型要看你研究的水面宽度。WFV 是 16 米分辨率、四波段蓝、绿、红、近红外幅宽大一条干流加两岸缓冲一景就能盖住如果是城区内河、支流这类窄水面16 米一个像元里混了岸线、建筑和树冠反演出来的指数基本没有物理意义这种情况得改用全色多光谱数据代价是幅宽窄覆盖同样范围要拼好几景正射和匀色的工作量翻倍。检索条件建议取值说明时间窗枯水期看黑臭夏秋季看水华避开强降水后三到五天的浑水期云量不大于 10%整景统计值水面区域仍需目视复核传感器宽水面用 WFV窄河道改用全色多光谱16 米对城区内河偏粗覆盖范围研究区外扩 2 到 3 公里给正射校正和水域掩膜留缓冲余量下载下来的压缩包里通常有三样东西影像文件、描述成像参数和定标系数的元数据文件、RPC 系数文件。这三个文件必须放在同一个目录下、保持原始文件名ENVI 打开影像时会自动去同级目录找元数据和 RPC改名或者拆开存放后面正射校正会直接报找不到 RPC 系数。2.2 辐射定标GUI 里那几个必须改的项辐射定标的意义是把无量纲的 DN 值换算成辐射亮度让像元值有物理量纲同时消除传感器响应差异。这一步没做好FLAASH 输入的就不是辐射亮度输出的大气校正结果会整体跑偏而且报错信息通常含糊属于典型的“黑匣子”环节。操作路径是在 ENVI 里打开影像工具栏搜索 rad点 Radiometric Calibration选中刚打开的影像然后按下表设置参数。参数设置理由Calibration TypeRadiance送进 FLAASH 的必须是辐射亮度不是表观反射率Output InterleaveBILFLAASH 对 BIL 组织方式更友好内存占用也更低Output Data TypeFloat用整型会把小数部分直接截掉Scale Factor1.0不做线性缩放保持浮点物理量参数里唯一需要动脑的是 Calibration Type。有人习惯直接选 Reflectance 一步到位省掉大气校正这在做植被指数时勉强能用但黑臭水体反演靠的就是可见光波段的细微差异表观反射率里还混着气溶胶散射两类水体在数值上会糊成一片。真正的坑在定标系数上。ENVI 能不能自动完成定标取决于元数据文件里有没有它认识的增益和偏移字段。国产卫星数据在不同 ENVI 版本下的识别情况不一致常见现象是工具跑完了、界面也提示成功但输出值明显不对。这时候要手动打开影像头文件File 菜单下的 Edit ENVI Header找到增益和偏移两个字段把数据提供方随数据发布的绝对定标系数按波段顺序填进去。增益的量纲一般是 DN 转辐射亮度填错量纲比不填更麻烦。2.3 定标结果自检三个数字判断有没有翻车定标跑完不要急着往下走先看一眼统计值。打开定标结果的统计面板重点看三件事数值量级、有没有大面积负值、波段之间的相对关系是否合理。如果整幅影像的最大值只有个位数或者十几基本可以判定定标系数没生效你拿到的还是 DN 值。如果出现大范围负值多半是偏移量填错了符号或者量纲。如果陆地区域的近红外波段反射值反而低于红波段那不是定标的问题是波段顺序搞错了——WFV 常规顺序是蓝、绿、红、近红外但不同分发批次偶尔会有差异打开影像后先看头文件里的波长字段确认一遍比事后返工便宜得多。还有一个容易被忽略的点定标输出的文件名和路径里不要带中文和空格。ENVI 本身能读但后面 FLAASH 调用外部辐射传输模型时路径里的非 ASCII 字符会让它悄无声息地失败日志里只留一行看不出所以然的报错。我一般会在项目根目录下建一个全英文的临时目录所有中间结果都往那儿放。3. 从 FLAASH 到 RPC 正射让像元值有物理意义让影像落到坐标上大气校正和正射校正是两个方向的修正前者修的是“值”把传感器收到的辐射亮度还原成地表反射率后者修的是“位置”把有倾斜、有地形位移的影像拉到统一的地图坐标上。这两步顺序不能反先做大气校正再做正射理由是 FLAASH 需要用到原始成像几何下的辐亮度统计重采样过的影像会引入插值误差。3.1 进 FLAASH 之前的三项确认在工具栏搜索 fla 打开 FLAASH 之前有三件事必须先确认否则后面全是白跑。第一输入影像必须是辐射亮度且已转为浮点型整型输入会让 FLAASH 报数据类型不匹配。第二头文件里要有完整的成像时间格林尼治时间和太阳几何参数这些是 FLAASH 计算太阳天顶角和方位角的依据缺了就只能在界面里手填手填的经纬度一旦偏了半度大气模型选得再对也没用。第三影像中心经纬度要从元数据里读不要凭印象填——一景 WFV 幅宽很大用研究区中心的坐标代替景中心边缘区域的大气校正会有系统性偏差。3.2 FLAASH 参数逐项拆解FLAASH 的参数面板看着吓人真正影响结果的就那么七八项其余保持默认即可。参数取值思路说明传感器类型列表里没有对应项就选未知多光谱类型国产卫星基本都要走自定义光谱响应这条路传感器高度按元数据填WFV 约 645 公里填错会直接影响大气路径长度计算像元大小16 米与原始分辨率一致不要在这里重采样地面平均高程用 DEM 统计研究区均值山区和平原差别很大不能一律填 0大气模型按纬度和季节选中纬度夏季或中纬度冬季研究区在亚热带夏季水华监测用热带模型更稳气溶胶模型城区用城市型郊区用乡村型黑臭水体多在城区选城市型气溶胶反演选无反演手动给初始能见度四波段数据缺短波红外反演算法跑不动初始能见度城区 20 到 40 公里有当日气象观测就用实测值气溶胶反演这一项值得多说两句。经典的双波段反演方法依赖短波红外波段来估算气溶胶光学厚度WFV 四个波段里没有短波红外硬开这个选项要么报错要么给出离谱的结果。所以常规做法是关掉反演、手动给能见度用当日地面气象站的水平能见度观测值没有观测就按城区 20 到 40 公里估。这个近似在清洁天气下误差可以接受在雾霾天就要谨慎解读反演结果。水汽反演同样建议关闭。它需要 1.1 微米附近的水汽吸收波段WFV 没有开着只会浪费计算时间。3.3 光谱响应函数那个光谱曲线文件怎么加载进入多光谱设置面板在滤波函数文件一栏点打开切到光谱库选中与数据对应的高分一号宽幅相机光谱响应曲线文件。这一步是整个 FLAASH 流程里最容易翻车的地方因为 GF-1 的四台宽幅相机响应曲线并不完全一致用错了相机的曲线输出反射率会整体偏移几个百分点植被光谱形状看着还行水体的细微差异已经被抹平了。确认相机编号的办法是从原始文件名或元数据里找 WFV1 到 WFV4 的标识别靠猜。我见过有人图省事所有数据统一用同一个曲线文件跑出来的 NDWI 在同一区域不同时相之间差异巨大排查了两天才发现问题出在这里。3.4 输出量纲反射率被放大了 10000 倍FLAASH 默认输出 16 位整型反射率乘以 10000 存储。这个细节必须记住因为它对后续运算的影响是分裂的比值型指数NDWI、NDVI里分子分母同时带这个系数会自然抵消结果不受影响但单波段阈值判断、波段差值、以及任何涉及绝对数值的运算都必须先把结果除以 10000否则你设定的阈值全部失效。另一个心理准备是输出值会出现负数和大于 1 的数。这是大气校正残差和噪声放大造成的正常现象尤其在水体这类低反射目标上更明显不要一看到负值就以为流程错了。3.5 RPC 正射参考影像、高程和那个必须改的忽略值正射校正用有理多项式系数模型配合参考影像和高程数据把影像纠正到垂直视图。工具栏搜索 RPC打开基于参考影像的正射校正工具输入大气校正结果参考影像选同区域已有地理编码的正射底图高程数据用公开的 30 米数字高程模型。重采样方法选三次卷积输出像元大小保持 16 米输出投影要和研究区水域矢量一致这一步没统一后面的掩膜会直接空掉。正射输出后影像四周会有一圈黑色填充值这些 0 值如果不处理会参与后续所有统计和波段运算把均值往下拉、把直方图在 0 处堆出一根尖峰。处理办法是双击图层打开元数据编辑添加一个数据忽略值项把值设为 0。设完之后再打开影像黑边会被自动忽略统计值立刻恢复正常。3.6 校正结果怎么验验证不用太复杂三个动作够了。把正射结果叠加到参考底图上放大到岸线位置看是否套合偏移超过两三个像元就得回去检查高程数据和控制点再看一眼水体的光谱曲线清洁水体的可见光反射率从蓝到红递减、近红外接近 0如果近红外反而最高说明波段对应关系或者定标环节出了问题最后扫一遍反射率的数值范围乘过 10000 之后大部分地表应该在 0 到 10000 之间大面积超出这个区间就要回头查定标。4. 水域掩膜与波段运算指数反演和阈值分割的落地写法前面三步做完你手里是一景位置正确、数值有物理意义的反射率影像。接下来的目标很明确把水体从陆地里抠出来在抠出来的范围内算出能区分黑臭、水华和正常水体的指数再切成等级。4.1 构建掩膜从矢量到掩膜文件工具栏搜索 mask打开构建掩膜工具输入正射校正结果再把研究区水域矢量加载进来。这里有个顺序问题先选影像再加载矢量工具会根据矢量的坐标系自动检查是否匹配不匹配会直接提示比事后发现掩膜为空要省事。输出掩膜文件建议存成英文路径格式保持默认。掩膜文件的本质是一张二值图矢量范围内为 1、范围外为 0后面应用掩膜时就是拿这张图去乘影像。4.2 应用掩膜两个参数决定输出质量打开应用掩膜工具输入文件选正射校正结果掩膜文件选刚生成的掩膜然后设置掩膜值。掩膜值指的是“被保留区域对应的像素值”因为构建时范围内是 1这里就填 1。输出选项里可以选把掩膜外区域设为 0 或者设为 NaN做统计和出图建议选 0转成其他格式时兼容性更好做波段运算建议选 NaN能避免 0 值混进比值计算。掩膜后的结果就是纯水体影像像元总数会大幅下降后续波段运算的速度会明显提升。4.3 波段运算公式除以零和量纲这两件事工具栏搜索 band 打开波段运算公式框里写表达式然后用变量映射把 b1、b2 这类占位符指到具体波段。归一化差异水体指数的写法如下# NDWI绿波段减近红外除以两者之和 # b1 映射绿波段b2 映射近红外波段 # 末尾的 (b1b2) ne 0 用来屏蔽分母为零的像元避免产生 NaN 和无穷值 ((float(b1) - float(b2)) / (float(b1) float(b2))) * ((b1 b2) ne 0)这段公式里有三个细节。float 强制转换不能省否则整型相减再相除小数部分全被丢掉指数会变成一堆 0 和 1。末尾的逻辑判断也不能省水体边缘、云影、掩膜边界上经常出现两个波段同时为 0 的像元不做保护就会产生无效值这些无效值在出图时会变成刺眼的白点或者黑斑属于典型的血泪经验。如果输入是 FLAASH 输出、没做过量纲还原比值型指数可以直接算因为 10000 这个系数上下抵消了。但如果你要算的是单波段阈值或者波段差值公式里就得先除以 10000这一点在写公式的时候就要想清楚别等结果出来发现阈值怎么调都不对。4.4 指数选择四个波段够不够用原始流程里提到 NDWI 和 MNDWI 两条路实际拿 WFV 数据做的时候只有 NDWI 走得通。MNDWI 需要短波红外波段参与运算WFV 的四个波段里没有这是传感器能力的硬边界不是参数调一调就能绕过去的。要用 MNDWI 就得换成带短波红外的数据源代价是分辨率、幅宽、重访周期全都要重新评估。水华的识别相对好办。藻类在水面堆积会显著抬高近红外反射率用归一化差异植被指数就能把它们和清洁水体分开# NDVI近红外减红除以两者之和 # b1 映射红波段b2 映射近红外波段 # 水华区域该值明显高于清洁水体清洁水体接近 0 或为负 ((float(b2) - float(b1)) / (float(b2) float(b1))) * ((b1 b2) ne 0)黑臭水体的识别是这套流程里最没有标准答案的部分。它的光谱特征是可见光波段反射率整体偏低、红绿波段差异收窄、光谱曲线趋于平缓常见做法是用红绿波段的差值或者比值构造一个差异度指标再用现场采样点或者高分影像目视解译点去标定阈值。这里必须提醒一句文献里给出的阈值不能直接套用不同区域、不同季节、不同大气条件下的数值差异很大同一套阈值在 A 时相好用换到 B 时相可能整片水体都被判成黑臭。我的做法是每个时相至少留二十个已知类型的采样点做阈值校准宁可多花半天也别交一张自己都不信的分级图。4.5 阈值分割从浮点指数到分级图波段运算出来的是浮点指数还要切成等级。ENVI 里有两种做法用密度分割工具交互式拖动阈值适合快速看效果或者在波段运算里直接写条件判断输出分级码适合批量处理。# 把 NDWI 浮点结果切成三级大于 0.3 判为稳定水体0 到 0.3 判为疑似水体小于 0 判为非水体 # b1 映射上一步算出的 NDWI 结果 (b1 gt 0.3) * 1 (b1 gt 0) * (b1 le 0.3) * 2 (b1 le 0) * 3这段公式输出的是分级码栅格1、2、3 分别对应三类。逻辑运算在 ENVI 波段运算里返回 0 或 1乘上类别码再相加就能把多个条件拼成一张分级图。写法看着别扭但比在 GIS 里做多次重分类要快得多。分级完就该进 ArcGIS 了。把分级栅格加载进去按值符号化配一套从深蓝到浅蓝再到黄绿的色带黑臭水体和水华单独给醒目的颜色。出图前记得做一次栅格转面然后用字段计算几何统计各类面积这是汇报时最常被问到的数字。5. 避坑与排查五个让反演结果整段作废的细节这一章的五条都是实际作业里反复出现的问题每一条的表现形式都不一样但共同点是在中间环节几乎不报错直到最后出图才发现结果不能用。5.1 定标系数没生效反射率整体失真现象FLAASH 跑完输出反射率大面积超出正常范围水体像元值明显偏高或者偏低直方图形状和预期完全对不上。原因ENVI 没有从元数据里读到增益和偏移辐射定标工具“成功”执行了但实际只做了类型转换没有做物理量换算。国产数据的元数据字段命名和 ENVI 的识别规则不完全一致这种情况相当常见。解决定标后立刻看统计值最大值只有个位数就说明系数没生效。手动编辑头文件按波段顺序填入数据提供方公布的绝对定标系数重新执行定标再重新走一遍 FLAASH。已经跑完大气校正的影像没法修补只能重来。5.2 光谱响应函数选错相机现象反射率数值整体偏移几个百分点肉眼看不出来但同一区域两个时相的指数对比时差异异常大植被光谱形状也略微变形。原因高分一号四台宽幅相机的响应曲线不同用了 A 相机的曲线去处理 B 相机的数据。这个错误隐蔽性极强因为流程不报错输出图看着也“像那么回事”。解决从原始文件名或元数据里确认相机编号选对应的曲线文件。已经处理过的数据要重跑大气校正因为响应函数是在辐射传输模型内部参与计算的事后没法用系数修正。5.3 正射黑边没设忽略值统计和运算全被拉偏现象水体指数均值接近 0直方图在 0 处出现一根异常尖峰分级图边缘出现成片的无效斑块。原因正射校正后影像四周填充了 0 值这些 0 值在波段运算里参与除法产生大量无效结果同时把整体统计量往下拉。构建掩膜时如果掩膜范围没完全覆盖影像黑边还会被保留下来。解决正射输出后立刻在元数据里添加数据忽略值并设为 0确认黑边消失再继续。已经生成的指数结果无法修补只能从正射那一步重来。另外掩膜前先检查矢量范围是否完整覆盖水体边缘缺一块会让整条江的下游消失。5.4 波段顺序搞错指数算出来是另一回事现象NDWI 结果里陆地是正值、水体是负值或者水体和水华分不开分级图看着像随机噪声。原因不同分发批次的波段排列顺序可能有差异公式里的 b1、b2 映射到了错误的波段。这类错误不会报错因为公式本身语法没问题只是物理含义完全变了。解决打开影像后先看头文件里的波长字段确认每个波段对应的中心波长再去做变量映射。映射完成后可以先在小范围上算一遍用水体和陆地的已知位置做一次正负号检查确认无误再跑全图。5.5 水域矢量与影像坐标系不一致掩膜结果为空现象应用掩膜后输出全黑或者只剩下零星几个像元矢量明明叠上去看着是对的。原因矢量文件是地理坐标系影像经过正射后是投影坐标系两者在软件里叠加快照看着重合实际坐标数值差了一个投影转换。另一种情况是矢量边界画得过紧贴着水陆交界线混合像元被整片切掉水体面积偏小一到两成。解决掩膜前统一投影把矢量转到影像的投影坐标系下并做一到两个像元的缓冲把水陆交界的混合像元保留进来。已经生成的掩膜没法补救重新构建一次比反复调试省时间。6. 精度怎么验、图怎么出一份能交差的成图规范反演结果出来只是半成品能不能交差取决于两件事精度站不站得住图面说不说得清。精度验证最朴素的办法是混淆矩阵。在研究区里随机撒点每一类至少三十个用同时期的高分影像或者现场记录判读真实类别和反演结果做交叉表算总体精度和卡帕系数。经验上NDWI 提取水体的总体精度做到 90% 以上不算难难的是黑臭和水华的细分能到 75% 到 85% 就算可用。撒点的时候记得避开云影和岸边混合像元那些点会把精度拉低但拉低的原因是数据本身不可判不是算法不行。面积统计用脚本比手工点选靠谱得多ArcGIS 的字段计算能直接出平方公里import arcpy arcpy.CheckOutExtension(Spatial) # 分级栅格转面用于统计各类水体面积 arcpy.RasterToPolygon_conversion(rD:\water\grade.tif, rD:\water\grade.shp, NO_SIMPLIFY, Value) # 新增面积字段按几何计算单位平方公里 arcpy.AddField_management(rD:\water\grade.shp, AREA_KM2, DOUBLE) arcpy.CalculateField_management(rD:\water\grade.shp, AREA_KM2, !shape.areasquarekilometers!, PYTHON3)注意面积计算的前提是数据框使用等面积投影用地理坐标系算出来的平方公里是错的这一点在出图前就要检查。低版本 ArcGIS 如果不认这个面积单位写法就改用几何对象的面积字段除以一百万。制图环节有几个固定动作分层设色要给黑臭、水华、正常水体三种状态留出明显的色相差别用同一色系的不同明度打印出来分不清图例里写清阈值区间不要只写“高”“低”版面加上比例尺、指北针和数据时相缺一样审图的人就会问。我一般还会在角落放一张小尺寸的原影像缩略图方便对照异常区域的位置。从那以后我每次跑这套流程都会在辐射定标和正射两步之后各强制走一遍统计值检查宁可多花十分钟也不愿意在出图前夜发现整条江的水体面积是零。希望帮到你。本文还有配套的精品资源点击获取 
长河君 · 2026/10/31 18:56:39 · 阅读 1.2k 高分一号WFV影像黑臭水体与水华遥感反演全流程
简介面向遥感与GIS应用学习者及环境监测从业者的实操文档围绕高分一号影像在重庆部分区域的水体监测展开目标是从数据获取到成果出图完整走通黑臭水体与水华反演流程兼顾ArcGIS与ENVI协同使用。压缩包仅含1个docx文件大小约13.95MB以图文步骤组织操作说明便于对照软件界面复现实验。目前已有560人学习适合具备基础遥感概念、希望补齐工程化处理链条的读者。内容覆盖数据检索与下载、辐射定标、基于FLAASH的大气校正、RPC正射校正与黑色背景值处理、水域掩膜裁剪以及用Band Math计算NDWI/MNDWI提取水体并在ArcGIS中完成符号化、制图与面积统计。读者可理解各环节参数设置依据掌握从原始影像到专题图的完整方法并积累可迁移到其他高分影像水体识别任务的实践思路。1. 一景 WFV 影像到一张水情分级图这套链路能解决什么一景 16 米分辨率的高分一号 WFV 影像摊在屏幕上水面只是灰绿一片哪里在发黑发臭、哪里是藻类堆积肉眼能感觉到不对却说不清位置、边界和面积。这套流程就是把这种“感觉不对”变成能交付的成果从检索下载开始经辐射定标、FLAASH 大气校正、RPC 正射校正用水域矢量把陆地抠掉再用 ENVI 的波段运算反演出黑臭水体和水华的空间分布最后进 ArcGIS 分级出图并统计面积。它适合做水环境遥感监测的从业者、要交毕业设计的学生也适合已经会用 ENVI 却每次卡在大气校正参数上的人。这不是一键工具而是一条有明确参数、也有明确边界的作业链路参数对了一景影像从下载到成图能压在一个工作日内参数错了后面每一步都在放大误差。跑通一次换区域、换时相都能复用。2. 数据获取与辐射定标云量、时间窗和定标参数怎么定死这一章把两件事放在一起讲——下数据和把 DN 值变成有物理意义的辐射亮度。原因很直接检索时选错时相或者定标系数根本没生效错误会一路传导到最后的反演图上而且在前面的环节几乎看不出来。2.1 检索下单时间窗、云量和传感器编号要一起定死在某陆地观测卫星数据服务平台上按通用卫星载荷检索界面看起来选项很多真正决定成败的就四个时间窗、云量、传感器、覆盖范围。时间窗的选法取决于你要看什么——黑臭水体优先选枯水期或非汛期这时水体流动性差、外源输入相对稳定光谱特征更容易和正常水体拉开差距水华则要选夏秋季藻类高发期最好是连续晴好天气之后三五天。如果一景影像想同时看两种现象春末或初秋是折中方案但要有心理准备两类目标的阈值都不太好标定。云量筛选是这里最容易想当然的地方。平台给的云量是整景统计值一景 WFV 幅宽很大整景云量 5% 完全可能刚好压在你那条江上。我的习惯是云量卡在 10% 以内下载后先目视扫一遍水面区域有薄云或云影就直接换时相别指望后面用算法修补。传感器选型要看你研究的水面宽度。WFV 是 16 米分辨率、四波段蓝、绿、红、近红外幅宽大一条干流加两岸缓冲一景就能盖住如果是城区内河、支流这类窄水面16 米一个像元里混了岸线、建筑和树冠反演出来的指数基本没有物理意义这种情况得改用全色多光谱数据代价是幅宽窄覆盖同样范围要拼好几景正射和匀色的工作量翻倍。检索条件建议取值说明时间窗枯水期看黑臭夏秋季看水华避开强降水后三到五天的浑水期云量不大于 10%整景统计值水面区域仍需目视复核传感器宽水面用 WFV窄河道改用全色多光谱16 米对城区内河偏粗覆盖范围研究区外扩 2 到 3 公里给正射校正和水域掩膜留缓冲余量下载下来的压缩包里通常有三样东西影像文件、描述成像参数和定标系数的元数据文件、RPC 系数文件。这三个文件必须放在同一个目录下、保持原始文件名ENVI 打开影像时会自动去同级目录找元数据和 RPC改名或者拆开存放后面正射校正会直接报找不到 RPC 系数。2.2 辐射定标GUI 里那几个必须改的项辐射定标的意义是把无量纲的 DN 值换算成辐射亮度让像元值有物理量纲同时消除传感器响应差异。这一步没做好FLAASH 输入的就不是辐射亮度输出的大气校正结果会整体跑偏而且报错信息通常含糊属于典型的“黑匣子”环节。操作路径是在 ENVI 里打开影像工具栏搜索 rad点 Radiometric Calibration选中刚打开的影像然后按下表设置参数。参数设置理由Calibration TypeRadiance送进 FLAASH 的必须是辐射亮度不是表观反射率Output InterleaveBILFLAASH 对 BIL 组织方式更友好内存占用也更低Output Data TypeFloat用整型会把小数部分直接截掉Scale Factor1.0不做线性缩放保持浮点物理量参数里唯一需要动脑的是 Calibration Type。有人习惯直接选 Reflectance 一步到位省掉大气校正这在做植被指数时勉强能用但黑臭水体反演靠的就是可见光波段的细微差异表观反射率里还混着气溶胶散射两类水体在数值上会糊成一片。真正的坑在定标系数上。ENVI 能不能自动完成定标取决于元数据文件里有没有它认识的增益和偏移字段。国产卫星数据在不同 ENVI 版本下的识别情况不一致常见现象是工具跑完了、界面也提示成功但输出值明显不对。这时候要手动打开影像头文件File 菜单下的 Edit ENVI Header找到增益和偏移两个字段把数据提供方随数据发布的绝对定标系数按波段顺序填进去。增益的量纲一般是 DN 转辐射亮度填错量纲比不填更麻烦。2.3 定标结果自检三个数字判断有没有翻车定标跑完不要急着往下走先看一眼统计值。打开定标结果的统计面板重点看三件事数值量级、有没有大面积负值、波段之间的相对关系是否合理。如果整幅影像的最大值只有个位数或者十几基本可以判定定标系数没生效你拿到的还是 DN 值。如果出现大范围负值多半是偏移量填错了符号或者量纲。如果陆地区域的近红外波段反射值反而低于红波段那不是定标的问题是波段顺序搞错了——WFV 常规顺序是蓝、绿、红、近红外但不同分发批次偶尔会有差异打开影像后先看头文件里的波长字段确认一遍比事后返工便宜得多。还有一个容易被忽略的点定标输出的文件名和路径里不要带中文和空格。ENVI 本身能读但后面 FLAASH 调用外部辐射传输模型时路径里的非 ASCII 字符会让它悄无声息地失败日志里只留一行看不出所以然的报错。我一般会在项目根目录下建一个全英文的临时目录所有中间结果都往那儿放。3. 从 FLAASH 到 RPC 正射让像元值有物理意义让影像落到坐标上大气校正和正射校正是两个方向的修正前者修的是“值”把传感器收到的辐射亮度还原成地表反射率后者修的是“位置”把有倾斜、有地形位移的影像拉到统一的地图坐标上。这两步顺序不能反先做大气校正再做正射理由是 FLAASH 需要用到原始成像几何下的辐亮度统计重采样过的影像会引入插值误差。3.1 进 FLAASH 之前的三项确认在工具栏搜索 fla 打开 FLAASH 之前有三件事必须先确认否则后面全是白跑。第一输入影像必须是辐射亮度且已转为浮点型整型输入会让 FLAASH 报数据类型不匹配。第二头文件里要有完整的成像时间格林尼治时间和太阳几何参数这些是 FLAASH 计算太阳天顶角和方位角的依据缺了就只能在界面里手填手填的经纬度一旦偏了半度大气模型选得再对也没用。第三影像中心经纬度要从元数据里读不要凭印象填——一景 WFV 幅宽很大用研究区中心的坐标代替景中心边缘区域的大气校正会有系统性偏差。3.2 FLAASH 参数逐项拆解FLAASH 的参数面板看着吓人真正影响结果的就那么七八项其余保持默认即可。参数取值思路说明传感器类型列表里没有对应项就选未知多光谱类型国产卫星基本都要走自定义光谱响应这条路传感器高度按元数据填WFV 约 645 公里填错会直接影响大气路径长度计算像元大小16 米与原始分辨率一致不要在这里重采样地面平均高程用 DEM 统计研究区均值山区和平原差别很大不能一律填 0大气模型按纬度和季节选中纬度夏季或中纬度冬季研究区在亚热带夏季水华监测用热带模型更稳气溶胶模型城区用城市型郊区用乡村型黑臭水体多在城区选城市型气溶胶反演选无反演手动给初始能见度四波段数据缺短波红外反演算法跑不动初始能见度城区 20 到 40 公里有当日气象观测就用实测值气溶胶反演这一项值得多说两句。经典的双波段反演方法依赖短波红外波段来估算气溶胶光学厚度WFV 四个波段里没有短波红外硬开这个选项要么报错要么给出离谱的结果。所以常规做法是关掉反演、手动给能见度用当日地面气象站的水平能见度观测值没有观测就按城区 20 到 40 公里估。这个近似在清洁天气下误差可以接受在雾霾天就要谨慎解读反演结果。水汽反演同样建议关闭。它需要 1.1 微米附近的水汽吸收波段WFV 没有开着只会浪费计算时间。3.3 光谱响应函数那个光谱曲线文件怎么加载进入多光谱设置面板在滤波函数文件一栏点打开切到光谱库选中与数据对应的高分一号宽幅相机光谱响应曲线文件。这一步是整个 FLAASH 流程里最容易翻车的地方因为 GF-1 的四台宽幅相机响应曲线并不完全一致用错了相机的曲线输出反射率会整体偏移几个百分点植被光谱形状看着还行水体的细微差异已经被抹平了。确认相机编号的办法是从原始文件名或元数据里找 WFV1 到 WFV4 的标识别靠猜。我见过有人图省事所有数据统一用同一个曲线文件跑出来的 NDWI 在同一区域不同时相之间差异巨大排查了两天才发现问题出在这里。3.4 输出量纲反射率被放大了 10000 倍FLAASH 默认输出 16 位整型反射率乘以 10000 存储。这个细节必须记住因为它对后续运算的影响是分裂的比值型指数NDWI、NDVI里分子分母同时带这个系数会自然抵消结果不受影响但单波段阈值判断、波段差值、以及任何涉及绝对数值的运算都必须先把结果除以 10000否则你设定的阈值全部失效。另一个心理准备是输出值会出现负数和大于 1 的数。这是大气校正残差和噪声放大造成的正常现象尤其在水体这类低反射目标上更明显不要一看到负值就以为流程错了。3.5 RPC 正射参考影像、高程和那个必须改的忽略值正射校正用有理多项式系数模型配合参考影像和高程数据把影像纠正到垂直视图。工具栏搜索 RPC打开基于参考影像的正射校正工具输入大气校正结果参考影像选同区域已有地理编码的正射底图高程数据用公开的 30 米数字高程模型。重采样方法选三次卷积输出像元大小保持 16 米输出投影要和研究区水域矢量一致这一步没统一后面的掩膜会直接空掉。正射输出后影像四周会有一圈黑色填充值这些 0 值如果不处理会参与后续所有统计和波段运算把均值往下拉、把直方图在 0 处堆出一根尖峰。处理办法是双击图层打开元数据编辑添加一个数据忽略值项把值设为 0。设完之后再打开影像黑边会被自动忽略统计值立刻恢复正常。3.6 校正结果怎么验验证不用太复杂三个动作够了。把正射结果叠加到参考底图上放大到岸线位置看是否套合偏移超过两三个像元就得回去检查高程数据和控制点再看一眼水体的光谱曲线清洁水体的可见光反射率从蓝到红递减、近红外接近 0如果近红外反而最高说明波段对应关系或者定标环节出了问题最后扫一遍反射率的数值范围乘过 10000 之后大部分地表应该在 0 到 10000 之间大面积超出这个区间就要回头查定标。4. 水域掩膜与波段运算指数反演和阈值分割的落地写法前面三步做完你手里是一景位置正确、数值有物理意义的反射率影像。接下来的目标很明确把水体从陆地里抠出来在抠出来的范围内算出能区分黑臭、水华和正常水体的指数再切成等级。4.1 构建掩膜从矢量到掩膜文件工具栏搜索 mask打开构建掩膜工具输入正射校正结果再把研究区水域矢量加载进来。这里有个顺序问题先选影像再加载矢量工具会根据矢量的坐标系自动检查是否匹配不匹配会直接提示比事后发现掩膜为空要省事。输出掩膜文件建议存成英文路径格式保持默认。掩膜文件的本质是一张二值图矢量范围内为 1、范围外为 0后面应用掩膜时就是拿这张图去乘影像。4.2 应用掩膜两个参数决定输出质量打开应用掩膜工具输入文件选正射校正结果掩膜文件选刚生成的掩膜然后设置掩膜值。掩膜值指的是“被保留区域对应的像素值”因为构建时范围内是 1这里就填 1。输出选项里可以选把掩膜外区域设为 0 或者设为 NaN做统计和出图建议选 0转成其他格式时兼容性更好做波段运算建议选 NaN能避免 0 值混进比值计算。掩膜后的结果就是纯水体影像像元总数会大幅下降后续波段运算的速度会明显提升。4.3 波段运算公式除以零和量纲这两件事工具栏搜索 band 打开波段运算公式框里写表达式然后用变量映射把 b1、b2 这类占位符指到具体波段。归一化差异水体指数的写法如下# NDWI绿波段减近红外除以两者之和 # b1 映射绿波段b2 映射近红外波段 # 末尾的 (b1b2) ne 0 用来屏蔽分母为零的像元避免产生 NaN 和无穷值 ((float(b1) - float(b2)) / (float(b1) float(b2))) * ((b1 b2) ne 0)这段公式里有三个细节。float 强制转换不能省否则整型相减再相除小数部分全被丢掉指数会变成一堆 0 和 1。末尾的逻辑判断也不能省水体边缘、云影、掩膜边界上经常出现两个波段同时为 0 的像元不做保护就会产生无效值这些无效值在出图时会变成刺眼的白点或者黑斑属于典型的血泪经验。如果输入是 FLAASH 输出、没做过量纲还原比值型指数可以直接算因为 10000 这个系数上下抵消了。但如果你要算的是单波段阈值或者波段差值公式里就得先除以 10000这一点在写公式的时候就要想清楚别等结果出来发现阈值怎么调都不对。4.4 指数选择四个波段够不够用原始流程里提到 NDWI 和 MNDWI 两条路实际拿 WFV 数据做的时候只有 NDWI 走得通。MNDWI 需要短波红外波段参与运算WFV 的四个波段里没有这是传感器能力的硬边界不是参数调一调就能绕过去的。要用 MNDWI 就得换成带短波红外的数据源代价是分辨率、幅宽、重访周期全都要重新评估。水华的识别相对好办。藻类在水面堆积会显著抬高近红外反射率用归一化差异植被指数就能把它们和清洁水体分开# NDVI近红外减红除以两者之和 # b1 映射红波段b2 映射近红外波段 # 水华区域该值明显高于清洁水体清洁水体接近 0 或为负 ((float(b2) - float(b1)) / (float(b2) float(b1))) * ((b1 b2) ne 0)黑臭水体的识别是这套流程里最没有标准答案的部分。它的光谱特征是可见光波段反射率整体偏低、红绿波段差异收窄、光谱曲线趋于平缓常见做法是用红绿波段的差值或者比值构造一个差异度指标再用现场采样点或者高分影像目视解译点去标定阈值。这里必须提醒一句文献里给出的阈值不能直接套用不同区域、不同季节、不同大气条件下的数值差异很大同一套阈值在 A 时相好用换到 B 时相可能整片水体都被判成黑臭。我的做法是每个时相至少留二十个已知类型的采样点做阈值校准宁可多花半天也别交一张自己都不信的分级图。4.5 阈值分割从浮点指数到分级图波段运算出来的是浮点指数还要切成等级。ENVI 里有两种做法用密度分割工具交互式拖动阈值适合快速看效果或者在波段运算里直接写条件判断输出分级码适合批量处理。# 把 NDWI 浮点结果切成三级大于 0.3 判为稳定水体0 到 0.3 判为疑似水体小于 0 判为非水体 # b1 映射上一步算出的 NDWI 结果 (b1 gt 0.3) * 1 (b1 gt 0) * (b1 le 0.3) * 2 (b1 le 0) * 3这段公式输出的是分级码栅格1、2、3 分别对应三类。逻辑运算在 ENVI 波段运算里返回 0 或 1乘上类别码再相加就能把多个条件拼成一张分级图。写法看着别扭但比在 GIS 里做多次重分类要快得多。分级完就该进 ArcGIS 了。把分级栅格加载进去按值符号化配一套从深蓝到浅蓝再到黄绿的色带黑臭水体和水华单独给醒目的颜色。出图前记得做一次栅格转面然后用字段计算几何统计各类面积这是汇报时最常被问到的数字。5. 避坑与排查五个让反演结果整段作废的细节这一章的五条都是实际作业里反复出现的问题每一条的表现形式都不一样但共同点是在中间环节几乎不报错直到最后出图才发现结果不能用。5.1 定标系数没生效反射率整体失真现象FLAASH 跑完输出反射率大面积超出正常范围水体像元值明显偏高或者偏低直方图形状和预期完全对不上。原因ENVI 没有从元数据里读到增益和偏移辐射定标工具“成功”执行了但实际只做了类型转换没有做物理量换算。国产数据的元数据字段命名和 ENVI 的识别规则不完全一致这种情况相当常见。解决定标后立刻看统计值最大值只有个位数就说明系数没生效。手动编辑头文件按波段顺序填入数据提供方公布的绝对定标系数重新执行定标再重新走一遍 FLAASH。已经跑完大气校正的影像没法修补只能重来。5.2 光谱响应函数选错相机现象反射率数值整体偏移几个百分点肉眼看不出来但同一区域两个时相的指数对比时差异异常大植被光谱形状也略微变形。原因高分一号四台宽幅相机的响应曲线不同用了 A 相机的曲线去处理 B 相机的数据。这个错误隐蔽性极强因为流程不报错输出图看着也“像那么回事”。解决从原始文件名或元数据里确认相机编号选对应的曲线文件。已经处理过的数据要重跑大气校正因为响应函数是在辐射传输模型内部参与计算的事后没法用系数修正。5.3 正射黑边没设忽略值统计和运算全被拉偏现象水体指数均值接近 0直方图在 0 处出现一根异常尖峰分级图边缘出现成片的无效斑块。原因正射校正后影像四周填充了 0 值这些 0 值在波段运算里参与除法产生大量无效结果同时把整体统计量往下拉。构建掩膜时如果掩膜范围没完全覆盖影像黑边还会被保留下来。解决正射输出后立刻在元数据里添加数据忽略值并设为 0确认黑边消失再继续。已经生成的指数结果无法修补只能从正射那一步重来。另外掩膜前先检查矢量范围是否完整覆盖水体边缘缺一块会让整条江的下游消失。5.4 波段顺序搞错指数算出来是另一回事现象NDWI 结果里陆地是正值、水体是负值或者水体和水华分不开分级图看着像随机噪声。原因不同分发批次的波段排列顺序可能有差异公式里的 b1、b2 映射到了错误的波段。这类错误不会报错因为公式本身语法没问题只是物理含义完全变了。解决打开影像后先看头文件里的波长字段确认每个波段对应的中心波长再去做变量映射。映射完成后可以先在小范围上算一遍用水体和陆地的已知位置做一次正负号检查确认无误再跑全图。5.5 水域矢量与影像坐标系不一致掩膜结果为空现象应用掩膜后输出全黑或者只剩下零星几个像元矢量明明叠上去看着是对的。原因矢量文件是地理坐标系影像经过正射后是投影坐标系两者在软件里叠加快照看着重合实际坐标数值差了一个投影转换。另一种情况是矢量边界画得过紧贴着水陆交界线混合像元被整片切掉水体面积偏小一到两成。解决掩膜前统一投影把矢量转到影像的投影坐标系下并做一到两个像元的缓冲把水陆交界的混合像元保留进来。已经生成的掩膜没法补救重新构建一次比反复调试省时间。6. 精度怎么验、图怎么出一份能交差的成图规范反演结果出来只是半成品能不能交差取决于两件事精度站不站得住图面说不说得清。精度验证最朴素的办法是混淆矩阵。在研究区里随机撒点每一类至少三十个用同时期的高分影像或者现场记录判读真实类别和反演结果做交叉表算总体精度和卡帕系数。经验上NDWI 提取水体的总体精度做到 90% 以上不算难难的是黑臭和水华的细分能到 75% 到 85% 就算可用。撒点的时候记得避开云影和岸边混合像元那些点会把精度拉低但拉低的原因是数据本身不可判不是算法不行。面积统计用脚本比手工点选靠谱得多ArcGIS 的字段计算能直接出平方公里import arcpy arcpy.CheckOutExtension(Spatial) # 分级栅格转面用于统计各类水体面积 arcpy.RasterToPolygon_conversion(rD:\water\grade.tif, rD:\water\grade.shp, NO_SIMPLIFY, Value) # 新增面积字段按几何计算单位平方公里 arcpy.AddField_management(rD:\water\grade.shp, AREA_KM2, DOUBLE) arcpy.CalculateField_management(rD:\water\grade.shp, AREA_KM2, !shape.areasquarekilometers!, PYTHON3)注意面积计算的前提是数据框使用等面积投影用地理坐标系算出来的平方公里是错的这一点在出图前就要检查。低版本 ArcGIS 如果不认这个面积单位写法就改用几何对象的面积字段除以一百万。制图环节有几个固定动作分层设色要给黑臭、水华、正常水体三种状态留出明显的色相差别用同一色系的不同明度打印出来分不清图例里写清阈值区间不要只写“高”“低”版面加上比例尺、指北针和数据时相缺一样审图的人就会问。我一般还会在角落放一张小尺寸的原影像缩略图方便对照异常区域的位置。从那以后我每次跑这套流程都会在辐射定标和正射两步之后各强制走一遍统计值检查宁可多花十分钟也不愿意在出图前夜发现整条江的水体面积是零。希望帮到你。本文还有配套的精品资源点击获取