Add the ability to control turbulence rate as well as magnitude.

Square the normalized direction acceleration for the y and z axes, so
that turbulence predominantly affects pitch.

Bind to the /environment/turbulence/magnitude-norm and
/environment/turbulence/rate-hz properties in FlightGear.

src/Environment/environment.cxx

@@ -126,7 +126,8 @@ FGEnvironment::FGEnvironment()
     wind_from_north_fps(0),
     wind_from_east_fps(0),
     wind_from_down_fps(0),
-    turbulence_norm(0)
+    turbulence_magnitude_norm(0),
+    turbulence_rate_hz(1)
 {
   _setup_tables();
   _recalc_density();
@@ -157,14 +158,15 @@ FGEnvironment::copy (const FGEnvironment &env)
     wind_from_north_fps = env.wind_from_north_fps;
     wind_from_east_fps = env.wind_from_east_fps;
     wind_from_down_fps = env.wind_from_down_fps;
-    turbulence_norm = env.turbulence_norm;
+    turbulence_magnitude_norm = env.turbulence_magnitude_norm;
+    turbulence_rate_hz = env.turbulence_rate_hz;
 }
 
 static inline bool
 maybe_copy_value (FGEnvironment * env, const SGPropertyNode * node,
                   const char * name, void (FGEnvironment::*setter)(double))
 {
-    const SGPropertyNode * child = node->getChild(name);
+    const SGPropertyNode * child = node->getNode(name);
                                 // fragile: depends on not being typed
                                 // as a number
     if (child != 0 && child->getStringValue()[0] != '\0') {
@@ -205,8 +207,11 @@ FGEnvironment::read (const SGPropertyNode * node)
     maybe_copy_value(this, node, "elevation-ft",
                      &FGEnvironment::set_elevation_ft);
 
-    maybe_copy_value(this, node, "turbulence-norm",
-                     &FGEnvironment::set_turbulence_norm);
+    maybe_copy_value(this, node, "turbulence/magnitude-norm",
+                     &FGEnvironment::set_turbulence_magnitude_norm);
+
+    maybe_copy_value(this, node, "turbulence/rate-hz",
+                     &FGEnvironment::set_turbulence_rate_hz);
 }
 
 
@@ -289,9 +294,15 @@ FGEnvironment::get_wind_from_down_fps () const
 }
 
 double
-FGEnvironment::get_turbulence_norm () const
+FGEnvironment::get_turbulence_magnitude_norm () const
 {
-  return turbulence_norm;
+  return turbulence_magnitude_norm;
+}
+
+double
+FGEnvironment::get_turbulence_rate_hz () const
+{
+  return turbulence_rate_hz;
 }
 
 double
@@ -394,9 +405,15 @@ FGEnvironment::set_wind_from_down_fps (double d)
 }
 
 void
-FGEnvironment::set_turbulence_norm (double t)
+FGEnvironment::set_turbulence_magnitude_norm (double t)
 {
-  turbulence_norm = t;
+  turbulence_magnitude_norm = t;
+}
+
+void
+FGEnvironment::set_turbulence_rate_hz (double r)
+{
+  turbulence_rate_hz = r;
 }
 
 void
@@ -588,9 +605,14 @@ interpolate (const FGEnvironment * env1, const FGEnvironment * env2,
                    env2->get_elevation_ft(),
                    fraction));
 
-    result->set_turbulence_norm
-        (do_interp(env1->get_turbulence_norm(),
-                   env2->get_turbulence_norm(),
+    result->set_turbulence_magnitude_norm
+        (do_interp(env1->get_turbulence_magnitude_norm(),
+                   env2->get_turbulence_magnitude_norm(),
+                   fraction));
+
+    result->set_turbulence_rate_hz
+        (do_interp(env1->get_turbulence_rate_hz(),
+                   env2->get_turbulence_rate_hz(),
                    fraction));
 }
 

src/Environment/environment.hxx

@@ -72,7 +72,8 @@ public:
   virtual double get_wind_from_east_fps () const;
   virtual double get_wind_from_down_fps () const;
 
-  virtual double get_turbulence_norm () const;
+  virtual double get_turbulence_magnitude_norm () const;
+  virtual double get_turbulence_rate_hz () const;
 
   virtual void set_visibility_m (double v);
 
@@ -89,7 +90,8 @@ public:
   virtual void set_wind_from_east_fps (double e);
   virtual void set_wind_from_down_fps (double d);
 
-  virtual void set_turbulence_norm (double t);
+  virtual void set_turbulence_magnitude_norm (double t);
+  virtual void set_turbulence_rate_hz (double t);
 
   virtual double get_elevation_ft () const;
   virtual void set_elevation_ft (double elevation_ft);
@@ -120,7 +122,8 @@ private:
   double pressure_inhg;
   double density_slugft3;
 
-  double turbulence_norm;
+  double turbulence_magnitude_norm;
+  double turbulence_rate_hz;
 
   double wind_from_heading_deg;
   double wind_speed_kt;

src/Environment/environment_mgr.cxx

@@ -106,10 +106,14 @@ FGEnvironmentMgr::bind ()
 	&FGEnvironment::get_wind_from_down_fps,
 	&FGEnvironment::set_wind_from_down_fps);
   fgSetArchivable("/environment/wind-from-down-fps");
