气溶胶参数差异条件下对地观测辐射图像的快速仿真方法

摘要

本发明公开了一种气溶胶参数差异条件下对地观测辐射图像的快速仿真方法，包括：步骤1：选定标准气溶胶廓线参数；步骤2：查表法反演场景中各像素单元的气溶胶参数；步骤3：查询各像元对应的真实地表高程信息；步骤4：计算各像元的消光系数等效地表高度；步骤5：二流近似方法求解辐射传输方程，得到各类辐射的辐射矩阵；步骤6：获得各像元的光谱特性信息；步骤7：辐射矩阵插值方法仿真对地观测辐射图像。本发明图像仿真速度快，与传统方法相比，相同平台下仿真计算速度可提高10‑20倍；该方法的通用性好：所提出的仿真方法不受场景中气溶胶参数设定方法的限制，既适用于仿真真实在轨传感器的观测图像，也适用于仿真想定观测场景的红外辐射图像。

说明书

技术领域

本发明涉及对地观测辐射图像的仿真技术领域，具体涉及一种气溶胶参数差异条件下对地观测辐射图像的快速仿真方法。

背景技术

地球大气中，除了包含氮气(N₂)、氧气(O₂)、二氧化碳(CO₂)、氩气(Ar)等气体成分外，还含有悬浮的尘埃、液滴、冰晶等固体或液体颗粒，这些颗粒统称为气溶胶(参见文献[1]张建奇，方小平.红外物理[M].西安：西安电子科技大学出版社.2004.6.)。受各类气体及气溶胶的散射及吸收作用的影响，红外辐射在大气中传输将产生衰减。另一方面，通过理论计算仿真对地观测辐射图像，并同天基平台传感器的遥感实际结果进行比对，进而反演实测图像中的气溶胶空间分布特性及成分组成，在气象观测、农业遥感及空间监视等领域有重要应用。

