CN110913306B - 一种实现阵列麦克风波束形成的方法 - Google Patents

一种实现阵列麦克风波束形成的方法 Download PDF

Info

Publication number
CN110913306B
CN110913306B CN201911213687.6A CN201911213687A CN110913306B CN 110913306 B CN110913306 B CN 110913306B CN 201911213687 A CN201911213687 A CN 201911213687A CN 110913306 B CN110913306 B CN 110913306B
Authority
CN
China
Prior art keywords
array
sub
frequency
frequency part
microphone
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
CN201911213687.6A
Other languages
English (en)
Other versions
CN110913306A (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.)
Beijing Philisense Electronics Co ltd
Original Assignee
Beijing Philisense Electronics 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 Beijing Philisense Electronics Co ltd filed Critical Beijing Philisense Electronics Co ltd
Priority to CN201911213687.6A priority Critical patent/CN110913306B/zh
Publication of CN110913306A publication Critical patent/CN110913306A/zh
Application granted granted Critical
Publication of CN110913306B publication Critical patent/CN110913306B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • HELECTRICITY
    • H04ELECTRIC COMMUNICATION TECHNIQUE
    • H04RLOUDSPEAKERS, MICROPHONES, GRAMOPHONE PICK-UPS OR LIKE ACOUSTIC ELECTROMECHANICAL TRANSDUCERS; DEAF-AID SETS; PUBLIC ADDRESS SYSTEMS
    • H04R1/00Details of transducers, loudspeakers or microphones
    • H04R1/20Arrangements for obtaining desired frequency or directional characteristics
    • H04R1/22Arrangements for obtaining desired frequency or directional characteristics for obtaining desired frequency characteristic only 
    • H04R1/26Spatial arrangements of separate transducers responsive to two or more frequency ranges
    • H04R1/265Spatial arrangements of separate transducers responsive to two or more frequency ranges of microphones

Landscapes

  • Health & Medical Sciences (AREA)
  • Otolaryngology (AREA)
  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Acoustics & Sound (AREA)
  • Signal Processing (AREA)
  • Circuit For Audible Band Transducer (AREA)

Abstract

本发明提供了一种实现阵列麦克风波束形成的方法,其特征在于,包括以下步骤,1)将需处理的声音范围分为高频部分a-b和低频部分b-c;2)选取比例常数R,使c×Rm=b,b×Rn=a;3)以[b×R0],[b×R1],[b×R2],.....,[b×Rn]分别为中心频率,计算每个中心频率对应的二元阵孔径;4)将阵列麦克风每个阵元分解成子带信号,对于其中高频部分的n+1个子带按照二阵元直接相加,低频部分的m个子带按照n+2个阵元加权相加;5)对步骤4)中所有子带按照汉明窗平滑合成,得到了在a-c频率范围内的空间响应。本发明的实现阵列麦克风波束形成的方法采用阵列麦克风孔径叠加和自适应波束形成的方法还原宽带信号,能够使用较少的麦克风单元和较小的阵列尺寸,具有很好的语音还原质量。

Description

