CN112149295B - 一种全球通用植被总初级生产力遥感指数估算方法 - Google Patents

一种全球通用植被总初级生产力遥感指数估算方法 Download PDF

Info

Publication number
CN112149295B
CN112149295B CN202010976999.9A CN202010976999A CN112149295B CN 112149295 B CN112149295 B CN 112149295B CN 202010976999 A CN202010976999 A CN 202010976999A CN 112149295 B CN112149295 B CN 112149295B
Authority
CN
China
Prior art keywords
vegetation
ref
reflectivity
remote sensing
rchl
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
CN202010976999.9A
Other languages
English (en)
Other versions
CN112149295A (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.)
Aerospace Information Research Institute of CAS
Original Assignee
Aerospace Information Research Institute of CAS
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 Aerospace Information Research Institute of CAS filed Critical Aerospace Information Research Institute of CAS
Priority to CN202010976999.9A priority Critical patent/CN112149295B/zh
Publication of CN112149295A publication Critical patent/CN112149295A/zh
Application granted granted Critical
Publication of CN112149295B publication Critical patent/CN112149295B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • YGENERAL 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
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02ATECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/10Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Investigating Or Analysing Materials By Optical Means (AREA)

Abstract

本发明公开一种全球通用植被总初级生产力遥感指数估算方法。该方法基于全球通用植被总初级生产力遥感指数模型,采用遥感数据对陆地植被GPP实现大范围的快速估算。拟合效果显著优于目前的NIRv指数获取的结果,而且该指数模型表现在不同植被类型和纬度之间有较好的一致性,可以准确追踪不同植被类型的物候变化,同时不易受雪的影响。

Description

