Friday, May 31, 2013

Efficiently Sample from Normal (Gaussian) Distribution

The following code snippet samples from a standard normal distribution using the Box-Muller transform, which avoids nastiness like cropping to a fixed number of standard-deviations or rejection sampling.  I have shamelessly cribbed this from Wikipedia (the previous link), but use it all the time, so I'm putting it up as a snippet.

void utils_sample_normal( float &x, float &y ){
    float u=drand48(), v=drand48();
    float lnu = sqrt( -2.0*log(u) );
    x = lnu*cos(2.0*M_PI*v);
    y = lnu*sin(2.0*M_PI*v);

No comments: