ARTICLE · INTELLIGENCE

战地情报 · 详情页

来自尧图项目组的一线实战观察与深度解析

遥感图像预处理全流程:辐射定标、大气校正与正射校正

遥感图像预处理全流程:辐射定标、大气校正与正射校正

刚接触遥感那会儿,我从数据平台下了一景多光谱影像,兴冲冲地拖进软件里一看,灰蒙蒙一片,地物轮廓依稀能辨,可数值完全不是那么回事——水体算出来的反射率比植被还高,两个时相的同一块农田颜色对不上,拼接的时候接缝像刀切的一样生硬。折腾大半天才搞明白,卫星下传的原始数据只是"半成品",传感器记录的不过是入瞳处的辐射能量,中间隔着大气散射、地形起伏、传感器响应差异、平台姿态偏差一大堆干扰项。遥感图像预处理要干的活,就是把这些干扰一层层剥掉,把像素值还原成能反映地表真实属性的物理量,再把每景影像放到正确的地理位置上。这套流程看着步骤不少,但顺序和参数都有讲究,踩过几次坑之后基本就形成固定套路了。下面我按实际处理顺序,把辐射定标、大气校正、几何校正、正射校正、图像融合、镶嵌裁剪这些环节拆开讲,每一步为什么做、参数怎么定、哪里容易出问题,尽量说透,给刚入门的朋友一条能直接照着走的路子。

1. 遥感预处理到底在解决什么问题

1.1 从传感器记录值到可用地表信息的落差

很多人第一次接触遥感数据,会默认"打开就能用",等到做植被指数或者分类才发现结果离谱。原因在于传感器上天以后,它测量的是到达镜头的那部分电磁波能量,这个能量在传播路径上已经被大气吸收散射改造过一遍,又被地形坡度朝向重新分配过一遍,最后还被传感器的响应曲线压缩进有限的量化位数里。换句话说,你拿到的DN值(Digital Number,数字量化值)是一个经过多重"滤镜"的混合产物,它既不是辐亮度,更不是反射率。

举个具体的例子,同一片针叶林在晴空和薄雾两种天气下成像,DN值能差出一大截,可林木本身的光谱特性并没有变。如果直接拿DN值做时序对比,结论必然是错的。预处理的第一步辐射定标,就是先把这个无量纲的整数还原成有物理单位的辐亮度;第二步大气校正,再把辐亮度里的气溶胶散射、水汽吸收剥离掉,得到地表反射率;几何相关的校正则负责把"图像坐标"变成"地理坐标"。这几件事的先后顺序不能乱,因为大气校正的输入必须是定标后的辐亮度,而正射校正又需要几何信息已经相对准确。

我见过有新手把大气校正放在几何校正之后做,结果重采样引入的混合像元让暗像元法估气溶胶的精度大打折扣。也有人图省事跳过辐射定标直接跑大气校正模块,软件虽然不报错,但出来的反射率量级明显不对。这些坑的共同点是:流程顺序的内在依赖关系没搞清楚,单看每个模块都能跑通,串起来就出问题。

1.2 一条标准处理流水线的结构设计

把上面这些依赖关系理清楚,一条通用的预处理流水线就出来了:辐射定标 → 大气校正 → 几何粗校正(如果数据自带几何信息这一步可省) → 正射校正 → 图像融合 → 镶嵌与裁剪。中间还会穿插一些辅助步骤,比如条带噪声去除、坏像元修复、云掩膜生成,这些属于"按需插入"的环节,不是每景数据都要做。

为什么把辐射相关的处理排在几何之前?除了前面说的依赖关系,还有一个实际考虑:大气校正模型通常假设像元是纯净的、未被重采样污染的,几何重采样会引入邻近像元混合,破坏这个假设。而几何校正放在后面,是因为它涉及重采样,会改变像元值,如果先做几何再做辐射校正,重采样后的像元值就不再是原始观测值,定标系数和大气参数就对不上了。

至于融合、镶嵌、裁剪,它们本质上是在"整理成品",放在最后能避免反复重采样。我个人的习惯是每一步都保留中间成果,定标后存一版、大气校正后存一版、正射后存一版,虽然占硬盘,但出问题时能快速定位是哪一环崩的,比从头重跑一遍省时间得多。

1.3 不同数据源的流程差异

流程框架是通用的,但不同数据源在细节上差别很大。中分辨率多光谱数据一般自带RPC(有理多项式系数)文件和定标系数,正射校正比较省心;高分辨率数据往往全色和多光谱分开下传,需要额外做融合;一些早期数据或者特殊载荷可能连几何粗校正都没做,得自己找控制点。

