基于不同开源DEM的地形校正模型应用对比研究

李 旭,王 志 慧,李 春 意,刘 杨,马 力

(1.河南理工大学测绘与国土信息工程学院,河南 焦作 454003;2.黄河水利委员会黄河水利科学研究院,河南 郑州 450003;3.水利部黄土高原水土保持重点实验室,河南 郑州 450003)

0 引言

地表反射率表征不同物体对太阳辐射的反射和吸收能力[1],受传感器、大气、地形和光照等多种因素影响。山区地形起伏变化改变了太阳、传感器和地表三者之间的位置关系,从而产生“同物异谱,异物同谱”现象[2],因此,地形校正是山区遥感影像预处理中不可缺少的一步。地形校正指通过一系列的变换,将斜坡上的像元辐射值转化为理想表面(通常取水平面)上的像元辐射值,以此减弱或消除因地形效应引起的辐射值变化,从而获取目标像元的真实光谱特征[3]。现有地形校正模型可分为三类[4]:①经验模型,应用简单、适用性强,且不需要额外的物理参数,但模型输出的结果不具备实际的物理意义,如波段比值法[5]和经验统计模型[6]等;②物理模型,通过使用大量物理参数获得精确的结果,但模型复杂且需要较多参数,如Proy校正模型[7]、Sandmeier校正模型[8]和ACTOR校正模型[9];③半经验模型,结合了经验模型和物理模型,利用DEM数据计算坡度、坡向和太阳入射角等地形参数并应用于校正模型中,计算简便且输出结果较准确,如C校正(C-correction)模型[6]、SCS+C校正(Sun-Canopy-Sensor+C-correction)模型[10]和Minnaert(MIN)校正模型[11]。

地形校正结果受DEM的分辨率与精度[12]、校正方法[13]、大气条件[14]等多因素的综合影响,国内外学者对不同校正模型的校正效果展开了对比分析。例如:Hantson等[15]对余弦校正、经验统计校正(Statistical-Empirical correction,SE)[16]、C模型、Minnaert模型和改进Minnaert模型进行对比分析,结果表明SE模型校正效果最好;Richter等[17]对C模型、改进Minnaert模型和Gamma模型进行对比实验,结果表明改进Minnaert模型在多数情况下校正效果优秀;王飒等[18]使用余弦模型、C模型、SCS+C模型和Minnaert模型进行地形校正,指出C模型的应用范围最广泛。随着AW3D30、SRTM、GDEM、MERIT[19]和ALOS等开源DEM的免费开放和广泛应用[20],基于开源DEM的半经验地形校正模型成为地形校正中常用的方法,然而针对不同开源DEM最合适的校正模型以及最佳的组合方式仍缺乏研究[21]。为此,本文以黄土丘陵区为例,基于Landsat8 OLI影像,利用3种开源DEM(AW3D、SRTM、GDEM)和4种地形校正模型(C模型、SCS+C模型、SE模型、MIN模型)开展地形校正,并通过构建地形校正效果评价指标,对由局地入射角引起的光谱异常、坡向引起的光谱异常和同地类光谱变异性的校正效果进行定量评估,以期为黄土丘陵区地形校正模型选择提供参考,对山区地表定量遥感应用具有重要意义。

1 数据与方法

1.1 研究区概况

研究区为陕西省榆林市绥德县和子洲县(图1),平均海拔为1 028 m,坡度在0.02°～50.7°之间,平均坡度为20.02°,草地、农田、裸地面积占比分别为72.3%、22.7%和4.95%。

图1 研究区概况

1.2 遥感数据

从美国地质勘探局(United States Geological Survey,USGS)下载2020年9月16日Landsat8 OLI地表反射率产品(Path=127,Row=34)并进行大气校正(图1a),空间分辨率30 m,包括可见光、近红外和短波红外范围的7个光谱波段,太阳高度角和太阳方位角分别为39.22°和149.13°。

