基于生成对抗网络的三维云构建方法

2022-09-28程文聪王志刚张文军史小康

计算机仿真 2022年8期

程文聪，王志刚，张文军，史小康

(北京航空气象研究所，北京 100085)

1 引言

云会影响飞行员的视觉判断，飞机在起降过程中要考虑低云大雾的影响，积雨云则可能会引发飞机颠簸和积冰，云的结构在一定程度上也反映了强对流风暴等天气系统的发展方向和强度，因此云的三维结构对于航空飞行安全和飞行任务实施有着巨大的影响，一直是民航和军事部门非常关注的气象要素之一。在飞行航线规划中应尽量避免飞行器在云层中穿行，实时航空气象服务保障过程中也应及时提醒机组人员注意大片密集云团，此外云对于红外、电视、激光等空中制导武器的制导性能也有较大影响，掌握云的三维结构对于航空气象服务保障具有重要的实际意义。

当前对于天气雷达未覆盖区域的云三维结构的获取主要有两种途径，一种是基于天气数值预报模式的方法，另一种是基于卫星云图产品反演的方法。天气数值预报模式是能够定量预测未来天气变化的一种技术手段，通过对天气数值预报模式的诊断分析可根据大气物理特征和规律获取未来数十小时各个标准等压面层的云量产品以及云底高、云量等云结构参数产品。但其在实际业务使用过程中存在一定不足和限制，云是目前公认的数值预报模式中重要的不确定因素之一，数值预报的计算是一个复杂的过程，一般只在固定时间点发布数小时后的预报产品，不能及时融合具有较好时效性的实时资料，因而对临近时间点的云参数预测准确率偏低。而且该方法对计算资源有较高的要求，难以在机场气象台实际航空气象业务工作中部署应用。卫星作为云探测的有效工具，其云图能够提供大范围、多时次的宏观和微观云信息，不仅可以检测行星尺度的天气系统，也能对单个对流云团进行监测，基于卫星云图反演的方法通过直接对卫星云图的处理可获取云顶高、云分类、云覆盖等结构信息，但难以获取云底高、云层厚度等重要的云结构参数，因而无法获取云三维结构的完整信息。

在云的三维可视化方面，当前有多种方法可以构建视觉效果较为逼真的三维云，如基于体过程的建模方法、基于分形的建模方法、基于文法的建模方法以及基于粒子系统的建模方法等。由云的物理形成原理可知，云是由无数的水滴和小冰晶随机运动产生的，而基于粒子系统的建模方法正是由基本粒子群按照一定的规则组成物体的形态，因此基于粒子系统生成的云具有丰富的三维结构，与自然环境中的云具有物理相似性，能够构建可视化及交互效果较为真实的三维云。

目前深度学习在图像分类、图像识别和语音识别等多个领域的应用取得了较好的结果，基于卷积神经网络(CNN)的方法在计算机视觉等领域有了广泛的应用，但在图像生成领域中使用欧式距离作为相似性度量的常规深度学习方法通常会倾向于产生模糊的图像。深度生成对抗网络是近年来人工智能领域有较多研究的一个方向，通过生成器和判别器的对抗，实现生成器参数和判别器参数的训练优化，从而建模从特定领域A向特定领域B的数据映射关系，实现产品从领域A向领域B的转换，生成对抗网络模型目前在图像生成领域有较广泛的应用。

本文提出一种基于生成对抗网络的三维云构建方法。在统一的模型下，利用深度学习的非线性映射能力，使用生成对抗网络建立卫星云图产品与云结构参数产品的映射模型，将卫星云图产品重构为云顶高、云底高和云量产品，再通过这些云结构参数构建基于粒子系统的三维云，从而为航空气象服务保障提供云信息资料支持。

2 数据资料

深度学习需要有大量数据作为训练的支撑，本文工作的目标是将卫星云图产品作为源数据，反演云顶高、云底高和云量这三个反映云结构的参数，进而通过引入粒子系统构建三维云。因此需要使用卫星云图产品作为源数据，对应时次的云顶高、云底高和云量产品作为训练目标。具体研究中选取了FY-4A(风云4A)气象卫星4km分辨率各通道产品作为卫星云图源数据，欧洲中期天气预报中心提供的0.25°×0.25°(经纬度)分辨率的ERA5再分析场产品作为云底高和云量产品的训练目标，国家气象卫星中心网站上公布的云顶高2级反演产品作为云顶高产品的训练目标。各产品的覆盖范围均选取为北纬10°～北纬50°、东经90°～东经130°。

