CN110223250A - 基于单应变换的sar几何校正方法 - Google Patents

Abstract

本发明提出了一种基于单应变换的SAR几何校正方法，用以解决现有技术难以对待校正SAR图像在与成像平面垂直的空间纵深上产生的几何形变进行校正的问题。实现步骤为：获取参考SAR图像和待校正SAR图像；获取参考SAR图像的匹配特征点集合和待校正SAR图像的匹配特征点集合；计算待校正SAR图像的匹配特征点集合到参考SAR图像的匹配特征点集合的单应变换矩阵；获取待校正SAR图像的几何校正结果。本发明通过单应变换对待校正SAR图像进行几何校正，避免了因平台姿态和天线伺服存在抖动给几何校正精度带来的影响，有效地提高了几何校正精度。可用于条带SAR模式下SAR图像的拼接与融合。

Description

基于单应变换的SAR几何校正方法

技术领域

本发明属于雷达技术领域，涉及一种SAR几何校正方法，具体涉及一种基于单应变换的SAR几何校正方法，可用于条带SAR模式下SAR图像的拼接与融合。

背景技术

SAR作为一种主动式对地观测系统，具有全天候、全天时、高分辨率的成像特点，能够在不同极化方式、不同波段以及不同视角条件下获得感兴趣地面区域的观测数据。近年来，它被广泛地应用到国防、经济、农业以及环境监测等众多领域中。随着SAR技术的不断发展，特别是在条带SAR模式下，通常需对每个合成孔径中录取的数据分别进行成像处理，再对所得到的多幅SAR图像进行综合分析。然而，所观测区域的地形起伏、载机飞行状态的不稳定等因素，往往会造成SAR图像存在一定的几何形变，从而对SAR图像拼接、融合等后续处理形成障碍。

目前，SAR图像几何校正的方法可分为两类，一类是利用SAR成像几何构型，并结合全球定位系统和数字高程模型，通过重采样的方法将多幅SAR图像在同一个地面参考系下重建，以此来完成几何校正。然而，此类方法在实施时，依赖SAR成像几何构型，而决定SAR成像几何构型的相关参数均由平台携带的惯性测量装置提供，当惯性测量装置的精度较低时，会导致所获取参数的误差较大，进一步导致SAR成像几何构型失配，造成几何校正精度大大降低。因此，另一类以图像处理的方式来进行几何校正的方法应运而生，它的基本思路是选取参考SAR图像和待校正SAR图像，对两者进行特征匹配，并建立两者之间的图像变换模型，再通过匹配特征估计变换模型参数，从而以图像变换的形式完成几何校正。由于整个过程为自适应处理，并不需要SAR成像几何构型的相关参数介入，从而有效地消除了测量误差的影响。然而，此类方法的关键问题在于对待校正SAR图像相对于参考SAR图像的几何形变进行定性定量的描述，进而确立变换模型，估计模型参数，实施图像变换，完成几何校正。因此，描述量的准确与否决定了几何校正精度的好坏。现有技术中通常将几何形变描述为平移量、旋转角和尺度因子，以此构建仿射变换模型，通过仿射变换实现对待校正SAR图像的几何校正。

例如申请公布号为CN106780309A，名称为“一种合成孔径雷达图像拼接方法”的中国专利申请，公开了一种基于仿射变换的SAR几何校正方法，该方法通过对两幅SAR图像进行特征匹配，并建立两者之间的仿射变换模型，通过最小二乘法从匹配特征中计算出仿射变换所用的平移量、旋转角和尺度因子，进一步对待校正图像进行仿射变换，有效地补偿了待校正SAR图像相对于参考SAR图像的平移量、旋转角和尺度因子，实现了SAR图像的几何校正。然而，在实际SAR对地观测时，由于平台姿态在横滚维上的抖动和天线伺服在俯仰维上的抖动，会导致对地观测时下视角的变化，从而造成所获取SAR图像在与成像平面垂直的空间纵深上产生几何形变，而该方法采用的仿射变换模型只能对SAR图像在成像平面上的几何形变进行校正，因此，当SAR图像在与成像平面垂直的空间纵深上存在较大的几何形变时，该方法并不适用。

发明内容

