CN103197293B - 一种机载干涉sar多路径误差的提取与补偿方法 - Google Patents

Abstract

本发明提供一种机载干涉SAR多路径误差的提取与补偿方法，能够提取并补偿由于多路径效应造成的干涉相位误差。第一步：根据多路径反射原理推导多路径误差的数学模型，得到多路径误差φ_err关于视角θ的函数；第二步：通过干涉定标和高程反演获取该区域的数字高程模型DEM，将数字高程模型DEM的平均值视为该参考区域的真实高程值；第三步：计算出真实视角θ与真实干涉相位φ_true，通过与干涉相位测量值φ_measure比较获得此区域的多路径误差φ_err；第四步：获得多路径误差的模型参数值；第五步：将回归计算出的模型参数值应用于该架次的所有观测数据进行多路径误差的补偿。

Description

一种机载干涉SAR多路径误差的提取与补偿方法

技术领域

本发明属于信号处理技术领域，尤其是一种机载干涉SAR多路径误差的提取与补偿方法。

背景技术

合成孔径雷达（Synthetic Aperture Radar,SAR）是一种主动式二维高分辨成像雷达，其在距离向（波束照射方向）通过发射大时间带宽积的线性调频信号，采用脉冲压缩技术来获取高分辨率，在方位向（平台运动方向，通常与距离向垂直）利用目标和雷达的相对运动形成的轨迹来构成一个合成孔径来取代庞大的阵列实孔径来获取高分辨率。干涉SAR（Interferometric SAR，InSAR）是指利用两部或多部不同位置处的SAR对同一场景进行观测，并通过数据的后处理获取场景的高程信息。由InSAR系统的回波数据获取场景数学高程模型(Digital Elevation Model，DEM)的主要步骤包括运动补偿、成像、配准、滤波、相位解缠、干涉定标、高程反演。

在机载双天线InSAR系统中，部分雷达回波经机身反射后再进入天线，与直接进入天线的回波叠加在一起，此现象称为多路径反射。由于天线位置与载机平台的几何关系不同，以及机身涂层的差异，各机载InSAR系统的多路径反射影响强弱不等，有些系统中的多路径反射影响可以忽略，有些系统中的多路径反射影响能够造成米级甚至十米级的高程震荡误差。

为了消除多路径反射的不利影响，通常采取前期抑制加后期补偿的措施。前期抑制是指合理地协调天线与载机平台的位置关系、喷刷吸波性能好的涂层等，后期补偿是指利用数据处理的方法来补偿雷达回波中由于多路径反射影响引入的干涉相位误差，即多路径误差。InSAR回波数据中包含多种误差，每一种误差都会对干涉相位有所贡献。虽然能够从产生机理出发，建立各种误差的数学模型，但并不能从总相位误差中分离出不同来源的误差并逐一补偿。因此，在处理多路径误差时，要保证其他类型的误差已基本补偿或相对较弱。在InSAR系统中，最常见的干涉相位误差是干涉相位噪声，干涉相位噪声的数值通常较小，而且可以通过滤波的方法补偿。因此，在进行多路径误差补偿前需要进行干涉相位滤波。

发明内容

本发明的目的是解决现有技术中存在的问题，提供一种机载干涉SAR多路径误差的提取与补偿方法，能够提取并补偿由于多路径效应造成的干涉相位误差。

本发明所提出的一种机载干涉SAR多路径误差的提取与补偿方法，包括以下步骤：

第一步：根据多路径反射原理推导多路径误差的数学模型，得到多路径误差φ_err关于视角θ的函数；

${\mathrm{\φ}}_{\mathrm{err}}=\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_{A,k^{\sin }}\left[\frac{2\mathrm{\π}}{\mathrm{\λ}}(x_k\sin \mathrm{\θ}-y_k\cos \mathrm{\θ}-\sqrt{x_k^2+y_k^2})\right]$

$-\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_{B,K}\sin \left\{\frac{2\mathrm{\π}}{\mathrm{\λ}}\left[\begin{array}{c}(x_k-b\cos \mathrm{\α})\sin \mathrm{\θ}-(y_k-b\sin \mathrm{\α})\cos \mathrm{\θ}\\ -\sqrt{{(x_k-b\cos \mathrm{\α})}^2+{(y_k-b\sin \mathrm{\α})}^2}\end{array}\right]\right\}$

式中λ为雷达波长，b为基线长度，α为基线倾角，反射点数目K、反射系数ε_A，k和ε_B，k、反射点位置(x_k,y_k)为模型参数；

