CN105487129A - 一种地空时域电磁数据高度校正方法 - Google Patents
一种地空时域电磁数据高度校正方法 Download PDFInfo
- Publication number
- CN105487129A CN105487129A CN201610005275.3A CN201610005275A CN105487129A CN 105487129 A CN105487129 A CN 105487129A CN 201610005275 A CN201610005275 A CN 201610005275A CN 105487129 A CN105487129 A CN 105487129A
- Authority
- CN
- China
- Prior art keywords
- magnetic
- data
- ground
- rate
- change
- 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
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V3/00—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation
- G01V3/38—Processing data, e.g. for analysis, for interpretation, for correction
Abstract
本发明涉及一种地空时域电磁数据高度校正方法。主要包括采用无人机或飞艇在空中进行飞行测量,先将测量的数据进行电磁噪声去除、双极性叠加取样;再应用Kriging插值、双三次B样条插值法,将飞行测量中各测点的磁场变化率进行均匀插值,得到二维均匀网格中的磁场变化率;应用逆延拓理论,将空中z=-h高度处测量的磁场变化率进行延拓计算,获得地面z=0处的磁场变化率数据,进而实现了飞行高度的校正。解决了LOTEM方法解释地空时域电磁数据中忽略高度的问题。可直接使用磁场变化率数据进行校正,无需将磁场变化率转化为磁场的进行校正,只要3条飞行测线就可进行高度校正,计算过程仅有插值和傅里叶变换,计算简单、计算速度快。
Description
技术领域
本发明涉及一种地球物理勘探中的航空物探方法,尤其是地空时域电磁探测数据高度校正方法。
背景技术
地空时域电磁法(Ground-airborneTime-domainelectromagneticsystem,简称GATEM)是通过接地长导线供以双极性脉冲电流激励大地,在空间形成大面积的电磁波,通过空中无人飞行器或飞艇上搭载的接收系统进行飞行测量,获取地下磁场变化率数据,通过电阻率成像就可发现地下异常地质体的存在,并确定异常体的电性结构和空间分布形态。地空时域电磁探测方法,不仅具有航空时域电磁法的探测高效等优势,还具有地面时域电磁法的信噪比高、勘探深度大的优势,适合在沼泽地带、滩涂地区、无人山区等特殊地区开展探测工作,尤其适用于我国地形复杂的山区资源探测,近几年成为国内外研究热点,有着广泛应用前景。
CN201410081433.4公开了电性源瞬变电磁地空探测方法,公开了电性源地空探测方法,采用全域视电阻率方法对观测数据进行处理,以及采用逆合成孔径方法对地空数据进行解释处理。
CN200910263449.6公开了一种时间域航空电磁高度校正方法,公开了利用x、z两个分量的感生电动势数据联合计算等效均匀半空间电阻率序列,并通过计算均匀半空间x、z两个分量的感生电动势理论幅值比作为查询数据,依次采用查询和局部样条插值获得对应时刻的等效均匀半空间电阻率,利用其正演计算均匀半空间模型的感生电动势,即为校正结果。
以上现有技术,一个公开了电性源瞬变电磁地空探测的方法和解释方法,但并未对高度进行校正;一个提供了一种航空电磁高度校正方法,虽然可以采用均匀半空间模型等效思想对地空数据进行高度校正,但是由于实际的地质情况是复杂的地质结构体,感生电动势在不同介质之间存在过渡过程,导致利用均匀半空间模型等效时,无法准确找到电阻率值,因此只有已知地下准确的电导值及其深度,才能正确校正高度。
发明内容
本发明的目的在于针对现有技术的不足,提供一种基于逆延拓理论的地空时域电磁数据高度校正方法。
本发明的目的是通过以下技术方案实现的:
一种地空时域电磁数据高度校正方法,包括以下步骤:
a、按照勘探的设计要求在测区布置接地长导线,在测区布设平行于长导线三条以上飞行航线;
b、采用接地长导线作为发射源,供以双极性脉冲电流,通过无人飞行器或飞艇搭载的接收系统飞行测量,获取地下磁场变化率数据;
c、将接收到的地下磁场数据进行双极性叠加、去噪及数据取样处理;
d、采用Kriging插值和双三次B样条插值法,将飞行测量中各测点的磁场变化率进行均匀插值,获得二维均匀网格中的磁场变化率;
e、采用逆延拓理论,将空中z=-h高度处测量的磁场变化率进行延拓计算,获得地面z=0处的磁场变化率数据,进而实现飞行高度校正。
步骤c所述的双极性叠加和去噪,包括以下步骤:
A、输入实测数据的单点衰减曲线;
B、双极性叠加;
C、加权均值滤波;
D、指数拟合-自适应卡尔曼滤波;
E、判断信噪比是否大于预设值,是;
F、输出滤波后结果。
步骤c所述的数据取样处理,包括以下步骤:
Ⅰ、输入滤波后的衰减曲线;
Ⅱ、确定取样起始点;
Ⅲ、计算数据窗宽;
Ⅳ、取样积分计算:
式中:Xi为数据窗宽内第i个数据,N为数据窗宽,X为取样结果值;
Ⅴ、判断数据是否结束,是;
Ⅵ、输出取样后结果。
步骤d所述的采用Kriging插值,是将3条飞行测线上各测点滤波取样后的磁场变化率插值到z=-h高度平面;在z=-h高度平面,再应用双三次B样条插值,分别对横轴、纵轴的各测点的磁场变化率进行插值,得到均匀网格中的磁场变化率。
步骤e所述的采用逆延拓理论推导空中z=-h高度处测量的磁场变化率和地面z=0的磁场变化率的变换关系式为:
式中:x,y为接收点所在的空间坐标,kx与ky是对应于x与y的波数域变量,为空中z=-h高度测量的磁场变化率在波数域中的表示形式,为地面z=0的磁场变化率。
有益效果:本发明与现有技术相比计算简单,无需已知地下准确的电阻率及其深度,就能够有效地实现地空时域电磁高度校正,解决了LOTEM(长偏移瞬变电磁法,LongOffsetTransientElectromagneticMethod)方法解释地空时域电磁数据中忽略高度的问题。(这种近似等效的方法,忽略接收高度的影响,导致计算的电阻率与真实电阻率存在较大偏差)。可直接使用磁场变化率数据进行校正,无需将磁场变化率转化为磁场的进行校正,只要至少3条飞行测线就可进行高度校正,计算过程仅有插值和傅里叶变换,计算简单、快速,完全符合地空时域电磁要求。
附图说明
图1为地空时域电磁数据高度校正方法测量布置图;
图2为地空时域电磁数据高度校正方法流程图;
图3为地空时域电磁实测数据单点衰减曲线滤波流程图;
图4为数据取样算法流程图;
图5为空中磁场变化率逆延拓方法示意图;
图6为本发明一个实施例的地空时域电磁数据校正效果图。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白,下面结合附图和实施例对本发明作进一步的详细说明。
地空时域电磁数据高度校正方法,包括以下步骤:
a、按照勘探的设计要求在测区布置接地长导线,在测区布设平行于长导线三条以上飞行航线;
b、采用接地长导线作为发射源,供以双极性脉冲电流,通过无人飞行器或飞艇搭载的接收系统飞行测量,获取地下磁场变化率数据;
c、将接收到的地下磁场数据进行双极性叠加、去噪及数据取样处理;
d、采用Kriging插值和双三次B样条插值法,将飞行测量中各测点的磁场变化率进行均匀插值,获得二维均匀网格中的磁场变化率;
e、采用逆延拓理论,将空中z=-h高度处测量的磁场变化率进行延拓计算,获得地面z=0处的磁场变化率数据,进而实现飞行高度校正。
步骤c所述的双极性叠加和去噪,包括以下步骤:
A、输入实测数据的单点衰减曲线;
B、双极性叠加;
C、加权均值滤波;
D、指数拟合-自适应卡尔曼滤波;
E、判断信噪比是否大于预设值,是;
F、输出滤波后结果。
步骤c所述的数据取样处理,包括以下步骤:
Ⅰ、输入滤波后的衰减曲线;
Ⅱ、确定取样起始点;
Ⅲ、计算数据窗宽;
Ⅳ、取样积分计算:
式中:Xi为数据窗宽内第i个数据,N为数据窗宽,X为取样结果值;
Ⅴ、判断数据是否结束,是;
Ⅵ、输出取样后结果。
步骤d所述的采用Kriging插值,是将3条飞行测线上各测点滤波取样后的磁场变化率插值到z=-h高度平面;在z=-h高度平面,再应用双三次B样条插值,分别对横轴、纵轴的各测点的磁场变化率进行插值,得到均匀网格中的磁场变化率。
步骤e所述的采用逆延拓理论推导空中z=-h高度处测量的磁场变化率和地面z=0的磁场变化率的变换关系式为:
式中:x,y为接收点所在的空间坐标,kx与ky是对应于x与y的波数域变量,为空中z=-h高度测量的磁场变化率在波数域中的表示形式,为地面z=0的磁场变化率。
实施例
参见图2结合图1所示,一种基于逆延拓理论的地空时域电磁数据高度校正方法,包括:
a、合理布置接地长导线的位置,在探测区域设计平行于长导线的至少三条飞行航线;
如图1所示的步骤1中,接地导线源长度为1-3千米,设置空中飞行测线需要与接地导线源平行,测线长度应不超过接地导线源长度、并且测线需位于收发距为0.2-1倍区间内。因为飞行测线上各测点的磁场变化率要插值成均匀网格点的磁场变化率,若飞行测线间距均匀则不需要测线间插值,校正精度更高。
b、采用接地长导线作为发射源供以双极性脉冲电流,通过空中无人飞行器或飞艇上搭载的接收系统进行飞行测量,获取地下磁场变化率数据;
如图1所示的步骤2中,接地导线源发射电流为不小于10A,结合实际飞行测量工区的地质条件,无人飞行器或飞艇飞行高度在保证飞行安全情况下为10-50米。
c、将空中实测数据进行双极性叠加、电磁噪声去除以及数据取样处理;
地空时域电磁系统在低空飞行测量过程中,飞行器受风向、气流、飞行不平稳、地形等影响产生运动噪声,以及工频噪声、天电干扰等电磁噪声,地空系统测量的电磁数据中还含有白噪声,导致其数据受到较为复杂的噪声干扰。为了提高数据质量,需要对实测的数据进行双极性叠加去除白噪声和随机噪声,再通过指数拟合-自适应卡尔曼滤除电磁噪声。
如图3所示的步骤3中地空实测电磁数据单点衰减曲线滤波,基于地空时域电磁信号的特点,采用双极性叠加、卡尔曼滤波的消噪方法,来压制数据中的白噪声以及其他电磁噪声。
进如图4所示的数据取样算法,先确定第一个取样道的起始点即从第一个采样点开始,由近似对数等间隔进行计算,数据窗宽为100.1(1.26)倍关系;取样积分计算公式为:
Xi为数据窗宽内第i个数据,N为数据窗宽,X为取样结果值,按此方法处理至数据结束为止。
d、应用Kriging插值、双三次B样条插值法,将飞行测量中各测点的磁场变化率进行均匀插值,获得二维均匀网格中的磁场变化率;
步骤4中,由于测量的参数为感应电动势,将感应电动势V变换为磁场变化率:
式中Vz(t)为测量感应电动势,变换后的磁场变换率,n为线圈匝数,S为接收线圈有效面积,μ0为真空中磁场导率。
进一步,如图5所示的步骤4中,由于飞行不平稳,各接收点高度存在细微偏差,应用Kriging插值,将3条飞行测线上各测点滤波取样后的磁场变化率插值到z=-h高度平面。在z=-h高度平面,平行于接地电性源水平线为x轴,垂直为y轴,二维傅里叶变换的前提条件是中的变量x和y应该是等间距的,而飞行测量中不同测线的测点可能不等间距,所以再应用双三次B样条插值,分别对x,y轴各测点的磁场变化率进行插值,得到均匀网格中的磁场变化率。
e、采用逆延拓理论,将空中z=-h高度处测量的磁场变化率进行延拓计算,获得地面z=0处的磁场变化率数据,进而实现飞行高度校正。
如图5所示的步骤5中,应用逆延拓理论,将均匀网格中的磁场变化率进行二维傅里叶变换,实现空间域变换到波数域:
再根据矢量拉普拉斯方程,推导空中磁场变化率和地面磁场变化率的变换关系,在任意取样时刻,空中z=-h高度处测量的磁场变化率和地面z=0的磁场变化率关系表达式为:
表达式(4)中:kx与ky是对应于x与y的波数域变量,为空中z=-h高度测量的磁场变化率在波数域中的表示形式,为地面z=0的磁场变化率在波数域中的表示形式。再在波数域进行逆傅里叶变换:
即可得到地面z=0空间域的磁场变化率。
图6为采用图2所示的本发明一个实施例的地空时域电磁数据校正结果,校正后的感应电动势曲线明显接近地面解析解结果,相对误差为12%,其中相对误差计算公式为(校正后的磁场变化率-理论模型在地面的磁场变化率)/理论模型在地面的磁场变化率*100%,充分验证了基于逆延拓理论的地空时域电磁高度校正方法的有效性,为地空时域电磁数据法野外高精度解释提供了新的思路和方法。
Claims (5)
1.一种地空时域电磁数据高度校正方法,其特征在于,包括以下步骤:
a、按照勘探的设计要求在测区布置接地长导线,在测区布设平行于长导线三条以上飞行航线;
b、采用接地长导线作为发射源,供以双极性脉冲电流,通过无人飞行器或飞艇搭载的接收系统飞行测量,获取地下磁场变化率数据;
c、将接收到的地下磁场数据进行双极性叠加、去噪及数据取样处理;
d、采用Kriging插值和双三次B样条插值法,将飞行测量中各测点的磁场变化率进行均匀插值,获得二维均匀网格中的磁场变化率;
e、采用逆延拓理论,将空中z=-h高度处测量的磁场变化率进行延拓计算,获得地面z=0处的磁场变化率数据,进而实现飞行高度校正。
2.按照权利要求1所述的地空时域电磁数据高度校正方法,其特征在于,步骤c所述的双极性叠加和去噪,包括以下步骤:
A、输入实测数据的单点衰减曲线;
B、双极性叠加;
C、加权均值滤波;
D、指数拟合-自适应卡尔曼滤波;
E、判断信噪比是否大于预设值,是;
F、输出滤波后结果。
3.按照权利要求1所述的地空时域电磁数据高度校正方法,其特征在于,步骤c所述的数据取样处理,包括以下步骤:
Ⅰ、输入滤波后的衰减曲线;
Ⅱ、确定取样起始点;
Ⅲ、计算数据窗宽;
Ⅳ、取样积分计算;
Ⅴ、判断数据是否结束,是;
Ⅵ、输出取样后结果。
4.按照权利要求1所述的地空时域电磁数据高度校正方法,其特征在于,步骤d所述的采用Kriging插值,是将3条飞行测线上各测点滤波取样后的磁场变化率插值到z=-h高度平面;在z=-h高度平面,再应用双三次B样条插值,分别对横轴、纵轴的各测点的磁场变化率进行插值,得到均匀网格中的磁场变化率。
5.按照权利要求1所述的地空时域电磁数据高度校正方法,其特征在于,
步骤e所述的采用逆延拓理论推导空中z=-h高度处测量的磁场变化率和地面z=0的磁场变化率的变换关系式为:
式中:x,y为接收点所在的空间坐标,kx与ky是对应于x与y的波数域变量,为空中z=-h高度测量的磁场变化率在波数域中的表示形式,为地面z=0的磁场变化率。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201610005275.3A CN105487129B (zh) | 2016-01-06 | 2016-01-06 | 一种地空时域电磁数据高度校正方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201610005275.3A CN105487129B (zh) | 2016-01-06 | 2016-01-06 | 一种地空时域电磁数据高度校正方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN105487129A true CN105487129A (zh) | 2016-04-13 |
CN105487129B CN105487129B (zh) | 2017-09-26 |
Family
ID=55674234
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201610005275.3A Active CN105487129B (zh) | 2016-01-06 | 2016-01-06 | 一种地空时域电磁数据高度校正方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN105487129B (zh) |
Cited By (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106772633A (zh) * | 2017-01-20 | 2017-05-31 | 吉林大学 | 一种具有姿态记录功能的地空tem接收系统及测量方法 |
CN108072906A (zh) * | 2016-11-18 | 2018-05-25 | 北京自动化控制设备研究所 | 一种分布式磁探测磁目标识别方法 |
CN109739263A (zh) * | 2019-01-25 | 2019-05-10 | 清华大学 | 一种基于磁信号延拓算法进行潜艇探测的探潜机导航方法 |
CN110794464A (zh) * | 2019-11-07 | 2020-02-14 | 中国地质科学院地球物理地球化学勘查研究所 | 一种去除时间域航空电磁数据背景场的方法 |
CN112730051A (zh) * | 2020-12-23 | 2021-04-30 | 西安交通大学 | 一种基于傅里叶微分变换的金属板应变测量方法及其系统 |
CN114814963A (zh) * | 2022-03-03 | 2022-07-29 | 吉林大学 | 一种基于b样条抽道的时间域电磁数据甚低频噪声压制方法 |
Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101636669A (zh) * | 2007-03-23 | 2010-01-27 | 阿克斯有限责任公司 | 地形校正系统 |
CN101710187A (zh) * | 2009-12-17 | 2010-05-19 | 成都理工大学 | 一种时间域航空电磁高度校正方法 |
CN102096113A (zh) * | 2010-12-03 | 2011-06-15 | 吉林大学 | 时间域地空电磁探测系统及标定方法 |
CN104360403A (zh) * | 2014-11-08 | 2015-02-18 | 长安大学 | 一种多辐射场源瞬变电磁法多分量全域视电阻率定义方法 |
WO2015088205A1 (ko) * | 2013-12-09 | 2015-06-18 | 한국지질자원연구원 | 3차원 항공 자력 탐사 시스템 및 이를 이용한 3차원 항공 자력 탐사 방법 |
-
2016
- 2016-01-06 CN CN201610005275.3A patent/CN105487129B/zh active Active
Patent Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101636669A (zh) * | 2007-03-23 | 2010-01-27 | 阿克斯有限责任公司 | 地形校正系统 |
CN101710187A (zh) * | 2009-12-17 | 2010-05-19 | 成都理工大学 | 一种时间域航空电磁高度校正方法 |
CN102096113A (zh) * | 2010-12-03 | 2011-06-15 | 吉林大学 | 时间域地空电磁探测系统及标定方法 |
WO2015088205A1 (ko) * | 2013-12-09 | 2015-06-18 | 한국지질자원연구원 | 3차원 항공 자력 탐사 시스템 및 이를 이용한 3차원 항공 자력 탐사 방법 |
CN104360403A (zh) * | 2014-11-08 | 2015-02-18 | 长安大学 | 一种多辐射场源瞬变电磁法多分量全域视电阻率定义方法 |
Non-Patent Citations (3)
Title |
---|
于长春 等: "航磁剖面异常高度改正方法研究", 《"九五"全国地质科技重要成果论文集》 * |
嵇艳鞠 等: "ATEM瞬变电磁系统在长春市活断层勘探中的应用", 《吉林大学学报(地球科学版)》 * |
嵇艳鞠 等: "无人飞艇长导线源时域地空电磁勘探系统及其应用", 《地球物理学报》 * |
Cited By (9)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN108072906A (zh) * | 2016-11-18 | 2018-05-25 | 北京自动化控制设备研究所 | 一种分布式磁探测磁目标识别方法 |
CN108072906B (zh) * | 2016-11-18 | 2020-01-10 | 北京自动化控制设备研究所 | 一种分布式磁探测磁目标识别方法 |
CN106772633A (zh) * | 2017-01-20 | 2017-05-31 | 吉林大学 | 一种具有姿态记录功能的地空tem接收系统及测量方法 |
CN109739263A (zh) * | 2019-01-25 | 2019-05-10 | 清华大学 | 一种基于磁信号延拓算法进行潜艇探测的探潜机导航方法 |
CN109739263B (zh) * | 2019-01-25 | 2020-06-30 | 清华大学 | 一种基于磁信号延拓算法进行潜艇探测的探潜机导航方法 |
CN110794464A (zh) * | 2019-11-07 | 2020-02-14 | 中国地质科学院地球物理地球化学勘查研究所 | 一种去除时间域航空电磁数据背景场的方法 |
CN112730051A (zh) * | 2020-12-23 | 2021-04-30 | 西安交通大学 | 一种基于傅里叶微分变换的金属板应变测量方法及其系统 |
CN112730051B (zh) * | 2020-12-23 | 2022-04-05 | 西安交通大学 | 一种基于傅里叶微分变换的金属板应变测量方法及其系统 |
CN114814963A (zh) * | 2022-03-03 | 2022-07-29 | 吉林大学 | 一种基于b样条抽道的时间域电磁数据甚低频噪声压制方法 |
Also Published As
Publication number | Publication date |
---|---|
CN105487129B (zh) | 2017-09-26 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN105487129A (zh) | 一种地空时域电磁数据高度校正方法 | |
Xue et al. | The application of ground-airborne TEM systems for underground cavity detection in China | |
Liu et al. | 3D inversion for multipulse airborne transient electromagnetic data | |
CN105353428B (zh) | 一种地面参考区磁场延拓的地空协同电磁数据校正方法 | |
Wu et al. | The development and applications of the semi-airborne electromagnetic system in China | |
CN105589108B (zh) | 基于不同约束条件的瞬变电磁快速三维反演方法 | |
CN112949134B (zh) | 基于非结构有限元方法的地-井瞬变电磁反演方法 | |
CN104237956A (zh) | 电性源瞬变电磁地空探测方法 | |
Jiu‐Long et al. | Synthetic aperture imaging in advanced detection of roadway using mine transient electromagnetic method | |
CN105785455A (zh) | 一种基于b样条插值的二维地面核磁共振反演方法 | |
CN106291725A (zh) | 一种快速反演地下地质体空间位置的方法 | |
CN104360403A (zh) | 一种多辐射场源瞬变电磁法多分量全域视电阻率定义方法 | |
Zhou et al. | Migration velocity analysis and prestack migration of common-transmitter GPR data | |
CN103675926A (zh) | 一种航空瞬变电磁资料的电导率-深度转换方法 | |
Lin et al. | The helicopter time-domain electromagnetic technology advances in China | |
Wu et al. | The development and applications of the helicopter-borne transient electromagnetic system CAS-HTEM | |
CN106405665A (zh) | 基于dbim的瞬变电磁电导率反演方法 | |
CN105372718B (zh) | 一种2d位场构建3d位场的方法 | |
Yanju et al. | Interpretation research on electrical source of time domain ground-airborne electromagnetic data | |
Liu* et al. | A review for cross-hole electromagnetic methods | |
Mao et al. | Reverse-time migration by combining Laplacian filtering with wavefield decomposition for high-resolution subsurface imaging | |
CN105549031A (zh) | 一种卫星信号的电离层传播时延的时域数值计算方法 | |
Wang et al. | Characterization of excitation effects and data interpretation of combined time-domain multi-waveform transmission currents | |
Zhang et al. | Inversion of airborne transient electromagnetic data based on reference point lateral constraint | |
Zhou et al. | Divergence of tipper real induction vector in tensor frequency-domain ground-airborne electromagnetic method |
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 |