2

I'm trying to implement connected component labeling in OpenCV using recursive algorithm. I'm not sure what I have implemented wrongly. the algorithm is like this

B is Binary Image, LB is Labelled Binary Image

procedure connected_component(B,LB)
{
    LB:=negate(B);
    label:=0;
    findComponents(LB,label);
    display(LB);
}

procedure findComponents(LB,label)
{
    for L:=0 to maxRow
        for P:= 0 to maxCol
            if LB[L,P] == -1 then
              {
               label:=label+1;
               search(LB,label,L,P);
              }
}

procedure search(LB,label,L,P)
{
    LB[L,P]:=label;;
    Nset:= neighbours(L,P);
      for each(L',P') in Nset
      {
        if(LB[L',P'] == -1) then
        search(LB,label,L',P');
      }
}

I have written the code in OpenCV as follows

#include<iostream>
#include<opencv2\opencv.hpp>
using namespace cv;
using namespace std;

void findComponents(Mat res, int label);
void search(Mat res, int label, int row, int col);

int main()
{
    Mat src = imread("D:/My Library/test/peppers.bmp",0);
    src.convertTo(src,CV_8S);
    Mat th = src.clone();
    threshold(src,th,128,255,CV_8S);
    Mat res = th.clone();

    for(int i=0;i<res.rows;i++)
        for(int j=0;j<res.cols;j++)
            res.at<signed char>(i,j) = 0 - th.at<signed char>(i,j);

    int label = 0;
    findComponents(res,label);


    waitKey(0);
    return 0;
}

void findComponents(Mat res, int label)
{
    for (int i = 1; i < res.rows - 1; i++)
    {
        for (int j = 1; j < res.cols - 1; j++)
        {
          if (res.at<signed char>(i, j) == -255)
          {
            label++;
            search(res, label, i, j);
          }
        }
    }
    imshow("CC Image", res);

}

void search(Mat res, int label, int row, int col)
{
    res.at<signed char>(row, col) = label;
    if (res.at<signed char>(row, col + 1) == -255) search(res, label, row, col + 1);
    if (res.at<signed char>(row + 1, col + 1) == -255) search(res, label, row+1, col + 1);
    if (res.at<signed char>(row + 1, col) == -255) search(res, label, row + 1, col);
    if (res.at<signed char>(row + 1, col - 1) == -255) search(res, label, row + 1, col - 1);
    else return;

}

The code is does not works. What have I made wrong in implementing the algorithm? I'm new to OpenCV.

Gopiraj
  • 647
  • 7
  • 19

1 Answers1

4

You have a few problems in your code. The most important is that you shouldn't use CV_8S matrices. Why?

  • They have values limited in range [-128, 127]
  • checking for values equal to -255 won't work correctly
  • you are limited to at most 127 connected components per image
  • threshold won't work as expected
  • maybe others...

I re-implemented your code to correct for these issues:

  • you should use CV_32S for your labels.
  • you should account for borders
  • you can use Mat_<Tp> for easy access, instead of .at<Tp>

Below is the code. I used applyCustomColorMap to better visualize results.

#include <opencv2/opencv.hpp>
#include <algorithm>
#include <vector>
#include <stack>
using namespace cv;

void search(Mat1i& LB, int label, int r, int c)
{
    LB(r, c) = label;

    // 4 connected
    if ((r - 1 > 0) && LB(r - 1, c) == -1)       { search(LB, label, r - 1, c    ); }
    if ((r + 1 < LB.rows) && LB(r + 1, c) == -1) { search(LB, label, r + 1, c    ); }
    if ((c - 1 > 0) && LB(r, c - 1) == -1)       { search(LB, label, r    , c - 1); }
    if ((c + 1 < LB.cols) && LB(r, c + 1) == -1) { search(LB, label, r    , c + 1); }

    // 8 connected
    if ((r - 1 > 0) && (c - 1 > 0) && LB(r - 1, c - 1) == -1)             { search(LB, label, r - 1, c - 1); }
    if ((r - 1 > 0) && (c + 1 < LB.cols) && LB(r - 1, c + 1) == -1)       { search(LB, label, r - 1, c + 1); }
    if ((r + 1 < LB.rows) && (c - 1 > 0) && LB(r + 1, c - 1) == -1)       { search(LB, label, r + 1, c - 1); }
    if ((r + 1 < LB.rows) && (c + 1 < LB.cols) && LB(r + 1, c + 1) == -1) { search(LB, label, r + 1, c + 1); }

}

