使用OpenCV的一个MSER/MSCR实现

#ifndef GUARD_cvmser_h
#define GUARD_cvmser_h
#include "cv.h"
typedef struct CvMSERParams
{
// delta, in the code, it compares (size_{i}-size_{i-delta})/size_{i-delta}
int delta;
// prune the area which bigger/smaller than max_area/min_area
int max_area;
int min_area;
// prune the area have simliar size to its children
float max_variation;
// trace back to cut off mser with diversity < min_diversity
float min_diversity;
/* the next few params for MSER of color image */
// for color image, the evolution steps
int max_evolution;
// the area threshold to cause re-initialize
double area_threshold;
// ignore too small margin
double min_margin;
// the aperture size for edge blur
int edge_blur_size;
}
CvMSERParams;
CvMSERParams cvMSERParams( int delta = 5, int min_area = 60, int max_area = 14400, float max_variation = .25, float min_diversity = .2, int max_evolution = 200, double area_threshold = 1.01, double min_margin = .003, int edge_blur_size = 5 );
void cvExtractMSER( CvArr* _img, CvArr* _mask, CvSeq** contours, CvMemStorage* storage, CvMSERParams params );
#endif

/* Redistribution and use in source and binary forms, with or
* without modification, are permitted provided that the following
* conditions are met:
* 	Redistributions of source code must retain the above
* 	copyright notice, this list of conditions and the following
* 	disclaimer.
* 	Redistributions in binary form must reproduce the above
* 	copyright notice, this list of conditions and the following
* 	disclaimer in the documentation and/or other materials
* 	provided with the distribution.
* 	The name of Contributor may not be used to endorse or
* 	promote products derived from this software without
* 	specific prior written permission.
*
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND
* CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES,
* INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF
* MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
* DISCLAIMED. IN NO EVENT SHALL THE CONTRIBUTORS BE LIABLE FOR ANY
* DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
* CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
* PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA,
* OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
* THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR
* TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
* OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY
* OF SUCH DAMAGE.
* Copyright© 2009, Liu Liu All rights reserved.
*
* OpenCV functions for MSER extraction
*
* 1. there are two different implementation of MSER, one for grey image, one for color image
* 2. the grey image algorithm is taken from: Linear Time Maximally Stable Extremal Regions;
*    the paper claims to be faster than union-find method;
*    it actually get 1.5~2m/s on my centrino L7200 1.2GHz laptop.
* 3. the color image algorithm is taken from: Maximally Stable Colour Regions for Recognition and Match;
*    it should be much slower than grey image method ( 3~4 times );
*    the chi_table.h file is taken directly from paper's source code which is distributed under GPL.
* 4. though the name is *contours*, the result actually is a list of point set.
*/
#include "cvmser.h"
#include "cxmisc.h"
typedef struct CvLinkedPoint
{
struct CvLinkedPoint* prev;
struct CvLinkedPoint* next;
CvPoint pt;
}
CvLinkedPoint;
// the history of region grown
typedef struct CvMSERGrowHistory
{
struct CvMSERGrowHistory* shortcut;
struct CvMSERGrowHistory* child;
int stable; // when it ever stabled before, record the size
int val;
int size;
}
CvMSERGrowHistory;
typedef struct CvMSERConnectedComp
{
CvLinkedPoint* head;
CvLinkedPoint* tail;
CvMSERGrowHistory* history;
unsigned long grey_level;
int size;
int dvar; // the derivative of last var
float var; // the current variation (most time is the variation of one-step back)
}
CvMSERConnectedComp;
// Linear Time MSER claims by using bsf can get performance gain, here is the implementation
// however it seems that will not do any good in real world test
#if 0
#ifdef __GNUC__
inline unsigned char _BitScanForward(unsigned long * Index, unsigned long Mask)
{
unsigned int EFlags = 0;
__asm__ ( "bsf %[Mask], %[Index];"
"pushf;"
"pop %[EFlags];"
: [Index]"=r"(*Index), [EFlags]"=r"(EFlags)
: [Mask]"r"(Mask) );
return !(EFlags & 0x40);
}
#define __INTRIN_ENABLED__
#elif _MSC_BUILD
#include <intrin.h>
#define __INTRIN_ENABLED__
#endif
#ifdef __INTRIN_ENABLED__
inline void _bitset(unsigned long * a, unsigned long b)
{
*a |= 1<<b;
}
inline void _bitreset(unsigned long * a, unsigned long b)
{
*a &= ~(1<<b);
}
#endif
#endif
CvMSERParams cvMSERParams( int delta, int min_area, int max_area, float max_variation, float min_diversity, int max_evolution, double area_threshold, double min_margin, int edge_blur_size )
{
CvMSERParams params;
params.delta = delta;
params.min_area = min_area;
params.max_area = max_area;
params.max_variation = max_variation;
params.min_diversity = min_diversity;
params.max_evolution = max_evolution;
params.area_threshold = area_threshold;
params.min_margin = min_margin;
params.edge_blur_size = edge_blur_size;
return params;
}
// clear the connected component in stack
CV_INLINE static void
icvInitMSERComp( CvMSERConnectedComp* comp )
{
comp->size = 0;
comp->var = 0;
comp->dvar = 1;
comp->history = NULL;
}
// add history of size to a connected component
CV_INLINE static void
icvMSERNewHistory( CvMSERConnectedComp* comp,
CvMSERGrowHistory* history )
{
history->child = history;
if ( NULL == comp->history )
{
history->shortcut = history;
history->stable = 0;
} else {
comp->history->child = history;
history->shortcut = comp->history->shortcut;
history->stable = comp->history->stable;
}
history->val = comp->grey_level;
history->size = comp->size;
comp->history = history;
}
// merging two connected component
CV_INLINE static void
icvMSERMergeComp( CvMSERConnectedComp* comp1,
CvMSERConnectedComp* comp2,
CvMSERConnectedComp* comp,
CvMSERGrowHistory* history )
{
CvLinkedPoint* head;
CvLinkedPoint* tail;
comp->grey_level = comp2->grey_level;
history->child = history;
// select the winner by size
if ( comp1->size >= comp2->size )
{
if ( NULL == comp1->history )
{
history->shortcut = history;
history->stable = 0;
} else {
comp1->history->child = history;
history->shortcut = comp1->history->shortcut;
history->stable = comp1->history->stable;
}
if ( NULL != comp2->history && comp2->history->stable > history->stable )
history->stable = comp2->history->stable;
history->val = comp1->grey_level;
history->size = comp1->size;
// put comp1 to history
comp->var = comp1->var;
comp->dvar = comp1->dvar;
if ( comp1->size > 0 && comp2->size > 0 )
{
comp1->tail->next = comp2->head;
comp2->head->prev = comp1->tail;
}
head = ( comp1->size > 0 ) ? comp1->head : comp2->head;
tail = ( comp2->size > 0 ) ? comp2->tail : comp1->tail;
// always made the newly added in the last of the pixel list (comp1 ... comp2)
} else {
if ( NULL == comp2->history )
{
history->shortcut = history;
history->stable = 0;
} else {
comp2->history->child = history;
history->shortcut = comp2->history->shortcut;
history->stable = comp2->history->stable;
}
if ( NULL != comp1->history && comp1->history->stable > history->stable )
history->stable = comp1->history->stable;
history->val = comp2->grey_level;
history->size = comp2->size;
// put comp2 to history
comp->var = comp2->var;
comp->dvar = comp2->dvar;
if ( comp1->size > 0 && comp2->size > 0 )
{
comp2->tail->next = comp1->head;
comp1->head->prev = comp2->tail;
}
head = ( comp2->size > 0 ) ? comp2->head : comp1->head;
tail = ( comp1->size > 0 ) ? comp1->tail : comp2->tail;
// always made the newly added in the last of the pixel list (comp2 ... comp1)
}
comp->head = head;
comp->tail = tail;
comp->history = history;
comp->size = comp1->size + comp2->size;
}
CV_INLINE static float
icvMSERVariationCalc( CvMSERConnectedComp* comp,
int delta )
{
CvMSERGrowHistory* history = comp->history;
int val = comp->grey_level;
if ( NULL != history )
{
CvMSERGrowHistory* shortcut = history->shortcut;
while ( shortcut != shortcut->shortcut && shortcut->val + delta > val )
shortcut = shortcut->shortcut;
CvMSERGrowHistory* child = shortcut->child;
while ( child != child->child && child->val + delta <= val )
{
shortcut = child;
child = child->child;
}
// get the position of history where the shortcut->val <= delta+val and shortcut->child->val >= delta+val
history->shortcut = shortcut;
return (float)(comp->size-shortcut->size)/(float)shortcut->size;
// here is a small modification of MSER where cal ||R_{i}-R_{i-delta}||/||R_{i-delta}||
// in standard MSER, cal ||R_{i+delta}-R_{i-delta}||/||R_{i}||
// my calculation is simpler and much easier to implement
}
return 1.;
}
CV_INLINE static bool
icvMSERStableCheck( CvMSERConnectedComp* comp,
CvMSERParams params )
{
// tricky part: it actually check the stablity of one-step back
if ( comp->history == NULL || comp->history->size <= params.min_area || comp->history->size >= params.max_area )
return 0;
float div = (float)(comp->history->size-comp->history->stable)/(float)comp->history->size;
float var = icvMSERVariationCalc( comp, params.delta );
int dvar = ( comp->var < var || comp->history->val + 1 < comp->grey_level );
int stable = ( dvar && !comp->dvar && comp->var < params.max_variation && div > params.min_diversity );
comp->var = var;
comp->dvar = dvar;
if ( stable )
comp->history->stable = comp->history->size;
return stable;
}
// add a pixel to the pixel list
CV_INLINE static void
icvAccumulateMSERComp( CvMSERConnectedComp* comp,
CvLinkedPoint* point )
{
if ( comp->size > 0 )
{
point->prev = comp->tail;
comp->tail->next = point;
point->next = NULL;
} else {
point->prev = NULL;
point->next = NULL;
comp->head = point;
}
comp->tail = point;
comp->size++;
}
// convert the point set to CvSeq
CV_INLINE static CvContour*
icvMSERToContour( CvMSERConnectedComp* comp,
CvMemStorage* storage )
{
CvSeq* _contour = cvCreateSeq( CV_SEQ_KIND_GENERIC|CV_32SC2, sizeof(CvContour), sizeof(CvPoint), storage );
CvContour* contour = (CvContour*)_contour;
cvSeqPushMulti( _contour, 0, comp->history->size );
CvLinkedPoint* lpt = comp->head;
for ( int i = 0; i < comp->history->size; i++ )
{
CvPoint* pt = CV_GET_SEQ_ELEM( CvPoint, _contour, i );
pt->x = lpt->pt.x;
pt->y = lpt->pt.y;
lpt = lpt->next;
}
cvBoundingRect( contour );
return contour;
}
// to preprocess src image to following format
// 32-bit image
// > 0 is available, < 0 is visited
// 17~19 bits is the direction
// 8~11 bits is the bucket it falls to (for BitScanForward)
// 0~8 bits is the color
static int*
icvPreprocessMSER_8UC1( CvMat* img,
int*** heap_cur,
CvMat* src,
CvMat* mask )
{
int srccpt = src->step-src->cols;
int cpt_1 = img->cols-src->cols-1;
int* imgptr = img->data.i;
int* startptr;
int level_size[256];
for ( int i = 0; i < 256; i++ )
level_size[i] = 0;
for ( int i = 0; i < src->cols+2; i++ )
{
*imgptr = -1;
imgptr++;
}
imgptr += cpt_1-1;
uchar* srcptr = src->data.ptr;
if ( mask )
{
startptr = 0;
uchar* maskptr = mask->data.ptr;
for ( int i = 0; i < src->rows; i++ )
{
*imgptr = -1;
imgptr++;
for ( int j = 0; j < src->cols; j++ )
{
if ( *maskptr )
{
if ( !startptr )
startptr = imgptr;
*srcptr = 0xff-*srcptr;
level_size[*srcptr]++;
*imgptr = ((*srcptr>>5)<<8)|(*srcptr);
} else {
*imgptr = -1;
}
imgptr++;
srcptr++;
maskptr++;
}
*imgptr = -1;
imgptr += cpt_1;
srcptr += srccpt;
maskptr += srccpt;
}
} else {
startptr = imgptr+img->cols+1;
for ( int i = 0; i < src->rows; i++ )
{
*imgptr = -1;
imgptr++;
for ( int j = 0; j < src->cols; j++ )
{
*srcptr = 0xff-*srcptr;
level_size[*srcptr]++;
*imgptr = ((*srcptr>>5)<<8)|(*srcptr);
imgptr++;
srcptr++;
}
*imgptr = -1;
imgptr += cpt_1;
srcptr += srccpt;
}
}
for ( int i = 0; i < src->cols+2; i++ )
{
*imgptr = -1;
imgptr++;
}
heap_cur[0][0] = 0;
for ( int i = 1; i < 256; i++ )
{
heap_cur[i] = heap_cur[i-1]+level_size[i-1]+1;
heap_cur[i][0] = 0;
}
return startptr;
}
static void
icvExtractMSER_8UC1_Pass( int* ioptr,
int* imgptr,
int*** heap_cur,
CvLinkedPoint* ptsptr,
CvMSERGrowHistory* histptr,
CvMSERConnectedComp* comptr,
int step,
int stepmask,
int stepgap,
CvMSERParams params,
int color,
CvSeq* contours,
CvMemStorage* storage )
{
comptr->grey_level = 256;
comptr++;
comptr->grey_level = (*imgptr)&0xff;
icvInitMSERComp( comptr );
*imgptr |= 0x80000000;
heap_cur += (*imgptr)&0xff;
int dir[] = { 1, step, -1, -step };
#ifdef __INTRIN_ENABLED__
unsigned long heapbit[] = { 0, 0, 0, 0, 0, 0, 0, 0 };
unsigned long* bit_cur = heapbit+(((*imgptr)&0x700)>>8);
#endif
for ( ; ; )
{
// take tour of all the 4 directions
while ( ((*imgptr)&0x70000) < 0x40000 )
{
// get the neighbor
int* imgptr_nbr = imgptr+dir[((*imgptr)&0x70000)>>16];
if ( *imgptr_nbr >= 0 ) // if the neighbor is not visited yet
{
*imgptr_nbr |= 0x80000000; // mark it as visited
if ( ((*imgptr_nbr)&0xff) < ((*imgptr)&0xff) )
{
// when the value of neighbor smaller than current
// push current to boundary heap and make the neighbor to be the current one
// create an empty comp
(*heap_cur)++;
**heap_cur = imgptr;
*imgptr += 0x10000;
heap_cur += ((*imgptr_nbr)&0xff)-((*imgptr)&0xff);
#ifdef __INTRIN_ENABLED__
_bitset( bit_cur, (*imgptr)&0x1f );
bit_cur += (((*imgptr_nbr)&0x700)-((*imgptr)&0x700))>>8;
#endif
imgptr = imgptr_nbr;
comptr++;
icvInitMSERComp( comptr );
comptr->grey_level = (*imgptr)&0xff;
continue;
} else {
// otherwise, push the neighbor to boundary heap
heap_cur[((*imgptr_nbr)&0xff)-((*imgptr)&0xff)]++;
*heap_cur[((*imgptr_nbr)&0xff)-((*imgptr)&0xff)] = imgptr_nbr;
#ifdef __INTRIN_ENABLED__
_bitset( bit_cur+((((*imgptr_nbr)&0x700)-((*imgptr)&0x700))>>8), (*imgptr_nbr)&0x1f );
#endif
}
}
*imgptr += 0x10000;
}
int i = imgptr-ioptr;
ptsptr->pt = cvPoint( i&stepmask, i>>stepgap );
// get the current location
icvAccumulateMSERComp( comptr, ptsptr );
ptsptr++;
// get the next pixel from boundary heap
if ( **heap_cur )
{
imgptr = **heap_cur;
(*heap_cur)--;
#ifdef __INTRIN_ENABLED__
if ( !**heap_cur )
_bitreset( bit_cur, (*imgptr)&0x1f );
#endif
} else {
#ifdef __INTRIN_ENABLED__
bool found_pixel = 0;
unsigned long pixel_val;
for ( int i = ((*imgptr)&0x700)>>8; i < 8; i++ )
{
if ( _BitScanForward( &pixel_val, *bit_cur ) )
{
found_pixel = 1;
pixel_val += i<<5;
heap_cur += pixel_val-((*imgptr)&0xff);
break;
}
bit_cur++;
}
if ( found_pixel )
#else
heap_cur++;
unsigned long pixel_val = 0;
for ( unsigned long i = ((*imgptr)&0xff)+1; i < 256; i++ )
{
if ( **heap_cur )
{
pixel_val = i;
break;
}
heap_cur++;
}
if ( pixel_val )
#endif
{
imgptr = **heap_cur;
(*heap_cur)--;
#ifdef __INTRIN_ENABLED__
if ( !**heap_cur )
_bitreset( bit_cur, pixel_val&0x1f );
#endif
if ( pixel_val < comptr[-1].grey_level )
{
// check the stablity and push a new history, increase the grey level
if ( icvMSERStableCheck( comptr, params ) )
{
CvContour* contour = icvMSERToContour( comptr, storage );
contour->color = color;
cvSeqPush( contours, &contour );
}
icvMSERNewHistory( comptr, histptr );
comptr[0].grey_level = pixel_val;
histptr++;
} else {
// keep merging top two comp in stack until the grey level >= pixel_val
for ( ; ; )
{
comptr--;
icvMSERMergeComp( comptr+1, comptr, comptr, histptr );
histptr++;
if ( pixel_val <= comptr[0].grey_level )
break;
if ( pixel_val < comptr[-1].grey_level )
{
// check the stablity here otherwise it wouldn't be an ER
if ( icvMSERStableCheck( comptr, params ) )
{
CvContour* contour = icvMSERToContour( comptr, storage );
contour->color = color;
cvSeqPush( contours, &contour );
}
icvMSERNewHistory( comptr, histptr );
comptr[0].grey_level = pixel_val;
histptr++;
break;
}
}
}
} else
break;
}
}
}
static void
icvExtractMSER_8UC1( CvMat* src,
CvMat* mask,
CvSeq* contours,
CvMemStorage* storage,
CvMSERParams params )
{
int step = 8;
int stepgap = 3;
while ( step < src->step+2 )
{
step <<= 1;
stepgap++;
}
int stepmask = step-1;
// to speedup the process, make the width to be 2^N
CvMat* img = cvCreateMat( src->rows+2, step, CV_32SC1 );
int* ioptr = img->data.i+step+1;
int* imgptr;
// pre-allocate boundary heap
int** heap = (int**)cvAlloc( (src->rows*src->cols+256)*sizeof(heap[0]) );
int** heap_start[256];
heap_start[0] = heap;
// pre-allocate linked point and grow history
CvLinkedPoint* pts = (CvLinkedPoint*)cvAlloc( src->rows*src->cols*sizeof(pts[0]) );
CvMSERGrowHistory* history = (CvMSERGrowHistory*)cvAlloc( src->rows*src->cols*sizeof(history[0]) );
CvMSERConnectedComp comp[257];
// darker to brighter (MSER-)
imgptr = icvPreprocessMSER_8UC1( img, heap_start, src, mask );
icvExtractMSER_8UC1_Pass( ioptr, imgptr, heap_start, pts, history, comp, step, stepmask, stepgap, params, -1, contours, storage );
// brighter to darker (MSER+)
imgptr = icvPreprocessMSER_8UC1( img, heap_start, src, mask );
icvExtractMSER_8UC1_Pass( ioptr, imgptr, heap_start, pts, history, comp, step, stepmask, stepgap, params, 1, contours, storage );
// clean up
cvFree( &history );
cvFree( &heap );
cvFree( &pts );
cvReleaseMat( &img );
}
struct CvMSCRNode;
typedef struct CvTempMSCR
{
CvMSCRNode* head;
CvMSCRNode* tail;
double m; // the margin used to prune area later
int size;
} CvTempMSCR;
typedef struct CvMSCRNode
{
CvMSCRNode* shortcut;
// to make the finding of root less painful
CvMSCRNode* prev;
CvMSCRNode* next;
// a point double-linked list
CvTempMSCR* tmsr;
// the temporary msr (set to NULL at every re-initialise)
CvTempMSCR* gmsr;
// the global msr (once set, never to NULL)
int index;
// the index of the node, at this point, it should be x at the first 16-bits, and y at the last 16-bits.
int rank;
int reinit;
int size, sizei;
double dt, di;
double s;
} CvMSCRNode;
typedef struct CvMSCREdge
{
double chi;
CvMSCRNode* left;
CvMSCRNode* right;
} CvMSCREdge;
#include "chi_table.h"
CV_INLINE static double
icvChisquaredDistance( uchar* x, uchar* y )
{
return (double)((x[0]-y[0])*(x[0]-y[0]))/(double)(x[0]+y[0]+1e-10)+
(double)((x[1]-y[1])*(x[1]-y[1]))/(double)(x[1]+y[1]+1e-10)+
(double)((x[2]-y[2])*(x[2]-y[2]))/(double)(x[2]+y[2]+1e-10);
}
CV_INLINE static void
icvInitMSCRNode( CvMSCRNode* node )
{
node->gmsr = node->tmsr = NULL;
node->reinit = 0xffff;
node->rank = 0;
node->sizei = node->size = 1;
node->prev = node->next = node->shortcut = node;
}
// the preprocess to get the edge list with proper gaussian blur
static int
icvPreprocessMSER_8UC3( CvMSCRNode* node,
CvMSCREdge* edge,
double* total,
CvMat* src,
CvMat* mask,
CvMat* dx,
CvMat* dy,
int Ne,
int edge_blur_size )
{
int srccpt = src->step-src->cols*3;
uchar* srcptr = src->data.ptr;
uchar* lastptr = src->data.ptr+3;
double* dxptr = dx->data.db;
for ( int i = 0; i < src->rows; i++ )
{
for ( int j = 0; j < src->cols-1; j++ )
{
*dxptr = icvChisquaredDistance( srcptr, lastptr );
dxptr++;
srcptr += 3;
lastptr += 3;
}
srcptr += srccpt+3;
lastptr += srccpt+3;
}
srcptr = src->data.ptr;
lastptr = src->data.ptr+src->step;
double* dyptr = dy->data.db;
for ( int i = 0; i < src->rows-1; i++ )
{
for ( int j = 0; j < src->cols; j++ )
{
*dyptr = icvChisquaredDistance( srcptr, lastptr );
dyptr++;
srcptr += 3;
lastptr += 3;
}
srcptr += srccpt;
lastptr += srccpt;
}
// get dx and dy and blur it
if ( edge_blur_size >= 1 )
{
cvSmooth( dx, dx, CV_GAUSSIAN, edge_blur_size, edge_blur_size );
cvSmooth( dy, dy, CV_GAUSSIAN, edge_blur_size, edge_blur_size );
}
dxptr = dx->data.db;
dyptr = dy->data.db;
// assian dx, dy to proper edge list and initialize mscr node
// the nasty code here intended to avoid extra loops
if ( mask )
{
Ne = 0;
int maskcpt = mask->step-mask->cols+1;
uchar* maskptr = mask->data.ptr;
CvMSCRNode* nodeptr = node;
icvInitMSCRNode( nodeptr );
nodeptr->index = 0;
*total += edge->chi = *dxptr;
if ( maskptr[0] && maskptr[1] )
{
edge->left = nodeptr;
edge->right = nodeptr+1;
edge++;
Ne++;
}
dxptr++;
nodeptr++;
maskptr++;
for ( int i = 1; i < src->cols-1; i++ )
{
icvInitMSCRNode( nodeptr );
nodeptr->index = i;
if ( maskptr[0] && maskptr[1] )
{
*total += edge->chi = *dxptr;
edge->left = nodeptr;
edge->right = nodeptr+1;
edge++;
Ne++;
}
dxptr++;
nodeptr++;
maskptr++;
}
icvInitMSCRNode( nodeptr );
nodeptr->index = src->cols-1;
nodeptr++;
maskptr += maskcpt;
for ( int i = 1; i < src->rows-1; i++ )
{
icvInitMSCRNode( nodeptr );
nodeptr->index = i<<16;
if ( maskptr[0] )
{
if ( maskptr[-mask->step] )
{
*total += edge->chi = *dyptr;
edge->left = nodeptr-src->cols;
edge->right = nodeptr;
edge++;
Ne++;
}
if ( maskptr[1] )
{
*total += edge->chi = *dxptr;
edge->left = nodeptr;
edge->right = nodeptr+1;
edge++;
Ne++;
}
}
dyptr++;
dxptr++;
nodeptr++;
maskptr++;
for ( int j = 1; j < src->cols-1; j++ )
{
icvInitMSCRNode( nodeptr );
nodeptr->index = (i<<16)|j;
if ( maskptr[0] )
{
if ( maskptr[-mask->step] )
{
*total += edge->chi = *dyptr;
edge->left = nodeptr-src->cols;
edge->right = nodeptr;
edge++;
Ne++;
}
if ( maskptr[1] )
{
*total += edge->chi = *dxptr;
edge->left = nodeptr;
edge->right = nodeptr+1;
edge++;
Ne++;
}
}
dyptr++;
dxptr++;
nodeptr++;
maskptr++;
}
icvInitMSCRNode( nodeptr );
nodeptr->index = (i<<16)|(src->cols-1);
if ( maskptr[0] && maskptr[-mask->step] )
{
*total += edge->chi = *dyptr;
edge->left = nodeptr-src->cols;
edge->right = nodeptr;
edge++;
Ne++;
}
dyptr++;
nodeptr++;
maskptr += maskcpt;
}
icvInitMSCRNode( nodeptr );
nodeptr->index = (src->rows-1)<<16;
if ( maskptr[0] )
{
if ( maskptr[1] )
{
*total += edge->chi = *dxptr;
edge->left = nodeptr;
edge->right = nodeptr+1;
edge++;
Ne++;
}
if ( maskptr[-mask->step] )
{
*total += edge->chi = *dyptr;
edge->left = nodeptr-src->cols;
edge->right = nodeptr;
edge++;
Ne++;
}
}
dxptr++;
dyptr++;
nodeptr++;
maskptr++;
for ( int i = 1; i < src->cols-1; i++ )
{
icvInitMSCRNode( nodeptr );
nodeptr->index = ((src->rows-1)<<16)|i;
if ( maskptr[0] )
{
if ( maskptr[1] )
{
*total += edge->chi = *dxptr;
edge->left = nodeptr;
edge->right = nodeptr+1;
edge++;
Ne++;
}
if ( maskptr[-mask->step] )
{
*total += edge->chi = *dyptr;
edge->left = nodeptr-src->cols;
edge->right = nodeptr;
edge++;
Ne++;
}
}
dxptr++;
dyptr++;
nodeptr++;
maskptr++;
}
icvInitMSCRNode( nodeptr );
nodeptr->index = ((src->rows-1)<<16)|(src->cols-1);
if ( maskptr[0] && maskptr[-mask->step] )
{
*total += edge->chi = *dyptr;
edge->left = nodeptr-src->cols;
edge->right = nodeptr;
Ne++;
}
} else {
CvMSCRNode* nodeptr = node;
icvInitMSCRNode( nodeptr );
nodeptr->index = 0;
*total += edge->chi = *dxptr;
dxptr++;
edge->left = nodeptr;
edge->right = nodeptr+1;
edge++;
nodeptr++;
for ( int i = 1; i < src->cols-1; i++ )
{
icvInitMSCRNode( nodeptr );
nodeptr->index = i;
*total += edge->chi = *dxptr;
dxptr++;
edge->left = nodeptr;
edge->right = nodeptr+1;
edge++;
nodeptr++;
}
icvInitMSCRNode( nodeptr );
nodeptr->index = src->cols-1;
nodeptr++;
for ( int i = 1; i < src->rows-1; i++ )
{
icvInitMSCRNode( nodeptr );
nodeptr->index = i<<16;
*total += edge->chi = *dyptr;
dyptr++;
edge->left = nodeptr-src->cols;
edge->right = nodeptr;
edge++;
*total += edge->chi = *dxptr;
dxptr++;
edge->left = nodeptr;
edge->right = nodeptr+1;
edge++;
nodeptr++;
for ( int j = 1; j < src->cols-1; j++ )
{
icvInitMSCRNode( nodeptr );
nodeptr->index = (i<<16)|j;
*total += edge->chi = *dyptr;
dyptr++;
edge->left = nodeptr-src->cols;
edge->right = nodeptr;
edge++;
*total += edge->chi = *dxptr;
dxptr++;
edge->left = nodeptr;
edge->right = nodeptr+1;
edge++;
nodeptr++;
}
icvInitMSCRNode( nodeptr );
nodeptr->index = (i<<16)|(src->cols-1);
*total += edge->chi = *dyptr;
dyptr++;
edge->left = nodeptr-src->cols;
edge->right = nodeptr;
edge++;
nodeptr++;
}
icvInitMSCRNode( nodeptr );
nodeptr->index = (src->rows-1)<<16;
*total += edge->chi = *dxptr;
dxptr++;
edge->left = nodeptr;
edge->right = nodeptr+1;
edge++;
*total += edge->chi = *dyptr;
dyptr++;
edge->left = nodeptr-src->cols;
edge->right = nodeptr;
edge++;
nodeptr++;
for ( int i = 1; i < src->cols-1; i++ )
{
icvInitMSCRNode( nodeptr );
nodeptr->index = ((src->rows-1)<<16)|i;
*total += edge->chi = *dxptr;
dxptr++;
edge->left = nodeptr;
edge->right = nodeptr+1;
edge++;
*total += edge->chi = *dyptr;
dyptr++;
edge->left = nodeptr-src->cols;
edge->right = nodeptr;
edge++;
nodeptr++;
}
icvInitMSCRNode( nodeptr );
nodeptr->index = ((src->rows-1)<<16)|(src->cols-1);
*total += edge->chi = *dyptr;
edge->left = nodeptr-src->cols;
edge->right = nodeptr;
}
return Ne;
}
#define cmp_mscr_edge(edge1, edge2) /
((edge1).chi < (edge2).chi)
static CV_IMPLEMENT_QSORT( icvQuickSortMSCREdge, CvMSCREdge, cmp_mscr_edge )
// to find the root of one region
CV_INLINE static CvMSCRNode*
icvFindMSCR( CvMSCRNode* x )
{
CvMSCRNode* prev = x;
CvMSCRNode* next;
for ( ; ; )
{
next = x->shortcut;
x->shortcut = prev;
if ( next == x ) break;
prev= x;
x = next;
}
CvMSCRNode* root = x;
for ( ; ; )
{
prev = x->shortcut;
x->shortcut = root;
if ( prev == x ) break;
x = prev;
}
return root;
}
// the stable mscr should be:
// bigger than min_area and smaller than max_area
// differ from its ancestor more than min_diversity
CV_INLINE static bool
icvMSCRStableCheck( CvMSCRNode* x,
CvMSERParams params )
{
if ( x->size <= params.min_area || x->size >= params.max_area )
return 0;
if ( x->gmsr == NULL )
return 1;
double div = (double)(x->size-x->gmsr->size)/(double)x->size;
return div > params.min_diversity;
}
static void
icvExtractMSER_8UC3( CvMat* src,
CvMat* mask,
CvSeq* contours,
CvMemStorage* storage,
CvMSERParams params )
{
CvMSCRNode* map = (CvMSCRNode*)cvAlloc( src->cols*src->rows*sizeof(map[0]) );
int Ne = src->cols*src->rows*2-src->cols-src->rows;
CvMSCREdge* edge = (CvMSCREdge*)cvAlloc( Ne*sizeof(edge[0]) );
CvTempMSCR* mscr = (CvTempMSCR*)cvAlloc( src->cols*src->rows*sizeof(mscr[0]) );
double emean = 0;
CvMat* dx = cvCreateMat( src->rows, src->cols-1, CV_64FC1 );
CvMat* dy = cvCreateMat( src->rows-1, src->cols, CV_64FC1 );
Ne = icvPreprocessMSER_8UC3( map, edge, &emean, src, mask, dx, dy, Ne, params.edge_blur_size );
emean = emean / (double)Ne;
icvQuickSortMSCREdge( edge, Ne, 0 );
CvMSCREdge* edge_ub = edge+Ne;
CvMSCREdge* edgeptr = edge;
CvTempMSCR* mscrptr = mscr;
// the evolution process
for ( int i = 0; i < params.max_evolution; i++ )
{
double k = (double)i/(double)params.max_evolution*(TABLE_SIZE-1);
int ti = floor(k);
double reminder = k-ti;
double thres = emean*(chitab3[ti]*(1-reminder)+chitab3[ti+1]*reminder);
// to process all the edges in the list that chi < thres
while ( edgeptr < edge_ub && edgeptr->chi < thres )
{
CvMSCRNode* lr = icvFindMSCR( edgeptr->left );
CvMSCRNode* rr = icvFindMSCR( edgeptr->right );
// get the region root (who is responsible)
if ( lr != rr )
{
// rank idea take from: N-tree Disjoint-Set Forests for Maximally Stable Extremal Regions
if ( rr->rank > lr->rank )
{
CvMSCRNode* tmp;
CV_SWAP( lr, rr, tmp );
} else if ( lr->rank == rr->rank ) {
// at the same rank, we will compare the size
if ( lr->size > rr->size )
{
CvMSCRNode* tmp;
CV_SWAP( lr, rr, tmp );
}
lr->rank++;
}
rr->shortcut = lr;
lr->size += rr->size;
// join rr to the end of list lr (lr is a endless double-linked list)
lr->prev->next = rr;
lr->prev = rr->prev;
rr->prev->next = lr;
rr->prev = lr;
// area threshold force to reinitialize
if ( lr->size > (lr->size-rr->size)*params.area_threshold )
{
lr->sizei = lr->size;
lr->reinit = i;
if ( lr->tmsr != NULL )
{
lr->tmsr->m = lr->dt-lr->di;
lr->tmsr = NULL;
}
lr->di = edgeptr->chi;
lr->s = 1e10;
}
lr->dt = edgeptr->chi;
if ( i > lr->reinit )
{
double s = (double)(lr->size-lr->sizei)/(lr->dt-lr->di);
if ( s < lr->s )
{
// skip the first one and check stablity
if ( i > lr->reinit+1 && icvMSCRStableCheck( lr, params ) )
{
if ( lr->tmsr == NULL )
{
lr->gmsr = lr->tmsr = mscrptr;
mscrptr++;
}
lr->tmsr->size = lr->size;
lr->tmsr->head = lr;
lr->tmsr->tail = lr->prev;
lr->tmsr->m = 0;
}
lr->s = s;
}
}
}
edgeptr++;
}
if ( edgeptr >= edge_ub )
break;
}
for ( CvTempMSCR* ptr = mscr; ptr < mscrptr; ptr++ )
// to prune area with margin less than min_margin
if ( ptr->m > params.min_margin )
{
CvSeq* _contour = cvCreateSeq( CV_SEQ_KIND_GENERIC|CV_32SC2, sizeof(CvContour), sizeof(CvPoint), storage );
cvSeqPushMulti( _contour, 0, ptr->size );
CvMSCRNode* lpt = ptr->head;
for ( int i = 0; i < ptr->size; i++ )
{
CvPoint* pt = CV_GET_SEQ_ELEM( CvPoint, _contour, i );
pt->x = (lpt->index)&0xffff;
pt->y = (lpt->index)>>16;
lpt = lpt->next;
}
CvContour* contour = (CvContour*)_contour;
cvBoundingRect( contour );
contour->color = 0;
cvSeqPush( contours, &contour );
}
cvReleaseMat( &dx );
cvReleaseMat( &dy );
cvFree( &mscr );
cvFree( &edge );
cvFree( &map );
}
void
cvExtractMSER( CvArr* _img,
CvArr* _mask,
CvSeq** _contours,
CvMemStorage* storage,
CvMSERParams params )
{
CvMat srchdr, *src = cvGetMat( _img, &srchdr );
CvMat maskhdr, *mask = _mask ? cvGetMat( _mask, &maskhdr ) : 0;
CvSeq* contours = 0;
CV_FUNCNAME( "cvExtractMSER" );
__BEGIN__;
CV_ASSERT(src != 0);
CV_ASSERT(CV_MAT_TYPE(src->type) == CV_8UC1 || CV_MAT_TYPE(src->type) == CV_8UC3);
CV_ASSERT(mask == 0 || (CV_ARE_SIZES_EQ(src, mask) && CV_MAT_TYPE(mask->type) == CV_8UC1));
CV_ASSERT(storage != 0);
contours = *_contours = cvCreateSeq( 0, sizeof(CvSeq), sizeof(CvSeq*), storage );
// choose different method for different image type
// for grey image, it is: Linear Time Maximally Stable Extremal Regions
// for color image, it is: Maximally Stable Colour Regions for Recognition and Matching
switch ( CV_MAT_TYPE(src->type) )
{
case CV_8UC1:
icvExtractMSER_8UC1( src, mask, contours, storage, params );
break;
case CV_8UC3:
icvExtractMSER_8UC3( src, mask, contours, storage, params );
break;
}
__END__;
}