一种实现阵列麦克风波束形成的方法
技术领域
本发明属于语音信号处理技术领域,具体涉及一种实现阵列麦克风波束形成的方法。
背景技术
在现场会议场景中,分割子带、对子带进行降噪、再对子带合并的谱减法的单通道降噪方法已应用广泛,但这种方法只是抑制了被噪声影响的频谱分量,造成了处理的语音信号失真而降低了语音质量和可理解性。人们在借鉴雷达、天线等阵列的应用后,设计出带有指向性的阵列麦克风。阵列麦克风有诸多优点,例如不需要鹅颈话筒因而更加美观,具有指向性因而有助于抑制会场噪声的干扰,抑制音响设备与发言设备之间可能产生的啸叫等等。阵列麦克风的语音增强算法一般采用波束形成噪音抑制方法,其思想是期望在阵列中心到目标语音间形成一个波束,而波束外的信号被抑制掉,然后对波束内渗入的噪音采用单道噪音抑制算法进行处理,得到最终的目标语音输出。
不过,现有的阵列麦克风噪声抑制技术,普遍存在算法复杂度高,计算量大或者阵列结构特殊、适应性不强的缺陷。主要原因在于,借鉴于雷达、天线等阵列的处理方法,必须适应非常宽泛的电磁波谱。而人耳能明确感知的声波主要在130Hz-8000Hz,语音信号主要在300Hz-3400Hz。虽然声音信号的频率跨度非常大,难以在阵列麦克风有限的阵列尺寸和阵元个数下,使其空间响应在整个频带达到理想的一致性。但也不能简单地认为采用固定或自适应波束成形算法就能实现阵列麦克风语音交互的良好效果,毕竟人耳接收的不是冰冷的电平,而是有很大的主观性。
因此,需要一种综合考虑阵列麦克风孔径和自适应波束形成的方法,结合人类语音信号频率范围,而可以使用较少的麦克风单元和较小的阵列尺寸,在语音的频率范围内实现一致的空间响应。查阅国内外文献,还不存在相关技术方法。
发明内容
本发明的发明目的是提供一种实现阵列麦克风波束形成的方法,通过对人耳接收的声音频率分为高频和低频两部分进行子带划分,采用阵列麦克风孔径叠加和自适应波束形成的方法还原宽带信号,能够使用较少的麦克风单元和较小的阵列尺寸,在语音的频率范围内实现一致的空间响应,并具有很好的语音还原质量。
本发明的具体技术方案是一种实现阵列麦克风波束形成的方法,其特征在于,包括以下步骤,
1)将需处理的声音范围分为高频部分a-b和低频部分b-c,其中a>b>c,a的取值范围为[4kHz,16kHz],b的取值范围为[850Hz,3.4kHz],c的取值范围为[20Hz,1kHz];
2)选取比例常数R,使c×Rm=b,b×Rn=a,其中,m、n为正整数;
3)以[b×R0],[b×R1],[b×R2],.....,[b×Rn]分别为中心频率,计算每个中心频率对应的二元阵孔径
Figure BDA0002298880480000021
得到包括n+2个阵元的阵列麦克风,其中最左侧阵元为参考阵元,其右侧n+1个阵元分别与参考阵元组合成n+1个子阵,n+1个子阵的孔径分别为d0,d1,d2,.....,dn
4)以[c×R0],[c×R1],[c×R2],.....,[c×Rm],[b×R1],[b×R2],.....,[b×Rn]为中心频率将阵列麦克风每个阵元分解成子带信号,其中[c×Rm]=[b×R0],对于其中高频部分[b×R0],[b×R1],[b×R2],.....,[b×Rn]的n+1个子带按照二阵元直接相加,低频部分[c×R0],[c×R1],[c×R2],.....,[c×Rm-1]的m个子带按照n+2个阵元加权相加;
5)对步骤4)中所有子带按照汉明窗平滑合成,得到了在a-c频率范围内的空间响应。
更进一步地,所述的步骤4)中低频部分的m个子带按照n+2个阵元加权相加的方法为LCMV方法。
更进一步地,所述的a、b、c分别为8000Hz,b为1543Hz,c为131Hz,所述的比例常数R为1.509,可以得到中心频率分别为131Hz、197Hz、298Hz、450Hz、678Hz、1022Hz、1543Hz、2328Hz、3513Hz、5302Hz和8000Hz,6个阵元中一个为参考阵元,其余阵元距参考阵元的孔径距离分别为2.13cm、3.21cm、4.84cm、7.30cm和11.02cm。
本发明的有益效果是1)本发明的实现阵列麦克风波束形成的方法通过对人耳接收的声音频率分为高频和低频两部分进行子带划分,由高频部分设定阵元数量和阵列麦克风孔径,对于高频部分子带按照二阵元直接相加的方法处理,对于低频子带采用LCMV(Linear Constraint Minimal Variance)方法波束形成,得到了在整个频率范围内的空间响应;2)以不同尺寸的子阵处理高频部分,以不同权矢量处理低频部分,通过两种方法的综合,扩大了波束形成的频率范围,在130Hz至8000Hz,超过60倍的频率范围内实现了一致性较好的空间响应;3)通过汉明窗平滑合成的方法,使子带之间平滑过渡,进一步提高了空间响应的一致性;4)在划分子带时,使各相邻子带之间部分重叠,然后利用汉明窗使重叠部分加权叠加,因此抵消了子带之间响应曲线的变化,使子带合成结果更加平滑,增强了空间响应的一致性;5)本发明以LCMV方法处理低频子带,降低了对阵列尺寸的需求。另外,嵌套的阵列结构使多个二阵元子阵共用一个参考阵元,进一步降低了阵列尺寸,也大大地减少了算法计算量。
本发明的实现阵列麦克风波束形成的方法能够使用较少的麦克风单元和较小的阵列尺寸,以更小的计算量,在人耳接收的声音频率范围内实现一致的空间响应,并具有很好的语音还原质量。
附图说明
图1本发明的实现阵列麦克风波束形成的方法实施例的高频部分子带划分图;
图2本发明的实现阵列麦克风波束形成的方法实施例的低频部分子带划分图;
图3二元阵及其空间响应示意图;
图4本发明的实现阵列麦克风波束形成的方法实施例的阵列麦克风结构示意图;
图5本发明的实现阵列麦克风波束形成的方法实施例的高频部分子带的汉明窗平滑合成示意图;
图6本发明的实现阵列麦克风波束形成的方法实施例的高频部分一个子带汉明窗合成前后比较示意图;
图7本发明的实现阵列麦克风波束形成的方法实施例的高频部分所有子带汉明窗合成示意图;
图8本发明的实现阵列麦克风波束形成的方法实施例的高频部分和低频部分所有子带汉明窗合成示意图;
具体实施方式
下面结构说明书附图对本发明的具体技术方案作进一步地描述。
本发明的一种实现阵列麦克风波束形成的方法,其特征在于,包括以下步骤,
1)将需处理的声音范围分为高频部分a-b和低频部分b-c,其中a>b>c,a的取值范围为[4kHz,16kHz],b的取值范围为[850Hz,3.4kHz],c的取值范围为[20Hz,1kHz]。对于a的取值范围,我们知道人耳的敏感频率范围上限约为4kHz,在高音质需求中可以扩展上限,但通常不有关会超过16kHz。对于b的取值范围,我们一般认为与我们阵列最大孔径dmax有关,考虑物理设备的限制和人类语音的中间范围,对应的阵列麦克风的尺寸一般在5cm至20cm范围内,根据340/(2*dmax)转化成对应的取值范围。对于c的取值范围,因为人耳的敏感频率范围下限约为1kHz,在高音质需求中可以扩展其下限,但极限情况也不会低于20Hz;
2)选取比例常数R,使c×Rm=b,b×Rn=a,其中,m、n为正整数;
3)以[b×R0],[b×R1],[b×R2],.....,[b×Rn]分别为中心频率,计算每个中心频率对应的二元阵孔径
Figure BDA0002298880480000051
得到包括n+2个阵元的阵列麦克风,其中最左侧阵元为参考阵元,其右侧n+1个阵元分别与参考阵元组合成n+1个子阵,n+1个子阵的孔径分别为d0,d1,d2,.....,dn
4)以[c×R0],[c×R1],[c×R2],.....,[c×Rm],[b×R1],[b×R2],.....,[b×Rn]为中心频率将阵列麦克风每个阵元分解成子带信号,其中[c×Rm]=[b×R0],对于其中高频部分[b×R0],[b×R1],[b×R2],.....,[b×Rn]的n+1个子带按照二阵元直接相加,低频部分[c×R0],[c×R1],[c×R2],.....,[c×Rm-1]的m个子带按照n+2个阵元加权相加;
5)对步骤4)中所有子带按照汉明窗平滑合成,得到了在a-c频率范围内的空间响应。
采用本发明的实现阵列麦克风波束形成的方法,制成了一个具体的阵列麦克风实施例。
如图1-2所示,分别确定a、b、c分别为8000Hz,b为1543Hz,c为131Hz后,将高频部分分为中心频率为8000Hz、5302Hz、3513Hz、2328Hz、1543Hz的8000-5302Hz、8000-3513Hz、5302-2328Hz、3513-1543Hz、2328-1022Hz共5个子带,将低频部分分为中心频率为1022Hz、678Hz、450Hz、298Hz、197Hz、131Hz的1543-678Hz、1022-450Hz、678-298Hz、450-197Hz、298-131Hz、197-131Hz共6个子带。在划分子带时,使高频部分的子带2328-1022Hz与低频部分的1543-678Hz有部分重合,这是为了在最后波束形成时有平滑的过渡。
在高频部分,我们使用二阵元阵列,这是最简单的阵列结构。如附图3所示,是一个孔径为d的二阵元阵列,将此两个阵元的信号直接相加即可对频率为
Figure BDA0002298880480000061
的声音信号形成如图3右侧的空间响应,其中340m/s是声速。该响应曲线符合我们对阵列麦克风指向性的要求:在保护方向上形成一个增益为1且宽度合适的主瓣,在抑制方向上形成增益为0的零点。
对于高频部分的子带,以合适孔径的阵列分别处理。首先选定比例常数R,以等比的方式划分子带。一般比例常数R不宜过大。但该比例常数R过小则子带很窄,需要划分更多子带,系统的复杂度和计算量会很大。本实施例中比例常数R=1.509。
高频部分第一个子带的中心频率为8000Hz,然后按照等比方式依次划分其它子带。随着频率下降,二阵元阵列的孔径也成比例上升,因此仅以二阵元阵列的方式处理前5个子带,其中心频率依次为8000Hz、5302Hz、3513Hz、2328Hz和1543Hz,第1子带没有高于8000Hz的部分,以8000Hz为截止频率。
高频部分5个子带对应的二阵元阵列孔径分别为各中心频率的半波长,因此5个二阵元阵列的孔径分别为2.13cm、3.21cm、4.84cm、7.30cm和11.02cm,得到如附图4所示的嵌套结构,它包含6个阵元,以最左侧阵元为参考阵元,其右侧5个阵元分别与之组合成5个子阵。
为了控制阵列尺寸,对低频部分的子带不再使用二阵元处理,而是利用附图4所示的6个阵元的信号经加权相加处理。
本发明实施例中采用LCMV(Linear Constraint Minimal Variance)方法,即线性约束最小方差方法来计算权矢量。LCMV的基本原理是将波束形成问题转换为一个约束最小化问题:
Wopt=argminw(|WX|2) (1)
subj
ect to:
W·d(fi,θi)=Gi (2)
其中(2)式表示一个约束条件,它指定频率为fi、入射角度为θi的信号的增益为Gi,可以指定多个这样的约束条件;(1)式表示在约束条件下求解使输出能量最小的权矢量,最优解称为最优权矢量Wopt,对应的空间响应称为最优空间响应。
对于本发明的方法实施例的阵列麦克风,可以指定3个约束条件:从90°入射信号的增益为1,从0°和180°入射信号的增益为0,在附图4的阵列结构下,对指定的频率可以求出最优权矢量。LCMV方法需要的阵列孔径明显更小,这正是本发明的方法在低频部分应用LCMV方法的原因。
如附图5所示,相邻子带有部分重叠,此重叠部分在子带合成时起到过渡作用,使合成结果更加平滑。一个子带内不同频率处信号的空间响应不同,图5中A行5条曲线展示了中心频率为3513Hz的第3子带在五个频率点上的空间响应,在中心频率处的响应曲线最理想,越往边界处变形越严重。图5中B行和C行分别展示了中心频率为5302Hz的第2子带及中心频率为2328Hz的第4子带与第3子带的重叠部分在相同频率点上的响应曲线。可以见到,若在5302Hz处取B曲线;在4408Hz处取A和B曲线的平均;在3513Hz处取A曲线;在2920Hz处取A和C曲线的平均;在2328Hz处取C曲线,则可以在这五个频率上得到较一致的空间响应。为了达到这一效果,一般应用汉明窗进行平滑合成。附图5中D展示了我们使用的汉明窗,将各子带做波束形成的结果按汉明窗相乘叠加,可以使空间响应在整个频率范围内达到很好的一致性。附图6展示了第3子带内各频率的响应曲线的变化,左侧是做汉明窗叠加之前的,其变化范围很大,右侧是做汉明窗叠加之后的,所有曲线都收敛到了理想响应附近,一致性显著提高。应用汉明窗平滑合成方法,高频部分5个子带的合成结果如附图7所示。
本发明的方法应用的实例中,总共将131Hz至8000Hz的频率范围划分为了11个子带,其中高频部分的5个子带按照二阵元直接相加的方法处理,低频部分的的6个子带按照六阵元加权相加的方法处理。对所有子带按照汉明窗平滑合成,得到了在整个频率范围内的空间响应,如附图8所示。
虽然本发明已经以较佳实施例公开如上,但实施例并不是用来限定本发明的。在不脱离本发明之精神和范围内,所做的任何等效变化或润饰,同样属于本发明之保护范围。因此本发明的保护范围应当以本申请的权利要求所界定的内容为标准。