-  fgTie("/environment/turbulence-norm", _environment,
-        &FGEnvironment::get_turbulence_norm,
-        &FGEnvironment::set_turbulence_norm);
-  fgSetArchivable("/environment/turbulence-norm");
+  fgTie("/environment/turbulence/magnitude-norm", _environment,
+        &FGEnvironment::get_turbulence_magnitude_norm,
+        &FGEnvironment::set_turbulence_magnitude_norm);
+  fgSetArchivable("/environment/turbulence/magnitude-norm");
+  fgTie("/environment/turbulence/rate_hz", _environment,
+        &FGEnvironment::get_turbulence_rate_hz,
+        &FGEnvironment::set_turbulence_rate_hz);
+  fgSetArchivable("/environment/turbulence/rate_hz");
 
   for (int i = 0; i < MAX_CLOUD_LAYERS; i++) {
     char buf[128];

src/FDM/JSBSim/FGAtmosphere.cpp

@@ -92,6 +92,7 @@ FGAtmosphere::FGAtmosphere(FGFDMExec* fdmex) : FGModel(fdmex)
   turbType = ttStandard;
 //   turbType = ttBerndt;
   TurbGain = 0.0;
+  TurbRate = 1.0;
   
   bind();
   Debug(0);
@@ -272,11 +273,13 @@ void FGAtmosphere::Turbulence(void)
                                 // away from the peaks
     MagnitudedAccelDt = ((MagnitudedAccelDt - Magnitude) /
                          (1 + fabs(Magnitude)));
-    MagnitudeAccel    += MagnitudedAccelDt*rate*State->Getdt();
+    MagnitudeAccel    += MagnitudedAccelDt*rate*TurbRate*State->Getdt();
     Magnitude         += MagnitudeAccel*rate*State->Getdt();
 
     vDirectiondAccelDt.Normalize();
-    vDirectionAccel += vDirectiondAccelDt*rate*State->Getdt();
+    vDirectiondAccelDt(eY) *= vDirectiondAccelDt(eY);
+    vDirectiondAccelDt(eZ) *= vDirectiondAccelDt(eZ);
+    vDirectionAccel += vDirectiondAccelDt*rate*TurbRate*State->Getdt();
     vDirectionAccel.Normalize();
     vDirection      += vDirectionAccel*rate*State->Getdt();
 
@@ -284,7 +287,7 @@ void FGAtmosphere::Turbulence(void)
     
                                 // Diminish turbulence within three wingspans
                                 // of the ground
-    vTurbulence = TurbGain*Magnitude * vDirection;
+    vTurbulence = TurbGain * Magnitude * vDirection;
     double HOverBMAC = Position->GetHOverBMAC();
     if (HOverBMAC < 3.0)
         vTurbulence *= (HOverBMAC / 3.0) * (HOverBMAC / 3.0);

src/FDM/JSBSim/FGAtmosphere.h

@@ -150,6 +150,7 @@ public:
   inline double GetWindPsi(void) const { return psiw; }
   
   inline void SetTurbGain(double tt) {TurbGain = tt;}
+  inline void SetTurbRate(double tt) {TurbRate = tt;}
   
   inline double GetTurbPQR(int idx) const {return vTurbPQR(idx);}
   inline FGColumnVector3& GetTurbPQR(void) {return vTurbPQR;}
@@ -176,6 +177,7 @@ private:
   
   double MagnitudedAccelDt, MagnitudeAccel, Magnitude;
   double TurbGain;
+  double TurbRate;
   FGColumnVector3 vDirectiondAccelDt;
   FGColumnVector3 vDirectionAccel;
   FGColumnVector3 vDirection;

src/FDM/JSBSim/JSBSim.cxx

@@ -182,7 +182,8 @@ FGJSBsim::FGJSBsim( double dt )
     temperature = fgGetNode("/environment/temperature-degc",true);
     pressure = fgGetNode("/environment/pressure-inhg",true);
     density = fgGetNode("/environment/density-slugft3",true);
-    turbulence = fgGetNode("environment/turbulence-norm",true);
+    turbulence_gain = fgGetNode("/environment/turbulence/magnitude-norm",true);
+    turbulence_rate = fgGetNode("/environment/turbulence/rate-hz",true);
     
     wind_from_north= fgGetNode("/environment/wind-from-north-fps",true);
     wind_from_east = fgGetNode("/environment/wind-from-east-fps" ,true);
@@ -201,6 +202,8 @@ FGJSBsim::~FGJSBsim(void) {
 // each subsequent iteration through the EOM
 
 void FGJSBsim::init() {
+
+    double tmp;
     
     SG_LOG( SG_FLIGHT, SG_INFO, "Starting and initializing JSBsim" );
 
@@ -216,9 +219,13 @@ void FGJSBsim::init() {
                   9.0/5.0*(temperature->getDoubleValue()+273.15) );
       Atmosphere->SetExPressure(pressure->getDoubleValue()*70.726566);
       Atmosphere->SetExDensity(density->getDoubleValue());
-      Atmosphere->SetTurbGain(turbulence->getDoubleValue() *
-                              turbulence->getDoubleValue() *
-                              100.0);
+
+      tmp = turbulence_gain->getDoubleValue();
+      Atmosphere->SetTurbGain(tmp * tmp * 100.0);
+
+      tmp = turbulence_rate->getDoubleValue();
+      Atmosphere->SetTurbRate(tmp);
+
     } else {
       Atmosphere->UseInternal();
     }
@@ -364,6 +371,7 @@ FGJSBsim::update( double dt ) {
 // Convert from the FGInterface struct to the JSBsim generic_ struct
 
 bool FGJSBsim::copy_to_JSBsim() {
+    double tmp;
     unsigned int i;
 
     // copy control positions into the JSBsim structure
@@ -412,9 +420,12 @@ bool FGJSBsim::copy_to_JSBsim() {
                   9.0/5.0*(temperature->getDoubleValue()+273.15) );
     Atmosphere->SetExPressure(pressure->getDoubleValue()*70.726566);
     Atmosphere->SetExDensity(density->getDoubleValue());
-    Atmosphere->SetTurbGain(turbulence->getDoubleValue() *
-                            turbulence->getDoubleValue() *
-                            100.0);
+
+    tmp = turbulence_gain->getDoubleValue();
+    Atmosphere->SetTurbGain(tmp * tmp * 100.0);
+
+    tmp = turbulence_rate->getDoubleValue();
+    Atmosphere->SetTurbRate(tmp);
 
     Atmosphere->SetWindNED( wind_from_north->getDoubleValue(),
                             wind_from_east->getDoubleValue(),

src/FDM/JSBSim/JSBSim.hxx

@@ -265,7 +265,8 @@ private:
     SGPropertyNode *temperature;
     SGPropertyNode *pressure;
     SGPropertyNode *density;
-    SGPropertyNode *turbulence;
+    SGPropertyNode *turbulence_gain;
+    SGPropertyNode *turbulence_rate;
     
     SGPropertyNode *wind_from_north;
     SGPropertyNode *wind_from_east;

src/Main/options.cxx

@@ -995,7 +995,7 @@ fgOptWind( const char *arg )
 static int
 fgOptTurbulence( const char *arg)
 {
-    fgDefaultWeatherValue("turbulence-norm", atof(arg));
+    fgDefaultWeatherValue("turbulence/magnitude-norm", atof(arg));
     return FG_OPTIONS_OK;
 }
 

src/Network/native_ctrls.cxx

@@ -180,7 +180,8 @@ void FGProps2NetCtrls( FGNetCtrls *net, bool net_byte_order ) {
 
     net->wind_speed_kt = fgGetDouble("/environment/wind-speed-kt");
     net->wind_dir_deg = fgGetDouble("/environment/wind-from-heading-deg");
-    net->turbulence_norm = fgGetDouble("/environment/turbulence-norm");
+    net->turbulence_norm =
+        fgGetDouble("/environment/turbulence/magnitude-norm");
 
     // cur_fdm_state->get_ground_elev_ft() is what we want ... this
     // reports the altitude of the aircraft.
@@ -338,7 +339,7 @@ void FGNetCtrls2Props( FGNetCtrls *net, bool net_byte_order ) {
     node = fgGetNode( "/environment", true );
     node->setDoubleValue( "wind-speed-kt", net->wind_speed_kt );
     node->setDoubleValue( "wind-from-heading-deg", net->wind_dir_deg );
-    node->setDoubleValue( "turbulence-norm", net->turbulence_norm );
+    node->setDoubleValue( "turbulence/magnitude-norm", net->turbulence_norm );
     node->setBoolValue( "magnetic-variation-deg", net->magvar );
 
     // ground elevation ???