1. 高分遥感数据预处理到底在解决什么问题
做高分遥感这一行,最磨人的往往不是算法调参,而是影像预处理。我接手的项目里,十有八九的原始数据都不能直接用:辐射值漂移、大气雾霾、几何偏移、条带噪声、波段错位……这些干扰要是不在预处理阶段处理干净,后面无论是目视解译、地物分类还是变化检测,结果都会走样。今天想把这些年踩过的坑整理下来,重点说说高分遥感数据预处理里的核心环节、参数选择和实操细节,给正要入门的同学一份能直接照着做的清单。
1.1 原始高分影像“脏”在哪里
高分影像听起来高大上,但传感器记录下来的原始数据,本质上只是“亮度值”的集合。以常见的多光谱影像为例,卫星在几百公里外的轨道上接收的是地物反射的太阳辐射,但因为太阳高度角、大气散射、传感器响应差异、地形起伏等因素,传感器输出的DN值(Digital Number)和地物真实的反射率之间,隔着一大堆误差。说到底,原始影像是一个“被污染的信号”,预处理的本质就是尽可能把污染剥离掉,还原地物真实的辐射特征和几何位置。
我们实际碰到的“脏”,大概可以分成四类。第一类是辐射失真,包括传感器各波段响应不一致带来的条带、噪声,以及太阳高度角变化带来的光照差异。第二类是大气影响,薄云、雾霾、气溶胶都会让影像灰蒙蒙的,地物边界模糊,反射率被“抬”高或者“压”低。第三类是几何变形,卫星姿态、轨道高度、地形起伏和地球自转都会造成像素位置偏移,高山区或者影像边缘尤其明显。第四类是波段之间的空间错位,多光谱相机不同波段探测器存在一定的安装夹角和扫描时差,同一地物在不同波段上可能会差出几个像素。这些“脏数据”如果不修,就算你反演算法再先进,输出结果也只是在噪声里做文章。
1.2 预处理结果对下游应用的影响
很多新手对预处理不够重视,总觉得这是“体力活”,真正考验功力的是后面的分类、检测、反演。但事实恰恰相反:预处理结果直接决定了整个应用链路的稳定性和精度。我举个很直观的例子,如果你要做建筑物轮廓提取,原始影像边缘地带因为几何误差,楼顶和道路的边界会错位几个像素,算法提取出来的轮廓就会“毛刺”非常重,后处理要花大量时间去平滑。再比如做植被覆盖度反演,如果大气校正没做好,红光波段的反射率被气溶胶干扰,计算出的植被指数会明显偏高,同一个区域、同一期影像,不同时间处理出来的结果可能都对不上。
从工程角度看,预处理还是在给项目“定基线”。我参与过好几个需要多期影像对比的变化检测项目,第一件事就是把不同期次的影像一一预处理到统一的坐标系、统一的辐射基准下。坐标基准不一致,两期影像的同一地物可能偏移好几个像素,变化检测的“伪变化”一大堆;辐射基准不一致,一个月的NDVI变化可能不是植被真实变化,而是两次大气校正参数不同造成的伪信号。所以预处理不是可有可无的“清灰”,它是整个高分遥感应用里最底层、最不能省的一环。
2. 预处理核心环节拆解与技术选型
高分影像预处理的标准化流程,业内基本有一套比较成熟的套路:辐射定标、大气校正、几何正射校正、融合、镶嵌、裁剪。每一步都有自己的目标,也都有不少细节需要抠。我按照实际执行顺序,把每个环节的原理和参数选择拆开讲清楚,这样你带着自己的数据走一遍流程,起码不会稀里糊涂。
2.1 辐射定标:把DN值变成物理量
第一步是辐射定标。传感器记录的DN值只是一个无量纲的整数,它和入瞳处的辐亮度之间是线性关系,通常表示为 L = Gain × DN + Offset。标定参数一般放在影像自带的元数据文件里,也可能是卫星地面系统定期更新的绝对定标系数。这一步的意义在于让不同时间、不同传感器获取的影像,在辐射尺度上有可比性。如果你只做单景目视解译,不定标也能看;但只要你做定量反演,或者多期影像对比,这步一定不能跳。
实际操作中,很多平台都提供了定标工具,比如常见遥感软件里的“Radiometric Calibration”,你只需要选择输入影像,软件会自动读取元数据里的定标参数。需要注意两点:一是定标后的单位,有的是辐亮度(W·m⁻²·sr⁻¹·μm⁻¹),有的是反射率,你要根据后续需求选择,如果后面要做大气校正,一般先定标到辐亮度;二是少数影像的元数据定标系数更新滞后,如果发现定标结果和已知地物反射率差得离谱,要去官方渠道核对最新系数。另外,热红外波段如果有,定的量纲和可见光不一样,不要用同一套流程盲目套用。
2.2 大气校正:去掉空气这层“毛玻璃”
大气校正是所有环节里概念最绕、参数最多的一个。通俗地说,大气就像一层“毛玻璃”,太阳光经过大气照射到地物,再反射回传感器,中间经过了吸收和散射,传感器看到的是“地物信号 + 大气信号 + 邻近地物散射信号”的混合体。大气校正的目的,就是把这层“毛玻璃”的影响去掉,获取地表真实反射率。
常见的大气校正方法分三类:统计方法、经验方法和辐射传输模型方法。统计方法比如内部平均法、对数残差法,适合没有气象数据的特殊情况,精度一般,主要用于矿物填图这类特定场景;经验方法如暗像元法,通过找影像中最暗的像元(通常假设反射率接近0)来估算大气路径辐射,操作简单,但对地物分布有严格要求,大面积水体多的时候效果会好一些;辐射传输模型方法如6S、MODTRAN以及专业遥感软件里的大气校正模块(如FLAASH、ATCOR),模拟大气对辐射的散射和吸收过程,输入参数包括太阳天顶角、能见度/气溶胶类型、大气水汽含量、传感器高度等,精度最高,也是高分影像预处理的主流选择。
我自己的习惯是:4波段多光谱影像,优先用模型法;如果现场实在没有大气参数,可以尝试暗像元法做快速校正。但要注意,模型法不是参数越多越好,反而容易“过拟合”。很多时候,默认的气溶胶模型(乡村型/城市型)、默认的初始能见度(比如40公里),就能把结果修得不错;真正影响最大的是你选择的波段中心波长和传感器响应函数,这个必须和影像元数据对应准确。如果软件支持自动从元数据读取,就尽量用自动读取,不要手动填错。
2.3 几何校正与正射校正:让像素落到正确的位置
几何校正解决的问题是“像素在哪”。高分影像原始数据带有RPC(Rational Polynomial Coefficients)文件,就是用来描述影像像元与地面三维坐标之间关系的数学模型。正射校正会结合RPC和数字高程模型(DEM),把影像从中心投影转换成近似正射投影,消除地形起伏引起的位移误差。对于平原地区,平坦地形下RPC正射校正效果就够;但对于山区、丘陵地带,没有DEM参与的正射校正,山脊和山谷的错位会非常严重,和矢量边界套不上。
实操中,几何校正分为“有控制点”和“无控制点”两条路线。无控制点的情况下,直接利用RPC和DEM做正射校正,精度大约在几个像素到十几像素,适合对精度要求不高的底图生产。有控制点的情况下,通过人工选点或者自动匹配控制点,再进行三角网纠正或多项式纠正,可以把误差控制到1~2个像素内。这里有一个容易被忽视的点:控制点不一定要多,但分布一定要均匀,尤其是影像边缘和地形起伏大的区域,控制点稀疏容易造成局部扭曲。我还遇到过因为DEM精度太差导致正射影像“糊掉”的情况,后来换了更高分辨率的DEM数据,问题才解决。所以不要迷信算法,基础数据的精度同样重要。
2.4 影像融合:兼顾分辨率与光谱信息
高分影像传感器通常同时提供一个全色波段(高空间分辨率,比如0.5米)和多个多光谱波段(低空间分辨率,比如2米)。融合的目的就是把全色的纹理细节和多光谱的色彩信息合在一起,生成高频细节丰富且光谱不失真的高分辨率多光谱影像。方法很多,常见的有Brovey变换、IHS变换、主成分分析(PCA)融合、Gram-Schmidt(GS)融合、Pansharpening等。早期经典方法Brovey和IHS算法简单、速度快,但容易产生光谱畸变,尤其在植被和水体区域,颜色会变得“发飘”;GS融合因为能更好地保持光谱信息,目前是主流选择。
我说一个实用判断:融合后影像的红、绿、蓝波段直方图,和原始多光谱影像直方图分布应该大致接近,不能出现明显的拉伸或者偏移。如果影像颜色出现“假彩”或者地物光谱曲线明显变形,多半是融合算法和参数没选对。另外,融合前必须保证全色和多光谱影像已经完成精确配准,如果两个影像错位超过半个像素,融合后会出现重影和边缘发亮。很多新手为了追求细节,喜欢把多光谱波段也做了锐化增强,结果造成彩色噪点,实际上融合后适量锐化就够了。
2.5 镶嵌与匀色:拼出干净的大范围影像
当你处理的不是一个景,而是覆盖某个区域的多景影像时,镶嵌和匀色就成了最麻烦的环节。镶嵌不只是简单“拼图”,它要解决三个问题:几何接边、色彩一致和羽化过渡。几何接边首先要求所有待镶嵌影像都完成正射校正,并在同一个投影坐标系下;色彩一致要求不同期次、不同侧摆角度下获取的影像,经过辐射定标和大气校正后,色调差异尽可能小。但即使做了统一预处理,因为太阳高度角、雾霾残留的影响,相邻影像之间仍可能存在肉眼可见的色差,这时就需要匀色处理。
匀色的思路主要有两种:一是直方图匹配,选定一幅参考影像,把其他影像的直方图统计特征匹配到参考影像上,速度快但可能会压掉部分真实信息;二是在重叠区计算线性回归或整体色彩校正系数,让两景影像的均值、方差统一。我建议如果重叠区覆盖面积大且包含丰富地物类型,优先用回归法匀色,这样能更好保持地物光谱的相对关系。镶嵌线上,“强制羽化”不是越宽越好,羽化太宽会在复杂地物上产生模糊重影,一般3~5个像素过渡就够,必要时把镶嵌线尽量避开房屋、山脊这种结构明显的区域。
2.6 裁剪与输出:定义产品边界
最后一步是裁剪与输出。很多项目只关心“把需要的地方裁出来”就行,但输出参数的规范往往决定了产品能否直接进入业务流程。首先是坐标系:一定要统一到项目要求的投影坐标系,比如UTM分带或者地方坐标系,不要软件默认是什么就导什么。其次是保存格式:一般用GeoTIFF,能带坐标信息和金字塔,方便后续发布和浏览;如果用于深度学习样本制作,还要注意位深和波段顺序。再有就是输出范围:边界外扩或者内缩多少、是否生成掩膜文件、是否包含无效值标记,这些都需要在项目前期定义清楚。
裁剪时有个容易踩的坑:使用ArcGIS/QGIS或者GDAL裁剪的时候,如果不设置“NoData”值,影像边缘会产生一圈黑色或者白色的无效值区域,叠加到底图上会很难看。所以裁剪后一定要检查无效值范围,把它显式标记为NoData或者设置成固定的透明色。另外,如果输出影像超过2~3GB,建议将金字塔文件和TFW世界文件一并生成,否则后期发布成地图服务的速度会很慢。
3. 实操记录:一景高分影像的完整预处理流程
光讲理论容易空,我带大家走一遍我最近处理的某区域高分多光谱影像的完整流程。这景影像是我从数据供应商那里拿到的原始L1级产品,包含一个全色波段和四个多光谱波段(蓝、绿、红、近红外),覆盖区域有平原、丘陵和一片城区,属于比较典型的“难啃”数据。整个流程我在常用遥感软件里完成了大约四十分钟,下面会把关键步骤和参数都写出来。
3.1 前期检查与数据准备
拿到数据以后,我做的第一件事不是急着开软件跑流程,而是先做“体检”。打开影像文件夹,确认包含了影像文件、RPC文件、元数据XML文件和影像质量说明。用软件打开元数据,查看成像日期、天底角、太阳高度角、传感器名称、波段中心波长等关键字段。同时快速浏览影像整体亮度和云量,估算一下无效边界范围。
特别要检查三样东西:一是波段顺序和数量,不同产品批次的波段排序可能不一样,我见过把红波段和近红外波段顺序搞反的训练样本,整个项目差点废掉;二是数据位深,16位还是8位,直接影响到后续直方图拉伸和保存格式;三是云量和云雪掩膜情况,如果覆盖区有大量厚云,最好直接放弃该景或者计划后续补景,预处理做得再漂亮也无法恢复被云遮挡的地表信息。完成这些检查后,我会在影像上给云区域画一个粗略的矢量范围,作为后续质检时的参考。
3.2 辐射定标与大气校正的参数设置
在软件中打开原始影像,先执行辐射定标。这里的关键参数是“定标类型”,我选择了辐亮度(Radiance),确认软件读取的定标增益和偏移系数正确,并勾选了“输出波段顺序与输入一致”。定标完成的数据会输出为一个新的IMG或者TIFF文件,理论上不需要再做其他处理。
紧接着做大气校正,我用的是模型法。参数设置如下:传感器类型选择“高分辨率多光谱”;成像日期和经纬度由元数据自动读取,只需要核对;飞行高度填写卫星轨道高度;地面平均海拔可以填覆盖区的平均DEM高程,我这次填了约300米;大气模型选“中纬度夏季”;气溶胶模型选“乡村型”,初始能见度填了40公里;水汽反演用“自动”,软件会根据近红外波段自动估算水汽柱含量。需要提醒的是,如果影像覆盖的是城市建成区,气溶胶模型选“城市型”可能更合适;如果影像日期是冬季或高纬度,大气模型要改成“中纬度冬季”或“亚极地冬季”。这些看起来很小的参数,对最终反射率的影响可达到5%~10%。
校正完成后,我会快速抽几个典型地物的光谱曲线出来看:清洁水体近红外反射率应该很低,植被的绿波段反射率应该高于蓝波段,裸土的红波段和近红外波段比较接近。如果这些现象都不符合,就要回头检查参数。
3.3 正射校正与控制点配准
大气校正完成后,下一步是正射校正。我使用的软件支持RPC正交化,选择输入影像,加载RPC文件,再指定DEM。这里我选了一个分辨率约为30米的开源DEM数据,因为覆盖区地形起伏不大,足够用了;如果是陡峭山区,最好找12米或者更高精度的DEM。输出投影我设定为UTM 50N、WGS84坐标系,像素大小保持原始分辨率(全色0.8米,多光谱3.2米),重采样方法选了三次卷积(Cubic Convolution),因为它能保留更好的纹理,但抗噪性能略差,如果原始影像噪点较多,改用双线性内插更安全。
正射校正的中间结果,我会把它和已有的矢量参考数据(比如过往的高精度影像或地形图)叠加,检查是否有系统性偏移。如果偏差超过三五个像素,就需要导入控制点做二次几何精纠正。我这次碰到的情况是影像整体向西北偏移了约4个像素,原因可能是原始RPC精度不足。我在影像的四个角落和中部均匀选了8个控制点,用参考影像上清晰的道路交叉口进行匹配,然后使用二次多项式纠正。纠正后的检查点误差控制在1个像素以内,效果比较理想。注意控制点最好不在云区、水域和阴影区选,这些地方特征不稳定,容易引入误差。
3.4 融合、镶嵌、裁剪的落地操作
正射校正完成之后,就需要把全色影像和多光谱影像融合起来。我选的是GS融合方法,因为它在保持光谱信息方面表现稳定。操作上,第一步将全色和多光谱影像导入同一个文件,软件会自动识别全色作为高分辨率输入;第二步设置融合参数,选择“从多光谱统计量生成模拟全色波段”,传感器类型保持默认;第三步输出前检查输出像素大小,我设成全色分辨率0.8米,这样融合后的影像同时具备0.8米的空间分辨率和4波段的光谱信息。融合完成之后,我习惯再做一个快速色彩对比,把融合影像和原始多光谱影像同时打开,在植被、水体、裸露地三处分别观察颜色变化。我这次融合出来后近红外波段的水体数值比原始多光谱略微抬高了一点,但还在可接受范围内,就没有返工。
然后是多景镶嵌。我手头一共三景影像,重叠率在10%~25%之间。先通过软件自动生成镶嵌线,再把明显穿过房屋和山脊、容易产生错位的镶嵌线手动调整到主干道或者农田边界上。匀色选择以中间那景影像为参考,使用直方图匹配。由于三景影像成像日期相差不到一个月,多数地物光谱差别不大,匀色后接边处基本看不出明显色差。最后按项目研究区范围做矩形裁剪,设置输出为GeoTIFF格式,位深16位,NoData值设为0,同时生成金字塔和世界文件,这样数据就算正式出活了。
3.5 结果质量检查与精度评价
预处理做完不代表万事大吉,质量检查是最后一道防线。我会从三个维度抽检。第一是几何精度,把预处理成果叠加到参考影像上,随机选20个左右的检查点统计平面位置误差。我这次的结果是平均误差在1.2个像素左右,最大误差2.8个像素,主要出现在丘陵阴影区,基本符合1:1万制图精度要求。第二是辐射精度,抽取光谱曲线和统计直方图,确认地物光谱形态合理,四个波段的最小值、最大值、均值没有出现极端异常。第三是目视效果,打开全色锐化后的融合影像,检查是否有融合重影、拉花、色彩断裂、黑边、云雾残留等问题。如果发现局部重影,需要用原始波段重新检查配准精度,定位问题出在融合前还是融合后。
另外,我会保留中间过程的全部文件,包括定标后影像、大气校正影像、正射校正影像、融合影像、镶嵌和裁剪结果,每一个都配上简单的说明文档。这样做的好处是,一旦发现问题,能确定是哪个环节造成的,不用全部从头跑一遍。而且项目复盘时,别人能按着你的产物链条快速跟上思路,在团队协作里特别重要。
4. 预处理中的高频问题与排查思路
即便流程走得很熟练,预处理过程中还是会碰到各种奇奇怪怪的问题。这里我把这几年高频出现的坑集中排一下,并给出排查路径和解法,你可以直接当速查表用。
4.1 影像偏色和去云不干净怎么办
影像偏色是最高频的问题之一。大气校正后整个影像“红通通”或者“蓝汪汪”,首先检查波段顺序是否错位,很多人把红蓝通道弄反了都不知道;其次检查大气模型和气溶胶模型选得是否合适,比如冬季影像用夏季大气模型,蓝波段吸收就会被高估,导致影像偏红。还有一个常见原因是传感器响应函数没有正确加载,软件在校正时用了一个近似光谱响应,和真实传感器差别较大,也会明显偏色。
云和雾处理不干净,常见于暗像元法和一些快速大气校正流程。我的建议是:如果覆盖区云影很重,不要期望大气校正能“去云”。大气校正只能减轻霾和薄云的影响,对厚云无能为力。薄云残留明显的区域,可以尝试多时相云检测和去云算法,或者干脆在后续应用阶段把云区掩膜掉。如果只是轻微雾感,可以降低气溶胶光学厚度(AOD)或者提高初始能见度试试,有些人习惯把能见度填到20公里以下,结果把整景影像整体抬得发白。
4.2 拼接缝明显与匀色失败
多景影像拼接出来后接缝两侧颜色泾渭分明,是最让人头疼的事。匀色失败的原因,首先要看是系统性问题还是局部问题。如果整景影像都和参考影像色调不一致,说明辐射定标或大气校正的输入参数存在差异,可能是两景影像的成像时间差太大、太阳高度角差异大,也可能是大气校正参数设置不一致。解决办法是重新对问题影像做一次大气校正,参数与参考影像保持一致,而不是在镶嵌阶段强行匀色。
如果只是某些地物区域颜色断层,比如农田区和城区过渡处特别明显,多半是匀色算法对异质地物不太适用。这时可以考虑把影像按地物类型分区域匀色,或者改用回归法动态计算重叠区校正系数。如果所有方法都不理想,还有一个“笨办法”:在镶嵌线上多拉几个手动调节点,把颜色过渡做得尽量柔和。虽然耗时,但有时候反而比自动匀色更自然。
4.3 几何定位偏差大怎么追查
正射校正后定位偏差较大,可能的原因很多。首先是DEM精度不够,山区尤其明显,DEM如果有几十米的垂直误差,正射影像的水平位移可能达到十几个像素;其次是RPC文件本身精度不够,这常见于一些宽幅影像边缘;再有就是投影坐标系选择错误,比如目标投影是UTM 49N,你输出成了UTM 50N,导致整景偏移非常夸张。
排查顺序建议先看整体偏移量是不是恒定常数,如果是,先检查投影坐标系和基准面是否一致。如果偏移方向和大小随地形变化,重点怀疑DEM。如果偏移在影像边缘特别大,中心区域正常,重点怀疑RPC模型。解决方法是先提高DEM分辨率,再做控制点精纠正,一般都能把误差降到可接受范围内。注意控制点选取时,野外的裸地或者水体边角往往GPS信号也不稳,优先选道路十字口、桥头、独立建筑物角点这类地物。
4.4 内存爆炸与批处理速度慢的优化思路
预处理大范围高分影像时,我最常被问到“为什么软件一处理大文件就卡死”。高分影像动辄几个GB,全色和多光谱叠加之后内存开销很大。几个实用的优化办法:一是分块处理,不要一次性把整景影像都放进内存,把影像按网格切块,逐块执行定标、大气校正和几何校正,最后再镶嵌回去,这是工业级产品生产线的常见做法;二是降低中间产物位深,如果不需要16位动态范围,可以保存为8位,能省一半空间和大量I/O时间;三是在融合和正射校正时提前建立影像金字塔,软件读取会更高效;四是不要同时打开多个大影像窗口,能关就关,肉眼看到的画面也会缓存大量资源。
另外,如果项目量大且重复性强,我建议用脚本或命令行批处理代替GUI点击。比如用GDAL系列的gdalwarp、gdal_merge.py 处理几何和镶嵌步骤,写一个循环批处理几十景影像,节省的时间非常可观。当然,脚本里的参数、坐标系、输出路径一定要配置清楚,我见过脚本人为把EPSG代码写错,一大批影像全部要重新处理的教训。
5. 预处理的一些个人经验
文章写到这里,聊点不太会在正式文档里出现的个人经验。第一,永远不要相信“原始数据已经做过初步处理,可以直接用”这句话。不同渠道拿到的产品级别不同,哪怕标称同一个级别,实际质量也可能天差地别。我养成了一个习惯:上手先做一次覆盖全区域的快速浏览,有异常早发现、早反馈,不要等预处理做完才发现整景数据废了。
第二,预处理过程中的每个参数都要记录。不是光记一个“大气校正完成”就完了,要把大气模型、气溶胶类型、能见度、DEM来源、控制点数量、误差值、重采样方法全部记录下来。后期如果结果被别人质疑,或者需要复现,这份元数据日志就是你最有力的凭证。我甚至见过一个项目因为处理日志缺失,导致整批成果被审计推翻、重新返工,那种痛苦大家务必提前避免。
第三,可以和算法同学多聊两句。有时候模型效果不好,你以为是自己分类网络调参的问题,实际上可能是预处理时用了过强的锐化和融合,导致影像纹理失真。预处理不一定是“做得越多越好”,不同的下游算法对预处理的敏感度不一样,最好在项目初期就和下游约定好影像产品的规范。比如做深度学习目标检测的,可能只需要辐射定标和几何校正,不需要过度的大气校正;做定量遥感的,大气校正则是核心步骤,一步都省不得。
最后再分享一个小技巧:正式批量处理前,一定拿一小块影像试跑完整流程,确认所有参数、输出格式、命名规则都没有问题,再铺开到全区域。我早年因为偷懒,直接对几十景影像跑了全流程,结果发现融合算法参数设置错了,只能全部重来,白白浪费了两天时间。现在不管多急,我都会先花二十分钟试跑一景,这个习惯帮我拦下了不少低级错误。