Created
May 1, 2012 08:59
-
-
Save edouardp/2566569 to your computer and use it in GitHub Desktop.
Transform Random and Halton Uniform Generators into Normal Distribution
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
class UniformGenerator: | |
def __init__(self): | |
self.uniform = random.uniform | |
def next(self): | |
return self.uniform(0,1) | |
def haltonterm(i, base=2): | |
h = 0 | |
fac = 1.0/base | |
while i != 0: | |
digit = i % base | |
h = h + digit*fac | |
i = (i-digit)/base | |
fac = fac/base | |
return h | |
class HaltonGenerator: | |
def __init__(self, index=1, base=2): | |
self.index = index | |
self.base = base | |
def next(self): | |
value = haltonterm(self.index, self.base) | |
self.index = self.index + 1 | |
return value | |
class HaltonGeneratorTwo: | |
def __init__(self, index=1, base_one=2, base_two=3): | |
self.index = index | |
self.base_one = base_one | |
self.base_two = base_two | |
self.cache = 0 | |
self.cache_full = False | |
def next(self): | |
if self.cache_full: | |
self.cache_full = False | |
return self.cache | |
else: | |
value = haltonterm(self.index, self.base_one) | |
self.cache = haltonterm(self.index, self.base_two) | |
self.cache_full = True | |
self.index = self.index + 1 | |
return value | |
class MarsagliaPolarGaussianGenerator: | |
def __init__(self, mu = 0, sigma = 1, uniform = UniformGenerator()): | |
self.mu = mu | |
self.sigma = sigma | |
self.uniform = uniform | |
self.cache = 0 | |
self.cache_full = False | |
def next(self): | |
if self.cache_full: | |
self.cache_full = False | |
return self.mu + self.sigma * self.cache | |
else: | |
while True: | |
u = self.uniform.next() * 2.0 - 1.0 | |
v = self.uniform.next() * 2.0 - 1.0 | |
s = u*u + v*v | |
if 0 < s < 1: | |
break; | |
self.cache = v * math.sqrt(-2.0 * math.log(s) / s) | |
self.cache_full = True | |
return self.mu + self.sigma * u * math.sqrt(-2.0 * math.log(s) / s ) |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment