基于特征点匹配的图像拼接matlab

比眉伴天荒 2023-01-22 05:54 36阅读 0赞

## 说明 ##

接上文《基于图像相位相关性进行相邻图像检测》，此处完成基于傅里叶-梅林变换和SIFT图像配准的学习笔记。  
先占个坑。。。待续。。。  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3dlaXhpbl80NjI1NzQ1OA_size_16_color_FFFFFF_t_70_pic_center]

### 5.18 ###

本阶段完成傅里叶-梅林变换实现图像配准的方法学习。  
图像配准方法主要分为三类：一种是灰度方法信息方法，另一种是基于特征的方法，可细分为特征点、直线段、边缘轮廓、特征结构以及矩不变统计特征等，还有一种就是基于变换域的方法，如相位相关、Walsh Transform等方法。傅里叶-梅林变换就是一种变换域的方法。  
简单地学习之后，大致了解了傅里叶-梅林变换实现图像配准的实现思路。  
关于对抗旋转：与图像的相位相关性去求解位移同理，该方法利用相位相关也实现了旋转角度的求解，但对于旋转来说复杂之处在于当两张图存在旋转关系时，两个坐标之间是相互耦合的，难以通过傅里叶反变换直接求解出旋转角。于是傅里叶梅林变换图像配准法中，采用了坐标系变换的技巧，将直角坐标系中的量投射到极坐标系中(ρ，θ)，此时，旋转角度作为一个单独变量是解耦的，再采取傅里叶反变换就能够求出相应的旋转角度。这个问题的解决是通过坐标变换将问题划归到与位移相同的解决思路上来。  
关于对抗缩放：利用了频域方法，即傅里叶尺度变换特性，见下式。  
![在这里插入图片描述][20210523161411765.png]

计算得到两张图像的相对尺度变化。

但从我个人设计实验结果来看，这种方法对于大范围的旋转和缩放表现欠佳，但对于大范围的稳定的图像之间变换矩阵的计算相对有效。

--------------------

### 5.19 ###

本阶段完成SIFT算法阅读和学习。

### 概况梳理（SIFT特征匹配技术讲义1） ###

#### 1.多尺度表示 ####

1.1 高斯滤波器  
高斯滤波器是一类根据高斯函数的形状来选择权值的线性平滑滤波器。高斯平滑滤波器对于抑制服从正态分布的噪声非常有效。  
![在这里插入图片描述][20210523161603954.png]

高斯分布参数Sigma决定了高斯函数的宽度。理论上，高斯分布在所有定义域上都有非负值，这就需要一个无限大的卷积核。实际上，仅需要取均值周围3倍标准差内的值，以外部份直接去掉即可。对于图像处理来说，常用二维零均值离散高斯函数作平滑滤波器，高斯函数的图形：  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3dlaXhpbl80NjI1NzQ1OA_size_16_color_FFFFFF_t_70]  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3dlaXhpbl80NjI1NzQ1OA_size_16_color_FFFFFF_t_70 1]  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3dlaXhpbl80NjI1NzQ1OA_size_16_color_FFFFFF_t_70 2]  
1.2 高斯滤波性质  
1）二维高斯函数具有旋转对称性，即滤波器在各个方向上的平滑程度是相同的．一般来说，一幅图像的边缘方向是事先不知道的，因此，在滤波前是无法确定一个方向上比另一方向上需要更多的平滑．旋转对称性意味着高斯平滑滤波器在后续边缘检测中不会偏向任一方向．

（2）高斯函数是单值函数．这表明，高斯滤波器用像素邻域的加权均值来代替该点的像素值，而每一邻域像素点权值是随该点与中心点的距离单调递减的．这一性质是很重要的，因为边缘是一种图像局部特征，如果平滑运算对离算子中心很远的像素点仍然有很大作用，则平滑运算会使图像失真．

