Friday, May 3, 2013

Hough Transformation - C++ Implementation


fence - example 1
  
accumulator
flower stem - example 2

The Hough Transformation is a great way to detect lines in an image and it is quite useful for a number of Computer Vision tasks. 

A good definition can be found at http://en.wikipedia.org/wiki/Hough_transform

This post will guide you through (a little) theory and focusses on the C++ implementation of the linear Hough Transformation.


Theory

Let's take an image (Fig 1) with two lines A and B. Obviously both lines are each made of its own set of pixels laying on a straight line. Now, one way or another we need to learn our software which pixels are on a straight line and, if so, to what line they belong to.


Fig 1
Fig 2

The trick of the Hough Transformation is to represent lines in a polar coordinate system (see Fig 2), or in this case, also called "Hough Space". 

You can read this wiki for more details but important in Fig 2 is that we imagine a line from the center (the blue line) which is perpendicular to the line we want to find (the green line) and we call the distance from the center r and θ(theta) is the angle to the X axis (the angle between the blue line and the X axis).

The transformation from any (x, y) pixel in Fig1 to the (r, θ) representation int Fig2 is the following:

r = x.cos(θ) + y.sin(θ)

Now, if you lean back and think a bit about Fig2. Every possible line can be represented by a unique r and θ. And further, every pixel on a given line will transform to the exact same r and θ for that line. And this is immediately the most important aspect of the whole Hough stuff. If you get that, you understand the transformation.

So, given the above; every line in polar space is defined by its unique coordinates (r, θ). Let play with this and say we decide to use  θ ∈ { 1 ... 180 } with a resolution of 1 degree. 1 to 180 allow us to make any line if we allow r to be negative. Check Fig2, if we use the same θ but make r negative that line A' would be at the oposite side (left) and parallel to A.

That means that for every point (x,y) we can compute 180 r values.
(x1,y1)(x2,y2)
r1_1 = x1.cos(1) + y1.sin(1)
r1_2 = x1.cos(2) + y1.sin(2)
r1_3 = x1.cos(3) + y1.sin(3)
...
r1_180 = x1.cos(180) + y1.sin(180)
r2_1 = x2.cos(1) + y2.sin(1)
r2_2 = x2.cos(2) + y2.sin(2)
r2_3 = x2.cos(3) + y2.sin(3)
...
r2_180 = x2.cos(180) + y2.sin(180)


So far, all fine, but to make it really useful we need to add one extra element and that is called an Accumulator. An Accumulator is nothing more than a 2 dimensional matrix with θmax columns and rmax rows. So each cell, or bin, represents a unique coordinate (r, θ) and thus one unique line in Hough Space. 

The algorithm is as follow:

  • Start with an accumulator with all 0's
  • For every point in your image
    • compute r for θ = 1 to 180
    • increment the Accumulator by 1 for the every found r, so increment bin(r, θ)

Fig 3 - Illustration 















Once done, the accumulator bins contain a count of pixels for every possible line. Go over the binss and only keep the ones that are over a certain threshold (for example, you can decide to only consider it a real line if it contains at least 80 pixels). 

The result is a set of coordinates (r, θ), one for every straight line in the image thats above the threshold. And these polar coordinates easily can be computed back to 2 points in the original image space, just solve the equation for x and y. 



         




C++ Implementation

The following code is a command line implementation that take an image path as input (png, jpg, ...) and displays a view with the detected lines drawn over the image. It also show a view of the accumulator, which can look pretty neat.

The procedure is as follow
  1. Load the color image
  2. Do Canny edge detection on it  so we have a BW image with lines and points representing the edges
  3. Do a hough transformation on the edge image
  4. Search the accumulator for all lines given a certain threshold
  5. Visualize all on the screen 
