青藏高原地区夏季易受降水影响而引发地质灾害，而利用气象卫星云图进行降水的观察与监测成为了灾害预警的一项必要的工作^[1]，在这个过程中，云检测工作成为了应用气象卫星资料进行定量工作的重要基础之一^[2]，它受下垫面、时间、地理位置等因素的影响非常大^[3]，在青藏高原地区尤其如此，国内外众多学者对云检测工作投入了诸多精力。针对不同类型的遥感卫星数据，针对性地提出了不同的云检测方式^[2]，总结来看迄今为止主要的云检测算法包括基于卫星通道数据的阈值法、基于统计学原理的聚类分析方法以及基于神经网络深度学习的方法^[4]。阈值法经历了从固定阈值、自适应阈值到多通道组合阈值的过程。李微等^[5]针对不同下垫面物体光谱特性及其差异，利用MODIS多通道阈值进行云检测，但是研究也存在不同区域选择对云检测结果的影响的问题。马芳等^[6]综合运用FY-2E数据多通道阈值组合以及红外分裂窗区进云检测，检测效果与可见光云区基本吻合，但也存在着地理位置的影响，并且受限于FY-2E卫星通道的情况，所谓的多通道法也仅仅只利用了两个通道。文雄飞等^[7]利用FY-2C数据提出基于云指数的云检测过程，较普通阈值法提高了检测精度。而随着高时间高光谱分辨率卫星的发展，Martinuzzi等^[8]结合地面信息与高分辨率卫星的特点进行云检测，并通过区域目标的特点进行修正；从统计学角度，Christodoulou等^[9]对卫星云图进行多纹理特征提取的方式进行云检测；金祎等^[10]提出了基于密度聚类支持向量机的云检测方法，较好地完成云检测的工作,但是存在检测过程复杂，预处理缓慢的问题；赵敏等^[11]利用SVM以及无监督分类的方式进行云分类，但同样存在分类效果和分类效率的问题。而近年来随着深度学习的火热，越来越多的学者开始利用深度学习方法进行云检测。裴亮等^[12]利用深度学习方法对资源三号卫星遥感影像进行云检测，效果上优于传统方式，检测需要构建训练集测试集以及网络的搭建训练等过程；董蕴雅等^[13]证明了深度学习方法应用到遥感技术领域的可行性；上述方法进行检测的过程中普遍基于单时相的卫星图像进行云检测，算法上较稳定且业务化简单，但是容易出现漏检与误检。近年来随着高时间分辨率卫星数据的增加，利用多时相方法进行云检测具有很大的优势，通过相隔较短的时间的卫星影像，能够有效去除地面静止的干扰信息^[14]。杨昌军等^[15]利用FY-2C时间序列数据进行云检测的研究，提供了如何去利用多时相卫星数据进行云检测的思路。但总的来说，利用静止卫星高时间分辨率的特点进行云检测的研究还较少。

在云检测的基础上，进一步地对云相态进行识别；针对云相态这一热力学过程，能够有效监测出各种天气系统的变化与发展，因此在青藏高原地区进行云相态的识别十分重要。

针对青藏高原地区高海拔特性，云检测过程将面临着高亮地表影响以及地面积雪的影响等诸多问题，传统的检测方法并不适用于该地区。FY-4A卫星是我国新一代静止气象卫星，本文利用其高时间分辨率，高光谱分辨率的特点，首先结合前人云检测经验提出一种多时相多通道阈值组合的云检测方法；其次通过选取15个时间样本进行云检测过程阈值的确定，同时选取另外时间的案例进行算法效果验证，对比云产品结果与传统单时相云检测方法均有较大的提升，证明了本文算法对于青藏高原地区的云检测有明显的效果；最后在云检测基础上进行云相态的识别，利用GPM降水资料验证了云相态识别的精度。

青藏高原是我国最大，世界海拔最高的高原，其南起喜马拉雅山脉，北至昆仑山、阿尔金山，西部为帕米尔高原，东及秦岭山脉并与黄土高原相接，东西长约2 800 km，南北宽约300~1 500 km，有“世界屋脊”之称，对于青藏高原地区进行气象研究与监测是一项重要的工作。青藏高原地区全年大部分降水集中于夏季，极易引发地质灾害^[1]；利用卫星气象资料进行灾害的监测预警需要大量的前期工作，包括云检测。而在青藏高原地区，其特殊的地理特性容易引起云检测的诸多问题，有必要对该地区针对性地进行云检测以及云相态的研究工作。

