[发明专利]基于深度玻尔兹曼机提取高层特征的高分辨率遥感图像飞机检测方法在审
申请/专利权人：
公开/公告号：CNA
发明/设计人：;;;
公开/公告日：
主分类号：
分类号：;;
搜索关键词：
【说明书】：
技术领域本发明属于计算机视觉算法研究领域，涉及一种基于深度玻尔兹曼机提取高层特征的高分辨率遥感图像飞机检测方法，可以在遥感图像数据库中准确、鲁棒地检测出图像中的飞机。背景技术遥感技术的飞速发展促使许多卫星和航空传感器可以提供具有高空间分辨率的光学图像，这些图像有广泛的应用，如灾害治理，土地规划，监控和交通规划。在这些应用中，自然或人造物体的自动检测是一个基础任务，而且已经吸引了越来越多的研究兴趣。光学遥感图像的飞机中包含的丰富的空间信息和详细的结构信息，为我们解决这一挑战性的任务提供了新的机遇。特征表示在遥感图像的飞机检测任务中起着重要的作用，然而，目前存在的特征表示在遥感图像分析中的表现并不好，主要问题在于仅仅利用基于像素的频谱信息提取特征是不够的，因为这样忽略了图像的纹理空间信息，因而并不能真正挖掘出图像中飞机更深度的结构信息，随着遥感技术的进步，高分辨率遥感卫星和航拍器使图像空间和结构信息的获取成为可能，许多方法开始使用底层特征(如SIFT,HOG,GLCM)或是中层特征(如BOW和PLSA)表示图块，虽然这些方法可以在一定程度上提高分类的效果和检测的准确率，但是还会面临一些问题，特别是这些底层特征能表示的图片的空间几何信息是有限的，不能直接用于描述图块的结构信息，中层特征是基于底层特征的统计信息来提取空间区域的结构信息，但是对于背景复杂的情况下也不具备强大的表述能力。发明内容要解决的技术问题为了避免现有技术的不足之处，本发明提出一种基于深度玻尔兹曼机提取高层特征的高分辨率遥感图像飞机检测方法，通过深度玻尔兹曼机提取高层特征以得到由底层和中层特征编码构成的空间和结构信息。深度玻尔兹曼机是一种具有深度结构的无监督的特征学习工具，并且已被证明在高层特征学习上具有非常令人满意的表现。技术方案一种基于深度玻尔兹曼机的高分辨率遥感图像飞机检测方法，其特征在于：利用方形滑动窗口在拍摄的遥感图像中提取图块，然后利用训练得到的飞机检测器对图块进行分类，并利用非极大值抑制的方法保留得分最高的窗口以解决不同尺度下所取窗口可能高度重叠的问题，从而得到最终的飞机检测结果；所述飞机检测器的训练步骤如下：步骤1：对遥感图像数据库中多幅遥感图像进行正负样本分类，将遥感图像中含有飞机信息的遥感图像作为正样本图像，未含有飞机信息的遥感图像作为负样本图像；步骤2：利用多尺度方形滑动窗口在每一幅遥感图像中提取图块，将正样本图像得到的图块作为正样本集中图块负样本图像得到的图块作为负样本集中图块所述滑动步长为所选取窗口大小的1/3；步骤3、提取正负样本图像中所有图块的底层特征：利用尺度不变特征转换(SIFT)提取每一个图块中的关键点，作为图块的底层特征；步骤4、获取正负样本图像中所有图块的中层特征：利用局部约束线性编码(LLC)算法对底层特征编码得到中层特征；首先利用K-means聚类方法将底层特征聚类以获取字典，用D＝[d1,d2,...,dN]表示获取的底层特征，CB＝[cb1,cb2,...,cbM]表示聚类后的字典，对于底层特征中的每一个元素dn,n∈[1,N]，在CB中取5个与dn特征距离最近的元素作为dn的局部字典LBn，然后利用LBn对dn进行编码，编码所用函数为其中为局部编码系数，根据得到Cn，其中Cn为dn基于全局字典CB的编码系数，Cn为含有5个非零元素的M×1维向量,最后通过求取所产生的全局编码系数中每一维的最大值产生图块的中层特征；步骤5、学习正负样本图像中所有图块的高层特征：利用一个三层的深度玻尔兹曼机从中层特征中得到每一个图块的结构和空间格局，它包含一个可视层V∈{0,1}M和两个隐含层其中H1和H2分别表示第一个隐含层和第二个隐含层的节点数目，整个单元{v,h1,h2}的能量函数为E(v,h1,h2,Θ)＝-vTW1h1-h1TW2h2，其中Θ＝{W1,W2}是模型的参数，W1和W2分别表示可视层到隐含层和隐含层到隐含层之间的对称交互项，通过最小化能量函数学习深度玻尔兹曼机模型中的参数，并将第三层的输出作为每个图块的高层特征；步骤6、训练分类器：基于提取的图块的高层信息，训练检测飞机的分类器,训练所用到的的得分函数为其中w1和b1分别是支持向量机的决策面和偏置项，为正样本集中的图块，为图块的高层特征，该得分函数大于阈值时，表示图块是飞机，否则表示图块不是飞机，所述阈值为其中取0或1，表示图块是否为飞机，P表示所有正样本图像中的图块总数，通过计算该得分函数判断图块是否包含飞机。所述多尺度方形滑动窗口大小为size∈{60,100,135}。所述遥感图像数据库采用Google Earth库。所述步骤3中的尺度不变特征转换SIFT采用论文B.Sirmacek,and C.Unsalan,“Urban-area and building detection using SIFT keypoints and graph theory,”IEEE Trans.Geosci.Remote Sens.,vol.47,no.4,pp.09.中的方法。所述步骤4中的局部约束线性编码LLC采用论文J.Wang,J.Yang,K.Yu,F.Lv,T.Huang,and Y.Gong,Locality-constrained linear coding for image classification,inProc.IEEE put.Vision Pattern Recognit.,2010.pp..中的方法。所述步骤5中的深度玻尔兹曼机采用论文R.Salakhutdinov,and G.E.Hinton,Deep boltzmann machines,inProc.Int.Conf.Artificial Intelligence and Statistics,2009.pp.448-455.中的方法。所述非极大值抑制采用论文G.Cheng,J.Han,L.Guo,X.Qian,P.Zhou,X.Yao,and X.Hu,“Object detection in remote sensing imagery using a discriminatively trained mixture model,”ISPRS J.Photogramm.Remote Sens.,vol.85,pp.32-43,2013.中的方法。有益效果本发明提出的一种基于深度玻尔兹曼机提取高层特征的高分辨率遥感图像飞机检测方法，首先将图片分成若干图块，然后利用利用尺度不变特征转换(SIFT)提取图块中的关键点，作为图块的底层特征，再利用局部约束线性编码(LLC)算法对底层特征编码得到中层特征，接着利用一个三层的深度玻尔兹曼机从中层特征中得到图块的高层特征，然后利用高层特征训练支持向量机分类器，最后用该分类器检测测试图片的飞机，取得更具准确性，鲁棒性的飞机检测结果。本发明利用深度玻尔兹曼机将中层特征映射到高层中以提高各种不同场景和飞机的表示效果，高层特征的提取分三个阶段：1.利用尺度不变特征转换(SIFT)得到图片的局部底层特征；2.利用局部约束线性编码(LLC)算法对底层特征编码得到中层特征；3.利用深度玻尔兹曼机从中层特征中学习得到表现力更强的高层特征。将得到的图块的高层特征训练支持向量机分类器，利用训练得到的分类器取得更具准确性，鲁棒性的飞机检测结果。附图说明图1：本发明方法的基本流程图图2：实验结果图图3：PR结果图具体实施方式现结合实施例、附图对本发明作进一步描述：用于实施的硬件环境是：Intel Pentium 2.93GHz CPU计算机、2.0GB内存，运行的软件环境是：Matlab R2011b和Windows XP。实验选取了Google Earth库中的所有图像作为测试数据，该数据库中包含120幅高分辨率遥感图像，是国际公开的用于测试高分辨率遥感图像飞机检测计算模型的数据库。
友情链接：交换友情链接需要网站权重大于2，网站收录10W以上，如符合条件，请联系QQ：。
行业网站：相关推荐：
400-周一至周五 9:00-18:00
服务热线：400-投诉建议：022-
扫一扫，微信关注高智网
高智&让创新无法想象2000万件&专利数据查看: 1132|回复: 18
在三维场景下，影像、高程数据结合快速查看
马上注册地信网，享受更多功能，学习更多知识，成就人生精彩！
才可以下载或查看，没有帐号？
& && & LocaSpace Viewer具备便捷的影像、高程、倾斜摄影数据阅读功能。用户能够快速地浏览、测量、分析和标注三维地理信息数据，可以添加tif图层、shp图层、kml图层等多种格式数据。实现三维场景的飞行浏览和多视角浏览，快捷的对地理信息数据进行格式转换。& && & 工具：
LocaSpace Viewer影像、高程、倾斜摄影数据等
& && & 步骤：
打开LocaSpace Viewer
开启界面.png (293.1 KB, 下载次数: 22)
15:11 上传
加载图层——单击 开始-&加载图层 ，即可弹出如图本地文件对话框，找到要加载图层数据的位置，选取后点击打开即可在右侧视窗中看到加载的图层。
打开图层.png (90.36 KB, 下载次数: 23)
17:38 上传
影像数据格式的批量转换——单击 开始-&影像拼接 ，弹出如图对话框，
影像批处理.png (33.19 KB, 下载次数: 20)
17:39 上传
点击添加数据，可以添加单个的影像数据。点击添加目录，可以同时添加整个文件夹的多个影像数据。点击移除选中，可以选择对话框中已添加的单个数据进行移除。点击移除全部，可以选择对话框中已添加的所有数据进行同时移除。在输出路径可以设置影像拼接后文件的存储位置。启用无效值过滤，可以设置特定的无效值进行过滤。
影像数据的拼接—单击 开始-&影像拼接 ，弹出如图对话框，点击添加数据，可以添加单个的影像数据。
影像拼接.png (32.78 KB, 下载次数: 24)
17:37 上传
& && &点击添加目录，可以同时添加整个文件夹的多个影像数据。& && &点击移除选中，可以选择对话框中已添加的单个数据进行移除。& && &点击移除全部，可以选择对话框中已添加的所有数据进行同时移除。& && &在输出路径可以设置影像拼接后文件的存储位置。& && &启用无效值过滤，可以设置特定的无效值进行过滤。
地形数据格式的批量转换——
& && & 通过单击 开始-&地形批处理 ，弹出如图对话框，
地形批处理.png (34.45 KB, 下载次数: 22)
17:37 上传
& && & 点击添加数据，可以添加单个的影像数据。& && & 点击添加目录，可以同时添加整个文件夹的多个影像数据。& && & 点击移除选中，可以选择对话框中已添加的单个数据进行移除。& && & 点击移除全部，可以选择对话框中已添加的所有数据进行同时移除。& && & 点击确定开始后，能够看到处理进度。
地形数据的拼接——
& && & 单击 开始-&地形拼接 ，弹出如图对话框，
地形拼接.png (32.13 KB, 下载次数: 22)
17:37 上传
& && & 点击添加数据，可以添加单个的影像数据。& && & 点击添加目录，可以同时添加整个文件夹的多个影像数据。& && & 点击移除选中，可以选择对话框中已添加的单个数据进行移除。& && & 点击移除全部，可以选择对话框中已添加的所有数据进行同时移除。& && & 在输出路径可以设置影像拼接后文件的存储位置。& && & 启用无效值过滤，可以设置特定的无效值进行过滤。效果如图
suox.png (469.24 KB, 下载次数: 20)
17:41 上传
楼主可以发个安装包给大家吗！谢谢
挺好的软件！
楼主可以发个安装包给大家吗！谢谢
百度搜索LocaSpace Viewer可以去软件官网免费下载，软件大概25兆，这边超过上传的限制。
楼主可以发个安装包给大家吗！谢谢
安装包太大，系统上传不了。搜locaspace viewer可以上官网下载使用。
这个球快不快的？三维球很多，是基于osg改的吗？
那个淹没模拟能不能设置的啊，看起来不错。下载来试试,谢谢楼主！
界面看起来还不错
那个倾斜摄影的数据是自己的还是软件上可以看到的？
挺不错的软件，赞一个
论坛推荐话题 /1
地信网全新手机版隆重上线了！手机“每日签到”3秒搞定！即刻起用手机访问地信网，体验强大功能！
热线: 400-
Powered by利用卫星遥感数据中的云影信息检测云高的方法
专利名称利用卫星遥感数据中的云影信息检测云高的方法
本发明涉及大气探测领域，特别涉及一种利用卫星遥感数据中的云影信息检测云高的方法。
云在天气系统的发展、降水的形成和大气辐射传输等物理过程中都扮演着极其重要的角色。如何准确、及时地获取云的信息，对于天气气候学研究、气象预报、人工影响天气以及国民经济和军事等诸多领域都十分重要。通常对云的气象观测主要有云状、云量和云高三项，其中云高是指检测站地平面至云底的垂直距离，通常以米为单位。
目前，对于云量和云状的观测主要还是靠目测，但云高测量的各项技术业已成熟，出现了许多用于测量云高的仪器。因此，云高测量方法大致可分为以下三类估测云高、仪器检测云高、利用卫星数据和空间几何反演获得云高。
(一)估测云高包括目测和经验公式估算。
主要是根据云的形状来估测云高，所以必须正确地判定云的形状，然后再根据云的结构、大小、亮度、颜色、移动速度等情况，结合当地经常出现的云高范围进行估计。
2.经验公式估算
除了目测云高的方法外，还可以利用经验公式估算云高，这主要针对某些特定云的估测，需要专业技术人员才能估算。
(二)仪器测云高。
1.云幕灯测云高
云幕灯是一种云高指示器，是最早用来测云高的仪器，主要由光源和测角器两部分组成，通过测得对应的仰角、光源和测高仪之间的距离是来计算云高。
2.旋转光束云高仪测云高
该系统由一个发射器、一个接收器和一台记录设备组成。通过测试，结合相应公式来计算云高。
3.气球测云高
主要利用一只充满氢气或氦气的橡胶气球，从地面匀速至云底所经历的时间来测量云高。
4.微脉冲激光雷达(Micropulse Lidar，MPL)测云高
MPL探测云高，主要是对激光回波数据进行采集，通过数据分析反演云层高度。由于云的种类繁多，形状千姿百态，每块云中云滴的尺度、浓度和含水量不同，具有强散射性，所以会引起后向散射的激光回波也不相同，所以由该仪器测量确定的云底始终是一个有争议的问题。
5.半导体激光云高仪
半导体激光云高仪是激光雷达的一种，主要用来测量云的高度。半导体激光云高仪是近年来发展起来的，采用固态的近红外激光二极管和数字信号处理技术，独特的光学镜头设计，小型性便携，能使实时在线监测等优点是其它仪器不可比拟的，所以成为目前研究的一个热点领域。
通过对以上云高探测方法的分析比较得出云幕灯和旋转光束云高仪采用的都是普通的光源，探测距离不高，只能用于低云探测，在降雨时测量不准确，但是，成本低、方法简单、容易操作。探空气球测云成本高，不容易操作，很少使用，但是能够同时探测到其它的一些气象参数。随着激光和电子技术的发展，测云激光雷达(激光测云仪)技术得到很大提高，很多国家都在花费大量的人力、物力进行研究，激光雷达是已成为测量云层高度的有效工具。目前，国内外对于云高的测量，主要采用半导体激光云高仪来测量。国外对云高仪的研制比较早，技术相对比较成熟，已经研发出很多测云的仪器；而国内在这方面起步较晚，主要还是靠进口的仪器，与国外还有一定技术差距。
(三)利用卫星数据信息和空间几何反演获得云层高度。
1.双站数字摄像测量云高方法
2005年，中国科学技术大学地球和空间科学学院与中国气象科学研究院发表双站数字摄像测量云高方法，其基本原理是利用双站测距法。不足之处在于该测试方法受观测场地和光学机械加工精度的限制，测试中需对测量系统(CCD等)进行精密的设计和严格的校准，否则会产生较大测试误差。
2.基于MODIS卫星数据的二氧化碳(CO2)薄片算法
MODIS的全称为中分辨率成像光谱仪，是搭载在Terra和Aqua卫星上的一个重要的传感器，MODIS是第一个拥有CO2薄片波段的高空间分辨率探测器。CO2薄片算法是使用被动遥感资料反演云顶性质的重要方法之一。CO2薄片技术推断云顶气压和有效云量或辐射，是利用15μm CO2吸收带的几个吸收不同的MODIS红外波段，这几个波段对于大气不同高度的层次都有不同的敏感度。CO2薄片技术可以反演出云顶气压和有效云量主要是中高层云，但对于云顶气压大于700hPa(接近地面)的云层，因为云信息的减少，所以不能再用这种方法。因此，CO2薄片技术对低层云的反演存在一些问题，对于云高检测有一定局限性。
3.双星立体观测云顶高度几何反演方法
2008年，北京大学物理学院和中国民用航空总局天津空中交通管理站气象台联合发表双星立体观测云顶高度几何反演方法。该方法存在的问题是云顶高度反演方法实验中采用的“立体图像对”是卫星立体观测模拟资料。该方法中的立体测量依赖两颗卫星对给定区域的同步观测，然后对卫星资料进行云顶高度的反演。但是，目前该方法中提到的两颗卫星同步观测是很难达到的，即使按预定的相同时间开始摄像，也很难保证卫星所覆盖的区域相同。因此，两颗卫星测量过程中的不同步，就会导致云顶高度的计算出现误差。
综上所述，目前常用的云高测量方法都存在一些局限性激光云高仪只能在测试区域中有限的地区进行测量，而且一台仪器同一时刻只能测天空中云的一个点的云高，存在测量值的局限性问题；双站数字摄像测量云高方法需要采用图像采集仪器进行本地测量，测量仪器调整过程比较复杂，仪器的调整会极大的影响检测结果的精度；基于MODIS遥感图像CO2薄片测云高的方法，对低云高度的检测有一定局限性；双星立体观测云顶高度几何反演方法对于卫星立体观测的模拟资料可以得到较好的检测结果，对于真实数据的检测精度还未证实；其他的云高测试方法也有一些局限性。因此，研究精度高、可靠性好、经济实用、自动化高的探测方法具有重要的科学意义和应用价值。
本发明的目的在于提供一种利用卫星遥感数据中的云影信息检测云高的方法，以解决目前常用的云高测量方法都存在局限性、精度低、可靠性差、成本高及自动化成度低等问题。
MODIS(中分辨率成像光谱仪)最大空间分辨率可达250米，扫描宽度2330公里，是当前世界上新一代“图谱合一”的光学遥感仪器。MODIS仪器与其它卫星相比，有以下特点和优势1)空间分辨率大幅提高。空间分辨率提高了一个量级，由NOAA气象卫星的千米级提高到了MODIS的百米级。2)时间分辨率有优势。一天可过境4次，对各种突发性、快速变化的自然灾害有更强的实时监测能力。3)光谱分辨率大大提高。具有36个波段，这种多通道观测大大增强了对地球复杂系统的观测能力和对地表类型的识别能力。4)MODIS数据全球免费接收。该数据可以从国内、外MODIS卫星数据地面接收站免费接收。
本发明以MODIS遥感图像处理技术为基础，提出平原区域云高测量的新方法。其应用条件是分辨率为250m的MODIS数据卫星遥感数据中的平原区域，使用的MODIS数据是经过常规纠正处理后的数据。例如选取250m的MODIS数据卫星遥感数据MOD02QKM.A0.005.5.hdf。
本发明的技术方案包括以下步骤
一、插值计算由于需要云和云阴影点的具体地理信息数据，包括经度、纬度、太阳天顶角、太阳方位角、卫星天顶角、卫星方位角数据，因此需要MODIS数据星历表MOD03的数据，即MOD03.A0.005.5.hdf。针对250m分辨率的MOD03数据中的某些数据量重叠的情况，需要从1000m分辨率的数据中重新计算获得250m分辨率对应的MOD03数据，才能获得云和云阴影上点的更加精确的信息。因此这里采用通过采用线性插值公式解决该问题。线性内插是假设在二个已知数据中的变化为线性关系，因此可由已知二点a和c，去计算内插点b的值，公式如下
其中a＜b＜c在上式的b点即是代表要内插的点，f(b)则是要计算的内插函数值。
具体实现步骤如下
(1)读取1000m MOD03中的地理信息数据，如经度。
(2)1000m MOD03数据格式是2030行，1354列；插值后的250m MOD03数据格式是8120行，5416列，其行和列大小分别为1000m数据的4倍。采用将1000mMOD03数据中每10行×1354列的数据通过公式(1)插值成为40行×5416列的数据。由于1000m MOD03数据有2030行，因此需要进行203次的插值运算。
(3)利用Matlab软件中的INTERP2函数可以实现线性内插，通过计算得到250mMOD03中的地理信息数据。
二、云和云阴影检测从原图像中手动选取具有云和云阴影的图像，选取图像在原图中的位置信息是起始点在3960行、4480列，即()，范围是40×40大小的图像范围，选取的波段是可见光波段2。该步骤主要完成从选取图像中检测云和云阴影图像的具体位置，选取数据是MOD02QKM。
目前，绝大多数的遥感图像云特征建模都是从图像的光谱特征出发，考虑云与其它地物在可见光、近红外和热红外波段上的差异。云在可见光和红外波段与植被、土壤以及水域等下垫面介质的反射率和辐射亮温值存在差异，云具有较高的反射率而具有低的亮温值等特点进行云与地物的区分。在MODIS数据的可见光波段0.645μm处，晴空数据地物一般具有较低的反射率，反射值在0～0.3之间，云则有高的反射率，一般高于0.3，是进行云检测处理的首选波段。1.375μm处于强的水汽吸收波段，地面的反射因水汽吸收难以到达传感器，导致反射率低，而高云上方的水汽对云的反射值削弱较少，因此该波段对高层卷云敏感，具有高的反射率，低云则具有低的反射率，通常低于0.02。此波段可以用来进行卷云的检测。8.6μm和11μm波段对云特性敏感。短波中有强的水汽吸收作用，导致晴空中8.6μm处的辐射值低于11μm处的辐射值；冰晶在8.6μm的吸收远小于11μm处的吸收，所以，冰晶组成的卷云在8.6μm处的辐射值高于11μm处的辐射值；而水滴组成的光学厚度较大的云则在8.6μm处的辐射值低于11μm处的辐射值。
采用多光谱云检测主要是利用以上各个波段对云的敏感特性，以及MODIS数据的各通道的波谱特性实现云与地物的区分。在多光谱云检测结果的基础上，对于可见光波段的图像，结合常用的图像检测方法如区域标定法，区域生长法，分水岭法，阈值分割法等将云与非云区域进行更有效的区分。对于云阴影的检测采取同云检测一样，首先利用云影在各个波段的敏感特性，主要利用1.2、0.94、0.87、0.66μm这四个波段反射率来检测云阴影；另外，利用3.9μm和4.0μm波段的亮度温度差大于零也可以完成云阴影的检测。
这里采用双阈值分割法在可见光波段2中实现云和云影图像的检测。图像阈值分割是一种广泛使用的图像分割技术，它利用了图像中要提取的目标与其背景在灰度特性上的差异，把图像视为具有不同灰度级的两类区域(目标和背景)的组合，选取一个合适的阈值，以确定图像中的每一个像素点应该属于目标还是背景区域，从而产生相应的二值图像。采用的双阈值分割的数学表达式为(2)式所示
通过对云和云影图像的光谱分析，选取云图像的阈值是T2，阴影图像的T1，f(i，j)是图像中像素的值。然后采用区域标定法，在双阈值检测结果的基础上，将云图像的位置和云阴影图像提取出来。
区域标定是给连接在一起的像素(连接成分)附上相同的标记，不同的连接成分附上不同的标记。区域标定在二值图像处理中占有非常重要的地位。通过该处理将各个连接成分区分开，然后就可以调查各个连接成分的形状特征。
区域标定也有许多方法，本发明采用了一种简单的方法。步骤如下
(1)扫描图像，遇到没加标记的目标像素(白像素)P时，添加一个新的标记(LABEL)。
(2)将与P连接在一起(即相同连接成分)的像素添加相同的标记。
(3)进一步给所有与加标记像素连接在一起的像素添加相同的标记。直到连接在一起的像素全部被添加标记。这样，一个连接成分就被添加了相同的标记。
(4)回到第1步，重新查找新的没加标记的像素，重复上述各个步骤。图像全部被扫描后，处理结束。
采用区域标定后，会检测出很多候选目标。通过设定面积参数，使得云区域和云阴影区域能够从候选目标中提取出来，同时去掉图像中噪声对检测结果的影响。
三、MODIS数据的云和云影匹配
根据云和云影检测的结果，进行云和云影的图像匹配。可以采用数字图像处理学中的图像匹配方法进行云和云影匹配，通过匹配得到云和云影之间的距离。
在云影图像中，通过云影检测，云和阴影区域已确定。每个云都有对应的影子图像，而且通常云阴影的图像要小于云的图像。基于模板匹配的原理，求得云与影间的距离的匹配方法为以云阴影图像为模板，搜索云图像中与阴影图像匹配的位置。图像匹配中，通过采用自相关算法来计算模板与搜索区域的相关程度，将具有最大相关度的位置定为匹配最佳位置。如S为云阴影图像，C为云图像，M×N是云阴影图像的大小，则相关系数Coef(coefficient自定义的参数)计算公式如(3)所示，(i，j)为云和云阴影图像的匹配相对位置矢量，i为水平相对位置矢量，j为垂直相对位置矢量。
其中各变量为(4)式和(5)式
通过公式(6)得到Coef最大时对应的最佳匹配相对位置矢量(i，j)
将所有检测到的云图像的像素信息代入公式(5)的C，所有检测到的云阴影图像的像素信息代入S，根据公式(3)可以得到多个相关系数Coef。利用公式(6)选取最大的相关系数Coef，并据此得到最佳匹配相对位置矢量(i，j)。
四、建立云高模型求得云层高度
根据云与云影之间的距离和卫星检测到的数据，通过云高模型的计算得到云层高度。由于空间中存在云，所以太阳的光线会被云完全或者部分遮挡，因此在地面上形成云的阴影图像。云阴影图像的位置与太阳相对于空间中云的位置有关，太阳-真云-云阴影共线；而当卫星观测地球时，空间中的云在MODIS探测器成像的位置与卫星的拍摄角度有关，即卫星-真云-投影云共线；因此MODIS数据中，云阴影图像和投影云图像都与空间中的云有着密切关系。可以根据云阴影图像、投影云图像和空间中的真云三者之间的空间位置关系建立模型，求得空间中的真云相对于地面的高度。
由于测试选择平原区域的云高测量，云和对应云阴影在同一平原区域，所以可以忽略地球是一个椭圆球体，并近似定义投影云和云阴影的地理信息在同一个平面上。
已知空间云上一点Point，太阳对其辐射形成的云阴影图像上对应的点是点Shadow，卫星对空间云中点Point投影后对应的点是点cloud。通过云和云影检测以及云和云影匹配处理，可以得到点shadow和点cloud之间的距离。从MOD03数据中读出点shadow和点cloud包含的地理信息，这里选取cloud点的位置信息是3971行、4502列，即()，则由最佳匹配相对位置(2，9)可知其对应的点shadow的坐标为()。根据这两点的位置从MOD03中读取点shadow的纬度lat1＝34.3604°、经度lon1＝-106.2759°，点cloud的纬度lat2＝34.3468°、经度lon2＝-106.2365°，将图像上云和云影像素之间的距离(单位像素)转换为实际地理上的距离Distance(单位公里)，所用计算公式为(7)-(9)
其中a＝lat1-lat2(8)
b＝lon1-lon2
经过公式计算，Distance＝3.9247(km)。由于测试选择平原地区的云和云阴影图像，因此可以忽略云和云阴影的高度差异，认为点shadow和点cloud是在同一平面上。从MOD03数据中读出点shadow的太阳方位角SunAZ＝134.905°＝2.3545rad、太阳天顶角SunZ＝14.954°＝0.261rad、点cloud的卫星方位角SensorAZ＝-77.2409°＝4.93508rad、卫星天顶角SensorZ＝41.2025°＝0.71912rad，将其结合云和云影距离等信息建立云高的空间几何模型，求得空间云中的一点Point相对于地面的高度Height。
令点Point为点1，其坐标为(x1，y1，z1)；点Shadow为点2，其坐标为(x2，y21，z2)；点Cloud为点3，其坐标为(x3，y3，z3)。
已知SunZ，SensorZ，因此A1，A2可知，即有公式(10)
SunZ＝A1，SensorZ＝A2(10)
已知SunAZ和Sensor AZ，因此B1，B2可求，即有公式(11)
B1＝π-SunAZ，B2＝SensorAZ-π(11)
令点1到点3的单位向量为
点1’的坐标为(x1’，y1’，z1’)，该向量在x、y、z轴的投影分别为
则有公式(12)-(15)
令点1到点2的向量为
该向量在x、y、z轴的投影分别为
且有公式(16)-(19)
由向量内积公式有公式(20)
的夹角，0≤α≤π
同理令点2到点3的单位向量为
点2’的坐标为(x2’，y2’，z2’)，该向量在x、y、z轴的投影分别为
则有公式(24)-(27)
令点2到点1的向量为
该向量在x、y、z轴的投影分别为
且有公式(28)-(31)
由向量内积公式有
的夹角，0≤β≤π
因此可得公式(36)
γ＝π-α-β(36)
令D13为点1到点3的距离，令D23为点2到点3的距离，令D12为点1到点2的距离，由三角形正弦定理有
且有D12＝‖r12‖(38)
点3的坐标3(x3，y3，z3)可用点1表示为
x3＝x1+D13sinA1cosB1(41)
y3＝y1+D13sinA1sinB1(42)
z3＝z1+D13cosA1 (43)
或者用点2表示为
x3＝x2+D23sinA2cosB2(44)
y3＝y2+D23sinA2sinB2(45)
z3＝z2+D23cosA2 (46)
理论上的云高Height的值计算为公式(47)
Height＝z1+D13cosA1＝z2+D23cosA2(47)
利用上述公式，将已知数据代入，可以求得云高。
五、云高计算和分析寻找符合云高检测方法测试条件的时间和地区，确定测试地点；通过本发明提出的云高检测方法，将检测云的MODIS空间数据信息与云和云影之间的距离代入云高模型，对该检测云的多个点进行云高计算，并将计算结果求均值，作为检测云的云层高度，采用一至四的步骤，对周围9点进行云高计算，最后求得的云的云层高度的均值。
本发明的有益效果在于可以利用MODIS遥感图像数据来实现满足检测条件的云高的检测，测试条件主要是MODIS卫星拍摄的平原区域中有云和云影的信息，从而完成对拍摄区域中单层云高的检测。该方法可以在不借助其它仪器的条件下，利用MODIS遥感图像数据的信息获得云层高度。与传统激光探测云高的方法相比，本发明中提出的方法可以实现从多点进行云高检测，进而提高云高检测结果的可靠性。该方法研究的最重要意义是可以完成对全球区域内满足测试条件的云进行高度测量，解决了云高测量区域局限性的问题。
图1是本发明实施例1选取出的云和云阴影图像；
图2是本发明实施例1检测得到的云和云阴影图像；
图3是本发明实施例1匹配的云和云阴影图像；
图4是本发明云高模型建立的理论基础图；
图5是本发明云高模型建立的几何模型；
图6是本发明云高的具体计算几何模型；
图7是本发明实施例2选取出的云和云阴影图像；
图8是本发明实施例2检测得到的云和云阴影图像；
图9是本发明实施例2匹配的云和云阴影图像。
具体实施例方式 实施例1
本发明利用卫星遥感数据中的云影信息检测云高的方法，该方法的应用条件是分辨率为250m的MODIS卫星遥感数据中的平原区域的数据。
本发明以MODIS遥感图像处理技术为基础，提出平原区域云高测量的新方法。其应用条件是分辨率为250m的MODIS数据卫星遥感数据中的平原区域，使用的MODIS数据是经过常规纠正处理后的数据。例如选取250m的MODIS数据卫星遥感数据MOD02QKM.A0.005.5.hdf。
本发明的技术方案包括以下步骤
一、插值计算由于需要云和云阴影点的具体地理信息数据，包括经度、纬度、太阳天顶角、太阳方位角、卫星天顶角、卫星方位角数据，因此需要MODIS数据星历表MOD03的数据，即MOD03.A0.005.5.hdf。针对250m分辨率的MOD03数据中的某些数据量重叠的情况，需要从1000m分辨率的数据中重新计算获得250m分辨率对应的MOD03数据，才能获得云和云阴影上点的更加精确的信息。因此这里采用通过采用线性插值公式解决该问题。线性内插是假设在二个已知数据中的变化为线性关系，因此可由已知二点a和c，去计算内插点b的值，公式如下
其中a＜b＜c在上式的b点即是代表要内插的点，f(b)则是要计算的内插函数值。
具体实现步骤如下
(1)读取1000m MOD03中的地理信息数据，如经度。
(2)1000m MOD03数据格式是2030行，1354列；插值后的250m MOD03数据格式是8120行，5416列，其行和列大小分别为1000m数据的4倍。采用将1000mMOD03数据中每10行×1354列的数据通过公式(1)插值成为40行×5416列的数据。由于1000m MOD03数据有2030行，因此需要进行203次的插值运算。
(3)利用Matlab软件中的INTERP2函数可以实现线性内插，通过计算得到250mMOD03中的地理信息数据。
二、云和云阴影检测从原图像中手动选取具有云和云阴影的图像，选取图像在原图中的位置信息是起始点在3960行、4480列，即()，范围是40×40大小的图像范围，选取的波段是可见光波段2。具体选取图像如图1所示。该步骤主要完成从选取图像中检测云和云阴影图像的具体位置，选取数据是MOD02QKM。
目前，绝大多数的遥感图像云特征建模都是从图像的光谱特征出发，考虑云与其它地物在可见光、近红外和热红外波段上的差异。云在可见光和红外波段与植被、土壤以及水域等下垫面介质的反射率和辐射亮温值存在差异，云具有较高的反射率而具有低的亮温值等特点进行云与地物的区分。在MODIS数据的可见光波段0.645μm处，晴空数据地物一般具有较低的反射率，反射值在0～0.3之间，云则有高的反射率，一般高于0.3，是进行云检测处理的首选波段。1.375μm处于强的水汽吸收波段，地面的反射因水汽吸收难以到达传感器，导致反射率低，而高云上方的水汽对云的反射值削弱较少，因此该波段对高层卷云敏感，具有高的反射率，低云则具有低的反射率，通常低于0.02。此波段可以用来进行卷云的检测。8.6μm和11μm波段对云特性敏感。短波中有强的水汽吸收作用，导致晴空中8.6μm处的辐射值低于11μm处的辐射值；冰晶在8.6μm的吸收远小于11μm处的吸收，所以，冰晶组成的卷云在8.6μm处的辐射值高于11μm处的辐射值；而水滴组成的光学厚度较大的云则在8.6μm处的辐射值低于11μm处的辐射值。
采用多光谱云检测主要是利用以上各个波段对云的敏感特性，以及MODIS数据的各通道的波谱特性实现云与地物的区分。在多光谱云检测结果的基础上，对于可见光波段的图像，结合常用的图像检测方法如区域标定法，区域生长法，分水岭法，阈值分割法等将云与非云区域进行更有效的区分。对于云阴影的检测采取同云检测一样，首先利用云影在各个波段的敏感特性，主要利用1.2、0.94、0.87、0.66μm这四个波段反射率来检测云阴影；另外，利用3.9μm和4.0μm波段的亮度温度差大于零也可以完成云阴影的检测。
这里采用双阈值分割法在可见光波段2中实现云和云影图像的检测。图像阈值分割是一种广泛使用的图像分割技术，它利用了图像中要提取的目标与其背景在灰度特性上的差异，把图像视为具有不同灰度级的两类区域(目标和背景)的组合，选取一个合适的阈值，以确定图像中的每一个像素点应该属于目标还是背景区域，从而产生相应的二值图像。采用的双阈值分割的数学表达式为(2)式所示
通过对云和云影图像的光谱分析，选取云图像的阈值是T2＝0.315μm，阴影图像的T1＝0.134μm，f(i，j)是图像中像素的值。然后采用区域标定法，在双阈值检测结果的基础上，将云图像的位置和云阴影图像提取出来。
区域标定是给连接在一起的像素(连接成分)附上相同的标记，不同的连接成分附上不同的标记。区域标定在二值图像处理中占有非常重要的地位。通过该处理将各个连接成分区分开，然后就可以调查各个连接成分的形状特征。
区域标定也有许多方法，本实施例采用了一种简单的方法。步骤如下
(1)扫描图像，遇到没加标记的目标像素(白像素)P时，添加一个新的标记(LABEL)。
(2)将与P连接在一起(即相同连接成分)的像素添加相同的标记。
(3)进一步给所有与加标记像素连接在一起的像素添加相同的标记。直到连接在一起的像素全部被添加标记。这样，一个连接成分就被添加了相同的标记。
(4)回到第1步，重新查找新的没加标记的像素，重复上述各个步骤。图像全部被扫描后，处理结束。
采用区域标定后，会检测出很多候选目标。通过设定面积参数，使得云区域和云阴影区域能够从候选目标中提取出来，同时去掉图像中噪声对检测结果的影响。经过云和云阴影检测后的图像如图2所示。其中灰色的像素代表背景，黑色的像素代表云阴影图像，白色的像素代表云图像。从图像可见，该方法可以较好地从原图像中检测出云和云影图像。
三、MODIS数据的云和云影匹配
根据云和云影检测的结果，进行云和云影的图像匹配。可以采用数字图像处理学中的图像匹配方法进行云和云影匹配，通过匹配得到云和云影之间的距离。
在云影图像中，通过云影检测，云和阴影区域已确定。每个云都有对应的影子图像，而且通常云阴影的图像要小于云的图像。基于模板匹配的原理，求得云与影间的距离的匹配方法为以云阴影图像为模板，搜索云图像中与阴影图像匹配的位置。图像匹配中，通过采用自相关算法来计算模板与搜索区域的相关程度，将具有最大相关度的位置定为匹配最佳位置。如S为云阴影图像，C为云图像，M×N是云阴影图像的大小，则相关系数Coef(coefficient自定义的参数)计算公式如(3)所示，(i，j)为云和云阴影图像的匹配相对位置矢量，i为水平相对位置矢量，j为垂直相对位置矢量。
其中各变量为(4)式和(5)式
通过公式(6)得到Coef最大时对应的最佳匹配相对位置矢量(i，j)
将所有检测到的云图像的像素信息代入C，所有检测到的云阴影图像的像素信息代入S，根据公式(3)可以得到多个相关系数Coef。利用公式(6)选取最大的相关系数Coef，并据此得到最佳匹配相对位置矢量(i，j)。通过计算，最大的相关系数是0.95，对应的(i，j)为(2，9)，即云和云影图像的对应匹配点水平方向相差9个像素，垂直方向相差2个像素。匹配效果如图3所示。
四、建立云高模型求得云层高度
根据云与云影之间的距离和卫星检测到的数据，通过云高模型的计算得到云层高度。图4为云高模型建立的理论基础图。由于空间中存在云，所以太阳的光线会被云完全或者部分遮挡，因此在地面上形成云的阴影图像。云阴影图像的位置与太阳相对于空间中云的位置有关，太阳-真云-云阴影共线；而当卫星观测地球时，空间中的云在MODIS探测器成像的位置与卫星的拍摄角度有关，即卫星-真云-投影云共线；因此MODIS数据中，云阴影图像和投影云图像都与空间中的云有着密切关系。可以根据云阴影图像、投影云图像和空间中的真云三者之间的空间位置关系建立模型，求得空间中的真云相对于地面的高度。
图5是云高模型建立的几何模型由于测试选择平原区域的云高测量，投影云和对应云阴影在同一平原区域，所以可以忽略地球是一个椭圆球体，并近似定义云和云阴影的地理信息在同一个平面上。
已知空间云上一点Point，太阳对其辐射形成的云阴影图像上对应的点是点Shadow，卫星对空间云中点Point投影后对应的点是点cloud。通过云和云影检测以及云和云影匹配处理，可以得到点shadow和点cloud之间的距离。从MOD03数据中读出点shadow和点cloud包含的地理信息，这里选取cloud点的位置信息是3971行、4502列，即()，则由最佳匹配相对位置(2，9)可知其对应的点shadow的坐标为()。根据这两点的位置从MOD03中读取点shadow的纬度lat1＝34.3604°、经度lon1＝-106.2759°，点cloud的纬度lat2＝34.3468°、经度lon2＝-106.2365°，将图像上云和云影像素之间的距离(单位像素)转换为实际地理上的距离Distance(单位公里)，所用计算公式为(7)-(9)
其中a＝lat1-lat2(8)
b＝lon1-lon2
经过公式计算，Distance＝3.9247(km)。由于测试选择平原地区的云和云阴影图像，因此可以忽略云和云阴影的高度差异，认为点shadow和点cloud是在同一平面上。从MOD03数据中读出点shadow的太阳方位角SunAZ＝134.905°＝2.3545rad、太阳天顶角SunZ＝14.954°＝0.261rad、点cloud的卫星方位角SensorAZ＝-77.2409°＝4.93508rad、卫星天顶角SensorZ＝41.2025°＝0.71912rad，将其结合云和云影距离等信息建立云高的空间几何模型，求得空间云中的一点Point相对于地面的高度Height。
为了便于表达云高模型，将图5的模型表示为图6。令点Point为点1，其坐标为(x1，y1，z1)；点Shadow为点2，其坐标为(x2，y2，z2)；点Cloud为点3，其坐标为(x3，y3，z3)。
已知SunZ，SensorZ，因此A1，A2可知，即有公式(10)
SunZ＝A1，SensorZ＝A2
已知SunAZ和Sensor AZ，因此B1，B2可求，即有公式(11)
B1＝π-SunAZ，B2＝SensorAZ-π(11)
令点1到点3的单位向量为
点1’的坐标为(x1’，y1’，z1’)，该向量在x、y、z轴的投影分别为
则有公式(12)-(15)
令点1到点2的向量为
该向量在x、y、z轴的投影分别为
且有公式(16)-(19)
由向量内积公式有公式(20)
的夹角，0≤α≤π
同理令点2到点3的单位向量为
点2’的坐标为(x2’，y2’，z2’)，该向量在x、y、z轴的投影分别为
则有公式(24)-(27)
令点2到点1的向量为
该向量在x、y、z轴的投影分别为
且有公式(28)-(31)
由向量内积公式有
的夹角，0≤β≤π
因此可得公式(36)
γ＝π-α-β(36)
令D13为点1到点3的距离，令D23为点2到点3的距离，令D12为点1到点2的距离，由三角形正弦定理有
且有D12＝‖r12‖
点3的坐标3(x3，y3，z3)可用点1表示为
x3＝x1+D13sinA1cosB1(41)
y3＝y1+D13sinA1sinB1(42)
z3＝z1+D13cosA1 (43)
或者用点2表示为
x3＝x2+D23sinA2cosB2(44)
y3＝y2+D23sinA2sinB2(45)
z3＝z2+D23cosA2 (46)
理论上的云高Height的值计算为公式(47)
Height＝z1+D13cosA1＝z2+D23cosA2(47)
利用上述公式，将已知数据代入，可以求得云高Height＝3.547(km)。
五、云高计算和分析寻找符合云高检测方法测试条件的时间和地区，确定测试地点；通过本发明提出的云高检测方法，将检测云的MODIS空间数据信息与云和云影之间的距离代入云高模型，对该检测云的多个点进行云高计算，并将计算结果求均值，作为检测云的云层高度，采用一至四的步骤，对周围9点进行云高计算，最后求得的云的云层高度的均值是3.554km。
本发明利用卫星遥感数据中的云影信息检测云高的方法，该方法的应用条件是分辨率为250m的MODIS卫星遥感数据中的平原区域的数据。
本发明以MODIS遥感图像处理技术为基础，提出平原区域云高测量的新方法。其应用条件是分辨率为250m的MODIS数据卫星遥感数据中的平原区域，使用的MODIS数据是经过常规纠正处理后的数据。例如选取250m的MODIS数据卫星遥感数据MOD02QKM.A0.005.5.hdf。
一、插值计算由于需要云和云阴影点的具体地理信息数据，包括经度、纬度、太阳天顶角、太阳方位角、卫星天顶角、卫星方位角数据，因此需要MODIS数据星历表MOD03的数据，即MOD03.A0.005.5.hdf。针对250m分辨率的MOD03数据中的某些数据量重叠的情况，需要从1000m分辨率的数据中重新计算获得250m分辨率对应的MOD03数据，才能获得云和云阴影上点的更加精确的信息。因此这里采用通过采用线性插值公式解决该问题。线性内插是假设在二个已知数据中的变化为线性关系，因此可由已知二点a和c，去计算内插点b的值，公式如下
其中a＜b＜c在上式的b点即是代表要内插的点，f(b)则是要计算的内插函数值。
具体实现步骤如下
(1)读取1000m MOD03中的地理信息数据，如经度。
(2)1000m MOD03数据格式是2030行，1354列；插值后的250m MOD03数据格式是8120行，5416列，其行和列大小分别为1000m数据的4倍。采用将1000mMOD03数据中每10行×1354列的数据通过公式(1)插值成为40行×5416列的数据。由于1000m MOD03数据有2030行，因此需要进行203次的插值运算。
(3)利用Matlab软件中的INTERP2函数可以实现线性内插，通过计算得到250mMOD03中的地理信息数据。
二、云和云阴影检测从原图像中手动选取具有云和云阴影的图像，选取图像在原图中的位置信息是起始点在1930行、2985列，即()，范围是30×30大小的图像范围，选取的波段是可见光波段2。具体选取图像如图7所示。该步骤主要完成从选取图像中检测云和云阴影图像的具体位置，选取数据是MOD02QKM。
目前，绝大多数的遥感图像云特征建模都是从图像的光谱特征出发，考虑云与其它地物在可见光、近红外和热红外波段上的差异。云在可见光和红外波段与植被、土壤以及水域等下垫面介质的反射率和辐射亮温值存在差异，云具有较高的反射率而具有低的亮温值等特点进行云与地物的区分。在MODIS数据的可见光波段0.645μm处，晴空数据地物一般具有较低的反射率，反射值在0～0.3之间，云则有高的反射率，一般高于0.3，是进行云检测处理的首选波段。1.375μm处于强的水汽吸收波段，地面的反射因水汽吸收难以到达传感器，导致反射率低，而高云上方的水汽对云的反射值削弱较少，因此该波段对高层卷云敏感，具有高的反射率，低云则具有低的反射率，通常低于0.02。此波段可以用来进行卷云的检测。8.6μm和11μm波段对云特性敏感。短波中有强的水汽吸收作用，导致晴空中8.6μm处的辐射值低于11μm处的辐射值；冰晶在8.6μm的吸收远小于11μm处的吸收，所以，冰晶组成的卷云在8.6μm处的辐射值高于11μm处的辐射值；而水滴组成的光学厚度较大的云则在8.6μm处的辐射值低于11μm处的辐射值。
采用多光谱云检测主要是利用以上各个波段对云的敏感特性，以及MODIS数据的各通道的波谱特性实现云与地物的区分。在多光谱云检测结果的基础上，对于可见光波段的图像，结合常用的图像检测方法如区域标定法，区域生长法，分水岭法，阈值分割法等将云与非云区域进行更有效的区分。对于云阴影的检测采取同云检测一样，首先利用云影在各个波段的敏感特性，主要利用1.2、0.94、0.87、0.66μm这四个波段反射率来检测云阴影；另外，利用3.9μm和4.0μm波段的亮度温度差大于零也可以完成云阴影的检测。
这里采用双阈值分割法在可见光波段2中实现云和云影图像的检测。图像阈值分割是一种广泛使用的图像分割技术，它利用了图像中要提取的目标与其背景在灰度特性上的差异，把图像视为具有不同灰度级的两类区域(目标和背景)的组合，选取一个合适的阈值，以确定图像中的每一个像素点应该属于目标还是背景区域，从而产生相应的二值图像。采用的双阈值分割的数学表达式为(2)式所示
通过对云和云影图像的光谱分析，选取云图像的阈值是T2＝0.382μm，阴影图像的T1＝0.136μm，f(i，j)是图像中像素的值。然后采用区域标定法，在双阈值检测结果的基础上，将云图像的位置和云阴影图像提取出来。
区域标定是给连接在一起的像素(连接成分)附上相同的标记，不同的连接成分附上不同的标记。区域标定在二值图像处理中占有非常重要的地位。通过该处理将各个连接成分区分开，然后就可以调查各个连接成分的形状特征。
区域标定也有许多方法，本发明采用了一种简单的方法。步骤如下
(1)扫描图像，遇到没加标记的目标像素(白像素)P时，添加一个新的标记(LABEL)。
(2)将与P连接在一起(即相同连接成分)的像素添加相同的标记。
(3)进一步给所有与加标记像素连接在一起的像素添加相同的标记。直到连接在一起的像素全部被添加标记。这样，一个连接成分就被添加了相同的标记。
(4)回到第1步，重新查找新的没加标记的像素，重复上述各个步骤。图像全部被扫描后，处理结束。
采用区域标定后，会检测出很多候选目标。通过设定面积参数，使得云区域和云阴影区域能够从候选目标中提取出来，同时去掉图像中噪声对检测结果的影响。经过云和云阴影检测后的图像如图8所示。其中灰色的像素代表背景，黑色的像素代表云阴影图像，白色的像素代表云图像。从图像可见，该方法可以较好地从原图像中检测出云和云影图像。
三、MODIS数据的云和云影匹配
根据云和云影检测的结果，进行云和云影的图像匹配。可以采用数字图像处理学中的图像匹配方法进行云和云影匹配，通过匹配得到云和云影之间的距离。
在云影图像中，通过云影检测，云和阴影区域已确定。每个云都有对应的影子图像，而且通常云阴影的图像要小于云的图像。基于模板匹配的原理，求得云与影间的距离的匹配方法为以云阴影图像为模板，搜索云图像中与阴影图像匹配的位置。图像匹配中，通过采用自相关算法来计算模板与搜索区域的相关程度，将具有最大相关度的位置定为匹配最佳位置。如S为云阴影图像，C为云图像，M×N是云阴影图像的大小，则相关系数Coef(coefficient自定义的参数)计算公式如(3)所示，(i，j)为云和云阴影图像的匹配相对位置矢量，i为水平相对位置矢量，j为垂直相对位置矢量。
其中各变量为(4)式和(5)式
通过公式(6)得到Coef最大时对应的最佳匹配相对位置矢量(i，j)
将所有检测到的云图像的像素信息代入C，所有检测到的云阴影图像的像素信息代入S，根据公式(3)可以得到多个相关系数Coef。利用公式(6)选取最大的相关系数Coef，并据此得到最佳匹配相对位置矢量(i，j)。通过计算，最大的相关系数是0.96，对应的(i，j)为(4，8)，即云和云影图像的对应匹配点水平方向相差8个像素，垂直方向相差4个像素。匹配效果如图9所示。
四、建立云高模型求得云层高度
根据云与云影之间的距离和卫星检测到的数据，通过云高模型的计算得到云层高度。图4为云高模型建立的理论基础图。由于空间中存在云，所以太阳的光线会被云完全或者部分遮挡，因此在地面上形成云的阴影图像。云阴影图像的位置与太阳相对于空间中云的位置有关，太阳-真云-云阴影共线；而当卫星观测地球时，空间中的云在MODIS探测器成像的位置与卫星的拍摄角度有关，即卫星-真云-投影云共线；因此MODIS数据中，云阴影图像和投影云图像都与空间中的云有着密切关系。可以根据云阴影图像、投影云图像和空间中的真云三者之间的空间位置关系建立模型，求得空间中的真云相对于地面的高度。
图5是云高模型建立的几何模型由于测试选择平原区域的云高测量，投影云和对应云阴影在同一平原区域，所以可以忽略地球是一个椭圆球体，并近似定义云和云阴影的地理信息在同一个平面上。
已知空间云上一点Point，太阳对其辐射形成的云阴影图像上对应的点是点Shadow，卫星对空间云中点Point投影后对应的点是点cloud。通过云和云影检测以及云和云影匹配处理，可以得到点shadow和点cloud之间的距离。从MOD03数据中读出点shadow和点cloud包含的地理信息，这里选取cloud点的位置信息是1942行、2998列，即()，则由最佳匹配相对位置(4，8)可知其对应的点shadow的坐标为()。根据这两点的位置从MOD03中读取点shadow的纬度lat1＝39.7052°、经度lon1＝-109.8255°，点cloud的纬度lat2＝39.6927°、经度lon2＝-109.8047°，将图像上云和云影像素之间的距离(单位像素)转换为实际地理上的距离Distance(单位公里)，所用计算公式为(7)-(9)
其中a＝lat1-lat2(8)
b＝lon1-lon2
经过公式计算，Distance＝2.2606(km)。由于测试选择平原地区的云和云阴影图像，因此可以忽略云和云阴影的高度差异，认为点shadow和点cloud是在同一平面上。从MOD03数据中读出点shadow的太阳方位角SunAZ＝137.1709°＝2.3941rad、太阳天顶角SunZ＝20.9475°＝0.3656rad、点cloud的卫星方位角SensorAZ＝-80.4722°＝4.8787rad、卫星天顶角SensorZ＝6.4925°＝0.1133rad，将其结合云和云影距离等信息建立云高的空间几何模型，求得空间云中的一点Point相对于地面的高度Height。
为了便于表达云高模型，将图5的模型表示为图6。令点Point为点1，其坐标为(x1，y1，z1)；点Shadow为点2，其坐标为(x2，y2，z2)；点Cloud为点3，其坐标为(x3，y3，z3)。
已知SunZ，SensorZ，因此A1，，A2可知，即有公式(10)
SunZ＝A1，SensorZ＝A2(10)
已知SunAZ和Sensor AZ，因此B1，B2可求，即有公式(11)
B1＝π-SunAZ，B2＝SensorAZ-π(11)
令点1到点3的单位向量为
点1’的坐标为(x1’，y1’，z1’)，该向量在x、y、z轴的投影分别为
则有公式(12)-(15)
令点1到点2的向量为
该向量在x、y、z轴的投影分别为
且有公式(16)-(19)
由向量内积公式有公式(20)
的夹角，0≤α≤π
同理令点2到点3的单位向量为
点2’的坐标为(x2’，y2’，z2’)，该向量在x、y、z轴的投影分别为
则有公式(24)-(27)
令点2到点1的向量为
该向量在x、y、z轴的投影分别为
且有公式(28)-(31)
由向量内积公式有
的夹角，0≤β≤π
因此可得公式(36)
γ＝π-α-β(36)
令D13为点1到点3的距离，令D23为点2到点3的距离，令D12为点1到点2的距离，由三角形正弦定理有
且有D12＝‖r12‖(38)
点3的坐标3(x3，y3，z3)可用点1表示为
x3＝x1+D13sinA1cosB1(41)
y3＝y1+D13sinA1sinB1(42)
z3＝z1+D13cosA1 (43)
或者用点2表示为
x3＝x2+D23sinA2cosB2(44)
y3＝y2+D23sinA2sinB2(45)
z3＝z2+D23cosA2 (46)
理论上的云高Height的值计算为公式(47)
Height＝z1+D13cosA1＝z2+D23cosA2(47)
利用上述公式，将已知数据代入，可以求得云高Height＝4.7181(km)。
五、云高计算和分析寻找符合云高检测方法测试条件的时间和地区，确定测试地点；通过本发明提出的云高检测方法，将检测云的MODIS空间数据信息与云和云影之间的距离代入云高模型，对该检测云的多个点进行云高计算，并将计算结果求均值，作为检测云的云层高度，采用一至四的步骤，对周围9点进行云高计算，最后求得的云的云层高度的均值是4.6736km。
1、一种利用卫星遥感数据中的云影信息检测云高的方法，该方法的应用条件是分辨率为250m的MODIS卫星遥感数据中的平原区域的数据，该方法包括以下四个部分
一、插值计算线性内插是由已知二点a和c，去计算内插点b的值，公式如下
其中a＜b＜c，上式的b点是代表要内插的点，f(b)则是要计算的内插函数值；
具体实现步骤如下
1)读取1000m MOD03中的地理信息数据；
2)1000m MOD03数据格式是2030行，1354列；插值后的250m MOD03数据格式是8120行，5416列，其行和列大小分别为1000m数据的4倍；采用将1000m MOD03数据中每10行×1354列的数据通过公式(1)插值成为40行x5416列的数据；由于1000m MOD03数据有2030行，因此需要进行203次的插值运算；
3)利用Matlab软件中的INTERP2函数实现线性内插，通过计算得到250m MOD03中的地理信息数据；
二、云和云阴影检测
从气象卫星提供的数据中选取具有云和云阴影的图像；
采用双阈值分割法在可见光波段2中实现云和云影图像的检测，采用的双阈值分割的数学表达式为(2)式所示
通过对云和云影图像的光谱分析，选取云图像的阈值是T2，阴影图像的T1，f(i，j)是图像中像素的值；然后采用区域标定法，在双阈值检测结果的基础上，将云图像的位置和云阴影图像提取出来；
区域标定法的步骤如下
(a)扫描图像，遇到白像素P时，添加一个新的标记LABEL；
(b)将与白像素P连接在一起，即相同连接成分的像素添加相同的标记；
(c)进一步给所有与加标记像素连接在一起的像素添加相同的标记，直到连接在一起的像素全部被添加标记；这样，一个连接成分就被添加了相同的标记；
(d)回到第(a)步，重新查找新的白像素，重复上述各个步骤；图像全部被扫描后，处理结束；
采用区域标定后，对检测出的候选目标通过设定面积参数，使得云区域和云阴影区域能够从候选目标中提取出来，同时去掉图像中噪声对检测结果的影响；
三、MODIS数据的云和云影匹配
根据云和云影检测的结果，采用数字图像处理学中的图像匹配方法进行云和云影匹配，通过匹配得到云和云影之间的距离；
具体步骤如下
设S为云阴影图像，C为云图像，MxN是云阴影图像的大小，则相关系数Coef计算公式如(3)所示，(i，j)为云和云阴影图像的匹配相对位置矢量，i为水平相对位置矢量，j为垂直相对位置矢量；
其中各变量为(4)式和(5)式
通过公式(6)得到Coef最大时对应的最佳匹配相对位置矢量(i，j)
根据上述公式进行云和云阴影的匹配，将所有检测到的云图像的像素信息代入公式(5)的C，所有检测到的云阴影图像的像素信息代入S，根据公式(3)可以得到多个相关系数Coef；利用公式(6)选取最大的相关系数Coef，并据此得到最佳匹配相对位置矢量(i，j)；
四、建立云高模型求得云层高度
根据云与云影之间的距离和卫星检测到的数据，通过云高模型的计算得到云层高度；根据云阴影图像、投影云图像和空间中的真云三者之间的空间位置关系建立模型，求得空间中的真云相对于地面的高度；
已知空间云上一点Point，太阳对其辐射形成的云阴影图像上对应的点是点Shadow，卫星对空间云中点Point投影后对应的点是点cloud；通过云和云影检测以及云和云影匹配处理，可以得到点shadow和点cloud之间的距离；从MOD03数据中读出点shadow和点cloud包含的地理信息点shadow的纬度lat1、经度lon1，点cloud的纬度lat2、经度lon2，将图像上云和云影像素之间的单位为像素的距离，转换为实际地理上的单位为公里的距离Distance，所用计算公式为(7)-(9)
其中a＝lat1-lat2 (8)
b＝lon1-lon2
从MOD03数据中读出点shadow的太阳方位角SunAZ、太阳天顶角SunZ、点cloud的卫星方位角SensorAZ、卫星天顶角SensorZ，将其结合云和云影距离等信息建立云高的空间几何模型，求得空间云中的一点Point相对于地面的高度Height；
令点Point为点1，其坐标为(x1，y1，z1)；点Shadow为点2，其坐标为(x2，y2，z2)；点Cloud为点3，其坐标为(x3，y3，z3)；
已知SunZ，SensorZ，因此A1，A2可知，即有公式(10)
SunZ＝A1，SensorZ＝A2
已知SunAZ和Sensor AZ，因此B1，B2可求，即有公式(11)
B1＝π-SunAZ，B2＝SensorAZ-π
令点1到点3的单位向量为
点1’的坐标为(x1’，y1’，z1’)，该向量在x、y、z轴的投影分别为
则有公式(12)-(15)
令点1到点2的向量为
该向量在x、y、z轴的投影分别为
且有公式(16)-(19)
由向量内积公式有公式(20)
的夹角，0≤α≤π，则有
同理令点2到点3的单位向量为
点2’的坐标为(x2’，y2’，z2’)，该向量在x、y、z轴的投影分别为
则有公式(24)-(27)
令点2到点1的向量为
该向量在x、y、z轴的投影分别为
且有公式(28)-(31)
由向量内积公式有
的夹角，0≤β≤π，则有
因此可得公式(36)
γ＝π-α-β(36)
令D13为点1到点3的距离，令D23为点2到点3的距离，令D12为点1到点2的距离，由三角形正弦定理有
且有D12＝‖r12‖
已知点3的坐标3(x3，y3，z3)可用点1表示为
x3＝x1+D13sinA1cosB1
y3＝y1+D13sinA1sinB1
z3＝z1+D13cosA1
或者用点2表示为
x3＝x2+D23sinA2cosB2
y3＝y2+D23sinA2sinB2
z3＝z2+D23cosA2
云高Height的值计算为公式(47)
Height＝z1+D13cosA1＝z2+D23cosA2(47)
利用上述公式，将已知数据代入，求得云高。
2、根据权利要求1所述的利用卫星遥感数据中的云影信息检测云高的方法，其特征在于该方法包括云高计算和分析，将检测云的MODIS空间数据信息与云和云影之间的距离代入云高模型，对该检测云的多个点进行云高计算，并将计算结果求均值。
本发明涉及一种利用卫星遥感数据中的云影信息检测云高的方法，属于大气探测领域。包括插值计算、云和云阴影检测、MODIS数据的云和云影匹配、建立云高模型求得云层高度和云高计算和分析五部分。本方法对遥感图像数据中云和云影的信息进行检测；结合数字图像技术进行云和云影匹配，得到云和云影的之间的距离；最后，建立云高模型求得云层高度。本方法可以同一时刻对一片云多个点进行高度计算，以此求得的高度值将极大的提高云高的检测精度。通过本方法可以对符合测量条件的任何地区进行云高测量，而且不用在当地设置测量仪器。
文档编号G01S17/95GKSQ
公开日日 申请日期日 优先权日日
发明者顾玲嘉, 任瑞治, 郭树旭, 爽 张, 王昊丰, 单江东 申请人:吉林大学