MSER最穩定極值區域源碼分析

最穩定極值區域介紹

如把灰度圖看成高低起伏的地形圖，其中灰度值看成海平面高度的話，MSER的作用就是在灰度圖中找到符合條件的坑窪。條件為坑的最小高度，坑的大小，坑的傾斜程度，坑中如果已有小坑時大坑與小坑的變化率。

上圖展示了幾種不同的坑窪，根據最小高度，大小，傾斜程度這些條件的不同，選擇的坑也就不同。

 上圖展示了最后一個條件，大坑套小坑的情況。根據條件的不同，選擇也不同。

 以上便是對坑的舉例，MSER主要流程就三部分組成：

　　　　1.預處理數據

　　　　2.遍歷灰度圖

　　　　3.判斷一個區域（坑窪）是否滿足條件

簡單來說，就如將水注入這個地形中。水遇到低處就往低處流，如果沒有低處了，水位就會一點點增長，直至淹沒整個地形。在之前預處理下數據，在水位提高時判斷下是否滿足條件。

預處理數據

 先說下流程中的主要部件，如下：

　　1.img，由原8位單通道灰度圖轉化的更容易遍歷和記錄數據的32位單通道圖。預處理內容為：

　　　　32位值記錄從這點是否探索過，探索過的方向，灰度值；圖大小也擴大了，最外添加了一個像素的完整一圈，值為-1可看作牆，寬度也改變為2的整數次方，用於加快運算。

　　　　如果由掩碼圖，如下：

　　2.heap，記錄坑窪邊界的堆棧，每個灰度值都有自己的堆棧。預處理內容為：

　　　　計算所有灰度值的個數，這樣提前就可以分配堆棧大小。例如知道了灰度2的像素由4個，就可以將灰度2的堆棧大小分配為5（多一個位標志位空）。

　　3.comp，記錄水坑數據的堆棧，有水位值（灰度值），面積（像素個數和像素位置）等。預處理內容為：

　　　　僅僅是分配內存，分配257個（0-255外多一個用作結束）

　　4.history，記錄水位抬高的歷史，就是一個小坑抬高水位后一點點變成大坑的歷史。預處理內容為：

　　　　僅僅是分配內存，大小為像素點個數（就是寬*高）。可以想成所有點都不同都可以形成歷史的最大個數。

遍歷灰度圖

在重復下整個簡單的過程：就如將水注入這個地形中。水遇到低處就往低處流，如果沒有低處了，水位就會一點點增長，直至淹沒整個地形。先說下主要部件：

　　1.img，由原8位單通道灰度圖轉化的更容易遍歷和記錄數據的32位單通道圖。遍歷時：

　　　　當前像素位置中有3位記錄方向（除了東南西北還有一個用來代表結束），逐個改變方向遍歷。還有最高1位記錄是否發現過了。根據方向遍歷相鄰像素，如果4個方向都探索過了，就從heap邊界中找到一個最小灰度的邊界，出棧來用作當前像素。最終將所有像素的4個方向都走完，也是所有像素都被發現了，遍歷就結束。

　　2.heap，記錄坑窪邊界的堆棧，每個灰度值都有自己的堆棧。遍歷時：

　　　　當水遇到低處時入棧當前位置為低處的邊界，當水遇到相等高度或高處時入棧那個邊界；當抬高水位時出棧被抬高到的邊界。

　　3.comp，記錄水坑數據的堆棧，有水位值（灰度值），面積（像素個數和像素位置）等。遍歷時：

　　　　當水位下降時新入棧，水位提高時出棧並可能與之前的合並。

　　4.history，記錄水位抬高的歷史，就是一個小坑抬高水位后一點點變成大坑的歷史。遍歷時：

　　　　history主要是記錄用來判斷最大穩定極值區域的數據，沒有遍歷的作用。主要記錄時刻有兩種：提高水位到邊界heap中的最小高度，提高水位到comp中上一項的高度。要記錄灰度值，像素數，快捷指針，孩子指針，穩定時的像素數。

下面舉例子，走下遍歷的流程（並不是依次就是一步，一些步驟合並了）（紅色為有變動位置，時間匆忙沒有仔細校准每個位置）：

 