特定观测几何关系条件下，通过求解辐射传输方程(参见文献[2]K.Stamnes,S.Tsay,and I.Laszlo.DISORT,a General-Purpose FORTRAN Program for Discrete-Ordinate-Method Radiative Transfer in Scattering and Emitting Layered Media:Documentation of Methodology[R/OL].ftp://climate.gsfc.nasa.gov/pub/wiscombe/Multiple_Scatt/.Mar.2000/Jul.2014)可以得到传输路径中某一位置的光谱辐出度E：

$\frac{\pm dE(\lambda ,v)}{dv}=E(\lambda ,v)-S(\lambda ,v)---\left(1\right)$

式中，表示光谱辐出度随传输路径的变化率；λ为波长；v为光学厚度；S为辐射源项。

图1所示为天基平台传感器接收辐射的基本构成示意图(参见文献[3]R.A.Schowengerdt.Remote Sensing:Models and Methods for Image Processing:Third Edition[M].USA:Academic Press.2007.)。如图1所示，传感器口面处的总光谱辐照度E_t由地表直接反射辐射E_r、大气路径散射辐射E_s、地表热辐射E_b及大气热辐射E_e四部分组成，即：

E_t(λ)＝E_r(λ)+E_s(λ)+E_b(λ)+E_e(λ)(2)

求解辐射传输方程式(1)后，可得各类辐射的计算表达式：

E_r(λ,θ_v,φ_v,θ_s,φ_s)＝E₀(λ)·τ_s(λ,θ_s)·ρ_r(λ,θ_v,θ_s,φ_v,φ_s)·τ_v(λ,θ_v)(3)

E_s(λ,θ_v,φ_v,θ_s,φ_s)＝E₀(λ)·ω(λ)·[1-τ_s(λ,θ_s)·τ_v(λ,θ_v)]·P(λ,θ_v,θ_s,φ_v,φ_s)(4)

E_b(λ,θ_v)＝τ_v(λ,θ_v)·ε_s(λ)·B(λ,T_b)(5)

E_e(λ,θ_v)＝[1-ω(λ)]·[1-τ_v(λ,θ_v)]·B(λ,T_a)(6)

式中，θ_v为观测天顶角；θ_s为太阳天顶角；为观测方位角；为太阳方位角；E₀为大气顶层的太阳光谱辐照；τ_s为入射路径的光谱透过率；τ_v为观测路径的光谱透过率；ρ_r为地表的双向光谱反射率；ω为大气的单次光谱反照率；P为大气散射相函数；ε_s为地表的光谱发射率；B为黑体光谱辐出度；T_b为地表温度；T_a为大气温度。

若辐射传输路径上的气体浓度分布不均匀，则可以采用Curtis-Godson近似(参见文献[4]A.R.Curtis.Discussion of a Statistical Model for Water Vapour Absorption[J].Quarterly Journal of the Royal Meteorological Society.1952,78:638–640和文献[5]W.L.Godson.The Evaluation of Infrared-Radiative Fluxes due to Atmospheric Water Vapour[J].Quarterly Journal of the Royal Meteorological Society.1953,79:367–379.)将非均匀大气等效为若干均匀子层的总和。如图2所示，以大气顶层为光学深度零点，依海拔高度将大气划分为均匀的n层。图中q_i为大气顶层至第i子层的光谱光学深度；v_i为大气顶层至第i子层的光谱光学厚度；h_i为第i子层的高度；T_i为第i子层的温度；k_i为第i子层的光谱消光系数；E⁺表示向上传输的光谱辐出度；E^-表示向下传输的光谱辐出度。

根据定义，消光系数k_i与光学深度q_i有如下关系：

$q_i\left(\lambda \right)=\underset{j=0}{\overset{i}{\Sigma }}{k}_j\left(\lambda \right)---\left(7\right)$

消光系数k_i与光学厚度v_i的关系定义如下：

$v_i\left(\lambda \right)=\underset{j=0}{\overset{i}{\Sigma }}[l_j·k_j\left(\lambda \right)]---\left(8\right)$

式中，l_j为第j子层内的传输路径长度；k_j为第j子层的光谱消光系数。

式(3)至式(6)中的大气透过率τ与光学厚度v为e指数关系：

τ＝exp(-v)(9)

故将式(8)代入式(9)可得入射路径光谱透过率τ_s及观测路径光谱透过率τ_v的计算式：

${\tau }_s(\lambda ,{\theta }_s)=\exp \{-\underset{i=0}{\overset{n}{\Sigma }}[l_i\left({\theta }_s\right)·k_i\left(\lambda \right)\left]\right\}---\left(10\right)$

${\tau }_v(\lambda ,{\theta }_v)=\exp \{-\underset{i=0}{\overset{n}{\Sigma }}[l_i\left({\theta }_v\right)·k_i\left(\lambda \right)\left]\right\}---\left(11\right)$

第i子层的光谱消光系数k_i由三部分组成：

k_i(λ)＝k_a(λ)+k_r(λ)+k_m(λ)(12)

式中，k_a为大气吸收消光系数，简称吸收系数，可采用逐线积分(参见文献[6]S.A.Clough,M.W.Shephard,E.J.Mlawer,et>r为大气分子的散射消光系数，简称分子散射系数，由Rayleigh散射计算方法(参见文献[9]A.A.Bodhaine,N.B.Wood,E.G.Dution>m为气溶胶的散射消光系数，简称气溶胶散射系数，由Mie散射计算方法(参见文献[10]A.Matzler.MATLAB>

大气散射系数定义为分子散射系数与气溶胶散射系数的总和，故式(4)与式(6)中的大气单次光谱反照率ω有定义式：

$\omega \left(\lambda \right)=\frac{k_r\left(\lambda \right)+k_m\left(\lambda \right)}{k_i\left(\lambda \right)}---\left(13\right)$

同理，由于大气散射包含气体分子散射与气溶胶散射两部分，式(4)中的大气散射相函数P与Rayleigh散射相函数P_r及Mie散射相函数P_m有关，即：

$P(\lambda ,{\theta }_v,{\theta }_s,{\phi }_v,{\phi }_s)=\frac{P_r({\theta }_v,{\theta }_s,{\phi }_v,{\phi }_s)·k_r\left(\lambda \right)+P_m(\lambda ,{\theta }_v,{\theta }_s,{\phi }_v,{\phi }_s)·k_m\left(\lambda \right)}{k_r\left(\lambda \right)+k_m\left(\lambda \right)}---\left(14\right)$

其中，Rayleigh散射相函数P_r的计算式为：

式中，为散射角，由立体角公式可知：

由于精确计算Mie散射相函数P_m过于耗时，故采用Heney-Greenstein函数近似计算Mie散射相函数P_m，即：

式中，g为非对称因子，由Mie散射计算方法计算；为散射角。

在已知大气气体成分组成及气溶胶的粒子浓度、尺度分布及复折射指数等物理参数的条件下，通过Rayleigh散射及Mie散射计算方法，容易求得大气各子层内的光学参数，如消光系数、散射系数及非对称因子等(参见文献[11]E.P.Shettle and R.W.Fenn.Models for the aerosols of the Lower Atmosphere and the Effects of Humidity Variations on Their Optical Properties[R].AFGL-TR-79-0214.Sep.1979.)。利用求得的光学参数及计算式(2)至(6)，可以计算传感器接收的总辐射。

与本发明相关的现有技术介绍如下：

1.1现有技术一的技术方案

对于某一确定观测场景，在确定了各条传输路径的大气气体成分组成及气溶胶参数后，采用大气辐射传输代码逐条路径求解辐射传输方程式(1)，解得仿真图像中各个像元的辐射强度，并最终仿真得到红外辐射图像。基于物理机理的红外场景仿真模型通常采用此方案(参见文献[12]W.M.Cornette,J.M.Alfred,and J.M.Goldspiel.Moderate Spectral Atmospheric Radiance and Transmittance Code(MOSART v2.1)User Reference Manual[R].MOSART-NRL-DOC-URM-001-v2.1-111206.Dec.2011.和文献[13]L.Labarre,K.Caillault,S.Fauqueux,et al.An Overiew of Matisse-v2.0[C].Optics in Atmospheric Propagation and Adaptive Systems XIII,Proceedings of SPIE,2010,7828:1-10.)。

1.2现有技术一的缺点

当仿真天基平台对地观测传感器的观测图像时，由于该类传感器条带的覆盖区域很广，观测图像包含的像元数目很多，逐个像素单元计算辐射图像的方法耗时长。且图像仿真过程中，受计算机内存容量的限制，能够完成的仿真图像场景尺寸十分有限。

2.1现有技术的技术方案二

文献[14](X.He and X.Xu.Fast Calculation of Scattered Radiance in Multispectral Imagery Simulation.SPIE Infrared Remote Sensing and Instruments XXIII,2015,9608(96080O):1-11.)提出一种基于辐射矩阵插值的对地观测辐射图像仿真方法。当场景中的大气参数设定相同时，大气散射辐射及地表反射辐射随观测天顶角与太阳天顶角均匀变化，此时只需要计算几个特定像元的各类辐射强度，并分别构成辐射矩阵，接着通过插值方法快速得到其他像元的各类辐射。插值方法的计算时间取决于选取的像元数目而非仿真场景的原始像元数目，故该方案的计算效率较方案一有较大提升，而计算精度略有下降。

2.2现有技术二的缺点

仿真中假定了场景中各个像元的大气条件相同，这一假设并不适用于图像条带覆盖几千公里的天基平台对地观测传感器。且该方案中尚未考虑地表高度及气溶胶气象视程对传感器 接收辐射强度的影响，仿真图像的纹理(辐射强弱变化)不能反映气溶胶的空间分布特性。

发明内容

本发明所要解决的技术问题为：针对现有技术以上方案的不足之处，本发明提出一种考虑气溶胶参数差异条件时，对地观测辐射图像的快速仿真方法。基于传感器像元的辐射强度与气溶胶浓度及传输距离的关系，采用消光系数等效地表高度法仿真天基平台传感器对地观测辐射图像，实现条带覆盖区域很广的传感器辐射图像快速仿真。

本发明采用的技术方案为：一种气溶胶参数差异条件下对地观测辐射图像的快速仿真方法，包括如下步骤：

步骤1：选定标准气溶胶廓线参数；

步骤2：查表法反演场景中各像素单元的气溶胶参数；

步骤3：查询各像元对应的真实地表高程信息；

步骤4：计算各像元的消光系数等效地表高度；

步骤5：二流近似方法求解辐射传输方程，得到各类辐射的辐射矩阵；

步骤6：获得各像元的光谱特性信息；

步骤7：辐射矩阵插值方法仿真对地观测辐射图像。

其中，步骤1中所述的选定标准气溶胶廓线参数具体步骤包括：选定气溶胶类型与气象视程，并计算随大气高度变化的消光系数及散射系数廓线。

其中，步骤2中所述的查表法反演场景中各像素单元的气溶胶参数具体步骤包括：计算不同类型气溶胶及气象视程条件下的像元色比，得到色比表格；在色比表格中查找与实测图像色比最接近的值，该色比值对应的气溶胶参数即为反演结果。

其中，步骤3中所述的查询各像元对应的真实地表高程信息具体步骤包括：通过几何关系计算各像元地表投影位置的经纬度坐标；通过查询地表高程信息图，获得像元对应的地表高度。

其中，步骤4中所述的计算各像元的消光系数等效地表高度具体步骤包括：计算各像元实际的光学厚度，依据步骤1中指定的气溶胶标准廓线求解等效地表高度，使得各像元的光学厚度保持不变。

其中，步骤5中所述的二流近似方法求解辐射传输方程得到各类辐射的辐射矩阵具体步骤包括：采用步骤1中选定的标准气溶胶廓线参数，利用二流近似方法计算选定像元的各类辐射；计算所得各类辐射的强度随太阳天顶角及观测天顶角变化，故构成辐射强度随太阳天顶角及观测天顶角变化的辐射矩阵。

其中，步骤6中所述的获得各像元的光谱特性信息具体步骤包括：根据观测几何关系计 算各像元的地表投影位置的经纬度坐标，通过查询地表覆盖物分类图可以确定地物类别；依据地物类别查询地物光谱信息数据库可以得到地物的光谱反射率及光谱发射率廓线。

其中，步骤7中所述的辐射矩阵插值方法仿真对地观测辐射图像具体步骤包括：计算各像元的太阳天顶角及观测天顶角，利用步骤4中的等效地表高度及步骤5中得到的辐射矩阵计算各像元的辐射。

与现有对地观测红外辐射图像仿真方法相比，本发明提出的对地观测红外辐射图像仿真方法具有以下优点：

(1)、本发明图像仿真速度快：通过提出消光系数等效地表高度的概念并计算得到等效地表高度图，可实现包含多类气溶胶及多种气象视程下对地观测场景红外辐射场景的图像仿真，与传统方法相比，相同平台下仿真计算速度可提高10-20倍；

(2)、本发明方法的通用性好：所提出的仿真方法不受场景中气溶胶参数设定方法的限制，既适用于仿真真实在轨传感器的观测图像，也适用于仿真想定观测场景的红外辐射图像。

附图说明

图1为天基平台传感器接收辐射的基本组成示意图，其中，图1(a)为地表直接反射，图1(b)为大气路径散射，图1(c)地表热辐射，图1(d)为大气热辐射；

图2为Curtis-Godson近似子层划分示意图；

图3为考虑各像元气溶胶参数差异时的对地观测辐射图像快速仿真流程图；

图4为消光系数等效地表高度的计算流程；

图5为消光系数快速计算方法验证举例；

图6为散射系数快速计算方法验证举例；

图7为气溶胶气象视程反演举例，其中，图7(a)为实测色比图像，图7(b)为反演所得气象视程图像；

图8为MODIS谱段7消光系数等效地表高度计算举例，其中，图8(a)为场景实际高程图像，图8(b)为消光系数等效地表高度图像；

图9为2015年1月19日，MODIS谱段7辐射图像仿真举例，其中，图9(a)为MODIS实测辐射图像，图9(b)为仿真辐射图像。

具体实施方式

下面结合附图以及具体实施方式进一步说明本发明。

本发明的技术方案框图如图3所示，基本技术原理如下。

步骤1：选定标准气溶胶廓线参数。

现有的大气辐射传输代码中，常用的对流层气溶胶模型有四类，分别为乡村型、城市型、 海洋型及对流层型。各类气溶胶的粒子组成、粒子尺度分布及复折射指数不同，故四类典型气溶胶模型的光谱消光系数k_i，光谱散射系数k_m及非对称因子g均不相同。

在实际工程应用中，确切知道气溶胶粒子的尺度分布非常困难，通常采用气象视程估算气溶胶的粒子浓度。根据气象视程的定义，若某个谱段内的大气光谱吸收系数很小，近似为零，则气象视程V与大气散射系数有如下关系：

$V=\frac{3.91}{k_r\left({\lambda }_0\right)+k_m\left({\lambda }_0\right)}---\left(18\right)$

式中，λ₀为选定谱段的波长，通常选取为λ₀＝0.61μm或0.55μm。

对于某类气溶胶，其光谱散射系数k_m与粒子浓度N有关，有估算式：

$k_{m-s}\left(\lambda \right)=\frac{N_s·k_{m-0}\left(\lambda \right)}{N_0}---\left(19\right)$

式中，N₀为标准条件下的气溶胶粒子密度；k_m-0为采用标准气溶胶粒子密度N₀计算所得的气溶胶光谱散射系数；N_s为实际的气溶胶粒子密度；k_m-s为估算得到的气溶胶光谱散射系数。

由式(18)及式(19)可知，气象视程V与气溶胶浓度N存在一定联系，可以采用气象视程的远近反映气溶胶浓度的大小，则气溶胶浓度的空间分布可以由像元气象视程图像表示。此外，可以为各类气溶胶预先计算某一确定气象视程条件下的光谱散射系数。对于其他气象视程的观测情形，可以采用式(19)中的比例关系近似计算气溶胶散射系数，节约计算时间。

步骤2：查表法反演场景中各像素单元的气溶胶参数。

天基平台对地观测传感器的条带覆盖范围很广，通过地面气象站获取各像元的气溶胶光学特性十分困难，故通常采用对地观测传感器的遥感数据反演气溶胶的光学参数。

气溶胶光学参数反演利用了两个谱段内实测图像的辐亮度比值，即色比。定义色比R为：

$R=\frac{{\mathrm{\Delta \lambda }}_1·E_t\left({\lambda }_1\right)}{{\mathrm{\Delta \lambda }}_2·E_t\left({\lambda }_2\right)}---\left(20\right)$

式中，E_t为传感器接收的总光谱辐照度，由式(2)计算；λ₁为谱段1的中心波长；Δλ₁为谱段1的谱宽；λ₂为谱段2的中心波长；Δλ₂为谱段2的谱宽。

对于位于可见光、近红外或短波红外光谱区间内的观测谱段，热辐射相较于散射辐射很弱，可以忽略。即忽略式(2)中的地表热辐射E_b及大气热辐射E_e，将传感器接收的光谱辐射照度E_t重新书写如下(为简化表达，省略几何参数θ_v、θ_s、及)：

E_t(λ)＝E₀(λ)·τ_s(λ)·ρ_r(λ)·τ_v(λ)+E₀(λ)·ω(λ)·[1-τ_s(λ)·τ_v(λ)]·P(λ)(21)

若视场中的气溶胶浓度增加，使得大气的透过率τ下降，即：

τ′(λ)＝τ(λ)-Δτ(λ)(22)

式中，τ′为气溶胶浓度增加后的大气透过率；Δτ为透过率的减少量。

将式(22)代入式(21)，光谱辐射照度E_t发生变化，有：

E_t′(λ)＝E_t(λ)+ΔE_t(λ)(23)

式中，E_t′为气溶胶浓度增加后的光谱辐射照度；ΔE_t为光谱辐射照度的变化量，有表达式：

ΔE_t(λ)＝E₀(λ)·Δτ(λ)·[ω(λ)·P(λ)-ρ_r(λ)](24)

在某一特定谱段内，场景中任意一像元的E₀,τ,ρ_r,ω及P在观测时间内保持不变。若瞬时视场中不存在云或雾，则可以认为该像元的地表反射率ρ_r大于大气反射率ωP。故随着气溶胶浓度的增加，谱段内的辐射照度减弱。

利用式(20)与式(23)可得气溶胶浓度变化后，像元色比R′的计算式：

$R^{\prime }=\frac{{\mathrm{\Delta \lambda }}_1·{E_t}^{\prime }\left({\lambda }_1\right)}{{\mathrm{\Delta \lambda }}_2·{E_t}^{\prime }\left({\lambda }_2\right)}=\frac{{\mathrm{\Delta \lambda }}_1·[E_t\left({\lambda }_1\right)+{\mathrm{\Delta E}}_t\left({\lambda }_1\right)]}{{\mathrm{\Delta \lambda }}_2·[E_t\left({\lambda }_2\right)+{\mathrm{\Delta E}}_t\left({\lambda }_2\right)]}---\left(25\right)$

由于大气顶层的太阳光谱辐照度E₀与气溶胶散射系数k_m在波长大于0.45μm的光谱区间内随波长的增长而减小。故选取两个位于该波长区间内的谱段，若λ₁<λ₂，则有E₀(λ₁)>E₀(λ₂),k_m(λ₁)>k_m(λ₂)。故谱段1的背景辐射变化量ΔE_t(λ₁)大于谱段2的背景辐射变化量ΔE_t(λ₂)。由于辐射照度的变化量ΔE_t通常为负值，故式(25)中的色比值降低。

若某瞬时视场中气溶胶的浓度很低，大气透过率τ趋近于1，则该像元的色比取得极值，且由两个谱段内的地表反射辐出度决定，即：

$R_{max}=\frac{E_0\left({\lambda }_1\right)·{\rho }_r\left({\lambda }_1\right)·{\mathrm{\Delta \lambda }}_1}{E_0\left({\lambda }_2\right)·{\rho }_r\left({\lambda }_2\right)·{\mathrm{\Delta \lambda }}_2}---\left(26\right)$

由以上分析可知，对于某一像元，当视场中几乎不存在气溶胶时，色比取得极大值；而随着视场中气溶胶浓度的增加，色比值逐渐减小。从理论上说，可以取任意两个位于大气透过窗口中的谱段，利用其实测数据反演气溶胶光学参数。

若采用MODIS反演得到的气溶胶光学参数数据作为仿真场景的气溶胶参数输入，则可以略过该步骤。略过步骤2不会影响后续的图像仿真。

步骤3：查询各像元对应的真实地表高程信息。

地球表面某一区域的平均海拔高度随观测时间变化比较缓慢，因此可以利用几何关系计算各像元地表投影的地理坐标(经纬度坐标)，通过查询地表高程信息图，获得各像元的地表高程信息。

步骤4：计算各像元的消光系数等效地表高度。

由式(2)至式(6)可知，若将地球大气视为整体，且观测几何关系确定，则传感器接收的 总光谱辐照度E_t仅受大气透过率τ的影响。由透过率的计算式(10)与式(11)可知，大气透过率不仅受气溶胶浓度的影响(式中的消光系数k_i)，还受传输路程长度的影响(式中的路径长度l_i)。因此，从数值计算结果而言，气溶胶浓度增加或传输路程长度增加对传输路径透过率的影响相同，所引起的辐射强度变化也理应相同。换言之，理论上存在换算关系，当气溶胶浓度增大引起传感器接收的总辐射减弱时，可以保持气溶胶的浓度不变，增加适当的传输距离使得传感器接收的总辐射产生相等数值的衰减。

气溶胶浓度随海拔高度z的变化可由e指数函数近似：

N(z)＝N(0)exp(-z/z₀)(27)

式中，N(z)为距离地表高度为z位置的气溶胶浓度；z₀为标高，是与气象视程有关的常数；N(0)为地表位置的气溶胶浓度。

由步骤1中的式(19)可知，气溶胶浓度与消光系数呈正比例关系。由式(27)可知，气溶胶浓度在低海拔地区的更高，故观测地表高度较低的区域，其观测路径的光学厚度更大，而光学厚度越大，透过率越低。因此，在选定某一组气溶胶参数后，以该组参数的消光系数垂直廓线为参考，计算各传输路径对应的等效地表高度，使得各传输路径的总光学厚度与利用步骤2中气溶胶参数求得的总光学厚度相同。消光系数等效地表高度的计算流程如图4所示。

图4中传感器至第i子层的大气光学厚度v_i由式(8)计算，实测光学厚度v为仿真场景的输入参数。如图3所示，大气顶层为h₀位置，故有消光系数等效高度h_e的计算式：

$h_e=h_{i-1}-\frac{v-v_i}{k_i·l_i}·(h_{i-1}-h_i)---\left(28\right)$

式中，h_i为第i子层至地表的高度，且h_i<h_i-1。

步骤5：二流近似方法求解辐射传输方程，得到各类辐射的辐射矩阵。

文献[14]中给出了球面平行大气条件下，二流近似求解式(1)的基本方法，并分别为地表热辐射、地表散射辐射、大气热辐射及大气散射辐射建立了辐射矩阵。由于各类型辐射的强度随观测天顶角及太阳天顶角变化平缓，故可以采用插值方法仿真背景图像。

应当指出，采用其他方法求解式(1)同样能够得到各辐射的辐射矩阵。

步骤6：获得各像元的光谱特性信息。

光谱特性信息主要指各像素单元的光谱反射率和光谱发射率，由地物的物理特性决定。与地表高程信息类似，可以利用各像元的地理坐标，通过查询地物覆盖物分类图及地物光谱特性数据库得到。

步骤7：辐射矩阵插值方法仿真对地观测辐射图像。

由步骤4可知，利用消光系数等效地表高度，在生成辐射矩阵的过程中能够继续采用单 一气溶胶参数廓线的假设。在生成背景图像的过程中，利用等效后的地表高度反映气溶胶的空间分布特性。

实施例：

下面以仿真MODIS实测图像为例对本发明作进一步说明。目标场景为2015年1月19日Aqua-MODIS于5:30a.m至5:35a.m的探测区域。

步骤1：利用式(19)中的比例关系快速计算气溶胶的光谱消光系数及光谱散射系数。

图5所示为相对湿度0.7的乡村型气溶胶，气象视程2km条件下，近似计算气象视程10km、21km及50km条件下的光谱消光系数，在MODIS谱段3(0.459-0.479μm)内的结果比较图。图中点划线为Mie散射计算方法的精确计算结果，虚线为式(19)的近似计算结果。

由图5可知，Mie散射的精确计算曲线与式(19)的近似计算曲线几乎重合。分析其他气溶胶类型、相对湿度及观测谱段的计算结果后，所得结论相同。

如图6所示为采用图5中相同参数计算得到的光谱散射系数比对图。与图5中的结论类似，Mie散射的精确计算廓线与式(19)的近似计算廓线几乎重合。在其他仿真参数条件下的计算结果验证了该结论。

由式(13)可知，光谱单次反照率由光谱散射系数与光谱消光系数将决定。故采用步骤1中的方法可以节约图像仿真时间。

步骤2：由步骤2中的理论推导可知，理论上可以采用任意两个谱段反演气溶胶的气象视程。

这里采用MODIS谱段1(0.620-0.670μm)与谱段3(0.459-0.479μm)的色比反演气象视程，其结果如图7所示。

图7(a)MODIS实测图像谱段1与谱段3的色比；图7(b)为图7(a)反演所得的气溶胶气象视程图像，色标单位为km。由于MODIS谱段1的中心波长小于MODIS谱段3，因此色比值小的区域气象视程大。

步骤4：由反演得到的气溶胶气象视程可以计算消光系数等效地表高度。由理论推导可知，由于消光系数在各个观测谱段内不同，因此，各谱段的等效地表高度图像应该存在差异。

图8所示为采用乡村型气溶胶，气象视程2km为标准气溶胶廓线，在MODIS谱段7(2.105-2.155μm)内计算得到的地表等效高度图像，色标单位为km。图8(a)所示为该场景区域的实际地表高程图像，图8(b)为该场景在谱段7内的等效地表高度图像。从图中可以看出，等效高度图不仅反映了真实的地理高程信息，还包含了场景中气溶胶的空间分布信息。由于随着海拔高度增加，气溶胶浓度降低，因此气象视程距离远的像元路径透过率高，而路径透过率高其对应的等效高度也将越高。

步骤7：利用等效地表高度图及辐射矩阵插值方法能够快速仿真对地观测辐射图像。图9所示为仿真辐射图像与MODIS谱段7的实测辐射图像对比结果，图像色标单位为W·m^-2·μm^-1·sr^-1。图9(a)为实测辐射图像；图9(b)为未考虑云层的仿真辐射图像。由比对结果可知，仿真辐射图像的纹理反映了气溶胶的空间分布特性。

本发明还可以采用的替代方案同样能完成发明目的：由式(7)、式(8)及式(9)可知，消光系数与光学深度、光学厚度以及透过率之间均存在换算关系，因此本项发明所提出的方法不限于采用消光系数、也可采用光学深度、光学厚度或透过率等参数计算对应的等效高度，实现红外场景的快速仿真。采用其他参数时，其仿真计算方法在本质上与本发明中的技术方案相同。