/*
**  Source file for mscr. (c) 2007 Per-Erik Forssen
**
**  This program is free software; you can redistribute it and/or modify
**  it under the terms of the GNU General Public License as published by
**  the Free Software Foundation; either version 2 of the License, or
**  (at your option) any later version.
**
**  See the file COPYING for details.
**
*/
#define TABLE_SIZE 400
static double chitab1[]={0,  9.86678e-06,  3.94675e-05,  8.88034e-05,
0.000157876,  0.000246689,  0.000355245,  0.000483549,
0.000631605,  0.00079942,  0.000987,  0.00119435,
0.00142149,  0.00166841,  0.00193513,  0.00222166,
0.00252802,  0.0028542,  0.00320024,  0.00356613,
0.0039519,  0.00435756,  0.00478312,  0.00522861,
0.00569404,  0.00617943,  0.00668479,  0.00721016,
0.00775554,  0.00832097,  0.00890646,  0.00951204,
0.0101377,  0.0107836,  0.0114495,  0.0121357,
0.0128421,  0.0135688,  0.0143157,  0.0150829,
0.0158704,  0.0166784,  0.0175067,  0.0183554,
0.0192247,  0.0201144,  0.0210247,  0.0219555,
0.022907,  0.0238791,  0.0248719,  0.0258854,
0.0269197,  0.0279749,  0.0290509,  0.0301478,
0.0312656,  0.0324045,  0.0335643,  0.0347453,
0.0359474,  0.0371708,  0.0384153,  0.0396812,
0.0409684,  0.042277,  0.043607,  0.0449586,
0.0463317,  0.0477265,  0.0491429,  0.0505811,
0.0520411,  0.053523,  0.0550268,  0.0565526,
0.0581005,  0.0596705,  0.0612627,  0.0628771,
0.0645139,  0.0661731,  0.0678548,  0.069559,
0.0712859,  0.0730354,  0.0748078,  0.076603,
0.0784211,  0.0802622,  0.0821265,  0.0840139,
0.0859246,  0.0878587,  0.0898161,  0.0917972,
0.0938018,  0.0958301,  0.0978823,  0.0999583,
0.102058,  0.104182,  0.106331,  0.108503,
0.1107,  0.112922,  0.115168,  0.117438,
0.119734,  0.122054,  0.1244,  0.12677,
0.129166,  0.131587,  0.134034,  0.136506,
0.139004,  0.141528,  0.144077,  0.146653,
0.149255,  0.151884,  0.154538,  0.15722,
0.159928,  0.162663,  0.165425,  0.168214,
0.17103,  0.173874,  0.176745,  0.179644,
0.182571,  0.185526,  0.188509,  0.19152,
0.19456,  0.197629,  0.200726,  0.203852,
0.207008,  0.210193,  0.213407,  0.21665,
0.219924,  0.223228,  0.226562,  0.229926,
0.23332,  0.236746,  0.240202,  0.24369,
0.247208,  0.250759,  0.254341,  0.257955,
0.261601,  0.265279,  0.26899,  0.272734,
0.276511,  0.28032,  0.284164,  0.288041,
0.291952,  0.295897,  0.299876,  0.303891,
0.307939,  0.312024,  0.316143,  0.320298,
0.324489,  0.328716,  0.33298,  0.33728,
0.341617,  0.345992,  0.350404,  0.354854,
0.359341,  0.363868,  0.368433,  0.373037,
0.377681,  0.382364,  0.387087,  0.39185,
0.396655,  0.4015,  0.406386,  0.411314,
0.416284,  0.421297,  0.426352,  0.43145,
0.436592,  0.441778,  0.447008,  0.452282,
0.457602,  0.462967,  0.468378,  0.473835,
0.479339,  0.48489,  0.490489,  0.496136,
0.501831,  0.507575,  0.513369,  0.519212,
0.525106,  0.531051,  0.537047,  0.543095,
0.549195,  0.555349,  0.561556,  0.567817,
0.574132,  0.580503,  0.586929,  0.593412,
0.599952,  0.606549,  0.613204,  0.619918,
0.626692,  0.633525,  0.64042,  0.647375,
0.654393,  0.661474,  0.668618,  0.675827,
0.6831,  0.690439,  0.697845,  0.705318,
0.712859,  0.720469,  0.728149,  0.7359,
0.743722,  0.751617,  0.759584,  0.767626,
0.775743,  0.783937,  0.792207,  0.800556,
0.808983,  0.817491,  0.82608,  0.834751,
0.843506,  0.852345,  0.86127,  0.870282,
0.879382,  0.888572,  0.897852,  0.907224,
0.916689,  0.926249,  0.935904,  0.945658,
0.95551,  0.965462,  0.975516,  0.985674,
0.995937,  1.00631,  1.01678,  1.02737,
1.03807,  1.04888,  1.05981,  1.07086,
1.08202,  1.09331,  1.10472,  1.11625,
1.12791,  1.13971,  1.15163,  1.16369,
1.17588,  1.18822,  1.20069,  1.21331,
1.22608,  1.239,  1.25207,  1.26529,
1.27868,  1.29222,  1.30593,  1.31981,
1.33386,  1.34808,  1.36248,  1.37707,
1.39184,  1.4068,  1.42195,  1.4373,
1.45286,  1.46862,  1.48459,  1.50078,
1.51719,  1.53383,  1.5507,  1.5678,
1.58515,  1.60275,  1.6206,  1.63872,
1.6571,  1.67576,  1.6947,  1.71393,
1.73346,  1.75329,  1.77344,  1.7939,
1.8147,  1.83585,  1.85734,  1.87919,
1.90142,  1.92403,  1.94703,  1.97044,
1.99427,  2.01854,  2.04325,  2.06843,
2.09409,  2.12024,  2.14691,  2.17411,
2.20186,  2.23018,  2.25909,  2.28862,
2.31879,  2.34962,  2.38115,  2.4134,
2.4464,  2.48019,  2.5148,  2.55027,
2.58663,  2.62393,  2.66223,  2.70155,
2.74197,  2.78353,  2.8263,  2.87034,
2.91574,  2.96256,  3.01089,  3.06083,
3.11249,  3.16597,  3.22139,  3.27891,
3.33867,  3.40085,  3.46563,  3.53323,
3.6039,  3.67791,  3.75558,  3.83726,
3.92338,  4.01442,  4.11095,  4.21364,
4.3233,  4.4409,  4.56764,  4.70499,
4.85482,  5.01955,  5.20235,  5.40752,
5.64113,  5.91206,  6.2341,  6.63044,
7.14465,  7.87491,  9.13601,  400 };
static double chitab3[]={0,  0.0150057,  0.0239478,  0.0315227,
0.0383427,  0.0446605,  0.0506115,  0.0562786,
0.0617174,  0.0669672,  0.0720573,  0.0770099,
0.081843,  0.0865705,  0.0912043,  0.0957541,
0.100228,  0.104633,  0.108976,  0.113261,
0.117493,  0.121676,  0.125814,  0.12991,
0.133967,  0.137987,  0.141974,  0.145929,
0.149853,  0.15375,  0.15762,  0.161466,
0.165287,  0.169087,  0.172866,  0.176625,
0.180365,  0.184088,  0.187794,  0.191483,
0.195158,  0.198819,  0.202466,  0.2061,
0.209722,  0.213332,  0.216932,  0.220521,
0.2241,  0.22767,  0.231231,  0.234783,
0.238328,  0.241865,  0.245395,  0.248918,
0.252435,  0.255947,  0.259452,  0.262952,
0.266448,  0.269939,  0.273425,  0.276908,
0.280386,  0.283862,  0.287334,  0.290803,
0.29427,  0.297734,  0.301197,  0.304657,
0.308115,  0.311573,  0.315028,  0.318483,
0.321937,  0.32539,  0.328843,  0.332296,
0.335749,  0.339201,  0.342654,  0.346108,
0.349562,  0.353017,  0.356473,  0.35993,
0.363389,  0.366849,  0.37031,  0.373774,
0.377239,  0.380706,  0.384176,  0.387648,
0.391123,  0.3946,  0.39808,  0.401563,
0.405049,  0.408539,  0.412032,  0.415528,
0.419028,  0.422531,  0.426039,  0.429551,
0.433066,  0.436586,  0.440111,  0.44364,
0.447173,  0.450712,  0.454255,  0.457803,
0.461356,  0.464915,  0.468479,  0.472049,
0.475624,  0.479205,  0.482792,  0.486384,
0.489983,  0.493588,  0.4972,  0.500818,
0.504442,  0.508073,  0.511711,  0.515356,
0.519008,  0.522667,  0.526334,  0.530008,
0.533689,  0.537378,  0.541075,  0.54478,
0.548492,  0.552213,  0.555942,  0.55968,
0.563425,  0.56718,  0.570943,  0.574715,
0.578497,  0.582287,  0.586086,  0.589895,
0.593713,  0.597541,  0.601379,  0.605227,
0.609084,  0.612952,  0.61683,  0.620718,
0.624617,  0.628526,  0.632447,  0.636378,
0.64032,  0.644274,  0.648239,  0.652215,
0.656203,  0.660203,  0.664215,  0.668238,
0.672274,  0.676323,  0.680384,  0.684457,
0.688543,  0.692643,  0.696755,  0.700881,
0.70502,  0.709172,  0.713339,  0.717519,
0.721714,  0.725922,  0.730145,  0.734383,
0.738636,  0.742903,  0.747185,  0.751483,
0.755796,  0.760125,  0.76447,  0.768831,
0.773208,  0.777601,  0.782011,  0.786438,
0.790882,  0.795343,  0.799821,  0.804318,
0.808831,  0.813363,  0.817913,  0.822482,
0.827069,  0.831676,  0.836301,  0.840946,
0.84561,  0.850295,  0.854999,  0.859724,
0.864469,  0.869235,  0.874022,  0.878831,
0.883661,  0.888513,  0.893387,  0.898284,
0.903204,  0.908146,  0.913112,  0.918101,
0.923114,  0.928152,  0.933214,  0.938301,
0.943413,  0.94855,  0.953713,  0.958903,
0.964119,  0.969361,  0.974631,  0.979929,
0.985254,  0.990608,  0.99599,  1.0014,
1.00684,  1.01231,  1.01781,  1.02335,
1.02891,  1.0345,  1.04013,  1.04579,
1.05148,  1.05721,  1.06296,  1.06876,
1.07459,  1.08045,  1.08635,  1.09228,
1.09826,  1.10427,  1.11032,  1.1164,
1.12253,  1.1287,  1.1349,  1.14115,
1.14744,  1.15377,  1.16015,  1.16656,
1.17303,  1.17954,  1.18609,  1.19269,
1.19934,  1.20603,  1.21278,  1.21958,
1.22642,  1.23332,  1.24027,  1.24727,
1.25433,  1.26144,  1.26861,  1.27584,
1.28312,  1.29047,  1.29787,  1.30534,
1.31287,  1.32046,  1.32812,  1.33585,
1.34364,  1.3515,  1.35943,  1.36744,
1.37551,  1.38367,  1.39189,  1.4002,
1.40859,  1.41705,  1.42561,  1.43424,
1.44296,  1.45177,  1.46068,  1.46967,
1.47876,  1.48795,  1.49723,  1.50662,
1.51611,  1.52571,  1.53541,  1.54523,
1.55517,  1.56522,  1.57539,  1.58568,
1.59611,  1.60666,  1.61735,  1.62817,
1.63914,  1.65025,  1.66152,  1.67293,
1.68451,  1.69625,  1.70815,  1.72023,
1.73249,  1.74494,  1.75757,  1.77041,
1.78344,  1.79669,  1.81016,  1.82385,
1.83777,  1.85194,  1.86635,  1.88103,
1.89598,  1.91121,  1.92674,  1.94257,
1.95871,  1.97519,  1.99201,  2.0092,
2.02676,  2.04471,  2.06309,  2.08189,
2.10115,  2.12089,  2.14114,  2.16192,
2.18326,  2.2052,  2.22777,  2.25101,
2.27496,  2.29966,  2.32518,  2.35156,
2.37886,  2.40717,  2.43655,  2.46709,
2.49889,  2.53206,  2.56673,  2.60305,
2.64117,  2.6813,  2.72367,  2.76854,
2.81623,  2.86714,  2.92173,  2.98059,
3.04446,  3.1143,  3.19135,  3.27731,
3.37455,  3.48653,  3.61862,  3.77982,
3.98692,  4.2776,  4.77167,  133.333 };

