[Work] How to draw anti-aliased circles and disks

How to draw anti-aliased circles and disks
Theory:
We will use the basic equation for a circle. Using two coordinates, X and Y, and the center of the circle located at Xc and Yc, and a radius R:

sqr(X-Xc) + sqr(Y-Yc) = sqr(R)

Notice that the value sqr(R) is always the same, so we can pre-calculate it. On each row of pixels (called a scanline), the value of sqr(Y-Yc) is also constant, because the Y value is constant. And just so, for each column of pixels the value of sqr(X - Xc) is constant.

Using this, we can precalculate these values, and then the only thing left to do is to compare the correct values in each pixel with simple addition. This results in very fast code.

How to do anti-aliasing
Anti-aliasing means that we must find "inbetween" values whenever we are "on" the edge of the circle. Suppose we want to draw a disk, with radius R:

if sqr(X-Xc) + sqr(Y-Yc) < sqr(R) we are on the disk, if not, we are outside of the disk.

An easy way of implementing anti-aliasing here is to assume that there are two circles. We have an inner circle R1 and an outer circle R2:

if sqr(X-Xc) + sqr(Y-Yc) < sqr(R1) we are on the disk,

if sqr(X-Xc) + sqr(Y-Yc) > sqr(R2) we are outside,

and in the other cases we have fallen inbetween both circles. This is the area where the color must change from background color to disk color.

First we find our local radius: R = sqrt(sqr(X-Xc) + sqr(Y-Yc)). This requires a square root calculation, but fortunately only in the relatively small area of the anti-aliasing zone.

Next, we find the factor ranging from 0 (outside) to 1 (inside) and use this one to scale our color:

Factor = (R2 - R) / (R2 - R1), which is 1 if R=R1 and 0 if R=R2

#include 

#include "array.h"
#include "math.h"

#ifndef A2D
typedef Array_2D A2D;
#endif

#ifndef SQR
#define SQR(a) ((a)*(a))
#endif

/***************************************************************************************************/
// Draw a disk on Bitmap. Bitmap must be a 256 color (pf8bit) palette bitmap,
// and parts outside the disk will get palette index 0, parts inside will get
// palette index 255, and in the antialiased area (feather), the pixels will
// get values inbetween.
// ***Parameters***
// Bitmap:
//   The bitmap to draw on
// CenterX, CenterY:
//   The center of the disk (float precision). Note that [0, 0] would be the
//   center of the first pixel. To draw in the exact middle of a 100x100 bitmap,
//   use CenterX = 49.5 and CenterY = 49.5
// Radius:
//   The radius of the drawn disk in pixels (float precision)
// Feather:
//   The feather area. Use 1 pixel for a 1-pixel antialiased area. Pixel centers
//   outside 'Radius + Feather / 2' become 0, pixel centers inside 'Radius - Feather/2'
//   become 255. Using a value of 0 will yield a bilevel image.
/***************************************************************************************************/
void draw_disk(float CenterX, float CenterY, float Radius, float Feather, A2D *img){

  int width, height;
  width=(*img).x_size();
  height=(*img).y_size();
  

  int x, y;
  int LX, RX, LY, RY;
  //  int  Fact;
  float RPF2, RMF2;
  float SqY, SqDist;
  std::vector SqX;

  // Determine some helpful values (singles)
  RPF2 = SQR(Radius + Feather/2);
  RMF2 = SQR(Radius - Feather/2);

  // Determine bounds:
  LX = (int)MAX(floor(CenterX - RPF2), 0);
  RX = (int)MIN(ceil (CenterX + RPF2), width - 1);
  LY = (int)MAX(floor(CenterY - RPF2), 0);
  RY = (int)MIN(ceil (CenterY + RPF2), height - 1);

  // Optimization run: find squares of X first
  SqX.resize(RX - LX + 1, 0.0f);

  for( x = LX ; x<= RX; x++){
    SqX[x - LX] = SQR(x - CenterX);
  }

  // Loop through Y values
  for(y = LY ; y<= RY; y++) {
    SqY = SQR(y - CenterY);
    // Loop through X values
    for( x= LX; x<=RX; x++){
      // determine squared distance from center for this pixel
      SqDist = SqY + SqX[x - LX];
      // inside inner circle? Most often..
      if( SqDist < RMF2){
 // inside the inner circle.. just give the scanline the
 (*img)(x,y)= 1.0; 
      }else{
 // inside outer circle?
 if( SqDist < RPF2){
   // We are inbetween the inner and outer bound, now
   // mix the color
//    Fact = int(((Radius - sqrt(SqDist)) * 2 / Feather)
//    * 127.5 + 127.5 +0.5);
//    (*img)(x, y)= float(MAX(0, MIN(Fact, 255)));
   (*img)(x,y)= MAX(0, ((Radius - sqrt(SqDist)) * 2 / Feather)*0.5+0.5);
 }else{
   (*img)(x, y) = 0.0;
 }
      }
    }
  }
}

