CN111123404A - 一种地震和直流电法巷道超前探测的数据融合方法 - Google Patents

一种地震和直流电法巷道超前探测的数据融合方法 Download PDF

Info

Publication number
CN111123404A
CN111123404A CN202010033559.XA CN202010033559A CN111123404A CN 111123404 A CN111123404 A CN 111123404A CN 202010033559 A CN202010033559 A CN 202010033559A CN 111123404 A CN111123404 A CN 111123404A
Authority
CN
China
Prior art keywords
resistivity
data
data set
velocity
grid node
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
Application number
CN202010033559.XA
Other languages
English (en)
Other versions
CN111123404B (zh
Inventor
张平松
邱实
李圣林
郭立全
胡雄武
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Anhui University of Science and Technology
Original Assignee
Anhui University of Science and Technology
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Anhui University of Science and Technology filed Critical Anhui University of Science and Technology
Priority to CN202010033559.XA priority Critical patent/CN111123404B/zh
Publication of CN111123404A publication Critical patent/CN111123404A/zh
Priority to ZA2020/02934A priority patent/ZA202002934B/en
Application granted granted Critical
Publication of CN111123404B publication Critical patent/CN111123404B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V11/00Prospecting or detecting by methods combining techniques covered by two or more of main groups G01V1/00 - G01V9/00

Landscapes

  • Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • General Physics & Mathematics (AREA)
  • Geophysics (AREA)
  • Geophysics And Detection Of Objects (AREA)

Abstract

本发明公开了一种地震和直流电法巷道超前探测的数据融合方法,首先进行地震反射波法速度的单一反演和直流电阻率法电阻率的单一反演,得到速度数据集和电阻率数据集;选择两种探测方法探测的公共区域为数据融合区域,对数据融合区域内数据进行归一化处理;对归一化处理后的数据集计算交叉梯度值并进行离散化;基于离散化的交叉梯度值和数据融合区域内的速度数据集和电阻率数据集计算速度更新量和电阻率更新量;基于速度更新量和电阻率更新量得到新的速度数据集和电阻率数据集,并基于新的速度数据集和电阻率数据集进行速度成像和电阻率成像。本发明的技术方案对地震反射波法与直流电阻率法超前探测的数据进行融合,提高了超前探测的探测精度。

Description

一种地震和直流电法巷道超前探测的数据融合方法
技术领域
本发明涉及工程地球物理勘探技术领域,更具体的说是涉及一种地震反射波法与直流电阻率法超前探测的数据融合方法。
背景技术
在巷道掘进过程中,探明前方的地质条件是十分必要的,这也是保证安全掘进的前提。长期以来,巷道掌子面前方的构造异常、采空区、岩溶空洞、陷落柱等不良地质体,是困扰巷道安全掘进的几大难题,如何更准确的探测前方不良地质体,是现在技术研究的热点。
目前巷道探测多使用地震反射波法、直流电阻率法、瞬变电磁法等物探方法,分别对单一物性数据进行处理,然后进行结果的相互验证。考虑到巷道有限的空间,单独使用上述某一种方法探测时,探测精度相对较低。
因此,如何提供一种高精度巷道超前探测的方法是本领域技术人员亟需解决的问题。
发明内容
有鉴于此,本发明提供了一种地震和直流电法巷道超前探测的数据融合方法,基于交叉梯度约束对地震反射波法与直流电阻率法超前探测的数据进行融合计算,提高了超前探测的探测精度。
为了实现上述目的,本发明采用如下技术方案:
一种地震和直流电法巷道超前探测的数据融合方法,具体包括以下步骤:
(1)首先进行地震反射波法速度的单一反演,得到速度数据集Dms;同时进行直流电阻率法电阻率的单一反演,得到电阻率数据集Dmr
(2)选取地震反射波法和直流电阻率法探测的公共区域为数据融合区域,得到数据融合区域内的速度数据集D'ms和电阻率数据集D'mr
(3)对数据融合区域内的速度数据集D'ms和电阻率数据集D'mr进行数据的归一化处理;
(4)对归一化处理后的数据集计算交叉梯度值t(x,z),对交叉梯度值t(x,z)进行离散化,得到离散化后的数据ti,j
(5)基于离散化后的数据ti,j和步骤(2)中数据融合区域内的速度数据集D'ms计算速度更新量V″i,j;基于离散化后的数据ti,j和步骤(2)中数据融合区域内的电阻率数据集D'mr计算电阻率更新量R″i,j
(6)基于速度更新量V″i,j构成新的速度数据集D″ms,并基于新的速度数据集D″ms进行速度成图;基于电阻率更新量R″i,j构成新的电阻率数据集D″mr,并基于新的电阻率数据集D″mr进行电阻率成图。
优选的,所述步骤(3)中,利用Min-max标准法对速度数据集D'ms和电阻率数据集D'mr进行数据归一化处理,具体包括:
设A为速度数据集或者电阻率数据集,minA和maxA分别为速度数据集或电阻率数据集的最大值和最小值,将A中每个原始数据x通过Min-max标准化映射成在区间[0,1]中的新数据x’,其公式为:新数据x’=(原始数据x-最小值minA)/(最大值maxA-最小值minA),所得新数据构成归一化后的速度数据集和电阻率数据集。
优选的,所述步骤(4)具体包括如下步骤:
a、基于二维速度模型与电阻率模型,定义两者之间的交叉梯度值为
Figure BDA0002365220910000021
式中
Figure BDA0002365220910000022
为梯度运算符,ms(x,z)表示地震模型速度数据,mr(x,z)表示电法模型电阻率数据,t(x,z)为定义的交叉梯度值;
b、采用中心差分法,对a中交叉梯度值t(x,z)进行数据离散化,得到:
Figure BDA0002365220910000023
式中i,j=1,2,...,N
其中,ms(i,j)表示当前网格节点速度数值,ms(i,j+1)表示当前网格节点沿垂直方向的网格节点数加1位置的速度数值,ms(i,j-1)表示当前网格节点沿垂直方向的网格节点数减1位置的速度数值,mr(i,j)表示当前网格节点电阻率数值,mr(i+1,j)表示当前网格节点沿水平方向的网格节点数加1位置的电阻率数值,mr(i-1,j)表示当前网格节点沿水平方向的网格节点数减1位置的电阻率数值;Δxj表示沿水平方向每次计算时网格节点序号相差的个数,即当前网格节点左右相邻的网格节点数值的差值,Δzi表示沿垂直方向每次计算时网格节点序号相差的个数,即当前网格节点上下相邻的网格节点数值的差值。
优选的,所述步骤(4)中,速度更新量
Figure BDA0002365220910000031
其中i,j=1,2,3……N;式中,αs为地震参数的权重因子,
Figure BDA0002365220910000032
表示数据融合区域的速度集数据;电阻率更新量
Figure BDA0002365220910000033
其中i,j=1,2,3……N;式中,βr为电阻率参数的权重因子,
Figure BDA0002365220910000034
表示数据融合区域的电阻率集数据。
优选的,在所述步骤(1)之前还包括:初始化速度模型和初始化电阻率模型。
经由上述的技术方案可知,与现有技术相比,本发明公开提供了一种地震和直流电法巷道超前探测的数据融合方法,首先进行地震反射波法速度的单一反演和直流电阻率法电阻率的单一反演,得到速度数据集和电阻率数据集;选择地震反射波法和直流电阻率法探测的公共区域为数据融合区域,对数据融合区域内的速度数据集和电阻率数据集的数据进行归一化处理;对归一化处理后的数据集计算交叉梯度值并进行离散化;基于离散化的交叉梯度值和数据融合区域内的速度数据集和电阻率数据集计算速度更新量和电阻率更新量;基于速度更新量和电阻率更新量得到新的速度数据集和电阻率数据集,并基于新的速度数据集进行速度成像,基于新的电阻率数据集进行电阻率成像。本发明的技术方案基于交叉梯度约束对地震反射波法与直流电阻率法超前探测的数据进行融合计算,提高了超前探测的探测精度。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案,下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本发明的实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据提供的附图获得其他的附图。
图1为本发明提供的地震和直流电法巷道超前探测的数据融合方法的流程图。
图2为本发明提供的交叉梯度离散化网格示意图。
图3为本发明提供的数据融合区域示意图。
图4a为本发明提供的模拟数据单一反演速度剖面图。
图4b为本发明提供的模拟数据数据融合后速度剖面图。
图5a为本发明提供的模拟数据单一反演电阻率剖面图。
图5b为本发明模拟数据数据融合后电阻率剖面图。
其中,图3中的1表示巷道掌子面;2表示巷道;3表示电法探测区域;4表示地震探测区域;5表示数据融合区域。
具体实施方式
下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅仅是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
参见图1,本发明实施例公开了一种地震和直流电法巷道超前探测的数据融合方法,具体步骤如下:
(1)参见图4a和图5a,首先进行地震反射波法速度的单一反演,得到速度数据集Dms;同时进行直流电阻率法电阻率的单一反演,得到电阻率数据集Dmr
(2)参见图3,选取地震反射波法和直流电阻率法探测的公共区域为数据融合区域,得到数据融合区域内的速度数据集D'ms和电阻率数据集D'mr
(3)对数据融合区域内的速度数据集D'ms和电阻率数据集D'mr进行数据的归一化处理;
(4)参见图2,对归一化处理后的数据集计算交叉梯度值t(x,z),对交叉梯度值t(x,z)进行离散化,得到离散化后的数据ti,j
(5)基于离散化后的数据ti,j和步骤(2)中数据融合区域内的速度数据集D'ms计算速度更新量V″i,j;基于离散化后的数据ti,j和步骤(2)中数据融合区域内的电阻率数据集D'mr计算电阻率更新量R″i,j
(6)参见图4b和图5b,基于速度更新量V″i,j构成新的速度数据集D″ms,并基于新的速度数据集D″ms进行速度成图;基于电阻率更新量R″i,j构成新的电阻率数据集D″mr,并基于新的电阻率数据集D″mr进行电阻率成图。
优选的,步骤(3)中,利用Min-max标准法对速度数据集D'ms和电阻率数据集D'mr进行数据归一化处理,具体包括:
设A为速度数据集或者电阻率数据集,minA和maxA分别为速度数据集或电阻率数据集的最大值和最小值,将A中每个原始数据x通过Min-max标准化映射成在区间[0,1]中的新数据x’,其公式为:新数据x’=(原始数据x-最小值minA)/(最大值maxA-最小值minA),所得新数据构成归一化后的速度数据集和电阻率数据集。
优选的,参见图2,步骤(4)具体包括如下步骤:
a、基于二维速度模型与电阻率模型,定义两者之间的交叉梯度值为
Figure BDA0002365220910000051
式中
Figure BDA0002365220910000052
为梯度运算符,ms(x,z)表示地震模型速度数据,mr(x,z)表示电法模型电阻率数据,t(x,z)为定义的交叉梯度值;
b、采用中心差分法,对a中交叉梯度值t(x,z)进行数据离散化,得到:
Figure BDA0002365220910000053
式中i,j=1,2,...,N
其中,ms(i,j)表示当前网格节点速度数值,ms(i,j+1)表示当前网格节点沿垂直方向的网格节点数加1位置的速度数值,ms(i,j-1)表示当前网格节点沿垂直方向的网格节点数减1位置的速度数值,mr(i,j)表示当前网格节点电阻率数值,mr(i+1,j)表示当前网格节点沿水平方向的网格节点数加1位置的电阻率数值,mr(i-1,j)表示当前网格节点沿水平方向的网格节点数减1位置的电阻率数值;Δxj表示沿水平方向每次计算时网格节点序号相差的个数,即当前网格节点左右相邻的网格节点数值的差值,Δzi表示沿垂直方向每次计算时网格节点序号相差的个数,即当前网格节点上下相邻的网格节点数值的差值。
优选的,步骤(4)中,速度更新量
Figure BDA0002365220910000054
其中i,j=1,2,3……N;式中,αs为地震参数的权重因子,
Figure BDA0002365220910000055
表示数据融合区域的速度集数据;电阻率更新量
Figure BDA0002365220910000061
其中i,j=1,2,3……N;式中,βr为电阻率参数的权重因子,
Figure BDA0002365220910000062
表示数据融合区域的电阻率集数据。
优选的,在步骤(1)之前还包括:初始化速度模型和初始化电阻率模型。
经由上述的实施例可知,与现有技术相比,本发明公开提供了一种地震和直流电法巷道超前探测的数据融合方法,首先进行地震反射波法速度的单一反演和直流电阻率法电阻率的单一反演,得到速度数据集和电阻率数据集;选择地震反射波法和直流电阻率法探测的公共区域为数据融合区域,对数据融合区域内的速度数据集和电阻率数据集的数据进行归一化处理;对归一化处理后的数据集计算交叉梯度值并进行离散化;基于离散化的交叉梯度值和数据融合区域内的速度数据集和电阻率数据集计算速度更新量和电阻率更新量;基于速度更新量和电阻率更新量得到新的速度数据集和电阻率数据集,并基于新的速度数据集进行速度成像,基于新的电阻率数据集进行电阻率成像。本发明的实施例基于交叉梯度约束对地震反射波法与直流电阻率法超前探测的数据进行融合计算,提高了超前探测的探测精度。
本说明书中各个实施例采用递进的方式描述,每个实施例重点说明的都是与其他实施例的不同之处,各个实施例之间相同相似部分互相参见即可。对于实施例公开的装置而言,由于其与实施例公开的方法相对应,所以描述的比较简单,相关之处参见方法部分说明即可。
对所公开的实施例的上述说明,使本领域专业技术人员能够实现或使用本发明。对这些实施例的多种修改对本领域的专业技术人员来说将是显而易见的,本文中所定义的一般原理可以在不脱离本发明的精神或范围的情况下,在其它实施例中实现。因此,本发明将不会被限制于本文所示的这些实施例,而是要符合与本文所公开的原理和新颖特点相一致的最宽的范围。

Claims (5)

1.一种地震和直流电法巷道超前探测的数据融合方法,其特征在于,包括以下步骤:
(1)首先进行地震反射波法速度的单一反演,得到速度数据集Dms;同时进行直流电阻率法电阻率的单一反演,得到电阻率数据集Dmr
(2)选取地震反射波法和直流电阻率法探测的公共区域为数据融合区域,得到数据融合区域内的速度数据集D'ms和电阻率数据集D'mr
(3)对数据融合区域内的速度数据集D'ms和电阻率数据集D'mr进行数据的归一化处理;
(4)对归一化处理后的数据集计算交叉梯度值t(x,z),对交叉梯度值t(x,z)进行离散化,得到离散化后的数据ti,j
(5)基于离散化后的数据ti,j和步骤(2)中数据融合区域内的速度数据集D'ms计算速度更新量V″i,j;基于离散化后的数据ti,j和步骤(2)中数据融合区域内的电阻率数据集D'mr计算电阻率更新量R″i,j
(6)基于速度更新量V″i,j构成新的速度数据集D″ms,并基于新的速度数据集D″ms进行速度成图;基于电阻率更新量R″i,j构成新的电阻率数据集D″mr,并基于新的电阻率数据集D″mr进行电阻率成图。
2.根据权利要求1所述的一种地震和直流电法巷道超前探测的数据融合方法,其特征在于,所述步骤(3)中,利用Min-max标准法对速度数据集D'ms和电阻率数据集D'mr进行数据归一化处理,具体包括:
设A为速度数据集或者电阻率数据集,minA和maxA分别为速度数据集或电阻率数据集的最大值和最小值,将A中每个原始数据x通过Min-max标准化映射成在区间[0,1]中的新数据x’,其公式为:新数据x’=(原始数据x-最小值minA)/(最大值maxA-最小值minA),所得新数据构成归一化后的速度数据集和电阻率数据集。
3.根据权利要求1所述的一种地震和直流电法巷道超前探测的数据融合方法,其特征在于,所述步骤(4)具体包括如下步骤:
a、基于二维速度模型与电阻率模型,定义两者之间的交叉梯度值为
Figure FDA0002365220900000011
式中
Figure FDA0002365220900000012
为梯度运算符,ms(x,z)表示地震模型速度数据,mr(x,z)表示电法模型电阻率数据,t(x,z)为定义的交叉梯度值;
b、采用中心差分法,对a中交叉梯度值t(x,z)进行数据离散化,得到:
Figure FDA0002365220900000021
其中,ms(i,j)表示当前网格节点速度数值,ms(i,j+1)表示当前网格节点沿垂直方向的网格节点数加1位置的速度数值,ms(i,j-1)表示当前网格节点沿垂直方向的网格节点数减1位置的速度数值,mr(i,j)表示当前网格节点电阻率数值,mr(i+1,j)表示当前网格节点沿水平方向的网格节点数加1位置的电阻率数值,mr(i-1,j)表示当前网格节点沿水平方向的网格节点数减1位置的电阻率数值;Δxj表示沿水平方向每次计算时网格节点序号相差的个数,即当前网格节点左右相邻的网格节点数值的差值,Δzi表示沿垂直方向每次计算时网格节点序号相差的个数,即当前网格节点上下相邻的网格节点数值的差值。
4.根据权利要求1所述的一种地震和直流电法巷道超前探测的数据融合方法,其特征在于,所述步骤(4)中,速度更新量
Figure FDA0002365220900000025
其中i,j=1,2,3……N;式中,αs为地震参数的权重因子,
Figure FDA0002365220900000022
表示数据融合区域的速度集数据;电阻率更新量
Figure FDA0002365220900000023
其中i,j=1,2,3……N;式中,βr为电阻率参数的权重因子,
Figure FDA0002365220900000024
表示数据融合区域的电阻率集数据。
5.根据权利要求1所述的一种地震和直流电法巷道超前探测的数据融合方法,其特征在于,在所述步骤(1)之前还包括:初始化速度模型和初始化电阻率模型。
CN202010033559.XA 2020-01-13 2020-01-13 一种地震和直流电法巷道超前探测的数据融合方法 Active CN111123404B (zh)

Priority Applications (2)

Application Number Priority Date Filing Date Title
CN202010033559.XA CN111123404B (zh) 2020-01-13 2020-01-13 一种地震和直流电法巷道超前探测的数据融合方法
ZA2020/02934A ZA202002934B (en) 2020-01-13 2020-05-20 Method for fusing data of seismic and dc method advanced detection for roadways

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010033559.XA CN111123404B (zh) 2020-01-13 2020-01-13 一种地震和直流电法巷道超前探测的数据融合方法

Publications (2)

Publication Number Publication Date
CN111123404A true CN111123404A (zh) 2020-05-08
CN111123404B CN111123404B (zh) 2023-01-13

Family

ID=70489246

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010033559.XA Active CN111123404B (zh) 2020-01-13 2020-01-13 一种地震和直流电法巷道超前探测的数据融合方法

Country Status (2)

Country Link
CN (1) CN111123404B (zh)
ZA (1) ZA202002934B (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111722292A (zh) * 2020-06-03 2020-09-29 山东大学 一种三维电阻率层析成像方法及系统

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105785476A (zh) * 2016-05-17 2016-07-20 浙江水利水电学院 一种高密度电阻率成像方法
CN107765337A (zh) * 2016-08-19 2018-03-06 中国石油化工股份有限公司 电法与地震同步联合反演方法及系统

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105785476A (zh) * 2016-05-17 2016-07-20 浙江水利水电学院 一种高密度电阻率成像方法
CN107765337A (zh) * 2016-08-19 2018-03-06 中国石油化工股份有限公司 电法与地震同步联合反演方法及系统

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
NINFA L. BENNINGTON 等: "Joint Inversion of Seismic and Magnetotelluric Data in the Parkfield Region of California Using the Normalized Cross-Gradient Constraint", 《PURE AND APPLIED GEOPHYSICS》 *
王俊 等: "交叉梯度理论及其在地球物理联合反演中的应用", 《地球物理学进展》 *
高级 等: "基于交叉梯度交替结构约束的二维地震走时与全通道直流电阻率联合反演", 《地球物理学报》 *

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111722292A (zh) * 2020-06-03 2020-09-29 山东大学 一种三维电阻率层析成像方法及系统
CN111722292B (zh) * 2020-06-03 2021-06-01 山东大学 一种三维电阻率层析成像方法及系统
WO2021243967A1 (zh) * 2020-06-03 2021-12-09 山东大学 一种三维电阻率层析成像方法及系统
US11960048B2 (en) 2020-06-03 2024-04-16 Shandong University Three-dimensional electrical resistivity tomography method and system

Also Published As

Publication number Publication date
ZA202002934B (en) 2021-03-31
CN111123404B (zh) 2023-01-13

Similar Documents

Publication Publication Date Title
CN112381937B (zh) 一种基于钻孔和复杂地质剖面的多源地质数据耦合建模方法
GB2532590B (en) Simulating fluid flow using a stairstepped grid to represent a geological fault
CN111221048B (zh) 基于跨孔电阻率ct多尺度反演的孤石边界识别与成像方法
KR20080023946A (ko) 오일러 해를 이용한 지하공동의 3차원 중력역산 방법 및이를 이용한 3차원 영상화 방법
CN110286416B (zh) 一种基于物性函数的快速二维密度反演方法
CN103592697B (zh) 隧道等地下工程先验形态信息约束的联合超前预报方法
CN108760825A (zh) 适应于城市环境的任意布孔的三维跨孔电阻率ct成像方法
CN102877828A (zh) 一种三维多井联合井地ct成像方法
CN111123404B (zh) 一种地震和直流电法巷道超前探测的数据融合方法
CN113433595B (zh) 基于自然电场隧道裂隙水的超前预报方法
CN103217715B (zh) 多尺度规则网格层析反演静校正方法
KR101157792B1 (ko) 구조화 지수를 이용한 3차원 지질구조 해석방법
Ji et al. An automated method to build 3D multi-scale geological models for engineering sedimentary layers with stratum lenses
CN117372639B (zh) 一种稀疏骨干剖面结合任意位置辅助面的三维地质-地球物理建模方法
CN112666612B (zh) 基于禁忌搜索的大地电磁二维反演方法
WO2024078134A1 (zh) 基于多参数约束和结构校正的随掘巷道全波形反演方法
CN114491774B (zh) 一种深部背斜构造及地层结构三维数值模型构建方法
CN112882124B (zh) 一种三维铀成矿有利部位的圈定方法
CN106970426B (zh) 一种基于测线差分与主成分分析的航空电磁数据调平方法
CN106405648B (zh) 绕射波的成像方法和装置
CN110687591B (zh) 基于先验数据的波形匹配确定煤层及围岩物性参数的方法
CN106405672A (zh) 一种tbm施工隧道三维激发极化超前预报干扰去除方法
CN104834005B (zh) 一种确定地质界面的方法
CN116522688B (zh) 一种井控多信息融合工程地质建模方法和装置
CN114706126B (zh) 一种地质认识约束下的地震地质体定量雕刻方法及系统

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