第二步：选取地形平坦、介质单一的区域作为参考区域，通过干涉定标和高程反演获取该区域的数字高程模型DEM，将数字高程模型DEM的平均值视为该参考区域的真实高程值；

第三步：根据第二步得到的参考区域的真实高程值，计算出真实视角θ与真实干涉相位φ_true，通过与干涉相位测量值φ_measure比较获得此区域的多路径误差φ_err；

第四步：根据第三步得到的参考区域的多路径误差φ_err与真实视角θ代入第一步中得到的数学模型做回归计算，获得多路径误差的模型参数值；

第五步：将回归计算出的模型参数值应用于该架次的所有观测数据进行多路径误差的补偿。

本发明方法能够提取并补偿由于多路径效应造成的干涉相位误差，即多路径误差，消除多路径误差对DEM反演结果的影响。

附图说明

图1是本发明机载干涉SAR多路径误差的提取与补偿方法的流程图；

图2是机载干涉SAR多路径反射几何关系图；

图3是机身反射点位置分布图；

图4是第一类反射点的多路径反射几何关系；

图5是第二类反射点的多路径反射几何关系；

图6是第三类反射点对多路径反射几何关系；

图7是某一距离向剖面上多路径误差补偿前后的高程比较，其中(a)为含多路径误差的高程值，(b)为补偿多路径误差后的高程值。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚明白，以下结合具体实施例，并参照附图，对本发明进一步详细说明。

如图1所示，本发明所提出的一种机载干涉SAR多路径误差的提取与补偿方法包括以下步骤：

步骤S1，从多路径误差的产生机理出发，推导出其数学模型，得到其是关于视角的函数；

多路径反射几何关系如图所示，当机身周围存在反射点时，天线收到的回波包括直接返回的和经机身反射进入的，其中经机身反射后进入天线的回波是造成多路径误差的根源。机身反射点会改变信号的传播路径和幅度，造成信号的延时和衰减。

天线A、B处的回波S_A、S_B可分别表示为

$S_A=\mathrm{aexp}\{{\mathrm{j\φ}}_{\mathrm{obj}}-j\frac{4\mathrm{\π}}{\mathrm{\λ}}{\mathrm{\ρ}}_A\}+\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{a\ϵ}}_{A,k^{\exp }}\{{\mathrm{j\φ}}_{\mathrm{obj}}-j\frac{2\mathrm{\π}}{\mathrm{\λ}}({\mathrm{\ρ}}_A+{\mathrm{\ρ}}_{m,k}+b_{A,k})\}---\left(1\right)$

$S_B=\mathrm{aexp}\{{\mathrm{j\φ}}_{\mathrm{obj}}-j\frac{2\mathrm{\π}}{\mathrm{\λ}}[{2\mathrm{\ρ}}_A+Q({\mathrm{\ρ}}_B-{\mathrm{\ρ}}_A)\left]\right\}+\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{a\ϵ}}_{B,k}\exp \{{\mathrm{j\φ}}_{\mathrm{obj}}-j\frac{2\mathrm{\π}}{\mathrm{\λ}}({\mathrm{\ρ}}_A+{\mathrm{\ρ}}_{m,k}+b_{B,k})\}---\left(2\right)$

其中，λ为雷达波长；ε_A,k、ε_B,k表示经由第k(k＝1…K)个反射点F_k进入天线A、天线B的回波的衰减系数，即多路径反射系数；ρ_A、ρ_B、ρ_m,k分别表示地物目标与天线A、天线B、反射点F_k之间的距离；b_A,k、b_B,k分别表示反射点F_k与天线天线A、天线B之间的距离；

表示地物目标对雷达回波的复反射系数；Q用于表征雷达的工作模式，标准收发模式时Q＝1，乒乓收发模式时Q＝2。

由式（1）、式（2）可得，天线A与天线B回波信号的干涉结果为

$S_A{S}_B^{*}=a^2\exp \{-j\frac{2\mathrm{Q\π}}{\mathrm{\λ}}({\mathrm{\ρ}}_A-{\mathrm{\ρ}}_B)\}\·\left[\begin{array}{c}1+\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_{A,k}\exp \{-j\frac{2\mathrm{\π}}{\mathrm{\λ}}({\mathrm{\ρ}}_{m,k}+b_{A,k}-{\mathrm{\ρ}}_A)\}\\ +\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_{B,k}\exp \{-j\frac{2\mathrm{\π}}{\mathrm{\λ}}({\mathrm{\ρ}}_B-{\mathrm{\ρ}}_{m,k}-b_{B,k})\}\end{array}\right]$