还有一类情况是数据本身已经做过部分处理,比如"正射校正产品""地表反射率产品",这时候你要先看产品说明,搞清楚它已经做了哪几步,避免重复处理。我就犯过把已经大气校正过的地表反射率产品再跑一遍大气校正的错误,结果整景影像的数值被压得极低,植被指数全乱。所以拿到数据的第一件事不是急着开软件,而是把元数据文件从头到尾读一遍,看清楚处理级别、定标参数、投影信息、云量这些关键字段。

2. 辐射定标:把DN值翻译成物理量

2.1 定标参数从哪来、怎么用

辐射定标的本质是一次线性变换,把整数DN映射成有物理意义的辐亮度。最常见的形式是:

L = gain × DN + offset

其中L是辐亮度,单位通常是 W/(m²·sr·μm),gain和offset就是定标系数,一般写在元数据文件里。有些数据给的是增益和偏移,有些给的是最大最小辐亮度配合量化范围,形式不同但本质一样:

L = (Lmax - Lmin) / (Qmax - Qmin) × (DN - Qmin) + Lmin

这里的Qmax、Qmin是量化位数对应的极值,比如12位数据就是4095和0。两种写法算出来的结果应该一致,如果不一致,多半是元数据里的系数版本对不上,或者你用的量化范围搞错了。

实操中最容易出问题的是"系数选错波段"和"单位没统一"。一景多光谱数据动辄七八个波段,每个波段的gain、offset都不同,手动填的时候错位一格,整条曲线就歪了。我的做法是写个小脚本把元数据解析成表,波段号和系数一一对应,批量处理时直接读表,不给人手出错的机会。单位方面,有些系数给出的辐亮度单位是 W/(m²·sr·nm),跟常用的 μm 差了一千倍,如果不换算,后续大气校正的输入量级就会离谱。

2.2 表观反射率的计算与太阳高度角修正

定标到辐亮度之后,很多应用其实还需要再往前一步,算表观反射率(也叫大气顶层反射率)。公式是:

ρ_TOA = π × L × d² / (ESUN × cosθ)

d是日地距离校正因子,ESUN是波段对应的太阳辐照度,θ是太阳天顶角。这一步的意义在于把不同时间、不同太阳角度下获取的辐亮度归一化,让数据之间具备可比性。

太阳高度角这个参数经常被忽略。一景影像成像时太阳高度角可能只有三四十度,如果不做 cosθ 修正,算出来的反射率会系统性偏高,而且高纬度地区偏差更明显。有些软件在定标模块里会自动读取元数据里的太阳角度并完成修正,有些则需要你手动确认,别默认它一定做了。

提示:定标完成后,挑几个典型地物(深水体、浓密植被、裸沙)看看辐亮度或反射率的数值范围是否合理。深水体在可见光波段的反射率通常很低,如果算出来比植被还高,基本可以断定定标环节出了问题。

2.3 定标环节的几个实操细节

第一,注意数据的量化位数和是否有符号。有些传感器用有符号整型存储,直接按无符号读会把负值读成很大的正数,定标结果面目全非。第二,部分数据存在"填充值"(比如0或者-9999表示无效像元),定标前要把这些值识别出来,否则它们参与统计会拉偏直方图。第三,定标系数有时会随传感器状态更新,同一颗卫星不同时期的数据系数可能不同,批量处理前务必逐景核对元数据,不要拿一套系数套所有数据。

我个人的习惯是定标后立刻做一次快速统计,输出每个波段的最小值、最大值、均值、标准差,存成一张表。这张表在后续排查问题时非常有用——如果某景数据某个波段的均值明显偏离同期其他景,要么是成像条件特殊,要么是定标系数错了,早发现早处理。

3. 大气校正:去掉那层"雾"

3.1 大气校正到底在修正什么

大气对成像的影响可以粗略分成两部分:散射和吸收。散射主要由气溶胶和空气分子引起,它给影像叠加了一层"亮度底噪",让暗地物变亮、对比度下降,这就是为什么原图看着发灰。吸收则主要由水汽、臭氧、二氧化碳等气体造成,在特定波段形成吸收谷,削弱信号。

大气校正的目标就是把这两部分影响从辐亮度里剥离,反演出地表反射率。听起来简单,做起来难,因为大气状态(气溶胶光学厚度、水汽含量)在成像瞬间并没有被完整测量,只能靠模型反演或者辅助数据估计。

常见的思路有三类:辐射传输模型法(如6S、MODTRAN这类模型)、暗像元法、以及基于影像本身的经验线性法。辐射传输模型精度高但需要输入大气参数,暗像元法假设影像中存在反射率很低的暗地物(深水体、浓密植被阴影),用它来估计气溶胶贡献,实现简单但依赖场景假设。

3.2 方法选型与实操参数设置

选哪种方法,取决于你的数据、精度要求和手头有多少辅助信息。做定量遥感,比如反演叶绿素含量、监测水质,那必须上辐射传输模型,而且尽量用同步的大气参数;如果只是做分类、变化检测这类定性或半定量应用,暗像元法往往够用,速度快、依赖少。