int findComponents(Mat1i& LB)
{
    int label = 0;

    for (int r = 0; r < LB.rows; ++r) {
        for (int c = 0; c < LB.cols; ++c) {
            if (LB(r, c) == -1) {
                ++label;
                search(LB, label, r, c);
            }
        }
    }
    return label;
}

int connected_components(const Mat1b& B, Mat1i& LB)
{
    // Foreground is > 0
    // Background is 0

    LB = Mat1i(B.rows, B.cols, 0);
    LB.setTo(-1, B > 0);

    // Foreground labels are initialized to -1
    // Background labels are initialized to 0

    return findComponents(LB);
}

void applyCustomColormap(const Mat1i& src, Mat3b& dst);

int main()
{
    // Load grayscale image
    Mat1b img = imread("path_to_image", IMREAD_GRAYSCALE);

    // Binarize the image
    Mat1b bin;
    threshold(img, bin, 127, 255, THRESH_BINARY);

    // Find labels
    Mat1i labels;   
    int n_labels = connected_components(bin, labels);

    // Show results
    Mat3b out;
    applyCustomColormap(labels, out);

    imshow("Labels", out);
    waitKey();

    return 0;
}


void applyCustomColormap(const Mat1i& src, Mat3b& dst)
{
    // Create JET colormap

    double m;
    minMaxLoc(src, nullptr, &m);
    m++;

    int n = ceil(m / 4);
    Mat1d u(n * 3 - 1, 1, double(1.0));

    for (int i = 1; i <= n; ++i) {
        u(i - 1) = double(i) / n;
        u((n * 3 - 1) - i) = double(i) / n;
    }

    std::vector<double> g(n * 3 - 1, 1);
    std::vector<double> r(n * 3 - 1, 1);
    std::vector<double> b(n * 3 - 1, 1);
    for (int i = 0; i < g.size(); ++i)
    {
        g[i] = ceil(double(n) / 2) - (int(m) % 4 == 1 ? 1 : 0) + i + 1;
        r[i] = g[i] + n;
        b[i] = g[i] - n;
    }

    g.erase(std::remove_if(g.begin(), g.end(), [m](double v){ return v > m; }), g.end());
    r.erase(std::remove_if(r.begin(), r.end(), [m](double v){ return v > m; }), r.end());
    b.erase(std::remove_if(b.begin(), b.end(), [](double v){ return v < 1.0; }), b.end());

    Mat1d cmap(m, 3, double(0.0));
    for (int i = 0; i < r.size(); ++i) { cmap(int(r[i]) - 1, 0) = u(i); }
    for (int i = 0; i < g.size(); ++i) { cmap(int(g[i]) - 1, 1) = u(i); }
    for (int i = 0; i < b.size(); ++i) { cmap(int(b[i]) - 1, 2) = u(u.rows - b.size() + i); }

    Mat3d cmap3 = cmap.reshape(3);

    Mat3b colormap;
    cmap3.convertTo(colormap, CV_8U, 255.0);

    // Apply color mapping
    dst = Mat3b(src.rows, src.cols, Vec3b(0, 0, 0));
    for (int r = 0; r < src.rows; ++r)
    {
        for (int c = 0; c < src.cols; ++c)
        {
            dst(r, c) = colormap(src(r, c));
        }
    }
}

Please take care that a recursive implementation is not a good idea for labeling:

  • it's quite slow
  • it may fail if you go too deep with recursion, i.e. your components are very big

