src/Instrumentation/vertical_speed_indicator.cxx

   1 // vertical_speed_indicator.cxx - a regular VSI.
   2 // Written by David Megginson, started 2002.
   3 //
   4 // Last change by E. van den Berg, 17.02.1013
   5 //
   6 // This file is in the Public Domain and comes with no warranty.
   7
   8 #ifdef HAVE_CONFIG_H
   9 #  include "config.h"
  10 #endif
  11
  12 #include <simgear/constants.h>
  13 #include <simgear/math/interpolater.hxx>
  14
  15 #include "vertical_speed_indicator.hxx"
  16 #include <Main/fg_props.hxx>
  17 #include <Main/util.hxx>
  18
  19 //** NOTE: do not change these values. If you change one of them the others need to be changed too */
  20 //** these values calibrate the VSI at SL. */
  21 #define Vol_casing 1.25e-4          //m3
  22 #define A_orifice 7.853982e-9       //m2
  23 #define Factor_cal 189.145628       //-
  24
  25 using std::string;
  26
  27 VerticalSpeedIndicator::VerticalSpeedIndicator ( SGPropertyNode *node )
  28     : _casing_pressure_Pa(101325),
  29       _name(node->getStringValue("name", "vertical-speed-indicator")),
  30       _num(node->getIntValue("number", 0)),
  31       _static_pressure(node->getStringValue("static-pressure", "/systems/static/pressure-inhg")),
  32       _static_temperature(node->getStringValue("static-temperature", "/environment/temperature-degc"))
  33 {
  34 }
  35
  36 VerticalSpeedIndicator::~VerticalSpeedIndicator ()
  37 {
  38 }
  39
  40 void
  41 VerticalSpeedIndicator::init ()
  42 {
  43     string branch;
  44     branch = "/instrumentation/" + _name;
  45
  46     SGPropertyNode *node = fgGetNode(branch.c_str(), _num, true );
  47     _serviceable_node = node->getChild("serviceable", 0, true);
  48     _pressure_node = fgGetNode(_static_pressure.c_str(), true);
  49     _temperature_node = fgGetNode(_static_temperature.c_str(), true);
  50     _speed_fpm_node = node->getChild("indicated-speed-fpm", 0, true);
  51     _speed_mps_node = node->getChild("indicated-speed-mps", 0, true);
  52     _speed_kts_node = node->getChild("indicated-speed-kts", 0, true);
  53     _speed_up_node = fgGetNode("/sim/speed-up", true);
  54
  55     reinit();
  56 }
  57
  58 void
  59 VerticalSpeedIndicator::reinit ()
  60 {
  61     // Initialize at ambient conditions
  62     double casing_pressure_inHg = _pressure_node->getDoubleValue();
  63     _casing_pressure_Pa =  casing_pressure_inHg * SG_INHG_TO_PA;
  64     double casing_temperature_C = _temperature_node->getDoubleValue();
  65     double casing_temperature_K = casing_temperature_C + 273.15;
  66     _casing_density_kgpm3 = _casing_pressure_Pa / (casing_temperature_K * SG_R_m2_p_s2_p_K);
  67     _casing_airmass_kg = _casing_density_kgpm3 * Vol_casing;
  68     _orifice_massflow_kgps = 0.0;
  69 }
  70
  71 void
  72 VerticalSpeedIndicator::update (double dt)
  73 {
  74     if (_serviceable_node->getBoolValue()) {
  75         double pressure_inHg = _pressure_node->getDoubleValue() ;
  76         double pressure_Pa = pressure_inHg * SG_INHG_TO_PA;
  77         double speed_up = _speed_up_node->getDoubleValue();
  78         double Fsign = 0.;
  79         double orifice_mach = 0.0;
  80         if( speed_up > 1 )
  81             dt *= speed_up;
  82
  83 // This is a thermodynamically correct model of a mechanical vertical speed indicator:
  84 // It represents an aneroid in a closed (constant volume) casing with the aneroid internal pressure = static pressure
  85 // The casing has an orifice to static pressure
  86 // the mass flow through the orifice is calculated using compressible aerodynamics (but adiabatic and of course a perfect gas)
  87 // using the pressure in the casing and static pressure
  88 //
  89 // sadly at very low flows (small VS) in conjunction with the fact discrete timesteps (dt) are used, a numerical instability is formed.
  90 // this is counteracted by setting the massflow 0 at very small pressure differentials
  91 // this causes a small funny jump of your VSI when passing through 0...cannot be helped!
  92 //
  93 // also note the calibration is only valid for 0ft, so at higher altitudes, the vertical speed is not correct, but would indicate as a real mechanical VSI.
  94 // Only use for conventional mechanical VSI-s. Dont use in an Air Data Computer.
  95 //
  96 // (...and it is supposed to lag!)
  97
  98         _casing_airmass_kg = _casing_airmass_kg - _orifice_massflow_kgps * dt;
  99         double new_density_kgpm3 = _casing_airmass_kg / Vol_casing;
 100         _casing_pressure_Pa = _casing_pressure_Pa * pow(new_density_kgpm3 / _casing_density_kgpm3 , SG_gamma);
 101         double casing_temperature_K = _casing_pressure_Pa / (new_density_kgpm3 * SG_R_m2_p_s2_p_K);
 102
 103         if( _casing_pressure_Pa - pressure_Pa > 0.0 ) {
 104             Fsign = 1.0;         //outflow, pos VS
 105         } else {
 106             Fsign = -1.0;        //inflow, neg VS
 107         }
 108
 109         if( fabs(_casing_pressure_Pa - pressure_Pa) < 0.01 ) {
 110             orifice_mach = 0.0;
 111         } else {
 112             orifice_mach = sqrt(fabs (2.0*SG_cp_m2_p_s2_p_K / (SG_gamma * SG_R_m2_p_s2_p_K) * ( pow(pressure_Pa / _casing_pressure_Pa ,(SG_gamma-1)/SG_gamma ) -1 ) ) );
 113         }
 114
 115         _orifice_massflow_kgps = Fsign * _casing_pressure_Pa / sqrt(casing_temperature_K) * sqrt(SG_gamma/SG_R_m2_p_s2_p_K) * orifice_mach * pow(1+(SG_gamma-1)/2*orifice_mach*orifice_mach,-(SG_gamma+1)/(2*(SG_gamma-1))) * A_orifice;
 116
 117         double vs_fpm = Fsign * sqrt( fabs( pressure_Pa - _casing_pressure_Pa ) ) * Factor_cal;
 118         double vs_kts = vs_fpm / 60 * SG_FPS_TO_KT;
 119         double vs_mps = vs_fpm / 60 * SG_FEET_TO_METER;
 120
 121         _speed_fpm_node
 122           ->setDoubleValue(vs_fpm);
 123         _speed_kts_node
 124           ->setDoubleValue(vs_kts);
 125         _speed_mps_node
 126           ->setDoubleValue(vs_mps);
 127
 128         _casing_density_kgpm3 = new_density_kgpm3;
 129
 130     }
 131 }
 132
 133 // end of vertical_speed_indicator.cxx