FY-4A卫星数据从2018年4月开始可以获得，卫星一级数据总共包含14个通道。本文主要使用的通道包括1、2、3、5、8、10、12、13、14等，主要的技术指标如表1所示。

本文所使用的数据是2018年7月FY-4A卫星15 min间隔的全圆盘数据，数据格式为HDF5。本文采用分辨率为4 km的数据。数据前期的预处理包括辐射定标，滤波去除噪声以及对所研究的青藏高原地区进行等经纬投影。

本文使用FY-4A时间序列数据，通过不同波段组合成不同的判断指标，结合多时相数据进行云检测；在云检测的基础上，再对云相态进行识别。算法的流程图如图1所示。其中，各项判断指标如表2所示。以下针对云检测部分与云相态识别分别做介绍。

高原地区云检测的难点就是其海拔高，且存在地面积雪容易引起误判。为进一步说明云雪误检问题，通过人工判别选择了一例包含有地面积雪，高层云以及低层云的卫星数据，分析云雪在各波段的数值情况。

在图2(a)中使用一条直线穿过卫星图像，并将低层云区、高层云区以及地面积雪三个部分标注出来；图2(b)是直线上像元的反射率变化，发现无论是云还是雪，各个波段的反射率都比地面反射率高，但是低层云区与地面积雪的反射率无论是数值还是变化趋势都非常接近；图2(c)是直线上亮温变化的情况，可以很明显地发现高层云区的亮温要明显低于地面积雪以及低层云区，但同样地面积雪与低层云区的亮温非常接近;因此若仅利用单时相亮温值或反射率进行云检测会产生一定程度上的误判。

而本文提出的多通道多时相的云检测算法，在多通道组合阈值的基础上引入了多时相的思路，能够有效避免积雪对于云检测的干扰。算法使用待检测时相卫星数据通道1、通道2、通道3、通道5（中心波长分别为0.46、0.64、0.86和1.61 μm）反射率和通道12（中心波长为11.0 μm）亮温值，以及前6个时相通道12的亮温值。各项检验对应一个判断指标，每通过一项检验则对应的计算指标I记作1，否则为0。各项检验的原理以及作用如下：

（1）阈值a，指标I1。前人的研究表明，在可见光0.64 µm处，云的反射率一般大于30%，而晴空等区域的发射率则处于0~30%之间，因此日间的云检测过程中，0.64 µm是一个首要考虑的波段。

（2）归一化云指数NDCI,阈值b，指标I2。研究表明，大气中的水汽特征能够在云的0.8~0.95 µm波谱特征反映出来，又由于在可见光0.64 µm波段与短波红外0.86 µm波段，云与地面各种物体的波谱特征反差明显，因此通过归一化云指数可以具体得到云的信息。云在0.64 µm的反射率一般较高，而在0.86 µm处由于水汽的影响，导致其反射率较低归一化云指数R_0.64-0.86如下所示：

其中：R_0.64为0.64 µm波段反射率，R_0.86为0.86 μm波段反射率。

（3）归一化雪指数NDSI，阈值c，指标I3。在可见光波段，云与雪的反射率特性非常相似，所以很容易混合影响判断；而在1.61 µm处，雪会吸收太阳辐射导致反射率降低，而云吸收相对少所以反射率高，因此通过归一化雪指数NDSI，将干扰性很大的雪可以提取出来。NDSI公式如下所示：

其中：R_1.61是1.61 μm波段反射率；

（4）R_0.46 & R_0.64，阈值d,e，指标I4。归一化过程中，部分轮廓线会被误判为云，因此需要用到0.46 µm紫光波段进行修正；同时，利用可见光波段去修正水体的误判；两个阈值通过具体的实验确定。