中上，要遍歷的灰度圖。為了方便觀看，上文提到周圍一圈的-1被去掉了。

左下，history是抬高水位的歷史。

中下，comp是水位數據。預先入棧一個256的灰度作為頂，用來抬高水位時判斷邊界值小還是上一個水位數據的灰度值小。

右下，heap是邊界堆棧，heap_start是每個灰度指向heap堆棧的指針。特殊說明下，heap是一個個堆棧連接在一起的一個數組，由於上面說的預處理過了，已經知道每個灰度的像素個數，所以提前指定了heap_start中每個灰度指向heap中的位置，指向0代表所在堆棧沒有數據。例如灰度2有4個像素，所以灰度3的指針從灰度2指針后5個后開始，4個是像素數，1個是代表空的0。

從A1位置開始，comp中入棧個灰度2的數據項，並將heap_cur當前指針設置為2灰度的指針。

探索A1右邊B1，標識為已發現。B1的值2沒有小於當前水位值2，作為邊界入棧。

探索A1下面的A2。值1小於當前水位2，將2入棧邊界棧，入棧水位數據1，調整邊界指針heap_cur為指向1的指針，當前像素為A2。

探索A2右邊B3與下邊A3，都沒有比當前水位1小，分別入棧所屬灰度的邊界棧。

A2所有方向都探索完，將A2加入當前水位數據comp中。

在邊界棧中找到最小灰度的一個值出棧（圖5里邊界里有灰度2的和灰度3的，從當前灰度1開始一點點加大所以找到了灰度2），出棧了A3。A3的灰度2，所以抬高水位。記錄歷史histroy，修改當前水位數據灰度為2，邊界指針heap_cur指向2灰度的堆棧。

探索A3周邊，發現B3，灰度3比當前大作為邊界入棧。

A3所有方向也都探索完，將A3加入當前水位數據comp中。

邊界中找到A1。由於A1灰度還是2，沒有提升水位。將A1作為當前像素。

剛剛的A1周圍也早就探索完了，將A1加入當前水位數據comp中。

又在邊界中找到了B1，並出棧作為當前像素。

B1右邊探索到了C1，加入灰度3的邊界棧。

這時，B1周圍已經探索完畢，將B1加入當前水位數據comp中。

B1被加入在邊界棧中從灰度2開始查找，找到灰度3中C1作為當前像素。然后記錄歷史history，提高當前水位數據comp的灰度值，設置heap_cur指針到灰度3的邊界棧。

從當前像素C1向下找到C2，C2灰度比當前低。將當前像素C1入棧邊界棧，新建灰度2的水位數據comp，邊界指針heap_cur指向灰度2，設置C2為當前指針。

探索C2下面最后一個像素C3，將C3加入邊界棧。

將C2加入水位數據comp中。

需要抬高水位了，從灰度3的邊界棧中出棧C3，發現灰度和上一個水位數據comp的灰度一樣，需要合並這兩個comp數據。添加歷史history，合並兩個comp數，改變邊界棧heap_cur到灰度3，設置C3為當前像素。

最后的C3,C1,B3,B2周圍都沒有可以探索的像素了，依次出棧加入水位數據。

至此所有9個像素都探索完畢了。

判斷一個區域（坑窪）是否滿足條件

先看下參數：

  int delta;      // 兩個區域間的灰度差
  int minArea;     // 區域最小像素數
  int maxArea;     // 區域最大像素數
  double maxVariation;  // 兩個區域的偏差
  double minDiversity;  // 當前區域與穩定區域的變化率

一個水坑的變化如下圖A，隨着水位的提高，面積由R_a變為R_b在到R_c，R_a為R_b的父區域；判斷極值區域的方法如圖B，在delta水位差間兩個區域面積是否滿足一定條件；還有一個判斷條件如圖C，如果已經有一個候選區域R_stable了，R_candidate是否可以作為一個極值區域，也就是大坑套小坑的情況。

maxVariation是上圖B的情況，值為下面的公式A；minDiversity是上圖C的情況，值為下面公式B：