本发明的目的在于克服上述现有技术存在的缺陷，提出了一种基于单应变换的SAR几何校正方法，用以解决现有技术中存在的几何校正精度较低的问题。

为实现上述目的，本发明采取的技术方案包括如下步骤：

(1)获取参考SAR图像I₁和待校正SAR图像I₂：

利用SAR对部分重叠的两个地面区域进行拍摄，并将对一个区域拍摄的图像作为参考SAR图像I₁，将对另一个区域拍摄的图像作为待校正SAR图像I₂；

(2)获取参考SAR图像I₁的匹配特征点集合C₁和待校正SAR图像I₂的匹配特征点集合C₂：

(2a)提取参考SAR图像I₁的N₁个特征点和每个特征点对应的特征描述向量，N₁≥4，得到参考SAR图像I₁的特征点集合Q₁和特征描述向量集合P₁，Q₁＝{q₁(i)|i＝1,2,…,N₁}，P₁＝{p₁(i)|i＝1,2,…,N₁}，同时提取待校正SAR图像I₂的N₂个特征点和每个特征点对应的特征描述向量，N₂≥4，得到待校正SAR图像I₂的特征点集合Q₂和特征描述向量集合P₂，Q₂＝{q₂(j)|j＝1,2,…,N₂}，P₂＝{p₂(j)|j＝1,2,…,N₂}，其中，q₁(i)为集合Q₁中第i个特征点，p₁(i)为集合P₁中第i个特征描述向量，q₂(j)为集合Q₂中第j个特征点，p₂(j)为集合P₂中第j个特征描述向量；

(2b)对集合P₁中的每一个特征描述向量p₁(i)与集合P₂中的每一个特征描述向量p₂(j)进行匹配，当min(||p₁(i)-p₂(j)||₂)＜D₀时，提取特征描述向量p₁(i)对应的特征点q₁(i)，构成参考SAR图像I₁的匹配特征点集合C₁，C₁＝{c₁(m)|m＝1,2,…,M}，同时提取特征描述向量p₂(j)对应的特征点q₂(j)，构成待校正SAR图像I₂的匹配特征点集合C₂，C₂＝{c₂(m)|m＝1,2,…,M}，其中，||(·)||₂为2范数运算，min(·)为取最小值运算，D₀为筛选阈值，c₁(m)为集合C₁中第m个匹配特征点，c₂(m)为集合C₂中第m个匹配特征点，M为匹配特征点的个数，M≤min(N₁,N₂)；

(3)计算集合C₂到集合C₁的单应变换矩阵H：

(3a)从集合C₂中随机抽取W个匹配特征点构成集合K，K＝{k(w)|w＝1,2,…,W}，同时从集合C₁中抽取对应的W个匹配特征点构成集合K′，K′＝{k′(w)|w＝1,2,…,W}，并计算集合K到K′的单应变换矩阵H′，其中，k(w)为集合K中第w个匹配特征点，k′(w)为集合K′中第w个匹配特征点，W≥4；

(3b)通过H′对集合C₂进行单应变换，得到集合C′₂：

C′₂＝H′·C₂＝{c′₂(m)|c＝1,2,…,M}，

其中，c′₂(m)为集合C′₂中第m个匹配特征点；

(3c)计算集合C′₂中每一个匹配特征点c′₂(m)和集合C₁中的对应匹配特征点c₁(m)的差的2范数||c′₂(m)-c₁(m)||₂，并构成集合E，E＝{e(m)|m＝1,2,…,M}，选取当e(m)＜E₀时集合C₁中的匹配特征点c₁(m)，构成集合同时选取当e(m)＜E₀时集合C₂中的匹配特征点c₂(m)，构成集合其中，e(m)为集合E中第m个2范数，E₀为筛选阈值，为集合中第n个匹配特征点，为集合中第n个匹配特征点，N为e(m)＜E₀时集合E中2范数的个数；

(3d)判断N＞Mσ是否成立，若是，将计算的集合到的单应变换矩阵作为所求单应变换矩阵H的第s次估计值γ(s)，并将N作为与γ(s)对应的第s次判决值ω(s)，执行步骤(3e)，否则，执行步骤(3a)，其中，σ为判决比例；