一种全球通用植被总初级生产力遥感指数估算方法
技术领域
本发明涉及全球植被生产力估算技术领域,具体涉及一种全球通用植被总初级生产力遥感指数估算方法。
背景技术
目前对植被总初级生产力(gross primary productivity,GPP)的估算主要采用光能利用率模型(light-use efficiency,LUE)和太阳诱导叶绿素荧光模型(sun-inducedchlorophyll fluorescence,SIF)两种方式。其中LUE模型将GPP建模为植被对光合有效辐射的吸收部分,光合有效辐射总量和,以及光能转化效率的乘积(即GPP=fPAR×PAR×ε)。而ε与植被类型和空气温湿度等参数存在复杂的关系,需要针对不同植被类型建立查找表来获取,导致计算过程十分复杂。同时由于空气温湿度等气象参数通常采用再分析模型获得,空间分辨率较低,且气象站点在全球的不均匀分布导致不同区域精度也存在差异,给GPP的计算带来复杂的不确定性。而基于SIF的GPP估算需要对植被的荧光信号进行测量,这种星载传感器目前还无法实现大范围空间的连续覆盖,并且SIF与GPP的关系在不同区域和不同植被类型中存在一定差异。
为了更方便对GPP进行估算,也有一些研究通过建立新型植被指数直接对GPP进行拟合,如叶绿素/胡萝卜素指数(chlorophyll/carotenoid index,CCI)通过中分辨率成像光谱仪(Moderate Resolution Imaging Spectroradiometer,MODIS)波段11和波段1的归一化比值(b11-b1)/(b11+b1)进行计算,并认为与常绿针叶林的光合作用物候具有很好的一致性;植被近红外指数(near-infrared reflectance of vegetation,NIRv)采用近红外反射率和归一化植被指数(normalized difference vegetation index,NDVI)的乘积获得,被认为可以更好反映GPP的动态。这些指数虽然在特定区域应用中效果可观,但与GPP的回归关系在不同地理位置和植被类型中存在差异,在研究中更多仅用于对GPP的动态进行反映,而非计算GPP的实际数值。估算GPP时需要针对具体区域和树种采用不同的计算公式,目前未见有全球通用回归模型的相关报道。
发明内容
为了解决上述问题,本发明公开一种全球通用植被总初级生产力遥感指数估算方法,包括以下步骤:
首先,提出全球通用植被总初级生产力遥感指数模型,模型的计算公式为:
GPP=A×CC×PAR×Ev+B,
其中,A=0.15,B=1.3,为回归模型的经验系数,该系数全球统一,不必随不同的植被类型和地理位置进行调整;CC代表植被覆盖度;PAR是总体下行短波辐射;Ev是植被对光合有效辐射的总体转化效率,为吸收效率和利用效率的乘积;
将CC通过归一化植被指数NDVI进行近似;吸收效率通过近红外波段反射率-NIR进行表征,利用效率通过叶绿素含量占总色素含量的比例-Rchl进行表征,
则该模型可以进一步推导为:
GPP=A×NDVI×PAR×NIR×Rchl+B;
利用MODIS反射率数据分别计算每日的NDVI,NIR和Rchl均值,并通过相应产品获取PAR;
基于所述全球通用植被总初级生产力遥感指数模型求得最终的植被总初级生产力值。
本发明的全球通用植被总初级生产力遥感指数估算方法中,优选为,所述Rchl可以根据以下两个波段的反射率Ref531和Ref667进行近似:
已知531nm为类胡萝卜素和花青素的共同吸收波段,叶绿素无吸收;667nm波段为叶绿素的吸收峰值波段,类胡萝卜素和花青素无吸收,
基于比值植被指数概念,某一色素含量与该色素的无吸收波段反射率和该色素显著吸收波段反射率的比值显著线性相关,则叶绿素含量Cchl可近似为
Cchl=a×Ref858/Ref531
类胡萝卜素和花青素含量Ccaro+anth可近似为
Ccaro+anth=b×Ref858/Ref667
其中a,b为该反射率比值和色素含量真值之间的线性拟合斜率,由于类胡萝卜素和花青素在531nm附近的吸收率与叶绿素在667nm附近的吸收率近似相等,可以近似认为a≈b,则
Rchl=Cchl/(Cchl+Ccaro+anth)=Ref531/(Ref531+Ref667),
Ref531采用MODIS的第11波段获取,Ref667采用MODIS的第13波段获取。
本发明的全球通用植被总初级生产力遥感指数估算方法中,优选为,利用MODIS反射率数据分别计算每日的NDVI,NIR和Rchl均值,并通过相应产品获取PAR,具体包括以下步骤:
获取MODIS数据陆地反射率产品MOD09GA和MYD09GA,海洋反射率产品MODOCGA和MYDOCGA,以及光合有效辐射产品MCD18;
在产品中分别获取模型所需的波段参数bx和PAR,bx代表经过辐射和大气校正的MODIS第x波段的反射率;
结合产品的质量控制波段对质量不好的数据进行剔除;
基于上述所得结果分别计算每日的NDVI,NIR和Rchl均值,
其中,
NDVI=(b1-b2)/(b1+b2);
NIR=b2;
Rchl=b11/(b11+b13),
其中,b1和b2来自陆地反射率产品MOD09GA和MYD09GA产品,b11和b13来自海洋反射率产品MODOCGA和MYDOCGA,PAR来自MCD18产品。
本发明的全球通用植被总初级生产力遥感指数估算方法所采用的模型源自光能利用率模型,具有切实物理解析意义。而且,所使用的参数全球统一,无需查找表和气温、湿度等气象参数。此外,全部输入均可通过遥感观测数据获取,计算简单快速,结合遥感数据可快速生成每日全球陆地GPP空间分布图。
附图说明
图1是基于广泛采用的叶片光学和冠层二向反射率模型模拟不同色素含量下,真实Rchl与基于光谱反射率估算Rchl的对应关系。
图2是基于本发明的全球通用植被总初级生产力遥感指数模型获得植被总初级生产力的流程图。
图3是不同植被类型模型预测GPP与Fluxnet GPP之间的Pearson相关系数分布情况。
图4是不同纬度区模型预测GPP与Fluxnet GPP之间的Pearson相关系数分布情况。
图5是典型站点Fluxnet GPP与模型估测GPP时序图。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白,下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、完整地描述,应当理解,此处所描述的具体实施例仅用以解释本发明,并不用于限定本发明。所描述的实施例仅仅是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其它实施例,都属于本发明保护的范围。
在本发明的描述中,需要说明的是,术语“上”、“下”、“垂直”“水平”等指示的方位或位置关系为基于附图所示的方位或位置关系,仅是为了便于描述本发明和简化描述,而不是指示或暗示所指的装置或元件必须具有特定的方位、以特定的方位构造和操作,因此不能理解为对本发明的限制。此外,术语“第一”、“第二”仅用于描述目的,而不能理解为指示或暗示相对重要性。
此外,在下文中描述了本发明的许多特定的细节,例如器件的结构、材料、尺寸、处理工艺和技术,以便更清楚地理解本发明。但正如本领域的技术人员能够理解的那样,可以不按照这些特定的细节来实现本发明。除非在下文中特别指出,器件中的各个部分可以由本领域的技术人员公知的材料构成,或者可以采用将来开发的具有类似功能的材料。
为了采用遥感数据对陆地植被GPP实现大范围的快速估算,本发明提出全球通用植被总初级生产力遥感指数模型,模型的计算公式为:
GPP=A×CC×PAR×Ev+B
其中,A=0.15,B=1.3,为回归模型的经验系数,该系数全球统一,不必随不同的植被类型和地理位置进行调整。
CC代表植被覆盖度,可以通过归一化植被指数NDVI进行近似;
Ev是植被对光合有效辐射的总体转化效率,为吸收效率和利用效率的乘积。其中,吸收效率通过近红外波段反射率-NIR进行表征,利用效率通过叶绿素占总色素(叶绿素、类胡萝卜素、花青素)含量的比例指数-Rchl进行表征。
该模型可以进一步推导为:
GPP=A×NDVI×PAR×NIR×Rchl+B
已知531nm为类胡萝卜素(Caro)和花青素(Anth)的共同吸收波段,叶绿素(chl)无吸收,667nm波段为叶绿素的吸收峰值波段,类胡萝卜素和花青素无吸收。基于比值植被指数概念,某一色素含量与该色素的无吸收波段反射率和该色素显著吸收波段反射率的比值显著线性相关。则叶绿素含量(Cchl)可近似为
Cchl=a×Ref858/Ref531
类胡萝卜素和花青素含量(Ccaro+anth)可近似为
Ccaro+anth=b×Ref858/Ref667
其中a,b为该反射率比值和色素含量真值之间的线性拟合斜率,由于类胡萝卜素和花青素在531nm附近的吸收率与叶绿素在667nm附近的吸收率近似相等,可以近似认为a≈b,则
Rchl=Cchl/(Cchl+Ccaro+anth)=Ref531/(Ref531+Ref667),
Ref531采用MODIS的第11波段获取,Ref667采用MODIS的第13波段获取。
图1是基于广泛采用的叶片光学和冠层二向反射率模型(PROSAIL-D模型)模拟不同色素含量下,真实Rchl与基于光谱反射率估算Rchl的对应关系。虚线代表y=x。可以看到本发明提出的基于反射率的Rchl估算方法可解释约99%的真实Rchl变化。
基于上述模型对植被总初级生产力进行计算,如图2所示,首先,获取MODIS数据陆地反射率产品MOD09GA和MYD09GA,海洋反射率产品MODOCGA和MYDOCGA,以及光合有效辐射产品MCD18。
然后,在产品中分别获取模型所需的波段参数bx,bx代表MODIS第x波段的反射率(经过辐射和大气校正)。其中,b1和b2来自MOD09GA和MYD09GA产品,b11和b13来自MODOCGA和MYDOCGA产品,PAR来自MCD18产品的总体下行短波辐射(total DSR at MODISoverpass)。
再结合产品的质量控制(QC)波段对质量不好的数据进行剔除。
然后,基于上述所得结果分别计算每日的NDVI,NIR和Rchl均值,其中
NDVI=(b1-b2)/(b1+b2);
NIR=b2;
Rchl=b11/(b11+b13)。
最后,根据全球通用植被总初级生产力遥感指数模型GPP=0.15×NDVI×PAR×NIR×Rchl+1.3求得最终的GPP值,单位为gC/m2
本发明的全球通用植被总初级生产力遥感指数估算方法所采用的模型源自光能利用率模型,具有切实物理解析意义。而且,所使用的参数全球统一,无需查找表和气温、湿度等气象参数。此外,全部输入均可通过遥感观测数据获取,计算简单快速,结合遥感数据可快速生成每日全球陆地GPP空间分布图。
采用Fluxnet2015数据集位于全球的163个不同植被类型通量站点的日GPP数据(有效点数>13.4万)对本发明的全球通用植被总初级生产力遥感指数模型进行验证r=0.70(p<10-200),显著优于仅使用目前拟合效果最好的NIRv指数获取的结果(r=0.59)
图3是不同植被类型模型预测GPP与Fluxnet GPP之间的Pearson相关系数分布情况:CRO-农田,DBF-落叶阔叶林,EBF-常绿阔叶林,ENF-常绿针叶林,GRA-草地,MF-混生林,OSH-灌丛,SAV-热带草原,WET-实地,WSA-稀树草原。图4是不同纬度区模型预测GPP与Fluxnet GPP之间的Pearson相关系数分布情况:HN-北半球高纬(60°N以北),L-低纬度地区(30°S-30°N),MN-北半球中纬(30°N-60°N),MS-南半球中纬(30°S-60°S)。如图3和图4所示,本发明的全球通用植被总初级生产力遥感指数模型表现在不同植被类型和纬度之间有较好的一致性。
图5是典型站点Fluxnet GPP与模型估测GPP(model GPP)时序图:(a)CA-NS2常绿针叶林,(b)AU-Wom常绿阔叶林,(c)DK-Sor落叶阔叶林,(d)CH-FRU草地。可以看出本模型可以准确追踪不同植被类型的物候变化。此外,可以看到该指数模型不像NDVI和CCI一样易受雪的影响。
以上所述,仅为本发明的具体实施方式,但本发明的保护范围并不局限于此,任何熟悉本技术领域的技术人员在本发明揭露的技术范围内,可轻易想到的变化或替换,都应涵盖在本发明的保护范围之内。

