Ellipse Detection using Hough Transform

using Hough Transform, how can I detect and get coordinates of (x0,y0) and "a" and "b" of an ellipse in 2D space?

This is ellipse01.bmp:

I = imread('ellipse01.bmp');
[m n] = size(I);
c=0;
for i=1:m
    for j=1:n
        if I(i,j)==1
        c=c+1;
        p(c,1)=i;
        p(c,2)=j;
        end
    end
end
Edges=transpose(p);
Size_Ellipse = size(Edges);
B = 1:ceil(Size_Ellipse(1)/2);
Acc = zeros(length(B),1);
a1=0;a2=0;b1=0;b2=0;
Ellipse_Minor=[];Ellipse_Major=[];Ellipse_X0 = [];Ellipse_Y0 = [];
Global_Threshold = ceil(Size_Ellipse(2)/6);%Used for Major Axis Comparison
Local_Threshold = ceil(Size_Ellipse(1)/25);%Used for Minor Axis Comparison
[Y,X]=find(Edges);
Limit=numel(Y);
Thresh = 150;
Para=[];

for Count_01 =1:(Limit-1)
  for Count_02 =(Count_01+1):Limit
    if ((Count_02>Limit) || (Count_01>Limit))
      continue
    end
    a1=Y(Count_01);b1=X(Count_01);
    a2=Y(Count_02);b2=X(Count_02);
    Dist_01 = (sqrt((a1-a2)^2+(b1-b2)^2));
    if (Dist_01 >Global_Threshold)
      Center_X0 = (b1+b2)/2;Center_Y0 = (a1+a2)/2;
      Major = Dist_01/2.0;Alpha = atan((a2-a1)/(b2-b1));
      if(Alpha == 0)
        for Count_03 = 1:Limit
          if( (Count_03 ~= Count_01) || (Count_03 ~= Count_02))
            a3=Y(Count_03);b3=X(Count_03);
            Dist_02 = (sqrt((a3 - Center_Y0)^2+(b3 - Center_X0)^2));
            if(Dist_02 > Local_Threshold)
              Cos_Tau = ((Major)^2 + (Dist_02)^2 - (a3-a2)^2 - (b3-b2)^2)/(2*Major*Dist_02);
              Sin_Tau = 1 - (Cos_Tau)^2;
              Minor_Temp = ((Major*Dist_02*Sin_Tau)^2)/(Major^2 - ((Dist_02*Cos_Tau)^2));
              if((Minor_Temp>1) && (Minor_Temp<B(end)))
                Acc(round(Minor_Temp)) = Acc(round(Minor_Temp))+1;
              end
            end
          end
        end
      end
      Minor = find(Acc == max(Acc(:)));
      if(Acc(Minor)>Thresh)
        Ellipse_Minor(end+1)=Minor(1);Ell开发者_StackOverflow社区ipse_Major(end+1)=Major;
        Ellipse_X0(end+1) = Center_X0;Ellipse_Y0(end+1) = Center_Y0;
        for Count = 1:numel(X)
          Para_X = ((X(Count)-Ellipse_X0(end))^2)/(Ellipse_Major(end)^2);
          Para_Y = ((Y(Count)-Ellipse_Y0(end))^2)/(Ellipse_Minor(end)^2);
          if (((Para_X + Para_Y)>=-2)&&((Para_X + Para_Y)<=2))
            Edges(X(Count),Y(Count))=0;
          end
        end
      end
      Acc = zeros(size(Acc));
    end
  end
end

---

Although this is an old question, perhaps what I found can help someone.

The main problem of using the normal Hough Transform to detect ellipses is the dimension of the accumulator, since we would need to vote for 5 variables (the equation is explained here):

There is a very nice algorithm where the accumulator can be a simple 1D array, for example, and that runs in . If you wanna see code, you can look at here (the image used to test was that posted above).

---

If you use circle for rough transform is given as rho = xcos(theta) + ysin(theta) For ellipse since it is

You could transform the equation as rho = axcos(theta) + bysin(theta) Although I am not sure if you use standard Hough Transform, for ellipse-like transforms, you could manipulate the first given function.

---

If your ellipse is as provided, being a true ellipse and not a noisy sample of points; the search for the two furthest points gives the ends of the major axis, the search for the two nearest points gives the ends of the minor axis, the intersection of these lines (you can check it's a right angle) occurs at the centre.

---

If you know the 'a' and 'b' of an ellipse then you can rescale the image by factor of a/b in one direction and look for circle. I am still thinking about what to do when a and b are unknown.

If you know that it is circle then use Hough transform for circles. Here is a sample code:

int accomulatorResolution  = 1;  // for each pixel
    int minDistBetweenCircles  = 10; // In pixels
    int cannyThresh            = 20;
    int accomulatorThresh      = 5*_accT+1;
    int minCircleRadius        = 0;
    int maxCircleRadius        = _maxR*10;
    cvClearMemStorage(storage);
    circles = cvHoughCircles( gryImage, storage,
                              CV_HOUGH_GRADIENT, accomulatorResolution, 
                              minDistBetweenCircles,
                              cannyThresh , accomulatorThresh,
                              minCircleRadius,maxCircleRadius );    
    // Draw circles
    for (int i = 0; i < circles->total; i++){
        float* p = (float*)cvGetSeqElem(circles,i);
        // Draw center
        cvCircle(dstImage, cvPoint(cvRound(p[0]),cvRound(p[1])),
                           1, CV_RGB(0,255,0), -1, 8, 0 );
        // Draw circle
        cvCircle(dstImage, cvPoint(cvRound(p[0]),cvRound(p[1])),
                           cvRound(p[2]),CV_RGB(255,0,0), 1, 8, 0 );
    }