（3）高斯函数的傅立叶变换频谱是单瓣的．正如下面所示，这一性质是高斯函数付立叶变换等于高斯函数本身这一事实的直接推论．图像常被不希望的高频信号所污染(噪声和细纹理)．而所希望的图像特征（如边缘），既含有低频分量，又含有高频分量．高斯函数付立叶变换的单瓣意味着平滑图像不会被不需要的高频信号所污染，同时保留了大部分所需信号．

（4）高斯滤波器宽度(决定着平滑程度)是由参数σ表征的，而且σ和平滑程度的关系是非常简单的．σ越大，高斯滤波器的频带就越宽，平滑程度就越好．通过调节平滑程度参数σ，可在图像特征过分模糊(过平滑)与平滑图像中由于噪声和细纹理所引起的过多的不希望突变量(欠平滑)之间取得折衷．

（5）由于高斯函数的可分离性，较大尺寸的高斯滤波器可以得以有效地实现．二维高斯函数卷积可以分两步来进行，首先将图像与一维高斯函数进行卷积，然后将卷积结果与方向垂直的相同一维高斯函数卷积．因此，二维高斯滤波的计算量随滤波模板宽度成线性增长而不是成平方增长．  
1.3 高斯滤波应用  
高斯滤波后图像被平滑的程度取决于标准差。它的输出是领域像素的加权平均，同时离中心越近的像素权重越高。因此，相对于均值滤波（mean filter）它的平滑效果更柔和，而且边缘保留的也更好。  
　　高斯滤波被用作为平滑滤波器的本质原因是因为它是一个低通滤波器，见下图。而且，大部分基于卷积平滑滤波器都是低通滤波器。  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3dlaXhpbl80NjI1NzQ1OA_size_16_color_FFFFFF_t_70 3]  
1.4 关于DoG（Difference-of-Gaussian）  
从特征检测说起。  
关于边缘检测，直观上认知，如果对于一个灰度图像中的边缘部分，其邻域上的变化是比较剧烈的，对应一个极值的概念，相应的对一阶偏导再进行求导，边缘部分的点二阶偏导趋于0。因此可以依靠两种方法检测边缘。  
与边缘检测类似地，如果想实现特征点（角点=斑点=兴趣点）的检测，事实上我们需要检测一个局部极值。  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3dlaXhpbl80NjI1NzQ1OA_size_16_color_FFFFFF_t_70 4]  
特征点（稳定点）：想法是首先找到图像中的一些“稳定点”，这样如果两幅图像有相同的“稳定点”，这样就能实现匹配。因此，SIFT/SURF算法的基础是稳定点。简单的思路是寻找极值点，当两个物体的大小比例不一样时，大图像的局部极值点在小图像的对应位置上有可能不是极值点。于是SIFT/SURF都采用图像金字塔的方法，每一个截面与原图像相似，这样两个金字塔中就有可能包含大小最近似的两个截面了。这样找到的特征点会比较多，经过一些处理后滤掉一些相对不稳定的点。  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3dlaXhpbl80NjI1NzQ1OA_size_16_color_FFFFFF_t_70 5]  
特征描述子：找到稳定点后，配准时我们需要知道一张图像中地某些点在另一张图像中的位置，但如何实现这种匹配呢？SIFT采用的方法是建立一个128维度的特征描述子去描述一个特征点，采用的方式是统计周围像素的梯度。这一方法是与我们人类的视觉认知是相符的，我们人眼感知图画时，对于高频的信息更为敏感，表现为梯度变化大的我们更容易关注。  
关键思路：  
·利用不同标准差的高斯滤波器对图像进行滤波，相邻标准差之间得到的结果做差，以此来对抗尺度缩放。我理解的是，单纯取极值的话，在两张存在缩放关系的图像中，可能并不能得到配对，因为极值邻域点占有的像素个数是不同的。但通过这种高斯滤波，也即加权平均的方式可以一定程度上实现将极值点“摊开”，原本占有一个像素点的极值点（假设不是噪声点），通过高斯滤波或许将其尺度放大，从而可以与相应尺度的图像实现配准。  
·利用梯度计算出一个主方向，并取主方向和与主方向垂直的方向建立坐标系来描述梯度，从而实现旋转不变性。我认为这是非常巧妙的，但同时也是非常耗费计算复杂度的。需要统计128维度的特征去构成一个特征描述符，特征的维  
数非常大，然后根据两张图象中特征点的特征描述符判断两个点是否是实际中的同一个点，此处采用的方法是计算两个特征描述符的相似度，采用欧氏距离等。  
·消除错配，关于这一点并没有做自习的了解，但是消除错配过程中存在随机性，这一点在我做拼接实现的时候是很直观的，几乎每一次运行都能得到不同的拼接结果。

