/****************************************************************************** 
 * 
 * box_muller_transformation 
 * 
 * A Box-Muller transformation picks two numbers using a uniformly-random 
 * random number generator, and then uses those two numbers to come up with a 
 * normally-distributed number. In essence, it is a normally-distributed 
 * random number generator (using the standard normal [IE - avg = 0, std = 1]) 
 * 
 * Returns - A normally distributed number 
 * 
 ******************************************************************************/ 
 
float box_muller_transformation() { 

  /* Two uniformly random numbers,  and a temporary for the magnitude. */ 
  float ur1, ur2, magnitude; 
    

  /* See below for the link to the full explanation of how and why this works. */        
  do { 
    ur1 = 2.0 * (rand() % 10000) / 10000.0 - 1.0; 
    ur2 = 2.0 * (rand() % 10000) / 10000.0 - 1.0; 
    magnitude = ur1 * ur1 + ur2 * ur2;  
  } while ( magnitude >= 1.0 ); 
 
  magnitude = sqrt( (-2.0 * log( magnitude ) ) / magnitude ); 
 
  /* It's MAGIC!  Actually, check out: 
   * http://www.wikipedia.org/wiki/Box-Muller_transformation 
   * 
   * This is basically algorithm 2. 
   */ 
 
  return ur1 * magnitude; 
}

/****************************************************************************** 
 * 
 * random_normal_value 
 * 
 * Generates a random value based on the normal distribution centered at Avg 
 * with standard deviation std.  NOTE: Technically, a normal distribution 
 * extends to infinity in both directions.  However, to make this plausible in 
 * game terms we cap the distribution at 3 standard-deviations in either 
 * direction.  IE - if a magic missile is supposed to do an average of 12 
 * damage, with standard-deviation of 2, it will never do less than 12 - 2*3 = 6, 
 * and never do more than 12 + 2*3 = 18.  
 * 
 * avg - The average value we should give out. 
 * std - The standard deviation of the distribution. 
 * 
 * Returns - a random value. 
 *  
 ******************************************************************************/ 
 
float random_normal_value(float avg, float std) { 
 
  /* Get a normally distributed random number.  This is based on the standard 
   * normal, so it will be between -3 and 3 99% of the time, the average will 
   * be 0, and the standard deviation is 1. 
   */ 
  float val = box_muller_transformation(); 
 
  /* We cap this value at +-3.0.  If we don't, an occaisonal call to this could 
   * feasibly be 10+, which makes it really hard to deal with the formula. 
   * (IE - imagine a magic missile that once every 1000 casts or so does 1000 
   * damage... 
   */ 
  while (val > 3.0 || val < -3.0) 
    val = box_muller_transformation(); 
 
  /* Now, we have a standard-normal number. But we want to translate this into 
   * a normal with a given average and standard deviation.  The logic of the 
   * following is ass follows - the avg term merely shifts the whole to be 
   * centered around our arbitrary average.  The std * val is tricky - this 
   * may not actually be statisticaly valid, but it works well.  The idea is 
   * to scale the deviation from the average based on the supplied standard 
   * deviation.  
   * 
   * Example:  We want avg = 50, std = 10.  We grab a normally distributed 
   * value of say, .5.  This translates into 50 + 10 * .5 == 55.  This is 
   * one-half standard deviation right of the average of 50.  The original 
   * is also one-half standard-deviation to the right of the average (of 0) 
   * See? It's happy, and shiny, and good. 
   */ 
  return avg + std * val; 
}