Since we focus on the implementation of step 3 & 4 (Hough) we are going to use some open source for step 1, 2 and 5. I decided to use OpenCV (http://opencv.org), its open source, popular, awesome and pretty easy to use. 

OpenCV itself has a Hough Transformation, but we are not going to use that, our goal is to get insight in the inner working of the transformation and implement it ourselves. Also I think OpenCV implemented a Probabilistic or Randomized Hough Transform, wich is a faster but different approach http://en.wikipedia.org/wiki/Randomized_Hough_transform'.



1, 2 - Load the image and do the edge detection


Let's work with the following image:

Original Image

1:       cv::Mat img_edge;  
2:       cv::Mat img_blur;  
3:    
4:       cv::Mat img_ori = cv::imread( file_path, 1 );  
5:       cv::blur( img_ori, img_blur, cv::Size(5,5) );  
6:       cv::Canny(img_blur, img_edge, 100, 150, 3);  
7:    
8:       int w = img_edge.cols;  
9:       int h = img_edge.rows;  

This is all OpenCV stuff, line 4 loads the image from path 'file_path'. Line 5 blurs it a bit and line 6 does the actual Edge detection. The result is a Black and White image 'img_edge' with the edges.

Canny Edge Detection















3 - The Hough Transformation

The following code does the actual transformation, img_data contains the edge data. The function loops over all pixels in the edge image and increments the accumulator at  the computed (r, θ).

1:       int Hough::Transform(unsigned char* img_data, int w, int h)  
2:       {  
3:            _img_w = w;  
4:            _img_h = h;  
5:    
6:            //Create the accu  
7:            double hough_h = ((sqrt(2.0) * (double)(h>w?h:w)) / 2.0);  
8:            _accu_h = hough_h * 2.0; // -r -> +r  
9:            _accu_w = 180;  
10:    
11:            _accu = (unsigned int*)calloc(_accu_h * _accu_w, sizeof(unsigned int));  
12:    
13:            double center_x = w/2;  
14:            double center_y = h/2;  
15:    
16:    
17:            for(int y=0;y<h;y++)  
18:            {  
19:                 for(int x=0;x<w;x++)  
20:                 {  
21:                      if( img_data[ (y*w) + x] > 250 )  
22:                      {  
23:                           for(int t=0;t<180;t++)  
24:                           {  
25:                                double r = ( ((double)x - center_x) * cos((double)t * DEG2RAD)) + (((double)y - center_y) * sin((double)t * DEG2RAD));  
26:                                _accu[ (int)((round(r + hough_h) * 180.0)) + t]++;  
27:                           }  
28:                      }  
29:                 }  
30:            }  
31:    
32:            return 0;  
33:       }  
34:    

Line 11 creates the empty accumulator, we have chosen to compute r for θ = {1..180} so the width contain 180 bins, the maximum height depends on the image size and is computed at line 7.

Line 17 to 30 loops over every point in the edge image and computes r for every 
θ at line 25, you should see r = x.cos(θ) + y.sin(θ) in that. Little extra math is added to work from the center of the image. The found accumulator bin is incremented at line 26.


4. Search the accumulator for all lines given a certain threshold

A few steps are done here; we loop over the accumulator and consider every bin which value is on or above the threshold. Than we check if that bin is a local maxima. If we decide the line at coordinates (r, θ) is a valid one we compute them back to two points in the image space.


1:       std::vector< std::pair< std::pair<int, int>, std::pair<int, int> > > Hough::GetLines(int threshold)  
2:       {  
3:            std::vector< std::pair< std::pair<int, int>, std::pair<int, int> > > lines;  
4:    
5:            if(_accu == 0)  
6:                 return lines;  
7:    
8:            for(int r=0;r<_accu_h;r++)  
9:            {  
10:                 for(int t=0;t<_accu_w;t++)  
11:                 {  
12:                      if((int)_accu[(r*_accu_w) + t] >= threshold)  
13:                      {  
14:                           //Is this point a local maxima (9x9)  
15:                           int max = _accu[(r*_accu_w) + t];  
16:                           for(int ly=-4;ly<=4;ly++)  
17:                           {  
18:                                for(int lx=-4;lx<=4;lx++)  
19:                                {  
20:                                     if( (ly+r>=0 && ly+r<_accu_h) && (lx+t>=0 && lx+t<_accu_w) )  
21:                                     {  
22:                                          if( (int)_accu[( (r+ly)*_accu_w) + (t+lx)] > max )  
23:                                          {  
24:                                               max = _accu[( (r+ly)*_accu_w) + (t+lx)];  
25:                                               ly = lx = 5;  
26:                                          }  
27:                                     }  
28:                                }  
29:                           }  
30:                           if(max > (int)_accu[(r*_accu_w) + t])  
31:                                continue;  
32:    
33:    
34:                           int x1, y1, x2, y2;  
35:                           x1 = y1 = x2 = y2 = 0;  
36:    
37:                           if(t >= 45 && t <= 135)  
38:                           {  
39:                                //y = (r - x cos(t)) / sin(t)  
40:                                x1 = 0;  
41:                                y1 = ((double)(r-(_accu_h/2)) - ((x1 - (_img_w/2) ) * cos(t * DEG2RAD))) / sin(t * DEG2RAD) + (_img_h / 2);  
42:                                x2 = _img_w - 0;  
43:                                y2 = ((double)(r-(_accu_h/2)) - ((x2 - (_img_w/2) ) * cos(t * DEG2RAD))) / sin(t * DEG2RAD) + (_img_h / 2);  
44:                           }  
45:                           else  
46:                           {  
47:                                //x = (r - y sin(t)) / cos(t);  
48:                                y1 = 0;  
49:                                x1 = ((double)(r-(_accu_h/2)) - ((y1 - (_img_h/2) ) * sin(t * DEG2RAD))) / cos(t * DEG2RAD) + (_img_w / 2);  
50:                                y2 = _img_h - 0;  
51:                                x2 = ((double)(r-(_accu_h/2)) - ((y2 - (_img_h/2) ) * sin(t * DEG2RAD))) / cos(t * DEG2RAD) + (_img_w / 2);  
52:                           }  
53:    
54:                           lines.push_back(std::pair< std::pair<int, int>, std::pair<int, int> >(std::pair<int, int>(x1,y1), std::pair<int, int>(x2,y2)));  
55:    
56:                      }  
57:                 }  
58:            }  
59:    
60:            std::cout << "lines: " << lines.size() << " " << threshold << std::endl;  
61:            return lines;  
62:       }  


We loop over the accumulator at line 8 - 10, we check if the value of a bin is on or above our threshold at line 12, if so we check if this point is a local maxima at line 14 - 32. If we have a valid line the (r, θ) coordinates are computed back to two points in image space at lines 34 - 52. The results are store in a vector at line 54




5. Visualize all on the screen

1:            //Search the accumulator  
2:            std::vector< std::pair< std::pair<int, int>, std::pair<int, int> > > lines = hough.GetLines(threshold);  
3:    
4:            //Draw the results  
5:            std::vector< std::pair< std::pair<int, int>, std::pair<int, int> > >::iterator it;  
6:            for(it=lines.begin();it!=lines.end();it++)  
7:            {  
8:                 cv::line(img_res, cv::Point(it->first.first, it->first.second), cv::Point(it->second.first, it->second.second), cv::Scalar( 0, 0, 255), 2, 8);  
9:            }  
10:    
11:            //Visualize all  
12:            int aw, ah, maxa;  
13:            aw = ah = maxa = 0;  
14:            const unsigned int* accu = hough.GetAccu(&aw, &ah);  
15:    
16:            for(int p=0;p<(ah*aw);p++)  
17:            {  
18:                 if(accu[p] > maxa)  
19:                      maxa = accu[p];  
20:            }  
21:            double contrast = 1.0;  
22:            double coef = 255.0 / (double)maxa * contrast;  
23:    
24:            cv::Mat img_accu(ah, aw, CV_8UC3);  
25:            for(int p=0;p<(ah*aw);p++)  
26:            {  
27:                 unsigned char c = (double)accu[p] * coef < 255.0 ? (double)accu[p] * coef : 255.0;  
28:                 img_accu.data[(p*3)+0] = 255;  
29:                 img_accu.data[(p*3)+1] = 255-c;  
30:                 img_accu.data[(p*3)+2] = 255-c;  
31:            }  
32:    
33:    
34:            cv::imshow(CW_IMG_ORIGINAL, img_res);  
35:            cv::imshow(CW_IMG_EDGE, img_edge);  
36:            cv::imshow(CW_ACCUMULATOR, img_accu);  
37:    


After we searched the accumulator at line 2 an OpenCV function is used to draw lines over the original image at line 5 - 9. Line 12 - 31 builds a visual representation of the accumulator. Finally line 34 - 36 displays all on screen.

The result with a threshold of 175 is the following:
accumulator
result



Sources

Full source code can be found at: https://github.com/brunokeymolen/hough

Make sure to install OpenCV (I used 2.4.4) prior to compile; http://opencv.org/




46 comments:

  1. Dear actually I am new to computer vision. I want to implement Hough Transform for line detection I tried your but it is not working properly. I am using visual studio 2010. Can you please help me in this regards??

    ReplyDelete
  2. it is giving error like ''Hough" is not class or name space
    "getopt" identifier not found

    ReplyDelete
  3. This comment has been removed by the author.

    ReplyDelete
  4. thank you so much. I will try it.
    Best regards

    ReplyDelete
  5. Hi! Thanks for the code, im running it on my VC2010.
    I have a question: It is possible to tweak the code in order to detect circles?
    Thanks!

    ReplyDelete
    Replies
    1. Yes it is; see for example, 4.2 in this ppt:
      http://ehm.kocaeli.edu.tr/duyuru/dosyalar/1058/Hough_Transform_SE.ppt

      In the Transform in hough.cpp it would look something like this:

      https://dl.dropboxusercontent.com/u/23324693/circle_hough.png

      Delete
    2. Thanks for your answer! now im getting this error:
      error: C2065: '_r' : undeclared identifier

      How should I initialize it?

      Delete
    3. Hi,

      It is the circle radius, set it to the radius you are looking for.

      There are a bit more changes to do to get it working (like draw a circle instead of a line based on the values in the accumulator).

      It's more meant to illustrate that the core of it is to replace the equations.

      Delete
    4. Thanks, im going to try to figure it out. If you are planning to implement this soon, please let me know, thanks for your work and time.

      Delete
    5. Can you any body share the code for hough circle transform with accumulator voting value for each pixel???Thanks @Bruno Keymolen @Bruno Keymolen

      Delete
  6. Hi,

    I tried to adapt this code to a c# project and have understanding problems. I have two questions:
    1) What does the following code geometrically mean?

    7: double hough_h = ((sqrt(2.0) * (double)(h>w?h:w)) / 2.0);
    8: _accu_h = hough_h * 2.0; // -r -> +r

    2) Could you please explain shortly how you determine the index here? Can't understand what you do here, honestly :(

    _accu[ (int)((round(r + hough_h) * 180.0)) + t]++;

    Thanks for sharing this approach!

    ReplyDelete
    Replies
    1. I am looking for it too. Is it possible to get an explanation?

      Delete
    2. it is like consider your image with height h and width w.
      then maximum distance will be designated by diagonal of the image.
      then r = sqrt(h*h + w*w)
      provided h>w, we can replace it
      r < sqrt (h*h + w*w)
      r < sqrt(2) h

      so this value has been taken for accumulator allocation.

      Delete
    3. How about the accu ? Why multiple by 180.0?

      Delete
    4. Markus;
      To multiple by 180.0 is required because it is the width of the accumulator. You can find the exact index of the r in the array by doing this.

      Delete
  7. Thank you for your brief explanation. Could you explain what are lx and ly variables and why are they in [-4, 4] range? This is in part 4 Search the accumulator for all lines given a certain threshold.

    ReplyDelete
    Replies
    1. We need a local maxima, otherwise there would be noise because other points in the accumulator close to the maxima would reach the threshold too.
      The -4,4 is an easy way to do this, the values are arbitrary.

      Delete
    2. Thanks for your quick reply. So basically I can use any value I want?

      Delete
    3. Ivana;
      That -4 to 4 statement is just a neighboor not more. Remember from the image processing in blurring.

      Delete
  8. Is it possible to explain this line of code. thanks

    double r = ( ((double)x - center_x) * cos((double)t * DEG2RAD)) + (((double)y - center_y) * sin((double)t * DEG2RAD));

    ReplyDelete
  9. It's this equation where r(adius) is computed for every t(heta) : r = x.cos(θ) + y.sin(θ)

    ReplyDelete
  10. hi! your review is really really awesome - it's even BETTER than textbook e.g. http://adaptiveart.eecs.umich.edu/2011/wp-content/uploads/2011/09/The-pocket-handbook-of-image-processing-algorithms-in-C.pdf

    implementing it in plain labview, i had to make one adjustment regarding to de-houging from polar space:
    http://asset-9.soupcdn.com/asset/12913/6779_99ec.jpeg

    ReplyDelete
  11. Why u multiple the _accu (round( r+hough_h)) by 180?

    ReplyDelete
    Replies
    1. Because he uses 1d array like 2d array and the 180 stands for width of 2d array.
      round( r+hough_h) stands for row number and t for column number.
      You can read more about this here: http://stackoverflow.com/questions/2151084/map-a-2d-array-onto-a-1d-array-c

      Delete
  12. Hi Bruno, very nice explanation and helped me a lot.
    Although, finding local maxima can be improved by checking also equals values, not only greater.
    It can be done by introduction new list that will contain ignored indexes of _accu array.
    Then, you need to add another condition, that checks if (int)_accu[( (r+ly)*_accu_w) + (t+lx)] equals max and if so, add the point to ignored array (of course if lx != 0 and ly != 0 and the point has not been added to ignored list before).
    To complete the task you have to modify this line (number 12):
    if((int)_accu[(r*_accu_w) + t] >= threshold) by adding another condition that checks if ignored list not contains ((r*_accu_w) + t) index.
    I hope you get it.
    I wrote this i C# so, but if you'd like to I can share or send it directly to you.

    ReplyDelete
    Replies
    1. Hi Maciej, thanks for the reply! If you have anything to share, please do, it can be helpful for other people.

      Delete
    2. Ok, no problem. Here is my implementation od GetCircles method in C#. Just FYI: A is 2d array containing (t, r) pairs and _image is a processed image.

      public List GetLines(int threshold)
      {
      List lines = new List();
      List ignored = new List();

      for(int theta = 0; theta < A.GetLength(0); theta++)
      {
      for(int r = 0; r < A.GetLength(1); r++)
      {
      if(A[theta, r] >= threshold && !ignored.Exists(s => s.Theta == theta && s.R == r))
      {
      //Is this point a local maxima (9x9)
      int max = A[theta, r];
      bool neighbour = false;
      for(int lx = -4; lx <= 4; lx++)
      {
      for(int ly = -4; ly <= 4; ly++)
      {
      if(((lx + theta) >= 0 && (lx + theta) < A.GetLength(0)) && ((ly + r) >= 0 && (ly + r) < A.GetLength(1)))
      {
      if(A[(theta + lx), (r + ly)] > max)
      {
      max = A[(theta + lx), (r + ly)];
      lx = ly = 5;
      neighbour = true;
      }
      else if (A[(theta + lx), (r + ly)] == max)
      {
      if (lx == 0 && ly == 0) continue;
      if (!ignored.Exists(s => s.Theta == (theta + lx) && s.R == (r + ly)))
      {
      ignored.Add(new PolarCoordinate(theta + lx, r + ly));
      }
      }
      }
      }
      }
      if(neighbour)
      continue;

      int x1, y1, x2, y2;

      if(theta >= 45 && theta <= 135)
      {
      //y = (r - x cos(t)) / sin(t)
      x1 = 0;
      y1 = (int) (((r - (A.GetLength(1)/2)) - ((x1 - (_image.Width/2) ) * MathUtil.CosinusTable[theta])) / MathUtil.SinusTable[theta] + (_image.Height/ 2));
      x2 = _image.Width - 0;
      y2 = (int) (((r - (A.GetLength(1)/2)) - ((x2 - (_image.Width/2) ) * MathUtil.CosinusTable[theta])) / MathUtil.SinusTable[theta] + (_image.Height/ 2));
      }
      else
      {
      //x = (r - y sin(t)) / cos(t);
      y1 = 0;
      x1 = (int) (((r - (A.GetLength(1)/2)) - ((y1 - (_image.Height/2) ) * MathUtil.SinusTable[theta])) / MathUtil.CosinusTable[theta] + (_image.Width / 2));
      y2 = _image.Height - 0;
      x2 = (int) (((r - (A.GetLength(1)/2)) - ((y2 - (_image.Height/2)) * MathUtil.SinusTable[theta])) / MathUtil.CosinusTable[theta] + (_image.Width / 2));
      }

      lines.Add(new LinePoints(new Point(x1, y1), new Point(x2, y2)));
      }
      }
      }

      return lines;
      }

      And below the LinePoints structure:
      public struct LinePoints
      {
      public LinePoints(Point p1, Point p2)
      : this()
      {
      P1 = p1;
      P2 = p2;
      }

      public Point P1 { get; set; }
      public Point P2 { get; set; }
      }

      Delete
  13. Thanks a ton. I understand it a lot easier than any other article around.

    ReplyDelete
  14. Hi Bruno! Could you, please, explain the condition "if(t >= 45 && t <= 135)" you put while calculating back x and y points

    ReplyDelete
  15. This comment has been removed by the author.

    ReplyDelete
  16. Thank you so much for sharing, i wanna ask you something.

    When i try to compile your code, i get error in "getopt" and "optarg" function. I could not find how i can fix it

    ReplyDelete
  17. I wanna ask one more question to you, in your code i also get "unistd.h" could not find error. How can i fix it.

    ReplyDelete
  18. This comment has been removed by the author.

    ReplyDelete
  19. can you share the code for hough circle to see the voting value or accumulator image?

    ReplyDelete
  20. Can you any body share the code for hough circle transform with accumulator voting value for each pixel???Thanks @Bruno Keymolen @Bruno Keymolen

    ReplyDelete
  21. Adnan, the code for a circle transform is here:

    https://github.com/brunokeymolen/hough/tree/master/tryout/hough_circle

    ReplyDelete
  22. Can someone please help me?
    I tried to implement this to find fences in my game and I'm getting strange results. My "images" are actually maps for a game, for which 2's and 3's mean the "pixel" represents a location that contains a piece of fence.
    Here is the code I came up with, based on your post, but to suit my needs:

    #define DEG2RAD (M_PI/180.0f)

    void detect_runs(int _world=-1) {
    int SW=0,EW=MAXWORLDS-1;
    if(_world>-1) { SW=_world; EW=_world; }
    for(int world=SW;world<=EW;world++) {
    int SIZ=nybWid[world];

    // std::vector< std::pair< std::pair, std::pair > > Hough::GetLines(int threshold)
    int threshold=6;

    //Create the accu
    double hough_h = ((sqrt(2.0) * (double)(SIZ)) / 2.0);
    int _accu_h = hough_h * 2.0; // -r -> +r
    int _accu_w = 180;
    unsigned int *_accu = (unsigned int*)calloc(_accu_h * _accu_w, sizeof(unsigned int));
    double center_x = SIZ/2;
    double center_y = SIZ/2;

    for(int y=0;y= threshold) {
    //Is this point a local maxima (9x9)
    int max = _accu[(r*_accu_w) + t];
    for(int ly=-4;ly<=4;ly++) {
    for(int lx=-4;lx<=4;lx++) {
    if( (ly+r>=0 && ly+r<_accu_h) && (lx+t>=0 && lx+t<_accu_w) ) {
    if( (int)_accu[( (r+ly)*_accu_w) + (t+lx)] > max ) {
    max = _accu[( (r+ly)*_accu_w) + (t+lx)];
    ly = lx = 5;
    }
    }
    }
    }
    if(max > (int)_accu[(r*_accu_w) + t]) continue;

    int x1, y1, x2, y2;
    x1 = y1 = x2 = y2 = 0;

    if(t >= 45 && t <= 135) {
    //y = (r - x cos(t)) / sin(t)
    x1 = 0;
    y1 = ((double)(r-(_accu_h/2)) - ((x1 - (SIZ/2) ) * cos(t * DEG2RAD))) / sin(t * DEG2RAD) + (SIZ / 2);
    x2 = SIZ - 0;
    y2 = ((double)(r-(_accu_h/2)) - ((x2 - (SIZ/2) ) * cos(t * DEG2RAD))) / sin(t * DEG2RAD) + (SIZ / 2);
    } else {
    //x = (r - y sin(t)) / cos(t);
    y1 = 0;
    x1 = ((double)(r-(_accu_h/2)) - ((y1 - (SIZ/2) ) * sin(t * DEG2RAD))) / cos(t * DEG2RAD) + (SIZ / 2);
    y2 = SIZ - 0;
    x2 = ((double)(r-(_accu_h/2)) - ((y2 - (SIZ/2) ) * sin(t * DEG2RAD))) / cos(t * DEG2RAD) + (SIZ / 2);
    }

    // lines.push_back(std::pair< std::pair, std::pair >(std::pair(x1,y1), std::pair(x2,y2)));
    printf("line: %d,%d - %d,%d\n",x1,y1,x2,y2);

    }
    }
    }
    free(_accu);
    }
    }

    (Output in next post.)

    ReplyDelete
    Replies
    1. There should be a single diagonal line starting around 1278, 1216 that is about 40 "pixels" or so, but here is the result I got: (as you can see, there are NEGATIVE co-ordinates and lots of zeroes. I can see where it sets x or y to 0 in the code but I cannot figure out why). PLEASE help :)

      line: 607,0 - 1742,2048
      line: 584,0 - 1766,2048
      line: 556,0 - 1787,2048
      line: 552,0 - 1782,2048
      line: 526,0 - 1806,2048
      line: 549,0 - 1780,2048
      line: 496,0 - 1826,2048
      line: 459,0 - 1841,2048
      line: 365,0 - 1908,2048
      line: 330,0 - 1930,2048
      line: 295,0 - 1953,2048
      line: 260,0 - 1978,2048
      line: 258,0 - 1977,2048
      line: 222,0 - 2002,2048
      line: 184,0 - 2028,2048
      line: 144,0 - 2054,2048
      line: 103,0 - 2080,2048
      line: 0,-60 - 2048,1987
      line: 0,-16 - 2048,1961
      line: 0,26 - 2048,1936
      line: 0,68 - 2048,1912
      line: 0,107 - 2048,1887
      line: 0,108 - 2048,1888
      line: 0,110 - 2048,1890
      line: 0,146 - 2048,1864
      line: 0,187 - 2048,1845
      line: 0,223 - 2048,1824
      line: 0,259 - 2048,1803
      line: 0,291 - 2048,1779
      line: 0,325 - 2048,1759
      line: 0,358 - 2048,1740
      line: 0,395 - 2048,1725
      line: 0,427 - 2048,1707
      line: 0,458 - 2048,1689
      line: 0,513 - 2048,1648
      line: 0,571 - 2048,1615

      Please note that the reason the variable SIZ is references for all "image" widths *and* heights is that my maps are always square.

      Delete
  23. Hi Paul, is this all of your code ?
    I don't see the actual transformation in it.
    see "int Hough::Transform(unsigned char* img_data, int w, int h)"

    ReplyDelete
    Replies
    1. This comment has been removed by the author.

      Delete
    2. Thanks Bruno,
      I just realized this site keeps cropping out portions of my code for some reason. It wasn't displaying everything. Let me try to post just that first part, where it scans the image/map and builds the accumulator:
      int threshold=6;

      //Create the accu
      double hough_h = ((sqrt(2.0) * (double)(SIZ)) / 2.0);
      int _accu_h = hough_h * 2.0; // -r -> +r
      int _accu_w = 180;
      unsigned int *_accu = (unsigned int*)calloc(_accu_h * _accu_w, sizeof(unsigned int));
      double center_x = SIZ/2;
      double center_y = SIZ/2;

      for(int y=0;y<SIZ;y++) {
      for(int x=0;x<SIZ;x++) {
      int nyb=get_nybINT(world,x,y);
      if(nyb==2||nyb==3) {
      for(int t=0;t<180;t++) {
      double r = ( ((double)x - center_x) * cos((double)t * DEG2RAD)) + (((double)y - center_y) * sin((double)t * DEG2RAD));
      _accu[ (int)((round(r + hough_h) * 180.0)) + t]++;
      }
      }
      }
      }

      Delete
    3. Actually, the more I look at this code, the more I realize it's probably not going to get me what I want. This is checking all angles/radii from the center of the "image" and finding the lines that have a certain threshold (number of points) within them. But I gather that it's finding lines, not line segments (hence all the vector points on the perimeter of the image). I could have potentially have a section of fence on my map that has the same angle as another section of fence that's not connected and way off in another area of the map and just happens to line up. I need to be able to find the end points of every line segment on my map and build an array. Consider a huge matrix of integers where 2's and 3's mean a fence runs through that location and anything else can be considered empty. This algorithm seems like a great start to be able to compute where all straight fence runs start and end but I can't quite figure out how to modify the code to achieve that. Fences can also change direction at points and new endpoints should be detected at those junctions.

      Also, from the output, it seems like it's finding too much. There is a single angled line segment in the test matrix I used. I didn't even have any extra "noise" for it to find.

      Delete
  24. Now this is in actual fact cooperative. It’s very openhanded of you to share this with us.
    tutorial on c++

    ReplyDelete