一种针对大视角差异光学与SAR遥感图像的配准方法

一种针对大视角差异光学与sar遥感图像的配准方法
技术领域：
：1.本发明属于多源遥感图像配准
技术领域：
：，具体涉及一种针对大视角差异光学与sar遥感图像的配准方法。
背景技术：
：：2.随着航天技术和多源观测的发展，多平台、多波段以及多光谱的图像数据量在不断的增加。远程成像传感器包括了光学传感器、合成孔径雷达(syntheticapertureradar,sar)和红外雷达，它们各有优势。作为一种主动微波遥感系统，sar不受天气和时间的影响，可以提供多视角、多波段以及多极化的高分辨率图像。光学传感器是一种被动传感器，在良好的地面条件下，可以获得丰富的灰度和纹理信息。由于不同的传感器优势互补，多源遥感图像配准已经成为了众多学者研究的热点。3.已经提出的大部分算法，例如os-sift(optical-to-sarscale-invariantfeaturetransform，os-sift)算法，虽然考虑到了光图与sar图中存在着不同类型的噪声影响并采用了不同的方式计算图像梯度，但依然无法提取到足够多的同名点，并且在建立描述子时，仍无法较好地处理非线性辐射畸变带来的影响，导致最后的配准结果并不理想。rift(radiation-variationinsensitivefeaturetransform，rift)算法，虽然利用了相位一致性代替图像强度进行特征点检测并采用了在最大索引图上描述特征，建立的描述子对于非线性辐射畸变具有很好的鲁棒性，但却忽略了sar图中存在的相干斑噪声影响，导致在sar图中检测到的关键点很多都分布在平滑区域，关键点的重复率较低。ass(adjacentself-similarity，ass)算法，虽然参考lee滤波器设计了权重函数来对相干斑噪声进行抑制并获得了较好的特征点检测结果，但建立的描述子不够鲁棒，最后得到的正确匹配点对数量较少。另外，已经提出的大部分算法都没有考虑大视角差异下的配准情况，当参考图与待配准图之间存在严重的几何畸变时，配准精度往往较差。4.专利cn202310420988.6(“一种针对大视角差异sar图像的配准方法”)提供了一种针对大视角差异sar图像的配准方法。该方法提出的三阶段配准框架虽然在sar图像上可以得到大量正确匹配点对并实现较高的配准精度，但在异源的光学和sar遥感图像数据集上却表现较差。该方法虽然可以解决光学和sar图像配准中因为大视角差异带来的几何畸变影响，但却没有考虑到异源图像间存在的噪声类型是不同的，需要采用不同的方式对噪声进行抑制，并且异源图像间还存在着严重的非线性辐射畸变影响。技术实现要素：5.为了解决相关技术中当参考图与待配准图之间存在严重的几何畸变时，配准精度较差的问题，本发明提供了一种针对大视角差异光学与sar遥感图像的配准方法。本发明要解决的技术问题通过以下技术方案实现：6.本发明提供一种针对大视角差异光学与sar遥感图像的配准方法，如图1，该方法包括：7.(1)利用os-sift算法、多尺度harris-affine算法以及多尺度mser(maximallystableextremalregions，mser)算法分别提取参考图像i1与待配准图像i2中的关键点，然后，建立os-sift描述子并进行匹配，最终得到了三种具有不同性质的匹配点对集合：未受视角变化影响的点对集合so、具有仿射不变性的角点点对集合sh，以及具有仿射不变性的区域点对集合sm；8.(2)根据步骤(1)得到的所有匹配点对，计算全局变换模型γ，利用全局变换模型γ对图像i2进行变换，得到图像在参考图像i1与待配准图像上利用risfm(radiation-insensitivestructuralfeaturematching，risfm)算法进行匹配，得到匹配点对集合sr；其中，risfm算法具体步骤为：首先建立参考图像i1与待配准图像的最小自相似图，通过极大值检测与非极大值抑制在最小自相似图上检测关键点；之后，基于二维log-garbor小波变换，计算参考图像i1与待配准图像在各个方向下的log-gabor响应，寻找最大响应的方向，并将方向的索引值取出建立最大索引图mim；然后，基于检测到的关键点以及最大索引图mim建立描述子；最后，利用b-nndr(bidirectionalnearestneighbordistanceratio，b-nndr)算法找到参考图像i1与待配准图像之间的采样集与一致集，再经过fsc(fastsampleconsensus，fsc)算法去除误匹配点对后，得到匹配点对集合sr；9.(3)将点对集合so、sh、sm与sr进行汇总后，利用los-flow(localoptical-to-sarflow，los-flow)算法进行精配，再利用lsor(length-andslope-basedoutlierremovalmethod，lsor)算法去除其中可能包含的误匹配点对，最终得到匹配点对集合sfinal，利用sfinal重新计算全局变换模型。10.在一些实施例中，步骤(1)包括：11.1a)未受视角变化影响的点对集合so：利用os-sift算法提取图像中的关键点并建立os-sift描述子，使用b-nndr算法找到参考图像i1与待配准图像i2之间的采样集与一致集，再经过fsc算法去除误匹配点对后，最终得到匹配点对集合so；12.1b)具有仿射不变性的角点点对集合sh：利用多尺度harris-affine算法提取图像中的关键点并得到每个关键点处的局部仿射变换矩阵，其中，提取图像中的关键点的步骤与os-sift算法相同；之后，在每个关键点处利用局部仿射变换矩阵对图像进行变换后建立os-sift描述子，使用b-nndr算法找到参考图像i1与待配准图像i2之间的采样集与一致集，再经过fsc算法去除误匹配点对后，最终得到匹配点对集合sh；13.1c)具有仿射不变性的区域点对集合sm：利用多尺度mser算法提取图像中的关键点，并得到每个关键点处的局部仿射变换矩阵，在每个关键点处利用局部仿射变换矩阵对图像进行变换后建立os-sift描述子，使用b-nndr算法找到参考图像i1与待配准图像i2之间的采样集与一致集，再经过fsc算法去除误匹配点对后，最终得到匹配点对集合sm。14.在一些实施例中，步骤1a)和1b)中利用os-sift算法提取图像中的关键点的具体步骤包括：15.2a)选取一组指数加权尺度因子[α0,α1,...,αn-1]，其中，初始值α0＝2，αi＝α0*ki，(i∈[1,n-1])，k＝21/3，n＝8；[0016]2b)根据α的尺度序列[α0,α1,...,αn-1]，对参考图像i1与待配准图像i2分别利用sobel算子与roewa(ratioofexponentiallyweightedaverages，roewa)算子计算图像在不同尺度因子下的水平梯度gx,α与垂直梯度gy,α，从而得到图像在不同尺度下的梯度幅值magα与梯度方向oriα：[0017][0018][0019]2c)根据已经求出的梯度计算图像中每个像素点处的harris矩阵ch(x,y,α)：[0020][0021]其中，表示标准差为的高斯核，*表示卷积；[0022]2d)利用ch(x,y,α)计算图像中每个像素点处的harris响应值rh(x,y,α)：[0023]rh(x,y,α)＝det(ch(x,y,α))-d·tr(ch(x,y,α))2；[0024]其中，d为0.04～0.06之间的数值，det表示行列式的值，tr表示行列式的迹；[0025]2e)在尺度空间图像的同一层中，将每个像素点的harris响应值与它周围的8邻域内的像素点的响应值，以及一个全局阈值dh分别作比较，若在邻域内中心点的响应值是最大的，并且大于全局阈值dh，那么该点为检测到的关键点，其中，全局阈值dh取0.85。[0026]在一些实施例中，步骤1a)、1b)和1c)中建立os-sift描述子的步骤具体包括：[0027]3a)利用梯度幅值magα与梯度方向oriα，以关键点为圆心，取半径为6α的圆形邻域，在圆形邻域内，将0～360°平均分为18份，每份代表20°的范围，横坐标表示梯度方向角，纵坐标表示梯度幅值，遍历该邻域内的所有像素点，对邻域内的每一个像素点，首先找到该像素点梯度方向角对应的直方图幅值柱，然后，将该像素点的梯度幅值累加在该点所在方向角的直方图幅值柱上，从而得到该关键点的主方向直方图，并对直方图进行平滑处理；α表示该关键点所在的图像层尺度；[0028]3b)将平滑后的直方图峰值对应的方向角，作为该关键点的主方向，大于峰值能量80％的柱对应的方向角作为辅方向；[0029]3c)以关键点为圆心，分别取半径为rmax＝8α、12α以及16α的圆形邻域，将每个圆形邻域以主方向角为参考进行旋转，使特征描述子具有旋转不变性，然后，将每个圆形邻域由内到外划分成一个圆形和两个圆环，同心圆半径由内到外分别为0.25·rmax，0.75·rmax和rmax，再将两个圆环分别在0～360°度范围内以45°为间隔等分的划分成8个子区域，加上中心的圆形则将每个邻域一共划分为了17个子区域，三个圆形邻域一共有51个子区域，在每个子区域上将0～360°平均分为8份，建立梯度方向直方图，将三个圆形邻域的一共51个子区域的直方图向量拼接起来并归一化就生成了408维的os-sift描述子。[0030]在一些实施例中，步骤1a)、1b)、1c)和(2)中利用b-nndr算法找到参考图像与待配准图像之间的采样集与一致集的步骤具体包括：[0031]4a)在关键点空间内，以特征向量之间的欧氏距离衡量两个关键点的远近程度，越近则越相似；对于参考图像上的一个关键点，通过计算关键点的描述子之间的欧氏距离，在待配准图像上寻找距离该关键点最近与次近的关键点，最近距离与次近距离分别用dmin与dnd表示；如果则此关键点与距离它最近的关键点互为正确匹配点对；遍历参考图像上的每一个关键点，当阈值distratio取0.9时得到采样集当阈值distratio取0.999时得到一致集[0032]4b)遍历待配准图像上的每一个关键点，在参考图像上寻找与该关键点互为正确匹配点对的对应点。当阈值distratio取0.9时得到采样集当阈值distratio取0.999时得到一致集[0033]4c)参考图像与待配准图像之间的采样集与一致集为ch和cl：[0034][0035]∪表示取并集。[0036]在一些实施例中，步骤1a)、1b)、1c)和(2)中利用fsc算法去除误匹配点对的步骤具体包括：[0037]5a)设置迭代次数为n，在第t次迭代过程中，从点对集合ch中随机挑选三对匹配点：t为1～n中的整数；[0038]5b)利用这三对点计算图像的变换模型θt；[0039]5c)利用得到的变换模型θt计算点对集合cl中的匹配点对ci的变换误差e(ci,θt)；[0040][0041][0042]其中，(xi,yi)表示匹配点对ci在参考图像上的关键点坐标，表示匹配点对ci在待配准图像上的对应点坐标；t((xi,yi),θt)表示利用变换模型θt对(xi,yi)变换后在待配准图像上的对应位置，e(ci,θt)表示在变换模型θt下，匹配点对ci的变换误差；[0043]5d)遍历点对集合cl中的每一对匹配点，将变换误差e＜3的所有匹配点对进行汇总，得到点对集合ct；[0044]5e)当算法迭代了n次后结束，此时得到了c1,c2,...,cn，一共n个点对集合，取出其中拥有最多点对数量的点对集合，作为fsc算法最终得到的匹配点对集合。[0045]在一些实施例中，步骤(2)中建立参考图像i1与待配准图像的最小自相似图，通过极大值检测与非极大值抑制在最小自相似图上检测关键点的步骤具体包括：[0046]9a)对参考图像i1进行裁剪，构建子图像：以图像i1的中心像素点为中心，建立一个大小为lsub×wsub的搜索框，其中，lsub＝l-10，wsub＝w-10，l与w为参考图像i1的长和宽；当搜索框不移动时对参考图像i1按照搜索框所在位置进行裁剪得到中心子图像块subic，再将搜索框分别向0°、45°、90°、135°和180°偏移一个像素后，对参考图像i1按照当前搜索框所在位置进行裁剪得到偏移子图像块subi1、subi2、subi3、subi4和subi5；[0047]9b)通过以下公式获得在各个方向上偏移一个像素后的偏移子图像块subiθ：[0048][0049]θ＝180(o-1)/no,o＝1,2,...,no；其中，θ表示偏移的方向，no表示将θ∈[0°,180°)分成no份；[0050]9c)计算中心子图像块subic上每个像素点(x,y)的权重值v(x,y)；[0051][0052][0053][0054]σl(x,y)＝ηl(x,y)-(μl(x,y))2；[0055]其中，μl、ηl和σl分别表示局部均值、均方值和方差值；(x,y)表示中心子图像块subic上的一个像素点；wl(x,y,rl)表示以(x,y)为中心，半径为rl的一个圆形邻域，当l取1时，r1等于2，当l取2时，r2等于4；nl表示圆形邻域wl内的像素点个数；[0056]9d)利用均值滤波器计算θ方向下的自相似图sθ：[0057][0058]其中，subic表示中心子图像块，subiθ表示向θ方向偏移一个像素后的偏移子图像块，v表示中心子图像块subic上每个像素点的权重值，w1表示半径为2的一个圆形邻域，表示对每个像素点取半径为2的圆形邻域进行均值滤波；[0059]9e)寻找最小自相似图smin：[0060][0061]9f)在最小自相似图smin上进行极大值检测以及非极大值抑制，寻找图像i1的关键点；[0062]9g)对待配准图像进行以上相同操作，寻找图像的关键点。[0063]在一些实施例中，步骤(2)中基于二维log-garbor小波变换，计算参考图像i1与待配准图像在各个方向下的log-gabor响应，寻找最大响应的方向，并将方向的索引值取出建立最大索引图mim的步骤具体包括：[0064]10a)选取4个不同的尺度s＝1,2,3,4，6个不同的方向o＝0°,30°,60°,90°,120°,150°，分别对应索引号1、2、3、4、5和6，建立二维偶对称log-gabor滤波器leven(x,y,s,o)和二维奇对称log-gabor滤波器lodd(x,y,s,o)；[0065]10b)计算图像上的像素点(x,y)在o方向、s尺度下的响应幅值aso(x,y)：[0066]eso(x,y)＝i(x,y)*leven(x,y,s,o)，[0067]oso(x,y)＝i(x,y)*lodd(x,y,s,o)，[0068]其中，*表示卷积；[0069][0070]10c)将o方向下所有不同尺度的响应幅值进行累加，得到像素点(x,y)在o方向下的响应值：[0071][0072]10d)取出最大响应对应的方向索引值作为点(x,y)处的值，遍历图像上每一个像素点，建立最大索引图mim。[0073]在一些实施例中，步骤(2)中基于检测到的关键点以及最大索引图mim建立描述子的步骤具体包括：[0074]11a)对于最大索引图mim上的每一个关键点，以其为中心，裁剪出一个大小为96×96的图像块，再将图像块划分成6×6个子网格，在每一个网格内，建立直方图，横坐标表示方向的索引号，取值范围为1～6，纵坐标表示该方向的索引号出现的次数，如此，可以将每个子网格转化为一个六维直方图向量，最后，将36个子区域的六维直方图向量拼接起来并进行归一化，就生成了该关键点的216维的描述子。[0075]在一些实施例中，步骤(3)中将点对集合so、sh、sm与sr进行汇总后，利用los-flow算法进行精配的步骤具体包括：[0076]12a)利用点对集合so、sh、sm计算得到的待配准图像i2到参考图像i1的变换模型γ，使用变换模型γ对待配准图像i2进行变换，得到变换后的图像[0077]对点对集合so、sh、sm中在待配准图像i2上的一个点(xs,ys)，利用变换模型γ计算点(xs,ys)在变换后的图像上的对应位置坐标[0078]将点对集合so、sh、sm中的所有点都进行以上相同的变换，从而得到新的点对集合中的所有点都进行以上相同的变换，从而得到新的点对集合[0079]12b)对于图像i1与图像分别利用sobel算子与roewa算子计算在尺度因子α＝2下的水平梯度gx与垂直梯度gy，从而得到图像的梯度幅值mag与梯度方向ori：[0080][0081][0082]对于图像i1与上的一个像素点，以该像素点为圆心，取半径为rmax＝24的圆形邻域，并将此圆形邻域由内到外划分成一个圆形和两个圆环，同心圆半径由内到外分别为0.25·rmax，0.75·rmax和rmax，再将两个圆环分别在0～360°度范围内以45°为间隔等分的划分成8个子区域，加上中心的圆形则将此邻域一共划分为了17个子区域，在每个子区域内，将0～360°平均分为8份，每份代表45°的范围，横坐标表示梯度方向角，纵坐标表示梯度幅值，遍历该子区域内的所有像素点，对区域内的每一个像素点，首先找到该像素点梯度方向角对应的直方图幅值柱，然后，将该像素点的梯度幅值累加在该点所在方向角的直方柱上，从而得到该点的梯度方向直方图，并将其转化为一个八维直方图向量，17个子区域的直方图向量拼接起来并进行归一化后，就生成了该像素点的136维的os-sift描述子；[0083]通过分别计算图像i1与上每个像素点的136维os-sift描述子，就形成了图像i1的os-sift描述子图i1_desc，以及图像的os-sift描述子图[0084]12c)待精配点对集合为和sr的并集，对于待精配点对集合中在图像i1上的一个点(xi,yi)，以点(xi,yi)为中心，在os-sift描述子图i1_desc上取出一个边长为2·rlf+1的局部正方形区域描述子图像i1_squ，对点(xi,yi)在图像上的对应点进行相同的操作，得到局部正方形区域描述子图像为其中，rlf取61；[0085]12d)将得到的i1_squ和代入损失函数e(w)，计算出光流矢量w，损失函数e(w)为：[0086][0087]其中，p＝(x,y)表示局部正方形区域描述子图像上的某个像素点，w(p)＝(u(p),v(p))表示p点的光流矢量，其中，u(p)表示p点在水平方向上的偏移量，v(p)表示p点在垂直方向上的偏移量；ε表示以p点为中心，加上其相邻8个点组成的区域，q表示该区域中不包含中心点p的某个点；参数η和α分别取0.001与0.03，参数t和d分别取0.1和0.6；上述(1)式为数据项，约束了沿着光流矢量w(p)的os-sift描述子要相互匹配；(2)式为小位移项，约束了在没有其他可用信息的情况下，光流矢量要尽可能的小；(3)式为平滑项，约束了相邻像素的光流矢量要相似；[0088]12e)计算图像上点(xit，yit)的新坐标(xit_new，yit_new)：[0089]xit_new＝xit+u(rlf+1,rlf+1)；[0090]yit_new＝yit+v(rlf+1,rlf+1)；[0091]12f)遍历待精配点对集合中的每一对匹配点后得到了更加准确的点对集合[0092]本发明具有如下有益技术效果：[0093]1)本发明深入分析了光学和sar遥感图像中包含的不同噪声影响，在第一阶段，采用了不同的方式计算光图与sar图的梯度，使得多尺度harris-affine算法与多尺度mser算法获得了更多的具有仿射不变性的匹配点对，通过计算全局变换矩阵对待配准图进行变换，初步消除了图像间几何畸变带来的影响。[0094]2)本发明针对异源遥感图像配准中非线性辐射畸变带来的影响，提出了risfm算法。算法利用最小自相似图检测关键点，极大的抑制了sar图像中相干斑噪声的影响，在光学图像和sar图像中获得了大量重复率高且性质稳定的关键点；算法基于二维log-gabor小波变换，求图像在多方向下的响应，通过取出最大响应的方向索引值，建立最大索引图，在最大索引图上建立描述子并进行匹配，通过利用结构信息而避免了图像间非线性辐射畸变带来的影响，最终获得了大量正确的匹配点对。[0095]3)本发明针对大视角差异光学和sar遥感图像配准问题，提出了三阶段配准框架，有效提升了图像配准的精度。[0096]以下将结合附图及实施例对本发明做进一步详细说明。附图说明[0097]图1为本发明实施例提供的针对大视角差异光学与sar遥感图像的配准方法的一个流程图；[0098]图2为本发明实施例提供的示例性的对参考图像和待配准图像进行配准的阶段示意图；[0099]图3为本发明实施例提供的示例性的os-sift描述子的图像；[0100]图4为本发明实施例提供的示例性的在参考图为光学图像，待配准图为sar图像时risfm算法的一个实现流程图；[0101]图5为本发明实施例提供的示例性的test1输入的参考图像；[0102]图6为本发明实施例提供的示例性的test1输入的待配准图像；[0103]图7为本发明实施例提供的示例性的test2输入的参考图像；[0104]图8为本发明实施例提供的示例性的test2输入的待配准图像；[0105]图9为本发明实施例提供的示例性的test1的配准后的棋盘图像；[0106]图10为本发明实施例提供的示例性的test1的正确匹配点对图；[0107]图11为本发明实施例提供的示例性的test2的配准后的棋盘图；[0108]图12为本发明实施例提供的示例性的test2的正确匹配点对图。具体实施方式[0109]下面结合具体实施例对本发明做进一步详细的描述，但本发明的实施方式不限于此。[0110]在本发明的描述中，术语“第一”、“第二”仅用于描述目的，而不能理解为指示或暗示相对重要性或者隐含指明所指示的技术特征的数量。由此，限定有“第一”、“第二”的特征可以明示或者隐含地包括一个或者更多个该特征。在本发明的描述中，“多个”的含义是两个或两个以上，除非另有明确具体的限定。[0111]在本说明书的描述中，参考术语“一个实施例”、“一些实施例”、“示例”、“具体示例”、或“一些示例”等的描述意指结合该实施例或示例描述的具体特征、结构、材料或者特点包含于本发明的至少一个实施例或示例中。在本说明书中，对上述术语的示意性表述不必须针对的是相同的实施例或示例。而且，描述的具体特征、结构、材料或者特点可以在任何的一个或多个实施例或示例中以合适的方式结合。此外，本领域的技术人员可以将本说明书中描述的不同实施例或示例进行接合和组合。[0112]尽管在此结合各实施例对本发明进行了描述，然而，在实施所要求保护的本发明过程中，本领域技术人员通过查看所述附图、公开内容、以及所附权利要求书，可理解并实现所述公开实施例的其他变化。在权利要求中，“包括”(comprising)一词不排除其他组成部分或步骤，“一”或“一个”不排除多个的情况。单个处理器或其他单元可以实现权利要求中列举的若干项功能。相互不同的从属权利要求中记载了某些措施，但这并不表示这些措施不能组合起来产生良好的效果。[0113]图1是本发明实施例提供的针对大视角差异光学与sar遥感图像的配准方法的一个流程图；图2为本发明实施例提供的示例性的对参考图像和待配准图像进行配准的阶段示意图；如图1和图2所示，本发明在对参考图像i1和对应的待配准图像i2进行配准时分为三个阶段：[0114](1)利用os-sift算法、多尺度harris-affine算法以及多尺度mser算法分别提取参考图像i1与待配准图像i2中的关键点，然后，建立os-sift描述子并进行匹配，最终得到了三种具有不同性质的匹配点对集合：未受视角变化影响的点对集合so、具有仿射不变性的角点点对集合sh，以及具有仿射不变性的区域点对集合sm；[0115]1a)未受视角变化影响的点对集合so：利用os-sift算法提取图像中的关键点并建立os-sift描述子，使用b-nndr算法找到参考图像i1与待配准图像i2之间的采样集与一致集，再经过fsc算法去除误匹配点对后，最终得到匹配点对集合so；[0116]1b)具有仿射不变性的角点点对集合sh：利用多尺度harris-affine算法提取图像中的关键点并得到每个关键点处的局部仿射变换矩阵，其中，提取图像中的关键点的步骤与os-sift算法相同；之后，在每个关键点处利用局部仿射变换矩阵对图像进行变换后建立os-sift描述子，使用b-nndr算法找到参考图像i1与待配准图像i2之间的采样集与一致集，再经过fsc算法去除误匹配点对后，最终得到匹配点对集合sh；[0117]1c)具有仿射不变性的区域点对集合sm：利用多尺度mser算法提取图像中的关键点，并得到每个关键点处的局部仿射变换矩阵，在每个关键点处利用局部仿射变换矩阵对图像进行变换后建立os-sift描述子，使用b-nndr算法找到参考图像i1与待配准图像i2之间的采样集与一致集，再经过fsc算法去除误匹配点对后，最终得到匹配点对集合sm。[0118](2)根据步骤(1)得到的所有匹配点对，计算全局变换模型γ，利用全局变换模型γ对图像i2进行变换，得到图像在参考图像i1与待配准图像上利用risfm算法进行匹配，得到匹配点对集合sr；其中，risfm算法具体步骤为：首先建立参考图像i1与待配准图像的最小自相似图，通过极大值检测与非极大值抑制在最小自相似图上检测关键点；之后，基于二维log-garbor小波变换，计算参考图像i1与待配准图像在各个方向下的log-gabor响应，寻找最大响应的方向，并将方向的索引值取出建立最大索引图mim；然后，基于检测到的关键点以及最大索引图mim建立描述子；最后，利用b-nndr算法找到参考图像i1与待配准图像之间的采样集与一致集，再经过fsc算法去除误匹配点对后，得到匹配点对集合sr；[0119](3)将点对集合so、sh、sm与sr进行汇总后，利用los-flow算法进行精配，再利用lsor算法去除其中可能包含的误匹配点对，最终得到匹配点对集合sfinal，利用sfinal重新计算全局变换模型。[0120]本发明中，步骤1a)和1b)中利用os-sift算法提取图像中的关键点的操作步骤如下：[0121]2a)选取一组指数加权尺度因子[α0,α1,...,αn-1]，其中，初始值α0＝2，αi＝α0*ki，(i∈[1,n-1])，k＝21/3，n＝8；[0122]2b)根据α的尺度序列[α0,α1,...,αn-1]，对参考图像i1与待配准图像i2分别利用sobel算子与roewa算子计算图像在不同尺度因子下的水平梯度gx,α与垂直梯度gy,α，从而得到图像在不同尺度下的梯度幅值magα与梯度方向oriα：[0123][0124][0125]2c)根据已经求出的梯度计算图像中每个像素点处的harris矩阵ch(x,y,α)：[0126][0127]其中，表示标准差为的高斯核，*表示卷积；[0128]2d)利用ch(x,y,α)计算图像中每个像素点处的harris响应值rh(x,y,α)：[0129]rh(x,y,α)＝det(ch(x,y,α))-d·tr(ch(x,y,α))2；[0130]其中，d为任意大小的参数，一般为0.04～0.06之间，det表示行列式的值，tr表示行列式的迹。[0131]2e)在尺度空间图像的同一层中，将每个像素点的harris响应值与它周围的8邻域内的像素点的响应值，以及一个全局阈值dh分别作比较，若在邻域内中心点的响应值是最大的，并且大于全局阈值dh，那么该点为检测到的关键点，其中，全局阈值dh取0.85。[0132]本发明中，步骤1a)、1b)和1c)中在找到的关键点处建立os-sift描述子的操作步骤如下：[0133]3a)利用梯度幅值magα与梯度方向oriα，以关键点为圆心，取半径为6α的圆形邻域，在圆形邻域内，将0～360°平均分为18份，每份代表20°的范围，横坐标表示梯度方向角，纵坐标表示梯度幅值，遍历该邻域内的所有像素点，对邻域内的每一个像素点，首先找到该像素点梯度方向角对应的直方图幅值柱，然后，将该像素点的梯度幅值累加在该点所在方向角的直方图幅值柱上，从而得到该关键点的主方向直方图，并对直方图进行平滑处理；[0134]3b)将平滑后的直方图峰值对应的方向角，作为该关键点的主方向，大于峰值能量80％的柱对应的方向角作为辅方向。[0135]3c)以关键点为圆心，分别取半径为rmax＝8α、12α以及16α的圆形邻域，将每个圆形邻域以主方向角为参考进行旋转，使特征描述子具有旋转不变性，然后，将每个圆形邻域由内到外划分成一个圆形和两个圆环，同心圆半径由内到外分别为0.25·rmax，0.75·rmax和rmax，再将两个圆环分别在0～360°度范围内以45°为间隔等分的划分成8个子区域，加上中心的圆形则将每个邻域一共划分为了17个子区域，三个圆形邻域一共有51个子区域，在每个子区域上将0～360°平均分为8份，建立梯度方向直方图，将三个圆形邻域的一共51个子区域的直方图向量拼接起来并归一化就生成了408维的os-sift描述子。示例性的，图3为生成单个圆形领域的描述子的图像。[0136]本发明中，步骤1a)、1b)、1c)和(2)中利用b-nndr算法找到参考图像与待配准图像之间的采样集与一致集的操作步骤如下：[0137]4a)在关键点空间内，以特征向量之间的欧氏距离衡量两个关键点的远近程度，越近则越相似；对于参考图像上的一个关键点，通过计算关键点的描述子之间的欧氏距离，在待配准图像上寻找距离该关键点最近与次近的关键点，最近距离与次近距离分别用dmin与dnd表示；如果则此关键点与距离它最近的关键点互为正确匹配点对；遍历参考图像上的每一个关键点，当阈值distratio取0.9时得到采样集当阈值distratio取0.999时得到一致集[0138]4b)遍历待配准图像上的每一个关键点，在参考图像上寻找与该关键点互为正确匹配点对的对应点。当阈值distratio取0.9时得到采样集当阈值distratio取0.999时得到一致集[0139]4c)参考图像与待配准图像之间的采样集与一致集为ch和cl：[0140][0141]∪表示取并集。[0142]本发明中，步骤1a)、1b)、1c)和(2)中利用fsc算法去除误匹配点对的操作步骤如下：[0143]5a)设置迭代次数为n，在第t次迭代过程中，从点对集合ch中随机挑选三对匹配点：t为1～n中的整数。[0144]5b)利用这三对点计算图像的变换模型θt。[0145]5c)利用得到的变换模型θt计算点对集合cl中的匹配点对ci的变换误差e(ci,θt)。[0146][0147][0148]其中，(xi,yi)表示匹配点对ci在参考图像上的关键点坐标，表示匹配点对ci在待配准图像上的对应点坐标；t((xi,yi),θt)表示利用变换模型θt对(xi,yi)变换后在待配准图像上的对应位置，e(ci,θt)表示在变换模型θt下，匹配点对ci的变换误差；[0149]5d)遍历点对集合cl中的每一对匹配点，将变换误差e＜3的所有匹配点对进行汇总，得到点对集合。[0150]5e)当算法迭代了n次后结束，此时得到了c1,c2,...,cn，一共n个点对集合，取出其中拥有最多点对数量的点对集合，作为fsc算法最终得到的匹配点对集合。[0151]本发明中，1b)中利用多尺度harris-affine算法计算每个关键点处的仿射变换矩阵的操作步骤如下：[0152]6a)对于参考图像i1上的一个关键点(xi,yi)，设置迭代次数k＝15并初始化形状自适应矩阵u(1)为单位矩阵e。[0153]6b)在第k次迭代过程中(此处，k为1～k中整数)，利用形状自适应矩阵u(k)对参考图像i1以及其上的关键点(xi,yi)进行变换，得到变换后的图像以及关键点坐标[0154][0155][0156]其中，t是利用形状自适应矩阵对图像或坐标进行变换的一种映射操作。[0157]6c)以为中心，在图像上取一个边长为4α的正方形区域w，其中，α表示该关键点所在的图像层尺度。[0158]6d)对参考图像i1与待配准图像i2分别利用sobel算子与roewa算子计算正方形区域w在尺度因子α下的水平梯度gx,α与垂直梯度gy,α，从而得到w的梯度幅值magα：[0159][0160]根据已经求出的水平和垂直梯度计算w中每个像素点处的harris矩阵ch(x,y,α)：[0161][0162]其中，表示标准差为的高斯核；[0163]利用ch(x,y,α)计算w中每个像素点处的harris响应值rh(x,y,α)：[0164]rh(x,y,α)＝det(ch(x,y,α))-d·tr(ch(x,y,α))2；其中，d为任意大小的参数，一般为0.04～0.06之间。[0165]将正方形区域w内拥有最大harris响应值的点作为新的关键点坐标[0166]6e)更新关键点在参考图像i1上的坐标：[0167][0168]6f)以为中心，在图像上重新选取一个边长为4α的正方形区域w_new，并计算新区域w_new内中心像素点的harris矩阵[0169]6g)更新形状自适应矩阵：[0170][0171]u(k+1)＝(μ(k))-1u(k)；[0172]之后，规范化更新后的形状自适应矩阵u(k+1)，使其最大特征值等于1。[0173]6h)计算第k次迭代时的收敛率：[0174][0175]其中，λmin(μ(k))与λmax(μ(k))分别表示矩阵μ(k)的最小特征值与最大特征值。[0176]6i)当ratio＜0.1时退出循环，并得到在关键点(xi,yi)处的仿射变换矩matrixi＝u(k)，否则，迭代次数加1并继续循环，利用更新后的形状自适应矩阵u(k+1)重新计算在关键点(xi,yi)处的收敛率ratio；若迭代了k次都无法满足ratio＜0.1的条件，则抛弃这个关键点。[0177]6j)对参考图像i1与待配准图像i2上的每一个关键点都进行相同的以上操作。[0178]本发明中，步骤1b)和1c)中在每个关键点处利用仿射变换矩阵对图像(例如，以参考图像i1为例)进行变换的操作步骤如下：[0179]7a)对于参考图像i1上的某个关键点(xi,yi)，利用仿射变换矩阵matrixi对参考图像i1的梯度幅值图magα、梯度方向图oriα以及关键点坐标(xi,yi)进行变换：[0180]magαt＝t(magα,matrixi)；[0181]oriαt＝t(oriα,matrixi)；[0182](xit,yit)＝t((xi,yi),matrixi)。[0183]本发明中，步骤1c)中利用多尺度mser算法提取图像中的关键点并计算每个关键点处的仿射变换矩阵的操作步骤如下：[0184]8a)选取一组尺度空间因子[α0,α1,...,αn-1]，其中，初始值α0＝2，αi＝α0*ki，(i∈[1,n-1])，k＝21/3，n＝4，即为高斯核的尺度，选取窗长w＝4α，利用高斯模糊核建立尺度空间；[0185]8b)对于每一层图像，如果该层图像是彩色图像，需要先将该彩色图像转换成灰度图像，之后，对于每一层图像，根据灰度值进行排序；预先为该层图像中的每个像素点分配一个节点，节点索引号即为该像素点对应的灰度值。根据像素点的排序结果，将它们逐一置入成分树中，置入顺序为每个像素点所对应的节点索引号；在置入的过程中，首先将该像素点置入，然后检查该像素点的四邻域位置，如果存在节点，则查找各自的根节点，把两个节点区域进行合并。当像素点全部置入成分树后，就获得了该层图像所对应的全部极值区域；其中，极值区域定义为：若某个区域内所有像素的灰度值都大于其边界像素的灰度值，则该区域定义为最大极值区域；若区域内所有像素的灰度值都小于其边界像素的灰度值，则定义为最小极值区域。[0186]8c)使用最大稳定判定条件来获取mser区域：如果q1,...,qi-1,qi,...为一系列相互包含的极值区域，即若极值区域qi*为最大稳定极值区域，当且仅当区域变化率q(i)＝|qi+δ-qi-δ|/|qi|在i*处取得局部极小值，其中，-表示取出两个区域的非公共部分，|·|表示区域内像素点的个数，下标i∈[0,255]表示灰度等级，δ表示微小的灰度变化量；[0187]8d)将不规则的极大稳定值区域近似拟合成一个椭圆区域：首先将极大稳定值区域的重心作为椭圆的中心，计算椭圆的中心：[0188]计算极大稳定值区域的几何零阶距和几何一阶距：[0189]m00＝∑ie(x,y)；[0190]m01＝∑yie(x,y)；[0191]m10＝∑xie(x,y)；[0192]其中m00、m01与m10分别为极大稳定值区域的几何零阶距和几何一阶距，ie(x,y)表示极大稳定值区域，从而可以得到椭圆的中心坐标，即mser算法检测到的关键点坐标(xc,yc)：[0193][0194][0195]计算极大稳定值区域的几何二阶距：[0196]其中，μ20＝∑(x-xc)2ie(x,y)，μ02＝∑(y-yc)2ie(x,y)，μ11＝∑(x-xc)(y-yc)ie(x,y)；[0197]计算几何二阶距的两个特征值：[0198][0199][0200]计算椭圆的长半轴w、短半轴l与长轴方向：[0201][0202][0203][0204]8e)利用椭圆的长半轴w、短半轴l与长轴方向计算关键点(xc,yc)处的仿射变换矩阵：[0205][0206]8f)将原始图像的灰度取反，重复上述操作；[0207]8g)对尺度空间中的每一层图像进行以上相同操作。[0208]本发明中，步骤(2)中建立参考图像i1与待配准图像的最小自相似图，通过极大值检测与非极大值抑制在最小自相似图上检测关键点的操作步骤如下：[0209]9a)对参考图像i1进行裁剪，构建子图像：以图像i1的中心像素点为中心，建立一个大小为lsub×wsub的搜索框，其中，lsub＝l-10，wsub＝w-10，l与w为参考图像i1的长和宽；当搜索框不移动时对参考图像i1按照搜索框所在位置进行裁剪得到中心子图像块subic，再将搜索框分别向0°、45°、90°、135°和180°偏移一个像素后，对参考图像i1按照当前搜索框所在位置进行裁剪得到偏移子图像块subi1、subi2、subi3、subi4和subi5；[0210]9b)通过以下公式获得在各个方向上偏移一个像素后的偏移子图像块subiθ：[0211][0212]θ＝180(o-1)/no,o＝1,2,...,no[0213]其中，θ表示偏移的方向，no表示将θ∈[0°,180°)分成no份，no可以根据实际需要设定。[0214]9c)计算中心子图像块subic上每个像素点(x,y)的权重值v(x,y)；[0215][0216][0217][0218]σl(x,y)＝ηl(x,y)-(μl(x,y))2；[0219]其中，μl、ηl和σl分别表示局部均值、均方值和方差值；(x,y)表示中心子图像块subic上的一个像素点；wl(x,y,rl)表示以(x,y)为中心，半径为rl的一个圆形邻域，当l取1时，r1等于2，当l取2时，r2等于4；nl表示圆形邻域wl内的像素点个数；[0220]9d)利用均值滤波器计算θ方向下的自相似图sθ：[0221][0222]其中，subic表示中心子图像块，subiθ表示向θ方向偏移一个像素后的偏移子图像块，v表示中心子图像块subic上每个像素点的权重值，w1表示半径为2的一个圆形邻域，表示对每个像素点取半径为2的圆形邻域进行均值滤波；[0223]9e)寻找最小自相似图smin：[0224][0225]9f)在最小自相似图smin上进行极大值检测以及非极大值抑制，寻找图像i1的关键点；[0226]9g)对待配准图像进行以上相同操作，寻找图像的关键点。[0227]本发明中，步骤(2)中基于二维log-garbor小波变换，计算参考图像i1与待配准图像在各个方向下的log-gabor响应，寻找最大响应的方向，并将方向的索引值取出建立最大索引图mim的操作步骤如下：[0228]10a)选取4个不同的尺度s＝1,2,3,4，6个不同的方向o＝0°,30°,60°,90°,120°,150°，分别对应索引号1、2、3、4、5和6，建立二维偶对称log-gabor滤波器leven(x,y,s,o)和二维奇对称log-gabor滤波器lodd(x,y,s,o)。[0229]10b)计算图像上的像素点(x,y)在o方向、s尺度下的响应幅值aso(x,y)：[0230]eso(x,y)＝i(x,y)*leven(x,y,s,o)，[0231]oso(x,y)＝i(x,y)*lodd(x,y,s,o)，[0232]其中，*表示卷积；[0233][0234]10c)将o方向下所有不同尺度的响应幅值进行累加，得到像素点(x,y)在o方向下的响应值：[0235][0236]10d)取出最大响应对应的方向索引值作为点(x,y)处的值，遍历图像上每一个像素点，建立最大索引图mim。[0237]本发明中，步骤(2)中基于检测到的关键点以及最大索引图mim建立描述子的操作步骤如下：[0238]11a)对于最大索引图mim上的每一个关键点，以其为中心，裁剪出一个大小为96×96的图像块，再将图像块划分成6×6个子网格，在每一个网格内，建立直方图，横坐标表示方向的索引号，取值范围为1～6，纵坐标表示该方向的索引号出现的次数，如此，可以将每个子网格转化为一个六维直方图向量，最后，将36个子区域的六维直方图向量拼接起来并进行归一化，就生成了该关键点的216维的描述子。[0239]示例性的，如图4所示，对于光学图像和粗矫正后sar图像，一方面通过建立最小自相似图并利用极大值检测与非极大值抑制检测关键点，另一方面通过计算多方向下的log-garbor响应并取出最大响应的方向索引值建立索引图；之后，基于检测到的关键点在索引图上建立描述子；最后，根据关键点的描述子对光学图像和粗矫正后sar图像进行匹配，得到大量正确匹配点对。[0240]本发明中，步骤(3)中将点对集合so、sh、sm与sr进行汇总后，利用los-flow算法进行精配的操作步骤如下：[0241]12a)利用点对集合so、sh、sm计算得到的待配准图像i2到参考图像i1的变换模型γ，使用变换模型γ对待配准图像i2进行变换，得到变换后的图像[0242]对点对集合so、sh、sm中在待配准图像i2上的一个点(xs,ys)，利用变换模型γ计算点(xs,ys)在变换后的图像上的对应位置坐标[0243]将点对集合so、sh、sm中的所有点都进行以上相同的变换，从而得到新的点对集合中的所有点都进行以上相同的变换，从而得到新的点对集合[0244]12b)对于图像i1与图像分别利用sobel算子与roewa算子计算在尺度因子α＝2下的水平梯度gx与垂直梯度gy，从而得到图像的梯度幅值mag与梯度方向ori：[0245][0246][0247]对于图像i1与上的一个像素点，以该像素点为圆心，取半径为rmax＝24的圆形邻域，并将此圆形邻域由内到外划分成一个圆形和两个圆环，同心圆半径由内到外分别为0.25·rmax，0.75·rmax和rmax，再将两个圆环分别在0～360°度范围内以45°为间隔等分的划分成8个子区域，加上中心的圆形则将此邻域一共划分为了17个子区域，在每个子区域内，将0～360°平均分为8份，每份代表45°的范围，横坐标表示梯度方向角，纵坐标表示梯度幅值，遍历该子区域内的所有像素点，对区域内的每一个像素点，首先找到该像素点梯度方向角对应的直方图幅值柱，然后，将该像素点的梯度幅值累加在该点所在方向角的直方柱上，从而得到该点的梯度方向直方图，并将其转化为一个八维直方图向量，17个子区域的直方图向量拼接起来并进行归一化后，就生成了该像素点的136维的os-sift描述子；[0248]通过分别计算图像i1与上每个像素点的136维os-sift描述子，就形成了图像i1的os-sift描述子图i1_desc，以及图像的os-sift描述子图[0249]12c)待精配点对集合为和sr的并集，对于待精配点对集合中在图像i1上的一个点(xi,yi)，以点(xi,yi)为中心，在os-sift描述子图i1_desc上取出一个边长为2·rlf+1的局部正方形区域描述子图像i1_squ，对点(xi,yi)在图像上的对应点(xit，yit)进行相同的操作，得到局部正方形区域描述子图像为其中，rlf取61；[0250]12d)将得到的i1_squ和代入损失函数e(w)，计算出光流矢量w，损失函数e(w)为：[0251][0252]其中，p＝(x,y)表示局部正方形区域描述子图像上的某个像素点，w(p)＝(u(p),v(p))表示p点的光流矢量，其中，u(p)表示p点在水平方向上的偏移量，v(p)表示p点在垂直方向上的偏移量；ε表示以p点为中心，加上其相邻8个点组成的区域，q表示该区域中不包含中心点p的某个点；参数η和α分别取0.001与0.03，参数t和d分别取0.1和0.6；上述(1)式为数据项，约束了沿着光流矢量w(p)的os-sift描述子要相互匹配；(2)式为小位移项，约束了在没有其他可用信息的情况下，光流矢量要尽可能的小；(3)式为平滑项，约束了相邻像素的光流矢量要相似；[0253]12e)计算图像上点(xit,yit)的新坐标(xit_new，yit_new)：[0254]xit_new＝xit+u(rlf+1,rlf+1)；[0255]yit_new＝yit+v(rlf+1,rlf+1)。[0256]12f)遍历待精配点对集合中的每一对匹配点后得到了更加准确的点对集合[0257]本发明中，步骤(3)中利用lsor算法去除可能包含的误匹配点对的操作步骤如下：[0258]13a)计算点对集合中每一对匹配点连线的长度以及斜率，对其求和后平均，得到平均长度dist_ave与平均斜率slope_ave。[0259]13b)对于点对集合中的每一对匹配点，计算连线的长度disti与斜率slopei，通过设置阈值thd和ths进行筛选，保留满足筛选条件的匹配点对，得到筛选条件如下：[0260]|disti-dist_ave|＜thd；[0261]|slopei-slope_ave|＜ths；其中，thd取0.1，ths取5°。[0262]13c)对于点对集合中在图像的每一个点，利用变换模型γ-1计算该点在原图像i2上的对应点，得到点对集合sfinal。[0263]以下通过仿真实验数据，对本发明实施例的技术效果进行进一步说明。[0264]选取两组公开的数据集(第一组数据集出自论文“adeeptranslation(gan)basedchangedetectionnetworkforopticalandsarremotesensingimages”，xinghuali等人，2021；第二组数据集出自论文“self-supervisedkeypointdetectionandcross-fusionmatchingnetworksformultimodalremotesensingimageregistration”，liangzhili等人，2022)，并为了进一步验证本发明提出的针对大视角差异光学和sar遥感图像的配准方法的性能，分别对两组公开数据集中的某一张图像进行变换，模拟大视角差异下的配准环境。[0265]下表1给出了本发明的方法在两组数据集上的配准结果，并与现有的配准方法affine-sift算法(简称asift，出自论文“asift：anewframeworkforfullyaffineinvariantimagecomparison”，siamjournalonimagingsciences，j.m.moreal等人，2009)以及os-sift算法(出自论文“os-sift:arobustsift-likealgorithmforhigh-resolutionoptical-to-sarimageregistrationinsuburbanareas”，ieeetransactionsongeoscienceandremotesensing，yumingxiang等人，2018)进行了对比，表1中的precision表示精确率，即最后得到的点对集合中正确点对所占的百分比；rmse表示均方根误差；test1输入的参考图与待配准图分别为图5和图6，test2输入的参考图与待配准图分别为图7和图8。图9和图11分别为test1的配准后的棋盘图和test2的配准后的棋盘图；图10和图12分别是test1的正确匹配点对图和test2的正确匹配点对图。[0266]表1本发明方法与现有方法的配准性能对比[0267][0268]从表1和上述结果图中可以看出，本发明提出的方法具有更好的配准性能。[0269]以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明，不能认定本发明的具体实施只局限于这些说明。对于本发明所属
技术领域：
：的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干简单推演或替换，都应当视为属于本发明的保护范围。当前第1页12当前第1页12
技术特征：
1.一种针对大视角差异光学与sar遥感图像的配准方法，其特征在于，包括：(1)利用os-sift算法、多尺度harris-affine算法以及多尺度mser算法分别提取参考图像i1与待配准图像i2中的关键点，然后，建立os-sift描述子并进行匹配，最终得到了三种具有不同性质的匹配点对集合：未受视角变化影响的点对集合s
o
、具有仿射不变性的角点点对集合s
h
，以及具有仿射不变性的区域点对集合s
m
；(2)根据步骤(1)得到的所有匹配点对，计算全局变换模型γ，利用全局变换模型γ对图像i2进行变换，得到图像在参考图像i1与待配准图像上利用risfm算法进行匹配，得到匹配点对集合s
r
；其中，risfm算法具体步骤为：首先建立参考图像i1与待配准图像的最小自相似图，通过极大值检测与非极大值抑制在最小自相似图上检测关键点；之后，基于二维log-garbor小波变换，计算参考图像i1与待配准图像在各个方向下的log-gabor响应，寻找最大响应的方向，并将方向的索引值取出建立最大索引图mim；然后，基于检测到的关键点以及最大索引图mim建立描述子；最后，利用b-nndr算法找到参考图像i1与待配准图像之间的采样集与一致集，再经过fsc算法去除误匹配点对后，得到匹配点对集合s
r
；(3)将点对集合s
o
、s
h
、s
m
与s
r
进行汇总后，利用los-flow算法进行精配，再利用lsor算法去除其中可能包含的误匹配点对，最终得到匹配点对集合s
final
，利用s
final
重新计算全局变换模型。2.根据权利要求1所述的针对大视角差异光学与sar遥感图像的配准方法，其特征在于，步骤(1)包括：1a)未受视角变化影响的点对集合s
o
：利用os-sift算法提取图像中的关键点并建立os-sift描述子，使用b-nndr算法找到参考图像i1与待配准图像i2之间的采样集与一致集，再经过fsc算法去除误匹配点对后，最终得到匹配点对集合s
o
；1b)具有仿射不变性的角点点对集合s
h
：利用多尺度harris-affine算法提取图像中的关键点并得到每个关键点处的局部仿射变换矩阵，其中，提取图像中的关键点的步骤与os-sift算法相同；之后，在每个关键点处利用局部仿射变换矩阵对图像进行变换后建立os-sift描述子，使用b-nndr算法找到参考图像i1与待配准图像i2之间的采样集与一致集，再经过fsc算法去除误匹配点对后，最终得到匹配点对集合s
h
；1c)具有仿射不变性的区域点对集合s
m
：利用多尺度mser算法提取图像中的关键点，并得到每个关键点处的局部仿射变换矩阵，在每个关键点处利用局部仿射变换矩阵对图像进行变换后建立os-sift描述子，使用b-nndr算法找到参考图像i1与待配准图像i2之间的采样集与一致集，再经过fsc算法去除误匹配点对后，最终得到匹配点对集合s
m
。3.根据权利要求2所述的针对大视角差异光学与sar遥感图像的配准方法，其特征在于，步骤1a)和1b)中利用os-sift算法提取图像中的关键点的具体步骤包括：2a)选取一组指数加权尺度因子[α0,α1,...,α
n-1
]，其中，初始值α0＝2，α
i
＝α0*k
i
，(i∈[1,n-1])，k＝2
1/3
，n＝8；2b)根据α的尺度序列[α0,α1,...,α
n-1
]，对参考图像i1与待配准图像i2分别利用sobel算子与roewa算子计算图像在不同尺度因子下的水平梯度g
x,α
与垂直梯度g
y,α
，从而得到图像在不同尺度下的梯度幅值mag
α
与梯度方向ori
α
：
2c)根据已经求出的梯度计算图像中每个像素点处的harris矩阵c
h
(x,y,α)：其中，表示标准差为的高斯核，*表示卷积；2d)利用c
h
(x,y,α)计算图像中每个像素点处的harris响应值r
h
(x,y,α)：r
h
(x,y,α)＝det(c
h
(x,y,α))-d
·
tr(c
h
(x,y,α))2；其中，d为0.04～0.06之间的数值，det表示行列式的值，tr表示行列式的迹；2e)在尺度空间图像的同一层中，将每个像素点的harris响应值与它周围的8邻域内的像素点的响应值，以及一个全局阈值d
h
分别作比较，若在邻域内中心点的响应值是最大的，并且大于全局阈值d
h
，那么该点为检测到的关键点，其中，全局阈值d
h
取0.85。4.根据权利要求3所述的针对大视角差异光学与sar遥感图像的配准方法，其特征在于，步骤1a)、1b)和1c)中建立os-sift描述子的步骤具体包括：3a)利用梯度幅值mag
α
与梯度方向ori
α
，以关键点为圆心，取半径为6α的圆形邻域，在圆形邻域内，将0～360
°
平均分为18份，每份代表20
°
的范围，横坐标表示梯度方向角，纵坐标表示梯度幅值，遍历该邻域内的所有像素点，对邻域内的每一个像素点，首先找到该像素点梯度方向角对应的直方图幅值柱，然后，将该像素点的梯度幅值累加在该点所在方向角的直方图幅值柱上，从而得到该关键点的主方向直方图，并对直方图进行平滑处理；3b)将平滑后的直方图峰值对应的方向角，作为该关键点的主方向，大于峰值能量80％的柱对应的方向角作为辅方向；3c)以关键点为圆心，分别取半径为r
max
＝8α、12α以及16α的圆形邻域，将每个圆形邻域以主方向角为参考进行旋转，使特征描述子具有旋转不变性，然后，将每个圆形邻域由内到外划分成一个圆形和两个圆环，同心圆半径由内到外分别为0.25
·
r
max
，0.75
·
r
max
和r
max
，再将两个圆环分别在0～360
°
度范围内以45
°
为间隔等分的划分成8个子区域，加上中心的圆形则将每个邻域一共划分为了17个子区域，三个圆形邻域一共有51个子区域，在每个子区域上将0～360
°
平均分为8份，建立梯度方向直方图，将三个圆形邻域的一共51个子区域的直方图向量拼接起来并归一化就生成了408维的os-sift描述子。5.根据权利要求1所述的针对大视角差异光学与sar遥感图像的配准方法，其特征在于，步骤1a)、1b)、1c)和(2)中利用b-nndr算法找到参考图像与待配准图像之间的采样集与一致集的步骤具体包括：4a)在关键点空间内，以特征向量之间的欧氏距离衡量两个关键点的远近程度，越近则越相似；对于参考图像上的一个关键点，通过计算关键点的描述子之间的欧氏距离，在待配准图像上寻找距离该关键点最近与次近的关键点，最近距离与次近距离分别用d
min
与d
nd
表
示；如果则此关键点与距离它最近的关键点互为正确匹配点对；遍历参考图像上的每一个关键点，当阈值distratio取0.9时得到采样集当阈值distratio取0.999时得到一致集4b)遍历待配准图像上的每一个关键点，在参考图像上寻找与该关键点互为正确匹配点对的对应点。当阈值distratio取0.9时得到采样集当阈值distratio取0.999时得到一致集4c)参考图像与待配准图像之间的采样集与一致集为c
h
和c
l
：：∪表示取并集。6.根据权利要求5所述的针对大视角差异光学与sar遥感图像的配准方法，其特征在于，步骤1a)、1b)、1c)和(2)中利用fsc算法去除误匹配点对的步骤具体包括：5a)设置迭代次数为n，在第t次迭代过程中，从点对集合c
h
中随机挑选三对匹配点：t为1～n中的整数；5b)利用这三对点计算图像的变换模型θ
t
；5c)利用得到的变换模型θ
t
计算点对集合c
l
中的匹配点对c
i
的变换误差e(c
i
,θ
t
)；)；其中，(x
i
,y
i
)表示匹配点对c
i
在参考图像上的关键点坐标，表示匹配点对c
i
在待配准图像上的对应点坐标；t((x
i
,y
i
),θ
t
)表示利用变换模型θ
t
对(x
i
,y
i
)变换后在待配准图像上的对应位置，e(c
i
,θ
t
)表示在变换模型θ
t
下，匹配点对c
i
的变换误差；5d)遍历点对集合c
l
中的每一对匹配点，将变换误差e＜3的所有匹配点对进行汇总，得到点对集合c
t
；5e)当算法迭代了n次后结束，此时得到了c1,c2,...,c
n
，一共n个点对集合，取出其中拥有最多点对数量的点对集合，作为fsc算法最终得到的匹配点对集合。7.根据权利要求1所述的针对大视角差异光学与sar遥感图像的配准方法，其特征在于，步骤(2)中建立参考图像i1与待配准图像的最小自相似图，通过极大值检测与非极大值抑制在最小自相似图上检测关键点的步骤具体包括：9a)对参考图像i1进行裁剪，构建子图像：以图像i1的中心像素点为中心，建立一个大小为l
sub
×
w
sub
的搜索框，其中，l
sub
＝l-10，w
sub
＝w-10，l与w为参考图像i1的长和宽；当搜索框不移动时对参考图像i1按照搜索框所在位置进行裁剪得到中心子图像块subi
c
，再将搜索框分别向0
°
、45
°
、90
°
、135
°
和180
°
偏移一个像素后，对参考图像i1按照当前搜索框所在位置进行裁剪得到偏移子图像块subi1、subi2、subi3、subi4和subi5；
9b)通过以下公式获得在各个方向上偏移一个像素后的偏移子图像块subi
θ
：θ＝180(o-1)/n
o
,o＝1,2,...,n
o
；其中，θ表示偏移的方向，n
o
表示将θ∈[0
°
,180
°
)分成n
o
份；9c)计算中心子图像块subi
c
上每个像素点(x,y)的权重值v(x,y)；上每个像素点(x,y)的权重值v(x,y)；上每个像素点(x,y)的权重值v(x,y)；σ
l
(x,y)＝η
l
(x,y)-(μ
l
(x,y))2；其中，μ
l
、η
l
和σ
l
分别表示局部均值、均方值和方差值；(x,y)表示中心子图像块subi
c
上的一个像素点；w
l
(x,y,r
l
)表示以(x,y)为中心，半径为r
l
的一个圆形邻域，当l取1时，r1等于2，当l取2时，r2等于4；n
l
表示圆形邻域w
l
内的像素点个数；9d)利用均值滤波器计算θ方向下的自相似图s
θ
：其中，subi
c
表示中心子图像块，subi
θ
表示向θ方向偏移一个像素后的偏移子图像块，v表示中心子图像块subi
c
上每个像素点的权重值，w1表示半径为2的一个圆形邻域，表示对每个像素点取半径为2的圆形邻域进行均值滤波；9e)寻找最小自相似图s
min
：9f)在最小自相似图s
min
上进行极大值检测以及非极大值抑制，寻找图像i1的关键点；9g)对待配准图像进行以上相同操作，寻找图像的关键点。8.根据权利要求1所述的针对大视角差异光学与sar遥感图像的配准方法，其特征在于，步骤(2)中基于二维log-garbor小波变换，计算参考图像i1与待配准图像在各个方向下的log-gabor响应，寻找最大响应的方向，并将方向的索引值取出建立最大索引图mim的步骤具体包括：10a)选取4个不同的尺度s＝1,2,3,4，6个不同的方向o＝0
°
,30
°
,60
°
,90
°
,120
°
,150
°
，分别对应索引号1、2、3、4、5和6，建立二维偶对称log-gabor滤波器l
even
(x,y,s,o)和二维奇
对称log-gabor滤波器l
odd
(x,y,s,o)；10b)计算图像上的像素点(x,y)在o方向、s尺度下的响应幅值a
so
(x,y)：e
so
(x,y)＝i(x,y)*l
even
(x,y,s,o)，o
so
(x,y)＝i(x,y)*l
odd
(x,y,s,o)，其中，*表示卷积；10c)将o方向下所有不同尺度的响应幅值进行累加，得到像素点(x,y)在o方向下的响应值：10d)取出最大响应对应的方向索引值作为点(x,y)处的值，遍历图像上每一个像素点，建立最大索引图mim。9.根据权利要求1所述的针对大视角差异光学与sar遥感图像的配准方法，其特征在于，步骤(2)中基于检测到的关键点以及最大索引图mim建立描述子的步骤具体包括：11a)对于最大索引图mim上的每一个关键点，以其为中心，裁剪出一个大小为96
×
96的图像块，再将图像块划分成6
×
6个子网格，在每一个网格内，建立直方图，横坐标表示方向的索引号，取值范围为1～6，纵坐标表示该方向的索引号出现的次数，如此，可以将每个子网格转化为一个六维直方图向量，最后，将36个子区域的六维直方图向量拼接起来并进行归一化，就生成了该关键点的216维的描述子。10.根据权利要求1所述的针对大视角差异光学与sar遥感图像的配准方法，其特征在于，步骤(3)中将点对集合s
o
、s
h
、s
m
与s
r
进行汇总后，利用los-flow算法进行精配的步骤具体包括：12a)利用点对集合s
o
、s
h
、s
m
计算得到的待配准图像i2到参考图像i1的变换模型γ，使用变换模型γ对待配准图像i2进行变换，得到变换后的图像对点对集合s
o
、s
h
、s
m
中在待配准图像i2上的一个点(x
s
,y
s
)，利用变换模型γ计算点(x
s
,y
s
)在变换后的图像上的对应位置坐标将点对集合s
o
、s
h
、s
m
中的所有点都进行以上相同的变换，从而得到新的点对集合中的所有点都进行以上相同的变换，从而得到新的点对集合12b)对于图像i1与图像分别利用sobel算子与roewa算子计算在尺度因子α＝2下的水平梯度g
x
与垂直梯度g
y
，从而得到图像的梯度幅值mag与梯度方向ori：，从而得到图像的梯度幅值mag与梯度方向ori：对于图像i1与上的一个像素点，以该像素点为圆心，取半径为r
max
＝24的圆形邻域，并将此圆形邻域由内到外划分成一个圆形和两个圆环，同心圆半径由内到外分别为0.25
·
r
max
，0.75
·
r
max
和r
max
，再将两个圆环分别在0～360
°
度范围内以45
°
为间隔等分的划分成8个子区域，加上中心的圆形则将此邻域一共划分为了17个子区域，在每个子区域内，将0～360
°
平均分为8份，每份代表45
°
的范围，横坐标表示梯度方向角，纵坐标表示梯度幅值，遍历该子区域内的所有像素点，对区域内的每一个像素点，首先找到该像素点梯度方向角对应的直方图幅值柱，然后，将该像素点的梯度幅值累加在该点所在方向角的直方柱上，从而得到该点的梯度方向直方图，并将其转化为一个八维直方图向量，17个子区域的直方图向量拼接起来并进行归一化后，就生成了该像素点的136维的os-sift描述子；通过分别计算图像i1与上每个像素点的136维os-sift描述子，就形成了图像i1的os-sift描述子图i1_desc，以及图像的os-sift描述子图12c)待精配点对集合为和s
r
的并集，对于待精配点对集合中在图像i1上的一个点(x
i
,y
i
)，以点(x
i
,y
i
)为中心，在os-sift描述子图i1_desc上取出一个边长为2
·
r
lf
+1的局部正方形区域描述子图像i1_squ，对点(x
i
,y
i
)在图像上的对应点(x
it
,y
it
)进行相同的操作，得到局部正方形区域描述子图像为其中，r
lf
取61；12d)将得到的i1_squ和代入损失函数e(w)，计算出光流矢量w，损失函数e(w)为：其中，p＝(x,y)表示局部正方形区域描述子图像上的某个像素点，w(p)＝(u(p),v(p))表示p点的光流矢量，其中，u(p)表示p点在水平方向上的偏移量，v(p)表示p点在垂直方向上的偏移量；ε表示以p点为中心，加上其相邻8个点组成的区域，q表示该区域中不包含中心点p的某个点；参数η和α分别取0.001与0.03，参数t和d分别取0.1和0.6；上述(1)式为数据项，约束了沿着光流矢量w(p)的os-sift描述子要相互匹配；(2)式为小位移项，约束了在没有其他可用信息的情况下，光流矢量要尽可能的小；(3)式为平滑项，约束了相邻像素的光流矢量要相似；12e)计算图像上点(x
it
,y
it
)的新坐标(x
it
_new，y
it
_new)：x
it
_new＝x
it
+u(r
lf
+1,r
lf
+1)；y
it
_new＝y
it
+v(r
lf
+1,r
lf
+1)；12f)遍历待精配点对集合中的每一对匹配点后得到了更加准确的点对集合

技术总结
本发明公开了一种针对大视角差异光学和SAR遥感图像的配准方法，该方法利用OS-SIFT算法、多尺度Harris-Affine算法以及多尺度MSER算法分别提取参考图与待配准图之间的关键点，建立描述子并匹配；接下来利用前一阶段得到的所有匹配点对计算全局变换模型，对待配准图进行变换后，在原参考图与变换后的待配准图上利用RISFM算法进行匹配；最后汇总前两阶段得到的所有匹配点对，精配后再进行误匹配点对去除，计算新的全局变换模型。本发明的方法深入考虑了SAR图像中存在的相干斑噪声影响以及异源遥感图像之间存在的几何畸变和非线性辐射畸变带来的影响，能够显著提升大视角差异下光学和SAR遥感图像的配准性能。学和SAR遥感图像的配准性能。学和SAR遥感图像的配准性能。


技术研发人员：王英华 刘俊 李欢 刘宏伟
受保护的技术使用者：西安电子科技大学
技术研发日：2023.07.06
技术公布日：2023/10/15