加入高程因子的航空高光谱影像大气辐射校正

图1 实验区数据源

Fig.1 Data of the study area

2 实验方法

鉴于本实验重点在于分析高程因子变化对航空高光谱成像的影响,因此研究思路为: 假定每一个高程点为水平均一朗伯体,然后选择相应的高程范围,分别解算出不同高程对应的大气辐射传输参量并对比分析。

2.1 基于MODTRAN辐射传输分量计算

在假定地表为水平均一朗伯体的情况下,传感器接收的辐射亮度可表示为^[14]

(1)

L = L_{path} + \frac{ρ_{t}}{1 - ρ_{t} S} F_{d} \cdot T

式中: L_path为大气上行辐射; S为大气半球反照率; F_d为下行总辐射; T为地面像元与传感器之间的大气透过率; ρ_t为地表反射率。

采用MODTRAN软件,选择中纬度夏季(middle latitude summer)大气模式、乡村(rural,visible=23 km)气溶胶模式; 设置传感器垂直观测。

选择3 500 m,4 000 m,4 500 m,5 000 m,5 500 m共计5个地面高程值。以高程3 500 m为例,在MODTRAN4.0软件中输入设置的反射率和传感器高度参数(如表2所示),分别模拟计算传感器在航飞高度(6 900 m)、地面高度(3 501 m)以及不同地表反射率数值(0,0.1,0.2)情况下接收的辐射亮度。然后将5组模拟计算的结果代入式(1),解算出每个波段对应的L_path,S,F_d和T。

表2 模拟输入参数

Tab.2 Input parameters

地面高程/m	地表反射率	传感器高度/m
3 500	0	6 900
	0.1	6 900
	0.2	6 900
	0.1	3 501
	0.2	3 501

新窗口打开| 下载CSV

参照以上参数设置,分别解算出其他高程值对应的辐射传输分量。进一步对5个高程值对应的大气辐射传输分量对比分析,以确定地形高程变化对传感器接收信号的影响。

2.2 加入高程因子的大气辐射校正

依据式(1),ρ_t可表示为^[15]

(2)

$ρ_{t} = \frac{L - L_{path}}{F_{d} \cdot T + S (L - L_{path})}$ 。

若经2.1节结果分析,高程变化对大气辐射传输分量有较大影响,则依据式(2)开展航空高光谱影像大气辐射校正须考虑到各项参数受高程变化造成的影响。

实验思路为依据一定数量高程值对应的大气辐射传输分量,通过函数拟合与插值计算,并结合高光谱像元对应的高程值,计算出任意高程值对应的各项大气辐射传输参数,然后完成航空高光谱影像的反射率反演计算。

基于IDL语言,在ENVI+IDL编程模型下实现以上过程。具体实现步骤如下:

1)首先,将计算的各项大气辐射传输分量依据标准高程值和波段存储构建为查找表,格式如图2所示。将各项参数按照对应的波段分开,逐波段将大气参数生成矩阵,每一行对应一个高程值,每一列对应一项大气参数。具体到本实验中,包括3 500 m,4 000 m,4 500 m,5 000 m,5 500 m共5行,以及L_path,T,S和F_d共4列。

图2

图2 查找表结构示意图

Fig.2 Structure of look-up table

2)读取DEM,与高光谱影像进行坐标范围和空间分辨率匹配,并获取每个高光谱像元对应的高程值,生成高程数组。具体IDL代码为

dem=envi_get_data(fid=fid_dem,dims=dims1,pos=0)

envi_enter_data,dem,r_fid=r_fid

resize_doit,fid=r_fid,dims=dims1,pos=0,/in_memory,interp=3,r_fid=r_fid1,rfact=rfact

elevation=envi_get_data(fid=r_fid1,dims=dims1,pos=0)

3)逐波段对各项大气参数进行回归计算,解算出回归系数,进一步与高程数组运算,解算出每个像元对应的大气辐射校正参数。具体IDL代码为

lut=envi_get_data(fid=fid_lut,dims=dims_lut,pos=pos_lut[k])

HeightArr=[3500,4000,4500,5000,5500]

Lup=reform(lut[0,*,*],5,1)

Trans=reform(lut[1,*,*],5,1)

Sr=reform(lut[2,*,*],5,1)

Fd=reform(lut[3,*,*],5,1)

L0=regress(heightArr,Lup,const=L1)

T0=regress(heightArr,Trans,const=t1)

S0=regress(heightArr,Sr,const=s1)

F0=regress(heightArr,Fd,const=f1)

L[mask]=L1+elevation[mask]*L0

