CN110045323B - 一种基于矩阵填充的互质阵稳健自适应波束形成算法 - Google Patents

一种基于矩阵填充的互质阵稳健自适应波束形成算法 Download PDF

Info

Publication number
CN110045323B
CN110045323B CN201910194101.XA CN201910194101A CN110045323B CN 110045323 B CN110045323 B CN 110045323B CN 201910194101 A CN201910194101 A CN 201910194101A CN 110045323 B CN110045323 B CN 110045323B
Authority
CN
China
Prior art keywords
vector
matrix
array
prime
covariance matrix
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.)
Expired - Fee Related
Application number
CN201910194101.XA
Other languages
English (en)
Other versions
CN110045323A (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.)
University of Electronic Science and Technology of China
Original Assignee
University of Electronic Science and Technology of China
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 University of Electronic Science and Technology of China filed Critical University of Electronic Science and Technology of China
Priority to CN201910194101.XA priority Critical patent/CN110045323B/zh
Publication of CN110045323A publication Critical patent/CN110045323A/zh
Application granted granted Critical
Publication of CN110045323B publication Critical patent/CN110045323B/zh
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S3/00—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received
    • G01S3/02—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received using radio waves
    • G01S3/14—Systems for determining direction or deviation from predetermined direction
    • G01S3/143—Systems for determining direction or deviation from predetermined direction by vectorial combination of signals derived from differently oriented antennae
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S3/00—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received
    • G01S3/78—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received using electromagnetic waves other than radio waves
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S3/00—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received
    • G01S3/80—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received using ultrasonic, sonic or infrasonic waves
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S3/00—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received
    • G01S3/80—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received using ultrasonic, sonic or infrasonic waves
    • G01S3/82—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received using ultrasonic, sonic or infrasonic waves with means for adjusting phase or compensating for time-lag errors
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S3/00—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received
    • G01S3/80—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received using ultrasonic, sonic or infrasonic waves
    • G01S3/86—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received using ultrasonic, sonic or infrasonic waves with means for eliminating undesired waves, e.g. disturbing noises

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Electromagnetism (AREA)
  • Radar Systems Or Details Thereof (AREA)
  • Radio Transmission System (AREA)

Abstract

本发明提出一种基于矩阵填充的互质阵稳健自适应波束形成算法，包括：计算接收数据的样本协方差矩阵；将样本协方差矩阵向量化得到向量，然后对向量进行去冗余和向量重排，得到完整互质阵差分优化阵的接收数据向量；将接收数据向量中波程差不连续的元素之间全部填0，得到向量，然后取向量的正半部分的信息，得到向量；将向量扩展成托普利兹矩阵；将托普利兹矩阵恢复得到填充协方差矩阵；在干扰信号角度区域内进行谱峰搜索得到各干扰信号到达角的估计；利用估计出的干扰信号的到达角、以及互质阵的物理阵列信息以及重构出互质阵物理阵列的干扰加噪声协方差矩阵；利用干扰加噪声协方差矩阵和期望信号导向矢量的估计计算自适应波束形成器的加权矢量。

Description

