2. 福建省龙岩市气象局,龙岩 364000;
3. 北京市气象局,北京 100089
2. Longyan Meteorological Office of Fujian Province, Longyan 364000;
3. Beijing Municipal Meteorological Service, Beijing 100089
云在天气系统发展、降水形成和大气辐射传输等物理过程中扮演着重要角色,准确、及时地获取云的信息 (云底高度、云量、云状),对气象、民航等诸多领域都有着十分重要的意义。然而,目前气象台站云的观测仍以人工目测为主,人为主观因素大,如何准确有效地进行云自动化观测是云探测领域的一项重要任务[1]。
目前常用的自动云底高度测量方法[2]有激光测云仪和红外辐射测云仪。激光云高仪采用主动测量方式,具有测量精度高等优点,但受气溶胶和云类型的影响,有一定比例漏测。红外辐射测云仪接收云底辐射亮温,根据辐射传输方程反演云底高度,能够昼夜观测,数据获取率高。但云底辐射亮温受到气溶胶、水汽和云体密实程度等因素影响,造成算法复杂,测量误差较大。如何提高云底高度测量能力是气象探测领域面临的重要课题。
随着数字摄像技术和立体视觉传感器的发展,尤其是双目成像视觉传感器以结构简单、使用方便、测量精度高等诸多优点被广泛应用。双目成像测距系统正是基于立体视角差测距原理,利用两个相对固定的摄像机,从不同位置采集到两幅图像,通过匹配特征点,得到相对视差,然后计算出目标物的三维信息[3]。
国内外采用立体成像技术测量云底高度已经有了一定的研究基础,文献[4]对人工观测、无线电探空仪、云幂仪、激光雷达、云雷达和立体成像仪等设备做了云底高度测量对比,得出基于立体成像方式测量云底高度具有精度最高、成本最小等优点。
1996年美国Sandia国家实验室Mark等[5]利用一对长基线的鱼眼相机在不同位置拍摄云体,利用光流匹配等算法获取云图同名点,计算云底高度,为满足大视角探测要求,该方法增加基线长度,造成基线过长很难应用于气象台站观测。谭涌波等[6]用一对基线长73 m的相机拍摄云体,利用模块匹配等算法寻找同名点计算云底高度,该方法很好地解决了短基线摄影测量问题,但由于天空背景复杂,通过亮度差来分离云和天空背景,匹配区域很难确定,所以数据获取率不高,采用高精度光学仪器对两台相机进行相对姿态角标定,标定过程复杂。
针对气象台站云底高度自动化观测的特殊性,一方面数据获取率要求高,另一方面要求探测精度高,安装区域小。本研究旨在解决自动获取同名点、数据质量控制和标校等关键技术问题,构成了采用较短基线 (60 m) 的双目成像系统,并与激光云底高度仪进行对比分析。
1 系统结构与框图双目成像云底高度测量系统采用一对相机置于可旋转的云台中,两个云台相距60 m,通过对云台的控制实现全天空扫描,获取全天空云图,并计算天顶方向50°范围内的云底高度。系统硬件结构主要包括图像采集传感器 (一对相机) 及图像处理和终端显示,系统结构如图 1所示。
|
|
| 图 1. 系统结构示意图 Fig 1. System structure diagram | |
图像采集传感器采用Sony公司的FCB-1010P机芯,视场角为54°,图像分辨率为44万像素,最低照度为0.1 lx。旋转云台采用一体化高速球机,球机旋转角范围,水平角为0°~360°,俯仰角为0°~93°,球机预置位重复精度小于0.05°。工作温度为-40 ℃~60 ℃。
系统功能包括:图像采集、处理、标校、全天空云图拼接、云底高度分层与计算以及数据存储和终端显示等 (如图 2所示)。先对相机进行内外方位元素标校,包括实验室标校和现场标校。图像采集和处理主要是完成同步采集两幅图像, 并进行图像处理自动获取同名点,计算云底高度。全天空云图拼接是利用旋转云台,全方位多仰角获取全天空云图,利用球面投影原理拼接图像,计算全天空云量。云底高度分层是计算同一幅图像中云底高度分布直方图,利用大津阈值分割方法等得到不同云层的高度。
|
|
| 图 2. 系统功能框图 Fig 2. The block diagram of system function | |
2 双目成像云底高度测量 2.1 双目测距原理
当透镜物距U远大于焦距f时,可以将相机成像模型等效为小孔成像模型[7-8]。在实际应用中,由于两个像平面很难保持在同一个水平面上,可以通过求取两像平面的相对姿态角 (横滚 (α)、俯仰 (β)、旋转 (γ)) 进行校正,校正之后基于小孔成像模型,双目成像测距原理示意图如图 3所示。
|
|
| 图 3. 双目成像测距原理示意图 Fig 3. The principle diagram of dual-camera stereovision distance measurement | |
设两相机之间的距离 (基线长) 为|O1O2|=d,A为被测目标物,Z为目标物相对像平面的垂直距离,P为垂足,F1,F2为焦点,f为焦距 (单位像素点数);A1,A2为A点在两像平面上的投影点坐标,分别为 (x1, y1),(x2, y2); F1F′1, F2F′2, AB垂直于A1A2。
根据相似三角形原理可知:|AB|=|O1O2|·


