您的位置:首页 > 其它

最佳缝合线拼接算法

2017-10-09 15:27 204 查看
最佳缝合线算法有助于消除鬼影,以得到较好的视觉效果,这里转载该博主的文章,方便自己查阅学习,若有不妥,请联系本人。

原文:http://blog.csdn.net/wd1603926823/article/details/49536691

理论根据《图像拼接的改进算法》,据说这个算法可以消除重叠部分运动物体的重影和鬼影问题,所以就编下试试看,反正之前编的那种很老的取平均值法融合、渐入渐出(基于距离)融合、以及改进的三角函数权重融合都只是适合静态图像融合 有重影和鬼影问题  不适合有运动物体的图像融合,所以还是要最佳缝合线算法:看论文上就四步 很清晰也很好懂  结果自己写的时候才发现看起来很容易的也许编起来没那么容易 之前 想得太简单了。

缝合线算法:

function D=bestlinefusion(A,B)

%《图像拼接的改进算法》最佳缝合线算法 图像融合

%先根据之前得到的H矩阵计算重叠区域Rect

[H,W,k]=size(A);

rdata1=-118;

L=W+1+rdata1;

R=W;

n=R-L+1;

%计算得到重叠区域的差值图像  其实我不懂计算差值图像干嘛  只要计算重叠区域的大小就好了 为什么还要计算差值图 后面又没用到

Rect=zeros(H,n);

for i=1:H

    for j=L:R

        Rect(i,j-L+1)=A(i,j)-B(i,j-L+1);

    end

end

Rect=uint8(Rect);%这句要不要呢?

%最终融合图的大小

rdata1=-118;

rdata2=3;

Y=2*W+rdata1+1;

D=zeros(H,Y);

%放路径的矩阵

path=zeros(H,n);

%放强度值 每条路径的强度值strength=color^2+geometry

color=zeros(1,n);

geometry=zeros(1,n);

strength1=zeros(1,n);

strength2=zeros(1,n);

%计算第一行即初始化的强度值

for j=L:R

    y=j-L+1;

    color(y)=A(1,j)-B(1,y);

    if(y==1)

        Bxdao=B(1,y+1)+2*B(2,y+1);

        Bydao=B(2,y)+2*B(2,y+1);

        Aydao=2*A(2,j-1)+A(2,j)+2*A(2,j+1);

        Axdao=A(1,j+1)+2*A(2,j+1)-A(1,j-1)-2*A(2,j-1);

        geometry(y)=(Axdao-Bxdao)*(Aydao-Bydao);

        strength1(y)=color(y)^2+geometry(y);

        path(1,y)=y;

        continue;

    end

    if(j==R)

        Axdao=A(1,j-1)-2*A(2,j-1);

        Aydao=2*A(2,j-1)+A(2,j);

        Bxdao=B(1,y+1)+2*B(2,y+1)-B(1,y-1)-2*B(2,y-1);

        Bydao=2*B(2,y-1)+B(2,y)+2*B(2,y+1);

        geometry(y)=(Axdao-Bxdao)*(Aydao-Bydao);

        strength1(y)=color(y)^2+geometry(y);

        path(1,y)=y;

        continue;

    end

    Axdao=A(1,j+1)+2*A(2,j+1)-A(1,j-1)-2*A(2,j-1);

    Bxdao=B(1,y+1)+2*B(2,y+1)-B(1,y-1)-2*B(2,y-1);

    Aydao=2*A(2,j-1)+A(2,j)+2*A(2,j+1);

    Bydao=2*B(2,y-1)+B(2,y)+2*B(2,y+1);

    geometry(y)=(Axdao-Bxdao)*(Aydao-Bydao);

    strength1(y)=color(y)^2+geometry(y);

    path(1,y)=y;

end

color=zeros(1,n);

geometry=zeros(1,n);

small=0;

%开始扩展 向下一行 从第二行到倒数第二行 最后一行单独拿出来 像第一行一样 因为它的结构差值geometry不好算