一种基于矩阵填充的互质阵稳健自适应波束形成算法
技术领域
本发明属于阵列信号处理算法领域，具体涉及一种基于矩阵填充的互质阵稳健自适应波束形成算法。
背景技术
自适应波束形成是阵列信号处理中的的核心技术之一，在诸如移动通信、雷达、医学成像、声呐探测等领域被广泛应用。自适应波束形成技术能够根据接收的训练序列或者通信信号，自适应地改变阵列天线的加权向量，从而达到抑制干扰、无失真接收期望信号的目的。自适应波束形成技术已经被理论和实践证明其对期望信号导向矢量误差、协方差矩阵误差和采样快拍次数不足等问题十分敏感。造成这些误差的因素主要包括阵元位置误差，阵元间互藕，阵列通道幅相误差，期望信号波达方向(direction ofarrival，DOA)估计误差，以及局部散射等。只要已知的先验信息与真实值有偏差，就会严重影响自适应波束形成器的性能。所以如何降低自适应波束形成器对误差的敏度，提高自适应波束形成算法的稳健性，是目前波束形成技术领域的热点。
近年来，一种新型稀疏阵——互质阵开始受到研究者们的关注。互质阵是通过在均匀线阵中抽掉一些特定的阵元来构成的。它能够以更少的阵元、更大的阵元间距，获得更大的阵列孔径与更多的自由度。和均匀线阵相比，互质阵能够显著提高阵列分辨率，并且降低系统成本与复杂度。与其他类型的稀疏阵相比，互质阵阵元间距较小的相邻阵元数更少，因此互质阵拥有更小的阵元互藕。在对互质阵接收数据进行处理时，一般由其高阶统计量得到相应的差分优化阵，从而获得虚拟阵列的更大自由度的优势。目前，互质阵的研究主要集中在DOA估计上，而在自适应波束形成方面的研究较少。不过，由于互质阵优良的特性，最近对基于互质阵的稳健自适应波束形成技术的研究也开始渐渐开展。目前互质阵自适应波束形成的文章的算法都只能利用互质阵差分优化阵的一部分连续阵元构成的ULA部分来处理虚拟阵列的接收信号，而互质阵的差分优化阵中是有许多孔洞的，因此，这些算法就无法利用到互质阵虚拟阵列全部的自由度。近年来，相关研究表明可以通过矩阵填充技术来对互质阵的差分优化阵进行插值，基于此理论，一些基于差分优化阵插值技术的互质阵DOA估计算法被提出。这些算法的核心思想是，通过矩阵填充技术，对互质阵的差分优化阵的非连续部分的孔洞进行插值，从而使整个差分优化阵形成一个连续的ULA，使其所有虚拟阵元的数据得到有效利用，从而提高了算法实际能利用上的互质阵的差分优化阵的自由度。与DOA估计类似，对稳健自适应波束形成算法来说，用于处理接收信号的阵列的孔径越大、自由度越高，干扰信号就能被更准确地定位和抑制，从而获得更好的波束形成器性能。基于此事实，可以考虑利用矩阵填充技术对互质阵的差分优化阵进行插值来提高稳健自适应波束形成算法的性能。
发明内容
鉴于以上所述现有技术的缺点，本发明的目的在于提供一种基于矩阵填充的互质阵稳健自适应波束形成算法，本发明利用矩阵填充对互质阵的差分优化阵进行插值，将其变成一个连续的ULA，使得更多的自由度能够得到利用，从而提高了自适应波束形成器的输出性能。
为实现上述目的及其他相关目的，本发明提供一种基于矩阵填充的互质阵稳健自适应波束形成算法，该互质阵包括第一ULA和第二ULA，所述第一ULA与所述第二ULA互质；所述算法包括以下步骤：
计算接收数据的样本协方差矩阵
Figure GDA0003083536180000021
对样本协方差矩阵
Figure GDA0003083536180000022
进行特征分解计算其最小的特征值γmin
将所述样本协方差矩阵
Figure GDA0003083536180000023
向量化得到向量z，然后对所述向量z进行去冗余和向量重排，得到完整互质阵差分优化阵的接收数据向量z1；
将所述接收数据向量z1中波程差不连续的元素之间全部填0，得到向量
Figure GDA0003083536180000024
然后取所述向量
Figure GDA0003083536180000025
的正半部分的信息，得到向量
Figure GDA0003083536180000026
将所述向量
Figure GDA0003083536180000027
扩展成托普利兹矩阵Rv∈C(2MN-N+1)×(2MN-N+1)，其中，2M为第一ULA的阵元数，N为第二ULA的阵元个数，M、N互为质数，且M<N；
将所述托普利兹矩阵Rv恢复得到填充协方差矩阵
Figure GDA0003083536180000028
对所述填充协方差矩阵
Figure GDA0003083536180000029
进行特征分解，利用填充协方差矩阵
Figure GDA00030835361800000210
的噪声子空间Un，构建MUSIC谱；在期望信号角度区域Θs内进行谱峰搜索得到期望信号的到达角
Figure GDA00030835361800000211
计算得到期望信号的导向矢量
Figure GDA00030835361800000212
在干扰信号角度区域
Figure GDA00030835361800000213
内进行谱峰搜索得到各干扰信号到达角的估计
Figure GDA00030835361800000214
利用估计出的干扰信号的到达角、以及互质阵的物理阵列信息以及重构出互质阵物理阵列的干扰加噪声协方差矩阵
Figure GDA00030835361800000215
利用所述干扰加噪声协方差矩阵
Figure GDA00030835361800000216
和期望信号导向矢量的估计
Figure GDA00030835361800000217
计算自适应波束形成器的加权矢量w。
可选地，利用K个采样快拍计算样本协方差矩阵：
Figure GDA00030835361800000218
其中：K为采样的次数，k为采样的序号，x(k)为互质阵的各个阵元接收到数据，xH(k)为x(k)的共轭转置。
可选地，对所述向量z进行去冗余和向量重排，具体包括：
对具有相同波程差的元素求和取平均，计算如下
Figure GDA0003083536180000031
其中，m表示差分优化阵的波程差，|·|代表求集合中的元素个数，
Figure GDA0003083536180000032
表示协方差矩阵
Figure GDA0003083536180000033
对应波程差为n1-n2处的值，<x>n1表示信号矢量x在阵元参考位置n1处的值，<x>n2表示信号矢量x在阵元参考位置n2处的值；
集合T(m)表示波程差为m的阵元位置二元组的集合
T(m)＝{(n1,n2)∈S2|n1-n2＝m}
将<z1>m按照波程差从小到大的顺序排列，得到对应的完整互质阵差分优化阵的接收数据向量z1∈C3MN+M-N。
可选地，所述将所述接收数据向量z1中波程差不连续的元素之间全部填0，得到向量
Figure GDA0003083536180000034
然后取所述向量
Figure GDA0003083536180000035
的正半部分的信息，得到向量
Figure GDA0003083536180000036
具体包括：
将所述接收数据向量z1中波程差不连续的元素之间全部填0，使其构成一个更高维的接收数据向量
Figure GDA0003083536180000037
即
Figure GDA0003083536180000038
其中，Sdiff代表差分优化阵阵元位置的数值集合，Sdiff＝{n1-n2|n1,n2∈S}，S＝{dn/d,n＝1,2,...,2M+N-1}，d表示半波长距离λ/2；
取所述接收数据向量
Figure GDA0003083536180000039
的正半部分的信息，得到维度为2MN-N+1的复向量
Figure GDA00030835361800000310
Figure GDA00030835361800000311
可选地，构造如下矩阵填充问题，来由托普利兹矩阵Rv恢复得到填充协方差矩阵
Figure GDA00030835361800000312
Figure GDA0003083536180000041
Figure GDA0003083536180000042
Figure GDA0003083536180000043
n1,n2∈S+＝{n|n∈S,n≥0}
其中，rank(·)代表求矩阵的秩，将秩最小化问题转化为核范数最小化问题，如下
Figure GDA0003083536180000044
Figure GDA0003083536180000045
Figure GDA0003083536180000046
n1,n2∈S+
其中，||·||*代表求矩阵的核范数。
可选地，所述各干扰信号到达角的估计
Figure GDA0003083536180000047
通过以下方法获得：
对所述填充协方差矩阵
Figure GDA0003083536180000048
进行特征分解，
Figure GDA0003083536180000049
其中，Us为
Figure GDA00030835361800000410
的信号子空间，Σs为Us中每个向量对应的特征值，Un为
Figure GDA00030835361800000411
的噪声子空间，Σn为Un中每个向量对应的特征值；
利用所述噪声子空间Un，在整个角度区域上构建MUSIC谱，
Figure GDA00030835361800000412
其中，θ∈[-90°,90°]表示搜索中的角度，d(θ)表示填充协方差矩阵
Figure GDA00030835361800000413
对应于角度θ处的导向矢量，d(θ)＝[1,ej2πdsinθ/λ,ej2π2dsinθ/λ...,ej2π(2MN-N)dsinθ/λ]T；
假定期望信号位于区间
Figure GDA00030835361800000414
其中
Figure GDA00030835361800000415
是假定的期望信号的到达角，干扰信号位于Θs的补区间
Figure GDA00030835361800000416
中；
在期望信号区间Θs内对MUSIC谱PMUSIC(θ)进行谱峰搜索，将区域内最大的谱峰作为期望信号的谱峰，并将其角度位置
Figure GDA00030835361800000417
记录为期望信号的到达角的估计；
期望信号的导向矢量的估计如下；
Figure GDA00030835361800000418
在干扰信号区间
Figure GDA00030835361800000419
内对MUSIC谱进行谱峰搜索来获得各干扰信号的到达角信息，选择其中最大的(L-1)个谱峰后，得到各干扰信号到达角的估计
Figure GDA00030835361800000420
可选地，所述干扰加噪声协方差矩阵
Figure GDA00030835361800000421
为，
Figure GDA0003083536180000051
其中，γmin表示样本协方差矩阵
Figure GDA0003083536180000052
的最小特征值，
Figure GDA0003083536180000053
为第l个信源的估计功率，I2M+N-1为维度为2MN+N-1的单位阵；
可选地，所述自适应波束形成器的加权矢量w：
Figure GDA0003083536180000054
如上所述，本发明的一种基于矩阵填充的互质阵稳健自适应波束形成算法，具有以下有益效果：
本发明通过矩阵填充技术，对互质阵的差分优化阵的非连续部分的孔洞进行插值，从而使整个差分优化阵形成一个连续的ULA，使其所有虚拟阵元的数据得到有效利用，提高了算法实际能利用上的互质阵的差分优化阵的自由度。而在自适应波束形成技术中，用于处理接收信号的阵列的自由度越高，干扰信号就能被更准确地定位和抑制，从而获得更好的波束形成性能。
附图说明
为了进一步阐述本发明所描述的内容，下面结合附图对本发明的具体实施方式作进一步详细的说明。应当理解，这些附图仅作为典型示例，而不应看作是对本发明的范围的限定。
图1为互质阵阵列结构；
图2为所测试波束形成器的输出SINR随着输入SNR变化的曲线；
图3为所测试波束形成器的输出SINR随采样快拍数变化的曲线；
图4为一种基于矩阵填充的互质阵稳健自适应波束形成算法的流程图。
具体实施方式
以下通过特定的具体实例说明本发明的实施方式，本领域技术人员可由本说明书所揭露的内容轻易地了解本发明的其他优点与功效。本发明还可以通过另外不同的具体实施方式加以实施或应用，本说明书中的各项细节也可以基于不同观点与应用，在没有背离本发明的精神下进行各种修饰或改变。需说明的是，在不冲突的情况下，以下实施例及实施例中的特征可以相互组合。
需要说明的是，以下实施例中所提供的图示仅以示意方式说明本发明的基本构想，遂图式中仅显示与本发明中有关的组件而非按照实际实施时的组件数目、形状及尺寸绘制，其实际实施时各组件的型态、数量及比例可为一种随意的改变，且其组件布局型态也可能更为复杂。
如图4所示，本发明提供一种基于矩阵填充的互质阵稳健自适应波束形成算法，该互质阵包括第一ULA和第二ULA，所述第一ULA与所述第二ULA互质；所述算法包括以下步骤：
S1计算接收数据的样本协方差矩阵
Figure GDA0003083536180000061
S2将所述样本协方差矩阵
Figure GDA0003083536180000062
向量化得到向量z，然后对所述向量z进行去冗余和向量重排，得到完整互质阵差分优化阵的接收数据向量z1；
S3将所述接收数据向量z1中波程差不连续的元素之间全部填0，得到向量
Figure GDA0003083536180000063
然后取所述向量
Figure GDA0003083536180000064
的正半部分的信息，得到向量
Figure GDA0003083536180000065
S4将所述向量
Figure GDA0003083536180000066
扩展成托普利兹矩阵Rv∈C(2MN-N+1)×(2MN-N+1)，其中，2M为第一ULA的阵元数，N为第二ULA的阵元个数，M、N互为质数，且M<N；
S5将所述托普利兹矩阵Rv恢复得到填充协方差矩阵
Figure GDA0003083536180000067
S6对所述填充协方差矩阵
Figure GDA0003083536180000068
进行特征分解，利用填充协方差矩阵
Figure GDA0003083536180000069
的噪声子空间Un，构建MUSIC谱；在期望信号角度区域Θs内进行谱峰搜索得到期望信号的到达角
Figure GDA00030835361800000610
计算得到期望信号的导向矢量
Figure GDA00030835361800000611
在干扰信号角度区域
Figure GDA00030835361800000612
内进行谱峰搜索得到各干扰信号到达角的估计
Figure GDA00030835361800000613
S7利用估计出的干扰信号的到达角、以及互质阵的物理阵列信息以及重构出互质阵物理阵列的干扰加噪声协方差矩阵
Figure GDA00030835361800000614
S8利用所述干扰加噪声协方差矩阵
Figure GDA00030835361800000615
和期望信号导向矢量的估计
Figure GDA00030835361800000616
计算自适应波束形成器的加权矢量w。
于一实施例中，在步骤S1中，利用K个采样快拍计算样本协方差矩阵：
Figure GDA00030835361800000617
其中：K为采样的次数，k为采样的序号，x(k)为互质阵的各个阵元接收到数据，xH(k)为x(k)的共轭转置。
于一实施例中，在步骤S2中，对所述向量z进行去冗余和向量重排，具体包括：
对具有相同波程差的元素求和取平均，计算如下：
Figure GDA0003083536180000071
其中，m表示差分优化阵的波程差，|·|代表求集合中的元素个数，<x>i代表一个信号x在阵元参考位置i处的值，
Figure GDA0003083536180000072
代表协方差矩阵R对应波程差为i-j处的值；
Figure GDA0003083536180000073
代表协方差矩阵
Figure GDA0003083536180000074
对应波程差为n1-n2处的值，<x>n1代表信号矢量x在阵元参考位置n1处的值，
Figure GDA0003083536180000075
代表信号矢量x在阵元参考位置n2处的值，
集合T(m)表示波程差为m的阵元位置二元组的集合
T(m)＝{(n1,n2)∈S2|n1-n2＝m}
将<z1>m按照波程差从小到大的顺序排列，得到对应的完整互质阵差分优化阵的接收数据向量z1∈C3MN+M-N。
于一实施例中，在步骤S3中，所述将所述接收数据向量z1中波程差不连续的元素之间全部填0，得到向量
Figure GDA0003083536180000076
然后取所述向量
Figure GDA0003083536180000077
的正半部分的信息，得到向量
Figure GDA0003083536180000078
具体包括：
将所述接收数据向量z1中波程差不连续的元素之间全部填0，使其构成一个更高维的接收数据向量
Figure GDA0003083536180000079
即
Figure GDA00030835361800000710
其中，Sdiff代表差分优化阵阵元位置的数值集合，Sdiff＝{n1-n2|n1,n2∈S}，S＝{dn/d,n＝1,2,...,2M+N-1}，d表示半波长距离λ/2；
取所述接收数据向量
Figure GDA00030835361800000711
的正半部分的信息，得到维度为2MN-N+1的复向量
Figure GDA00030835361800000712
Figure GDA00030835361800000713
于一实施例中，在步骤S5中，构造如下矩阵填充问题，来由托普利兹矩阵Rv恢复得到填充协方差矩阵
Figure GDA00030835361800000714
Figure GDA0003083536180000081
Figure GDA0003083536180000082
Figure GDA0003083536180000083
n1,n2∈S+＝{n|n∈S,n≥0}
其中，rank(·)代表求矩阵的秩，将秩最小化问题转化为核范数最小化问题，如下
Figure GDA0003083536180000084
Figure GDA0003083536180000085
Figure GDA0003083536180000086
n1,n2∈S+
其中，||·||*代表求矩阵的核范数。
于一实施例中，在步骤S6中，所述各干扰信号到达角的估计
Figure GDA0003083536180000087
通过以下方法获得：
对所述填充协方差矩阵
Figure GDA0003083536180000088
进行特征分解，
Figure GDA0003083536180000089
其中，Us为
Figure GDA00030835361800000810
的信号子空间，Σs为Us中每个向量对应的特征值，Un为
Figure GDA00030835361800000811
的噪声子空间，Σn为Un中每个向量对应的特征值；
利用所述噪声子空间Un，在整个角度区域上构建MUSIC谱，
Figure GDA00030835361800000812
其中，θ∈[-90°,90°]表示搜索中的角度，d(θ)表示填充协方差矩阵
Figure GDA00030835361800000813
对应于角度θ处的导向矢量，d(θ)＝[1，ej2πdsinθλ，ej2π2dsinθλθ/λ...，ej2π(2MN-N)dsinθλ]T；
假定期望信号位于区间
Figure GDA00030835361800000815
其中
Figure GDA00030835361800000816
是假定的期望信号的到达角，干扰信号位于Θs的补区间
Figure GDA00030835361800000817
中；
在期望信号区间Θs内对MUSIC谱PMUSIC(θ)进行谱峰搜索，将区域内最大的谱峰作为期望信号的谱峰，并将其角度位置
Figure GDA00030835361800000818
记录为期望信号的到达角的估计；
期望信号的导向矢量的估计如下；
Figure GDA0003083536180000091
在干扰信号区间
Figure GDA0003083536180000092
内对MUSIC谱进行谱峰搜索来获得各干扰信号的到达角信息，选择其中最大的(L-1)个谱峰后，得到各干扰信号到达角的估计
Figure GDA0003083536180000093
于一实施例中，在步骤S7中，所述干扰加噪声协方差矩阵
Figure GDA0003083536180000094
为，
Figure GDA0003083536180000095
其中，γmin表示样本协方差矩阵
Figure GDA0003083536180000096
的最小特征值，
Figure GDA0003083536180000097
为第l个信源的估计功率，I2M+N-1为维度为2MN+N-1的单位阵；
于一实施例中，在步骤S8中，所述自适应波束形成器的加权矢量w：
Figure GDA0003083536180000098
本发明通过矩阵填充技术，对互质阵的差分优化阵的非连续部分的孔洞进行插值，从而使整个差分优化阵形成一个连续的ULA，使其所有虚拟阵元的数据得到有效利用，提高了算法实际能利用上的互质阵的差分优化阵的自由度。而在自适应波束形成技术中，用于处理接收信号的阵列的自由度越高，干扰信号就能被更准确地定位和抑制，从而获得更好的波束形成性能。
下面以具体的数据来对本发明进行说明，
一种基于矩阵填充的互质阵稳健自适应波束形成算法，包括以下步骤：
步骤11.设置天线阵列：
天线阵列为一互质阵，该互质阵由一对阵元间距互质的ULA(规则线性阵列)构造而成。确定一对互质的数M＝3、N＝5，M<N。第一个ULA由2M＝6个阵元间距为Nd＝5d的阵元构成，第二个ULA由N＝5个阵元间距为Md＝3d的阵元构成，其中d代表半波长距离λ/2，本实施方式信号波长λ＝0.375m，所以d＝0.1875m。构成的互质阵最终一共有(2M+N-1＝10)个阵元。现有L＝4个远场窄带信号入射到互质阵上，其DOA分别为{θ1＝0°,θ2＝-30°,θ3＝30°,θ4＝45°}，期望信号的DOA为θ1，其余DOA为干扰信号的到达角。
步骤12：计算接收数据的样本协方差矩阵
Figure GDA0003083536180000099
利用500个采样快拍计算样本协方差矩阵如下：
Figure GDA00030835361800000910
其中：k为采样的序号，x(k)为互质阵的各个阵元接收到数据按阵元的顺序排成一列。
接下来对样本协方差矩阵
Figure GDA0003083536180000101
进行特征分解计算其最小的特征值γmin，以备后续计算。
步骤13：将样本协方差矩阵
Figure GDA0003083536180000102
向量化得到向量z，
Figure GDA0003083536180000103
对向量z根据对应的波程差进行去冗余和向量重排，得到完整互质阵差分优化阵的接收数据向量z1。具体操作如下：
对于相同波程差的元素进行去冗余的时候是对这些元素求和取平均，即
Figure GDA0003083536180000104
将<z1>m按照波程差从小到大的顺序排列，便可得到对应的完整互质阵差分优化阵的接收数据向量z1∈C43。
步骤14：将接收数据向量z1中波程差不连续的元素之间全部填0，使其构成一个更高维的接收数据向量
Figure GDA0003083536180000105
即
Figure GDA0003083536180000106
取接收数据向量
Figure GDA0003083536180000107
的正半部分的信息，得到向量
Figure GDA0003083536180000108
Figure GDA0003083536180000109
步骤15：根据信号协方差矩阵共轭对称的特性，将向量
Figure GDA00030835361800001010
扩展成托普利兹矩阵Rv∈C26×26，
Figure GDA00030835361800001011
步骤16：由于接收数据协方差矩阵的低秩特性，托普利兹矩阵Rv中的0元素可以通过构造如下核范数最小化问题来得到有效恢复，从而得到填充协方差矩阵
Figure GDA00030835361800001012
Figure GDA00030835361800001013
Figure GDA00030835361800001014
Figure GDA00030835361800001015
n1,n2∈S+
核范数最小化问题是一个凸优化问题，通过半正定规划可以有效地进行求解，可以方便地利用凸优化工具包CVX直接对其求解。
步骤17：对填充协方差矩阵
Figure GDA0003083536180000111
进行特征分解，利用其噪声子空间Un，构建MUSIC谱。在期望信号角度区域Θs内进行谱峰搜索得到期望信号的到达角
Figure GDA0003083536180000112
计算得到期望信号的导向矢量
Figure GDA0003083536180000113
在干扰信号角度区域
Figure GDA0003083536180000114
内进行谱峰搜索得到各干扰信号到达角的估计
Figure GDA0003083536180000115
具体步骤如下：
步骤1701：对填充协方差矩阵
Figure GDA0003083536180000116
进行特征分解，
Figure GDA0003083536180000117
其中，Us由
Figure GDA0003083536180000118
的主要特征向量组成，为
Figure GDA0003083536180000119
的信号子空间，Σs为Us中每个向量对应的特征值；Un由
Figure GDA00030835361800001110
的小特征向量组成，为
Figure GDA00030835361800001111
的噪声子空间，Σn为Un中每个向量对应的特征值。
步骤1702：利用噪声子空间Un，在整个角度区域上构建MUSIC谱，
Figure GDA00030835361800001112
其中，θ∈[-90°,90°]表示搜索中的角度，d(θ)表示填充协方差矩阵
Figure GDA00030835361800001113
对应于角度θ处的导向矢量，d(θ)＝[1,ejπsinθ,ej2πsinθ,...,ej25πsinθ]T。
步骤1703：在本实施例中，Θs＝[-5°,5°]，干扰信号位于Θs的补区间
Figure GDA00030835361800001114
中。
首先，在期望信号区间Θs内对MUSIC谱PMUSIC(θ)进行谱峰搜索，将区域内最大的谱峰作为期望信号的谱峰，并将其角度位置
Figure GDA00030835361800001115
记录为期望信号的到达角的估计；在输入信噪比较低时，在Θs内往往是难以搜索到期望信号的谱峰的，此时用假定的期望信号角度作为期望信号到达角的估计，计算期望信号的导向矢量如下。
Figure GDA00030835361800001116
接下来，在干扰信号区间
Figure GDA00030835361800001117
内对MUSIC谱进行谱峰搜索来获得各干扰信号的到达角信息。选择其中最大的3个谱峰后，可得到各干扰信号到达角的估计
Figure GDA00030835361800001118
步骤18：令
Figure GDA00030835361800001119
代表估计的各信源DOA的集合。接下来，利用θm可以构造出如下的最小二乘问题
Figure GDA00030835361800001120
subject to p(θm)＞0,
其中，
Figure GDA00030835361800001121
代表集合θm上的功率分布，
Figure GDA00030835361800001122
为第i个信源的估计功率，是优化问题求解的目标。diag{·}代表向量的对角化操作。
Figure GDA0003083536180000121
为θm对应的阵列流型。噪声功率的估计
Figure GDA0003083536180000122
取为样本协方差矩阵
Figure GDA0003083536180000123
的最小特征值。实际上，在求解上述优化问题时，可以忽略不等式约束，从而得到如下的最小二乘闭式解
p(θm)＝(GHG)-1GHr
其中，
Figure GDA0003083536180000124
步骤19：利用估计出的干扰信号的到达角
Figure GDA0003083536180000125
以及互质阵的物理阵列信息，重构出互质阵物理阵列的干扰加噪声协方差矩阵
Figure GDA0003083536180000126
Figure GDA0003083536180000127
步骤110：利用干扰加噪声协方差矩阵
Figure GDA0003083536180000128
和期望信号导向矢量的估计
Figure GDA0003083536180000129
计算自适应波束形成器的加权矢量w：
Figure GDA00030835361800001210
为了进一步验证上述实施方式的有效性，设计了一个仿真实验，在该仿真实验中考虑期望信号导向矢量被准确获知。在仿真实验中，加性噪声为复高斯零均值白噪声，并且在每个阵元处拥有相同的方差。在样本数据中总是含有期望信号成分。干扰信号在每个传感器处的干噪比INR都等于30dB。在分析输出信干噪比随快拍数变化的实验中，输入信噪比固定为SNR＝20dB；在分析输出信干噪比随输入信噪比变化的实验中，快拍数固定为K＝500。为了获取各实验中的每个数据点，进行了500次蒙特卡洛实验。
图2描绘了所测试波束形成器的输出SINR随着输入SNR变化的曲线，图3为所测试波束形成器的输出SINR随采样快拍数变化的曲线。可以看出，所提出算法性能依然优于其他对比算法，并且收敛速度要快于其他对比算法。
上述实施例仅例示性说明本发明的原理及其功效，而非用于限制本发明。任何熟悉此技术的人士皆可在不违背本发明的精神及范畴下，对上述实施例进行修饰或改变。因此，举凡所属技术领域中具有通常知识者在未脱离本发明所揭示的精神与技术思想下所完成的一切等效修饰或改变，仍应由本发明的权利要求所涵盖。