### 实验 ###

对二十张图像进行拼接，预先是不知道相邻关系和拼接顺序的。  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3dlaXhpbl80NjI1NzQ1OA_size_16_color_FFFFFF_t_70 6]  
拼接效果如图  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3dlaXhpbl80NjI1NzQ1OA_size_16_color_FFFFFF_t_70 7]

### 部分重要代码 ###

function [ssp] = ssp(I1,I2)
    % 加黑边以便于相位相似度计算
    if sum(size(I1)) ~= sum(size(I2))
        Iex = zeros(size(I1));
        [x,y,~] = size(I2);
        Iex(1:x,1:y,:) = I2;
        Iex = uint8(Iex);
        I2 = Iex;
    end
    % 转灰度
    I1gray = rgb2gray(I1);
    I2gray = rgb2gray(I2);
    % 傅里叶变换
    FI1 = fft2(I1gray);
    FI2 = fft2(I2gray);  
    % 取相位差做傅里叶反变换
    ssp = ifft2((FI1.*conj(FI2))./sqrt((FI1.*conj(FI1)).*(FI2.*conj(FI2))));  %频率域反变换到空间域
    ssp = max(max(ssp));
    end

function [IterationChoose] = nextChoose(I,maintainingList,Q)
    % 功能说明:迭代过程中，根据相位相似性检测下一张需要拼接的图像
    % input: 
    % I:目前拼接完成的图像
    % maintainingList:剩余的未进行拼接的图像列表
    % 
    % output:
    % IterationChoose.nextPair:下一个拼接的图像
    % IterationChoose.maintainingList:下一步迭代中，剩余的未进行拼接的图像
    
    
    SSP = 0; % 相位相似性
    nextPair = 0;
    for i = 1:length(maintainingList)
        j = maintainingList(i);
        I2 = readimage(Q,j);
        tempSSP = ssp(I,I2);
        if tempSSP > SSP
            SSP = tempSSP;
            nextPair = j;
        end
    end
    maintainingList(maintainingList == nextPair) = [];% 删除下一个进行拼接的图像
    
    IterationChoose.nextPair = nextPair;
    IterationChoose.maintainingList = maintainingList;
    end

function [sspPair] = init(Q)
    % 为第一次拼接做选择;
    % 选择相位相似度最高的一对图像之一
    % Q:图片文件夹
    % sspPair(1):相位相似度最高的图像对之一
    % sspPair(2):相位相似度最高的图像对之二
    % sspPair(3):相位相似度ssp
    
    sspPair = [0 0 0];% 前两列表示相似对，第三列表示ssp
    for i = 1:19
        I1 = readimage(Q,i);
        for j = i+1:20
            I2 = readimage(Q,j);
            temp = ssp(I1,I2);
            if temp > sspPair(3)
                sspPair(1)=i;
                sspPair(2)=j;
                sspPair(3) = temp;
            end
        end
    end
    fprintf('相位相似度最高的两张图像为%d和%d\n相似度为%d\n',sspPair(1),sspPair(2),sspPair(3));
    end

