DSPRelated.com
Forums

gaussian white noise generation

Started by sita January 18, 2006
Gordon Sande wrote:

   ...

> Ask google about Box-Muller transformation. Turns uniforms into > Gaussians exactly. There ae a couple related transforms that > are also possible. Then there are the various schemes put forward > by Marsaglia.
Thanks for contributing to my education. I think "Gaussian exactly" is a bit strong, though. Few computers can exactly represent the tails. :-) From http://mathworld.wolfram.com/Box-MullerTransformation.html: Eric W. Weisstein, Box-Muller Transformation. A transformation which transforms from a two-dimensional continuous uniform distribution to a two-dimensional bivariate normal distribution (or complex normal distribution). If x_1 and x_2 are uniformly and independently distributed between 0 and 1, then z_1 and z_2 as defined below have a normal distribution with mean mu==0 and variance sigma^2==1. z_1 = sqrt(-2lnx_1)cos(2pix_2) z_2 = sqrt(-2lnx_1)sin(2pix_2). Is seems that it needs floating point, two calls to random(), a log, and a sine (or cosine). In some cases, half a dozen or so calls to an integer random() and as many adds may be preferable. Whatever, it's real neat. I'll remember it. Jerry -- Engineering is the art of making what you want from things you can get. �����������������������������������������������������������������������
Naebad wrote:
> "Jerry Avins" <jya@ieee.org> wrote in message > news:wdqdnQYS9OboxVPeRVn-qg@rcn.net... > >>Naebad wrote: > > >>What is the DC level of a single average? Assuming that the RNG is >>itself zero mean, you can omit the subtraction step. >> >>Jerry >>-- >>Engineering is the art of making what you want from things you can get. >> > > > I think it's 6 but it's so long since I did this...
It depends on the range of the averaged random variables. If they run from -1 to 1, then the mean is zero if 0 t0 1, than half the number averaged. Jerry -- Engineering is the art of making what you want from things you can get. &#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;
On 2006-01-18 16:31:48 -0400, Jerry Avins <jya@ieee.org> said:

> Gordon Sande wrote: > > ... > >> Ask google about Box-Muller transformation. Turns uniforms into >> Gaussians exactly. There ae a couple related transforms that >> are also possible. Then there are the various schemes put forward >> by Marsaglia. > > Thanks for contributing to my education.
> I think "Gaussian exactly" is a bit strong, though.
The Box-Muller transformation is exact, unlike the piecewise polynomial that arises from summing 12 uniforms and subtracting 6. It is zero for absolute values above 6.
> Few computers can exactly represent the tails. :-)
If the "uniforms" are not exact because they are discrete then the resulting "Gaussians" will be discrete as well. The discreteness is not caused by the transformation. If you choose to keep only the Gaussian tails than you are doing an experiment that is very prone to expose the discreteness of the uniforms. If you keep the middle of the Gaussians then you will not notice the discreteness. The issue is quite general as it also applies to exponential variables that might be used for low rate Poisson processes with a long time between events. Etc. Etc. If you want to chase the Gaussian tails then you had better get some multiple precision software and start with a very high precision uniform such as the Mesenne Twister or one of the other high order recurrence generators.
Gordon Sande wrote:
> On 2006-01-18 16:31:48 -0400, Jerry Avins <jya@ieee.org> said: > >> Gordon Sande wrote: >> >> ... >> >>> Ask google about Box-Muller transformation. Turns uniforms into >>> Gaussians exactly. There ae a couple related transforms that >>> are also possible. Then there are the various schemes put forward >>> by Marsaglia. >> >> >> Thanks for contributing to my education. > > > >> I think "Gaussian exactly" is a bit strong, though. > > > The Box-Muller transformation is exact, unlike the piecewise > polynomial that arises from summing 12 uniforms and subtracting 6. > It is zero for absolute values above 6. > >> Few computers can exactly represent the tails. :-) > > > If the "uniforms" are not exact because they are discrete then > the resulting "Gaussians" will be discrete as well. The discreteness > is not caused by the transformation. If you choose to keep only > the Gaussian tails than you are doing an experiment that is > very prone to expose the discreteness of the uniforms. If you > keep the middle of the Gaussians then you will not notice the > discreteness. The issue is quite general as it also applies to > exponential variables that might be used for low rate Poisson > processes with a long time between events. Etc. Etc. > > If you want to chase the Gaussian tails then you had better get > some multiple precision software and start with a very high precision > uniform such as the Mesenne Twister or one of the other high order > recurrence generators.
Of course you're right. I tend to bridle at the term "exact". Someone who asked for a hole in the exact center of something, refusing to give a tolerance because any error is too large, sensitized me. (Told that the best my setup could do was .2 mil, he asked what a mil was. When told, he settled for 1 mil, the precision of my three-jaw chuck.) Jerry -- Engineering is the art of making what you want from things you can get. &#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;&#4294967295;
>>>>> "Gordon" == Gordon Sande <g.sande@worldnet.att.net> writes:
Gordon> Ask google about Box-Muller transformation. Turns uniforms into Gordon> Gaussians exactly. There ae a couple related transforms that Gordon> are also possible. Then there are the various schemes put forward Gordon> by Marsaglia. I thought I'd add a few of the other relatively simple methods for generating Gaussians. 1. Box-Muller using trig functions already given. 2. Box-Muller using rejection technique: a. u1, u2 = uniform random variable in [-1, 1] b. w = u1^2 + u2^2 c. Repeat a-b until w < 1. d. s = sqrt(-2*log(w)/w) e. x1 = s*u1, x2 = s*u2 This is faster than 1 if your trig functions are slow. 3. Ahrens, NA algorithm. a. u = random [0,1] b. B = 0 if u < 1/2, otherwise B = 1 c. Generate an exponential variate E. (One method is -log(u2), u2 being uniform [0,1], but there are many others.) d. S = E + E. e. Generate a Cauchy variate, C. (One method is tan(pi*(u2-1/2)), where u2 is uniform [0,1], but it's possible to use the bits of u above for this.) f. x = sqrt(S/(1+C*C)). g. Y = C*X. h. If B = 0, return X, otherwise return -X. On the next call, we can just return Y. 4. Ratio of uniforms. Knuth describes this. a. u = uniform [0, 1], u2 = uniform [0,1] b. x = sqrt(8/e)*(u2-0.5)/u c. xs = x^2 d. if xs < 1/2 - u*(4*exp(-1/4)) return e. if (1.4 + 4*exp(-1.35)/u) <= xs <= -4*log(u) then return x. f. Otherwise repeat. 5. Marsaglia's Ziggurat method. This is a bit complex, so I won't describe it here. Of these, I've found Marsaglia's method to be the fastest, but it does require storage of 127 numbers for a table. Ray