CN111289976B - 阵列3-d成像检测***以及成像方法 - Google Patents

阵列3-d成像检测***以及成像方法 Download PDF

Info

Publication number
CN111289976B
CN111289976B CN202010168238.0A CN202010168238A CN111289976B CN 111289976 B CN111289976 B CN 111289976B CN 202010168238 A CN202010168238 A CN 202010168238A CN 111289976 B CN111289976 B CN 111289976B
Authority
CN
China
Prior art keywords
frequency
target
imaging
array
distance
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
CN202010168238.0A
Other languages
English (en)
Other versions
CN111289976A (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.)
Suzhou Weimo Electronic Information Technology Co.,Ltd.
Original Assignee
Suzhou Weimo Electronic Information Technology 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 Suzhou Weimo Electronic Information Technology Co ltd filed Critical Suzhou Weimo Electronic Information Technology Co ltd
Priority to CN202010168238.0A priority Critical patent/CN111289976B/zh
Publication of CN111289976A publication Critical patent/CN111289976A/zh
Application granted granted Critical
Publication of CN111289976B publication Critical patent/CN111289976B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO 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/00Systems 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/88Radar or analogous systems specially adapted for specific applications
    • G01S13/89Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO 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/00Systems 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/02Systems using reflection of radio waves, e.g. primary radar systems; Analogous systems
    • G01S13/06Systems determining position data of a target
    • G01S13/08Systems for measuring distance only
    • G01S13/10Systems for measuring distance only using transmission of interrupted, pulse modulated waves
    • G01S13/18Systems for measuring distance only using transmission of interrupted, pulse modulated waves wherein range gates are used

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Physics & Mathematics (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • General Physics & Mathematics (AREA)
  • Electromagnetism (AREA)
  • Radar Systems Or Details Thereof (AREA)

Abstract

本发明的阵列3‑D成像检测***以及成像方法,其中成像检测***包括照射天线、接收天线阵、发射机、接收机、显示***。本发明阵列3‑D成像方法,包括:步骤一,确定***参数;步骤二,确定天线单元的间距;步骤三,确定阵列形状;步骤四,确定FFT运算长度;步骤五,确定***采样频率;步骤六,确定调频连续波照射源的调频斜率;步骤七,提供阵列3‑D成像检测;步骤八,提供数字信号处理***;步骤九,计算不同距离切面上的像;步骤十,精确测距。本发明提供了一种高精度、高分辨率的成像检测***以及成像方法,优化选择合适的阵列参数,以较小的阵列规模实现宽视角、高精度的目标探测与成像识别。

Description

阵列3-D成像检测***以及成像方法
技术领域
本发明涉及目标探测与成像识别领域,特别涉及一种阵列3-D成像检测***以及成像方法。
背景技术
在目标探测与成像识别领域,现有的技术主要分为以下三类:
第一类为基于可见光、红外、激光的光学探测与成像***,如摄像头、夜视仪、激光雷达等,该类***往往会受到使用条件的限制,在复杂天气情况或有遮挡的情况下会造成***失效,此外,强激光还会危害人身安全,这使得该类***在使用上具有一定的局限性。
第二类为基于电离辐射的探测与成像***,如X射线、计算机断层扫描(CT)等***,该类设备能够对遮挡目标进行探测与成像,但往往设备笨重、体积庞大,且高剂量照射对人体具有较大的危害。
第三类为基于非电离辐射的探测与成像***,如雷达、B超、核磁共振、微波CT等,该类设备往往工作频率越高探测精度越高,并且往往设备笨重、价格昂贵。
综上所述,在目标探测与成像识别领域,现有技术在成本和普适性等方面存在不足,特别在微波成像方面,需要开发低成本、高可行性、高精度的探测与成像识别技术。
发明内容
为了解决现有成像技术存在的问题,本发明的目的在于,提供一种既能够提供高精度的目标探测与成像识别,又结构简单、成本低、非接触式的阵列3-D成像检测***以及成像方法。
本发明的阵列3-D成像方法,包括:
步骤一、确定***基本参数,***基本参数包括工作频率、探测空域远界、探测视角范围、距离分辨率、角度分辨率、目标最大轴向速度;
步骤二、确定天线单元的基准间距,根据探测视角确定天线单元的最大间距,对于最远处的目标,***探测视角范围为β,则天线单元最大间距的计算公式为:
Figure BDA0002408231520000021
其中,λ为波长;
步骤三、确定阵列形状及规模,根据分辨率选择阵列的等效焦距或配置探测远界来压缩阵列孔径,当天线单元的最大间隔确定后,根据横向距离成像分辨率或角度分辨率确定阵列规模,若满足有效分辨最远处两个横向间距为D的目标,则接收阵列的横向最小单元数为:
Figure BDA0002408231520000022
其中,γ为透镜焦距F与目标最远轴向距离Rmax的比值,即:γ=F/Rmax。γ<Rmin/Rmax,其中Rmin为最小检测距离,虚拟透镜的焦距值为:
F=γRmax
阵列形状为矩形、圆形或椭圆形,阵列为均匀阵列或非均匀阵列;
步骤四、确定FFT运算长度;
当探测远界为Rmax,距离分辨率为ΔR时,对每个距离单元进行频率抽取,选择最小FFT运算长度Nfft为:
Figure BDA0002408231520000023
步骤五、确定***采样频率,最小采样频率为:
Figure BDA0002408231520000031
其中,v为目标的最大轴向运动速度;
步骤六、确定调频连续波照射源的调频斜率,最小的线性调频斜率为:
Figure BDA0002408231520000032
步骤七、提供阵列3-D成像检测***;
步骤八、提供数字信号处理***;
步骤九、计算不同距离切面上的像,采用数字成像算法计算不同距离切面上的像,形成一系列成像切片;
步骤十、精确测距。
本发明的阵列3-D成像方法,其中,所述计算不同距离切面上的像,采用数字成像算法计算不同距离切面上的像,形成一系列成像切片包括:
分步骤一:计算天线单元相移,建立空间直角坐标系,建立天线阵的法线方向为-z轴的空间直角坐标系,坐标原点位于天线阵中心,每个天线单元的坐标为(Xmn,Ymn,0),其中m,n为天线单元的横向与纵向编号,计算天线单元对应虚拟透镜的移相量:
Figure BDA0002408231520000033
其中RL为天线阵列的有效半径,dmn为天线单元距离天线阵中心的距离,其计算公式为:
Figure BDA0002408231520000034
分步骤二:计算通道信号的频谱,对每个通道的时域信号进行快速傅里叶变换,求出其频谱分布,
Figure BDA0002408231520000035
分步骤三:恒虚警检测与测频,选取一个通道的FFT数据进行目标检测,检测门限采用单元平均法进行计算:
Figure BDA0002408231520000041
其中,NCFAR不大于Nfft,当频谱幅度超过上述门限时,确认该频率通道存在目标,从而确定目标对应的频率;
分步骤四:频率抽取,根据所述分步骤三确定存在目标的频率通道,依次串行抽取相应的频率通道FFT数据送后级进行处理;
分步骤五:计算目标距离并离散化目标所在球面区域,根据频率抽取对应的频率值,计算目标的距离:
最大频偏为fs/2,对应着目标最远距离,当目标出现在频率为f的位置时,对应的目标距离为:
Figure BDA0002408231520000042
目标在半径为R的球面上,对此球面区域按p行q列进行剖分和离散化,剖分单元的立体角应满足关系式:
Figure BDA0002408231520000043
进而求出离散化球面区域的坐标值(xu,yu,zu);
分步骤六:计算像场区域并离散化,将目标所在球面区域映射到像场区域,获得离散化的像场区域坐标,成像区域所在的位置坐标为:
Figure BDA0002408231520000044
Figure BDA0002408231520000045
Figure BDA0002408231520000046
从而获得成像区域每个剖分节点的坐标(xpq,ypq,zpq),其中p、q为剖分节点的横向与纵向编号;
分步骤七:计算传播相移,计算天线单元到像场网格节点的传播路程:
Figure BDA0002408231520000051
计算传播相移:
Figure BDA0002408231520000052
分步骤八:计算像场值,像场计算公式如下:
Figure BDA0002408231520000053
其中,Pmn为抽取出的频率为f的频谱数据,Γmn,pq为传播系数,其计算公式如下:
Figure BDA0002408231520000054
在简化条件下可直接选取Γmn,pq=1;
分步骤九:像场反演与坐标测量,
将三维立体像进行坐标变换,使得像与真实目标的尺寸、方向一致,坐标为(xv,yv,zv)的像,其反演公式为:
Figure BDA0002408231520000055
Figure BDA0002408231520000056
Figure BDA0002408231520000057
得出目标坐标后求得目标的角度信息;
分步骤十:判断是否需要精细分辨目标,当需要更精细的分辨目标细节时,仅需要减小目标所在球面的离散间隔,执行分步骤五;
分步骤十一:判断程序是否结束,“是”则终止程序,“否”则重新执行分步骤二。
本发明的阵列3-D成像方法,其中,NCFAR=Nfft/10。
本发明的阵列3-D成像方法,其中,分步骤一、分步骤二之间还包括:
距离选择步骤,当阵列3-D成像检测***采用脉冲体制时,采用不同的距离选择波门对每个通道的时域信号进行距离滤波,利用常规门函数确定距离波门:
Figure BDA0002408231520000061
其中t0为距离波门的中心坐标,τ为距离波门的宽度,经距离波门滤波后的时域信号为:
E(t)=E0(t)·g(t)。
本发明的阵列3-D成像方法,其中,当需要更精确的测量目标坐标时,根据精度要求选择更长的FFT运算长度,测距精度为:
Figure BDA0002408231520000062
采用三角形线性调频方案时,分别计算正向和负向调频斜率时的目标频率fu、fd,目标的精确距离为:
Figure BDA0002408231520000063
其中,c=3×108m/s为光速。
本发明的阵列3-D成像检测***,包括照射天线、接收天线阵、发射机、接收机、本地振荡器、A/D、信号处理机、显示***,其中:
所述照射天线用于射频信号的发射;
发射机用于对本振信号进行放大,与照射天线连接;
接收天线阵用于接收目标反射信号;
接收机用于对接收天线阵的射频信号的放大、混频、滤波,将射频信号变频为基带信号或中频信号,与接收天线阵的天线一一对应连接;
本地振荡器用于产生基准射频信号;
A/D用于将模拟信号转换为数字信号,与接收机连接;
数字信号处理***用于信号检测、信号合成与目标成像;
显示***用于三维立体成像的显示以及人机交互,所述阵列3-D成像检测***由本发明的阵列3-D成像方法确定。
本发明的阵列3-D成像检测***,其中,数字信号处理***包括距离选择模块、FFT模块、CFAR模块、目标检测/测频模块、频率抽取模块、数字透镜成像***、图像处理***,其中:
距离选择模块用于采用不同的距离选择波门对每个通道的时域信号进行距离滤波;
FFT模块用于每个接收通道信号的快速傅里叶变换,将时域信号变换为频域信号;
CFAR模块用于采用单元平均法生成检测门限并进行过门限目标检测;
目标检测/测频模块用于检测不同的目标并根据FFT数据进行频率测量,以确定每个目标对应的频率值;
频率抽取模块用于在测频模块的控制下,抽取出存在目标的频域数据送后级电路或根据用户实际需求按顺序依次抽取频域数据送后级电路;
数字透镜成像***用于对接收到的频域抽取数据进行成像处理,顺序生成不同频率通道的目标的三维像;
图像处理***用于对数字透镜成像***生成的三维像相进行整理、反演变换,转换为人眼可识别的图像。
本发明的阵列3-D成像检测***,其中,照射天线采用宽波束天线,接收天线阵包括多个天线单元,发射机包括放大电路、驱动电路和功率放大器,接收机包括低噪声放大器、混频器、滤波器,本地振荡器包括直接数字频率合成器或锁相环;A/D包括模数转换芯片及板卡;信号处理***为搭载嵌入式***的数字板卡。
本发明的技术方案提供了一种阵列目标探测与成像识别一体化解决方案,克服了现有微波焦平面成像***、微波摄像头成像精度不高的缺点,能够利用有限的阵列孔径实现高精度的目标成像识别。
本发明的技术方案提供了一种基于连续波/脉冲微波照射、通过距离选择和角度高分辨成像实现了类似计算机断层扫描(CT)的目标探测与三维成像方案,从而能够更准确的分辨目标形状。
本发明的技术方案提供了一种高精度、高分辨率的阵列综合方法,综合考虑目标距离、探测视角、成像分辨率等因素,优化选择合适的阵列参数,以较小的阵列规模实现宽视角、高精度的成像识别。
附图说明
图1为本发明的阵列3-D成像检测***的组成框图;
图2为数字信号处理***的组成框图;
图3为本发明的阵列3-D成像方法的示意图;
图4为采用数字成像算法计算不同距离切面上的像的流程图;
图5为本发明的阵列3-D成像检测***的接收天线阵列的可选形状示意图;
图6为金属物体以及利用本发明的技术方案对上述金属物体探测成像的效果图。
具体实施方式
本发明的目的是提供一种目标探测与成像识别一体化解决方案,在实现目标探测的同时提供高精度目标成像,实现了高性能的一体化目标探测和三维成像识别。
如图1、图2、图3、图4、图5、图6所示,本发明的阵列3-D成像方法,包括:
步骤一、确定***基本参数,***基本参数包括工作频率、探测空域远界、探测视角范围、距离分辨率、角度分辨率、目标最大轴向速度;
步骤二、确定天线单元的基准间距,根据探测视角确定天线单元的最大间距,对于最远处的目标,***探测视角范围为β,则天线单元最大间距的计算公式为:
Figure BDA0002408231520000091
其中,λ为波长;
步骤三、确定阵列形状及规模,根据分辨率选择阵列的等效焦距或配置探测远界来压缩阵列孔径,当天线单元的最大间隔确定后,根据横向距离成像分辨率或角度分辨率确定阵列规模,若满足有效分辨最远处两个横向间距为D的目标,则接收阵列的横向最小单元数为:
Figure BDA0002408231520000092
其中,γ为透镜焦距F与目标最远轴向距离Rmax的比值,即:γ=F/Rmax。γ<Rmin/Rmax,其中Rmin为最小检测距离,虚拟透镜的焦距值为:
F=γRmax
阵列形状为矩形、圆形或椭圆形,阵列为均匀阵列或非均匀阵列;
步骤四、确定FFT运算长度;
当探测远界为Rmax,距离分辨率为ΔR时,对每个距离单元进行频率抽取,选择最小FFT运算长度Nfft为:
Figure BDA0002408231520000093
步骤五、确定***采样频率,最小采样频率为:
Figure BDA0002408231520000094
其中,v为目标的最大轴向运动速度;
步骤六、确定调频连续波照射源的调频斜率,最小的线性调频斜率为:
Figure BDA0002408231520000101
步骤七、提供阵列3-D成像检测***;
步骤八、提供数字信号处理***;
步骤九、计算不同距离切面上的像,采用数字成像算法计算不同距离切面上的像,形成一系列成像切片;
步骤十、精确测距。
本发明的阵列3-D成像方法,其中,所述计算不同距离切面上的像,采用数字成像算法计算不同距离切面上的像,形成一系列成像切片包括:
分步骤一:计算天线单元相移,建立空间直角坐标系,建立天线阵的法线方向为-z轴的空间直角坐标系,坐标原点位于天线阵中心,每个天线单元的坐标为(Xmn,Ymn,0),其中m,n为天线单元的横向与纵向编号,计算天线单元对应虚拟透镜的移相量:
Figure BDA0002408231520000102
其中RL为天线阵列的有效半径,dmn为天线单元距离天线阵中心的距离,其计算公式为:
Figure BDA0002408231520000103
分步骤二:计算通道信号的频谱,对每个通道的时域信号进行快速傅里叶变换,求出其频谱分布,
Figure BDA0002408231520000104
分步骤三:恒虚警检测与测频,选取一个通道的FFT数据进行目标检测,检测门限采用单元平均法进行计算:
Figure BDA0002408231520000111
其中,NCFAR不大于Nfft,当频谱幅度超过上述门限时,确认该频率通道存在目标,从而确定目标对应的频率;
分步骤四:频率抽取,根据所述分步骤三确定存在目标的频率通道,依次串行抽取相应的频率通道FFT数据送后级进行处理;
分步骤五:计算目标距离并离散化目标所在球面区域,根据频率抽取对应的频率值,计算目标的距离:
最大频偏为fs/2,对应着目标最远距离,当目标出现在频率为f的位置时,对应的目标距离为:
Figure BDA0002408231520000112
目标在半径为R的球面上,对此球面区域按p行q列进行剖分和离散化,剖分单元的立体角应满足关系式:
Figure BDA0002408231520000113
进而求出离散化球面区域的坐标值(xu,yu,zu);
分步骤六:计算像场区域并离散化,将目标所在球面区域映射到像场区域,获得离散化的像场区域坐标,成像区域所在的位置坐标为:
Figure BDA0002408231520000114
Figure BDA0002408231520000115
Figure BDA0002408231520000116
从而获得成像区域每个剖分节点的坐标(xpq,ypq,zpq),其中p、q为剖分节点的横向与纵向编号;
分步骤七:计算传播相移,计算天线单元到像场网格节点的传播路程:
Figure BDA0002408231520000121
计算传播相移:
Figure BDA0002408231520000122
分步骤八:计算像场值,像场计算公式如下:
Figure BDA0002408231520000123
其中,Pmn为抽取出的频率为f的频谱数据,Γmn,pq为传播系数,其计算公式如下:
Figure BDA0002408231520000124
在简化条件下可直接选取Γmn,pq=1;
分步骤九:像场反演与坐标测量,
将三维立体像进行坐标变换,使得像与真实目标的尺寸、方向一致,坐标为(xv,yv,zv)的像,其反演公式为:
Figure BDA0002408231520000125
Figure BDA0002408231520000126
Figure BDA0002408231520000127
得出目标坐标后求得目标的角度信息;
分步骤十:判断是否需要精细分辨目标,当需要更精细的分辨目标细节时,仅需要减小目标所在球面的离散间隔,执行分步骤五;
分步骤十一:判断程序是否结束,“是”则终止程序,“否”则重新执行分步骤二。
本发明的阵列3-D成像方法,其中,NCFAR=Nfft/10。
本发明的阵列3-D成像方法,其中,分步骤一、分步骤二之间还包括:
距离选择步骤,当阵列3-D成像检测***采用脉冲体制时,采用不同的距离选择波门对每个通道的时域信号进行距离滤波,利用常规门函数确定距离波门:
Figure BDA0002408231520000131
其中t0为距离波门的中心坐标,τ为距离波门的宽度,经距离波门滤波后的时域信号为:
E(t)=E0(t)·g(t)。
本发明的阵列3-D成像方法,其中,当需要更精确的测量目标坐标时,根据精度要求选择更长的FFT运算长度,测距精度为:
Figure BDA0002408231520000132
采用三角形线性调频方案时,分别计算正向和负向调频斜率时的目标频率fu、fd,目标的精确距离为:
Figure BDA0002408231520000133
其中,c=3×108m/s为光速。
本发明的阵列3-D成像检测***,包括照射天线、接收天线阵、发射机、接收机、本地振荡器、A/D、信号处理机、显示***,其中:
所述照射天线用于射频信号的发射;
发射机用于对本振信号进行放大,与照射天线连接;
接收天线阵用于接收目标反射信号;
接收机用于对接收天线阵的射频信号的放大、混频、滤波,将射频信号变频为基带信号或中频信号,与接收天线阵的天线一一对应连接;
本地振荡器用于产生基准射频信号;
A/D用于将模拟信号转换为数字信号,与接收机连接;
数字信号处理***用于信号检测、信号合成与目标成像;
显示***用于三维立体成像的显示以及人机交互,所述阵列3-D成像检测***由本发明的阵列3-D成像方法确定。
本发明的阵列3-D成像检测***,其中,数字信号处理***包括距离选择模块、FFT模块、CFAR模块、目标检测/测频模块、频率抽取模块、数字透镜成像***、图像处理***,其中:
距离选择模块用于采用不同的距离选择波门对每个通道的时域信号进行距离滤波;
FFT模块用于每个接收通道信号的快速傅里叶变换,将时域信号变换为频域信号;
CFAR模块用于采用单元平均法生成检测门限并进行过门限目标检测;
目标检测/测频模块用于检测不同的目标并根据FFT数据进行频率测量,以确定每个目标对应的频率值;
频率抽取模块用于在测频模块的控制下,抽取出存在目标的频域数据送后级电路或根据用户实际需求按顺序依次抽取频域数据送后级电路;
数字透镜成像***用于对接收到的频域抽取数据进行成像处理,顺序生成不同频率通道的目标的三维像;
图像处理***用于对数字透镜成像***生成的三维像相进行整理、反演变换,转换为人眼可识别的图像。
本发明的阵列3-D成像检测***,其中,照射天线采用宽波束天线,接收天线阵包括多个天线单元,发射机包括放大电路、驱动电路和功率放大器,接收机包括低噪声放大器、混频器、滤波器,本地振荡器包括直接数字频率合成器或锁相环;A/D包括模数转换芯片及板卡;信号处理***为搭载嵌入式***的数字板卡。
本发明的技术方案提供了一种阵列目标探测与成像识别一体化解决方案,克服了现有微波焦平面成像***、微波摄像头成像精度不高的缺点,能够利用有限的阵列孔径实现高精度的目标成像识别。
本发明的技术方案提供了一种基于连续波/脉冲微波照射、通过距离选择和角度高分辨成像实现了类似计算机断层扫描(CT)的目标探测与三维成像方案,从而能够更准确的分辨目标形状。
本发明的技术方案提供了一种高精度、高分辨率的阵列综合方法,综合考虑目标距离、探测视角、成像分辨率等因素,优化选择合适的阵列参数,以较小的阵列规模实现宽视角、高精度的成像识别。
本发明的技术方案具体包括以下步骤:
步骤一,确定***基本参数。***基本参数包含工作频率、探测空域远界、探测视角范围、距离分辨率、角度分辨率、目标最大轴向速度等。***基本参数主要取决于用户需求,由***实际工作时面临的工作场景来确定相关参数。
步骤二,确定天线单元的基准间距、阵列形状。天线单元的间距主要取决于探测视角,选择阵列形状时可综合考虑探测视角、阵列成本等。
天线单元的间距主要取决于探测视角范围,对于最远处的目标,若***探测视角范围为β,则天线单元的最大间隔为:
Figure BDA0002408231520000151
其中λ为波长,可以采用等间距阵列或不等间距的非均匀阵列,推荐使用非均匀阵列,能够改善探测性能并能够对副像有效抑制。
步骤三,确定阵列规模。分辨率是影响阵列规模和孔径的主要因素,在分辨率一定的情况下,可通过合理选择阵列的等效焦距或配置合适的探测远界来压缩阵列孔径。
当天线单元的最大间隔确定后,可以根据横向距离成像分辨率或角度分辨率确定阵列规模,若要求能够有效分辨最远处两个横向间距为D的目标,则接收阵列的横向最小单元数为:
Figure BDA0002408231520000161
其中γ为透镜焦距F与目标最远轴向距离Rmax的比值,即:γ=F/Rmax。一般情况下要求γ<Rmin/Rmax,其中Rmin为最小检测距离,适当减小γ可以获得更小的Nmin值,但实际上通过减小最远探测距离Rmax来降低阵列横向单元数目更有效。在某些对成像分辨率要求不高的应用中,可以设置较大的D值,此时***仍然能够具有较好的目标角度分辨能力。
确定γ值后,虚拟透镜的焦距值为:
F=γRmax
在确定阵列最大单元间距Δ和最小横向单元数目后,进而可确定阵列形状,可选择的阵列形状有矩形、圆形甚至椭圆形等几种形式。在阵列间距设置方面,优先选择不等间距的非均匀阵列。可供选择的阵列形形状见附图5。
步骤四,确定FFT运算长度。FFT运算主要用来实现距离分辨,与***的探测远界和距离分辨率有关。
当探测远界为Rmax,距离分辨率为ΔR时,需要对每个距离单元进行频率抽取,因而可以选择最小FFT运算长度Nfft为:
Figure BDA0002408231520000162
除了采用FFT运算来进行快速测频外,还可以采用滤波器组来进行频率选择与距离分辨,但会耗费大量硬件资源。
步骤五,确定***采样频率。选择A/D采样频率需要综合考虑探测远界、距离分辨率、目标最大轴向速度、***硬件资源等因素。
在保证距离分辨精度的情况下,需要确保一个完整FFT运算长度内,目标轴向运动距离不超过距离分辨率ΔR,从而有最小采样频率为:
Figure BDA0002408231520000163
其中v为目标的最大轴向运动速度。
步骤六,确定调频连续波照射源的调频斜率。调频范围对应着目标的距离范围,在A/D采样频率确定的情况下,最大频偏不得大于采样率的一半,在一个FFT运算长度内需要扫描覆盖最大频偏。
最小的线性调频斜率为:
Figure BDA0002408231520000171
步骤七,提供阵列3-D成像检测***。可选择线性调频连续波测距+阵列成像或脉冲测距+阵列成像的***体制,***既具有测距能力,又具有角度高分辨成像能力,从而实现空间三维成像,达到类似CT的断层扫描成像效果。
本发明的阵列3-D成像检测***采用一发多收+数字透镜成像的技术体制,主要由照射天线、接收天线阵、发射机、接收机、本地振荡器、A/D、信号处理机、显示***、电源等模块组成,见附图1,其中照射天线完成射频信号的发射,由接收天线阵列的法线方向对目标进行照射,天线形式可以采用宽波束天线。接收天线阵完成目标反射信号的接收,天线单元间距及天线阵的规模和形状按步骤二和步骤三的方法进行提供,天线单元可根据实际需要灵活选择,推荐优先使用微带贴片天线阵列。发射机对本振信号进行放大,由其末级功率放大器形成大功率微波射频信号,其主要由放大电路、驱动电路和功率放大器组成。接收机完成天线接收射频信号的放大、混频、滤波,将射频信号变频为基带信号或中频信号,其主要由低噪声放大器、混频器、滤波器等组成。本地振荡器主要用来产生基准射频信号,本地振荡器主要由直接数字频率合成器或锁相环等组成。A/D完成模拟信号到数字信号的转换,其中A/D采样模块可集成到信号处理***中。数字信号处理***完成信号检测、信号合成与目标成像等功能,其组成形式为搭载嵌入式***的数字板卡,软件和算法是其技术核心。显示***完成三维立体成像的显示以及人机交互。电源用来提供***工作所需的各种电压电流。
在某些对测距无要求或不需要精确分辨目标细节的应用场景中,可以采用点频连续波进行照射或采用脉冲调制波形进行照射,对应的本振模块需要工作在点频模式,发射机需要增加脉冲调制功能。相对连续波线性调频体制,采用脉冲测距体制,***的硬件构造复杂、成本高,需要根据实际需求情况来优选***方案。
步骤八,提供数字信号处理***。数字信号处理***是整个***的灵魂和核心,需要完成目标检测、距离测量、角度测量、三维成像等一系列功能。
数字信号处理***主要由距离选择、FFT模块、CFAR模块、目标检测/测频模块、频率抽取模块、数字透镜成像***、图像处理***等功能模块组成,见附图2,硬件主要由FPGA、DSP等大规模集成电路构成。距离选择模块仅在脉冲测距方案中使用,在时域用不同的距离波门对回波信号进行选择和过滤。FFT模块完成每个接收通道信号的快速傅里叶变换,将时域信号变换为频域信号。CFAR模块完成恒虚警检测,采用单元平均法生成检测门限并进行过门限检测。目标检测/测频模块检测不同的目标并根据FFT数据进行频率测量,确定每个目标对应的频率值。频率抽取模块可在测频模块的控制下,抽取出存在目标的频域数据送后级电路,也可以根据用户实际需求按顺序依次抽取频域数据送后级电路。数字成像***对接收到的频域抽取数据进行成像处理,顺序生成不同频率通道的目标的三维像。图像处理***对数字式透镜成像***生成的三维像相进行整理、反演变换,转换为人眼可识别的图像。
步骤九,计算不同距离切面上的像。采用数字成像算法计算不同距离切面上的像,形成一系列成像切片,实现类似CT的成像效果。
算法组成框图见附图3,算法流程见附图4,成像效果图见附图6。
图6中左边分别是十字形、带圈十字形的金属物,相互之间有一定的距离,经过本发明的阵列3-D成像检测***探测成像后,得出右边的两个空间像分布,进而可分辨出两个物体的形状、相对位置关系等。
分步骤一:计算天线单元相移。建立空间直角坐标系,例如,建立天线阵的法线方向为-z轴的空间直角坐标系,坐标原点位于天线阵中心。每个天线单元的坐标为(Xmn,Ymn,0),其中m,n为天线单元的横向与纵向编号。
计算天线单元对应虚拟透镜的移相量:
Figure BDA0002408231520000191
其中RL为天线阵列的有效半径,dmn为天线单元距离天线阵中心的距离,其计算公式为:
Figure BDA0002408231520000192
分步骤二(可选):距离选择。采用不同的距离选择波门对每个通道的时域信号进行距离滤波。
距离波门通常可以选择常规门函数:
Figure BDA0002408231520000193
其中t0为距离波门的中心坐标,τ为距离波门的宽度。经距离波门滤波后的时域信号为:
E(t)=E0(t)·g(t)
分步骤三:计算通道信号的频谱。对每个通道的时域信号进行快速傅里叶变换,求出其频谱分布。
Figure BDA0002408231520000194
分步骤四:恒虚警检测与测频。选取一个通道的FFT数据进行目标检测,检测门限采用单元平均法进行计算:
Figure BDA0002408231520000195
其中NCFAR不大于Nfft,通常可选取为NCFAR=Nfft/10。当频谱幅度超过上述门限时认为该频率通道存在目标,从而可确定目标对应的频率。
分步骤五:频率抽取。根据分步骤四确定存在目标的频率通道,依次串行抽取相应的频率通道FFT数据送后级进行处理。
分步骤六:计算目标距离并离散化目标所在球面区域。根据频率抽取对应的频率值,计算目标的距离:
最大频偏为fs/2,对应着目标最远距离。当目标出现在频率为f的位置时,对应的目标距离为:
Figure BDA0002408231520000201
目标在半径为R的球面上,对此球面区域按p行q列进行剖分和离散化,为保证能够分辨距离为D的两个目标,剖分单元的立体角应满足关系式:
Figure BDA0002408231520000202
进而求出离散化球面区域的坐标值(xu,yu,zu)。
分步骤七:计算像场区域并离散化。将目标所在球面区域映射到像场区域,获得离散化的像场区域坐标。
成像区域所在的位置坐标为:
Figure BDA0002408231520000203
Figure BDA0002408231520000204
Figure BDA0002408231520000205
从而获得成像区域每个剖分节点的坐标(xpq,ypq,zpq),其中p,q为剖分节点的横向与纵向编号。
分步骤八:计算传播相移。先计算天线单元到像场网格节点的传播路程:
Figure BDA0002408231520000211
然后计算传播相移:
Figure BDA0002408231520000212
分步骤九:计算像场值。像场计算公式如下:
Figure BDA0002408231520000213
其中,Pmn为抽取出的频率为f的频谱数据,Γmn,pq为传播系数,其计算公式如下:
Figure BDA0002408231520000214
在简化条件下可直接选取Γmn,pq=1。
分步骤十:像场反演与坐标测量。采用数字成像算法合成的像与真实目标存在差异,需要进行像场反演来获得与目标尺寸、方向一致的像。
像场反演的本质是坐标变换,是分步骤七的逆运算,将三维立体像进行坐标变换,使得像与真实目标的尺寸、方向一致。坐标为(xv,yv,zv)的像,其反演公式为:
Figure BDA0002408231520000221
Figure BDA0002408231520000222
Figure BDA0002408231520000223
上式也是像场坐标变换为目标坐标的计算公式,得出目标坐标后即可求得目标的角度信息。
分步骤十一:判断是否需要精细分辨目标。当需要更精细的分辨目标细节时,仅需要减小目标所在球面的离散间隔,重新执行分步骤六即可。
分步骤十二:判断程序是否结束。“是”则终止程序,“否”则重新执行分步骤二。
步骤十,精确测距。
当需要更精确的测量目标坐标时,可以通过选择更长的FFT运算长度,对某一通道的信号进行FFT变换,当数据长度不足时可以进行补零,通过扩展FFT长度,可获得更高的频率分辨力和测距精度。
为降低多普勒频移对测距的影响,可以采用三角形线性调频方案,分别计算正向和负向调频斜率时的目标频率fu、fd,目标的精确距离为:
Figure BDA0002408231520000224
其中,c=3×108m/s为光速。
以上所述仅是本发明的优选实施方式,应当指出,对于本技术领域的普通技术人员来说,在不脱离本发明原理的前提下,还可以做出若干改进和润饰,这些改进和润饰也应视为本发明的保护范围。

Claims (7)

1.一种阵列3-D成像方法,其特征在于,包括:
步骤一、确定***基本参数,***基本参数包括工作频率、探测空域远界、探测视角范围、距离分辨率、角度分辨率、目标最大轴向速度;
步骤二、确定天线单元的基准间距,根据探测视角确定天线单元的最大间距,对于最远处的目标,***探测视角范围为β,则天线单元最大间距的计算公式为:
Figure FDA0003369223480000011
其中,λ为波长;
步骤三、确定阵列形状及规模,根据分辨率选择阵列的等效焦距或配置探测远界来压缩阵列孔径,当天线单元的最大间隔确定后,根据横向距离成像分辨率或角度分辨率确定阵列规模,若满足有效分辨最远处两个横向间距为D的目标,则接收阵列的横向最小单元数为:
Figure FDA0003369223480000012
其中,γ为透镜焦距F与目标最远轴向距离Rmax的比值,即:γ=F/Rmax,γ<Rmin/Rmax,其中Rmin为最小检测距离,虚拟透镜的焦距值为:
F=γRmax
阵列形状为矩形、圆形或椭圆形,阵列为均匀阵列或非均匀阵列;
步骤四、确定FFT运算长度;
当探测远界为Rmax,距离分辨率为ΔR时,对每个距离单元进行频率抽取,选择最小FFT运算长度Nfft为:
Figure FDA0003369223480000021
步骤五、确定***采样频率,最小采样频率为:
Figure FDA0003369223480000022
其中,v为目标的最大轴向运动速度;
步骤六、确定调频连续波照射源的调频斜率,最小的线性调频斜率为:
Figure FDA0003369223480000023
步骤七、提供阵列3-D成像检测***;
步骤八、提供数字信号处理***;
步骤九、计算不同距离切面上的像,采用数字成像算法计算不同距离切面上的像,形成一系列成像切片;
步骤十、精确测距;
其中,所述计算不同距离切面上的像,采用数字成像算法计算不同距离切面上的像,形成一系列成像切片包括:
分步骤一:计算天线单元相移,建立空间直角坐标系,建立天线阵的法线方向为-z轴的空间直角坐标系,坐标原点位于天线阵中心,每个天线单元的坐标为(Xmn,Ymn,0),其中m,n为天线单元的横向与纵向编号,计算天线单元对应虚拟透镜的移相量:
Figure FDA0003369223480000024
其中RL为天线阵列的有效半径,dmn为天线单元距离天线阵中心的距离,其计算公式为:
Figure FDA0003369223480000025
分步骤二:计算通道信号的频谱,对每个通道的时域信号进行快速傅里叶变换,求出其频谱分布,
Figure FDA0003369223480000031
分步骤三:恒虚警检测与测频,选取一个通道的FFT数据进行目标检测,检测门限采用单元平均法进行计算:
Figure FDA0003369223480000032
其中,NCFAR不大于Nfft,当频谱幅度超过上述门限时,确认该频率通道存在目标,从而确定目标对应的频率;
分步骤四:频率抽取,根据所述分步骤三确定存在目标的频率通道,依次串行抽取相应的频率通道FFT数据送后级进行处理;
分步骤五:计算目标距离并离散化目标所在球面区域,根据频率抽取对应的频率值,计算目标的距离:
最大频偏为fs/2,对应着目标最远距离,当目标出现在频率为f的位置时,对应的目标距离为:
Figure FDA0003369223480000033
目标在半径为R的球面上,对此球面区域按p行q列进行剖分和离散化,剖分单元的立体角应满足关系式:
Figure FDA0003369223480000034
进而求出离散化球面区域的坐标值(xu,yu,zu);
分步骤六:计算像场区域并离散化,将目标所在球面区域映射到像场区域,获得离散化的像场区域坐标,成像区域所在的位置坐标为:
Figure FDA0003369223480000035
Figure FDA0003369223480000041
Figure FDA0003369223480000042
从而获得成像区域每个剖分节点的坐标(xpq,ypq,zpq),其中p、q为剖分节点的横向与纵向编号;
分步骤七:计算传播相移,计算天线单元到像场网格节点的传播路程:
Figure FDA0003369223480000043
计算传播相移:
Figure FDA0003369223480000044
分步骤八:计算像场值,像场计算公式如下:
Figure FDA0003369223480000045
其中,Pmn为抽取出的频率为f的频谱数据,Γmn,pq为传播系数,其计算公式如下:
Figure FDA0003369223480000046
在简化条件下可直接选取Γmn,pq=1;
分步骤九:像场反演与坐标测量,
将三维立体像进行坐标变换,使得像与真实目标的尺寸、方向一致,坐标为(xv,yv,zv)的像,其反演公式为:
Figure FDA0003369223480000047
Figure FDA0003369223480000048
Figure FDA0003369223480000051
得出目标坐标后求得目标的角度信息;
分步骤十:判断是否需要精细分辨目标,当需要更精细的分辨目标细节时,仅需要减小目标所在球面的离散间隔,执行分步骤五;
分步骤十一:判断程序是否结束,″是″则终止程序,″否″则重新执行分步骤二。
2.根据权利要求1所述的阵列3-D成像方法,其特征在于,NCFAR=Nfft/10。
3.根据权利要求1所述的阵列3-D成像方法,其特征在于,分步骤一、分步骤二之间还包括:
距离选择步骤,当阵列3-D成像检测***采用脉冲体制时,采用不同的距离选择波门对每个通道的时域信号进行距离滤波,利用常规门函数确定距离波门:
Figure FDA0003369223480000052
其中t0为距离波门的中心坐标,τ为距离波门的宽度,经距离波门滤波后的时域信号为:
E(t)=E0(t)·g(t)。
4.根据权利要求1所述的阵列3-D成像方法,其特征在于,当需要更精确的测量目标坐标时,根据精度要求选择更长的FFT运算长度,测距精度为:
Figure FDA0003369223480000053
采用三角形线性调频方案时,分别计算正向和负向调频斜率时的目标频率fu、fd,目标的精确距离为:
Figure FDA0003369223480000061
其中,c=3×108m/s为光速。
5.一种阵列3-D成像检测***,其特征在于,包括照射天线、接收天线阵、发射机、接收机、本地振荡器、A/D、信号处理机、显示***,其中:
所述照射天线用于射频信号的发射;
发射机用于对本振信号进行放大,与照射天线连接;
接收天线阵用于接收目标反射信号;
接收机用于对接收天线阵的射频信号的放大、混频、滤波,将射频信号变频为基带信号或中频信号,与接收天线阵的天线一一对应连接;
本地振荡器用于产生基准射频信号;
A/D用于将模拟信号转换为数字信号,与接收机连接;
数字信号处理***用于信号检测、信号合成与目标成像;
显示***用于三维立体成像的显示以及人机交互,所述阵列3-D成像检测***由如权利要求1-4任一项的阵列3-D成像方法确定。
6.根据权利要求5所述的阵列3-D成像检测***,其特征在于,数字信号处理***包括距离选择模块、FFT模块、CFAR模块、目标检测/测频模块、频率抽取模块、数字透镜成像***、图像处理***,其中:
距离选择模块用于采用不同的距离选择波门对每个通道的时域信号进行距离滤波;
FFT模块用于每个接收通道信号的快速傅里叶变换,将时域信号变换为频域信号;
CFAR模块用于采用单元平均法生成检测门限并进行过门限目标检测;
目标检测/测频模块用于检测不同的目标并根据FFT数据进行频率测量,以确定每个目标对应的频率值;
频率抽取模块用于在测频模块的控制下,抽取出存在目标的频域数据送后级电路或根据用户实际需求按顺序依次抽取频域数据送后级电路;
数字透镜成像***用于对接收到的频域抽取数据进行成像处理,顺序生成不同频率通道的目标的三维像;
图像处理***用于对数字透镜成像***生成的三维像相进行整理、反演变换,转换为人眼可识别的图像。
7.根据权利要求5所述的阵列3-D成像检测***,其特征在于,照射天线采用宽波束天线,接收天线阵包括多个天线单元,发射机包括放大电路、驱动电路和功率放大器,接收机包括低噪声放大器、混频器、滤波器,本地振荡器包括直接数字频率合成器或锁相环;A/D包括模数转换芯片及板卡;信号处理***为搭载嵌入式***的数字板卡。
CN202010168238.0A 2020-03-11 2020-03-11 阵列3-d成像检测***以及成像方法 Active CN111289976B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010168238.0A CN111289976B (zh) 2020-03-11 2020-03-11 阵列3-d成像检测***以及成像方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010168238.0A CN111289976B (zh) 2020-03-11 2020-03-11 阵列3-d成像检测***以及成像方法

Publications (2)

Publication Number Publication Date
CN111289976A CN111289976A (zh) 2020-06-16
CN111289976B true CN111289976B (zh) 2022-02-01

Family

ID=71030994

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010168238.0A Active CN111289976B (zh) 2020-03-11 2020-03-11 阵列3-d成像检测***以及成像方法

Country Status (1)

Country Link
CN (1) CN111289976B (zh)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112612024A (zh) * 2020-12-13 2021-04-06 张艺恒 微波阵列快速成像方法
CN113933834B (zh) * 2021-10-13 2022-07-29 苏州威陌电子信息科技有限公司 圆柱扫描微波成像方法
CN113917461B (zh) * 2021-10-21 2022-10-28 苏州威陌电子信息科技有限公司 一种mimo雷达成像方法及***
CN113917465B (zh) * 2021-10-21 2022-07-26 苏州威陌电子信息科技有限公司 一种sar雷达成像方法及***
CN117687107B (zh) * 2024-01-26 2024-05-07 浙江华视智检科技有限公司 一种毫米波成像方法和相关装置

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103630905A (zh) * 2013-08-29 2014-03-12 中国科学院电子学研究所 阵列天线sar极坐标交叠子孔径成像方法
CN103969837A (zh) * 2014-05-26 2014-08-06 中国科学技术大学 一种高分辨率的集成成像立体显示方法及装置
CN109061638A (zh) * 2018-06-02 2018-12-21 张继龙 相控阵近距离数字成像方法

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8184043B2 (en) * 2010-03-12 2012-05-22 The Boeing Company Super-resolution imaging radar
JP2012168157A (ja) * 2011-02-11 2012-09-06 National Univ Corp Shizuoka Univ 車載用のマルチビーム方式レーダ装置、マルチビーム方式レーダ方法およびマルチビーム方式レーダプログラム

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103630905A (zh) * 2013-08-29 2014-03-12 中国科学院电子学研究所 阵列天线sar极坐标交叠子孔径成像方法
CN103969837A (zh) * 2014-05-26 2014-08-06 中国科学技术大学 一种高分辨率的集成成像立体显示方法及装置
CN109061638A (zh) * 2018-06-02 2018-12-21 张继龙 相控阵近距离数字成像方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
Self-supporting design of a time-encoded aperture, gamma-neutron imaging system;Xiuzuo Liang et al.;《Nuclear Inst. and Methods in Physics Research》;20191019;正文第1-8页 *
微波透镜成像技术在目标识别中的应用研究;张继龙 等;《中国电子科学研究院学报》;20111231;第639-642页 *
针孔式点衍射干涉仪的无镜成像方法;黄磊 等;《光学学报》;20170331;第0312002-1 - 0312002-11页 *

Also Published As

Publication number Publication date
CN111289976A (zh) 2020-06-16

Similar Documents

Publication Publication Date Title
CN111289976B (zh) 阵列3-d成像检测***以及成像方法
CN110988862B (zh) 一种基于极近距离毫米波雷达感知方法及***
Smith et al. Robust through-the-wall radar image classification using a target-model alignment procedure
EP1798570A2 (en) Stand off microwave imaging system and method
CN110794471B (zh) 一种毫米波稀疏阵列远程监视成像方法及***
US11391836B2 (en) Liveliness detection using radar
Sakamoto et al. Fast imaging method for security systems using ultrawideband radar
CN110045367B (zh) 柱面阵列天线目标体三维成像的装置
Soldovieri et al. Sparse tomographic inverse scattering approach for through-the-wall radar imaging
Ford et al. Use of a plane-wave synthesis technique to obtain target RCS from near-field measurements, with selective feature extraction capability
KR20200065827A (ko) 레이더 영상 재구성 기반 객체 추적 장치 및 방법
Bleh W-Band FMCW MIMO radar demonstrator system for 3D imaging.
KR20230090853A (ko) 합성개구레이다 영상으로부터 물체의 지표면 상의 전력 레벨 반사도를 산출하는 방법
RU2524399C1 (ru) Способ обнаружения малоразмерных подвижных объектов
CN109884622B (zh) 柱面阵列天线三维成像的方法
Agarwal et al. Non-invasive concealed weapon detection and identification using V band millimeter wave imaging radar system
Shen et al. Forward scatter shadow ratio: Concept and its application in shadow profile retrieval
CN116520318A (zh) 毫米波成像实时校准方法、装置、计算机设备和存储介质
CN115469306A (zh) 手持式毫米波成像***及成像方法
RU2522853C1 (ru) Способ и устройство обнаружения и идентификации предметов, спрятанных под одеждой на теле человека
KR102080332B1 (ko) 전자파 반사단면적 측정 및 영상화 방법 및 시스템
Hayashi et al. Iterative data clustering algorithm of Doppler-associated RPM imaging for UWB human body imaging radar
DeMartinis et al. A 100 ghz polarimetric compact radar range for scale-model radar cross section measurements
Kong et al. Three-dimensional human imaging for through-the-wall radar
McMakin et al. Millimeter-wave imaging for concealed weapon detection

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
TA01 Transfer of patent application right
TA01 Transfer of patent application right

Effective date of registration: 20210923

Address after: Room 1506, building 1, GANGLONG business building, 99 Xiahe Road, Kunshan Development Zone, Suzhou, Jiangsu 215335

Applicant after: Suzhou Weimo Electronic Information Technology Co.,Ltd.

Address before: Room 501, gate 4, 15 / F, courtyard 11, anningzhuang Road, Haidian District, Beijing 100085

Applicant before: Song Yuhua

GR01 Patent grant
GR01 Patent grant