Claims (2)

1.一种实现阵列麦克风波束形成的方法,其特征在于,包括以下步骤,
1)将需处理的声音范围分为高频部分a-b和低频部分b-c,其中a>b>c,a的取值范围为[4kHz,16kHz],b的取值范围为[850Hz,3.4kHz],c的取值范围为[20Hz,1kHz];
2)选取比例常数R=1.509,使c×Rm=b,b×Rn=a,其中,m、n为正整数;
3)以[b×R0],[b×R1],[b×R2],.....,[b×Rn]分别为中心频率,计算每个中心频率对应的二元阵孔径
Figure FDA0003047284960000011
得到包括n+2个阵元的阵列麦克风,其中最左侧阵元为参考阵元,其右侧n+1个阵元分别与参考阵元组合成n+1个子阵,n+1个子阵的孔径分别为d0,d1,d2,.....,dn
4)以[c×R0],[c×R1],[c×R2],.....,[c×Rm],[b×R1],[b×R2],.....,[b×Rn]为中心频率将阵列麦克风每个阵元分解成子带信号,其中[c×Rm]=[b×R0],对于其中高频部分[b×R0],[b×R1],[b×R2],.....,[b×Rn]的n+1个子带按照二阵元直接相加,低频部分[c×R0],[c×R1],[c×R2],.....,[c×Rm-1]的m个子带按照n+2个阵元加权相加,低频部分的m个子带按照n+2个阵元加权相加的方法为线性约束最小方差(LCMV)方法;
5)对步骤4)中所有子带按照汉明窗平滑合成,得到了在a-c频率范围内的空间响应。
2.根据权利要求1所述的一种实现阵列麦克风波束形成的方法,其特征在于,所述的a为8000Hz,b为1543Hz,c为131Hz,所述的比例常数R为1.509,可以得到中心频率分别为131Hz、197Hz、298Hz、450Hz、678Hz、1022Hz、1543Hz、2328Hz、3513Hz、5302Hz和8000Hz,6个阵元中一个为参考阵元,其余阵元距参考阵元的孔径距离分别为2.13cm、3.21cm、4.84cm、7.30cm和11.02cm。
CN201911213687.6A 2019-12-02 2019-12-02 一种实现阵列麦克风波束形成的方法 Active CN110913306B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201911213687.6A CN110913306B (zh) 2019-12-02 2019-12-02 一种实现阵列麦克风波束形成的方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201911213687.6A CN110913306B (zh) 2019-12-02 2019-12-02 一种实现阵列麦克风波束形成的方法