下面是在條件判斷時兩個有用的部件（其他沒有任何作用）：

　　3.comp，記錄水坑數據的堆棧，有水位值（灰度值），面積（像素個數和像素位置）等。條件判斷時：

　　　　有個history指向當前區域的歷史的指針，用來查找當前區域之前的變化歷史；var用來記錄上次計算的variation；div用來記錄上次計算的diversity。（var與div用來確保坑越來越穩定，如果與上次的值比較發散了則不滿足條件）

　　4.history，記錄水位抬高的歷史，就是一個小坑抬高水位后一點點變成大坑的歷史。條件判斷時：

　　　　每一個歷史項都有指向孩子歷史的指針child，與指向相差delta灰度歷史的快捷指針shortcut，還有上次穩定時的像素數stable，最后就是那個歷史時刻的灰度值val與像素數size。（快捷指針是用來加速計算的，在歷史里一個一個向前找也能找到，但總沒有直接在上次找到的位置前后找更快吧：））

源碼

 基本結構：

typedef struct LinkedPoint
    {
        struct LinkedPoint* prev;
        struct LinkedPoint* next;
        Point pt;
    }
    LinkedPoint;

    // the history of region grown
    typedef struct MSERGrowHistory
    {
        // 快捷路徑，是指向以前歷史的指針。因為不是一個一個連接的，所以不是parent。算法中是記錄灰度差為delta的歷史的指針。
        // 例如：當前是灰度是10，delta=3，這個指針就指向灰度為7時候的歷史
        struct MSERGrowHistory* shortcut;
        // 指向更新歷史的指針，就是從這個歷史繁衍的新歷史，所以叫孩子
        struct MSERGrowHistory* child;
        // 大於零代表穩定，值是穩定是的像素數。這個值在不停的繼承
        int stable; // when it ever stabled before, record the size
        // 灰度值
        int val;
        // 像素數
        int size;
    }
    MSERGrowHistory;

    typedef struct MSERConnectedComp
    {
        // 像素點鏈的頭
        LinkedPoint* head;
        // 像素點鏈的尾
        LinkedPoint* tail;
        // 區域上次的增長歷史，可以通過找個歷史找到之前的記錄
        MSERGrowHistory* 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)
    }
    MSERConnectedComp;

struct MSERParams
    {
        MSERParams(int _delta, int _minArea, int _maxArea, double _maxVariation,
            double _minDiversity, int _maxEvolution, double _areaThreshold,
            double _minMargin, int _edgeBlurSize)
            : delta(_delta), minArea(_minArea), maxArea(_maxArea), maxVariation(_maxVariation),
            minDiversity(_minDiversity), maxEvolution(_maxEvolution), areaThreshold(_areaThreshold),
            minMargin(_minMargin), edgeBlurSize(_edgeBlurSize)
        {}

        // MSER使用
        int delta;                        // 兩個區域間的灰度差
        int minArea;                    // 區域最小像素數
        int maxArea;                    // 區域最大像素數
        double maxVariation;        // 兩個區域的偏差
        double minDiversity;        // 當前區域與穩定區域的變化率
        // MSCR使用
        int maxEvolution;
        double areaThreshold;
        double minMargin;
        int edgeBlurSize;
    };

 預處理：

// 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
    /** @brief 將所給原單通道灰度圖和掩碼圖 預處理為一張方便遍歷與記錄數據的32位單通道圖像；並且根據像素灰度值分配邊緣棧。
    * 32位格式如下：
    * > 0 可用，< 0 已經被訪問
    * 17~19位用於記錄下一個要探索的方向，5個值
    * 8~11位 用於優化的二值搜索
    * 0~8位用於記錄灰度值
    *@param heap_cur 邊緣棧
    *@param src 原單通道灰度圖
    *@param mask 掩碼圖
    */
    static int* preprocessMSER_8UC1(CvMat* img,
        int*** heap_cur,
        CvMat* src,
        CvMat* mask)
    {
        // 數據有效內容是在img中，由一圈-1包圍着，靠左的區域。也就是被一圈-1的牆包圍着。

        // 原始數據跳轉到下一行的偏移量。
        int srccpt = src->step - src->cols;

        // 跳轉到下一行的偏移量，最后減一是因為，例如：xoooxxx，o是有效數據，x是擴充出來的。偏移量應該是3，就是ooo最
        // 右邊的xxx個數。為了計算，就需要減去ooo最左面的一個x。
        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;

        // 設置第一行為-1
        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++)
            {
                // 最左面設置為-1
                *imgptr = -1;
                imgptr++;
                for (int j = 0; j < src->cols; j++)
                {
                    if (*maskptr)
                    {
                        if (!startptr)
                            startptr = imgptr;

                        // 灰度值取反!!!!! !!!!! !!!!! !!!!!
                        *srcptr = 0xff - *srcptr;

                        // 所在灰度值個數自增
                        level_size[*srcptr]++;

                        // 寫入0~8位，8~13位用作BitScanForward
                        *imgptr = ((*srcptr >> 5) << 8) | (*srcptr);
                    }
                    else {
                        // 標為-1，就是當作一個已經被發現的位置，和外圍-1牆的原理一樣
                        *imgptr = -1;
                    }
                    imgptr++;
                    srcptr++;
                    maskptr++;
                }

                // 最右面設置為-1
                *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;
            }
        }

        // 設置最后一行為-1
        for (int i = 0; i < src->cols + 2; i++)
        {
            *imgptr = -1;
            imgptr++;
        }

        // 確定每個灰度在邊界堆中的指針位置。0代表沒有值。
        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 extractMSER_8UC1_Pass(int* ioptr,
        int* imgptr,
        int*** heap_cur,                            // 邊界棧的堆，里面是每一個灰度的棧
        LinkedPoint* ptsptr,
        MSERGrowHistory* histptr,
        MSERConnectedComp* comptr,
        int step,
        int stepmask,
        int stepgap,
        MSERParams params,
        int color,
        CvSeq* contours,
        CvMemStorage* storage)
    {
        // ER棧第一項為結束的標識項，值為大於255的256
        comptr->grey_level = 256;

        // 將當前位置值入棧，並初始化
        comptr++;
        comptr->grey_level = (*imgptr) & 0xff;
        initMSERComp(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
                        // 如果相鄰像素的灰度小於當前像素，將當前像素加入邊界棧堆，並把相鄰像素設置為當前像素，並新建ER棧項
                        // 將當前加入邊界棧堆
                        (*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;

                        // 新建ER棧項，並設置灰度為當前像素灰度
                        comptr++;
                        initMSERComp(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 imsk = (int)(imgptr - ioptr);

            // 記錄x&y，
            ptsptr->pt = cvPoint(imsk&stepmask, imsk >> stepgap);
            // get the current location
            accumulateMSERComp(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)
                    {
                        // 剛從邊界獲得灰度如果小於上一個MSER組件灰度值，需要提高當前水位到邊界的灰度值
                        // check the stablity and push a new history, increase the grey level
                        if (MSERStableCheck(comptr, params))
                        {
                            CvContour* contour = MSERToContour(comptr, storage);
                            contour->color = color;
                            cvSeqPush(contours, &contour);
                        }

                        // 由於水位要有變化了，添加一個歷史
                        MSERNewHistory(comptr, histptr);

                        // 提高水位到邊界的水位
                        comptr[0].grey_level = pixel_val;

                        // 指向下一個未使用歷史空間
                        histptr++;
                    }
                    else {
                        // 剛從邊界獲得灰度如果不小於上一個MSER組件灰度值，其實就是和上一個灰度值一樣。
                        // 例如：當前水位2，上一個水位3，從邊界出棧的水位為3.

                        // keep merging top two comp in stack until the grey level >= pixel_val
                        for (;;)
                        {
                            // 合並MSER組件，里面也隨帶完成了一個歷史
                            comptr--;
                            MSERMergeComp(comptr + 1, comptr, comptr, histptr);
                            histptr++;

                            if (pixel_val <= comptr[0].grey_level)
                                break;

                            // 到這里，等於comptr[0].grey_level < pixel_val，也是當前像素的灰度與MSER組件的不一致，要提高MSER組件灰度
                            if (pixel_val < comptr[-1].grey_level)
                            {
                                // 其實就是comptr[0].grey_level < pixel_val < comptr[-1].grey_level
                                // 當前灰度大於當前MSER灰度小於上一個MSER組件灰度。同上面的代碼情況一樣。
                                // check the stablity here otherwise it wouldn't be an ER
                                if (MSERStableCheck(comptr, params))
                                {
                                    CvContour* contour = MSERToContour(comptr, storage);
                                    contour->color = color;
                                    cvSeqPush(contours, &contour);
                                }
                                MSERNewHistory(comptr, histptr);
                                comptr[0].grey_level = pixel_val;
                                histptr++;
                                break;
                            }
                        }
                    }
                }
                else
                    break;
            }
        }
    }

    /** @brief 通過8UC1類型的圖像提取MSER
    *@param mask 掩碼
    *@param contours 輪廓結果
    *@param storage 輪廓內存空間
    *@param params 參數
    */
    static void extractMSER_8UC1(CvMat* src,
        CvMat* mask,
        CvSeq* contours,
        CvMemStorage* storage,
        MSERParams params)
    {
        // 為了加速計算，將每行數據大小擴展為大於原大小的第一個2的整指數。
        // 這樣在后面計算y時，只要右移stepgap就算除以2^stepgap了
        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;                                        // 用於指向mser遍歷的當前像素（所有數據）

        // 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
        // 預分配連接像素點，用於將區域中的像素連接起來，大小就為所有像素個數
        LinkedPoint* pts = (LinkedPoint*)cvAlloc(src->rows*src->cols * sizeof(pts[0]));
        // 預分配增長歷史，用於記錄區域在太高水位后的父子關系，最大個數為所有像素個數。
        MSERGrowHistory* history = (MSERGrowHistory*)cvAlloc(src->rows*src->cols * sizeof(history[0]));
        // 預分配區域，用於記錄每個區域的數據，大小為所有灰度值+1個超大灰度值代表頂
        MSERConnectedComp comp[257];

        // darker to brighter (MSER-)
        // 提取mser亮區域（preprocessMSER_8UC1中將灰度值取反）
        imgptr = preprocessMSER_8UC1(img, heap_start, src, mask);
        extractMSER_8UC1_Pass(ioptr, imgptr, heap_start, pts, history, comp, step, stepmask, stepgap, params, -1, contours, storage);
        // brighter to darker (MSER+)
        // 提取mser暗區域
        imgptr = preprocessMSER_8UC1(img, heap_start, src, mask);
        extractMSER_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);
    }

