一种实现浑浊水体区域MODIS影像大气校正的方法

技术领域

本发明涉及一种实现浑浊水体区域MODIS(Moderate-Resolution Imaging  Spectroradiometer)影像大气校正的方法，属于水色定量遥感理论与应用技 术领域，适用于水色定量遥感的理论方法和应用技术研究。

背景技术

标准MODIS大气校正算法称之为“清洁”水体大气校正算法。该算法的 基本假设在于研究区域内存在“清洁”水体(叶绿素a浓度＜0.5-1.0μg/l， 悬浮浓度几乎0)，则MODIS近红外两个波段的离水反射率近似为零，但是 由于大气散射作用，使得传感器探测到的信号不为零，因此可以将这部分信 号归根为大气程辐射的贡献，并用于大气信息的提取。标准MODIS大气校正 算法就是利用“清洁”水体的这种光学性质进行大气校正的。因此，标准 MODIS大气校正算法成败的关键在于是否能够在研究区域内找到一种为“清 洁”水体的像元。

众所周知，中国海岸带海域，尤其是黄河口、长江口和珠江口等三角洲 地带的水体，黄河、长江和珠江携带大量的泥沙入海，以致入海口很大面积 的水域含沙量非常之高。研究结果表明，在近红外波段，水体的离水反射率 与悬浮泥沙浓度成正比，哪怕水体中仅含微量的悬浮泥沙，都将引起近红外 波段离水辐亮度成倍的增长。这就意味着标准MODIS大气校正算法在这些区 域将会失败。要想在浑浊的II类水体中实现大气精校正，则需要更多的先 验知识参与大气校正过程。

发明内容

本发明的目的在于提供一种实现浑浊水体区域MODIS影像大气校正的方 法，以解决现有技术中标准MODIS大气校正算法在浑浊水体区域无法适用的 缺陷，在浑浊的II类水体中实现大气精校正。

本发明一种实现浑浊水体区域MODIS影像大气校正的方法，具体所采用 的技术方案，包括如下步骤：

步骤(1)：读入531、551、667和678nm波段实测的标准离水反射率之 间的经验关系、Rayleigh散射查找表、模型所涉及波段的MODIS影像、同步 大气压强、风速和臭氧厚度数据。其中标准离水反射率之间的关系如下：

[ρ_w]_N(λ₂)≈a₁[ρ_w]_N(λ₁)

[ρ_w]_N(λ₄)≈a₂[ρ_w]_N(λ₃)

式中，x_m、x_n和y_m为经验常数，可以通过回归计算得到。

步骤(2)：在步骤(1)涉及的气象数据的基础上，利用Rayleigh散射 查找表算法，计算MODIS各波段Rayleigh散射贡献的反射率，并加以订正。

ρ_t(λ)＝ρ_r(λ)+ρ_a(λ)+ρ_ra(λ)+t(λ)[ρ_w]_N(λ)

ρ_at(λ)＝ρ_path-ρ_r(λ)＝ρ_a(λ)+ρ_ra(λ)+t(λ)ρ_n(λ)

式中，ρ_r(λ)为Rayleigh散射贡献的反射率；ρ_a(λ)为气溶胶散射贡献 的反射率；ρ_at(λ)为经过Rayleigh散射校正后的表观反射率；ρ_ra(λ)为气溶 胶和大气分子多次散射贡献的反射率，该项可以解释光子在大气分子和气溶 胶颗粒之间多次穿梭的结果；ρ_n(λ)为待求的表观反射率；t(λ)为漫射透 过系数，其表达式如下：

$t\left(\lambda \right)=\exp \{-[0.5{\tau }_r\left(\lambda \right)+{\tau }_{\mathrm{oz}}\left(\lambda \right)\left](\frac{1}{\cos {\theta }_{\mathrm{sat}}}+\frac{1}{\cos {\theta }_{\mathrm{sol}}})\right\}$

式中，θ_sat为卫星天顶角；θ_sol为太阳高度角；τ_r(λ)为Rayleigh散射 光学厚度。Rayleigh散射光学厚度与气压之间存在如下关系：

${\tau }_r\left(\lambda \right)=\frac{P_0}{P}{\tau }_{r,0}\left(\lambda \right)$