(3e)判断s≥S是否成立，若是，将前S次得到的估计值γ(s)构成估计值的集合Γ，Γ＝{γ(s)|s＝1,2,…,S}，并将前S次得到的判决值ω(s)构成判决值集合Ω，Ω＝{ω(s)|s＝1,2,…,S}，执行步骤(3f)，否则，执行步骤(3a)，其中，S为判决次数，S≥2；

(3f)在集合Ω选取ω(s)的最大值，并在集合Γ选取与ω(s)对应的γ(s)，作为单应变换矩阵H；

(4)获取待校正SAR图像I₂的几何校正结果：

通过H对待校正SAR图像I₂进行单应变换I′₂＝H·I₂，实现对待校正SAR图像I₂的几何校正，得到校正SAR图像I′₂。

本发明与现有的技术相比，具有以下优点：

本发明通过计算待校正SAR图像的匹配特征点集合到参考SAR图像的匹配特征点集合的单应变换矩阵，并通过该单应变换矩阵，对待校正SAR图像进行单应变换，完成了对待校正SAR图像相对于参考SAR图像几何形变的校正，解决了现有技术难以对待校正SAR图像在与成像平面垂直的空间纵深上产生的几何形变进行校正的问题，避免了因平台姿态和天线伺服存在抖动给几何校正精度带来的影响，有效地提高了几何校正精度。

附图说明

图1是本发明实现流程图；

图2是本发明参考SAR图像和待校正SAR图像；

图3是本发明参考SAR图像和待校正SAR图像的特征匹配结果图；

图4是本发明与现有技术几何校正结果的对比图。

具体实施方式

下面结合附图和具体实施例，对本发明作进一步详细描述。

参照图1，本发明包括如下步骤：

步骤1)获取参考SAR图像I₁和待校正SAR图像I₂：

利用机载Ka波段SAR对部分重叠的两个地面区域进行拍摄，以保证所得两幅区域拍摄的图像具有相同的部分，从而保证特征匹配能够顺利进行，并将对一个区域拍摄的0.3米分辨率图像作为如图2(a)所示的参考SAR图像I₁，将对另一个区域拍摄的0.3米分辨率图像作为如图2(b)所示的待校正SAR图像I₂，两幅SAR图像的重叠部分由白色虚线框中给出。

步骤2)获取参考SAR图像I₁的匹配特征点集合C₁和待校正SAR图像I₂的匹配特征点集合C₂：

(2a)提取参考SAR图像I₁的N₁个特征点和每个特征点对应的特征描述向量，N₁≥4，得到参考SAR图像I₁的特征点集合Q₁和特征描述向量集合P₁，Q₁＝{q₁(i)|i＝1,2,…,N₁}，P₁＝{p₁(i)|i＝1,2,…,N₁}，同时提取待校正SAR图像I₂的N₂个特征点和每个特征点对应的特征描述向量，N₂≥4，得到待校正SAR图像I₂的特征点集合Q₂和特征描述向量集合P₂，Q₂＝{q₂(j)|j＝1,2,…,N₂}，P₂＝{p₂(j)|j＝1,2,…,N₂}，其中，N₁取值为6707，N₂取值为4281，q₁(i)为集合Q₁中第i个特征点，p₁(i)为集合P₁中第i个特征描述向量，q₂(j)为集合Q₂中第j个特征点，p₂(j)为集合P₂中第j个特征描述向量；

本实施例在对特征点和特征描述向量进行提取时采用的SURF方法，其利用海森矩阵行列式值的正负作为判定条件进行特征点提取，利用二维离散小波变换响应在积分图像上进行特征描述向量的提取，相比于SIFT、HOG等方法，SURF方法在保证了所提取特征尺度不变性的前提下，利用积分图像提高运算效率，是一种稳健、高效的图像特征提取方法；

(2b)对集合P₁中的每一个特征描述向量p₁(i)与集合P₂中的每一个特征描述向量p₂(j)进行匹配，当min(||p₁(i)-p₂(j)||₂)＜D₀时，提取特征描述向量p₁(i)对应的特征点q₁(i)，构成参考SAR图像I₁的匹配特征点集合C₁，C₁＝{c₁(m)|m＝1,2,…,M}，同时提取特征描述向量p₂(j)对应的特征点q₂(j)，构成待校正SAR图像I₂的匹配特征点集合C₂，C₂＝{c₂(m)|m＝1,2,…,M}，其中，||(·)||₂为2范数运算，min(·)为取最小值运算，D₀为筛选阈值，取值为4，c₁(m)为集合C₁中第m个匹配特征点，c₂(m)为集合C₂中第m个匹配特征点，M为匹配特征点的个数，取值为61，M≤min(N₁,N₂)；