$+a^2\underset{u=1}{\overset{K}{\mathrm{\Σ}}}\underset{v=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_{A,u}{\mathrm{\ϵ}}_{B,v}\exp \{-j\frac{2\mathrm{\π}}{\mathrm{\λ}}[b_{A,u}-b_{B,v}+(Q-1)({\mathrm{\ρ}}_A-{\mathrm{\ρ}}_B)\left]\right\}$

$\≈a^2\exp \{-j\frac{2\mathrm{Q\π}}{\mathrm{\λ}}({\mathrm{\ρ}}_A-{\mathrm{\ρ}}_B)\}\·\left[\begin{array}{c}1+\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_{A,k}\exp \{-j\frac{2\mathrm{\π}}{\mathrm{\λ}}({\mathrm{\ρ}}_{m,k}+b_{A,k}-{\mathrm{\ρ}}_A)\}\\ +\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_{B,k}\exp \{-j\frac{2\mathrm{\π}}{\mathrm{\λ}}({\mathrm{\ρ}}_B-{\mathrm{\ρ}}_{m,k}-b_{B,k})\}\end{array}\right]---\left(3\right)$

$\≈a^2\exp \{-j\frac{2\mathrm{Q\π}}{\mathrm{\λ}}({\mathrm{\ρ}}_A-{\mathrm{\ρ}}_B)\}\·\exp \left\{\begin{array}{c}-j\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_{A,k}\sin \left[\frac{2\mathrm{\π}}{\mathrm{\λ}}({\mathrm{\ρ}}_{m,k}+b_{A,k}-{\mathrm{\ρ}}_A)\right]\\ -j\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_{B,k}\sin \left[\frac{2\mathrm{\π}}{\mathrm{\λ}}({\mathrm{\ρ}}_B-{\mathrm{\ρ}}_{m,k}-b_{B,k})\right]\end{array}\right\}$

在式（3）中，

表示没有误差时的理想干涉相位，最后一个指数项即为误差项，故多路径反射造成的相位误差可表示为

${\mathrm{\φ}}_{\mathrm{err}}=-\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_{A,k}\sin \left[\frac{2\mathrm{\π}}{\mathrm{\λ}}({\mathrm{\ρ}}_{m,k}+b_{A,k}-{\mathrm{\ρ}}_A)\right]-\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_{B,k}\sin \left[\frac{2\mathrm{\π}}{\mathrm{\λ}}({\mathrm{\ρ}}_B-{\mathrm{\ρ}}_{m,k}-b_{B,k})\right]---\left(4\right)$

在式（4）中，ρ_m,k、ρ_A、ρ_B等变量不互相独立，并且均是随视角变化的量，因此可将诸变量统一为视角，从而获得更加简洁的表达形式。在将式（4）中的诸变量统一为视角的过程中，需要依据反射点的位置，对几何关系进行分析。反射点的位置可分三种情况：水平线以下（F₁点）；水平线以上、基线以下（F₂点）；基线以上（F₃点）。下面针对反射点的位置情况，对多路径反射的几何关系进行具体推导。

（1）反射点位于水平线以下

该类反射点处的多路径反射几何关系如图所示，b为基线长度，α为基线倾角，α_A,k、α_B,k分别表示矢量

的辐角（逆时针方向为正值）。则有α_A,k＜0，

α_B,k＜0， $\∠{\mathrm{QBF}}_k={-\mathrm{\α}}_{B,k}-(\frac{\mathrm{\π}}{2}-\mathrm{\θ})>0,$ 故有

$\left.\begin{array}{c}{\mathrm{\ρ}}_A-{\mathrm{\ρ}}_{m,k}=\left|\mathrm{AP}\right|=b_{A,k}\cos \∠{\mathrm{PAF}}_k=b_{A,k}\sin (\mathrm{\θ}-{\mathrm{\α}}_{A,k})\\ {\mathrm{\ρ}}_B-{\mathrm{\ρ}}_{m,k}=\left|\mathrm{BQ}\right|=b_{B,k}\cos \∠{\mathrm{QBF}}_k=b_{B,k}\sin (\mathrm{\θ}-{\mathrm{\α}}_{B,k})\end{array}\right\}---\left(5\right)$

（2）反射点位于水平线以上、基线以下

该类反射点处的多路径反射几何关系如图所示，则有α_A,k＞0，

α_B，k＜0， $\∠{\mathrm{QBF}}_k=-{\mathrm{\α}}_{B,k}-(\frac{\mathrm{\π}}{2}-\mathrm{\θ})>0,$ 故有