I suggest to use another algorithm. Here is an implementation of (almost) your algorithm in iterative form. I strongly recommend this one over yours. It can be trivially modified to output the points for each connected component as vector<vector<Point>>, just like cv::findContours would do:

int connected_components2(const Mat1b& img, Mat1i& labels)
{
    Mat1b src = img > 0;
    labels = Mat1i(img.rows, img.cols, 0);

    int label = 0;
    int w = src.cols;
    int h = src.rows;
    int i;

    cv::Point point;
    for (int y = 0; y<h; y++)
    {
        for (int x = 0; x<w; x++)
        {
            if ((src(y, x)) > 0)   // Seed found
            {
                std::stack<int, std::vector<int>> stack2;
                i = x + y*w;
                stack2.push(i);

                // Current component
                std::vector<cv::Point> comp;

                while (!stack2.empty())
                {
                    i = stack2.top();
                    stack2.pop();

                    int x2 = i%w;
                    int y2 = i / w;

                    src(y2, x2) = 0;

                    point.x = x2;
                    point.y = y2;
                    comp.push_back(point);

                    // 4 connected
                    if (x2 > 0 && (src(y2, x2 - 1) != 0))
                    {
                        stack2.push(i - 1);
                        src(y2, x2 - 1) = 0;
                    }
                    if (y2 > 0 && (src(y2 - 1, x2) != 0))
                    {
                        stack2.push(i - w);
                        src(y2 - 1, x2) = 0;
                    }
                    if (y2 < h - 1 && (src(y2 + 1, x2) != 0))
                    {
                        stack2.push(i + w);
                        src(y2 + 1, x2) = 0;
                    }
                    if (x2 < w - 1 && (src(y2, x2 + 1) != 0))
                    {
                        stack2.push(i + 1);
                        src(y2, x2 + 1) = 0;
                    }

                    // 8 connected
                    if (x2 > 0 && y2 > 0 && (src(y2 - 1, x2 - 1) != 0))
                    {
                        stack2.push(i - w - 1);
                        src(y2 - 1, x2 - 1) = 0;
                    }
                    if (x2 > 0 && y2 < h - 1 && (src(y2 + 1, x2 - 1) != 0))
                    {
                        stack2.push(i + w - 1);
                        src(y2 + 1, x2 - 1) = 0;
                    }
                    if (x2 < w - 1 && y2>0 && (src(y2 - 1, x2 + 1) != 0))
                    {
                        stack2.push(i - w + 1);
                        src(y2 - 1, x2 + 1) = 0;
                    }
                    if (x2 < w - 1 && y2 < h - 1 && (src(y2 + 1, x2 + 1) != 0))
                    {
                        stack2.push(i + w + 1);
                        src(y2 + 1, x2 + 1) = 0;
                    }
                }

                ++label;
                for (int k = 0; k <comp.size(); ++k)
                {
                    labels(comp[k]) = label;
                }
            }
        }
    }

    return label;
}
Community
  • 1
  • 1
Miki
  • 40,887
  • 13
  • 123
  • 202
  • In the algorithm I've given, We don't have to scan all the neighbors right? In the search function, we don't have to call the recursion for LB(r - 1, c), LB(r, c - 1), LB(r - 1, c - 1), LB(r - 1, c + 1) as we are scanning from left to right and top down. upon commenting those lines I don't get stack overflow error while using connected_components function. But yeah it is slightly slower compared to connected_components2 for bigger images. – Gopiraj Mar 23 '16 at 03:52
  • No, you need to search also in those directions. Otherwise you need to recover label equivalency using a second pass. That would be a [two-pass algorithm](https://en.wikipedia.org/wiki/Connected-component_labeling#Two-pass) – Miki Mar 23 '16 at 11:58
  • 1
    If you need to go really fast, you could start reading [this](https://www.lri.fr/~cabaret/documents/SIPS2014.pdf) – Miki Mar 23 '16 at 12:03