简介:围绕高分一号影像的重庆水体监测实操文档,面向遥感与GIS初学者、环境监测技术人员及需要复现水环境反演流程的读者,重点解决从卫星数据获取、预处理到黑臭水体与水华识别出图的完整链路问题。压缩包内仅含1个docx文件,约13.95MB,以图文操作记录形式组织,便于按章节查阅ENVI与ArcGIS的协同处理思路。内容覆盖辐射定标、FLAASH大气校正、RPC正射校正、水域掩膜裁剪、Band Math水体指数反演,以及ArcGIS符号化制图与面积统计等环节,可帮助读者理解NDWI、MNDWI在真实区域影像上的应用方式,并积累预处理参数设置和流程排错经验。目前已有560人学习,适合用作课程设计、项目实训或水体遥感监测入门的参考案例。
1. 从一张"看不出问题"的影像说起
做水体监测的人常常遇到同一个困境:现场采样数据很扎实,但一个断面只能代表几十米范围,想说明整条河、整个湖库的状况,很难有说服力;反过来,处理大片水域时如果把真彩色影像直接存成一套图,又很容易变成凭眼睛猜的过程,深一点就说黑臭,绿一点就说是水华,评审时一个问题就能把你问住。这套方案要解决的就是这件事:用高分一号 WFV 影像作为数据源,在 ENVI 里把辐射定标、大气校正、正射校正、水体掩膜裁剪一路走完,再用波段运算把水体范围、黑臭程度、水华范围变成可量化的图层,最后转到 ArcGIS 做出分级图和面积统计表。它适合遥感、GIS 方向的毕业设计和课程设计,也适合环境监测岗位想搭一套可复用流程的人。整条处理链路在 ENVI 和 ArcGIS 里完成,参数可复现、结果可复查,这才是能交付的方案。
2. 需求与数据源界定:三个可计算量加一套 WFV 影像
反演这件事最容易出问题的地方不在软件操作,而在动手之前没有把"要算出什么"写清楚。很多方案写到一半卡住,是因为一开始把目标定成了"识别黑臭水体",而这是一个定性描述,软件里没有对应的按钮。合格的做法是把它拆成能用波段运算表达的量,再倒推需要什么样的影像和辅助数据。
2.1 把监测目标拆成三个可计算的量
这套流程最终要产出三个图层。第一个是水体范围,用归一化水体指数把水体与建成区、植被、裸土分开,它决定了后面所有统计的分母。第二个是黑臭程度分级,黑臭水体的光谱特征通常表现为可见光波段整体反射率偏低、绿红波段相对抬升、近红外明显下压,所以用红波段与近红外波段的归一化差值做一个辅助指数,再按阈值切成若干等级。第三个是水华范围,藻类在近红外波段的反射率显著高于红波段,用归一化植被指数的高值区配合绿波段反射率抬升来圈定。
三个量必须共用同一套预处理结果,否则水体范围和水华范围对不上,图做出来自己都解释不通。所以流程设计上只有一个预处理出口,三套波段运算都从这个出口分叉出去,这是整条链路的结构基础。
提示:分级阈值一定要说明是初始值还是率定值。没有实测点位支撑时,阈值只能算经验区间,报告里要写清楚它的来源和不确定性,这是答辩时最常被追问的一点。
2.2 影像与辅助数据清单
高分一号 WFV 相机是这套方案的主力数据。它提供四个多光谱波段,空间分辨率 16 米,单台相机幅宽较大,四台组合后重访能力强,对内陆水体的季节变化监测够用。波段顺序通常按蓝、绿、红、近红外排列,这个顺序直接决定了后面波段运算表达式里 b1 到 b4 分别代表什么,弄反了指数就全错。
| 数据项 | 用途 | 关键参数 |
|---|---|---|
| GF-1 WFV 多光谱影像 | 反演主数据 | 4 波段,16 m 分辨率,含 RPC 参数文件 |
| 光谱响应函数文件 | FLAASH 大气校正的滤波函数输入 | 按相机编号选取,如 wfv3 |
| 数字高程数据 | 正射校正与地形阴影判断 | 与影像覆盖范围匹配 |
| 水域矢量边界 | 制作水体掩膜 | 面状要素,坐标系与影像一致 |
| 实测点位或历史水质资料 | 阈值率定与精度验证 | 点位坐标加水质分级 |
这里要强调一个常被忽略的约束:WFV 只有可见光和近红外波段,没有短波红外。这意味着在很多教材和论文里被广泛推荐的改进型水体指数在这套数据上根本算不出来,硬套公式只会得到一条报错或者一堆无意义的数值。选指数之前先看传感器波段配置,这一步比记公式重要得多。
2.3 影像筛选与获取参数
数据获取走陆地观测卫星数据服务平台的通用载荷检索入口,不按专题产品检索。检索条件一般设四个:时间范围、空间范围、云覆盖量、传感器类型。时间范围建议覆盖同一个水文年的枯水期和丰水期各一景,这样能看出水华随季节的变化;云覆盖量控制在 10% 以内,山区多云,放宽到 20% 也常有云影残留,云影会被误判成水体,处理起来很麻烦。
空间范围按水域矢量边界的外接矩形略微外扩即可,范围开太大只会拖慢后续每一步的处理速度。传感器类型勾选 WFV,注意记录下相机的编号,因为大气校正要用的光谱响应文件是按相机编号区分的,下载时漏记,到 FLAASH 那一步就得回头翻订单信息。
检索结果里通常会有同一区域的多景影像,挑选原则是优先选成像时间接近正午、太阳高度角较大的景,阴影短,水体边界清楚。下载完成后按"传感器_成像日期_轨道号"的规则重命名,这套命名习惯看着不起眼,等到处理十几景数据时能省下大量找文件的时间。
3. 预处理流水线:从辐射定标到水体掩膜的六道工序
预处理是这条链路里最耗时也最容易返工的部分。它的目标很明确:把原始 DN 值变成有物理意义的地表反射率,把影像摆到正确的地理位置上,再把研究区之外的部分裁掉。三道处理环环相扣,前一步参数设错,后一步要么报错要么结果失真,而且失真往往不会立刻暴露,等到反演出图才发现数值全挤在一个窄区间里,那就只能从头再来。
3.1 辐射定标:把 DN 值换成辐射亮度
在 ENVI 工具栏搜索 rad,打开 Radiometric Calibration 工具,选中刚打开的原始影像。参数面板里定标类型选辐射亮度值,输出数据类型设为浮点型,输出文件名按"影像名_rad"命名,方便后续步骤识别。
这一步的核心是读取影像自带的定标系数,把无量纲的 DN 值乘以增益再加上偏移,转换成具有物理意义的辐射亮度。参数本身不多,但有两个细节值得留意:一是输出必须是浮点型,整型会把小数部分截断,后面大气校正的输入精度就不够了;二是定标前后的数值范围应该有明显变化,如果两者几乎一样,多半是定标系数没被正确读取,需要检查元数据文件是否完整。
定标结果可以用统计工具快速看一眼,水体在辐射亮度图上仍然是低值区,但整体数值应该落在合理的量级区间内。如果出现大面积饱和或者负值,先别急着往下走,回头核对输入影像是否选错波段组合。
3.2 FLAASH 大气校正:光谱响应文件是关键输入
大气校正在工具栏搜索 fla,打开 FLAASH Atmospheric Correction。输入辐射定标后的影像,设置输出路径与文件名,然后按下面的思路逐项配置参数。
传感器类型这一栏里通常找不到高分一号的预设选项,这是正常现象,处理方式是选通用多光谱类型,再手动挂载光谱响应函数。打开多光谱设置面板,在滤波函数文件处点击打开,选择光谱库,找到对应相机编号的光谱响应曲线文件载入。这个文件的作用是告诉算法每个波段的响应范围,没有它,大气校正的波段匹配就是错的,结果会整体偏移。
成像时间、传感器高度、像元大小、地面平均海拔这几项从影像元数据和区域地形资料里取。大气模型按成像季节和纬度带选择,水汽反演这一项要关闭,因为 WFV 缺少用于水汽反演的吸收波段,强行开启只会得到无效值。气溶胶模型按区域类型选城市型或乡村型,初始能见度给一个中等偏上的经验值,后续可以通过结果里水体反射率曲线的形状做一次粗判。
运行完成后,最直观的检查方式是看植被像元的波谱曲线:健康植被在红波段应有明显吸收谷,在近红外有明显反射峰。如果曲线平滑得像一条直线,或者近红外反而低于红波段,说明大气校正没有生效,要回头查响应文件和大气模型。
3.3 RPC 正射校正与黑边取值处理
正射校正用来消除地形起伏和传感器姿态造成的几何变形。工具栏搜索 RPC,打开基于参考影像的正射校正工具,输入大气校正结果,参考影像选一景已经做过正射处理的同区域影像,高程数据选对应的数字高程模型,像元大小保持 16 米,重采样方法用三次卷积,山区地形起伏大时这个选择对边界精度更友好。
输出结果往往会带一圈黑色背景,这不是错误,是旋转校正后画布扩展留下的无效值。处理方式是双击图层打开元数据面板,通过添加元数据项找到无效数据值这一项,把它设为 0。设置完成后重新加载影像,黑边就会变成透明,后续做掩膜和统计时不会被算进有效像元。
这一步还有一个容易忽略的检查点:校正后的影像叠到参考底图上,水体边界应当基本贴合,如果出现整体平移或者明显拉伸,多半是 RPC 参数文件缺失或者参考影像的坐标系不匹配。这种几何偏差在单景图上不显眼,一旦和水域矢量叠加裁剪,就会切掉一部分水体或者带进大量陆地像元。
3.4 水体掩膜制作与裁剪
裁剪的目的是把计算量限制在研究水域内,同时排除山体阴影、深色建筑屋顶这类容易和水体混淆的地物。工具栏搜索 mask,先打开构建掩膜工具,输入正射校正结果,载入水域矢量数据,设置输出掩膜文件的路径和名称。构建完成后会得到一个二值掩膜,水体区域为有效值,其余为无效值。
接着打开应用掩膜工具,输入文件选正射校正结果,掩膜文件选刚生成的掩膜,输出参数里注意把掩膜值的方向设对,让水体保留、非水体被剔除。输出文件建议命名为 water.dat,后续波段运算直接引用这个名字,链路清晰。
裁剪后建议做一次像元数统计,把水域面积粗算一遍,和矢量边界给出的面积做对比。两者差距在几个百分点以内属于正常,差距过大说明掩膜方向设反了或者坐标系不一致。这个检查只要一分钟,却能挡掉后面一大类莫名其妙的错误。
4. 反演与出图:从波段运算到 ArcGIS 分级制图
预处理出口只有一个,反演则要分三路走。水体范围、黑臭程度、水华范围分别用不同的波段组合表达,每一路都要有明确的表达式、阈值和验证方式。做完这三层,再把栅格转进 ArcGIS 做符号化和统计,成果才算完整。
4.1 Band Math 表达式与波段映射
工具栏搜索 band,打开波段运算工具,输入表达式后点击变量映射,把 b1 到 b4 依次对应到裁剪结果文件的四个波段。表达式按下面的写法组织,浮点转换和除零保护都要带上。
# ENVI Band Math 表达式,波段映射:b1=蓝 b2=绿 b3=红 b4=近红外 # 1) 归一化水体指数:水体在绿波段高、近红外低,结果落在 -1 ~ 1 NDWI = "(float(b2) - float(b4)) / (float(b2) + float(b4) + 0.000001)" # 2) 黑臭辅助指数:红波段相对近红外的抬升程度,值越高越接近黑臭特征 BOI = "(float(b3) - float(b4)) / (float(b3) + float(b4) + 0.000001)" # 3) 水华识别指数:藻类近红外反射强,红波段吸收深 NDVI = "(float(b4) - float(b3)) / (float(b4) + float(b3) + 0.000001)"三行表达式共用一套映射关系,所以映射只在第一次做,后面两次可以直接复用。float 转换是为了避免整型相减产生的截断误差,分母加一个极小值是为了防止近红外与绿波段同时为 0 时出现除零,虽然在实际影像里极少发生,但一旦出现就是整幅图出现无效值,排查起来很费时间。
三个指数的取值范围都在 -1 到 1 之间。水体指数大于 0 通常判为水体,具体阈值要结合本景影像的直方图看双峰分布的位置;黑臭辅助指数的高值区要叠加水体掩膜才有意义,否则会把裸土也圈进来;水华识别指数的高值区同样需要限制在水体范围内,不然岸边植被会被误判。
运算结果保存为独立的栅格文件,命名里带上指数名称,不要用默认的运算结果名,否则后面做分级时根本分不清哪个文件对应哪个指数。
4.2 黑臭分级阈值与率定方法
黑臭程度不是一个可以直接测量的物理量,它是多个水质指标的综合表现,遥感能捕捉的是它的光谱表征。比较稳妥的做法是采用双指标约束:先用黑臭辅助指数切出候选区,再用可见光波段的整体反射率做二次筛选,两个条件同时满足才判为黑臭。
分级一般设三档:疑似黑臭、轻度黑臭、重度黑臭。初始阈值可以参考同类区域已发表的研究区间,但一定要在报告里注明这是初始值。有实测点位时,用点位坐标提取对应像元的指数值,做一次简单的阈值寻优,让分级结果和实测等级的一致率尽量高。没有实测数据时,退而求其次的做法是选取影像上已知问题河段作为验证样本,看分级结果是否落在合理位置。
分级结果建议转成整型分类栅格,用 1、2、3 表示三个等级,0 表示非黑臭水体。整型分类栅格在 ArcGIS 里符号化最方便,做面积统计也不会有小数误差。
4.3 水华识别与范围提取
水华识别的难点在于区分水面藻团和岸边植被。两者在归一化植被指数上都可能是高值,区别在于空间位置和形状:水华贴着水体内部成片或成条带分布,岸边植被则沿着水陆边界呈条带状。
处理方式是在水华识别指数高值区的基础上,再叠加一层距岸距离约束,把紧贴水陆边界的像元排除掉,剩下的才是水华候选区。这个约束可以用 ArcGIS 的缓冲区分析实现,也可以在 ENVI 里对水体掩膜做一次腐蚀运算近似实现。
水华范围提取后同样要做一次目视核对,把分类结果和真彩色影像叠在一起看,重点检查湖湾、入河口这些容易聚集藻类的区域是否被漏掉。漏检往往出现在薄云或者云影边缘,这类区域在预处理阶段就该标记出来,反演时直接排除,比事后修补省事。
4.4 ArcGIS 分级制图与面积统计
把分类栅格和指数栅格导入 ArcGIS 后,先统一坐标系,再叠加水域矢量边界和行政边界作为底图要素。分类栅格用分级设色符号化,三个黑臭等级用冷暖递进的色带区分,水华范围用醒目的强调色叠加,避免和黑臭等级撞色。
面积统计用栅格统计的方式最直接,按分类值统计像元数,乘以单个像元的面积(16 米乘 16 米,即 256 平方米)就得到各类面积。统计结果导出成表格,按水域分段或者按行政区汇总,这张表是报告里最有分量的成果之一。
# 按分类值统计像元数并换算面积(示意,可在 ArcGIS 字段计算或外部脚本中实现) pixel_area = 16 * 16 # 单像元面积,单位平方米 for cls, count in class_counts.items(): area_m2 = count * pixel_area # 该类别的总面积,单位平方米 area_km2 = area_m2 / 1_000_000 # 换算成平方公里,便于报告引用 print(cls, count, round(area_km2, 3))出图规范上有两个细节值得坚持:一是导出分辨率不低于 300 dpi,二是图面必须包含指北针、比例尺、图例和坐标系说明。评审看一张图的时间通常不超过半分钟,缺了这些要素,再准确的结果也会被打折扣。
5. 避坑与常见问题:五个最容易翻车的环节
这条链路走通一次不难,难的是换一景影像还能稳定复现。下面这五个问题在实操中出现的频率最高,每一个都对应过真实返工。
5.1 照搬短波红外指数导致计算失败
现象是波段运算执行后结果全为无效值,或者数值分布完全不符合水体的光谱规律。原因很直接:高分一号 WFV 只有四个可见光近红外波段,没有短波红外通道,凡是依赖短波红外的指数在这套数据上都无从计算。解决办法是回到传感器波段配置表,只用蓝、绿、红、近红外四个波段能表达的指数,水体识别用绿波段与近红外的组合,需要增强对比时可以把红波段也纳入构造比值型指数。选指数之前先查波段,这个习惯能省掉大量试错。
5.2 大气校正未挂载光谱响应文件
现象是大气校正跑完后,植被波谱曲线异常平滑,红波段吸收谷和近红外反射峰都不明显,水体反射率整体偏高。原因是 FLAASH 的传感器列表里没有高分一号预设,如果直接按默认传感器运行,算法用的波段响应范围和实际数据对不上。解决办法是在多光谱设置里手动载入对应相机编号的光谱响应函数文件,载入后重新运行。判断是否生效的快速办法是对比校正前后的同一像元波谱曲线,形状应该更接近典型地物特征。
5.3 正射校正后的黑边未剔除
现象是裁剪后的水域范围比矢量边界大出一圈,面积统计结果明显偏高,分类图上出现大片零值区域被判为某一等级。原因是正射校正旋转画布产生的无效值默认为 0,而 0 在波段运算里会被当成真实反射率参与计算。解决办法是在元数据里把无效数据值设为 0,让软件识别这些像元为无效;如果已经进入反演阶段才发现,可以在波段运算表达式里加一层判断,把四个波段同时为 0 的像元输出为无效值。这个坑的隐蔽性很强,因为出图时黑边和深色水体肉眼很难区分。
5.4 水体掩膜方向设反或坐标系不一致
现象是裁剪结果只剩陆地,或者水域边界和水体完全错位。原因通常有两个:应用掩膜时保留与剔除的方向设反了,或者矢量数据和影像的坐标系不一致,导致叠加时发生偏移。解决办法是先检查两者的坐标系定义,必要时对矢量做投影转换;掩膜方向可以在输出后用一个简单办法验证,统计结果的有效像元数应该和矢量面积除以像元面积得到的数量级接近,差一个数量级基本就是方向问题。
5.5 直接套用文献阈值
现象是分级结果要么全是同一等级,要么等级分布和影像上的实际情况明显不符。原因是不同区域、不同季节、不同大气条件下,同一指数的数值分布差异很大,文献里的阈值只对原文的研究区有效。解决办法是把阈值当作初始值,用本景影像的指数直方图找双峰之间的谷底位置作为参考,有实测点位时用点位提取值做一次率定。报告里要写清楚阈值的确定过程,这一步写扎实了,答辩时关于精度的问题会好回答很多。
6. 验收与答辩:用一条可复现的链路说话
这套方案的验收标准不适合按"图好不好看"来定,而应该按"换个人能不能复现"来定。我在整理交付材料时通常要求三样东西齐全:一份参数记录表、一套中间结果、一张对照图。
参数记录表按处理步骤逐条列出,每一项写清楚工具名称、关键参数取值和取值依据。比如辐射定标的输出类型为什么选浮点、大气校正的大气模型为什么选这一档、正射校正的重采样为什么用三次卷积。这张表看起来繁琐,但它是答辩时最有力的材料,因为大部分提问都集中在"为什么这么设"上,有据可查比临场解释可靠得多。
中间结果要保留辐射定标、大气校正、正射校正、掩膜裁剪四个节点各一份文件。它们的作用是出错时能快速定位问题出在哪一环,而不是从头重跑。评审如果对某个数值有疑问,可以顺着中间结果往回查,这种可追溯性本身就是方案质量的一部分。
对照图建议做两张:一张是反演分类结果叠在真彩色影像上,一张是分类结果叠在水域矢量边界上。前者用来看空间位置对不对,后者用来看范围有没有溢出。两张图放在一起,判别结果的可信度基本就呈现出来了。
精度验证部分,如果有实测点位,用混淆矩阵给出总体精度和黑臭等级的一致率;如果没有,至少给出抽样目视核查的结果,说明抽了多少个点、其中多少个判读一致。数字不必追求漂亮,把方法和局限写清楚,比给一个没有依据的高精度更有说服力。
答辩时被问得最多的问题往往是阈值怎么来的、精度怎么保证的、这套流程能不能推广到别的区域。前两个问题靠参数记录表和验证材料回答,第三个问题的答案是:把水体掩膜和阈值率定这两步替换掉,其余流程可以整段复用,这正是把流程做成链路而不是做成单次操作的价值所在。
从那以后我每次交付这类反演成果,都会强制自己先跑一遍参数记录表,再动笔写报告——顺序反过来的那次,我在答辩现场被追问大气模型的选择依据,翻遍文件也找不到当时的设置,只能含糊带过。希望这套结构能帮到你,少走一次这样的弯路。
本文还有配套的精品资源,点击获取