OpenCV学习(21) Grabcut算法详解

grab cut算法是graph cut算法的改进。在理解grab cut算之前，应该学习一下graph cut算法的概念及实现方式。

我搜集了一些graph cut资料：http://yunpan.cn/QGDVdBXwkXutH

     grab cut算法详细描述见资料中的pdf文件：“GrabCut” — Interactive Foreground Extraction using Iterated Graph Cuts

     grab cut算法是一种基于图论的图像分割方法，首先要定义一个Gibbs能量函数，然后求解这个函数的min-cut，这个min-cut就是前景背景的分割像素集合。

1. 能量函数的定义 

      在grab cut算法中，能量函数定义为：

      其中U函数部分表示能量函数的区域数据项，V函数表示能量函数的光滑项(边界项)。

     我们使用混合多高斯模型D(x)表示某个像素属于前景或背景的概率，这里K=5，表示第i个单高斯函数对概率贡献的权重系数，所以有， 为第i单高斯函数。为第i个单高斯函数的均值，为第i个单高斯函数的协方差。

单高斯函数g公式为：

 

区域数据项U函数为：

这里n表示图像中的第n个像素。【注：这儿的，是对前面的D(x)取负对数】：

grab cut算法的输入图像是RGB 3通道的图像,对于输入图像，我们用两个混合多高斯模型来分别表示前景和背景。

 

光滑项V函数为：

C是相邻颜色对的集合，是一个常量值50，

 

这些公式直接理解有点困难，下面我们结合程序看看grabcut算法中，如何计算能量公式，以及如何分段。

1、首先定义混合多高斯模型

class GMM
{
public:
    static const int componentsCount = 5; //对应K=5

    GMM( Mat& _model );
    double operator()( const Vec3d color ) const;
    double operator()( int ci, const Vec3d color ) const;
    int whichComponent( const Vec3d color ) const;

    void initLearning();
    void addSample( int ci, const Vec3d color );
    void endLearning();

private:
    void calcInverseCovAndDeterm( int ci );
    Mat model;
    double* coefs; //权重系数
    double* mean; //均值
    double* cov;  //协方差

    double inverseCovs[componentsCount][3][3];
    double covDeterms[componentsCount];

    double sums[componentsCount][3]; //所有样本bgr三个颜色分量的和，用来计算权重系数
    double prods[componentsCount][3][3]; //所有样本bgr颜色的行列式值，用来计算协方差
    int sampleCounts[componentsCount]; //每个单高斯函数的样本数
    int totalSampleCount;  //样本总数

};

     从上面混合多高斯公式可以知道，只要确定了三个参数：权重系数、均值、协方差，就可以根据当前像素点的bgr值确定当前像素属于前景和背景的概率D(x)，所以在GMM类中，我们定义三个指针，分别表示权重系数，均值和协方差。因为当前像素用bgr值表示，所以均值其实为3个double数，再加上K=5(5个单高斯函数组成的多高斯混合函数)，总共15双精度值，而权重系数则为5个双精度值，cov公共3*3*5=45个双精度值。

double* coefs; //权重系数
double* mean; //均值
double* cov; //协方差
在GMM的构造函数中，我们会创建一个1维的矩阵，总共65个双精度数，权重系数指向矩阵数据头，均值指向第6个双精度数，协方差指向第21个双精度数。

    _model.create( 1, modelSize*componentsCount, CV_64FC1 );
    _model.setTo(Scalar(0));

coefs = model.ptr<double>(0);
mean = coefs + componentsCount;
cov = mean + 3*componentsCount;

2.初始化GMM变量

我们定义了两个变量

    GMM bgdGMM( bgdModel ), fgdGMM( fgdModel );

分别表示和前景的混合多高斯模型。

      首先我们根据选定的四边形框来初始化mask图像，四边形框外的像素是背景，值为GC_BGD ，四边形内的像素可能是前景，值为GC_PR_FGD。