所以双目成像测距公式可表示为
|
(1) |
由于CCD获取的云图像具有纹理贫乏、信噪比低等特征,对其进行图像匹配之前,必须进行图像增强和滤波等预处理,以满足图像特征提取和匹配的要求。现有的图像增强算法主要有直方图均衡化、对比度限制自适应直方图均衡化、Wallis滤波等,本文采用直方图均衡化的方法进行图像增强。
直方图均衡化就是将一幅图像的灰度级归一化分布在0≤r≤1范围内,以累积分布函数作为变换函数,产生一幅具有均匀概率密度的灰度级分布图像。扩展像素取值的动态范围,提高匹配正确率。
累积分布函数为
|
(2) |
式 (2) 中, Pr (r) 为当图像灰度值为r时的概率密度。
2.3 特征提取与匹配基于地基拍摄的云图,受投影角度、光照效应等影响较大,很难直接进行区域匹配。角点是图像中具有高曲率的局部特征点,对光照和几何畸变有很好的鲁棒性 (表征该系统对某特性或参数摄动的不敏感性),被广泛用于特征匹配。Harris角点检测算子相比于其他算子具有最高的稳定性[9-10],本文采用Harris角点检测器对图像进行角点检测,然后利用极线约束消除两图中不一致的角点,再从一个图中每个角点所在的区域,按一定半径在另一幅图像中寻找对应角点所在的区域,并计算两区域的相关系数,将相关系数最大的角点作为特征匹配点。利用角点作为特征匹配点,进行相关区域匹配可以减少计算量,提高匹配正确率,其不足之处在于检测精度只能达到像素级,为了提高测量精度,本文在Harris角点检测基础上再进行亚像素级角点检测。
2.3.1 亚像素角点检测现有的亚像素级角点检测算法很多,其中利用梯度向量的方法简单直接[11-13],该方法是基于对向量正交性的观测而实现的,在以角点为中心点的某个领域内,假设亚像素角点为q,则到其邻域点p的向量与p点处的图像梯度正交,由于噪声等原因造成正交值不为零存在误差为ε,可以通过遍历所有正交值求取累积和最小的值此时的q点位置则为亚像素级角点位置。亚像素角点检测示意图如图 4所示,带虚线箭头方向为梯度方向。
|
|
| 图 4. 亚像素点检测示意图 Fig 4. The schematic diagram of subpixel detection | |
计算表达式[14]为εi=▽HiT·(q-pi),其中,▽HiT表示在q点的一个邻域点pi处的图像梯度,q的值可以通过求解Harris角点的某个区域内所有pi点的联立方程
|
(3) |
使εi最小化得到。然后通过多次迭代,找到低于某个阈值点的q值作为新的亚像素级角点位置。
2.3.2 极线约束极线约束[15]是唯一不依赖于场景的约束方法被广泛用到图像匹配中,由对极几何关系可知,第1幅图像上的每一个点U,在第2幅图像上的匹配点U′一定在其对应的对极线L′上,但是由于噪声的存在一般情况下并不一定落在对极线上而是在对极线的附近通过极线约束匹配误差,可以提高匹配精度,获取最终匹配点。
2.3.3 归一化互相关图像匹配互相关图像匹配方法具有很高的准确性和自适应性,但是该方法找到的是匹配区域以真实位置为中心的平缓峰值的位置,往往找不到准确的尖峰位置,为了克服这一缺点需要对图像做增强或边缘处理,实际上当两幅图像高度相关时,其相关性集中在图像的边缘部分。
归一化互相关计算公式[6]为
|
(4) |
式 (4) 中,R (i, j) 为两特征角点对应区域的归一化相关系数,T为角点对应M×N个像点的模板,Si,j为待搜索区域上的角点对应的区域。
该方法最大的缺点是旋转和缩放对计算相关系数影响很大,所以要求两相机相对姿态具有很好的一致性,尽量减少旋转和缩放的影响。
通过旋转变换、角点检测,图像匹配之后的同名点图如图 5所示,图中十字箭头表示同名点,与箭头的连接线表示两幅图中同名点的矢量位移。
|
|
| 图 5. 2011年5月10日10:32 (北京时,下同) 同一时刻两幅云图同名点 (同名点数为1713) (a) 经过图像匹配之后的同名点图,(b) 经过旋转变换之后的同名点图 Fig 5. The drawing of the closest points (the number: 1713) in two images at 1032 BT 10 May in 2011 (a) the closest points map through image matching, (b) the closest points map through rotation transformation | |
从图 5可以看出, 同名点集中在云的边缘部分,角点和同名点的获取概率主要受云图纹理清晰程度以及云图中包含云区域大小的影响较大,当云图纹理清晰时同名点获取概率较高,相反获取率较低,当云底平整无纹理时无法找到同名点。
3 相机现场标校技术 3.1 基于星星相对位置的相机姿态角计算由于在实际安装过程中很难保证两相机具有一致的姿态角,假设相对姿态角 (横滚 (α)、俯仰 (β)、旋转 (γ)),根据立体像对相对定向原理可以求出两相机的相对姿态角,但是由于云体边缘的模糊很难精确地找到同名点再加上相机内方位元素的影响,通过共面方程计算相对姿态角时会产生较大误差,有时无法收敛。本文通过选择特殊目标物星星,对相机进行现场标校,星星可近似认为是无穷远处目标物,当两个相机姿态角一致时从理论上在两个相机中星星成像位置是相同的。根据实测星星在图像中成像的位置差,计算两相机的相对姿态角。
在相机的内方位元素和畸变系数已知的情况下,由于星星为无穷远点,则两相机的相对位置关系等效为一个相机以焦点为中心点自身进行横滚 (α)、俯仰 (β)、旋转 (γ) 变换 (如图 6所示)。以焦点为中心建立坐标系,水平方向为x轴方向,则设水平相机中星星像点三维坐标为 (x1, y1, -f),根据旋转公式[3]得到旋转之后带有姿态角的相机中星星像点三维坐标为
|
(5) |
|
|
| 图 6. 相机姿态投影示意图 Fig 6. The diagram of camera attitude projection | |
其中, 旋转矩阵
|
根据式 (5),可以通过3个相互独立的星星坐标位置求取旋转矩阵R,设水平相机中3个星星坐标矩阵为

