1-- Title: Show Redshifts of Galaxies
2
3function get_distance(obj_pos, ref_pos)
4   -- Returns distance Earth-object in Mpc.
5   distance = ref_pos:distanceto(obj_pos) * km2Mpc
6   return distance
7end
8
9function get_z(distance)
10   -- Returns redshift z.
11   -- Hubble constant = 73.2 (km/s)/Mpc, WMAP 3years, best
12   -- http://map.gsfc.nasa.gov/m_mm/tp_links.html
13   local H0 = 73.2
14   local c = 299792.458
15   z0 = (H0 * distance)/c
16   z = z0*(1 - z0/4)/(1 - z0/2)^2
17   -- d <= 2 * c/H0
18   return z
19end
20
21function rgb2hex(r,g,b)
22   -- Converts color code from RGB to Hex
23   local hex = string.format("#%.2x%.2x%.2x", r,g,b)
24   return hex
25end
26
27function get_color(z_rel)
28   -- Returns Hex color code from z_rel
29   if z_rel > 1 then
30       z_rel = 1
31   end
32   --local green = 255 * (1 - z_rel)/(1 + math.sqrt(z_rel))
33   local green = 255 * (1 - z_rel)^2
34   local red   = 255
35   local blue  = 255* math.sqrt(green/255)
36   hex_color   = rgb2hex(red, green, blue)
37   return hex_color
38end
39
40function mark_galaxies(sel_pos)
41   zz = z_max(sel_pos)
42   for dso in celestia:dsos() do
43      if dso:type() == "galaxy" then
44         dso_pos = dso:getposition()
45         local d = get_distance(dso_pos,sel_pos)
46         local z_rel = get_z(d)/zz
47         local hex_color = get_color(z_rel)
48         dso:mark( hex_color, "disk", 1, 1 )
49      end
50   end
51end
52
53function z_max(ref_pos)
54   -- determine maximal redshift in catalog wrto ref_pos
55   z_old = 0
56   for dso in celestia:dsos() do
57      if dso:type() == "galaxy" then
58         dso_pos = dso:getposition()
59         local d = get_distance(dso_pos, ref_pos)
60         local z = get_z(d)
61         if z > z_old then
62             z_max = z
63             z_old = z_max
64             dsomax = dso
65         end
66      end
67   end
68   return z_max
69end
70
71----------
72-- main --
73----------
74celestia:unmarkall()
75
76celestia:show("markers")
77km2Mpc = 1/3.08568025e19
78MW = celestia:find("Milky Way")
79MW_pos = MW:getposition()
80--
81-- select and specially mark Milky Way
82--
83celestia:select(MW)
84celestia:mark(MW)
85--
86-- color encode all other galaxies according to their redshift
87-- relative to Milky Way
88--
89mark_galaxies(MW_pos)
90--
91-- move observer to a distance of 1000 Mpc from Milky Way
92--
93observer = celestia:getobserver()
94observer:gotodistance(MW, 1000/km2Mpc,5)
95
96sel_old = MW
97
98while true do
99   sel = celestia:getselection()
100   --
101   -- specially mark possible new selection, unmark the old one
102   --
103   if sel ~= sel_old then
104       sel:unmark()
105       sel_old:unmark()
106       sel:mark("#00FF00","disk",10,1)
107       sel_old = sel
108   end
109
110   sel_pos = sel:getposition()
111
112   -- obs_pos = observer:getposition()
113
114   if sel:type() == "galaxy" then
115      local d = get_distance(MW_pos,sel_pos)
116      local z = get_z(d)
117      celestia:print(sel:name()..": "..string.format("redshift z = %5.3f,  distance = %5.2f  Mpc", z,d).."\nmax. redshift: "..dsomax:name(),5,-1,-1,0,6)
118   end
119   wait(0)
120end