以辐射传输模型法为例,实操中需要关注几个参数:

参数作用获取方式
气溶胶光学厚度决定散射强度影像反演或外部数据插值
水汽含量影响近红外吸收影像反演或再分析资料
大气模式描述温压廓线按纬度和季节选标准模式
传感器光谱响应匹配波段随数据附带的响应函数文件

气溶胶光学厚度是误差最大的来源。没有同步观测时,很多软件会用影像自身的暗像元反演,这时候"暗像元"选得准不准直接决定结果。城市区域找深水体相对容易,干旱区就麻烦,有时整景找不到足够暗的目标,反演结果会偏。

还有一个容易踩的坑是光谱响应函数。不同传感器的波段响应曲线不同,如果大气校正时用了默认的通用响应函数,而不是该传感器真实的响应函数,在吸收波段附近会有明显偏差。这个细节在入门资料里很少提,但做定量分析时影响不小。

3.3 校正效果的快速验证

大气校正做完,怎么判断靠不靠谱?我常用三个土办法。一是看植被光谱曲线,健康植被应该在红光波段有吸收谷、近红外有反射峰,如果校正后曲线还是平的一坨,说明没校好。二是看深水体的反射率,清洁水体在近红外波段应该接近零,如果数值明显偏高,说明气溶胶去除不干净。三是和已知地物的实测光谱对比,哪怕是粗略对比,也能发现量级上的错误。

注意:大气校正后的反射率理论上应该落在0到1之间,出现负值是正常现象(模型反演的数值噪声导致),但如果大面积出现负值或者超过1的值,说明输入参数或者量纲有问题,别急着往下走。

4. 几何校正与正射校正:让像素归位

4.1 几何粗校正与精校正的分工

几何校正解决的是"影像上的像素位置和真实地理坐标对不上"的问题。造成偏差的原因包括平台姿态抖动、地球曲率、地形起伏、传感器扫描方式等。处理上一般分两步:粗校正和精校正。

粗校正利用数据自带的轨道参数、姿态参数做系统性偏差改正,很多数据产品在分发前已经做过这一步,你拿到手时影像大致已经落在正确的地理范围内,只是还存在局部残余误差。精校正则要用地面控制点(GCP)或者参考影像来消除这些残余误差,常用的模型有多项式、有理函数模型等。

多项式校正适合地形平坦、范围不大的场景,一二阶多项式就能压得比较准;地形起伏大的区域,多项式就不够用了,因为它是平面模型,处理不了高程引起的位移。这时候就得上正射校正。

4.2 正射校正与DEM的匹配问题

正射校正的核心是引入高程信息,把每个像元投影到真实地表位置。它需要两样东西:传感器几何模型(RPC参数或者严格物理模型)和高程数据(DEM)。

DEM的选择直接影响校正精度。用低分辨率DEM去校正高分辨率影像,山区的位移改正会不够精细;用过高分辨率DEM又可能引入不必要的噪声,而且处理速度慢。经验上,DEM分辨率与影像分辨率保持在同一量级或略粗即可,比如米级影像配十几米到三十米的DEM,基本够用。

另一个常被忽视的点是DEM的基准面。有些DEM是大地高,有些是正高,如果和传感器模型要求的基准不一致,校正后会有系统性高程偏移,在平地上看不出来,到山区就很明显。我遇到过整景影像在山区整体"漂移"几十米的情况,查了半天才发现是DEM基准没对齐。

RPC参数的质量也值得关注。多数商业数据附带的RPC是经过检校的,精度较好;一些免费数据或者早期数据的RPC可能需要用户自己优化。如果发现正射后影像在部分区域仍有明显偏移,可以考虑用少量控制点做RPC优化,再重新校正。

4.3 重采样方法的选择逻辑

几何校正和正射校正都涉及重采样,也就是把原始像元值重新分配到新的网格上。常用方法有最近邻、双线性、三次卷积。最近邻不改变原始像元值,适合分类图或者需要保持光谱纯净度的场景;双线性平滑但会轻微模糊;三次卷积在保持细节和光谱之间平衡较好,是多数定量应用的默认选择。

这里有个取舍:重采样越"平滑",光谱保真度往往越差。做光谱分析时我倾向最近邻或者双线性,尽量避免三次卷积带来的像元混合;做目视解译或者出图,三次卷积的视觉效果更好。没有绝对的对错,看你的下游应用要什么。

5. 图像融合、镶嵌与裁剪:把数据整理成可用成品

5.1 全色锐化的算法选择

高分辨率数据常常是全色波段分辨率高、多光谱波段分辨率低,融合的目的就是把全色的空间细节注入多光谱,得到既有高分辨率又有颜色信息的影像。常用算法有Brovey变换、IHS变换、主成分法、Gram-Schmidt法、以及各种基于调制传递函数的锐化方法。