Claims (3)

1.一种全球通用植被总初级生产力遥感指数估算方法,其特征在于,
包括以下步骤:
提出全球通用植被总初级生产力遥感指数模型,模型的计算公式为:
GPP=A×CC×PAR×Ev+B,
其中,A=0.15,B=1.3,为回归模型的经验系数,该系数全球统一,不必随不同的植被类型和地理位置进行调整;CC代表植被覆盖度;PAR是总体下行短波辐射;Ev是植被对光合有效辐射的总体转化效率,为吸收效率和利用效率的乘积;
将CC通过归一化植被指数NDVI进行近似;吸收效率通过近红外波段反射率-NIR进行表征,利用效率通过叶绿素含量占总色素含量的比例-Rchl进行表征,
则该模型可以进一步推导为:
GPP=A×NDVI×PAR×NIR×Rchl+B;
利用MODIS反射率数据分别计算每日的NDVI,NIR和Rchl均值,并通过相应产品获取PAR;
基于所述全球通用植被总初级生产力遥感指数模型求得最终的植被总初级生产力值。
2.根据权利要求1所述的全球通用植被总初级生产力遥感指数估算方法,其特征在于,
所述Rchl可以根据以下两个波段的反射率Ref531和Ref667进行近似:
已知531nm为类胡萝卜素和花青素的共同吸收波段,叶绿素无吸收;667nm波段为叶绿素的吸收峰值波段,类胡萝卜素和花青素无吸收,
基于比值植被指数概念,某一色素含量与该色素的无吸收波段反射率和该色素显著吸收波段反射率的比值显著线性相关,则叶绿素含量Cchl可近似为
Cchl=a×Ref858/Ref531
类胡萝卜素和花青素含量Ccaro+anth可近似为
Ccaro+anth=b×Ref858/Ref667
其中a,b为该反射率比值和色素含量真值之间的线性拟合斜率,由于类胡萝卜素和花青素在531nm附近的吸收率与叶绿素在667nm附近的吸收率近似相等,可以近似认为a≈b,则
Rchl=Cchl/(Cchl+Ccaro+anth)=Ref531/(Ref531+Ref667),
Ref531采用MODIS的第11波段获取,Ref667采用MODIS的第13波段获取。
3.根据权利要求2所述的全球通用植被总初级生产力遥感指数估算方法,其特征在于,
利用MODIS反射率数据分别计算每日的NDVI,NIR和Rchl均值,并通过相应产品获取PAR,具体包括以下步骤:
获取MODIS数据陆地反射率产品MOD09GA和MYD09GA,海洋反射率产品MODOCGA和MYDOCGA,以及光合有效辐射产品MCD18;
在产品中分别获取模型所需的波段参数bx和PAR,bx代表经过辐射和大气校正的MODIS第x波段的反射率;
结合产品的质量控制波段对质量不好的数据进行剔除;
基于上述所得结果分别计算每日的NDVI,NIR和Rchl均值,
其中,
NDVI=(b1-b2)/(b1+b2);
NIR=b2;
Rchl=b11/(b11+b13),
其中,b1和b2来自陆地反射率产品MOD09GA和MYD09GA产品,b11和b13来自海洋反射率产品MODOCGA和MYDOCGA,PAR来自MCD18产品。
CN202010976999.9A 2020-09-17 2020-09-17 一种全球通用植被总初级生产力遥感指数估算方法 Active CN112149295B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010976999.9A CN112149295B (zh) 2020-09-17 2020-09-17 一种全球通用植被总初级生产力遥感指数估算方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010976999.9A CN112149295B (zh) 2020-09-17 2020-09-17 一种全球通用植被总初级生产力遥感指数估算方法