参考SAR图像I₁的匹配特征点集合C₁和构成待校正SAR图像I₂的匹配特征点集合C₂如图3所示，左半侧为参考SAR图像I₁，右半侧为待校正SAR图像I₂，参考SAR图像I₁的匹配特征点集合C₁中的每个匹配特征点c₁(m)在左半侧图像上标注，待校正SAR图像I₂的匹配特征点集合C₂中的每个匹配特征点c₂(m)在右半侧图像上标注，横跨两幅图像的白色直线段为反映了集合C₁中的每个匹配特征点c₁(m)和集合C₂中的每个匹配特征点c₂(m)的一一对应关系。

步骤3)计算集合C₂到集合C₁的单应变换矩阵H：

(3a)从集合C₂中随机抽取W个匹配特征点构成集合K，K＝{k(w)|w＝1,2,…,W}，同时从集合C₁中抽取对应的W个匹配特征点构成集合K′，K′＝{k′(w)|w＝1,2,…,W}，并计算集合K到K′的单应变换矩阵H′，其中，k(w)为集合K中第w个匹配特征点，k′(w)为集合K′中第w个匹配特征点，W取值为4，W≥4，

为对行向量按行重塑为方阵的运算，(·)^-1为矩阵求逆运算，(·)^T为矩阵转置运算，A＝[A₁,A₂,…,A_w,…,A_W]^T，a＝[a₁,a₂,…,a_w,…,a_W]^T，且：

x_w为集合K中第w个匹配特征点k(w)的横坐标，y_w为集合K中第w个匹配特征点k(w)的纵坐标，x′_w为集合K′中第w个匹配特征点k′(w)的横坐标，y′_w为集合K′中第w个匹配特征点k′(w)的纵坐标；

(3b)通过H′对集合C₂进行单应变换，得到集合C′₂：

C′₂＝H′·C₂＝{c′₂(m)|c＝1,2,…,M}，

其中，c′₂(m)为集合C′₂中第m个匹配特征点；

(3c)计算集合C′₂中每一个匹配特征点c′₂(m)和集合C₁中的对应匹配特征点c₁(m)的差的2范数||c′₂(m)-c₁(m)||₂，并构成集合E，E＝{e(m)|m＝1,2,…,M}，选取当e(m)＜E₀时集合C₁中的匹配特征点c₁(m)，构成集合同时选取当e(m)＜E₀时集合C₂中的匹配特征点c₂(m)，构成集合其中，e(m)为集合E中第m个2范数，E₀为筛选阈值，取值为2，为集合中第n个匹配特征点，为集合中第n个匹配特征点，N为e(m)＜E₀时集合E中2范数的个数；

(3d)判断N＞Mσ是否成立，若是，将计算的集合到的单应变换矩阵作为所求单应变换矩阵H的第s次估计值γ(s)，并将N作为与γ(s)对应的第s次判决值ω(s)，执行(3e)，否则，执行(3a)，其中，σ为判决比例，取值为0.2，

为对行向量按行重塑为方阵的运算，(·)^-1为矩阵求逆运算，(·)^T为矩阵转置运算，B＝[B₁,B₂,…,B_n,…,B_N]^T，b＝[b₁,b₂,…,b_n,…,b_N]^T，且：

x_n为集合中第n个匹配特征点的横坐标，y_n为集合中第n个匹配特征点的纵坐标，x′_n为集合中第n个匹配特征点的横坐标，y′_n为集合中第n个匹配特征点的纵坐标；

(3e)判断s≥S是否成立，若是，将前S次得到的估计值γ(s)构成估计值的集合Γ，Γ＝{γ(s)|s＝1,2,…,S}，并将前S次得到的判决值ω(s)构成判决值集合Ω，Ω＝{ω(s)|s＝1,2,…,S}，执行步骤(3f)，否则，执行步骤(3a)，其中，S为判决次数，取值为64，S≥2；

