1 // environment.cxx -- routines to model the natural environment
3 // Written by David Megginson, started February 2002.
5 // Copyright (C) 2002 David Megginson - david@megginson.com
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.
36 #include <simgear/constants.h>
37 #include <simgear/debug/logstream.hxx>
38 #include <simgear/math/interpolater.hxx>
40 #include <Main/fg_props.hxx>
41 #include "environment.hxx"
44 FGEnvironment::FGEnvironment()
45 : _temperature_degc_table(new SGInterpTable),
46 _pressure_inhg_table(new SGInterpTable),
49 temperature_sea_level_degc(15),
51 dewpoint_sea_level_degc(5), // guess
53 pressure_sea_level_inhg(29.92),
55 wind_from_heading_deg(0),
57 wind_from_north_fps(0),
58 wind_from_east_fps(0),
64 FGEnvironment::FGEnvironment (const FGEnvironment &env)
65 : _temperature_degc_table(new SGInterpTable),
66 _pressure_inhg_table(new SGInterpTable),
67 elevation_ft(env.elevation_ft),
68 visibility_m(env.visibility_m),
69 temperature_sea_level_degc(env.temperature_sea_level_degc),
70 dewpoint_sea_level_degc(env.dewpoint_sea_level_degc),
71 pressure_sea_level_inhg(env.pressure_sea_level_inhg),
72 wind_from_heading_deg(env.wind_from_heading_deg),
73 wind_speed_kt(env.wind_speed_kt),
74 wind_from_north_fps(env.wind_from_north_fps),
75 wind_from_east_fps(env.wind_from_east_fps),
76 wind_from_down_fps(env.wind_from_down_fps)
81 FGEnvironment::~FGEnvironment()
83 delete _temperature_degc_table;
84 delete _pressure_inhg_table;
89 FGEnvironment::get_visibility_m () const
95 FGEnvironment::get_temperature_sea_level_degc () const
97 return temperature_sea_level_degc;
101 FGEnvironment::get_temperature_degc () const
103 return temperature_degc;
107 FGEnvironment::get_dewpoint_sea_level_degc () const
109 return dewpoint_sea_level_degc;
113 FGEnvironment::get_dewpoint_degc () const
115 return dewpoint_degc;
119 FGEnvironment::get_pressure_sea_level_inhg () const
121 return pressure_sea_level_inhg;
125 FGEnvironment::get_pressure_inhg () const
127 return pressure_inhg;
131 FGEnvironment::get_density_slugft3 () const
133 return density_slugft3;
137 FGEnvironment::get_wind_from_heading_deg () const
139 return wind_from_heading_deg;
143 FGEnvironment::get_wind_speed_kt () const
145 return wind_speed_kt;
149 FGEnvironment::get_wind_from_north_fps () const
151 return wind_from_north_fps;
155 FGEnvironment::get_wind_from_east_fps () const
157 return wind_from_east_fps;
161 FGEnvironment::get_wind_from_down_fps () const
163 return wind_from_down_fps;
167 FGEnvironment::get_elevation_ft () const
173 FGEnvironment::set_visibility_m (double v)
179 FGEnvironment::set_temperature_sea_level_degc (double t)
181 temperature_sea_level_degc = t;
182 _recalc_alt_temperature();
187 FGEnvironment::set_temperature_degc (double t)
189 temperature_degc = t;
190 _recalc_sl_temperature();
195 FGEnvironment::set_dewpoint_sea_level_degc (double t)
197 dewpoint_sea_level_degc = t;
198 _recalc_alt_dewpoint();
203 FGEnvironment::set_dewpoint_degc (double t)
206 _recalc_sl_dewpoint();
211 FGEnvironment::set_pressure_sea_level_inhg (double p)
213 pressure_sea_level_inhg = p;
214 _recalc_alt_pressure();
219 FGEnvironment::set_pressure_inhg (double p)
222 _recalc_sl_pressure();
227 FGEnvironment::set_wind_from_heading_deg (double h)
229 wind_from_heading_deg = h;
234 FGEnvironment::set_wind_speed_kt (double s)
241 FGEnvironment::set_wind_from_north_fps (double n)
243 wind_from_north_fps = n;
248 FGEnvironment::set_wind_from_east_fps (double e)
250 wind_from_east_fps = e;
255 FGEnvironment::set_wind_from_down_fps (double d)
257 wind_from_down_fps = d;
262 FGEnvironment::set_elevation_ft (double e)
265 _recalc_alt_temperature();
266 _recalc_alt_dewpoint();
267 _recalc_alt_pressure();
273 // Copied from YASim Atmosphere.cxx, with m converted to ft, degK
274 // converted to degC, Pa converted to inHG, and kg/m^3 converted to
275 // slug/ft^3; they were then converted to deltas from the sea-level
276 // defaults (approx. 15degC, 29.92inHG, and 0.00237slugs/ft^3).
278 // Original comment from YASim:
280 // Copied from McCormick, who got it from "The ARDC Model Atmosphere"
281 // Note that there's an error in the text in the first entry,
282 // McCormick lists 299.16/101325/1.22500, but those don't agree with
283 // R=287. I chose to correct the temperature to 288.20, since 79F is
284 // pretty hot for a "standard" atmosphere.
286 // Elevation (ft), temperature factor (degK), pressure factor (inHG)
287 static double atmosphere_data[][3] = {
289 2952.76, 0.98, 0.898,
290 5905.51, 0.96, 0.804,
291 8858.27, 0.94, 0.719,
292 11811.02, 0.92, 0.641,
293 14763.78, 0.90, 0.570,
294 17716.54, 0.88, 0.506,
295 20669.29, 0.86, 0.447,
296 23622.05, 0.84, 0.394,
297 26574.80, 0.82, 0.347,
298 29527.56, 0.80, 0.304,
299 32480.31, 0.78, 0.266,
300 35433.07, 0.76, 0.231,
301 38385.83, 0.75, 0.201,
302 41338.58, 0.75, 0.174,
303 44291.34, 0.75, 0.151,
304 47244.09, 0.75, 0.131,
305 50196.85, 0.75, 0.114,
306 53149.61, 0.75, 0.099,
307 56102.36, 0.75, 0.086,
308 59055.12, 0.75, 0.075,
309 62007.87, 0.75, 0.065,
314 FGEnvironment::_setup_tables ()
316 for (int i = 0; atmosphere_data[i][0] != -1; i++) {
317 _temperature_degc_table->addEntry(atmosphere_data[i][0],
318 atmosphere_data[i][1]);
319 _pressure_inhg_table->addEntry(atmosphere_data[i][0],
320 atmosphere_data[i][2]);
325 FGEnvironment::_recalc_hdgspd ()
329 if (wind_from_east_fps == 0) {
330 angle_rad = (wind_from_north_fps >= 0 ? SGD_PI/2 : -SGD_PI/2);
332 angle_rad = atan(wind_from_north_fps/wind_from_east_fps);
334 wind_from_heading_deg = angle_rad * SGD_RADIANS_TO_DEGREES;
335 if (wind_from_east_fps >= 0)
336 wind_from_heading_deg = 90 - wind_from_heading_deg;
338 wind_from_heading_deg = 270 - wind_from_heading_deg;
340 wind_speed_kt = fabs(wind_from_east_fps
341 * SG_METER_TO_NM * SG_FEET_TO_METER * 3600);
343 wind_speed_kt = (wind_from_north_fps / sin(angle_rad))
344 * SG_METER_TO_NM * SG_FEET_TO_METER * 3600;
348 FGEnvironment::_recalc_ne ()
351 wind_speed_kt * SG_NM_TO_METER * SG_METER_TO_FEET * (1.0/3600);
353 wind_from_north_fps = speed_fps *
354 cos(wind_from_heading_deg * SGD_DEGREES_TO_RADIANS);
355 wind_from_east_fps = speed_fps *
356 sin(wind_from_heading_deg * SGD_DEGREES_TO_RADIANS);
360 FGEnvironment::_recalc_sl_temperature ()
362 // If we're in the stratosphere, leave sea-level temp alone
363 if (elevation_ft < 38000) {
364 temperature_sea_level_degc =
365 (temperature_degc + 273.15)
366 /_temperature_degc_table->interpolate(elevation_ft)
372 FGEnvironment::_recalc_alt_temperature ()
374 if (elevation_ft < 38000) {
376 (temperature_sea_level_degc + 273.15) *
377 _temperature_degc_table->interpolate(elevation_ft) - 273.15;
379 temperature_degc = -56.49; // Stratosphere is constant
384 FGEnvironment::_recalc_sl_dewpoint ()
387 // FIXME: this will work only for low
389 dewpoint_sea_level_degc = dewpoint_degc + (elevation_ft * .0002);
390 if (dewpoint_sea_level_degc > temperature_sea_level_degc)
391 dewpoint_sea_level_degc = temperature_sea_level_degc;
395 FGEnvironment::_recalc_alt_dewpoint ()
398 // FIXME: this will work only for low
400 dewpoint_degc = dewpoint_sea_level_degc + (elevation_ft * .0002);
401 if (dewpoint_degc > temperature_degc)
402 dewpoint_degc = temperature_degc;
406 FGEnvironment::_recalc_sl_pressure ()
408 pressure_sea_level_inhg =
409 pressure_inhg / _pressure_inhg_table->interpolate(elevation_ft);
413 FGEnvironment::_recalc_alt_pressure ()
416 pressure_sea_level_inhg * _pressure_inhg_table->interpolate(elevation_ft);
420 FGEnvironment::_recalc_density ()
422 double pressure_psf = pressure_inhg * 70.7487;
424 // adjust for humidity
425 // calculations taken from USA Today (oops!) at
426 // http://www.usatoday.com/weather/basics/density-calculations.htm
427 double temperature_degk = temperature_degc + 273.15;
428 double pressure_mb = pressure_inhg * 33.86;
429 double vapor_pressure_mb =
430 6.11 * pow(10.0, 7.5 * dewpoint_degc / (237.7 + dewpoint_degc));
431 double virtual_temperature_degk = temperature_degk / (1 - (vapor_pressure_mb / pressure_mb) * (1.0 - 0.622));
432 double virtual_temperature_degr = virtual_temperature_degk * 1.8;
434 density_slugft3 = pressure_psf / (virtual_temperature_degr * 1718);
437 // end of environment.cxx