Publications (2)

Publication Number Publication Date
CN112149295A CN112149295A (zh) 2020-12-29
CN112149295B true CN112149295B (zh) 2023-07-18

Family

ID=73892305

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010976999.9A Active CN112149295B (zh) 2020-09-17 2020-09-17 一种全球通用植被总初级生产力遥感指数估算方法

Country Status (1)

Country Link
CN (1) CN112149295B (zh)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115952702B (zh) * 2022-08-30 2023-07-07 中国气象科学研究院 一种基于forcchn模型和遥感数据的森林nep的计算方法
CN115358095B (zh) * 2022-10-19 2023-03-17 中国科学院、水利部成都山地灾害与环境研究所 山地高空间分辨率植被总初级生产力估算方法
CN117314184B (zh) * 2023-09-19 2024-04-30 沈阳农业大学 年总初级生产力时空变化结果的确定方法、装置及设备

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103886213A (zh) * 2014-04-03 2014-06-25 北京大学 一种农作物总初级生产力的遥感估算方法及系统
CN110276304A (zh) * 2019-06-25 2019-09-24 北京师范大学 基于降尺度的高分辨率植被生产力遥感估算方法
CN110443504A (zh) * 2019-07-30 2019-11-12 中国环境科学研究院 一种植被净初级生产力的遥感估算方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
MX2011013334A (es) * 2009-06-11 2012-04-11 Pa Llc Indices de vegetacion para medir densidad de microcultivo y crecimiento en multicapas.

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103886213A (zh) * 2014-04-03 2014-06-25 北京大学 一种农作物总初级生产力的遥感估算方法及系统
CN110276304A (zh) * 2019-06-25 2019-09-24 北京师范大学 基于降尺度的高分辨率植被生产力遥感估算方法
CN110443504A (zh) * 2019-07-30 2019-11-12 中国环境科学研究院 一种植被净初级生产力的遥感估算方法