式中，τ_r，0(λ)为标准大气压下的Rayleigh散射的光学厚度，可从 Rayleigh散射查找表中获取；P₀为标准大气压，其值取1013.25mb；P为像 元所在位置的实际大气压强。臭氧光学厚度可以用下式表示：

${\tau }_{\mathrm{oz}}\left(\lambda \right)={\tau }_{\mathrm{oz},0}\left(\lambda \right)\frac{\mathrm{DU}}{1000}$

式中，DU为单位臭氧浓度。

步骤(3)：利用Rayleigh散射修正后的MODIS影像数据，构建531、551、 667和678nm波段的反射率之间的经验关系。

Rayleigh散射修正后的MODIS影像在531、551、667和678nm波段的反射 率之间的经验关系形式如下：

ρ_at(λ₂)＝x_mρ_at(λ₁)+y_m

ρ_at(λ₄)＝x_nρ_at(λ₃)

步骤(4)：在步骤(3)的基础上，利用531、551、667和678nm波段反射 率的经验关系，结合Angstrom气溶胶指数模型，计算MODIS各波段的气溶胶 散射和气溶胶-Rayleigh多次散射贡献的反射率。

将步骤(1)所示的经验关系，代入ρ_at(λ)的表达式，并化简可得：

${\rho }_{\mathrm{at}}\left({\lambda }_2\right)=a_1\frac{t\left({\lambda }_2\right)}{t\left({\lambda }_1\right)}{\rho }_{\mathrm{at}}\left({\lambda }_1\right)+{\rho }_{\mathrm{as}}\left({\lambda }_2\right)-a_1\frac{t\left({\lambda }_2\right)}{t\left({\lambda }_1\right)}{\rho }_{\mathrm{as}}\left({\lambda }_1\right)$

${\rho }_{\mathrm{at}}\left({\lambda }_4\right)=a_2\frac{t\left({\lambda }_4\right)}{t\left({\lambda }_3\right)}{\rho }_{\mathrm{at}}\left({\lambda }_3\right)+{\rho }_{\mathrm{as}}\left({\lambda }_4\right)-a_2\frac{t\left({\lambda }_4\right)}{t\left({\lambda }_3\right)}{\rho }_{\mathrm{as}}\left({\lambda }_3\right)$

Angstrom气溶胶指数模型ε(λ_s，λ₁)如下：

$ϵ(\lambda ,{\lambda }_1)\approx \frac{{\rho }_{\mathrm{as}}\left(\lambda \right)}{{\rho }_{\mathrm{as}}\left({\lambda }_1\right)}\approx \mathrm{Exp}\left[n({\lambda }_1-\lambda )\right]$

式中，n为Angstrom指数；λ₁为参考波段。将ε(λ_s，λ₁)代入本步骤 中的ρ_t(λ₂)和ρ_t(λ₄)的表达式可得：

$R_{\mathrm{at}}\left({\lambda }_2\right)=\frac{a_{1}t\left({\lambda }_2\right)}{t\left({\lambda }_1\right)}{R}_{\mathrm{at}}\left({\lambda }_1\right)+[ϵ({\lambda }_2,{\lambda }_1)-\frac{a_{1}t\left({\lambda }_2\right)}{t\left({\lambda }_1\right)}ϵ({\lambda }_1,{\lambda }_1)]R_{\mathrm{as}}\left({\lambda }_1\right)$

$R_{\mathrm{at}}\left({\lambda }_4\right)=\frac{a_{2}t\left({\lambda }_4\right)}{t\left({\lambda }_3\right)}{R}_{\mathrm{at}}\left({\lambda }_3\right)+[ϵ({\lambda }_4,{\lambda }_1)-\frac{a_{2}t\left({\lambda }_4\right)}{t\left({\lambda }_3\right)}ϵ({\lambda }_3,{\lambda }_1)]R_{\mathrm{as}}\left({\lambda }_1\right)$

式中，ρ_as(λ)为ρ_ra(λ)和ρ_a(λ)的和。将步骤(3)中的经验关系代入上 式可得参考波段的ρ_as(λ)和n的计算方程组：

$n=\frac{\ln x_n}{{\lambda }_3-{\lambda }_4}$

${\rho }_{\mathrm{as}}\left({\lambda }_1\right)=\frac{y_m}{ϵ({\lambda }_2,{\lambda }_1)-x_m}$