$\left.\begin{array}{c}{\mathrm{\ρ}}_A-{\mathrm{\ρ}}_{m,k}=\left|\mathrm{AP}\right|=b_{A,k}\cos \∠{\mathrm{PAF}}_k=b_{A,k}\sin (\mathrm{\θ}-{\mathrm{\α}}_{A,k})\\ {\mathrm{\ρ}}_B-{\mathrm{\ρ}}_{m,k}=\left|\mathrm{BQ}\right|=b_{B,k}\cos \∠{\mathrm{QBF}}_k=b_{B,k}\sin (\mathrm{\θ}-{\mathrm{\α}}_{B,k})\end{array}\right\}---\left(6\right)$

（3）反射点位于基线以上

该类反射点处的多路径反射几何关系如图所示，则有α_A,k＞0，

$\∠{\mathrm{PAF}}_k=\mathrm{\π}-(\frac{\mathrm{\π}}{2}-\mathrm{\θ})-{\mathrm{\α}}_{A,k}>0,{\mathrm{\α}}_{B,k}>0,\∠{\mathrm{QBF}}_k={\mathrm{\α}}_{B,k}-(\frac{\mathrm{\π}}{2}+\mathrm{\θ})>0,$ 故有

$\left.\begin{array}{c}{\mathrm{\ρ}}_A-{\mathrm{\ρ}}_{m,k}=-\left|\mathrm{AP}\right|={-b}_{A,k}\cos \∠{\mathrm{PAF}}_k=b_{A,k}\sin (\mathrm{\θ}-{\mathrm{\α}}_{A,k})\\ {\mathrm{\ρ}}_B-{\mathrm{\ρ}}_{m,k}=-\left|\mathrm{BQ}\right|={-b}_{B,k}\cos \∠{\mathrm{QBF}}_k=b_{B,k}\sin (\mathrm{\θ}-{\mathrm{\α}}_{B,k})\end{array}\right\}---\left(7\right)$

综上，由式（5）-（7）可知，无论反射点位置如何分布，ρ_A-ρ_m,k和ρ_B-ρ_m,k的表达形式是一致的，即

$\left.\begin{array}{c}{\mathrm{\ρ}}_A-{\mathrm{\ρ}}_{m,k}=b_{A,k}\sin (\mathrm{\θ}-{\mathrm{\α}}_{A,k})\\ {\mathrm{\ρ}}_B-{\mathrm{\ρ}}_{m,k}=b_{B,k}\sin (\mathrm{\θ}-{\mathrm{\α}}_{B,k})\end{array}\right\}---\left(8\right)$

将式（8）代入式（4）可得

${\mathrm{\φ}}_{\mathrm{err}}=\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_{B,k}\sin \left[\frac{2\mathrm{\π}}{\mathrm{\λ}}{b}_{B,k}(1-\sin (\mathrm{\θ}-{\mathrm{\α}}_{B,k}))\right]-\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_{A,k}\sin \left[\frac{2\mathrm{\π}}{\mathrm{\λ}}{b}_{A,k}(1-\sin (\mathrm{\θ}-{\mathrm{\α}}_{A,k}))\right]---\left(9\right)$

在式（9）中，反射点的位置由该点与两天线之间的相对几何关系给出，即由变量b_A,k、α_A,k、b_B,k、α_B,k表示，事实上此四个变量之间有耦合，可以做进一步的推导简化。

以天线A为原点，以水平线右向为X轴正向，竖直向上为Y轴正向，建立平面直角坐标系，设反射点F_k的位置坐标为(x_k,y_k)，则有

$\left.\begin{array}{c}{\sin \mathrm{\α}}_{A,k}=\frac{y_k}{b_{A,k}},{\sin \mathrm{\α}}_{B,k}=\frac{y_k-b\sin \mathrm{\α}}{b_{B,k}},{\cos \mathrm{\α}}_{A,k}=\frac{x_k}{b_{A,k}},{\cos \mathrm{\α}}_{B,k}=\frac{x_k-b\cos \mathrm{\α}}{b_{B,k}}\\ b_{A,k}=\sqrt{{x_k}^2+{y_k}^2},b_{B,k}=\sqrt{{(x_k-b\cos \mathrm{\α})}^2+{(y_k-b\sin \mathrm{\α})}^2}\end{array}\right\}---\left(10\right)$

将（10）代入（9），可得

${\mathrm{\φ}}_{\mathrm{err}}=\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_{A,k}\sin \left[\frac{2\mathrm{\π}}{\mathrm{\λ}}(x_k\sin \mathrm{\θ}-y_k\cos \mathrm{\θ}-\sqrt{x_k^2+y_k^2})\right]$