for i=2:H-1

    %先把下一行的强度值全部计算出来 到时候需要比较哪三个就拿出哪三个

    for j=L:R

        x=i;

        y=j-L+1;

        color(y)=A(i,j)-B(x,y);

        if(y==1)

            Axdao=2*A(i-1,j+1)+A(i,j+1)+2*A(i+1,j+1)-2*A(i-1,j-1)-A(i,j-1)-2*A(i+1,j-1);

            Bxdao=2*B(x-1,y+1)+B(x,y+1)+2*B(x+1,y+1);

            Aydao=-2*A(i-1,j-1)-A(i-1,j)-2*A(i-1,j+1)+2*A(i+1,j-1)+A(i+1,j)+2*A(i+1,j+1);

            Bydao=-B(x-1,y)-2*B(x-1,y+1)+B(x+1,y)+2*B(x+1,y+1);

            geometry(y)=(Axdao-Bxdao)*(Aydao-Bydao);

            strength2(y)=color(y)^2+geometry(y);

            continue;

        end

        if(j==R)

            Axdao=-2*A(i-1,j-1)-A(i,j-1)-2*A(i+1,j-1);

            Bxdao=2*B(x-1,y+1)+B(x,y+1)+2*B(x+1,y+1)-2*B(x-1,y-1)-B(x,y-1)-2*B(x+1,y-1);

            Aydao=-2*A(i-1,j-1)-A(i-1,j)+2*A(i+1,j-1)+A(i+1,j);

            Bydao=-2*B(x-1,y-1)-B(x-1,y)-2*B(x-1,y+1)+2*B(x+1,y-1)+B(x+1,y)+2*B(x+1,y+1);

            geometry(y)=(Axdao-Bxdao)*(Aydao-Bydao);

            strength2(y)=color(y)^2+geometry(y);

           continue;

        end

        Axdao=2*A(i-1,j+1)+A(i,j+1)+2*A(i+1,j+1)-2*A(i-1,j-1)-A(i,j-1)-2*A(i+1,j-1);

        Bxdao=2*B(x-1,y+1)+B(x,y+1)+2*B(x+1,y+1)-2*B(x-1,y-1)-B(x,y-1)-2*B(x+1,y-1);

        Aydao=-2*A(i-1,j-1)-A(i-1,j)-2*A(i-1,j+1)+2*A(i+1,j-1)+A(i+1,j)+2*A(i+1,j+1);

        Bydao=-2*B(x-1,y-1)-B(x-1,y)-2*B(x-1,y+1)+2*B(x+1,y-1)+B(x+1,y)+2*B(x+1,y+1);

        geometry(y)=(Axdao-Bxdao)*(Aydao-Bydao);

        strength2(y)=color(y)^2+geometry(y);

    end

    for j=1:n

        if(path(i-1,j)==1)

           if(strength2(1)<strength2(2))

              strength1(j)=strength1(j)+strength2(1);

              path(i,j)=1;

           else

              strength1(j)=strength1(j)+strength2(2);

              path(i,j)=2;

           end

        else

            if(path(i-1,j)==n)

              if(strength2(n-1)<strength2(n))

                strength1(j)=strength1(j)+strength2(n-1);

                path(i,j)=n-1;

              else

                strength1(j)=strength1(j)+strength2(n);

                path(i,j)=n;

              end

            else

                if(strength2(path(i-1,j)-1)<strength2(path(i-1,j)))

                    if(strength2(path(i-1,j)-1)<strength2(path(i-1,j)+1))

                        small=strength2(path(i-1,j)-1);

                        path(i,j)=path(i-1,j)-1;

                    else

                        small=strength2(path(i-1,j)+1);

                        path(i,j)=path(i-1,j)+1;

                    end

                else

                    if(strength2(path(i-1,j))<strength2(path(i-1,j)+1))

                        small=strength2(path(i-1,j));

                        path(i,j)=path(i-1,j);

                    else

                        small=strength2(path(i-1,j)+1);

                        path(i,j)=path(i-1,j)+1;

                    end

                end

                strength1(j)=strength1(j)+small;

            end

        end

        small=0;

    end

    strength2=zeros(1,
4000
n);

    color=zeros(1,n);

    geometry=zeros(1,n);

end

%单独计算最后一行

i=H;

for j=L:R

        x=i;

        y=j-L+1;

        color(y)=A(i,j)-B(x,y);

        if(y==1)

            Axdao=2*A(i-1,j+1)+A(i,j+1)-2*A(i-1,j-1)-A(i,j-1);

            Aydao=-2*A(i-1,j-1)-A(i-1,j)-2*A(i-1,j+1);

            Bxdao=2*B(x-1,y+1)+B(x,y+1);

            Bydao=-B(x-1,y)-2*B(x-1,y+1);

            continue;

        end

        if(j==R)

            Bxdao=2*B(x-1,y+1)+B(x,y+1)-2*B(x-1,y-1)-B(x,y-1);

            Bydao=-2*B(x-1,y-1)-B(x-1,y)-2*B(x-1,y+1);

            Axdao=-2*A(i-1,j-1)-A(i,j-1);

            Aydao=-2*A(i-1,j-1)-A(i-1,j);

            continue;

        end 

        Axdao=2*A(i-1,j+1)+A(i,j+1)-2*A(i-1,j-1)-A(i,j-1);

        Bxdao=2*B(x-1,y+1)+B(x,y+1)-2*B(x-1,y-1)-B(x,y-1);

        Aydao=-2*A(i-1,j-1)-A(i-1,j)-2*A(i-1,j+1);

        Bydao=-2*B(x-1,y-1)-B(x-1,y)-2*B(x-1,y+1);

        geometry(y)=(Axdao-Bxdao)*(Aydao-Bydao);

        strength2(y)=color(y)^2+geometry(y);

end

