金 健,卢 远,华 璀,何 文
(广西师范学院 遥感与GIS开放实验室,广西 南宁 530001)
随着珠江-西江经济带的发展上升为国家战略,粤桂交界地区的经济社会发展、工业化过程、生物质燃烧以及城市生活等都引起各界广泛关注.而在迅速发展起来的同时,大气污染也逐渐成为敏感话题.其中PM2.5对城市可见度、人类身体健康、城市环境等的危害尤其显著,所以对其定量遥感研究已引起决策者们的重点关注.利用遥感数据对气溶胶特性的反演研究已经开展30多年[1].在美国、欧洲等发达地区,利用观测数据反演和研究地面颗粒物浓度时空分布特征的工作已经大范围开展.而目前国内相关方面研究主要集中在北京、长江三角洲等发达地区[2-4].
本文以2005、2007、2009、2011年为研究时段,研究范围东经:107°43′~ 115°24 ′,北纬23°35 ′~21°20′N.遥感影像采用国际科学数据提供的MOD021KM Terra卫星1B 1KM分辨率的对地观测影像.并利用ARONET气溶胶观测网站上公布的Bac_Giang(北江)站、MuKdahan(穆达汉)站和Hong_kong_poly(香港)站L2.0级550nm陆地上气溶胶光学厚度监测值对反演结果进行精度验证[5].站点位置及研究范围见图1.
图1 ARONET观测站点位置分布与研究范围图Fig.1 The location of the observation of ARONET and the scope of the study
气溶胶光学厚度本身是一个无量纲的正值,一般从0到1计数,值越大,表明气溶胶光学厚度越大.本次研究利用暗像元法对研究区的气溶胶光学厚度进行反演.其基本原理是利用低反射率的地表在近红外(2.13~3.8μm)与红(0.47μm)、蓝(0.65μm)波段反射率具有很好的线性相关性(ρ0.47=ρ*2.1/4;ρ0.65=ρ*2.1/2;ρ:为地面反射率,ρ*为传感器在大气外界观测到的表观反射率)的性质来推算红、蓝波段的表观反射率[6-7].
在假定气溶胶类型、大气模式的基础上,根据相应方位角、天顶角等关系,利用IDL语言调用6S模型计算得出适应的查找表.最后通过输入地表反射率、相关的几何参数、波段表观反射率来查找气溶胶的光学厚度.
首先假设地表是均匀朗伯表面,大气沿垂直方向均匀变化,表观反射率Q(μv)表达为:
式(1)中μs为太阳天顶角的余弦,μv为观测天顶角的余弦;Qo(μv)为观测方向的路径辐射项;r为朗伯体地表反射率,S为大气下界的半球反射率;T为大气透过率;μsFo为大气层顶与太阳光垂直方向的辐射通量密度.
利用入射太阳辐射项Foμs对式(1)进行归一化可得大气顶部反射率ρTOA的表达式为:
其中,ρo为大气路径辐射项等效反射率;ρs为地表二向反射率,当地表为朗伯体时为r;Φ为相对方位角.从(2)中可以明显看出,ρs、S、ρo和T(μs)·T(μv)等参数是未知量.S、ρo和T(μs)·T(μv)参数本身是气溶胶函数,ρs的估算直接影响地表反射噪声的去除[8].为了更准确的反演气溶胶光学厚度,需要合理设定气溶胶模型,以提供单次散射反照率值和气溶胶相函数.在确定了气溶胶模式和地表反射率后,根据上述公式,可以从表观反射率反演得到气溶胶光学厚度.
步骤1 数据预处理;将影像进行辐射定标、双眼皮去除、几何校正、去云处理等.由于广西气候潮湿,常常为多云天气,获得的卫星影像受云的干扰尤为严重,这将极大影响对影像数据的使用.本次的研究采用改进后的云检测算法,该算法根据云检测原理中不同波段对云的敏感特性,以及MODIS数据各通道的波普特性,结合大气窗口和云的辐射传输特点[9],分别选择了B1、B6、B8、B26、B29、B31波段,从可见光到热红外6个通道数据进行云检测[10].
步骤2 暗像目标的选取;根据Kaufman等专家学者的研究,当波长2.1μm处的像元表观反射率ρ*2.1满足大于0.01小于0.35时,认为可以用作暗像目标,但对于浓密植被,0.47、0.66μm波长处的反射率之间的关系与植被的茂密程度有关[10].散射角由下述公式获得:
(3)式中θs为太阳天顶角;θv为观测天顶角;Φ为相对方位角.
用受气溶胶影响较小的短波红外来计算归一化植被指数:
式(4)中ρm1.24、ρm2.1为1.24μm、2.1μm波段处的表观反射率.
步骤3 计算地表反射率;通过ρ2.1和ρ0.47、ρ 0.66与植被指数、散射角与的关系求得0.47μm、0.66μm处的地表反射率[9]:
(5)、(6)式中a、b根据散射角和NDVISWIR 确定.
步骤4 选择气溶胶类型;由于桂东南地处与广东邻近的内陆地区,南面临北部湾,受地理环境以及自然条件的影响,所以分别选择城市型、大陆型、海洋型进行反演,再与ARONET网站公布的数据进行比较分析,选取最适合的气溶胶类型.
步骤5 生成相应查找表;首先确定几何参数(太阳天顶角、太阳方位角、卫星天顶角、卫星方位角)、大气模式、气溶胶类型参数、目标高度参数、传感器高度参数、光谱参数、地表反射率类型等.根据这七类参数建立适合的查找表.本次研究通过IDL语言调用6S辐射模型,经过循环输入参数,调整各类参数计算的步长,从而可以相应提高最终反演精度[8].本次研究的参数设置为:
气溶胶光学厚度范围(taer55=0.1~1.2,step=0.2);
太阳天顶角(asol=0~60,step=15);
卫星天顶角(avis=0~60,step=15);
卫星和太阳相对方位角(phiv=0~180,step=30);
大气模式(iastm=2,为中纬度夏季、iastm=3,为中纬度冬季大气);
光谱参数(iwave=45为绿光波段);
地表反射类型(inhomo=0,为均匀表面)
根据上述所设置的各类参数自动生成相应的查找表,在气溶胶光学厚度计算中输入经过预处理的地表反射率后,自动进行逐像元查找,最终生成反演结果.
将2009年的反演结果与AERONET气溶胶监测网站Bac_Giang站、MuKdahan站和Hongkong_poly站2009年公布的观测数据进行精度检验.选取以地面站点为中心的反演结果(0.55μm处光学厚度各月份平均值)作为待验证值,用卫星过境前后半小时地基观测数据的各月份平均值作为标准值,进行对比[11-13],验证结果见表1.
表1 2009年气溶胶光学厚度反演值与ARONET观测数据一元线性拟合统计结果Tab.1 Statistic results of the linear relation of aerosol optical thickness and ARONET data in 2009
可看出大陆型和海洋型气溶胶模式的反演结果与AERONET网站公布的观测结果有一定的相似性.主要是因为大陆型气溶胶70%由沙尘性粒子组成,城市型气溶胶主要有水溶性、煤烟性粒子组成,而海洋型气溶胶95%由海洋性粒子组成,说明大陆型气溶胶还是贴近乡村型大气模式.
反演结果中城市型气溶胶模式得出的结果与观测值整体偏差较大,这是由于桂东南地区的城市化水平还较低,城市型气溶胶中多含水溶性以及煤烟性粒子,这与桂东南实际地表气溶胶模式相差甚远.海洋型气溶胶中海洋性粒子占95%,水溶性粒子占5%,这与桂东南地区的自然环境也不相符合,而大陆型反演结果的值更接近观测值,说明桂东南地区还是以沙尘型粒子为主.
故本次研究为了更贴合实际情况,选用大陆型气溶胶模式进行反演.
通过上述反演步骤对2005、2007、2009、2011年桂东南地区气溶胶光学厚度进行反演,并将反演结果进行分析对比.选择每年各月份中研究区云量少于10%的影像生成年合成数据,气溶胶类型选用大陆型,最终反演结果见图2.
对比上述反演结果可看出,广西东南部的气溶胶光学厚度高值区在2007年有南移的现象,并在2009年和2011年呈现向周边扩散的趋势.从反演结果可以看出,由2005年至2011年梧州、肇庆两市的气溶胶光学厚度值有明显上升趋势,玉林、云浮以及湛江也逐渐成为较高值区.这是由于广西高层对桂东南地区的重视,桂东南地区掀起新一轮的经济发展大潮.受北部湾发展浪潮以及粤桂合作特别实验区的梧州-肇庆交界地区发展的影响,桂东南逐渐形成为珠三角地区资金、技术、产业直接向西部转移、打开大西南和东盟市场的重要通道节点,使得这一区域内的经济飞速发展,并推动了周边各区域经济的迅速增长,促进城市化进程,加快了工业化.桂东南地区也将会是广东省产业转移最好的目的地,一个桂东南地区经济的发展,也将会带动并影响整个广西经济的发展,但同时也加大各种染料的消耗、废气排放等.
本文使用经典暗像元法对桂东南地区及周边区域进行气溶胶光学厚度的反演,经过初步精度验证,其精度能满足一般研究等时机应用,对研究桂东南地区的气溶胶时空分布有很大的意义.但遥感气溶胶光学厚度定量反演表征的是在无云大气铅直方向上的积分,与实地地面定点观测数据必定存在一定的误差[14].基于暗像元法的原理分析,由于复杂地表为地表反射率推算的盲区,所以对城区反演效果较差.下一步研究可加入对下垫面类型的模拟,补充使用耦合算法,使得反演结果更贴近实际情况.
图2 2005年(a)、2007年(b)、2009年(c)、2011年(d)桂东南地区气溶胶光学厚度反演结果Fig.2 The inversion results of aerosol optical thickness in the southeast area of Guangxi in 2005年(a)、2007年(b)、2009年(c)、2011年(d)
[1]陈良富,李莘莘,陶金花,等.气溶胶遥感定量反演研究与应用[M].北京:科学出版社,2011:14-17.
[2]王静,杨复沫,王鼎益,等.北京市MODIS气溶胶光学厚度和PM2.5质量浓度的特征及其相关性[J].中国科学院研究生院学报,2010,27(1):10-16.
[3]宋巍巍,管东生.利用TM影像反演广州市气溶胶光学厚度空间分布[J].环境科学学报,2008,28(8):1638-1645.
[4]肖钟湧,江洪,陈健,等.利用MODIS遥感数据反演广州市气溶胶光学厚度[J].中国环境科学,2010,30(5):577-584.
[5]李晓静,张鹏,张里阳.MODIS气溶胶光学厚度反演算法技术报告[C].4版.北京:国家卫星气象中心出版社,2006.
[6]Kaufman Y J.Algorithm for antomatic atmosphieric corrections to visible and near-IR satelliteimagery Int J Rem[J].Sens,1988,9:1357-1381.
[7]Kanfman Y J,Wald A E,Remer L A.The MODIS 2.1μnChannel-Correlation with visible reflectance for use in remote sensing of aerosol[J].IEEE Trans Geosci Remote Sens,1997,35:1286-1298.
[8]彭妮娜,易维宁,方勇华.400~1000nm波段反演气溶胶光学厚度的暗像元法[J].红外与激光工程,2008,5(37):878-883.
[9]李晓静,张鹏,张兴赢,等.中国区域MODIS陆上气溶胶光学厚度产品检验[J].应用气象学报,2009,2(20):148-149.
[10]何全军,曹静,黄江,等.基于多光谱综合的MODIS数据云检测研究[J].国土资源遥感,2006,69(3):19-22.
[11]任佳,王振会,孙林,等.基于MODIS数据反演江浙皖地区气溶胶光学厚度[J].环境科学与技术,2010,8(33):168-171.
[12]郑卓云,陈良富,郑君瑜,等.高分辨率气溶胶光学厚度在珠三角及香港地区区域颗粒物监测中的应用研究[J].环境科学学报,2011,31(6):1154-1161.
[13]李冬辉.中国东部海域气溶胶光学特性及其季节变化和地理分布特征研究[D].青岛:信息科学与工程学院,2012.
[14]曹琦,杨金梅,淑杰.气溶胶光学厚度反演查找表法误差初步解析[J].安徽农业科学,2010,38(12):6356-6358.