步骤(5)：在步骤(2)、(4)的基础上计算大气的程辐射贡献的发射率， 并对MODIS影像进行大气校正，计算公式如下：

${\left[{\rho }_w\right]}_N\left(\lambda \right)=\frac{{\rho }_{\mathrm{at}}\left(\lambda \right)-ϵ(\lambda ,{\lambda }_1){\rho }_{\mathrm{as}}\left({\lambda }_1\right)}{t\left(\lambda \right)}.$

本发明一种实现浑浊水体区域MODI S影像大气校正的方法，其优点及功 效在于：本发明以波段之间的经验关系为基础，结合查找表Rayleigh散射计 算方法，研究与探讨适用于混杂II水体的MODIS大气校正算法。该算法克服 了浑浊I I类水体大气校正难的问题，有助于提高浑浊水体的大气校正精度， 进而有利于推动水色遥感的发展，具有重要的理论和应用意义。

附图说明

图1为本发明的气溶胶散射贡献的反射率。

图2为本发明的各波段大气校正误差。

具体实施方式

为了更好的说明本发明涉及的一种实现浑浊水体区域MODIS影像大气 校正的方法，利用ASD野外地物光谱仪分别测量得到黄河口、长江口和太湖 的光谱数据。2003年10月27日和28日，在太湖进行表光光学量测量的同 时，收集了MODIS影像数据各一景，用于算法分析计算和精度评估。本发明 涉及一种实现浑浊水体区域MODIS影像大气校正的方法，具体实现步骤如 下：

步骤(1)：读入531、551、667和678nm波段实测的标准离水反射率之 间的经验关系、Rayleigh散射查找表、模型所涉及波段的MODIS影像、同步 大气压强、风速和臭氧厚度数据。

臭氧厚度数据、风速和大气压强同步气象数据可以从NASA的官方网站 ftp://oceans.gsfc.nasa.gov下载得到。531、551、667和678nm波段实测的 标准离水反射率之间的经验关系如下：

[ρ_w]_N(555nm)≈0.9812[ρ_w]_N(53lnm)，R²＝0.9970

[ρ_w]_N(678nm)≈0.9812[ρ_w]_N(667nm)，R²＝0.9955

步骤(2)：在步骤(1)涉及的气象数据的基础上，利用Rayleigh散射 查找表算法，计算MODIS各波段Rayleigh散射贡献的反射率，并加以订正。

步骤(3)：利用Rayleigh散射修正后的2003年10月27日和28日两景MODIS 影像数据，构建531、551、667和678nm波段的反射率之间的经验关系。

2003年10月27日遥感影像的经验关系：

ρ_at(555nm)＝1.1681ρ_at(531nm)-0.0043，R²＝0.9795

ρ_at(678nm)＝0.9808ρ_at(667nm)，R²＝0.9934

2003年10月28日遥感影像的经验关系：

ρ_at(555nm)＝1.1106ρ_at(531nm)-0.0012，R²＝0.9936

ρ_at(678nm)＝0.9835ρ_at(667nm)，R²＝0.9880

步骤(4)：在步骤(3)的基础上，利用531、551、667和678nm波段反 射率之间的经验关系，结合Angstrom气溶胶指数模型，计算MODIS各波段的 气溶胶散射和气溶胶-Rayleigh多次散射贡献的反射率。

将步骤(3)所示的经验算法，即可求得R_as(531nm)和n的值。将上述 的R_as(531nm)和n的值代入Angstrom气溶胶指数模型，便可得到412、443、 488、531、551、667、678、745和870nm波段的气溶胶散射贡献的反射率， 如图1所示。

步骤(5)：在步骤(2)、(4)的基础上计算大气的程辐射贡献的反射率， 并对MODIS影像进行大气校正。

为了检验ACES大气校正算法的稳定性，本发明分别收集了2003年10 月27日和28日的MODIS影像及其同步标准离水反射率数据。利用同步的标 准离水反射率与经验算法的结果进行比对分析，并将该比对结果作为该经验 算法的误差。图2给出了ACES算法的误差状况。由图2可知，ACES算法在 412(11.96％)、443(13.50％)、488(14.94％)、531(3.56％)、555(4.19％)和 667nm(18.99％)波段具有较好的遥感反射率反演精度，但是在678、748和 870nm，该算法的误差较大(＞30.53％)。