1)FY-4A卫星通道产品

FY-4A地球静止气象卫星是我国第二代静止轨道气象卫星的首发星，搭载了多通道扫描成像仪、干涉式大气垂直探测仪、闪电成像仪、空间环境监视仪等多种载荷，本文的工作中选择FY-4A搭载的多通道扫描成像仪通道产品作为源数据。多通道扫描成像仪的主要任务是对地球表面和云物理状态参数进行高频次、高精度、多光谱定量遥感，直接为天气分析和预报、短期气候预测以及环境和灾害监测服务。观测波段覆盖可见光、近红外、短波红外、中波红外和长波红外，既可观测到大尺度天气系统的全貌，又可观测到中、小尺度天气系统的迅速演变过程。多通道扫描成像仪共配备14个通道，包括7个可见光-近红外通道和7个红外通道。14个通道中500m地面分辨率通道1个，1km通道2个，2km通道4个，4km通道 7个，全圆盘观测时间间隔为15分钟。

本文工作中为统一各通道产品分辨率，且不影响工作的一般性，选取4km分辨率的通道产品作为卫星云图源数据，则在考察范围 (北纬10°～北纬50°、东经90°～东经130°) 内单一通道产品分辨率为1000×1000；为配合数值模式产品的时效，云图通道产品也选择逐小时的整点产品。目前在国家卫星气象中心网站上该资料的起始时间为2018年3月12日。

2)欧洲中期天气预报中心第五代数值模式再分析场(ERA5)云底高和云量产品

云底高指云层底所处的高度，云量指云遮蔽天空视野的成数。欧洲中期天气预报中心发布的数值模式产品是实际业务中使用较为广泛的产品之一，各个气象业务部门一般均将其作为天气分析和预报的重要资料源之一。ERA5再分析数据是其发展的第5代全球大气再分析数据，在大气物理模型的基础上同化了30余颗卫星的观测资料以及全球多种观探测资料，具有较高的可信度和准确度。通过欧洲中期天气预报中心网站可下载1979年至今的历史ERA5再分析场产品，本文中使用0.25°×0.25°(经纬度)空间分辨率的逐小时ERA5云底高产品和总云量产品作为云底高和总云量训练目标产品，则在所考察范围 (北纬10°～北纬50°、东经90°～东经130°) 内产品分辨率为160×160。目前在欧洲数值中心网站上该资料的起始时间为1979年1月1日。

3)国家卫星气象中心FY-4卫星反演的云顶高产品

云顶高指云层顶所处的高度，欧洲中期天气预报中心的ERA5再分析产品中没有提供云顶高产品，而国家卫星气象中心利用FY-4A卫星的2个红外通道和1个CO2吸收通道，结合数值预报资料，通过最优估计的迭代计算，生成了全圆盘的云顶高度2级反演产品并公开发布，本文将其作为云顶高训练的目标产品。该产品分辨率为4km，则在所考察范围 (北纬10°～北纬50°、东经90°～东经130°) 内产品分辨率为1000×1000，为配合云底高产品和云量产品共同构建三维云结构，将其降采样为与云底高和云量产品相同的160×160分辨率。目前在国家卫星气象中心网站上该资料的起始时间为2019年1月18日。

3 基于生成对抗网络的云参数获取

采用基于生成对抗网络的模型获取云结构参数产品。产品生成过程如图1所示。

图1 云参数产品生成过程

该问题可形式化为如下表述：

(1)

为生成云参数产品，首先需要收集同时次、同区域的历史云图通道产品和目标云参数产品组成训练数据集，继而构建深度生成对抗网络模型，经训练提取训练数据集中云图通道产品与云参数产品间的对应信息，最后将待推算时次的云图通道产品作为模型输入，模型的输出即为对应时次的云参数产品。

图2 云参数产品生成对抗网络架构