Also Published As

Publication number Publication date
CN112149295A (zh) 2020-12-29

Similar Documents

Publication Publication Date Title
CN112149295B (zh) 一种全球通用植被总初级生产力遥感指数估算方法
CN110222475B (zh) 一种基于无人机多光谱遥感反演冬小麦植株含水率的方法
Migliavacca et al. Using digital repeat photography and eddy covariance data to model grassland phenology and photosynthetic CO2 uptake
CN110865040A (zh) 天空地一体化高光谱水质监测分析方法
Gamon et al. A mobile tram system for systematic sampling of ecosystem optical properties
Martin et al. A generalizable method for remote sensing of canopy nitrogen across a wide range of forest ecosystems
US10663445B2 (en) Method and system for identifying plant species based on hyperspectral data
Wu et al. The potential of the satellite derived green chlorophyll index for estimating midday light use efficiency in maize, coniferous forest and grassland
Rautiainen et al. Coupling forest canopy and understory reflectance in the Arctic latitudes of Finland
CN113554232A (zh) 一种农作物产量预测方法及系统
US20210247297A1 (en) Systems and Methods for Converting Satellite Images to Surface Reflectance Using Scene Statistics
Matese et al. Mapping of vine vigor by UAV and anthocyanin content by a non-destructive fluorescence technique
Harris et al. Retrieval of the photochemical reflectance index for assessing xanthophyll cycle activity: a comparison of near-surface optical sensors
Huete et al. Monitoring photosynthesis from space
CN114821349A (zh) 顾及谐波模型系数和物候参数的森林生物量估算方法
CN113673490A (zh) 一种物候期自适应的作物生理参数遥感估测方法及系统
Cho Hyperspectral remote sensing of biochemical and biophysical parameters: the derivate red-edge" double-peak feature", a nuisance or an opportunity?
Yule et al. Mapping nutrient concentration in pasture using hyperspectral imaging
CN115759524A (zh) 一种基于遥感影像植被指数的土壤生产力等级识别方法
Chen et al. Leaf chlorophyll contents dominates the seasonal dynamics of SIF/GPP ratio: Evidence from continuous measurements in a maize field
Páscoa In situ visible and near-infrared spectroscopy applied to vineyards as a tool for precision viticulture
Barela et al. Applications of proximal remote sensing in agriculture: A review
Yuke Greenness index from phenocams performs well in linking climatic factors and monitoring grass phenology in a temperate prairie ecosystem
Aragon et al. Revisiting the spatial scale effects on remotely sensed evaporation
US20230237622A1 (en) Model generation and application for removing atmospheric effects in imagery

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