Publications (2)

Publication Number Publication Date
CN110913306A CN110913306A (zh) 2020-03-24
CN110913306B true CN110913306B (zh) 2021-07-02

Family

ID=69821250

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201911213687.6A Active CN110913306B (zh) 2019-12-02 2019-12-02 一种实现阵列麦克风波束形成的方法

Country Status (1)

Country Link
CN (1) CN110913306B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112565981B (zh) * 2021-02-20 2021-06-15 深圳市智听科技有限公司 啸叫抑制方法、装置、助听器及存储介质

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105223544A (zh) * 2015-08-26 2016-01-06 南京信息工程大学 近场线性约束最小方差自适应加权频率不变波束形成方法
US10210882B1 (en) * 2018-06-25 2019-02-19 Biamp Systems, LLC Microphone array with automated adaptive beam tracking

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US10448152B2 (en) * 2015-09-21 2019-10-15 Northeastern University Systems and methods for monitoring and classifying marine animals based on acoustic signals
CN106782590B (zh) * 2016-12-14 2020-10-09 南京信息工程大学 基于混响环境下麦克风阵列波束形成方法
CN110428851B (zh) * 2019-08-21 2022-02-18 浙江大华技术股份有限公司 基于麦克风阵列的波束形成方法和装置、存储介质

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105223544A (zh) * 2015-08-26 2016-01-06 南京信息工程大学 近场线性约束最小方差自适应加权频率不变波束形成方法
US10210882B1 (en) * 2018-06-25 2019-02-19 Biamp Systems, LLC Microphone array with automated adaptive beam tracking

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
FIRST ORDER DIFFERENTIAL MICROPHONE ARRAYS FOR AUTOMOTIVE APPLICATIONS;Markus Buck等;《7th international workshop on acoustic echo and noise control darmstadt university of technology》;20010901;全文 *

