src/Environment/environment_ctrl.cxx

   1 // environment_ctrl.cxx -- manager for natural environment information.
   2 //
   3 // Written by David Megginson, started February 2002.
   4 //
   5 // Copyright (C) 2002  David Megginson - david@megginson.com
   6 //
   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.
  11 //
  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.
  16 //
  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., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.
  20 //
  21 // $Id$
  22
  23 #ifdef HAVE_CONFIG_H
  24 #  include "config.h"
  25 #endif
  26
  27 #include <algorithm>
  28
  29 #include <simgear/debug/logstream.hxx>
  30 #include <simgear/structure/commands.hxx>
  31 #include <simgear/structure/exception.hxx>
  32
  33 #include <Airports/simple.hxx>
  34 #include <Main/fg_props.hxx>
  35 #include <Main/util.hxx>
  36
  37 #include "atmosphere.hxx"
  38 #include "fgmetar.hxx"
  39 #include "environment_ctrl.hxx"
  40
  41 using std::sort;
  42
  43 class AirportWithMetar : public FGAirport::AirportFilter {
  44 public:
  45         virtual bool passAirport(FGAirport* aApt) const {
  46                 return aApt->getMetar();
  47         }
  48
  49   // permit heliports and seaports too
  50   virtual FGPositioned::Type maxType() const
  51   { return FGPositioned::SEAPORT; }
  52 };
  53
  54 static AirportWithMetar airportWithMetarFilter;
  55 \f
  56 ////////////////////////////////////////////////////////////////////////
  57 // Implementation of FGEnvironmentCtrl abstract base class.
  58 ////////////////////////////////////////////////////////////////////////
  59
  60 FGEnvironmentCtrl::FGEnvironmentCtrl ()
  61         : _environment(0),
  62         _lon_deg(0),
  63         _lat_deg(0),
  64         _elev_ft(0)
  65 {
  66 }
  67
  68 FGEnvironmentCtrl::~FGEnvironmentCtrl ()
  69 {
  70 }
  71
  72 void
  73 FGEnvironmentCtrl::setEnvironment (FGEnvironment * environment)
  74 {
  75         _environment = environment;
  76 }
  77
  78 void
  79 FGEnvironmentCtrl::setLongitudeDeg (double lon_deg)
  80 {
  81         _lon_deg = lon_deg;
  82 }
  83
  84 void
  85 FGEnvironmentCtrl::setLatitudeDeg (double lat_deg)
  86 {
  87         _lat_deg = lat_deg;
  88 }
  89
  90 void
  91 FGEnvironmentCtrl::setElevationFt (double elev_ft)
  92 {
  93         _elev_ft = elev_ft;
  94 }
  95
  96 void
  97 FGEnvironmentCtrl::setPosition (double lon_deg, double lat_deg, double elev_ft)
  98 {
  99         _lon_deg = lon_deg;
 100         _lat_deg = lat_deg;
 101         _elev_ft = elev_ft;
 102 }
 103
 104
 105 \f
 106 ////////////////////////////////////////////////////////////////////////
 107 // Implementation of FGInterpolateEnvironmentCtrl.
 108 ////////////////////////////////////////////////////////////////////////
 109
 110
 111 FGInterpolateEnvironmentCtrl::FGInterpolateEnvironmentCtrl ()
 112 {
 113         altitude_n = fgGetNode("/position/altitude-ft", true);
 114         altitude_agl_n = fgGetNode("/position/altitude-agl-ft", true);
 115         boundary_transition_n = fgGetNode("/environment/config/boundary-transition-ft", false );
 116         boundary_n = fgGetNode("/environment/config/boundary", true );
 117         aloft_n = fgGetNode("/environment/config/aloft", true );
 118 }
 119
 120 FGInterpolateEnvironmentCtrl::~FGInterpolateEnvironmentCtrl ()
 121 {
 122         unsigned int i;
 123         for (i = 0; i < _boundary_table.size(); i++)
 124                 delete _boundary_table[i];
 125         for (i = 0; i < _aloft_table.size(); i++)
 126                 delete _aloft_table[i];
 127 }
 128
 129
 130
 131 void
 132 FGInterpolateEnvironmentCtrl::init ()
 133 {
 134         read_table( boundary_n, _boundary_table);
 135         read_table( aloft_n, _aloft_table);
 136 }
 137
 138 void
 139 FGInterpolateEnvironmentCtrl::reinit ()
 140 {
 141         init();
 142 }
 143
 144 void
 145 FGInterpolateEnvironmentCtrl::read_table (const SGPropertyNode * node, vector<bucket *> &table)
 146 {
 147         double last_altitude_ft = 0.0;
 148         double sort_required = false;
 149         size_t i;
 150
 151         for (i = 0; i < (size_t)node->nChildren(); i++) {
 152                 const SGPropertyNode * child = node->getChild(i);
 153                 if ( strcmp(child->getName(), "entry") == 0
 154                  && child->getStringValue("elevation-ft", "")[0] != '\0'
 155                  && ( child->getDoubleValue("elevation-ft") > 0.1 || i == 0 ) )
 156         {
 157                         bucket * b;
 158                         if( i < table.size() ) {
 159                                 // recycle existing bucket
 160                                 b = table[i];
 161                         } else {
 162                                 // more nodes than buckets in table, add a new one
 163                                 b = new bucket;
 164                                 table.push_back(b);
 165                         }
 166                         if (i > 0)
 167                                 b->environment.copy(table[i-1]->environment);
 168                         b->environment.read(child);
 169                         b->altitude_ft = b->environment.get_elevation_ft();
 170
 171                         // check, if altitudes are in ascending order
 172                         if( b->altitude_ft < last_altitude_ft )
 173                                 sort_required = true;
 174                         last_altitude_ft = b->altitude_ft;
 175                 }
 176         }
 177         // remove leftover buckets
 178         while( table.size() > i ) {
 179                 bucket * b = *(table.end() - 1);
 180                 delete b;
 181                 table.pop_back();
 182         }
 183
 184         if( sort_required )
 185                 sort(table.begin(), table.end(), bucket::lessThan);
 186 }
 187
 188 void
 189 FGInterpolateEnvironmentCtrl::update (double delta_time_sec)
 190 {
 191         double altitude_ft = altitude_n->getDoubleValue();
 192         double altitude_agl_ft = altitude_agl_n->getDoubleValue();
 193         double boundary_transition =
 194                 boundary_transition_n == NULL ? 500 : boundary_transition_n->getDoubleValue();
 195
 196         int length = _boundary_table.size();
 197
 198         if (length > 0) {
 199                                                                 // boundary table
 200                 double boundary_limit = _boundary_table[length-1]->altitude_ft;
 201                 if (boundary_limit >= altitude_agl_ft) {
 202                         do_interpolate(_boundary_table, altitude_agl_ft, _environment);
 203                         return;
 204                 } else if ((boundary_limit + boundary_transition) >= altitude_agl_ft) {
 205                         //TODO: this is 500ft above the top altitude of boundary layer
 206                         //shouldn't this be +/-250 ft off of the top altitude?
 207                                                                 // both tables
 208                         do_interpolate(_boundary_table, altitude_agl_ft, &env1);
 209                         do_interpolate(_aloft_table, altitude_ft, &env2);
 210                         double fraction =
 211                                 (altitude_agl_ft - boundary_limit) / boundary_transition;
 212                         interpolate(&env1, &env2, fraction, _environment);
 213                         return;
 214                 }
 215         }
 216                                                                 // aloft table
 217         do_interpolate(_aloft_table, altitude_ft, _environment);
 218 }
 219
 220 void
 221 FGInterpolateEnvironmentCtrl::do_interpolate (vector<bucket *> &table, double altitude_ft, FGEnvironment * environment)
 222 {
 223         int length = table.size();
 224         if (length == 0)
 225                 return;
 226
 227                                                                 // Boundary conditions
 228         if ((length == 1) || (table[0]->altitude_ft >= altitude_ft)) {
 229                 environment->copy(table[0]->environment);
 230                 return;
 231         } else if (table[length-1]->altitude_ft <= altitude_ft) {
 232                 environment->copy(table[length-1]->environment);
 233                 return;
 234         }
 235                                                                 // Search the interpolation table
 236         for (int i = 0; i < length - 1; i++) {
 237                 if ((i == length - 1) || (table[i]->altitude_ft <= altitude_ft)) {
 238                                 FGEnvironment * env1 = &(table[i]->environment);
 239                                 FGEnvironment * env2 = &(table[i+1]->environment);
 240                                 double fraction;
 241                                 if (table[i]->altitude_ft == table[i+1]->altitude_ft)
 242                                         fraction = 1.0;
 243                                 else
 244                                         fraction =
 245                                                 ((altitude_ft - table[i]->altitude_ft) /
 246                                                  (table[i+1]->altitude_ft - table[i]->altitude_ft));
 247                                 interpolate(env1, env2, fraction, environment);
 248
 249                                 return;
 250                 }
 251         }
 252 }
 253
 254 bool
 255 FGInterpolateEnvironmentCtrl::bucket::operator< (const bucket &b) const
 256 {
 257         return (altitude_ft < b.altitude_ft);
 258 }
 259
 260 bool
 261 FGInterpolateEnvironmentCtrl::bucket::lessThan(bucket *a, bucket *b)
 262 {
 263         return (a->altitude_ft) < (b->altitude_ft);
 264 }
 265
 266 \f
 267 ////////////////////////////////////////////////////////////////////////
 268 // Implementation of FGMetarCtrl.
 269 ////////////////////////////////////////////////////////////////////////
 270
 271 FGMetarCtrl::FGMetarCtrl( SGSubsystem * environmentCtrl )
 272         :
 273         metar_valid(false),
 274         setup_winds_aloft(true),
 275         wind_interpolation_required(true),
 276         // Interpolation constant definitions.
 277         EnvironmentUpdatePeriodSec( 0.2 ),
 278         MaxWindChangeKtsSec( 0.2 ),
 279         MaxVisChangePercentSec( 0.05 ),
 280         MaxPressureChangeInHgSec( 0.0033 ),
 281         MaxCloudAltitudeChangeFtSec( 20.0 ),
 282         MaxCloudThicknessChangeFtSec( 50.0 ),
 283         MaxCloudInterpolationHeightFt( 5000.0 ),
 284         MaxCloudInterpolationDeltaFt( 4000.0 ),
 285         _environmentCtrl(environmentCtrl)
 286 {
 287         windModulator = new FGBasicWindModulator();
 288
 289         metar_base_n = fgGetNode( "/environment/metar", true );
 290         station_id_n = metar_base_n->getNode("station-id", true );
 291         station_elevation_n = metar_base_n->getNode("station-elevation-ft", true );
 292         min_visibility_n = metar_base_n->getNode("min-visibility-m", true );
 293         max_visibility_n = metar_base_n->getNode("max-visibility-m", true );
 294         base_wind_range_from_n = metar_base_n->getNode("base-wind-range-from", true );
 295         base_wind_range_to_n = metar_base_n->getNode("base-wind-range-to", true );
 296         base_wind_speed_n = metar_base_n->getNode("base-wind-speed-kt", true );
 297         base_wind_dir_n = metar_base_n->getNode("base-wind-dir-deg", true );
 298         gust_wind_speed_n = metar_base_n->getNode("gust-wind-speed-kt", true );
 299         temperature_n = metar_base_n->getNode("temperature-degc", true );
 300         dewpoint_n = metar_base_n->getNode("dewpoint-degc", true );
 301         humidity_n = metar_base_n->getNode("rel-humidity-norm", true );
 302         pressure_n = metar_base_n->getNode("pressure-inhg", true );
 303         clouds_n = metar_base_n->getNode("clouds", true );
 304         rain_n = metar_base_n->getNode("rain-norm", true );
 305         hail_n = metar_base_n->getNode("hail-norm", true );
 306         snow_n = metar_base_n->getNode("snow-norm", true );
 307         snow_cover_n = metar_base_n->getNode("snow-cover", true );
 308         magnetic_variation_n = fgGetNode( "/environment/magnetic-variation-deg", true );
 309         ground_elevation_n = fgGetNode( "/position/ground-elev-m", true );
 310         longitude_n = fgGetNode( "/position/longitude-deg", true );
 311         latitude_n = fgGetNode( "/position/latitude-deg", true );
 312         environment_clouds_n = fgGetNode("/environment/clouds");
 313
 314         boundary_wind_speed_n = fgGetNode("/environment/config/boundary/entry/wind-speed-kt");
 315         boundary_wind_from_heading_n = fgGetNode("/environment/config/boundary/entry/wind-from-heading-deg");
 316         boundary_visibility_n = fgGetNode("/environment/config/boundary/entry/visibility-m");
 317         boundary_sea_level_pressure_n = fgGetNode("/environment/config/boundary/entry/pressure-sea-level-inhg");
 318 }
 319
 320 FGMetarCtrl::~FGMetarCtrl ()
 321 {
 322 }
 323
 324 void FGMetarCtrl::bind ()
 325 {
 326         fgTie("/environment/metar/valid", this, &FGMetarCtrl::get_valid );
 327         fgTie("/environment/params/metar-updates-environment", this, &FGMetarCtrl::get_enabled, &FGMetarCtrl::set_enabled );
 328         fgTie("/environment/params/metar-updates-winds-aloft", this, &FGMetarCtrl::get_setup_winds_aloft, &FGMetarCtrl::set_setup_winds_aloft );
 329 }
 330
 331 void FGMetarCtrl::unbind ()
 332 {
 333         fgUntie("/environment/metar/valid");
 334         fgUntie("/environment/params/metar-updates-environment");
 335         fgUntie("/environment/params/metar-updates-winds-aloft");
 336 }
 337
 338 // use a "command" to set station temp at station elevation
 339 static void set_temp_at_altitude( float temp_degc, float altitude_ft ) {
 340         SGPropertyNode args;
 341         SGPropertyNode *node = args.getNode("temp-degc", 0, true);
 342         node->setFloatValue( temp_degc );
 343         node = args.getNode("altitude-ft", 0, true);
 344         node->setFloatValue( altitude_ft );
 345         globals->get_commands()->execute("set-outside-air-temp-degc", &args);
 346 }
 347
 348 static void set_dewpoint_at_altitude( float dewpoint_degc, float altitude_ft ) {
 349         SGPropertyNode args;
 350         SGPropertyNode *node = args.getNode("dewpoint-degc", 0, true);
 351         node->setFloatValue( dewpoint_degc );
 352         node = args.getNode("altitude-ft", 0, true);
 353         node->setFloatValue( altitude_ft );
 354         globals->get_commands()->execute("set-dewpoint-temp-degc", &args);
 355 }
 356
 357 /*
 358  Setup the wind nodes for a branch in the /environment/config/<branchName>/entry nodes
 359
 360  Output properties:
 361  wind-from-heading-deg
 362  wind-speed-kt
 363  turbulence/magnitude-norm
 364
 365  Input properties:
 366  wind-heading-change-deg       how many degrees does the wind direction change at this level
 367  wind-speed-change-rel         relative change of wind speed at this level
 368  turbulence/factor             factor for the calculated turbulence magnitude at this level
 369  */
 370 static void setupWindBranch( string branchName, double dir, double speed, double gust )
 371 {
 372         SGPropertyNode_ptr branch = fgGetNode("/environment/config", true)->getNode(branchName,true);
 373         vector<SGPropertyNode_ptr> entries = branch->getChildren("entry");
 374         for ( vector<SGPropertyNode_ptr>::iterator it = entries.begin(); it != entries.end(); it++) {
 375
 376                 // change wind direction as configured
 377                 double layer_dir = dir + (*it)->getDoubleValue("wind-heading-change-deg", 0.0 );
 378                 if( layer_dir >= 360.0 ) layer_dir -= 360.0;
 379                 if( layer_dir < 0.0 ) layer_dir += 360.0;
 380                 (*it)->setDoubleValue("wind-from-heading-deg", layer_dir);
 381
 382                 double layer_speed = speed*(1 + (*it)->getDoubleValue("wind-speed-change-rel", 0.0 ));
 383                 (*it)->setDoubleValue("wind-speed-kt", layer_speed );
 384
 385                 // add some turbulence
 386                 SGPropertyNode_ptr turbulence = (*it)->getNode("turbulence",true);
 387
 388                 double turbulence_norm = speed/50;
 389                 if( gust > speed ) {
 390                         turbulence_norm += (gust-speed)/25;
 391                 }
 392                 if( turbulence_norm > 1.0 ) turbulence_norm = 1.0;
 393
 394                 turbulence_norm *= turbulence->getDoubleValue("factor", 0.0 );
 395                 turbulence->setDoubleValue( "magnitude-norm", turbulence_norm );
 396         }
 397 }
 398
 399 static void setupWind( bool setup_aloft, double dir, double speed, double gust )
 400 {
 401         setupWindBranch( "boundary", dir, speed, gust );
 402         if( setup_aloft )
 403                 setupWindBranch( "aloft", dir, speed, gust );
 404 }
 405
 406 double FGMetarCtrl::interpolate_val(double currentval, double requiredval, double dt)
 407 {
 408         double dval = EnvironmentUpdatePeriodSec * dt;
 409
 410         if (fabs(currentval - requiredval) < dval) return requiredval;
 411         if (currentval < requiredval) return (currentval + dval);
 412         if (currentval > requiredval) return (currentval - dval);
 413         return requiredval;
 414 }
 415
 416 void
 417 FGMetarCtrl::init ()
 418 {
 419         first_update = true;
 420         wind_interpolation_required = true;
 421 }
 422
 423 void
 424 FGMetarCtrl::reinit ()
 425 {
 426         init();
 427 }
 428
 429 static inline double convert_to_360( double d )
 430 {
 431         if( d < 0.0 ) return d + 360.0;
 432         if( d >= 360.0 ) return d - 360.0;
 433         return d;
 434 }
 435
 436 static inline double convert_to_180( double d )
 437 {
 438         return d > 180.0 ? d - 360.0 : d;
 439 }
 440
 441 // Return the sea level pressure for a metar observation, in inHg.
 442 // This is different from QNH because it accounts for the current
 443 // temperature at the observation point.
 444 // metarPressure in inHg
 445 // fieldHt in ft
 446 // fieldTemp in C
 447
 448 static double reducePressureSl(double metarPressure, double fieldHt,
 449                                double fieldTemp)
 450 {
 451     double elev = fieldHt * SG_FEET_TO_METER;
 452     double fieldPressure
 453         = FGAtmo::fieldPressure(elev, metarPressure * atmodel::inHg);
 454     double slPressure = P_layer(0, elev, fieldPressure,
 455                                 fieldTemp + atmodel::freezing,
 456                                 atmodel::ISA::lam0);
 457     return slPressure / atmodel::inHg;
 458 }
 459
 460 void
 461 FGMetarCtrl::update(double dt)
 462 {
 463         if( dt <= 0 || !metar_valid ||!enabled)
 464                 return;
 465
 466         windModulator->update(dt);
 467         // Interpolate the current configuration closer to the actual METAR
 468
 469         bool reinit_required = false;
 470         bool layer_rebuild_required = false;
 471         double station_elevation_ft = station_elevation_n->getDoubleValue();
 472
 473         if (first_update) {
 474                 double dir = base_wind_dir_n->getDoubleValue()+magnetic_variation_n->getDoubleValue();
 475                 double speed = base_wind_speed_n->getDoubleValue();
 476                 double gust = gust_wind_speed_n->getDoubleValue();
 477                 setupWind(setup_winds_aloft, dir, speed, gust);
 478
 479                 double metarvis = min_visibility_n->getDoubleValue();
 480                 fgDefaultWeatherValue("visibility-m", metarvis);
 481
 482                 double metarpressure = pressure_n->getDoubleValue();
 483                 fgDefaultWeatherValue("pressure-sea-level-inhg",
 484                                       reducePressureSl(metarpressure,
 485                                                        station_elevation_ft,
 486                                                        temperature_n->getDoubleValue()));
 487
 488                 // We haven't already loaded a METAR, so apply it immediately.
 489                 vector<SGPropertyNode_ptr> layers = clouds_n->getChildren("layer");
 490                 vector<SGPropertyNode_ptr>::const_iterator layer;
 491                 vector<SGPropertyNode_ptr>::const_iterator layers_end = layers.end();
 492
 493                 int i;
 494                 for (i = 0, layer = layers.begin(); layer != layers_end; ++layer, i++) {
 495                         SGPropertyNode *target = environment_clouds_n->getChild("layer", i, true);
 496
 497                         target->setStringValue("coverage",
 498                                         (*layer)->getStringValue("coverage", "clear"));
 499                         target->setDoubleValue("elevation-ft",
 500                                         (*layer)->getDoubleValue("elevation-ft"));
 501                         target->setDoubleValue("thickness-ft",
 502                                         (*layer)->getDoubleValue("thickness-ft"));
 503                         target->setDoubleValue("span-m", 40000.0);
 504                 }
 505
 506                 first_update = false;
 507                 reinit_required = true;
 508                 layer_rebuild_required = true;
 509
 510         } else {
 511                 if( wind_interpolation_required ) {
 512                         // Generate interpolated values between the METAR and the current
 513                         // configuration.
 514
 515                         // Pick up the METAR wind values and convert them into a vector.
 516                         double metar[2];
 517                         double metar_speed = base_wind_speed_n->getDoubleValue();
 518                         double metar_heading = base_wind_dir_n->getDoubleValue()+magnetic_variation_n->getDoubleValue();
 519
 520                         metar[0] = metar_speed * sin(metar_heading * SG_DEGREES_TO_RADIANS );
 521                         metar[1] = metar_speed * cos(metar_heading * SG_DEGREES_TO_RADIANS);
 522
 523                         // Convert the current wind values and convert them into a vector
 524                         double current[2];
 525                         double speed = boundary_wind_speed_n->getDoubleValue();
 526                         double dir_from = boundary_wind_from_heading_n->getDoubleValue();;
 527
 528                         current[0] = speed * sin(dir_from * SG_DEGREES_TO_RADIANS );
 529                         current[1] = speed * cos(dir_from * SG_DEGREES_TO_RADIANS );
 530
 531                         // Determine the maximum component-wise value that the wind can change.
 532                         // First we determine the fraction in the X and Y component, then
 533                         // factor by the maximum wind change.
 534                         double x = fabs(current[0] - metar[0]);
 535                         double y = fabs(current[1] - metar[1]);
 536
 537                         // only interpolate if we have a difference
 538                         if (x + y > 0.01 ) {
 539                                 double dx = x / (x + y);
 540                                 double dy = 1 - dx;
 541
 542                                 double maxdx = dx * MaxWindChangeKtsSec;
 543                                 double maxdy = dy * MaxWindChangeKtsSec;
 544
 545                                 // Interpolate each component separately.
 546                                 current[0] = interpolate_val(current[0], metar[0], maxdx);
 547                                 current[1] = interpolate_val(current[1], metar[1], maxdy);
 548
 549                                 // Now convert back to polar coordinates.
 550                                 if ((fabs(current[0]) > 0.1) || (fabs(current[1]) > 0.1)) {
 551                                         // Some real wind to convert back from. Work out the speed
 552                                         // and direction value in degrees.
 553                                         speed = sqrt((current[0] * current[0]) + (current[1] * current[1]));
 554                                         dir_from = (atan2(current[0], current[1]) * SG_RADIANS_TO_DEGREES );
 555
 556                                         // Normalize the direction.
 557                                         if (dir_from < 0.0)
 558                                                 dir_from += 360.0;
 559
 560                                         SG_LOG( SG_GENERAL, SG_DEBUG, "Wind : " << dir_from << "@" << speed);
 561                                 } else {
 562                                         // Special case where there is no wind (otherwise atan2 barfs)
 563                                         speed = 0.0;
 564                                 }
 565                                 double gust = gust_wind_speed_n->getDoubleValue();
 566                                 setupWind(setup_winds_aloft, dir_from, speed, gust);
 567                                 reinit_required = true;
 568                         } else {
 569                                 wind_interpolation_required = false;
 570                         }
 571                 } else { // if(wind_interpolation_required)
 572                         // interpolation of wind vector is finished, apply wind
 573                         // variations and gusts for the boundary layer only
 574
 575
 576                         bool wind_modulated = false;
 577
 578                         // start with the main wind direction
 579                         double wind_dir = base_wind_dir_n->getDoubleValue()+magnetic_variation_n->getDoubleValue();
 580                         double min = convert_to_180(base_wind_range_from_n->getDoubleValue()+magnetic_variation_n->getDoubleValue());
 581                         double max = convert_to_180(base_wind_range_to_n->getDoubleValue()+magnetic_variation_n->getDoubleValue());
 582                         if( max > min ) {
 583                                 // if variable winds configured, modulate the wind direction
 584                                 double f = windModulator->get_direction_offset_norm();
 585                                 wind_dir = min+(max-min)*f;
 586                                 double old = convert_to_180(boundary_wind_from_heading_n->getDoubleValue());
 587                                 wind_dir = convert_to_360(fgGetLowPass(old, wind_dir, dt ));
 588                                 wind_modulated = true;
 589                         }
 590
 591                         // start with main wind speed
 592                         double wind_speed = base_wind_speed_n->getDoubleValue();
 593                         max = gust_wind_speed_n->getDoubleValue();
 594                         if( max > wind_speed ) {
 595                                 // if gusts are configured, modulate wind magnitude
 596                                 double f = windModulator->get_magnitude_factor_norm();
 597                                 wind_speed = wind_speed+(max-wind_speed)*f;
 598                                 wind_speed = fgGetLowPass(boundary_wind_speed_n->getDoubleValue(), wind_speed, dt );
 599                                 wind_modulated = true;
 600                         }
 601                         if( wind_modulated ) {
 602                                 setupWind(false, wind_dir, wind_speed, max);
 603                                 reinit_required = true;
 604                         }
 605                 }
 606
 607                 // Now handle the visibility. We convert both visibility values
 608                 // to X-values, then interpolate from there, then back to real values.
 609                 // The length_scale is fixed to 1000m, so the visibility changes by
 610                 // by MaxVisChangePercentSec or 1000m X MaxVisChangePercentSec,
 611                 // whichever is more.
 612                 double vis = boundary_visibility_n->getDoubleValue();;
 613                 double metarvis = min_visibility_n->getDoubleValue();
 614                 if( vis != metarvis ) {
 615                         double currentxval = log(1000.0 + vis);
 616                         double metarxval = log(1000.0 + metarvis);
 617
 618                         currentxval = interpolate_val(currentxval, metarxval, MaxVisChangePercentSec);
 619
 620                         // Now convert back from an X-value to a straightforward visibility.
 621                         vis = exp(currentxval) - 1000.0;
 622                         fgDefaultWeatherValue("visibility-m", vis);
 623                         reinit_required = true;
 624                 }
 625
 626                 double pressure = boundary_sea_level_pressure_n->getDoubleValue();
 627                 double metarpressure = pressure_n->getDoubleValue();
 628                 double newpressure = reducePressureSl(metarpressure,
 629                                                       station_elevation_ft,
 630                                                       temperature_n->getDoubleValue());
 631                 if( pressure != newpressure ) {
 632                         pressure = interpolate_val( pressure, newpressure, MaxPressureChangeInHgSec );
 633                         fgDefaultWeatherValue("pressure-sea-level-inhg", pressure);
 634                         reinit_required = true;
 635                 }
 636
 637                 // Set the cloud layers by interpolating over the METAR versions.
 638                 vector<SGPropertyNode_ptr> layers = clouds_n->getChildren("layer");
 639                 vector<SGPropertyNode_ptr>::const_iterator layer;
 640                 vector<SGPropertyNode_ptr>::const_iterator layers_end = layers.end();
 641
 642                 double aircraft_alt = fgGetDouble("/position/altitude-ft");
 643                 int i;
 644
 645                 for (i = 0, layer = layers.begin(); layer != layers_end; ++layer, i++) {
 646                         SGPropertyNode *target = environment_clouds_n->getChild("layer", i, true);
 647
 648                         // In the case of clouds, we want to avoid writing if nothing has
 649                         // changed, as these properties are tied to the renderer and will
 650                         // cause the clouds to be updated, reseting the texture locations.
 651
 652                         // We don't interpolate the coverage values as no-matter how we
 653                         // do it, it will be quite a sudden change of texture. Better to
 654                         // have a single change than four or five.
 655                         const char *coverage = (*layer)->getStringValue("coverage", "clear");
 656                         SGPropertyNode *cov = target->getNode("coverage", true);
 657                         if (strcmp(cov->getStringValue(), coverage) != 0) {
 658                                 cov->setStringValue(coverage);
 659                                 layer_rebuild_required = true;
 660                         }
 661
 662                         double required_alt = (*layer)->getDoubleValue("elevation-ft");
 663                         double current_alt = target->getDoubleValue("elevation-ft");
 664                         double required_thickness = (*layer)->getDoubleValue("thickness-ft");
 665                         SGPropertyNode *thickness = target->getNode("thickness-ft", true);
 666
 667                         if (current_alt < -9000 || required_alt < -9000 ||
 668                                 fabs(aircraft_alt - required_alt) > MaxCloudInterpolationHeightFt ||
 669                                 fabs(current_alt - required_alt) > MaxCloudInterpolationDeltaFt) {
 670                                 // We don't interpolate any layers that are
 671                                 //  - too far above us to be visible
 672                                 //  - too far below us to be visible
 673                                 //  - with too large a difference to make interpolation sensible
 674                                 //  - to or from -9999 (used as a placeholder)
 675                                 //  - any values that are too high above us,
 676                                 if (current_alt != required_alt)
 677                                         target->setDoubleValue("elevation-ft", required_alt);
 678
 679                                 if (thickness->getDoubleValue() != required_thickness)
 680                                         thickness->setDoubleValue(required_thickness);
 681
 682                         } else {
 683                                 // Interpolate the other values in the usual way
 684                                 if (current_alt != required_alt) {
 685                                         current_alt = interpolate_val(current_alt, required_alt, MaxCloudAltitudeChangeFtSec);
 686                                         target->setDoubleValue("elevation-ft", current_alt);
 687                                 }
 688
 689                                 double current_thickness = thickness->getDoubleValue();
 690
 691                                 if (current_thickness != required_thickness) {
 692                                         current_thickness = interpolate_val(current_thickness,
 693                                                                                                  required_thickness,
 694                                                                                                  MaxCloudThicknessChangeFtSec);
 695                                         thickness->setDoubleValue(current_thickness);
 696                                 }
 697                         }
 698                 }
 699         }
 700         set_temp_at_altitude(temperature_n->getDoubleValue(), station_elevation_ft);
 701         set_dewpoint_at_altitude(dewpoint_n->getDoubleValue(), station_elevation_ft);
 702         //TODO: check if temperature/dewpoint have changed. This requires reinit.
 703
 704         // Force an update of the 3D clouds
 705         if( layer_rebuild_required )
 706                 fgSetInt("/environment/rebuild-layers", 1 );
 707
 708         // Reinitializing of the environment controller required
 709         if( reinit_required )
 710                 _environmentCtrl->reinit();
 711 }
 712
 713 const char * FGMetarCtrl::get_metar(void) const
 714 {
 715         return metar.c_str();
 716 }
 717
 718 static const char *coverage_string[] = { "clear", "few", "scattered", "broken", "overcast" };
 719 static const double thickness_value[] = { 0, 65, 600, 750, 1000 };
 720
 721 void FGMetarCtrl::set_metar( const char * metar_string )
 722 {
 723         int i;
 724
 725         metar = metar_string;
 726
 727         SGSharedPtr<FGMetar> m;
 728         try {
 729                 m = new FGMetar( metar_string );
 730         }
 731         catch( sg_io_exception ) {
 732                 fprintf( stderr, "can't get metar: %s\n", metar_string );
 733                 metar_valid = false;
 734                 return;
 735         }
 736
 737         wind_interpolation_required = true;
 738
 739         min_visibility_n->setDoubleValue( m->getMinVisibility().getVisibility_m() );
 740         max_visibility_n->setDoubleValue( m->getMaxVisibility().getVisibility_m() );
 741
 742         const SGMetarVisibility *dirvis = m->getDirVisibility();
 743         for (i = 0; i < 8; i++, dirvis++) {
 744                 SGPropertyNode *vis = metar_base_n->getChild("visibility", i, true);
 745                 double v = dirvis->getVisibility_m();
 746
 747                 vis->setDoubleValue("min-m", v);
 748                 vis->setDoubleValue("max-m", v);
 749         }
 750
 751         base_wind_dir_n->setIntValue( m->getWindDir() );
 752         base_wind_range_from_n->setIntValue( m->getWindRangeFrom() );
 753         base_wind_range_to_n->setIntValue( m->getWindRangeTo() );
 754         base_wind_speed_n->setDoubleValue( m->getWindSpeed_kt() );
 755         gust_wind_speed_n->setDoubleValue( m->getGustSpeed_kt() );
 756         temperature_n->setDoubleValue( m->getTemperature_C() );
 757         dewpoint_n->setDoubleValue( m->getDewpoint_C() );
 758         humidity_n->setDoubleValue( m->getRelHumidity() );
 759         pressure_n->setDoubleValue( m->getPressure_inHg() );
 760
 761
 762         // get station elevation to compute cloud base
 763         double station_elevation_ft = 0;
 764         {
 765                 // 1. check the id given in the metar
 766                 FGAirport* a = FGAirport::findByIdent(m->getId());
 767
 768                 // 2. if unknown, find closest airport with metar to current position
 769                 if( a == NULL ) {
 770                         SGGeod pos = SGGeod::fromDeg(longitude_n->getDoubleValue(), latitude_n->getDoubleValue());
 771                         a = FGAirport::findClosest(pos, 10000.0, &airportWithMetarFilter);
 772                 }
 773
 774                 // 3. otherwise use ground elevation
 775                 if( a != NULL ) {
 776                         station_elevation_ft = a->getElevation();
 777                         station_id_n->setStringValue( a->ident());
 778                 } else {
 779                         station_elevation_ft = ground_elevation_n->getDoubleValue() * SG_METER_TO_FEET;
 780                         station_id_n->setStringValue( m->getId());
 781                 }
 782         }
 783
 784         station_elevation_n->setDoubleValue( station_elevation_ft );
 785
 786         vector<SGMetarCloud> cv = m->getClouds();
 787         vector<SGMetarCloud>::const_iterator cloud, cloud_end = cv.end();
 788
 789         int layer_cnt = environment_clouds_n->getChildren("layer").size();
 790         for (i = 0, cloud = cv.begin(); i < layer_cnt; i++) {
 791
 792
 793                 const char *coverage = "clear";
 794                 double elevation = -9999.0;
 795                 double thickness = 0.0;
 796                 const double span = 40000.0;
 797
 798                 if (cloud != cloud_end) {
 799                         int c = cloud->getCoverage();
 800                         coverage = coverage_string[c];
 801                         elevation = cloud->getAltitude_ft() + station_elevation_ft;
 802                         thickness = thickness_value[c];
 803                         ++cloud;
 804                 }
 805
 806                 SGPropertyNode *layer = clouds_n->getChild("layer", i, true );
 807
 808                 // if the coverage has changed, a rebuild of the layer is needed
 809                 if( strcmp(layer->getStringValue("coverage"), coverage ) ) {
 810                         layer->setStringValue("coverage", coverage);
 811                 }
 812                 layer->setDoubleValue("elevation-ft", elevation);
 813                 layer->setDoubleValue("thickness-ft", thickness);
 814                 layer->setDoubleValue("span-m", span);
 815         }
 816
 817         rain_n->setDoubleValue(m->getRain());
 818         hail_n->setDoubleValue(m->getHail());
 819         snow_n->setDoubleValue(m->getSnow());
 820         snow_cover_n->setBoolValue(m->getSnowCover());
 821         metar_valid = true;
 822 }
 823
 824 #if defined(ENABLE_THREADS)
 825 /**
 826  * This class represents the thread of execution responsible for
 827  * fetching the metar data.
 828  */
 829 class MetarThread : public OpenThreads::Thread {
 830 public:
 831         MetarThread( FGMetarFetcher * f ) : metar_fetcher(f) {}
 832         ~MetarThread() {}
 833
 834         /**
 835          * Fetche the metar data from the NOAA.
 836          */
 837         void run();
 838
 839 private:
 840         FGMetarFetcher * metar_fetcher;
 841 };
 842
 843 void MetarThread::run()
 844 {
 845         for( ;; ) {
 846                 string airport_id = metar_fetcher->request_queue.pop();
 847
 848                 if( airport_id.size() == 0 )
 849                         break;
 850
 851                 if( metar_fetcher->_error_count > 3 ) {
 852                         SG_LOG( SG_GENERAL, SG_WARN, "Too many erros fetching METAR, thread stopped permanently.");
 853                         break;
 854                 }
 855
 856                 metar_fetcher->fetch( airport_id );
 857         }
 858 }
 859 #endif
 860
 861 FGMetarFetcher::FGMetarFetcher() :
 862 #if defined(ENABLE_THREADS)
 863         metar_thread(NULL),
 864 #endif
 865         fetch_timer(0.0),
 866         search_timer(0.0),
 867         error_timer(0.0),
 868         _stale_count(0),
 869         _error_count(0),
 870         enabled(false)
 871 {
 872         longitude_n = fgGetNode( "/position/longitude-deg", true );
 873         latitude_n  = fgGetNode( "/position/latitude-deg", true );
 874         enable_n    = fgGetNode( "/environment/params/real-world-weather-fetch", true );
 875
 876         proxy_host_n = fgGetNode("/sim/presets/proxy/host", true);
 877         proxy_port_n = fgGetNode("/sim/presets/proxy/port", true);
 878         proxy_auth_n = fgGetNode("/sim/presets/proxy/authentication", true);
 879         max_age_n    = fgGetNode("/environment/params/metar-max-age-min", true);
 880
 881         output_n         = fgGetNode("/environment/metar/data", true );
 882 #if defined(ENABLE_THREADS)
 883         metar_thread = new MetarThread(this);
 884 // FIXME: do we really need setProcessorAffinity()?
 885 //      metar_thread->setProcessorAffinity(1);
 886         metar_thread->start();
 887 #endif // ENABLE_THREADS
 888 }
 889
 890
 891 FGMetarFetcher::~FGMetarFetcher()
 892 {
 893 #if defined(ENABLE_THREADS)
 894         request_queue.push("");
 895         metar_thread->join();
 896         delete metar_thread;
 897 #endif // ENABLE_THREADS
 898 }
 899
 900 void FGMetarFetcher::init ()
 901 {
 902         fetch_timer = 0.0;
 903         search_timer = 0.0;
 904         error_timer = 0.0;
 905         _stale_count = 0;
 906         _error_count = 0;
 907         current_airport_id.clear();
 908         /* Torsten Dreyer:
 909            hack to stop startup.nas complaining if metar arrives after nasal-dir-initialized
 910            is fired. Immediately fetch and wait for the METAR before continuing. This gets the
 911            /environment/metar/xxx properties filled before nasal-dir is initialized.
 912            Maybe the runway selection should happen here to make startup.nas obsolete?
 913         */
 914         const char * startup_airport = fgGetString("/sim/startup/options/airport");
 915         if( *startup_airport ) {
 916                 FGAirport * a = FGAirport::getByIdent( startup_airport );
 917                 if( a ) {
 918                         SGGeod pos = SGGeod::fromDeg(a->getLongitude(), a->getLatitude());
 919                         a = FGAirport::findClosest(pos, 10000.0, &airportWithMetarFilter);
 920                         current_airport_id = a->getId();
 921                         fetch( current_airport_id );
 922                 }
 923         }
 924 }
 925
 926 void FGMetarFetcher::reinit ()
 927 {
 928         init();
 929 }
 930
 931 /* search for closest airport with metar every xx seconds */
 932 static const int search_interval_sec = 60;
 933
 934 /* fetch metar for airport, even if airport has not changed every xx seconds */
 935 static const int fetch_interval_sec = 900;
 936
 937 /* reset error counter after xxx seconds */
 938 static const int error_timer_sec = 3;
 939
 940 void FGMetarFetcher::update (double delta_time_sec)
 941 {
 942         fetch_timer -= delta_time_sec;
 943         search_timer -= delta_time_sec;
 944         error_timer -= delta_time_sec;
 945
 946         if( error_timer <= 0.0 ) {
 947                 error_timer = error_timer_sec;
 948                 _error_count = 0;
 949         }
 950
 951         if( enable_n->getBoolValue() == false ) {
 952                 enabled = false;
 953                 return;
 954         }
 955
 956         // we were just enabled, reset all timers to
 957         // trigger immediate metar fetch
 958         if( !enabled ) {
 959                 search_timer = 0.0;
 960                 fetch_timer = 0.0;
 961                 error_timer = error_timer_sec;
 962                 enabled = true;
 963         }
 964
 965         FGAirport * a = NULL;
 966
 967         if( search_timer <= 0.0 ) {
 968                 // search timer expired, search closest airport with metar
 969                 SGGeod pos = SGGeod::fromDeg(longitude_n->getDoubleValue(), latitude_n->getDoubleValue());
 970                 a = FGAirport::findClosest(pos, 10000.0, &airportWithMetarFilter);
 971                 search_timer = search_interval_sec;
 972         }
 973
 974         if( a == NULL )
 975                 return;
 976
 977
 978         if( a->ident() != current_airport_id || fetch_timer <= 0 ) {
 979                 // fetch timer expired or airport has changed, schedule a fetch
 980                 current_airport_id = a->ident();
 981                 fetch_timer = fetch_interval_sec;
 982 #if defined(ENABLE_THREADS)
 983                 // push this airport id into the queue for the worker thread
 984                 request_queue.push( current_airport_id );
 985 #else
 986                 // if there is no worker thread, immediately fetch the data
 987                 fetch( current_airport_id );
 988 #endif
 989         }
 990 }
 991
 992 void FGMetarFetcher::fetch( const string & id )
 993 {
 994         if( enable_n->getBoolValue() == false )
 995                 return;
 996
 997         SGSharedPtr<FGMetar> result = NULL;
 998
 999         // fetch current metar data
1000         try {
1001                 string host = proxy_host_n->getStringValue();
1002                 string auth = proxy_auth_n->getStringValue();
1003                 string port = proxy_port_n->getStringValue();
1004
1005                 result = new FGMetar( id, host, port, auth);
1006
1007                 long max_age = max_age_n->getLongValue();
1008                 long age = result->getAge_min();
1009
1010                 if (max_age && age > max_age) {
1011                         SG_LOG( SG_GENERAL, SG_WARN, "METAR data too old (" << age << " min).");
1012                         if (++_stale_count > 10) {
1013                                 _error_count = 1000;
1014                                 throw sg_io_exception("More than 10 stale METAR messages in a row." " Check your system time!");
1015                         }
1016                 } else {
1017                         _stale_count = 0;
1018                 }
1019
1020         } catch (const sg_io_exception& e) {
1021                 SG_LOG( SG_GENERAL, SG_WARN, "Error fetching live weather data: " << e.getFormattedMessage().c_str() );
1022                 result = NULL;
1023                 // remove METAR flag from the airport
1024                 FGAirport * a = FGAirport::findByIdent( id );
1025                 if( a ) a->setMetar( false );
1026                 // immediately schedule a new search
1027                 search_timer = 0.0;
1028         }
1029
1030         // write the metar to the property node, the rest is done by the methods tied to this property
1031         // don't write the metar data, if real-weather-fetch has been disabled in the meantime
1032         if( result != NULL && enable_n->getBoolValue() == true )
1033                 output_n->setStringValue( result->getData() );
1034 }
1035
1036 // end of environment_ctrl.cxx
1037