Claims (7)

1.一种基于矩阵填充的互质阵稳健自适应波束形成算法，其特征在于，该互质阵包括第一ULA和第二ULA，所述第一ULA与所述第二ULA互质；所述算法包括以下步骤：
计算接收数据的样本协方差矩阵
Figure FDA0003370289740000011
将所述样本协方差矩阵
Figure FDA0003370289740000012
向量化得到向量z，然后对所述向量z进行去冗余和向量重排，得到完整互质阵差分优化阵的接收数据向量z1；
将所述接收数据向量z1中波程差不连续的元素之间全部填0，得到向量
Figure FDA0003370289740000013
然后取所述向量
Figure FDA0003370289740000014
的正半部分的信息，得到向量
Figure FDA0003370289740000015
将所述向量
Figure FDA0003370289740000016
扩展成托普利兹矩阵Rv∈C(2MN-N+1)×(2MN-N+1)，其中，2M为第一ULA的阵元数，N为第二ULA的阵元个数，M、N互为质数，且M<N；
将所述托普利兹矩阵Rv恢复得到填充协方差矩阵
Figure FDA0003370289740000017
对所述填充协方差矩阵
Figure FDA0003370289740000018
进行特征分解，利用填充协方差矩阵
Figure FDA0003370289740000019
的噪声子空间Un，构建MUSIC谱；在期望信号角度区域Θs内进行谱峰搜索得到期望信号的到达角
Figure FDA00033702897400000110
计算得到期望信号的导向矢量
Figure FDA00033702897400000111
在干扰信号角度区域
Figure FDA00033702897400000112
内进行谱峰搜索得到各干扰信号到达角的估计
Figure FDA00033702897400000113
利用估计出的干扰信号的到达角、以及互质阵的物理阵列信息以及重构出互质阵物理阵列的干扰加噪声协方差矩阵
Figure FDA00033702897400000114
利用所述干扰加噪声协方差矩阵
Figure FDA00033702897400000115
和期望信号导向矢量的估计
Figure FDA00033702897400000116
计算自适应波束形成器的加权矢量w；
所述干扰加噪声协方差矩阵
Figure FDA00033702897400000117
为，
Figure FDA00033702897400000118
其中，γmin表示样本协方差矩阵
Figure FDA00033702897400000119
的最小特征值，
Figure FDA00033702897400000120
为第l个信源的估计功率，I2M+N-1为维度为2MN+N-1的单位阵。
2.根据权利要求1所述一种基于矩阵填充的互质阵稳健自适应波束形成算法，其特征在于，利用K个采样快拍计算样本协方差矩阵：
Figure FDA00033702897400000121
其中：K为采样的次数，k为采样的序号，x(k)为互质阵的各个阵元接收到数据，xH(k)为x(k)的共轭转置。
3.根据权利要求2所述的一种基于矩阵填充的互质阵稳健自适应波束形成算法，其特征在于，对所述向量z进行去冗余和向量重排，具体包括：
对具有相同波程差的元素求和取平均，计算如下
Figure FDA0003370289740000021
其中，m表示差分优化阵的波程差，|·|代表求集合中的元素个数，
Figure FDA0003370289740000022
表示协方差矩阵
Figure FDA0003370289740000023
对应波程差为n1-n2处的值，
Figure FDA0003370289740000024
表示信号矢量x在阵元参考位置n1处的值，
Figure FDA0003370289740000025
表示信号矢量x在阵元参考位置n2处的值；
集合T(m)表示波程差为m的阵元位置二元组的集合
T(m)＝{(n1,n2)∈S2|n1-n2＝m}
将<z1>m按照波程差从小到大的顺序排列，得到对应的完整互质阵差分优化阵的接收数据向量z1∈C3MN+M-N。
4.根据权利要求3所述的一种基于矩阵填充的互质阵稳健自适应波束形成算法，其特征在于，所述将所述接收数据向量z1中波程差不连续的元素之间全部填0，得到向量
Figure FDA0003370289740000026
然后取所述向量
Figure FDA0003370289740000027
的正半部分的信息，得到向量
Figure FDA0003370289740000028
具体包括：
将所述接收数据向量z1中波程差不连续的元素之间全部填0，使其构成一个更高维的接收数据向量
Figure FDA0003370289740000029
即
Figure FDA00033702897400000210
其中，Sdiff代表差分优化阵阵元位置的数值集合，Sdiff＝{n1-n2|n1,n2∈S}，S＝{dn/d,n＝1,2,...,2M+N-1}，d表示半波长距离λ/2；
取所述接收数据向量
Figure FDA00033702897400000211
的正半部分的信息，得到维度为2MN-N+1的复向量
Figure FDA00033702897400000212
5.根据权利要求4所述的一种基于矩阵填充的互质阵稳健自适应波束形成算法，其特征在于，构造如下矩阵填充问题，来由托普利兹矩阵Rv恢复得到填充协方差矩阵
Figure FDA00033702897400000213
Figure FDA0003370289740000031
Figure FDA0003370289740000032
Figure FDA0003370289740000033
n1,n2∈S+＝{n|n∈S,n≥0}
其中，rank(·)代表求矩阵的秩，将秩最小化问题转化为核范数最小化问题，如下
Figure FDA0003370289740000034
Figure FDA0003370289740000035
Figure FDA0003370289740000036
n1,n2∈S+
其中，||·||*代表求矩阵的核范数。
6.根据权利要求5所述的一种基于矩阵填充的互质阵稳健自适应波束形成算法，其特征在于，所述各干扰信号到达角的估计
Figure FDA0003370289740000037
通过以下方法获得：
对所述填充协方差矩阵
Figure FDA0003370289740000038
进行特征分解，
Figure FDA0003370289740000039
其中，Us为
Figure FDA00033702897400000310
的信号子空间，Σs为Us中每个向量对应的特征值，Un为
Figure FDA00033702897400000311
的噪声子空间，Σn为Un中每个向量对应的特征值；
利用所述噪声子空间Un，在整个角度区域上构建MUSIC谱，
Figure FDA00033702897400000312
其中，θ∈[-90°,90°]表示搜索中的角度，d(θ)表示填充协方差矩阵
Figure FDA00033702897400000313
对应于角度θ处的导向矢量，d(θ)＝[1,ej2πdsinθ/λ,ej2π2dsinθ/λ...,ej2π(2MN-N)dsinθ/λ]T；
假定期望信号位于区间
Figure FDA00033702897400000314
其中
Figure FDA00033702897400000315
是假定的期望信号的到达角，干扰信号位于Θs的补区间
Figure FDA00033702897400000316
中；
在期望信号区间Θs内对MUSIC谱PMUSIC(θ)进行谱峰搜索，将区域内最大的谱峰作为期望信号的谱峰，并将其角度位置
Figure FDA00033702897400000317
记录为期望信号的到达角的估计；
期望信号的导向矢量的估计如下；
Figure FDA00033702897400000318
在干扰信号区间
Figure FDA0003370289740000041
内对MUSIC谱进行谱峰搜索来获得各干扰信号的到达角信息，选择其中最大的(L-1)个谱峰后，得到各干扰信号到达角的估计
Figure FDA0003370289740000042
7.根据权利要求1所述的一种基于矩阵填充的互质阵稳健自适应波束形成算法，其特征在于，所述自适应波束形成器的加权矢量w：
Figure FDA0003370289740000043
CN201910194101.XA 2019-03-14 2019-03-14 一种基于矩阵填充的互质阵稳健自适应波束形成算法 Expired - Fee Related CN110045323B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910194101.XA CN110045323B (zh) 2019-03-14 2019-03-14 一种基于矩阵填充的互质阵稳健自适应波束形成算法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910194101.XA CN110045323B (zh) 2019-03-14 2019-03-14 一种基于矩阵填充的互质阵稳健自适应波束形成算法

