[发明专利]基于高精度数值风场资料的涡旋中心自动识别方法

摘要

本发明公开了一种基于高精度数值风场资料的涡旋中心自动识别方法，包括：1)依据850pha数值风场的矢量特征，利用衡量风向杂乱度的特征提取出包含轴对称型风场的感兴趣区域；2)利用圆形数据理论对理想的轴对称型风场进行特征描述和分类；3)相对于复杂的轴对称风场依据步骤2)的结论和构建变形系数，初步检索出涡旋候选中心点集；4)改进Ward聚类算法得到属于同一涡旋系统的中心点集；5)利用步骤3)的候选中心点集的涡度对涡旋系统分类，并精确定位全局的涡旋中心。本发明可自动识别和定位多种类型的轴对称风场中心，能识别多尺度的涡旋系统及复杂风场中带有一定形变的涡旋系统；从而实现涡旋系统的准确、完整自动识别和涡旋中心定位。

主权项

1.一种基于高精度数值风场资料的涡旋中心自动识别方法，其特征在于，包括以下步骤：步骤一、依据850pha数值风场的矢量特征，利用衡量风向杂乱度的特征提取出包含轴对称型风场的感兴趣区域；步骤如下：1‑1)设输入的风场矩阵为W_M,N，其中，M,N分别为矩阵的行数和列数；设风场矩阵W_M,N中的两个矢量的方向夹角为Λ，构建一个统计特征即：在L×L窗口内，分别计算风场矩阵W_M,N中的网格点(i,j)和该网格点的相邻点(i，j+1)与(i‑1,j)的方向夹角Λ，并累积求和；设风矢量的角度分量为θ，则两个矢量的方向夹角为Λ(θ₁,θ₂)，统计特征的计算公式如下：式(1)中，i，j，m和n分别是风场矩阵W_M,N中网格的东西方向和南北方向的索引；m＝r+2，…，M‑r；n＝r+1，…，N‑(r+1)，L为窗口的尺度，且为奇数；1‑2)对得出的统计特征进行归一化得到归一化的统计特征对输入的给定分辨率的高精度网格风场数据，计算出统计特征的最大值，即对进行归一化得到计算公式为最后，遍历计算W_M,N，得到归一化的特征矩阵即为：1‑3)利用Bayesian决策对参数L进行优化，和对公式(1)得出的η^L进行阈值选取，以获得输入风场矩阵的最优感兴趣区；从给定分辨率的风场数据中搜集涡旋样本V，等分为三组测试样本，包括：V_A，V_B和V_C，其中，风场数据包含不同季节、时间和海拔高度；然后，对杂乱风样本M的三个子样本：M_A，M_B和M_C和一致风样本S的三个子样本：S_A，S_B和S_C分别计算当L＝5,7,9,11和13的五类训练子特征：和当L取上述给定某一尺度值时，分别计算对于类别V_A，M_A和S_A的条件概率密度：和根据最小错误率的原则，获得的最优阈值和相应的最小错误率同样的过程应用于[V_B,M_B,S_B]和[V_C,M_C,S_C]，获得相应的最小错误率和比较不同的L取值下的最小错误率的平均值，求得对于给定分辨率的风场数据的最优L取值；基于最优L值，确定的最优阈值，该阈值由与杂乱区分界的最大值和与一致风向区域分界的最小值组成；步骤二、利用圆形数据理论，对理想的轴对称型风场，包括：顺时针旋转场、逆时针旋转场、辐合场、辐散场和鞍形场进行特征描述和分类；步骤如下：2‑1)对于两个圆形数据变量α和β，它们之间的一种有向的方向差Q(α,β)，定义为：定义α到β的一种走向为：其中，α和β的取值范围各为[‑180°,180°]；α到β的走向α→β为：在极坐标系下从α到β的走向；2‑2)结论1：设θ_k是以理想轴对称风场中心为中心的单位圆弧上k位置处风的角度，当0°＜Q(θ_k,θ_k+ε)＜180°,ε＞0时，则θ_k→θ_k+ε:counterclockwise；对于所有的0°≤k＜360°，结论1都成立时，判定θ_k为第一类轴对称场，即涡旋场或辐合场或辐散场；结论2：设θ_k是以理想轴对称风场中心为中心的单位圆弧上k位置处风的角度，当‑180°＜Q(θ_k,θ_k+ε)＜0°,ε＞0时，则θ_k→θ_k+ε:clockwise；对于所有的0°≤k＜360°，结论2都成立时，判定θ_k为第二类轴对称场，即鞍形场；2‑3)设是轴对称型风场，对于所有的0°≤k＜360°，当时，和θ_k同属一类轴对称场；如果θ_k为逆时针涡旋场，当δ＝90°时，为辐合场；当δ＝180°时，为顺时针涡旋场；当δ＝270°时，为辐散场；如果θ_k为鞍形场，δ取0°～360°任意值，为鞍形场；步骤三、相对于理想轴对称风场，实际的涡旋场是复杂的轴对称风场，依据步骤二的结论和构建变形系数，初步检索出涡旋候选中心点集；步骤如下：3‑1)对步骤一计算出的感兴趣区域的每个格点风向值进行角度量化：依次为：东风、东北风、北风、西北风、西风、西南风、南风和东南风八个方向；3‑2)以步骤一计算出的感兴趣区域为当前检索范围，依次对每个以当前格点为中心的d×d窗口下分析域的风向进行特征提取，判定当前格点是否为涡旋中心点；从中心点向X轴方向的射线开始，以45°角为步长，依次引出八条射线对分析域进行八等分；3‑3)涡旋场的特征提取；对于当前分析域的第i个子区域，所有格点数为S_i、具有相同角度α的连通格点数为则该子区域的特征值定义为c_i，因此当前分析域的特征向量描述为C＝{c₁,c₂,…c₈}；C的提取规则为：IF存在α使得c_i＝α；ELSE停止并检索下一格点；END其中，ρ是阈值变量，默认值是0.4；3‑4)将C简化为设c_i和c_j是两个相邻子区域的特征值；即，其中，如果i＜8,j＝i+1,如果i＝8,j＝1；简化规则为：IFc_i＝c_j去掉c_jEND3‑5)剔除不合理的设q是的长度，和是的相邻元素，其中，如果i＜8,j＝i+1,如果i＝8,j＝1；则规则为：IF q＜4或存在c_i使得停止并检索下一格点；END3‑6)设为中方向角所对应的空间子象限的索引向量；对通过步骤3‑5)的和映射为单位圆的内接多边形；规则是：以或中元素的值作为单位圆的内接多边形的顶点；或映射的19种圆内接多边形如表1所示；表13‑7)对步骤3‑6)中的所有多边形构建形变指数λ，去除形变过度的候选中心，并获得确定下的涡旋的形变指数，如表1所示；λ的构建方法如下：设圆形数据向量A＝{α₁,α₂,...,α_n}，其中，α_i≠α_j，当0°＜Q(α_i,α_i+1)＜180°且0°＜Q(α_n,α₁)＜180°时，构建形变指数λ；λ是圆形数据向量A的合成向量长度、圆形数据向量A均匀波动性和圆形数据向量A的空缺率的线性组合；(1)合成向量长度：r_or_o的取值范围为[0，1]；r_o越接近1，α₁,α₂,...,α_n越聚焦于它们的均值；(2)均匀波动性：δ_o圆形数据向量A的有向方向差为Q(A),则Q(A)的标准差为：同样地，表1中所有圆内接多边形的δ值中的最大值为1.73，得出标准化后的δ_o：(3)空缺率：ρ_o设圆形数据向量A的维度为n,空缺率ρ_o的计算如下：最后，将上述三种特征线性组合，得到一个综合的形变指数λ：λ＝ar_o+bδ_o+cρ_o                                               (9)式(9)中，0＜a＜1,0＜b＜1,0＜c＜1且a+b+c＝1，a＝b＝0.25,c＝0.5时，所有圆内接多边形的形变指数λ值如表1所示；3‑8)排除涡旋的伪中心点从风场数据中挑选涡旋样本估计λ的值；首先，每个样本的和分别由d＝5,7,9的分析域计算得到；利用步骤1‑3)中的Bayesian决策得到最优λ值，从而得到：或停止并检索下一格点；END步骤四、改进Ward聚类算法，对属于同一个涡旋系统的中心点集聚类；步骤如下：4‑1)对步骤三使用d＝5，7和9对感兴区域进行多尺度窗口的搜索候选中心点集P，其中P包含了多个涡旋系统的候选中心点集；利用Ward聚类算法将P中的各点聚类到各自的涡旋系统得到P_k，k＝1,2,…,K；K为类别数；4‑2)对于类别数K进行自适应选取；步骤4‑1)中，每步类合并计算的最小离差平方和为D_i，计算D_i在最大离散曲率T_i处对应的类别数为最终认定的K；T_i的计算方法如下：T_i＝1‑|cosω_i|                        (10)式(10)中，ω_i＝angle(p_i‑1,p_i,p_i+1)是线段[p_i‑1,p_i]和[p_i,p_i+1]的夹角，而K＝argmax_i{T_i}；步骤五、利用步骤三的候选中心点集的涡度对涡旋系统分类，并精确定位全局的涡旋中心；5‑1)得到了每类涡旋系统的候选中心点集P_k后，采用数学形态学中的膨胀算法对P_k进行膨胀处理，以完全覆盖涡旋中心区域，方法如下：式(11)中，是膨胀操作符，B是一个3×3的结构元素；5‑2)将涡旋类别判定为顺时针方向或逆时针方向；引进二维连续风场的垂直涡度ζ：将式(12)转换为离散风向场的垂直涡度式(13)中，θ_i,j是单位矢量的角，i和j是水平方向和竖直方向的索引；令l₁或l₂记为Z_k元素的正或负涡度的个数；l₁或l₂定义为如下：式(14)中，#{·}表示集合{·}的基数；如果l₁＞l₂，当前第k个涡旋判定为逆时针涡旋，否则判定为顺时针涡旋；5‑3)根据涡度值确定涡旋中心参考点；在涡旋中心区域Z_k，由公式(13)计算每个格点的垂直涡度将垂直涡度绝对值最大的检测点作为涡旋中心的参考点；5‑4)根据零风速约束精确定位涡旋中心；在到涡旋中心参考点距离小于ε的区域内，ε＝3，将网格划分为两个三角网格，在三角网格内通过线性插值算法计算出零风速点坐标位置，连接两个零风速点即可得到零等值线段；由此得到风速东西分量U和南北分量V在一个网格内的零等值线段，同一网格上U、V的等值线段的交点，即为风速零点；若整个局部区域存在多个风速零点，则选择距离气旋中心参考点最近的点作为气旋中心；若整个局部区域不存在风速零点，则对风速分量U、V分别进行线性插值，求解风速矢量U、V的强度极小值，将气旋中心点校正到风速矢量U、V极小值点位置。