|
(6) |
某时刻拍摄到的星星位置坐标如表 1所示,利用式 (5)、式 (6) 计算出两相机的相对姿态角度差,横滚 (α)、俯仰 (β) 和旋转 (γ) 分别为0.328°,-3.33°和-2.339°。
|
|
表 1 某一时刻的星星位置 Table 1 The star positions |
3.2 基于相机姿态角的云图校正
在实际云底高度测量中先对带有姿态角相机做姿态校正之后再计算云底高度,根据式 (6) 可以得到云底高度测量修正公式
|
(7) |
式 (7) 中,dLx为由姿态角引起的像点在x轴分量上的位移量:
|
(8) |
对某一时刻采集到的云图做姿态校正,校正前后的云底高度直方图分布如图 7所示,横轴是每对同名点利用式 (1) 计算出来的云底高度,纵轴是在每个云底高度区间内出现同名点的个数,从图 7可以看出,姿态校正前后云底高度直方图的变化情况,姿态校正之前直方图谱宽较宽,具有多个峰值,校正之后云底高度直方图有明显改进。通过大量试验数据统计校正之后的云底高度直方图基本服从正态分布,可以通过平均的方法消除误差提高测量精度。
|
|
| 图 7. 2011年5月10日10:32云图经过姿态角校正前 (a)、校正后 (b) 的云底高度直方图分布 Fig 7. The histogram of cloud base height before (a) and after (b) attitude angle correction at 1032 BT 10 May 2011 | |
4 试验数据对比与误差分析 4.1 数据对比
该系统与安装在北京市观象台的维萨拉公司生产的CL31激光云高仪进行对比试验[16-17],2011年5月1日—6月30日,采集样本总数为17568个,包括无云天和有云天。双成像云高仪和激光云高仪同时获取到的总样本云底高度数据为488个,双成像云高仪相对于激光云高仪的误差统计如表 2所示,根据标准偏差可以看出,云底高度越高其标准偏差越大,与双成像云高仪测量精度随着云底高度增加而降低相一致。云底高度分布如图 8所示,其中云底高度样本序列按照维萨拉激光云高仪观测到的云底高度从高到低排序,云底高度散点图如图 9所示。
|
|
表 2 双成像云高仪相对于激光云高仪的误差统计 (单位:m) Table 2 The error of cloud base height between dual-camera stereovision and CL31 (unit:m) |
|
|
| 图 8. 2011年5月1日—6月30日双成像云高仪与维萨拉激光云高仪云底高度对比分布 Fig 8. Comparison of cloud base height by dual-camera stereovision and CL31 from 1 May to 30 June in 2011 | |
|
|
| 图 9. 双成像云高仪与维萨拉激光云高仪云底高度散点图 Fig 9. Scattered plot of cloud base height by dual-camera stereovision and CL31 | |
从图 8可以看出, 云底高度在5000 m以下,双成像云高仪略低于维萨拉激光云高仪测到的云底高度,云底高度误差平均值为143 m,标准偏差为435 m。从图 9可以看出,5000 m以下两种设备同时探测到的云底高度数据具有很高的相关性。存在较少偏差较大的点可能原因是双成像探测的是一定视场角范围内的平均云底高度,而激光云高仪探测的是单点云底高度,在空间上测量区域可能没有很好的对应上。
当云底高度大于5000 m时,从图 8可以看出,双成像数据波动比较大,云底高度误差平均值为46.5 m,标准偏差为1109.0 m。主要原因为仪器测量精度随着云底高度增加而降低,在5000 m左右时双成像云底高度数据存在一次跳变过程,平均高度略高于激光云高仪,与5000 m以内的云底高度趋势相反随着云底高度的增加而下降,主要原因是相机标校后仍存在一个小的Δx的位移量,以及相机内方位元素、相对姿态角和高空风等影响造成一定的系统误差。
根据式 (1) 得到带有系统误差的实测云底高度
由于采集卡同时采集两幅图像时产生延时,使云底高度测量产生误差,假设两相机同步采集时间差为dt,风向为θ,风速为v,本系统双成像基线方向为东西走向,根据云底高度测量式 (1),可以推导出由风速带来的云底高度测量误差为
|
(9) |
两幅图像采集延时可控制在30 ms以内,利用同一相机间隔4 s拍摄两张云图可以计算出云图中的每个同名点对应的云迹风,某时刻实测云迹风部分数据如表 3所示。
|
|
表 3 2011年5月21日某时刻实测云迹风部分数据 Table 3 The data of cloud motion winds on 21 May 2011 |
根据式 (9) 可以计算出由两相机采集云图时间延时造成最大云底高度测量误差平均值为±21.0 m左右。
4.2.2 相机成像分辨率对云底高度测量的影响在视场角一定的情况下图像分辨率与焦距长度成线性关系,所以提高图像分辨率能够提高双成像云底高度测量精度。
4.2.3 姿态角的精度分析姿态角对测量精度的影响是个复杂的过程,包括旋转、俯仰和横滚,在旋转为零的情况下,相对于基线方向的横滚对云底高度测量精度起主要作用,根据式 (8) 可以计算出由横滚带来云底高度测量相对误差η,
|
(10) |
当实际云底高度为5000 m,焦距f=690,基线长为60 m,横滚角为0.328°时,根据式 (10) 云底高度测量平均相对误差η=0.67,实测云底高度为3371 m,与实际云底高度相差较大。
4.2.4 云体对云底高度测量精度的影响由于云体是个弥散体,其底部的形状不规则,很难准确测量云底的高度。
角点主要集中在云体的边缘部分,所测的云底高度值实际上代表云底边缘部分的高度。
当云块不在相机天顶方向时,由于测量角度的原因,测量对象主要是云体高度,该系统中相机视场角相对较小,可以近似认为是天顶方向的云块。
4.2.5 安装位置的影响安装时将激光云高仪放在两相机基线中心位置左右与两相机的相对距离约25 m,由于相机有一定的视场角,测量范围包含了激光云高仪的测量“点”。获取的测量数据可近似地认为同时、同点,两者由于位置不同带来的云底高度测量误差相比于测量精度带来的误差较小,后期会做更深入的研究。
5 结论双目成像云底高度测量方法采用的是直接摄影测量方法,克服被动反演精度不高等不足,提高了测量精度。通过对双目成像云底高度测量方法的研究和影响测量精度因素的分析得到以下结论:
1) 通过图像处理、角点检测、区域相关匹配和基于小孔成像的云底高度测量方法能够实现部分云底高度的自动测量,通过星星相对位置关系对两相机进行现场标校,能够很好地解决两相机的相对姿态问题。通过与维萨拉激光云高仪云底高度测量对比分析得出双目成像云底高度测量系统具有较高的测量精度。
2) 通过对影响测量精度主要原因进行分析得出:提高相机分辨率、对相机内外方位元素进行标校、对两相机相对姿态角进行校正以及控制采样同步时间等,可以减小双目成像云底高度测量误差,提高测量精度。
3) 图像纹理清晰程度是影响云底高度数据获取率的主要原因,采用可见光图像传感器获取云图时受到光照和视程障碍等影响,只能在白天和云图纹理清晰的条件下观测,数据获取率较低。为了减少光照等影响,提高数据获取率,可采用红外图像传感器构成双目成像云底高度测量系统,实现云底高度测量。
| [1] | 张春光, 张玉钧, 韩道文, 等. 测云技术研究进展. 光散射学报, 2007, 19, (4): 388–394. |
| [2] | 高太长, 刘磊, 赵世军, 等. 全天空测云技术现状及进展. 应用气象学报, 2010, 21, (1): 101–109. DOI:10.11898/1001-7313.20100114 |
| [3] | 张剑清, 潘励, 王树根. 摄影测量学. 武汉: 武汉大学出版社, 2003. |
| [4] | Gabrlela S, Emmanuel P B, Annln G. Cloud mapping from the ground:Use of photogrammetric methods. Photogrammetric Engineering & Remote Sensing, 2002, 68, (9): 941–951. |
| [5] | Mark C A, Philip K Jr. The computation of cloud-base height from paired whole-sky imaging cameras. Machine Vision and Applications, 1995, 9, (4): 160–165. |
| [6] | 谭涌波, 陶善昌, 吕伟涛, 等. 双站数字摄像测量云底高度. 应用气象学报, 2005, 16, (5): 629–637. DOI:10.11898/1001-7313.20050509 |
| [7] | 白晓宁.空中目标的双目立体视觉测距.西安:西安电子科技大学, 2004. |
| [8] | 王建华, 韩红艳, 王春平, 等. CCD双目立体视觉测量系统的理论研究. 电光与控制, 2007, 14, (4): 95–116. |
| [9] | Harris C, Stephens M. A Combined Comer and Edge Detector. Proc 4th Alvey Vision Conference, 1988: 147–151. |
| [10] | Tissainayagam Suter. Assessing the performance of corner detectors for point feature tracking applications. Image and Vision Computing, 2004, 22, (8): 663–679. DOI:10.1016/j.imavis.2004.02.001 |
| [11] | Lucchese L, Mitra S K. Using Saddle Points for Subpixel Feature Detection in Camera Calibration Targets. Proc of the 2002 Asia Pacific Conference on Circuits and Systems, 2002: 191–195. |
| [12] | Chen Dazhi, Zhang Guangjun. A New Sub-Pixel Detector for X-Corners in Camera Calibration Targets. International Conference in Central Europe on Computer Graphics and Visualization—WSCG, 2005: 97–100. |
| [13] | Wang Sheyang, Song Shenmin, Qiang Wenyi, et al. Subpixel corner detection using spatial moment. Acta Automatica Sinica, 2005, 9, 31, (5): 714–719. |
| [14] | 梁志敏, 高洪明. 摄像机标校中亚像素级角点检测算法. 焊接学报, 2006, 27, (2): 102–104. |
| [15] | 邓志燕, 陈炽坤. 利用外极线约束的图像匹配新算法. 工程图学学报, 2009, 30, (5): 104–107. |
| [16] | 胡树贞, 马舒庆, 陶法, 等. 地基双波段测云系统及其对比试验. 应用气象学报, 2012, 23, (4): 441–449. DOI:10.11898/1001-7313.20120407 |
| [17] | 杨俊, 吕伟涛, 马颖, 等. 基于自适应阈值的地基云自动检测方法. 应用气象学报, 2009, 20, (6): 713–721. DOI:10.11898/1001-7313.20090609 |
2013, 24 (3): 323-331



