CN115326019A - 一种基于光学畸变加权的外方位元素转移算法 - Google Patents

一种基于光学畸变加权的外方位元素转移算法 Download PDF

Info

Publication number
CN115326019A
CN115326019A CN202210896678.7A CN202210896678A CN115326019A CN 115326019 A CN115326019 A CN 115326019A CN 202210896678 A CN202210896678 A CN 202210896678A CN 115326019 A CN115326019 A CN 115326019A
Authority
CN
China
Prior art keywords
data
optical distortion
algorithm
station
external orientation
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
CN202210896678.7A
Other languages
English (en)
Other versions
CN115326019B (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.)
NINGBO SPECIAL EQUIPMENT INSPECTION CENTER
Original Assignee
NINGBO SPECIAL EQUIPMENT INSPECTION CENTER
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 NINGBO SPECIAL EQUIPMENT INSPECTION CENTER filed Critical NINGBO SPECIAL EQUIPMENT INSPECTION CENTER
Priority to CN202210896678.7A priority Critical patent/CN115326019B/zh
Publication of CN115326019A publication Critical patent/CN115326019A/zh
Application granted granted Critical
Publication of CN115326019B publication Critical patent/CN115326019B/zh
Priority to JP2023122769A priority patent/JP7427128B1/ja
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01CMEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C11/00Photogrammetry or videogrammetry, e.g. stereogrammetry; Photographic surveying
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01CMEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C11/00Photogrammetry or videogrammetry, e.g. stereogrammetry; Photographic surveying
    • G01C11/04Interpretation of pictures
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/80Analysis of captured images to determine intrinsic or extrinsic camera parameters, i.e. camera calibration

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Multimedia (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Theoretical Computer Science (AREA)
  • Length Measuring Devices By Optical Means (AREA)

Abstract

本发明公开了一种基于光学畸变加权的外方位元素转移算法,涉及近景摄影测量基础及图像处理领域,包括以下步骤:对移站基准的数据进行分组;对分组后的数据进行共面检测,分为合格数据组和废弃数据组;对合格数据组进行基于光学畸变加权的外方位元素转移计算。本发明的算法旨在高效、科学地利用移站过程中可供使用的冗余观测数据,达到提升移站计算结果精度的目的,从而为最终测量结果的精度提供保障。

Description

一种基于光学畸变加权的外方位元素转移算法
技术领域
本发明涉及近景摄影测量基础及图像处理领域,具体地指一种基于光学畸变加权的外方位元素转移算法。
背景技术
对于摄影测量系统,最新的技术为无标定摄影测量,该技术采用两台前站相机、一套后站双目相机以及与之配套的室内标定场组成。在实验室内,通过标定场建立一个整体的坐标系,实现后方双目的外方位元素标定,以及前方相机一个机位外放元素及对应定位关键点的标定后,即可脱离标定场,调整前方相机以合适距离、角度拍摄被测物体,通过后站双目相机重新测定各前站相机移站关键点的物方坐标,实现前方相机外方位元素的迁移计算,从而实现无标定的远距离测量。这一系统工作流程可主要区分为两个主要部分:
(1)标定阶段:这一阶段的主要工作位标定后站双目相机的外方位元素,测定前站相机一个机位的外方位元素以及与之相对应的定位关键点;
(2)测量阶段:这一阶段将由两部前站相机拍摄待测物体,并由后方相机拍摄前站相机,通过对各前站相机定位关键点的再测定,结合移站算法,实现前站相机外方位元素的计算,从而实现以此为基础,不依赖现场控制点,完成测量目标点物方坐标的解算。
由于测量误差的不可避免以及随机性特征,各用于移站操作的观测点的可信度不可避免地存在着差距,在平差过程中采用等权算法对各测量点的数据加以利用显然存在着理论缺陷,不利于测量结果精度的最大化。
发明内容
有鉴于此,本发明提供一种基于光学畸变加权的外方位元素转移算法,用于现有最新无标定摄影测量。本发明的算法旨在高效、科学地利用移站过程中可供使用的冗余观测数据,达到提升移站计算结果精度的目的,从而为最终测量结果的精度提供保障。
为了实现上述目的,本发明采用如下技术方案:
一种基于光学畸变加权的外方位元素转移算法,包括以下步骤:
对移站基准的数据进行分组;
对分组后的数据进行共面检测,分为合格数据组和废弃数据组;
对合格数据组进行基于光学畸变加权的外方位元素转移计算。
可选的,对分组后的数据进行共面检测采用光学畸变模型进行数据筛选。
可选的,对移站基准的数据进行分组,具体为,设可用控制点的数量为n,单次匹配所用的点数为k,则可用的数组数量为
Figure BDA0003769118740000021
可选的,各点光学畸变导致的位移大小SG与该点数据的可信度TP的关系为:
TP∝SG-1
可选的,各点在平差公式中所占的权重Q与各点因光学畸变导致的位移大小SG的关系为:
Q∝SG-1
可选的,对合格数据组进行基于光学畸变加权的外方位元素转移计算,各数组内所包含点在各照片的像平面坐标系内的因光学畸变导致的位移大小之和SGH为:
Figure BDA0003769118740000031
其中:z=1,2,...,j,k1、k2为径向畸变参数,q1、q2为偏心畸变参数,s1、s2为薄棱镜畸变参数,δ为薄棱镜畸变参数。
经由上述的技术方案可知,与现有技术相比,本发明公开提供了一种基于光学畸变加权的外方位元素转移算法,高效、科学地利用移站过程中可供使用的冗余观测数据,达到提升移站计算结果精度的目的,从而为最终测量结果的精度提供保障。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案,下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本发明的实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据提供的附图获得其他的附图。
图1为本发明的流程示意图;
图2为本发明的A站点编号图;
图3为本发明的B站点编号图;
图4为本发明的B站点编号图;
图5为本发明的左侧像片图;
图6为本发明的右侧像片图;
图7为本发明的前站同名点示意图;
图8为本发明的前站同名点示意图。
具体实施方式
下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅仅是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
本发明实施例公开了一种基于光学畸变加权的外方位元素转移算法,如图1-8所示,首先对可作为移站基准的数据进行分组,经过共面检测后,针对各分组数据,依光学畸变模型进行评估,并依据阈值筛选堪用的数据,随后经过权重计算、参数归一化后参与旋转矩阵与平移矩阵的平差优化,最终输出移站数据。
具体过程包括以下步骤如下:
随机数据分组模式:设可用控制点的数量为n,单次匹配所用的点数为k(k≥4),则可用的数组数量为
Figure BDA0003769118740000041
其中,非共面数组的数量为
Figure BDA0003769118740000042
若用于解算的照片数量为m,各数组各特征点在各解算照片像平面坐标系内的坐标分别为:
Figure BDA0003769118740000043
与之相对应的,外方位元素转移的旋转矩阵分别为R1,R2,...,Rj,其中:
Figure BDA0003769118740000044
与各旋转矩阵所对应的平移矩阵为t1,t2,...,tj,其中:
Figure BDA0003769118740000051
基于光学畸变模型的数组加权方法:尽管在实际工程测量中,往往因为无法获取到准确的光学畸变参数而无法完成对光学畸变的精准修正,但利用建议标定手段得到的光学畸变参数仍然可以反映像平面坐标系内各点因光学畸变导致的位移大小。
显然,各点因光学畸变导致的位移大小SG与该点数据的可信度TP的关系为:
TP∝SG-1. (4)
因此,各点在平差公式中所占的权重Q与各点因光学畸变导致的位移大小SG的关系为:
Q∝SG-1. (5)
依据光学畸变计算模型可知,各数组内所包含点在各照片的像平面坐标系内的因光学畸变导致的位移大小之和SGH为:
Figure BDA0003769118740000052
结合公式(5),各数组的权重Qz为:
Figure BDA0003769118740000053
为了计算过程的简明,还应对个数组的权重进行归一化处理:
Figure BDA0003769118740000061
依据公式(2)、(3)中的数据,以及加权算法,移站算法中的旋转矩阵R的加权平差算法为:
Figure BDA0003769118740000062
R对应的平移矩阵t的加权平差算法为:
Figure BDA0003769118740000063
设旋转矩阵R为
Figure BDA0003769118740000064
依据经典摄影测量算法可知:
Figure BDA0003769118740000065
Figure BDA0003769118740000066
Figure BDA0003769118740000067
b1=cosωsinκ;
b2=cosωcosκ;
b3=-sinω;
Figure BDA0003769118740000071
Figure BDA0003769118740000072
Figure BDA0003769118740000073
因此:
Figure BDA0003769118740000074
ω2=arcsin(-b3) (12)
Figure BDA0003769118740000075
依据上述算法,即可完成基于光学畸变加权的外方位元素转移计算。本设计所述的基于光学畸变加权的外方位元素转移计算,可以通过以下方式进行验证:
试验前利用三坐标光笔测量仪对标定场内的标定板进行测量,由于标定场范围较大,三坐标光笔测量仪无法一次覆盖整个场地,所以测量时设了A、B两个站,再进行数据拼接。其中A1-A46为在测量站A测量的点,点的编号如图2所示,B61-B88为在测量站B测量的点,点的编号如图3所示。详细的试验点物方坐标数据见表A2(坐标皆为标定板的中心点)。为方便以后设站问题,设计了移站算法,可利用同名点完成两个测量站间的数据拼接。
为将B站的数据转换为A站的,在B站测量时选取了部分A站已经测量过的同名点,如表1所示。
表1A、B站同名点物方坐标
Figure BDA0003769118740000076
Figure BDA0003769118740000081
利用B62、B63、B64、B79这4个同名点通过SVD方法进行数据转换,B站数据转换后见附录表A3。将转换后的同名点数据与A站数据对比,得到转换后的数据误差,如表2所示。
表2同名点物方坐标转换误差
试验点编号 X(mm) Y(mm) Z(mm)
B62(A39) -0.075 -0.040 -0.331
B63(A38) 0.090 0.089 -0.200
B64(A40) -0.095 0.085 -0.292
B70(A43) -0.143 -0.002 0.499
B72(A44) -0.348 0.021 0.564
B74(A45) -0.300 0.038 0.592
B78(A42) 0.045 -0.050 0.619
B79(A41) 0.079 -0.134 0.823
B84(A46) 0.317 0.046 1.205
根据转换结果可知,转换后的数据误差大部分在0.5(mm)之内,精度较高。
根据试验流程,首先需要标定后站相机的外方位元素,由章节2可知,测量影像的外方位元素分别为摄影中心物方坐标X、Y、Z以及像片的空间姿态角度ψ、ω、κ。设前站在左边为A位置,右边为B位置,以A位置解算为例,由于双目系统还未完成,试验是通过移动后方相机拍摄了图4、图5左右两张像片,进行后站的外方位元素求解。
试验选点如图6所示,其中点1、2、3、4为控制点,点a、b为指针点,这些点与A、B站测量点之间的对应关系已在图中标注。
试验点的物方坐标和解算的像平面坐标分别如表3、表4所示。
表3后站试验点物方坐标
编号 X(mm) Y(mm) Z(mm)
1 -429.538 464.726 -2381.242
2 329.613 383.497 -1765.016
3 -1060.380 -201.241 -1706.850
4 306.310 -236.008 -1838.063
a -1051.724 34.540 -1678.448
b 333.265 -91.484 -1822.235
表4后站试验点像平面坐标值
编号 x-left(mm) y-left(mm) x-right(mm) y-right(mm)
1 4273.297 3718.121 2663.973 3645.044
2 5989.812 3794.625 5472.393 3700.001
3 1930.194 1967.729 1486.014 1467.790
4 5967.369 1838.329 5457.942 1824.827
a 1923.625 2660.054 1475.577 2322.504
b 6041.010 2292.288 5513.994 2265.284
根据上述表中的数据,解算外方位元素,并根据第二章所述加权算法标定f,得到的结果如表5所示。
表5后站外方位元素和f
Figure BDA0003769118740000091
Figure BDA0003769118740000101
要求解前站在A、B两个位置的旋转和平移变换,需要两个位置的同名点坐标,同名点示意图如图5-6所示。
在已知后站外方位元素的前提下,通过前方交会,可将同名点的二维平面坐标还原为三维空间中的物方坐标,在A位置,前站上四个同名点的物方坐标如表6所示。
表6A位置同名点物方坐标
编号 X(mm) Y(mm) Z(mm)
1 -282.095 -6.613 -1667.238
2 -52.385 -12.083 -1645.973
3 -46.361 257.496 -1630.095
4 -227.105 273.466 -1847.402
同理,可以求解前站在B位置上述四个点的物方坐标,如表7所示。
表7B位置同名点物方坐标
编号 X(mm) Y(mm) Z(mm)
1 888.364 -16.476 -1920.513
2 1103.374 -11.845 -2004.304
3 1104.716 257.912 -1989.135
4 844.241 265.881 -2106.483
测量阶段需要利用标定场直接解算前站在A、B两个位置的外方位元素,再通过外方位元素转移算法,利用A位置的外方位元素解算出B位置的外方位元素。
将前站置于标定场中,在A、B两个位置拍摄像片。并根据三坐标光笔测量仪测量的坐标,进行试验点的选取,如图8所示,控制点为点1、2、3、4(分别对应A站的点26、31、22、23),指针点为点a、b(分别对应A站的点25、24)。
前站试验点的物方坐标和解算的像平面坐标分别如表8、表9所示。
表8前站试验点物方坐标
编号 X(mm) Y(mm) Z(mm)
1 -26.583 177.359 -4428.766
2 201.399 172.540 -4449.419
3 84.745 691.320 -4500.980
4 101.169 582.998 -4515.902
a 92.647 352.595 -4544.429
b 76.107 479.750 -4527.530
表9前站试验点像平面坐标值
编号 x-left(mm) y-left(mm) x-right(mm) y-right(mm)
1 3817.843 2304.901 4127.341 2358.629
2 4895.568 2302.425 5250.930 2321.857
3 4267.687 4674.857 4816.816 4857.287
4 4347.981 4171.364 4908.215 4325.848
a 4323.138 3106.053 4897.102 3190.698
b 4240.926 3692.840 4803.197 3814.568
根据上述表中的数据,解算前站的外方位元素,解算结果如表10所示。
表10前站外方位元素
参数 A位置 B位置
X(mm) -144.042 940.878
Y(mm) 107.178 106.926
Z(mm) 1817.681 -2113.814
ψ(rad) 0.088 -0.379
ω(rad) 6.356 0.071
κ(rad) 29.821 -20.401
解算出前站在左右两个位置的同名点物方坐标后,可基于A位置的外方位元素解算出B位置的外方位元素,其结果如表11所示。
表11B位置外方位元素
参数 B位置
X(mm) 939.039
Y(mm) 102.309
Z(mm) -2118.769
ψ(rad) -0.373
ω(rad) 0.066
κ(rad) -1.552
将用摄影测量方法解算的B位置外方位元素,与通过外方位元转移算法素解算的B位置外方位元素进行对比,如表12所示。
表12前站B位置外方位元素对比
Figure BDA0003769118740000121
将在标定场利用摄影测量方法解算的B位置外方位元素作为真值,根据对比结果可知,按照本项目提出的预标定算法解算的外方位元素与摄影测量方法解算的外方位元素相差很小(κ的单位为弧度,-1.552与-20.401相差了6π,并非误差大),满足精度要求。
本说明书中各个实施例采用递进的方式描述,每个实施例重点说明的都是与其他实施例的不同之处,各个实施例之间相同相似部分互相参见即可。对于实施例公开的装置而言,由于其与实施例公开的方法相对应,所以描述的比较简单,相关之处参见方法部分说明即可。
对所公开的实施例的上述说明,使本领域专业技术人员能够实现或使用本发明。对这些实施例的多种修改对本领域的专业技术人员来说将是显而易见的,本文中所定义的一般原理可以在不脱离本发明的精神或范围的情况下,在其它实施例中实现。因此,本发明将不会被限制于本文所示的这些实施例,而是要符合与本文所公开的原理和新颖特点相一致的最宽的范围。

Claims (6)

1.一种基于光学畸变加权的外方位元素转移算法,其特征在于,包括以下步骤:
对移站基准的数据进行分组;
对分组后的数据进行共面检测,分为合格数据组和废弃数据组;
对合格数据组进行基于光学畸变加权的外方位元素转移计算。
2.根据权利要求1所述的一种基于光学畸变加权的外方位元素转移算法,其特征在于,对分组后的数据进行共面检测采用光学畸变模型进行数据筛选。
3.根据权利要求1所述的一种基于光学畸变加权的外方位元素转移算法,其特征在于,对移站基准的数据进行分组,具体为,设可用控制点的数量为n,单次匹配所用的点数为k,则可用的数组数量为
Figure FDA0003769118730000011
4.根据权利要求1所述的一种基于光学畸变加权的外方位元素转移算法,其特征在于,各点光学畸变导致的位移大小SG与该点数据的可信度TP的关系为:
TP∝SG-1
5.根据权利要求1所述的一种基于光学畸变加权的外方位元素转移算法,其特征在于,各点在平差公式中所占的权重Q与各点因光学畸变导致的位移大小SG的关系为:
Q∝SG-1
6.根据权利要求1所述的一种基于光学畸变加权的外方位元素转移算法,其特征在于,对合格数据组进行基于光学畸变加权的外方位元素转移计算,各数组内所包含点在各照片的像平面坐标系内的因光学畸变导致的位移大小之和SGH为:
Figure FDA0003769118730000021
其中:z=1,2,...,j,k1、k2为径向畸变参数,q1、q2为偏心畸变参数,s1、s2为薄棱镜畸变参数,δ为薄棱镜畸变参数。
CN202210896678.7A 2022-07-28 2022-07-28 一种基于光学畸变加权的外方位元素转移算法 Active CN115326019B (zh)

Priority Applications (2)

Application Number Priority Date Filing Date Title
CN202210896678.7A CN115326019B (zh) 2022-07-28 2022-07-28 一种基于光学畸变加权的外方位元素转移算法
JP2023122769A JP7427128B1 (ja) 2022-07-28 2023-07-27 光学歪み重み付けに基づく外部標定要素の転移アルゴリズム

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210896678.7A CN115326019B (zh) 2022-07-28 2022-07-28 一种基于光学畸变加权的外方位元素转移算法

Publications (2)

Publication Number Publication Date
CN115326019A true CN115326019A (zh) 2022-11-11
CN115326019B CN115326019B (zh) 2023-06-06

Family

ID=83918682

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210896678.7A Active CN115326019B (zh) 2022-07-28 2022-07-28 一种基于光学畸变加权的外方位元素转移算法

Country Status (2)

Country Link
JP (1) JP7427128B1 (zh)
CN (1) CN115326019B (zh)

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108088383A (zh) * 2017-12-05 2018-05-29 武汉理工大学 一种应用于起重机械的摄影测量算法
CN110567437A (zh) * 2019-06-28 2019-12-13 武汉理工大学 港口起重机械摄影测量方法
CN112857328A (zh) * 2021-03-30 2021-05-28 宁波市特种设备检验研究院 一种无标定摄影测量方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2013079854A (ja) 2011-10-03 2013-05-02 Topcon Corp 3次元計測システム及び3次元計測方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108088383A (zh) * 2017-12-05 2018-05-29 武汉理工大学 一种应用于起重机械的摄影测量算法
CN110567437A (zh) * 2019-06-28 2019-12-13 武汉理工大学 港口起重机械摄影测量方法
CN112857328A (zh) * 2021-03-30 2021-05-28 宁波市特种设备检验研究院 一种无标定摄影测量方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
单杰: "光束法平差简史与概要", 《武汉大学学报(信息科学版)》 *

Also Published As

Publication number Publication date
JP7427128B1 (ja) 2024-02-02
CN115326019B (zh) 2023-06-06
JP2024020165A (ja) 2024-02-14

Similar Documents

Publication Publication Date Title
CN108765328B (zh) 一种高精度多特征平面模板及其畸变优化和标定方法
CN105716542B (zh) 一种基于柔性特征点的三维数据拼接方法
CN110285827B (zh) 一种距离约束的摄影测量高精度目标定位方法
CN110099267A (zh) 梯形校正系统、方法以及投影仪
CN105424058B (zh) 基于摄影测量技术的数码相机投影中心位置精确标定方法
CN105931222A (zh) 用低精度二维平面靶标实现高精度相机标定的方法
WO2018201677A1 (zh) 基于光束平差的远心镜头三维成像系统的标定方法及装置
CN109272555B (zh) 一种rgb-d相机的外部参数获得及标定方法
CN110443879B (zh) 一种基于神经网络的透视误差补偿方法
CN111707187B (zh) 一种大型零件的测量方法及系统
WO2020199439A1 (zh) 基于单双目混合测量的三维点云计算方法
CN114331924B (zh) 大工件多相机视觉测量方法
CN110966935A (zh) 基于标志点的偏折测量系统一体化几何标定方法
CN112229323B (zh) 基于手机单目视觉的棋盘格合作目标的六自由度测量方法及其应用
CN110108203B (zh) 一种基于摄影测量技术的丝线位置测量方法及系统
CN109544642B (zh) 一种基于n型靶标的tdi-ccd相机参数标定方法
CN105405135B (zh) 基于基本配置点的两步法摄影物点、像点自动匹配方法
Yuan et al. A precise calibration method for line scan cameras
CN114926538A (zh) 单目激光散斑投影系统的外参数标定方法和装置
CN113888641A (zh) 一种基于机器视觉和深度学习的立木胸径测量方法
CN112308926B (zh) 一种无公共视场的相机外参标定方法
CN107976146B (zh) 一种线阵ccd相机的自标定方法及测量方法
CN113963067A (zh) 一种采用小靶标对大视场视觉传感器进行标定的标定方法
CN110490941B (zh) 一种基于法向量的远心镜头外参数标定方法
CN115326019A (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