如图2所示，生成对抗网络模型需要训练两个不同的网络，分别为生成网络G和辨别网络D。选用U-Net结构再结合一个下采样层构成生成网络G。U-Net结构是一种增加了跳跃连接的编码-解码器网络，在图像分割领域有较广泛的使用，文中模型使用的U-Net为64层结构。由于卫星云图产品与云参数产品具有不同的分辨率，因此在U-Net结构后增加下采样层将输出产品调整到目标云参数产品的分辨率，下采样层采用双三次插值方法将输出结果的分辨率降为目标云参数产品的分辨率。辨别网络D为一个64层的卷积分类网络，通过计算其输入为真实云参数产品的概率，辨别输入的是真实云参数产品还是模型生成的云参数产品。在模型训练过程中，判别器尝试正确辨别真实云参数产品与生成的云参数产品，生成器则尝试生成尽可能真实的云参数产品使辨别器无法辨别真伪。

令(，)为辨别器准确识别真伪云参数产品的概率，(，)为由卫星云图产品生成云参数产品的函数(为随机噪声)，则对于辨别器，目标函数为找到如下使最大化的模型参数：

argmax((，))+log(1-(，(，))))

(2)

对于生成器，目标函数为找到如下最优化参数

argmax((，(，))

(3)

模型具体实现中，使用二元交叉熵作为损失度量。即对于辨别器而言，损失函数如下

((，)，1)+((，(，))，0)

(4)

其中

(5)

对于生成器而言，文献[14]的研究表明将生成对抗损失与传统的损失函数结合将得到更好的结果。因此令的损失函数为生成对抗损失与1损失的按比例合成，具体损失函数如下

=((，(，))，1)+|-(，)|

(6)

其中和为生成对抗损失和1损失的比例系数。

训练过程以迭代方式执行。首先对辨别器进行训练，对批量输入真实云参数产品和生成器生成的云参数产品，利用损失通过反向传播的方式更新辨别器的参数；继而冻结辨别器的参数，对生成器和辨别器批量输入卫星云图通道产品和对应的真实云参数，计算生成器损失，再通过反向传播的方式更新生成器的参数。重复上述过程直至生成器和辨别器的能力达到平衡。

训练优化方法为使用mini-batch SGD的ADAM方法，学习率为0.0002，优化器动量参数=05，=0999。在推理测试过程中，使用生成器生成目标云参数产品，生成器网络架构与训练阶段的生成器网络模型相同，参数为训练后的生成器模型参数。

4 基于粒子系统的三维云构建

基于深度学习计算得到的云顶高和云底高参数产品，结合特定经纬度网格所限定的平面边框，可构建该网格区域的云立体轮廓，再利用粒子系统根据云量参数生成一定量的云粒子对立体轮廓进行填充，即可得到可视化及交互效果逼真的三维云体，云立体轮廓如图3所示。

图3 地理网格单元上的云立体轮廓示意图

传统的基于粒子系统的三维云构建方法一般强调追求视觉效果而没有使用真实的云结构参数，本文在云立体轮廓中按如下方式填充粒子系统以使生成的三维云与实际空域中的云形态保持基本一致：

1)为每一个地理网格单元生成一个粒子系统，设定该网格单元上粒子系统的轮廓为与云立体轮廓相同的立方体。

2)每一个粒子系统中的粒子数量由式(7)确定，使得粒子的密度与云量线性相关。

=(-)

(7)

其中为粒子数量，为调整系数，为云量，为云顶高，为云底高。

3)在粒子系统中，云粒子需要指定纹理材质以进行渲染，本文使用的粒子纹理材质如图4所示。为获取较好的云粒子视觉效果，设定该纹理材质的亮度按高斯分布从中心向边缘连续递减。

图4 云粒子渲染材质

5 实验分析与演示软件

通过实验分析基于深度学习生成的云顶高、云底高以及云量参数的质量，并展示基于云参数的三维云演示验证软件。

设定云参数质量判别方法为平均绝对误差，在实验中仅计算真实云参数产品和生成云参数产品中不同时为0的格点。对于尺寸为×的生成云参数产品数据以及真实云参数产品数据，相关计算公式如式(8)所示

(8)

其中为样本数据纬度方向的格点数，为样本数据经度方向的格点数。

实验环境配置如下，CPU：Intel Xeon 4116×2、内存：64GB、GPU：NVIDIA GTX 1080ti 11GB×2、硬盘：512GB SSD、操作系统：Ubuntu 18.04、深度学习平台：Pytorch 1.0.1。