for j=1:n

        if(path(i-1,j)==1)

           if(strength2(1)<strength2(2))

              strength1(j)=strength1(j)+strength2(1);

              path(i,j)=1;

           else

              strength1(j)=strength1(j)+strength2(2);

              path(i,j)=2;

           end

        else

            if(path(i-1,j)==n)

              if(strength2(n-1)<strength2(n))

                strength1(j)=strength1(j)+strength2(n-1);

                path(i,j)=n-1;

              else

                strength1(j)=strength1(j)+strength2(n);

                path(i,j)=n;

              end

            else

                if(strength2(path(i-1,j)-1)<strength2(path(i-1,j)))

                    if(strength2(path(i-1,j)-1)<strength2(path(i-1,j)+1))

                        small=strength2(path(i-1,j)-1);

                        path(i,j)=path(i-1,j)-1;

                    else

                        small=strength2(path(i-1,j)+1);

                        path(i,j)=path(i-1,j)+1;

                    end

                else

                    if(strength2(path(i-1,j))<strength2(path(i-1,j)+1))

                        small=strength2(path(i-1,j));

                        path(i,j)=path(i-1,j);

                    else

                        small=strength2(path(i-1,j)+1);

                        path(i,j)=path(i-1,j)+1;

                    end

                end

                strength1(j)=strength1(j)+small;

            end

        end

        small=0;

end        

%比较strength1里放的每条路径的强度值的总和 谁最小 就选path中对应的那一列的路径

[minzhi,minth]=min(strength1);

mypath=path(:,minth);

%mypath放的就是最佳缝合线选出的路径 这条路径坐标是参考图A 右边是目标图B

for i=1:H

    for j=1:mypath(i)+L-1

        D(i,j,1)=A(i,j,1);

        D(i,j,2)=A(i,j,2);

        D(i,j,3)=A(i,j,3);

    end

    for j=mypath(i)+L-1:Y-1

         x=i;

         y=j-L+1;

         D(i,j,1)=B(x,y,1);

         D(i,j,2)=B(x,y,2);

         D(i,j,3)=B(x,y,3);

    end

end

 D=uint8(D);

还是以那两幅图为例:





结果:

差值图像Rect:


可是最终得到的融合后的图像D怎么是下面这个鬼样子:

这就是选出来的最佳缝合线??戳瞎我的狗眼啊!这哪里最佳了??!!

貌似理论部分没编错啊 哪里有问题呢  还是运行最佳缝合线本来就应该是这样?为了验证下  在另一篇也讲最佳缝合线算法的论文里 《全景图像拼接关键技术研究》下载了两张图片,看找出来最佳缝合线是否和这个作者的近似或者一样

原图:

通过放进OPENCV求匹配对 有32对匹配对 然后将这些匹配对输出给MATLAB计算H变换矩阵 我用了30对坐标 得到H矩阵 再进行H矩阵变换
后:

然后用下渐入渐出融合 结果:

用下改进的三角函数权重的融合
结果

果然改进的肉眼看不出来效果  但它们两个都有重影 比较严重的重影!用下刚刚编的最佳缝合线算法看看 这是差值图像

然后这是最佳缝合线得到的D:

可以看到用最佳缝合线后
没有重影了  那棵树那里以及那棵大树旁边两棵树没有重影了  可以和渐入渐出以及改进的那个对比下   这个已经没有重影了   原来最佳缝合线还是有用的喔   

我想把这幅图中找到的最佳缝合线画出来 ,其实我放大后看得到那条最佳缝,我想想怎么画出来:

L=W+1+rdata1;

imshow(D)

hold on;

for i=1:137

     D(i,L+mypath(i),1)=255;

     D(i,L+mypath(i),2)=255;

     D(i,L+mypath(i),3)=255;

end

hold off;


这条白线就是用上面的最佳缝合线程序找到的最佳缝合线  感觉好像一道白色的闪电 上面看不太出来  用黑色画下:


这样就看得清了
这就是找到的最佳缝合线  的确是避开了那些重影物体  我是这样来写的 其实就是和上面的H矩阵不一样 这个计算出来平移距离是76  竖直方向其实有5个像素的平移量  但我懒得改上面的了   直接只改了rdata1=-76 所以和那个论文中找的最佳缝合线有点点差别  不过没多大关系  因为我指望着像《图像拼接的改进算法》里说的那样找到最佳缝合线后还要进行多分辨率拼接的   就像之前我编渐入渐出时初步得到H矩阵变换后的结果也是有点没对齐的,但是融合时用渐入渐出后它就完完全全对齐了   这个也是一样的   我猜想后续的多分辨率拼接也可以做到这样的效果
  虽然现在得到的最佳缝合线有点没对齐的样子  因为我没用竖直方向的平移量rdata2=5,其实加上一样的   但太麻烦了因为我再看一次程序 然后重新计算 我下次再去改!现在这样就行了 竖直方向没对齐的我就依靠后续的多分辨率算法好了   但当竖直方向的平移量rdata2很大时以及有旋转关系时 也要重新计算  不能只算一个rdata1   绝对不能这样  会影响后续多分辨率的效果!
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息