研究选用3种常用的开源DEM数据集:日本宇宙航空研究开发机构(JAXA)发布的AW3D30 DEM数据集(http://www.eorc.jaxa.jp/alos/en/aw3d30/index.htm),其高程基准和参考椭球面分别为EGM96高程和WGS-84椭球;美国航空航天局(NASA)发布的ASTER GDEMV3和SRTM-1 DEM数据集(http://earthexplorer.usgs.gov/),高程基准和参考椭球均为EGM96高程和WGS84椭球,空间分辨率均为30 m。以Landsat8 OLI影像为基准,利用ENVI软件对3种DEM数据分别进行几何校正处理,几何配准误差控制在0.5个像元以内。另外,地表覆盖数据采用2020年全球地表覆盖GlobeLand30[22](http://www.globallandcover.com/)。

1.3 研究方法

1.3.1 技术路线 首先,通过分析遥感影像地表反射率与基于3种开源DEM计算得到的太阳局地入射角余弦值(cosi)之间的相关性,筛选出精度最高的开源DEM作为基准DEM;其次,利用3种开源DEM和4种地形校正模型(C、SCS+C、SE、MIN)构建不同组合对研究区遥感影像进行地形校正;最后,通过构建地形校正评价指标对不同DEM和校正模型组合的校正效果进行定量评价,筛选出适用于黄土丘陵区遥感影像地表反射率地形校正的最佳组合方式(图2)。

图2 技术路线

1.3.2 地形校正模型 地形校正的前提是计算太阳局地入射角i的余弦值(cosi)(式(1))[23],据此运行C模型(式(2))[6]、SCS+C模型(式(3))[10]、SE模型(式(4))[16]和MIN模型(式(5))[11]。

cosi=cosφZcosθS-sinφZsinθScos(φA-θA)

(1)

式中:φZ和φA分别为太阳天顶角和太阳方位角,θS和θA分别为坡度和坡向,由DEM经空间计算得到。

(2)

(3)

(4)

(5)

1.3.3 DEM刻画局地入射角精度评价 本文通过对原始地表反射率和不同开源DEM的cosi值进行皮尔逊相关性分析,相关性最高的DEM即为能最精准刻画局地入射角、最准确反映地表反射率随地形变化规律的DEM数据集,并将该数据集作为评价地形校正效果的基准DEM。

1.3.4 光谱地形校正效果评价

1)目视效果对比法。目视效果对比是评价地形校正效果常用的方法之一,通过观察影像校正前后的整体色彩差异和阴影区域信息恢复判断基于不同开源DEM和不同校正模型的效果[24]。

2)局地入射角引起光谱异常校正效果评价。由于植被对近红外波段的变化更敏感,本文以草地为例,选用近红外波段进行比较。首先,基于地表覆盖分类图,在草地类型区域随机选取10 000个像元,分别计算地形校正前和校正后的地表反射率与基准DEM计算的局地太阳入射角余弦值(cosi)之间的皮尔逊相关性系数R,其绝对值越小,说明局地入射角对地表反射率的影响越小[13]。校正前后R的相对变化率RΔ(式(6))越大,表明对局地入射角所引起的光谱异常的校正效果越好。

RΔ=(|R|-|Rcorr|)/|R|×100%

(6)

式中:R、Rcorr分别为校正前、后地表反射率与cosi的相关系数。

3)同地类光谱变异性校正效果评价。基于地表覆盖分类图,在草地类型区域随机选取10 000个像元,分别统计校正前和校正后地表反射率的四分位数,然后使用上下四分位数之差(IQR)[4]评价草地类型的光谱变异性,其值越小,说明光谱变异越小,该方法可消除由影像光谱异常值所引起的地形校正效果偏差。校正前后IQR的相对变化率QΔ(式(7))越大,表明对同地类光谱稳定性校正效果越好。

QΔ=(Q-Qcorr)/Q×100%

(7)

式中:Q、Qcorr分别为校正前、后的IQR。