$-\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_{B,k}\sin \left\{\frac{2\mathrm{\π}}{\mathrm{\λ}}\left[\begin{array}{c}(x_k-b\cos \mathrm{\α})\sin \mathrm{\θ}-(y_k-b\sin \mathrm{\α})\cos \mathrm{\θ}\\ -\sqrt{{(x_k-b\cos \mathrm{\α})}^2+{(y_k-b\sin \mathrm{\α})}^2}\end{array}\right\}\right\}---\left(11\right)$

式（11）即为多路径相位误差的数学模型的最终形式。在模型中，多路径相位误差φ_err是关于视角θ的函数，反射点数目K、反射系数ε_A,k和ε_B,k、反射点位置(x_k,y_k)为模型参数。

步骤S2，选取地形平坦、介质单一的区域作为参考区域，通过干涉定标和高程反演获取该区域的DEM，将其平均值视为其真实高程；

步骤S3，根据第二步得到的参考区域的高程值，计算出真实的视角θ与干涉相位φ_true，通过与相位测量值φ_measure比较获得此区域的多路径相位误差φ_err；

$\mathrm{\θ}=\arccos \left(\frac{H-h}{{\mathrm{\ρ}}_A}\right)---\left(12\right)$

${\mathrm{\φ}}_{\mathrm{true}}=-\frac{2\mathrm{Q\π}}{\mathrm{\λ}}({\mathrm{\ρ}}_A-{\mathrm{\ρ}}_B)$

$=-\frac{2\mathrm{Q\π}}{\mathrm{\λ}}\left[\begin{array}{c}{\mathrm{\ρ}}_A-\sqrt{{\mathrm{\ρ}}_A^2+b^2-{2\mathrm{\ρ}}_A\·b\·\cos (90-\mathrm{\θ}+\mathrm{\α})}\end{array}\right]---\left(13\right)$

$=-\frac{2\mathrm{Q\π}}{\mathrm{\λ}}\left[\begin{array}{c}{\mathrm{\ρ}}_A-\sqrt{{\mathrm{\ρ}}_A^2+b^2-{2\mathrm{\ρ}}_A\·b\·\sin (\mathrm{\θ}-\mathrm{\α})}\end{array}\right]$

步骤S4，将该参考区域的相位误差与视角数据代入数学模型做回归计算，即可获得模型参数值；

步骤S5，将回归计算出的模型参数值应用于该架次的所有观测数据进行多路径误差的补偿，获取补偿多路径误差前后DEM的距离向剖面，经过多路径误差补偿后，DEM的震荡误差大大降低。

以上所述的具体实施例，对本发明的目的、技术方案和有益效果进行了进一步详细说明，所应理解的是，以上所述仅为本发明的具体实施例而已，并不用于限制本发明，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (1)

1.一种机载干涉SAR多路径误差的提取与补偿方法，首先根据多路径反射原理推导多路径误差的数学模型，其特征在于，该方法还包括以下步骤：

第一步：得到多路径误差φ_err关于视角θ的函数；

$\begin{array}{c}{\mathrm{\φ}}_{\mathrm{err}}=\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_{A,k}\sin \left[\frac{2\mathrm{\π}}{\mathrm{\λ}}(x_k\sin \mathrm{\θ}-y_k\cos \mathrm{\θ}-\sqrt{x_k^2+y_k^2})\right]\\ -\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_{B,k}\sin \left\{\begin{array}{c}\frac{2\mathrm{\π}}{\mathrm{\λ}}\left[\begin{array}{c}(x_k-b\cos \mathrm{\α})\sin \mathrm{\θ}-(y_k-b\sin \mathrm{\α})\cos \mathrm{\θ}\\ -\sqrt{{(x_k-b\cos \mathrm{\α})}^2+{(y_k-b\sin \mathrm{\α})}^2}\end{array}\right]\end{array}\right\}\end{array}$

式中λ为雷达波长，b为基线长度，α为基线倾角，反射点数目K、反射系数ε_A,k和ε_B,k、反射点位置(x_k,y_k)为模型参数；

第二步：选取地形平坦、介质单一的区域作为参考区域，通过干涉定标和高程反演获取该区域的数字高程模型DEM，将数字高程模型DEM的平均值视为该参考区域的真实高程值；

第三步：根据第二步得到的参考区域的真实高程值，计算出真实视角θ与真实干涉相位φ_true，通过与干涉相位测量值φ_measure比较获得此区域的多路径误差φ_err；

第四步：根据第三步得到的参考区域的多路径误差φ_err与真实视角θ代入第一步中得到的数学模型做回归计算，获得多路径误差的模型参数值；

第五步：将回归计算出的模型参数值应用于该架次的所有观测数据进行多路径误差的补偿。

Images