图像拼接_图像拼接python - CSDN
  • 图像拼接

    2020-03-31 19:55:14
    图像拼接算法一般基于如下思路: 对每张图像进行特征点提取; 特征点匹配; 特征点筛选 计算投影矩阵; 进行投影变换 进行图像简单拼接; 对重叠部分进行融合; 特征点提取 选择不同的特征点提取算法,会对...

    图像拼接算法一般基于如下思路:

    1. 对每张图像进行特征点提取;
    2. 特征点匹配;
    3. 特征点筛选
    4. 计算投影矩阵;
    5. 进行投影变换
    6. 进行图像简单拼接;
    7. 对重叠部分进行融合;

    特征点提取

    选择不同的特征点提取算法,会对拼接的结果产生影响。

    分别测试了sift,surf,orb对特征点进行提取,在其他步骤一样的情况下,最终的拼接效果sift≈surf>orb

    特征点匹配

    特征点筛选

    排序筛选:根据匹配度进行排序,选择前特定数量或者特定百分比的特征点作为优秀特征点。

    使用orb特征,当分别选择前30%和前50%特征点的最终效果。

    阈值筛选:当匹配度大于一定阈值,认为该点为优秀特征点。

    findhomography计算单应性矩阵(3*3)

    单应性矩阵(Homography)是一个从一张图像到另一张图像映射关系的转换矩阵(3*3)。

     

    以上图为例,点(x1,y1)到点(x2,y2)的单应性变换为:

    相比于使用getPerspectiveTransform根据4个点来获得投影变换矩阵,findHomography使用优化算法从一系列点中筛选出最优解。

    Mat homo = findHomography(imagePoints1, imagePoints2, CV_RANSAC);

     

    计算配准图的四个顶点坐标

    在执行投影变换之前,计算投影变换后四个顶点在拼接后图像的位置,进而得出存放拼接后图像的尺寸大小。

     

    按照单应性变换公式,求解顶点变换后的坐标位置。

    例如,以左上角顶点(0,0)为例,变换后的坐标为:

        v2[0] = 0;
        v2[1] = 0;
        v2[2] = 1;
        V2 = Mat(3, 1, CV_64FC1, v2);  
        V1 = Mat(3, 1, CV_64FC1, v1);  
        V1 = H * V2;
        x = v1[0] / v1[2];
        y = v1[1] / v1[2];

    执行投影变换

    对右边图像执行投影变换。

    图像简单拼接

    将左边图像覆盖叠加在右边图像上。

    左边为surf的拼接效果,右边为orb的拼接效果

    简单拼接完成后就可以看出最终效果了,简单拼接错位严重的话,说明特征点提取步骤就存在问题,融合完成后效果也不会好到哪里去。

    重叠部分融合

    左边为surf的融合效果,右边为orb的融合效果,右边存在明显的重影现象,栏杆区域重叠度不好。

    图像融合对重叠部分采用加权融合的方式,在重叠部分由前一幅图像慢慢过渡到第二幅图像,非重叠部分直接采用原图像。

    重叠部分的权值公式:

    其中x为当前处理点到重叠区域左边界距离。

    参考:

    https://www.cnblogs.com/skyfsm/p/7411961.html

    https://blog.csdn.net/weixin_42717395/article/details/85768313

     

    opencv中使用stitcher将图像拼接算法封装起来,便于调用

    https://blog.csdn.net/zhaocj/article/details/78989548

     

    展开全文
  • 图像拼接之MATLAB实现

    2018-11-24 22:31:44
    图像拼接是一项应用广泛的图像处理技术。根据特征点的相互匹配,可以将多张小视角的图像拼接成为一张大视角的图像,在广角照片合成、卫星照片处理、医学图像处理等领域都有应用。早期的图像拼接主要是运用像素值匹配...

    转自http://www.cnblogs.com/naive/p/3579610.html

    背景介绍

    图像拼接是一项应用广泛的图像处理技术。根据特征点的相互匹配,可以将多张小视角的图像拼接成为一张大视角的图像,在广角照片合成、卫星照片处理、医学图像处理等领域都有应用。早期的图像拼接主要是运用像素值匹配的方法。后来,人们分别在两幅图像中寻找拐点、边缘等稳定的特征,用特征匹配的方法拼接图像。本实验根据Matthew Brown (2005) 描述的方法,实现多张生活照的拼接。

     

    特征点捕捉 (Interest Point Detection)

    首先,拍摄两张场景有重合的照片。为了保证有足够多的公共特征点,照片的重合度应该保证在30%以上。将两张照片转换为灰度图像,对图像做σ=1的高斯模糊。在Matthew的文章中,他建立了一个图像金字塔,在不同尺度寻找Harris关键点。考虑到将要拼接的照片视野尺寸接近,故简化此步骤,仅在原图提取特征点。

    接下来用sobel算子计算图像在x、y两个方向亮度的梯度,用σ=1.5的高斯函数对梯度做平滑处理,减小噪点对亮度的影响。很容易发现,若我们求一小块区域内亮度的累加值,在图像变化平缓的区域上下左右移动窗口累加值的变化并不明显;在物体的边缘,沿着边缘方向的变化也不明显;而在关键点附近,轻微的移动窗口都会强烈改变亮度的累加值,如图1所示。

     

    图1 http://www.cse.psu.edu/~rcollins/CSE486/lecture06.pdf

    亮度的变化值可以用下面的公式计算得到:

            (1)

    其中,w(x, y) 是高斯函数的权重,I(x, y)是该点亮度的梯度。

    在计算时,上面的公式又可以近似为如下:

            (2)

    通过比较矩阵的特征值l1和l2,我们可以判断该点所处的状态。若l1>>l2或者l2<<l1,表示该点位于纵向或者横向的边缘;若l1和l2近似且值很小,表示该点位于平滑区域;若l1和l2近似但值很大,表示该点位于关键点。根据Harris and Stephens (1988) 的介绍,我们并不需要直接计算两个特征值,用R = Det(H)/Tr(H)2的值就可以反映两个特征值的比值,这样可以减少运算量。我们保留R > 2的点。除此之外,每个点的R和周围8邻域像素的R值比较,仅保留局部R值最大的点。最后,去除图片边界附近的关键点。

    至此,我们在两幅图片分别得到了一组关键点,如图2所示。

     

    图2 Harris Corner

     

    自适应非极大值抑制 (Adaptive Non-Maximal Suppression)

    由于上一步得到的关键点很多,直接计算会导致很大的运算量,也会增加误差。接下去就要去除其中绝大部分的关键点,仅保留一些特征明显点,且让关键点在整幅图像内分布均匀。Matthew发明了adaptive non-maximal suppression (ANMS) 方法来择优选取特定数量的关键点。

    ANMS的思想是有一个半径r,初始值为无限远。当r不断减小时,保留在半径r以内其它关键点R值均小于中心点R值的关键点,将其加入队列。队列内的关键点数达到预设值后停止搜索。

     

    Xi是上一步得到的关键点的2维坐标,G是所有关键点的集合,c=0.9。

    实际计算时,我们将上述过程相反。这里我设定每幅图像各提取500个关键点。首先找出整幅图片R值最大的关键点Rmax,加入队列,并且得到Rmax*0.9的值。遍历所有关键点,若该关键点xi的Ri> Rmax*0.9, 该点的半径设为无限远;若该关键点xi的Ri< Rmax*0.9,计算该点到离它最近的Rj>0.9R的点xi,记录两点间的距离ri。最后将所有r排序,找出r最大的500个点,如图3所示。

     

    图3 Harris corner after ANMS

     

    关键点的描述 (Feature Descriptor)

    关键点的描述方法有很多种,包括局部梯度描述、尺度不变特征变换 (SIFT、SUFT) 等等。因为生活照的旋转角度通常不超过15°,所以这里不考虑关键点的旋转不变性。

    对图像做适度的高斯模糊,以关键点为中心,取40x40像素的区域。将该区域降采样至8x8的大小,生成一个64维的向量。对向量做归一化处理。每个关键点都用一个64维的向量表示,于是每幅图像分别得到了一个500x64的特征矩阵。

     

    关键点的匹配

    首先,从两幅图片的500个特征点中筛选出配对的点。筛选的方法是先计算500个特征点两两之间的欧氏距离,按照距离由小到大排序。通常情况下选择距离最小的一对特征向量配对。Lowe(2004)认为,仅仅观察最小距离并不能有效筛选配对特征点,而用最小的距离和第二小的距离的比值可以很好的进行筛选。如图4所示, 使用距离的比值能够获得更高的true positive, 同时控制较低的false positive。我使用的阈值是r1/r2<0.5。经过筛选后的配对特征点 如图5所示

     

    图 4. 配对正确率和配对方法、阈值选择的关系

     

    图 5. 筛选后的配对特征点

    关键点的匹配使用Random Sample Consensus (RANSAC) 算法。以一幅图像为基准,每次从中随机选择8个点,在另一幅图像中找出配对的8个点。用8对点计算得到一个homography,将基准图中剩余的特征点按照homography变换投影到另一幅图像,统计配对点的个数。

    重复上述步骤2000次,得到准确配对最多的一个homography。至此,两幅图像的投影变换关系已经找到。

     

    新图像的合成

    在做图像投影前,要先新建一个空白画布。比较投影后两幅图像的2维坐标的上下左右边界,选取各个方向边界的最大值作为新图像的尺寸。同时,计算得到两幅图像的交叉区域。

    在两幅图像的交叉区域,按照cross dissolve的方法制作两块如图6所示的蒙版,3个通道的像素值再次区间内递减(递升)。

     

    效果展示

    下面展示几张照片拼接的效果图。

     

    图 7. 拼接完成的新图像

     

    图 8. 以左边照片为基准拼接

     

     

     

    附Matlab代码:

    function [output_image] = image_stitching(input_A, input_B)
    % -------------------------------------------------------------------------
    % 1. Load both images, convert to double and to grayscale.
    % 2. Detect feature points in both images.
    % 3. Extract fixed-size patches around every keypoint in both images, and
    % form descriptors simply by "flattening" the pixel values in each patch to
    % one-dimensional vectors.
    % 4. Compute distances between every descriptor in one image and every descriptor in the other image. 
    % 5. Select putative matches based on the matrix of pairwise descriptor
    % distances obtained above. 
    % 6. Run RANSAC to estimate (1) an affine transformation and (2) a
    % homography mapping one image onto the other. 
    % 7. Warp one image onto the other using the estimated transformation.
    % 8. Create a new image big enough to hold the panorama and composite the
    % two images into it. 
    % 
    % Input:
    % input_A - filename of warped image
    % input_B - filename of unwarped image
    % Output:
    % output_image - combined new image
    % 
    % Reference:
    % [1] C.G. Harris and M.J. Stephens, A combined corner and edge detector, 1988.
    % [2] Matthew Brown, Multi-Image Matching using Multi-Scale Oriented Patches.
    % 
    % zhyh8341@gmail.com
    
    % -------------------------------------------------------------------------
    
    % READ IMAGE, GET SIZE INFORMATION
    image_A = imread(input_A);
    image_B = imread(input_B);
    [height_wrap, width_wrap,~] = size(image_A);
    [height_unwrap, width_unwrap,~] = size(image_B);
    
    % CONVERT TO GRAY SCALE
    gray_A = im2double(rgb2gray(image_A));
    gray_B = im2double(rgb2gray(image_B));
    
    
    % FIND HARRIS CORNERS IN BOTH IMAGE
    [x_A, y_A, v_A] = harris(gray_A, 2, 0.0, 2);
    [x_B, y_B, v_B] = harris(gray_B, 2, 0.0, 2);
    
    % ADAPTIVE NON-MAXIMAL SUPPRESSION (ANMS)
    ncorners = 500;
    [x_A, y_A, ~] = ada_nonmax_suppression(x_A, y_A, v_A, ncorners);
    [x_B, y_B, ~] = ada_nonmax_suppression(x_B, y_B, v_B, ncorners);
    
    % EXTRACT FEATURE DESCRIPTORS
    sigma = 7;
    [des_A] = getFeatureDescriptor(gray_A, x_A, y_A, sigma);
    [des_B] = getFeatureDescriptor(gray_B, x_B, y_B, sigma);
    
    % IMPLEMENT FEATURE MATCHING
    dist = dist2(des_A,des_B);
    [ord_dist, index] = sort(dist, 2);
    % THE RATIO OF FIRST AND SECOND DISTANCE IS A BETTER CRETIA THAN DIRECTLY
    % USING THE DISTANCE. RATIO LESS THAN .5 GIVES AN ACCEPTABLE ERROR RATE.
    ratio = ord_dist(:,1)./ord_dist(:,2);
    threshold = 0.5;
    idx = ratio<threshold;
    
    x_A = x_A(idx);
    y_A = y_A(idx);
    x_B = x_B(index(idx,1));
    y_B = y_B(index(idx,1));
    npoints = length(x_A);
    
    
    % USE 4-POINT RANSAC TO COMPUTE A ROBUST HOMOGRAPHY ESTIMATE
    % KEEP THE FIRST IMAGE UNWARPED, WARP THE SECOND TO THE FIRST
    matcher_A = [y_A, x_A, ones(npoints,1)]'; %!!! previous x is y and y is x,
    matcher_B = [y_B, x_B, ones(npoints,1)]'; %!!! so switch x and y here.
    [hh, ~] = ransacfithomography(matcher_B, matcher_A, npoints, 10);
    
    % s = load('matcher.mat');
    % matcher_A = s.matcher(1:3,:);
    % matcher_B = s.matcher(4:6,:);
    % npoints = 60;
    % [hh, inliers] = ransacfithomography(matcher_B, matcher_A, npoints, 10);
    
    
    % USE INVERSE WARP METHOD
    % DETERMINE THE SIZE OF THE WHOLE IMAGE
    [newH, newW, newX, newY, xB, yB] = getNewSize(hh, height_wrap, width_wrap, height_unwrap, width_unwrap);
    
    [X,Y] = meshgrid(1:width_wrap,1:height_wrap);
    [XX,YY] = meshgrid(newX:newX+newW-1, newY:newY+newH-1);
    AA = ones(3,newH*newW);
    AA(1,:) = reshape(XX,1,newH*newW);
    AA(2,:) = reshape(YY,1,newH*newW);
    
    AA = hh*AA;
    XX = reshape(AA(1,:)./AA(3,:), newH, newW);
    YY = reshape(AA(2,:)./AA(3,:), newH, newW);
    
    % INTERPOLATION, WARP IMAGE A INTO NEW IMAGE
    newImage(:,:,1) = interp2(X, Y, double(image_A(:,:,1)), XX, YY);
    newImage(:,:,2) = interp2(X, Y, double(image_A(:,:,2)), XX, YY);
    newImage(:,:,3) = interp2(X, Y, double(image_A(:,:,3)), XX, YY);
    
    % BLEND IMAGE BY CROSS DISSOLVE
    [newImage] = blend(newImage, image_B, xB, yB);
    
    % DISPLAY IMAGE MOSIAC
    imshow(uint8(newImage));
    
    

    ------------------------------- other functions -------------------------

    function [xp, yp, value] = harris(input_image, sigma,thd, r)
    % Detect harris corner 
    % Input:
    % sigma - standard deviation of smoothing Gaussian
    % r - radius of region considered in non-maximal suppression
    % Output:
    % xp - x coordinates of harris corner points
    % yp - y coordinates of harris corner points
    % value - values of R at harris corner points
    
    % CONVERT RGB IMAGE TO GRAY-SCALE, AND BLUR WITH G1 KERNEL
    g1 = fspecial('gaussian', 7, 1);
    gray_image = imfilter(input_image, g1);
    
    % FILTER INPUT IMAGE WITH SOBEL KERNEL TO GET GRADIENT ON X AND Y
    % ORIENTATION RESPECTIVELY
    h = fspecial('sobel');
    Ix = imfilter(gray_image,h,'replicate','same');
    Iy = imfilter(gray_image,h','replicate','same');
    
    % GENERATE GAUSSIAN FILTER OF SIZE 6*SIGMA (± 3SIGMA) AND OF MINIMUM SIZE 1x1
    g = fspecial('gaussian',fix(6*sigma), sigma);
    
    Ix2 = imfilter(Ix.^2, g, 'same').*(sigma^2); 
    Iy2 = imfilter(Iy.^2, g, 'same').*(sigma^2);
    Ixy = imfilter(Ix.*Iy, g, 'same').*(sigma^2);
    
    % HARRIS CORNER MEASURE
    R = (Ix2.*Iy2 - Ixy.^2)./(Ix2 + Iy2 + eps); 
    % ANOTHER MEASUREMENT, USUALLY k IS BETWEEN 0.04 ~ 0.06
    % response = (Ix2.*Iy2 - Ixy.^2) - k*(Ix2 + Iy2).^2;
    
    % GET RID OF CORNERS WHICH IS CLOSE TO BORDER
    R([1:20, end-20:end], :) = 0;
    R(:,[1:20,end-20:end]) = 0;
    
    % SUPRESS NON-MAX 
    d = 2*r+1; 
    localmax = ordfilt2(R,d^2,true(d)); 
    R = R.*(and(R==localmax, R>thd));
    
    % RETURN X AND Y COORDINATES 
    [xp,yp,value] = find(R);
    
    function [newx, newy, newvalue] = ada_nonmax_suppression(xp, yp, value, n)
    % Adaptive non-maximun suppression 
    % For each Harris Corner point, the minimum suppression radius is the
    % minimum distance from that point to a different point with a higher 
    % corner strength. 
    % Input:
    % xp,yp - coordinates of harris corner points
    % value - strength of suppression
    % n - number of interesting points
    % Output:
    % newx, newy - new x and y coordinates after adaptive non-maximun suppression
    % value - strength of suppression after adaptive non-maximun suppression
    
    % ALLOCATE MEMORY
    % newx = zeros(n,1);
    % newy = zeros(n,1);
    % newvalue = zeros(n,1);
    
    if(length(xp) < n)
    newx = xp;
    newy = yp;
    newvalue = value;
    return;
    end
    
    radius = zeros(n,1);
    c = .9;
    maxvalue = max(value)*c;
    for i=1:length(xp)
    if(value(i)>maxvalue)
    radius(i) = 99999999;
    continue;
    else
    dist = (xp-xp(i)).^2 + (yp-yp(i)).^2;
    dist((value*c) < value(i)) = [];
    radius(i) = sqrt(min(dist));
    end
    end
    
    [~, index] = sort(radius,'descend');
    index = index(1:n);
    
    newx = xp(index);
    newy = yp(index);
    newvalue = value(index);
    
    function n2 = dist2(x, c)
    % DIST2	Calculates squared distance between two sets of points.
    % Adapted from Netlab neural network software:
    % http://www.ncrg.aston.ac.uk/netlab/index.php
    %
    %	Description
    %	D = DIST2(X, C) takes two matrices of vectors and calculates the
    %	squared Euclidean distance between them. Both matrices must be of
    %	the same column dimension. If X has M rows and N columns, and C has
    %	L rows and N columns, then the result has M rows and L columns. The
    %	I, Jth entry is the squared distance from the Ith row of X to the
    %	Jth row of C.
    %
    %
    %	Copyright (c) Ian T Nabney (1996-2001)
    
    [ndata, dimx] = size(x);
    [ncentres, dimc] = size(c);
    if dimx ~= dimc
    error('Data dimension does not match dimension of centres')
    end
    
    n2 = (ones(ncentres, 1) * sum((x.^2)', 1))' + ...
    ones(ndata, 1) * sum((c.^2)',1) - ...
    2.*(x*(c'));
    
    % Rounding errors occasionally cause negative entries in n2
    if any(any(n2<0))
    n2(n2<0) = 0;
    end
    
    function [descriptors] = getFeatureDescriptor(input_image, xp, yp, sigma)
    % Extract non-rotation invariant feature descriptors
    % Input:
    % input_image - input gray-scale image
    % xx - x coordinates of potential feature points
    % yy - y coordinates of potential feature points
    % output:
    % descriptors - array of descriptors
    
    % FIRST BLUR WITH GAUSSIAN KERNEL
    g = fspecial('gaussian', 5, sigma);
    blurred_image = imfilter(input_image, g, 'replicate','same');
    
    % THEN TAKE A 40x40 PIXEL WINDOW AND DOWNSAMPLE TO 8x8 PATCH
    npoints = length(xp);
    descriptors = zeros(npoints,64);
    
    for i = 1:npoints
    %pA = imresize( blurred_image(xp(i)-20:xp(i)+19, yp(i)-20:yp(i)+19), .2);
    patch = blurred_image(xp(i)-20:xp(i)+19, yp(i)-20:yp(i)+19);
    patch = imresize(patch, .2);
    descriptors(i,:) = reshape((patch - mean2(patch))./std2(patch), 1, 64); 
    end
    
    function [hh] = getHomographyMatrix(point_ref, point_src, npoints)
    % Use corresponding points in both images to recover the parameters of the transformation 
    % Input:
    % x_ref, x_src --- x coordinates of point correspondences
    % y_ref, y_src --- y coordinates of point correspondences
    % Output:
    % h --- matrix of transformation
    
    % NUMBER OF POINT CORRESPONDENCES
    x_ref = point_ref(1,:)';
    y_ref = point_ref(2,:)';
    x_src = point_src(1,:)';
    y_src = point_src(2,:)';
    
    % COEFFICIENTS ON THE RIGHT SIDE OF LINEAR EQUATIONS
    A = zeros(npoints*2,8);
    A(1:2:end,1:3) = [x_ref, y_ref, ones(npoints,1)];
    A(2:2:end,4:6) = [x_ref, y_ref, ones(npoints,1)];
    A(1:2:end,7:8) = [-x_ref.*x_src, -y_ref.*x_src];
    A(2:2:end,7:8) = [-x_ref.*y_src, -y_ref.*y_src];
    
    % COEFFICIENT ON THE LEFT SIDE OF LINEAR EQUATIONS
    B = [x_src, y_src];
    B = reshape(B',npoints*2,1);
    
    % SOLVE LINEAR EQUATIONS
    h = A\B;
    
    hh = [h(1),h(2),h(3);h(4),h(5),h(6);h(7),h(8),1];
    
    function [hh, inliers] = ransacfithomography(ref_P, dst_P, npoints, threshold);
    % 4-point RANSAC fitting
    % Input:
    % matcher_A - match points from image A, a matrix of 3xN, the third row is 1
    % matcher_B - match points from image B, a matrix of 3xN, the third row is 1
    % thd - distance threshold
    % npoints - number of samples
    % 
    % 1. Randomly select minimal subset of points
    % 2. Hypothesize a model
    % 3. Computer error function
    % 4. Select points consistent with model
    % 5. Repeat hypothesize-and-verify loop
    % 
    % Yihua Zhao 02-01-2014
    % zhyh8341@gmail.com
    
    ninlier = 0;
    fpoints = 8; %number of fitting points
    for i=1:2000
    rd = randi([1 npoints],1,fpoints);
    pR = ref_P(:,rd);
    pD = dst_P(:,rd);
    h = getHomographyMatrix(pR,pD,fpoints);
    rref_P = h*ref_P;
    rref_P(1,:) = rref_P(1,:)./rref_P(3,:);
    rref_P(2,:) = rref_P(2,:)./rref_P(3,:);
    error = (rref_P(1,:) - dst_P(1,:)).^2 + (rref_P(2,:) - dst_P(2,:)).^2;
    n = nnz(error<threshold);
    if(n >= npoints*.95)
    hh=h;
    inliers = find(error<threshold);
    pause();
    break;
    elseif(n>ninlier)
    ninlier = n;
    hh=h;
    inliers = find(error<threshold);
    end 
    end
    
    function [newH, newW, x1, y1, x2, y2] = getNewSize(transform, h2, w2, h1, w1)
    % Calculate the size of new mosaic
    % Input:
    % transform - homography matrix
    % h1 - height of the unwarped image
    % w1 - width of the unwarped image
    % h2 - height of the warped image
    % w2 - height of the warped image
    % Output:
    % newH - height of the new image
    % newW - width of the new image
    % x1 - x coordate of lefttop corner of new image
    % y1 - y coordate of lefttop corner of new image
    % x2 - x coordate of lefttop corner of unwarped image
    % y2 - y coordate of lefttop corner of unwarped image
    % 
    % Yihua Zhao 02-02-2014
    % zhyh8341@gmail.com
    %
    
    % CREATE MESH-GRID FOR THE WARPED IMAGE
    [X,Y] = meshgrid(1:w2,1:h2);
    AA = ones(3,h2*w2);
    AA(1,:) = reshape(X,1,h2*w2);
    AA(2,:) = reshape(Y,1,h2*w2);
    
    % DETERMINE THE FOUR CORNER OF NEW IMAGE
    newAA = transform\AA;
    new_left = fix(min([1,min(newAA(1,:)./newAA(3,:))]));
    new_right = fix(max([w1,max(newAA(1,:)./newAA(3,:))]));
    new_top = fix(min([1,min(newAA(2,:)./newAA(3,:))]));
    new_bottom = fix(max([h1,max(newAA(2,:)./newAA(3,:))]));
    
    newH = new_bottom - new_top + 1;
    newW = new_right - new_left + 1;
    x1 = new_left;
    y1 = new_top;
    x2 = 2 - new_left;
    y2 = 2 - new_top;
    
    
    function [newImage] = blend(warped_image, unwarped_image, x, y)
    % Blend two image by using cross dissolve
    % Input:
    % warped_image - original image
    % unwarped_image - the other image
    % x - x coordinate of the lefttop corner of unwarped image
    % y - y coordinate of the lefttop corner of unwarped image
    % Output:
    % newImage
    % 
    % Yihua Zhao 02-02-2014
    % zhyh8341@gmail.com
    %
    
    
    % MAKE MASKS FOR BOTH IMAGES 
    warped_image(isnan(warped_image))=0;
    maskA = (warped_image(:,:,1)>0 |warped_image(:,:,2)>0 | warped_image(:,:,3)>0);
    newImage = zeros(size(warped_image));
    newImage(y:y+size(unwarped_image,1)-1, x: x+size(unwarped_image,2)-1,:) = unwarped_image;
    mask = (newImage(:,:,1)>0 | newImage(:,:,2)>0 | newImage(:,:,3)>0);
    mask = and(maskA, mask);
    
    % GET THE OVERLAID REGION
    [~,col] = find(mask);
    left = min(col);
    right = max(col);
    mask = ones(size(mask));
    if( x<2)
    mask(:,left:right) = repmat(linspace(0,1,right-left+1),size(mask,1),1);
    else
    mask(:,left:right) = repmat(linspace(1,0,right-left+1),size(mask,1),1);
    end
    
    % BLEND EACH CHANNEL
    warped_image(:,:,1) = warped_image(:,:,1).*mask;
    warped_image(:,:,2) = warped_image(:,:,2).*mask;
    warped_image(:,:,3) = warped_image(:,:,3).*mask;
    
    % REVERSE THE ALPHA VALUE
    if( x<2)
    mask(:,left:right) = repmat(linspace(1,0,right-left+1),size(mask,1),1);
    else
    mask(:,left:right) = repmat(linspace(0,1,right-left+1),size(mask,1),1);
    end
    newImage(:,:,1) = newImage(:,:,1).*mask;
    newImage(:,:,2) = newImage(:,:,2).*mask;
    newImage(:,:,3) = newImage(:,:,3).*mask;
    
    newImage(:,:,1) = warped_image(:,:,1) + newImage(:,:,1);
    newImage(:,:,2) = warped_image(:,:,2) + newImage(:,:,2);
    newImage(:,:,3) = warped_image(:,:,3) + newImage(:,:,3);
    
    展开全文
  • matlab全景图像拼接技术
  • 图像拼接技术包括三大部分:特征点提取与匹配、图像配准、图像融合。 1、基于SRUF 的特征点的提取与匹配 为了使拼接具有良好的精度和鲁棒性,同时又使其具有较好的实时性,本实验采用SURF 算法完成图像序列特征点...
    图像的拼接技术包括三大部分:特征点提取与匹配、图像配准、图像融合。
    1、基于SRUF 的特征点的提取与匹配
    为了使拼接具有良好的精度和鲁棒性,同时又使其具有较好的实时性,本实验采用SURF 算法完成图像序列特征点的提取。
    SURF 算法又称快速鲁棒特征,借鉴了SIFT 中简化近似的思想,将DoH 中的高斯二阶微分模板进行了近似简化,使得模板对图像的滤波只需要进行几个简单的加减法运算,并且这种运算与滤波模板的尺寸无关。实验证明,SURF 算法较SIFT 在运算速度上要快3 倍左右,综合性能要优于SIFT 算法。
    SURF 特征点提取与描述主要包含4 个步骤:
    1)检测尺度空间极值。
    2)精炼特征点位置。
    3)计算特征点的描述信息。

    4)生成描述特征点的特征向量。

    SURF的匹配算法是通过计算两个特征点藐视算子之间的欧式距离得到的,即找出与特征点描述符pi欧式距离最近和次近的两个邻居特征点描述符qi'和qi'',然后计算pi与qi'以及pi与qi'' 两组描述符之间欧式距离的比值r。如果比值r小于规定阈值规则视为匹配成功,(pi,qi')点对则为图像序列中的一对匹配点,否则匹配失败。这种匹配方法简便快捷,但是会差生误匹配。在图像配准模块会使用RANSAC进行误匹配点的筛选。

    2、图像配准

    图像配准是一种确定待拼接图像间的重叠区域以及重叠位置的技术,它是整个图像拼接的核心。本节采用的是基于特征点的图像配准方法,即通过匹配点对构建图像序列之间的变换矩阵,从而完成全景图像的拼接。

    变换矩阵H求解是图像配准的核心,其求解的算法流程如下。

    1)检测每幅图像中特征点。

    2)计算特征点之间的匹配。

    3)计算图像间变换矩阵的初始值。

    4)迭代精炼H变换矩阵。

    5)引导匹配。用估计的H去定义对极线附近的搜索区域,进一步确定特征点的对应。

    6)重复迭代4)和5)直到对应点的数目稳定为止。

    设图像序列之间的变换为投影变换

    可用4组最佳匹配计算出H矩阵的8 个自由度参数hi=( i=0,1,...,7),并以此作为初始值。

    为了提高图像配准的精度,本节采用RANSAC算法对图像变换矩阵进行求解与精炼,达到了较好的图像拼接效果。RANSAC算法的基本原理可通过图9-1来描述。


    图9-1:RANSAC算法

    RANSAC算法的思想简单而巧妙:首先随机地选择两个点,这两个点确定了一条直线,并且称在这条直线的一定范围内的点为这条直线的支撑。这样的随机选择重复数次,然后,具有最大支撑集的直线被确认为是样本点集的拟合。在拟合的误差距离范围内的点被认为是内点,它们构成一致集,反之则为外点。根据算法描述,可以很快判断,如果只有少量外点,那么随机选取的包含外点的初始点集确定的直线不会获得很大的支撑,值得注意的是,过大比例的外点将导致RANSAC算法失败。在直线拟合的例子中,由点集确定直线至少需要两个点;而对于透视变换,这样的最小集合需要有4个点。

    图9-1中蓝色点属于内点(正确点),而第红色点属于外点(偏移点)。此时用最小二乘法拟合这组数据,实际的最佳拟合直线是那条穿越了最多点数的蓝色实线。

    3、图像合成

    根据图像间变换矩阵H,可以对相应图像进行变换以确定图像间的重叠区域,并将待融和图像映射到到一幅新的空白图像中形成拼接图。需要注意的是,由于普通的相机在拍摄照片时会自动选取曝光参数,这会使输入图像间存在亮度差异,导致拼接后的图像缝合线两端出现明显的明暗变化。因此,在融和过程中需要对缝合线进行处理。进行图像拼接缝合线处理的方法有很多种,如颜色插值和多分辨率样条技术等,本节采用了快速简单的加权平滑算法处理拼接缝问题。该算法的主要思想是:图像重叠区域中像素点的灰度值Pixel 由两幅图像中对应点的灰度值Pixel_L和_R加权平均得到,即Pixel=k×Pixel_L+(1-k)× Pixel_R,其中k是可调因子。图9-2为加权平滑算法图。


    图9-2:加权平滑算法

    通常情况下0<k<1,即在重叠区域中,沿图像1向图像2的方向,k由1渐变为0,从而实现重叠区域的平滑拼接。为使图像重叠区域中的点与两幅图像建立更大的相关性,令k=d1/(d1+d2),其中d1,d2分别表示重叠区域中的点到两幅图像重叠区域的左边界和右边界的距离。即使用公式Pixel = ×Pixel_L + ×Pixel_R进行缝合线处理。



    展开全文
  • 图像拼接技术

    2018-09-03 20:13:36
    图像拼接在实际的应用场景很广,比如无人机航拍,遥感图像等等,图像拼接是进一步做图像理解基础步骤,拼接效果的好坏直接影响接下来的工作,所以一个好的图像拼接算法非常重要。 再举一个身边的例子吧,你用你的...

    转载来源:http://www.cnblogs.com/skyfsm/p/7411961.html

    图像拼接在实际的应用场景很广,比如无人机航拍,遥感图像等等,图像拼接是进一步做图像理解基础步骤,拼接效果的好坏直接影响接下来的工作,所以一个好的图像拼接算法非常重要。

    再举一个身边的例子吧,你用你的手机对某一场景拍照,但是你没有办法一次将所有你要拍的景物全部拍下来,所以你对该场景从左往右依次拍了好几张图,来把你要拍的所有景物记录下来。那么我们能不能把这些图像拼接成一个大图呢?我们利用opencv就可以做到图像拼接的效果!

    比如我们有对这两张图进行拼接。

    从上面两张图可以看出,这两张图有比较多的重叠部分,这也是拼接的基本要求。

    那么要实现图像拼接需要那几步呢?简单来说有以下几步:

    1. 对每幅图进行特征点提取
    2. 对对特征点进行匹配
    3. 进行图像配准
    4. 把图像拷贝到另一幅图像的特定位置
    5. 对重叠边界进行特殊处理

    好吧,那就开始正式实现图像配准。

    第一步就是特征点提取。现在CV领域有很多特征点的定义,比如sift、surf、harris角点、ORB都是很有名的特征因子,都可以用来做图像拼接的工作,他们各有优势。本文将使用ORB和SURF进行图像拼接,用其他方法进行拼接也是类似的。

    基于SURF的图像拼接

    用SIFT算法来实现图像拼接是很常用的方法,但是因为SIFT计算量很大,所以在速度要求很高的场合下不再适用。所以,它的改进方法SURF因为在速度方面有了明显的提高(速度是SIFT的3倍),所以在图像拼接领域还是大有作为。虽说SURF精确度和稳定性不及SIFT,但是其综合能力还是优越一些。下面将详细介绍拼接的主要步骤。

    1.特征点提取和匹配

    特征点提取和匹配的方法我在上一篇文章《OpenCV探索之路(二十三):特征检测和特征匹配方法汇总》中做了详细的介绍,在这里直接使用上文所总结的SURF特征提取和特征匹配的方法。

    //提取特征点    
    SurfFeatureDetector Detector(2000);  
    vector<KeyPoint> keyPoint1, keyPoint2;
    Detector.detect(image1, keyPoint1);
    Detector.detect(image2, keyPoint2);
    
    //特征点描述,为下边的特征点匹配做准备    
    SurfDescriptorExtractor Descriptor;
    Mat imageDesc1, imageDesc2;
    Descriptor.compute(image1, keyPoint1, imageDesc1);
    Descriptor.compute(image2, keyPoint2, imageDesc2);
    
    FlannBasedMatcher matcher;
    vector<vector<DMatch> > matchePoints;
    vector<DMatch> GoodMatchePoints;
    
    vector<Mat> train_desc(1, imageDesc1);
    matcher.add(train_desc);
    matcher.train();
    
    matcher.knnMatch(imageDesc2, matchePoints, 2);
    cout << "total match points: " << matchePoints.size() << endl;
    
    // Lowe's algorithm,获取优秀匹配点
    for (int i = 0; i < matchePoints.size(); i++)
    {
        if (matchePoints[i][0].distance < 0.4 * matchePoints[i][1].distance)
        {
            GoodMatchePoints.push_back(matchePoints[i][0]);
        }
    }
    
    Mat first_match;
    drawMatches(image02, keyPoint2, image01, keyPoint1, GoodMatchePoints, first_match);
    imshow("first_match ", first_match);

    2.图像配准

    这样子我们就可以得到了两幅待拼接图的匹配点集,接下来我们进行图像的配准,即将两张图像转换为同一坐标下,这里我们需要使用findHomography函数来求得变换矩阵。但是需要注意的是,findHomography函数所要用到的点集是Point2f类型的,所有我们需要对我们刚得到的点集GoodMatchePoints再做一次处理,使其转换为Point2f类型的点集。

    vector<Point2f> imagePoints1, imagePoints2;
    
    for (int i = 0; i<GoodMatchePoints.size(); i++)
    {
        imagePoints2.push_back(keyPoint2[GoodMatchePoints[i].queryIdx].pt);
        imagePoints1.push_back(keyPoint1[GoodMatchePoints[i].trainIdx].pt);
    }

    这样子,我们就可以拿着imagePoints1, imagePoints2去求变换矩阵了,并且实现图像配准。值得注意的是findHomography函数的参数中我们选泽了CV_RANSAC,这表明我们选择RANSAC算法继续筛选可靠地匹配点,这使得匹配点解更为精确。

    //获取图像1到图像2的投影映射矩阵 尺寸为3*3  
    Mat homo = findHomography(imagePoints1, imagePoints2, CV_RANSAC);
    ////也可以使用getPerspectiveTransform方法获得透视变换矩阵,不过要求只能有4个点,效果稍差  
    //Mat   homo=getPerspectiveTransform(imagePoints1,imagePoints2);  
    cout << "变换矩阵为:\n" << homo << endl << endl; //输出映射矩阵     
    
    //图像配准  
    Mat imageTransform1, imageTransform2;
    warpPerspective(image01, imageTransform1, homo, Size(MAX(corners.right_top.x, corners.right_bottom.x), image02.rows));
    //warpPerspective(image01, imageTransform2, adjustMat*homo, Size(image02.cols*1.3, image02.rows*1.8));
    imshow("直接经过透视矩阵变换", imageTransform1);
    imwrite("trans1.jpg", imageTransform1);

    3. 图像拷贝

    拷贝的思路很简单,就是将左图直接拷贝到配准图上就可以了。

    //创建拼接后的图,需提前计算图的大小
    int dst_width = imageTransform1.cols;  //取最右点的长度为拼接图的长度
    int dst_height = image02.rows;
    
    Mat dst(dst_height, dst_width, CV_8UC3);
    dst.setTo(0);
    
    imageTransform1.copyTo(dst(Rect(0, 0, imageTransform1.cols, imageTransform1.rows)));
    image02.copyTo(dst(Rect(0, 0, image02.cols, image02.rows)));
    
    imshow("b_dst", dst);

    4.图像融合(去裂缝处理)

    从上图可以看出,两图的拼接并不自然,原因就在于拼接图的交界处,两图因为光照色泽的原因使得两图交界处的过渡很糟糕,所以需要特定的处理解决这种不自然。这里的处理思路是加权融合,在重叠部分由前一幅图像慢慢过渡到第二幅图像,即将图像的重叠区域的像素值按一定的权值相加合成新的图像。

    
    //优化两图的连接处,使得拼接自然
    void OptimizeSeam(Mat& img1, Mat& trans, Mat& dst)
    {
        int start = MIN(corners.left_top.x, corners.left_bottom.x);//开始位置,即重叠区域的左边界  
    
        double processWidth = img1.cols - start;//重叠区域的宽度  
        int rows = dst.rows;
        int cols = img1.cols; //注意,是列数*通道数
        double alpha = 1;//img1中像素的权重  
        for (int i = 0; i < rows; i++)
        {
            uchar* p = img1.ptr<uchar>(i);  //获取第i行的首地址
            uchar* t = trans.ptr<uchar>(i);
            uchar* d = dst.ptr<uchar>(i);
            for (int j = start; j < cols; j++)
            {
                //如果遇到图像trans中无像素的黑点,则完全拷贝img1中的数据
                if (t[j * 3] == 0 && t[j * 3 + 1] == 0 && t[j * 3 + 2] == 0)
                {
                    alpha = 1;
                }
                else
                {
                    //img1中像素的权重,与当前处理点距重叠区域左边界的距离成正比,实验证明,这种方法确实好  
                    alpha = (processWidth - (j - start)) / processWidth;
                }
    
                d[j * 3] = p[j * 3] * alpha + t[j * 3] * (1 - alpha);
                d[j * 3 + 1] = p[j * 3 + 1] * alpha + t[j * 3 + 1] * (1 - alpha);
                d[j * 3 + 2] = p[j * 3 + 2] * alpha + t[j * 3 + 2] * (1 - alpha);
    
            }
        }
    
    }

    多尝试几张,验证拼接效果

    测试一

    测试二

    测试三

    最后给出完整的SURF算法实现的拼接代码。

    #include "highgui/highgui.hpp"    
    #include "opencv2/nonfree/nonfree.hpp"    
    #include "opencv2/legacy/legacy.hpp"   
    #include <iostream>  
    
    using namespace cv;
    using namespace std;
    
    void OptimizeSeam(Mat& img1, Mat& trans, Mat& dst);
    
    typedef struct
    {
        Point2f left_top;
        Point2f left_bottom;
        Point2f right_top;
        Point2f right_bottom;
    }four_corners_t;
    
    four_corners_t corners;
    
    void CalcCorners(const Mat& H, const Mat& src)
    {
        double v2[] = { 0, 0, 1 };//左上角
        double v1[3];//变换后的坐标值
        Mat V2 = Mat(3, 1, CV_64FC1, v2);  //列向量
        Mat V1 = Mat(3, 1, CV_64FC1, v1);  //列向量
    
        V1 = H * V2;
        //左上角(0,0,1)
        cout << "V2: " << V2 << endl;
        cout << "V1: " << V1 << endl;
        corners.left_top.x = v1[0] / v1[2];
        corners.left_top.y = v1[1] / v1[2];
    
        //左下角(0,src.rows,1)
        v2[0] = 0;
        v2[1] = src.rows;
        v2[2] = 1;
        V2 = Mat(3, 1, CV_64FC1, v2);  //列向量
        V1 = Mat(3, 1, CV_64FC1, v1);  //列向量
        V1 = H * V2;
        corners.left_bottom.x = v1[0] / v1[2];
        corners.left_bottom.y = v1[1] / v1[2];
    
        //右上角(src.cols,0,1)
        v2[0] = src.cols;
        v2[1] = 0;
        v2[2] = 1;
        V2 = Mat(3, 1, CV_64FC1, v2);  //列向量
        V1 = Mat(3, 1, CV_64FC1, v1);  //列向量
        V1 = H * V2;
        corners.right_top.x = v1[0] / v1[2];
        corners.right_top.y = v1[1] / v1[2];
    
        //右下角(src.cols,src.rows,1)
        v2[0] = src.cols;
        v2[1] = src.rows;
        v2[2] = 1;
        V2 = Mat(3, 1, CV_64FC1, v2);  //列向量
        V1 = Mat(3, 1, CV_64FC1, v1);  //列向量
        V1 = H * V2;
        corners.right_bottom.x = v1[0] / v1[2];
        corners.right_bottom.y = v1[1] / v1[2];
    
    }
    
    int main(int argc, char *argv[])
    {
        Mat image01 = imread("g5.jpg", 1);    //右图
        Mat image02 = imread("g4.jpg", 1);    //左图
        imshow("p2", image01);
        imshow("p1", image02);
    
        //灰度图转换  
        Mat image1, image2;
        cvtColor(image01, image1, CV_RGB2GRAY);
        cvtColor(image02, image2, CV_RGB2GRAY);
    
    
        //提取特征点    
        SurfFeatureDetector Detector(2000);  
        vector<KeyPoint> keyPoint1, keyPoint2;
        Detector.detect(image1, keyPoint1);
        Detector.detect(image2, keyPoint2);
    
        //特征点描述,为下边的特征点匹配做准备    
        SurfDescriptorExtractor Descriptor;
        Mat imageDesc1, imageDesc2;
        Descriptor.compute(image1, keyPoint1, imageDesc1);
        Descriptor.compute(image2, keyPoint2, imageDesc2);
    
        FlannBasedMatcher matcher;
        vector<vector<DMatch> > matchePoints;
        vector<DMatch> GoodMatchePoints;
    
        vector<Mat> train_desc(1, imageDesc1);
        matcher.add(train_desc);
        matcher.train();
    
        matcher.knnMatch(imageDesc2, matchePoints, 2);
        cout << "total match points: " << matchePoints.size() << endl;
    
        // Lowe's algorithm,获取优秀匹配点
        for (int i = 0; i < matchePoints.size(); i++)
        {
            if (matchePoints[i][0].distance < 0.4 * matchePoints[i][1].distance)
            {
                GoodMatchePoints.push_back(matchePoints[i][0]);
            }
        }
    
        Mat first_match;
        drawMatches(image02, keyPoint2, image01, keyPoint1, GoodMatchePoints, first_match);
        imshow("first_match ", first_match);
    
        vector<Point2f> imagePoints1, imagePoints2;
    
        for (int i = 0; i<GoodMatchePoints.size(); i++)
        {
            imagePoints2.push_back(keyPoint2[GoodMatchePoints[i].queryIdx].pt);
            imagePoints1.push_back(keyPoint1[GoodMatchePoints[i].trainIdx].pt);
        }
    
    
    
        //获取图像1到图像2的投影映射矩阵 尺寸为3*3  
        Mat homo = findHomography(imagePoints1, imagePoints2, CV_RANSAC);
        ////也可以使用getPerspectiveTransform方法获得透视变换矩阵,不过要求只能有4个点,效果稍差  
        //Mat   homo=getPerspectiveTransform(imagePoints1,imagePoints2);  
        cout << "变换矩阵为:\n" << homo << endl << endl; //输出映射矩阵      
    
       //计算配准图的四个顶点坐标
        CalcCorners(homo, image01);
        cout << "left_top:" << corners.left_top << endl;
        cout << "left_bottom:" << corners.left_bottom << endl;
        cout << "right_top:" << corners.right_top << endl;
        cout << "right_bottom:" << corners.right_bottom << endl;
    
        //图像配准  
        Mat imageTransform1, imageTransform2;
        warpPerspective(image01, imageTransform1, homo, Size(MAX(corners.right_top.x, corners.right_bottom.x), image02.rows));
        //warpPerspective(image01, imageTransform2, adjustMat*homo, Size(image02.cols*1.3, image02.rows*1.8));
        imshow("直接经过透视矩阵变换", imageTransform1);
        imwrite("trans1.jpg", imageTransform1);
    
    
        //创建拼接后的图,需提前计算图的大小
        int dst_width = imageTransform1.cols;  //取最右点的长度为拼接图的长度
        int dst_height = image02.rows;
    
        Mat dst(dst_height, dst_width, CV_8UC3);
        dst.setTo(0);
    
        imageTransform1.copyTo(dst(Rect(0, 0, imageTransform1.cols, imageTransform1.rows)));
        image02.copyTo(dst(Rect(0, 0, image02.cols, image02.rows)));
    
        imshow("b_dst", dst);
    
    
        OptimizeSeam(image02, imageTransform1, dst);
    
    
        imshow("dst", dst);
        imwrite("dst.jpg", dst);
    
        waitKey();
    
        return 0;
    }
    
    
    //优化两图的连接处,使得拼接自然
    void OptimizeSeam(Mat& img1, Mat& trans, Mat& dst)
    {
        int start = MIN(corners.left_top.x, corners.left_bottom.x);//开始位置,即重叠区域的左边界  
    
        double processWidth = img1.cols - start;//重叠区域的宽度  
        int rows = dst.rows;
        int cols = img1.cols; //注意,是列数*通道数
        double alpha = 1;//img1中像素的权重  
        for (int i = 0; i < rows; i++)
        {
            uchar* p = img1.ptr<uchar>(i);  //获取第i行的首地址
            uchar* t = trans.ptr<uchar>(i);
            uchar* d = dst.ptr<uchar>(i);
            for (int j = start; j < cols; j++)
            {
                //如果遇到图像trans中无像素的黑点,则完全拷贝img1中的数据
                if (t[j * 3] == 0 && t[j * 3 + 1] == 0 && t[j * 3 + 2] == 0)
                {
                    alpha = 1;
                }
                else
                {
                    //img1中像素的权重,与当前处理点距重叠区域左边界的距离成正比,实验证明,这种方法确实好  
                    alpha = (processWidth - (j - start)) / processWidth;
                }
    
                d[j * 3] = p[j * 3] * alpha + t[j * 3] * (1 - alpha);
                d[j * 3 + 1] = p[j * 3 + 1] * alpha + t[j * 3 + 1] * (1 - alpha);
                d[j * 3 + 2] = p[j * 3 + 2] * alpha + t[j * 3 + 2] * (1 - alpha);
    
            }
        }
    
    }

    基于ORB的图像拼接

    利用ORB进行图像拼接的思路跟上面的思路基本一样,只是特征提取和特征点匹配的方式略有差异罢了。这里就不再详细介绍思路了,直接贴代码看效果。

    #include "highgui/highgui.hpp"    
    #include "opencv2/nonfree/nonfree.hpp"    
    #include "opencv2/legacy/legacy.hpp"   
    #include <iostream>  
    
    using namespace cv;
    using namespace std;
    
    void OptimizeSeam(Mat& img1, Mat& trans, Mat& dst);
    
    typedef struct
    {
        Point2f left_top;
        Point2f left_bottom;
        Point2f right_top;
        Point2f right_bottom;
    }four_corners_t;
    
    four_corners_t corners;
    
    void CalcCorners(const Mat& H, const Mat& src)
    {
        double v2[] = { 0, 0, 1 };//左上角
        double v1[3];//变换后的坐标值
        Mat V2 = Mat(3, 1, CV_64FC1, v2);  //列向量
        Mat V1 = Mat(3, 1, CV_64FC1, v1);  //列向量
    
        V1 = H * V2;
        //左上角(0,0,1)
        cout << "V2: " << V2 << endl;
        cout << "V1: " << V1 << endl;
        corners.left_top.x = v1[0] / v1[2];
        corners.left_top.y = v1[1] / v1[2];
    
        //左下角(0,src.rows,1)
        v2[0] = 0;
        v2[1] = src.rows;
        v2[2] = 1;
        V2 = Mat(3, 1, CV_64FC1, v2);  //列向量
        V1 = Mat(3, 1, CV_64FC1, v1);  //列向量
        V1 = H * V2;
        corners.left_bottom.x = v1[0] / v1[2];
        corners.left_bottom.y = v1[1] / v1[2];
    
        //右上角(src.cols,0,1)
        v2[0] = src.cols;
        v2[1] = 0;
        v2[2] = 1;
        V2 = Mat(3, 1, CV_64FC1, v2);  //列向量
        V1 = Mat(3, 1, CV_64FC1, v1);  //列向量
        V1 = H * V2;
        corners.right_top.x = v1[0] / v1[2];
        corners.right_top.y = v1[1] / v1[2];
    
        //右下角(src.cols,src.rows,1)
        v2[0] = src.cols;
        v2[1] = src.rows;
        v2[2] = 1;
        V2 = Mat(3, 1, CV_64FC1, v2);  //列向量
        V1 = Mat(3, 1, CV_64FC1, v1);  //列向量
        V1 = H * V2;
        corners.right_bottom.x = v1[0] / v1[2];
        corners.right_bottom.y = v1[1] / v1[2];
    
    }
    
    int main(int argc, char *argv[])
    {
        Mat image01 = imread("t1.jpg", 1);    //右图
        Mat image02 = imread("t2.jpg", 1);    //左图
        imshow("p2", image01);
        imshow("p1", image02);
    
        //灰度图转换  
        Mat image1, image2;
        cvtColor(image01, image1, CV_RGB2GRAY);
        cvtColor(image02, image2, CV_RGB2GRAY);
    
    
        //提取特征点    
        OrbFeatureDetector  surfDetector(3000);  
        vector<KeyPoint> keyPoint1, keyPoint2;
        surfDetector.detect(image1, keyPoint1);
        surfDetector.detect(image2, keyPoint2);
    
        //特征点描述,为下边的特征点匹配做准备    
        OrbDescriptorExtractor  SurfDescriptor;
        Mat imageDesc1, imageDesc2;
        SurfDescriptor.compute(image1, keyPoint1, imageDesc1);
        SurfDescriptor.compute(image2, keyPoint2, imageDesc2);
    
        flann::Index flannIndex(imageDesc1, flann::LshIndexParams(12, 20, 2), cvflann::FLANN_DIST_HAMMING);
    
        vector<DMatch> GoodMatchePoints;
    
        Mat macthIndex(imageDesc2.rows, 2, CV_32SC1), matchDistance(imageDesc2.rows, 2, CV_32FC1);
        flannIndex.knnSearch(imageDesc2, macthIndex, matchDistance, 2, flann::SearchParams());
    
        // Lowe's algorithm,获取优秀匹配点
        for (int i = 0; i < matchDistance.rows; i++)
        {
            if (matchDistance.at<float>(i, 0) < 0.4 * matchDistance.at<float>(i, 1))
            {
                DMatch dmatches(i, macthIndex.at<int>(i, 0), matchDistance.at<float>(i, 0));
                GoodMatchePoints.push_back(dmatches);
            }
        }
    
        Mat first_match;
        drawMatches(image02, keyPoint2, image01, keyPoint1, GoodMatchePoints, first_match);
        imshow("first_match ", first_match);
    
        vector<Point2f> imagePoints1, imagePoints2;
    
        for (int i = 0; i<GoodMatchePoints.size(); i++)
        {
            imagePoints2.push_back(keyPoint2[GoodMatchePoints[i].queryIdx].pt);
            imagePoints1.push_back(keyPoint1[GoodMatchePoints[i].trainIdx].pt);
        }
    
    
    
        //获取图像1到图像2的投影映射矩阵 尺寸为3*3  
        Mat homo = findHomography(imagePoints1, imagePoints2, CV_RANSAC);
        ////也可以使用getPerspectiveTransform方法获得透视变换矩阵,不过要求只能有4个点,效果稍差  
        //Mat   homo=getPerspectiveTransform(imagePoints1,imagePoints2);  
        cout << "变换矩阵为:\n" << homo << endl << endl; //输出映射矩阵      
    
                                                    //计算配准图的四个顶点坐标
        CalcCorners(homo, image01);
        cout << "left_top:" << corners.left_top << endl;
        cout << "left_bottom:" << corners.left_bottom << endl;
        cout << "right_top:" << corners.right_top << endl;
        cout << "right_bottom:" << corners.right_bottom << endl;
    
        //图像配准  
        Mat imageTransform1, imageTransform2;
        warpPerspective(image01, imageTransform1, homo, Size(MAX(corners.right_top.x, corners.right_bottom.x), image02.rows));
        //warpPerspective(image01, imageTransform2, adjustMat*homo, Size(image02.cols*1.3, image02.rows*1.8));
        imshow("直接经过透视矩阵变换", imageTransform1);
        imwrite("trans1.jpg", imageTransform1);
    
    
        //创建拼接后的图,需提前计算图的大小
        int dst_width = imageTransform1.cols;  //取最右点的长度为拼接图的长度
        int dst_height = image02.rows;
    
        Mat dst(dst_height, dst_width, CV_8UC3);
        dst.setTo(0);
    
        imageTransform1.copyTo(dst(Rect(0, 0, imageTransform1.cols, imageTransform1.rows)));
        image02.copyTo(dst(Rect(0, 0, image02.cols, image02.rows)));
    
        imshow("b_dst", dst);
    
    
        OptimizeSeam(image02, imageTransform1, dst);
    
    
        imshow("dst", dst);
        imwrite("dst.jpg", dst);
    
        waitKey();
    
        return 0;
    }
    
    
    //优化两图的连接处,使得拼接自然
    void OptimizeSeam(Mat& img1, Mat& trans, Mat& dst)
    {
        int start = MIN(corners.left_top.x, corners.left_bottom.x);//开始位置,即重叠区域的左边界  
    
        double processWidth = img1.cols - start;//重叠区域的宽度  
        int rows = dst.rows;
        int cols = img1.cols; //注意,是列数*通道数
        double alpha = 1;//img1中像素的权重  
        for (int i = 0; i < rows; i++)
        {
            uchar* p = img1.ptr<uchar>(i);  //获取第i行的首地址
            uchar* t = trans.ptr<uchar>(i);
            uchar* d = dst.ptr<uchar>(i);
            for (int j = start; j < cols; j++)
            {
                //如果遇到图像trans中无像素的黑点,则完全拷贝img1中的数据
                if (t[j * 3] == 0 && t[j * 3 + 1] == 0 && t[j * 3 + 2] == 0)
                {
                    alpha = 1;
                }
                else
                {
                    //img1中像素的权重,与当前处理点距重叠区域左边界的距离成正比,实验证明,这种方法确实好  
                    alpha = (processWidth - (j - start)) / processWidth;
                }
    
                d[j * 3] = p[j * 3] * alpha + t[j * 3] * (1 - alpha);
                d[j * 3 + 1] = p[j * 3 + 1] * alpha + t[j * 3 + 1] * (1 - alpha);
                d[j * 3 + 2] = p[j * 3 + 2] * alpha + t[j * 3 + 2] * (1 - alpha);
    
            }
        }
    
    }

    看一看拼接效果,我觉得还是不错的。

    看一下这一组图片,这组图片产生了鬼影,为什么?因为两幅图中的人物走动了啊!所以要做图像拼接,尽量保证使用的是静态图片,不要加入一些动态因素干扰拼接。

    opencv自带的拼接算法stitch

    opencv其实自己就有实现图像拼接的算法,当然效果也是相当好的,但是因为其实现很复杂,而且代码量很庞大,其实在一些小应用下的拼接有点杀鸡用牛刀的感觉。最近在阅读sticth源码时,发现其中有几个很有意思的地方。

    1.opencv stitch选择的特征检测方式

    一直很好奇opencv stitch算法到底选用了哪个算法作为其特征检测方式,是ORB,SIFT还是SURF?读源码终于看到答案。

    #ifdef HAVE_OPENCV_NONFREE
            stitcher.setFeaturesFinder(new detail::SurfFeaturesFinder());
    #else
            stitcher.setFeaturesFinder(new detail::OrbFeaturesFinder());
    #endif

    在源码createDefault函数中(默认设置),第一选择是SURF,第二选择才是ORB(没有NONFREE模块才选),所以既然大牛们这么选择,必然是经过综合考虑的,所以应该SURF算法在图像拼接有着更优秀的效果。

    2.opencv stitch获取匹配点的方式

    以下代码是opencv stitch源码中的特征点提取部分,作者使用了两次特征点提取的思路:先对图一进行特征点提取和筛选匹配(1->2),再对图二进行特征点的提取和匹配(2->1),这跟我们平时的一次提取的思路不同,这种二次提取的思路可以保证更多的匹配点被选中,匹配点越多,findHomography求出的变换越准确。这个思路值得借鉴。

    
    matches_info.matches.clear();
    
    Ptr<flann::IndexParams> indexParams = new flann::KDTreeIndexParams();
    Ptr<flann::SearchParams> searchParams = new flann::SearchParams();
    
    if (features2.descriptors.depth() == CV_8U)
    {
        indexParams->setAlgorithm(cvflann::FLANN_INDEX_LSH);
        searchParams->setAlgorithm(cvflann::FLANN_INDEX_LSH);
    }
    
    FlannBasedMatcher matcher(indexParams, searchParams);
    vector< vector<DMatch> > pair_matches;
    MatchesSet matches;
    
    // Find 1->2 matches
    matcher.knnMatch(features1.descriptors, features2.descriptors, pair_matches, 2);
    for (size_t i = 0; i < pair_matches.size(); ++i)
    {
        if (pair_matches[i].size() < 2)
            continue;
        const DMatch& m0 = pair_matches[i][0];
        const DMatch& m1 = pair_matches[i][1];
        if (m0.distance < (1.f - match_conf_) * m1.distance)
        {
            matches_info.matches.push_back(m0);
            matches.insert(make_pair(m0.queryIdx, m0.trainIdx));
        }
    }
    LOG("\n1->2 matches: " << matches_info.matches.size() << endl);
    
    // Find 2->1 matches
    pair_matches.clear();
    matcher.knnMatch(features2.descriptors, features1.descriptors, pair_matches, 2);
    for (size_t i = 0; i < pair_matches.size(); ++i)
    {
        if (pair_matches[i].size() < 2)
            continue;
        const DMatch& m0 = pair_matches[i][0];
        const DMatch& m1 = pair_matches[i][1];
        if (m0.distance < (1.f - match_conf_) * m1.distance)
            if (matches.find(make_pair(m0.trainIdx, m0.queryIdx)) == matches.end())
                matches_info.matches.push_back(DMatch(m0.trainIdx, m0.queryIdx, m0.distance));
    }
    LOG("1->2 & 2->1 matches: " << matches_info.matches.size() << endl);

    这里我仿照opencv源码二次提取特征点的思路对我原有拼接代码进行改写,实验证明获取的匹配点确实较一次提取要多。

    //提取特征点    
    SiftFeatureDetector Detector(1000);  // 海塞矩阵阈值,在这里调整精度,值越大点越少,越精准 
    vector<KeyPoint> keyPoint1, keyPoint2;
    Detector.detect(image1, keyPoint1);
    Detector.detect(image2, keyPoint2);
    
    //特征点描述,为下边的特征点匹配做准备    
    SiftDescriptorExtractor Descriptor;
    Mat imageDesc1, imageDesc2;
    Descriptor.compute(image1, keyPoint1, imageDesc1);
    Descriptor.compute(image2, keyPoint2, imageDesc2);
    
    FlannBasedMatcher matcher;
    vector<vector<DMatch> > matchePoints;
    vector<DMatch> GoodMatchePoints;
    
    MatchesSet matches;
    
    vector<Mat> train_desc(1, imageDesc1);
    matcher.add(train_desc);
    matcher.train();
    
    matcher.knnMatch(imageDesc2, matchePoints, 2);
    
    // Lowe's algorithm,获取优秀匹配点
    for (int i = 0; i < matchePoints.size(); i++)
    {
        if (matchePoints[i][0].distance < 0.4 * matchePoints[i][1].distance)
        {
            GoodMatchePoints.push_back(matchePoints[i][0]);
            matches.insert(make_pair(matchePoints[i][0].queryIdx, matchePoints[i][0].trainIdx));
        }
    }
    cout<<"\n1->2 matches: " << GoodMatchePoints.size() << endl;
    
    #if 1
    
    FlannBasedMatcher matcher2;
    matchePoints.clear();
    vector<Mat> train_desc2(1, imageDesc2);
    matcher2.add(train_desc2);
    matcher2.train();
    
    matcher2.knnMatch(imageDesc1, matchePoints, 2);
    // Lowe's algorithm,获取优秀匹配点
    for (int i = 0; i < matchePoints.size(); i++)
    {
        if (matchePoints[i][0].distance < 0.4 * matchePoints[i][1].distance)
        {
            if (matches.find(make_pair(matchePoints[i][0].trainIdx, matchePoints[i][0].queryIdx)) == matches.end())
            {
                GoodMatchePoints.push_back(DMatch(matchePoints[i][0].trainIdx, matchePoints[i][0].queryIdx, matchePoints[i][0].distance));
            }
            
        }
    }
    cout<<"1->2 & 2->1 matches: " << GoodMatchePoints.size() << endl;
    #endif

    最后再看一下opencv stitch的拼接效果吧~速度虽然比较慢,但是效果还是很好的。

    #include <iostream>
    #include <opencv2/core/core.hpp>
    #include <opencv2/highgui/highgui.hpp>
    #include <opencv2/imgproc/imgproc.hpp>
    #include <opencv2/stitching/stitcher.hpp>
    using namespace std;
    using namespace cv;
    bool try_use_gpu = false;
    vector<Mat> imgs;
    string result_name = "dst1.jpg";
    int main(int argc, char * argv[])
    {
        Mat img1 = imread("34.jpg");
        Mat img2 = imread("35.jpg");
    
        imshow("p1", img1);
        imshow("p2", img2);
    
        if (img1.empty() || img2.empty())
        {
            cout << "Can't read image" << endl;
            return -1;
        }
        imgs.push_back(img1);
        imgs.push_back(img2);
    
    
        Stitcher stitcher = Stitcher::createDefault(try_use_gpu);
        // 使用stitch函数进行拼接
        Mat pano;
        Stitcher::Status status = stitcher.stitch(imgs, pano);
        if (status != Stitcher::OK)
        {
            cout << "Can't stitch images, error code = " << int(status) << endl;
            return -1;
        }
        imwrite(result_name, pano);
        Mat pano2 = pano.clone();
        // 显示源图像,和结果图像
        imshow("全景图像", pano);
        if (waitKey() == 27)
            return 0;
    }

     

    展开全文
  • 翻译:图像拼接

    2019-10-24 21:24:42
    图像拼接1. 简介1.1 特征点提取1.2 特征点匹配1.3 单应估计1.4 图像合成1.5 结果2.源码2.1 特征点提取2.2 特征点匹配2.3 单应估计2.4 小结3. 例子4. 备注5. 结论6. 引用7.参考 1. 简介 图像拼接在无人机航拍、遥感...

    0 简介

    图像拼接是将同一场景的多个重叠图像拼接成较大的图像的一种方法,在医学成像、计算机视觉、卫星数据、军事目标自动识别等领域具有重要意义。图像拼接的输出是两个输入图像的并集。通常用到五个步骤:

    Created with Raphaël 2.2.0输入图像特征提取Feature Extraction图像配准 Image Registration用RANSAC计算单应矩阵H变形和融合 ImageWraping and Blending输出图像

    特征提取 Feature Extraction:在所有输入图像中检测特征点
    图像配准 Image Registration:建立了图像之间的几何对应关系,使它们可以在一个共同的参照系中进行变换、比较和分析。
    大致可以分为以下几个类

    1. 直接使用图像的像素值的算法,例如,correlation methods
    2. 在频域处理的算法,例如,基于快速傅里叶变换(FFT-based)方法;
    3. 低水平特征的算法low level features,通常用到边缘和角点,例如,基于特征的方法,
    4. 高水平特征的算法high-level features,通常用到图像物体重叠部分,特征关系,例如,图论方法(Graph-theoretic methods)

    图像变形 Warping:
    图像变形是指将其中一幅图像的图像重投影,并将图像放置在更大的画布上。
    图像融合 Blending
    图像融合是通过改变边界附近的图像灰度级,去除这些缝隙,创建混合图像,从而在图像之间实现平滑过渡。混合模式(Blend modes)用于将两层融合到一起。

    特征点提取

    特征是要匹配的两个输入图像中的元素,它们是在图像块的内部。这些图像块是图像中的像素组。对输入图像进行Patch匹配。具体解释如下: 如下图所示,fig1和fig2给出了一个很好的patch匹配,因为fig2中有一个patch看起来和fig1中的patch非常相似。当我们考虑到fig3和fig4时,这里的patch并不匹配,因为fig4中有很多类似的patch,它们看起来与fig3中的patch很相似。由于像素强度很相近,所以无法进行精确的特征匹配,

    为了给图像对提供更好的特征匹配,采用角点匹配,进行定量测量。角点是很好的匹配特性。在视点变化时,角点特征是稳定的。此外,角点的邻域具有强度突变。利用角点检测算法对图像进行角点检测。角点检测算法有Harris角点检测算法、SIFT特征点检测算法((Scale Invariant Feature Transform),FAST算法角点检测算法,SURF特征点检测算法(Speeded-up robust feature)

    Harris角点检测算法

    Harris算法是一种基于Moravec算法的点特征提取算法。1988年C. Harris 和 M.J Stephens设计了一种图像局部检测窗口。通过在不同的方向上移动少量窗口,可以确定强度的平均变化。我们可以通过观察小窗口内的强度值很容易地识别角点。在移动窗口时,平坦区域在所有方向上均不会显示强度的变化。边缘区域在沿边缘方向强度不会发生变化。对于角点,则在各个方向上产生显著强度变化。Harris角点探测器给出了一种检测平坦区域、边缘和角点的数学方法。Harris检测的特征较多,具有旋转不变性和尺度变异性。位移[u,v][u, v]下的强度变化:E(u,v)=x,yw(x,y)[I(x+u,y+v)I(x,y)]2E(u,v)=∑_{x,y}w(x,y)[I(x+u,y+v)−I(x,y)]^2其中,w(x,y)w(x,y)是窗口函数,I(x+u,y+v)I(x+u,y+v)是移动后的强度,I(x,y)I(x,y)是单个像素位置的强度。

    Harris角点检测算法如下:

    1. 对图像中的每个像素点(x,y)(x,y)计算自相关矩阵MM(autocorrelation matrix M):
      M=x,y[Ix2IxIyIxIyIy2]M=\sum_{x,y} \begin{bmatrix}I_{x}^{2} & I_{x}I_{y}\\ I_{x}I_{y} & I_{y}^{2}\end{bmatrix}其中Ix,IyI_{x},I_{y}I(x,y)I(x,y)的偏导数
    2. 对图像中的每个像素点做高斯滤波,获得新的矩阵MM,离散二维零均值高斯函数为Gauss=exp(u2+v2)/2δ2Gauss = exp(-u^2+v^2)/2\delta^2
    3. 计算每个像素点(x,y)的角点度量,得到R=Det(M)ktrace(M)R=Det(M)-k*trace(M)kk 的范围是0.04k0.060.04≤k≤0.06
    4. 选择局部最大值点。Harris方法认为特征点与局部最大兴趣点的像素值对应。
    5. 设置阈值T,检测角点。如果 RR 的局部最大值高于阈值TT,那么此点为角点。

    SIFT角点检测算法

    SIFT算法是尺度不变的特征点检测算法,可用于识别其他图像中的相似目标。SIFT的图像特征表示为关键点描述符(key-point-descriptors)。在检查图像匹配时,将两组关键点描述符作为输入提供给最近邻搜索(Nearest Neighbor Search,NNS),并生成一个紧密匹配的关键点描述符(matching key-point-descriptors)。

    Blog_1-s2.0-S1047320315002059-gr5.jpg
    SIFT的计算分为四个阶段:

    1. 尺度空间构造(Scale-space construction)
    2. 尺度空间极值检测(Scale-space extrema detection)
    3. 关键点定位(key-point localization)
    4. 方向分配(orientation assignment)和关键点描述符定义(defining key-point descriptors)

    第一阶段确定潜在的兴趣点。它利用高斯函数的差分(difference of Gaussian function,DOG)搜索所有尺度和图像位置。第一阶段中发现的所有兴趣点的location和scale是确定的。根据关键点的稳定性来选择关键点。一个稳定的关键点能够抵抗图像失真。在方向分配环节,SIFT算法计算稳定关键点周围梯度的方向。根据局部图像梯度方向,为每个关键点分配一个或多个方向。对于一组输入帧,SIFT提取特征。图像匹配使用Best Bin First(BBF)算法来估计输入帧之间的初始匹配点。为了去除不属于重叠区域的不需要的角,使用RANSAC算法。它删除图像对中的错误匹配。通过定义帧的大小、长度和宽度来实现帧的重投影。最后进行拼接,得到最终的输出拼接图像。在拼接时,检查场景每帧中的每个像素是否属于扭曲的第二帧。如果是,则为该像素分配来自第一帧的对应像素的值。SIFT算法既具有旋转不变性,又具有尺度不变性。SIFT非常适合于高分辨率图像中的目标检测。它是一种鲁棒的图像比较算法,虽然速度较慢。SIFT算法的运行时间很大,因为比较两幅图像需要更多的时间。

    FAST 算法

    FAST是Trajkovic和Hedley在1998年创建的角点检测算法。对于FAST,角点的检测优于边缘检测,因为角点有二维强度变化,容易从邻近点中区分出来。适用于实时图像处理应用程序。

    FAST角点探测器应该满足以下要求

    1. 检测到的位置要一致,对噪声变化不敏感,对同一场景的多幅图像不能移动。
    2. 准确;检测到的角点应该尽可能接近正确的位置。
    3. 速度;角落探测器应该足够快。

    原理:首先围绕一个候选角点选择16个像素点。如果其中有n(n一般为12)个连续的像素都比候选角点加上一个阈值要高,或者比候选角点减去一个阈值要低,那么此点即为角点(如图4所示)

    为了加快FAST算法的速度,通常会使用角点响应函数( corner response function, CRF)。该函数根据局部邻域的图像强度给出角点强度的数值。
    对图像进行CRF计算,并将CRF的局部最大值作为角点,采用多网格(multi-grid)技术提高了算法的计算速度,并对检测到的假角点进行了抑制。FAST是一种精确、快速的算法,具有良好的定位(位置精度)和较高的点可靠性。FAST的角点检测的算法难点在于最佳阈值的选择。

    SURF算法

    Speed-up Robust Feature(SURF)角点探测器采用三个特征检测步骤;检测(Detection)、描述(Description)、匹配(Matching),SURF通过考虑被检测点的质量,加快了位移的检测过程。它更注重加快匹配步骤。使用Hessian矩阵和低维描述符来显著提高匹配速度。SURF在计算机视觉社区中得到了广泛的应用。SURF在不变特征定位上十分有效和鲁棒

    图像配准

    在特征点被检测出来之后,我们需要以某种方式将它们关联起来,可以通过NCC或者SDD(Sum of Squared Difference)方法来确定其对应关系。

    归一化互相关(normalized cross correlation,NCC)

    互相关的工作原理是分析第一幅图像中每个点周围的像素窗口,并将它们与第二幅图像中每个点周围的像素窗口关联起来。将双向相关性最大的点作为对应的对。

    基于图像强度值计算在两个图像中的每个位移(shifts)的“窗口”之间的相似性

    NCC(u)=i[I1(xi)I1ˉ][I2(xi+u)I2ˉ]i[I1(xi)I1ˉ]2[I2(xi+u)I2ˉ]2NCC(u)=\frac{\sum_i[I_1(x_i)-\bar{I_1}][I_2(x_i+u)-\bar{I_2}] }{\sqrt{\sum_i[I_1(x_i)-\bar{I_1}]^2[I_2(x_i+u)-\bar{I_2}]^2} }
    其中,I1ˉ,I2ˉ\bar{I_1},\bar{I_2}是窗口的平均值图像
    I1ˉ=1NiI1(xi)\bar{I_1}=\frac{1}{N}\sum _i I_1(x_i)
    I2ˉ=1NiI2(xi+u)\bar{I_2}=\frac{1}{N}\sum _i I_2(x_i+u)
    I1(x,y)I_1(x,y)I2(x,y)I_2(x,y)分别是两张图片。xi=(xi,yi)x_i=(x_i,y_i) 是窗口的像素坐标,u=(u,v)u=(u,v) 是通过NCC系数计算出的位移或偏移。NCC系数的范围为[1,1][-1,1]。 NCC峰值相对应的位移参数表示两个图像之间的几何变换。此方法的优点是计算简单,但是速度特别慢。此外,此类算法要求源图像之间必须有显著的重叠。

    互信息(Mutual Information, MI)

    互信息测量基于两个图像之间共享信息数量的相似性。

    两个图像I1(X,Y)I_1(X,Y)I2(X,Y)I_2(X,Y)之间的MI以熵表示:

    MI(I1,I2)=E(I1)+E(I2)E(I1,I2)MI(I_1,I_2)=E(I_1)+E(I_2)−E(I_1,I_2)
    其中,E(I1)E(I_1)E(I2)E(I_2)分别是I1(x,y)I_1(x,y)I2(x,y)I_2(x,y)的熵。E(I1,I2)E(I_1,I_2)表示两个图像之间的联合熵。
    E(I1)=gpI1(g)log(pI1(g))E(I_1)=−∑_gp_{I1}(g)log(p_{I1}(g))
    ggI1(x,y)I_1(x,y)可能的灰度值,pI1(g)p_{I1}(g)gg的概率分布函数
    E(I1,I2)=g,hpI1,I2(g,h)log(pI1,I2(g,h))E(I1,I2)=−∑_{g,h}p_{I_1,I_2}(g,h)log(p_{I_1,I_2}(g,h))

    然而,从图中我们可以看到,许多点被错误地关联在一起。
    1093303-20170822154344449-509144423.png

    计算单应矩阵

    单应矩阵估计是图像拼接的第三步。在单应矩阵估计中,不属于重叠区域的不需要的角被删除。采用RANSAC算法进行单应。

    随机样本一致算法RANSAC(random sample consensus)

    RANSAC算法从可能含有异常值的观测数据集中拟合数学模型,是一种鲁棒参数估计的迭代方法。该算法是不确定性的,因为它只在一定的概率下产生一个合理的结果,随着执行更多的迭代,这个概率会增加。RANSAC算法用于在存在大量可用数据外行的情况下以鲁棒的方式拟合模型。RANSAC算法在计算机视觉中有许多应用。
    20190610124718901.png

    RANSAC原理

    从数据集中随机选取一组数据并认为是有效数据(内点)来确定待定参数模型,以此模型测试数据集中的所有数据,满足该模型的数据成为内点,反之为外点(通常为噪声、错误测量或不正确数据的点),迭代执行,直到某一个参数模型得到的内点数最大,则该模型为最优模型。
    考虑如下假设:

    1. 参数可以从N个数据项中估计。
    2. 可用的数据项总共是M。
    3. 随机选择的数据项成为好模型的一部分的概率为PgP_g
    4. 如果存在一个很好的拟合,那么算法在没有找到一个很好的拟合的情况下退出的概率是PfailP_{fail}

    RANSAC步骤

    1. 随机选取N个数据(3个点对)
    2. 估计参数x(计算变换矩阵H)
    3. 根于使用者设定的阈值,找到M中合适该模型向量x的的数据对总数量K( 计算每个匹配点经过变换矩阵后到对应匹配点的距离,根据预先设定的阈值将匹配点集合分为内点和外点,如果内点足够多,则H足够合理,用所有内点重新估计H)。
    4. 如果符合的数量K足够大,则接受该模型并退出
    5. 重复1-4步骤 L次
    6. 到这一步退出

    K有多大取决于我们认为属于合适结构的数据的百分比以及图像中有多少结构。如果存在多个结构,则在成功拟合后,删除拟合数据并重做RANSAC。

    迭代次数L可以用如下公式计算:
    Pfail=LP_{fail} = L连续失败的概率
    Pfail=()LP_{fail}=(给定试验失败的概率)L
    Pfail=(1)LP_{fail}=(1 - 给定试验成功的概率)L
    Pfail=(1()N)LP_{fail}=(1-(随机数据项符合模型的概率)N)L
    Pfail=(1(Pg)N)LP_{fail}=(1-(Pg)^N)^L
    L=log(Pfail)/log(1(Pg)N)L = log(P_{fail})/log(1-(Pg)N)

    优点:可以robust地估计模型参数
    缺点:迭代次数无上限,设置的迭代次数会影响算法时间复杂度和精确程度,并且需要预设阈值

    在执行RANSAC之后,我们只能在图像中看到正确的匹配,因为RANSAC找到了一个与大多数点相关的单应矩阵,并将不正确的匹配作为异常值丢弃

    单应矩阵(Homography)

    有了两组相关点,接下来就需要建立两组点的转换关系,也就是图像变换关系。单应性是两个空间之间的映射,常用于表示同一场景的两个图像之间的对应关系,可以匹配大部分相关的特征点,并且能实现图像投影,使一张图通过投影和另一张图实现大面积的重合。

    设2个图像的匹配点分别是X=[x,y]TX=[x,y]^T,X=[x,y]TX'=[x',y']^T,则必须满足公式:
    X=HXX'=HX且由于两向量共线,所以X×HX=0X'\times HX = 0其中,HH 为8参数的变换矩阵,可知四点确定一个H
    (xy1)=(h11h12h13h21h22h23h31h321)(xy1)\begin{pmatrix}x' \\y'\\1 \end{pmatrix} =\begin{pmatrix} h_{11} & h_{12} & h_{13}\\ h_{21} & h_{22} & h_{23}\\ h_{31} & h_{32} & 1 \end{pmatrix}\begin{pmatrix}x\\y\\1\\\end{pmatrix}

    h=(h11:h12:h13:h21:h22:h23:h31:h32:h33)Th =(h11:h12:h13:h21:h22:h23:h31:h32:h33)T则有
    Bh=0Bh=0N个点对给出2N个线性约束。
    minhBh2h=1\underset{h}{min}║Bh║^2 ,║h║ = 1
    用RANSAC方法估算H:

    1. 首先检测两边图像的角点
    2. 在角点之间应用方差归一化相关,收集相关性足够高的对,形成一组候选匹配。
    3. 选择四个点,计算H
    4. 选择与单应性一致的配对。如果对于某些阈值:Dist(Hp、q) <ε,则点对(p, q)被认为与单应性H一致
    5. 重复34步,直到足够多的点对满足H
    6. 使用所有满足条件的点对,通过公式重新计算H

    图像变形和融合

    最后一步是将所有输入图像变形并融合到一个符合的输出图像中。基本上,我们可以简单地将所有输入的图像变形到一个平面上,这个平面名为复合全景平面。

    图像变形步骤

    1. 首先计算每个输入图像的变形图像坐标范围,得到输出图像大小,可以很容易地通过映射每个源图像的四个角并且计算坐标(x,y)的最小值和最大值确定输出图像的大小。最后,需要计算指定参考图像原点相对于输出全景图的偏移量的偏移量x_offset和偏移量y_offset。
    2. 下一步是使用上面所述的反向变形,将每个输入图像的像素映射到参考图像定义的平面上,分别执行点的正向变形和反向变形。

    1093303-20170822154403277-379911836.png
    平滑过渡(transition smoothing)图像融合方法包括 羽化(feathering), 金字塔(pyramid), 梯度(gradient)

    图形融合

    最后一步是在重叠区域融合像素颜色,以避免接缝。最简单的可用形式是使用羽化(feathering),它使用加权平均颜色值融合重叠的像素。我们通常使用alpha因子,通常称为alpha通道,它在中心像素处的值为1,在与边界像素线性递减后变为0。当输出拼接图像中至少有两幅重叠图像时,我们将使用如下的alpha值来计算其中一个像素处的颜色:
    假设两个图像 I1,I2I_1,I_2,在输出图像中重叠;每个像素点(x,y)(x,y)在图像Ii(x,y)=(αiR,αiG,αiB,αj,)I_i(x,y)=(\alpha iR,\alpha iG,\alpha iB,\alpha j,),其中(R,G,B)是像素的颜色值,我们将在缝合后的输出图像中计算(x, y)的像素值:
    [(α1R,α1G,α1B,α1)+(α2R,α2G,α2B,α2)]/(α1+α2) [(α1R, α1G, α1B, α1) + (α2R, α2G, α2B, α2)]/(α1+α2).

    1093303-20170822154425449-1866547572.png
    7.png

    参考

    1. OpenCV探索之路(二十四)图像拼接和图像融合技术

    2. Debabrata Ghosh,Naima Kaabouch. A survey on image mosaicing techniques[J]. Journal of Visual Communication and Image Representation,2016,34.地址

    3. 图像拼接综述

    展开全文
  • 图像拼接算法(zz)

    2018-07-25 20:32:54
    图像拼接算法原理 1: http://planckscale.info/?p=7 http://planckscale.info/?p=84 http://blog.csdn.net/xiaolizi399/article/details/41931809 文章来源(planckscale.info) 360°全景拼接技术简介   核心...
  • 转载自:OpenCV探索之路(二十四)图像拼接图像融合技术 测试环境为Vs2015+OpenCV2.4.9 一、OpenCV自带的拼接算法stitch 这部分在opencv自带的例程中也有相应的sample,可以参照。 #include &lt;...
  • from: ... For the latest version of the code, which may contain the latest enhancements and corrections, please download the latest ver...
  • 图像拼接算法及实现

    2014-10-30 18:14:55
    图像拼接解决的问题一般式,通过对齐一系列空间重叠的图像,构成一个无缝的、高清晰的图像,它具有比单个图像更高的分辨率和更大的视野。  早期的图像拼接研究一直用于照相绘图学,主要是对大量
  • 答案是当然有,图像拼接的第二春13年开始的,这里再简单总结和展望一下。一个(全局)单应性对齐+柱面或球面投影+光束平差+多频带融合为核心的老一代拼接算法以BROWN大神03’ICCV和07’IJCV的AutoStitch AutoStitch为...
  • CV——图像拼接

    2019-03-24 15:31:00
    何谓图像拼接,顾名思义,就是把两张具有重叠区域的图片拼接在一起,形成一个更完整的图片,如下图: 我们把第一个图片和第二个图片拼接起来,形成了第三个图片。 图像拼接步骤 检测特征点 特征点匹配 确定对应...
  • python opencv 图像拼接

    2018-08-21 10:28:50
    高级图像拼接也叫作基于特征匹配的图像拼接拼接时消去两幅图像相同的部分,实现拼接合成全景图。 具有相同尺寸的图A和图B含有相同的部分与不同的部分,如图所示:   用基于特征的图像拼接实现后: 设图像...
  • python最全的图像拼接

    2019-10-30 15:02:06
    搞这个图片拼接真的是搞了很久,尝试了很多种方法,现在都在这里列举一下,与大家分享一下,相互激励一下彼此吧!!! 首先说明一下: 我的测试图片文件夹是: test_img,一共六张图片 与py文件同级 (小姐姐是...
  • MATLAB 程序 实现图像拼接,自动匹配特征点 将图片拼接成长图,MATLAB 程序 实现图像拼接,自动匹配特征点 将图片拼接成长图MATLAB 程序 实现图像拼接,自动匹配特征点 将图片拼接成长图
  • 图像拼接(十三):OpenCV单应变换模型拼接多幅图像 这篇博客中实现了用单应变换模型拼接多幅图像图像拼接使用的OpenCV库函数warpPerspective()。因为这个函数只有在右侧图像变换到左侧图像时才能完整显示,所以...
  • matlab图像拼接的四种方法 1、直接拼接, 2、亮度调整后拼接, 3、按距离比例融合, 4、亮度调整后按距离比例融合 流程: 1。读入左,右图,并取出重合部分,并转化为亮度图 2。分别把每点的亮度值相加,得到一个...
  • 图像拼接,就是将多张图片拼接成一个大图,类似手机里的全景照片。按照某大神的思路,可以自己按照如下步骤实现: 对每幅图进行特征点提取 对对特征点进行匹配 进行图像配准 把图像拷贝到另一幅图像的特定位置 ...
  • 图像拼接是利用连续帧图像生成全景图或更高分辨率的图像,通常图像拼接技术需要消除图像拼接部分的缝隙间隔,因些需要进行重叠区域匹配修复。通常图像拼接分为以下三个步骤: ⑴特征点检测。图像拼接操作对序列图像...
  • 1、单张图片拼接 # 图片拼接 from PIL import Image # pil paste可以进行图片拼接 import cv2 import numpy as np path="F:/out/"+str(0)+".jpg" img_out=cv2.imread(path) num=5 for i in ...
1 2 3 4 5 ... 20
收藏数 83,206
精华内容 33,282
关键字:

图像拼接