基于多普勒天气雷达资料的阵风锋自动识别方法技术

本发明专利技术公开了基于多普勒天气雷达资料的阵风锋自动识别方法，所述方法包括以下步骤：1)依据阵风锋的雷达表现特征利用局部二值化算法提取出弱窄带回波疑似区域；2)对弱窄带回波疑似区域进行分割、连接和筛选处理，得到弱脊状带对应的骨架图像；3)由当前时刻和前一时刻的两幅低仰角雷达图像得到光流场，将步骤2)得到的弱脊状带对应的骨架图像中前后时刻匹配的骨架拟定为疑似阵风锋，根据该疑似阵风锋的位置和速度与风暴单体的位置和速度的关系及该疑似阵风锋的走向与速度的关系识别出阵风锋。本方法实现了阵风锋的自动检测，对灾害进行及时的预警，减少了经济损失和人员伤亡；并通过实验验证了本方法的有效性。

【技术实现步骤摘要】

本专利技术涉及气象学领域，特别涉及一种基于多普勒天气雷达资料的阵风锋自动识别方法。
技术介绍
天气雷达是对强对流天气进行监测和预警的主要工具之一。天气雷达发射脉冲形式的电磁波，当电磁波遇到降水物质(雪花、雨滴和冰雹等)，大部分能量继续前进，而有一小部分能量被降水物质向四面八方散射，向后散射的能量被雷达接收。在对流单体的成熟阶段，冷性下沉气流作为一股冷空气，在近地面的底层向外扩展，与单体运动前方的暖湿气流交汇形成锋面。这种现象反映在雷达反射率图上，通常可以观测到一条若隐若现、粗细不等、强度较弱且取值不定的弱脊带状区域，气象上称之为阵风锋或雷暴的出流边界。阵风锋常引起气压突变、风速风向剧烈变化以及气温骤降等强烈天气现象，会导农作物倒伏、大树或树枝折断、威胁飞机的起降，严重影响人民大众的生命财产安全。从图像角度看，做阵风锋区域的横剖面，其反射率值的分布呈脊状，但由于区域自身反射率值普遍较低，致使由区域两侧过渡到区域外围的反射率值的梯度变化微弱，且时常出现部分区域段混于大面积的弱回波区域之中的现象，致使原本就不太强的梯度变化在这部分区域段基本消失。这就使常规的基于边缘和区域图像分割方法用于检测这种阵风锋区域时奏效甚微,另外，弱回波区域中通常含有脊状非阵风锋带的干扰，这些干扰甚至与真正的阵风锋相邻或相交。使识别阵风锋方法的空识率与击中率之间的矛盾愈发突出。在基于多普勒雷达数据的阵风锋识别算法中，Delanoy和Troxel[1]较早提出了基于模糊理论的函数模板相关法,该方法将模板尺度固定为17×7，为使检测算法与阵风锋的方向无关，该模板需围绕其中心多次旋转。OsamaAlkhouli和VictorDeBrunner[2]应用熵函数模板法自动识别边界层辐合线，且方法无需旋转模板。郑佳峰等[3]结合速度场与强度场数据检测阵风锋。在速度场数据中，检测辐合带；在反射率强度场中，该方法主要使用的是双向梯度法寻找弱窄带回波，两者结合确定阵风锋。专利技术人在实现本专利技术的过程中，发现现有技术中至少存在以下缺点和不足：首先，人为观测费时费力，会影响到预报的时效性。文献[1]只适合检测宽度不大于3与周边区域存在较大梯度的阵风锋。文献[2]应用熵函数模板法提高了检测速度，但仅对典型的阵风锋有效，对于部分阵风锋区域融入背景的情况，会使检测结果出现断裂。文献[3]的方法在速度场数据中出现无效速度的区域时无法发挥；在反射率强度场时，双向梯度法限定窄带回波与母体回波距离。至今为止，在对阵风锋检测的文献中尚未见到阵风锋与其他弱窄带回波碰撞产生交叠现象的解决方案，也未见到在阵风锋的检测中出现断裂现象的处理方法。[参考文献][1]DelanoyRL,TroxelSW.Amachineintelligentgustfrontalgorithmfordopplerweatherradars[C].ContributionstotheAmericanMeteorologicalSociety’s26thInternationalConferenceonRadarMeteorology.1993:9.[2]AlkhouliO,DeBrunnerV.Gustfrontdetectioninweatherradarimagesbyentropymatchedfunctionaltemplate[C].ImageProcessing,2005.ICIP2005.IEEEInternationalConferenceon.IEEE,2005,1:I-645-8.[3]郑佳锋,张杰,朱克云,等.阵风锋自动识别与预警[J].应用气象学报,2013,24(1):117.[4]OjalaT,M,HarwoodD.Acomparativestudyoftexturemeasureswithclassificationbasedonfeatureddistributions[J].Patternrecognition,1996,29(1):51-59.[5]G.Two-framemotionestimationbasedonpolynomialexpansion[M].ImageAnalysis.SpringerBerlinHeidelberg,2003:363-370.
技术实现思路
本专利技术公开了一种基于多普勒天气雷达资料的阵风锋自动识别方法，可以解决的技术问题包括：自动识别不同宽度的阵风锋；能识别与周边区域的梯度较小或部分域融入背景的阵风锋；能分离与其他弱窄带回波碰撞产生交叠的阵风锋；对阵风锋检测中出现断裂现象进行有效连接的处理；实现准确、完整自动识别阵风锋的目的。为了解决上述技术问题，本专利技术提出的一种基于多普勒天气雷达资料的阵风锋自动识别方法，包括以下步骤：步骤一、依据阵风锋的雷达表现特征，利用局部二值化算法提取出弱窄带回波疑似区域；步骤如下：1-1)设一幅低仰角雷达图像的大小为N×N，对中间[N-(2n+1)×(2n+1)]×[N-(2n+1)×(2n+1)]区域内的像素点pij进行其是否属于弱脊状区域的判断；1-2)求出区域[N-(2n+1)×(2n+1)]×[N-(2n+1)×(2n+1)]内大于等于35dBZ的连通区域ωi，i＝1,2，…s，计算该连通区域ωi的面积si及该连通区域ωi外包矩形长度li，对于面积si小于S1或外包矩形长度li小于L1的连通区域进行标记，得到ωj′，j＝1,2，…m，m≤s；1-3)若像素点pij的反射率值f(i,j)∈x1或f(i,j)∈x2且f(i,j)∈ωj′，其中x1＝[5,35)dBZ，x2＝[35,40)dBZ，执行1-4)，否则，执行1-5)1-4)在以该像素点pij为中心的区域根据式(1)做卷积运算，得到卷积运算结果g1(i,j)和g2(i,j)，当f(i,j)≥g1(i,j)且f(i,j)>g2(i,j)，则认为像素点pij属于弱窄带回波疑似区域Ω，并将该像素点pij置为前景，否则，置为背景；1-5)将像素点pij置为背景，至此将低仰角雷达图像转化为一张二值图；1-6)计算步骤1-5)所形成的二值图中各个连通区域的面积和外接矩形长度，将连通区域面积小于S2或外接矩形长度小于L2的连通区域置为背景；从而提取出弱窄带回波疑似区域；步骤二、对弱窄带回波疑似区域进行分割、连接和筛选处理，得到弱脊状带对应的骨架图像；步骤如下：2-1)对步骤一提取出的弱窄带回波疑似区域的轮廓进行除毛刺处理，再进行细化得到区域的骨架图像A；2-2)在上述骨架图像A中的骨架交叉点和折点处断开骨架，得到骨架图像B，包括：依照折点的特点通过计算骨架上某点两侧切线夹角识别出折点，即沿骨架行进，若由某点到其前侧第n个点的向量与该点到其后侧第n个点的向量的夹角小于135度，则认为该点为折点；利用局部二值模式算子检测出端点和交叉点，即针对步骤2-1)得到的骨架图像，骨架图像中取值为1的点p为可能的端点或交叉点，以该点p为中心，考察其3×3区域及5×5区域边界的取值分布；其中，5×5区域边界点若取值为1，但在5×5范围内与区域中心不连通，则将该点置为0，分别从上述3×3区域及5×5区域的左上角开始沿逆时针方向形成8位01链码和16位01链码来描述所本文档来自技高网...

