面向对象的多时相遥感影像标准化方法研究

胥海威，吉柏锋，张晓磊

（1.河南工程学院 土木工程学院，河南 郑州 451191；2.武汉理工大学 土木工程与建筑学院，湖北 武汉 430070；3.河南省科学院地理研究所，河南 郑州 450001）

通过遥感影像可了解地球表面的分布与变化[1]。不同地物的变化周期不同，搭载可见光传感器的对地观测卫星由于时间分辨率适中，常被用以观察与监测地球环境的变化[2]，如基于多时相卫星影像的地表覆被变化监测，可为研究环境变化[3]、城市变迁[4]、智能农业[5]和灾害评估[6]提供帮助。地表覆被变化监测的前提是多时相遥感影像间具有辐射亮度稳定性和光谱连续性[7]，但由于成像条件的差异性，目前多时相遥感影像普遍存在辐射畸变[8-9]，为解决该问题，遥感应用前的辐射校正是必不可少的预处理步骤。国内外学者提出许多辐射校正理论与方法[9-10]，可归结为绝对辐射校正[11]与相对辐射校正[12]两类，绝对辐射校正能提供高精度地表反射率，但需要很多野外实地获取的大气参数，且获取参数时间需与影像成像时间一致，实际操作难度颇高[13]；相对辐射校正是在同一区域多时相影像同波段灰度值存在线性关系的假设上提出的，主要分为基于像元对和基于分布两类[14]，前者通过在多时相遥感影像中选择共同点的方式建立辐射校正方程[15]，计算校正系数，存在一定的主观性，对几何配准精度要求较高[15]，后者则通过分析多时相遥感影像分布特征，建立影像间的辐射校正方程，消除一定的主观性，性能不依赖几何配准精度[16]。辐射校正对保留地物真实变化信息，有效去除其他因素造成的差异具有积极作用[17]，但目前的辐射校正方法不仅无法避免主观性对辐射校正结果的影响，而且逐波段纠正遥感影像，不可避免地破坏了地物的光谱曲线。

在前人研究的基础上，本文提出了面向对象的影像标准化方法。该方法的处理对象并非像元，而是提取的地物对象。参照eCognition 软件的思想，通过预分类基准影像数据提取地物图斑，并以地物类型为单位，根据地物的光谱特征，通过设置聚类指数来选择参与计算转换函数的像元数量，计算不同地物类型的转换函数；再对其余时相遥感影像进行处理，从而实现对多时相遥感影像辐射差异的标准化。该方法基本避免了训练样本选择的主观性，且待处理多时相影像的地物类型图斑亮度被逐波段转换为基准影像的图斑亮度，其地物光谱特征能得到较好保护。本文利用4期Landsat8 OLI遥感影像进行验证，结果表明该方法能在保持地物光谱特征的基础上，有效消除时相不同对遥感影像造成的辐射畸变，也可用于地物变化监测拓宽该领域的数据选择范围。

1 研究方法

1.1 基本假设

多时相遥感影像的辐射差异是由影像成像环境不同造成的，其中大气环境是主要因素。本文提出的方法基于以下假定：①多时相遥感影像中同一地物的光谱差异是由于成像时间不同，太阳辐射能量通过大气时衰减程度差异造成的；②遥感影像内同类地物反射率相近并服从正态分布。

1.2 影像预处理

遥感影像预处理主要包括两部分：辐射定标，将图像的亮度灰度值转换为实际物理意义的大气顶层辐射亮度或反射率；大气校正，消除大气和光照等因素对地物反射的影响，再将定标值还原为地表真实信息，并能高保真地恢复地物波谱信息。

影像由于成像时间不同，行列数难免参差不齐，通过对多期影像进行几何校正来消除位置原因造成的影像差异。由于遥感影像用途众多，用户关心的遥感信息各异，因此在影像标准化处理前，对遥感影像进行分类并根据用户关注程度设置不同的权重是必要的，如国土部门的应用可按照土地利用分类标准，林业部门的应用可按照林地的分类体系等，而与主题无关的影像内容，权重可设置为零。分类精度由遥感影像的分辨率和分类方法决定，本文不做讨论。本文采用支持向量机对基准遥感影像进行水体、植被、人工建筑和裸地等感兴趣信息的提取。

1.3 标准化方法

根据大气传输方程[1]，第k个波段进入方位为(θ,φ)传感器的辐射能量为：

式中，Ik为卫星接收到的总辐射；τk为大气层透过率；Ratk↑为大气上行辐射；Rsk↑为太阳能量被大气多次后向散射到达传感器的部分。

