精华内容
下载资源
问答
  • 立体匹配NCC算法
    千次阅读
    2015-04-02 09:49:54

    FROM:http://blog.csdn.net/tulun/article/details/6388759

    NCC算法(Normal Cross Correlation),具体原理见相关图像处理书籍。

    该程序是opencv中文论坛的牛人贡献的,感谢他的工作。

     

    (程序所需图片可以在网上找如http://vision.middlebury.edu/stereo/data/scenes2003/)

     

    #include <iostream>
    #include <stdio.h>
    #include <stdlib.h>
    #include <cv.h>
    #include <cxcore.h>
    #include <highgui.h>
    #include <math.h>
    #include <ctime>

    using namespace std;
    template<class T> class Image
    {
    private:
        IplImage* imgp;

    public:
        Image(IplImage* img=0){imgp=img;}
        ~Image(){imgp=0;}
        void operator=(IplImage* img){imgp=img;}
        inline T* operator[](const int rowIndx)
        {
            return((T*)(imgp->imageData+rowIndx*imgp->widthStep));
        }

    };

    typedef struct
    {
        unsigned char b,g,r;
    }RgbPixel;

    typedef struct
    {
        float b,g,r;
    }RgbPixelFloat;

    typedef Image<RgbPixel> RgbImage;
    typedef Image<RgbPixelFloat> RgbImageFloat;
    typedef Image<unsigned char> BwImage;
    typedef Image<float> BwImageFloat;


    void displayImageProperty(IplImage* image)
    {
        cout<<"-------Image Properties--------"<<endl;
        cout<<"Image width="<<image->width<<endl;
        cout<<"Image height="<<image->height<<endl;
        cout<<"Image depth="<<image->depth<<endl;
        cout<<"Image nSize="<<image->nSize<<endl;
        cout<<"Image nChannels="<<image->nChannels<<endl;

        char* origin;
        char* dataOrder;
        if (image->origin==0)
        {
            origin="Top-left";
        }
        else
        {
            origin="Below-left";//image->origin=1
        }

        cout<<"Image origin="<<origin<<endl;

        if (image->dataOrder==0)
        {
            dataOrder="Order_Pixel(Interleaved)";
        }
        else
        {
            dataOrder="Order_Plane";//image->dataOrder=1
        }
        cout<<"Image dataOrder="<<dataOrder<<endl;
        cout<<"Image widthStep="<<image->widthStep<<" Bytes"<<endl;
    }

    // display an image in a new window with title to be given.
    void displayImageNewWindow(char* title,CvArr* img)
    {
        cvNamedWindow(title, CV_WINDOW_AUTOSIZE );
        cvShowImage(title,img);
    }


    int getMaxMin(double value[],int valueSize, int maxmin)
    {
        int pos=0;
        int i=0;
        double max1=-1;//?-999999;
        double min1=999999;
          
        if (maxmin==1)
        {
            //find max
            for (i=0;i<valueSize;i++)
            {
                //find the index with the max ncc;
                if (value[i]>max1)
                {
                    pos=i;
                    max1=value[i];
                }
             }
        }
          
        if (maxmin==0)
        {
            //find min
            for (i=0;i<valueSize;i++)
            {
                //find the index with the max ncc;
                if (value[i]<min1)
                {
                    pos=i;
                    min1=value[i];
                }
             }
        }
        return pos;
    }

    IplImage* generateDisparityImage(IplImage* greyLeftImg32,IplImage* greyRightImg32,int windowSize,int DSR)
    {
        int offset=floor((double)windowSize/2);
        int height=greyLeftImg32->height;
        int width=greyLeftImg32->width;
        double* localNCC=new double[DSR];

        int x=0, y=0,d=0,m=0;
        int N=windowSize;           

        IplImage* leftWinImg=cvCreateImage(cvSize(N,N),32,1);//mySubImage(greyLeftImg32,cvRect(0,0,N,N));
        IplImage* rightWinImg=cvCreateImage(cvSize(N,N),32,1);;//mySubImage(greyRightImg32,cvRect(0,0,N,N));
        IplImage* disparity=cvCreateImage(cvSize(width,height),8,1);//or IPL_DEPTH_8U
        BwImage imgA(disparity);
          
        for (y=0;y<height;y++)
        {
            for (x=0;x<width;x++)
            {
                imgA[y][x]=0;
            }
        }
          
        CvScalar s1;
        CvScalar s2;
        for (y=0;y<height-N;y++)
        {
            //height-N
            for (x=0;x<width-N;x++)
            {
                //width-N
                //getWindow(i,j,leftim,wl,N);
                cvSetImageROI(greyLeftImg32, cvRect(x,y,N,N));
                s1=cvAvg(greyLeftImg32,NULL);
                cvSubS(greyLeftImg32,s1,leftWinImg,NULL);//zero-means
                cvNormalize(leftWinImg,leftWinImg,1,0,CV_L2,NULL);//0变成1
                d=0;
                
                //initialise localNCC
                for (m=0;m<DSR;m++)
                {
                    localNCC[m]=0;
                }
                
                do{
                    if (x-d>=0)
                    {
                        cvSetImageROI(greyRightImg32, cvRect(x-d,y,N,N));
                        s2=cvAvg(greyRightImg32,NULL);
                        cvSubS(greyRightImg32,s2,rightWinImg,NULL);//zero-means
                        cvNormalize(rightWinImg,rightWinImg,1,0,CV_L2,NULL);//0变成1
                    }
                    else
                    {
                        break;
                    }
                    localNCC[d]=cvDotProduct(leftWinImg,rightWinImg);
                    cvResetImageROI(greyRightImg32);
                    d++;
                 }while(d<=DSR);
                
                 //to find the best d and store
                 imgA[y+offset][x+offset]=getMaxMin(localNCC,DSR,1)*16;
                 cvResetImageROI(greyLeftImg32);
              }//x
              if (y%10==0)
                  cout<<"row="<<y<<" of "<<height<<endl;
           }//y
          
           cvReleaseImage(&leftWinImg);
           cvReleaseImage(&rightWinImg);
              
           return disparity;
    }

    int main (int argc, char * const argv[])
    {
        // insert code here...
        cout << "Stereo Normalized Cross Correlation"<<endl;
          
        //**********image input*********************//
          

        char* filename1="im0.ppm";//im2_cone.png
        IplImage* greyLeftImg= cvLoadImage(filename1,0);
        char* filename2="im1.ppm";
        IplImage* greyRightImg= cvLoadImage(filename2,0);
          
        if (greyLeftImg==NULL)
        {
            cout << "No valid image input."<<endl;
            //char c=getchar();
            return 1;
        }
        else
        {
            displayImageProperty(greyLeftImg);
        }
          
        if (greyRightImg==NULL)
        {
            cout << "No valid image input."<<endl;
            //char c=getchar();
            return 1;
        }

        int width=greyLeftImg->width;
        int height=greyLeftImg->height;

        /****************8U to 32F**********************/
        IplImage* greyLeftImg32=cvCreateImage(cvSize(width,height),32,1);//IPL_DEPTH_32F
        IplImage* greyRightImg32=cvCreateImage(cvSize(width,height),32,1);
        cvConvertScale(greyLeftImg, greyLeftImg32, 1/255.);
        cvConvertScale(greyRightImg, greyRightImg32, 1/255.);//1/255. equals to 1/255.0
          
        //-------------Computing stereo matching----------------
        time_t tstart, tend;
        tstart = time(0);
        int windowSize=11,DSR=20;//Disparity Search Range
        IplImage* disparity32=generateDisparityImage(greyLeftImg32,greyRightImg32,windowSize,DSR);
        tend = time(0);
        cout << "It took "<< difftime(tend, tstart) <<" second(s)."<< endl;

        displayImageNewWindow("Dispairty Image",disparity32);
        displayImageNewWindow("Left Image",greyLeftImg32);
        displayImageNewWindow("Right Image",greyRightImg32);
        //cvSaveImage("D:/OpenCV_stuff/SampleImages/disparity.jpg",disparity32);

        //********destroy window************/
        cvWaitKey(0);
        cvReleaseImage(&greyLeftImg32);
        cvReleaseImage(&greyRightImg32);
        cvReleaseImage(&greyLeftImg);
        cvReleaseImage(&greyRightImg);
        cvReleaseImage(&disparity32);
        cvDestroyWindow("Left Image");
        cvDestroyWindow("Right Image");
        cvDestroyWindow("Dispairty Image");
        return 0;
    }
    更多相关内容
  • 为此,本文提出一种改进的NCC立体匹配算法,通过引入积分图像和平方积分图像,将矩形窗口区域像素求和运算转化为四个像素点值的简单相加减,同时剔除基准图像中无法匹配区域以减小搜索范围,使计算复杂度得到简化,...
  • 三个基础性的立体匹配算法,直接下载下来就可以运行,学立体匹配必须掌握的算法,只要5个币,真的是白菜价。 希望能够帮助到大家,下载后绝对不会后悔的。
  • 立体匹配NCC算法代码

    2015-09-21 22:24:14
    在VS上实现的立体匹配NCC算法,用C++编写。很适合初学者学习,注意配置opencv
  • 前言一、立体视觉原理1.1 几何原理1.2 立体匹配的步骤1.3 窗口匹配方法匹配代价1.4 归一化互相关(NCC)二、实验内容2.1 实验目的与要求2.2 代码实现2.3 实验结果与分析三、总结 前言 立体匹配主要是通过找出每对...

    前言

    立体匹配主要是通过找出每对图像间的对应关系,根据三角测量原理,得到视差图;在获得了视差信息后,根据投影模型很容易地可以得到原始图像的深度信息和三维信息。立体视觉匹配在立体视觉研究中是比较核心的问题。立体视觉应用:车导航,3D场景重建等。本文主要介绍窗口匹配方法NCC视差匹配

    一、立体视觉原理

    1.1 几何原理

    一个多视图成像的特殊例子是主体视觉(或者立体成像),即使用两台只有水平(向一侧) 偏移的照相机观测同一场景。当照相机的位置如上设置,两幅图像具有相同的图像平面,图像的行是垂直对齐的,那么称图像对是经过矫正的。
    假设两幅图像经过了矫正,那么对应点的寻找限制在图像的同一行上。一旦找到对应点,由于深度是和偏移成正比的,那么深度(Z坐标)可以直接由水平偏移来计算。那么Z该如何计算呢?先看下图:
    在这里插入图片描述
    通过几何知识运算,得到 T + x l − x r T = Z − f Z \frac{T+x_{l}-x_{r}}{T}=\frac{Z-f}{Z} TT+xlxr=ZZf进一步化简,可得 Z = f T x r − x l Z=f\frac{T}{x_{r}-x_{l}} Z=fxrxlT其中,f是经过矫正图像的焦距,T是两个照相机中心之间的距离, x l 和 x r x_{l}和x_{r} xlxr是左右两幅图像中对应点的x坐标。分开照相机中心的距离称为基线。

    1.2 立体匹配的步骤

    1)匹配代价计算

    计算匹配代价,即计算参考图像上每个像素点 I R ( P ) IR(P) IR(P) ,以所有视差可能性去匹配目标图像上对应点 I T ( p d ) IT(pd) IT(pd)的代价值,因此计算得到的代价值可以存储在一个 h ∗ w ∗ d ( M A X ) h*w*d(MAX) hwd(MAX)的三维数组中,通常称这个三维数组为视差空间图。匹配代价时立体匹配的基础,设计抗噪声干扰、对光照变化不敏感的匹配代价,能提高立体匹配的精度。因此,匹配代价的设计在全局算法和局部算法中都是研究的重点。

    2)代价聚合:

    通常全局算法不需要代价聚合,而局部算法需要通过求和、求均值或其他方法对一个支持窗口内的匹配代价进行聚合而得到参考图像上一点p在视差d处的累积代价 C A ( p , d ) CA(p,d) CA(p,d),这一过程称为代价聚合。通过匹配代价聚合,可以降低异常点的影响,提高信噪比进而提高匹配精度。代价聚合策略通常是局部匹配算法的核心,策略的好坏直接关系到最终视差图(Disparity maps)的质量。

    3)视差计算

    局部立体匹配算法的思想,在支持窗口内聚合完匹配代价后,获取视差的过程通常采用‘胜者为王’策略(WTA,Winner Take All),即在视差搜索范围内选择累积代价最优的点作为对应匹配点,与之对应的视差即为所求的视差。即P点的视差为 d = a r g m i n C A ( p , d ) d=arg minCA(p,d) d=argminCA(p,d)但是这种方法可能导致噪声大,因为如果直接用点的像素值坐标匹配,风险太大,找到的匹配点不一定是正确的,这时可以采用窗口匹配方法,关于这种方法在下面会做详细介绍。

    4)后处理:

    一般的,分别以左右两图为参考图像,完成上述三个步骤后可以得到左右两幅视差图像。但所得的视差图还存在一些问题,如遮挡点视差不准确、噪声点、误匹配点等存在,因此还需要对视差图进行优化,采用进一步执行后处理步骤对视差图进行修正。常用的方法有插值(Interpolation)、亚像素增强(Subpixel Enhancement)、精细化(Refinement)、图像滤波(Image Filtering)等操作。

    1.3 窗口匹配方法

    在上面提到过,如果使用点匹配方法的话,可能导致噪声大,因为直接找最匹配的点误差太大,如下图,比如我们要在(T)图中找(R)图标示的眼睛位置匹配点,由于两边眼睛的像素值相近,很可能找错,造成误差。
    在这里插入图片描述
    这时,我们就可以考虑改进的方法:窗口匹配方法。它的主要思路是:以当前点为中心,切出一个图像块,每个图像块依次做对比,找到最接近的图像块。这时就不只是拿单个点做对比了,可以有效减小误差。举例见下图:
    在这里插入图片描述

    匹配代价

    相似性测度函数用于度量参考图像中的匹配基元和目标图像中的匹配基元的相似性,即判断参考图像和目标图像中两点为对应匹配点的可能性,也称为匹配代价,这里记为 C ( x , y , d ) C(x,y,d) C(x,y,d),表示像素点 ( x , y ) (x,y) (x,y)在视差为d情况下的匹配误差。
    下面介绍最常见的三种匹配代价:
    (1)绝对差值和:(Sum of Absolute Differences,SAD C ( x , y , d ) = ∑ x ε S ∣ I R ( x , y ) − I T ( x + d , y ) ∣ C(x,y,d)=\sum_{x\varepsilon S}^{}\left | I_{R}(x,y)-I_{T}(x+d,y)\right | C(x,y,d)=xεSIR(x,y)IT(x+d,y)(2)差值平方和(Sum of squared Differences,SSD C ( x , y , d ) = ∑ x ε S ( I R ( x , y ) − I T ( x + d , y ) ) 2 C(x,y,d)=\sum_{x\varepsilon S}^{}(I_{R}(x,y)-I_{T}(x+d,y))^{2} C(x,y,d)=xεS(IR(x,y)IT(x+d,y))2(3)截断绝对差值和(Sum of Truncated Absolute Differences,STAD C ( x , y , d ) = ∑ x ε S m i n { ∣ I R ( x , y ) − I T ( x + d , y ) ∣ , T } C(x,y,d)=\sum_{x\varepsilon S}^{}min\left \{\left | I_{R}(x,y)-I_{T}(x+d,y)\right |,T\right \} C(x,y,d)=xεSmin{IR(x,y)IT(x+d,y),T}另外常用的匹配代价还有归一化互相关NCC(Normalized Cross Correlation),在下面做详细介绍。

    1.4 归一化互相关(NCC)

    原理:

    对于原始的图像内任意一个像素点 ( p x , p y ) (p_{x},p_{y}) (px,py)构建一个n x n的邻域作为匹配窗口。然后对于目标相素位置 ( p x + d , p y ) (p_{x}+d,p_{y}) (px+d,py)同样构建一个n x n大小的匹配窗口,对两个窗口进行相似度度量,注意这里的d有一个取值范围。对于两幅图像来说,在进行NCC计算之前要对图像处理,也就是将两帧图像校正到水平位置,即光心处于同一水平线上,此时极线是水平的,否则匹配过程只能在倾斜的极线方向上完成,这将消耗更多的计算资源。

    NCC计算公式:
    在这里插入图片描述
    其中 N C C ( p , d ) NCC(p,d) NCC(p,d)得到的值得范围将在[−1,1]之间。

    • W p W_{p} Wp为之前提到的匹配窗口。
    • I 1 ( x , y ) I_{1}(x,y) I1(x,y)为原始图像的像素值。
    • I 1 ˉ ( p x , p y ) \bar{I_{1}}(p_{x},p_{y}) I1ˉ(px,py)为原始窗口内像素的均值。
    • I 2 ( x + d , y ) I_{2}(x+d,y) I2(x+d,y)为原始图像在目标图像上对应点位置在x方向上偏移d后的像素值。
    • I 2 ˉ ( p x + d , p y ) \bar{I_{2}}(p_{x}+d,p_{y}) I2ˉ(px+d,py)为目标图像匹配窗口像素均值。

    若NCC=−1,则表示两个匹配窗口完全不相关,相反,若NCC=1时,表示两个匹配窗口相关程度非常高。

    匹配流程:

    • 采集图像:通过标定好的双目相机采集图像,当然也可以用两个单目相机来组合成双目相机。

    • 极线校正:校正的目的是使两帧图像极线处于水平方向,或者说是使两帧图像的光心处于同一水平线上。通过校正极线可以方便后续的NCC操作。

      由标定得到的内参中畸变信息中可以对图像去除畸变。
      通过校正函数校正以后得到相机的矫正变换R和新的投影矩阵P,接下来是要对左右视图进行去畸变,并得到重映射矩阵。

    • 特征匹配:这里便是我们利用NCC做匹配的步骤啦,匹配方法如上所述,右视图中与左视图待测像素同一水平线上相关性最高的即为最优匹配。完成匹配后,我们需要记录其视差d,即待测像素水平方向xl与匹配像素水平方向xr之间的差值 d = x r − x l d=x_{r}-x_{l} d=xrxl,最终我们可以得到一个与原始图像尺寸相同的视差图D。

    • 深度恢复:通过上述匹配结果得到的视差图D,我们可以很简单的利用相似三角形反推出以左视图为参考系的深度图。计算原理同上述1.1中介绍的几何原理(这里不做重述),即通过公式: Z = f T x r − x l Z=f\frac{T}{x_{r}-x_{l}} Z=fxrxlT可以简单得到以左视图为参考系的深度图了。

    二、实验内容

    2.1 实验目的与要求

    • 实现NCC 视差匹配方法,即给定左右两张视图,根据NCC计算视差图
    • 分析不同窗口值对匹配结果的影响,重点考查那些点(或者哪些类型的点)在不同窗口大小下的匹配精度影响

    2.2 代码实现

    # -*- coding: utf-8
    
    from pylab import *
    from numpy import *
    
    np.seterr(divide='ignore', invalid='ignore')
    from PIL import Image
    import numpy as np
    import math
    from scipy import ndimage
    
    
    def plane_sweep_ncc(im_l, im_r, start, steps, wid):
        """ 使用归一化的互相关计算视差图像 该函数返回每个像素的最佳视差"""
        m, n = im_l.shape
        # 保存不同求和值的数组
        mean_l = np.zeros((m, n))
        mean_r = np.zeros((m, n))
        s = np.zeros((m, n))
        s_l = np.zeros((m, n))
        s_r = np.zeros((m, n))
        # 保存深度平面的数组
        dmaps = np.zeros((m, n, steps))
        # 计算图像块的平均值
        ndimage.filters.uniform_filter(im_l, wid, mean_l)
        ndimage.filters.uniform_filter(im_r, wid, mean_r)
        # 归一化图像
        norm_l = im_l - mean_l
        norm_r = im_r - mean_r
        # 尝试不同的视差
        for displ in range(steps):
            # 将左边图像移动到右边,计算加和
            ndimage.filters.uniform_filter(np.roll(norm_l, -displ - start) * norm_r, wid, s)  # 和归一化
            ndimage.filters.uniform_filter(np.roll(norm_l, -displ - start) * np.roll(norm_l, -displ - start), wid, s_l)
            ndimage.filters.uniform_filter(norm_r * norm_r, wid, s_r)  # 和反归一化
            # 保存 ncc 的分数
            dmaps[:, :, displ] = s / np.sqrt(s_l * s_r)
        # 为每个像素选取最佳深度
        return np.argmax(dmaps, axis=2)
    
    
    def plane_sweep_gauss(im_l, im_r, start, steps, wid):
        """ 使用带有高斯加权周边的归一化互相关计算视差图像 """
        m, n = im_l.shape
        # 保存不同加和的数组
        mean_l = np.zeros((m, n))
        mean_r = np.zeros((m, n))
        s = np.zeros((m, n))
        s_l = np.zeros((m, n))
        s_r = np.zeros((m, n))
        # 保存深度平面的数组
        dmaps = np.zeros((m, n, steps))
        # 计算平均值
        ndimage.filters.gaussian_filter(im_l, wid, 0, mean_l)
        ndimage.filters.gaussian_filter(im_r, wid, 0, mean_r)
        # 归一化图像
        norm_l = im_l - mean_l
        norm_r = im_r - mean_r
        # 尝试不同的视差
        for displ in range(steps):
            # 将左边图像移动到右边,计算加和
            ndimage.filters.gaussian_filter(np.roll(norm_l, -displ - start) * norm_r, wid, 0, s)  # 和归一化
            ndimage.filters.gaussian_filter(np.roll(norm_l, -displ - start) * np.roll(norm_l, -displ - start), wid, 0, s_l)
            ndimage.filters.gaussian_filter(norm_r * norm_r, wid, 0, s_r)  # 和反归一化
            # 保存 ncc 的分数
            dmaps[:, :, displ] = s / np.sqrt(s_l * s_r)
        # 为每个像素选取最佳深度
        return np.argmax(dmaps, axis=2)
    
    
    """载入图像,并使用该函数计算偏移图"""
    im_l = np.array(Image.open('../data/data_six/p1.png').convert('L'), 'f')
    im_r = np.array(Image.open('../data/data_six/p2.png').convert('L'), 'f')
    figure()
    gray()
    subplot(221)
    imshow(im_l)
    subplot(222)
    imshow(im_r)
    
    # 开始偏移,并设置步长
    steps = 12
    start = 4
    # ncc 的宽度
    wid = 13
    res = plane_sweep_ncc(im_l, im_r, start, steps, wid)
    wid2 = 5
    gauss = plane_sweep_gauss(im_l, im_r, start, steps, wid2)
    subplot(223)
    imshow(res)
    subplot(224)
    imshow(gauss)
    show()
    import scipy.misc
    
    scipy.misc.imsave('depth.png', res)
    

    2.3 实验结果与分析

    图片来源:
    下载链接为: 视差匹配实验图片集.

    NCC版本与高斯版本的视差图:
    在这里插入图片描述

    分析:如上图所示,在标准版本中设置wid=9,高斯版本中设置wid=3。上面一行图像所示为图像对,左下方图像是标准的NCC扫平面法重建的视差图,右下方是高斯版本的视差图。可以看到,与标准版本相比,高斯版本具有较少的噪声,但缺少很多细节信息。

    NCC视差匹配及不同窗口值对匹配结果的影响:
    为了实现一次性输出不同wid值下的视差图,修改代码:

    wid = np.array([1,3,5,7,9,11])
    res=[]
    for i in range(0, 6):
        res.append(plane_sweep_ncc(im_l, im_r, start, steps, wid[i]))
    a = np.array([231,232,233,234,235,236])
    color_d = ['ncc wid=1', 'ncc wid=3', 'ncc wid=5','ncc wid=7','ncc wid=9','ncc wid=11']
    for i in range(0, 6):
        subplot(a[i])
        imshow(res[i])
        title(color_d[i])
        axis('off')
    

    原图:
    在这里插入图片描述NCC视差匹配(不同窗口值):
    分别为wid=1、3、5、7、9、11的视差图:
    在这里插入图片描述
    分析:

    • 如上图所示,首先对单张NCC视差匹配的结果做分析,比如观察当wid=9时的结果(如下图):
      在这里插入图片描述
      可以看到离相机越近的物体亮度越高,比如靠近相机的笔,它们整体呈现的更亮,而后面的背景则更暗。这是因为前景的位移比背景的位移更多,越靠近摄像机的目标,它在左右视图位移的距离会更多,这是由近大远小导致。结合原理来分析,我们可以观察公式: Z = f T x r − x l Z=f\frac{T}{x_{r}-x_{l}} Z=fxrxlT发现 Z 和 x r − x l Z和x_{r}-x_{l} Zxrxl呈反比关系,因此视差越大,Z越小,即目标越近,视差越大,深度越小。
    • 观察不同窗口值(wid)的运行结果,初步观察可得知:随着wid值的增大,图像深度信息的图像的轮廓逐渐明显,可以看到离照相机最近的目标是笔,随着wid增大,笔的轮廓逐渐清晰,但较远处的物体,比如墙壁以及画布,则越来越难以辨认。
    • 进一步观察结果,这次拿出wid=3与wid=11的结果作对比,如下图:
      在这里插入图片描述 可以看到wid越小,能得到更多的细节,比如见wid=3时的结果图,虽然图像整体比较模糊,但是可以看到我们能够得到更多笔的细节信息,从原理分析,当wid值较小,进行滤波时考虑的因素少,只受某像素点自身影响以及周围少量像素点的影响,所以保留下来更多的自身特征。但同时噪声也大,见wid=3时的背景墙,上面出现了许多噪声点。
    • 而当窗口值增大时,鲁棒性增强,噪声明显减少,但我们可以发现图像的细节也在逐渐减少,可物体的轮廓逐渐清晰,尤其是近景目标,如wid=11的视差图,我们能够更清晰地看出笔的轮廓,但是由于细节减少,远处的物体就比较不易判断。在wid值较大情况下,NCC匹配下充分考虑了周围像素点的影响,能更好的观察出视差。

    三、总结

    实验内容总结:
        经过本次实验,主要实现了NCC视差匹配,发现越靠近摄像机的目标,视差越大,因此视差图中近景目标的图像更亮。对于不同wid值,wid越小,能得到更多的细节,但是噪声大;wid越大,具有更好的鲁棒性,但是得到的细节少。

    实验过程总结:
    ①问题:在使用高斯滤波器做视差匹配时报错:
    在这里插入图片描述解决办法:改变sqrt()函数调用的库,即将原本的从math库调用改为从numpy库调用,修改代码:

    dmaps[:, :, displ] = s / np.sqrt(s_l * s_r)
    

    ②问题:输出结果模糊不清,无法看到轮廓。
    解决办法:调整图片的像素值,而且注意左右两张图拍摄时位移不能太大。

    展开全文
  • 立体匹配NCC算法代码

    千次阅读 2015-09-21 22:25:35
    NCC算法是立体匹配较为经典的一个算法,其定义为: 其中Wp为以像素p=(px,py)为中心的匹配窗口。I1为左图匹配窗口内的像素,I2为右图匹配窗口内的像素。 下面是用C++编写的算法:/******************************...
                   NCC算法是立体匹配较为经典的一个算法,其定义为:

    其中Wp为以像素p=(px,py)为中心的匹配窗口。I1为左图匹配窗口内的像素,I2为右图匹配窗口内的像素。


    下面是用C++编写的算法:

    <pre name="code" class="cpp">/*********************************************** 
     *
     *   left_sq存放了左图窗口内像素与均值差值的平方
     *   right_sq存放了右图窗口内像素与均值差值的平方
     *   left_avg存放了左图窗口内像素的均值
     *   right_avg存放了右图窗口内像素的均值
     *
     **************************************************/
    void compute_sq(IplImage* left_img, IplImage* right_img, float *left_sq, float *right_sq, float *left_avg, float *right_avg)
    {
    	//图像的高度和宽度
    	int height = left_img->height;
    	int width = left_img->width;
    	//窗口半径,为奇数
    	int N = 5;
    	//图像匹配的起始行和终止行
    	int line_start = N;
    	int line_end = height-N;
    	//图像需要视差搜索的起始列和终止列
    	int row_start = N;
    	int row_end = width-N;
    	int addr = 0;
    	float temp_l = 0, temp_r = 0, suml = 0, sumr = 0;
    
    
    	for (int j = line_start; j < line_end; j++)
    	{
    		for (int i = row_start; i < row_end; i++)
    		{
    			suml = 0.0, sumr = 0.0;
    			temp_l = 0.0; temp_r = 0.0;
    			for (int m = j - N; m <= j + N; m++)
    			{
    				for (int n = i - N; n <= i + N; n++)
    				{
    					suml += ((uchar*)(left_img->imageData + m*left_img->widthStep))[n];
    					//cout << "l_px:" << (int)((uchar*)(left_img->imageData + m*left_img->widthStep))[n] << endl;
    					sumr += ((uchar*)(right_img->imageData + m*right_img->widthStep))[n];
    					//cout << "r_px:" << (int)((uchar*)(right_img->imageData + m*right_img->widthStep))[n]<<endl;
    				}
    			}
    			addr = j*width + i;
    			left_avg[addr] = suml / pow((2 * N + 1), 2);
    			right_avg[addr] = sumr / pow((2 * N + 1), 2);
    			//cout << "l_avg:" << (float)left_avg[addr]<<endl;
    			//cout << "r_avg:" << (float)right_avg[addr]<<endl;
    			for (int m = j - N; m <= j + N; m++)
    			{
    				for (int n = i - N; n <= i + N; n++)
    				{
    					temp_l += pow((((uchar*)(left_img->imageData + m*left_img->widthStep))[n] - left_avg[addr]), 2);
    					temp_r += pow((((uchar*)(right_img->imageData + m*right_img->widthStep))[n] - right_avg[addr]), 2);
    				}
    
    				
    			}
    			left_sq[addr] = temp_l;
    			right_sq[addr] = temp_r;
    			//cout << "l_sq:" << (float)left_sq[addr] << endl;
    			//cout << "r_sq:" << (float)right_sq[addr] << endl;
    		}
    		
    
    	}
    
    }
    void compute_DP(IplImage* left_img, IplImage* right_img, IplImage* depth_img, float *left_sq, float *right_sq, float *left_avg, float *right_avg)
    	{
    	    //图像的高度和宽度
    		int height = left_img->height;
    		int width = left_img->width;
    		//窗口半径,为奇数
    		int N = 5; 
    		//搜索的视差值范围
    		int maxD = 10;
    		int minD = 1;
    		//图像匹配的起始行和终止行
    		int line_start = N;
    		int line_end = height-N;
    		//图像需要视差搜索的起始列和终止列
    		int row_start = N;
    		int row_end = width-N-maxD;
    
    		int addr;
    		float max_tmp, cov_tmp;
    		//视差
    		int vd = 0;
    
    		for (int j = line_start; j < line_end; j++)
    		{
    			for (int i = row_start; i < row_end; i++)
    			{
    				for (int d = minD; d <= maxD; d++)
    				{
    					cov_tmp = 0.0;
    					addr = j*width + i;
    					for (int m = j - N; m <= j + N; m++)
    					{
    						for (int n = i - N; n <= i + N; n++)
    						{
    							//cout << "l_px:" << (int)((uchar*)(left_img->imageData + m*left_img->widthStep))[n] << endl;
    							//cout << "l-avg:" << (float)(((uchar*)(left_img->imageData + m*left_img->widthStep))[n] - left_avg[addr]) << endl;
    							//cout << "r_px:" << (int)((uchar*)(right_img->imageData + m*right_img->widthStep))[n + d] << endl;
    							//cout << "r-avg:" << (float)(((uchar*)(right_img->imageData + m*right_img->widthStep))[n + d] - right_avg[addr+d])<<endl;
    							cov_tmp += (((uchar*)(left_img->imageData + m*left_img->widthStep))[n]-left_avg[addr])*(((uchar*)(right_img->imageData + m*right_img->widthStep))[n + d]-right_avg[addr+d]);
    						}
    					}
    					//cout << "a:" << cov_tmp << endl;		
    					cov_tmp = (float)(cov_tmp / sqrt(left_sq[addr] * right_sq[addr+d]));
    					//cout << "cov_tmp:" << cov_tmp << endl;
    					if (d == minD)
    					{
    						max_tmp = cov_tmp;
    						vd = d;
    					}
    					else
    					{
    						if (cov_tmp > max_tmp)
    						{
    							max_tmp = cov_tmp;
    							vd = d;
    						}
    					}
    				}
    				vd = vd * 20;
    				if (vd > 255)
    				{
    					vd = 255;
    				}
    				else if (vd < 0)
    				{
    					vd = 0;
    				}
    				//cout << "d:" << vd << endl;
    				((uchar*)(depth_img->imageData + j*depth_img->widthStep))[i] = vd;
    			}
    		}	
    }
    void NCC_Match(IplImage* left_img, IplImage* right_img, IplImage* depth_img)
    {
    	float *left_sq, *right_sq, *left_avg, *right_avg;
    	int img_size = left_img->height*left_img->width;
    
    	left_sq = (float *)malloc(img_size * sizeof(float));
    	right_sq = (float *)malloc(img_size * sizeof(float));
    	left_avg = (float *)malloc(img_size * sizeof(float));
    	right_avg = (float *)malloc(img_size * sizeof(float));
    
    	memset(left_sq, 0, img_size);
    	memset(right_sq, 0, img_size);
    
    	compute_sq(left_img, right_img, left_sq, right_sq,left_avg,right_avg);
    	compute_DP(left_img, right_img, depth_img, left_sq, right_sq,left_avg,right_avg);
    
    	free(left_sq);
    	free(right_sq);
    
    	return;
    }

     
    
     
    

          在配置好opencv的VS环境下,输出标准的左右图运行即可,注意标准图是不需要做矫正的图,网上可以下到几张经典的图片。

          运行结果:

       

        效果只能说非常的一般,只适合初学者入门学习。

       转载请注明出处:http://blog.csdn.net/lvlitc/article/details/46431245

    展开全文
  • 文章目录实验原理视差图计算归一化互相关双目立体匹配的步骤实验内容代码实现实验过程实验总结遇到的问题 实验原理 视差图计算 深度信息可以通过计算1幅图像和其它图像的特征位置的像素差获得。视差图和深度图很像,...
  • 6SSD算法,SAD算法,NCC算法,DP算法,DP_5算法,census算法
  • 立体匹配1.1 概述1.2 主要立体匹配算法分类1.3 立体匹配的基本步骤2. 归一化互相关(NCC)视差匹配法2.1 原理2.2 匹配流程:3. 不同窗口值对匹配结果的影响 实验要求: 从理论角度,分析以窗口代价计算视差的原理 ...

    1. 立体匹配

    1.1 概述

    立体匹配是立体视觉研究中的关键部分。其目标是在两个或多个视点中匹配相应像素点,计算视差。通过建立一个能量代价函数,对其最小化来估计像素点的视差,求得深度
    在这里插入图片描述
    点P和Q,映射到左相机OR像面上的同一点p≡q,只要找到p和q在右相机OT像面上的对应点就可以通过三角计算估计深度找到对应点的过程,即立体匹配

    为了找到对应点,需要增加约束,最常用的是极线约束

    在这里插入图片描述
    P和Q映射到左相机QR像面上的同一点p≡q,直线pq上的点对应点一定位于右像面的直线p’q’上,p’q’即为直线pq的极线,这就是极线约束

    接下来就可以根据视差估计深度,然后通过Graph cuts算法给每一个像素点分配视差从而得到深度图,不再详细说明。

    1.2 主要立体匹配算法分类

    1)根据采用图像表示的基元不同,立体匹配算法分为:
    A、区域立体匹配算法:可获取稠密视差图。缺点:受图像的仿射畸变和辐射畸变影响较大;像素点约束窗口的大小与形状选择比较困难,选择过大,在深度不连续处,视差图中会出现过度平滑现象;选择过小,对像素点的约束比较少,图像信息没有得到充分利用,容易产生误匹配。
    B、基于特征的立体匹配算法:可获得稀疏的视差图,经差值估计可获得稠密视差图。可提取点、线、面等局部特征,也可提取多边形和图像结构等全局特征。缺点:特征提取易受遮挡、光线、重复纹理等影响较大;差值估计计算量大。
    C、基于相位立体匹配算法:假定在图像对应点中,其频率范围内,其局部相位是相等的,在频率范围内进行视差估计。
    2)依据采用最优化理论方法的不同,立体匹配算法可以分为:
    A、局部的立体匹配算法
    B、全局的立体匹配算法

    还有立体匹配算法介绍的更详细内容:参考博客

    1.3 立体匹配的基本步骤

    立体匹配过程:

    1. 匹配代价计算: 一般是通过计算左右两图对应像素3个通道的灰度值差来决定匹配代价的,常用的就是基于像素点匹配代价计算,一般有AD,SD,TAD什么的,基于区域的匹配代价计算一般有SAD,SSD, STAD之类的。匹配代价计算会生成一个disparity space image,也就是DSI。这个DSI是一个三维的空间,也就是每一个视差,得到一张代价图。假如视差范围是0~16,则会得到17幅代价图。视差搜索范围就是MiddleBurry网站上的stereopair值,也就是说在视差范围内(比如0-16)内搜索匹配代价,得到17张匹配代价图,然后找到匹配代价最小的对应的视差值就是此像素对应的视差
    2. 代价聚合:其实就是一个滤波的过程,对每一幅代价图进行聚合,最简单的就是采用boxfilter。第一步代价计算只是得到了图像上所有孤立像素的视差值,但是这些时差值都是孤立的,引入了过多噪声,比如一片区域的视差值都是10,可是引入噪声后就会导致这一片的视差值都不一样,那么就需要一个滤波的过程,也就是我们所说的局部立体匹配方法,即采用窗口卷积达到局部滤波的目的
    3. 计算视差:常用的方法就是WTA算法(局部),对于图像中的同一个点,选出17幅代价图中匹配代价最小的那张图,该幅图对应的视差值就选取为最终的视差。或者在全局立体匹配中采用能量函数的方法,分为数据项和平滑项,数据项其实就是代价计算,平滑项就是代价聚合,只不过窗口大小是整幅图像,也可以试试如果把平滑项前面的系数lamda设为0,那么得到的结果和单纯代价计算的局部立体匹配是一样的。
    4. 视差精化:也就是对得到的视差进行优化的过程,如:左右一致性检测、区域投票等;这步其实是很多立体匹配的遮羞布,比如用遮挡处理,中值滤波,左右一致性检测等,都能使最后的是视差图提升1%左右,它是很多论文的遮羞布。但是不可否认的是,立体匹配最关键的步骤仍然是代价计算和代价聚合步骤。
      在立体匹配方法中,基于全局和局部的算法有些区别。不过基本步骤都差不多。有些时候,基于局部的算法,第一步和第二步是合并在一起进行的,基于全局的算法,会跳过第二步。

    2. 归一化互相关(NCC)视差匹配法

    2.1 原理

    对于原始的图像内任意一个像素点 ( p x , p y ) (px,py) (px,py)构建一个 n × n n×n n×n的邻域作为匹配窗口。然后对于目标相素位置 ( p x + d , p y ) (px+d,py) (px+d,py)同样构建一个 n × n n×n n×n大小的匹配窗口,对两个窗口进行相似度度量,注意这里的dd dd有一个取值范围。对于两幅图像来说,在进行 N C C NCC NCC计算之前要对图像处理,也就是将两帧图像校正到水平位置,即光心处于同一水平线上,此时极线是水平的,否则匹配过程只能在倾斜的极线方向上完成,这将消耗更多的计算资源。

    N C C NCC NCC计算公式如下:
    N C C ( p , d ) = ∑ ( x , y ) ∈ W p ( I 1 ( x , y ) − I 1 ‾ ( P x , P y ) ) ⋅ ( I 1 ( x + d , y ) − I 2 ‾ ( P x + d , P y ) ) ∑ ( x , y ) ∈ W p ( I 1 ( x , y ) − I 1 ‾ ( P x , P y ) ) 2 ⋅ ∑ ( x , y ) ∈ W p ( I 2 ( x + d , y ) − I 2 ‾ ( P x + d , P y ) ) 2 NCC(p,d)=\frac{\sum_{(x,y)\in W_{p}}(I_{1}(x,y)-\overline{I_{1}}(P_{x},P_{y}))\cdot (I_{1}(x+d,y)-\overline{I_{2}}(P_{x}+d,P_{y}))}{\sum_{(x,y)\in W_{p}}(I_{1}(x,y)-\overline{I_{1}}(P_{x},P_{y}))^2\cdot \sum_{(x,y)\in W_{p}}(I_{2}(x+d,y)-\overline{I_{2}}(P_{x}+d,P_{y}))^2} NCC(p,d)=(x,y)Wp(I1(x,y)I1(Px,Py))2(x,y)Wp(I2(x+d,y)I2(Px+d,Py))2(x,y)Wp(I1(x,y)I1(Px,Py))(I1(x+d,y)I2(Px+d,Py))

    其中 N C C ( p , d ) NCC(p,d) NCC(p,d)得到的值的范围将在[−1,1]之间。
    W p Wp Wp为之前提到的匹配窗口
    I 1 ( x , y ) I_{1}(x,y) I1(x,y)原始图像的像素值
    I 1 ‾ ( p x , p y ) \overline{I_{1}}(p_x,p_y) I1(px,py)原始窗口内像素的均值
    I 2 ( x + d , y ) I_{2}(x+d,y) I2(x+d,y)原始图像在目标图像上对应点位置在 x x x方向上偏移 d d d后的像素值
    I 2 ‾ ( p x + d , p y ) \overline{I_{2}}(p_x+d, p_y) I2(px+d,py)目标图像匹配窗口像素均值

    N C C = − 1 NCC=−1 NCC=1,则表示两个匹配窗口完全不相关,相反,若 N C C = 1 NCC=1 NCC=1时,表示两个匹配窗口相关程度非常高。

    2.2 匹配流程

    1. 采集图像:通过标定好的双目相机采集图像,当然也可以用两个单目相机来组合成双目相机。

    2. 极线校正:校正的目的是使两帧图像极线处于水平方向,或者说是使两帧图像的光心处于同一水平线上。通过校正极线可以方便后续的 N C C NCC NCC操作。
      1)由标定得到的内参中畸变信息中可以对图像去除畸变
      2)通过校正函数校正以后得到相机的矫正变换R和新的投影矩阵P,接下来是要对左右视图进行去畸变,并得到重映射矩阵

    3. 特征匹配:这里便是我们利用 N C C NCC NCC做匹配的步骤啦,匹配方法如上所述,右视图中与左视图待测像素同一水平线上相关性最高的即为最优匹配。完成匹配后,我们需要记录其视差 d d d,即待测像素水平方向 x l xl xl与匹配像素水平方向 x r xr xr之间的差值 d = x r − x l d=xr−xl d=xrxl,最终我们可以得到一个与原始图像尺寸相同的视差图 D D D

    4. 深度恢复:通过上述匹配结果得到的视差图 D D D,我们可以很简单的利用相似三角形反推出以左视图为参考系的深度图。计算原理如下图所示:
      在这里插入图片描述

    如图, T x Tx Tx为双目相机基线, f f f为相机焦距,这些可以通过相机标定步骤得到。而 x r − x l xr−xl xrxl就是视差 d d d
    通过公式 z = f × T x d z=\frac{f×Tx}{d} z=df×Tx可以很简单地得到以左视图为参考系的深度图了。
    至此,我们便完成了双目立体匹配。倘若只是用于图像识别,那么到步骤3时已经可以结束了。

    2.3 代码实现

    NCCfaster.py

    import numpy as np
    import cv2
    from PCV.localdescriptors import sift
    
    im1 = 'C://Users//Garfield//PycharmProjects//untitled//NCC-master//im2.ppm'
    im2 = 'C://Users//Garfield//PycharmProjects//untitled//NCC-master//im6.ppm'
    
    img1 = cv2.imread(im1, cv2.CV_8UC1)
    img2 = cv2.imread(im2, cv2.CV_8UC1)
    rows, cols = img1.shape
    print(img1.shape)
    
    def translaton(image, shape):
        step = round((shape[0]-1)/2)
        print(step)
        shifted = []
        for i in range(0, step+1):
            for j in range(0, step+1):
                if i==0 and j==0:
                    M1 = np.float32([[1, 0, i], [0, 1, j]])
                    shifted.append(cv2.warpAffine(image, M1, (image.shape[1], image.shape[0])))
                elif i==0 and j!=0:
                    M1 = np.float32([[1, 0, i], [0, 1, j]])
                    M2 = np.float32([[1, 0, i], [0, 1, -j]])
                    shifted.append(cv2.warpAffine(image, M1, (image.shape[1], image.shape[0])))
                    shifted.append(cv2.warpAffine(image, M2, (image.shape[1], image.shape[0])))
                elif i!=0 and j==0:
                    M1 = np.float32([[1, 0, i], [0, 1, j]])
                    M2 = np.float32([[1, 0, -i], [0, 1, j]])
                    shifted.append(cv2.warpAffine(image, M1, (image.shape[1], image.shape[0])))
                    shifted.append(cv2.warpAffine(image, M2, (image.shape[1], image.shape[0])))
                else:
                    M1 = np.float32([[1, 0, i], [0, 1, j]])
                    M2 = np.float32([[1, 0, -i], [0, 1, -j]])
                    M3 = np.float32([[1, 0, -i], [0, 1, j]])
                    M4 = np.float32([[1, 0, i], [0, 1, -j]])
                    shifted .append(cv2.warpAffine(image, M1, (image.shape[1], image.shape[0])))
                    shifted.append(cv2.warpAffine(image, M2, (image.shape[1], image.shape[0])))
                    shifted.append(cv2.warpAffine(image, M3, (image.shape[1], image.shape[0])))
                    shifted.append(cv2.warpAffine(image, M4, (image.shape[1], image.shape[0])))
    
        print(len(shifted))
        return np.array(shifted)
    
    #I(x,y)-avg(I(x,y))
    def img_sub_avg(img_shifted, avg_img):
        len, height, width = img1_shifted.shape
        tmp_ncc1 = np.zeros([len, height, width])
        for i in range(len):
            tmp_ncc1[i] = img_shifted[i] - avg_img
        print(tmp_ncc1)
        return tmp_ncc1
    
    def NCC(img1_sub_avg,img2_sub_avg, threshold, max_d):
        #设立阈值
        len, height, width = img1_sub_avg.shape
        thershould_shifted = np.zeros([len, height, width])
        ncc_max = np.zeros([height, width])
        ncc_d = np.zeros([height, width])
        for j in range(3, max_d):
            tmp_ncc1 = np.zeros([height, width])
            tmp_ncc2 = np.zeros([height, width])
            tmp_ncc3 = np.zeros([height, width])
            for k in range(len):
                M1 = np.float32([[1, 0, -j - 1], [0, 1, 0]])
                thershould_shifted[k] = cv2.warpAffine(img1_sub_avg[k], M1, (img1_sub_avg.shape[2], img1_sub_avg.shape[1]))
            for i in range(len):
                tmp_ncc1 += (img2_sub_avg[i])*(thershould_shifted[i])
                tmp_ncc2 += pow(img2_sub_avg[i], 2)
                tmp_ncc3 += pow(thershould_shifted[i], 2)
    
            tmp_ncc2 = tmp_ncc2*tmp_ncc3
            tmp_ncc2 = np.sqrt(tmp_ncc2)
            tmp_ncc4 = tmp_ncc1/tmp_ncc2
            for m in range(height):
                for n in range(width):
                    if tmp_ncc4[m, n] > ncc_max[m ,n] and tmp_ncc4[m, n] > threshold:
                        ncc_max[m, n] = tmp_ncc4[m, n]
                        ncc_d[m , n] = j
        for i in ncc_d:
            print(i)
        return ncc_max, ncc_d
    
    if __name__ == "__main__":
        disparity = np.zeros([rows, cols])
        NCC_value = np.zeros([rows, cols])
        deeps = np.zeros([rows, cols])
        # 用3*3卷积核做均值滤波
        avg_img1 = cv2.blur(img1, (7, 7))
        avg_img2 = cv2.blur(img2, (7, 7))
        fimg1 = img1.astype(np.float32)
        fimg2 = img2.astype(np.float32)
        avg_img1 = avg_img1.astype(np.float32)
        avg_img2  = avg_img2.astype(np.float32)
        img1_shifted = translaton(fimg1, [7, 7])
        img2_shifted = translaton(fimg2, [7, 7])
        img1_sub_avg = img_sub_avg(img1_shifted, avg_img1)
        img2_sub_avg = img_sub_avg(img2_shifted, avg_img2)
        ncc_max, ncc_d = NCC(img1_sub_avg,img2_sub_avg, threshold = 0.5, max_d = 64)
        print(img1_shifted.shape)
        disp = cv2.normalize(ncc_d, ncc_d, alpha=0, beta=255, norm_type=cv2.NORM_MINMAX,
                              dtype=cv2.CV_8U)
        cv2.imshow("left", img1)
        cv2.imshow("right", img2)
        cv2.imshow("depth", disp)
        cv2.waitKey(0)  # 等待按键按下
        cv2.destroyAllWindows()#清除所有窗口
    

    原始图像:
    在这里插入图片描述

    运行结果:
    在这里插入图片描述
    在这里插入图片描述
    代码2:

    # -*- coding: utf-8 -*-
    import scipy.misc
    from PIL import Image
    from pylab import *
    import cv2
    from numpy import *
    from numpy.ma import array
    from scipy.ndimage import filters
    def plane_sweep_ncc(im_l,im_r,start,steps,wid):
        """ 使用归一化的互相关计算视差图像 """
        m,n = im_l.shape
        # 保存不同求和值的数组
        mean_l = zeros((m,n))
        mean_r = zeros((m,n))
        s = zeros((m,n))
        s_l = zeros((m,n))
        s_r = zeros((m,n))
        # 保存深度平面的数组
        dmaps = zeros((m,n,steps))
        # 计算图像块的平均值
        filters.uniform_filter(im_l,wid,mean_l)
        filters.uniform_filter(im_r,wid,mean_r)
        # 归一化图像
        norm_l = im_l - mean_l
        norm_r = im_r - mean_r
        # 尝试不同的视差
        for displ in range(steps):
            # 将左边图像移动到右边,计算加和
            filters.uniform_filter(np.roll(norm_l, -displ - start) * norm_r, wid, s) # 和归一化
            filters.uniform_filter(np.roll(norm_l, -displ - start) * np.roll(norm_l, -displ - start), wid, s_l)
            filters.uniform_filter(norm_r*norm_r,wid,s_r) # 和反归一化
            # 保存 ncc 的分数
            dmaps[:,:,displ] = s / sqrt(s_l * s_r)
            # 为每个像素选取最佳深度
        return np.argmax(dmaps, axis=2)
    
    def plane_sweep_gauss(im_l,im_r,start,steps,wid):
     """ 使用带有高斯加权周边的归一化互相关计算视差图像 """
     m,n = im_l.shape
     # 保存不同加和的数组
     mean_l = zeros((m,n))
     mean_r = zeros((m,n))
     s = zeros((m,n))
     s_l = zeros((m,n))
     s_r = zeros((m,n))
     # 保存深度平面的数组
     dmaps = zeros((m,n,steps))
     # 计算平均值
     filters.gaussian_filter(im_l,wid,0,mean_l)
     filters.gaussian_filter(im_r,wid,0,mean_r)
     # 归一化图像
     norm_l = im_l - mean_l
     norm_r = im_r - mean_r
     # 尝试不同的视差
     for displ in range(steps):
         # 将左边图像移动到右边,计算加和
         filters.gaussian_filter(np.roll(norm_l, -displ - start) * norm_r, wid, 0, s) # 和归一化
         filters.gaussian_filter(np.roll(norm_l, -displ - start) * np.roll(norm_l, -displ - start), wid, 0, s_l)
         filters.gaussian_filter(norm_r*norm_r,wid,0,s_r) # 和反归一化
         # 保存 ncc 的分数
         dmaps[:,:,displ] = s / np.sqrt(s_l * s_r)
     # 为每个像素选取最佳深度
     return np.argmax(dmaps, axis=2)
    
    im_l = array(Image.open(r'C://Users//Garfield//Desktop//towelmatch//jidian//1.jpg').convert('L'), 'f')
    im_r = array(Image.open(r'C://Users//Garfield//Desktop//towelmatch//jidian//2.jpg').convert('L'),'f')
    # 开始偏移,并设置步长
    steps = 12
    start = 4
    # ncc 的宽度
    wid = 2000
    res = plane_sweep_ncc(im_l,im_r,start,steps,wid)
    
    imsave('C://Users//Garfield//PycharmProjects//untitled//NCC-master//depth2000.png',res)
    show()
    
    

    此代码可以修改窗口值的大小,具体结果与分析放在下一小节

    3. 不同窗口值对匹配结果的影响

    原图:
    在这里插入图片描述

    以下分别是当窗口值设定为15,20,50,100,200,500,700,1000,2000的运行结果:
    在这里插入图片描述
    在这里插入图片描述

    总的来看,窗口值设置的越大,区域化分得越为明显,而设置的越小噪声越大,区块划分的越不明显。同时,从本次实验结果图来看,wid设置为700时,结果最为理想。
    另外可以看到
    在这里插入图片描述
    瓶身光照最为强烈的地方会与瓶身其他其他地方像素值有较大区分,以及窗口大小设置越大,则会更多的提取出有光照的像素点并连结成块,区分出亮部暗部。

    4. 实验遇到的问题与解决

    问题:

    • ‘NoneType’ object has no attribute ‘shape’

    求解cols,rows = img.shape时出现上述报错,是图片路径不对,程序无法获取到图片。

    解决办法:
    右键下图目录灰色处
    在这里插入图片描述
    —>copypath,再在原本的读取图片路径代码中粘贴
    注意格式改成C://Users//…
    参考博客:立体匹配_数据结构与算法
    参考博客:立体匹配过程
    参考博客:双目立体匹配算法–归一化互相关(NCC)详解和代码实现(python)
    NCC的更详细相关
    延伸学习:真实场景的双目立体匹配(Stereo Matching)获取深度图详解
    原理详解:双目视觉(三)立体匹配算法

    展开全文
  • MATLAB关于立体匹配,实现方法SAD SSD NCC CENSUS,适合初学者。
  • 6SSD算法,SAD算法,NCC算法,DP算法,DP_5算法,census算法
  • 用于立体匹配的三种算法,显示视差图,matlab代码
  • Stereo Match 立体匹配算法,包括SSD SAD SSD NCC NSSD Census 等基本算法,。用C/C++语言编写,适合初学者学习!
  • 双目立体匹配的原理类似于人类的视觉原理:分别单独使用左眼和单独使用右眼去观察一个物体,该物体在两个视角中的相对位置会存在一定的偏移,人类视觉系统可以根据该偏移去判断该物体的距离远近。 同样的,在双目...
  • NCC和SSDA算法的图像匹配实现

    热门讨论 2012-11-05 17:10:43
    NCC和SSDA算法的图像匹配实现,图片的读取用opencv实现,算法是纯C++代码。
  • Opencv——立体匹配

    千次阅读 2019-11-26 17:16:31
    Opencv——立体匹配 立体匹配 立体匹配,就是匹配两个不同的摄像头视图的3D点,只有在两个摄像头的重叠视图内的可视区域才能被计算,这就是为什么要使摄像头靠近前向平行了。 立体匹配的目的是通过匹配得到视差。 ...
  • 计算机视觉7—立体视觉之NCC算法

    千次阅读 2020-04-26 17:33:15
    1.立体视觉 一个多视图成像的特殊例子是立体视觉(或者立体成像),即使用两台只有水平(向一侧)偏移的照相机观测同一场景。当照相机的位置如上设置,两幅图像具有相同的图像平面,图像的行是垂直对齐的,那么称...
  • 老外写的NCC匹配算法,备注十分详细,通俗易懂,对新接触到双目视觉立体匹配的学生很有参考价值,还有SAD算法, SSD算法我也从不同老外的要到了,如果需要的话,也可以去我的主页下载
  • 归一化相关性,normalization cross-correlation,因此简称NCC,下文... NCC,顾名思义,就是用于归一化待匹配目标之间的相关程度,注意这里比较的是原始像素。通过在待匹配像素位置p(px,py)构建3*3邻域匹配窗口
  • 计算SSD 与NCC用于匹配 完整matlab程序及图片 计算SSD 与NCC用于匹配 完整matlab程序及图片 计算SSD 与NCC用于匹配 完整matlab程序及图片
  • 在传统的NCC算法上采用卷积加速
  • 对于原始的图像内任意一个像素点(px,py)(p_x,p_y)(px​,py​)构建一个n×nn\times nn×n的邻域作为匹配窗口。然后对于目标相素位置(px+d,py)(p_x+d, p_y)(px​+d,py​)同样构建一个n×nn\times nn×n大小的匹配窗口...
  • 文章目录1.NCC匹配介绍1.1原理1.2匹配流程2.代码3.结果和分析3.1实验结果展示3.2小结4.遇到的问题及解决方法 1.NCC匹配介绍 1.1原理 对于原始的图像内任意一个像素点(px,py) 构建一个n×n的邻域作为匹配窗口。然后...
  • NCC,normalization cross-correlation(归一化相关性),用于归一化待匹配目标之间的相关程度。 对于原始的图像内任意一个像素点(px,py)(px,py)(px,py)构建一个n×nn×nn×n 的邻域作为匹配窗口。然后对于目标相素...
  • OpenCv中实现了三种立体匹配算法:

    千次阅读 2016-03-23 16:39:06
    OpenCv中实现了三种立体匹配算法: BM算法 SGBM算法 Stereo Processing by Semiglobal Matching and Mutual Information GC算法 算法文献:Realistic CG Stereo Image Dataset with Ground Truth ...
  • 老外写的SSD匹配算法,备注十分详细,通俗易懂,对新接触到双目视觉立体匹配的学生很有参考价值,还有ncc算法, SAD算法我也从不同老外的要到了,如果需要的话,也可以去我的主页下载
  • 4.SGM 根据论文Accurate and Efficient Stereo Processing by Semi-Global Matching and Mutual Information写的双目立体匹配代码,matlab编写方便阅读,带测试图片,注意算法只实现了4个方向,即左右,右左,上下,...
  • NCC视差匹配方法

    2020-04-26 17:17:28
    文章目录一、实验原理1.1 窗口代价计算视差的...立体匹配算法的原理: 就是找出两张图像的对应关系,根据三角测量原理,得到视差图;在获得了视差信息后,根据投影模型很容易地可以得到原始图像的深度信息和三维信息...

空空如也

空空如也

1 2 3 4 5 ... 20
收藏数 458
精华内容 183
关键字:

ncc立体匹配