CN112799141B - 一种快速二维核磁共振测井信号处理及t1t2谱反演方法 - Google Patents

一种快速二维核磁共振测井信号处理及t1t2谱反演方法 Download PDF

Info

Publication number
CN112799141B
CN112799141B CN202110010570.9A CN202110010570A CN112799141B CN 112799141 B CN112799141 B CN 112799141B CN 202110010570 A CN202110010570 A CN 202110010570A CN 112799141 B CN112799141 B CN 112799141B
Authority
CN
China
Prior art keywords
echo signals
echo
filtering
signal
wavelet
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.)
Active
Application number
CN202110010570.9A
Other languages
English (en)
Other versions
CN112799141A (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.)
China National Petroleum Corp
China Petroleum Logging Co Ltd
Original Assignee
China National Petroleum Corp
China Petroleum Logging Co Ltd
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 China National Petroleum Corp, China Petroleum Logging Co Ltd filed Critical China National Petroleum Corp
Priority to CN202110010570.9A priority Critical patent/CN112799141B/zh
Publication of CN112799141A publication Critical patent/CN112799141A/zh
Application granted granted Critical
Publication of CN112799141B publication Critical patent/CN112799141B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

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
    • 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/18—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation specially adapted for well-logging
    • G01V3/32—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation specially adapted for well-logging operating with electron or nuclear magnetic resonance
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
    • 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/30—Assessment of water resources

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Mathematical Physics (AREA)
  • Geophysics (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • Remote Sensing (AREA)
  • Data Mining & Analysis (AREA)
  • Computational Mathematics (AREA)
  • Geology (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Environmental & Geological Engineering (AREA)
  • Pure & Applied Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computing Systems (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Algebra (AREA)
  • High Energy & Nuclear Physics (AREA)
  • Magnetic Resonance Imaging Apparatus (AREA)

Abstract

本发明提供一种快速二维核磁共振测井信号处理及T1T2谱反演方法，包括：S1，从每个组分的回波信号中按照设定规则抽取部分回波信号，得到抽取的回波信号；S2，根据核磁共振采集模式并利用抽取的回波信号构建系数矩阵A；S3，对系数矩阵A进行T1T2谱求解；S4，输出T1T2谱。本发明方法针对核磁共振回波信号特点，通过信号抽取，实现了大矩阵的压缩，缩短了反演时间。

Description

一种快速二维核磁共振测井信号处理及T1T2谱反演方法
技术领域
本发明属于地球物理测井领域，涉及一种快速二维核磁共振测井信号处理及T1T2谱反演方法。
背景技术
一维核磁共振测井在油气识别和孔隙结构评价中已经得到广泛应用，但在实际应用中也存在一系列不足，比如多相流体存在时流体识别率偏低等问题凸显。由此发展了二维核磁共振测井方法，通过二维核磁信号处理得到二维谱来实现对于多相流体的识别。
采用不同等待时间进行回波信号测量，对回波信号构建方程式，通过反演可以得到T1T2谱，对应理论图谱或岩心实验图谱可以完成流体的准确识别。
目前二维核磁测井资料的处理技术并不成熟，较一维核磁回波信号相比，二维核磁回波信号数据量成倍增长，在处理效率上很低，无法投入正常生产应用。
发明内容
本发明的主要目的为实现一种快速二维核磁共振测井信号处理及T1T2谱反演方法，解决二维核磁反演时间长、效率低问题。
本发明是通过以下技术方案来实现：
一种快速二维核磁共振测井信号处理及T1T2谱反演方法，包括：
S1，从每个组分的回波信号中按照设定规则抽取部分回波信号，得到抽取的回波信号；
S2，根据核磁共振采集模式并利用抽取的回波信号构建系数矩阵A；
S3，对系数矩阵A进行T1T2谱求解；
S4，输出T1T2谱。
优选的，S1之前，对回波信号进行滤波处理，S1中从每个组分的滤波处理后的回波信号中按照设定规则抽取部分回波信号，得到抽取的回波信号。
进一步的，所述滤波为小波阀值滤波。
进一步的，小波阀值滤波采用小波软阀值滤波，并采用heursure阀值选择规则。
进一步的，滤波处理具体为：
S01，对回波信号进行小波分解，得到不同分解层数的小波，进行小波滤波计算，得到不同分解层数下的小波滤波结果；
S02，计算不同分解层数下小波滤波结果的小波定量评价指标IDC；
S03，输出最大IDC值对应的小波滤波结果作为滤波处理后的回波信号。
优选的，S1中，对于每个组分，设定规则为：
(1)若回波个数ECHONUM<N，则保留全部回波信号；
(2)若回波个数ECHONUM>＝N，保留前N个点的回波信号，对于N<ECHONUM<＝N’部分的回波信号，保留回波数为奇数的回波信号，对于ECHONUM>N’部分的回波信号，保留ECHONUM/4为整数的回波信号。
进一步的，N为50，N’为150。
优选的，S2中，系数矩阵A是n×m矩阵。
进一步的，S3具体为：
S31，对系数矩阵A奇异值分解，得到两个正交矩阵U,S和对角阵V；其中U是n×m列正交矩阵，S是m×m阶正交矩阵，
Figure BDA0002884157450000021
是m阶对角矩阵，其中
Figure BDA0002884157450000022
Figure BDA0002884157450000023
是奇异值；
AX＝Y
A＝UVST
Figure BDA0002884157450000031
式中，X为向量型式的方程解，Y是抽取的回波信号；
S32，计算回波信号信噪比SNR，X在
Figure BDA0002884157450000032
处进行信号截断，保留小于
Figure BDA0002884157450000033
的奇异值，初始解表示为
Figure BDA0002884157450000034
S33，设定迭代计算最大次数inter，设定迭代步进因子alpha，精度要求eprecison；令迭代计算次数k＝1；
S34，进行第k次迭代计算：对X进行非负判断，保留X(i)＞＝0对应的系数矩阵A中的列，得到新的系数矩阵A；同时令X中小于零的项改为0.001，如果||AX-Y||2＜epricison，转至S36；如果A为空，转至S36；如果k＝inter，转至S36；
S35，计算矩阵M＝A*(A)T，计算新的X＝(M+alpha*eye(mRow))\Y，其中eye(mRow)表示阶数为mRow的单位矩阵，令k＝k+1，转S34；
S36，输出X。
进一步的，设定迭代计算最大次数inter＝500，设定迭代步进因子alpha＝0.01，精度要求eprecison＝0.0001。
与现有技术相比，本发明具有以下有益的技术效果：
本发明方法针对核磁共振回波信号特点，通过信号抽取，实现了大矩阵的压缩，缩短了反演时间。应用本发明，能够快速的反演得到核磁共振T1T2谱，为二维核磁共振流体识别打下基础，满足生产作业需求。
进一步的，现有的核磁共振测井信号强度弱，信噪比低，在低信噪比条件下直接进行T1T2谱反演的处理结果失真，无法准确识别油气。针对核磁共振测井信号信噪比低的问题，本发明采用小波阀值自动去噪方法提高了信号信噪比，最大限度克服仪器自身信噪比低带来的难题。
进一步的，一维核磁采用的多指数反演算法已经无法满足现有的二维核磁处理在精度和效率上的需求，需要改进算法。本发明改进现有SVD一维核磁反演算法，提高了反演效率和精度，实现技术突破，提升二维核磁T1T2谱反演效果。
附图说明
图1基本流程图；
图2 T1T2模型；
图3回波信号图；
图4 T1T2谱反演结果；
图5实际井资料处理图。
具体实施方式
下面结合具体的实施例对本发明做进一步的详细说明，所述是对本发明的解释而不是限定。
如图1所示，本发明所述的快速二维核磁共振测井信号处理及T1T2谱反演方法，包括如下步骤：
步骤1，对核磁共振的回波信号进行小波阀值自动滤波处理，得到小波滤波后的回波信号；
步骤2，对小波滤波后的回波信号进行信号抽取，即从每个组分小波滤波后的回波信号中按照设定规则抽取部分回波信号，得到抽取的回波信号；
步骤3，根据核磁共振采集模式并利用抽取的回波信号构建系数矩阵A；
步骤4，采用改进SVD方法对系数矩阵A进行T1T2谱求解；
步骤5，输出T1T2谱；
步骤1包含：
步骤11：选择小波软阀值滤波，采用heursure阀值选择规则，对回波信号进行小波分解，小波分解层数分别设置为N1～N2层，进行小波滤波计算，得到不同分解层数下的小波滤波结果。
步骤12：计算小波滤波效果评价指标，包括均方误差RMSE，信噪比SNR，平滑度指标r以及小波定量评价指标IDC，IDC＝SNR/(RMSE*S)。
步骤13：选择IDC值最大时的分解层数，作为小波滤波最佳效果。
步骤14：输出IDC值最大时的分解层数下的小波滤波结果作为小波滤波后的回波信号。
步骤2包含：
针对回波信号ECHO特点进行信号抽取，给每个组分回波信号进行编号，按照设定原则进行信号抽取。
(1)回波个数ECHONUM<50的组分保留全部回波信号。
(2)当回波个数ECHONUM>＝50时，保留前50个点回波信号。50<ECHONUM<＝150部分，保留回波数为奇数的回波信号，即ECHO51，ECHO53，…ECHO149。ECHONUM>150部分，保留ECHONUM/4为整数的回波信号，即ECHO152，ECHO156，…ECHOECHONUM/4。
步骤4包含：
步骤41：对系数矩阵A奇异值分解，得到两个正交矩阵U,S和对角阵V，其中A是n×m矩阵，U是n×m列正交矩阵，S是m×m阶正交矩阵，
Figure BDA0002884157450000051
是m阶对角矩阵，其中
Figure BDA0002884157450000052
Figure BDA0002884157450000053
是奇异值。
AX＝Y
A＝UVST
Figure BDA0002884157450000061
式中，A是根据核磁共振采集模式构建的系数矩阵，X为向量型式的方程解，Y是抽取的回波信号。
步骤42：计算回波信号信噪比SNR，在
Figure BDA0002884157450000062
处进行信号截断，保留小于
Figure BDA0002884157450000063
的奇异值，初始解可以表示为：
Figure BDA0002884157450000064
步骤43：设定迭代计算最大次数inter，设定迭代步进因子alpha，精度要求eprecison；令迭代计算次数k＝1；
步骤44：进行第k次迭代计算，对X进行非负判断，保留方程解中X(i)＞＝0时对应的系数矩阵A中的列，得到新的系数矩阵A；同时令X中小于零的项改为0.001，如果||AX-Y||2＜epricison，转步骤46；如果A为空，转步骤46；如果k＝inter，转步骤46。
步骤45：计算矩阵M＝A*(A)T，计算新的X＝(M+alpha*eye(mRow))\Y，其中eye(mRow)表示阶数为mRow的单位矩阵，令k＝k+1，转步骤44。
步骤46：输出X。
实施例
本发明为一种快速二维核磁共振信号处理及T1T2谱反演方法，包括如下步骤：
步骤一，按照图2设定正演模型，设计采集观测模式信息：
等待时间TW＝[8028，2675，10，30，100，300 80 1000]；回波间隔TE＝[1.2，1.2，0.6，0.6，0.6，0.6 0.6 0.6]；回波个数NE＝[400，400，10，10，10，1020,80]；根据设定的正演模型构建回波信号，往回波信号中加入噪声信号，信噪比为30，得到原始回波信号，如图3所示。
步骤二，对原始回波信号进行小波阀值自动滤波处理。选择小波软阀值滤波，采用heursure阀值选择规则，小波分解层数分别设置为4～11层，进行小波滤波计算，得到不同分解层数下的小波滤波结果。计算不同分解层数小波滤波效果评价指标，包括均方误差RMSE，信噪比SNR，平滑度指标r以及小波定量评价指标IDC，IDC＝SNR/(RMSE*S)。选择IDC值最大时的分解层数，作为小波滤波最佳效果。输出小波滤波结果。
步骤三：针对回波信号ECHO特点对小波滤波后的回波信号进行信号抽取，给每个组分回波信号进行编号，按照固定原则进行信号抽取。
(1)回波个数ECHONUM<50的组分保留全部回波信号。
(2)当ECHONUM>＝50时，保留前50个点回波。
50<ECHONUM<＝150部分，保留回波数为奇数部分,即ECHO51，ECHO53，…ECHO149。ECHONUM>150部分，保留ECHONUM/4为整数部分,即ECHO152，ECHO156，…ECHOECHONUM/4。
步骤四：根据采集模式构建系数矩阵A
步骤五：设定迭代计算最大次数inter＝500，设定迭代步进因子alpha＝0.01，精度要求eprecison＝0.0001，进行T1T2谱反演得到输出结果，如图4。
步骤六：绘图输出T1T2谱。
该次单点反演处理耗时0.9s，实际地层下采用0.1m采样间隔，五个采样点做一次T1T2处理，100m井资料处理耗时仅需3分钟，满足实际井资料处理需求。采用该方法实现了实际井资料的处理，见图5。
本发明并不局限于上述实施例，在本发明公开的技术方案的基础上，本领域的技术人员根据所公开的技术内容，不需要创造性的劳动就可以对其中的一些技术特征作出一些替换和变形，这些替换和变形均在本发明的保护范围内。

Claims (7)

1.一种快速二维核磁共振测井信号处理及T1T2谱反演方法，其特征在于，包括：
S1，从每个组分的回波信号中按照设定规则抽取部分回波信号，得到抽取的回波信号；
S2，根据核磁共振采集模式并利用抽取的回波信号构建系数矩阵A；
S3，对系数矩阵A进行T1T2谱求解；
S4，输出T1T2谱；
S1中，对于每个组分，设定规则为：
(1)若回波个数ECHONUM<N，则保留全部回波信号；
(2)若回波个数ECHONUM>＝N，保留前N个点的回波信号，对于N<ECHONUM<＝N’部分的回波信号，保留回波数为奇数的回波信号，对于ECHONUM>N’部分的回波信号，保留ECHONUM/4为整数的回波信号；
S2中，系数矩阵A是n×m矩阵；
S3具体为：
S31，对系数矩阵A奇异值分解，得到两个正交矩阵U,S和对角阵V；其中U是n×m列正交矩阵，S是m×m阶正交矩阵，
Figure QLYQS_1
是m阶对角矩阵，其中
Figure QLYQS_2
Figure QLYQS_3
是奇异值；
AX＝Y
A＝UVST
Figure QLYQS_4
式中，X为向量型式的方程解，Y是抽取的回波信号；
S32，计算回波信号信噪比SNR，X在
Figure QLYQS_5
处进行信号截断，保留小于
Figure QLYQS_6
的奇异值，初始解表示为
Figure QLYQS_7
S33，设定迭代计算最大次数inter，设定迭代步进因子alpha，精度要求eprecison；令迭代计算次数k＝1；
S34，进行第k次迭代计算：对X进行非负判断，保留X(i)＞＝0对应的系数矩阵A中的列，得到新的系数矩阵A；同时令X中小于零的项改为0.001，如果||AX-Y||2＜epricison，转至S36；如果A为空，转至S36；如果k＝inter，转至S36；
S35，计算矩阵M＝A*(A)T，计算新的X＝(M+alpha*eye(mRow))\Y，其中eye(mRow)表示阶数为mRow的单位矩阵，令k＝k+1，转S34；
S36，输出X。
2.根据权利要求1所述的快速二维核磁共振测井信号处理及T1T2谱反演方法，其特征在于，S1之前，对回波信号进行滤波处理，S1中从每个组分的滤波处理后的回波信号中按照设定规则抽取部分回波信号，得到抽取的回波信号。
3.根据权利要求2所述的快速二维核磁共振测井信号处理及T1T2谱反演方法，其特征在于，所述滤波为小波阈值滤波。
4.根据权利要求3所述的快速二维核磁共振测井信号处理及T1T2谱反演方法，其特征在于，小波阈值滤波采用小波软阈值滤波，并采用heursure阈值选择规则。
5.根据权利要求4所述的快速二维核磁共振测井信号处理及T1T2谱反演方法，其特征在于，滤波处理具体为：
S01，对回波信号进行小波分解，得到不同分解层数的小波，进行小波滤波计算，得到不同分解层数下的小波滤波结果；
S02，计算不同分解层数下小波滤波结果的小波定量评价指标IDC；
S03，输出最大IDC值对应的小波滤波结果作为滤波处理后的回波信号。
6.根据权利要求1所述的快速二维核磁共振测井信号处理及T1T2谱反演方法，其特征在于，N为50，N’为150。
7.根据权利要求1所述的快速二维核磁共振测井信号处理及T1T2谱反演方法，其特征在于，设定迭代计算最大次数inter＝500，设定迭代步进因子alpha＝0.01，精度要求eprecison＝0.0001。
CN202110010570.9A 2021-01-05 2021-01-05 一种快速二维核磁共振测井信号处理及t1t2谱反演方法 Active CN112799141B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110010570.9A CN112799141B (zh) 2021-01-05 2021-01-05 一种快速二维核磁共振测井信号处理及t1t2谱反演方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110010570.9A CN112799141B (zh) 2021-01-05 2021-01-05 一种快速二维核磁共振测井信号处理及t1t2谱反演方法

Publications (2)

Publication Number Publication Date
CN112799141A CN112799141A (zh) 2021-05-14
CN112799141B true CN112799141B (zh) 2023-05-26

Family

ID=75808364

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110010570.9A Active CN112799141B (zh) 2021-01-05 2021-01-05 一种快速二维核磁共振测井信号处理及t1t2谱反演方法

Country Status (1)

Country Link
CN (1) CN112799141B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115201925B (zh) * 2022-05-26 2024-08-20 长江大学 一种基于核磁共振二维谱识别方法、装置、设备及介质

Family Cites Families (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102565865B (zh) * 2011-12-07 2014-02-12 中国石油大学(北京) 降噪核磁共振测井回波信号的获得方法及装置
CN103116148B (zh) * 2013-01-30 2015-04-01 上海理工大学 一种核磁共振二维谱反演的方法
CN106156505B (zh) * 2016-07-05 2019-07-23 中国科学技术大学 一种基于正交匹配追踪算法的核磁共振t2谱反演方法
CN110320227B (zh) * 2018-03-29 2022-08-05 中国石油化工股份有限公司 一种二维核磁共振d-t2谱反演方法及装置
CN110454153B (zh) * 2019-08-20 2022-11-04 中国海洋石油集团有限公司 一种核磁共振测井弛豫反演方法
CN111827968B (zh) * 2020-07-15 2023-08-18 长江大学 一种基于核磁共振测井的储层非均质性评价方法及装置

Also Published As

Publication number Publication date
CN112799141A (zh) 2021-05-14

Similar Documents

Publication Publication Date Title
CN111340046A (zh) 基于特征金字塔网络和通道注意力的视觉显著性检测方法
CN103871058B (zh) 基于压缩采样矩阵分解的红外小目标检测方法
CN103336305B (zh) 一种基于灰色理论划分致密砂岩储层岩石物理相的方法
CN112180433B (zh) 地震初至波拾取方法及装置
CN106491131A (zh) 一种磁共振的动态成像方法和装置
CN110490219A (zh) 一种基于纹理约束的U-net网络进行地震数据重建的方法
CN109376589A (zh) 基于卷积核筛选ssd网络的rov形变目标与小目标识别方法
CN106156505A (zh) 一种基于正交匹配追踪算法的核磁共振t2谱反演方法
CN105044794B (zh) 一种核磁共振回波数据的压缩方法及装置
CN108345033A (zh) 一种微地震信号时频域初至检测方法
CN111580181A (zh) 一种基于多场多特征信息融合的导水陷落柱识别方法
CN112799141B (zh) 一种快速二维核磁共振测井信号处理及t1t2谱反演方法
CN112991483B (zh) 一种非局部低秩约束的自校准并行磁共振成像重构方法
CN114138919B (zh) 一种基于非局部注意力卷积神经网络的地震数据重建方法
CN116776734A (zh) 基于物理约束神经网络的地震速度反演方法
CN111461988A (zh) 一种基于多任务学习的地震速度模型超分辨率技术
CN117114235A (zh) 一种河流生态输水后植被要素的动态变化评价方法
CN108596383B (zh) 储层分类的方法及装置
CN114969994B (zh) 基于单测点多向数据融合的rv减速器故障诊断方法及系统
CN115393733A (zh) 一种基于深度学习的水体自动识别方法及系统
CN108416816B (zh) 多维核磁共振测井数据的压缩处理方法和装置
CN114942473A (zh) 基于注意力门神经网络的叠前地震速度反演方法
CN113869289A (zh) 基于熵的多通道舰船辐射噪声特征提取方法
WO2020224012A1 (zh) 一种基于lsqr-rsvd的二维核磁共振快速反演算法
CN109064477A (zh) 用改进的U-Net检测细胞核边缘的方法

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