CN114004259A - 基于改进社区合并的雷达信号密度峰值聚类方法 - Google Patents
基于改进社区合并的雷达信号密度峰值聚类方法 Download PDFInfo
- Publication number
- CN114004259A CN114004259A CN202111297714.XA CN202111297714A CN114004259A CN 114004259 A CN114004259 A CN 114004259A CN 202111297714 A CN202111297714 A CN 202111297714A CN 114004259 A CN114004259 A CN 114004259A
- Authority
- CN
- China
- Prior art keywords
- cluster
- distance
- pulse
- doa
- clusters
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000000034 method Methods 0.000 title claims abstract description 39
- 239000011159 matrix material Substances 0.000 claims abstract description 40
- 238000004088 simulation Methods 0.000 claims abstract description 10
- 238000010606 normalization Methods 0.000 claims description 13
- 230000001174 ascending effect Effects 0.000 claims description 6
- 230000000007 visual effect Effects 0.000 claims description 6
- 238000004364 calculation method Methods 0.000 claims description 4
- 238000012163 sequencing technique Methods 0.000 claims description 4
- 230000008859 change Effects 0.000 claims description 3
- 238000012935 Averaging Methods 0.000 claims description 2
- 230000004075 alteration Effects 0.000 claims description 2
- 230000004927 fusion Effects 0.000 claims description 2
- 230000000630 rising effect Effects 0.000 claims description 2
- 238000012216 screening Methods 0.000 abstract 1
- 238000004422 calculation algorithm Methods 0.000 description 17
- 230000005855 radiation Effects 0.000 description 7
- 238000001514 detection method Methods 0.000 description 3
- 238000010586 diagram Methods 0.000 description 3
- 238000005516 engineering process Methods 0.000 description 3
- 238000005457 optimization Methods 0.000 description 3
- 230000008901 benefit Effects 0.000 description 2
- 230000007547 defect Effects 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000007781 pre-processing Methods 0.000 description 2
- 238000012545 processing Methods 0.000 description 2
- 241000544061 Cuculus canorus Species 0.000 description 1
- 241000255581 Drosophila <fruit fly, genus> Species 0.000 description 1
- 230000003044 adaptive effect Effects 0.000 description 1
- 230000002776 aggregation Effects 0.000 description 1
- 238000004220 aggregation Methods 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 238000006243 chemical reaction Methods 0.000 description 1
- 238000000354 decomposition reaction Methods 0.000 description 1
- 238000009792 diffusion process Methods 0.000 description 1
- 239000006185 dispersion Substances 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000000605 extraction Methods 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 238000013507 mapping Methods 0.000 description 1
- 238000005295 random walk Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000010187 selection method Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
- G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
- G01S7/41—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section
- G01S7/418—Theoretical aspects
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/22—Matching criteria, e.g. proximity measures
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/23—Clustering techniques
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
- G06F2218/08—Feature extraction
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation
Landscapes
- Engineering & Computer Science (AREA)
- Data Mining & Analysis (AREA)
- Theoretical Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Evolutionary Biology (AREA)
- Evolutionary Computation (AREA)
- Bioinformatics & Computational Biology (AREA)
- General Engineering & Computer Science (AREA)
- Artificial Intelligence (AREA)
- Life Sciences & Earth Sciences (AREA)
- Computer Networks & Wireless Communication (AREA)
- Radar, Positioning & Navigation (AREA)
- Remote Sensing (AREA)
- Radar Systems Or Details Thereof (AREA)
Abstract
本发明公开了一种基于改进社区合并的雷达信号密度峰值聚类方法，主要解决现有技术初始聚类误差较大和在密度峰值聚类中通过人为确定参数实用性差的问题。其方案为：仿真生成多种雷达脉冲信号的交错脉冲描述字，并根据其生成邻接矩阵；由邻接矩阵和欧氏距离生成相似度并进行社区合并，得到初聚类标签，用其相同标签组成初始簇；依次计算簇间距离和基尼系数，确定截止距离，并根据这两个参数计算簇权值指数和，确定聚类待选中心簇；从待选中心簇中筛选聚类中心簇，将非中心簇合并到距离最近的中心簇，得到最终聚类标签。本发明减少了初聚类造成的误差，能自适应地确定截止距离和聚类中心，保证分选的准确率，可用于多体制雷达脉冲信号的分选中。
Description
技术领域
本发明属于信号处理技术领域，特别涉及一种雷达信号密度峰值聚类方法，可应用于对多体制雷达的信号分选。
背景技术
多体制雷达的使用越来越普遍，这使得侦察接收机对未知雷达辐射源信号分选的难度系数变高。在电子侦察领域，未知雷达辐射源的信号分选是核心技术。对未知雷达辐射源信号的有效分选，对雷达辐射源信号进行下一步分析有重要意义。为了增强雷达信号的抗干扰能力，雷达辐射源信号具有调制方式复杂多变，特征参数的动态范围大以及组合方式多的特点，传统的信号分选技术的已不能满足识别率的要求。因此，深入研究雷达信号分选技术，探索信号分选的新方法，具有迫切的必要和重大的价值。
基于聚类的信号分选算法在辐射源分选领域得到了广泛的应用，一系列针对聚类参数确定、数据簇形状识别、算法速度提升的改进方法被相继提出。在复杂电磁环境中，利用聚类算法对先验知识无要求的优势可以有效地对未知雷达辐射源预分选。
K.Chi等人在其发表的论文“基于复杂网络的多功能雷达信号分选”(系IEEE信号处理快报，第一卷，28，第91-95页，2021年)中提出了一种基于复杂网络的多功能雷达信号分选方法。其通过改进有限穿越可视图算法将时间序列转换为复杂网络，得到邻接矩阵。然后改进标签传播算法，得到初始聚类标签，最后使用密度峰值聚类，对子网络进行合并聚类。该方法存在以下不足：
1、初始聚类形成的子网络个数和正确簇的个数有较大差距，影响后期合并聚类效果，使待合并簇内部的错误节点增加，聚类的总准确率下降；
2、在密度峰值聚类中需要人工设置截止距离和聚类中心数目，设置的截止距离不能统一，易出现过小或过大，造成分选的准确率下降；由于聚类中心的数目在实际分选中是未知的，这种通过人为确定方式缺乏实用性。
发明内容
本发明的目的在于针对上述现有技术存在的不足，提出一种基于改进社区合并的雷达信号密度峰值聚类方法，以根据节点间的距离和时间关系，增加复杂网络社区检测初聚类的正确簇个数，并以密度峰值聚类为基础，自适应地确定截止距离和聚类中心，保证分选的准确率，且更符合实际应用情况。
为实现上述目的，本发明基于社区合并的雷达信号密度峰值聚类方法，其实现方案包括如下：
1)通过仿真生成大小为2400×4的交错脉冲描述字PDW，其含有2400个脉冲信号，每一个脉冲信号的四维特征分别是频率RF、脉宽PW、到达时间TOA和到达角DOA；
2)利用脉冲描述字PDW生成邻接矩阵：
2a)将交错脉冲描述字PDW归一化预处理，将每一维特征归一化到0～1范围内；
2b)设滑动窗长度W为20，有限穿越个数N为6，滑动步长为1，将所有脉冲信号的频率RF、脉宽PW、到达角DOA这三个按时间排序的特征序列[rf1,…,rfj,…,rf2400]T、[pw1,…,pwj,…,pw2400]T和[doa1,…,doaj,…,doa2400]T，分别通过有限穿越可视图方法得到三个邻接矩阵A1,A2,A3，并将三个矩阵取并集得到融合邻接矩阵A＝A1∪A2∪A3，其中j＝1,2,…,2400；
3)把每一个脉冲视为一个节点，并将单个节点作为一个社区，设置社区的初始标签为序列号，计算所有节点间的欧氏距离dist，由邻接矩阵A和欧氏距离dist计算节点间的相似度s；
4)获得初聚类的节点标签：
4a)根据节点间相似度s计算社区间平均相似度，将平均相似度最大的两个社区进行合并，将平均相似度最大的社区标签作为合并后社区内节点的标签；
4b)重复4a)，直到所有节点合并成为c个社区，得到初聚类的节点标签；
5)根据初聚类的节点得到二次聚类标签：
5a)用初聚类中的相同标签组成一个初始簇，由节点间的欧氏距离得到各初始簇间的距离，并将其按升序排列，取前1％～5％个簇间距离组成待选簇间距离，由待选簇间距离计算基尼指数，选择该指数最小时的距离作为截止距离dc；
5b)根据簇间距离和截止距离dc计算每一个簇的局部密度和距离，再将这二者的归一化值对应相乘得到c个簇权值，对这些簇权值进行降序排列，结合簇权值和特征维数计算每个簇权值指数和Dk，再根据肘部法则寻找簇权值指数和Dk的上升率变缓的拐点的横坐标值，将其作为初始聚类中心的数目pmax；
5c)对于排序后簇权值的前pmax个候选中心簇，排除簇间距离小于截止距离dc的簇，将剩余的簇定为中心簇，最后将非中心簇合并到距离最近的中心簇，得到二次聚类的标签，完成对雷达脉冲信号的聚类。
本发明具有如下优点：
第一、采用节点间相似度来进行社区合并，相对于传统社区检测合并算法中利用边介数和模块度来计算，相似度的计算复杂度降低，且相似度结合了节点间欧式距离和邻接矩阵，将数据的节点特征和邻居关系联合考虑，更全面地表现节点间的关系，进而保证正确簇的个数尽量多，减少子网络个数和正确簇的个数的差距。
第二、采用基尼指数，选择使基尼指数最小的截止距离，可以令数据不确定性降到最低，使与截止距离相关的聚类中心的选取的结果更可靠，可避免因人工选取的不适宜而导致的中心个数选取错误，保证分选的准确率。
第三、采用自适应的聚类中心个数选取方法，通过获取簇权值指数和拐点来确定聚类中心个数，避免了现有技术因为缺少先验知识导致无法设定聚类中心个数的问题，使其更符合实际应用情况。
附图说明
图1是本发明的实现总流程图；
图2是本发明中仿真生成交错脉冲描述字的子流程图；
图3是本发明中生成邻接矩阵的子流程图。
具体实施方式
以下结合附图对本发明的实施例进行详细说明。
参照图1，本发明的基于改进社区合并的雷达信号密度峰值聚类方法，其实现步骤如下：
步骤1：生成雷达信号交错的脉冲描述字PDW数据集。
雷达信号交错的脉冲描述字PDW数据集，可以从雷达侦察机接收获取，也可以通过仿真生成，本实例通过仿真生成大小为2400×4的交错脉冲描述字PDW，其实现如下：
参照图2，本步骤的具体实现如下：
1.1)设置仿真条件：
选取四部雷达，并设置其脉冲信号个数及调制方式和参数：
第一部雷达的脉冲信号为500个，脉冲重复周期PRIrd1通过滑变调制为500us～1000us，频率RFrd1通过脉组捷变调制为8000MHz～8400MHz，脉宽PWrd1固定为200us，到达角DOArd1为35°～38°；
第二部雷达的脉冲信号的到达角DOArd2为36°～39°，且分为3种模式：
第1种模式的脉冲信号有200个，脉冲重复周期PRIrd2通过抖动调制为450us～550us，频率RFrd2是参差调制的7900MHz～8300MHz，脉宽PWrd2通过捷变调制为90us～110us；
第2种模式的脉冲信号有300个，脉冲重复周期PRIrd2通过脉组捷变为500us～600us，频率RFrd2通过参差调制为8300MHz～8700MHz，脉宽PWrd2为固定120us；
第3种模式的脉冲信号有200个，脉冲重复周期PRIrd2通过参差调制为600us～700us，频率RFrd2通过脉组捷变调制为8800MHz～9100MHz，脉宽PWrd2通过参差调制为130us～150us；
第三部雷达的脉冲信号为500个，脉冲重复周期PRIrd3通过脉组捷变调制为650us～800us，频率RFrd3通过参差调制为8200MHz～8600MHz，脉宽PWrd3为固定的100us，到达角DOArd3为37°～40°；
第四部雷达脉冲信号的到达角DOArd4为38°～41°，且分为三种模式：
第一种模式的脉冲信号有200个，脉冲重复周期PRIrd4通过抖动调制为700us～800us，频率RFPRIrd4通过参差调制为8400MHz～8800MHz，脉宽PWrd4通过捷变调制为140us～160us；
第二种模式的脉冲信号有200个，脉冲重复周期PRIrd4通过脉组捷变调制为800us～900us，频率RFrd4通过参差调制为8800MHz～9200MHz，脉宽PWrd4通过参差调制为160us～180us；
第三种模式的脉冲信号有300个，脉冲重复周期PRIrd4通过参差调制为900us～1000us，频率RFrd4通过脉组捷变调制为9200MHz～9400MHz，脉宽PWrd4通过参差调制为180us～200us；
以上仿真条件如表1所示：
表1四部雷达的脉冲信号个数及调制方式和参数
1.2)根据上述参数，得到每一部雷达脉冲信号的脉冲描述字PDWrdi,i＝1,2,3,4，即：
第一部雷达脉冲信号的脉冲描述字：PDWrd1＝{RFrd1,PWrd1,TOArd1,DOArd1}，
第二部雷达脉冲信号的脉冲描述字：PDWrd2＝{RFrd2,PWrd2,TOArd2,DOArd2},
第三部雷达脉冲信号的脉冲描述字：PDWrd3＝{RFrd3,PWrd3,TOArd3,DOArd3}，
第四部雷达脉冲信号的脉冲描述字：PDWrd4＝{RFrd4,PWrd4,TOArd4,DOArd4}；
其中，RFrdi＝{RFrdi(1),...,RFrdi(li)},
PWrdi＝{PWrdi(1),...,PWrdi(li)}，
TOArdi＝{TOArdi(1),...,TOArdi(li)},
DOArdi＝{DOArdi(1),...,DOArdi(li)}，
li表示第i个雷达的脉冲个数；
1.3)将所有脉冲信号的到达时间[TOArd1,TOArd2,TOArd3,TOArd4]T按升序排列，获得排序后的到达时间TOA＝[TOA(1),…,TOA(j),…,TOA(2400)]T，其中：j＝1,2,…,2400；
1.4)将四部雷达脉冲信号的脉冲描述字的频率、脉宽和到达角分别合并为：
[RFrd1,RFrd2,RFrd3,RFrd4]T，
[PWrd1,PWrd2,PWrd3,PWrd4]T，
[DOArd1,DOArd2,DOArd3,DOArd4]T；
1.5)按排序后脉冲信号到达时间TOA的顺序排序，得到排序后的脉冲信号的频率RF、脉宽PW和到达角DOA：
RF＝[RF(1),…,RF(j),…,RF(2400)]T，
PW＝[PW(1),…,PW(j),…,PW(2400)]T，
DOA＝[DOA(1),…,DOA(j),…,DOA(2400)]T，
1.6)将排序后的频率RF、脉宽PW、到达时间TOA和到达角DOA共四个序列按行方向合并，得到交错的脉冲描述字PDW为：
步骤2：利用脉冲描述字PDW生成邻接矩阵。
参照图3，本步骤的实现如下：
2.1)对交错脉冲描述字PDW进行归一化预处理：
数据的归一化预处理一般是将数据映射到指定的范围，用于去除不同维度数据的量纲以及量纲单位。常见的映射范围有[0,1]和[-1,1]，常见的归一化方法就是区间放缩法，本实例采用但不限于离差标准化Min-Max的归一化方法，即将每一维特征频率RF、脉宽PW、到达时间TOA和到达角DOA归一化到[0,1]范围内，计算公式如下：
其中，rfj、pwj、toaj和doaj分别表示每一列特征序列归一化后的第j个值，RF(j)、PW(j)、TOA(j)和DOA(j)分别是交错的脉冲描述字PDW中的每一列特征序列的第j个值，min(·)和max(·)表示·的最小值和最大值；
由归一化的频率特征值rfj构成频率RF的归一化序列[rf1,…,rfj,…,rf2400]T，
由归一化的脉宽特征值pwj构成脉宽PW的归一化序列[pw1,…,pwj,…,pw2400]T,
由归一化的到达角特征值doaj构成到达角DOA的归一化序列[doa1,…,doaj,…,doa2400]T；
2.2)生成邻接矩阵：
时间序列的邻接矩阵生成，就是将时间序列转化为复杂网络，网络的表示是邻接矩阵，常见的生成邻接矩阵的方法有可视图算法，基于频域复杂网络分解的局部特征提取的方法，基于随机图方法等，本实例采用但不限于有限穿越可视图方法生成邻接矩阵，其实现如下：
2.2.1)取交错脉冲描述字PDW归一化后的3维特征
的第D列，其中D＝1,2,3，初始化设置一个2400×2400的零矩阵作为统计矩阵count，使用宽为W＝20的滑窗沿着第D维特征序列，以步长为1的间隔开始滑动；
2.2.2)滑动到第1个滑窗内时，判断滑窗内的第k点是否满足如下条件：
若满足该条件，则令统计矩阵元素count(m,n)增加1，
否则，count(m,n)不增加，
其中，m≥1，m＜k＜n，n＝2,...,W，x(m)、x(n)和x(k)分别表示第D列特征的第m个、第n个和第k个值；
2.2.3)重复2.2.2)，对点m到点n间所有的点都进行判断，得到统计矩阵元素count(m,n)＝h，h为点m到点n间符合上述判断条件的点数；
2.2.4)滑动到第g个滑窗内时，判断滑窗内的k点是否满足如下条件：
若满足该条件，则将统计矩阵元素count(m,g+W-1)增加1，
否则，count(m,g+W-1)不增加，
其中，m≥g，m＜k＜g+W-1，x(m)、x(g+W-1)和x(k)分别表示第D列特征的第m个、第g+W-1个和第k个值；
2.2.5)重复2.2.4)，对点m到点g+W-1间所有的点都进行判断，得到统计矩阵元素count(m,g+W-1)＝u，u为点g+W-1到点间符合上述判断条件的点数；
2.2.6)重复2.2.5)，滑动滑窗从g＝2直至g＝2400-W+1滑窗结束，用第D维特征序列的点m和点n间符合判断条件2.2.2)或判断条件2.2.4)的点数构成统计矩阵元素count(m,n)，m,n＝1,2,…,2400,
根据统计矩阵count和有限穿越数，定义第D维特征的点m和点n的邻接关系
其中N＝6为有限穿越数，
表示第D维特征的邻接矩阵AD的第m行第n列的元素，当
时，表示点m和点n是邻接的，且
2.2.7)重复2.2.1)至2.2.6)，直到完成D＝1,2,3维特征的全部判完，根据D维特征序列的邻接关系
得到三个对称的邻接矩阵A1、A2和A3表示如下：
步骤3：计算所有节点间的欧氏距离dist，由邻接矩阵A和欧氏距离dist计算节点间的相似度s。
3.1)计算所有节点间的欧氏距离dist。
描述两个点的距离，可以有多种方法，例如欧氏距离、马氏距离、余弦距离、汉明距离、曼哈顿距离等,本实例采用但不限于使用欧氏距离描述两点间距离，其公式如下：
其中，dist(m,n)表示节点m和节点n之间的欧氏距离，
rfm和rfn分别属于归一化频率[rf1,…,rfj,…,rf2400]T中的第m和第n个值，
pwm和pwn分别属于归一化脉宽[pw1,…,pwj,…,pw2400]T中的第m和第n个值，
pwm和pwn分别属于归一化到达角[doa1,…,doaj,…,doa2400]T中的第m和第n个值；
3.2)由邻接矩阵A和欧氏距离dist计算节点间的相似度s。
相似度越大，表示两个元素越相似；相似度越小，表示两个元素越不相似。可以使用多种方法如余弦相似度、皮尔逊相似度、Jaccard相似度等来描述两点之间的相似度，本实例采用但不限于使用基于邻接矩阵和欧氏距离描述节点间的相似度，其公式如下：
其中，s(m,n)表示节点m和节点n之间的相似度，dist(m,n)表示节点m和节点n之间的欧氏距离，当两节点的距离dist(m,n)越大，则相似度越小，当两节点不相邻时，即A(m,n)＝0时，相似度为0。
步骤4：获得初聚类的节点标签。
常见的社区聚类算法包括GN算法、Newman快速算法、鲁文社区检测算法、随机游走算法等，本实例采用但不限于使用基于改进GN算法进行初始聚类获得节点的初聚类标签，其实现如下：
4.1)根据节点间相似度s计算社区间平均相似度s_ave(cn,cm)：
其中cm和cn分别表示第cm和第cn个社区序号，|cm|和|cn|分别表示属于第cm和第cn个社区的节点个数；
4.2)将平均相似度矩阵中最大值s_ave(cmmax,cnmax)＝max(s_ave)对应的社区cmmax和社区cnmax进行合并，将社区cnmax内所有节点的标签也合并成与社区cmmax内节点相同的标签；
4.3)重复4.1)和4.2)，直到所有节点合并成为c个社区，得到初聚类的节点标签。
步骤5：获得截止距离。
现今为止，已有多种针对截止距离这一参数进行改进的密度峰值聚类算法，例如：结合K近邻的概念重新定义截止距离，利用果蝇优化算法得到最优截止距离，还有基于布谷鸟优化、基于热扩散、基于自适应聚合策略优化的密度峰值聚类算法等，本实例采用基于最小不确定度获取截止距离的方法，其实现如下：
5.1)用初聚类中的相同标签组成一个初始簇，将分别属于簇ci和簇cj的两节点间的欧式距离dist(m,n)求和再平均得到各初始簇的距离：
其中,d_averaged(ci,cj)表示初始簇ci和簇cj的之间的距离，dist(m,n)表示节点m和节点n之间的欧氏距离，|ci|和|cj|分别表示属于初始聚类得到的簇ci和簇cj的节点个数，ci,cj＝1,2,...,c，c为初始聚类得到的社区个数；
5.2)将初始簇的距离阵d_averaged中的所有值按升序排列，取前1％～5％个簇间距离组成待选簇间距离；
5.3)由待选簇间距离计算基尼系数：
5.3.1)计算待选截止距离dcnd对应的局部密度ρnd和距离δnd：
其中，nd＝1,2,...,n，n是前1％～5％个簇间距离的个数，dcnd是这些待选簇间距离的第nd个，ρnd(ci)和距离δnd(ci)分别表示由待选截至距离dcnd计算得到的第ci个簇的局部密度和距离，d_averaged(ci,cj)表示初始簇ci与簇cj的之间的距离，d_averaged(ci,:)表示初始簇ci与其他所有簇的之间的距离；
5.3.2)将第nd个待选簇间距离dcnd所对应的局部密度ρnd和距离δnd归一化，得到归一化后的局部密度ρ'nd和距离δ'nd：
其中min(ρnd)和max(ρnd)分别表示局部密度ρnd的最小值和最大值，min(δnd)和max(δnd)分别表示距离δnd的最小值和最大值；
5.3.3)根据归一化的局部密度ρ'nd和距离δ'nd计算基尼系数：
其中γ(ci)＝ρ'nd(ci)·δ'nd 2(ci)表示归一化的局部密度ρ'nd和距离δ'nd两者融合的第ci个值，
为所有融合值γ(ci)之和，基尼系数越小，则数据的不确定性越小。
5.4)选择基尼系数最小时的距离作为截止距离dc。
步骤6：获得待选中心簇。
密度峰值聚类原算法通过“决策图”人工判断类中心点，针对聚类中心的选择，可以引入K近邻计算全局参数和局部密度，判断密度可达点为聚类中心，或者通过模糊聚类找局部聚类中心，再或者根据决策函数进行中心点的自动选择等，本实例采用寻找定义的簇权值指数和的拐点的方法来获取待选的中心簇，其实现如下：
6.1)根据簇间距离和截止距离dc计算每一个簇的局部密度和距离，计算公式如下：
其中，ρ(ci)和δ(ci)分别表示选定的截至距离dc对应的第ci个簇的局部密度和距离，d_averaged(ci,cj)表示初始簇ci与簇cj的之间的距离，d_averaged(ci,:)表示初始簇ci与其他所有簇的之间的距离；
6.2)将每一个相同簇归一化的局部密度ρ和δ相乘得到c个簇权值，公式如下：
其中，r(ci)为第ci个簇权值，ρ(ci)和δ(ci)分别表示第ci个簇的局部密度和距离，将这些簇权值r进行降序排列获取排序后的簇权值R；
6.3)结合簇权值和特征维数计算每个簇权值指数和Dk，公式如下：
其中，Dk(ci)表示第ci个簇的权值指数，R(ci)和R(cj)分别为第ci和第cj个归一化的簇权值，c为初始簇的个数，d为3表示特征的维数；
6.4)根据肘部法则寻找簇权值指数和Dk的上升率变缓的拐点的横坐标值，并将其作为初始聚类中心的数目pmax：
6.4.1)选择簇权值的前4％个，计算每一个簇权值指数和Dk的梯度差ddk(p)：
ddk(p)＝(Dk(p+1)-Dk(p))-(Dk(p)-Dk(p-1))
其中p＝2,…,NP-1，NP为选择的簇权值个数，Dk(p)为第p个簇权值指数和；
6.4.2)记录梯度差ddk最大时的簇pmax，该pmax即为簇权值指数和Dk上升率变缓的拐点的横坐标值；
6.5)将pmax作为待选中心簇的个数，取排序后的归一化簇权值R的前pmax个值R(1),R(2),…R(pmax)对应的簇作为待选中心簇。
步骤7：获得二次聚类标签。
7.1)对于R(1),R(2),…R(pmax)这pmax个簇权值对应的待选中心簇，首先排除簇权值R(2),…R(pmax)对应簇中与最大簇权值R(1)对应簇的簇间距离小于截止距离dc的簇，然后排除剩余的待选中心簇中与次最大簇权值对应簇的簇间距离小于截止距离dc的簇，以此类推，将所有簇间距离小于截止距离，且簇权值较小的簇都排除后，剩余的待选簇即为中心簇；
7.2)将非中心簇合并到簇间距离最小的中心簇，得到二次聚类的标签，完成对雷达脉冲信号的聚类。
以上描述仅是本发明的一个具体实例，并未构成对本发明的任何限制，显然对于本领域的专业人员来说，在了解本发明的内容和原理后，都可能在不背离发明原理、结构的情况下，进行形式和细节上的修改和改变，但是这些基于本发明思想的修正和改变仍在本发明的权利要求保护范围之内。
Claims (10)
1.一种基于社区合并的雷达信号密度峰值聚类方法，其特征在于，包括：
1)通过仿真生成大小为2400×4的交错脉冲描述字PDW，其含有2400个脉冲信号，每一个脉冲信号的四维特征分别是频率RF、脉宽PW、到达时间TOA和到达角DOA；
2)利用脉冲描述字PDW生成邻接矩阵：
2a)将交错脉冲描述字PDW归一化预处理，将每一维特征归一化到0～1范围内；
2b)设滑动窗长度W为20，有限穿越个数N为6，滑动步长为1，将所有脉冲信号的频率RF、脉宽PW、到达角DOA这三个按时间排序的特征序列[rf1,…,rfj,…,rf2400]T、[pw1,…,pwj,…,pw2400]T和[doa1,…,doaj,…,doa2400]T，分别通过有限穿越可视图方法得到三个邻接矩阵A1,A2,A3，并将三个矩阵取并集得到融合邻接矩阵A＝A1∪A2∪A3，其中j＝1,2,…,2400；
3)把每一个脉冲视为一个节点，并将单个节点作为一个社区，设置社区的初始标签为序列号，计算所有节点间的欧氏距离dist，由邻接矩阵A和欧氏距离dist计算节点间的相似度s；
4)获得初聚类的节点标签：
4a)根据节点间相似度s计算社区间平均相似度，将平均相似度最大的两个社区进行合并，将平均相似度最大的社区标签作为合并后社区内节点的标签；
4b)重复4a)，直到所有节点合并成为c个社区，得到初聚类的节点标签；
5)根据初聚类的节点得到二次聚类标签：
5a)用初聚类中的相同标签组成一个初始簇，由节点间的欧氏距离得到各初始簇间的距离，并将其按升序排列，取前1％～5％个簇间距离组成待选簇间距离，由待选簇间距离计算基尼系数，选择该系数最小时的距离作为截止距离dc；
5b)根据簇间距离和截止距离dc计算每一个簇的局部密度和距离，再将这二者对应相乘得到c个簇权值，对这些簇权值进行降序排列，结合簇权值和特征维数计算每个簇权值指数和Dk，再根据肘部法则寻找簇权值指数和Dk的上升率变缓的拐点的横坐标值，将其作为初始聚类中心的数目pmax；
5c)对于排序后簇权值的前pmax个候选中心簇，排除簇间距离小于截止距离dc的簇，将剩余的簇定为中心簇，最后将非中心簇合并到距离最近的中心簇，得到二次聚类的标签，完成对雷达脉冲信号的聚类。
2.根据权利要求1所述的方法，其特征在于，1)中通过仿真生成大小为2400×4的交错脉冲描述字PDW，实现如下：
1.1)设置仿真条件：
选取四部雷达，并设置其脉冲信号的调制方式和参数：
第一部雷达的脉冲信号为500个，脉冲重复周期PRIrd1通过滑变调制为500us～1000us,频率RFrd1通过脉组捷变调制为8000MHz～8400MHz，脉宽PWrd1固定为200us，到达角DOArd1为35°～38°；
第二部雷达的脉冲信号的到达角DOArd2为36°～39°，且分为3种模式，第一种模式的脉冲信号有200个，脉冲重复周期PRIrd2通过抖动调制为450us～550us，频率RFrd2是参差调制的7900MHz～8300MHz，脉宽PWrd2通过捷变调制为90us～110us；第二种模式的脉冲信号有300个，脉冲重复周期PRIrd2通过脉组捷变为500us～600us，频率RFrd2通过参差调制为8300MHz～8700MHz，脉宽PWrd2为固定120us；第三种模式的脉冲信号有200个，脉冲重复周期PRIrd2通过参差调制为600us～700us，频率RFrd2通过脉组捷变调制为8800MHz～9100MHz，脉宽PWrd2通过参差调制为130us～150us；
第三部雷达的脉冲信号为500个，脉冲重复周期PRIrd3通过脉组捷变调制为650us～800us，频率RFrd3通过参差调制为8200MHz～8600MHz，脉宽PWrd3为固定的100us，到达角DOArd3为37°～40°；
第四部雷达脉冲信号的到达角DOArd4为38°～41°，且分为3种模式，第一种模式的脉冲信号有200个，脉冲重复周期PRIrd4通过抖动调制为700us～800us，频率RFPRIrd4通过参差调制为8400MHz～8800MHz，脉宽PWrd4通过捷变调制为140us～160us；第二种模式的脉冲信号有200个，脉冲重复周期PRIrd4通过脉组捷变调制为800us～900us，频率RFrd4通过参差调制为8800MHz～9200MHz，脉宽PWrd4通过参差调制为160us～180us；第三种模式的脉冲信号有300个，脉冲重复周期PRIrd4通过参差调制为900us～1000us，频率RFrd4通过脉组捷变调制为9200MHz～9400MHz，脉宽PWrd4通过参差调制为180us～200us。
1.2)根据上述参数，得到每一部雷达脉冲信号的脉冲描述字PDWrdi,i＝1,2,3,4为：
PDWrd1＝{RFrd1,PWrd1,TOArd1,DOArd1}，
PDWrd2＝{RFrd2,PWrd2,TOArd2,DOArd2},
PDWrd3＝{RFrd3,PWrd3,TOArd3,DOArd3}，
PDWrd4＝{RFrd4,PWrd4,TOArd4,DOArd4}；
其中，RFrdi＝{RFrdi(1),...,RFrdi(li)},PWrdi＝{PWrdi(1),...,PWrdi(li)}，TOArdi＝{TOArdi(1),...,TOArdi(li)},DOArdi＝{DOArdi(1),...,DOArdi(li)}，li表示第i个雷达的脉冲个数；
1.3)将所有脉冲信号的到达时间[TOArd1,TOArd2,TOArd3,TOArd4]T按升序排列，获得排序后的到达时间TOA＝[TOA(1),…,TOA(j),…,TOA(2400)]T，其中j＝1,2,…,2400；
1.4)将四部雷达脉冲信号的脉冲描述字的频率、脉宽和到达角分别合并为：
[RFrd1,RFrd2,RFrd3,RFrd4]T，
[PWrd1,PWrd2,PWrd3,PWrd4]T，
[DOArd1,DOArd2,DOArd3,DOArd4]T；
1.5)按排序后脉冲信号到达时间TOA的顺序排序，得到排序后的脉冲信号的频率RF、脉宽PW和到达角DOA：
RF＝[RF(1),…,RF(j),…,RF(2400)]T，
PW＝[PW(1),…,PW(j),…,PW(2400)]T，
DOA＝[DOA(1),…,DOA(j),…,DOA(2400)]T，
1.6)将排序后的频率RF、脉宽PW、到达时间TOA和到达角DOA共四个序列按行方向合并，得到交错的脉冲描述字PDW为：
3.根据权利要求1所述的方法，其特征在于2b)中通过有限穿越可视图方法得到三个邻接矩阵A1,A2,A3，实现表示如下：
2b1)取交错脉冲描述字PDW归一化后的3维特征
的第D列，其中D＝1,2,3，初始化设置一个2400×2400的零矩阵作为统计矩阵count，使用宽为W＝20的滑窗沿着第D维特征序列，以步长为1的间隔开始滑动；
2b2)滑动到第1个滑窗内时，判断滑窗内的第k点是否满足如下条件：
若满足该条件，则令统计矩阵元素count(m,n)增加1，否则，count(m,n)不增加，其中，m≥1，m＜k＜n，n＝2,...,W，x(m)、x(n)和x(k)分别表示第D列特征的第m个、第n个和第k个值；
2b3)重复2b2)，对点m到点n间所有的点都进行判断，得到统计矩阵元素count(m,n)＝h，h为点m到点n间符合上述判断条件的点数；
2b4)滑动到第g个滑窗内时，判断滑窗内的k点是否满足如下条件：
若满足该条件，则令统计矩阵元素count(m,g+W-1)增加1，否则，count(m,g+W-1)不增加，其中，m≥g，m＜k＜g+W-1，x(m)、x(g+W-1)和x(k)分别表示第D列特征的第m个、第g+W-1个和第k个值；
2b5)重复2b4)，对点m到点g+W-1间所有的点都进行判断，得到统计矩阵元素count(m,g+W-1)＝u，u为点g+W-1到点间符合上述判断条件的点数；
2b6)重复2b5)，滑动滑窗从g＝2直至g＝2400-W+1滑窗结束，用第D维特征序列的点m和点n间符合判断条件2b2)或判断条件2b4)的点数构成统计矩阵元素count(m,n)，m,n＝1,2,…,2400,根据统计矩阵count和有限穿越数，定义第D维特征的点m和点n的邻接关系
其中N＝6为有限穿越数，
表示第D维特征的邻接矩阵AD的第m行第n列的元素，当
时，表示点m和点n是邻接的，且
2b7)重复2b1)至2b6)，直到完成D＝1,2,3维特征的全部判完，根据D维特征序列的邻接关系
即可得到三个对称的邻接矩阵A1、A2和A3表示如下：
4.根据权利要求1所述的方法，其特征在于3)中计算所有节点间的欧氏距离dist，公式如下：
其中，dist(m,n)表示节点m和节点n之间的欧氏距离，
rfm和rfn分别属于归一化频率[rf1,…,rfj,…,rf2400]T中的第m和第n个值，
pwm和pwn分别属于归一化脉宽[pw1,…,pwj,…,pw2400]T中的第m和第n个值，
pwm和pwn分别属于归一化到达角[doa1,…,doaj,…,doa2400]T中的第m和第n个值。
5.根据权利要求1所述的方法，其特征在于3)中由邻接矩阵A和欧氏距离dist计算节点间的相似度s，公式如下：
其中，s(m,n)表示节点m和节点n之间的相似度，dist(m,n)表示节点m和节点n之间的欧氏距离，当两节点的距离dist(m,n)越大，则相似度越小，当两节点不相邻时，即A(m,n)＝0时，相似度为0。
6.根据权利要求1所述的方法，其特征在于5a)中由节点间的欧氏距离得到各初始簇间的距离，计算公式如下：
将分别属于簇ci和簇cj的两节点间的距离dist(m,n)求和再平均得到的各初始簇的距离如下：
其中,d_averaged(ci,cj)表示初始簇ci和簇cj的之间的距离，dist(m,n)表示节点m和节点n之间的欧氏距离，|ci|和|cj|分别表示属于初始聚类得到的簇ci和簇cj的节点个数，ci,cj＝1,2,...,c，c为初始聚类得到的社区个数。
7.根据权利要求1所述的方法，其特征在于5a)中由待选簇间距离计算基尼系数，计算如下：
5a1)计算待选截止距离dcnd对应的局部密度ρnd和距离δnd：
其中，nd＝1,2,...,n，n是前1％～5％个簇间距离的个数，dcnd是这些待选簇间距离的第nd个，ρnd(ci)和距离δnd(ci)分别表示由待选截至距离dcnd计算得到的第ci个簇的局部密度和距离，d_averaged(ci,cj)表示初始簇ci与簇cj的之间的距离，d_averaged(ci,:)表示初始簇ci与其他所有簇的之间的距离；
5a2)将第nd个待选簇间距离dcnd所对应的局部密度ρnd和距离δnd归一化，得到归一化后的局部密度ρ'nd和距离δ'nd：
其中min(ρnd)和max(ρnd)分别表示局部密度ρnd的最小值和最大值，min(δnd)和max(δnd)分别表示距离δnd的最小值和最大值；
5a3)根据归一化的局部密度ρ'nd和距离δ'nd计算基尼系数：
其中γ(ci)＝ρ'nd(ci)·δ'nd 2(ci)表示归一化的局部密度ρ'nd和距离δ'nd两者融合的第ci个值，
为所有融合值γ(ci)之和，基尼系数越小，则数据的不确定性越小。
8.根据权利要求1所述的方法，其特征在于5b)中根据簇间距离和截止距离dc计算每一个簇的局部密度和距离，计算公式如下：
其中，ρ(ci)和δ(ci)分别表示选定的截至距离dc对应的第ci个簇的局部密度和距离，d_averaged(ci,cj)表示初始簇ci与簇cj的之间的距离，d_averaged(ci,:)表示初始簇ci与其他所有簇的之间的距离。
9.根据权利要求1所述的方法，其特征在于5b)中结合簇权值和特征维数计算每个簇权值指数和Dk，公式如下：
其中，Dk(ci)表示第ci个簇的权值指数，R(ci)和R(cj)分别为第ci和第cj个归一化的簇权值，c为初始簇的个数，d为3表示特征的维数。
10.根据权利要求1所述的方法，其特征在于5b)中根据肘部法则寻找簇权值指数和Dk的上升率变缓的拐点的横坐标值，实现如下：
5b1)选择簇权值的前4％个，计算每一个簇权值指数和Dk的梯度差ddk(p)：
ddk(p)＝(Dk(p+1)-Dk(p))-(Dk(p)-Dk(p-1))
其中p＝2,…,NP-1，NP为选择的簇权值个数，Dk(p)为第p个簇权值指数和；
5b2)记录梯度差ddk最大时的簇pmax，该pmax即为簇权值指数和Dk上升率变缓的拐点的横坐标值。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202111297714.XA CN114004259B (zh) | 2021-11-04 | 2021-11-04 | 基于改进社区合并的雷达信号密度峰值聚类方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202111297714.XA CN114004259B (zh) | 2021-11-04 | 2021-11-04 | 基于改进社区合并的雷达信号密度峰值聚类方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN114004259A true CN114004259A (zh) | 2022-02-01 |
| CN114004259B CN114004259B (zh) | 2024-09-27 |
Family
ID=79927031
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202111297714.XA Active CN114004259B (zh) | 2021-11-04 | 2021-11-04 | 基于改进社区合并的雷达信号密度峰值聚类方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN114004259B (zh) |
Cited By (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115659162A (zh) * | 2022-09-15 | 2023-01-31 | 云南财经大学 | 雷达辐射源信号脉内特征提取方法、系统及设备 |
| CN116361678A (zh) * | 2023-05-26 | 2023-06-30 | 西南石油大学 | 基于图增强结构的准周期时间序列分割方法及终端 |
| CN116489741A (zh) * | 2023-06-19 | 2023-07-25 | 南昌大学 | 一种无线传感网络跨层节能方法及系统 |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108596230A (zh) * | 2018-04-13 | 2018-09-28 | 济南大学 | 一种参数自适应的密度峰值聚类方法 |
| WO2019216469A1 (ko) * | 2018-05-11 | 2019-11-14 | 서울대학교 산학협력단 | 자동차 레이더 시스템에서 탐지된 타겟들의 클러스터링 방법 및 이를 위한 장치 |
| CN112036074A (zh) * | 2020-07-27 | 2020-12-04 | 西安电子科技大学 | 一种高脉冲密度环境下雷达信号分选方法及系统 |
-
2021
- 2021-11-04 CN CN202111297714.XA patent/CN114004259B/zh active Active
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108596230A (zh) * | 2018-04-13 | 2018-09-28 | 济南大学 | 一种参数自适应的密度峰值聚类方法 |
| WO2019216469A1 (ko) * | 2018-05-11 | 2019-11-14 | 서울대학교 산학협력단 | 자동차 레이더 시스템에서 탐지된 타겟들의 클러스터링 방법 및 이를 위한 장치 |
| CN112036074A (zh) * | 2020-07-27 | 2020-12-04 | 西安电子科技大学 | 一种高脉冲密度环境下雷达信号分选方法及系统 |
Non-Patent Citations (2)
| Title |
|---|
| PING LANG等: "Subspace Decomposition Based Adaptive Density Peak Clustering for Radar Signals Sorting", 《IEEE SIGNAL PROCESSING LETTERS》, 31 December 2021 (2021-12-31) * |
| 叶雨嘉: "对复杂体制雷达的信号分选技术研究", 《CNKI》, 1 June 2023 (2023-06-01) * |
Cited By (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115659162A (zh) * | 2022-09-15 | 2023-01-31 | 云南财经大学 | 雷达辐射源信号脉内特征提取方法、系统及设备 |
| CN115659162B (zh) * | 2022-09-15 | 2023-10-03 | 云南财经大学 | 雷达辐射源信号脉内特征提取方法、系统及设备 |
| CN116361678A (zh) * | 2023-05-26 | 2023-06-30 | 西南石油大学 | 基于图增强结构的准周期时间序列分割方法及终端 |
| CN116361678B (zh) * | 2023-05-26 | 2023-08-25 | 西南石油大学 | 基于图增强结构的准周期时间序列分割方法及终端 |
| CN116489741A (zh) * | 2023-06-19 | 2023-07-25 | 南昌大学 | 一种无线传感网络跨层节能方法及系统 |
| CN116489741B (zh) * | 2023-06-19 | 2023-09-01 | 南昌大学 | 一种无线传感网络跨层节能方法及系统 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN114004259B (zh) | 2024-09-27 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN114004259A (zh) | 基于改进社区合并的雷达信号密度峰值聚类方法 | |
| CN113378632B (zh) | 一种基于伪标签优化的无监督域适应行人重识别方法 | |
| Briggs et al. | Federated learning with hierarchical clustering of local updates to improve training on non-IID data | |
| CN107885778B (zh) | 一种基于动态临近点谱聚类的个性化推荐方法 | |
| CN113326731A (zh) | 一种基于动量网络指导的跨域行人重识别算法 | |
| Li et al. | Simultaneous localized feature selection and model detection for Gaussian mixtures | |
| CN114386466B (zh) | 一种用于脉冲星搜寻中候选体信号挖掘的并行的混合聚类方法 | |
| CN105405118A (zh) | 基于量子衍生混合蛙跳的水下声纳图像目标检测方法 | |
| CN113917527A (zh) | 基于多重量子神经网络进行含气性检测的方法 | |
| CN115965318A (zh) | 一种基于变中心演化聚类的物流中心选址方法 | |
| Wu et al. | An improved k-means algorithm based on density normalization | |
| CN115510959A (zh) | 基于自然最近邻和多簇合并的密度峰值聚类方法 | |
| CN113111774B (zh) | 一种基于主动增量式微调的雷达信号调制方式识别方法 | |
| CN112182287B (zh) | 一种基于时空视觉词组和分层匹配的视频拷贝检测方法 | |
| Tao et al. | Radar emitter signal recognition method based on improved collaborative semi-supervised learning | |
| CN116486594A (zh) | 一种基于k-邻近算法的车联网车辆id连续性追踪确定方法 | |
| CN111046914B (zh) | 基于动态构图的半监督分类方法 | |
| CN118194179B (zh) | 带噪音数据分类方法、系统及社交网络用户分类方法 | |
| Cui et al. | Radar signal sorting based on adaptive SOFM and Coyote optimization | |
| CN116881750B (zh) | 一种轨迹聚类方法 | |
| Lu et al. | A recursive sorting method based on multi-station TDOA vector | |
| Gong et al. | A point symmetry-based automatic clustering approach using differential evolution | |
| Gao et al. | Seafloor Sediment Classification using BP Neural Network | |
| YASAR et al. | FN-DBSCAN-GM: A parameter free and robust version of DBSCAN algorithm | |
| Alhanjouri et al. | Integrating bat algorithm to inverse weighted k-means |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |