Efficient Ellipse Drawing - Part 1

This is the first in a multi-part series on drawing ellipses, pixel by pixel. Part 1 (this part) covers the basics of shape rendering, using lines as an example. Part 2 covers drawing circles. Following parts will cover ellipses.

My thesis is that any programmer can learn the math and tricks to write shape-drawing code. That's the point of this series -- to show you how to write your own code. If you just want code to copy-n-paste, go to the bottom of this post.

WARNING: Blogger.com likes to eat my code. I've reviewed this a few times, but there's no promise that it hasn't suddenly gotten hungry again. If something looks screwy, drop a comment and I'll take a look.

Background

I needed to implement efficient ellipse-drawing code in the tile editor I'm making, and I searched the net for some sample code but was thwarted. I vaguely remembered the algorithm, but I knew there was some math involved and some tricky stuff to make it fast, so I wanted to get it right. I looked at some sample code, converted it, but it didn't work. I wound up having to re-derive the equations by hand.

Part of the problem is a general problem with all sample code -- you have a specific purpose which may or may not match the sample. In my case, I was working in C# and had to clip the ellipse to the small tile. Plus little stuff, like how you specify the ellipse. A corner-to-corner definition (a bounding rect) can define an ellipse whose center is a pixel, the point between four pixels, or possibly even a point on the midpoint of an edge between two pixels. How do you adapt ellipse-drawing code to those situations?

The answer, really, is that you have to understand the algorithm itself. If you just want a library to do it for you, you can do that. In this case, that's not what I needed. I didn't want to spend an hour or two finding, downloading, installing, and integrating some third-party library -- or even just a single function -- to get something that I could do by myself in the same amount of time. What I wound up coding isn't actually an "efficient ellipse-drawing algorithm" -- I used a subfunction to calculate the error value for each pixel, with no strength-reduction accumulation code to make it faster. For what I was doing, that was fast enough, and the extra work to optimize the function wasn't necessary.

Goal

My thesis is this: line and ellipse (and circle) drawing is simple enough to learn.

Once you understand the basics, implementing your own code is straightforward.

Drawing a Line

Let's say you want to draw a line on a bitmap. If it's a straight vertical or horizontal line, the code is simple enough:

for (int x=0; x<length; ++x)
{
PlotPixel( x, y );
}

When you implement this algorithm in your own code, changes are it's gonna change to something like this:

int data = bitmap->GetDataPointer();
int offset = x + y * bitmap->GetRowStride();
for (int i=0; i<length; ++i)
data[offset+i] = pixelColor;

This is an example of what I meant above when I mentioned rolling your own. Most likely, any time you take an algorithm like this, you're going to change it enough that it won't look like the sample code you're using as a reference. That's why it's important to understand the algorithm. Straight lines are pretty simple, and you might not see much point here, but ellipses are much more complex.

Let's say you want to draw a diagonal line. That's also pretty easy.

for (int x=0; x<length; ++x)
PlotPixel( x, x );

That's kinda cheating, tho. When you adapt that bit of sample code, more likely you'll do something like this:

for (int i=0; i<length; ++i)
{
x = x1 + i;
y = y1 + i;
PlotPixel( x, y );
}

But that just draws a diagonal line. What if you want to draw a horizontal line, or a vertical one, or a diagonal that's at an angle other than 45º? Usually, you have a routine like "draw a line between these two points," and you have no idea where those two points are relative to each other. The complications pile up.

Let me start with changing that function to one that will draw a 45º line in any direction, including vertical and horizontal lines. Say that x2 and x1 are forty pixels apart. I'm going to draw forty pixels; one at each X coordinate. This means I start at X1, draw a pixel, move one pixel to the right, and repeat. This code handles going left or right with a loop:
for (x = x1; x!=x2; x+=inc)
The trick there is that 'inc' variable. If x2 is to the right (greater than) x1, then the increment is positive (+1). If x2 is to the left (less than) x1, then the increment is negative. I actually do the same thing with y.

Another option is to check to see if x2 is to the right of x1; if it's not, then swap the two points.

Anyway, this function will draw a line between two points, either flat (horizontal or vertical) or at a 45º diagonal.

