Skip to content

Instantly share code, notes, and snippets.

@FloooD
Created February 6, 2012 20:35
Show Gist options
  • Save FloooD/1754661 to your computer and use it in GitHub Desktop.
Save FloooD/1754661 to your computer and use it in GitHub Desktop.
circle antialiased by using integration to calculate pixel coverage. gamma corrected for srgb.
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#define PADDING 4
#define GAMMA_CORRECTED
#define SUBPIXEL_RATIO 0.5 //1 = full subpx aa, 0 = grayscale
double indef_integ(double r, double x)
{
if (fabs(x) > fabs(r)) {
printf("indef_integ: |%f| > |%f|\n", x, r);
exit(1);
}
return 0.5*(x*sqrt(r*r-x*x)+r*r*asin(x/r));
}
double def_integ(double r, double left, double right)
{
return indef_integ(r, right) - indef_integ(r, left);
}
int in(double r, double x, double y)
{
return x*x + y*y < r*r;
}
double *circle(double rad, int r) //r = half of length of square
{
int y, x;
double *c = malloc(4*r*r*sizeof(double));
#define L(a,b,c,d) (((a)<<3)|((b)<<2)|((c)<<1)|(d))
for (y = r-1; y >= 0; y--) for (x = 0; x < r && x <=y; x++) {
switch (L(in(rad,x,y),in(rad,x+1,y),in(rad,x,y+1),in(rad,x+1,y+1))) {
case L(0,0,0,0):
c[(y+r)*2*r+(x+r)] = 0.0;
break;
case L(1,1,1,1):
c[(y+r)*2*r+(x+r)] = 1.0;
break;
case L(1,1,0,0):
c[(y+r)*2*r+(x+r)] = def_integ(rad, x, x+1) - y;
break;
case L(1,0,0,0): {
double d = sqrt(rad*rad-y*y) - x;
c[(y+r)*2*r+(x+r)] = def_integ(rad, x, x+d) - y*d;
break;
}
case L(1,1,1,0):
c[(y+r)*2*r+(x+r)] = def_integ(rad, x, x+1) - y - c[(y+r+1)*2*r+(x+r)];
break;
default:
c[(y+r)*2*r+(x+r)] = 0.0;
break;
}
}
#undef L
for (y = 0; y < r; y++) for (x = r-1; x >= 0 && x > y; x--)
c[(y+r)*2*r+(x+r)] = c[(x+r)*2*r+(y+r)];
for (y = r; y < 2*r; y++) for (x = 0; x < r; x++)
c[y*2*r+x] = c[y*2*r+(2*r-1-x)];
for (y = 0; y < r; y++) for (x = 0; x < 2*r; x++)
c[y*2*r+x] = c[(2*r-1-y)*2*r+x];
return c;
}
#ifdef GAMMA_CORRECTED
int d2i(double d)
{
if (d<=0.0)
return 0;
else if (d>=1.0)
return 255;
else if (d>0.0031308)
return (int)(pow(d, 1.0/2.4)*269.025-13.525);
else
return (int)(3294.6*d+0.5);
}
#else
int d2i(double d)
{
return (int)(255*d+0.5);
}
#endif
void pgm(double *a, int w, int h)
{
int y, x;
printf("P2\n"
"%d %d 255\n\n", w, h);
for (y = 0; y < h; y++) {
for (x = 0; x < w; x++)
printf("%3d ", d2i(a[y*w+x]));
printf("\n");
}
}
void ppm(double *a, int w, int h)
{
int y, x;
printf("P3\n"
"%d %d 255\n\n", w, h);
for (y = 0; y < h; y++) {
for (x = 0; x < w; x++)
printf("%3d %3d %3d ", d2i(a[3*(y*w+x)]), d2i(a[3*(y*w+x)+1]), d2i(a[3*(y*w+x)+2]));
printf("\n");
}
}
int main(int argc, char *argv[])
{
if (argc != 2)
printf("Usage: %s radius > output.pbm\n", argv[0]);
else {
double rad = atof(argv[1]);
int r = ceil(fabs(rad)) + PADDING;
double *a = circle(rad, r);
#ifndef SUBPIXEL_RATIO
pgm(a, 2*r, 2*r);
free(a);
#else
double *b = circle(3*rad, 3*r);
double *c = malloc(3*4*r*r*sizeof(double));
int y, x, s;
for (y = 0; y < 2*r; y++) for (x = 0; x < 2*r; x++) for (s = 0; s < 3; s++)
c[3*(y*2*r+x)+s] =
SUBPIXEL_RATIO*(b[3*y*6*r+3*x+s]+b[(3*y+1)*6*r+3*x+s]+b[(3*y+2)*6*r+3*x+s])/3.0
+(1-SUBPIXEL_RATIO)*a[y*2*r+x];
ppm(c, 2*r, 2*r);
free(b);
free(c);
#endif
}
return 0;
}
all:
gcc -O3 -lm -Wall -pedantic -std=c99 -o circ circ.c
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment