基于干涉条纹方向的多基线高程干涉相位估计方法与流程

本发明属于无线通信
技术领域：
，具体涉及一种基于干涉条纹方向的多基线高程干涉相位估计方法。
背景技术：
：随着无线通信技术的飞速发展，现代雷达系统也得到了很好的发展。其中，雷达数据处理与雷达信号处理都属于现代雷达系统中的重要组成部分。地形高程干涉相位是雷达系统存在垂直航迹基线时，由高程起伏导致的雷达回波信号到不同雷达接收天线的相位差。干涉条纹是因高程干涉相位存在2π跳变而产生的周期性分布，具有明显的方向性。系统噪声和去相干使得干涉相位图中存在大量的随机相位，影响雷达信号处理的后续流程，因此需要估计无噪声的高程干涉相位。在干涉测高和地面运动目标检测处理中，都需要估计得到无噪声的高程干涉相位信息用于后续处理，高程干涉相位的估计精度严重影响处理性能。目前，高程干涉相位的估计方法主要有两大类，第一类是基于观测场景先验信息，利用数字地形高程模型和雷达系统参数来估计高程干涉相位；第二类是基于雷达观测数据估计高程干涉相位，利用空间平滑去除部分相位噪声。但是，现有的高程干涉相位估计方法存在不同程度的缺陷，其一是数字地形高程模型和系统参数信息的精度影响高程干涉相位的估计精度，其二是空间平滑滤波会损失边缘细节信息。为了提高干涉相位估计精度，降低估计误差，需要一种新的适应干涉条纹方向变化的高程干涉相位估计方法。目前，电子科技大学在其申请的专利文献“一种自适应的高精度干涉sar相位估计方法”(专利申请号201110099514.3，申请公布号cn102270341b)中公开了一种干涉相位估计方法。该方法通过对协方差矩阵特征分解得到信号子空间和噪声子空间，利用music算法中信号子空间和噪声子空间的正交性构造空间谱函数，通过谱峰搜索的方式得到干涉相位的估计结果。该方法利用现代谱估计技术可以提高干涉相位的估计精度。但是，该方法在利用谱峰搜索的方式估计干涉相位时，需要大量的重复计算，不利于工程实现；并且当干涉相位受噪声污染时，通过对协方差矩阵特征分解得到信号子空间和噪声子空间正交性降低，信号子空间向噪声子空间扩散，不适用于噪声干扰背景下稳健估计干涉相位。桂林电子科技大学在其申请的专利文献“一种基于频率补偿的多基线sar干涉相位估计方法”(专利申请号201410251346.9，申请公布号cn104062658a)中公开了一种多基线sar干涉相位估计方法。该方法通过对最短基线干涉相位展开，来提取复干涉图的每一像元随基线变化的粗略估计频率，再将粗估计的频率补偿给抽样的复干涉信号，通过对频率补偿后的信号做傅里叶变换后，利用谱峰搜索的方式得到降频抽样复干涉信号频率，再与粗估计频率叠加得到复干涉图的每一像元随基线变化的频率，进而由最短基线干涉相位估计得到长基线的干涉相位。该方法利用频率补偿技术和基线结构关系可以降低长基线干涉相位的估计误差。但是，该方法存在的不足之处是，对最短基线干涉相位做相位展开处理，干涉相位估计精度受限于相位展开精度，并且当干涉相位受噪声污染时，相位展开误差增大，难以避免误差传递效应，不适用于噪声干扰背景下稳健估计干涉相位。技术实现要素：为了解决现有技术中存在的上述问题，本发明提供了一种基于干涉条纹方向的多基线高程干涉相位估计方法。本发明要解决的技术问题通过以下技术方案实现：一种基于干涉条纹方向的多基线高程干涉相位估计方法，包括：获取多基线雷达复干涉图像数据，以得到张量矩阵；根据所述张量矩阵计算梯度平方张量矩阵，并根据所述梯度平方张量矩阵估计干涉条纹方向；沿所述干涉条纹方向选取样本数据，并对所述样本数据进行处理，得到当前位置索引高程干涉相位的估计值；遍历所述数据张量中的全部位置索引，得到多基线高程干涉相位的估计结果。在本发明的一个实施例中，根据所述张量矩阵计算梯度平方张量矩阵，并根据所述梯度平方张量矩阵估计干涉条纹方向，包括：在所述张量矩阵中一定大小的局部窗口内估计水平方向的梯度向量和垂直方向的梯度向量；根据所述水平方向的梯度向量和所述垂直方向的梯度向量计算梯度平方张量矩阵；对所述梯度平方张量矩阵进行空间平滑，得到平滑后的梯度平方张量矩阵；对所述平滑后的梯度平方张量矩阵进行特征分解，得到所述梯度平方张量矩阵的最小特征值对应的特征张量；根据所述梯度平方张量矩阵的最小特征值对应的特征张量计算干涉条纹方向。在本发明的一个实施例中，所述水平方向的梯度向量和所述垂直方向的梯度向量的计算公式为：其中，gx(m,n)表示水平方向的梯度向量，gy(m,n)表示垂直方向的梯度向量，m表示数据张量的垂直方向位置索引，n表示数据张量的水平方向位置索引，angle()表示取相位角度操作，s表示数据张量，⊙表示哈达玛积，()*表示取共轭操作。在本发明的一个实施例中，所述梯度平方张量矩阵的计算公式为：其中，g(m,n)表示梯度平方张量矩阵，m表示数据张量的垂直方向位置索引，n表示数据张量的水平方向位置索引，()t表示矩阵转置操作，·表示矩阵相乘。在本发明的一个实施例中，所述干涉条纹方向表示为：其中，θ(m,n)表示干涉条纹方向的角度，arctg()表示反正切函数，v2表示梯度平方张量矩阵的最小特征值对应的特征张量，v2(1)表示取向量v2中的第一个元素，v2(2)表示取向量v2中的第二个元素。在本发明的一个实施例中，沿所述干涉条纹方向选取样本数据，并对所述样本数据进行处理，得到当前位置索引高程干涉相位的估计值，包括：以当前位置索引为中心沿干涉条纹方向选取样本数据，得到样本数据矢量集合；根据所述样本数据矢量集合计算样本协方差矩阵；对所述样本协方差矩阵进行特征分解，得到协方差矩阵的最大特征值对应的特征向量；对所述协方差矩阵的最大特征值对应的特征向量取相位操作，得到当前位置索引高程干涉相位的估计值。在本发明的一个实施例中，以当前位置索引为中心沿干涉条纹方向选取样本数据，得到样本数据矢量集合，包括：以当前位置索引为中心点，设置一定大小的角度搜索窗；计算所述角度搜索窗内的每个点相对于所述中心点的空间几何角度；找出所述空间几何角度中与当前位置对应的干涉条纹方向角度最靠近的点作为下一个中心位置；根据所述干涉条纹方向角度计算当前位置的样本数据矢量；重复以上步骤，直到当前索引位置超出预设样本数据搜索范围，得到样本数据矢量集合。在本发明的一个实施例中，所述当前位置的样本数据矢量的计算公式为：其中，表示根据干涉条纹方向角度计算得到的当前位置的样本数据矢量，cos()表示余弦函数，θ(m,n)表示干涉条纹方向的角度，m表示数据张量的垂直方向位置索引，n表示数据张量的水平方向位置索引，φ(k,r)表示空间几何角度，k表示角度搜索窗的垂直位置索引，r表示角度搜索窗的水平位置索引。在本发明的一个实施例中，所述样本协方差矩阵的计算公式为：其中，rh表示包含高程干涉相位关系的样本协方差矩阵，n表示需要求和的样本总数，l表示样本集合中的某一个样本数据矢量的索引序号，表示样本集合中序号为l的样本数据矢量，∈表示数学中的属于符号，ωs表示样本集合。在本发明的一个实施例中，所述当前位置索引高程干涉相位的估计值的表达式为：其中，表示当前位置索引高程干涉相位的估计值，a1表示协方差矩阵的最大特征值对应的特征向量。本发明的有益效果：1、本发明通过采用梯度平方张量矩阵方法估计干涉条纹的方向，使得同一走向但方向相反的梯度矢量不会相互抵消，还可以相互增强，避免了干涉相位梯度在空间邻域平滑过程中相互抵消问题，使得估计干涉条纹方向更具稳健性；克服了现有技术中以干涉相位梯度估计干涉条纹方向时受相位噪声影响导致的估计准确性降低的问题，适用于噪声干扰背景下稳健估计干涉相位；2、本发明采用干涉条纹方向自适应搜索局部范围内走向相同的样本数据矢量来估计样本协方差矩阵，克服了现有技术基于规则结构选取的样本数据带来的边缘细节信息损失问题，且能够适应干涉条纹方向变化，具有更高的估计精度。以下将结合附图及实施例对本发明做进一步详细说明。附图说明图1是本发明实施例提供的一种基于干涉条纹方向的多基线高程干涉相位估计方法的流程示意图；图2是本发明实施例提供的沿干涉条纹方向选取样本数据的示意图；图3是本发明实施例提供的另一种基于干涉条纹方向的多基线高程干涉相位估计方法的流程示意图；图4是本发明实施例提供的多基线含噪声高程干涉相位图；图5是本发明实施例提供的多基线高程干涉估计结果图；图6是本发明实施例提供的多基线高程干涉相位估计误差图。具体实施方式下面结合具体实施例对本发明做进一步详细的描述，但本发明的实施方式不限于此。实施例一请参见图1，图1是本发明实施例提供的一种基于干涉条纹方向的多基线高程干涉相位估计方法的流程示意图，包括以下步骤：步骤1：获取多基线雷达复干涉图像数据，以得到张量矩阵。具体地，从雷达数据文件中读取多基线复干涉图像数据，并将每一组基线相同像素单元的数据排列成数据失量；然后遍历每一像素，构成多基线雷达复干涉图像数据张量矩阵。步骤2：根据张量矩阵计算梯度平方张量矩阵，并根据梯度平方张量矩阵估计干涉条纹方向，具体包括：2a)在张量矩阵中一定大小的局部窗口内估计水平方向的梯度向量和垂直方向的梯度向量。具体地，本实施例在3×3的局部窗内估计水平和垂直方向的梯度向量，估算公式为：其中，gx(m,n)表示水平方向的梯度向量，gy(m,n)表示垂直方向的梯度向量，m表示数据张量的垂直方向位置索引，n表示数据张量的水平方向位置索引，angle()表示取相位角度操作，s表示数据张量，⊙表示哈达玛积，()*表示取共轭操作。2b)根据水平方向的梯度向量和垂直方向的梯度向量计算梯度平方张量矩阵。在本实施例中，梯度平方张量矩阵的计算公式为：其中，g(m,n)表示梯度平方张量矩阵，()t表示矩阵转置操作，·表示矩阵相乘。2c)对梯度平方张量矩阵进行空间平滑，得到平滑后的梯度平方张量矩阵。具体地，本实施例按照下式对梯度平方张量矩阵进行空间平滑：其中，k表示以当前位置索引为中心的k邻域，∑()表示求和符号。2d)对平滑后的梯度平方张量矩阵进行特征分解，得到梯度平方张量矩阵的最小特征值对应的特征张量。具体地，对平滑后的梯度平方张量矩阵进行特征分解的过程可以表示为：其中，v1表示梯度平方张量矩阵的最大特征值对应的特征张量，v2表示梯度平方张量矩阵的最小特征值对应的特征张量，λ1表是梯度平方张量矩阵的最大特征值，λ2表示梯度平方张量矩阵的最小特征值。2e)根据梯度平方张量矩阵的最小特征值对应的特征张量计算干涉条纹方向。在本实施例中，干涉条纹方向用角度表示，具体为：其中，θ(m,n)表示干涉条纹方向角度，arctg()表示反正切函数，v2表示梯度平方张量矩阵的最小特征值对应的特征张量，v2(1)表示取向量v2中的第一个元素，v2(2)表示取向量v2中的第二个元素。进一步地，本实施例还包括重复步骤2a)到2e)，遍历张量矩阵中的全部位置索引，得到对应的干涉条纹方向角度矩阵。步骤3：沿干涉条纹方向选取样本数据，并对样本数据进行处理，得到当前位置索引高程干涉相位的估计值，具体包括：3a)以当前位置索引为中心沿干涉条纹方向选取样本数据，得到样本数据矢量集合。进一步地，首先，以当前位置索引为中心点，设置一定大小的角度搜索窗。在本实施例中，设置以当前位置索引为中心的p×p的范围作为角度搜索窗，并设置q×q的窗作为样本数据搜索范围。其次，计算角度搜索窗内的每个点相对于中心点的空间几何角度。具体地，计算公式如下：其中，φ(k,r)表示空间几何角度，水平向右表示0度，垂直向上表示90度，k表示当前中心位置的垂直位置索引，r表示当前中心位置的水平位置索引。然后，找出所述空间几何角度中与当前位置对应的干涉条纹方向角度最靠近的点作为下一个中心位置。接着，根据干涉条纹方向角度计算当前位置的样本数据矢量，计算公式为：其中，表示根据干涉条纹方向角度计算得到的当前位置的样本数据矢量，cos()表示余弦函数，θ(m,n)表示干涉条纹方向角度，m表示数据张量的垂直方向位置索引，n表示数据张量的水平方向位置索引，φ(k,r)表示空间几何角度，k表示角度搜索窗的垂直位置索引，r表示角度搜索窗的水平位置索引。最后，重复以上步骤，直到当前索引位置超出预设样本数据搜索范围q×q，得到样本数据矢量集合。为进一步说明本实施例提供的沿干涉条纹方向选取样本方法，下面结合附图进行详细说明。请参见图2，图2是本发明实施例提供的沿干涉条纹方向选取样本数据的示意图，其中，每一个方块表示sar图像中的一个像素单元，用数字1标注的方块表示当前中心位置，以该方块为中心的虚线框表示角度搜索窗，此处以大小为3×3的窗为例；箭头表示当前中心位置的干涉条纹方向；用数字2标注的方块表示以上述步骤3a)中确定的下一个中心位置；同理，用数字3标注的方块表示以数字2标注的方块为中心，按上述步骤3a)确定的下一个中心位置；以此类推直到搜索范围超过预设的样本搜索窗，也即样本数据搜索范围。在每一次方向搜索结束后，计算当前位置的样本数据矢量，并将所得样本数据矢量组合成一个样本数据矢量集合，该过程即为沿干涉条纹方向选取样本。3b)根据样本数据矢量集合计算样本协方差矩阵，计算公式如下：其中，rh表示包含高程干涉相位关系的样本协方差矩阵，n表示需要求和的样本总数，l表示样本集合中的某一个样本数据矢量的索引序号，表示样本集合中序号为l的样本数据矢量，∈表示数学中的属于符号，ωs表示样本集合。3c)对样本协方差矩阵进行特征分解，得到协方差矩阵的最大特征值对应的特征向量。具体地，对样本协方差矩阵进行特征分解的过程可以表示为：其中，a1a2…al表示分解之后得到的特征向量，l表示基线数量，λ1λ2…λl表示与特征向量对应的特征值，且λ1＞λ2＞……＞λl，则a1为协方差矩阵的最大特征值λ1对应的特征向量。3d)对协方差矩阵的最大特征值对应的特征向量取相位操作，得到当前位置索引高程干涉相位的估计值，具体表示为：其中，表示当前位置索引高程干涉相位的估计值。步骤4：遍历数据张量中的全部位置索引，重复步骤3，得到多基线高程干涉相位的估计结果。请参见图3，图3是本发明实施例提供的另一种基于干涉条纹方向的多基线高程干涉相位估计方法的流程示意图，本实施例将雷达接收的多通道回波数据排列成张量，利用梯度平方张量矩阵来估计图像中每一点的干涉条纹方向，然后在一定的局部范围内沿干涉条纹方向确定滤波样本，再以所选的样本估计协方差矩阵，最后对样本协方差矩阵做特征值分解，得到的最大特征值对应的特征向量即为多基线高程干涉相位构成的复干涉相位矢量，对其取相位操作得到高程干涉相位的估计值。其中，由于干涉相位存在2π跳变，导致干涉相位梯度在干涉条纹边缘处180度翻转。本实施例通过采用梯度平方张量矩阵方法估计干涉条纹的方向，使得同一走向但方向相反的梯度矢量不会相互抵消，还可以相互增强，避免了干涉相位梯度在空间邻域平滑过程中相互抵消问题，使得估计干涉条纹方向时更具稳健性；克服了现有技术中以干涉相位梯度估计干涉条纹方向时受相位噪声影响导致的估计准确性降低的问题，适用于噪声干扰背景下稳健估计干涉相位。此外，本实施例采用干涉条纹方向自适应搜索局部范围内走向相同的样本数据矢量来估计样本协方差矩阵，克服了现有技术基于规则结构选取的样本数据带来的边缘细节信息损失问题，且能够适应干涉条纹方向变化，具有更高的估计精度。实施例二下面通过仿真实验对本发明的效果进行进一步说明。1、仿真条件本实施例的仿真实验中仿真多基线复干涉sar图像数据，中心频率为9.6吉赫兹，带宽10兆赫兹，接收机采样率15兆赫兹，三组垂直基线高度分别为50米，100米，150米，场景大小为512×512点，地形场景为我国内蒙古南部山区精度为30×30米的dem地形图。2、仿真内容及结果分析仿真内容：本实施例以上述实施例一提供的基于干涉条纹方向的多基线高程干涉相位估计方法分别对多基线构型下的地形高程干涉相位仿真，首先基于理论模型并考虑随机散射和系统去相干等因素仿真得到含噪声的原始sar图像，通过干涉处理得到多基线复干涉图像数据，然后利用上述实施例一的方法估计地形高程干涉相位。结果分析：请参见图4、图5、图6，其中，图4、图5、图6的横坐标均表示图像的方位单元，纵坐标均表示图像的距离单元。具体地，图4是本发明实施例提供的多基线含噪声高程干涉相位图，其从上到下分别对应基线长度50米、100米和150米。从图4可看出高程干涉相位受噪声污染，干涉条纹密集的区域因噪声污染变得模糊，难以从中得到真实的高程干涉相位；图5是本发明实施例提供的多基线高程干涉估计结果图，其从上到下分别对应基线长度50米、100米和150米。从图5可看出高程干涉相位图中毛刺明显减少，高程干涉相位很平滑，干涉条纹结构清晰，证明了本发明方法的有效性。图6是本发明实施例提供的多基线高程干涉相位估计误差图，其从上到下分别对应基线长度50米、100米和150米。从图6可看出本发明方法估计得到的高程干涉相位误差较小，误差图中绝大部分值都接近零，证明了本发明方法的估计精度较高。为了进一步说明本发明方法估计高程干涉相位的有效性，下面对图6中的估计误差做数学统计分析，得到的误差均值和方差如表1所示。表1中的误差均值是指对图6中的高程干涉相位估计误差做统计平均得到的结果，方差是指对图6中的高程干涉相位估计误差做统计方差后得到的结果。表1高程干涉相位估计误差均值和方差统计结果估计误差统计基线长度50米基线长度100米基线长度150米均值0.017-0.0040.031方差0.1150.2160.327从仿真实验结果可以看出，采用本发明方法进行高程干涉相位估计，能够在噪声背景中估计得到高程干涉相位。本发明采用样本协方差矩阵的方式估计多基线高程干涉相位，能够保留干涉相位的空间导向结构；采样基于干涉条纹方向的方式选取样本，可以有效减小处理过程中的干涉相位细节损，提高干涉相位估计精度高。采用本发明方法进行多基线高程干涉相位估计，能够保留多基线之间的空间导向关系，在噪声背景中估计得到高程干涉相位，具有很好的稳健性，为雷达信号处理的后续过程提供有力的技术支持。以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明，不能认定本发明的具体实施只局限于这些说明。对于本发明所属
技术领域：
的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干简单推演或替换，都应当视为属于本发明的保护范围。当前第1页12