T[mask]=T1+elevation[mask]*T0

S[mask]=S1+elevation[mask]*S0

F[mask]=F1+elevation[mask]*F0

4)最后,逐波段、逐像元读取高光谱影像辐射亮度数据,结合大气参数,并基于式(2)进行反射率计算。具体IDL代码为

fork=0,num_bands-1do begin

band=envi_get_data(fid=fid_image,dims=dims_image,pos=pos_image[k])

image[mask,k]=(image[mask,k]-L[mask])/(F[mask]*T[mask]+S[mask]*(image[mask,k]-L[mask]))

endfor

image=reform(image,num_cols,num_rows,num_bands)

3 实验结果与分析

3.1 大气辐射传输分量分析

将3 500 m,4 000 m,4 500 m,5 000 m,5 500 m共5个高程值分别对应的L_path,T,S和F_d进行对比分析。图3为传感器接收的5个不同高程的大气上行辐射L_path计算结果对比。

图3

图3 不同高程对应的大气上行辐射

Fig.3 Path radiation corresponding to different elevation

图3中可以看出,高程越低,与传感器之间的距离越大,则L_path值越大。尤其在可见光波段范围内,5个不同高程(对应不同传输距离值)的L_path值差异较大,这主要与该波段范围内大气散射强烈有关。经统计,可见光波段范围内,3 500 m与5 500 m对应的L_path值差异为0.231.23 μWcm^-2·sr^-1·nm^-1,由于CASI影像的像元DN值单位为1 000 DN=1 μWcm^-2·sr^-1·nm^-1,则L_path差异造成的信号DN值差异可达2301 230之间,而且CASI高光谱影像大部分像元辐射亮度值(DN值)在2 0008 000之间,因此这种差异影响较大。在波长大于800 nm的范围内,L_path的差异逐渐变小,数值在0.050.14 μWcm^-2·sr^-1·nm^-1之间,影响相对较小。由此可见,航空传感器在获取数据过程中,每个地面像元对应的高程不同,其L_path值均存在一定差异,且这种差异不可忽略。

图4为传感器至5个不同高程间的大气透过率T计算结果对比。

图4

图4 不同高程至传感器的大气透过率

Fig.4 Atmospheric transmittance of different elevation to sensor

图4中可见,不同高程对应的T均不相同,高程越高,对传感器距离值越小,则T越高,但除第32波段(938.4 nm)和第33波段(957.5 nm)外,其他波段的T差别在10^-2数量级。而CASI数据的第32波段(938.4 nm)和第33波段(957.5 nm)处的T较低,主要是由于水汽在940 nm波长具有吸收作用。

图5为不同高程对应的大气半球反照率S的结果对比。

图5

图5 不同高程对应的大气半球反照率

Fig.5 Atmospheric hemisphere albedo corresponding to different elevation

图5中可见,首先,S的数值在短波波段范围偏高,但数值也低于0.25,在近红外波段范围内数值低于0.05。另外,高程越高,对应的S数值越小。

图6为不同高程接收的下行总辐射F_d的结果对比。

图6

图6 不同高程接收的下行总辐射

Fig.6 Downward radiance received at different elevation

图6中可见,高程越低,接收的F_d值越大。经统计,3 500 m与5 500 m对应的大部分波段下行辐射数值差异为12 μWcm^-2·sr^-1·nm^-1,只有第32波段(938.4 nm)和第33波段(957.5 nm),由于水汽影响大气透过率差异较大,从而造成不同高程对应的F_d差异更大,分别为3.1 μWcm^-2·sr^-1·nm^-1和3.8 μWcm^-2·sr^-1·nm^-1。因此,大部分波段下行辐射的DN值差异可达1 0002 000,假定地物反射率为0.3,大气透过率为0.95,则此部分对传感器信号造成的DN值差异约在285570范围,其影响亦不可忽略。

综合以上分析,不同高程对应的L_path,T,S和F_d均存在差异,且不可忽略。因此要实现精确的山地航空高光谱影像大气辐射校正,必须考虑加入高程因子的影响。

3.2 反射率计算结果对比

本文所选实验区为高海拔无人区,因此缺乏地面实测光谱数据为参考验证。为分析加入高程因子校正后的效果,将解算的反射率结果与FLAASH大气校正处理后的结果进行对比分析。FLAASH是遥感图像反射率反演的常用方法,基于MODTRAN辐射传输模型开发^[16,17]。在FLAASH大气校正过程中,同样选择中纬度夏季(middle latitude summer)大气模式、乡村(rural,visible=23 km)气溶胶模式,高程参数选择测区平均高程为4 500 m。