#include "cvmser.h"
#include <iostream>
#include "cv.h"
#include "highgui.h"
static CvScalar colors[] =
{
{{0,0,255}},
{{0,128,255}},
{{0,255,255}},
{{0,255,0}},
{{255,128,0}},
{{255,255,0}},
{{255,0,0}},
{{255,0,255}},
{{255,255,255}},
{{196,255,255}},
{{255,255,196}}
};

static uchar bcolors[][3] =
{
{0,0,255},
{0,128,255},
{0,255,255},
{0,255,0},
{255,128,0},
{255,255,0},
{255,0,0},
{255,0,255},
{255,255,255}
};
int main( int argc, char** argv )
{
assert( argc == 2 );
IplImage* img = cvLoadImage( argv[1], CV_LOAD_IMAGE_GRAYSCALE );
IplImage* rsp = cvLoadImage( argv[1], CV_LOAD_IMAGE_COLOR );
CvSeq* contours;
CvMemStorage* storage= cvCreateMemStorage();
IplImage* hsv = cvCreateImage( cvGetSize( rsp ), IPL_DEPTH_8U, 3 );
cvCvtColor( rsp, hsv, CV_BGR2YCrCb );
/*
IplImage* mask = cvCreateImage( cvGetSize( img ), IPL_DEPTH_8U, 1 );
cvZero( mask );
uchar* maskptr = (uchar*)mask->imageData+50*mask->width;
for ( int i = 0; i < mask->width*(mask->height-100); i++ )
{
*maskptr = 1;
maskptr++;
}
*/
double t = cvGetTickCount();
cvExtractMSER( hsv, NULL, &contours, storage, cvMSERParams( 5, 60, cvRound(.2*img->width*img->height), .25, .2 ) );
t = cvGetTickCount() - t;
printf( "MSER extracted %d in %g ms./n", contours->total, t/((double)cvGetTickFrequency()*1000.) );
uchar* rsptr = (uchar*)rsp->imageData;
// draw mser with different color
for ( int i = contours->total-1; i >= 0; i-- )
{
CvSeq* r = *(CvSeq**)cvGetSeqElem( contours, i );
for ( int j = 0; j < r->total; j++ )
{
CvPoint* pt = CV_GET_SEQ_ELEM( CvPoint, r, j );
rsptr[pt->x*3+pt->y*rsp->widthStep] = bcolors[i%9][2];
rsptr[pt->x*3+1+pt->y*rsp->widthStep] = bcolors[i%9][1];
rsptr[pt->x*3+2+pt->y*rsp->widthStep] = bcolors[i%9][0];
}
}
// find ellipse ( it seems cvfitellipse2 have error or sth?
for ( int i = 0; i < contours->total; i++ )
{
CvContour* r = *(CvContour**)cvGetSeqElem( contours, i );
CvBox2D box = cvFitEllipse2( r );
/*
if ( r->color > 0 )
cvEllipseBox( rsp, box, colors[9], 2 );
else
cvEllipseBox( rsp, box, colors[10], 2 );
*/
}
cvSaveImage( "rsp.png", rsp );
cvNamedWindow( "original", 1 );
cvShowImage( "original", img );

cvNamedWindow( "response", 1 );
cvShowImage( "response", rsp );
cvWaitKey(0);
cvDestroyWindow( "original" );
cvDestroyWindow( "response" );

}