C code from Marsaglia & Tsang’s “A simple method for generating gamma variables” (2000), with typos corrected http://dl.acm.org/citation.cfm?doid=358407.358414

rgama.c

#include <math.h>

extern float RNOR; // Normal random variable
extern float UNI;  // Uniform random variable

float rgama(float a) {
  float d, c, x, v, u;
  d = a - 1. / 3.;
  c = 1. / sqrt(9. * d);
  for (;;) {
    do {
      x = RNOR;
      v = 1. + c * x;
    } while (v <= 0.);
    v = v * v * v;
    u = UNI;
    if (u < 1. - 0.0331 * (x * x) * (x * x)) {
      return (d * v);
    }
    if (log(u) < 0.5 * x * x + d * (1. - v + log(v))) {
      return (d * v);
    }
  }
}