(3f)在集合Ω选取ω(s)的最大值，并在集合Γ选取与ω(s)对应的γ(s)，作为单应变换矩阵H；

理想情况下，只需任意抽取集合C₂中W个匹配特征点和集合C₁中与之对应的W个匹配特征点，即可唯一得到集合C₂到集合C₁的单应变换矩阵H，然而，实际情况中，SAR图像噪声的存在往往会影响步骤2特征匹配的正确率，造成所获取的集合C₁和集合C₂中的匹配特征点并不一定全部正确，如图3白色虚线框中所示，为一组存在问题的匹配特征点，因此，通过步骤3来降低集合C₁和集合C₂中存在问题的匹配特征点的对几何校正精度的影响。

步骤4)获取待校正SAR图像I₂的几何校正结果：

通过H对待校正SAR图像I₂进行单应变换I′₂＝H·I₂，实现对待校正SAR图像I₂的几何校正，得到校正SAR图像I′₂。

通过六个元素所唯一确定的仿射变换矩阵，只能描述两幅图像之间的平移量、旋转角和尺度因子，即平面上的三种几何形变，而对于非平面的几何形变，往往通过单应变换矩阵进行描述。单应变换矩阵在仿射变换矩阵的基础上，通过增加两个元素，引入了在与图像平面垂直的空间纵深上的几何形变的描述，因此，适用于解决现有技术中难以对待校正SAR图像在与成像平面垂直的空间纵深上产生的几何形变进行校正的问题，避免了因平台姿态和天线伺服存在抖动给几何校正精度带来的影响，有效地提高了几何校正精度。

为了直观地反映几何校正的效果，将校正SAR图像I′₂与参考SAR图像I₁的重叠部分做灰度融合处理，通过融合的精度来体现几何校正的精度，图4(a)给出了本发明的处理结果，对比图4(b)给出的现有技术的处理结果，特别是白色虚线框中的部分，不难看出，相比于现有技术，本发明较好地校正了待校正SAR图像I′₂相对于参考SAR图像I₁的几何形变。

Claims (3)

1.一种基于单应变换的SAR几何校正方法，其特征在于，包括如下步骤：

(1)获取参考SAR图像I₁和待校正SAR图像I₂：

利用SAR对部分重叠的两个地面区域进行拍摄，并将对一个区域拍摄的图像作为参考SAR图像I₁，将对另一个区域拍摄的图像作为待校正SAR图像I₂；

(2)获取参考SAR图像I₁的匹配特征点集合C₁和待校正SAR图像I₂的匹配特征点集合C₂：

(2a)提取参考SAR图像I₁的N₁个特征点和每个特征点对应的特征描述向量，N₁≥4，得到参考SAR图像I₁的特征点集合Q₁和特征描述向量集合P₁，Q₁＝{q₁(i)|i＝1,2,…,N₁}，P₁＝{p₁(i)|i＝1,2,…,N₁}，同时提取待校正SAR图像I₂的N₂个特征点和每个特征点对应的特征描述向量，N₂≥4，得到待校正SAR图像I₂的特征点集合Q₂和特征描述向量集合P₂，Q₂＝{q₂(j)|j＝1,2,…,N₂}，P₂＝{p₂(j)|j＝1,2,…,N₂}，其中，q₁(i)为集合Q₁中第i个特征点，p₁(i)为集合P₁中第i个特征描述向量，q₂(j)为集合Q₂中第j个特征点，p₂(j)为集合P₂中第j个特征描述向量；

(2b)对集合P₁中的每一个特征描述向量p₁(i)与集合P₂中的每一个特征描述向量p₂(j)进行匹配，当min(||p₁(i)-p₂(j)||₂)＜D₀时，提取特征描述向量p₁(i)对应的特征点q₁(i)，构成参考SAR图像I₁的匹配特征点集合C₁，C₁＝{c₁(m)|m＝1,2,…,M}，同时提取特征描述向量p₂(j)对应的特征点q₂(j)，构成待校正SAR图像I₂的匹配特征点集合C₂，C₂＝{c₂(m)|m＝1,2,…,M}，其中，||(·)||₂为2范数运算，min(·)为取最小值运算，D₀为筛选阈值，c₁(m)为集合C₁中第m个匹配特征点，c₂(m)为集合C₂中第m个匹配特征点，M为匹配特征点的个数，M≤min(N₁,N₂)；