/*
  Initialize mask using rectangular.
  设置mask的初始值，四边形框内圈定的像素值为GC_PR_FGD
*/
void gGrabCut::initMaskWithRect( Mat& mask, Size imgSize, Rect rect )
{
    mask.create( imgSize, CV_8UC1 );
    mask.setTo( GC_BGD );

    rect.x = max(0, rect.x);
    rect.y = max(0, rect.y);
    rect.width = min(rect.width, imgSize.width-rect.x);
    rect.height = min(rect.height, imgSize.height-rect.y);

    (mask(rect)).setTo( Scalar(GC_PR_FGD) );
}

      之后，我们会根据mask图像，读入样本数据。前景GMM的样本数据放在变量fgdSamples中，背景GMM的样本数据放入变量bgdSamples中。fgdSamples和bgdSamples中存放得都是一些bgr颜色值。

Mat bgdLabels, fgdLabels;
vector<Vec3f> bgdSamples, fgdSamples;
Point p;
for( p.y = 0; p.y < img.rows; p.y++ )
{
    for( p.x = 0; p.x < img.cols; p.x++ )
    {
        if( mask.at<uchar>(p) == GC_BGD || mask.at<uchar>(p) == GC_PR_BGD )
            bgdSamples.push_back( (Vec3f)img.at<Vec3b>(p) );
        else // GC_FGD | GC_PR_FGD
            fgdSamples.push_back( (Vec3f)img.at<Vec3b>(p) );
    }
}
     之后我们会根据kmeans聚类算法，计算得到当前像素属于前景或背景混合多高斯变量中的第几个单高斯函数，结果放在bgdSamples, fgdSamples中，值为0-4。

    //一行是一个数据样本，3列是b，g，r三个属性
    Mat _bgdSamples( (int)bgdSamples.size(), 3, CV_32FC1, &bgdSamples[0][0] );
    //GMM::componentsCount聚类的个数
    //KMEANS_PP_CENTERS是采用Arthur & Vassilvitskii (2007) k-means++: The Advantages of Careful Seeding获取初始化种子点
    kmeans( _bgdSamples, GMM::componentsCount, bgdLabels,
            TermCriteria( CV_TERMCRIT_ITER, kMeansItCount, 0.0), 0, kMeansType );
    Mat _fgdSamples( (int)fgdSamples.size(), 3, CV_32FC1, &fgdSamples[0][0] );
    kmeans( _fgdSamples, GMM::componentsCount, fgdLabels,
            TermCriteria( CV_TERMCRIT_ITER, kMeansItCount, 0.0), 0, kMeansType );

下面我们计算得到均值、协方差和权重系数：

      权重系数为属于某个单高斯函数的采样像素数量除以所有采样像素的数量。每个单高斯函数的均值为所有属于该函数的采样像素颜色和除以属于该函数的颜色采样数量。协方差的公式比较复杂，大家看看下面代码中c[0]-c[9]的计算就ok了。注意其中的函数calcInverseCovAndDeterm用来计算协方差矩阵的行列式值以及逆矩阵，这些值在计算能量公式的数据项函数U时候使用。

//计算得到均值、协方差以及权重系数
void GMM::endLearning()
{
    const double variance = 0.01;
    for( int ci = 0; ci < componentsCount; ci++ )
    {
        int n = sampleCounts[ci];
        if( n == 0 )
            coefs[ci] = 0;
        else
        {
            coefs[ci] = (double)n/totalSampleCount;

            double* m = mean + 3*ci;
            m[0] = sums[ci][0]/n; m[1] = sums[ci][1]/n; m[2] = sums[ci][2]/n;

            double* c = cov + 9*ci;
            c[0] = prods[ci][0][0]/n - m[0]*m[0]; c[1] = prods[ci][0][1]/n - m[0]*m[1]; c[2] = prods[ci][0][2]/n - m[0]*m[2];
            c[3] = prods[ci][1][0]/n - m[1]*m[0]; c[4] = prods[ci][1][1]/n - m[1]*m[1]; c[5] = prods[ci][1][2]/n - m[1]*m[2];
            c[6] = prods[ci][2][0]/n - m[2]*m[0]; c[7] = prods[ci][2][1]/n - m[2]*m[1]; c[8] = prods[ci][2][2]/n - m[2]*m[2];

            double dtrm = c[0]*(c[4]*c[8]-c[5]*c[7]) - c[1]*(c[3]*c[8]-c[5]*c[6]) + c[2]*(c[3]*c[7]-c[4]*c[6]);
            if( dtrm <= std::numeric_limits<double>::epsilon() )
            {
               // Adds the white noise to avoid singular covariance matrix.
                c[0] += variance;
                c[4] += variance;
                c[8] += variance;
            }

            calcInverseCovAndDeterm(ci);
        }
    }
}

