CN105956526B - 基于多尺度排列熵的低信噪比微震事件辨识方法 - Google Patents
基于多尺度排列熵的低信噪比微震事件辨识方法 Download PDFInfo
- Publication number
- CN105956526B CN105956526B CN201610257713.5A CN201610257713A CN105956526B CN 105956526 B CN105956526 B CN 105956526B CN 201610257713 A CN201610257713 A CN 201610257713A CN 105956526 B CN105956526 B CN 105956526B
- Authority
- CN
- China
- Prior art keywords
- microseismic
- scale
- signal
- permutation entropy
- data
- 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.)
- Expired - Fee Related
Links
- 238000000034 method Methods 0.000 title claims abstract description 26
- 238000012843 least square support vector machine Methods 0.000 claims abstract description 12
- 238000012549 training Methods 0.000 claims abstract description 11
- 238000012545 processing Methods 0.000 claims abstract description 7
- 238000012706 support-vector machine Methods 0.000 claims description 14
- 238000010606 normalization Methods 0.000 claims description 3
- 238000002790 cross-validation Methods 0.000 claims description 2
- 239000013307 optical fiber Substances 0.000 claims description 2
- 230000009286 beneficial effect Effects 0.000 abstract 1
- 238000010586 diagram Methods 0.000 description 5
- 239000011435 rock Substances 0.000 description 4
- 239000003245 coal Substances 0.000 description 3
- 238000004364 calculation method Methods 0.000 description 2
- 238000005065 mining Methods 0.000 description 2
- 238000012544 monitoring process Methods 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 230000001174 ascending effect Effects 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 238000003745 diagnosis Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 238000000605 extraction Methods 0.000 description 1
- 230000002349 favourable effect Effects 0.000 description 1
- 230000010365 information processing Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000035772 mutation Effects 0.000 description 1
- 230000008092 positive effect Effects 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
- G06F2218/12—Classification; Matching
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/288—Event detection in seismic signals, e.g. microseismics
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
- G06F18/214—Generating training patterns; Bootstrap methods, e.g. bagging or boosting
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Data Mining & Analysis (AREA)
- Artificial Intelligence (AREA)
- General Engineering & Computer Science (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Remote Sensing (AREA)
- Environmental & Geological Engineering (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Business, Economics & Management (AREA)
- Geophysics (AREA)
- General Life Sciences & Earth Sciences (AREA)
- Bioinformatics & Computational Biology (AREA)
- Geology (AREA)
- Evolutionary Biology (AREA)
- Evolutionary Computation (AREA)
- Acoustics & Sound (AREA)
- Emergency Management (AREA)
- Signal Processing (AREA)
- Geophysics And Detection Of Objects (AREA)
- Investigating Or Analyzing Materials By The Use Of Ultrasonic Waves (AREA)
Abstract
本发明公开了一种基于多尺度排列熵的低信噪比微震事件辨识方法,它是在计算机上实现的,首先将采集的微震数据进行处理转换为时间序列数据,然后粗粒化时间序列得到多尺度时间序列,计算多尺度时间序列排列熵;基于多尺度排列熵训练最小二乘支持向量机,最后利用训练好的LS‑SVM对待判别信号进行辨识。本发明将多尺度和排列熵结合应用于微震信号可以更加有效地分析信号特征,准确表达微震信号多个维度的波形特征,有利于区分微震事件及噪声事件。应用多尺度排列熵分别提取微震信号及噪声信号的特征数据,使用LS‑SVM对特征数据进行训练,得到LS‑SVM向量机,实现了对低信噪比微震事件及噪声事件的准确分类。
Description
技术领域
本发明涉及一种低信噪比微震事件辨识方法,具体的说,是一种使用多尺度排列熵表示微震波形特点,进而通过支持向量机判别该波形是否属于微震事件的方法,属于信息处理技术领域。
背景技术
在煤炭开采过程中,煤岩体受压破裂产生微地震信号,一次破裂产生的微地震信号称为一个微震事件,收集尽可能多的微震事件并对其进行定位及能量计算,可以直观反映围岩采动应力分布情况,进而揭示出应力场的时空演化规律,对冲击地压等煤岩动力灾害进行灾变预警具有重要意义。一般而言,微震事件持续时间很短,单个微震事件仅有几十毫秒,通常情况下数据采集现场的外部噪声干扰比较严重,有时甚至把微震事件淹没在噪声中,因此如何在连续的信号采集过程中准确辨识出微震事件成为微震技术应用的关键环节。
STA/LTA(Short Time Average/Long Time Average)法是最为常用的微震事件识别方法,该方法基于时间域信号能量变化构建特征函数,利用STA与LTA的比值识别震相初至,该方法实现简单,运算速度快,但会把干扰信号误判为微震事件,算法抗噪性能较差;AIC(Akaike Information Criterion,AIC)法是另外一种常用的微震事件辨识方法,该法依据微震波到达前后的波形数据统计差别,给出了判别震相初至的AIC准则,AIC函数在震相初至到来时局部最小值峰值尖锐,准确度较高,但AIC函数受时窗影响较大,不同的时窗下其最小值出现的位置也不一样,因此仅适用于已知震相初至的大致位置情形下;目前,微震事件辨识主要是在上述算法软件的辅助下进行人工判断,处理数据的效率低,且辨识标准不统一,极易出现误判。
发明内容
针对前述方法的局限性,本发明将排列熵引入微震事件辨识领域,根据信号的特点研究了尺度因子对排列熵的影响,给出了多尺度排列熵(Multi-scale PermutationEntropy,MPE)的概念,提出了基于MPE的微震信号特征提取方法,并与最小二乘支持向量机(Least Squares Support Vector Machine,LS-SVM)理论相结合,建立了低信噪比微震事件辨识模型及方法。
排列熵是衡量一维时间序列空间复杂度的平均熵参数,具有计算简单快捷、抗噪声能力强等特点,是一种新的动力学突变检测方法,能够很好的反映一维信号时间序列的微小变化,已在医学研究、故障诊断、信号处理等领域广泛应用,并已取得显著的应用效果。为便于叙述,先给出排列熵的定义如下:
对于单一时间序列{xi,i=1,2,…,n},对其中的任意一个元素进行相空间重构,得:
X(i)=[x(i),x(i+τ),…,x(i+(m-1)τ)], (1)
式中m和τ分别为嵌入维数和延迟时间。将X(i)的m个重构分量按照升序进行排列,得:
x(i+(j1-1)τ)≤x(i+(j2-1)τ)≤…≤x(i+(jm-1)τ). (2)
如果存在x(i+(ji1-1)τ)=x(i+(ji2-1)τ),此时就按j值大小来排序,即当ji1<ji2时,有:
x(i+(ji1-1)τ)≤x(i+(ji2-1)τ), (3)
所以,任意一个向量X(i)都可以得到一组符号序列
P(l)=[j1,j2,…,jm]. (4)
式中l=1,2,…,q,且q≤m!,m个不同的符号序列[j1,j2,…,jm]共有m!种排列,符号序列P(l)是其中的一种。如果每一种符号序列的出现的概率为P1,P2,…,Pq,按照Shannon信息熵的形式,时间序列x(i)的q种不同符号序列的排列熵可定义为:
当Pj=1/m!时,Hp(m)达到最大值ln(m!),将Hp(m)进行归一化处理,得
0≤Hp=Hp(m)/ln(m!)≤1. (6)
基于上述解释,本发明采用的技术方案如下:
一种基于多尺度排列熵的低信噪比微震事件辨识方法,其特征在于,它是在计算机上实现的,具体包括以下步骤:
步骤1:采集微震数据
在回采巷道布设拾震器,拾震器之间距离d=50-200米,设置拾震器采集频率为1kHz,使用光纤把各拾震器采集的微震信号传输到数据集中器后,数据集中器再把微震信号发送给计算机进行处理;
步骤2:计算机读取微震数据
当微震产生时,触发拾震器记录微震信号,并将微震信号通过数据集中器传输给计算机,计算机接收后转换为时间序列数据,记为时序序列Xn={x1,x2,…,xn},取n介于[2500,4500]之间,然后进入步骤3;
步骤3:粗粒化时间序列Xn
分别取尺度因子q=1,2,…,10,按照式(7)粗粒化时序序列Xn,则有:
步骤4:计算时序序列Xn的多尺度排列熵
则信号Xn的多尺度排列熵可记为:
Hmpe是一个10维的向量,用来表达微震信号的特征;
步骤5:基于多尺度排列熵训练最小二乘支持向量机
选择100组微震信号及100组噪声信号,重复步骤2-4,分别计算这200组信号的多尺度排列熵,分别标记100组微震信号的多尺度排列熵对应的判别结果为1,分别标记100组噪声信号的多尺度排列熵对应的判别结果为0,构成样本数据集S,使用LS-SVM算法对样本数据集S进行训练,训练LS-SVM时使用10-折交叉验证法确定正规化参数及核函数参数,得到一个能够区分微震事件与噪声事件的支持向量机SV;
步骤6:利用训练好的LS-SVM对待判别信号进行辨识
重复步骤2-4,计算待判别信号Yn={y1,y2,…,yn}的多尺度排列熵,取n介于[2500,4500]之间,得到一个10维的向量Hmpe,将Hmpe输入到支持向量机SV中,得到支持向量机输出值;
当支持向量机输出值为0时,表示待判别信号Ym不是微震事件;
当支持向量机输出值为1时,表示待判别信号Ym是微震事件。
本发明的积极效果是:排列熵是衡量信号复杂度的平均熵参数,具有计算简单、抗噪声能力强等特点,将多尺度和排列熵结合应用于微震信号可以更加有效地分析信号特征,准确表达微震信号多个维度的波形特征,有利于区分微震事件及噪声事件。应用多尺度排列熵分别提取微震信号及噪声信号的特征数据,使用LS-SVM对特征数据进行训练,得到LS-SVM向量机,实现了对低信噪比微震事件及噪声事件的准确分类。
附图说明
附图1为本发明微震监测系统布置图;
附图2为本发明中多尺度排列熵处理过程图;
附图3为本发明中LS-SVM训练与判断过程图;
附图4为本发明实施例中微震信号波形图;
附图5为本发明实施例中待判别的信号波形图。
具体实施方式
下面结合附图和具体实施例进一步说明本发明的技术方案。
实施例,一种基于多尺度排列熵的低信噪比微震事件辨识方法,包括以下步骤:
步骤1:按图1所示布置好微震监测系统,拾震器之间的距离d=100米,信号集中器采集微震数据并传输给计算机。
步骤2:计算机接收信号后转换为时间序列数据,记为时序序列,读取时序序列Xn={x1,x2,…,xn},n=3000,时序序列如图4所示,然后进入步骤3。
步骤3:分别取q=1,2,…,10,按照式(7)计算得到10个时间序列,分别记为X1,X2,…,X10。
步骤4:按照式(6)分别计算10个时间序列X1,X2,…,X10的排列熵H1,H2,…,H10,计算时,令嵌入维数m=4,延迟时间为1,排列熵组成一个向量代表时序序列Xn的多尺度排列熵,记为Hmpe=[H1,H2,…,H10];时序序列Xn的多尺度排列熵为Hmpe=[0.8966,0.6788,0.4927,0.3933,0.3309,0.2915,0.2574,0.2276,0.2063,0.1891];
上述步骤2-4的处理过程见图2。
步骤5:选择200组信号,其中微震信号及噪声信号各100组,且每组信号的长度为3000个采样点,按照步骤2-4分别对200组信号求取它们的多尺度排列熵,构成表1所示的样本数据,应用LS-SVM对表1中数据进行训练,得到能够区分微震事件与噪声信号的支持向量机Sv;
表1样本数据
注:“结果”栏内1表示该信号是微震事件;0表示该信号为噪声。
步骤6:对于待判别信号Yn={y1,y2,…,yn},n=3000,如图5所示,重复步骤2-4,计算其多尺度排列熵Hmpe=[0.9854,0.7625,0.5016,0.4102,0.3455,0.2816,0.2655,0.2318,0.2234,0.1798],将其作为支持向量机Sv的输入,支持向量机Sv输出判别结果为0,这表明待判别信号Yn不是微震事件。
上述步骤5-6的过程见图3.
以上所述为本发明最佳实施方式的举例,其中未详细述及的部分均为本领域普通技术人员的公知常识。本发明的保护范围以权利要求的内容为准,任何基于本发明的技术启示而进行的等效变换,也在本发明的保护范围之内。
Claims (1)
1.一种基于多尺度排列熵的低信噪比微震事件辨识方法,其特征在于,它是在计算机上实现的,具体包括以下步骤:
步骤1:采集微震数据
在回采巷道布设拾震器,拾震器之间距离d=50-200米,设置拾震器采集频率为1kHz,使用光纤把各拾震器采集的微震信号传输到数据集中器后,数据集中器再把微震信号发送给计算机进行处理;
步骤2:计算机读取微震数据
当微震产生时,触发拾震器记录微震信号,并将微震信号通过数据集中器传输给计算机,计算机接收后转换为时序序列数据,记为时序序列Xn={x1,x2,…,xn},取n介于[2500,4500]之间,然后进入步骤3;
步骤3:粗粒化时序序列Xn
分别取尺度因子q=1,2,…,10,按照下式粗粒化时序序列Xn,则有:
步骤4:计算时序序列Xn的多尺度排列熵
则Xn的多尺度排列熵记为:
Hmpe是一个10维的向量,用来表达Xn的特征;
步骤5:基于多尺度排列熵训练最小二乘支持向量机
选择100组微震信号及100组噪声信号,重复步骤3-4,分别计算这200组信号的多尺度排列熵,分别标记100组微震信号的多尺度排列熵对应的判别结果为1,分别标记100组噪声信号的多尺度排列熵对应的判别结果为0,构成样本数据集S,使用LS-SVM算法对样本数据集S进行训练,训练LS-SVM时使用10-折交叉验证法确定正规化参数及核函数参数,得到一个能够区分微震事件与噪声事件的支持向量机SV;
步骤6:利用训练好的LS-SVM对待判别信号进行辨识
按照步骤3-4记载的计算微震信号多尺度排列熵办法,计算待判别信号Yn={y1,y2,…,yn}的多尺度排列熵,取n介于[2500,4500]之间,得到一个10维的向量Hmpe,将Hmpe输入到支持向量机SV中,得到支持向量机输出值;
当支持向量机输出值为0时,表示待判别信号Yn不是微震事件;
当支持向量机输出值为1时,表示待判别信号Yn是微震事件。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201610257713.5A CN105956526B (zh) | 2016-04-22 | 2016-04-22 | 基于多尺度排列熵的低信噪比微震事件辨识方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201610257713.5A CN105956526B (zh) | 2016-04-22 | 2016-04-22 | 基于多尺度排列熵的低信噪比微震事件辨识方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN105956526A CN105956526A (zh) | 2016-09-21 |
CN105956526B true CN105956526B (zh) | 2020-01-24 |
Family
ID=56915804
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201610257713.5A Expired - Fee Related CN105956526B (zh) | 2016-04-22 | 2016-04-22 | 基于多尺度排列熵的低信噪比微震事件辨识方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN105956526B (zh) |
Families Citing this family (15)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
AU2016426667A1 (en) * | 2016-10-20 | 2019-03-07 | Landmark Graphics Corporation | Classifying well data using a support vector machine |
CN106896407B (zh) * | 2017-03-28 | 2018-07-13 | 吉林大学 | 一种基于近似负熵的微地震信号初至拾取方法 |
CN107422373B (zh) * | 2017-04-21 | 2018-05-18 | 山东科技大学 | 一种基于特征漂移的煤矿冲击地压灾害在线预警方法 |
CN107505652B (zh) * | 2017-07-26 | 2018-11-13 | 山东科技大学 | 一种基于能量分布特征的矿山微震信号辨识方法 |
CN107944469B (zh) * | 2017-11-03 | 2021-12-21 | 桂林电子科技大学 | 一种基于重排st的低信噪比微震事件辨识方法 |
CN107885940A (zh) * | 2017-11-10 | 2018-04-06 | 吉林大学 | 一种用于分布式光纤振动传感系统的信号特征提取方法 |
CN107944613A (zh) * | 2017-11-17 | 2018-04-20 | 湖南大学 | 一种基于符号网络的灾变预测方法 |
CN108362488A (zh) * | 2018-02-28 | 2018-08-03 | 河海大学 | 基于mpe与svm的oltc机械故障诊断方法 |
CN108983051B (zh) * | 2018-07-25 | 2021-02-09 | 武汉科技大学 | 基于同步挤压小波变换的局部放电类型识别方法 |
CN109143371A (zh) * | 2018-08-01 | 2019-01-04 | 中国石油化工股份有限公司江汉油田分公司物探研究院 | 一种地震数据的噪声去除方法及装置 |
CN110886599B (zh) * | 2018-09-07 | 2021-09-17 | 中国石油化工股份有限公司 | 基于破裂速度的非压裂事件识别方法及系统 |
WO2020220416A1 (zh) * | 2019-04-28 | 2020-11-05 | 山东科技大学 | 一种基于深度学习的微震信号分类辨识方法 |
CN110441654B (zh) * | 2019-07-29 | 2021-05-04 | 三峡大学 | 基于改进ewt和cmpe的配电网电能质量扰动检测方法 |
CN110333078B (zh) * | 2019-08-21 | 2021-04-20 | 佛山科学技术学院 | 一种滚动轴承退化状态阶段确定方法 |
CN114325823B (zh) * | 2021-12-30 | 2022-07-15 | 广西大学 | 一种基于nb-iot的岩体破裂失稳微震信号无线监测方法与装置 |
Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN104849050A (zh) * | 2015-06-02 | 2015-08-19 | 安徽工业大学 | 一种基于复合多尺度排列熵的滚动轴承故障诊断方法 |
CN104897403A (zh) * | 2015-06-24 | 2015-09-09 | 北京航空航天大学 | 一种基于排列熵和流形改进动态时间规整的自适应故障诊断方法 |
CN105004498A (zh) * | 2015-07-09 | 2015-10-28 | 西安理工大学 | 一种水电机组的振动故障诊断方法 |
Family Cites Families (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN103616719B (zh) * | 2013-12-05 | 2017-02-08 | 中国矿业大学(北京) | 具有噪声识别和自适应放大功能的微地震采集装置和方法 |
-
2016
- 2016-04-22 CN CN201610257713.5A patent/CN105956526B/zh not_active Expired - Fee Related
Patent Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN104849050A (zh) * | 2015-06-02 | 2015-08-19 | 安徽工业大学 | 一种基于复合多尺度排列熵的滚动轴承故障诊断方法 |
CN104897403A (zh) * | 2015-06-24 | 2015-09-09 | 北京航空航天大学 | 一种基于排列熵和流形改进动态时间规整的自适应故障诊断方法 |
CN105004498A (zh) * | 2015-07-09 | 2015-10-28 | 西安理工大学 | 一种水电机组的振动故障诊断方法 |
Non-Patent Citations (3)
Title |
---|
Simultaneous seismic data denoising and reconstruction via multichannel singular spectrum analysis;Vicente Oropeza 等;《GEOPHYSICS》;20110428;第76卷(第3期);第25-32页 * |
基于非线性信号分析的滚动轴承状态监测诊断研究;刘永斌;《中国博士学位论文全文数据库工程科技Ⅱ辑》;20110915(第9期);第5页 * |
多尺度排列熵及其在滚动轴承故障诊断中的应用;郑近德 等;《中国机械工程》;20131015;第24卷(第19期);第2642-2644页 * |
Also Published As
Publication number | Publication date |
---|---|
CN105956526A (zh) | 2016-09-21 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN105956526B (zh) | 基于多尺度排列熵的低信噪比微震事件辨识方法 | |
CN109974782B (zh) | 基于大数据敏感特征优化选取的设备故障预警方法及系统 | |
CN105044590B (zh) | 一种单测量节点模拟电路故障诊断方法 | |
CN108830328B (zh) | 融合空间知识的微震信号smote识别方法及监测系统 | |
CN105956623A (zh) | 一种基于模糊熵的癫痫脑电信号分类方法 | |
CN104198898A (zh) | 基于脉冲序列分析的局部放电发展过程诊断方法 | |
CN101666677B (zh) | 一种用于提取低速重载设备微弱故障特征信息的方法 | |
CN106382981B (zh) | 一种单站次声波信号识别提取方法 | |
CN110196448B (zh) | 一种滑坡次声信号识别方法 | |
CN103792000A (zh) | 一种基于稀疏表示的信号中瞬态成分检测方法及装置 | |
CN105067101A (zh) | 振源识别的基于振动信号的基音频率特征的提取方法 | |
CN103267652B (zh) | 一种智能早期设备故障在线诊断方法 | |
CN112364296A (zh) | 一种基于深度学习的p波到时自动拾取方法 | |
CN110398775B (zh) | 隧道突涌水灾害微震事件信号波动初至拾取方法及系统 | |
CN111275003B (zh) | 基于类最优高斯核多分类支持向量机的微震信号识别方法 | |
CN105139869A (zh) | 一种基于区间差分特征的婴儿啼哭检测方法 | |
CN112215307B (zh) | 一种应用机器学习自动检测地震仪器信号异常的方法 | |
CN116881712A (zh) | 一种混凝土大坝活动裂纹电磁脉冲信号识别方法 | |
CN104156339B (zh) | 一种利用二次排列熵识别周期微弱脉冲信号的方法 | |
CN113848253B (zh) | 仿真变电站主变压器基底渗水状态声发射监测方法及装置 | |
CN110346032A (zh) | 一种基于恒虚警与过零率相结合的φ-otdr振动信号端点检测方法 | |
CN116465975A (zh) | 基于声发射波形识别的断层滑移预警方法及系统 | |
CN115130748A (zh) | 一种建筑电气系统故障诊断方法 | |
CN104102834A (zh) | 录音地点的识别方法 | |
CN108984997B (zh) | 一种基于熵值的离心泵故障诊断方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant | ||
CF01 | Termination of patent right due to non-payment of annual fee | ||
CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20200124 Termination date: 20210422 |