CN104732532A - 一种遥感卫星多光谱图像配准方法 - Google Patents
一种遥感卫星多光谱图像配准方法 Download PDFInfo
- Publication number
- CN104732532A CN104732532A CN201510106747.XA CN201510106747A CN104732532A CN 104732532 A CN104732532 A CN 104732532A CN 201510106747 A CN201510106747 A CN 201510106747A CN 104732532 A CN104732532 A CN 104732532A
- Authority
- CN
- China
- Prior art keywords
- image
- point
- registration
- phi
- matching
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Landscapes
- Image Analysis (AREA)
- Image Processing (AREA)
Abstract
本发明一种遥感卫星多光谱图像配准方法,该方法包括以下步骤:①根据多光谱探测器的距离得到粗配准图像;②基于基准图像对待配准图像进行直方图映射;③分块划分上述两幅图像,针对每对图像块:提取每对图像块的角点、近红外图像的配准还需要提取各块的主方向相位一致特征、基于相位相关和曲面拟合得到各块的正反向亚像素匹配点集;④三角网格划分基准图像和待配准图像,计算变换系数并对待配准图像进行重采样得到配准图像。与现有技术相比,本发明方法对不同反射率地物在不同谱段上表现出的灰度差异大而导致的配准难问题,同时有效解决了大幅面遥感图像不同区域存在不同程度畸变的问题。
Description
技术领域
本发明涉及一种遥感卫星多光谱图像配准方法,属于遥感图像处理技术领域,特别适用于遥感卫星多光谱图像的配准。
背景技术
随着用户对遥感数据需求的增长以及遥感技术的发展,多光谱相机在遥感卫星上得到广泛应用。多光谱相机具有较宽谱段的探测能力,通常从蓝谱段到近红外谱段,一般包括蓝、绿、红、近红外等谱段,例如,CBERS-1的四个谱段分别为0.45-0.52、0.52-0.59、0.63-0.69、0.77-0.89(单位:μm),QuickBird的四个谱段分别为0.45-0.52、0.52-0.60、0.63-0.69、0.76-0.90(单位:μm),SPOT5的四个谱段分别为0.50-0.59、0.61-0.68、0.79-0.89、1.58-1.70(单位:μm)。
通过多光谱多个谱段图像的配准能够实现多光谱图像的融合,合成较高质量的彩色图像,进一步通过对彩色图像的分析,可以为资源普查、灾害监测等提供定量应用的数据。然而,不同谱段图像间的微小位置差异就能导致合成的彩色图像出现“重叠影”,严重影响图像质量,为提高多光谱图像的融合质量,在处理过程中需要各谱段间图像的严格配准。常用的图像配准方法主要包括基于图像灰度与基于图像特征的配准方法,其中,基于图像灰度的方法直接利用图像的灰度信息,建立两幅图像之间的相似性度量,但对于灰度差异大的两幅图像配准,该方法失去价值;基于图像特征的方法,通过提取基准图像和待配准图像的一些共同特征,建立特征之间的对应关系,但大多数的角点、边缘特征提取主要依赖于图像灰度信息,导致基准图像和待配准图像提取的特征很难一致。
由于地面物质对不同波段光的反射率不同,不同波段的图像灰度差别也较大,特别是蓝、绿、红谱段与近红外谱段间的灰度差异更大。而现有的基于灰度与基于特征的配准方法对灰度信息比较敏感,导致配准难度加大,同时由于大气、相机、卫星***等的影响导致不同谱段间存在一定的畸变。
发明内容
本发明解决的技术问题是:克服现有技术的不足,提供一种遥感卫星多光谱图像配准方法,用于解决遥感卫星多光谱谱段间的大幅面图像配准,特别解决了不同反射率地物在近红外谱段与蓝、绿、红谱段图像中表现出灰度差异大而导致的配准难问题。
本发明的技术方案是:一种遥感卫星多光谱图像配准方法,步骤如下:
1)根据多光谱相机蓝、绿、红、近红外四个谱段探测器间的距离,对获取的四个谱段图像分别进行裁剪,得到粗配准的多光谱图像,具体步骤如下:
多光谱相机焦平面上的蓝谱段探测器到绿、红、近红外谱段探测器的距离分别为d1、d2、d3,则从绿、红、近红外谱段图像的顶部分别裁剪δ1、δ2、δ3个像素,从蓝、绿、红谱段图像的底部分别裁剪δ3、δ2、δ1个像素,其中, τ为多光谱相机的单个探测器尺寸,表示下取整,裁剪后的图像即为粗配准的多光谱图像;
2)将粗配准的红谱段图像作为基准图像f(x,y),粗配准的蓝、绿、近红外图像中的任意一幅图像作为待配准图像,建立基准图像到待配准图像的直方图映射关系,以建立的直方图映射关系对待配准图像进行映射,得到映射图像g(x,y),所述f(x,y)与g(x,y)的大小为M×N个像素;M、N为正整数;
3)以大小为m×m个像素的图像块分别划分f(x,y)与g(x,y),得到的每对图像块表示为fi与gi,共得到对图像块,其中,m为正整数,表示上取整;针对每一对图像块fi与gi进行如下步骤:
3.1)分别提取图像块fi与gi的Harris角点,图像块fi的所有角点构成fi的正向角点集;
3.2)当待配准图像为近红外谱段图像时,计算图像块fi的梯度主方向Φi;构造方向为Φi的二维Log-Gabor滤波器对fi与gi进行相位一致特征检测,分别得到fi与gi的主方向相位一致特征PC(fi)与PC(gi);
3.3)将图像块fi的正向角点集中的每个角点作为待匹配点,在图像块gi中基于相位相关法和二元二次曲面拟合计算对应的亚像素匹配点,则得到在图像块gi中的正向亚像素匹配角点集;正向角点集和正向亚像素匹配角点集记为正向匹配角点对集Γi;
3.4)计算从gi提取的每个角点与正向亚像素匹配角点集中所有角点的最小距离,如果最小距离小于设定阈值则剔除gi中的该角点,否则保留gi中的该角点,所有符合上述条件的角点构成gi中的反向角点集;
3.5)将gi中的反向角点集中的每个角点作为待匹配点,在图像块fi中基于相位相关法和二元二次曲面拟合计算对应的亚像素匹配点,则得到在fi中的反向亚像素匹配角点集;反向角点集和反向亚像素匹配角点集构成反向匹配角点对集Γ'i;正反向匹配角点对集表示为Γi+Γ'i;
3.6)将图像块fi划分若干个大小为K×K个像素的块区域,K为正整数;遍历所有块区域,判断块区域内是否有匹配的角点对,即是否有集合Γi+Γ'i中的元素,若没有则将该块区域的中心点作为构造点,所有的构造点构成构造点集;
3.7)将构造点集中的每个构造点作为待匹配点,在图像块gi中基于相位相关法和二元二次曲面拟合计算对应的亚像素匹配点,则得到在gi中的亚像素匹配点集;构造点集及其亚像素匹配点集构成匹配构造点对集Γ”i;完整的匹配点集合记为Γi+Γ'i+Γ”i;
3.8)基于RANSAC法剔除误匹配点对,误匹配点集记为Γ”'i,则图像块fi与gi的匹配点对集
3.9)重复执行步骤3.1)~步骤3.8),直到f(x,y)和g(x,y)的对图像块都处理完毕,得到所有图像块的匹配点对集,即整幅图像的匹配点对集Γtotal;
4)在整幅图像的匹配点对集Γtotal的基础上,基于Delaunay三角网格划分基准图像和待配准图像,分别计算基准图像和待配准图像中对应三角网格的变换系数,将变换系数分别作用于待配准图像对应的三角网格内部的每个点,得到配准图像;
5)重复步骤2)~步骤4),直到蓝、绿、近红外谱段图像都完成配准,分别得到它们的配准图像。
所述步骤3.2)中计算图像块fi的梯度主方向Φi,通过以下步骤实现:
3.2.1)将[0,π)区间划分为n个方向,记为集合n为正整数;
3.2.2)计算图像块fi中每个像素的梯度幅值mag(x,y)和梯度方向θ(x,y);
其中,Δx=fi(x+1,y)-fi(x-1,y),Δy=fi(x,y+1)-fi(x,y-1);当Δx=0时,θ(x,y)=π/2,θ(x,y)的范围为[0,π);
3.2.3)根据划分的n个方向将步骤3.2.2)中得到的该像素的梯度幅值进行分解,具体步骤如下:
计算d=θ(x,y)×n/π,则该像素的梯度方向介于第个方向和第个方向之间,mod表示求余运算;将该像素的梯度幅值分解到第个方向和第个方向上,得到在第个方向的梯度幅值为在第个方向的梯度幅值为
3.2.4)对图像块fi中的所有像素重复步骤3.2.2)、步骤3.2.3),将所有像素的梯度幅值分解到n个方向上;
3.2.5)分别计算n个方向上的幅值和,最大的幅值和记为MAGmax,对应的方向记为属于 与相邻的两个方向表示为 对应的幅值和分别为MAGL、MAGR,利用MAGmax、MAGL、MAGR及其对应的方向ΦL、ΦR拟合一元二次曲线,则曲线最大峰值所对应的方向即为图像块fi的梯度主方向Φi。
所述步骤3.3)、步骤3.5)、步骤3.7)中任一项基于相位相关法和二元二次曲面拟合计算对应的亚像素匹配点,通过以下步骤实现:
(a)当待配准图像是蓝谱段或绿谱段图像时,以待匹配点为中心,在图像块fi、gi中分别截取大小为m1×m1个像素的块区域;当待配准图像是近红外谱段图像时,以待匹配点为中心,在主方向相位一致特征PC(fi)与PC(gi)中分别截取大小为m1×m1个像素的块区域,m1取128或256;
(b)基于相位相关计算两个块区域的相位相关系数,最大相位相关系数对应的平移量即为两个块区域的平移量,得到像素级匹配点;
(c)以像素级匹配点为中心的局部区域内的相位相关系数拟合二元二次曲面,所述局部区域的大小为3×3个像素;计算曲面最大值所对应的位置即为亚像素匹配点。
本发明与现有技术相比有益效果为:
(1)本发明提出一整套提取较均匀分布的角点、计算基准图像与待配准图像间的正反向匹配角点对集合、在无角点区域构造点并计算匹配点对集合、划分三角网格并计算变换系数等步骤的自动配准技术方案。
(2)本发明提出基于相位一致特征相位相关的配准方法,首先构造主方向的滤波器提取相位一致特征,然后通过相位相关和局部区域曲面拟合获得鲁棒的亚像素匹配精度,解决了近红外与蓝、绿、红谱段图像灰度差异大带来的配准难问题,同时仅提取主方向的相位一致特征减少了计算量。
(3)本发明提出从基准图像到待配准图像的正向匹配角点集、与从待配准图像到基准图像的反向匹配角点集的双向策略、以及在无明显特征区域构造点并计算匹配点对,可控制整幅图像划分的三角网格密度,有效地解决了由于地形、相机等导致的大幅面图像的不同区域存在不同程度畸变的问题。
附图说明
图1是本发明遥感卫星多光谱图像配准方法流程图;
图2a是本发明多光谱相机焦平面上的蓝谱段探测器到绿、红、近红外谱段探测器的距离,图2b是本发明对四个谱段图像进行裁剪获得粗配准图像的示意图;
图3是本发明有着显著灰度差异的红外与近红外谱段图像的角点匹配结果,其中图3a是红外谱段图像的结果,图3b是近红外谱段图像的结果;
图4是基于本发明方法的红外与近红外谱段图像的匹配结果及划分的三角网格,其中图4a是红外谱段图像的结果,图4b是近红外谱段图像的结果。
具体实施方式
一种遥感卫星多光谱图像配准方法,具体步骤如图1所示,该方法由以下步骤实现:
1、根据多光谱相机蓝、绿、红、近红外四个谱段探测器间的距离,对获取的四个谱段图像分别进行裁剪,得到粗配准的多光谱图像。
对多光谱相机采用分光棱镜获取多光谱图像的方式,四个谱段间没有偏移量,不需要对输入的四个谱段图像进行裁剪,获取的四个谱段图像即为粗配准的多光谱图像;对多光谱相机采用多色CCD探测器集成在一个器件上获取多光谱图像的方式,如图2a所示,焦平面上的蓝谱段探测器到绿、红、近红外谱段探测器的距离分别为d1、d2、d3,由于四个谱段探测器先后依次对地面同一区域成像,导致了谱段间图像的偏移量,如图2b所示,则从绿、红、近红外谱段图像的顶部分别裁剪δ1、δ2、δ3个像素,从蓝、绿、红谱段图像的底部分别裁剪δ3、δ2、δ1个像素,其中,τ为多光谱相机的单个探测器尺寸,表示下取整,裁剪后的图像即为粗配准的多光谱图像。
2、将粗配准的红谱段图像作为基准图像f(x,y),粗配准的蓝、绿、近红外图像中的任意一幅图像作为待配准图像,建立基准图像到待配准图像的直方图映射关系,以建立的直方图映射关系对待配准图像进行映射,得到映射图像g(x,y),所述f(x,y)与g(x,y)的大小为M×N个像素;M、N为正整数。
所述建立基准图像到待配准图像的直方图映射关系,其主要目的是减小两幅图像间的灰度差异,具体步骤包括:
(1)计算基准图像的直方图概率密度函数,作为期望直方图,计算期望直方图各级灰度值对应的累积概率密度。
(2)计算待配准图像的直方图概率密度函数,作为原始直方图,计算原始直方图各级灰度值对应的累积概率密度。
(3)基于直方图匹配建立原始直方图与期望直方图的直方图映射关系,将直方图映射关系作用到待配准图像,得到映射图像。
3、以大小为m×m个像素的图像块分别划分f(x,y)和g(x,y),得到的每对图像块表示为fi和gi,共得到对图像块,其中,表示上取整,m一般取值500~2000的整数。针对每一对图像块fi与gi进行如下步骤:
3.1、分别提取图像块fi与gi的Harris角点,图像块fi的所有角点构成fi的正向角点集。
Harris算子(Harris C.“A combined corner and edge detector”,1988)只用到图像灰度的一阶差分,具有计算简单、提取的角点特征均匀合理、可以定量提取特征点,并且算法稳定性好,对图像噪声等具有鲁棒性,特征点的数量可通过修改阈值来控制,因此比较适合大幅面的卫星图像的配准。
3.2、当待配准图像为近红外谱段图像时,计算图像块fi的梯度主方向Φi;构造方向为Φi的二维Log-Gabor滤波器,基于该滤波器对fi与gi进行相位一致特征检测,分别得到fi与gi的主方向相位一致特征PC(fi)与PC(gi)。
常用的空域检测算法主要采用亮度梯度值来表征边缘的强度,因此,检测出的边缘对图像亮度和对比度的变化非常敏感。1987年,Morrone等人在研究马赫带现象时发现,人类视觉感知的图像特征出现在图像Fourier谐波分量叠合最大的相位处,而特征的类型由相位的值决定。据此他提出了相位一致的理论和计算方法,并通过实验和理论证明了相位一致与人类视觉***对图像特征的认知相符合。
Kovesi提出了改进相位一致计算方法(参见:Peter Kovesi,“Image FeaturesFrom Phase Congruency”),但是需要多个方向(如6、8个方向)的滤波器与图像进行卷积,使得特征提取的计算量增加,对于大幅面的遥感图像提取相位一致特征耗费的内存大、时间长。本发明在Kovesi改进方法的基础上,只利用方向为Φi的二维Log-Gabor滤波器计算相位一致特征,极大地减少了计算量。
所述图像块fi的梯度主方向Φi的计算步骤如下:
(1)将[0,π)区间划分为n个方向,记为集合n为正整数,一般取值4~10。
(2)计算图像块fi中每个像素的梯度幅值mag(x,y)和梯度方向θ(x,y)。
θ(x,y)=tg-1|Δy/Δx|
其中,Δx=fi(x+1,y)-fi(x-1,y),Δy=fi(x,y+1)-fi(x,y-1);当Δx=0时,θ(x,y)=π/2,θ(x,y)的范围为[0,π)。
(3)根据划分的n个方向将步骤(2)中的得到的该像素的梯度幅值进行分解,具体步骤如下:
计算d=θ(x,y)×n/π,则该像素的梯度方向介于第个方向和第个方向之间,表示下取整,mod表示求余运算;将该像素的梯度幅值分解到第个方向和第个方向上,得到在第个方向的梯度幅值为在第个方向的梯度幅值为
(4)对图像块fi中的所有像素重复步骤(2)与(3),将所有像素的梯度幅值分解到n个方向上。
(5)分别计算n个方向上的幅值和,最大的幅值和记为MAGmax,对应的方向记为属于 与相邻的两个方向表示为 对应的幅值和分别为MAGL、MAGR,利用MAGmax、MAGL、MAGR及其对应的方向ΦL、ΦR拟合一元二次曲线,曲线方程为MAG=r·Φ2+s·Φ+t,MAG、Φ表示幅值和方向大小,r、s、t表示拟合得到的曲线系数,该曲线的最大峰值所对应的方向即为图像块fi的梯度主方向Φi。
所述构造方向为Φi的二维Log-Gabor滤波器对fi与gi进行相位一致特征检测,步骤如下:
(1)构造方向为Φi的二维Log-Gabor滤波器。
为减少计算量,仅利用计算出的主方向角Φi构造方向为Φi的二维Log-Gabor滤波器G(ω,θ),该滤波器定义在频域的极坐标上,表示形式为:
其中G(ω)为径向成分,ω0为中心频率,σω用于确定径向带宽;G(θ)为角度成分,Φi为滤波器的方向角,σθ是角度方向高斯函数的标准差,用于确定方向带宽。
可以看出,二维Log-Gabor滤波器在对数坐标系下为高斯函数,没有直流分量,并且包含更多的高频成分。本发明只利用计算出的主方向角建立滤波器,较传统相位一致特征提取的多个方向滤波器参与计算,在最大可能保留信息的同时,能够极大地减少计算量。
(2)利用构造的主方向滤波器对fi与gi进行相位一致特征检测。
Kovesi提出的改进相位一致表示为(参见:Peter Kovesi,“Image FeaturesFrom Phase Congruency”):
式中,Φi表示二维滤波器的方向,即计算得到的梯度主方向,n表示滤波器的尺度;为相位一致的权值函数;为相位偏离函数;为经过计算估计出的在方向Φi上的噪声干扰阈值;为在方向Φi、尺度n的幅值。
首先,将频域的二维滤波器G(ω,θ)变换到实域,则二维滤波器包括偶对称滤波器和奇对称滤波器两部分,它们与图像块fi或gi的卷积表示为:
计算幅值:
幅值与相位偏离函数相乘,分别得到:
其中:
根据公式(1)得到fi与gi在方向Φi上的主方向相位一致特征PC(fi)与PC(gi)。
3.3、将图像块fi的正向角点集中的每个角点作为待匹配点,在图像块gi中基于相位相关法和二元二次曲面拟合计算对应的亚像素匹配点,则得到图像块gi中的正向亚像素匹配角点集;正向角点集和正向亚像素匹配角点集记为正向匹配角点对集Γi。
特别是由于蓝、绿、红谱段与近红外谱段在灰度上有很大的差异,导致两幅图像提取的角点不一致,从而在不一致角点的基础上进行匹配得不到正确的匹配点,因此,直接用角点匹配得不到正确的匹配结果。图3为红谱段图像(图3a)和近红外谱段图像(图3b)的基于区域相关的角点匹配结果,图像大小1000×1000个像素,因此对于灰度差异较大的两幅图像基于区域相关的方法很难获得正确的匹配结果。
本发明提出的提取两幅图像中的Harris角点,不是建立基准图像的角点与待配准图像的角点的匹配关系,而是通过相位相关建立基准图像的角点到待配准图像的点、以及从待配准图像的角点到基准图像的点的双向匹配策略,一方面有效避免了两幅灰度差异大的图像中提取的角点不一致,同时充分利用了两幅图像的角点信息。1975年,Kuglin和Hines首次提出了基于傅里叶变换的相位相关算法用于检测两幅图像的平移关系(参见:The phase correlation imagealignment method),假设f2(x,y)为图像f1(x,y)在x和y方向分别平移x0和y0后的图像,即:
f2(x,y)=f1(x-x0,y-y0)
令F1(u,v)和F2(u,v)分别为f1(x,y)和f2(x,y)的傅里叶变换,由上式得:
图像f1(x,y)和f2(x,y)的互功率谱为:
式中,F1 *(u,v)为F1(u,v)的复共轭。互功率谱的傅里叶逆变换为一个二维脉冲函数δ(x-x0,y-y0),相位相关算法就是找到该脉冲函数的位置来确定平移参数。
所述基于相位相关法和二元二次曲面拟合计算对应的亚像素匹配点,主要步骤包括:
(1)当待配准图像是蓝谱段或绿谱段图像时,以待匹配点为中心,在图像块fi、gi中分别截取大小为m1×m1个像素的块区域;当待配准图像是近红外谱段图像时,以待匹配点为中心,在相位一致特征图PC(fi)与PC(gi)中分别截取大小为m1×m1个像素的块区域,m1取128或256;
(2)基于相位相关计算两个块区域的相位相关系数,最大相位相关系数对应的平移量即为两个块区域的平移量,得到像素级匹配点;
(3)以像素级匹配点为中心的局部区域内的相位相关系数拟合二元二次曲面,所述局部区域的大小为3×3个像素,曲面方程为z=ax2+by2+cxy+dx+ey+f,(x,y)表示点的坐标,z表示点(x,y)对应的相位相关系数,a~f为拟合得到的曲面系数,计算曲面最大值所对应的位置即为亚像素匹配点。
3.4、计算从gi提取的每个角点与正向亚像素匹配角点集中所有角点的最小距离,如果最小距离小于设定阈值则剔除gi中的该角点,设定阈值一般取5~10个像素,否则保留gi中的该角点,所有符合上述条件的角点构成gi中的反向角点集。
3.5、将gi中的反向角点集中的每个角点作为待匹配点,在图像块fi中基于相位相关法和二元二次曲面拟合计算其对应的亚像素匹配点,则得到在fi中的反向亚像素匹配角点集,反向角点集和反向亚像素匹配角点集构成反向匹配角点对集Γi;正反向匹配角点对集表示为Γi+Γ'i。
3.6、将图像块fi划分若干个大小为K×K个像素的块区域,K一般为200~500范围的正整数,遍历所有块区域,判断块区域内是否有匹配的角点对,即是否有集合Γi+Γ'i中的元素,若在区域内没有匹配的角点对则将该块区域的中心点作为构造点,否则继续遍历各个块区域,直到所有的块区域都处理完毕,则得到的所有构造点构成构造点集。
3.7、将构造点集中的每个点作为待匹配点,在图像块gi中基于相位相关法和二元二次曲面拟合计算其对应的亚像素匹配点,则得到在gi中的亚像素匹配点集,构造点集及其亚像素匹配点集构成匹配构造点对集Γ”i。因此,完整的匹配点集合记为Γi+Γ'i+Γ”i。
3.8、基于RANSAC法剔除误匹配点对,误匹配点集记为Γ”'i,则图像块fi与gi的匹配点对集Γi total=Γi+Γ'i+Γ”i-Γ”'i。
RANSAC法剔除误匹配点对的具体步骤如下:
(1)将图像块fi与gi都划分为大小为r×r个像素的块,r取值2048~4196的整数;
(2)针对每一块区域,统计最终匹配点对集合Γi total的匹配点对落入该区域的所有匹配点对;
(3)任取3对匹配点对计算从fi到gi的变换参数:
(4)对落入该区域的匹配点对集合中寻找所有符合距离小于阈值的点对,将它们作为内点,并记录满足条件的内点数量,如果内点数量大于事先设定的阈值则保留,否则舍弃;
(5)对步骤(3)和(4)重复L次,记录每一次的内点数量,迭代次数L的计算方法为:
其中,R为计算约束模型所需的最小匹配点对数量,在这里取3对匹配点对计算变换参数;ζ是期望达到的概率。
(6)选取对应内点数最大的变换参数,对应的内点作为最终的内点,也就是正确的匹配点对,不符合距离小于阈值的点对作为误匹配点对,予以剔除。
3.9、重复执行步骤3.1~3.8,直到f(x,y)和g(x,y)的对图像块都处理完毕,得到所有图像块的匹配点对集,即整幅图像的匹配点对集Γtotal。
如图4所示是本发明方法的红外谱段图像与近红外谱段图像的最终匹配点对示意图,图4a是红外谱段图像的结果,图4b是近红外谱段图像的结果。
4、在整幅图像的匹配点对集Γtotal的基础上,基于Delaunay三角网格划分基准图像和待配准图像,分别计算基准图像和待配准图像中对应三角网格的变换系数,将变换系数分别作用于待配准图像对应的三角网格内部的每个点,得到配准图像。在基准图像和待配准图像中划分的Delaunay三角网格分别如图4a和图4b所示。
5、重复步骤2~4,直到蓝、绿、近红外谱段图像都完成配准,分别得到它们的配准图像。
本发明未详细说明部分属本领域技术人员公知常识。
Claims (3)
1.一种遥感卫星多光谱图像配准方法,其特征在于步骤如下:
1)根据多光谱相机蓝、绿、红、近红外四个谱段探测器间的距离,对获取的四个谱段图像分别进行裁剪,得到粗配准的多光谱图像,具体步骤如下:
多光谱相机焦平面上的蓝谱段探测器到绿、红、近红外谱段探测器的距离分别为d1、d2、d3,则从绿、红、近红外谱段图像的顶部分别裁剪δ1、δ2、δ3个像素,从蓝、绿、红谱段图像的底部分别裁剪δ3、δ2、δ1个像素,其中, τ为多光谱相机的单个探测器尺寸,表示下取整,裁剪后的图像即为粗配准的多光谱图像;
2)将粗配准的红谱段图像作为基准图像f(x,y),粗配准的蓝、绿、近红外图像中的任意一幅图像作为待配准图像,建立基准图像到待配准图像的直方图映射关系,以建立的直方图映射关系对待配准图像进行映射,得到映射图像g(x,y),所述f(x,y)与g(x,y)的大小为M×N个像素;M、N为正整数;
3)以大小为m×m个像素的图像块分别划分f(x,y)与g(x,y),得到的每对图像块表示为fi与gi,共得到对图像块,其中,m为正整数,i=1,2,…,表示上取整;针对每一对图像块fi与gi进行如下步骤:
3.1)分别提取图像块fi与gi的Harris角点,图像块fi的所有角点构成fi的正向角点集;
3.2)当待配准图像为近红外谱段图像时,计算图像块fi的梯度主方向Φi;构造方向为Φi的二维Log-Gabor滤波器对fi与gi进行相位一致特征检测,分别得到fi与gi的主方向相位一致特征PC(fi)与PC(gi);
3.3)将图像块fi的正向角点集中的每个角点作为待匹配点,在图像块gi中基于相位相关法和二元二次曲面拟合计算对应的亚像素匹配点,则得到在图像块gi中的正向亚像素匹配角点集;正向角点集和正向亚像素匹配角点集记为正向匹配角点对集Γi;
3.4)计算从gi提取的每个角点与正向亚像素匹配角点集中所有角点的最小距离,如果最小距离小于设定阈值则剔除gi中的该角点,否则保留gi中的该角点,所有符合上述条件的角点构成gi中的反向角点集;
3.5)将gi中的反向角点集中的每个角点作为待匹配点,在图像块fi中基于相位相关法和二元二次曲面拟合计算对应的亚像素匹配点,则得到在fi中的反向亚像素匹配角点集;反向角点集和反向亚像素匹配角点集构成反向匹配角点对集Γ′i;正反向匹配角点对集表示为Γi+Γ′i;
3.6)将图像块fi划分若干个大小为K×K个像素的块区域,K为正整数;遍历所有块区域,判断块区域内是否有匹配的角点对,即是否有集合Γi+Γ′i中的元素,若没有则将该块区域的中心点作为构造点,所有的构造点构成构造点集;
3.7)将构造点集中的每个构造点作为待匹配点,在图像块gi中基于相位相关法和二元二次曲面拟合计算对应的亚像素匹配点,则得到在gi中的亚像素匹配点集;构造点集及其亚像素匹配点集构成匹配构造点对集Γ″i;完整的匹配点集合记为Γi+Γ′i+Γ″i;
3.8)基于RANSAC法剔除误匹配点对,误匹配点集记为Γ″′i,则图像块fi与gi的匹配点对集
3.9)重复执行步骤3.1)~步骤3.8),直到f(x,y)和g(x,y)的对图像块都处理完毕,得到所有图像块的匹配点对集,即整幅图像的匹配点对集Γtotal;
4)在整幅图像的匹配点对集Γtotal的基础上,基于Delaunay三角网格划分基准图像和待配准图像,分别计算基准图像和待配准图像中对应三角网格的变换系数,将变换系数分别作用于待配准图像对应的三角网格内部的每个点,得到配准图像;
5)重复步骤2)~步骤4),直到蓝、绿、近红外谱段图像都完成配准,分别得到它们的配准图像。
2.根据权利要求1所述的一种遥感卫星多光谱图像配准方法,其特征在于:所述步骤3.2)中计算图像块fi的梯度主方向Φi,通过以下步骤实现:
3.2.1)将[0,π)区间划分为n个方向,记为集合n为正整数;
3.2.2)计算图像块fi中每个像素的梯度幅值mag(x,y)和梯度方向θ(x,y);
其中,Δx=fi(x+1,y)-fi(x-1,y),Δy=fi(x,y+1)-fi(x,y-1);当Δx=0时,θ(x,y)=π/2,θ(x,y)的范围为[0,π);
3.2.3)根据划分的n个方向将步骤3.2.2)中得到的该像素的梯度幅值进行分解,具体步骤如下:
计算d=θ(x,y)×n/π,则该像素的梯度方向介于第个方向和第个方向之间,mod表示求余运算;将该像素的梯度幅值分解到第个方向和第个方向上,得到在第个方向的梯度幅值为在第个方向的梯度幅值为
3.2.4)对图像块fi中的所有像素重复步骤3.2.2)、步骤3.2.3),将所有像素的梯度幅值分解到n个方向上;
3.2.5)分别计算n个方向上的幅值和,最大的幅值和记为MAGmax,对应的方向记为属于与相邻的两个方向表示为 对应的幅值和分别为MAGL、MAGR,利用MAGmax、MAGL、MAGR及其对应的方向ΦL、ΦR拟合一元二次曲线,则曲线最大峰值所对应的方向即为图像块fi的梯度主方向Φi。
3.根据权利要求1所述的一种遥感卫星多光谱图像配准方法,其特征在于:所述步骤3.3)、步骤3.5)、步骤3.7)中任一项基于相位相关法和二元二次曲面拟合计算对应的亚像素匹配点,通过以下步骤实现:
(a)当待配准图像是蓝谱段或绿谱段图像时,以待匹配点为中心,在图像块fi、gi中分别截取大小为m1×m1个像素的块区域;当待配准图像是近红外谱段图像时,以待匹配点为中心,在主方向相位一致特征PC(fi)与PC(gi)中分别截取大小为m1×m1个像素的块区域,m1取128或256;
(b)基于相位相关计算两个块区域的相位相关系数,最大相位相关系数对应的平移量即为两个块区域的平移量,得到像素级匹配点;
(c)以像素级匹配点为中心的局部区域内的相位相关系数拟合二元二次曲面,所述局部区域的大小为3×3个像素;计算曲面最大值所对应的位置即为亚像素匹配点。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201510106747.XA CN104732532B (zh) | 2015-03-11 | 2015-03-11 | 一种遥感卫星多光谱图像配准方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201510106747.XA CN104732532B (zh) | 2015-03-11 | 2015-03-11 | 一种遥感卫星多光谱图像配准方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN104732532A true CN104732532A (zh) | 2015-06-24 |
CN104732532B CN104732532B (zh) | 2017-05-31 |
Family
ID=53456400
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201510106747.XA Expired - Fee Related CN104732532B (zh) | 2015-03-11 | 2015-03-11 | 一种遥感卫星多光谱图像配准方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN104732532B (zh) |
Cited By (13)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN105069743A (zh) * | 2015-07-28 | 2015-11-18 | 中国科学院长春光学精密机械与物理研究所 | 探测器拼接实时图像配准的方法 |
CN105261014A (zh) * | 2015-09-30 | 2016-01-20 | 西南交通大学 | 一种多传感器遥感影像匹配方法 |
CN106447663A (zh) * | 2016-09-30 | 2017-02-22 | 深圳市莫廷影像技术有限公司 | 一种去除重影的眼科oct图像高清配准方法及装置 |
CN108090869A (zh) * | 2017-11-27 | 2018-05-29 | 中国空间技术研究院 | 一种基于面阵cmos光学相机的星上超分辨率重建方法 |
CN108369740A (zh) * | 2015-12-10 | 2018-08-03 | 凯杰有限公司 | 用于将一个数字图像的至少一部分与另一个数字图像的至少一部分对准的方法 |
CN109523585A (zh) * | 2018-11-19 | 2019-03-26 | 武汉大学 | 一种基于方向相位一致性的多源遥感影像特征匹配方法 |
CN110599404A (zh) * | 2019-09-24 | 2019-12-20 | 陕西晟思智能测控有限公司 | 一种电路板显微图像拼接方法、装置、信息数据处理终端 |
CN111640142A (zh) * | 2019-12-25 | 2020-09-08 | 珠海大横琴科技发展有限公司 | 一种遥感图像多特征匹配方法、装置及电子设备 |
CN112598717A (zh) * | 2020-12-14 | 2021-04-02 | 珠海欧比特宇航科技股份有限公司 | 高光谱卫星影像全谱段配准方法及介质 |
CN112907486A (zh) * | 2021-03-18 | 2021-06-04 | 国家海洋信息中心 | 一种基于深度学习和色彩映射的遥感影像调色方法 |
CN113343747A (zh) * | 2021-03-30 | 2021-09-03 | 西南电子技术研究所(中国电子科技集团公司第十研究所) | 多模态图像鲁棒匹配vns的方法 |
CN114820739A (zh) * | 2022-07-01 | 2022-07-29 | 浙江工商大学 | 一种面向多光谱相机的图像快速配准方法及装置 |
CN118229749A (zh) * | 2024-05-24 | 2024-06-21 | 航天宏图信息技术股份有限公司 | 全月数字正射影像配准方法、装置、电子设备和存储介质 |
Citations (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101109818A (zh) * | 2006-07-20 | 2008-01-23 | 中国科学院自动化研究所 | 遥感影像高精度控制点自动选择方法 |
CN101751673A (zh) * | 2009-12-24 | 2010-06-23 | 中国资源卫星应用中心 | 基于相位一致特征的多光谱图像配准检测修正方法 |
CN102136142A (zh) * | 2011-03-16 | 2011-07-27 | 内蒙古科技大学 | 基于自适应三角形网格的非刚性医学图像配准方法 |
CN102982543A (zh) * | 2012-11-20 | 2013-03-20 | 北京航空航天大学深圳研究院 | 一种多源遥感图像配准方法 |
KR101412727B1 (ko) * | 2013-11-15 | 2014-07-01 | 동국대학교 산학협력단 | 얼굴 인식 장치 및 방법 |
CN104112278A (zh) * | 2014-08-01 | 2014-10-22 | 西安电子科技大学 | 基于协方差的多光谱图像实时配准方法 |
CN104240230A (zh) * | 2014-06-30 | 2014-12-24 | 华南理工大学 | 一种提高相位相关算法匹配精度的方法 |
-
2015
- 2015-03-11 CN CN201510106747.XA patent/CN104732532B/zh not_active Expired - Fee Related
Patent Citations (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101109818A (zh) * | 2006-07-20 | 2008-01-23 | 中国科学院自动化研究所 | 遥感影像高精度控制点自动选择方法 |
CN101751673A (zh) * | 2009-12-24 | 2010-06-23 | 中国资源卫星应用中心 | 基于相位一致特征的多光谱图像配准检测修正方法 |
CN102136142A (zh) * | 2011-03-16 | 2011-07-27 | 内蒙古科技大学 | 基于自适应三角形网格的非刚性医学图像配准方法 |
CN102982543A (zh) * | 2012-11-20 | 2013-03-20 | 北京航空航天大学深圳研究院 | 一种多源遥感图像配准方法 |
KR101412727B1 (ko) * | 2013-11-15 | 2014-07-01 | 동국대학교 산학협력단 | 얼굴 인식 장치 및 방법 |
CN104240230A (zh) * | 2014-06-30 | 2014-12-24 | 华南理工大学 | 一种提高相位相关算法匹配精度的方法 |
CN104112278A (zh) * | 2014-08-01 | 2014-10-22 | 西安电子科技大学 | 基于协方差的多光谱图像实时配准方法 |
Cited By (19)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN105069743B (zh) * | 2015-07-28 | 2018-06-26 | 中国科学院长春光学精密机械与物理研究所 | 探测器拼接实时图像配准的方法 |
CN105069743A (zh) * | 2015-07-28 | 2015-11-18 | 中国科学院长春光学精密机械与物理研究所 | 探测器拼接实时图像配准的方法 |
CN105261014A (zh) * | 2015-09-30 | 2016-01-20 | 西南交通大学 | 一种多传感器遥感影像匹配方法 |
CN105261014B (zh) * | 2015-09-30 | 2018-01-12 | 西南交通大学 | 一种多传感器遥感影像匹配方法 |
CN108369740A (zh) * | 2015-12-10 | 2018-08-03 | 凯杰有限公司 | 用于将一个数字图像的至少一部分与另一个数字图像的至少一部分对准的方法 |
CN106447663A (zh) * | 2016-09-30 | 2017-02-22 | 深圳市莫廷影像技术有限公司 | 一种去除重影的眼科oct图像高清配准方法及装置 |
CN108090869A (zh) * | 2017-11-27 | 2018-05-29 | 中国空间技术研究院 | 一种基于面阵cmos光学相机的星上超分辨率重建方法 |
CN109523585B (zh) * | 2018-11-19 | 2021-10-22 | 武汉大学 | 一种基于方向相位一致性的多源遥感影像特征匹配方法 |
CN109523585A (zh) * | 2018-11-19 | 2019-03-26 | 武汉大学 | 一种基于方向相位一致性的多源遥感影像特征匹配方法 |
CN110599404A (zh) * | 2019-09-24 | 2019-12-20 | 陕西晟思智能测控有限公司 | 一种电路板显微图像拼接方法、装置、信息数据处理终端 |
CN111640142A (zh) * | 2019-12-25 | 2020-09-08 | 珠海大横琴科技发展有限公司 | 一种遥感图像多特征匹配方法、装置及电子设备 |
CN112598717A (zh) * | 2020-12-14 | 2021-04-02 | 珠海欧比特宇航科技股份有限公司 | 高光谱卫星影像全谱段配准方法及介质 |
CN112598717B (zh) * | 2020-12-14 | 2024-05-17 | 珠海欧比特卫星大数据有限公司 | 高光谱卫星影像全谱段配准方法及介质 |
CN112907486A (zh) * | 2021-03-18 | 2021-06-04 | 国家海洋信息中心 | 一种基于深度学习和色彩映射的遥感影像调色方法 |
CN112907486B (zh) * | 2021-03-18 | 2022-12-09 | 国家海洋信息中心 | 一种基于深度学习和色彩映射的遥感影像调色方法 |
CN113343747A (zh) * | 2021-03-30 | 2021-09-03 | 西南电子技术研究所(中国电子科技集团公司第十研究所) | 多模态图像鲁棒匹配vns的方法 |
CN114820739A (zh) * | 2022-07-01 | 2022-07-29 | 浙江工商大学 | 一种面向多光谱相机的图像快速配准方法及装置 |
CN114820739B (zh) * | 2022-07-01 | 2022-10-11 | 浙江工商大学 | 一种面向多光谱相机的图像快速配准方法及装置 |
CN118229749A (zh) * | 2024-05-24 | 2024-06-21 | 航天宏图信息技术股份有限公司 | 全月数字正射影像配准方法、装置、电子设备和存储介质 |
Also Published As
Publication number | Publication date |
---|---|
CN104732532B (zh) | 2017-05-31 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN104732532A (zh) | 一种遥感卫星多光谱图像配准方法 | |
CN104318548B (zh) | 一种基于空间稀疏度和sift特征提取的快速图像配准实现方法 | |
Ye et al. | Hopc: A novel similarity metric based on geometric structural properties for multi-modal remote sensing image matching | |
CN104700404B (zh) | 一种果实定位识别方法 | |
Chen et al. | SAR and multispectral image fusion using generalized IHS transform based on à trous wavelet and EMD decompositions | |
CN102436652B (zh) | 一种多源遥感图像自动配准方法 | |
WO2018024030A1 (zh) | 一种夜视红外图像中的基于显著性的道路目标提取方法 | |
CN102073874B (zh) | 附加几何约束的航天三线阵ccd相机多影像立体匹配方法 | |
CN102800074B (zh) | 基于轮廓波变换的sar图像变化检测差异图生成方法 | |
CN103198483A (zh) | 基于边缘和光谱反射率曲线的多时相遥感图像配准方法 | |
CN103679714A (zh) | 一种基于梯度互相关的光学和sar图像自动配准方法 | |
CN102208109A (zh) | X射线图像和激光图像的异源图像配准方法 | |
CN102609918A (zh) | 基于图像特征配准的航空多光谱遥感影像几何精校正方法 | |
Fytsilis et al. | A methodology for near real-time change detection between Unmanned Aerial Vehicle and wide area satellite images | |
CN104361590A (zh) | 一种控制点自适应分布的高分辨率遥感影像配准方法 | |
Zhu et al. | Robust registration of aerial images and LiDAR data using spatial constraints and Gabor structural features | |
CN103116881A (zh) | 基于PCA与Shearlet变换的遥感图像融合方法 | |
CN107689051A (zh) | 一种基于变化因子的多时相sar影像变化检测方法 | |
CN105956544A (zh) | 一种基于结构指数特征的遥感影像道路交叉口提取的方法 | |
CN108446637B (zh) | 基于立体图模型的sar图像变化检测方法 | |
CN105513060A (zh) | 一种视觉感知启发的高分辨率遥感图像分割方法 | |
Li et al. | Automatic Road Extraction from High-Resolution Remote Sensing Image Based on Bat Model and Mutual Information Matching. | |
CN106910178A (zh) | 一种基于色调统计特性分类的多角度sar图像融合方法 | |
CN110874833A (zh) | 基于超图匹配的sar图像变化检测方法 | |
CN106355576A (zh) | 基于mrf图像分割算法的sar图像配准方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant | ||
CF01 | Termination of patent right due to non-payment of annual fee | ||
CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20170531 Termination date: 20190311 |