/*   Generate and plot a Gaussian Distribution   */
         
#declare Numbarands = 10000;       // Generate this number of Gaussianly distributed values
#declare Hyte = Numbarands/15;     // A scaling value for plotting the results
#declare Cell = array[12]          // Sorting cells to measure the distribution
#declare i = 1;                    // initialize the cells
#while (i<=11 )
  #declare Cell[i] = 0;
  #declare i = i + 1; 
#end
#warning "\n" 
 
#declare i =1;                    // a loop to generate the Gaussian values 
#while ( i<= Numbarands)
  #declare sd1 = seed(i*2123);    // a seed to generate a uniformly distributed value in the range 0<x<1
  #declare sd2 = seed(3212*i);    // a seed to generate a second 
  #declare rx = rand(sd1); #declare ry = rand(sd2);   // generate the uniformly distributed values

   // This the Boxmueller algorithm  
      #declare Arg1   = -2*log(rx);
      #declare Arg2   = sin(2*pi*ry);
      #declare Gnumba = sqrt(Arg1)*Arg2/3;

  #warning str(Gnumba,9,3)        // print and sort the values
  #switch ( Gnumba )
    #range(-1.1,-0.9)  #declare Cell[ 1] = Cell[ 1]+1;   #break
    #range(-0.9,-0.7)  #declare Cell[ 2] = Cell[ 2]+1;   #break 
    #range(-0.7,-0.5)  #declare Cell[ 3] = Cell[ 3]+1;   #break
    #range(-0.5,-0.3)  #declare Cell[ 4] = Cell[ 4]+1;   #break
    #range(-0.3,-0.1)  #declare Cell[ 5] = Cell[ 5]+1;   #break
    #range(-0.1, 0.1)  #declare Cell[ 6] = Cell[ 6]+1;   #break 
    #range( 0.1, 0.3)  #declare Cell[ 7] = Cell[ 7]+1;   #break 
    #range( 0.3, 0.5)  #declare Cell[ 8] = Cell[ 8]+1;   #break
    #range( 0.5, 0.7)  #declare Cell[ 9] = Cell[ 9]+1;   #break
    #range( 0.7, 0.9)  #declare Cell[10] = Cell[10]+1;   #break
    #range( 0.9, 1.1)  #declare Cell[11] = Cell[11]+1;   #break
  #end
  #declare j =1;                 
    #while ( j<=11)  
      #warning str(Cell[j],9,0)
      #declare j = j + 1;     
    #end
    #warning "\n"
  #declare i = i + 1;
#end 
  
#declare j=1;               // generate and visualize a set of boxes whose heights represent the distribution
#while ( j<=10) 
  box{ <-0.15,0,-0.15>,<0.15,5*Cell[j]/Hyte,0.15> pigment { rgb<1,1,0> } finish {ambient 1 }
          translate <-j+5,-10,1>}
  #declare j=j+1;
#end   

camera { location -z*120 up y right x*1.33 look_at 0 angle 15}