下面我们开始计算能量公式中光滑性函数V：

（注：=1)，

const double gamma = 50;
const double lambda = 9*gamma;
const double beta = calcBeta( img );

Mat leftW, upleftW, upW, uprightW;
calcNWeights( img, leftW, upleftW, upW, uprightW, beta, gamma );

     我们会在beta函数计算公式中的beta值，其中4*img.cols*img.rows - 3*img.cols - 3*img.rows + 2为邻接距离的数量。

/*
  计算光滑性函数中的beta值
  Calculate beta - parameter of GrabCut algorithm.
  beta = 1/(2*avg(sqr(||color[i] - color[j]||)))
*/
static double calcBeta( const Mat& img )
{
    double beta = 0;
    for( int y = 0; y < img.rows; y++ )
    {
        for( int x = 0; x < img.cols; x++ )
        {
            Vec3d color = img.at<Vec3b>(y,x);
            if( x>0 ) // left
            {
                Vec3d diff = color - (Vec3d)img.at<Vec3b>(y,x-1);
                beta += diff.dot(diff);
            }
            if( y>0 && x>0 ) // upleft
            {
                Vec3d diff = color - (Vec3d)img.at<Vec3b>(y-1,x-1);
                beta += diff.dot(diff);
            }
            if( y>0 ) // up
            {
                Vec3d diff = color - (Vec3d)img.at<Vec3b>(y-1,x);
                beta += diff.dot(diff);
            }
            if( y>0 && x<img.cols-1) // upright
            {
                Vec3d diff = color - (Vec3d)img.at<Vec3b>(y-1,x+1);
                beta += diff.dot(diff);
            }
        }
    }
    if( beta <= std::numeric_limits<double>::epsilon() )
        beta = 0;
    else //除以邻接距离的数量
        beta = 1.f / (2 * beta/(4*img.cols*img.rows - 3*img.cols - 3*img.rows + 2) );

    return beta;
}

      我们通过caclNWeights函数计算非终端顶点的权重值，计算公式依据V函数，权重结果放在四个矩阵leftW, upleftW, upW, uprightW中，最后，我们根据像素和权重值构建图，并用max-flow算法解得min-cut，求解的结果放在mask图像中，前景部分的值为GC_PR_FGD，背景部分的值为GC_PR_BGD

 

Mat leftW, upleftW, upW, uprightW;
calcNWeights( img, leftW, upleftW, upW, uprightW, beta, gamma );

for( int i = 0; i < iterCount; i++ )
{
    GCGraph<double> graph;
    assignGMMsComponents( img, mask, bgdGMM, fgdGMM, compIdxs );
    learnGMMs( img, mask, compIdxs, bgdGMM, fgdGMM );
    constructGCGraph(img, mask, bgdGMM, fgdGMM, lambda, leftW, upleftW, upW, uprightW, graph );
    estimateSegmentation( graph, mask );
}

时间： 2024-09-20 09:36:10

OpenCV学习(21) Grabcut算法详解的相关文章

Android 布局学习之——Layout(布局)详解二（常见布局和布局参数）

  [Android布局学习系列]   1.Android 布局学习之--Layout(布局)详解一   2.Android 布局学习之--Layout(布局)详解二(常见布局和布局参数)   3.Android 布局学习之--LinearLayout的layout_weight属性   4.Android 布局学习之--LinearLayout属性baselineAligned的作用及baseline      Layout Parameters(布局参数):            在XML文

机器学习（二）--- 分类算法详解

感觉狼厂有些把机器学习和数据挖掘神话了,机器学习.数据挖掘的能力其实是有边界的.机器学习.数据挖掘永远是给大公司的业务锦上添花的东西,它可以帮助公司赚更多的钱,却不能帮助公司在与其他公司的竞争中取得领先优势,所以小公司招聘数据挖掘/机器学习不是为了装逼就是在自寻死路.可是相比Java和C++语言开发来说,机器学习/数据挖掘确实是新一些老人占的坑少一些,而且可以经常接触一些新的东西.还是赶紧再次抓住机会集中的再总结一下吧,不能再拖拖拉拉了.  其实数据挖掘的主要任务是分类.聚类.关联分析.预测.时

