How can I detect an ellipse in a photo?

Basically, you want to fit a shape to a set of points with outliers. One common algorithm to do this is RANSAC (random sample consensus). The basic outline of this algorithm is:

1. Select N points at random (where N is the minimum number of points required for fitting the shape, i.e. 2 for a line, 3 for a circle and so on)
2. Fit the shape to these points
3. Repeat 1&2, pick the "best" shape (where "best" usually means closest to a randomly sampled test point - I'll use a simpler criterion below.)
4. Select the points close to that shape and fit an improved estimate to those points

Let's try this. (Or a slight modification of it, anyway.)

First, I perform standard edge extraction:

img = Import["http://i.stack.imgur.com/H63BK.jpg"];
edges = DeleteBorderComponents[EdgeDetect[img,5]];
edgePos = N[PixelValuePositions[edges, 1]];

Next, I define the shape I want to fit - in this case, a conic:

conic[{x_, y_}] :=  a[1] + a[2]*x + a[3]*y + a[4]*x*y + (1 + a[5])*x^2 + (1 - a[5])*y^2
conicNormal[{x_, y_}] = D[conic[{x, y}], {{x, y}}]
isValidEllipse[solution_] := Abs[a[4] /. solution[[2]]] < 0.1

isValidEllipse is a simple ad-hoc criterion to filter out completely deformed conics.

Normally, we'd need 5 points to fit a conic. The image contains about 50% outliers, so the chance to randomly select 5 points on the ellipse we're looking for is about $0.5^5$, which is not good. We'd need a lot of samples to be reasonably sure that at least one 5-tuple of points contains only non-outlier points.

But fortunately, the location of the edge is not the only information we have: We also know the local orientation at each pixel, which should be perpendicular to the conic's normal vector.

orientation = GradientOrientationFilter[img, 3];

Now 3 points and their orientations give us 6 equations, so we have an overdetermined equation system:

fitEllipse[samplePoints_] := Quiet[FindMinimum[
   Total[
    Flatten[
     {
      (conic /@ samplePoints)^2,
      With[{\[Alpha] = PixelValue[orientation, #]},
          (conicNormal[#].{Cos[\[Alpha]], Sin[\[Alpha]]})^2] & /@ samplePoints
      }]], Array[a, 5]]]

This function returns not only the best-fit conic's coefficients, but also the residual error, which is a cheap way to compare between different fitted conics. The assumption is: If randomly sampled 3 points are all parts of the ellipse, the residual error will be low. If the points don't belong to one conic, the residual error will be much higher.

potentialSolutions = 
  SortBy[Select[Table[fitEllipse[RandomSample[edgePos, 3]], {n, 100}],
     isValidEllipse], First];

result = potentialSolutions[[1]];

(There is room for improvement here: the ellipse you're looking for might not be contained in these 100 samples, or it might not be the one with the lowest residual error. A smarter algorithm would take e.g. 50 samples, take the best one or two of these and count the number of nearby points. If enough points are close to it, keep that ellipse, otherwise, take 50 new samples. But that's just a straightforward programming exercise left to the reader ;-) .)

Next step: Find points nearest to the conic. I've tried to calculate the geometric distance exactly (using FindMinimumValue) but that's very slow. A simple, fast approximation is to find N points on the ellipse and simply use Nearest to estimate the distance:

ellipsePts = 
  Select[Flatten[
    Table[{i, y} /. 
      Solve[conic[{x, y}] == 0 /. result[[2]] /. x -> i, y], {i, 0, w,
       10}], 1], Total[Abs[Im[#]]] == 0 &];
nf = Nearest[ellipsePts];

Now, we can calculate the distance of every edge pixel to this conic and pick the ones that are closer than some threshold:

distances = Norm[nf[#][[1]] - #] & /@ edgePos;
closestPoints = Pick[edgePos, # < 10 & /@ distances];

...and improve the ellipse estimate based on these points:

result = fitEllipse[closestPoints];

Repeat the last steps until convergence. (Possible improvements: You could try to reduce the distance threshold in each iteration. Or you could add a "weight" to the curve fitting and give points closer to the current estimate a higher weight when calculating the next estimate.)

The result (after 5 iterations) looks like this:

Show[
 img,
 ContourPlot[
  Evaluate[conic[{x, y}] == 0 /. result[[2]]], {x, 0, w}, {y, 0, h}, 
  ContourStyle -> {Red, Thick}],
 ListPlot[closestPoints]]

Disclaimer: The result will not always look like this, as the algorithm is randomized and the algorithm can get stuck in local minima.

Too late and not so general as @Nikie but I will still like to answer in order to highlight a method involving generalized Eigen values! Lets take the image.

img = Import["http://i.stack.imgur.com/H63BK.jpg"];
(* Detecting ellipsoidal edges and removing the straight line component *)
i = ImageSubtract @@ {EdgeDetect[img, 10],EdgeDetect[img, 10, "StraightEdges" -> .2]};
{wid, hig} = ImageDimensions[i];
(* Taking control of the two straight line segments *)
lines = ImageLines[EdgeDetect[img, 10, 0.025], .15,"Segmented" -> True];
(* Find the mid point of longest line and its hight *)
clip = Ceiling@Abs@Mean[Last@Flatten@Differences@
          First@Sort[#,(EuclideanDistance @@@ #1 >EuclideanDistance @@@ #2)&]&/@lines];
(*Trim the bottom content using above clip and extract the points*)
padim = ImagePad[ImageTake[i, clip], {{0, 0}, {hig - clip, 0}}, Black];
list = Position[Transpose@ImageData[padim, "Bit", DataReversed -> True], 1];

Lets see the above steps of image processing results clockwise starting from top left.

image = Show[img, Graphics[{Thick, Blue, Line /@ lines}],ImageSize -> 300];
GraphicsGrid[{{img, ColorNegate@i}, {image, ColorNegate@padim}},ImageSize -> 800,
Frame -> All, Spacings -> 0]

At this point we can use EigenSystem efficiently as very simplistically discussed in here!

{X, Y} = Transpose[list];
(* design matrix *)
d = Transpose@{X^2, X Y, Y^2, X , Y, ConstantArray[1, Length@X]};
s = Transpose@d . d;
(* constraints *)
c = ConstantArray[0, {6, 6}];
c[[1, 3]] = -2;
c[[2, 2]] = 1;
c[[3, 1]] = -2;
{val, vec} = Eigensystem[{N[s, 200], c}];
(* parameter for a general conic*)
{a, b, c, d, e, f} = vec // Last;

Now we are ready to check out the fitting.

(* Visualize fit *)
ptEig = Show[img, 
ContourPlot[
a x^2 + b x y + c y^2 + d x + e y + f == 0, {x, 0, 1385}, {y, 0,1023},
ContourStyle -> {Thickness[0.004], Blue}, 
PlotPoints -> 200, MaxRecursion -> 4, PlotRange -> All, 
Frame -> None, Background -> None], 
Graphics[{PointSize[0.001], Red, Point[#] & /@ list}], 
ImageSize -> 400]

This didn't work as well as I'd hoped, but it's quite interesting to see what went wrong.

i = Import["http://dailycoffeenews.com/wp-content/uploads/2011/09/disposable-coffee-cup.jpeg"];

j = ColorConvert[i, "GrayScale"];
k = Thinning[
      DeleteBorderComponents[
       DeleteSmallComponents[
         Opening[
           Dilation[
             EdgeDetect[j], 1], 1]]]

which gives us a trace of the cup:

(Aargh, I've lost the rim. Pressing on anyway:)

pts = PixelValuePositions[k, "Max"];

left = First@Sort[pts, #1[[1]] < #2[[1]] &];
right = Last@Sort[pts, #1[[1]] < #2[[1]] &];
top = First@Sort[pts, #1[[2]] > #2[[2]] &];
centre = Mean[{left, right}];
rim = Disk[centre, {First[centre - left], Last[top - centre]}];

Show[
 i,
 Graphics[{Opacity[0.5], Brown, rim}]]

It's close, but not close enough. And now I can see where I went wrong - I've assumed the cup is aligned with the camera and geometrically accurate, but I don't think it is. (And what exactly is the rim anyway?) I wonder whether this method would work if the cup was aligned.