(3)计算集合C₂到集合C₁的单应变换矩阵H：

(3a)从集合C₂中随机抽取W个匹配特征点构成集合K，K＝{k(w)|w＝1,2,…,W}，同时从集合C₁中抽取对应的W个匹配特征点构成集合K′，K′＝{k′(w)|w＝1,2,…,W}，并计算集合K到K′的单应变换矩阵H′，其中，k(w)为集合K中第w个匹配特征点，k′(w)为集合K′中第w个匹配特征点，W≥4；

(3b)通过H′对集合C₂进行单应变换，得到集合C′₂：

C′₂＝H′·C₂＝{c′₂(m)|c＝1,2,…,M}，

其中，c′₂(m)为集合C′₂中第m个匹配特征点；

(3c)计算集合C′₂中每一个匹配特征点c′₂(m)和集合C₁中的对应匹配特征点c₁(m)的差的2范数||c′₂(m)-c₁(m)||₂，并构成集合E，E＝{e(m)|m＝1,2,…,M}，选取当e(m)＜E₀时集合C₁中的匹配特征点c₁(m)，构成集合同时选取当e(m)＜E₀时集合C₂中的匹配特征点c₂(m)，构成集合其中，e(m)为集合E中第m个2范数，E₀为筛选阈值，为集合中第n个匹配特征点，为集合中第n个匹配特征点，N为e(m)＜E₀时集合E中2范数的个数；

(3d)判断N＞Mσ是否成立，若是，将计算的集合到的单应变换矩阵作为所求单应变换矩阵H的第s次估计值γ(s)，并将N作为与γ(s)对应的第s次判决值ω(s)，执行步骤(3e)，否则，执行步骤(3a)，其中，σ为判决比例；

(3e)判断s≥S是否成立，若是，将前S次得到的估计值γ(s)构成估计值的集合Γ，Γ＝{γ(s)|s＝1,2,…,S}，并将前S次得到的判决值ω(s)构成判决值集合Ω，Ω＝{ω(s)|s＝1,2,…,S}，执行步骤(3f)，否则，执行步骤(3a)，其中，S为判决次数，S≥2；

(3f)在集合Ω选取ω(s)的最大值，并在集合Γ选取与ω(s)对应的γ(s)，作为单应变换矩阵H；

(4)获取待校正SAR图像I₂的几何校正结果：

通过H对待校正SAR图像I₂进行单应变换I′₂＝H·I₂，实现对待校正SAR图像I₂的几何校正，得到校正SAR图像I′₂。

2.根据权利要求1所述的基于单应变换的SAR几何校正方法，其特征在于，步骤(3a)中所述的集合K到K′的单应变换矩阵H′，其计算公式为：

其中，为对行向量按行重塑为方阵的运算，(·)^-1为矩阵求逆运算，(·)^T为矩阵转置运算，A＝[A₁,A₂,…,A_w,…,A_W]^T，a＝[a₁,a₂,…,a_w,…,a_W]^T，且：

x_w为集合K中第w个匹配特征点k(w)的横坐标，y_w为集合K中第w个匹配特征点k(w)的纵坐标，x′_w为集合K′中第w个匹配特征点k′(w)的横坐标，y′_w为集合K′中第w个匹配特征点k′(w)的纵坐标。

3.根据权利要求1所述的基于单应变换的SAR几何校正方法，其特征在于，步骤(3d)中所述的集合到的单应变换矩阵其计算公式为：

其中，为对行向量按行重塑为方阵的运算，(·)^-1为矩阵求逆运算，(·)^T为矩阵转置运算，B＝[B₁,B₂,…,B_n,…,B_N]^T，b＝[b₁,b₂,…,b_n,…,b_N]^T，且：

x_n为集合中第n个匹配特征点的横坐标，y_n为集合中第n个匹配特征点的纵坐标，x′_n为集合中第n个匹配特征点的横坐标，y′_n为集合中第n个匹配特征点的纵坐标。