DrawLine45or90( int x1, int y1, int x2, int y2 )
{
int dx = x2 - x1;
int dy = y2 - y1;
int stepx = (dx < 0) ? -1 : 1;
int stepy = (dy < 0) ? -1 : 1;
int len = dx * stepx;
if (dx==0)
{
stepx = 0;
len = dy * stepy;
}
else if (dy==0)
stepy = 0;
int x = x1;
int y = y1;
for (int i=0; i<len; +=i)
{
PlotPixel( x, y );
x += stepx;
y += stepy;
}
}

On modern machines, conditional code is expensive. Whether you're working on a PC, a modern console, a mobile device, or writing code for a GPU, if statements should be avoided in inner loops. That's why I do this business with stepx and stepy. Those increments also mean you don't need separate implementations for vertical lines, horizontal lines, and diagonals -- this loop does them all. Once you've got len, stepx, and stepy set up, just let it run.

Next: What if you want to draw a line at a different angle, something other than 0° or 45°? That's really the heart of a line-drawing algorithm -- handling lines between two arbitrary points.

From here on (except for the last sample, which is the final working algo), I'll assume that you're drawing a line where dx is greater than dy, ie a mostly-horizontal line.

Here's a simple implementation:

DrawLine( int x1, int y1, int x2, int y2 )
{
int dx = x2 - x1;
int dy = y2 - y1;
int stepx = (dx<0) ? -1 : 1;
int length = dx * stepx;
int x = x1;
int y = y1;
for (int i=0; i<length; ++i)
{
PlotPixel( x, y );
x += stepx;
y += stepx * (dy / dx);
}
}

This doesn't work, though. If dy is less than dx, then you wind up with a horizontal line (image on the left). If dy is greater than dx, then the line doesn't end at your target point -- it climbs at a constant rate and doesn't match what you want (image on the right).





The y value needs to be calculated in the loop:

DrawLine( int x1, int y1, int x2, int y2 )
{
int dx = x2 - x1;
int dy = y2 - y1;
int stepx = (dx<0) ? -1 : 1;
int length = dx * stepx;
int x = x1;
for (int i=0; i<length; ++i)
{
int y = y1 + stepx * (i * dy) / dx;
PlotPixel( x, y );
x += stepx;
}
}

The next trick is to get rid of that division in the inner loop. The way we do that is by conditionally incrementing y. This is why we assume a mostly-horizontal line. Since Y will change by either 0 or 1 each time through the loop, we can always change X but only sometimes change Y. For a finished implementation here, you will probably want two functions -- one for steep (mostly-vertical) lines and a separate function for shallow (mostly-horizontal) lines. The code samples will still be able to handle cases where point one is above, below, to the left, or to the right of point two.

OK, so to get rid of the division, we just keep an accumulator, adding dy every time, and check to see if it's greater than dx. We're going to have to watch the signs, tho.

DrawHorizontalishLine( int x1, int y1, int x2, int y2 )
{
int dx = x2 - x1;
int dy = y2 - y1;
int stepx = (dx<0) ? -1 : 1;
int stepy = (dy<0) ? -1 : 1;
int length = dx * stepx;
int x = x1;
int y = y1;
int yAccumulator = 0;
int yIncrement = dy * stepy;
int yTest = dx * stepx;
for (int i=0; i<length; ++i)
{
PlotPixel( x, y );
x += stepx;
yAccumulator += yIncrement;
if (yAccumulator >= yTest)
{
yAccumulator -= yTest;
y += stepy;
}
}
}

Final Algorithm

DrawLine( int x1, int y1, int x2, int y2 )
{
int dx = x2 - x1;
int dy = y2 - y1;
int stepx = (dx<0) ? -1 : 1;
int stepy = (dy<0) ? -1 : 1;
int lenx = dx * stepx;
int leny = dy * stepy;
int x = x1;
int y = y1;
int accumulator = 0;
if (lenx >= leny)
{
int increment = leny; // note the substitution here
int test = lenx;
for (int i=0; i<lenx; ++i)
{
PlotPixel( x, y );
x += stepx;
accumulator += increment;
if (accumulator >= test)
{
accumulator -= test;
y += stepy;
}
}
}
else
{
int increment = lenx;
int test = leny;
for (int i=0; i<leny; ++i)
{
PlotPixel( x, y );
y += stepy;
accumulator += increment;
if (accumulator >= test)
{
accumulator -= test;
x += stepx;
}
}
}
}

Further Parts

Part 2 - Drawing a circle
Part 3 - pending; circle complications
Part 4 - ellipses