Also Published As

Publication number Publication date
CN110913306A (zh) 2020-03-24

Similar Documents

Publication Publication Date Title
EP2884763B1 (en) A headset and a method for audio signal processing
KR101519768B1 (ko) 멀티-마이크로폰 어레이를 가지고 노이즈를 제거하기 위한 방법, 장치, 및 시스템
US9749745B2 (en) Low noise differential microphone arrays
Doclo et al. Acoustic beamforming for hearing aid applications
Chen et al. On the design and implementation of linear differential microphone arrays
KR100978823B1 (ko) 확장된 위상차 값을 이용하여 노이즈 구별을 개선하기 위한방법 및 장치
KR100978827B1 (ko) 감쇄 인자를 이용하여 노이즈 구별을 개선하기 위한 방법및 장치
KR101060301B1 (ko) 센서 어레이에서 장치 또는 신호의 부정합을 조정하기 위한 방법 및 장치
CN105590631B (zh) 信号处理的方法及装置
US7756278B2 (en) Ultra-directional microphones
EP1994788B1 (en) Noise-reducing directional microphone array
US7619563B2 (en) Beam former using phase difference enhancement
US8861745B2 (en) Wind noise mitigation
WO2007025033A2 (en) Method and system for enhancing regional sensitivity noise discrimination
WO2010043998A1 (en) Microphone system and method of operating the same
WO2007025232A2 (en) System and method for improving time domain processed sensor signal output
US20180308503A1 (en) Real-time single-channel speech enhancement in noisy and time-varying environments
CN110913306B (zh) 一种实现阵列麦克风波束形成的方法
Ward et al. Broadband microphone arrays for speech acquisition
As’ad et al. Robust minimum variance distortionless response beamformer based on target activity detection in binaural hearing aid applications
Ahmed et al. Comparison of LMS and NLMS algorithm with the using of 4 Linear Microphone Array for Speech Enhancement.
Arora et al. Comparison of speech intelligibility parameter in cochlear implants by spatial filtering and coherence function methods
Chau et al. A subband beamformer on an ultra low-power miniature DSP platform
CN113329288B (zh) 一种基于陷波技术的蓝牙耳机降噪方法
Nguyen et al. A Study Of Dual Microphone Array For Speech Enhancement In Noisy Environment

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