How does one draw a circle? You wipe the dust off your Foley &
van Dam book and look up the circle drawing algorithm, or simply use
Cairo. But I had to explain to someone how to draw circles and going
directly to the Bresenham method wasn't reasonable. I started with a
basic implementation and we discussed ways to optimize it
in `#gimp`

. The code here was written a few days ago, and is
posted for anyone who searches for it.

Circles are invisible. When we say we want to draw a circle on the
screen, we want to set pixels (light up equally sized rectangular bits
of the screen) at integer coordinates to represent an approximation of
the circle. Let's start off by writing some code which implements
a **main**(), handles creating an image, setting pixels in
it and writing it to a Portable graymap
(`.pgm`

) file (which can be viewed with programs like Eye of
GNOME, GIMP, ImageMagick's `display`

, etc.). With that out of
the way, we can work on **draw_circle**() alone.

#include <stdio.h> #include <stdlib.h> #include <string.h> #include <math.h> typedef struct { size_t width; size_t height; unsigned char *data; } Image; static Image * image_new (size_t width, size_t height) { Image *image; image = malloc (sizeof *image); image->width = width; image->height = height; image->data = malloc (width * height); return image; } static void image_free (Image *image) { free (image->data); free (image); } static void image_fill (Image *image, unsigned char value) { memset (image->data, value, image->width * image->height); } /** * image_set_pixel: * * Sets a pixel passed in signed (x, y) coordinates, where (0,0) is at * the center of the image. **/ static void image_set_pixel (Image *image, ssize_t x, ssize_t y, unsigned char value) { size_t tx, ty; unsigned char *p; tx = (image->width / 2) + x; ty = (image->height / 2) + y; p = image->data + (ty * image->width) + tx; *p = value; } static void image_save (const Image *image, const char *filename) { FILE *out; out = fopen (filename, "wb"); if (!out) return; fprintf (out, "P5\n"); fprintf (out, "%zu %zu\n", image->width, image->height); fprintf (out, "255\n"); fwrite (image->data, 1, image->width * image->height, out); fclose (out); } static void draw_circle (Image *image, int radius, unsigned char value); int main (int argc, char *argv[]) { Image *image; image = image_new (600, 600); image_fill (image, 0xff); draw_circle (image, 200, 0); image_save (image, "circle.pgm"); image_free (image); return 0; }

**main**() basically creates a 600×600 pixels
image, fills it with white, and asks **draw_circle**() to
draw a black circle of radius 200 pixels centered in the image. Note
that **image_set_pixel**() takes signed coordinates, where
(x=0, y=0) is at the center of the image. So for an image created with
`image_new (400, 400)`

, calling
`image_set_pixel (image, 0, 0, 0x80)`

will set the pixel at
(200, 200) to 0x80.

So let's implement various versions of **draw_circle**()
now. The well-known equation of a circle is \(x^2+y^2=r^2\). The
variables are shown relative to each other in the figure below:

This type of an equation is called an *implicit equation*. You
can substitute values for \(x\), \(y\) and \(r\) to see if the equation
holds true. If it does, then the point \((x,y)\) is on the circle. So
the implicit equation lets you determine if a point is on the circle or
not.

In the figure above, you can see that for every line parallel to the \(y\)-axis through some \(x \in (-r,r)\), the circle intersects it in two points. If we scan \(x\) from \(-r\) to \(r\) and plot the two intersection points \((x,y_1)\) and \((x,y_2)\) in every vertical scanline, at the end of the scan we'll have drawn a circle. For this, we have to determine \(y_1\) and \(y_2\) for every \(x\). The implicit equation is no good for it, because it is only useful in checking values and not for determining them. So let's rearrange terms to find \(y\):

$$ \begin{eqnarray*} x^2+y^2 & = & r^2 \\ y^2 & = & r^2 - x^2 \\ y & = & \sqrt{r^2 - x^2} \end{eqnarray*} $$

The square-root returns a positive and negative value for \(y\), so
\(y_1=y\) and \(y_2=-y\), and the points at which the circle intersects
a scanline through \(x\) are \((x,y)\) and \((x,-y)\). So let's write
a **draw_circle**() based on this:

static void draw_circle (Image *image, int radius, unsigned char value) { int x, y; for (x = -radius; x <= radius; x++) { y = (int) sqrt ((double) (radius * radius) - (x * x)); image_set_pixel (image, x, y, value); image_set_pixel (image, x, -y, value); } }

Download circle1.c. This program generates the following image:

This looks like a circle, but it is discontinuous near the \(x\)-axis. At those places, for every new scanline \(x_{n+1}=x_{n}+1\), \(y\) changes by more than \(1\) (\(y_{n+1}-y_n > 1\)).. sometimes far more than \(1\). Notice that there are discontinuities between angles 0 and 45 degrees, and no discontinuities between angles 45 and 90 degrees in the drawn circle.

In a circle, all points are at the same distance from the center. In
the code above, we only compute a semi-circle and mirror it by plotting
at \((x,y)\) and \((x,-y)\). In the same way, we can only compute points
between the angles 45 degrees and 90 degrees and mirror it 8 ways by
plotting \((x,y)\), \((x,-y)\), \((-x,y)\), \((-x,-y)\), \((y,x)\),
\((y,-x)\), \((-y,x)\) and \((-y,-x)\). It plots 8 arcs, each of 45
degrees. So let's modify
**draw_circle**() to do this:

