一种基于自适应物方面元的最小二乘匹配方法与流程

(x0,y0)；
[0014]
步骤2.2，过所述物方面元主点p(x
c
,y
c
,z
c
)，建立平面p；
[0015]
步骤2.3，分别计算9个参考影像窗口r1(x0,y0)，r2(x0,y0)，
……
，r9(x0,y0)中各像素投影到所述平面p上得到的物点坐标，得到9个初始物方面元，分别记为p1(x,y,z)，p2(x,y,z)，
……
，p9(x,y,z)；
[0016]
步骤2.4，由空间投影共线条件方程分别计算9个初始物方面元投影到搜索影像上的像点坐标，得到9个搜索影像窗口，分别记为s1(x,y)，s2(x,y)，
……
，s9(x,y)；
[0017]
步骤2.5，分别计算影像窗口r1与s1，r2与s2，
……
，r9与s9的归一化互相关系数，确定初始匹配影像窗口；
[0018]
步骤2.6，以所述初始影像窗口为基础，扩展影像窗口尺寸，并计算参考影像窗口与搜索影像窗口间的归一化互相关系数，确定匹配影像窗口和匹配物方面元尺寸；
[0019]
步骤2.7，所述匹配影像窗口内逐像素建立以辐射畸变和几何畸变为未知数的基于影像灰度的误差方程；
[0020]
步骤2.8，计算所述误差方程，通过最小二乘法求解辐射畸变与几何畸变改正值；
[0021]
步骤2.9，畸变参数迭代优化，当所有改正值小于预设改正值阈值时，停止迭代，取匹配物方面元主点投影到参考影像与各搜索影像的投影点为最佳匹配点位；否则，计算改正值改正畸变参数，重复步骤2.7和2.8，并将迭代次数计算值+1，当迭代次数大于预设次数阈值仍未出现所有改正值小于预设改正值阈值时，判断匹配失败。
[0022]
优选的，步骤2.1中，待匹配像素(x
0c
,y
0c
)分别位于9个参考影像窗口的第1行第1列、第1行第(n+1)/2列、第1行第n列、第(n+1)/2行第1列、第(n+1)/2行第(n+1)/2列、第(n+1)/2行第n列、第n行第1列、第n行第(n+1)/2列和第n行第n列。
[0023]
优选的，步骤2.2中，在物空间坐标系中，以(a,b,c)为法向量，过所述物方面元主点p(x
c
,y
c
,z
c
)，建立平面p，法向量(a,b,c)的方向角(α,β)为预设值。
[0024]
优选的，所述步骤2.5中，分别计算影像窗口r1与s1，r2与s2，
……
，r9与s9的归一化互相关系数，取相关系数最接近于1的一对r
i
与s
i
(i＝1,2,
…
,9)作为初始匹配影像窗口，记为r0(x0,y0)；
[0025]
所述步骤2.6中，以初始匹配影像窗口r0(x0,y0)为基础，逐行或列增加影像窗口像素数目，得到待匹配影像窗口，记为r
k
(x0,y0)(k＝1,2,
…
)，并根据影像与物方面元间的投影关系计算搜索影像像点坐标s
k
(x
i
,y
i
)，直到待匹配影像窗口r
k
与搜索影像像点坐标s
k
间灰度值相关系数最接近1，以此时影像窗口的行列数作为匹配物方面元尺寸，得到匹配影像窗口。
[0026]
优选的，所述步骤2.6中，影像窗口的扩展方法如下：
[0027]
r1：向下增加一行并向右增加一列；
[0028]
r2：向下增加一行；
[0029]
r3：向下增加一行并向左增加一列；
[0030]
r4：向右增加一列；
[0031]
r5：下、上方向各增加一行，右、左方向各增加一列；
[0032]
r6：向左增加一列；
[0033]
r7：向上增加一行并向右增加一列；
[0034]
r8：向上增加一行；
[0035]
r9：向上增加一行并向左增加一列。
[0036]
优选的，所述步骤2.7中，所述误差方程表示如下：
[0037][0038]
其中，v为投影误差，g0(x0,y0)、g(x,y)分别为像点(x0,y0)与(x,y)的灰度值，h0、h1为辐射畸变系数，dh0和dh1分别为h0和h1的改正值，(dx,dy)为几何畸变改正值。
[0039]
本发明中提供的一个或多个技术方案，至少具有如下技术效果或优点：
[0040]
在发明中，首先利用初始匹配像对通过空间前方交会与共线方程投影关系得到初始物方三维点坐标p(x
c
,y
c
,z
c
)，作为物方面元主点；然后以物方面元主点为基准，根据物方面元投影到待匹配影像得到的影像窗口中各像素相关性，确定匹配物方面元尺寸和待匹配像素在影像窗口中的位置，根据空间几何投影关系建立误差方程，通过最小二乘迭代优化求得辐射畸变和几何畸变改正值，搜索最佳匹配点位。即本发明利用待匹配像素周围影像像素灰度值的相关性选择一个最相似的匹配影像窗口(影像窗口尺寸不固定，初始匹配像素在匹配影像窗口的位置不固定)建立误差方程，能够提高匹配精度和匹配成功率。
附图说明
[0041]
图1为本发明实施例提供的一种基于自适应物方面元的最小二乘匹配方法的流程图；
[0042]
图2为本发明实施例中主点位置计算相关系数(ncc)对比图；
[0043]
图3为本发明实施例中影像窗口尺寸计算相关系数对比图。
具体实施方式
[0044]
为了更好的理解上述技术方案，下面将结合说明书附图以及具体的实施方式对上述技术方案进行详细的说明。
[0045]
本实施例提供了一种基于自适应物方面元的最小二乘匹配方法，主要包括以下步骤：
[0046]
步骤1，利用初始匹配像对，通过空间前方交会与共线方程投影关系得到初始物方三维点坐标p(x
c
,y
c
,z
c
)，作为物方面元主点。
[0047]
步骤2，以物方面元主点为基准，根据物方面元投影到待匹配影像得到的影像窗口中各像素相关性，确定匹配物方面元尺寸和待匹配像素在影像窗口中的位置，根据空间几何投影关系建立误差方程，通过最小二乘迭代优化求得辐射畸变和几何畸变改正值，搜索最佳匹配点位。
[0048]
下面对本发明做进一步的说明。
[0049]
本发明技术方案可采用计算机软件技术实现自动运行流程。本发明的实施例是对两对城市区域低空遥感影像进行匹配，参照图1，本发明实施例流程步骤如下：
[0050]
步骤1，利用初始匹配像对，通过空间前方交会与共线方程投影关系得到初始物方三维点坐标p(x
c
,y
c
,z
c
)，作为物方面元主点。
[0051]
具体的，首先输入影像相关参数，采用前方交会方法(公式(1))计算匹配点对应的
初始物方三维点坐标p(x
c
,y
c
,z
c
)。
[0052]
x
c
＝x
s0
+λ0u
0c
＝x
s
+λu
c
[0053]
y
c
＝y
s0
+λ0v
0c
＝y
s
+λv
c
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(1)
[0054]
z
c
＝z
s0
+λ0w
0c
＝z
s
+λw
c
[0055]
其中，(x
s0
,y
s0
,z
s0
)和(x
s
,y
s
,z
s
)分别为参考影像和搜索影像的内方位元素，λ0和λ为投影系数，(u
0c
,v
0c
,w
0c
)和(u
c
,v
c
,w
c
)为分别为待匹配像点(x
0c
,y
0c
)和(x
c
,y
c
)的像空间辅助坐标。
[0056]
步骤2，以物方面元主点为基准，根据物方面元投影到待匹配影像得到的影像窗口中各像素相关性，确定匹配物方面元尺寸和待匹配像素在影像窗口中的位置，根据空间几何投影关系建立误差方程，通过最小二乘迭代优化求得辐射畸变和几何畸变改正值，搜索最佳匹配点位。
[0057]
即基于待匹配像素及其周围影像像素的相关性，确定最优匹配窗口，根据空间几何投影关系及影像像素灰度建立误差方程，通过迭代求解辐射畸变和几何畸变改正值对匹配像素坐标进行优化。
[0058]
例如，可设阈值为10-5
，即当求解的每个畸变参数的改正值均小于10-5
时停止迭代，否则，采用求得的改正值改正相关参数继续迭代。记录迭代次数，当迭代次数超过300次时，迭代失败，该点无效。
[0059]
实施例的步骤2中匹配优化过程具体实现步骤如下：
[0060]
步骤2.1，以参考影像上待匹配像点(x
0c
,y
0c
)为基准，从所述参考影像上取出9个待匹配像素周围n
×
n个像素大小的参考影像窗口，分别记为r1(x0,y0)，r2(x0,y0)，
……
，r9(x0,y0)。
[0061]
待匹配像素(x
0c
,y
0c
)分别位于9个参考影像窗口的第1行第1列、第1行第(n+1)/2列、第1行第n列、第(n+1)/2行第1列、第(n+1)/2行第(n+1)/2列、第(n+1)/2行第n列、第n行第1列、第n行第(n+1)/2列和第n行第n列。
[0062]
其中，n通常取值为7或5。例如，以参考影像上待匹配像点(x
0c
,y
0c
)为主点，从参考影像上取出9个主点周围7
×
7个像素大小的影像窗口，分别记为r1(x0,y0)，r2(x0,y0)，
……
，r9(x0,y0)，参照图2，主点(x
0c
,y
0c
)分别位于9个影像窗口的第1行第1列(西北)、第1行第4列(北)、第1行第7列(东北)、第4行第1列(西)、第4行第4列(中心)、第4行第7列(东)、7行第1列(西南)、第7行第4列(南)和第7行第7列(东南)。
[0063]
步骤2.2，过物方面元主点p(x
c
,y
c
,z
c
)在物空间坐标系中建立一张平面p，平面p的法向量为(a,b,c)，那么平面p的方程即为：
[0064]
a(x-x
c
)+b(y-y
c
)+c(z-z
c
)＝0
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(2)
[0065]
其中，法向量(a,b,c)可由其方向角(α,β)(实施例中初值设为(0
°
,90
°
))表示，即：
[0066][0067]
步骤2.3，分别计算9个参考影像窗口r1(x0,y0)，r2(x0,y0)，
……
，r9(x0,y0)中各像素投影到所述平面p上得到的物点坐标，得到9个初始物方面元，分别记为p1(x,y,z)，p2(x,y,z)，
……
，p9(x,y,z)。
[0068]
具体的，联立物方面元方程与投影方程(公式(4))：
[0069][0070]
其中，λ为投影系数；(u0,v0,w0)为r
i
(x0,y0)中的像点在像空间辅助坐标系中的坐标。
[0071]
步骤2.4，由空间投影共线条件方程分别计算9个初始物方面元投影到搜索影像上的像点坐标，得到9个搜索影像窗口，分别记为s1(x,y)，s2(x,y)，
……
，s9(x,y)。
[0072]
具体的，由共线条件方程(公式(5))分别计算物方面元p
i
上物点p
i
(x,y,z)投影到搜索影像上的像点坐标，得到9个搜索影像窗口s
i
(x,y)(i＝1,2,
…
,9)；
[0073][0074]
其中，f为搜索影像主距；为搜索影像旋转矩阵。
[0075]
步骤2.5，分别计算影像窗口r1与s1，r2与s2，
……
，r9与s9对应像素灰度值的归一化互相关系数(公式(6))，确定初始匹配影像窗口。
[0076]
具体的，参照图2，取相关系数最接近于1的一对r
i
与s
i
(i＝1,2,
…
,9)作为初始匹配影像窗口，记为r0(x0,y0)；
[0077][0078]
其中，g0(x0,y0)、g(x,y)分别为参考影像窗口中像点(x0,y0)与搜索影像窗口中像点(x,y)的灰度值。
[0079]
即在初始9个匹配窗口中选择一对最相似的初始影像窗口，确定初始匹配主点在影像窗口中的位置。
[0080]
步骤2.6，以初始影像窗口为基础，扩展影像窗口尺寸，并计算参考影像窗口与搜索影像窗口间的归一化互相关系数，确定匹配影像窗口和匹配物方面元尺寸。
[0081]
具体的，以初始匹配影像窗口r0(x0,y0)为基础，参照图3，逐行或列增加影像窗口像素数目，得到待匹配影像窗口，记为r
k
(x0,y0)(k＝1,2,
…
)，并根据影像与物方面元间的投影关系计算搜索影像像点坐标s
k
(x
i
,y
i
)，同时计算对应像素灰度值的归一化互相关系数(公式(6))，直到待匹配影像窗口r
k
与搜索影像像点坐标s
k
间灰度值相关系数最大(最接近1)，以此时影像窗口的行列数作为匹配物方面元尺寸，得到匹配影像窗口。
[0082]
其中，影像窗口扩展方法如下：
[0083]
r1：向南(下)增加一行并向东(右)增加一列；
[0084]
r2：向南(下)增加一行；
[0085]
r3：向南(下)增加一行并向西(左)增加一列；
[0086]
r4：向东(右)增加一列；
[0087]
r5：南(下)、北(上)方向各增加一行，东(右)、西(左)方向各增加一列；
[0088]
r6：向西(左)增加一列；
[0089]
r7：向北(上)增加一行并向东(右)增加一列；
[0090]
r8：向北(上)增加一行；
[0091]
r9：向北(上)增加一行并向西(左)增加一列。
[0092]
即在以n*n大小影像窗口为基础，逐行/列扩大影像窗口尺寸，知道待匹配影像窗口相似度最高，这样对于弱纹理区域可以采用大尺寸影像窗口匹配增加信息量，对于纹理复杂区域，可以采用小尺寸影像窗口提高匹配效率。
[0093]
步骤2.7，匹配影像窗口内逐像素建立以辐射畸变和几何畸变为未知数的基于影像灰度的误差方程(公式(7))。
[0094][0095]
其中，v为投影误差，g0(x0,y0)、g(x,y)分别为像点(x0,y0)与(x,y)的灰度值，h0、h1为辐射畸变系数，dh0和dh1分别为辐射畸变系数h0和h
1i
的改正值，(dx
i
,dy
i
)分别为像点坐标(x,y)几何畸变改正值。
[0096]
步骤2.8，将公式(1)-(5)带入公式(7)，计算误差方程，通过最小二乘法求解辐射畸变与几何畸变改正值。
[0097]
步骤2.9，畸变参数迭代优化，当所有改正值小于预设改正值阈值时，停止迭代，取匹配物方面元主点投影到参考影像与各搜索影像的投影点为最佳匹配点位；否则，计算改正值改正畸变参数，重复步骤2.7和2.8，并将迭代次数计算值+1，当迭代次数大于预设次数阈值仍未出现所有改正值小于预设改正值阈值时，判断匹配失败。
[0098]
具体实施时可根据精度要求预设改正值阈值和次数阈值，例如，设改正值阈值为10-5
，即当所有改正值小于10-5
时，停止迭代，如果改正值不满足阈值要求则将迭代次数+1并重复步骤2.7-2.8。
[0099]
匹配完成后的将结果记录为txt文档，文档内容包括：初始匹配点三维点坐标、改正后三维点坐标、参考影像像点坐标、搜索影像像点坐标、匹配迭代次数、是否迭代成功等。
[0100]
本发明首先以初始待匹配点为主点中心选取9个固定、小尺寸的影像窗口；然后通过投影关系和相关系数确定主点位置与影像窗口及投影物方面元尺寸，最后基于影像窗口中对应像素灰度值建立误差方程，进行辐射畸变和几何畸变改正参数的迭代优化，搜索最佳匹配点位。相对传统的基于物方面元的最小二乘匹配，本发明自适应地选取有效特征参与匹配，能够降低深度不连续对匹配结果的干扰并提高弱纹理区域匹配成功率，提高遥感影像匹配的精度和稳定性。与此同时，本发明在影像的弱纹理区域通过影像窗口扩展增加了匹配像素数目，提高了弱纹理区域匹配的成功率。
[0101]
以下通过仿真实验来验证本发明的有效性：
[0102]
仿真实验采用两对各两张航向重叠的城市区域(西北大学校区，广东省阳江市城区)真彩色低空遥感影像，并得到准确的内外方位元素，相邻影像间的重叠度均大于80％，标准影像大小分别为3888
×
2592和8206
×
6078。
[0103]
评价指标：对匹配成功率和匹配精度两个方面进行评价。
[0104]
(1)匹配成功率：统计比较采用基于多通道的最小二乘匹配方法和传统的最小二乘匹配方法分别对均匀分布在两张影像重叠部分46337个(西北大学)和165352个(阳江)初始像对匹配优化成功的概率。
[0105]
(2)匹配精度：统计比较匹配点前方交会得到的三维点相对高程精度，相对高程中误差越小匹配精度越高。
[0106]
仿真结果：根据评价指标，仿真实验的实验数据如表1所示：
[0107]
表1：实验结果统计表
[0108][0109]
从表1的实验结果可以看出，对比同一数据影像匹配结果，自适应物方面元的匹配成功率更高，误差更小。这是因为本发明利用影像本身的特征、纹理信息自适应选取相关性更高的影像像素参与匹配，突出了影像中的角点、边缘等特征纹理信息，同时通过面元扩展的方法增加了弱纹理区域影像窗口内的像素数目，初始匹配相关性越好，匹配精度和稳定性也越好，因而匹配结果优于传统的基于灰度影像的最小二乘匹配方法。
[0110]
本发明实施例提供的一种基于自适应物方面元的最小二乘匹配方法至少包括如下技术效果：
[0111]
(1)通过影像灰度相关性自适应确定待匹配像点在影像窗口的位置，突出匹配影像窗口中角点、边缘等特征纹理信息，提高了匹配精度。
[0112]
(2)通过面元扩展的方法增加了弱纹理区域参与匹配像素数目，丰富了影像匹配窗口的信息量，提高了影像匹配的精度和稳定性。
[0113]
最后所应说明的是，以上具体实施方式仅用以说明本发明的技术方案而非限制，尽管参照实例对本发明进行了详细说明，本领域的普通技术人员应当理解，可以对本发明的技术方案进行修改或者等同替换，而不脱离本发明技术方案的精神和范围，其均应涵盖在本发明的权利要求范围当中。