Publications (2)

Publication Number Publication Date
CN110045323A CN110045323A (zh) 2019-07-23
CN110045323B true CN110045323B (zh) 2022-02-18

Family

ID=67274691

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910194101.XA Expired - Fee Related CN110045323B (zh) 2019-03-14 2019-03-14 一种基于矩阵填充的互质阵稳健自适应波束形成算法

Country Status (1)

Country Link
CN (1) CN110045323B (zh)

Families Citing this family (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110554352A (zh) * 2019-09-11 2019-12-10 哈尔滨工业大学 基于vgg16神经网络的航天测控系统干扰源波达方向估计方法
CN111239714B (zh) * 2019-09-18 2022-06-07 中国人民解放军海军工程大学 柔性阵波束形成稳健化实现方法
CN112711013B (zh) * 2020-12-14 2022-10-21 中国船舶重工集团公司第七一五研究所 一种基于分块矩阵的快速自适应波束形成方法
CN113687395B (zh) * 2021-08-25 2023-04-18 电子科技大学 基于最大熵功率谱的卫星导航信号导向矢量估计方法
CN113779795B (zh) * 2021-09-13 2023-06-30 中国科学院声学研究所 一种阵列设计方法及装置
CN114019446B (zh) * 2021-10-19 2024-04-12 南京航空航天大学 一种基于去噪核范数最小化的互质相干信源估计方法
CN114994651B (zh) * 2022-05-18 2024-02-06 电子科技大学 一种存在通道幅相误差的双基地互质mimo雷达目标定位方法
CN114994595A (zh) * 2022-08-03 2022-09-02 四川太赫兹通信有限公司 波达方向获取方法、太赫兹相控阵波束控制方法和系统
CN116500625B (zh) * 2023-06-29 2023-10-20 天津知海科技有限公司 恢复成像方法、装置、系统、电子设备及可读存储介质

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8289211B1 (en) * 2010-04-14 2012-10-16 Lockheed Martin Corporation Method and apparatus for ranging
CN106569171A (zh) * 2016-11-08 2017-04-19 西安电子科技大学 基于双层混合阵列的波达方向角估计方法
CN108614234A (zh) * 2018-05-15 2018-10-02 浙江大学 基于多采样快拍互质阵列接收信号快速傅里叶逆变换的波达方向估计方法
CN108710758A (zh) * 2018-05-18 2018-10-26 电子科技大学 基于嵌套阵和协方差矩阵重构的自适应波束形成算法
CN108983168A (zh) * 2018-04-27 2018-12-11 常熟理工学院 基于孔径补全的压缩感知mimo雷达目标探测方法

Family Cites Families (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
TWI457585B (zh) * 2012-12-11 2014-10-21 Univ Nat Chiao Tung 角度偵測方法及其裝置
US10054666B2 (en) * 2014-09-26 2018-08-21 The United States Of America, As Represented By The Secretary Of The Navy Sparse space-time adaptive array architecture
CN105403856B (zh) * 2015-10-30 2017-10-24 西安电子科技大学 基于嵌套式最小冗余阵列的波达方向估计方法
WO2018045594A1 (zh) * 2016-09-12 2018-03-15 深圳大学 一种基于互质脉冲重复间隔的空时自适应处理方法及装置
CN108828551B (zh) * 2018-08-28 2022-02-22 中国人民解放军空军工程大学 一种基于压缩感知的灵活mimo雷达混合目标doa估计方法
CN109407045B (zh) * 2018-10-10 2020-05-22 苏州大学 一种非均匀传感器阵列宽带信号波达方向估计方法

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8289211B1 (en) * 2010-04-14 2012-10-16 Lockheed Martin Corporation Method and apparatus for ranging
CN106569171A (zh) * 2016-11-08 2017-04-19 西安电子科技大学 基于双层混合阵列的波达方向角估计方法
CN108983168A (zh) * 2018-04-27 2018-12-11 常熟理工学院 基于孔径补全的压缩感知mimo雷达目标探测方法
CN108614234A (zh) * 2018-05-15 2018-10-02 浙江大学 基于多采样快拍互质阵列接收信号快速傅里叶逆变换的波达方向估计方法
CN108710758A (zh) * 2018-05-18 2018-10-26 电子科技大学 基于嵌套阵和协方差矩阵重构的自适应波束形成算法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
"基于干扰矩阵重构的稳健自适应波束形成算法";杨金金 等;《信息工程大学学报》;20140630;第15卷(第3期);正文第1-2节 *
"基于矩阵填充的互质阵列欠定DOA估计方法";吴晨曦 等;《工程科学与技术》;20170930;第49卷(第5期);正文第2.1-2.3节 *

Also Published As

Publication number Publication date
CN110045323A (zh) 2019-07-23

Similar Documents

Publication Publication Date Title
CN110045323B (zh) 一种基于矩阵填充的互质阵稳健自适应波束形成算法
CN106980106B (zh) 阵元互耦下的稀疏doa估计方法
CN110113085B (zh) 一种基于协方差矩阵重构的波束形成方法及系统
CN109298383B (zh) 一种基于变分贝叶斯推断的互质阵波达方向角估计方法
WO2018094565A1 (zh) 脉冲噪声下的波束成形方法及装置
CN108375763B (zh) 一种应用于多声源环境的分频定位方法
CN107450047B (zh) 嵌套阵下基于未知互耦信息的压缩感知doa估计方法
CN111239678B (zh) 一种基于l型阵列的二维doa估计方法
CN107576931B (zh) 一种基于协方差低维度迭代稀疏重构的相关/相干信号波达方向估计方法
CN112379327A (zh) 一种基于秩损估计的二维doa估计与互耦校正方法
CN107290732B (zh) 一种量子大爆炸的单基地mimo雷达测向方法
CN111693947A (zh) 基于互质阵列doa估计的改进music方法
CN110727915B (zh) 一种基于数据相关约束的鲁棒自适应波束形成方法
Yan et al. MUSIC-like direction of arrival estimation based on virtual array transformation
CN115236584A (zh) 基于深度学习的米波雷达低仰角估计方法
Dong et al. Enhanced bnc approach for noncircular signals direction finding with sparse array in the scenario of impulsive noise
Anughna et al. Antenna reconfiguration based DOA estimation for AWGN channel in MIMO applications
CN107135026B (zh) 未知互耦存在时基于矩阵重构的稳健波束形成方法
CN113391266B (zh) 基于非圆多嵌套阵降维子空间数据融合的直接定位方法
CN106877918B (zh) 互耦条件下稳健自适应波束形成方法
CN114184999B (zh) 一种互耦小孔径阵列的生成式模型处理方法
CN115421098A (zh) 嵌套面阵下降维求根music的二维doa估计方法
CN114755628A (zh) 非均匀噪声下声矢量传感器阵列波达方向估计方法
CN112698263A (zh) 一种基于正交传播算子的单基地互质mimo阵列doa估计算法
CN113281698A (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
CF01 Termination of patent right due to non-payment of annual fee
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20220218