全天空成像仪云量计算方法的改进
2014-08-13周文君牛生杰许潇锋
周文君,牛生杰,许潇锋
(1.南京信息工程大学中国气象局大气物理与大气环境重点开放实验室,江苏南京210044;2.江苏省盐城市气象局,江苏盐城224000)
0 引言
云是表征天气、气候特点的重要因素之一,也是大气动力、热力、水分输送过程综合作用的外在表现。云对地气系统辐射平衡影响很大,其辐射强迫对地球气候变化有着非常重要的作用(陈洪滨,1997;Long and Thomas,2000;Calbo and Pages,2002;Long,2004)。因此,对云的探测极为重要。目前,已经有一系列的探测手段对云的各个要素进行观测分析,例如:飞机穿云观测云中微物理特征,卫星遥感为云的分布及微物理特征量的研究提供了重要依据,以及地面对云的目测或使用地基测云仪测量云量、云高等(陈渭民等,1997;宋庆利,2003;Martins et al.,2003)。
随着科技的不断进步,地基测云仪的技术也得到了快速发展。美国大气辐射研究机构为了定量研究云量对辐射的影响,与国家海洋大气科研所以及气候能源工作实验室合作研发了半球天空成像仪(hemispheric sky imager,HSI)。该仪器可以在白天拍摄全天空160°范围内的图像。Long and Deluisi(1998)详细地介绍了 HSI仪器的观测原理和计算方法,并对HSI的云量观测资料进行了分析,结果表明:HSI仪器的观测结果比卫星和雷达的观测结果更精确。之后,在美国YANKEE公司的合作之下,HSI被发展为TSI-440(total sky imager 440,全天成像仪)。该仪器在图像处理技术及构建方面均取得了进一步发展,是目前使用较为广泛的地基测云仪。
TSI-440仪器可以自动进行全天空云量的持续性观测,时空分辨率较高,得到的云量计算结果较人工目测、卫星探测的结果更精确(Hodges,1998)。Long et al.(2001)使用独立的分析方法将TSI仪器的观测结果与WSI仪器(whole sky imager)观测结果进行对比分析,发现两者相关性较好。
此外,Sabburg and Long(2004)和 Sabburg and Parisi(2008)使用TSI-440仪器观测云量,并根据云量的变化,分析云对紫外辐射增强的作用。Calbo and Sabburg(2008)采用TSI以及全天空相机同时拍摄天空,对两张图片中的同一块云体采用统计分析、傅里叶变换等获取图像信息,对晴空、波状云、卷云、层状云和积状云这五种情况进行识别,准确率达到了76%。Kassianov et al.(2005)等使用 TSI和 HSI两台仪器,根据双基站测量云高的原理进行试验分析,得到云高计算公式,并且将分析结果与微脉冲激光雷达实测云高进行对比,结果表明:TSI仪器对于平均云底高度低于3.5 km的单层云体的云高测量是可行的。按照地面观测自动化的发展要求,国内也逐步开展地基测云仪的利用和研制,霍娟和吕达仁(2006)利用地基测云仪分析云量对辐射变化的影响,并以此反演气溶胶的变化特征(霍娟和吕仁达,2005;霍娟,2007)。傅德胜和王新芝(1995)、孙学金等(2009a,2009b)、孙晓刚等(2008)使用红外测云仪提取云图特征,进行云状分类。高太长等(2010)、秦健等(2011)总结了目前地基测云仪的研究现状,详细分析了全天空云图的获取、镜头保护、云点识别、云量计算等多项测云技术问题。
虽然,TSI测云仪已经得到了广泛应用,计算结果也较为精确,但是在不同的天空情况下,仍然存在一定的云量计算误差。Pfister and McKenzie(2003)通过分析TSI云量和地表辐射量变化,指出:如果太阳被薄云遮掩但依然清晰可见时,很可能出现辐射增强的现象;此外,在观测过程中,该仪器对太阳周边云量的探测以及太阳是否被云体遮掩的判断准确度还需要进一步提高。
综上所述,TSI的云量计算存在一定误差,还需进一步研究。本文主要利用太湖地区2008年5—10月的TSI-440观测资料展开云量计算误差的分析。首先分析不同能见度条件下天空图像的成像特征,然后分析阴天图像、复杂天空(多云)图像的红蓝比值分布情况及云量计算误差,并对减小误差的方法进行了初步研究。
1 TSI基本工作原理及资料来源
TSI仪器主要的部件有:照相机、带有加热装置的半球镜面、镜面上方的遮光带以及下方的电子设备系统(图1)。该仪器可以在白天自动进行全天空云量的持续监测,其工作原理是:通过仪器上方的照相机垂直向下拍摄带有加热装置的半球镜面,得到当时天空所呈现的图像,并将拍摄得到的图像自动存储到计算机上进行云量计算和处理。
大气、云粒子对可见光具有不同的散射原理:大气分子对可见光的散射与λ4成反比,如果天空为晴空,对蓝光波段的散射远远大于对红光波段的散射,因此晴空呈现蓝色;而云粒子对可见光的散射,在不同的波段散射的程度是相当的,因此云体呈现白色。根据不同的散射原理,可以从图像信息值(RGB值)中区分云与蓝天。自然界中的可见光都是按照一定的比例混合而得到,任意一种颜色都可以用红(R)、绿(G)、蓝(B)三色混合而成,图像中的每个像素点均含有RGB三个亮度值,而RGB值在一定程度上反映了红(中心波段为650 nm)、绿(中心波段为570 nm)、蓝(中心波段为450 nm)这三个波段的辐射强度。因此,首先提取图像中的信息值(RGB值),计算每个像素点的红蓝比值(定义为:R/B),然后设定合理的阈值来判别像元点是否为云点,这种方法被称为阈值识别法(Morris,2005)。
图1 全天空成像仪TSI-440(a)及其组成部分(b)Fig.1 (a)Total sky imager 440 and(b)its components
本文所用资料取自南京信息工程大学太湖观测站(120.31°E,31.58°N)。该站自 2008 年 4 月—2009年4月在太湖地区进行自动观测,具备多种气象观测仪器,如:全天空成像仪(TSI-440)、微脉冲激光雷达、微波辐射计、多通道旋转光度计等,为云量及辐射的研究提供了全面的数据资料。本文为了研究TSI云量计算误差,采用了2008年5—10月共计6个月的TSI资料以及无锡站的常规观测资料。TSI资料包括:原始天空图像、处理后图像以及相应的云量计算结果。图像分辨率为352×288,是24位全彩色图像,图像采集率为1 min/张,所拍摄的天空范围是可视张角160°,其工作时间平均为07:00—17:00(北京时间,下同),仪器的具体工作时间与太阳高度角有关,当太阳高度角大于3°时,仪器自动开始获取图片。
2 TSI资料分析
通常人工观测时将天空划分为十等分来估算云量大小,文中的云量均按照0~1计算。根据分析图像特征发现:当天空完全被掩蔽,但云体仍显明亮时,这种云较高,如:厚薄不均的蔽光高层云或透光高层云。反之,若云呈现灰白或灰黑色,天空阴沉,这种云则较低,如:蔽光层积云、雨层云等。文中所定义的薄云是高云(云高大于5 km)且云层较薄,例如毛卷云、钩卷云、卷层云等,在图像中此类云体与非云大气的雾霾粒子难以区分,是云量计算中的难点。而不透光云体是指云层较厚,可以完全遮挡太阳直射光的中低层云,成像时可能出现阴影区域。
2.1 不同能见度条件下的图像分析
一般情况下,当气溶胶、雾霾不严重时,人眼所见的晴空应为蓝色,同样相机成像也是如此。但是,如果能见度较低,天空将呈现灰白色。图2是4种不同能见度条件下的晴空图像及仪器自带的处理结果,可以看出:随着能见度的减小,原始天空图像的晴空区域逐渐由深蓝色变为蓝色进而变为灰白色,而仪器自带的处理结果中云量区域随着能见度的减小而增大。主要是由于气溶胶前向散射作用使得太阳周边像素点饱和,能见度越小,受太阳影响的区域越大。
能见度是影响大气成像的重要因素,图2可以说明大气成像与能见度密切相关。能见度不同直接影响成像效果及后期的分析处理工作,因此首先分析不同能见度条件下图像的成像特征。图3中给出了3种不同能见度条件下晴空图像红蓝比值的分布情况,可以看出:随着能见度的减小,红蓝比值有明显增大的趋势,同时,比值为1(R=B=255的饱和像素点)的像元数随着能见度的减小而增加。基于能见度对天空成像的影响,本文中所选择的天空图像均为能见度大于5 km的图像。
同时,从晴空图像中也可以看出:在晴空天气情况下,由于气溶胶前向散射作用,会产生一定的云量计算误差。针对晴空图像所产生的云量误差,Long(2010)提出了云量改进方案。该方案是根据天空其他区域的云量分布情况及变化特征来估算太阳圈及近地平区域的云量。经与目测结果对比,发现该方案可以有效地减小晴空天气情况下全天空成像仪的云量计算误差。
图225 km(a、b)、15 km(c、d)、8 km(e、f)、3 km(g、h)能见度条件下的晴空图像(a、e、c、g)及仪器自带(b、f、d、h)的处理结果(蓝色、白色和灰色区域分别代表晴空、不透光云和薄云;黑色条带表示遮光带,遮光带上圆点表示太阳所在位置,其中黄点代表太阳未被云体掩蔽,白点表示太阳被云体掩蔽)Fig.2 (a,e,c,g)The clear sky images and(b,f,d,h)the corresponding TSI processed images at different visibilities(In the processed images,the blue regions represent cloudless and the write and gray regions stand for opaque cloud and thin cloud respectively.The black strip is shading belt.The dot on the shading belt denotes the position of the sun.The yellow dot represents sunny days without cloud and the white dot denotes cloudy days) a,b.25 km;c,d.15 km;e,f.8 km;g,h.3 km
图3 不同能见度条件下晴空图像的红蓝比值分布情况Fig.3 The red blue ratio distributions of the clear sky images at different visibilities
2.2 阴天图像红蓝比值分析及云量计算
不同天空状况,云的成像特征并不相同,例如:阴天图像比晴空图像明显偏暗。这主要是因为阴天里太阳直射光被云完全遮挡,图像的亮度整体偏小。TSI仪器在进行图像处理时,对所有的图像均采用统一阈值,一般情况下,红蓝比值小于0.7的为晴空点,比值大于0.8的为云点,介于0.7~0.8之间的为薄云点,此阈值对晴空图像或亮度较高的图像比较适用,而对于亮度偏暗的图像会造成一定的云量计算误差。图4中是2008年10月4日13:10的阴天图像及仪器自带的处理结果,地面观测资料显示当天全天云量均为10成,即云量为1。而仪器的结果为不透光云量0.332,薄云量0.441,晴空区域0.227,这与地面观测资料及目测识别图像的结果有较大差异。
文中选择了500幅阴天图像做RGB值直方图统计分析,RGB值取值范围在0~255之间,将其分为51档(以5为间隔单位),统计每幅图像天空范围像素点的RGB值分布情况。结果发现:图像中的红色像元值(R)多集中于125~240之间,绿色像元值(G)多集中于150~245之间,蓝色像元值(B)多集中于185~255之间。而晴空图像的R值多集中于80~200之间,G值多集中于150~220之间,B值多集中于230~255之间,太阳周边的像素点均为R=B=G=255的饱和像素点。由此可见:阴天图像的R值比晴空图像的R值偏大,而蓝色像元值B值相对较小,且阴天图像中的饱和像素点较少,说明太阳完全被云体遮掩。
图5是图4阴天图像对应的红蓝比值分布情况,可以看出:此时像元点红蓝比值多集中于0.65~0.9之间,且在红蓝比值为1处并无像元点,即无饱和像素点。说明此时太阳已被完全遮掩。容易看出:若对该图采用仪器的统一阈值(0.7;0.8),将会造成较大的云量计算误差。
图4 2008年10月4日13:10的阴天图像(a)、仪器自带程序的处理后图像(b)及改进阈值后的图像(c)Fig.4 (a)The overcast sky image,(b)the corresponding TSI processed image and(c)the image after setting the new threshold at 13:10 BST on October 4,2008
图5 2008年10月4日13:10的阴天图像的红蓝比值的分布Fig.5 The red blue ratio distribution of overcast sky images at 13:10 BST on October 4,2008
选择200幅阴天图像进行红蓝比值概率分布的分析,采用直方图统计分析法,将天空像素点的红蓝比值范围0.3~1.0分为35档(以0.02为间隔单位),首先对每幅阴天图像计算每个单位档上的像素点个数及概率分布值,然后将这200组数据每个单位档做平均值计算,得到像素点的平均概率分布。同时,分别计算这200组数据与平均概率分布之间的标准偏差,结果发现:所有数据的标准偏差均在0.03~0.06之间,标准偏差较小。说明像素点的平均概率分布可以反映阴天图像红蓝比值的分布情况。此外,为了进一步检验平均概率分布的代表性,另选取择了200幅阴天图像,分析红蓝比值的分布情况,并计算每幅图像的红蓝比值与平均概率分布的相关系数,共计200个相关系数。其中,22.7%的相关系数大于0.8,45.5%的相关系数介于0.7~0.8之间,还有27.2%的相关系数介于0.6~0.7之间,即相关系数大于0.6的数据占了整个样本数量的95%。其余一部分数据与平均概率分布的相关性不是非常好,通过分析这部分数据的红蓝比值分布发现:比值多数集中于0.8~1.0之间,观察图像资料发现这部分资料均为天空辐亮度较大的图像资料。
图6 阴天图像红蓝比值的平均概率分布Fig.6 The mean probability distribution of red blue ratio on cloudy days
图6给出了阴天图像红蓝比值的平均概率分布,可以看出,阴天图像比值集中于0.62~1.0之间,峰值位于0.7~0.8之间。因此,可以判定红蓝比值大于0.62的像素点均为云点,比值介于0.3~0.62的像素点均为晴空点,不透光云点阈值设定为0.66(即比值介于0.66~1.0像素点为云点),而0.62~0.66之间的像素点可能是薄云点,也可能是图像中的阴影区域,暂时无法定论,需要根据不同的天空图像进行分析。根据上述设定的阈值,重新计算后得到了图4a中阴天图像对应的计算结果及改进阈值之后的图像(图4c),其计算结果为:不透光云量1,总云量1,与实际观测结果一致。
由于阴天时辐亮度较小,使得图像中阴影区域被误识别为薄云或晴空,表1详细地给出了多个阴天图像的原始结果及改进阈值之后的结果,可以看出:原始结果高估了薄云量,误将不透光云量以薄云量计算,且总云量也存在一定的计算误差;改进阈值后得到的结果,薄云量的计算误差减小,不透光云量及总云量更接近实际观测资料。可见,针对阴天图像重新设定阈值可以有效地减小因辐亮度不同而产生的薄云计算误差,提高了云量计算的精度。
表1 2008年10月阴天图像的原始云量结果及改进阈值后的云量结果Table 1 The raw result and new result after setting the new threshold of cloud cover in the overcast sky image in October 2008
2.3 复杂天空图像分析
本文所指的复杂天空图像是云量大于0.2的多云图像,在图像分析过程中发现,有较多的图像成像时出现阴影区域,这些区域在图像处理过程中被作为薄云计算,造成了云量计算误差,此类情况与阴天时出现误差的情况类似。图7a是一幅蓝天白云分界清晰的天空图像,从仪器自带的处理结果中(图7b)可以看出:天空中云体较厚的部分(图中灰色区域)在图像处理过程中被作为薄云计算;而太阳未被云体遮掩,因此周边像素点饱和的区域,被识别为不透光云。图7c显示了该图像的红蓝比值分布情况,可以看出:图中晴空区域比值最小(蓝色区域);云层较高、较明亮的云体比值偏大,如太阳周边云体(红色区域);而云体较厚、较暗的部分红蓝比值相对较小,如图中的黄色区域(对应图7b的灰色区域)。由于此类图像的云量计算误差与阴天情况类似,因此重新设定的阈值与阴天图像相同,比值为0.3~0.62的像素点为晴空点,0.66~1.0的像素点为不透光云点,而介于0.62~0.66之间的像素点为薄云点。
图7 原始天空图像(a)、仪器自带的处理图像(b)及原始天空图像与仪器自带的处理图像比值(c)Fig.7 (a)The raw image,(b)the corresponding TSI processed image and(c)the red blue ratio distribution image
图8是2008年10月16日14:07的天空图像、仪器自带的处理图像及改进阈值后得到的图像。从原始图像中,可以清楚地分辨出红色框内区域为不透光云体,但由于云层较低,云体较厚,呈现灰色,在原始处理过程中被作为薄云计算,原始计算结果为:薄云量0.267,不透光云量0.239,总云量为0.506。通过目测识别可知红色框内区域为不透光云体,被误识别为薄云体。利用新阈值(晴空0.62,云0.66)对图像处理后,得的计算结果为:薄云量0.042,不透光云量为0.602。地面观测资料显示:总云量为6成,即0.6,低云量为0.6,更改阈值之后的计算结果与实际观测资料更为接近。由此可见,采用新阈值减小了仪器薄云量的计算误差,同时也提高了总云量的计算精度。
图82008年10月16日14:07的原始天空图像(a)、处理后图像(b)和改进阈值后的图像(c)Fig.8 (a)The raw image,(b)the corresponding processed image and(c)the simulated image after setting the new threshold at 14:07 BST on October 16,2008
结合太湖观测站的TSI资料及台站观测资料,选出每日08:00、11:00及14:00的观测资料共计450组数据,将不同天气下的云量进行误差统计分析。由于本文所采用的新阈值是对多云及阴天情况,因此对云量为0.3~1.0的8种情况进行误差计算。图9给出了不同情况下仪器自带处理结果、改进阈值后的结果分别与目测云量的绝对误差情况,可以看出,仪器自带的处理结果在云量大于等于0.7时,绝对误差大于0.1,其中阴天时误差达到最大。改进阈值之后得到的计算结果在多云及阴天时,绝对误差均小于0.08,说明重新设定阈值可以有效地减小多云及阴天情况下的云量计算误差。
3 结论与讨论
利用全天空成像仪代替人工观测云量,在多数情况下效果良好。但是,由于天空成像特征与天气变化情况、能见度等密切相关,如随着能见度的减小,晴空图像的红蓝比值逐渐增大,所以导致全天空成像仪自带程序在复杂天空状况下云量计算出现较大误差。
图9 不同天气情况下的云量绝对误差Fig.9 The absolute differences under different weather conditions
在分析太湖观测站TSI资料的过程中发现:该仪器在晴空、阴天及云体较厚的情况下易产生云量计算误差。针对阴天情况,文中对多幅天空图像进行红蓝比值统计分析,结果发现:比值主要集中于0.62~1.0之间。根据图像的成像特征,将晴空点阈值设定为0.62(<0.62判为晴空),云点阈值设定为0.66(>0.66判为云),采用新阈值得到的云量结果与常规气象观测结果更为接近。同时,对天空图像中存在阴影区域的其他图像也使用新阈值,得到的云量结果比原结果更准确,减小了因辐亮度不同而造成的云量计算误差。
在图像分析过程中,还发现图像的红蓝比值分布情况与太阳高度角、太阳方位角也密切相关。因此,要解决以上这些问题,应根据不同的天气情况、像素点的位置、太阳中心点位置等信息设定不同的阈值,建立一套完整的云点阈值识别法。另外,也可以尝试运用其他计算途径来获取云量信息,例如二维傅里叶变换、统计分析法或图像处理技术中的边缘提取法等进行云量的计算,进一步改善计算结果。
陈洪滨.1997.关于云和云天大气对太阳辐射的吸收异常[J].大气科学,21(6):750-757.
陈渭民,边多,郁凡.1997.由卫星资料估算晴空太阳直接辐射和散射辐射[J].南京气象学院学报,20(3):326-333.
傅德胜,王新芝.1995.云图纹理特征的抽取与云的自动分类[J].南京气象学院学报,18(4):530-535.
高太长,刘磊,赵世军,等.2010.全天空测云技术现状及发展[J].应用气象学报,21(1):100-109.
霍娟.2007.地基全天空成像系统云与气溶反演及其应用研究[D].北京:中国科学研究院.
霍娟,吕仁达.2005.利用全天空数字图像对北京上空云况分布特征的试验分析[J].气象科学,25(3):238-243.
霍娟,吕达仁.2006.晴空与有云大气辐射分布的数值模拟及其对全天空图像云识别的应用[J].气象学报,64(1):31-38.
秦健,刘磊,高太长,等.2011.基于统计方法的地基器测云天识别[J].气象科学,31(5):604-612.
宋庆利.2003.利用卫星资料研究云对地面净辐射的影响[D].南京:南京气象学院.
孙晓刚,孙学金,牛珍聪,等.2008.全天空云图获取的一种方式及算法实现[J].气象科学,28(3):338-341.
孙学金,刘磊,高太长,等.2009a.基于模糊纹理光谱的全天空红外图像云分类[J].应用气象学报,20(2):157-163.
孙学金,刘磊,高太长,等.2009b.基于LBP算法的全天空红外云图的分类研究[J].大气科学学报,32(4):490-497.
Calbo J,Pages D.2002.Empirical studies of cloud effects on UV radiation:A Review[J].J Geophys Res,43(2):1-28.
Calbo J,Sabburg J M.2008.Feature extraction from whole sky groundbased images for cloud-type recognition [J].Amer Meteor Soc,25:3-14.
Hodges G.1998.A Comparison of fractional cloud cover determination by a sky image with visual estimation by trained national weather service observers[EB/OL].[2012-10-10].http://www.noaa.gov/.
Kassianov E,Long C N,Christy J.2005.Cloud-base-height estimation from paired ground-based hemispherical observations[J].J Appl Meteor,44(8):1221-1233.
Long C N.2004.The shortwave clear-sky detection and fitting algorithm:Algorithm operational details and explanations[R]//Atmospheric radiation measurement program.Tech Rep ARM TR-004.
Long C N.2010.Correcting for circumsolar and near-horizon errors in sky cover retrievals from sky images[J].J Atmos Soc,4(1):45-52.
Long C N,Deluisi J J.1998.Development of an automated hemispheric sky imager for cloud fraction retrievals[C]//10th symposium on meteorological observations and instrumentation.Arizona:Amer Meteor Soc:171-174.
Long C N,Thomas P.2000.Identification of clear skies from broadband dynamometer measurements and calculation of down-welling shortwave cloud effects[J].J Geophys Res,105(12):15609-15626.
Long C N,Slater D W,Tooman T.2001.Total sky imager(TSI)model 880 status and testing results[R]//Atmospheric radiation measurement program.Tech Rep ARM TR-004.
Martins F R,Souza M P,Pereira E B.2003.Comparative study of satellite and ground techniques for cloud cover determination[J].Space Res,32(11):2275-2280.
Morris V R.2005.Total sky imager hand book[EB/OL].[2012-10-10].http://www.arm.gov/.
Pfister G,McKenzie R L.2003.Cloud coverage based on all-sky imaging and its impact on surface solar irradiance[J].Amer Meteor Soc,42:1421-1434.
Sabburg J M,Long C N.2004.Improved sky imaging for studies of enhanced UV irradiance[J].Atmos Chem Phys,4:2543-2552.
Sabburg J M,Parisi A V.2008.Cloud observation for the statistical evaluation of the UV index at Toowoomba,Australia[J].Int J Biometeorol,52:159-166.