static void draw_circle (Image *image, int radius, unsigned char value) { int x, y; int l; l = (int) radius * cos (M_PI / 4); for (x = 0; x <= l; x++) { y = (int) sqrt ((double) (radius * radius) - (x * x)); image_set_pixel (image, x, y, value); image_set_pixel (image, x, -y, value); image_set_pixel (image, -x, y, value); image_set_pixel (image, -x, -y, value); image_set_pixel (image, y, x, value); image_set_pixel (image, y, -x, value); image_set_pixel (image, -y, x, value); image_set_pixel (image, -y, -x, value); } }

Download circle2.c. This program generates the following image:

So we're done, right? If code to draw a circle is all we need, then we are done. On the other hand, if you want to optimize it a bit so it performs acceptably on very slow processors, keep reading.

At this point I posted a link to the code
in `#gimp`

. `pippin`

noted that filling the disk
(the area inside the circle) can be much faster than drawing the
circle. The square root computation is not needed as you're merely
checking every pixel's coordinates to see if it lies inside the disk or
outside. This just needs testing against the implicit equation of the
disk \(x^2+y^2 \leq r^2\). Note how it is very similar to the implicit
equation of the circle which just checks the boundary of the disk. As
the implicit equation is checked for every pixel, filling the disk can
be done in a ridiculously parallel manner. Here's a modified version
of **draw_circle**() which plots a filled disk:

static void draw_circle (Image *image, int radius, unsigned char value) { int x, y; for (y = -radius; y <= radius; y++) for (x = -radius; x <= radius; x++) if ((x * x) + (y * y) <= (radius * radius)) image_set_pixel (image, x, y, value); }

Download circle3.c. This program generates the following image:

Ideally you'd rewrite such code to inline
the **image_set_pixel**() as the pixels are filled in an
ordered manner. A good C compiler should also eliminate redundant
subexpressions and move invariants out of loops, but some embedded
compilers are not good. So you'd want to add temporaries in such
cases:

static void draw_circle (Image *image, int radius, unsigned char value) { unsigned char *buffer; size_t pad; int r2; int x, y; buffer = image->data + (((image->height / 2) - radius) * image->width); pad = (image->width / 2) - radius; r2 = radius * radius; for (y = -radius; y <= radius; y++) { int y2; unsigned char *p; y2 = y * y; p = buffer + pad; for (x = -radius; x <= radius; x++) { int x2 = x * x; if (x2 + y2 <= r2) *p = value; p++; } buffer += image->width; } }

Download circle4.c which generates the same
filled disk. Then the conversation in `#gimp`

dropped away
that night. I couldn't think of a way to avoid
that **sqrt**(), and resigned I went to bed.

The next morning I woke up to the news that Steve Jobs had died. Reading the Slashdot story, I saw a comment about rounded rectangles. The linked article mentioned how Bill Atkinson had used the arithmetic progression \(n^2=1+3+5+7+\cdots+(2n-1)\). Of course! We don't need exact square roots, only the nearest integer root. And between the angles 45 and 90 degrees where we work, \(y_{n}-y_{n+1} \leq 1\). So we can estimate whether \(y_{n+1}\) is a decrement by comparing \(y_{n+1}^2\) with \(y_n^2-(2y_n-1)\). Awesome, eh? This progression \(n^2=1+3+5+7+\cdots+(2n-1)\) was in the Knuth book, but somehow this and the circle case both didn't take their clothes off and jump into the same pool to get me all excited. Rather than give up that something can't be solved further, it's better to think that there's always a solution. It was there from before we were born just waiting to be found.

Also, let's see what happens to \(y=\sqrt{r^2-x^2}\) when \(x_{n+1}=x_n+1\) (the next vertical scanline).

$$ \begin{eqnarray*} y_{n+1}^2 & = & r^2 - x_{n+1}^2 \\ & = & r^2 - (x_n+1)^2 \\ & = & r^2 - (x^2_n+2x_n+1) \\ & = & r^2 - x_n^2 - 2x_n - 1 \\ & = & y_n^2 - 2x_n - 1. \end{eqnarray*} $$

So let's update **draw_circle**() to implement this:

static void draw_circle (Image *image, int radius, unsigned char value) { int x, y; int l; int r2, y2; int y2_new; int ty; /* cos pi/4 = 185363 / 2^18 (approx) */ l = (radius * 185363) >> 18; /* At x=0, y=radius */ y = radius; r2 = y2 = y * y; ty = (2 * y) - 1; y2_new = r2 + 3; for (x = 0; x <= l; x++) { y2_new -= (2 * x) - 3; if ((y2 - y2_new) >= ty) { y2 -= ty; y -= 1; ty -= 2; } image_set_pixel (image, x, y, value); image_set_pixel (image, x, -y, value); image_set_pixel (image, -x, y, value); image_set_pixel (image, -x, -y, value); image_set_pixel (image, y, x, value); image_set_pixel (image, y, -x, value); image_set_pixel (image, -y, x, value); image_set_pixel (image, -y, -x, value); } }

Download circle5.c which generates the same
circle. Notice how only additions, subtractions and shifts are enough in
the loop. **image_set_pixel**() can be inlined here too,
but I'm not going to do it.

So how do we draw rectangles with rounded corners? We currently draw the circle as 8 arcs, right? With the appropriate border radius, we simply translate and draw 2 arcs each at the 4 corners of the rectangle. :)