技术领域
[0001] 本
发明涉及防灾减灾技术领域，更具体的说是涉及一种地基干涉虚拟孔径形变监测雷达系统及工作方法。
背景技术
[0002] 目前，随着经济与科学技术的不断发展，我国的道路交通建设也进入了高速发展阶段。新建高速公路以及高速
铁路里程数居世界前列，城市轨道交通建设也迎来了飞速发展。隧道是埋置于
地层内的工程
建筑物，在高速公路、高速铁路、城市地铁建设总里程占有很大比重。在隧道施工和运营过程中，围岩的快速收敛或岩土的快速异常形变会导致隧道塌方，造成安全事故，对人们的生命和财产安全构成极大威胁。为研究隧道形变的触发机理并实现快速准确的预测预报，迫切需要有效的隧道监测仪器与设备。为了有效避免隧道坍塌造成伤害，需要对隧道壁的形变进行实时监测，有多种传统监测方法已用于隧道监测领域，例如收敛计、
水准仪、测距仪、全站仪等。虽然这些方法能够测量出隧道监测点处的位移量，但需要耗费大量的人
力物力与时间，无法实现自动化测量。
[0003] 近些年来，众多专家学者和工程技术人员针对隧道形变监测提出了一系列新的技术与方法，主要方法有：(1)光纤应变测量方法：使用黏结剂将光纤黏附在隧道内壁上，基于光的反射及干涉原理，利用光的某种变化量与应变的关系来计算位移形变量。该方法可实现位移测量与通信的一体化设计，但在实际使用中存在
温度与应变交叉敏感问题；另外，光纤非常脆弱，在施工过程中容易造成损坏。(2)微
机电系统传感器方法：该方法将内含
微机电系统的柔性线缆式传感设备安装在隧道壁上，对隧道周边的地层进行形变监测。但是作为一种
接触式测量手段，同样存在安装繁琐、不方便移动的问题。(3)近景摄影测量方法：该方法通过计算机技术和数字影像处理技术将二维数字影像恢复三维坐标，对比不同时刻的三维坐标变化得出结构的形变情况。该方法成本较低，能快速获得形变信息和实时画面，但测量
精度较低，易受环境因素影响。(4)三维激光扫描测量方法：该方法利用激光测距原理，记录被测物体表面大量密集点的三维坐标等信息。它的优点是可以进行多点测量，获取百万点以上的点
云数据绘制物体轮廓。但三维激光
扫描仪价格昂贵，在隧道中大规模使用成本较高；另外，激光易受周围环境因素的影响，隧道中的粉尘和光照会降低系统的测量精度，甚至完全失效。
[0004] 在上述监测测量方法的
基础上，还提出了雷达监测方法，雷达监测方法具有抗环境因素干扰强、可全天候工作、测量精度高和全自动非接触等优点，目前用于形变监测领域的雷达主要有两种类型：(1)地基真实孔径雷达，以澳大利亚GROUND PROBE公司的SSR(Slope Stability Radar)雷达为代表。(2)地基
合成孔径雷达，以意大利IDS公司的IBIS(Image By Interferometric Survey)和Hydra-U(HYper Definition Radar-Underground)雷达为代表。这些雷达系统具有全天候、自动化、高精度、非接触的优点，主要应用于远距离边坡形变监测。这些雷达系统整体造价高昂，体积大而笨重，安装使用不方便，不具有便携性，不适用于狭小的隧道场景。
[0005] 因此，如何实现适应于狭小隧道场景的监测，且安装简单、体积小的雷达监测系统是本领域技术人员亟需解决的问题。
发明内容
[0006] 有鉴于此，本发明提供了一种地基干涉虚拟孔径形变监测雷达系统及工作方法，由三
脚架、二维转台、雷达子系统和数据显示子系统构成，雷达子系统放置在三脚架和二维转台上，雷达子系统实现特定区域扫描，二维转台实现全域扫描。在隧道内壁上安装无源
角发射器。工作方法为，雷达子系统采用分布在矩形四条边上的均匀稀疏线阵向隧道内壁
辐射毫米波频段的宽带调频连续波
信号并接收回波信号，得到初始雷达回
波数据；经过第一次虚拟孔径处理即MIMO处理后，实现孔径扩展；执行距离向傅里叶变换；或者再执行第二次虚拟孔径处理即预测外推处理，进一步实现孔径扩展；经
相位补偿校正和数字波束形成处理后得到三维雷达复图像，从主副图像中获取包含差分相位的干涉图，对干涉图解缠校正后，可得到隧道内壁形变情况并传输到数据显示子系统。
[0007] 为了实现上述目的，本发明采用如下技术方案：
[0008] 一种地基干涉虚拟孔径形变监测雷达系统，所述雷达系统包括三脚架、二维转台、雷达子系统和数据显示子系统；所述二维转台设置于所述三脚架上，所述雷达子系统固定放置于所述二维转台上，实现所述雷达子系统随所述二维转台进行机械转动，所述雷达子系统与所述数据显示子系统无线通讯连接；所述雷达子系统发射雷达信号至所述隧道内壁，并接收所述隧道内壁的监测点反射的回波信号，并生成所述内壁的三维雷达复图像，获取所述内壁的形变信息。
[0009] 优选的，所述雷达子系统包括信号发生单元、控制
电路、若干功率
放大器、发射阵列、接收阵列、若干
低噪声放大器、若干
混频器、若干中频
滤波器、若干
模数转换器、
信号处理电路和无线传输模
块；所述发射阵列连接所述
功率放大器，所述功率放大器连接所述控制电路和所述信号发生单元，所述控制电路连接所述信号发生单元和所述信号处理电路；所述接收阵列连接所述低噪声放大器，所述低噪声放大器连接所述混频器，所述混频器连接所述中频滤波器，所述中频滤波器连接所述模数转换器，所述模数转换器连接所述信号处理电路；所述信号处理电路连接所述无线传输模块；所述信号发生单元连接所述混频器。
[0010] 优选的，所述发射阵列包括若干个发射天线阵元，所述功率放大器有若干个，每一个所述发射天线阵元均连接一个所述功率放大器；所述接收阵列包括若干个接收天线阵元，每一个所述接收天线阵元均连接由一个所述低噪声放大器、所述混频器、所述中频滤波器和所述
数模转换器构成的信号混频传输电路。
[0011] 优选的，所述发射阵列和所述接收阵列采用稀疏线阵收发分置的矩形布阵方法，在所述稀疏线阵的左右两条边上分别布置一个所述发射阵列，在所述稀疏线阵的上下两条边上分别布置一个所述接收阵列；每个所述发射天线阵元由多个微带贴片通过并联
馈线相连构成，贴片等间隔横向排列；每个所述接收天线阵元由多个所述微带贴片通过所述并联馈线相连构成，贴片等间隔纵向排列；每个所述发射天线阵元与每个所述接收天线阵元均采用超低副瓣技术对所述微带贴片构成的阵列进行幅度加权处理。
[0012] 优选的，所述地基干涉虚拟孔径雷达系统还包括无源角
反射器，所述无源角反射器安装在所述隧道内壁上，反射所述雷达系统发射的所述雷达信号。
[0013] 优选的，所述雷达子系统的所述发射阵列以及所述发射天线阵元、所述接收阵列以及所述接收天线阵元的排列参数的确定方法如下：
[0014] 步骤11：根据所述隧道内场景，计算远场条件下的所述稀疏线阵孔径的最大值；
[0015] 步骤12：根据所述稀疏线阵视场角，确定所述发射天线阵元和所述接收天线阵元的波束方向图；确定所述稀疏线阵中每个所述发射天线阵元和所述接收天线阵元的所述微带贴片个数及幅度加权系数，并得到超低副瓣的天线阵元方向图；
[0016] 步骤13：在所述视场角内，从半
波长开始逐渐增加所述稀疏线阵的阵元间距，直到栅瓣位于所述天线阵元方向图的第一副瓣内，从而确定所述稀疏线阵的阵元间距和所述阵元数量。
[0017] 一种地基干涉虚拟孔径形变监测雷达系统，所述的雷达子系统的所述信号发生单元产生毫米波频段的宽带调频连续波信号，通过由所述控制电路控制
开关的所述功率放大器到达所述发射阵列中的所述发射天线阵元；所述控制电路对所述功率放大器的打开和关闭进行控制，令每个时刻所述发射阵列中只有一个所述发射天线阵元发送所述宽带调频连续波信号；所述接收阵列中的所述接收天线阵元接收所述隧道内壁或所述无源角反射器反射的回波信号，并传输至所述低噪声放大器，再经过所述混频器、所述中频滤波器、所述模数转换器后传输至所述信号处理电路，所述信号处理电路对所述回波信号进行处理得到隧道内壁形变信息，并通过所述无线传输模块传输至所述数据显示子系统；所述接收阵列中的所有所述接收天线阵元同时接收所述回波信号；所述控制电路对所述功率放大器的控制，使得所述雷达子系统以时分复用的方式进行工作，保证MIMO信号
正交性。
[0018] 一种地基干涉虚拟孔径形变监测雷达工作方法，在所述信号处理电路中对所述回波信号进行处理的步骤如下：
[0019] 步骤21：对所述回波信号采用数据立方体形式存储，并进行MIMO处理，获得第一次虚拟孔径处理结果，实现第一次孔径扩展；获得方位向数据、
俯仰向数据和时域
采样数据；
[0020] 步骤22：对所述第一次虚拟孔径处理结果进行傅里叶变换，变换到频域，得到的频域采样点序号对应距离信息，实现距离向聚焦，获得距离向数据；
[0021] 步骤23：选取所述内壁到所述雷达子系统的最小距离与最大距离之间对应的频域采样点序号，并结合所述方位向数据、所述俯仰向数据和步骤22得到的数据，在方位向和俯仰向进行线性预测外推处理，进一步实现孔径扩展，得到第二次虚拟孔径处理结果，即获得更多虚拟阵元数据，第二次扩展阵列孔径；
[0022] 步骤24：对所述步骤22得到的所述第一次虚拟孔径处理结果，使用所述二维数字波束形成处理，实现所述方位向和所述俯仰向聚焦，得到高
分辨率的所述三维雷达复图像；或者对所述步骤23得到的所述第二次虚拟孔径处理结果，使用所述二维数字波束形成处理，实现所述方位向和所述俯仰向聚焦，得到更高分辨率的所述三维雷达复图像；
[0023] 步骤25：从所述三维雷达复图像得到雷达相位干涉图像；对所述三维雷达复图像采用球
坐标系表示，每个坐标包括距离、方位角和俯仰角，将不同时刻获取的所述三维雷达复图像与参考时刻雷达复图像结合，得到所述雷达相位干涉图像；
[0024] 步骤26：对干涉相位进行相位解缠与大气校正，得到实际干涉相位
[0025] 步骤27：根据所述实际干涉相位与位移Δr(r，θ，φ)的对应关系：
[0026]
[0027] 得到所述隧道内壁所述监测点的形变信息。
[0028] 优选的，所述步骤23中，利用基于
自回归模型进行线性预测，采用预测滤波器和预测误差滤波器，根据所述预测滤波器和所述预测误差滤波器的误差，调节所述预测滤波器的权值；所述线性预测外推的处理过程为：
[0029] 步骤231：所述步骤21得到虚拟面阵B，阵元个数为M1×N1，所述步骤22得到变换到频域的数据矩阵YB，所述YB的维数为M1×N1，所述YB的第(m，n)元素为yB(m，n)；所述线性预测采用二维线性预测方法，二维矩阵为Y，所述维数为M1×N1，所述Y的所述第(m，n)元素为y(m，n)；
[0030] 设定(K1，K2)阶线性所述预测滤波器和所述预测误差滤波器，所述预测误差滤波器系数为a(k1，k2)，0≤k1≤K1，0≤k2≤K2，a(0，0)＝1；
[0031] 步骤232：将所述预测误差滤波器系数转换成矢量形式：a＝[a(0)，a(1)，…，a(K2)]为行矢量，其中第k个元素a(k)＝[a(0，k)，a(1，k)，…，a(K1，k)]也为行矢量；
[0032] 步骤233：所述预测误差滤波器的M1×N1二维观测数据为y(m，n)，线性预测误差为e(m，n)：
[0033]
[0034] 所述预测误差滤波器输出的所述线性预测误差表示为矢量形式：
[0035] e(m，n)＝az(m，n) (3)
[0036] 其中，z(m，n)＝[y(m，n)，…，y(m-K1，n)，…，y(m，n-K2)，…，y(m-K1，n-K2)]T为列矢量；
[0037] 步骤234：根据所述M1×N1观测数据采用二维协方差方法求解线性预测的最小二乘解，构造总体均方误差函数为：
[0038]
[0039] 其中，Rzz为((K1+1)(K2+1))×((K1+1)(K2+1))大小的协方差矩阵；
[0040] 步骤235：利用所述公式(4)计算得到所述预测误差滤波器的系数a(k1，k2)，采用所述预测滤波器进行预测外推得到预测数据y(m，n)，公式为：
[0041]
[0042] 其中，1≤m≤M1/2+M2/2，1≤n≤N1/2+N2/2；由所述二维矩阵Y和所述预测数据，构成输出矩阵YP，所述维数为(M2/2)×(N2/2)，所述输出矩阵的元素为yp(m，n)＝y(m+M1/2，n+N1/2)，1≤m≤M2/2，1≤n≤N2/2；
[0043] 步骤236：令y(m，n)＝yB(m，n)；通过所述步骤233至所述步骤235进行预测外推，得到预测矩阵YC1为：yc1(m，n)＝yp(m，n)；
[0044] 步骤237：令y(m，n)＝yB(M1+1-m，n)，通过所述步骤233至所述步骤235进行所述预测外推，得到所述预测矩阵YC2为：yc2(m，n)＝yp(M2+1-m，n)；
[0045] 步骤238：令y(m，n)＝yB(M1+1-m，N1+1-n)，通过所述步骤233至所述步骤235进行所述预测外推，得到所述预测矩阵YC3为：yc3(m，n)＝yp(M2+1-m，N2+1-n)；
[0046] 步骤239：令y(m，n)＝yB(m，N1+1-n)，通过所述步骤233至所述步骤235进行所述预测外推，得到所述预测矩阵YC4为：yc4(m，n)＝yp(m，N2+1-n)；
[0047] 步骤2310：所述预测外推后的所述虚拟面阵记为C，阵元个数为M2×N2，对应的所述数据矩阵为YC，YC的所述维数为M2×N2。
[0048]
[0049] 优选的，步骤25的具体实现过程为：
[0050] 步骤251：将所述三维雷达复图像采用球坐标系，每个坐标值用距离、方位角、俯仰角三个值表示，坐标值I(r，θ，φ)为：
[0051]
[0052] 上式中 表示坐标原点到第l个监测点的距离，L为总的监测定个数；w(m，n)(θ，φ)为权系数；c为
电磁波的传播速度；f0为所述宽带调频连续波信号的起始
频率；T为所述宽带调频连续波信号的脉冲宽度，B为信号带宽，α＝B/T为调频率；
[0053] 步骤252：将不同时刻获取的所述三维雷达复图像与所述雷达系统初始获得的参考时刻雷达复图像结合，所述参考时刻雷达复图像为所述雷达系统安装后获得的第一幅所述雷达复图像，是所述对到未发生形变时的图像，得到所述雷达相位干涉图像，所述雷达相位干涉图像中每个
像素点的所述干涉相位值为
[0054]
[0055] 其中*表示复共轭，ang(·)表示求相角运算。
[0056] 优选的，接收阵列接收回波
信号传输至低噪声放大器，经过放大后传输至混频器，在混频器内与信号发生单元产生的宽带调频连续波信号进行混频，得到中频信号，中频信号与安装无源角反射器的监测点距离有关。
[0057] 优选的，所述数据显示子系统为加固
平板电脑终端。
[0058] 经由上述的技术方案可知，与
现有技术相比，本发明公开提供了一种地基干涉虚拟孔径形变监测雷达系统及工作方法，包括安装在隧道内壁上的无源角反射器和，放置于隧道内的雷达系统，雷达系统包括三脚架、设置于三脚架上的二维转台、固定放置于二维转台上的雷达子系统，雷达子系统可随二维转台机械转动，使得发射的雷达波
覆盖隧道内监测范围，雷达子系统中用于发射宽带调频连续波信号的发射阵列和接收无源角反射器反射回波信号的接收阵列构成稀疏线阵，发射天线阵元和接收天线阵元由微带贴片通过并联馈线构成，并采用超低副瓣技术进行幅度加权处理；对于接收的回波信号在雷达系统的信号处理电路中进行第一次虚拟孔径处理即MIMO处理、距离向傅里叶变换、第二次虚拟孔径处理即预测外推处理、数字波束形成和相位干涉测量处理，得到隧道内壁的三维雷达复图像以及相位干涉图，根据相位干涉图得到隧道内安装无源角反射器的监测点的位移情况信息，从而判断隧道内是否发生形变。
附图说明
[0059] 为了更清楚地说明本发明
实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据提供的附图获得其他的附图。
[0060] 图1附图为本发明提供的地基干涉虚拟孔径雷达系统整体结构示意图；
[0061] 图2附图为本发明提供的雷达子系统结构
框图；
[0062] 图3附图为本发明提供的雷达天线阵列结构示意图；
[0063] 图4附图为本发明提供的发射天线阵元结构示意图；
[0064] 图5附图为本发明提供的接收天线阵元结构示意图；
[0065] 图6附图为本发明提供的物理天线阵元和虚拟面阵阵元坐标表示图；
[0066] 图7附图为本发明提供的虚拟阵元线性预测外推示意图；
[0067] 图8附图为本发明提供的三维成像分辨率示意图；
[0068] 图9附图为本发明提供的雷达数据立方体处理示意图；
[0069] 图10附图为本发明提供的地基干涉虚拟孔径雷达系统工作
流程图。
[0070] 1-隧道，11-无源角反射器，21-三脚架，22-二维转台，23-雷达子系统，231-信号发生单元，322-功率放大器，233-发射阵列，2331-发射天线阵元，234-接收阵列，2341-接收天线阵元，235-低噪声放大器，236-混频器，237-中频滤波器，238-模数转换器，239-控制电路，2310-信号处理电路，2311-无线传输模块，24-数据显示子系统。
具体实施方式
[0071] 下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
[0072] 本发明实施例公开了一种地基干涉虚拟孔径形变监测雷达系统及工作方法，雷达系统包括三脚架21、二维转台22、雷达子系统23和数据显示子系统24；二维转台22设置于三脚架21上，雷达子系统23固定放置于二维转台22上，实现雷达子系统23随二维转台22进行机械转动，雷达子系统23与数据显示子系统24无线通讯连接；雷达子系统23发射雷达信号至隧道1内壁，并接收隧道1内壁的监测点反射的回波信号，并生成内壁的三维雷达复图像，获取内壁的形变信息。
[0073] 为了进一步优化上述技术方案，雷达子系统23包括信号发生单元231、控制电路239、若干功率放大器232、发射阵列233、接收阵列234、若干低噪声放大器235、若干混频器
236、若干中频滤波器237、若干模数转换器238、信号处理电路2310和无线传输模块2311；发射阵列233连接功率放大器232，功率放大器232连接控制电路239和信号发生单元231，控制电路239连接信号发生单元231和信号处理电路2310；接收阵列234连接低噪声放大器235，低噪声放大器235连接混频器236，混频器236连接中频滤波器237，中频滤波器237连接模数转换器238，模数转换器238连接信号处理电路2310；信号处理电路2310连接无线传输模块
2311；信号发生单元231连接混频器236。
[0074] 为了进一步优化上述技术方案，发射阵列233包括若干个发射天线阵元2331，功率放大器232有若干个，每一个发射天线阵元2331均连接一个功率放大器232；接收阵列234包括若干个接收天线阵元2341，每一个接收天线阵元2341均连接由一个低噪声放大器235、混频器236、中频滤波器237、
数模转换器238构成的信号混频传输电路。
[0075] 为了进一步优化上述技术方案，发射阵列233和接收阵列234采用稀疏线阵收发分置的矩形布阵方法，在稀疏线阵的左右两条边上分别布置一个发射阵列233，在稀疏线阵的上下两条边上分别布置一个接收阵列234；每个发射天线阵元2331由多个微带贴片通过并联馈线相连构成，贴片等间隔横向排列；每个接收天线阵元2341由多个微带贴片通过并联馈线相连构成，贴片等间隔纵向排列；每个发射天线阵元2331与每个接收天线阵元2341均采用超低副瓣技术对微带贴片构成的阵列进行幅度加权处理。
[0076] 为了进一步优化上述技术方案，地基干涉虚拟孔径雷达系统还包括无源角反射器11，无源角反射器11安装在隧道1内壁上，反射雷达子系统23发射的雷达信号。无源角反射器11的作用是增强反射回波，可在隧道关键区域安装角反射器，提高监测效果。
[0077] 为了进一步优化上述技术方案，数据显示子系统24为加固平板电脑终端。
[0078] 为了进一步优化上述技术方案，雷达子系统23的发射阵列233以及发射天线阵元2331、接收阵列234以及接收天线阵元2341的排列参数的确定方法如下：
[0079] S11：根据隧道1内的场景，计算远场条件下的稀疏线阵的最大值；
[0080] S12：根据稀疏线阵视场角，确定发射天线阵元2331和接收天线阵元2341的波束方向图；确定系数线阵中每个发射天线阵元2331和接收天线阵元2341的微带贴片个数及幅度加权系数，并得到超低副瓣的天线阵元方向图；其中天线阵元的波束宽度可以确定稀疏线阵的视场角范围；
[0081] S13：在视场角内，从半波长开始逐渐增加稀疏线阵的阵元间距，直到栅瓣位于天线阵元方向图的第一副瓣内，从而确定稀疏线阵的阵元间距和阵元数量。
[0082] 一种地基干涉虚拟孔径形变监测雷达系统及工作方法，雷达子系统23的信号发生单元231产生毫米波频段的宽带调频连续波信号，通过由控制电路239控制开关的功率放大器232到达发射阵列233中的发射天线阵元2331；控制电路239对功率放大器232的打开和关闭进行控制，令每个时刻发射阵列233中只有一个发射天线阵元2331发送宽带调频连续波信号；接收阵列234中的接收天线阵元2341接收隧道1内壁或无源角反射器11反射的回波信号，并传输至低噪声放大器235，再经过混频器236、中频滤波器237、模数转换器238后传输至信号处理电路2310，信号处理电路2310对回波信号进行处理得到隧道内壁形变信息，并通过无线传输模块2311传输至数据显示子系统24；接收阵列234中的所有接收天线阵元2341同时接收回波信号；控制电路239对功率放大器232的控制，使得雷达子系统23以时分复用的方式进行工作，保证MIMO信号正交性。
[0083] 为了进一步优化上述技术方案，在信号处理电路2310中对回波信号进行处理的步骤如下：
[0084] S21：对回波信号采用数据立方体形式存储，并进行MIMO处理，获得第一次虚拟孔径处理结果，得到虚拟阵元数据，实现第一次孔径扩展，并且获得方位向数据、俯仰向数据和时域采样数据；
[0085] S22：对时域采样数据进行傅里叶变换，变换到频域，得到的频域采样点序号对应距离信息，实现距离向聚焦，获得距离向数据；
[0086] S23：选取内壁到雷达子系统23的最小距离与最大距离之间对应的频域采样点序号，并结合方位向数据、俯仰向数据和对S22中的距离向数据在方位向和俯仰向进行线性预测外推处理，得到第二次虚拟孔径处理结果，获得更多的虚拟阵元数据，实现第二次孔径扩展；
[0087] 利用基于自回归模型进行线性预测，采用预测滤波器和预测误差滤波器，根据预测滤波器和预测误差滤波器的误差，调节预测滤波器的权值；线性预测外推的处理过程为：
[0088] S231：S21得到虚拟面阵B，阵元个数为M1×N1，S22得到变换到频域的数据矩阵YB，YB的维数为M1×N1，YB的第(m，n)元素为yB(m，n)；线性预测采用二维线性预测方法，二维矩阵为Y，维数为M1×N1，Y的第(m，n)元素为y(m，n)；
[0089] 设定(K1，K2)阶线性预测滤波器和预测误差滤波器，预测误差滤波器系数为a(k1，k2)，0≤k1≤K1，0≤k2≤K2，a(0，0)＝1；
[0090] S232：将预测误差滤波器系数转换成矢量形式：a＝[a(0)，a(1)，…，a(K2)]为行矢量，其中第k个元素a(k)＝[a(0，k)，a(1，k)，…，a(K1，k)]也为行矢量；
[0091] S233：预测误差滤波器的M1×N1二维观测数据为y(m，n)，线性预测误差为e(m，n)：
[0092]
[0093] 预测误差滤波器输出的线性预测误差表示为矢量形式：
[0094] e(m，n)＝az(m，n) (2)
[0095] 其中，z(m，n)＝[y(m，n)，…，y(m-K1，n)，…，y(m，n-K2)，…，y(m-K1，n-K2)]T为列矢量；
[0096] S234：根据M1×N1观测数据采用二维协方差方法求解线性预测的最小二乘解，构造总体均方误差函数为
[0097]
[0098] 其中，Rzz为((K1+1)(K2+1))×((K1+1)(K2+1))大小的协方差矩阵；
[0099] S235：利用公式(4)计算得到预测误差滤波器的系数a(k1，k2)，采用预测滤波器进行预测外推得到预测数据y(m，n)，公式为：
[0100]
[0101] 其中，1≤m≤M1/2+M2/2，1≤n≤N1/2+N2/2；由二维矩阵Y和预测数据，构成输出矩阵YP，维数为(M2/2)×(N2/2)，输出矩阵的元素为yp(m，n)＝y(m+M1/2，n+N1/2)，1≤m≤M2/2，1≤n≤N2/2；
[0102] S236：令y(m，n)＝yB(m，n)；通过S233至S235进行预测外推，得到预测矩阵YC1为：yc1(m，n)＝yp(m，n)；
[0103] S237：令y(m，n)＝yB(M1+1-m，n)，通过S233至S235进行预测外推，得到预测矩阵YC2为：yc2(m，n)＝yp(M2+1-m，n)；
[0104] S238：令y(m，n)＝yB(M1+1-m，N1+1-n)，通过S233至S235进行预测外推，得到预测矩阵YC3为：yc3(m，n)＝yp(M2+1-m，N2+1-n)；
[0105] S239：令y(m，n)＝yB(m，N1+1-n)，通过S233至S235进行预测外推，得到预测矩阵YC4为：yc4(m，n)＝yp(m，N2+1-n)；
[0106] S2310：预测外推后的虚拟面阵记为C，阵元个数为M2×N2，对应的数据矩阵为YC，YC的维数为M2×N2，
[0107]
[0108] S24：对S22或S23中的虚拟阵元数据使用二维数字波束形成处理，实现方位向和俯仰向聚焦，得到高分辨率的三维雷达复图像；
[0109] S25：从三维雷达复图像得到雷达相位干涉图像，获得雷达相位干涉图像的干涉相位；对三维雷达复图像采用球坐标系表示，每个坐标值包括距离、方位角和俯仰角，将不同时刻获取的三维雷达复图像与参考时刻雷达复图像结合，得到雷达相位干涉图像；具体步骤如下：
[0110] S251：将三维雷达复图像采用球坐标系，坐标值I(r，θ，φ)为：
[0111]
[0112] 上式中 表示坐标原点到第l个监测点的距离，w(m，n)(θ，φ)为权系数；c为电磁波的传播速度，f0为宽带调频连续波信号的起始频率，T为宽带调频连续波信号的脉冲宽度，B为信号带宽，α＝B/T为调频率；
[0113] S252：将不同时刻获取的三维雷达复图像与所述雷达系统初始获得的参考时刻雷达复图像结合，得到雷达相位干涉图像，雷达相位干涉图像中每个像素点的值，即干涉相位为
[0114]
[0115] 其中*表示复共轭，ang(·)表示求相角运算；
[0116] S26：对干涉相位进行相位解缠与大气校正，得到实际干涉相位
[0117] S27：根据实际干涉相位和位移Δr(r，θ，φ)的对应关系：
[0118]
[0119] 得到隧道内壁监测点的形变信息。
[0120] 实施例
[0121] 图1为本发明的系统整体结构图。由隧道1、地基干涉虚拟孔径雷达系统2构成；在监测区域的隧道1壁上可安装若干个无源角反射器11作为主要的监测点；地基干涉虚拟孔径雷达系统2由三脚架21、二维转台22、雷达子系统23和数据显示子系统24构成；数据显示子系统24包括加固平板电脑终端。雷达子系统23放置在固定在三脚架21上的二维转台22上，通过机械转动令雷达子系统21的波束范围覆盖隧道1内的所有监测区域。雷达子系统21通过发射阵列233向隧道1发射毫米波频段的宽带调频连续波信号，隧道1内壁及无源角反射器11反射信号，雷达子系统的接收阵列234接收到隧道1内壁以及无源角反射器11反射的回波信号并进行处理，得到雷达复图像和相位干涉图。通过相位信息计算出隧道1内壁的形变量，将数据通过无线传输模块2311传送到数据显示子系统24，技术人员可通过加固平板电脑终端查看监测点形变数据，分析隧道结构的
稳定性。
[0122] 图2为本发明的雷达子系统结构框图。雷达子系统23由信号发生单元231、功率放大器232、发射阵列233、接收阵列234、低噪声放大器235、混频器236、中频滤波器237、模数转换器238、控制电路239、信号处理电路2311和无线传输模块2311构成；发射阵列233由多个发射天线阵元2331组成，每个发射天线阵元又由多个微带贴片构成；接收阵列234由多个接收天线阵元2341组成，每个接收天线阵元又由多个微带贴片构成。信号发生单元231用于产生毫米波频段的宽带调频连续波信号，信号发生单元231产生的宽带调频连续波信号通过由控制电路239控制的功率放大器232到达发射阵列233中的发射天线阵元2331。
[0123] 接收阵列234中的接收天线阵元2341接收隧道1内壁反射的回波信号，依次经过低噪声放大器235、混频器236、中频滤波器237、模数转换器238后送入到信号处理电路2310。信号处理电路2310将得到的初始雷达回波数据进行处理，包括第一次虚拟孔径处理即MIMO处理、距离向傅里叶变换、第二次虚拟孔径处理即预测外推处理、数字波束形成和相位干涉测量等，得到隧道内壁的三维雷达复图像以及相位干涉图。无线传输模块2311将这些数据发送给数据显示子系统24。在雷达工作过程中，控制电路239对功率放大器232的打开和关闭进行控制。使得每个时刻发射阵列233中只有一个发射天线阵元2331发送毫米波频段的宽带调频连续波信号，接收阵列234中的所有接收天线阵元2341同时接收回波信号，采用这种方式使雷达子系统23以时分复用的方式进行工作。
[0124] 雷达子系统23发射信号为宽带调频连续波信号，可由数学公式表示：
[0125]
[0126] 其中t∈(0，T)；f0为初始频率；α＝B/T为调频率，为信号带宽B的比值和脉冲宽度T的比值；θ0为初始相位。
[0127] 雷达子系统使用发射天线照射监测点，使用接收天线接收监测点反射回波信号，假设收发总时延为τ，则接收天线收到的信号为：
[0128]
[0129] 其中τ＝τT+τR＝(rT+rR)/c；c为电磁波传播速度；rT为发射天线到目标的距离；rR为目标到接收天线之间的距离。
[0130] 将回波信号在混频器中与发射信号混频，得到中频信号为：
[0131]
[0132] 由于f0＞＞B，以及T＞＞τ，故上式可以改写为：
[0133] s3(t)＝cos(2πατt+2πf0τ) (12)
[0134] 将上式改写为复数形式表示：
[0135] s4(t)＝exp(j2πατt)·exp(j2πf0τ) (13)
[0136] 进行傅里叶变换可得：
[0137]
[0138] 其中
频谱峰值对应的频率为：fm＝ατ＝α(rT+rR)/c，由此频谱峰值的频率可以确定监测点与雷达子系统发射天线与接收天线之间的距离。中频信号的相位为：θ＝2πf0τ，根据此相位的变化量可计算得到监测点的位移量。
[0139] 图3为本发明的雷达天线阵列结构示意图。本发明采用一种稀疏线阵收发分置的矩形布阵方法，即在矩形的左右两条边上分别布置2个发射阵列233，在矩形的上下两条边上分别布置2个接收阵列234；或者在矩形的左右两条边上分别布置2个接收阵列234，在矩形的上下两条边上分别布置2个发射阵列233，发射阵列233与接收阵列234都是稀疏分布的均匀线阵。
[0140] 图6所示为物理天线阵元和虚拟面阵阵元坐标表示图，假设以坐标原点为参考点，用 表示第m个发射天线阵元的
位置坐标，用 表示第n个接收天线阵元的位置坐标，用表示等效虚拟接收阵元(m，n)的位置坐标，其表达式为：
[0141]
[0142] M发N收的MIMO阵列可以等效为1发MN收的虚拟阵列，虚拟阵列的等效虚拟发射阵元位于坐标原点，等效虚拟接收阵元的坐标为一对发射天线阵元和接收天线阵元的坐标之和。
[0143] 发射天线阵元由多个微带贴片通过并联馈线相连，贴片等间隔横向排列，每个发射天线阵元的结构都相同；接收天线阵元由多个微带贴片通过并联馈线相连，贴片等间隔纵向排列，每个接收天线阵元的结构都相同；相邻微带贴片之间的间隔均为λ/2，λ为波长。
[0144] 在本发明中，发射天线阵元、接收天线阵元采用微带贴片阵列，发射阵列与接收阵列采用稀疏线阵，从而可有效减少天线阵元之间的耦合与串扰，提升雷达系统的工作性能。
[0145] 合成阵列方向图是天线阵元(微带贴片阵列)方向图和稀疏线阵方向图的乘积。稀疏线阵中发射阵列的阵元间隔Dx＞λ/2，，和接收阵列的阵元间隔Dy＞λ/2，，稀疏线阵的方向图会产生栅瓣。天线阵元方向图为每个单元的天线阵元的方向图。合成阵列方向图为等效虚拟面阵的方向图，其孔径越大，波束越窄，为减小栅瓣对合成阵列方向图的影响，在本发明中，采用如下解决方法：天线阵元采用超低副瓣技术，对微带贴片阵列进行幅度加权处理，例如采用切比
雪夫权；改变稀疏线阵的阵元间距和天线阵元的贴片数量，使稀疏线阵的栅瓣位于视场角范围之外，落在天线阵元的副瓣内，从而减少对雷达系统的干扰。
[0146] 阵列排布需要遵循以下规则：
[0147] (1)天线阵元的微带贴片个数为2的整数次幂，相邻贴片的间距均为λ/2。贴片个数决定雷达天线视场角的大小，个数越多，天线方向图主瓣宽度越小，天线视场角也越小。
[0148] (2)稀疏线阵的阵元间距Dx＞λ/2，Dy＞λ/2，间距越大，栅瓣的个数越多。为保证合成阵列方向图中栅瓣对应位置的幅度尽可能小，需要使栅瓣在整个天线视场角范围内均落在天线阵元方向图的副瓣中。
[0149] (3)由于隧道空间狭小，物理天线孔径需要满足远场条件。物理天线孔径由稀疏线阵的阵元间距和阵元数量决定，如果增大阵元间距则阵元数量减少，阵列的栅瓣增加，越容易落在天线阵元方向图主瓣以内，从而造成较大干扰；如果减小阵元间距则阵元数量增加，一方面需要更多的发射接收通道，另一方面收发串扰也会增加。
[0150] 综上，确定阵列参数的过程如下：首先根据隧道场景，计算远场条件下的物理天线(即稀疏线阵)孔径的最大值；然后确定所需的天线阵(即稀疏线阵)的视场角范围，即天线阵元的波束宽度，确定天线阵元的微带贴片个数及幅度加权系数，得到天线阵元方向图；在天线阵的视场角内，从λ/2开始逐渐增加稀疏线阵的阵元间距，直到栅瓣位于天线阵元方向图的第一副瓣内，从而确定稀疏线阵的阵元间距、阵元数量和角度分辨率。
[0151] 图7为本发明虚拟阵元线性预测外推示意图。A为原始物理天线阵(2个发射稀疏线阵，2个接收稀疏线阵)，B为将原始物理天线阵通过第一次虚拟孔径处理即MIMO处理后得到的等效虚拟面阵，C为将先前得到的等效虚拟面阵通过第二次虚拟孔径处理即预测外推得到的等效虚拟面阵，进一步扩展了虚拟阵元个数和天线孔径。
[0152] 本发明提出一种二维线性预测外推方法，来扩展二维面阵。该方法基于自回归模型的线性预测技术，根据预测滤波器和预测误差滤波器的误差，计算预测滤波器的权值，实现对虚拟阵元快拍数据的预测和估计。空间中等间隔分布的各个阵元可以对应于线性预测滤波器中的各个延迟
节点，各节点的时间延迟对应于
空域中两相邻阵元之间的间距。
[0153] 第一次虚拟孔径处理后，得到虚拟面阵B，阵元个数为M1×N1，记变换到频域的数据矩阵为YB，YB的维数为M1×N1；记YB的第(m，n)个元素为yB(m，n)。将虚拟面阵B向4个方向(右上，左上，左下，右下)扩展，得到第二次虚拟孔径处理结果。第二次虚拟孔径处理的过程如下：
[0154] 步骤1，确定二维线性预测方法：假设二维数据矩阵为Y，维数为M1×N1，第(m，n)元素为y(m，n)。定义(K1，K2)阶线性预测滤波器和预测误差滤波器，滤波器系数为a(k1，k2)，0≤k1≤K1，0≤k2≤K2，a(0，0)＝1。
[0155] 将预测误差滤波器系数写成矢量形式：a＝[a(0)，a(1)，…，a(K2)]为行矢量，其中a(k)＝[a(0，k)，a(1，k)，…，a(K1，k)]为行矢量。
[0156] 定义线性预测误差e(m，n)：
[0157]
[0158] 写成矢量形式，上式变为：
[0159] e(m，n)＝az(m，n) (17)
[0160] 其中：z(m，n)＝[y(m，n)，…，y(m-K1，n)，…，y(m，n-K2)，…，y(m-K1，n-K2)]T为列矢量。
[0161] 构造总体均方误差函数为：
[0162]
[0163] Rzz为((K1+1)(K2+1))×((K1+1)(K2+1))大小的协方差矩阵。
[0164] 利用公式(18)可计算得到预测误差滤波器的系数a(k1，k2)。使用线性预测滤波器，进行预测外推，公式为：
[0165]
[0166] 其中，1≤m≤M1/2+M2/2，1≤n≤N1/2+N2/2。由原始数据和预测数据，构造一个新的输出矩阵YP，维数为(M2/2)×(N2/2)。矩阵的元素为yp(m，n)＝y(m+M1/2，n+N1/2)，这里的1≤m≤M2/2，1≤n≤N2/2。
[0167] 步骤2：令y(m，n)＝yB(m，n)；使用第一步所述方法进行预测外推，预测后的矩阵YC1为：yc1(m，n)＝yp(m，n)；
[0168] 步骤3：y(m，n)＝yB(M1+1-m，n)；使用第一步所述方法进行预测外推，预测后的矩阵YC2为：yc2(m，n)＝yp(M2+1-m，n)；
[0169] 步骤4：y(m，n)＝yB(M1+1-m，N1+1-n)；使用第一步所述方法进行预测外推，预测后的矩阵yC3为：yc3(m，n)＝yp(M2+1-m，N2+1-n)；
[0170] 步骤5：y(m，n)＝yB(m，N1+1-n)；使用第一步所述方法进行预测外推，预测后的矩阵yC4为：yC4(m，n)＝yp(m，N2+1-n)；
[0171] 步骤6：线性预测外推后的虚拟面阵记为C，阵元个数为M2×N2，对应数据矩阵为YC，YC的维数为M2×N2，可以由下式得到：
[0172]
[0173] 图8为本发明的三维成像分辨率示意图。假设电磁波的传播速度为c，调频连续波信号的起始频率为f0，带宽为B，波长为λ＝c/f0，则距离分辨率为δR＝c/2B。可以得出雷达距离向分辨率取决于信号带宽，带宽越宽，分辨率越高，本发明采用宽带调频连续波信号，雷达距离向分辨率可达到厘米级。假设天线阵的真实孔径为Drx×Dry，目标处于远场时的最小距离为Rmin＝max{2(Drx)2/λ，2(Dry)2/λ}，经过两次虚拟孔径扩展处理后的天线孔径为Dvx×Dvy，方位向孔径扩展倍数kx＝Dvx/Drx，俯仰向孔径扩展倍数为ky＝Dvy/Dry；对应的方位向角分辨率为δθ＝0.886·λ/Dvx＝51°·λ/Dvx，俯仰向角分辨率为δφ＝0.886·λ/Dvy＝51°·λ/Dvy。R1为俯仰分辨单元，R2为方位分辨单元，R3为距离分辨单元。可以得出雷达角分辨率与天线孔径参数有关，提高角分辨率最直接的方法是增加雷达天线尺寸，但是由于雷达体积、重量等物理因素的制约，无法安装大的物理孔径天线。在本发明中，采用两次虚拟孔径处理可以将天线孔径扩展至原物理尺寸的3～4倍，方位向和俯仰向角分辨率可达到1～2°。
[0174] 图9为本发明的雷达数据立方体处理示意图。本发明中采用数据立方体形式存储雷达数据，数据立方体的三个维度分别对应方位向数据、俯仰向数据和时域采样数据，数据个数分别为方位向虚拟接收阵元数量、俯仰向虚拟接收阵元数量和时域采样点数。通过MIMO处理，雷达增加了方位向和俯仰向的虚拟接收阵元数量，实现天线阵列孔径的第一次扩展。
[0175] 本发明采用宽带调频连续波，雷达收发混频后的中频信号频率与距离有关。因此，针对每一个时域采样数据序列，进行一维傅里叶变换后的结果与距离之间存在对应关系，由此可得到距离向的雷达数据。
[0176] 由于雷达到达监测目标区域的距离处于一定范围以内，可以只使用对应距离范围内的数据进行运算，即隧道内部到雷达的最小距离与最大距离之间，对应的频域采样点序号，实现雷达距离向数据的聚焦，由于只对有限采样点数进行运算，所以大大减少了后续信号处理过程中的运算量，提升系统的实时性。
[0177] 对第一次天线孔径扩展后的方位向数据和俯仰向数据进行线性预测外推，通过计算获取更多虚拟阵元数据。虚拟天线阵列在四个方向上得到扩展，增加了方位向和俯仰向的虚拟接收阵元数量，实现天线阵列孔径的第二次扩展。
[0178] 通过虚拟孔径处理技术，实现了孔径扩展，提高了雷达的角度分辨率，同时还减小了雷达信道数量与物理阵元个数，降低了系统的成本。
[0179] 图10为本发明的系统工作流程图。首先，需要在隧道1内安装无源角反射器11，布置地基干涉虚拟孔径雷达系统2。
[0180] 系统布置完成后，启动雷达子系统23，信号发生单元231在控制电路239控制下产生毫米波频段的宽带调频连续波；控制电路239在同一时刻启动一个功率放大器232，在提升信号发射功率的同时使雷达发射阵列233以时分复用的形式发射信号。雷达波束覆盖监测区域，接收阵列234同时接收回波信号，经过低噪声放大器235、混频器236、中频滤波器237和模数转换器238后进入信号处理电路2310。
[0181] 本发明对雷达原始数据即原始回波信号的处理过程如下：
[0182] 1)采用数据立方体形式存储雷达原始数据并进行MIMO处理，数据立方体的三个维度分别对应方位向数据、俯仰向数据和时域采样数据，得到第一次扩展阵列孔径。
[0183] 2)对每个时域采样数据序列进行傅里叶变换，得到距离向数据。根据雷达到达监测区域的距离选择合适数据范围内的数据进行计算，可以有效减少数据量，提升信号处理速度。
[0184] 或者，再次对方位向数据和俯仰向数据进行线性预测外推，通过计算获取更多虚拟阵元数据，使虚拟天线阵列(即虚拟面阵)在四个方向上得到扩展，得到第二次扩展天线孔径。
[0185] 3)一次或两次虚拟孔径处理完成后，使用二维数字波束形成处理实现雷达成像聚焦，获取雷达对监测区域的三维雷达复图像；
[0186] 三维雷达复图像采用球坐标系，每个坐标用距离、方位角、俯仰角三个值表示，坐标值I(r，θ，φ)为：
[0187]
[0188] 上式中 表示坐标原点到第l个监测点的距离，L为总的监测定个数；w(m，n)(θ，φ)为权系数；c为电磁波的传播速度；f0为所述宽带调频连续波信号的起始频率；T为所述宽带调频连续波信号的脉冲宽度，B为信号带宽，α＝B/T为调频率。
[0189] 4)将不同时刻获取的雷达复图像与参考时刻雷达复图像结合，得到雷达相位干涉图像，图像中每个像素点的干涉相位值 为：
[0190]
[0191] 其中*表示复共轭，ang(·)表示求相角运算。
[0192] 5)识别雷达复图像中监测点的坐标值，获取监测点在相位干涉图像中的数据，由于在实际处理中求得的实际干涉相位 为 在[-π，π]区间的主值，因此需要对干涉相位数据进行相位解缠与大气相位校正处理，得到实际干涉相位。
[0193] 6)最终计算出监测点在雷达视线方向的位移量：
[0194]
[0195] 从而获得隧道内壁监测点形变信息。
[0196] 传统的MIMO雷达回波信号为：
[0197]
[0198] 上式中 为第m个发射天线阵元到第l个监测点的时延， 为第l个监测点到第n个接收阵元的时延。为得到目标的高精度图像，需要补偿信号相位项相位中的 与发射天线阵元、接收天线阵元以及目标位置都有关，即存在距离和
方向耦合问题。传统的雷达采用后向投影(Back Projection，BP)成像
算法，将成像区域划分为一系列像素点，求解每一个像素点到不同发射阵元和接收阵元组合的时延相位项，补偿之后相干
叠加即可得到该像素点的聚焦成像，循环遍历所有像素点完成成像。这种串行处理方式导致运算效率低下。
[0199] 在本发明中，采用虚拟孔径处理使的天线阵的物理尺寸很小，使得监测点处于远场范围，到达天线的电磁波视为平面波，其中虚拟孔径处理可以进行一次或两次，式(24)变为：
[0200]
[0201] 上式中kl是波数并满足kl＝2π/λ[sin(θl)cos(φl)，sin(θl)sin(φl)]，上式中的第1项反映监测点与天线阵坐标原点的距离；第2项反映监测点的回波相位，用于位移测量；第3项与监测点来波方向以及发射天线阵元坐标和接收天线阵元坐标有关，与经过虚拟孔径处理得到的虚拟面阵导引矢量相同。由于解决了距离与方向的耦合问题，本发明中将传统BP算法中的延迟相加过程转化为应用二维数字波束形成处理中的数字波束扫描过程，大大简化了算法的运算量。
[0202] 本
说明书中各个实施例采用递进的方式描述，每个实施例重点说明的都是与其他实施例的不同之处，各个实施例之间相同相似部分互相参见即可。对于实施例公开的装置而言，由于其与实施例公开的方法相对应，所以描述的比较简单，相关之处参见方法部分说明即可。
[0203] 对所公开的实施例的上述说明，使本领域专业技术人员能够实现或使用本发明。对这些实施例的多种
修改对本领域的专业技术人员来说将是显而易见的，本文中所定义的一般原理可以在不脱离本发明的精神或范围的情况下，在其它实施例中实现。因此，本发明将不会被限制于本文所示的这些实施例，而是要符合与本文所公开的原理和新颖特点相一致的最宽的范围。