Floyd求最短路径算法详解

倘若我们要在计算机上建立一个交通咨询系统则可以采用图的结构来表示实际的交通网络.其实现最基本的功能,求出任意两点间的最短路径, 求最短路径的经典方法有很多种,最常用的便是迪杰斯特拉算法和佛洛依德(Floyd)算法,这篇文章就着重介绍Floyd算法. 求两点之间的最短路径无外乎有两种情况,一种就是从一点直接到另一点,另一种就是从一点经过n个节点后再到另一个节点,比如说要从A到B,则有两种情况就是A直接到B,或者是从A经过N个节点后再到B,所以,我们假设Dis(AB)为节点A到节点B的最短路径的距离

javascript快速排序算法详解_javascript技巧

"快速排序"的思想很简单,整个排序过程只需要三步: (1)在数据集之中,找一个基准点 (2)建立两个数组,分别存储左边和右边的数组 (3)利用递归进行下次比较 看一个demo:http://jsdo.it/norahiko/oxIy/fullscreen(网页打开可能较慢,慢慢等待吧) <script type="text/javascript"> function quickSort(arr){ if(arr.length<=1){ return

Javascript实现快速排序（Quicksort）的算法详解_javascript技巧

目前,最常见的排序算法大概有七八种,其中"快速排序"(Quicksort)使用得最广泛,速度也较快.它是图灵奖得主C. A. R. Hoare(1934--)于1960时提出来的. "快速排序"的思想很简单,整个排序过程只需要三步: (1)在数据集之中,选择一个元素作为"基准"(pivot). (2)所有小于"基准"的元素,都移到"基准"的左边:所有大于"基准"的元素,都移到"

OpenCV学习(20) grabcut分割算法

      在OpenCV中,实现了grabcut分割算法,该算法可以方便的分割出前景图像,操作简单,而且分割的效果很好.算法的原理参见papaer:"GrabCut" - Interactive Foreground Extraction using Iterated Graph Cuts 比如下面的一副图,我们只要选定一个四边形框,把框中的图像作为grabcut的一个输入参数,表示该框中的像素可能属于前景,但框外的部分一定属于背景. 然后调用grabcut函数,就可以分割出城堡来.

SSH框架网上商城项目第21战之详解易宝支付的流程_java

这一节我们先写一个简单点的Demo来测试易宝支付的流程,熟悉这个流程后,再做实际的开发,因为是一个Demo,所以我没有考虑一些设计模式的东西,就是直接实现支付功能.实现支付功能需要易宝给我们提供的API.那么问题来了,使用第三方支付平台最主要的一件事就是获取该平台的API,我们首先得获取他们的API以及开发文档,然后才可以做进一步的开发. 1. 获取易宝的API 获取API的第一步,要在易宝上注册一个账号,这个账号是商家的账号,后面买家付款后,会将钱款存入该账号中,然后商家自己提取到银行卡,易宝

Php学习之Apache服务器详解

Php学习之服务器--Apache服务器详解 Iis服务器:主要是服务于微软,基于运行Microsoft windows的互联网基本服务 Lighttpd服务器:开源软件,针对高性能,底内存开销,cpu占用底,成熟度还不够 Apache服务器的介绍 1. 抓包软件:httpwatch.rar,了解发送和接受数据包 2. Apache服务器的安装 1. 查看windows中已经安装的服务,确定原先没有安装apache 2. 安装apache服务器: 3. 测试是否真的成功.在浏览器中输入http:

php 二分查找法算法详解

一.概念:二分查找又称折半查找,优点是比较次数少,查找速度快,平均性能好:其缺点是要求待查表为有序表,且插入删除困难.因此,折半查找方法适用于不经常变动而查找频繁的有序列表.首先,假设表中元素是按升序排列,将表中间位置记录的关键字与查找关键字比较,如果两者相等,则查找成功:否则利用中间位置记录将表分成前.后两个子表,如果中间位置记录的关键字大于查找关键字,则进一步查找前一子表,否则进一步查找后一子表.重复以上过程,直到找到满足条件的记录,使查找成功,或直到子表不存在为止,此时查找不成功. 二.代