4)坡向引起光谱异常校正效果评价。基于地表覆盖分类图和坡向数据,在草地类型的阴坡区域(坡向=太阳方位角+180°±10°)和阳坡区域(坡向=太阳方位角±10°)分别随机选取5 000个像元,然后计算阴坡(阳坡)校正前、后地表反射率之差D,差值越接近0,表明坡向对地表反射率的影响越小。校正前、后D的相对变化率DΔ(式(8))越大,表明对坡向引起光谱异常的校正效果越好。

DΔ=(|D|-|Dcorr|)/|D|×100%

(8)

式中:D、Dcorr分别为校正前、后地表反射率的差值。

5)综合评价方法。为客观评价基于不同开源DEM和不同校正模型的地形校正效果,将上述3种校正效果评价指标进行算术平均,值最高的开源DEM和地形校正模型组合则为最佳地形校正方法。

2 结果与分析

2.1 不同开源DEM刻画局地入射角精度对比

从图1a选取局部区域(红色框线范围),对原始的假彩色合成影像与AW3D、SRTM和GDEM的cosi进行视觉对比(图3a),可以看出 AW3D的纹理信息更精细,与影像地表反射率明暗起伏更一致,SRTM次之,GDEM最模糊且缺失大量细节信息;提取图3a中黄线上近红外反射率和不同DEM的cosi,分析二者随距离的变化规律,可明显看出,AW3D的cosi变化趋势与近红外波段反射率最吻合(图3b);对局部区域可见光—近红外波段的地表反射率和cosi的相关系数进行统计可知,AW3D的相关性在各波段均最高,SRTM次之,GDEM最低(图3c)。综上,AW3D刻画局地入射角的精度最高,表征地表反射率对地形变化响应规律的能力最强,与Liu等[25]的结论一致,因此以AW3D作为基准DEM,对后续不同校正模型的校正效果进行评价。

图3 局部影像与3种DEM的cosi对比

2.2 目视效果对比

选取地形校正前后假彩色合成影像差异明显区域(图4绿色方框区域),可以看出,原始影像在未校正前凹凸不平,具有强烈的地形效应,经过地形校正后影像变得平坦,立体感降低,阴影区域的亮度提升,光谱信息得到恢复。对比4种校正模型可知,C模型、SCS+C模型和SE模型的整体亮度与校正前基本一致,而MIN模型的亮度发生明显改变,表明C、SCS+C和SE模型的校正效果较理想,MIN模型的校正效果较差;对比3种DEM发现,使用AW3D校正后的阴影区域信息恢复效果最好,图像的整体平整度最高,阳坡和阴坡的差异显著降低,而使用GDEM校正后阴影区域信息恢复效果最差,图像的立体感最强,坡度起伏明显,SRTM介于二者之间。

2.3 局地入射角引起光谱异常校正效果对比

影像校正前后近红外波段反射率与cosi的拟合结果如图5所示,可以看出原始反射率随局地入射角变化明显,拟合斜率较大,表明原始影像存在明显的地形效应。经过地形校正后, C模型、SCS+C模型和SE模型的校正结果接近,而MIN模型的校正效果偏差。在3种DEM中,基于AW3D校正后的散点拟合线斜率接近0,SRTM次之,而使用GDEM的校正效果不明显,拟合线倾斜严重。

图5 影像校正前后近红外波段拟合结果对比

由图6可知,基于AW3D校正的RΔ最大,其次为SRTM,GDEM的RΔ最小;4种校正模型中,SE模型和不同DEM组合的RΔ均最大,上述结果表明,AW3D结合SE模型对局地入射角引起的光谱异常的校正效果最好,减幅达98.2%。

图6 研究区校正前后各波段地表反射率与cosi的相关系数相对变化率对比

2.4 同地类光谱变异性校正效果对比

由图7可知,校正后的Q都有一定降低,3种DEM的校正效果排名为AW3D>SRTM>GDEM;4种校正模型中,SE模型的QΔ最大,SCS+C和C模型次之, MIN模型的QΔ最小,表明AW3D结合SE模型对同地类光谱变异性校正效果最好,减幅达97.4%。