Brovey简单快速,但容易让颜色失真,植被区域偏色明显;IHS在色彩保持上比Brovey好一些,但对光谱的扭曲依然存在;Gram-Schmidt和基于MTF的方法在光谱保真度上表现更好,代价是计算量大、参数多。做定量应用时我更倾向光谱保真度高的方法,宁可牺牲一点锐度;纯粹做底图展示,Brovey这类轻量方法完全够用。

融合效果的评价不能只看"清不清晰",还要看光谱有没有被改坏。一个简单的检验方法是:融合前后分别算同一区域的植被指数,如果两者差异很大,说明融合把光谱扭曲了,需要换算法或者调参数。

5.2 镶嵌的接缝与辐射归一化

多景影像拼接成一整幅时,最大的问题是接缝。相邻影像成像时间不同、太阳角度不同、大气条件不同,即使都做过大气校正,接缝处仍可能有明显的亮度差。处理办法一般分两步:先做相对辐射归一化,把各景的直方图往同一基准上靠;再在接缝处做羽化或者加权过渡。

相对辐射归一化常用直方图匹配,选一景作为参考,其他景向它靠拢。参考景的选择有讲究,最好选成像条件居中、云量少、覆盖地物类型全的一景,否则匹配后的结果会整体偏亮或偏暗。羽化过渡的宽度也要合适,太窄接缝还在,太宽会把两侧地物的真实差异抹平。

5.3 裁剪与掩膜的实际用途

裁剪看起来是最没技术含量的一步,但用好了能省很多事。按研究区矢量边界裁剪,可以去掉大量无关区域,减少后续计算量;按云掩膜裁剪,可以剔除云污染像元;按坡度、高程掩膜,可以只保留感兴趣的地形范围。

需要注意裁剪时的像元对齐问题。如果裁剪边界和影像网格不对齐,软件会做重采样,引入额外的值改变。我的做法是尽量让裁剪边界落在整数像元上,或者干脆先裁剪再校正,避免多次重采样叠加。

6. 常见问题与排查实录

6.1 典型问题速查表

预处理过程中遇到的问题,大部分有规律可循,我把踩过的和见过的整理成一张表:

现象可能原因排查方向
影像整体发灰、对比度低未做大气校正或校正不足检查气溶胶参数、暗像元选取
反射率量级明显偏差定标系数用错或单位未换算核对元数据、检查量纲
校正后山区影像整体偏移DEM基准不一致或分辨率过低检查DEM高程基准、更换DEM
拼接处出现明显亮暗条带未做相对辐射归一化加入直方图匹配、调整羽化宽度
融合后植被颜色异常融合算法光谱保真度差换用GS或MTF方法、调整参数
时序分析同一地物数值波动大太阳角度未修正、大气校正不稳定检查太阳高度角、统一处理参数

这张表不是穷举,但覆盖了入门阶段八成以上的困惑。遇到问题时先对号入座,能快速缩小排查范围。

6.2 几条踩坑心得

第一条,别迷信默认参数。软件里的默认大气模式、气溶胶类型都是"通用值",套到你的研究区未必合适。花十分钟看看研究区的气候类型、海拔、季节,手动调一下参数,结果可能差很多。

第二条,保留中间产品和处理日志。每一步的参数、输入输出文件、软件版本都记下来,出问题时能复现,跟别人交流时也能说清楚。我吃过没记录的亏,同一景数据两次处理结果不一样,查了半天才发现是某次手滑改了个参数没记住。

第三条,小范围先试跑。整景数据处理动辄几十分钟,参数错了重来很浪费时间。先裁一小块做试验,确认流程和参数没问题再上全图,这个习惯能省下大量等待时间。

第四条,光谱曲线是最诚实的裁判。不管中间用了什么高级算法,最后拿典型地物的光谱曲线和已知标准比一比,曲线形状对不对、吸收峰位置准不准,一眼就能看出问题。这个土办法比任何精度报告都直观。

6.3 关于流程取舍的个人体会

最后说一个观念上的问题:预处理不是步骤越多越好。每一步处理都会引入新的误差,也会消耗时间。如果下游应用只是做目视解译或者粗略分类,辐射定标加几何校正可能就够了,大气校正和融合未必需要;如果是做定量反演,那该做的步骤一个都不能少,而且参数要抠得细。

我现在的习惯是先明确下游任务需要什么级别的数据,再倒推该做哪几步。数据级别够用就停手,别为了"流程完整"而多做无用功。遥感预处理这门手艺,真正的门槛不在会点哪些按钮,而在于知道哪一步可以省、哪一步必须做、哪一步做过头了反而坏事。这个判断力,只能靠一次次实际项目慢慢磨出来。

RELATED READING

延伸阅读

更多一线实战笔记与深度复盘,助您持续精进