【技术保护点】
一种基于多普勒天气雷达资料的阵风锋自动识别方法，其特征在于，包括以下步骤：步骤一、依据阵风锋的雷达表现特征，利用局部二值化算法提取出弱窄带回波疑似区域；步骤如下：1‑1)设一幅低仰角雷达图像的大小为N×N，对中间[N‑(2n+1)×(2n+1)]×[N‑(2n+1)×(2n+1)]区域内的像素点pij进行其是否属于弱脊状区域的判断；1‑2)求出区域[N‑(2n+1)×(2n+1)]×[N‑(2n+1)×(2n+1)]内大于等于35dBZ的连通区域ωi，i＝1,2，…s，计算该连通区域ωi的面积si及该连通区域ωi外包矩形长度li，对于面积si小于S1或外包矩形长度li小于L1的连通区域进行标记，得到ωj′，j＝1,2，…m，m≤s；1‑3)若像素点pij的反射率值f(i,j)∈x1或f(i,j)∈x2且f(i,j)∈ωj′，其中x1＝[5,35)dBZ，x2＝[35,40)dBZ，执行1‑4)，否则，执行1‑5)1‑4)在以该像素点pij为中心的区域根据式(1)做卷积运算，得到卷积运算结果g1(i,j)和g2(i,j)，当f(i,j)≥g1(i,j)且f(i,j)>g2(i,j)，则认为像素点pij属于弱窄带回波疑似区域Ω，并将该像素点pij置为前景，否则，置为背景；g1(i,j)=Σs=-ns=nΣt=-nt=nw(s,t)f(x+s,y+t)w(s,t)=1n2g2(i,j)=Σs=-(2n+1)s=2n+1Σt=-(2n+1)t=2n+1w(s,t)f(x+s,y+t)w(s,t)=1(2n+1)2---(1)]]>1‑5)将像素点pij置为背景，至此将低仰角雷达图像转化为一张二值图；1‑6)计算步骤1‑5)所形成的二值图中各个连通区域的面积和外接矩形长度，将连通区域面积小于S2或外接矩形长度小于L2的连通区域置为背景；从而提取出弱窄带回波疑似区域；步骤二、对弱窄带回波疑似区域进行分割、连接和筛选处理，得到弱脊状带对应的骨架图像；步骤如下：2‑1)对步骤一提取出的弱窄带回波疑似区域的轮廓进行除毛刺处理，再进行细化得到区域的骨架图像A；2‑2)在上述骨架图像A中的骨架交叉点和折点处断开骨架，得到骨架图像B，包括：依照折点的特点通过计算骨架上某点两侧切线夹角识别出折点，即沿骨架行进，若由某点到其前侧第n个点的向量与该点到其后侧第n个点的向量的夹角小于135度，则认为该点为折点；利用局部二值模式算子检测出端点和交叉点，即针对步骤2‑1)得到的骨架图像，骨架图像中取值为1的点p为可能的端点或交叉点，以该点p为中心，考察其3×3区域及5×5区域边界的取值分布；其中，5×5区域边界点若取值为1，但在5×5范围内与区域中心不连通，则将该点置为0，分别从上述3×3区域及5×5区域的左上角开始沿逆时针方向形成8位01链码和16位01链码来描述所述两个区域边界的取值分布；然后分别沿两个区域边界循环一周，记录依次取值变化的次数n3(p)和n5(p)，若n3(p)＝2，则点p为端点；若n3(p)≥6或n5(p)≥6，则点p为交叉点；在上述折点和交叉点处断开骨架，从而得到骨架图像B；2‑3)通过端点匹配的方法将骨架图像B中任意两段不连通的曲线连接；根据阵风锋走向平缓的特点，设：曲线li的端点A和曲线lj(j≠i)的端点B的匹配条件如下：式(2)中，L＝30km，Φ1＝Φ2＝Φ3＝‑0.7，li和lj表示曲线i和j的长度，指由端点A沿曲线i在端点A处切线的方向指向远离该曲线的向量，指由端点B沿曲线j在端点B处切线的方向指向远离该曲线的向量；当端点A仅与一个端点符合匹配条件时，则该端点为端点B，连接端点A和端点B；当端点A与多个端点符合匹配条件时，则将其中值最小的端点作为端点B，连接端点A和端点B，其中，端点C是与端点A符合匹配条件的任意一端点，指由端点C沿其所在曲线在端点C处切线的方向指向远离该曲线的向量；从而形成骨架图像C；2‑4)判断低仰角雷达图像是否存在单侧脊区域，若存在，则剔除单侧脊区域的骨架；对骨架图像C中的每条曲线分别进行PCA处理，至少得到第二主成分方向的单位向量e2；当距离k分别为2,3,4,5km时，设变量αk的初值为0，对每条曲线上m个点中每一点pi，i＝1,2...m,求得所述曲线一侧的位置qi，qi＝pi+k e2对于每一点pi计算：αk=αkf(pi)>f(qi)αk=αk+1f(pi)=f(qi)αk=αk+2f(pi)<f(qi),i=1,2...m---(3)]]>式(3)中，f(x)表示点x的反射率值；设，所述曲线一侧的加权比例当距离k分别为‑2,‑3,‑4,‑5km时，设，所述...