图7 研究区校正前后草地地表反射率的IQR相对变化率对比

2.5 坡向引起光谱异常校正效果对比

由图8可以看出,地形校正前近红外波段反射率变化随坡向变化明显,两者呈高度相关,阳坡和阴坡处的反射率受地形效应影响最明显,阳坡反射率变大,阴坡反射率变小,其他坡向受地形影响较小,与Holben等[26]的结论一致。地形校正后,反射率随坡向变化的程度减弱,阳坡反射率降低,阴坡反射率得到补偿,反射率随坡向分布更均匀。3种DEM中,AW3D校正后的反射率分布最均匀,校正效果最好,SRTM次之,GDEM的校正效果最差,阳坡和阴坡均出现不完全校正现象;4种校正模型中,MIN模型出现过校正现象,其他模型的校正效果基本一致。综上,AW3D和SE模型组合对坡向引起的光谱异常校正效果最好,减幅达38.8%。

图8 草地校正前后近红外波段反射率随坡向的变化

由草地校正前后DΔ的对比(图9)可知,校正效果排名为AW3D>SRTM>GDEM,基于AW3D校正的DΔ最大,说明校正后阴坡和阳坡的差异几乎完全消除;4种校正模型均减小了阴坡、阳坡反射率的差异,排名为C>SE>SCS+C>MIN。

图9 研究区校正前后草地地表反射率之差的相对变化率对比

2.6 综合评价

由表1可知,基于AW3D的SE模型对地形引起的光谱异常减幅最高(78.16%),基于GDEM的MIN模型减幅最低(20.25%)。此外,AW3D和不同地形校正模型结合的校正效果均较好,这与Wu等[27]的结论一致,DEM的质量直接决定了地形校正后的效果,精度越高的DEM校正效果越好,AW3D的质量最高,所以校正效果最好。在4种校正模型中,SE模型和不同DEM结合的校正效果均较好,与Sola等[4]的结论一致,SE模型由于引入了经验参数和平均反射率,可以更好地消除地形效应。

表1 不同开源DEM和地形校正方法的评价指标平均值

3 结论

遥感影像的地表反射率受地形起伏的强烈扰动,消除地形效应至关重要。本文以黄土丘陵区为研究区,基于3种DEM(AW3D、SRTM、GDEM)和4种校正模型(C、SCS+C、SE、MIN),开展不同组合试验对研究区进行地形校正,并构建地形校正效果评价指标对不同组合校正效果进行综合评判。结果表明,AW3D刻画局地入射角的精度最高,能真实反映卫星观测地表反射率对地形变化的响应规律;AW3D和不同地形校正模型组合的校正效果均较好,SRTM次之,GDEM最差;与其他校正模型相比,SE模型和不同开源DEM组合的校正效果均较好;AW3D和SE模型组合为适用于黄土丘陵区卫星地表反射率地形校正的最佳组合。

不同地形校正模型消除地形效应的能力一直存在争议,需要多种评价标准来综合评价校正模型的校正效果[13],但不同评价指标的物理意义不同,量纲也不一致,不能直接对所有指标进行综合计算,无法体现校正方法的综合校正效果。针对该问题,本文采用相对变化率的形式对3项指标统一单位(%),通过求算术平均值实现综合评价。然而,本文仍存在如下不足:①黄土丘陵沟壑区具有塬、墚、峁等特殊地貌,本研究仅在该区域进行试验,缺乏不同复杂山区地形综合对比分析,未来将考虑加入对不同山区地貌的研究分析;②遥感影像地表反射率受地表覆盖类型、植被参数、土壤背景、地形特征等多种因素综合影响,本文随机选取样本点代表研究区所有可能的植被参数、土壤背景、地形特征,导致实验过程中产生不确定性,未来需要考虑更多的影响因素,进行更全面的验证。