（5）多时相检测，阈值f，指标I5。多时相检测是核心部分。对于固定区域来说，如果某处上空是晴空，那么在较短的一段时间内该处上空会保持晴空状态，并且温度值也是保持在一定的变化范围内；而如果某个位置由晴空转为有云，那么该区域卫星数据的亮温值就会超过这个范围。因此利用该原理，使用待检测时刻通道12，以及前6个时相通道12的亮温值。计算公式如下:

其中：∆BT为亮温变化范围，BTmax为前6个时刻以及待检测时刻通道12亮温值BT₁₁的最大值；如果亮温变化∆BT超过一定范围，那么就认为该区域为有云像元。

因此综上所述，利用FY-4A数据进行云检测的算法为：

（6）算法阈值确定。在前文算法分析的基础上，利用统计实验的方式确定算法阈值。本文从7月至8月按时间间隔采集15个时间样本进行实验，样本选择了不同地面积雪情况进行实验。通过综合对比调整得到最终的阈值如表3所示。

云相态识别中，通常利用冰云与水运在8~10 µm波段对辐射的吸收基本相同，而在10~13 µm波段上对辐射的吸收却有较大的差异性^[16]。本文使用FY-4A卫星数据，利用0.64 µm可见光波段以及1.61 µm近红外波段，结合10.8 µm与12.0 µm长波红外波段进行云相态检测。冰云在较低的温度下，可以认为是由冰粒子组成的；分析冰，水两相态云在所选择波段的数值变化特点。

如图3(a)所示，沿着纬度方向取一条穿过地面，对流云团直线，对这条线上像元在各波段数值进行分析。如图3(b)所示，在45~70像元区间内,即一对流云团中部，0.64 µm波段反射率在0.6~0.75之间，而1.61 µm波段反射率则处于0.1~0.2之间，两者最大相差约4倍，因此可以利用NDSI计算式可以很容易地将冰云识别出来；另外，冰云区域的10.8 µm亮温值有明显的降低，低于230 K，这也说明云层的厚度大，云顶高，属于典型的对流云特征^[17]；对于水云部分，其在0.64 µm和1.61 µm波段的反射率都较小，因此NDSI值较小；同时10.8 µm与12.0 µm波段的亮温差较大，可以综合利用两者进行识别。此外，结合数值曲线不仅能够帮助寻找不同相态云的波段特征，还能够初步定义选择阈值，根据算法分析部分的波段数值曲线以及15例样本实验的选择调整，最终确定冰，水两相态云的检测阈值如表4所示。

其余未能通过两项判断而剩下的云区认为是冰云与水云的混合态。

为了进一步验证本文算法阈值的适用性，另外选择了世界时间2018年8月29日8：00的青藏高原地区的FY-4A数据进行云检测以及云相态识别的实验，用于验证既定阈值的准确性。同时采用对比实验的方式验证结果准确性，将本文算法的检测结果与传统单时相通道阈值法检测以及国家气象卫星中心提供的云检测产品进行对比，结果如图3所示。

图4(a)是世界时间2018年8月29日08：00（北京时间16：00）FY-4A卫星的4 km分辨率通道3、2、1通道进行RGB合成的彩色图像，查阅青藏地区的天气情况发现，当日拉萨、日喀则等地区发生了对流天气，产生了中雨到小雨不等的降水情况，通过FY-4A卫星可以完整地观测到整个过程。从图4(a)可以发现，整个青藏地区东部以及南部地区上空覆盖有巨大的云区，包括亮温较低，颜色白而发亮的对流云团，以及颜色较浅的带状卷云。由于检测的时间正处于夏季，青藏地区的地面积雪较少，因此无云区域能够清晰地看到地表的情况。但是，虽然处于夏季，青藏地区仍然存在较为零散的积雪分布，部分无云地区的积雪，青藏地区的海拔过高云雪的差异性不明显，因此会对云检测的工作产生很大的影响。本文提出的多时相多通道云检测算法结合待检测时刻之前的数据，对地面上静止不动的积雪做到有效的排除。结合张永宏等^[18]提出的基于NDSI的青藏高原积雪判断方法，将该案例中的积雪提取出来，如图5所示，可以发现青藏高原西北部与新疆地区交界处，存在一长条状积雪，同时其他无云覆盖区域也存在零散的地面积雪的情况，用黄色虚线框标注到图4中。而这些积雪从图4(a)的彩色云图上来看并没有特别明显的差异性，这些积雪情况对于云检测是一个较大的影响。

