第一作者简介: 徐大琦(1979-),博士,高级工程师。主要从事光学卫星图像模拟、数据处理和遥感应用等研究,担任我国HJ-1A/1B/1C,CBERS-03/04等陆地观测卫星地面系统数据模拟分系统的主任设计师。Email: xudaqi@spacechina.com。
卫星图像模拟是在卫星发射之前,采用计算机模拟方法模拟图像的波段特征、空间几何特征、辐射特征、星历数据和格式编排的一项技术。为研究陆地卫星图像模拟技术,回顾了我国过去6 a自主开发卫星遥感图像模拟系统的发展过程,介绍了图像模拟系统设计及关键技术的实现情况。目前该系统具备的模拟波段包括可见光、近红外到热红外波段,模拟的空间分辨率在300~3 m。在模拟过程中,采用遥感辐射传输模型实现光谱特征模拟; 采用PROSPECT+SAIL模型模拟植被覆盖区的光谱,采用波谱库数据模拟非植被区的光谱; 基于对大气辐射传输过程的线性分解,建立了大气辐射传输过程查找表(LUT),在保证一定模拟精度的前提下显著提高了模拟计算的速度; 采用考虑地形起伏的高精度几何定位模型,逐像元计算出卫星观测视线与地球表层的交点,实现了几何信息的精确模拟; 最后,采用卫星发射后的观测数据和实验场测量数据验证了该技术的模拟精度。
Satellite image simulation technology aims at simulating the band features, space geometric features, radiation characteristics, ephemeris data and format of the satellite image by the method of computer simulation before the launching of satellites. To study the Landsat image simulation technology,this paper makes a review on the development of China’s own satellite remote sensing image simulation system over the past six years,and describes the design of China’s image simulation system as well as the key technology. At present, the simulated bands of the system include the bands from visible light, near infrared to thermal infrared band, with the simulated spatial resolution between 3 m to 300 m. In the process of simulation, the authors used the remote sensing radiative transfer model to simulate the spectrum characteristics, employed the PROSPECT+SAIL models to simulate the spectrum of the areas covered by vegetation, and adopted spectral library to simulate the spectrum of the non-vegetation area. Based on linear decomposition of the atmospheric radiative transfer process, the authors set up a look-up table (LUT)of the atmospheric radiative transfer process so as to improve the speed of simulation calculation on the premise of guaranteeing simulation precision significantly. In order to simulate the precise geometry information, the authors used high precision geometric positioning model on the basis of considering the topographic relief, calculated the intersection between the line of sight for satellite observation and the Earth’s surface pixel by pixel. Finally, The authors used the observation data after the launching of the satellite and the field measured data in the experimental field to verify the simulation precision of the image simulation technology described in this paper.
20世纪90年代后期, 随着遥感技术的发展, 逐渐出现了模拟遥感系统的方法和软件。1998年, Banon等[1]以高分辨率SPOT图像为基础模拟了CBERS波段4图像; 2001年, Borner等[2]发表了关于超光谱遥感成像系统模拟方面的研究论文, 介绍了模拟软件SENSOR及其被欧空局ESA项目APEX采用的情况。2002年, 美国国家宇航局(NASA)的Stennis 空间中心(SSC)开发了模拟传感器的算法和软件产品— — ART, 该软件是一个集成计算机算法和数学模型的遥感系统模拟工具箱, 通过对现有高光谱和高空间分辨率数据的特性进行提前定义, 生成模拟的多光谱图像; 另外, NASA还提出了基于Monte Carlo的三维辐射传输方程的高光谱遥感模拟系统[3]。2002 年, 美国纽约数字成像和遥感实验室的Brown等 [4, 5]开发了DIRSIG(the digital imaging and remote sensing image generation model)遥感成像模拟系统模型, 采用图像合成技术模拟了可见光— 热红外波段的多光谱或高光谱图像, 比较完整地考虑了地表、大气(MODTRAN和 FASCODE)与传感器3个方面的影响, 但地表场景模拟方法的物理意义还不甚明确; 2003年, Boggione 等组合了TM的多波段模拟全色波段[6]; 2009年和2013年, 法国图卢兹一个试验室介绍了一个新的基于MODTRAN的辐射传输模型DART, 利用高几何分辨率的机载数据, 对地表三维结构与大气影响进行了综合考虑和处理[7, 8]。此外, 美国陆军夜视电子传感器管理局研制了用于光电系统性能分析的FLIR92软件[9]; 以色列的光电工业有限公司遥感系统开发部开发了用于可见光光电侦查成像系统性能预估与评价的EOAPAS软件[10]; Multigen-Paradigm 公司推出的商业产品Vega系列模拟模块, 使用少量参数的简化传感器模型满足了实时动态显示的目的[10]。
上述成像模拟方法可大体归纳为“ 硬” 和“ 软” 2种, 两者的区别主要体现在对地表场景辐射的模拟方面。“ 硬” 方法是通过物理硬件实现成像模拟; 而“ 软” 方法则利用计算机及相关模型和数据实现成像模拟, 大致可分为如下几种类型:
1)图像到图像的成像模拟。即图像合成法, 利用已有的图像数据经变换后得到模拟图像, 例如用SPOT数据模拟CEBRS图像、用TM多波段数据合成全色波段图像等。
2)基于虚拟现实技术、光线追踪技术和纹理映射技术的成像模拟。通过建立三维实体的几何模型, 给定实体的纹理和材质, 映射出实体的组分温度和发射率; 根据给定的观测角度进行光线追踪, 并利用热传导等物理方程求出实体的热辐射, 例如Vega和DIRSIG地表场景热辐射模拟等。
3)利用相关辅助数据和遥感物理模型的正向成像模拟。利用已有数据提供的先验知识, 结合相关遥感物理模型模拟出遥感图像, 例如 DART, SENSOR和ART等。
4)基于数学方法的成像模拟。利用Monte Carlo等数学方法直接模拟整个成像过程, 不考虑明确的物理意义。对大气场景和传感器成像系统的模拟都是采用成熟的大气辐射传输模型和传感器成像几何物理模型, 没有太大差别。
为了快速发展我国的航空航天遥感事业, 研发具有我国自主知识产权的遥感成像模拟系统具有重要的应用价值和现实意义。本文从辐射特性和几何特性2个方面介绍了我国陆地观测卫星图像模拟技术的设计与关键技术的实现情况; 并通过一组卫星与地面实测数据对辐射模拟精度进行了验证; 最后对模拟技术未来发展的趋势进行了分析和展望。
笔者有幸参加了我国自2007年以来多颗国产民用卫星的成像模拟技术研究工作, 针对不同卫星的特点, 开发了不同侧重点的卫星图像模拟模块, 并建立了一套具有我国自主知识产权的卫星遥感图像模拟技术体系。在过去6 a的研究过程中, 大致经历了4个阶段:
1)以注重格式模拟为主的图图模拟阶段。格式模拟是把各种卫星数据按照约定的格式编排在一起。陆地观测卫星下传的数据中除包含有效载荷的观测数据外, 还包括卫星有效载荷的状态参数、卫星平台的姿态参数、卫星的轨道参数、系统的时间参数和定标系统参数等数据。格式模拟数据能够在卫星发射前替代真实的卫星数据, 检验卫星数据处理系统各分系统数据接口的正确性。
2)以注重辐射模拟为主的库图模拟阶段。遥感科学过去40 a的发展中, 对电磁波在地表大气中传播机理的研究一直是个热点问题。库图模拟是以测量得到的各种地表参数作为输入, 使用遥感模型定量地计算地表-大气-传感器的辐射传输过程, 模拟得到具有物理意义的卫星图像。这种模拟得到的图像既可以作为卫星数据处理算法调试的测试数据, 又可以支持各种卫星遥感数据定量处理算法的预先研究。
3)以注重几何模拟为主的传感器模拟阶段。随着我国国产卫星的几何分辨率和几何定位精度的不断提高, 图像几何特性模拟的要求也在显著提高。影响卫星图像几何特性的因素主要有相机的内方位元素、相机在卫星上的安装精度、卫星平台的飞行轨道与姿态测量精度以及地表起伏的三维特性等。为了满足新的需求, 本文对与几何特性相关的计算模块进行了重新设计和研发, 包括卫星轨道模拟模块、偏流角模拟模块、考虑高程起伏的视线几何定位的迭代计算模块等。
4)以注重通用性为主的模拟平台研究阶段。目前我国在已有的针对具体卫星、具体传感器的专用的图像模拟系统基础上, 正在开发通用的模拟平台。专用的模拟系统只能针对某个特定的传感器进行模拟, 无法进行灵活扩展; 正在开发的通用模拟平台则采用了完全不同的技术, 能够在模型定义的有效区间内任意模拟地表和大气的光学特性, 具备对多种传感器通用的模拟能力。
自2007年以来, 笔者开发了一套完整的全链路模拟的技术框架和软件平台。该平台曾参与过我国HJ-1A/1B, ZY-3, ZY-1-02C, CBERS-03, CBERS-04和SJ-9A等7颗民用陆地观测卫星地面系统的研制过程, 已被验证为可靠的软件平台。
卫星图像的形成是一个非常复杂的过程, 包括成像几何的重采样过程、地表和大气的辐射传输过程、传感器的光电转换和数字化过程等。在模拟技术框架中, 主要的技术环节包括: ①卫星的轨道和姿态模拟; ②地表光学特性模拟; ③大气辐射传输模拟; ④传感器模拟; ⑤格式模拟。实现该技术路线有2种途径: ①将模拟技术路线分为库图模拟和图图模拟2类; ②将技术路线分为辐射模拟和几何模拟2类。全链路模拟技术的框架见图1。
绕过地表和大气的辐射传输模拟, 使用已有的波段相近、分辨率相近的卫星或者航空遥感图像模拟出目标卫星图像的技术途径是图图模拟; 而使用地表/大气辐射传输模型计算得到目标卫星图像的技术途径是库图模拟。辐射模拟主要是保证目标卫星图像波段的辐射信息的正确生成; 而几何模拟是通过模拟卫星的轨道和姿态获得传感器的外方位元素, 结合传感器的内方位元素, 计算出每个像元在成像时刻的地理坐标和观测角度。
笔者采用了2种技术路线来实现地表辐射特性的模拟: ①基于地物波谱库和地表分类图联合的方法模拟; ②基于叶片尺度及冠层尺度的辐射传输模型模拟。
第一种技术路线的实现过程是: 根据分类图读取当前像元的地表分类代码, 按照分类代码检索地物波谱库, 得到当前像元的波谱数据; 使用相机的波段响应函数卷积光谱数据得到当前像元的反射率; 重复以上步骤, 得到整景图像的反射率。第二种技术路线是将设定的叶片参数(叶绿素含量、干物质含量和含水量等)输入到PROSPECT模型中, 模拟得到叶片的反射率和透过率; 将观测的天顶角、方位角和太阳的天顶角、方位角, 当前像元的叶面积指数(LAI)和叶片的反射率及透过率输入到冠层模型SAIL模型中, 模拟得到当前点的反射率; 重复上述步骤, 得到整景图像的反射率(图2)。
国产卫星的大部分工作波段设在可见光和近红外波段, 而HJ-1卫星设有中红外和热红外波段。目前, 能满足可见光/近红外/热红外大气辐射传输模拟的大气辐射模式只有MODTRAN, 但其运算量巨大。通过多次运行MODTRAN得到查找表(look-up table, LUT); 利用近似的大气辐射传输模型实现像元级大气效应叠加。具体流程如图3所示。
1)大气辐射传输解析模型建立。大气辐射传输解析模型基于现有大气辐射传输理论, 在已有解析模型基础上, 充分考虑红外大气辐射传输的特性, 重新对大气辐射的贡献进行分解, 并构建相应的非线性解析方程。根据不同大气廓线(H2O, CO2, O3及温度等)和观测几何参数等, 利用MODTRAN4进行非线性解析模型的分析和验证。模型中还充分考虑了邻近像元对目标像元交叉辐射的贡献, 使该模型更为精细地表达了完整的大气辐射传输过程。
2)LUT建立与搜索。大气参数组合主要选取大气模式、气溶胶光学厚度、水汽含量和大气温度等; 观测几何参数主要有太阳天顶角、卫星天顶角及二者相对于目标像元的方位角。查表过程是从上至下的子表检索过程和自下至上的插值过程: ①根据输入参量值的大小, 在LUT中找到与之最接近的相邻子表; ②由查到的最底端子表逐级向上插值, 最终获得输入参量对应结果; 加入波段响应函数, 最终得到积分后的结果。波段响应函数可采用数学方法模拟得到[11]; 采用线性插值方法对各级子表进行插值计算, 最后进行波段响应函数积分。
3)交叉辐射贡献分析。由大气散射引起的相邻像元间的交叉辐射显著, 其作用可用点扩散函数(PSF)或者在频率域用调制传递函数(modulation transfer function, MTF)来描述。交叉辐射贡献主要受目标像元和邻近像元反射率和发射率差异、温度差异及大气条件的影响, 其计算需要将大气辐射传输解析模型与近似PSF模型结合起来, 共同确定交叉辐射实际影响半径, 从而确定交叉辐射贡献大小。
传感器成像几何特性模拟的技术流程如图4所示。
像元定位过程就是根据成像时的卫星轨道参数和传感器的姿态参数建立构像方程, 实现图像坐标到地心坐标的转换。引进地心直角坐标系既便于与传感器坐标系直接进行三维空间线性变换, 又能结合大地测量知识和地图投影方程将图像坐标转换为经纬度坐标, 进而转换为任意一种地图投影坐标。
在低分辨率卫星遥感中, 可以假设地表为理想椭球面; 但对于高分辨率遥感图像的模拟, 必须考虑地形的起伏, 需要假设地球表面的高程是已知的, 采用插值、迭代的方法计算(图5)。
取地球半径与DEM最大值之和为半径的圆球面作为初始高程面, 投影光线与该球面的交点A作为迭代初始点。根据A点坐标, 在DEM中插值计算对应的高程值, 得到B点。比较A点与B点的距离, 若A与B的距离小于某一阈值, 则满足收敛要求; 否则在投影光线OA上计算与B点高程值相等的C点, 将C点作为进一步迭代的初始点, 按照上面的方法继续迭代, 直到F点成为满足收敛要求的交点时为止。
至此, 完成了探测器像元位置到景物入射光在地心直角坐标系中的对应关系模拟。把地心直角坐标转换到地理坐标、地图投影坐标乃至图图模拟的大气层顶(top of atmosphere, TOA)图像的像元坐标, 为后继的传感器成像模拟流程中的重采样提供基础。
成像系统入瞳处的信号经过成像系统的退化作用后到达探测器, 需要考虑光学系统、探测器和采集电路对光谱和空间响应的过程。这些因素决定了最终模拟成像的光谱分辨率、空间分辨率和辐射分辨率等重要参数。
像元的空间响应可通过点扩散函数进行描述, 而通常情况下, 对于光学成像系统都采用MTF或PSF进行计算。
成像系统的MTF包括镜头光学系统、像增强、光纤面板、探测器及电子学系统等各个部件的MTF, 这些函数乘积就是综合的系统MTF, 即
MTFsys= MTFoptMTFelcMTF∏ MTFfoMTFCCD , (1)
式中: MTFsys为成像系统的MTF; MTFopt为镜头光学系统的MTF; MTFelc为电子学系统的MTF; MTF∏ 为像增强的MTF; MTFfo为光纤面板的MTF; MTFCCD为CCD探测器的MTF。
如果将光学成像系统的总体响应看做是一个数学变换, 光学成像系统的作用就是导致图像发生退化和质量下降的主因。本研究采用式(2)对未发生退化图像的像元辐照度L(x, y)进行变换, 得到经过成像系统作用后的退化图像, 即模拟生成图像
IFinal(x, y)= {FFT-1{[FFT(L(x, y)+σ p )]u, vMTFsys}+σ N}G+B , (2)
式中: (x, y)为空间坐标; (u, v)为频率坐标; FFT和FFT-1分别为傅立叶正、反变换; σ p和σ N分别为光子噪声和预置增益噪声; G为信号增益; B为直流偏置常数。传感器成像特性模拟的技术流程如图6所示。
除了相机获取的图像数据之外, 遥感卫星还有轨道、姿态信息, 热控部件、定标部件等辅助部件的工作状态信息以及系统时间等信息。在带宽有限的条件下, 为了满足各种信息下传的效率、正确性和保密要求等, 需要对各项星上数据进行格式设计和编排。而定标、卫星轨道和姿态、系统时间等信息是与遥感数据地面处理密切相关的, 因此图像模拟也要把这些相关的辅助数据模拟出来, 才能保证模拟数据的完整性, 以便用这样的模拟数据测试地面系统其他的功能模块。在图像模拟中, 笔者采用STK软件来模拟卫星的轨道、姿态和系统时间等参数, 并使用这些参数进行传感器几何特性模拟, 按照数据格式的接口文件说明将配套的模拟图像和辅助数据编排到符合格式要求的数据文件中, 完成卫星数据格式模拟。
针对本文的卫星图像模拟方案, 笔者设计了一个验证其模拟精度的实验, 采用敦煌辐射定标场附近的地面实测数据进行大气辐射传输模拟, 将得到的入瞳辐亮度与卫星实测的辐亮度进行对比。
地面实测光谱数据包括4种地物类型: 低反射率的水体、中等反射率的定标场和高反射率的棉花场地和干河床。测量仪器为ASD光谱仪和SVC光谱仪。将地面测量得到的光谱数据与HJ-1B的CCD相机的通道响应函数积分, 得到卫星4个波段对应地表的反射率值。图7为地面实测的几种地物的光谱曲线。
![]() | 图7 地面实测光谱与HJ-1B CCD对应的地表反射率Fig.7 Spectrum measured in field and corresponding surface reflectance of HJ-1B CCD |
HJ-1B CCD 4个波段对应的4种地物地表反射率见表1。
![]() | 表1 HJ-1B CCD 4个波段对应的地物反射率 Tab.1 Reflectances of corresponding ground features in four bands of HJ-1B CCD |
2009年8月21日同步测量时的气象视距为43.150 km, 大气柱水汽含量约1.58 g/cm2, 太阳天顶角为30.352° , 太阳方位角为153.776° ; 卫星观测天顶角为17.4° (干河床)。
2009年8月25日同步测量时的气象视距为45.201 km, 大气柱水汽含量约1.58 g/cm2, 太阳天顶角为31.307° , 太阳方位角为156.491° ; 卫星观测天顶角为13.02° (寿昌海水体)。
将上述参数代入卫星图像模拟系统中, 模拟的结果与卫星观测的结果对比见表2和图8。
![]() | 表2 HJ-1B CCD 4个波段对应的辐射亮度与模拟亮度结果对比 Tab.2 Comparison between corresponding radiance in four bands of HJ-1B CCD and simulation results(W· m-2· sr-1· μ m-1) |
![]() | 图8 HJ-1B CCD 4个波段对应的辐射亮度与模拟结果对比Fig.8 Comparison between corresponding radiance in four bands of HJ-1B CCD and simulation results |
辐射模拟的平均误差大约为10%; 模拟的最大误差为14.59%, 最小误差为5.76%。B1, B2和B3的模拟结果偏高, 而B4的模拟结果偏低。
1)本文回顾了过去6 a中笔者参与我国卫星图像模拟技术研究的历程, 较系统地介绍了所开发卫星图像模拟系统的主要技术框架及突破的关键技术。研究表明, 经过多年的积累, 我国建立的卫星图像模拟系统已经具备了地表-大气-传感器完整成像链路的模拟能力, 模拟的精度也在多次陆地观测卫星地面系统建设中得到了检验。
2)我国陆地观测卫星已进入密集发射期, 未来即将发射的新型传感器对卫星图像模拟提出了新的要求, 即需要更深入研究适应更高空间分辨率、更高光谱分辨率卫星载荷特点, 并具有更高几何精度和辐射精度的图像模拟技术。
1)大数据量的快速模拟。模拟图像的空间分辨率提高之后, 相同空间范围内的数据量会增大很多倍。如何提高现有的模拟速度, 在规定的时间内模拟出一轨数据, 是保证模拟技术能否满足地面系统建设调试的重要条件。现有的查找表算法是提高速度的有效途径, 将来还应该考虑如何利用通用图形处理器(GPGPU, general purpose GPU)硬件来提高图像模拟的运算速度。
2)地表-大气联合模拟。目前, 地表和大气的辐射传输过程被看作2个独立的过程, 是在模拟系统中分别使用2个独立的模块模拟计算的。实际上, 地表-大气之间的辐射传输过程是一个互相影响的迭代过程。在高空间分辨率情况下, 邻近像元通过大气的多次散射会影响到其周边的像元。这种高精度的模拟需要将地表-大气统一考虑, 尤其要模拟光子在大气和地表的多次反弹过程。
3)统一场景下的多波段图像模拟。使用统一场景实现多波段联合模拟, 其意义在于真正通过模拟实现出多个波段之间的逻辑关联, 赋予模拟图像以真实的物理意义, 促进多波段联合反演算法的开发和测试。在高分遥感技术发展中, 多波段和多分辨率传感器会同时存在, 实现统一场景的多波段模拟, 可以将图像模拟的应用范围进一步扩大到反演算法的测试中, 探索如何更好地利用高分数据开发多波段联合处理算法。
The authors have declared that no competing interests exist.
[1] |
|
[2] |
|
[3] |
|
[4] |
|
[5] |
|
[6] |
|
[7] |
|
[8] |
|
[9] |
|
[10] |
|
[11] |
|
[12] |
|
[13] |
|
[14] |
|
[15] |
|
[16] |
|
[17] |
|
[18] |
|
[19] |
|
[20] |
|
[21] |
|
[22] |
|
[23] |
|
[24] |
|
[25] |
|