From 3066899c079b56ca6ca3c5943dcf823a635c754d Mon Sep 17 00:00:00 2001 From: cecilkorik Date: Tue, 27 Apr 2010 03:43:05 +0000 Subject: [PATCH] initial import --HG-- branch : stargen --- rnd.py | 22 +++++++ stargen.py | 176 +++++++++++++++++++++++++++++++++++++++++++++++++ stargen.py.bak | 174 ++++++++++++++++++++++++++++++++++++++++++++++++ 3 files changed, 372 insertions(+) create mode 100644 rnd.py create mode 100644 stargen.py create mode 100644 stargen.py.bak diff --git a/rnd.py b/rnd.py new file mode 100644 index 0000000..885416f --- /dev/null +++ b/rnd.py @@ -0,0 +1,22 @@ +import random as random_lib +import math + +def random(): + return random_lib.random() + +def gaussian(mean=0.0, stddev=1.0): + while True: + x1 = (2.0 * random()) - 1.0 + x2 = (2.0 * random()) - 1.0 + w = (x1 * x1) + (x2 * x2) + if w < 1.0: + break + + w = math.sqrt((-2.0 * math.log(w)) / w) + + return mean + (x1 * w * stddev) + + + +def randrange(min, max): + return min + (random() * (max - min)) \ No newline at end of file diff --git a/stargen.py b/stargen.py new file mode 100644 index 0000000..a2008f3 --- /dev/null +++ b/stargen.py @@ -0,0 +1,176 @@ +from rnd import * +import math + + +systemsize = 6e9 +stars = (1.0, 0.5) +planets = (8.0, 4.0) + +star_spectra = ['O', 'B', 'A', 'F', 'G', 'K', 'M'] +star_lumin = [(500000.0, 30000.0), (30000.0, 25.0), (25.0, 5.0), (5.0, 1.5), (1.5, 0.6), (0.6, 0.008), (0.08, 0.0008)] +star_masses = [(40.0, 16.0), (16.0, 2.1), (2.1, 1.4), (1.4, 1.04), (1.04, 0.8), (0.8, 0.45), (0.45, 0.1)] +star_radii = [(22.0, 8.0), (7.0, 1.8), (1.8, 1.4), (1.4, 1.15), (1.15, 0.96), (0.96, 0.7), (0.7, 0.25)] +star_prob = [0.0, 0.05, 0.35, 0.6, 0.8, 0.97, 1.0] + +num_stars = round(gaussian(*stars)) +if num_stars < 1.0: + num_stars = 1.0 +num_stars = int(num_stars) + +num_planets = round(gaussian(*planets)) +if num_planets < 0.0: + num_planets = 0.0 +num_planets = int(num_planets) + +stars = [] +x = y = z = 0.0 + +while len(stars) < num_stars: + r = random() + for i, e in enumerate(star_prob): + if r <= e: + si = i + break + r = random() + slumin = (star_lumin[si][1] * 0.8) + (star_masses[si][0] * 1.2 * r) + r2 = r + (random() * 0.3) - 0.15 + smass = (star_masses[si][1] * 0.8) + (star_masses[si][0] * 1.2 * r) + r2 = r + (random() * 0.3) - 0.15 + sradius = (star_radii[si][1] * 0.8) + (star_radii[si][0] * 1.2 * r) + print "Star lumin %s mass %s class %s" % (slumin, smass, star_spectra[si]) + + stars += [(slumin, smass, star_spectra[si], sradius, x, y, z)] + + # first star is always at 0.0, 0.0, 0.0, next star gets random coords + x = gaussian(0.0, systemsize) + y = gaussian(0.0, systemsize) + if len(stars) > 1: + z = gaussian(0.0, systemsize / 10e3) + +4 + +planets = [] +sys_size = math.sqrt(sys_lumin) + +while len(planets) < num_planets: + def dist(x, y, z): + return math.sqrt(x**2 + y**2 + z**2) + + x = gaussian(0.0, systemsize * sys_size) + y = gaussian(0.0, systemsize * sys_size) + z = gaussian(0.0, systemsize * sys_size / 10e3) + + reject = False + for star in stars: + slumin, smass, sclass, sradius, sx, sy, sz + d = dist(x-sx,y-sy,z-sz) + + deathzone = sradius + (sradius * 0.25 * math.sqrt(slumin)) + if d < deathzone: + reject = True + + if reject: + continue + + r = random() + rocky = False + gas_giant = False + if r < 0.6: + # rock planet + rocky = True + massmin = 1.0e19 + massmax = 8.8e25 + mass = 0.0 + while mass < massmin or mass > massmax: + massexp = gaussian(24.5, 1.25) + mass = 10.0 ** massexp + density = 0.0 + while density < 0.75: + density = gaussian(4.5, 1.0) + mag = -1.0 + while mag < 0.0: + mag = gaussian(density, 3.0) + + else: + # gas giant + gas_giant = True + massmin = 1e25 + massmax = 4e28 + mass = 0.0 + while mass < massmin or mass > massmax: + massexp = gaussian(26.5, 1.0) + mass = 10.0 ** massexp + density = 0.0 + while density < 0.3: + density = gaussian(1.5, 0.75) + mag = -1.0 + while mag < 2.0: + mag = gaussian(7.0, 3.0) + + + volume = mass / density / 1e12 + radius = pow(volume * 0.75 / math.pi, 1.0/3.0) + circsurfarea = math.pi * radius**2 + surfarea = circsurfarea * 4.0 + insol = 0.0 + for star in stars: + slumin, smass, sclass, sradius, sx, sy, sz + d = dist(x-sx,y-sy,z-sz) + + spower = 3.8e26 * slumin + distsphere = 4.0 * math.pi * d**2 + insol += (circsurfarea / distsphere) * spower + + areainsol = insol / surfarea / 1e6 + if areainsol > 75000: + continue + + atmo = None + while atmo == None or atmo > density: + if gas_giant: + atmo_exp = randrange(-3, 0.5) + else: + atmo_exp = randrange(-9.2, -0.9) + atmo = 10 ** atmo_exp + + albedo = -1.0 + while albedo < 0.0 or albedo > 1.0: + albedo = gaussian(0.35, 0.05) + + emissivity = 1.0 - albedo + emissivity = randrange(0.00001, 1.0) + + emiss_diff = (((10.0 - pow(atmo, 1.0/4.0)) / 10.0) ** 100.0) / 2.0 + emissivity_min = emissivity - (emissivity * emiss_diff) + emissivity_max = emissivity + (emissivity * emiss_diff) + + avg_temp = pow(((1.0 - albedo) * areainsol) / (4.0 * emissivity * 5.67e-8), 1.0/4.0) + min_temp = pow(((1.0 - albedo) * areainsol) / (4.0 * emissivity_min * 5.67e-8), 1.0/4.0) + max_temp = pow(((1.0 - albedo) * areainsol) / (4.0 * emissivity_max * 5.67e-8), 1.0/4.0) + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/stargen.py.bak b/stargen.py.bak new file mode 100644 index 0000000..69cee92 --- /dev/null +++ b/stargen.py.bak @@ -0,0 +1,174 @@ +from rnd import * +import math + + +systemsize = 6e9 +stars = (1.0, 0.5) +planets = (8.0, 4.0) + +star_spectra = ['O', 'B', 'A', 'F', 'G', 'K', 'M'] +star_lumin = [(500000.0, 30000.0), (30000.0, 25.0), (25.0, 5.0), (5.0, 1.5), (1.5, 0.6), (0.6, 0.008), (0.08, 0.0008)] +star_masses = [(40.0, 16.0), (16.0, 2.1), (2.1, 1.4), (1.4, 1.04), (1.04, 0.8), (0.8, 0.45), (0.45, 0.1)] +star_radii = [(22.0, 8.0), (7.0, 1.8), (1.8, 1.4), (1.4, 1.15), (1.15, 0.96), (0.96, 0.7), (0.7, 0.25)] +star_prob = [0.0, 0.05, 0.35, 0.6, 0.8, 0.97, 1.0] + +num_stars = round(gaussian(*stars)) +if num_stars < 1.0: + num_stars = 1.0 +num_stars = int(num_stars) + +num_planets = round(gaussian(*planets)) +if num_planets < 0.0: + num_planets = 0.0 +num_planets = int(num_planets) + +stars = [] +x = y = z = 0.0 + +while len(stars) < num_stars: + r = random() + for i, e in enumerate(star_prob): + if r <= e: + si = i + break + r = random() + slumin = (star_lumin[si][1] * 0.8) + (star_masses[si][0] * 1.2 * r) + r2 = r + (random() * 0.3) - 0.15 + smass = (star_masses[si][1] * 0.8) + (star_masses[si][0] * 1.2 * r) + r2 = r + (random() * 0.3) - 0.15 + sradius = (star_radii[si][1] * 0.8) + (star_radii[si][0] * 1.2 * r) + print "Star lumin %s mass %s class %s" % (slumin, smass, star_spectra[si]) + + stars += [(slumin, smass, star_spectra[si], sradius, x, y, z)] + + # first star is always at 0.0, 0.0, 0.0, next star gets random coords + x = gaussian(0.0, systemsize) + y = gaussian(0.0, systemsize) + if len(stars) > 1: + z = gaussian(0.0, systemsize / 10e3) + +planets = [] +sys_size = math.sqrt(sys_lumin) + +while len(planets) < num_planets: + def dist(x, y, z): + return math.sqrt(x**2 + y**2 + z**2) + + x = gaussian(0.0, systemsize * sys_size) + y = gaussian(0.0, systemsize * sys_size) + z = gaussian(0.0, systemsize * sys_size / 10e3) + + reject = False + for star in stars: + slumin, smass, sclass, sradius, sx, sy, sz + d = dist(x-sx,y-sy,z-sz) + + deathzone = sradius + (sradius * 0.25 * math.sqrt(slumin)) + if d < deathzone: + reject = True + + if reject: + continue + + r = random() + rocky = False + gas_giant = False + if r < 0.6: + # rock planet + rocky = True + massmin = 1.0e19 + massmax = 8.8e25 + mass = 0.0 + while mass < massmin or mass > massmax: + massexp = gaussian(24.5, 1.25) + mass = 10.0 ** massexp + density = 0.0 + while density < 0.75: + density = gaussian(4.5, 1.0) + mag = -1.0 + while mag < 0.0: + mag = gaussian(density, 3.0) + + else: + # gas giant + gas_giant = True + massmin = 1e25 + massmax = 4e28 + mass = 0.0 + while mass < massmin or mass > massmax: + massexp = gaussian(26.5, 1.0) + mass = 10.0 ** massexp + density = 0.0 + while density < 0.3: + density = gaussian(1.5, 0.75) + mag = -1.0 + while mag < 2.0: + mag = gaussian(7.0, 3.0) + + + volume = mass / density / 1e12 + radius = pow(volume * 0.75 / math.pi, 1.0/3.0) + circsurfarea = math.pi * radius**2 + surfarea = circsurfarea * 4.0 + insol = 0.0 + for star in stars: + slumin, smass, sclass, sradius, sx, sy, sz + d = dist(x-sx,y-sy,z-sz) + + spower = 3.8e26 * slumin + distsphere = 4.0 * math.pi * d**2 + insol += (circsurfarea / distsphere) * spower + + areainsol = insol / surfarea / 1e6 + if areainsol > 75000: + continue + + atmo = None + while atmo == None or atmo > density: + if gas_giant: + atmo_exp = randrange(-3, 0.5) + else: + atmo_exp = randrange(-9.2, -0.9) + atmo = 10 ** atmo_exp + + albedo = -1.0 + while albedo < 0.0 or albedo > 1.0: + albedo = gaussian(0.35, 0.05) + + emissivity = 1.0 - albedo + emissivity = randrange(0.00001, 1.0) + + emiss_diff = (((10.0 - pow(atmo, 1.0/4.0)) / 10.0) ** 100.0) / 2.0 + emissivity_min = emissivity - (emissivity * emiss_diff) + emissivity_max = emissivity + (emissivity * emiss_diff) + + avg_temp = pow(((1.0 - albedo) * areainsol) / (4.0 * emissivity * 5.67e-8), 1.0/4.0) + min_temp = pow(((1.0 - albedo) * areainsol) / (4.0 * emissivity_min * 5.67e-8), 1.0/4.0) + max_temp = pow(((1.0 - albedo) * areainsol) / (4.0 * emissivity_max * 5.67e-8), 1.0/4.0) + + + + + + + + + + + + + + + + + + + + + + + + + + +