在进行云检测之前，首先需要对所需要用到的FY-4A卫星各通道数据进行辐射定标以及几何校正的过程，从而使得图像数据能够转化为具有物理意义的反射率值以及亮温值。在这之后进行云检测的实验，通过对选定时刻的卫星数据，在青藏地区部分进行前文所介绍的多通道多时相云检测算法流程，最终结果如图4(b)所示。结合图4(b)与图4可以发现，本文所提出的算法能够非常有效地将地面积雪排除在外，同时对东南部大量云区也能够很好地完成检测。

图4(c)是从中国国家卫星气象中心下载的相对应时刻的FY-4A云检测产品数据，本文算法与该云产品的结果基本一致，但云产品的对于地面积雪的排除并没有很好的效果，仍然存在部分误检测的点,尤其是积雪的边缘区域有较为明显的误差，并且有较多零散的误检测的点；此外，在东部以及南部的大量有云区域，云产品上来看也存在较为明显的漏检以及误检的情况，有云与无云区域交界误检漏检问题较多。总体来说云检测产品存在以下两个问题：第一，存在较多的漏检，有部分很明显的亮白色块状云没有被检测出来；第二，云产品效果存在较为严重的误检情况，云产品将青藏高原地区中部云区绝大多数的地方都检测为云，而结合图4(a)可以发现，中部存在较多的无云区域。

图4(d)是利用传统的单时相通道阈值法对同一时刻数据进行实验的结果。传统通道阈值法仅利用单时相的数据进行检测，对比图4(d)与与本文算法的结果可以很明显地发现，传统的通道阈值法还是存在很明显的漏检以及误检的情况：第一，同样是对无云区域积雪情况的误检测—相对于图4(c)的云产品来说，传统的通道阈值法的检测效果更差，将一整块积雪区域误检为云区；第二，在东部以及南部有云区域存在较为严重的漏检情况。传统方法能够对绝大多数的云团起到很好的检测效果，但是对于有云与无云交界处，特别是一些颜色较浅，亮温较高的薄卷云无法有效完整地检测出来。

分析可知本文算法在云检测目视效果上的优越性，为进一步说明检测精度，利用误检率α，漏检率β以及准确率γ公式进行计算，公式如下：

其中：m₁表示将云像元检测为云像元的个数；m₂表示将云像元检测为非云像元的个数；m₃表示将非云像元检测为云像元的个数；m₄表示将非云像元检测为非云像元的个数。

针对检测时刻数据分别使用传统单时相检测法，云产品数据以及本文算法进行精度验证。选取统一区域，通过人工标注云区的方式，得到本文多时相多通道云检测方法与云产品及传统单时相方法的检测精度如表5所示。

误检率值主要表征将非云像元判断为云像元的情况。从表5可以看出，传统单时相法检测，其误检率为12.8%，要优于云产品的误检率22.5%，从效果图中也不难发现云产品中将大量的非云像元判断为云像元；而本文的算法误检率仅为7.2%，要优于其他两者。漏检率值主要表征对所有云像元的判断情况；由于传统单时相法很明显漏掉了大量的云区没有检测出来，因此它的漏检率是最高的；而云产品在检测效果图中可以发现存在部分白色块状云没有检测出来，导致它的漏检率也较高；本文算法的漏检率仍然是最优。准确率值表征判断为云像元的结果的精确性，此时传统方法与云产品准确率较为接近，但是本文算法的准确率是最高的，达到了94.4%。因此，本文提出的多通道多时相算法，无论是目视效果，还是漏检率、误检率以及准确率，在青藏地区的云检测效果都要优于传统方法以及云产品的效果。

在云检测的基础上进行云相态的识别。图6是云相态识别的结果，对比云相态检测结果与卫星云图可以发现，检测结果将云图中高亮白色的对流云团都检测为冰云，其余大部分云区被判断为水云。

GPM降水资料如图7所示，利用余弦相似度公式对图6与图7进行相似度计算以验证云相态识别的准确性，公式如下：

