1 /**************************************************************************
4 * Written 1997 by Durk Talsma, started October, 1997. For the flight gear
7 * This program is free software; you can redistribute it and/or
8 * modify it under the terms of the GNU General Public License as
9 * published by the Free Software Foundation; either version 2 of the
10 * License, or (at your option) any later version.
12 * This program is distributed in the hope that it will be useful, but
13 * WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
15 * General Public License for more details.
17 * You should have received a copy of the GNU General Public License
18 * along with this program; if not, write to the Free Software
19 * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
22 * (Log is kept at end of this file)
23 **************************************************************************/
37 #include <Debug/fg_debug.h>
38 #include <Include/fg_constants.h>
39 #include <Time/fg_time.hxx>
40 #include <Time/light.hxx>
43 #include "planets.hxx"
48 void fgCalculatePlanet(struct OrbElements planet, struct OrbElements theSun,
49 fgTIME t, int idx, struct CelestialCoord *result)
51 // struct CelestialCoord result;
53 fgSUNPOS SolarPosition;
55 double eccAnom, r, v, ecl, actTime, R, s, ir, Nr, B, FV, ring_magn,
56 xv, yv, xh, yh, zh, xg, yg, zg, xe, ye, ze;
58 actTime = fgCalcActTime(t);
59 /*calculate the angle between ecliptic and equatorial coordinate system */
60 ecl = DEG_TO_RAD * (23.4393 - 3.563E-7 * actTime);
62 /* calculate the eccentric anomaly */
63 eccAnom = fgCalcEccAnom(planet.M, planet.e);
65 /* calculate the planets distance (r) and true anomaly (v) */
66 xv = planet.a * (cos(eccAnom) - planet.e);
67 yv = planet.a * (sqrt(1.0 - planet.e*planet.e) * sin(eccAnom));
69 r = sqrt ( xv*xv + yv*yv);
71 /* calculate the planets position in 3-dimensional space */
72 xh = r * (cos (planet.N) * cos (v + planet.w) -
73 sin (planet.N) * sin (v + planet.w) * cos (planet.i));
74 yh = r * (sin (planet.N) * cos (v + planet.w) +
75 cos (planet.N) * sin (v + planet.w) * cos (planet.i));
76 zh = r * ( sin(v+planet.w) * sin(planet.i));
78 /* calculate the ecliptic longitude and latitude */
80 xg = xh + solarPosition.xs;
81 yg = yh + solarPosition.ys;
85 ye = yg * cos(ecl) - zg * sin(ecl);
86 ze = yg * sin(ecl) + zg * cos(ecl);
88 result->RightAscension = atan2(ye,xe);
89 result->Declination = atan2(ze, sqrt(xe*xe + ye*ye));
91 /* Let's calculate the brightness of the planet */
92 R = sqrt ( xg*xg + yg*yg + zg*zg);
93 s = SolarPosition.dist;
94 FV = acos( (r*r + R*R - s*s) / (2*r*R));
95 FV *= 57.29578; /* convert radians to degrees */
99 result->magnitude = -0.36 + 5*log10( r*R ) + 0.027 * FV + 2.2E-13 * pow(FV, 6);
102 result->magnitude = -4.34 + 5*log10( r*R ) + 0.013 * FV + 4.2E-07 * pow(FV,3);
105 result->magnitude = -1.51 + 5*log10( r*R ) + 0.016 * FV;
107 case 5: /* Jupiter */
108 result->magnitude = -9.25 + 5*log10( r*R ) + 0.014 * FV;
113 Nr = 2.9585076 + 6.6672E-7*actTime;
115 B = asin (sin (result->Declination) * cos (ir) -
116 cos (result->Declination) * sin (ir) *
117 sin (result->RightAscension - Nr));
118 ring_magn = -2.6 * sin (fabs(B)) + 1.2 * pow(sin(B),2);
119 result->magnitude = -9.0 + 5*log10( r*R ) + 0.044 * FV + ring_magn;
122 result->magnitude = -7.15 + 5*log10( r*R) + 0.001 * FV;
124 case 8: /* Neptune */
125 result->magnitude = -6.90 + 5*log10 (r*R) + 0.001 *FV;
128 fgPrintf( FG_ASTRO, FG_ALERT, "index %d out of range !!!!\n", idx);
130 fgPrintf( FG_ASTRO, FG_DEBUG,
131 " Planet found at %f (ra), %f (dec)\n",
132 result->RightAscension, result->Declination);
133 fgPrintf( FG_ASTRO, FG_DEBUG,
134 " Geocentric dist %f\n"
135 " Heliocentric dist %f\n"
136 " Distance to the sun %f\n"
138 " Brightness %f\n", R, r, s, FV, result->magnitude);
144 void fgPlanetsInit( void )
148 struct CelestialCoord pltPos;
151 l = &cur_light_params;
152 /* if the display list was already built during a previous init,
156 xglDeleteLists(planets, 1);
159 planets = xglGenLists(1);
160 xglNewList( planets, GL_COMPILE );
161 xglBegin( GL_POINTS );
164 /* Add the planets to all four display lists */
165 for ( i = 2; i < 9; i++ ) {
166 fgSolarSystemUpdate(&(pltOrbElements[i]), cur_time_params);
167 fgCalculatePlanet(pltOrbElements[i], pltOrbElements[0],
168 cur_time_params, i, &pltPos);
170 magnitude = (0.0 - pltPos.magnitude) / 5.0 + 1.0;
172 /* scale magnitudes again so they look ok */
174 magnitude * 0.7 + (((FG_STAR_LEVELS - 1) - i) * 0.1); */
176 /* the following statement could be made a little more sopisticated
177 for the moment: Stick to this one */
178 magnitude = magnitude * 0.7 + (3 * 0.1);
179 if ( magnitude > 1.0 ) { magnitude = 1.0; }
180 if ( magnitude < 0.0 ) { magnitude = 0.0; }
182 /* Add planets to the display list, based on sun_angle and current
183 magnitude It's pretty experimental... */
185 if ((double) (l->sun_angle - FG_PI_2) >
186 ((magnitude - 1.0) * -20 * DEG_TO_RAD))
188 xglColor3f (magnitude, magnitude, magnitude);
189 printf ("Sun Angle to Horizon (in Rads) = %f\n",
190 (double) (l->sun_angle - FG_PI_2));
191 printf ("Transformed Magnitude is :%f %f\n",
192 magnitude, (magnitude - 1.0) * -20 * DEG_TO_RAD);
194 xglVertex3f (50000.0 * cos (pltPos.RightAscension) *
195 cos (pltPos.Declination),
196 50000.0 * sin (pltPos.RightAscension) *
197 cos (pltPos.Declination),
198 50000.0 * sin (pltPos.Declination));
207 void fgPlanetsRender( void ) {
208 xglCallList(planets);
213 /* Revision 1.6 1998/06/27 16:51:54 curt
214 /* In fgCalculatePlanet() pass a pointer to a structure to be modified, rather
215 /* than returning the entire structure.
217 * Revision 1.5 1998/04/28 01:19:01 curt
218 * Type-ified fgTIME and fgVIEW
220 * Revision 1.4 1998/04/26 05:10:01 curt
221 * "struct fgLIGHT" -> "fgLIGHT" because fgLIGHT is typedef'd.
223 * Revision 1.3 1998/04/25 22:06:25 curt
224 * Edited cvs log messages in source files ... bad bad bad!
226 * Revision 1.2 1998/04/24 00:45:02 curt
227 * Wrapped "#include <config.h>" in "#ifdef HAVE_CONFIG_H"
228 * Fixed a bug when generating sky colors.
230 * Revision 1.1 1998/04/22 13:21:31 curt
231 * C++ - ifing the code a bit.
233 * Revision 1.9 1998/04/18 04:13:57 curt
234 * Moved fg_debug.c to it's own library.
236 * Revision 1.8 1998/04/03 21:52:50 curt
237 * Converting to Gnu autoconf system.
239 * Revision 1.7 1998/02/23 19:07:55 curt
240 * Incorporated Durk's Astro/ tweaks. Includes unifying the sun position
241 * calculation code between sun display, and other FG sections that use this
242 * for things like lighting.
244 * Revision 1.6 1998/02/12 21:59:36 curt
245 * Incorporated code changes contributed by Charlie Hotchkiss
246 * <chotchkiss@namg.us.anritsu.com>
248 * Revision 1.5 1998/02/03 23:20:12 curt
249 * Lots of little tweaks to fix various consistency problems discovered by
250 * Solaris' CC. Fixed a bug in fg_debug.c with how the fgPrintf() wrapper
251 * passed arguments along to the real printf(). Also incorporated HUD changes
252 * by Michele America.
254 * Revision 1.4 1998/02/02 20:53:23 curt
257 * Revision 1.3 1998/01/27 00:47:47 curt
258 * Incorporated Paul Bleisch's <pbleisch@acm.org> new debug message
259 * system and commandline/config file processing code.
261 * Revision 1.2 1998/01/19 19:26:59 curt
262 * Merged in make system changes from Bob Kuehne <rpk@sgi.com>
263 * This should simplify things tremendously.
265 * Revision 1.1 1998/01/07 03:16:18 curt
266 * Moved from .../Src/Scenery/ to .../Src/Astro/
268 * Revision 1.4 1997/12/30 20:47:52 curt
269 * Integrated new event manager with subsystem initializations.
271 * Revision 1.3 1997/12/30 16:36:52 curt
272 * Merged in Durk's changes ...
274 * Revision 1.2 1997/12/12 21:41:29 curt
275 * More light/material property tweaking ... still a ways off.
277 * Revision 1.1 1997/10/25 03:16:10 curt
278 * Initial revision of code contributed by Durk Talsma.