由于研究区影像为可见光近红外遥感，地表发射的热辐射与大气热辐射位于热红外波段，并不能被传感器接收，且与太阳辐射相比，地表与大气辐射位于可见光的电磁辐射可忽略，因此式（1）可略写为：

Landsat系列卫星属于太阳同步轨道卫星，不同时相影像的θ、φ、θs、Ek均不变，相同地物的εk相同，由于成像时间不同、大气状况略有差异、太阳辐射在大气中的折射路径发生改变导致Rsk↓、τk()θ,φ、τ(θs)、ρbk、Rsk↑(θ,φ)变化，因此标准化方法在上述因子中选择一组数据的辐射传输路径s作为基准，其余时相的传输路径通过相对辐射校正方法标准化为s，从而达到影像辐射标准化的目的。

影像预处理和分类后，设基准影像为f(i,j)，其余时相影像为ft(i,j)，i=1,2,…,M，j=1,2,…,N，所有影像均有K个波段，影像内有P类典型地物，则标准化方法的实质是在每个波段内找出变换函数gtkp，使得多时相影像ftkp，变换为：

式中，htkp为标准化后的其余时相影像，用于与基准影像进行比较，针对不同时相、不同波段、不同地物类型，函数形式相同，仅系数存在差异；t为时相，t=1,2,…,T；k为波段序号，k=1,2,…,K。

由于用户对影像数据不同地物的关注度差异，本文赋予P类典型地物的不同权重Wp，并根据最小二乘原理，给出单时相单波段的辐射畸变校正指数，即

式中，p为地物类型，p=1,2,…,P；i、j分别为图像的行数和列数。

计算过程中t、k均确定，根据曲线拟合情况采用二次至多次多项式回归模型计算相应的gtkp表达式。针对多时相遥感影像，需逐波段计算使Ltk最小的gtkp，再结合其余影像的统计特征，寻找最恰当的变换函数。

标准化方法的具体流程为：①从多时相遥感影像中选择一幅基准影像f，其余时相影像为ft；②对影像进行预处理，以基准影像为标准进行关心地物类型的影像信息提取，得到结果分类影像fc；③将fc与其他时相影像叠加，得到不同时相、不同波段、不同地物类型的图斑对象集合ftkp，并与基准影像叠加得到图斑对象集合fkp；④对ftkp进行聚类分析，根据聚类阈值a确定参与标准化计算的像元数目x；⑤根据对不同地物类型的关注程度设置权重Wp，逐时相、逐波段采用多项式拟合方法计算使Ltk最小的gtkp；⑥利用gtkp对影像进行变换。

1.4 噪声与地物变化点处理

由于自然界相邻地物的强相关性，不同地物类型的变化较缓慢，自然地物很少出现类似人工地物界线分明的现象，在遥感影像中则表现为影像灰度值的连续性和可微性。通过预分类过程，标准影像的分类结果与其余时相影像叠置后形成大小不一数目繁多的Q个地类图斑，在第q（q=1,2,…,Q）个图斑内地物反射率服从正态分布，包含m个像元，考虑到噪声，未选择所有像元参与标准化处理，而是对图斑内m个像元灰度值进行聚类分析得到n个类别，再按从多到少的顺序对类别进行排序与选择，计算选择像元的累积数目x。x/m定义为进行标准化处理的像元比例，若超过预设的聚类阈值α，则剩余像元不再选择，即不再参与影像的标准化处理，这部分像元归为噪声点或地物变化点。由于未参与标准化计算，因此像元灰度值与其余像元区别明显，简单处理后即可获取其范围。

2 实验与结果分析

本文采用4期Landsat8 OLI遥感影像进行实验，影像条带号为123，行号为38，日期分别为2017-10-30、2017-12-17、2018-03-23、2018-04-08，经纬度范围为112.9°～115.4°E、30.6°～32.8°N，位于河南省信阳市与湖北省交界处附近，见图1。本文选择2017-10-30的遥感影像作为基准影像，采用支持向量机进行分类，经过精度评价后得到分类精度为92.6%，获取实验区范围内分类图。由于本次实验对4类典型地物的关注度相同，因此4类地物的权重均设为0.25。将分类结果分别与其余3期影像进行叠置，分别按照聚类阈值a为85%、90%、95%逐类别进行标准化处理，得到后面3期遥感影像的标准化图像（图2）。

图1 研究区遥感影像图

图2 遥感影像标准化后图像

通过对比标准化前后的3期影像各波段的统计量，本文计算得到3个时相遥感影像标准化前后的辐射畸变指数Ltk（表1）；再绘制标准化后3期遥感影像典型地物的光谱曲线，得到4类地物在不同时相的光谱曲线，并与基准遥感影像的光谱曲线进行对比（图3）。

