grabcut.cpp

00001 /*M///////////////////////////////////////////////////////////////////////////////////////
00002 //
00003 //  IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
00004 //
00005 //  By downloading, copying, installing or using the software you agree to this license.
00006 //  If you do not agree to this license, do not download, install,
00007 //  copy or use the software.
00008 //
00009 //
00010 //                        Intel License Agreement
00011 //                For Open Source Computer Vision Library
00012 //
00013 // Copyright (C) 2000, Intel Corporation, all rights reserved.
00014 // Third party copyrights are property of their respective owners.
00015 //
00016 // Redistribution and use in source and binary forms, with or without modification,
00017 // are permitted provided that the following conditions are met:
00018 //
00019 //   * Redistribution's of source code must retain the above copyright notice,
00020 //     this list of conditions and the following disclaimer.
00021 //
00022 //   * Redistribution's in binary form must reproduce the above copyright notice,
00023 //     this list of conditions and the following disclaimer in the documentation
00024 //     and/or other materials provided with the distribution.
00025 //
00026 //   * The name of Intel Corporation may not be used to endorse or promote products
00027 //     derived from this software without specific prior written permission.
00028 //
00029 // This software is provided by the copyright holders and contributors "as is" and
00030 // any express or implied warranties, including, but not limited to, the implied
00031 // warranties of merchantability and fitness for a particular purpose are disclaimed.
00032 // In no event shall the Intel Corporation or contributors be liable for any direct,
00033 // indirect, incidental, special, exemplary, or consequential damages
00034 // (including, but not limited to, procurement of substitute goods or services;
00035 // loss of use, data, or profits; or business interruption) however caused
00036 // and on any theory of liability, whether in contract, strict liability,
00037 // or tort (including negligence or otherwise) arising in any way out of
00038 // the use of this software, even if advised of the possibility of such damage.
00039 //
00040 //M*/
00041 
00042 #include "precomp.hpp"
00043 #include "gcgraph.hpp"
00044 #include <limits>
00045 
00046 using namespace cv;
00047 
00048 /*
00049 This is implementation of image segmentation algorithm GrabCut described in
00050 "GrabCut — Interactive Foreground Extraction using Iterated Graph Cuts".
00051 Carsten Rother, Vladimir Kolmogorov, Andrew Blake.
00052  */
00053 
00054 /*
00055  GMM - Gaussian Mixture Model
00056 */
00057 class GMM
00058 {
00059 public:
00060     static const int componentsCount = 5;
00061 
00062     GMM( Mat& _model );
00063     double operator()( const Vec3d color ) const;
00064     double operator()( int ci, const Vec3d color ) const;
00065     int whichComponent( const Vec3d color ) const;
00066 
00067     void initLearning();
00068     void addSample( int ci, const Vec3d color );
00069     void endLearning();
00070 
00071 private:
00072     void calcInverseCovAndDeterm( int ci );
00073     Mat model;
00074     double* coefs;
00075     double* mean;
00076     double* cov;
00077 
00078     double inverseCovs[componentsCount][3][3];
00079     double covDeterms[componentsCount];
00080 
00081     double sums[componentsCount][3];
00082     double prods[componentsCount][3][3];
00083     int sampleCounts[componentsCount];
00084     int totalSampleCount;
00085 };
00086 
00087 GMM::GMM( Mat& _model )
00088 {
00089     const int modelSize = 3/*mean*/ + 9/*covariance*/ + 1/*component weight*/;
00090     if( _model.empty() )
00091     {
00092         _model.create( 1, modelSize*componentsCount, CV_64FC1 );
00093         _model.setTo(Scalar (0));
00094     }
00095     else if( (_model.type() != CV_64FC1) || (_model.rows != 1) || (_model.cols != modelSize*componentsCount) )
00096         CV_Error( CV_StsBadArg, "_model must have CV_64FC1 type, rows == 1 and cols == 13*componentsCount" );
00097 
00098     model = _model;
00099 
00100     coefs = model.ptr<double>(0);
00101     mean = coefs + componentsCount;
00102     cov = mean + 3*componentsCount;
00103 
00104     for( int ci = 0; ci < componentsCount; ci++ )
00105         if( coefs[ci] > 0 )
00106              calcInverseCovAndDeterm( ci );
00107 }
00108 
00109 double GMM::operator()( const Vec3d color ) const
00110 {
00111     double res = 0;
00112     for( int ci = 0; ci < componentsCount; ci++ )
00113         res += coefs[ci] * (*this)(ci, color );
00114     return res;
00115 }
00116 
00117 double GMM::operator()( int ci, const Vec3d color ) const
00118 {
00119     double res = 0;
00120     if( coefs[ci] > 0 )
00121     {
00122         CV_Assert( covDeterms[ci] > std::numeric_limits<double>::epsilon() );
00123         Vec3d diff = color;
00124         double* m = mean + 3*ci;
00125         diff[0] -= m[0]; diff[1] -= m[1]; diff[2] -= m[2];
00126         double mult = diff[0]*(diff[0]*inverseCovs[ci][0][0] + diff[1]*inverseCovs[ci][1][0] + diff[2]*inverseCovs[ci][2][0])
00127                    + diff[1]*(diff[0]*inverseCovs[ci][0][1] + diff[1]*inverseCovs[ci][1][1] + diff[2]*inverseCovs[ci][2][1])
00128                    + diff[2]*(diff[0]*inverseCovs[ci][0][2] + diff[1]*inverseCovs[ci][1][2] + diff[2]*inverseCovs[ci][2][2]);
00129         res = 1.0f/sqrt(covDeterms[ci]) * exp(-0.5f*mult);
00130     }
00131     return res;
00132 }
00133 
00134 int GMM::whichComponent( const Vec3d color ) const
00135 {
00136     int k = 0;
00137     double max = 0;
00138 
00139     for( int ci = 0; ci < componentsCount; ci++ )
00140     {
00141         double p = (*this)( ci, color );
00142         if( p > max )
00143         {
00144             k = ci;
00145             max = p;
00146         }
00147     }
00148     return k;
00149 }
00150 
00151 void GMM::initLearning()
00152 {
00153     for( int ci = 0; ci < componentsCount; ci++)
00154     {
00155         sums[ci][0] = sums[ci][1] = sums[ci][2] = 0;
00156         prods[ci][0][0] = prods[ci][0][1] = prods[ci][0][2] = 0;
00157         prods[ci][1][0] = prods[ci][1][1] = prods[ci][1][2] = 0;
00158         prods[ci][2][0] = prods[ci][2][1] = prods[ci][2][2] = 0;
00159         sampleCounts[ci] = 0;
00160     }
00161     totalSampleCount = 0;
00162 }
00163 
00164 void GMM::addSample( int ci, const Vec3d color )
00165 {
00166     sums[ci][0] += color[0]; sums[ci][1] += color[1]; sums[ci][2] += color[2];
00167     prods[ci][0][0] += color[0]*color[0]; prods[ci][0][1] += color[0]*color[1]; prods[ci][0][2] += color[0]*color[2];
00168     prods[ci][1][0] += color[1]*color[0]; prods[ci][1][1] += color[1]*color[1]; prods[ci][1][2] += color[1]*color[2];
00169     prods[ci][2][0] += color[2]*color[0]; prods[ci][2][1] += color[2]*color[1]; prods[ci][2][2] += color[2]*color[2];
00170     sampleCounts[ci]++;
00171     totalSampleCount++;
00172 }
00173 
00174 void GMM::endLearning()
00175 {
00176     const double variance = 0.01;
00177     for( int ci = 0; ci < componentsCount; ci++ )
00178     {
00179         int n = sampleCounts[ci];
00180         if( n == 0 )
00181             coefs[ci] = 0;
00182         else
00183         {
00184             coefs[ci] = (double)n/totalSampleCount;
00185 
00186             double* m = mean + 3*ci;
00187             m[0] = sums[ci][0]/n; m[1] = sums[ci][1]/n; m[2] = sums[ci][2]/n;
00188 
00189             double* c = cov + 9*ci;
00190             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];
00191             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];
00192             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];
00193 
00194             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]);
00195             if( dtrm <= std::numeric_limits<double>::epsilon() )
00196             {
00197                 // Adds the white noise to avoid singular covariance matrix.
00198                 c[0] += variance;
00199                 c[4] += variance;
00200                 c[8] += variance;
00201             }
00202 
00203             calcInverseCovAndDeterm(ci);
00204         }
00205     }
00206 }
00207 
00208 void GMM::calcInverseCovAndDeterm( int ci )
00209 {
00210     if( coefs[ci] > 0 )
00211     {
00212         double *c = cov + 9*ci;
00213         double dtrm =
00214               covDeterms[ci] = 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]);
00215 
00216         CV_Assert( dtrm > std::numeric_limits<double>::epsilon() );
00217         inverseCovs[ci][0][0] =  (c[4]*c[8] - c[5]*c[7]) / dtrm;
00218         inverseCovs[ci][1][0] = -(c[3]*c[8] - c[5]*c[6]) / dtrm;
00219         inverseCovs[ci][2][0] =  (c[3]*c[7] - c[4]*c[6]) / dtrm;
00220         inverseCovs[ci][0][1] = -(c[1]*c[8] - c[2]*c[7]) / dtrm;
00221         inverseCovs[ci][1][1] =  (c[0]*c[8] - c[2]*c[6]) / dtrm;
00222         inverseCovs[ci][2][1] = -(c[0]*c[7] - c[1]*c[6]) / dtrm;
00223         inverseCovs[ci][0][2] =  (c[1]*c[5] - c[2]*c[4]) / dtrm;
00224         inverseCovs[ci][1][2] = -(c[0]*c[5] - c[2]*c[3]) / dtrm;
00225         inverseCovs[ci][2][2] =  (c[0]*c[4] - c[1]*c[3]) / dtrm;
00226     }
00227 }
00228 
00229 /*
00230   Calculate beta - parameter of GrabCut algorithm.
00231   beta = 1/(2*avg(sqr(||color[i] - color[j]||)))
00232 */
00233 static double calcBeta( const Mat& img )
00234 {
00235     double beta = 0;
00236     for( int y = 0; y < img.rows; y++ )
00237     {
00238         for( int x = 0; x < img.cols; x++ )
00239         {
00240             Vec3d color = img.at<Vec3b>(y,x);
00241             if( x>0 ) // left
00242             {
00243                 Vec3d diff = color - (Vec3d)img.at<Vec3b>(y,x-1);
00244                 beta += diff.dot(diff);
00245             }
00246             if( y>0 && x>0 ) // upleft
00247             {
00248                 Vec3d diff = color - (Vec3d)img.at<Vec3b>(y-1,x-1);
00249                 beta += diff.dot(diff);
00250             }
00251             if( y>0 ) // up
00252             {
00253                 Vec3d diff = color - (Vec3d)img.at<Vec3b>(y-1,x);
00254                 beta += diff.dot(diff);
00255             }
00256             if( y>0 && x<img.cols-1) // upright
00257             {
00258                 Vec3d diff = color - (Vec3d)img.at<Vec3b>(y-1,x+1);
00259                 beta += diff.dot(diff);
00260             }
00261         }
00262     }
00263     if( beta <= std::numeric_limits<double>::epsilon() )
00264         beta = 0;
00265     else
00266         beta = 1.f / (2 * beta/(4*img.cols*img.rows - 3*img.cols - 3*img.rows + 2) );
00267 
00268     return beta;
00269 }
00270 
00271 /*
00272   Calculate weights of noterminal vertices of graph.
00273   beta and gamma - parameters of GrabCut algorithm.
00274  */
00275 static void calcNWeights( const Mat& img, Mat& leftW, Mat& upleftW, Mat& upW, Mat& uprightW, double beta, double gamma )
00276 {
00277     const double gammaDivSqrt2 = gamma / std::sqrt(2.0f);
00278     leftW.create( img.rows, img.cols, CV_64FC1 );
00279     upleftW.create( img.rows, img.cols, CV_64FC1 );
00280     upW.create( img.rows, img.cols, CV_64FC1 );
00281     uprightW.create( img.rows, img.cols, CV_64FC1 );
00282     for( int y = 0; y < img.rows; y++ )
00283     {
00284         for( int x = 0; x < img.cols; x++ )
00285         {
00286             Vec3d color = img.at<Vec3b>(y,x);
00287             if( x-1>=0 ) // left
00288             {
00289                 Vec3d diff = color - (Vec3d)img.at<Vec3b>(y,x-1);
00290                 leftW.at<double>(y,x) = gamma * exp(-beta*diff.dot(diff));
00291             }
00292             else
00293                 leftW.at<double>(y,x) = 0;
00294             if( x-1>=0 && y-1>=0 ) // upleft
00295             {
00296                 Vec3d diff = color - (Vec3d)img.at<Vec3b>(y-1,x-1);
00297                 upleftW.at<double>(y,x) = gammaDivSqrt2 * exp(-beta*diff.dot(diff));
00298             }
00299             else
00300                 upleftW.at<double>(y,x) = 0;
00301             if( y-1>=0 ) // up
00302             {
00303                 Vec3d diff = color - (Vec3d)img.at<Vec3b>(y-1,x);
00304                 upW.at<double>(y,x) = gamma * exp(-beta*diff.dot(diff));
00305             }
00306             else
00307                 upW.at<double>(y,x) = 0;
00308             if( x+1<img.cols && y-1>=0 ) // upright
00309             {
00310                 Vec3d diff = color - (Vec3d)img.at<Vec3b>(y-1,x+1);
00311                 uprightW.at<double>(y,x) = gammaDivSqrt2 * exp(-beta*diff.dot(diff));
00312             }
00313             else
00314                 uprightW.at<double>(y,x) = 0;
00315         }
00316     }
00317 }
00318 
00319 /*
00320   Check size, type and element values of mask matrix.
00321  */
00322 static void checkMask( const Mat& img, const Mat& mask )
00323 {
00324     if( mask.empty() )
00325         CV_Error( CV_StsBadArg, "mask is empty" );
00326     if( mask.type() != CV_8UC1 )
00327         CV_Error( CV_StsBadArg, "mask must have CV_8UC1 type" );
00328     if( mask.cols != img.cols || mask.rows != img.rows )
00329         CV_Error( CV_StsBadArg, "mask must have as many rows and cols as img" );
00330     for( int y = 0; y < mask.rows; y++ )
00331     {
00332         for( int x = 0; x < mask.cols; x++ )
00333         {
00334             uchar val = mask.at<uchar>(y,x);
00335             if( val!=GC_BGD && val!=GC_FGD && val!=GC_PR_BGD && val!=GC_PR_FGD )
00336                 CV_Error( CV_StsBadArg, "mask element value must be equel"
00337                     "GC_BGD or GC_FGD or GC_PR_BGD or GC_PR_FGD" );
00338         }
00339     }
00340 }
00341 
00342 /*
00343   Initialize mask using rectangular.
00344 */
00345 static void initMaskWithRect( Mat& mask, Size imgSize, Rect rect )
00346 {
00347     mask.create( imgSize, CV_8UC1 );
00348     mask.setTo( GC_BGD );
00349 
00350     rect.x = std::max(0, rect.x);
00351     rect.y = std::max(0, rect.y);
00352     rect.width = std::min(rect.width, imgSize.width-rect.x);
00353     rect.height = std::min(rect.height, imgSize.height-rect.y);
00354 
00355     (mask(rect)).setTo( Scalar (GC_PR_FGD) );
00356 }
00357 
00358 /*
00359   Initialize GMM background and foreground models using kmeans algorithm.
00360 */
00361 static void initGMMs( const Mat& img, const Mat& mask, GMM& bgdGMM, GMM& fgdGMM )
00362 {
00363     const int kMeansItCount = 10;
00364     const int kMeansType = KMEANS_PP_CENTERS;
00365 
00366     Mat bgdLabels, fgdLabels;
00367     std::vector<Vec3f> bgdSamples, fgdSamples;
00368     Point  p;
00369     for( p.y = 0; p.y < img.rows; p.y++ )
00370     {
00371         for( p.x = 0; p.x < img.cols; p.x++ )
00372         {
00373             if( mask.at<uchar>(p) == GC_BGD || mask.at<uchar>(p) == GC_PR_BGD )
00374                 bgdSamples.push_back( (Vec3f)img.at<Vec3b>(p) );
00375             else // GC_FGD | GC_PR_FGD
00376                 fgdSamples.push_back( (Vec3f)img.at<Vec3b>(p) );
00377         }
00378     }
00379     CV_Assert( !bgdSamples.empty() && !fgdSamples.empty() );
00380     Mat _bgdSamples( (int)bgdSamples.size(), 3, CV_32FC1, &bgdSamples[0][0] );
00381     kmeans( _bgdSamples, GMM::componentsCount, bgdLabels,
00382             TermCriteria( CV_TERMCRIT_ITER, kMeansItCount, 0.0), 0, kMeansType );
00383     Mat _fgdSamples( (int)fgdSamples.size(), 3, CV_32FC1, &fgdSamples[0][0] );
00384     kmeans( _fgdSamples, GMM::componentsCount, fgdLabels,
00385             TermCriteria( CV_TERMCRIT_ITER, kMeansItCount, 0.0), 0, kMeansType );
00386 
00387     bgdGMM.initLearning();
00388     for( int i = 0; i < (int)bgdSamples.size(); i++ )
00389         bgdGMM.addSample( bgdLabels.at<int>(i,0), bgdSamples[i] );
00390     bgdGMM.endLearning();
00391 
00392     fgdGMM.initLearning();
00393     for( int i = 0; i < (int)fgdSamples.size(); i++ )
00394         fgdGMM.addSample( fgdLabels.at<int>(i,0), fgdSamples[i] );
00395     fgdGMM.endLearning();
00396 }
00397 
00398 /*
00399   Assign GMMs components for each pixel.
00400 */
00401 static void assignGMMsComponents( const Mat& img, const Mat& mask, const GMM& bgdGMM, const GMM& fgdGMM, Mat& compIdxs )
00402 {
00403     Point  p;
00404     for( p.y = 0; p.y < img.rows; p.y++ )
00405     {
00406         for( p.x = 0; p.x < img.cols; p.x++ )
00407         {
00408             Vec3d color = img.at<Vec3b>(p);
00409             compIdxs.at<int>(p) = mask.at<uchar>(p) == GC_BGD || mask.at<uchar>(p) == GC_PR_BGD ?
00410                 bgdGMM.whichComponent(color) : fgdGMM.whichComponent(color);
00411         }
00412     }
00413 }
00414 
00415 /*
00416   Learn GMMs parameters.
00417 */
00418 static void learnGMMs( const Mat& img, const Mat& mask, const Mat& compIdxs, GMM& bgdGMM, GMM& fgdGMM )
00419 {
00420     bgdGMM.initLearning();
00421     fgdGMM.initLearning();
00422     Point  p;
00423     for( int ci = 0; ci < GMM::componentsCount; ci++ )
00424     {
00425         for( p.y = 0; p.y < img.rows; p.y++ )
00426         {
00427             for( p.x = 0; p.x < img.cols; p.x++ )
00428             {
00429                 if( compIdxs.at<int>(p) == ci )
00430                 {
00431                     if( mask.at<uchar>(p) == GC_BGD || mask.at<uchar>(p) == GC_PR_BGD )
00432                         bgdGMM.addSample( ci, img.at<Vec3b>(p) );
00433                     else
00434                         fgdGMM.addSample( ci, img.at<Vec3b>(p) );
00435                 }
00436             }
00437         }
00438     }
00439     bgdGMM.endLearning();
00440     fgdGMM.endLearning();
00441 }
00442 
00443 /*
00444   Construct GCGraph
00445 */
00446 static void constructGCGraph( const Mat& img, const Mat& mask, const GMM& bgdGMM, const GMM& fgdGMM, double lambda,
00447                        const Mat& leftW, const Mat& upleftW, const Mat& upW, const Mat& uprightW,
00448                        GCGraph<double>& graph )
00449 {
00450     int vtxCount = img.cols*img.rows,
00451         edgeCount = 2*(4*img.cols*img.rows - 3*(img.cols + img.rows) + 2);
00452     graph.create(vtxCount, edgeCount);
00453     Point  p;
00454     for( p.y = 0; p.y < img.rows; p.y++ )
00455     {
00456         for( p.x = 0; p.x < img.cols; p.x++)
00457         {
00458             // add node
00459             int vtxIdx = graph.addVtx();
00460             Vec3b color = img.at<Vec3b>(p);
00461 
00462             // set t-weights
00463             double fromSource, toSink;
00464             if( mask.at<uchar>(p) == GC_PR_BGD || mask.at<uchar>(p) == GC_PR_FGD )
00465             {
00466                 fromSource = -log( bgdGMM(color) );
00467                 toSink = -log( fgdGMM(color) );
00468             }
00469             else if( mask.at<uchar>(p) == GC_BGD )
00470             {
00471                 fromSource = 0;
00472                 toSink = lambda;
00473             }
00474             else // GC_FGD
00475             {
00476                 fromSource = lambda;
00477                 toSink = 0;
00478             }
00479             graph.addTermWeights( vtxIdx, fromSource, toSink );
00480 
00481             // set n-weights
00482             if( p.x>0 )
00483             {
00484                 double w = leftW.at<double>(p);
00485                 graph.addEdges( vtxIdx, vtxIdx-1, w, w );
00486             }
00487             if( p.x>0 && p.y>0 )
00488             {
00489                 double w = upleftW.at<double>(p);
00490                 graph.addEdges( vtxIdx, vtxIdx-img.cols-1, w, w );
00491             }
00492             if( p.y>0 )
00493             {
00494                 double w = upW.at<double>(p);
00495                 graph.addEdges( vtxIdx, vtxIdx-img.cols, w, w );
00496             }
00497             if( p.x<img.cols-1 && p.y>0 )
00498             {
00499                 double w = uprightW.at<double>(p);
00500                 graph.addEdges( vtxIdx, vtxIdx-img.cols+1, w, w );
00501             }
00502         }
00503     }
00504 }
00505 
00506 /*
00507   Estimate segmentation using MaxFlow algorithm
00508 */
00509 static void estimateSegmentation( GCGraph<double>& graph, Mat& mask )
00510 {
00511     graph.maxFlow();
00512     Point  p;
00513     for( p.y = 0; p.y < mask.rows; p.y++ )
00514     {
00515         for( p.x = 0; p.x < mask.cols; p.x++ )
00516         {
00517             if( mask.at<uchar>(p) == GC_PR_BGD || mask.at<uchar>(p) == GC_PR_FGD )
00518             {
00519                 if( graph.inSourceSegment( p.y*mask.cols+p.x /*vertex index*/ ) )
00520                     mask.at<uchar>(p) = GC_PR_FGD;
00521                 else
00522                     mask.at<uchar>(p) = GC_PR_BGD;
00523             }
00524         }
00525     }
00526 }
00527 
00528 void cv::grabCut( InputArray _img, InputOutputArray _mask, Rect rect,
00529                   InputOutputArray _bgdModel, InputOutputArray _fgdModel,
00530                   int iterCount, int mode )
00531 {
00532     Mat img = _img.getMat();
00533     Mat& mask = _mask.getMatRef();
00534     Mat& bgdModel = _bgdModel.getMatRef();
00535     Mat& fgdModel = _fgdModel.getMatRef();
00536 
00537     if( img.empty() )
00538         CV_Error( CV_StsBadArg, "image is empty" );
00539     if( img.type() != CV_8UC3 )
00540         CV_Error( CV_StsBadArg, "image must have CV_8UC3 type" );
00541 
00542     GMM bgdGMM( bgdModel ), fgdGMM( fgdModel );
00543     Mat compIdxs( img.size(), CV_32SC1 );
00544 
00545     if( mode == GC_INIT_WITH_RECT || mode == GC_INIT_WITH_MASK )
00546     {
00547         if( mode == GC_INIT_WITH_RECT )
00548             initMaskWithRect( mask, img.size(), rect );
00549         else // flag == GC_INIT_WITH_MASK
00550             checkMask( img, mask );
00551         initGMMs( img, mask, bgdGMM, fgdGMM );
00552     }
00553 
00554     if( iterCount <= 0)
00555         return;
00556 
00557     if( mode == GC_EVAL )
00558         checkMask( img, mask );
00559 
00560     const double gamma = 50;
00561     const double lambda = 9*gamma;
00562     const double beta = calcBeta( img );
00563 
00564     Mat leftW, upleftW, upW, uprightW;
00565     calcNWeights( img, leftW, upleftW, upW, uprightW, beta, gamma );
00566 
00567     for( int i = 0; i < iterCount; i++ )
00568     {
00569         GCGraph<double> graph;
00570         assignGMMsComponents( img, mask, bgdGMM, fgdGMM, compIdxs );
00571         learnGMMs( img, mask, compIdxs, bgdGMM, fgdGMM );
00572         constructGCGraph(img, mask, bgdGMM, fgdGMM, lambda, leftW, upleftW, upW, uprightW, graph );
00573         estimateSegmentation( graph, mask );
00574     }
00575 }
00576