/***************************************************************************************************/
// Draw a circle on Bitmap. Bitmap must be a 256 color (pf8bit) palette bitmap,
// and parts outside the circle will get palette index 0, parts inside will get
// palette index 255, and in the antialiased area (feather), the pixels will
// get values inbetween.
// ***Parameters***
// Bitmap:
//   The bitmap to draw on
// CenterX, CenterY:
//   The center of the circle (float precision). Note that [0, 0] would be the
//   center of the first pixel. To draw in the exact middle of a 100x100 bitmap,
//   use CenterX = 49.5 and CenterY = 49.5
// Radius:
//   The radius of the drawn circle in pixels (float precision)
// LineWidth
//   The line width of the drawn circle in pixels (float precision)
// Feather:
//   The feather area. Use 1 pixel for a 1-pixel antialiased area. Pixel centers
//   outside 'Radius + Feather / 2' become 0, pixel centers inside 'Radius - Feather/2'
//   become 255. Using a value of 0 will yield a bilevel image. Note that Feather
//   must be equal or smaller than LineWidth (or it will be adjusted internally)
/***************************************************************************************************/

void draw_circle(float CenterX, float CenterY, float Radius,float  LineWidth, float Feather, A2D *img){

  int width, height;
  width=(*img).x_size();
  height=(*img).y_size();

  int x, y;
  int   LX, RX, LY, RY;
  float ROPF2, ROMF2, RIPF2, RIMF2;
  float OutRad, InRad;

  //  P: PByteArray;
  float SqY, SqDist;
  
  std::vector SqX;

  // Determine some helpful values (singles)
  OutRad = Radius + LineWidth/2;
  InRad  = Radius - LineWidth/2;
  ROPF2 = SQR(OutRad + Feather/2);
  ROMF2 = SQR(OutRad - Feather/2);
  RIPF2 = SQR(InRad  + Feather/2);
  RIMF2 = SQR(InRad  - Feather/2);

  // Determine bounds:
  LX = (int) MAX(floor(CenterX - ROPF2), 0);
  RX = (int) MIN(ceil (CenterX + ROPF2), width - 1);
  LY = (int) MAX(floor(CenterY - ROPF2), 0);
  RY = (int) MIN(ceil (CenterY + ROPF2), height - 1);

  // Checks
  if(Feather > LineWidth)
    Feather = LineWidth;

  // Optimization run: find squares of X first
  SqX.resize(RX - LX + 1, 0.0f);

  for( x = LX; x<=RX; x++)
    SqX[x - LX] = SQR(x - CenterX);

  // Loop through Y values
  for(y= LY; y<= RY; y++){
    SqY = SQR(y - CenterY);
    // Loop through X values
    for( x = LX; x<= RX; x++){

      // determine squared distance from center for this pixel
      SqDist = SqY + SqX[x - LX];

      // now first check if we're completely inside (most often)
      if (SqDist < RIMF2){
        // We're on the disk inside everything
        (*img)(x,y) = 0;
      } else {
        // completely outside?
 if( SqDist < ROPF2){
          // inside outer line - feather?
   if (SqDist < ROMF2){
     // check if we're in inside feather area
     if( SqDist < RIPF2){
       // We are in the feather area of inner line, now mix the color
       //       Fact = int(((sqrt(SqDist) - InRad) * 2 / Feather) * 127.5 + 127.5+0.5);
       //(*img)(x,y) = MAX(0, MIN(Fact, 255)); // just in case limit to [0, 255]
       (*img)(x,y)=MAX(0,(sqrt(SqDist) - InRad) * 2 / Feather*0.5+0.5);
     } else{
       // on the line
       (*img)(x,y) = 1.0;
     }
   }else{
     // We are in the feather area of outer line, now mix the color
     // Fact = int(((OutRad - sqrt(SqDist)) * 2 / Feather) * 127.5 + 127.5+0.5);
     //(*img)(x,y) = MAX(0, MIN(Fact, 255)); // just in case limit to [0, 255]
     (*img)(x,y)=MAX(0, (OutRad - sqrt(SqDist)) * 2 / Feather);
   }
 } else{
   // outside everything
   (*img)(x,y) = 0;
 }
      
      }
    }
  }


}


/***************************************/

int main(int argc, char *argv[] )
{
  

  int width, height;
  width=1200;
  height=1000;
  A2D img(width, height);
  

  draw_disk(width/2, height/2, 400, 10.0, &img);
  write_A2D(img, "disk");

  draw_circle(width/2+0.5, height/2+0.5, 450,4, 1, &img);

  write_A2D(img, "circle");

  cout << "All done!" << endl;
}