其中：A，B分别为空间中的两个向量，θ为向量A，B的夹角。cosθ值在[-1,1]范围内。再将其与1做差再取绝对值的方式，将值归一化到[0,1]之间。其值越接近于1，说明两者之间的相似度越大，进而也说明云相态识别与GPM真值越接近。

为计算相似度，本文将图6与图7的降水区域对应一个矩阵，图中每一个降水的像元，对应于矩阵中相应位置为1，否则为0，这样就得到两个降水像元矩阵。利用式（6）对图6与图7两个降水像元矩阵的每一列进行相似度计算，再对列相似度求均值，得到整体的相似度为0.883，说明云相态识别结果与GPM降水结果的相似度高。

CALIPSO(The Cloud-Aerosol Lidar and Infrared Pathfinder Satellite Observation)卫星数据在气溶胶应用研究已经取得了很大的进展^[20，结合CALIPSO卫星所提供的二级VFM（垂直特征掩码，Vertical Feature Mask）产品数据，其能够准确地观测一条廓线上不同海拔高度的云类型分布。为进一步验证本文云相态识别的准确性，在使用GPM降水资料后，再利用CALIPSO卫星云相态观测结果进行验证。

根据卫星的运动轨迹可以将之绘制到云相态识别结果上，如图8(a)所示。图8(b)则分别将CALIPSO卫星实际观测结果与本文云相态的识别结果做对比。从图8(b)可以发现，本文对于云相态的检测结果与CALIPSO卫星的观测结果在26°~30° N之间存在一定的差别，其余区域的云相态是比较吻合的，因此也证明了本文算法在云相态检测上的准确性。

利用静止卫星高时间分辨率的特点进行云检测的研究还较少，本文正是借助我国新一代静止气象卫星FY-4A的全圆盘资料，利用其高光谱分辨率，高时间分辨率的特点，提出了一种可用于青藏高原地区的多时相多通道阈值组合的云检测方法，经过多次对比试验，算法可良好地适用于复杂高原地表下垫面的云检测工作；在云检测的基础上进一步对云相态特征进行研究，以此为基础分析了一次青藏高原地区降水的案例。

本文所提出的青藏高原地区多时相多通道阈值组合云检测与云相态识别方法，能够快速有效地适用于云检测与云相态识别业务场景，充分利用了我国新一代静止气象卫星FY-4A的性能优势。另外需要指出的是，阈值选择的效果受季节，区域的影响较大，因此需要根据不同的时间地域信息进行适当调整；此外受限于静止气象卫星空间分辨率的限制，对于一些细节上的检测还存在着不足，未来可以结合高空间分辨率的极轨卫星资料如FY-3号系列卫星与本文的研究内容相互补充，这将是后续工作的突破点。

（1）对比中国国家气象中心云检测产品以及传统单时相云检测方法，本文提出的方法漏检率为7.2%，误检率为5.6%，准确率为94.4%，均高于其余两种检测方式，充分显示出了本文多时相云检测的优越性。

（2）本文结合前人云相态检测工作，利用FY-4A卫星资料，对青藏高原地区冰云，水云以及地面区域在0.64、1.6、11以及12 μm等4个波段数值进行分析，利用反射率、NDSI以及红外通道亮温差确定云相态检测流程以及检测阈值，实验证明算法能够有效地检测出云顶的相态。

（3）实际案例分析中，对世界时间2018年8月29日08：00（北京时间16：00）的青藏高原地区进行云检测以及云相态识别，并且进一步结合美国国家宇航局的GPM降水资料以及CALIPSO卫星VFM云相态监测结果验证云相态识别的准确性，发现冰云的分布与GPM降水的分布较为一致，相似度达到0.883；检测出的云相态的分布与CALIPSO卫星的监测结果也较吻合，验证了云相态识别的准确性。实际案例中，冰云中心最高日累计降水量达35 mm以上，而由于青藏高原地区一年中绝大部分的降水都集中在了夏季，这样规模的日降水量持续一段时间极易引发诸如滑坡泥石流之类的地质灾害，利用本文研究内容可在一定程度上为青藏高原地区的降水监测、地质灾害预警提供辅助措施。