选择实验区中变质岩地表、苔藓植被、冰雪以及岩浆岩地表共4类地物(如图7所示),对反演后的反射率进行对比。

图7

图7 典型特征地物

Fig.7 Distribution of selected ground objects

对比图7中所选4类地物的FLAASH校正后与加入高程因子校正后的反射率,结果如图8所示。

图8

DOI:10.1175/1520-0450(2003)042<0113:SVOSRF>2.0.CO;2 URL [本文引用: 1]

图8 各类地物反射率光谱对比

Fig.8 Reflectance comparison of different methods

从图8可以看出,FLAASH大气校正结果与加入高程因子校正后的结果相比,其反射率谱形基本一致,主要差别在于短波波段,FLAASH校正计算的结果甚至出现负值,无疑是错误的。本文认为如第3.1节所述,不同高程对应的L_path,T,S和F_d之间均存在差异,且波长越短差异越大,而FLAASH校正在假定地表为均一朗伯体的前提下,参考单一高程值进行校正,从而导致反射率计算结果存在误差,且在短波波段尤为明显。因此,对于山区航空高光谱影像而言必须加入高程因子的影响,才能实现更为精确的大气辐射校正。

4 结论

针对地形高程因子,在假定每一个高程点为水平朗伯体的前提下,计算分析了地形起伏导致的辐射传输距离变化对航空高光谱遥感成像的影响,并实现了加入高程因子的大气辐射校正。主要结论如下:

1)基于MODTRAN的计算结果表明,地形高程的变化对航空高光谱影像的成像过程存在重要影响。

2)加入高程因子的大气辐射校正结果与FLAASH大气校正结果相比,两者在反射率曲线的谱形方面总体一致,但反射率数值存在差异,尤其是在短波波段FLAASH大气校正结果中甚至出现负值,无疑是错误的。

3)无论是基于MODTRAN计算的大气辐射校正参数(大气上行辐射、大气透过率),还是FLAASH大气校正结果,都可以发现高程因子变化的影响在短波波段更为显著,这主要与短波波段的大气散射辐射强烈有关。

高程因子的变化对山地航空高光谱影像成像过程的影响不可忽略,要实现精确的大气辐射校正就必须消除其影响。未来还需要从地形高程和地形角度2个方面共同开展研究,才能进一步实现高精度的航空影像地形辐射校正。

参考文献

原文顺序

文献年度倒序

文中引用次数倒序

被引期刊影响因子

[1]

Schaaf

, Li

, Strahler

Topographic effects on bidirectional and hemispherical reflectances calculated with a geometric-optical canopy model

[J]. IEEE Transactions on Geoscience and Remote Sensing, 1994,32(6):1186-1193.

DOI:10.1109/36.338367 URL [本文引用: 1]

[2]

Oliphant A

, Spronken-Smith

R A

,Sturman

A P

,et al.

Spatial variability of surface radiation fluxes in mountainous terrain

[J]. Journal of Applied Meteorology, 2003,42(1):113-128.

[3]

Wen J

, Zhao X

, Liu

, et al.

An improved land-surface albedo algorithm with DEM in rugged terrain

[J]. IEEE Geoscience and Remote Sensing Letters, 2014,11(4):883-887.

DOI:10.1109/LGRS.2013.2280696 URL [本文引用: 1]

[4]

李爱农, 边金虎, 张正健 , 等.

山地遥感主要研究进展、发展机遇与挑战

[J]. 遥感学报, 2016,20(5):1199-1215.

Li A

, Bian J

, Zhang Z

, et al.

Progresses,opportunities,and challenges of mountain remote sensing research

[J]. Journal of Remote Sensing, 2016,20(5):1199-1215.

[5]

宫鹏

遥感科学与技术中的一些前沿问题

[J]. 遥感学报, 2009,13(1):13-23.

Gong

Some essential questions in remote sensing science and technology

[J]. Journal of Remote Sensing, 2009,13(1):13-23.

[6]

王少楠, 李爱农 .

地形辐射校正模型研究进展

[J]. 国土资源遥感, 2012,24(2):1-6.doi: 10.6046/gtzyyg.2012.02.01.

从经验模型、物理模型及半经验模型3个方面综述了地形辐射校正模型研究的主要进展,并从模型的输入参数、前提假设条件和评价方法3个方面讨论了现有模型存在的一些问题。最后,对地形辐射校正模型的研究前景进行了展望,认为可以从多源或多时相数据、图像增强或信息填充等方面恢复阴影区域信息,提高模型的校正精度,考虑引入新的数理统计方法进行精确的定量评价,并建议地形辐射校正模型研究应该更加重视阴影区域光谱信息的恢复,尤其是在地形复杂的山区。