條件判斷和生成結果：

// clear the connected component in stack
    static void
        initMSERComp(MSERConnectedComp* comp)
    {
        comp->size = 0;
        comp->var = 0;
        comp->dvar = 1;
        comp->history = NULL;
    }

    // add history of size to a connected component
    static void
        /** @brief 通過當前ER項構建一個對應的歷史，也就是說找個ER項要准備改變了
        */
        MSERNewHistory(MSERConnectedComp* comp, MSERGrowHistory* history)
    {
        // 初始時將下一條歷史設置為自己
        history->child = history;
        if (NULL == comp->history)
        {
            // 從來沒有歷史過，將快捷路徑也設置為自己，穩定的像素數為0
            history->shortcut = history;
            history->stable = 0;
        }
        else {
            // 有歷史，將當前歷史設置為上一個歷史的下個歷史
            comp->history->child = history;

            // 快捷路徑與穩定值繼承至上一個歷史
            history->shortcut = comp->history->shortcut;
            history->stable = comp->history->stable;
        }

        // 記錄這時的ER項的灰度值與像素數
        history->val = comp->grey_level;
        history->size = comp->size;

        // 設置ER項的歷史為找個最新的歷史
        comp->history = history;
    }

    // merging two connected component
    static void
        MSERMergeComp(MSERConnectedComp* comp1,
            MSERConnectedComp* comp2,
            MSERConnectedComp* comp,
            MSERGrowHistory* history)
    {
        LinkedPoint* head;
        LinkedPoint* 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;
            }

            // 如果組件2有stable，並且大於1的，則stable使用2的值
            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;

            // 如果組件1和2都有像素點，將兩個鏈按照1->2連接在一起
            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;

        // 新ER的像素數量是兩個ER項的和
        comp->size = comp1->size + comp2->size;
    }

    /** @brief 通過delta計算給定ER項的偏差
    */
    static float MSERVariationCalc(MSERConnectedComp* comp, int delta)
    {
        MSERGrowHistory* history = comp->history;
        int val = comp->grey_level;
        if (NULL != history)
        {
            // 從快捷路徑開始往回找歷史，找到灰度差大於delta的歷史
            MSERGrowHistory* shortcut = history->shortcut;
            while (shortcut != shortcut->shortcut && shortcut->val + delta > val)
                shortcut = shortcut->shortcut;

            // 由於快捷路徑是直接跳過一些歷史的，要找到最准確的歷史還要從以前歷史往當前找
            MSERGrowHistory* 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;

            // 返回(R-R(-delta)) / (R-delta)
            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
        }

        // 沒有歷史，結果為1。也就是沒有-delta對應的值。
        // 如果按照(R-R(-delta)) / R(-delta) = 1公式推導:
        // R = 2R(-delta)
        // 就面積來說，怎么兩倍這種關系都比較奇怪，因為是xy兩個維度的，每個維度提高sqrt(2)倍
        return 1.;
    }

    /** @brief 檢查是否為最穩定極值區域
    */
    static bool MSERStableCheck(MSERConnectedComp* comp, MSERParams params)
    {
        // 檢查就是要確定水位的底是否是穩定的
        // tricky part: it actually check the stablity of one-step back
        // 穩定區域都是由比較而來的，不能沒有上一個歷史。
        if (comp->history == NULL || comp->history->size <= params.minArea || comp->history->size >= params.maxArea)
            return 0;

        // diversity : (R(-1) - R(stable)) / R(-1)
        // 使用水位的底與穩定時大小做比較
        float div = (float)(comp->history->size - comp->history->stable) / (float)comp->history->size;

        // variation
        float var = MSERVariationCalc(comp, params.delta);

        // 現在的variation要大於以前的variation，就是以前的更穩定
        // 灰度值差是否大於1
        int dvar = (comp->var < var || (unsigned long)(comp->history->val + 1) < comp->grey_level);
        int stable = (dvar && !comp->dvar && comp->var < params.maxVariation && div > params.minDiversity);
        comp->var = var;
        comp->dvar = dvar;
        if (stable)
            // 如果穩定的話，穩定值就是像素數
            comp->history->stable = comp->history->size;
        return stable != 0;
    }

    // add a pixel to the pixel list
    /** @brief 添加像素到給定的MSER項中
    */
    static void accumulateMSERComp(MSERConnectedComp* comp, LinkedPoint* 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
    static CvContour* MSERToContour(MSERConnectedComp* 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);
        LinkedPoint* 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;
    }

p.s. 以上代碼有點長了：（

應用 

下面對一些圖片做實驗，測試下mser的檢出能力。

    // 加載圖像
    Mat srcColor = imread("");

    //創建MSER類  
    MSER ms(4                    // delta
        , 60                        // min area
        , 1600                    // max area
        , 0.05f                    // max variation 
        , 0.4f                        // min diversity 
        );                            // edge blur size 

    // 轉換為灰度圖
    Mat srcGray;
    cvtColor(srcColor, srcGray, CV_BGR2GRAY);

    //用於組塊區域的像素點集  
    vector<vector<Point>> regions;
    ms(srcGray, regions, Mat());

    for (int i = 0; i < regions.size(); i++)
    {
        //用連線繪制  
        ////polylines(srcGray, regions[i], true, Scalar(0, 0, 255));
        //用橢圓形繪制  
        ellipse(srcColor, fitEllipse(regions[i]), Scalar(0, 0, 255), 1);
    }

可以看出mser對旋轉和不同大小的字符都有一定的檢出能力，但如果想對不同灰度變化也有能力，應該修改源碼來適應了。

下次介紹mscr，用於在彩色圖種查找穩定區域。

 

轉載請注明出處，謝謝~