function [panoramaCut] = edgeCut(panorama)
    % 为了图像拼接完成后不存在黑边，需要进行剪切
    % 分别在两个方向上找“黑边”，并且删除
    
    Ig = rgb2gray(panorama);
    blackInx = sum(Ig)==0;% 按列求和
    if ~isempty(blackInx)
        Inx = 1:size(Ig,2);
        Inx(blackInx)=[];
        panoramaCut = panorama(:,Inx,:);
    end
    
    % 另外一个方向上去黑边
    Ig = rgb2gray(panoramaCut);
    Ig = Ig';
    blackInx = sum(Ig)==0;% 按行求和
    if ~isempty(blackInx)
        Inx = 1:size(Ig,2);
        Inx(blackInx)=[];
        panoramaCut = panoramaCut(Inx,:,:);
    end
    end

详情见资源页：  
（将于2021.6.14上传）

参考  
一个发现图像里卷积的作用----平移算子  
https://www.cnblogs.com/yjphhw/p/10790876.html  
Matlab对图像进行傅里叶变换实例  
https://blog.csdn.net/weixin\_43637490/article/details/89196212  
基于傅里叶梅林变换的SURF算法  
https://www.docin.com/p-1721109734.html  
《An FFT-Based Technique for Translation,Rotation, and Scale-Invariant Image Registration》  
《Distinctive Image Features from Scale-Invariant Keypoints》  
图像变化与其频谱变化的关系及应用分析  
https://www.zhihu.com/search?type=content&q=%E5%82%85%E9%87%8C%E5%8F%B6%E6%A2%85%E6%9E%97%E5%8F%98%E6%8D%A2  
SIFT和SURF的比较  
https://blog.csdn.net/wgx571859177/article/details/79522982  
SIFT/SURF算法的深入剖析——谈SIFT的精妙与不足  
https://blog.csdn.net/cy513/article/details/4414352  
高斯滤波简介  
https://www.cnblogs.com/qiqibaby/p/5289977.html  
DOG算子  
https://blog.csdn.net/songzitea/article/details/8831909  
DoG和LoG算子  
https://zhuanlan.zhihu.com/p/49447503

[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3dlaXhpbl80NjI1NzQ1OA_size_16_color_FFFFFF_t_70_pic_center]: /images/20221020/adbcefa56094439f9412e361f873c918.png
[20210523161411765.png]: /images/20221020/c08f62ec7b304098a494506434c8b507.png
[20210523161603954.png]: /images/20221020/40f86d5f63e64355b16523bde01ce661.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3dlaXhpbl80NjI1NzQ1OA_size_16_color_FFFFFF_t_70]: /images/20221020/3564783ffd4e4e9ab00e337094c407fc.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3dlaXhpbl80NjI1NzQ1OA_size_16_color_FFFFFF_t_70 1]: /images/20221020/f911b0ddc878410fb50cba65ebdea14f.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3dlaXhpbl80NjI1NzQ1OA_size_16_color_FFFFFF_t_70 2]: /images/20221020/979c71b9875e443ca6b37283c1cceb15.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3dlaXhpbl80NjI1NzQ1OA_size_16_color_FFFFFF_t_70 3]: /images/20221020/ca555d38f18a4db7aa59eb702c2ea29e.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3dlaXhpbl80NjI1NzQ1OA_size_16_color_FFFFFF_t_70 4]: /images/20221020/9946715b77204b75b4be88faa905eb67.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3dlaXhpbl80NjI1NzQ1OA_size_16_color_FFFFFF_t_70 5]: /images/20221020/25b69dbcb607451eb80d2ed359616572.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3dlaXhpbl80NjI1NzQ1OA_size_16_color_FFFFFF_t_70 6]: /images/20221020/e3e1c8c0107d47c4a1b91d66573ea638.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3dlaXhpbl80NjI1NzQ1OA_size_16_color_FFFFFF_t_70 7]: /images/20221020/d703d20da33b4f6da6fd6ba438e2c01e.png