Compare DDA/Bresenham; midpoint ellipse

Comparison of DDA and Bresenham's line algorithms

Midpoint ellipse algorithm

Derivation

For an ellipse centred at the origin with semi-axes rx and ry:

f(x, y) = ry²x² + rx²y² − rx²ry²
f < 0 inside,  f = 0 on the ellipse,  f > 0 outside

By four-way symmetry only the first quadrant is computed. It is divided into two regions at the point where the slope is −1. The slope is dy/dx = −(ry²x)/(rx²y), so region 1 is where 2ry²x < 2rx²y (slope magnitude < 1) and region 2 is where 2ry²x ≥ 2rx²y.

Region 1: x is incremented by 1. The next y is y_k or y_k − 1; the midpoint is (x_k+1, y_k − ½).

p1_k = f(x_k+1, y_k−½) = ry²(x_k+1)² + rx²(y_k−½)² − rx²ry²
p1_k+1 − p1_k = 2ry²(x_k+1) + ry² + rx²[(y_k+1 − ½)² − (y_k − ½)²]
p1_k+1 = p1_k + 2ry²x_k+1 + ry²                    if p1_k < 0  (y same)
p1_k+1 = p1_k + 2ry²x_k+1 − 2rx²y_k+1 + ry²        if p1_k ≥ 0  (y−1)
p1_0   = f(1, ry−½) = ry² − rx²ry + rx²/4

Region 2: y is decremented by 1. The next x is x_k or x_k + 1; the midpoint is (x_k+½, y_k−1).

p2_k = f(x_k+½, y_k−1) = ry²(x_k+½)² + rx²(y_k−1)² − rx²ry²
p2_k+1 − p2_k = −2rx²(y_k − 1) + rx² + ry²[(x_k+1 + ½)² − (x_k + ½)²]
p2_k+1 = p2_k − 2rx²y_k+1 + rx²                    if p2_k > 0  (x same)
p2_k+1 = p2_k + 2ry²x_k+1 − 2rx²y_k+1 + rx²        if p2_k ≤ 0  (x+1)
p2_0 uses the last point (x0, y0) of region 1:
p2_0   = ry²(x0+½)² + rx²(y0−1)² − rx²ry²

Algorithm

Input rx, ry and centre (xc, yc).
Region 1 (slope magnitude < 1):
  x = 0, y = ry
  p1 = ry² − rx²·ry + rx²/4
  while 2·ry²·x < 2·rx²·y:
      plot the 4 symmetric points of (x, y)
      x = x + 1
      if p1 < 0:
          p1 = p1 + 2·ry²·x + ry²
      else:
          y = y − 1
          p1 = p1 + 2·ry²·x − 2·rx²·y + ry²
Region 2 (slope magnitude ≥ 1), start from last (x, y):
  p2 = ry²(x + ½)² + rx²(y − 1)² − rx²·ry²
  while y ≥ 0:
      plot the 4 symmetric points of (x, y)
      y = y − 1
      if p2 > 0:
          p2 = p2 − 2·rx²·y + rx²
      else:
          x = x + 1
          p2 = p2 + 2·ry²·x − 2·rx²·y + rx²
Each plotted point is shifted by (xc, yc): (xc±x, yc±y).