本文所使用的FY-4A卫星通道产品以及三个目标云参数产品的资料覆盖时间存在一定区别。FY-4A卫星通道产品的资料起始时间为2018年3月12日；目标云底高和云量数据来源于欧洲ERA5数值模式再分析产品，资料起始时间为1979年1月1日；目标云顶高数据来源于FY-4A卫星2级反演数据，资料起始时间为2019年1月18日。为最大化利用训练数据，使用上述生成对抗网络模型针对三个云参数分别独立进行建模训练。云底高和云量产品使用的输入数据为FY-4A全部14个通道的产品数据，训练数据时间区间为2018年3月12日至2019年6月31日的逐小时数据，由于夜晚可见光通道无数据，仅使用白天10个时次(北京时08时～17时)的数据，共计4455组；由于云顶高参数在理论上可由卫星红外通道产品计算，因此使用的输入数据为FY-4A的7个红外通道产品，训练数据时间区间为2019年1月18日至2019年6月31日全天的逐小时数据，共计3789组。

模型均经历200个训练周期(epoch)的训练，云顶高产品的训练时长为3小时16分，云底高产品的训练时长为3小时29分，云量产品的训练时长为3小时23分。

5.1 结果统计分析

参与统计分析的测试数据为2019年7月1日至2019年7月31日的逐小时数据。由于夜晚可见光通道无数据，因而仅取白天10个时次(北京时08时～17时)共302个时次的数据(存在缺测时次)。表1为云顶高产品、云底高产品和云量产品的平均MAE误差以及测试数据中最大、最小MAE误差统计。统计结果显示，测试数据中云顶高产品的平均MAE误差为1180米，云底高产品的平均MAE误差为1185米，根据航空气象保障的一般经验，在云高范围0米～20000米的尺度下，1000米左右的云底高平均误差、2000米以下的云顶高平均误差产品基本可用。测试数据中云量产品的平均MAE误差为1.78，在实际航空气象服务保障过程中一般云量误差低于2即认为基本准确。由统计结果可知基于本文方法根据卫星云图通道产品生成的云顶高、云底高和云量产品基本能够满足重构三维云的使用需要。

表1 生成的云参数产品结果统计

5.2 个例分析

2019年7月22日12时～15时(北京时)4个连续时次计算生成的云顶高、云底高以及云量参数产品与目标云参数产品和可见光云图的对比如图5～图7所示。图中可发现4个时次的计算云顶高产品、云底高产品以及云量产品均与可见光云图有相似的外形分布，表明从云图通道产品中提取到了相关有效信息；由热度图颜色对比可知，4个时次的计算生成产品与再分析云底高产品、云量产品以及卫星反演的云顶高产品在值域范围上保持了基本一致，说明深度学习网络准确建模了相应云参数产品的特征，完成了从卫星云图通道产品数据域向云参数产品数据域的转换。个例分析表明本文所提方法能够较准确地计算云参数产品。

图5 2019年7月22日12时～15时(北京时)云顶高计算结果对比

图6 2019年7月22日12时～15时云底高计算结果对比

图7 2019年7月22日12时～15时云量计算结果对比

5.3 三维云演示软件

为展示本文所提方法构建三维云的可行性和实用性，基于Unity3D平台开发了三维云演示软件，图8为基于云顶高、云底高以及云量参数构建的2019年7月22日12时北纬10°～北纬50°、东经90°～东经130°空域的三维云立体轮廓示例(使用透明度表示云量)。

图9为使用粒子系统构建的该时次三维云示例。在两个三维云示例下方叠加了对应时次可见光卫星云图产品进行对比，可以观察到构建出的三维云立体结构与实际云图保持了较高程度的一致性。

图8 云立体轮廓示例

图9 三维粒子云示例

图10 根据三维云信息进行航线规划的示例

图10演示了根据三维云信息进行飞行航线规划的一个操作示例。在图中右下方平面云图上指定航线后(红色线)，可在图中左下部分生成航线的云剖面图，通过在航线云剖面图中设定航线飞行高度(蓝色线)可规避航线中的大片云团。从而使三维云信息参与到航线规划过程中，图11展示了在第一视角下云环境中飞行的效果示例。通过点击右下方的平面云图，可以在右上方高度表上标识出云底高及云顶高，云量信息由高度表中白色区域的亮度标识，从而让使用人员了解感兴趣位置的云结构参数信息。由于该软件主要用于演示三维云信息对飞行任务规划以及航空气象服务保障的支持能力，因而与实际航线规划过程不完全吻合，仅为演示示例使用。