专利名称：环境一号hj-1a/b星ccd数据的云检测方法
技术领域：
本发明涉及卫星遥感数据云检测领域技术领域，涉及一种环境一号卫星CXD数据的云检测方法，尤其涉及一种基于地表反射率数据库的环境一号HJ-1A/B星CCD数据的云检测方法。
背景技术：
云在时间和空间分布上有极大的不确定性，其高度、厚度、组成种类以及太阳高度角、方位角的变化导致其卫星云图特征千变万化。在利用航天遥感手段获取地球空间信息的过程中，云是光信号传播的严重障碍，在很大程度上影响遥感信息获取的质量，从而降低了数据的利用率。云遮挡问题是遥感数据普遍存在的现象，也是卫星遥感数据定量化应用中的一大难题。准确的分离云和晴空是反演大气和地表各种参数的必要前提。因此采用实 时有效的云检测方法，剔除遥感影像上云覆盖区域，是缓解遥感图像海量数据对传输通道的压力，增加遥感数据利用率的一个重要途径，也是遥感定量化应用前，影像预处理中的一个重要环节。目前，关于遥感影像云检测的技术很多，主要分为3类阈值法、聚类分析法和人工神经网络法。阈值法是这三种方法当中相对较成熟且易于实现的方法。其主要思想是对需要分析的像元利用不同通道组合的亮温、亮温差、反射率与所设置的阈值进行比对，从而判断出像元是否为云污染像元。其中阈值的获取有基于影像自身直方图分析的方法，有采用滑动窗和嵌套窗的方法，也有对不同时相的影像资料进行统计分析，最终确定阈值的。目前阈值法主要有 ISCCP (The International Satellite Cloud Climatology Project)方法、APP0LL0 (The AVHRR Processing Scheme Over Cloud Land And Ocean)方法、CLAVR(The NOAA Cloud Advanced Very High Resolution Radiometer)方法等。ISCCP 算法假定观测辐射只来自晴空和云两种大气状况，与这两种大气状况相联系的辐射值其变化范围并不相互重叠。检测方法中使用了 0. 6 μ m可见光波段和Ilym红外窗区波段的资料，将像元和晴空的观测辐射值进行比较，当像元的辐射值明显有别于晴空像元时，即某个像元与晴空辐射值的差值明显大于晴空辐射值本身的变化时，该像元被判定为云点。这种方法使云检测的误差最小。但算法依赖于阈值，当像元只有部分云覆盖时则会发生误判现象。APP0LL0算法使用AVHRR五个全分辨率探测通道资料，对有云像元的判定须满足如下条件①通道2、1的比值范围为
像元的反射率高于所设定的阈值或温度低于所设定的阈值；③通道4、5其亮温差比所设定的阈值大对位于海洋上的像元，还要判断其空间均一性是否比所设定的阈值大。CLAVR算法同APP0LL0方法一样是采用了一系列判识阈值，不同的是，其判识单位为2*2的像元矩阵，算法根据下垫面性质和观测时间的不同，可分为四类白天海洋、白天陆地、夜间海洋和夜间陆地。聚类分析法在20世纪70年代被提出，80年代起该方法被用于研究云分类问题。通过将影像分割成为一系列像元阵列，采用聚类分析法获取影像中各像元所代表的云或地物类型。当像元阵内各种地类观测值存在明显差别时，检测结果可能较为合理。但由于该方法的结果较依赖于卫星观测角内观测对象的值域特征，从而在一定程度上影响了该方法的使用。聚类方法主要有直方图聚类、自适应阈值聚类、动态阈值聚类、动态聚类等。人工神经网络法通过对多个样本的学习，获取样本的知识并将其分布存贮于网络中，从而达到对云的理解和识别。该方法是模拟人脑神经系统的工作原理，具有很强的自组织自适应学习能力，但其结果对于训练资料的选取非常敏感。
随着环境一号卫星的发射成功，其卫星影像凭借获取方便、监测区域大、分辨率较高的优点，已在生态破坏、环境污染和灾后评估等领域得到广泛应用。然而，由于其波段设置较少、空间分辨率比较高等特点，已有的云检测技术都在不同程度上存在一定的局限性和不稳定性。而对于同一地区，不同年份同时期的地表反射率变化较�。缃乇矸瓷渎适菀氲絑 检测算法中，不仅有利于提闻Z 检测的效率，还可提闻Z 检测的精度。
发明内容
(一)要解决的技术问题本发明要解决的技术问题是，针对现有技术的不足，提供一种环境一号HJ-1A/B星CXD数据的云检测方法，可以稳定的，且快速、方便地对环境一号HJ-1A/B星CXD遥感影像进行云识别的检测方法。(二)技术方案本发明提供一种环境一号HJ-1A/B星CXD数据的云检测方法，包括如下步骤S1 利用M0D09Q1数据建立地表反射率数据库；S2 :对环境一号HJ-1A/B星的CXD数据进行辐射定标的预处理，获取CCD数据各波段数据的表观反射率；S3 :将S2中环境一号HJ-1A/B星的CCD数据各波段数据的表观反射率数据与从SI中地表反射率数据库中选取成像日期内的地表反射率数据进行比较进行云像元识别，并生成第一层云掩膜；S4 :利用环境一号HJ-1A/B星的CCD数据的第3波段的表观反射率数据值域范围进行云像元识别，并生成第二层云掩膜；S5 :利用环境一号HJ-1A/B星的CCD数据的第4波段与第3波段表观反射率数据比值的值域范围进行云像元识别，并生成第三层云掩膜；S6 :根据步骤S3-S5的云像元识别结果，得到环境一号HJ-1A/B星的CXD数据云掩膜。更好地，SI进一步包括如下步骤S1. I :获取一段时间内的M0D09Q1数据；S1. 2 将所述M0D09Q1数据中同时期的多个M0D09Q1影像的第I波段地表反射率数据和第2波段地表反射率数据分别进行拼接，形成一个第I波段地表反射率数据和一个第2波地表反射率段数据，并设置数据输出格式为GeoTiff，一个所述M0D09Q1影像包含多个第I波段地表反射率数据和多个第2波段地表反射率数据；S1. 3 :重复SI. 2，获取在SI. 2中设定的一段时间内的第I波段和第2波段地表反射率数据，并将所述第I波段和所述第2波段地表反射率数据分别组成两个数据库，所述数据库中存储的数据均为GeoTiff格式数据。更好地，S2进一步包括如下步骤S2. I :从环境一号HJ-1A/B星的C⑶数据的辅助数据文件中获取定标系数，并将所述的CCD数据的第3波段和4波段的像元亮度DN值转换为表观辐亮度；S2. 2 :从环境一号HJ-1A/B星的CCD数据的辅助数据文件中获取太阳天顶角，将S2. I中的表观辐亮度归一化为第3波段和4波段的表观反射率。更好地，S3进一步包括如下步骤S3. I :从环境一号HJ-1A/B星的C⑶数据的辅助数据文件中获取成像日期；S3. 2 :根据S3. I中的成像日期，从SI中的地表反射率数据库中选择最接近C⑶数据成像日期的地表反射率数据；S3. 3 :将环境一号HJ-1A/B星CXD数据的第3波段和4波段表观反射率数据分别与S3. 2中成像日期内的地表反射率数据的第I波段和2波段地表反射率数据做差，若CXD数据中像元所对应的这两个差值至少有一个大于O. 2，则该像元被判定为云像元，赋值为1，其余无云像元赋值为O，生成第一层云掩膜。更好地，在S4中，若环境一号HJ-1A/B星CXD数据中像元的第3波段表观反射率值大于O. 22，则该像元被判定为云像元，赋值为I，其余无云像元赋值为0，生成第二层云掩膜。更好地，在S5中，若环境一号HJ-1A/B星CXD数据中像元的第4波段表观反射率与第3波段表观反射率的比值大于O. 87，且小于I. 15，则该像元被判定为云像元，赋值为I，其余无云像元赋值为0，生成第三层云掩膜。 更好地，S6进一步包括如下步骤S6. I :将S3-S5中所述第一、二、三层云掩膜进行相加求和，生成新的云掩膜；S6. 2 :将所述新的云掩膜中，值为I的像元赋值为0，表示无云像元；值为2的像元保留原值不变，表示可能为云的像元；值为3的像元赋值为1，表示云像
J Li ο更好地，所述环境一号HJ-1A/B星C⑶数据为经过几何校正的第1、2、3、4四个波段的GeoTiff格式数据。更好地，所述M0D09Q1数据为M0D09L3级250m的8日合成格网数据，包含MODIS传感器的第I和第2波段地表反射率及经纬度信息，数据格式为HDF。更好地，所述环境一号卫星CXD数据的辅助数据文件为xml格式文件，其包括定标系数、成像时间辅助信息和观测几何数据。(三)有益效果本发明的方法利用M0D09Q1数据构建长时间序列、大空间区域的地表反射率库，通过选择合适的地表反射率数据与环境一号HJ-1A/B星的CCD数据进行比较，再结合CCD数据本身特征波段阈值的选�。芄豢焖佟⒎奖�、稳定地实现云检测，得到云掩膜，为进一步反演大气和地表各种参数提供必要的前提。
图I为本发明环境一号HJ-1A/B星CXD数据的云检测方法流程图。
具体实施例方式下面结合附图和实施例，对本发明的具体实施方式
作进一步详细描述。以下实施例用于说明本发明，但不用来限制本发明的范围。在本发明中的卫星选用环境一号HJ-1A/B星(也称环境减灾卫星)全称中国环境与灾害监测预报小卫星星座，是中国专门用于环境与灾害监测预报的卫星。主要用于对生态环境和灾害进行大范围、全天候动态监测，及时反映生态环境和灾害发生、发展过程，对生态环境和灾害发展变化趋势进行预测，对灾情进行快速评估，为紧急求援、灾后救助和重建工作提供科学依据。其搭载的C⑶相机像元分辨率为30m，包括3个可见光波段——第 I 波段(O. 43-0. 52μπι)、第 2 波段(O. 52-0. 60μπι)、第 3 波段(O. 63-0. 69 μ m)以及一个近红外波段——第4波段(O. 76-0. 90 μ m)。
MODIS的全称为中分辨率成像光谱仪,是美国EOS (Earth Observation System)计划系列卫星上的最主要的传感器。M0D09是陆地2级标准数据数据，内容为表面反射，空间分辨率250m，白天每日数据。GeoTiff是包含地理信息的一种Tiff格式文件。M0D098日合成的地表反射率数据中，选取8天中观测角较低，且经过云和气溶胶筛选过后，可信度最高的地表反射率值作为8日合成数据的像元值。本发明的方法主要基于250m空间分辨率的M0D09地表反射率数据，建立地表反射率库，通过比较环境一号HJ-1A/B星的CXD表观反射率数据与相应的M0D09地表反射率数据的差异，再结合CXD数据本身特征波段阈值的选�。詈笫迪忠缓盼佬荂CD数据的云检测。从而，本发明能够快速、方便、稳定地实现对环境一号HJ-1A/B星CXD遥感影像进行云检测。 图I为本发明环境一号HJ-1A/B星CXD数据的云检测方法流程图，如图I所示，环境一号HJ-1A/B星CCD数据的云检测方法包括如下步骤SI :从MODIS网站获取M0D09Q1数据，并利用M0D09Q1构建一段时间内(例如从2009年I月至2011年12月共三年)的空间分辨率为250m的地表反射率数据库，数据库中存储的数据均为GeoTifT格式数据；S2 :对环境一号HJ-1A/B星的C⑶数据进行辐射定标的预处理，以获取CXD数据各波段数据的表观反射率；S3 :将S2中所述环境一号HJ-1A/B星的CXD数据的表观反射率数据与从所述S I中地表反射率数据库中选取成像日期地表反射率数据进行比较进行云像元识别，并生成第
一层云掩膜；S4 :利用环境一号HJ-1A/B星的CXD数据的第3波段的表观反射率数据值域范围，进行云像元识别，若环境一号HJ-1A/B星的CXD数据中像元的第3波段表观反射率值大于
O.22，则该像元被判定为云像元，赋值为1，其余无云像元赋值为0，生成第二层云掩膜；S5 :利用环境一号HJ-1A/B星的CXD数据第4波段与第3波段表观反射率数据比值的值域范围，进行云像元识别，若环境一号HJ-1A/B星的CXD数据中像元的第4波段表观反射率与第3波段表观反射率的比值大于O. 87，且小于I. 15，则该像元被判定为云像元，赋值为1，其余无云像元赋值为0，生成第三层云掩膜；S6 :根据S3-S5的云像元识别结果，得到最终的环境一号HJ-1A/B星的CXD数据的
云掩膜。其中，SI进一步包括如下步骤SI. I :从MODIS网站获取一段时间内(例如从2009年I月至2011年12月共三年)的M0D09Q1数据(即M0D09L3级250m的8日合成格网数据)；SI. 2 :利用 MODIS 重投影工具(MRT，MODIS Reprojection Tool)将 M0D09Q1 数据中每8天的多个M0D09Q1影像的第I波段和第2波段数据分别进行拼接，行成一个第I波段数据和一个第2波段数据，并重新投影成与环境一号HJ-1A/B星的CCD数据相同的投影和坐标系，并设定其输出格式为GeoTiff 个所述M0D09Q1影像包含多个第I波段地表反射率数据和多个第2波段地表反射率数据SI. 3 :重复SI. 2，获取在SI. I中设定的一段时间内的第I波段和第2波段地表反射率数据，并将所述第I波段和所述第2波段地表反射率数据分别组成两个数据库，所述数据库中存储的数据均为GeoTifT格式数据。其中，S2进一步包括如下步骤S2. I :从环境一号HJ-1A/B星的C⑶数据的辅助数据文件中获取定标系数，并利用IDL计算机可视化语言编程实现将其第3波段和第4波段的像元亮度DN值转换为表观辐亮度L = DN/a+b (I)其中，L为表观福売度，a和b为定标系数；S2. 2 :从环境一号HJ-1A/B星的CXD数据的辅助数据文件中的观测几何数据中获取太阳天顶角，并利用IDL语言编程实现将S2. I中的表观辐亮度归一化为第3波段和第4 波段的表观反射率
TtL(2)P =---
h} COS^其中，P为表观反射率，Ea为大气层顶太阳辐照度，93为太阳天顶角。其中，S3进一步包括如下步骤S3. I :从环境一号HJ-1A/B星的CXD数据的辅助数据中的成像时间辅助信息中获取遥感影像的成像日期；S3. 2 :根据S3. I中成像日期，从SI中地表反射率数据库中选择最接近(XD数据成像日期的M0D09地表反射率数据，若所述该CCD数据所覆盖范围对应在所选取的M0D09地表反射数据中有超过10%的值为无效值，则开始搜索CCD数据成像时间前后15天内的其他地表反射率数据，若仍未选取出合适的地表反射率数据，则继续搜索不同年份同时期的其他地表反射率数据；S3. 3 :利用IDL语言编程，将环境一号HJ-1A/B星的CXD数据的第3波段和第4波段表观反射率数据，分别与地理位置匹配的、所述成像日期内的地表反射率数据的第I波段和第2波段地表反射率数据做差，若CXD数据中像元所对应的这两个差值至少有一个大于O. 2，则该像元被判定为云像元，赋值为1，其余无云像元赋值为0，生成第一层云掩膜。其中，S6进一步包括如下步骤S6. I :利用IDL语言编程，将步骤S3-S5中生成的第一、二、三层云掩膜进行相加求和，生成新的云掩膜；S6. 2 :利用IDL语言编程，将生成的新云掩膜中，值为I的像元赋值为0，表示无云像元；值为2的像元保留原值不变，表示可能为云的像元；值为3的像元赋值为1，表示云像元，得到最终的云掩膜，数据格式为GeoTiff。上述环境一号HJ-1A/B星的CXD数据为经过几何校正的第1、2、3、4四个波段的GeoTiff格式数据；所述环境一号HJ-1A/B星的CXD数据的辅助数据文件为xml格式文件，其为包括定标系数、成像时间辅助信息和观测几何数据；所述M0D09L3级250m的8日合成格网数据为M0D09Q1数据，包含MODIS传感器的第I和第2波段地表反射率及经纬度信息，数据格式为HDF。以上实施方式仅用于说明本发明，而并非对本发明的限制，有关技术领域的普通技术人员，在不脱离本发明的精神和范围的情况下，还可以做出各种变化和变型，因此所有等同的技术方案也属于本发明的范畴，本发明的专利保护范围应由权利要 求限定。
权利要求
1.一种环境一号HJ-1A/B星CXD数据的云检测方法，其特征在于，包括如下步骤 51:利用M0D09Q1数据建立地表反射率数据库； 52:对环境一号HJ-1A/B星的CCD数据进行辐射定标的预处理，获取CCD数据各波段数据的表观反射率； 53:将S2中环境一号HJ-1A/B星的CCD数据各波段数据的表观反射率数据与从SI中地表反射率数据库中选取成像日期内的地表反射率数据进行比较进行云像元识别，并生成第一层云掩膜； 54:利用环境一号HJ-1A/B星的CCD数据的第3波段的表观反射率数据值域范围进行云像元识别，并生成第二层云掩膜； 55:利用环境一号HJ-1A/B星的CCD数据的第4波段与第3波段表观反射率数据比值的值域范围进行云像元识别，并生成第三层云掩膜； 56:根据步骤S3-S5的云像元识别结果，得到环境一号HJ-1A/B星的CXD数据云掩膜。
2.如权利要求I所述的环境一号HJ-1A/B星CXD数据的云检测方法，其特征在于，SI进一步包括如下步骤 SI. I :获取一段时间内的M0D09Q1数据； SI. 2 :将所述M0D09Q1数据中同时期的多个M0D09Q1影像的第I波段地表反射率数据和第2波段地表反射率数据分别进行拼接，形成一个第I波段地表反射率数据和一个第2波地表反射率段数据，并设置数据输出格式为GeoTiff，一个所述M0D09Q1影像包含多个第I波段地表反射率数据和多个第2波段地表反射率数据； 51.3 :重复SI. 2，获取在SI. 2中设定的一段时间内的第I波段和第2波段地表反射率数据，并将所述第I波段和所述第2波段地表反射率数据分别组成两个数据库，所述数据库中存储的数据均为GeoTiff格式数据。
3.如权利要求I所述的环境一号HJ-1A/B星C⑶数据的云检测方法，其特征在于，S2进一步包括如下步骤 52.I :从环境一号HJ-1A/B星的CCD数据的辅助数据文件中获取定标系数，并将所述的CCD数据的第3波段和4波段的像元亮度DN值转换为表观辐亮度； 52.2 :从环境一号HJ-1A/B星的CXD数据的辅助数据文件中获取太阳天顶角，将S2. I中的表观辐亮度归一化为第3波段和4波段的表观反射率。
4.如权利要求I所述的环境一号HJ-1A/B星C⑶数据的云检测方法，其特征在于，S3进一步包括如下步骤 53.I :从环境一号HJ-1A/B星的CXD数据的辅助数据文件中获取成像日期； S3. 2 :根据S3. I中的成像日期，从SI中的地表反射率数据库中选择最接近(XD数据成像日期的地表反射率数据； S3. 3 :将环境一号HJ-1A/B星CXD数据的第3波段和4波段表观反射率数据分别与S3. 2中成像日期内的地表反射率数据的第I波段和2波段地表反射率数据做差，若CXD数据中像元所对应的这两个差值至少有一个大于O. 2，则该像元被判定为云像元，赋值为1，其余无云像元赋值为0，生成第一层云掩膜。
5.如权利要求I所述的环境一号HJ-1A/B星CXD数据的云检测方法，其特征在于，在S4中，若环境一号HJ-1A/B星(XD数据中像兀的第3波段表观反射率值大于O. 22,则该像元被判定为云像元，赋值为1，其余无云像元赋值为O，生成第二层云掩膜。
6.如权利要求I所述的环境一号HJ-1A/B星CXD数据的云检测方法，其特征在于，在S5中，若环境一号HJ-1A/B星CXD数据中像元的第4波段表观反射率与第3波段表观反射率的比值大于O. 87，且小于I. 15，则该像元被判定为云像元，赋值为1，其余无云像元赋值为O,生成第三层云掩膜。
7.如权利要求I所述的环境一号HJ-1A/B星C⑶数据的云检测方法，其特征在于，S6进一步包括如下步骤 S6. I :将S3-S5中所述第一、二、三层云掩膜进行相加求和，生成新的云掩膜； S6. 2 :将所述新的云掩膜中，值为I的像元赋值为O，表示无云像元；值为2的像元保留原值不变，表示可能为云的像元；值为3的像元赋值为1，表示云像元。
8.如权利要求1-7所述的环境一号HJ-1A/B星CXD数据的云检测方法，其特征在于，所述环境一号HJ-1A/B星CXD数据为经过几何校正的第1、2、3、4四个波段的GeoTiff格式数据。
9.如权利要求2所述的环境一号HJ-1A/B星CXD数据的云检测方法，其特征在于，所述M0D09Q1数据为M0D09L3级250m的8日合成格网数据，包含MODIS传感器的第I和第2波段地表反射率及经纬度信息，数据格式为HDF。
10.如权利要求3所述的环境一号HJ-1A/B星CXD数据的云检测方法，其特征在于，所述环境一号卫星CCD数据的辅助数据文件为xml格式文件，其包括定标系数、成像时间辅助信息和观测几何数据。
全文摘要
一种环境一号HJ-1A/B星CCD数据的云检测方法，其方法为利用MOD09Q1数据建立地表反射率数据库；获取HJ-1A/B星的CCD数据各波段数据的表观反射率；将HJ-1A/B星的CCD数据各波段数据的表观反射率数据与地表反射率数据库中选取成像日期内的地表反射率数据进行比较进行云像元识别，生成第一层云掩膜；并利用HJ-1A/B星的CCD数据的第3波段的表观反射率数据值域范围进行云像元识别，生成第二层云掩膜；然后将HJ-1A/B星的CCD数据的第4波段与第3波段表观反射率数据比值的值域范围进行云像元识别，生成第三层云掩膜；最后根据步骤S3-S5的云像元识别结果，得到云掩膜。
文档编号G01S17/95GK102854513SQ20121032657
公开日2013年1月2日 申请日期2012年9月5日 优先权日2012年9月5日
发明者王桥, 陈良富 申请人:环境保护部卫星环境应用中心, 中国科学院遥感应用研究所
专业数控铣床产品供应商
20年生产经验，实力工厂，精工品质，质量有保障