表1 标准化前后地物辐射差异

图3 标准化前后光谱曲线对比图

结合理论论述和实验成果图表发现：①Ltk作为最重要指标，可反映标准化处理前后，其余时相影像与基准影像的差异；②标准化处理后的多时相遥感影像内都存在细碎的杂色图斑，但随着聚类指数从85%增至95%，标准化后的影像零碎图斑明显减少（图2 a、2d、2g）；③横向对比图2 可知，无论设定何种聚类指数，在地物类型边缘地带都存在部分杂色图斑，对比原始影像发现这些杂色图斑是由于地物本身发生变化造成的，地物类型发生变化，光谱曲线随之发生改变，标准化方法灵敏捕捉到变化，在标准化后的影像中标识出来，因此标准化方法对地物的变化比较敏感，多期遥感影像的标准化可进行对地观测任务；④影像预处理的Flaash大气校正方法会使部分影像的部分波段地物反射率变为负值，在影像上虽无法体现，但在定量计算时会对模型造成影响，其原因在于Flaash 大气校正方法采用MODTRAN 辐射传输模型模拟成像中的大气过程，属于全球化模型，无法与实际情况完全吻合，且许多大气参数通过影像估算，当影像中有强吸收或高反射的地物时，就会出现反射率为负值或超过最大值的情况，标准化方法将其最小值设置为0（表1）；⑤标准化方法虽将异质像元标识出来，但并未对其进行删减处理，因此未出现像元标准差增加的情况，对比标准化方法处理前后的遥感影像统计值可知，标准化后3 期影像的反射率最值范围相同，均值小幅降低，Band 1～4标准差减少，Band 5～7标准差不变或增加（表1），表示影像灰度数据的集中度并没有明显变化；⑥实验采用3 期Landsat 影像共21 个波段，其中13 个波段的Ltk下降，总幅度为60.23，8 个波段的Ltk上升，总幅度为32.7，Ltk整体呈下降趋势，说明标准化方法能降低辐射畸变，但仍有继续提升的空间，Ltk上升的最大值为2017-12-17 的Band 5，达到了11.67，反推其原因，应是参与计算的gtkp的样本数目过多，由于循环次数有限，数据收敛速度太慢，未计算得到最恰当的gtkp；⑦2018-04-08 遥感影像有云存在，且云覆盖的地物类型较多，但标准化处理后各类地物都将云覆盖范围准确分离出来，可见较厚的云层破坏了地物的光谱连续性，标准化方法对此较敏感，因此标准化方法对去除云层影响同样有帮助；⑧聚类指数a的设置与光谱曲线特征无关，由于标准化方法将地物反射率按地物类型逐波段换算为基准影像的地物反射率，因此可以保证影像地物光谱的反射特征，且在实验中得到了验证，由图3 可知，标准化前后的地物光谱曲线虽有些轻微差异，但光谱曲线特征未发生明显变化，且4类地物的光谱曲线特征也与光谱库中典型地物的光谱曲线特征一致，可见标准化能较好地保持地物本身的光谱特性。

3 结语

针对多时相遥感影像由时相引起的辐射传输路径差异最终导致的遥感影像辐射畸变，本文提出了面向图斑对象的影像标准化方法，采用同一地区相同传感器4 期遥感影像进行实验，并对比分析了实验结果。结果表明：①本文提出的标准化方法能有效降低多时相影像辐射畸变，通过计算变换函数，将其余时相的遥感影像以基准影像为标准进行辐射畸变校正；②该方法可将多时相遥感影像的地物反射率以基准影像为标准进行统一，使波段间地物对象的光谱特征得到较好保持；③该方法将聚类后的数据全部参与运算，从而避免了样本选择的主观性；④辐射校正方法以像元为单位进行校正，易出现椒盐效应，标准化方法以图斑为单元，可避免这种现象，且更容易与地理实体对照；⑤该方法的缺点在于标准化结果较依赖基准影像的质量，必须选择成像环境较好的遥感影像作为基准影像；⑥该方法可根据不同应用需要，设置不同地物类型的聚类阈值、影响各地物类别的异质像元点的数目，与地物变化的检测精度有效对应。

后续研究计划包括两个方面：①gtkp的计算仍有不足，后续着重研究如何使用有限的像元数据快速计算得到有效的gtkp，尝试采用人工神经网络算法进行改进；②多源遥感辐射畸变的成因比多时相遥感影像更加复杂，但理论上同样可以采用标准化方法进行辐射校正，拟对多源遥感影像的标准化方法进行研究。