Wang S

, Li A

The progress in the study of topographic radiometric correction models

[J]. Remote Sensing for Land and Resources, 2012,24(2):1-6.doi: 10.6046/gtzyyg.2012.02.01.

[7]

黄博, 徐丽华 .

基于改进型Minnaert地形校正模型的应用研究

[J]. 遥感技术与应用, 2012,27(2):183-189.

<p>在复杂的山区,地表接受的太阳辐射受到地形起伏的影响严重,干扰提取遥感专题信息精度,地形校正成为提高定量化遥感数据处理精度的关键技术之一。在目前已有的多种地形校正模型中,陆灯盛改进型Minnaert模型不仅所需基础数据获取容易,并且解决了仅依靠单一的Minnaert系数k造成的校正精度不高的难题,是山区遥感信息处理过程中地形校正的有效方法之一。选取浙江省青田县的ETM影像和对应的DEM数据,分别利用余弦模型、原始的Minnaert模型和改进型Minnaert模型对遥感影像进行地形校正。研究结果表明改进型Minnaert模型能有效地去除地形的影响,同时又能保持地物原有的纹理特征和光谱特性,具有较强的实用性。基于图像统计分析,与余弦模型及原始的Minnaert模型校正结果相比,改进型Minnaert模型地形校正效果明显增强。</p>

Huang

, Xu L

Applied research of topographic correction based on the improved Minnaert model

[J]. Remote Sensing Technology and Application, 2012,27(2):183-189.

[8]

Lu D

, Ge H

, He S

, et al.

Pixel-based Minnaert correction method for reducing topographic effects on a Landsat 7 ETM+ image

[J]. Photogrammetric Engineering and Remote Sensing, 2008,74(11):1343-1350.

DOI:10.14358/PERS.74.11.1343 URL [本文引用: 1]

[9]

Wen J

, Liu Q

, Liu

, et al.

Parametrized BRDF for atmospheric and topographic correction and albedo estimation in Jiangxi rugged terrain,China

[J]. International Journal of Remote Sensing, 2009,30(11):2875-2896.

DOI:10.1080/01431160802558618 URL [本文引用: 1]

[10]

Li A

, Wang Q

, Bian J

, et al.

An improved physics-based model for topographic correction of Landsat TM images

[J]. Remote Sensing, 2015,7(5):6296-6319.

DOI:10.3390/rs70506296 URL [本文引用: 1]

[11]

Richter

, Kellenberger

, Kaufmann

Comparison of topographic correction methods

[J]. Remote Sensing, 2009,1(3):184-196.

DOI:10.3390/rs1030184 URL [本文引用: 1]

[12]

Hantson

, Chuvieco

Evaluation of different topographic correction methods for Landsat imagery

[J]. International Journal of Applied Earth Observation and Geoinformation, 2011,13(5):691-700.

DOI:10.1016/j.jag.2011.05.001 URL [本文引用: 1]

[13]

Balthazar

, Vanacker

, Lambin E

Evaluation and parameterization of ATCOR3 topographic correction method for forest cover mapping in mountain areas

[J]. International Journal of Applied Earth Observation and Geoinformation, 2012,18(1):436-450.

DOI:10.1016/j.jag.2012.03.010 URL [本文引用: 1]

[14]

童庆禧, 张兵, 郑兰芬 . 高光谱遥感——原理、技术与应用[M]. 北京: 高等教育出版社, 2006: 86-105.

Tong Q

, Zhang

, Zheng L

Hyperspectral Remote Sensing:Principle,Technology and Applications[M]. Beijing: Higher Education Press, 2006: 86-105.

[15]

Adler-Golden S

, Matthew M

, Bernstein L

, et al.

Atmospheric correction for short-wave spectral imagery based on MODTRAN4

[J]. Proceedings of SPIE,the International Society for Optical Engineering, 1999,3753:61-69.

[16]

Kaufman Y

Atmospheric effect on spectral signature-measurements and corrections

[J]. IEEE Transactions on Geoscience and Remote Sensing, 1988,26(4):441-450.

DOI:10.1109/36.3048 URL [本文引用: 1]

[17]

Anderson G

, Felde G

, Hoke M

, et al.

MODTRAN4-based atmospheric correction algorithm:FLAASH (fast line-of-sight atmospheric analysis of spectral hypercubes)

[J]. Proceedings of SPIE,the International Society for Optical Engineering, 2002,4725:65-71.