CN117491985B - 一种波达方向及速度联合估计方法及相关装置 - Google Patents
一种波达方向及速度联合估计方法及相关装置 Download PDFInfo
- Publication number
- CN117491985B CN117491985B CN202311831328.3A CN202311831328A CN117491985B CN 117491985 B CN117491985 B CN 117491985B CN 202311831328 A CN202311831328 A CN 202311831328A CN 117491985 B CN117491985 B CN 117491985B
- Authority
- CN
- China
- Prior art keywords
- space
- signal
- time
- matrix
- target
- 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
Links
- 238000000034 method Methods 0.000 title claims abstract description 92
- 238000012545 processing Methods 0.000 claims abstract description 67
- 238000009499 grossing Methods 0.000 claims abstract description 30
- 239000011159 matrix material Substances 0.000 claims description 182
- 239000013598 vector Substances 0.000 claims description 46
- 230000008569 process Effects 0.000 claims description 29
- 238000001514 detection method Methods 0.000 claims description 27
- 230000000903 blocking effect Effects 0.000 claims description 19
- 238000010276 construction Methods 0.000 claims description 15
- 238000000354 decomposition reaction Methods 0.000 claims description 14
- 238000004364 calculation method Methods 0.000 claims description 10
- 230000008859 change Effects 0.000 claims description 7
- 229940088594 vitamin Drugs 0.000 claims description 3
- 229930003231 vitamin Natural products 0.000 claims description 3
- 235000013343 vitamin Nutrition 0.000 claims description 3
- 239000011782 vitamin Substances 0.000 claims description 3
- 150000003722 vitamin derivatives Chemical class 0.000 claims description 3
- 238000000638 solvent extraction Methods 0.000 claims 2
- 230000008447 perception Effects 0.000 abstract description 2
- 230000002829 reductive effect Effects 0.000 description 9
- 230000000694 effects Effects 0.000 description 8
- 230000033001 locomotion Effects 0.000 description 6
- 238000004458 analytical method Methods 0.000 description 5
- 238000005070 sampling Methods 0.000 description 5
- 238000012360 testing method Methods 0.000 description 5
- 230000001427 coherent effect Effects 0.000 description 4
- 238000004891 communication Methods 0.000 description 4
- 238000002592 echocardiography Methods 0.000 description 4
- 239000000284 extract Substances 0.000 description 4
- 230000006870 function Effects 0.000 description 4
- 238000005259 measurement Methods 0.000 description 4
- 238000003491 array Methods 0.000 description 3
- 238000005516 engineering process Methods 0.000 description 3
- 238000001914 filtration Methods 0.000 description 3
- 230000004044 response Effects 0.000 description 3
- 230000009466 transformation Effects 0.000 description 3
- 238000012935 Averaging Methods 0.000 description 2
- 238000013459 approach Methods 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 230000002708 enhancing effect Effects 0.000 description 2
- 230000007613 environmental effect Effects 0.000 description 2
- 238000003384 imaging method Methods 0.000 description 2
- 230000001965 increasing effect Effects 0.000 description 2
- 230000000670 limiting effect Effects 0.000 description 2
- 238000001208 nuclear magnetic resonance pulse sequence Methods 0.000 description 2
- 230000003287 optical effect Effects 0.000 description 2
- 230000002123 temporal effect Effects 0.000 description 2
- 230000004304 visual acuity Effects 0.000 description 2
- 238000009825 accumulation Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000006835 compression Effects 0.000 description 1
- 238000007906 compression Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 230000001419 dependent effect Effects 0.000 description 1
- 230000001747 exhibiting effect Effects 0.000 description 1
- 238000000605 extraction Methods 0.000 description 1
- 230000004438 eyesight Effects 0.000 description 1
- 230000002349 favourable effect Effects 0.000 description 1
- 238000009432 framing Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 239000013307 optical fiber Substances 0.000 description 1
- 238000005457 optimization Methods 0.000 description 1
- 230000036961 partial effect Effects 0.000 description 1
- 238000007781 pre-processing Methods 0.000 description 1
- 238000003672 processing method Methods 0.000 description 1
- 230000006798 recombination Effects 0.000 description 1
- 238000005215 recombination Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 230000000717 retained effect Effects 0.000 description 1
- 230000035945 sensitivity Effects 0.000 description 1
- 239000007787 solid Substances 0.000 description 1
- 238000001228 spectrum Methods 0.000 description 1
- 238000010183 spectrum analysis Methods 0.000 description 1
- 230000003068 static effect Effects 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- 230000002195 synergetic effect Effects 0.000 description 1
- 230000000007 visual effect Effects 0.000 description 1
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
- G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
- G01S13/02—Systems using reflection of radio waves, e.g. primary radar systems; Analogous systems
- G01S13/50—Systems of measurement based on relative movement of target
- G01S13/58—Velocity or trajectory determination systems; Sense-of-movement determination systems
-
- 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
- G01S3/782—Systems for determining direction or deviation from predetermined direction
-
- 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
- G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
- G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
- G01S7/41—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section
- G01S7/415—Identification of targets based on measurements of movement associated with the target
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/10—Pre-processing; Data cleansing
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
- G06F18/213—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods
- G06F18/2131—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods based on a transform domain processing, e.g. wavelet transform
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
- G06F18/213—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods
- G06F18/2132—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods based on discrimination criteria, e.g. discriminant analysis
- G06F18/21322—Rendering the within-class scatter matrix non-singular
- G06F18/21328—Rendering the within-class scatter matrix non-singular involving subspace restrictions, e.g. nullspace techniques
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2123/00—Data types
- G06F2123/02—Data types in the time domain, e.g. time-series data
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
- G06F2218/02—Preprocessing
- G06F2218/04—Denoising
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
- G06F2218/08—Feature extraction
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Data Mining & Analysis (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Radar, Positioning & Navigation (AREA)
- Remote Sensing (AREA)
- Artificial Intelligence (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Life Sciences & Earth Sciences (AREA)
- Evolutionary Biology (AREA)
- Evolutionary Computation (AREA)
- General Engineering & Computer Science (AREA)
- Computer Networks & Wireless Communication (AREA)
- Electromagnetism (AREA)
- Mathematical Physics (AREA)
- Radar Systems Or Details Thereof (AREA)
Abstract
本申请公开了一种波达方向及速度联合估计方法及相关装置，该方法包括：获取雷达回波信号对应的目标距离门，基于目标距离门确定空时二维信号；利用多级正交投影方式构建空时二维信号对应的匹配滤波器集合；对空时二维信号进行平滑处理以得到空时信号，并基于匹配滤波器集合以及空时信号，确定空时二维信号对应的信源数；基于匹配滤波器集合和所述信源数确定多普勒信号子空间；通过多普勒信号子空间求解空时二维信号，得到速度和波达方向。本申请通过获取并处理雷达回波信号来分辨在相同距离门下的不同目标，通过复层次信号处理流程，提取目标的速度信息，增强目标波达方向的分辨率，满足自动驾驶中对于复杂环境的高精度感知需求。
Description
技术领域
本申请涉及信号处理领域，特别涉及一种波达方向及速度联合估计方法及相关装置。
背景技术
在自动驾驶技术进程中，毫米波雷达因其全天候全天时的探测能力成为关键传感器之一。毫米波雷达提供的点云数据为电动车辆提供了所需的关键环境信息，这包含了目标的距离、速度、方位角度以及俯仰角度等参数。然而，目标角度测量的准确性尤其重要，但却受到了现有技术所使用的天线孔径引起的波束宽度的限制。传统测角方法（例如快速傅里叶变换FFT）无法区分波束内的多个对象，这导致了在复杂环境下的探测分辨率不足。而基于子空间分解的MUSIC算法虽然能够提供较高的角度分辨率，但是其对于计算资源的需求较高，算法复杂度大，难以在实际工程中实时实现，对于需要快速并精确处理雷达数据的高级自动驾驶应用来说是存在局限性。
因而现有技术还有待改进和提高。
发明内容
本申请要解决的技术问题在于，针对现有技术的不足，提供一种波达方向及速度联合估计方法及相关装置。
为了解决上述技术问题，本申请实施例第一方面提供了一种波达方向及速度联合估计方法，其中，所述的波达方向及速度联合估计方法具体包括：
获取雷达回波信号对应的目标距离门，基于目标距离门确定空时二维信号；
利用多级正交投影方式构建空时二维信号对应的匹配滤波器集合；
对所述空时二维信号进行平滑处理以得到空时信号，并基于所述匹配滤波器集合以及所述空时信号，确定所述空时二维信号对应的信源数；
基于所述匹配滤波器集合和所述信源数确定多普勒信号子空间；
通过多普勒信号子空间求解所述空时二维信号，得到速度和波达方向。
根据上述技术手段，通过获取雷达回波信号并确定空时二维信号，利用多级正交投影构建匹配滤波器集合，并通过平滑处理与信源数确定拓展了波达方向估计的精确度，从而通过多普勒信号子空间求解得到速度和波达方向，实现了对复杂环境中目标的高精度波达方向与速度估计。
在一个实现方式中，所述获取雷达回波信号对应的目标距离门，基于目标距离门确定空时二维信号具体包括：
对雷达回波信号进行傅里叶处理和恒虚警率检测，得到目标距离门；
提取所述目标距离门对应的回波脉冲，并将提取到的所有回波脉冲排列为空时二维信号。
根据上述技术手段，通过对雷达回波信号的傅里叶处理和恒虚警率检测，得到目标距离门，并在此基础上利用回波脉冲构建空时二维信号，增强了对目标信号的提取效率与精确性。
在一个实现方式中，所述利用多级正交投影方式构建空时二维信号对应的匹配滤波器集合具体包括：
将空时二维信号中的位于最后的阵元的接收信息作为参考信息，将所述空时二维信号作为观测信息；
基于所述参考信息和所述观测信息的互相关构建所述空时二维信号对应的匹配滤波器集合。
根据上述技术手段，通过从空时二维信号中构建参考信息和观测信息，并利用互相关技术，构建适应性强的匹配滤波器集合，优化空时信号的处理。
在一个实现方式中，所述基于所述参考信息和所述观测信息的互相关构建所述空时二维信号对应的匹配滤波器集合具体包括：
将所述参考信息作为目标参考信息，将观测信息作为目标观测信息；
提取所述目标参考信息和目标观测信息的互相关，并基于所述互相关构建一匹配滤波器；
通过所述匹配滤波器对目标观测信号进行投影，得到更新参考信号；
构造与所述匹配滤波器正交的阻塞矩阵，并基于所述阻塞矩阵和所述目标观测信号确定更新参考信号；
将所述更新参考信号作为目标参考信号，将所述更新参考信号作为目标参考信号；
重新执行提取所述目标参考信息和目标观测信息的互相关，并基于所述互相关构建一匹配滤波器的步骤，直至获取到预设数量个匹配滤波器，以得到对应的匹配滤波器集合。
根据上述技术手段，通过不断提取更新参考信号，实现匹配滤波器的精细化，从而提升信号处理的分辨率和鲁棒性。
在一个实现方式中，所述对所述空时二维信号进行平滑处理以得到空时信号，并基于所述匹配滤波器集合以及所述空时信号，确定所述空时二维信号对应的信源数具体包括：
计算所述空时二维信号的第一协方差矩阵，并基于所述第一协方差矩阵计算空时二维信号前向空时平滑预设次数所得到的时空信息对应的第二协方差矩阵；
对所述第二协方差矩阵进行分块处理，并对分块后的第二协方差矩阵进行特征分解，以得到特征向量矩阵；
基于所述特征向量矩阵和所述匹配滤波器集合对所述第二协方差矩阵进行重构，以得到重构协方差矩阵；
基于重构协方差矩阵通过盖尔圆估计器，确定所述空时二维信号对应的信源数。
根据上述技术手段，通过计算空时信号的协方差矩阵并进行分块特征分解，以及基于盖尔圆估计器确定信源数，对空时信号性质的深入分析，为波达方向估计提供了更精确的信号源数量信息，提高处理信号的信息利用率和准确性。
在一个实现方式中，通过所述通过多普勒信号子空间求解所述空时二维信号，得到速度和波达方向具体包括：
将所述空时二维信号中的空时导向矢量构成为第一空时矩阵和第二空时矩阵，其中，所述第一空时矩阵和所述第二空时矩阵间滑动一个空域导向矢量长度，第二空时矩阵等于第一空时矩阵和旋转算子的乘积；
将所述多普勒信号子空间沿列方向滑动，以构造出第一信号子空间矩阵和第二信号子空间矩阵，其中，第一信号子空间矩阵等于第一空时矩阵和满秩非奇异阵的乘积，第二信号子空间矩阵等于第二空时矩阵和满秩非奇异阵的乘积；
根据第一信号子空间矩阵和第二信号子空间矩阵，通过最小二乘准则求解旋转算子以得到速度，并基于所述速度确定波达方向。
根据上述技术手段，通过构建信号子空间矩阵，运用最小二乘准则精确求解旋转算子，获得目标速度，并在此基础上确定波达方向，利用信号子空间信息，提高波达方向估计的精确度和实际应用效率。
在一个实现方式中，所述基于所述速度确定波达方向具体包括：
对多普勒信号子空间进行重构得到包含旋转算子的角度信号子空间；
通过对所述角度信号子空进行划分构建旋转算子对应的近似矩阵，并通过最小二乘准则求解所述近似矩阵以得到波达方向。
根据上述技术手段，通过构建旋转算子近似矩阵以及最小二乘准则求解波达方向，实现了旋转算子计算的优化，使得对波达方向的估计更加精确，强化波达方向估计在多目标情况下的分辨能力。
本申请实施例第二方面提供了一种波达方向及速度联合估计装置，其中，所述装置包括：
信号获取模块，获取雷达回波信号对应的目标距离门，基于目标距离门确定空时二维信号；
第一信号处理模块，利用多级正交投影方式构建空时二维信号对应的匹配滤波器集合；
第二信号处理模块，对所述空时二维信号进行平滑处理以得到空时信号，并基于所述匹配滤波器集合以及所述空时信号，确定所述空时二维信号对应的信源数；
信号空间构建模块，基于所述匹配滤波器集合和所述信源数确定多普勒信号子空间；
目标计算模块，通过多普勒信号子空间求解所述空时二维信号，得到速度和波达方向。
本申请实施例第三方面提供了一种计算机可读存储介质，所述计算机可读存储介质存储有一个或者多个程序，所述一个或者多个程序可被一个或者多个处理器执行，以实现如上任一所述的波达方向及速度联合估计方法中的步骤。
本申请实施例第四方面提供了一种终端设备，其包括：处理器和存储器；
所述存储器上存储有可被所述处理器执行的计算机可读程序；
所述处理器执行所述计算机可读程序时实现如上任一所述的波达方向及速度联合估计方法中的步骤。
有益效果：与现有技术相比，本申请首先对雷达接收的信号通过快时间维度的快速傅里叶变换（FFT）处理，然后运用恒虚警率检测（CFAR）提取目标距离信息，在保持处理速度的同时提高了目标检测的精确度；在获取目标距离门后，进一步从数据立方体中提取出与目标对应位置的单快拍空时信号，从而减少了数据的维度，简化后续计算；利用空间与时间维度的等效性，在空域-慢时间域上进行联合子空间估计，增加目标检测的自由度，允许本发明超越传统限制，提高角度估计的精度；最后，通过旋转不变技术（ESPRIT）算法联合估计每一个目标的方位角和对应的速度，避免了特征值分解，从而减轻了计算负担，并加快了估计过程。
本申请无需进行繁杂的协方差矩阵特征值分解，简化了波达方向和多普勒频率的联合估计处理流程，减少了时间复杂度，从而提升了处理速度，实现了快速波达方向估计，使得算法更加适用于在对实时性要求高的自动驾驶场景进行工程实现。此外，本申请通过增加目标检测的自由度，提高了波达方向角度和速度估计的精度，为自动驾驶车辆提供了更加可靠的环境感知数据，从而增强了车辆在复杂环境中的导航能力。
附图说明
为了更清楚地说明本申请实施例中的技术方案，下面将对实施例描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本申请的一些实施例，对于本领域普通技术人员而言，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
图1为本申请提供的波达方向及速度联合估计方法的流程图。
图2为本申请提供的波达方向及速度联合估计方法的四目标信号的速度估计效果图。
图3为本申请提供的波达方向及速度联合估计方法的四目标信号的角度估计效果图。
图4为本申请提供的波达方向及速度联合估计方法的测试场景下点云成像效果图。
图5为本申请提供的波达方向及速度联合估计装置的结构原理图。
图6为本申请提供的终端设备的结构原理图。
具体实施方式
本申请提供一种波达方向及速度联合估计方法及相关装置，为使本申请的目的、技术方案及效果更加清楚、明确，以下参照附图并举实施例对本申请进一步详细说明。应当理解，此处所描述的具体实施例仅用以解释本申请，并不用于限定本申请。
本技术领域技术人员可以理解，除非特意声明，这里使用的单数形式“一”、“一个”、“所述”和“该”也可包括复数形式。应该进一步理解的是，本申请的说明书中使用的措辞“包括”是指存在所述特征、整数、步骤、操作、元件和/或组件，但是并不排除存在或添加一个或多个其他特征、整数、步骤、操作、元件、组件和/或它们的组。应该理解，当我们称元件被“连接”或“耦接”到另一元件时，它可以直接连接或耦接到其他元件，或者也可以存在中间元件。此外，这里使用的“连接”或“耦接”可以包括无线连接或无线耦接。这里使用的措辞“和/或”包括一个或更多个相关联的列出项的全部或任一单元和全部组合。
本技术领域技术人员可以理解，除非另外定义，这里使用的所有术语(包括技术术语和科学术语)，具有与本申请所属领域中的普通技术人员的一般理解相同的意义。还应该理解的是，诸如通用字典中定义的那些术语，应该被理解为具有与现有技术的上下文中的意义一致的意义，并且除非像这里一样被特定定义，否则不会用理想化或过于正式的含义来解释。
应理解，本实施例中各步骤的序号和大小并不意味着执行顺序的先后，各过程的执行顺序以其功能和内在逻辑确定，而不应对本申请实施例的实施过程构成任何限定。
经过研究发现，在自动驾驶技术的不断演进中，精确快速地获取环境信息是确保操作安全和提升导航能力的核心环节。雷达传感器（特别是毫米波雷达），因其全天候全天时工作的优点，在自动驾驶系统中占有不可替代的地位，它能够在各种复杂气候条件下提供稳定可靠的目标探测功能。毫米波雷达通过发射和接收电磁波脉冲，从回波中提取目标的距离、速度、方位角度及俯仰角度等信息，这对于识别和定位周边物体至关重要。
然而，尽管毫米波雷达系统相对于光学传感器和激光雷达在能见度不良环境下显示出较强的适应性，其角度分辨能力却受限于天线的孔径大小，这导致了在探测同一距离门内不同目标时容易出现分辨率不足的问题。传统的角度估计方法，如快速傅里叶变换（FFT），在硬件简单和效率方面虽具有优势，但其性能受限于系统波束宽度，无法提供足够的角度分辨力以满足复杂场景分析的需求。尽管基于子空间分解的MUSIC算法理论上能够突破这些限制，但其对于协方差矩阵的特征值分解过程计算复杂且耗时，这在实际运用中降低了其工程可行性，特别是与自动驾驶中的实时性需求相冲突。
为了解决上述问题，在本申请实施例中，获取雷达回波信号对应的目标距离门，基于目标距离门确定空时二维信号；利用多级正交投影方式构建空时二维信号对应的匹配滤波器集合；对所述空时二维信号进行平滑处理以得到空时信号，并基于所述匹配滤波器集合以及所述空时信号，确定所述空时二维信号对应的信源数；基于所述匹配滤波器集合和所述信源数确定所述多普勒信号子空间；通过多普勒信号子空间求解所述空时二维信号，得到速度和波达方向。本申请通过获取并处理雷达回波信号来分辨在相同距离门下的不同目标，通过复层次信号处理流程，提取目标的速度信息，增强目标波达方向的分辨率，满足自动驾驶中对于复杂环境的高精度感知需求。
下面结合附图，通过对实施例的描述，对申请内容作进一步说明。
本实施例提供了一种波达方向及速度联合估计方法，如附图1所示，所述波达方向及速度联合估计方法包括：
S10、获取雷达回波信号对应的目标距离门，基于目标距离门确定空时二维信号。
具体地，所述目标距离门指的是雷达系统在距离维度上对回波信号进行处理后识别出来的特定的距离区间，这个区间内包含对应目标的反射信号。所述目标距离门类似于在雷达图像上用门（窗口）把目标的信号框起来的过程。通过设置距离门，可以提取雷达确认的目标所在的距离单元的信息，而忽略其他距离单元的信息。其中，在一些场景中，雷达可以同时扫描到多个目标，这时，根据回波信号的不同，可以确定多个对应的目标距离门。每一个目标距离门都对应着一个特定的距离区间，距离区间内包含了对应目标的反射信号。换句话说，多目标环境下，可形成多个目标距离门，每个距离门都对应一个目标的空时二维信号。
获取目标距离门的过程中，雷达接收到的信号是时域信号，首先通过傅里叶变换将其转换为频域信号，该过程可以将时间序列信号转换为频率成分和相应幅度的表示。然后，通过恒虚警率检测（CA-CFAR），在所得的频域信号中，根据预设的门限判断哪些频率成分代表了真正的目标信号，哪些仅仅是噪声或杂波。这需要设置保护单元和参考单元，计算这些参考单元的平均值作为功率估计，据此确定门限。最终，通过与门限比较，如果某个距离单元的信号幅度超过门限值，所述单元被标记为包含目标信号的距离门，从而在距离域中筛选出可能包含目标的范围。
基于目标距离门确定空时二维信号的过程中，首先提取所获得目标距离门内的回波信号，所述回波信号反映了处于对应距离门内的目标的存在和特性，可以反映出目标的距离，但尚未包含关于目标方位和速度的信息。将提取出的所有目标回波信号，按照脉冲序列和阵列天线的排列，整合成一个空时二维信号。这里，“空时二维”指的是结果数据包括空间（阵列所在的维度）上的变化，也包括时间（脉冲序列）上的变化。空时二维信号指的是经由雷达系统各个时间和空间采样点接收到的回波信号所形成的信号矩阵，该信号矩阵可以视为二维数据阵列，其中一维是时间维度（多个脉冲），另一维是空间维度（阵列天线元素）。所述信号矩阵通常包含了反映目标信息和环境噪声的丰富数据，所述空时二维信号被用作后续步骤的输入，为波达方向估计提供信息。
借助于上述步骤，从雷达接收的大量回波信号中筛选出与目标相关的信息，降低了后续处理的数据量和复杂度。同时确定的目标距离门可以提高目标检测的准确性和信噪比。此外，构造的空时二维信号为后续的信源数估计、匹配滤波器的构建以及波达方向和速度估计准备了必需的数据结构，实现了数据的有效利用，增强了波达方向估计的准确度和鲁棒性。
在本申请的具体实现中，所述获取雷达回波信号对应的目标距离门，基于目标距离门确定空时二维信号具体包括：
S11、对雷达回波信号进行傅里叶处理和恒虚警率检测，得到目标距离门。
S12、提取所述目标距离门对应的回波脉冲，并将提取到的所有回波脉冲排列为空时二维信号。
具体地，在步骤S11中，对雷达回波信号进行傅里叶处理，即执行快速傅里叶变换（FFT），将时域中的回波信号转换到频域中，得到的结果表征了目标对象在距离维度上的分布信息。随后，实施恒虚警率（Constant False Alarm Rate, CFAR）检测。在CFAR检测中，根据雷达系统环境噪声的统计特性设置一个动态门限，该门限用于区分目标返回信号和噪声，以保持虚警率恒定。其中，CFAR算法通过在设定的保护单元和参考单元中估计杂波水平，再以此计算出检测门限，以此来确定目标是否存在，进而得到目标存在的径向距离信息，以得到目标距离门。
在本申请的一个具体实现中，对雷达回波信号进行傅里叶处理和恒虚警率检测，得到目标距离门的实现过程中雷达回波信号可以采用使用均匀线阵模型的毫米波雷达天线获取雷达回波信号，其中天线数目为，天线之间的间隔为(为雷达信号的波长)，个相干脉冲。通过对雷达回波信号的数据立方的距离维度进行相干累积、傅里叶处理，然后采用恒虚警率检测获得目标的距离信息和距离门。基于此，目标距离门的获取过程可以包括：
首先，使用快速傅里叶变换对时域进行相干累积。
具体地，对第个脉冲中第个天线的（为2的整数次幂）个快拍信号进行点快速傅里叶变换进入频域（距离域）：
，，
其中，在上的进行单位均值恒虚警率检测（Cell Averaging CFAR，CA-CFAR）。
对第个待检测单元设置一定数量的保护单元，对两端个参考单元求平均得到杂波功率：
；
其中。设置待测单元的门限阈值为：
；
其中为门限因子，为设置的虚警率：
；
遍历上所有的待检测单元通过比对检测单元与门限阈值的大小进行目标检测。根据奈曼-皮尔逊准则，当检测单元的幅值时，表示目标存在于该检测单元并将处于的幅值置为一；反之，则表示该检测单元不存在目标，并将其幅值置零。
在步骤S12中，在获取到目标距离门后，提取目标距离门对应的回波脉冲。其中，对经过FFT处理且已通过CFAR门限的回波信号进行筛选，从多个相邻的距离单元中提取出含有目标反射信号的距离单元。每一个距离单元包含来自该目标距离门的回波信息，而每个雷达脉冲都为该距离单元提供一个时间样本。通过将所有脉冲对应的距离单元的回波信息提取出来并在时间维度上排列，得到空时二维信号。所述空时二维信号既含有空间域的信息（因为不同的天线在空间上是分开的），也含有时间域的信息（因为脉冲是随时间序列发射的）。换句话说，空时二维信号融合了空间和时间上的信息，为后续步骤提供更全面的数据。
在本申请的一个具体实现中，使用如下方式对回波数据建模：
从一维距离像中提取目标距离门的脉冲阵元截面，其中目标距离门中第k个脉冲N天线阵列回波数据建模如下：
，，
其中，为高斯白噪声；为空域导向矢量，包含了个目标回波：
，；
；
此外，是角度为的目标空域导向矢量，而为阵元间隔。
将提取到的回波脉冲排列为空时二维列向量的构造过程如下：
将全部脉冲，排列成的空时二维列向量：
，
其中，表示Kronecker积，为时域导向矢量：
，；
。
通过执行上述S11和S12步骤，将雷达回波信号转化为包含空间和时间信息的二维信号，提高了信号的可利用性，并降低了信号处理的复杂度，加速了后续信息提取的过程。
此外，利用傅里叶处理和恒虚警率检测协同作用，使得在获取距离门的基础上能够有效地压制杂波和噪声，提取出具有潜在目标信息的空时二维信号，不仅可以保持较低的虚警率，而且能适应不同的信号环境和背景噪声水平。通过排列为空时二维信号，使得后续的信号处理能够充分利用空间和时间维度上的信息，增加了复杂的动态目标分析的准确性，如多目标、快速移动目标的检测与参数估计。
S20、利用多级正交投影方式构建空时二维信号对应的匹配滤波器集合。
具体地，多级正交投影方式是指一种计算机视觉和信号处理中的技术，用于从所观测到的信号中提取特定的特征或分量。在本申请的具体实现中，多级正交投影方式用于构造能够匹配特定信号分量的滤波器集合，进而从空时二维信号中识别和分离出单独信源的具体特征。
其中匹配滤波器集合可以包含多个滤波器或多个滤波器，每一个滤波器都被设计为对特定特性的信号进行响应。比如说，在应对多目标场景时，可以为每个预期的单独信源设立一个对应的滤波器。因此，不同数量的目标会导致产生相应数量的匹配滤波器。
进一步，通过这些匹配滤波器集合，可以对信号进行多种阻滞操作。如根据任务需要，可以设置不同的目标距离门，进而滤除距离目标过远或过近的信号。
再者，这种方法可以应用于处理动态变化的信号。例如，当信号源呈现动态移动时，匹配滤波器可以针对多普勒频率移动进行设定。
这种多级正交投影方式的方法，提高了信号处理的效率和准确度，优化了信号识别和处理的性能。无论是处理静态信号还是动态信号，无论目标数量多少，无论信号特性如何，都可以通过设定适当的匹配滤波器集合来完成任务，展现出非常高的适用性和灵活性。
在步骤S20中，所述构建过程可以包括：
首先从多个雷达阵元接收到的空时二维信号中，选取一个阵元的接收信息作为参考信息；将上述参考信息与其他阵元接收到的信号进行互相关处理，得到各个阵元信号间的相关度量；基于互相关的结果，构建第一个匹配滤波器，该滤波器对应于空时信号中的主要信号分量；利用第一个匹配滤波器对信号进行投影操作，得到更新后的参考信号，该信号将用于下一个匹配滤波器的构建；为了去掉已经提取出的信号分量，构建一个与前述匹配滤波器正交的阻塞矩阵；在第一个滤波器基础上，迭代使用更新的信号与阻塞矩阵，构建一系列的正交匹配滤波器，直至获得足够数量的滤波器以涵盖所需的信号特征。
通过步骤S20提出的利用多级正交投影方式构建空时二维信号对应的匹配滤波器集合的方法，不仅可以提高雷达回波信号处理的准确性和鲁棒性，而且能适应动态变化的信号环境，并且在保证高精度波达方向估计的同时，优化整体的计算效率。
在本申请的具体实现中，所述利用多级正交投影方式构建空时二维信号对应的匹配滤波器集合具体包括：
S21、将空时二维信号中的位于最后的阵元的接收信息作为参考信息，将所述空时二维信号作为观测信息；
S22、基于所述参考信息和所述观测信息的互相关构建所述空时二维信号对应的匹配滤波器集合。
具体地，在步骤S21中，空时二维信号的最后一个阵元接收到的信号被选定为“参考信息”，而整个空时二维信号则作为“观测信息”。在信号处理中，通常，所选取的参考信息是希望通过处理手段，与其他信号进行比较或联合分析的信号。
而在步骤S22中，通过计算参考信息与观测信息之间的互相关性，构建针对空时二维信号的匹配滤波器集合。互相关是一种表征两个信号之间相似性的度量，通过评估两个信号在多个时间偏移量下的相似度来完成。在本申请的实现方式中，互相关不仅确定匹配滤波器如何根据参考信息调整自己的参数，以便最大化响应目标信号，同时也用于后续的多级正交投影工作。
通过利用互相关结果，可以构建第一个匹配滤波器，将观测信息投影到此滤波器上，获得对策的信号组分。进而，可以构建一个与匹配滤波器集合中的第一个滤波器正交的阻塞矩阵，即阻断该信号组分，从而使得新的参考信号仅包含之前未被第一个滤波器捕获的组分。这个过程可以递归地进行多次，每次都通过计算新的互相关来确定新的匹配滤波器，以此构造出正交的匹配滤波器集合。
在本申请的具体实现中，所述基于所述参考信息和所述观测信息的互相关构建所述空时二维信号对应的匹配滤波器集合具体包括：
S221、将所述参考信息作为目标参考信息，将观测信息作为目标观测信息；
S222、提取所述目标参考信息和目标观测信息的互相关，并基于所述互相关构建一匹配滤波器；
S223、通过所述匹配滤波器对目标观测信号进行投影，得到更新参考信号；
S224、构造与所述匹配滤波器正交的阻塞矩阵，并基于所述阻塞矩阵和所述目标观测信号确定更新参考信号；
S225、将所述更新参考信号作为目标参考信号，将所述更新参考信号作为目标参考信号；
S226、重新执行提取所述目标参考信息和目标观测信息的互相关，并基于所述互相关构建一匹配滤波器的步骤，直至获取到预设数量个匹配滤波器，以得到对应的匹配滤波器集合。
具体地，首先在S221步骤中选择最后一个阵元的数据作为目标参考信息，并以此为基准与整个观测信息（即空时二维信号）进行比较，提取出信号特征；其次在S222步骤中通过计算目标参考信息与目标观测信息之间的互相关，构建出一个初步的匹配滤波器，用于揭示了两者之间的相似性，高互相关指明部分观测信息与参考信息极为相似；再次在S223步骤中应用构建的匹配滤波器对目标观测信息进行滤波，即进行投影操作，从而在多信号环境中单独提取出与参考信息最匹配的信号成分，并生成更新参考信号；随后在S224步骤中构建一个与现有匹配滤波器正交的阻塞矩阵，阻塞矩阵能够抑制或“阻塞”已通过现有滤波器提取出来的信号成分，使得观测信息中剩余的成分被保留下来，得到更新参考信号；进而在S225步骤中将通过阻塞矩阵和目标观测信息计算得到的更新参考信号，作为新的目标参考信息，为下一轮的匹配滤波器构建提供基准；最后在S226步骤中重复执行S221至S225直到获取预设数量的匹配滤波器，通过递归式的处理方法允许系统逐步提取信号中的多个目标成分，并构建一套正交的匹配滤波器集合。
通过上述步骤S221至S225，对目标信号的一致性和区分度进行准确的建模，有助于同一时间从多目标环境中分辨和识别信号；通过正交的匹配滤波器集合从众多信号成分中分离出目标信号，从而减少对噪声和干扰的敏感度。
在本申请的一个具体实现方案中，使用多级正交投影的方法构建空时二维信号对应的匹配滤波器集合的过程如下：
首先，最后一个阵元的接收信号作为参考信号：
，
参考信号中包含了目标信号、杂波以及噪声。
其次，全部阵元接收信号作为观测信号，且参考信号与观测信号的互相关：
，
其中，代表期望操作，代表取共轭操作。
基于噪声信号的随机性，噪声子空间之间为非相干叠加，信号子空间为相干叠加的特征，定义第一个匹配滤波器：
，
其中表示范数，通过可以投影得到新的参考信号：
，
式中的代表共轭转置操作。
再次，构造与正交的阻塞矩阵，过滤相关信号成分：
；
将观测信号与阻塞矩阵相乘后输出为新的观测信号：
；
通过类似的方法得到第二个滤波器：
；
通过滤波器可以进一步得到新的参考信号：
；
经历次上述的迭代过程，得到个滤波器集合：
（1）
S30、对所述空时二维信号进行平滑处理以得到空时信号，并基于所述匹配滤波器集合以及所述空时信号，确定所述空时二维信号对应的信源数。
具体地，平滑处理通常涉及将原始空时信号按时间或空间分解为多个相重叠或不相重叠的子阵列，并在这些子阵列上做处理，从而提高信号处理的统计稳定性，减少因有限采样而造成的信号估计误差。平滑处理的目的是增加信号处理中独立采样点的数量，以增强目标信号并提升多目标分辨能力。
此外，信源数指的是雷达通过回波信号实际上检测到的目标数量。
在步骤S30中，首先对空时二维信号进行平滑处理，其中，平滑处理过程可以包括但不限于如下操作：对原始信号分块、在每块内部平均、应用信号窗函数等，从而减少随机噪声的影响，同时提高信号的相关性，为进一步的信号处理创造有利条件。
其次，经过平滑后的空时信号，具有提升后的信号特性，如更高的信噪比和更稳定的统计特性。随后，通过已经构建的匹配滤波器集合，对平滑后的空时信号进行处理，进一步消除非目标信号的影响，并且突出目标的特征。最后，基于处理过的空时信号，可以通过分析信号协方差矩阵、应用模型选择准则或者利用信号的统计特性等多种方法确定空时二维信号对应的信源数。
上述步骤通过平滑处理，有效地减少噪声和干扰，增加信号的连续性和一致性，从而为后续处理提供更准确的基础；此外，在平滑处理和匹配滤波的基础上，信源数的估计更接近真实情况，避免了过度估计或者低估信源数的风险。
在本申请的具体实现中，所述对所述空时二维信号进行平滑处理以得到空时信号，并基于所述匹配滤波器集合以及所述空时信号，确定所述空时二维信号对应的信源数具体包括：
S31、计算所述空时二维信号的第一协方差矩阵，并基于所述第一协方差矩阵计算空时二维信号前向空时平滑预设次数所得到的时空信息对应的第二协方差矩阵；
S32、对所述第二协方差矩阵进行分块处理，并对分块后的第二协方差矩阵进行特征分解，以得到特征向量矩阵；
S33、基于所述特征向量矩阵和所述匹配滤波器集合对所述第二协方差矩阵进行重构，以得到重构协方差矩阵；
S34、基于重构协方差矩阵通过盖尔圆估计器，确定所述空时二维信号对应的信源数。
具体地，所述第一协方差矩阵是从空时二维信号计算得到的统计矩阵，提供了信号的时间维度和空间维度内各元素之间的相互关系（即互协方差）的信息。
所述第二协方差矩阵是基于第一协方差矩阵经过一定次数的前向空时平滑处理后得到的矩阵。前向空时平滑是指将信号沿着时间或空间维度进行分割和重组的操作，通常用于增强信号特征，提高后续处理的稳定性和可靠性。
所述特征向量矩阵是指通过对第二协方差矩阵进行特征分解得到的矩阵，它的列向量（特征向量）表征了信号的主要模态，可以用来识别信号中的统计规律或特征。
所述重构协方差矩阵是利用特征向量矩阵和匹配滤波器集合，经过变换重构得到的协方差矩阵，所述协方差矩阵表示信号的统计特性，并用于信源数估计。
所述盖尔圆估计器是一种用于信源数估计的图形化工具，通过分析信号特性或分布来确定信号中包含的独立信源个数。
在S31步骤中，首先基于空时二维信号计算出第一协方差矩阵。对协方差矩阵进行多次前向空时平滑处理，以提高信号特征的可识别度。每次平滑处理可以包括信号的部分加权与再组合，从而减小随机噪声的效应。
在S32步骤中特征分解用于揭示信号内部的结构特性，即哪些是信号的主要部分，哪些可能是噪声或较不重要的部分。
在S33步骤中基于得到的特征向量矩阵，再结合匹配滤波器集合相关信息，对第二协方差矩阵进行重构，从而扣除噪声的影响，并更准确地表示信号的时间与空间特性。
在S34步骤中，借助盖尔圆估计器来分析重构后的协方差矩阵，从而确定信号中的独立信源数量。盖尔圆估计器通过圆的个数、大小和相对位置，从而确定空时二维信号对应的信源数。
上述步骤首先通过前向空时平滑处理，在不丢失重要信号信息的前提下，有效抑制噪声，从而提高信噪比；其次，特征向量矩阵和重构协方差矩阵使得信号中各个目标的特征更加突出，便于后续步骤的准确识别和处理；再次，通过盖尔圆估计器对信号特征的图形化分析，提供了一种直观便于理解的方法来确定信号中存在的独立信源数量。
在本申请的一个具体实现方式中，所述空时平滑处理空时信号确定信源数的具体过程包括：
首先，对于N个阵元、K个脉冲的单快拍空时二维信号，先计算协方差矩阵：
；
再进行次前向空时平滑得到第二协方差矩阵：
，
其中代表平滑后的第二协方差矩阵大小。
其次，对第二协方差矩阵进行矩阵分块处理：
，
其中代表最后一个元素，表示最后一列除以外的列向量，对进行特征值分解得到特征向量矩阵：
，
其中，代表特征值构成主对角元素的对角矩阵，其中特征值排列由大到小；表示由特征向量构成的矩阵；构造一个酉矩阵并对协方差矩阵进行变换得到重构协方差矩阵：
；
；
对于重构后的协方差矩阵，第个盖尔圆定义如下：
，，
其中重构协方差矩阵中第i行第i列的元素为第个盖尔圆的圆心，为第个盖尔圆的半径。通过以下判别式估计出接收信号所含信源数：
，，
其中，为权重系数，一般设置为，其中为常数。盖尔圆按照半径由大小排列，当时，该接收信号的信源数为。
S40、基于所述匹配滤波器集合和所述信源数确定多普勒信号子空间。
通常来说，多普勒效应在雷达信号处理中表现为目标运动引起的频率变化，其信息主要集中在接收到的雷达信号的频率差异中。多普勒信号子空间是指通过对信号频谱的分析获得的包含目标多普勒频率信息的部分信号空间。多普勒信号子空间内的维度与信号中独立目标的数量（即信源数）有直接关联。
在前述的步骤中，根据空时二维信号已经构建匹配滤波器集合，并估计出信号对应的信源数，其次通过匹配滤波器集合捕获的信号特征来构造多普勒信号子空间。由于匹配滤波器设计以响应特定的信号特征（如来自特定方向的回波信号），所述匹配滤波器集合直接揭示了存在于接收信号中的不同多普勒成分。在获得匹配滤波器集合后，根据信源数，选择相应数量的滤波器构成子空间，举例来说，如果通过前面步骤得出的信源数为N，则从匹配滤波器集合中选取与N相应的滤波器来定义此子空间。
在步骤S40中，多普勒信号子空间的建立有助于突出目标运动特性，为波达方向估计和速度测量提供更准确的信息；其次，通过选择与信源数相符合的匹配滤波器，减少了信号处理的复杂度，避免了处理无关的信号分量，从而提高了算法的运算效率。
S50、通过多普勒信号子空间求解所述空时二维信号，得到速度和波达方向。
具体地，多普勒信号子空间是指一个封闭的信号处理区域（或子集），在该区域内的信号包含了反映目标运动速度的多普勒频率信息。该子空间由前述步骤构建匹配滤波器集合和信源数决定，能够有效地分离并定量化目标运动的速度信息。
经过前面步骤的信号预处理和多普勒信号子空间构建之后，当前步骤S50通过解析在多普勒信号子空间中表现出来的频率变化进一步萃取目标的速度信息；解析过程通常涉及到对子空间中信号的频谱分析，例如使用快速傅里叶变换（FFT）或者其他频域转换算法，以捕捉频域上具体的频率分量，从而推导出目标的速度信息；同时，根据频率分量与波达方向之间的关系，求解波达方向。
其中，波达方向是由多普勒频率分量和波长之间的几何关系确定的，所述几何关系取决于雷达回波信号对应的目标相对于雷达天线阵列的方向，并且每个方向的波达方向与一个特定的多普勒频率相关联。
通过利用构建的多普勒信号子空间和对空时二维信号的深入分析，得到精确的速度信息，从而求解波达方向。
在本申请的具体实现中，所述通过多普勒信号子空间求解所述空时二维信号，得到速度和波达方向具体包括：
S51、将所述空时二维信号中的空时导向矢量构成为第一空时矩阵和第二空时矩阵，其中，所述第一空时矩阵和所述第二空时矩阵间滑动一个空域导向矢量长度，第二空时矩阵等于第一空时矩阵和旋转算子的乘积；
S52、将所述多普勒信号子空间沿列方向滑动，以构造出第一信号子空间矩阵和第二信号子空间矩阵，其中，第一信号子空间矩阵等于第一空时矩阵和满秩非奇异阵的乘积，第二信号子空间矩阵等于第二空时矩阵和满秩非奇异阵的乘积；
S53、根据第一信号子空间矩阵和第二信号子空间矩阵，通过最小二乘准则求解旋转算子以得到速度，并基于所述速度确定波达方向。
具体地，空时导向矢量是表示信号在特定时间和空间点上的方向性质的向量，用于描述信号在空间和时间上的传播；第一空时矩阵和第二空时矩阵表示在特定时间窗口内的多个空时导向矢量的集合，每个矩阵包含了对应于不同时间点上雷达天线阵列的响应，第一、第二空时矩阵用于时空域上的信号表述，反映了信号随着时间演变的特性；旋转算子是一个描述多普勒频移导致的信号相位变化的表述，用于估计目标速度所导致的信号相位差异；最小二乘用于求解由多个观测数据构成的方程组，目的是找到可使误差平方和最小的参数估计值。
在步骤S51中，将空时二维信号中的空时导向矢量构造成两个不同的矩阵，第一空时矩阵反映了信号的时空分布，第二空时矩阵通过与第一空时矩阵相差一个空间单元（空域导向矢量长度）以及多普勒频移效应（旋转算子）的乘积构成。
在步骤S52中，多普勒信号子空间通过对第一空时矩阵和第二空时矩阵的操作，所述操作可以是列方向滑动，从而构造出两个信号子空间矩阵。这两个信号子空间矩阵包含了信号的多普勒频移信息，对应于信号在时间上的变化。具体地，在本申请中，第一信号子空间矩阵等于第一空时矩阵和满秩非奇异阵的乘积，第二信号子空间矩阵等于第二空时矩阵和满秩非奇异阵的乘积。
在S53步骤中，利用最小二乘准则对两个信号子空间矩阵进行分析，求解得到旋转算子。旋转算子反映了信号因速度变化而发生的频率偏移，可用于精确估计目标的速度。基于所得速度进一步求解波达方向，确立目标的具体方向。
上述步骤S51至S53，通过准确求解旋转算子，使得目标速度的估计达到高精度，从而实现更精准的动态目标跟踪，此外通过与速度估计相关的信息，可对波达方向进行精确求解。
在本申请的具体实现中，所述基于所述速度确定波达方向具体包括：
S531、对多普勒信号子空间进行重构得到包含旋转算子的角度信号子空间；
S532、通过对所述角度信号子空间进行划分构建旋转算子对应的近似矩阵，并通过最小二乘准则求解所述近似矩阵以得到波达方向。
具体地，角度信号子空间用于从空时信号中提取角度信息，确定目标速度后，进一步对多普勒信号子空间进行操作，将其转换为一个新的子空间。
在S531步骤中，首先，对多普勒信号子空间执行数学上的重构操作，生成角度信号子空间。通过这个重构过程，将速度信息与角度信息分离，角度信号子空间便主要包含了波达方向的信息。
在S532步骤中，进一步分析和处理角度信号子空间中的信息，首先进行结构上的划分，构建与旋转算子相关联的近似矩阵。随后，应用数学上的最小二乘准则对这个近似矩阵进行求解，最终获得波达方向的估计值。
通过对多普勒信号子空间的重构和近似矩阵的构建，准确地定位目标的波达方向；此外，使用近似矩阵和最小二乘法简化了波达方向求解的复杂性，使得波达方向的获取变得更为直接和高效，并且采用近似矩阵和最小二乘法不仅简化了处理流程，还减少了计算成本和时间，提升了算法的实际应用价值。
在本申请的一个具体实现方式中，旋转算子的求解和速度估计步骤和波达方向的确定步骤如下：
利用旋转不变子空间的思想，求解空时导向矢量中包含的旋转算子，进而得到目标角度与速度估计。
首先，将空时导向矢量构造乘成第一空时矩阵和第二空时矩阵，
其中，两个矩阵之间滑动了一个空域导向矢量长度，脉冲间的旋转算子对应目标的多普勒频率：
（2）
（3）
由上式可推出，旋转算子亦可通过滑动信号子空间获得。
其次，根据阵列的平移不变性，将多普勒信号子空间沿列方向滑动一个阵列长度构造出两个子矩阵：
；
；
其中，为第一信号子空间矩阵，为第二信号子空间矩阵，为唯一的维满秩非奇异阵。
由式1和式2可得：
，
其中，与相似，具有相同的特征值，根据最小二乘准则：
；
对进行特征值分解：
；
对特征值对角阵进行反演即可得到目标的速度：
，，
其中表示取复数的相位角，表示取主对角线的第个元素。
同理，交换空时导向矢量的Kronecker积顺序构造时空导向矢量，并将其滑动一个时域导向矢量长度生成与旋转算子相关的两个子矩阵和：
；
；
重构后的信号表示为：
，
其中代表包含不同顺序Kronecker积的导向矢量位置变化信息的置换矩阵。
对信号子空间进行重构得到包含旋转算子的角度信号子空间：
；
从划分两个子矩阵和：
；
；
构造旋转算子造的近似矩阵为：
；
由式3可以知，的特征向量为，其中为置换矩阵。
通过共享的特征向量固定置换矩阵，使得与特征值及其特征向量的排序一致，进而推导出的特征值：
；
对特征值对角阵进行反演即可得到目标的角度：
，；
经过上述步骤，可以同时估算出目标的速度及对应的角度，且通过共享特征空间，节省了参数匹配的开销。
在一个具体地实现方式中，使用86均匀线阵对四目标信号进行参数估计，每个脉冲采样1024个快拍，累积64个脉冲，四个信源目标的速度分别为 [-5.3,-5,0,5]，角度分别为[-5,0,3,3.5]。本申请所提联合估计方法的效果如附图2、附图3所示，四个信源目标均能被准确估计出来。
综上所述，本实施例提供了一种波达方向及速度联合估计方法，该方法包括：获取雷达回波信号对应的目标距离门，基于目标距离门确定空时二维信号；利用多级正交投影方式构建空时二维信号对应的匹配滤波器集合；对所述空时二维信号进行平滑处理以得到空时信号，并基于所述匹配滤波器集合以及所述空时信号，确定所述空时二维信号对应的信源数；基于所述匹配滤波器集合和所述信源数确定所述多普勒信号子空间；通过多普勒信号子空间求解所述空时二维信号，得到速度和波达方向。本申请通过利用多级正交投影和盖尔圆估计器提升信号源数的精确估计，以及应用最小二乘法确保速度和波达方向的精确测量，整体提高了目标位置信息的准确性；通过平滑处理手段减小了噪声对信号处理的干扰，匹配滤波器集合进一步增强了目标信号相较于噪声信号的分辨能力；通过旋转算子和最小二乘法简化波达方向估计的数学处理，减少了算法的运行时间，提升了处理效率和实时性；通过造的多普勒信号子空间使得算法对于动态的目标（特别是高速运动目标）有更好的响应能力。
在本申请的一个具体实现方案中，采用四级联AWR2243雷达，可视为86阵元的均匀线阵，相邻阵元间距为雷达信号波长的一半。参数设计为：调频斜率设置为30MHz/us，相干脉冲数为64，脉冲重复频率为1.5MHz；ADC采样精度为16位、采样速率为8000ksps，每个脉冲采样1024个快拍。
以停车场场景为测试场景，以测试场景中的行人为动目标，对测试场景进行雷达回波采集，首先对接收到的信号回波使用快速傅里叶变换进行距离维压缩；接着使用步长为3、虚警率为0.1%的恒虚警率检测进行目标检测，提取维度大小为86*64的空时二维信号；接着对空时二维信号进行32个快拍的空时平滑处理，得到平滑后含有32个快拍的信号，维度为55*33*32的；将其第一第二维度合并，降维成维度为1815*32的空时信号，对使用多级正交投影计算特征空间，把得到的特征空间代入盖尔圆估计器估计信源数；然后使用二维旋转不变子空间算法进行波达方向、速度联合估计；最后通过距离、波达方向和速度生成点云图。
点云成像结果如图4所示，“x”点云代表运动物体，实心圆点云代表静止物体。从测试结果可以看出，本申请借助于波达方向及速度联合估计方法所生成的点云成像质量达到预期。
本实施例提供了一种波达方向及速度联合估计装置，如图5所示，所述装置包括：
信号获取模块100，获取雷达回波信号对应的目标距离门，基于目标距离门确定空时二维信号；
第一信号处理模块200，利用多级正交投影方式构建空时二维信号对应的匹配滤波器集合；
第二信号处理模块300，对所述空时二维信号进行平滑处理以得到空时信号，并基于所述匹配滤波器集合以及所述空时信号，确定所述空时二维信号对应的信源数；
信号空间构建模块400，基于所述匹配滤波器集合和所述信源数确定多普勒信号子空间；
目标计算模块500，通过多普勒信号子空间求解所述空时二维信号，得到速度和波达方向。
基于上述波达方向及速度联合估计方法，本实施例提供了一种计算机可读存储介质，所述计算机可读存储介质存储有一个或者多个程序，所述一个或者多个程序可被一个或者多个处理器执行，以实现如上述实施例所述的波达方向及速度联合估计方法中的步骤。
基于上述波达方向及速度联合估计方法，本申请还提供了一种终端设备，如图6所示，其包括至少一个处理器（processor）20；显示屏21；以及存储器（memory）22，还可以包括通信接口（Communications Interface）23和总线24。其中，处理器20、显示屏21、存储器22和通信接口23可以通过总线24完成相互间的通信。显示屏21设置为显示初始设置模式中预设的用户引导界面。通信接口23可以传输信息。处理器20可以调用存储器22中的逻辑指令，以执行上述实施例中的方法。
此外，上述的存储器22中的逻辑指令可以通过软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在一个计算机可读取存储介质中。
存储器22作为一种计算机可读存储介质，可设置为存储软件程序、计算机可执行程序，如本公开实施例中的方法对应的程序指令或模块。处理器20通过运行存储在存储器22中的软件程序、指令或模块，从而执行功能应用以及数据处理，即实现上述实施例中的方法。
存储器22可包括存储程序区和存储数据区，其中，存储程序区可存储操作系统、至少一个功能所需的应用程序；存储数据区可存储根据终端设备的使用所创建的数据等。此外，存储器22可以包括高速随机存取存储器，还可以包括非易失性存储器。例如，U盘、移动硬盘、只读存储器(Read-Only Memory，ROM)、随机存取存储器(Random Access Memory，RAM)、磁碟或者光盘等多种可以存储程序代码的介质，也可以是暂态存储介质。
此外，上述存储介质以及终端设备中的多条指令处理器加载并执行的具体过程在上述方法中已经详细说明，在这里就不再一一陈述。
最后应说明的是：以上实施例仅用以说明本申请的技术方案，而非对其限制；尽管参照前述实施例对本申请进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本申请各实施例技术方案的精神和范围。
Claims (7)
1.一种波达方向及速度联合估计方法，其特征在于，所述的波达方向及速度联合估计方法具体包括：
获取雷达回波信号对应的目标距离门，基于目标距离门确定空时二维信号；
利用多级正交投影方式构建空时二维信号对应的匹配滤波器集合；
对所述空时二维信号进行平滑处理以得到空时信号，并基于所述匹配滤波器集合以及所述空时信号，确定所述空时二维信号对应的信源数；
基于所述匹配滤波器集合和所述信源数确定多普勒信号子空间；
通过多普勒信号子空间求解所述空时二维信号，得到速度和波达方向；
其中，所述空时二维信号对应的信源数的具体过程包括：
对于N个阵元、K个脉冲的单快拍空时二维信号，先计算协方差矩阵：，再进行次前向空时平滑得到第二协方差矩阵：，其中代表平滑后的第二协方差矩阵大小；
对第二协方差矩阵进行矩阵分块处理，得到，其中，代表最后一个元素，表示最后一列除以外的列向量，对进行特征值分解得到特征向量矩阵，，其中，代表特征值构成主对角元素的对角矩阵，其中特征值排列由大到小；表示由特征向量构成的矩阵；构造一个酉矩阵并对协方差矩阵进行变换得到重构协方差矩阵，，；
对于重构后的协方差矩阵，第个盖尔圆定义如下：，其中，重构协方差矩阵中第行第列的元素为第个盖尔圆的圆心，为第个盖尔圆的半径；
通过以下判别式估计出接收信号所含信源数：，其中，为权重系数，一般设置为，其中为常数；盖尔圆按照半径由大小排列，当时，该接收信号的信源数为；
旋转算子的求解和速度估计步骤和波达方向的确定步骤如下：
将空时导向矢量构造乘成第一空时矩阵和第二空时矩阵，其中，两个矩阵之间滑动了一个空域导向矢量长度，脉冲间的旋转算子对应目标的多普勒频率：
；
；
根据阵列的平移不变性，将多普勒信号子空间沿列方向滑动一个阵列长度构造出两个子矩阵：
；
；
；
其中，为第一信号子空间矩阵，为第二信号子空间矩阵，为唯一的维满秩非奇异阵，与相似，具有相同的特征值；
根据最小二乘准则：，对进行特征值分解：
；
对特征值对角阵进行反演即可得到目标的速度：
，
其中表示取复数的相位角，表示取主对角线的第个元素；
同理，交换空时导向矢量的Kronecker积顺序构造时空导向矢量，并将其滑动一个时域导向矢量长度生成与旋转算子相关的两个子矩阵和：
，
，
重构后的信号表示为：
，
其中，代表包含不同顺序Kronecker积的导向矢量位置变化信息的置换矩阵；
对信号子空间进行重构得到包含旋转算子的角度信号子空间：
，
从划分两个子矩阵和：
；
；
构造旋转算子造的近似矩阵为：
，
其中，的特征向量为，为置换矩阵；
通过共享的特征向量固定置换矩阵，使得与特征值及其特征向量的排序一致，得到的特征值：
，
对特征值对角阵进行反演即可得到目标的角度：
。
2.根据权利要求1所述的波达方向及速度联合估计方法，其特征在于，所述获取雷达回波信号对应的目标距离门，基于目标距离门确定空时二维信号具体包括：
对雷达回波信号进行傅里叶处理和恒虚警率检测，得到目标距离门；
提取所述目标距离门对应的回波脉冲，并将提取到的所有回波脉冲排列为空时二维信号。
3.根据权利要求1所述的波达方向及速度联合估计方法，其特征在于，所述利用多级正交投影方式构建空时二维信号对应的匹配滤波器集合具体包括：
将空时二维信号中的位于最后的阵元的接收信息作为参考信息，将所述空时二维信号作为观测信息；
基于所述参考信息和所述观测信息的互相关构建所述空时二维信号对应的匹配滤波器集合。
4.根据权利要求3所述的波达方向及速度联合估计方法，其特征在于，所述基于所述参考信息和所述观测信息的互相关构建所述空时二维信号对应的匹配滤波器集合具体包括：
将所述参考信息作为目标参考信息，将观测信息作为目标观测信息；
提取所述目标参考信息和目标观测信息的互相关，并基于所述互相关构建一匹配滤波器；
通过所述匹配滤波器对目标观测信号进行投影，得到更新参考信号；
构造与所述匹配滤波器正交的阻塞矩阵，并基于所述阻塞矩阵和所述目标观测信号确定更新参考信号；
将所述更新参考信号作为目标参考信号，将所述更新参考信号作为目标参考信号；
重新执行提取所述目标参考信息和目标观测信息的互相关，并基于所述互相关构建一匹配滤波器的步骤，直至获取到预设数量个匹配滤波器，以得到对应的匹配滤波器集合。
5.一种波达方向及速度联合估计装置，其特征在于，所述装置包括：
信号获取模块，获取雷达回波信号对应的目标距离门，基于目标距离门确定空时二维信号；
第一信号处理模块，利用多级正交投影方式构建空时二维信号对应的匹配滤波器集合；
第二信号处理模块，对所述空时二维信号进行平滑处理以得到空时信号，并基于所述匹配滤波器集合以及所述空时信号，确定所述空时二维信号对应的信源数；
信号空间构建模块，基于所述匹配滤波器集合和所述信源数确定多普勒信号子空间；
目标计算模块，通过多普勒信号子空间求解所述空时二维信号，得到速度和波达方向；
其中，所述空时二维信号对应的信源数的具体过程包括：
对于N个阵元、K个脉冲的单快拍空时二维信号，先计算协方差矩阵：，再进行次前向空时平滑得到第二协方差矩阵：，其中代表平滑后的第二协方差矩阵大小；
对第二协方差矩阵进行矩阵分块处理，得到，其中，代表最后一个元素，表示最后一列除以外的列向量，对进行特征值分解得到特征向量矩阵，，其中，代表特征值构成主对角元素的对角矩阵，其中特征值排列由大到小；表示由特征向量构成的矩阵；构造一个酉矩阵并对协方差矩阵进行变换得到重构协方差矩阵，
，；
对于重构后的协方差矩阵，第个盖尔圆定义如下：，其中，重构协方差矩阵中第行第列的元素为第个盖尔圆的圆心，为第个盖尔圆的半径；
通过以下判别式估计出接收信号所含信源数：，其中，为权重系数，一般设置为，其中为常数；盖尔圆按照半径由大小排列，当时，该接收信号的信源数为；
旋转算子的求解和速度估计步骤和波达方向的确定步骤如下：
将空时导向矢量构造乘成第一空时矩阵和第二空时矩阵，其中，两个矩阵之间滑动了一个空域导向矢量长度，脉冲间的旋转算子对应目标的多普勒频率：
;
；
根据阵列的平移不变性，将多普勒信号子空间沿列方向滑动一个阵列长度构造出两个子矩阵：
;
;
;
其中，为第一信号子空间矩阵，为第二信号子空间矩阵，为唯一的维满秩非奇异阵，与相似，具有相同的特征值；
根据最小二乘准则：，对进行特征值分解：
；
对特征值对角阵进行反演即可得到目标的速度：
,
其中表示取复数的相位角，表示取主对角线的第个元素；
同理，交换空时导向矢量的Kronecker积顺序构造时空导向矢量，并将其滑动一个时域导向矢量长度生成与旋转算子相关的两个子矩阵和：
,
,
重构后的信号表示为：
，
其中，代表包含不同顺序Kronecker积的导向矢量位置变化信息的置换矩阵；
对信号子空间进行重构得到包含旋转算子的角度信号子空间：
,
从划分两个子矩阵和：
;
;
构造旋转算子造的近似矩阵为：
,
其中，的特征向量为，为置换矩阵；
通过共享的特征向量固定置换矩阵，使得与特征值及其特征向量的排序一致，得到的特征值：
,
对特征值对角阵进行反演即可得到目标的角度：
。
6.一种计算机可读存储介质，其特征在于，所述计算机可读存储介质存储有一个或者多个程序，所述一个或者多个程序可被一个或者多个处理器执行，以实现如权利要求1-4任意一项所述的波达方向及速度联合估计方法中的步骤。
7.一种终端设备，其特征在于，包括：处理器和存储器；
所述存储器上存储有可被所述处理器执行的计算机可读程序；
所述处理器执行所述计算机可读程序时实现如权利要求1-4任意一项所述的波达方向及速度联合估计方法中的步骤。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202311831328.3A CN117491985B (zh) | 2023-12-28 | 2023-12-28 | 一种波达方向及速度联合估计方法及相关装置 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202311831328.3A CN117491985B (zh) | 2023-12-28 | 2023-12-28 | 一种波达方向及速度联合估计方法及相关装置 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN117491985A CN117491985A (zh) | 2024-02-02 |
| CN117491985B true CN117491985B (zh) | 2024-05-10 |
Family
ID=89671183
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202311831328.3A Active CN117491985B (zh) | 2023-12-28 | 2023-12-28 | 一种波达方向及速度联合估计方法及相关装置 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN117491985B (zh) |
Citations (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101644760A (zh) * | 2009-08-27 | 2010-02-10 | 北京理工大学 | 一种适用于高分辨阵列的快速鲁棒的信源个数检测方法 |
| CN103383449A (zh) * | 2013-07-14 | 2013-11-06 | 西安电子科技大学 | 基于esprit算法的机载雷达近程杂波抑制方法 |
| CN110806557A (zh) * | 2019-11-22 | 2020-02-18 | 中国人民解放军国防科技大学 | 一种冲击噪声背景下单基地mimo雷达doa估计方法 |
| CN112379333A (zh) * | 2020-11-02 | 2021-02-19 | 武汉大学 | 基于空时维正交投影滤波的高频雷达海杂波抑制方法 |
| CN113406586A (zh) * | 2021-04-26 | 2021-09-17 | 北京理工大学 | 基于约束张量分解的mimo雷达二维波达方向估计方法 |
| CN114609595A (zh) * | 2022-03-03 | 2022-06-10 | 西安电子科技大学 | 一种频分正交mimo雷达信号处理方法 |
| CN115327505A (zh) * | 2022-10-12 | 2022-11-11 | 深圳大学 | 一种多径环境下低空目标角度估计方法及相关设备 |
| CN115514389A (zh) * | 2022-09-16 | 2022-12-23 | 西北工业大学 | 一种同步直扩信号的源数估计方法 |
-
2023
- 2023-12-28 CN CN202311831328.3A patent/CN117491985B/zh active Active
Patent Citations (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101644760A (zh) * | 2009-08-27 | 2010-02-10 | 北京理工大学 | 一种适用于高分辨阵列的快速鲁棒的信源个数检测方法 |
| CN103383449A (zh) * | 2013-07-14 | 2013-11-06 | 西安电子科技大学 | 基于esprit算法的机载雷达近程杂波抑制方法 |
| CN110806557A (zh) * | 2019-11-22 | 2020-02-18 | 中国人民解放军国防科技大学 | 一种冲击噪声背景下单基地mimo雷达doa估计方法 |
| CN112379333A (zh) * | 2020-11-02 | 2021-02-19 | 武汉大学 | 基于空时维正交投影滤波的高频雷达海杂波抑制方法 |
| CN113406586A (zh) * | 2021-04-26 | 2021-09-17 | 北京理工大学 | 基于约束张量分解的mimo雷达二维波达方向估计方法 |
| CN114609595A (zh) * | 2022-03-03 | 2022-06-10 | 西安电子科技大学 | 一种频分正交mimo雷达信号处理方法 |
| CN115514389A (zh) * | 2022-09-16 | 2022-12-23 | 西北工业大学 | 一种同步直扩信号的源数估计方法 |
| CN115327505A (zh) * | 2022-10-12 | 2022-11-11 | 深圳大学 | 一种多径环境下低空目标角度估计方法及相关设备 |
Non-Patent Citations (4)
| Title |
|---|
| 一种低复杂度的信号子空间拟合的新方法;黄磊等;电子学报;20050625(第06期);23-27 * |
| 基于二维CP-GTD模型的全极化ISAR超分辨成像;代大海等;自然科学进展;20071015(第10期);121-130 * |
| 基于实值信号子空间的虚拟阵列解相干算法;李欣等;哈尔滨工程大学学报;20120828(第09期);70-75 * |
| 高效相干源DOA估计的ESPRIT新方法;蔡会甫等;计算机应用研究;20100115(第01期);320-322 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN117491985A (zh) | 2024-02-02 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Potter et al. | Attributed scattering centers for SAR ATR | |
| Conte et al. | Distributed target detection in compound-Gaussian noise with Rao and Wald tests | |
| CN111736131A (zh) | 一种剔除一比特信号谐波虚假目标的方法及相关组件 | |
| CN109669182B (zh) | 无源双基地sar动/静目标联合稀疏成像方法 | |
| Wang et al. | Adaptive imaging for forward-looking ground penetrating radar | |
| CN109324322A (zh) | 一种基于被动相控阵天线的测向与目标识别方法 | |
| CN110231617A (zh) | 目标障碍物方位探测方法、装置、车载雷达及存储介质 | |
| Tao et al. | A knowledge aided SPICE space time adaptive processing method for airborne radar with conformal array | |
| Li et al. | Dimension-reduced space-time adaptive clutter suppression algorithm based on lower-rank approximation to weight matrix in airborne radar | |
| CN113504522A (zh) | 一种基于发射天线随机切换的时空解耦与超分辨测角方法 | |
| CN107576947A (zh) | 基于时间平滑的l型阵对相干信源二维波达方向估计方法 | |
| CN108020834B (zh) | 基于改进edpca的运动目标检测方法、装置及电子设备 | |
| Xu et al. | Joint Doppler and DOA estimation using (Ultra-) Wideband FMCW signals | |
| CN106772275A (zh) | 一种基于Alpha稳定分布的低秩矩阵恢复检测方法 | |
| CN115508799A (zh) | 一种基于矩空间的分布式被动雷达目标检测方法 | |
| Ali Khan et al. | Direction of arrival estimation using adaptive directional time-frequency distributions | |
| Yang et al. | IMIA: Interference mitigation via iterative approaches for automotive radar | |
| Nickel | Radar target parameter estimation with array antennas | |
| CN117491985B (zh) | 一种波达方向及速度联合估计方法及相关装置 | |
| Pauciullo et al. | Detection of partially coherent scatterers in multidimensional SAR tomography: A theoretical study | |
| Thammakhoune et al. | Moving target imaging for synthetic aperture radar via RPCA | |
| Anitori et al. | Array aperture extrapolation using sparse reconstruction | |
| CN111308436B (zh) | 基于体积相关函数的雷达空时自适应处理方法及装置 | |
| WO2021196165A1 (zh) | 频率分析方法、装置及雷达 | |
| Dehkordi et al. | Region of interest based adaptive high resolution parameter estimation with applications in automotive radar |
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 |