【技术特征摘要】
1.一种基于多普勒天气雷达资料的阵风锋自动识别方法，其特征在于，包括以下步骤：步骤一、依据阵风锋的雷达表现特征，利用局部二值化算法提取出弱窄带回波疑似区域；步骤如下：1-1)设一幅低仰角雷达图像的大小为N×N，对中间[N-(2n+1)×(2n+1)]×[N-(2n+1)×(2n+1)]区域内的像素点pij进行其是否属于弱脊状区域的判断；1-2)求出区域[N-(2n+1)×(2n+1)]×[N-(2n+1)×(2n+1)]内大于等于35dBZ的连通区域ωi，i＝1,2，…s，计算该连通区域ωi的面积si及该连通区域ωi外包矩形长度li，对于面积si小于S1或外包矩形长度li小于L1的连通区域进行标记，得到ωj′，j＝1,2，…m，m≤s；1-3)若像素点pij的反射率值f(i,j)∈x1或f(i,j)∈x2且f(i,j)∈ωj′，其中x1＝[5,35)dBZ，x2＝[35,40)dBZ，执行1-4)，否则，执行1-5)1-4)在以该像素点pij为中心的区域根据式(1)做卷积运算，得到卷积运算结果g1(i,j)和g2(i,j)，当f(i,j)≥g1(i,j)且f(i,j)>g2(i,j)，则认为像素点pij属于弱窄带回波疑似区域Ω，并将该像素点pij置为前景，否则，置为背景；g1(i,j)=Σs=-ns=nΣt=-nt=nw(s,t)f(x+s,y+t)w(s,t)=1n2g2(i,j)=Σs=-(2n+1)s=2n+1Σt=-(2n+1)t=2n+1w(s,t)f(x+s,y+t)w(s,t)=1(2n+1)2---(1)]]>1-5)将像素点pij置为背景，至此将低仰角雷达图像转化为一张二值图；1-6)计算步骤1-5)所形成的二值图中各个连通区域的面积和外接矩形长度，将连通区域面积小于S2或外接矩形长度小于L2的连通区域置为背景；从而提取出弱窄带回波疑似区域；步骤二、对弱窄带回波疑似区域进行分割、连接和筛选处理，得到弱脊状带对应的骨架图像；步骤如下：2-1)对步骤一提取出的弱窄带回波疑似区域的轮廓进行除毛刺处理，再进行细化得到区域的骨架图像A；2-2)在上述骨架图像A中的骨架交叉点和折点处断开骨架，得到骨架图像B，包括：依照折点的特点通过计算骨架上某点两侧切线夹角识别出折点，即沿骨架行进，若由某点到其前侧第n个点的向量与该点到其后侧第n个点的向量的夹角小于135度，则认为该点为折点；利用局部二值模式算子检测出端点和交叉点，即针对步骤2-1)得到的骨架图像，骨架图像中取值为1的点p为可能的端点或交叉点，以该点p为中心，考察其3×3区域及5×5区域边界的取值分布；其中，5×5区域边界点若取值为1，但在5×5范围内与区域中心不连通，则将该点置为0，分别从上述3×3区域及5×5区域的左上角开始沿逆时针方向形成8位01链码和16位01链码来描述所述两个区域边界的取值分布；然后分别沿两个区域边界循环一周，记录依次取值变化的次数n3(p)和n5(p)，若n3(p)＝2，则点p为端点；若n3(p)≥6或n5(p)≥6，则点p为交叉点；在上述折点和交叉点处断开骨架，从而得到骨架图像B；2-3)通过端点匹配的方法将骨架图像B中任意两段不连通的曲线连接；根据阵风锋走向平缓的特点，设：曲线li的端点A和曲线lj(j≠i)的端点B的匹配条件如下：式(2)中，L＝30km，Φ1＝Φ2＝Φ3＝-0.7，li和lj表示曲线i和j的长度，指由端点A沿曲线i在端点A处切线的方向指向远离该曲线的向量，指由端点B沿曲线j在端点B处切线的方向指向远离该曲线的向量；当端点A仅与一个端点符合匹配条件时，则该端点为端点B，连接端点A和端点B；当端点A与多个端点符合匹配条件时，则将其中值最小的端点作为端点B，连接端点A和端点B，其中，端点C是与端点A符合匹配条件的任意一端点，指由端点C沿其所在曲线在端点C处切线的方向指向远离该曲线的向量；从而形成骨架图像C；2-4)判断低仰角雷达图像是否存在单侧脊区域，若存在，则剔除单侧脊区域的骨架；对骨架图像C中的每条曲线分别进行PCA处理，至少得到第二主成分方向的单位向量e2；当距离k分...

【专利技术属性】
技术研发人员：王萍，袁悦，
申请(专利权)人：天津大学，
类型：发明
国别省市：天津;12