src/FDM/JSBSim/models/atmosphere/FGWinds.h

   1 /*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   2
   3  Header:       FGWinds.h
   4  Author:       Jon Berndt, Andreas Gaeb, David Culp
   5  Date started: 5/2011
   6
   7  ------------- Copyright (C) 2011  Jon S. Berndt (jon@jsbsim.org) -------------
   8
   9  This program is free software; you can redistribute it and/or modify it under
  10  the terms of the GNU Lesser General Public License as published by the Free Software
  11  Foundation; either version 2 of the License, or (at your option) any later
  12  version.
  13
  14  This program is distributed in the hope that it will be useful, but WITHOUT
  15  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
  16  FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public License for more
  17  details.
  18
  19  You should have received a copy of the GNU Lesser General Public License along with
  20  this program; if not, write to the Free Software Foundation, Inc., 59 Temple
  21  Place - Suite 330, Boston, MA  02111-1307, USA.
  22
  23  Further information about the GNU Lesser General Public License can also be found on
  24  the world wide web at http://www.gnu.org.
  25
  26 HISTORY
  27 --------------------------------------------------------------------------------
  28 5/2011   JSB   Created
  29
  30 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  31 SENTRY
  32 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
  33
  34 #ifndef FGWINDS_H
  35 #define FGWINDS_H
  36
  37 /*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  38 INCLUDES
  39 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
  40
  41 #include "models/FGModel.h"
  42 #include "math/FGColumnVector3.h"
  43 #include "math/FGMatrix33.h"
  44 #include "math/FGTable.h"
  45
  46 /*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  47 DEFINITIONS
  48 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
  49
  50 #define ID_WINDS "$Id: FGWinds.h,v 1.9 2011/11/19 14:14:57 bcoconni Exp $"
  51
  52 /*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  53 FORWARD DECLARATIONS
  54 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
  55
  56 namespace JSBSim {
  57
  58 /*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  59 CLASS DOCUMENTATION
  60 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
  61
  62 /** Models atmospheric disturbances: winds, gusts, turbulence, downbursts, etc.
  63
  64     Various turbulence models are available. They are specified
  65     via the property <tt>atmosphere/turb-type</tt>. The following models are
  66     available:
  67     - 0: ttNone (turbulence disabled)
  68     - 1: ttStandard
  69     - 2: ttCulp
  70     - 3: ttMilspec (Dryden spectrum)
  71     - 4: ttTustin (Dryden spectrum)
  72
  73     The Milspec and Tustin models are described in the Yeager report cited below.
  74     They both use a Dryden spectrum model whose parameters (scale lengths and intensities)
  75     are modelled according to MIL-F-8785C. Parameters are modelled differently
  76     for altitudes below 1000ft and above 2000ft, for altitudes in between they
  77     are interpolated linearly.
  78
  79     The two models differ in the implementation of the transfer functions
  80     described in the milspec.
  81
  82     To use one of these two models, set <tt>atmosphere/turb-type</tt> to 4 resp. 5,
  83     and specify values for <tt>atmosphere/turbulence/milspec/windspeed_at_20ft_AGL-fps<tt>
  84     and <tt>atmosphere/turbulence/milspec/severity<tt> (the latter corresponds to
  85     the probability of exceedence curves from Fig.&nbsp;7 of the milspec, allowable
  86     range is 0 (disabled) to 7). <tt>atmosphere/psiw-rad</tt> is respected as well;
  87     note that you have to specify a positive wind magnitude to prevent psiw from
  88     being reset to zero.
  89
  90     Reference values (cf. figures 7 and 9 from the milspec):
  91     <table>
  92       <tr><td><b>Intensity</b></td>
  93           <td><b><tt>windspeed_at_20ft_AGL-fps</tt></b></td>
  94           <td><b><tt>severity</tt></b></td></tr>
  95       <tr><td>light</td>
  96           <td>25 (15 knots)</td>
  97           <td>3</td></tr>
  98       <tr><td>moderate</td>
  99           <td>50 (30 knots)</td>
 100           <td>4</td></tr>
 101       <tr><td>severe</td>
 102           <td>75 (45 knots)</td>
 103           <td>6</td></tr>
 104     </table>
 105
 106     @see Yeager, Jessie C.: "Implementation and Testing of Turbulence Models for
 107          the F18-HARV" (<a
 108          href="http://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/19980028448_1998081596.pdf">
 109          pdf</a>), NASA CR-1998-206937, 1998
 110
 111     @see MIL-F-8785C: Military Specification: Flying Qualities of Piloted Aircraft
 112
 113 */
 114
 115 /*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 116 CLASS DECLARATION
 117 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
 118
 119 class FGWinds : public FGModel {
 120 public:
 121
 122   /// Constructor
 123   FGWinds(FGFDMExec*);
 124   /// Destructor
 125   ~FGWinds();
 126   /** Runs the winds model; called by the Executive
 127       Can pass in a value indicating if the executive is directing the simulation to Hold.
 128       @param Holding if true, the executive has been directed to hold the sim from
 129                      advancing time. Some models may ignore this flag, such as the Input
 130                      model, which may need to be active to listen on a socket for the
 131                      "Resume" command to be given.
 132       @return false if no error */
 133   bool Run(bool Holding);
 134   bool InitModel(void);
 135   enum tType {ttNone, ttStandard, ttCulp, ttMilspec, ttTustin} turbType;
 136
 137   // TOTAL WIND access functions (wind + gust + turbulence)
 138
 139   /// Retrieves the total wind components in NED frame.
 140   virtual const FGColumnVector3& GetTotalWindNED(void) const { return vTotalWindNED; }
 141
 142   /// Retrieves a total wind component in NED frame.
 143   virtual double GetTotalWindNED(int idx) const {return vTotalWindNED(idx);}
 144
 145   // WIND access functions
 146
 147   /// Sets the wind components in NED frame.
 148   virtual void SetWindNED(double wN, double wE, double wD) { vWindNED(1)=wN; vWindNED(2)=wE; vWindNED(3)=wD;}
 149
 150   /// Sets a wind component in NED frame.
 151   virtual void SetWindNED(int idx, double wind) { vWindNED(idx)=wind;}
 152
 153   /// Sets the wind components in NED frame.
 154   virtual void SetWindNED(const FGColumnVector3& wind) { vWindNED=wind; }
 155
 156   /// Retrieves the wind components in NED frame.
 157   virtual const FGColumnVector3& GetWindNED(void) const { return vWindNED; }
 158
 159   /// Retrieves a wind component in NED frame.
 160   virtual double GetWindNED(int idx) const {return vWindNED(idx);}
 161
 162   /** Retrieves the direction that the wind is coming from.
 163       The direction is defined as north=0 and increases counterclockwise.
 164       The wind heading is returned in radians.*/
 165   virtual double GetWindPsi(void) const { return psiw; }
 166
 167   /** Sets the direction that the wind is coming from.
 168       The direction is defined as north=0 and increases counterclockwise to 2*pi (radians). The
 169       vertical component of wind is assumed to be zero - and is forcibly set to zero. This function
 170       sets the vWindNED vector components based on the supplied direction. The magnitude of
 171       the wind set in the vector is preserved (assuming the vertical component is non-zero).
 172       @param dir wind direction in the horizontal plane, in radians.*/
 173   virtual void SetWindPsi(double dir);
 174
 175   virtual void SetWindspeed(double speed);
 176
 177   virtual double GetWindspeed(void) const;
 178
 179   // GUST access functions
 180
 181   /// Sets a gust component in NED frame.
 182   virtual void SetGustNED(int idx, double gust) { vGustNED(idx)=gust;}
 183
 184   /// Sets a turbulence component in NED frame.
 185   virtual void SetTurbNED(int idx, double turb) { vTurbulenceNED(idx)=turb;}
 186
 187   /// Sets the gust components in NED frame.
 188   virtual void SetGustNED(double gN, double gE, double gD) { vGustNED(eNorth)=gN; vGustNED(eEast)=gE; vGustNED(eDown)=gD;}
 189
 190   /// Retrieves a gust component in NED frame.
 191   virtual double GetGustNED(int idx) const {return vGustNED(idx);}
 192
 193   /// Retrieves a turbulence component in NED frame.
 194   virtual double GetTurbNED(int idx) const {return vTurbulenceNED(idx);}
 195
 196   /// Retrieves the gust components in NED frame.
 197   virtual const FGColumnVector3& GetGustNED(void) const {return vGustNED;}
 198
 199   /** Turbulence models available: ttNone, ttStandard, ttBerndt, ttCulp, ttMilspec, ttTustin */
 200   virtual void   SetTurbType(tType tt) {turbType = tt;}
 201   virtual tType  GetTurbType() const {return turbType;}
 202
 203   virtual void   SetTurbGain(double tg) {TurbGain = tg;}
 204   virtual double GetTurbGain() const {return TurbGain;}
 205
 206   virtual void   SetTurbRate(double tr) {TurbRate = tr;}
 207   virtual double GetTurbRate() const {return TurbRate;}
 208
 209   virtual void   SetRhythmicity(double r) {Rhythmicity=r;}
 210   virtual double GetRhythmicity() const {return Rhythmicity;}
 211
 212   virtual double GetTurbPQR(int idx) const {return vTurbPQR(idx);}
 213   virtual double GetTurbMagnitude(void) const {return Magnitude;}
 214   virtual const FGColumnVector3& GetTurbDirection(void) const {return vDirection;}
 215   virtual const FGColumnVector3& GetTurbPQR(void) const {return vTurbPQR;}
 216
 217   virtual void   SetWindspeed20ft(double ws) { windspeed_at_20ft = ws;}
 218   virtual double GetWindspeed20ft() const { return windspeed_at_20ft;}
 219
 220   /// allowable range: 0-7, 3=light, 4=moderate, 6=severe turbulence
 221   virtual void   SetProbabilityOfExceedence( int idx) {probability_of_exceedence_index = idx;}
 222   virtual int    GetProbabilityOfExceedence() const { return probability_of_exceedence_index;}
 223
 224   // Stores data defining a 1 - cosine gust profile that builds up, holds steady
 225   // and fades out over specified durations.
 226   struct OneMinusCosineProfile {
 227     bool Running;           ///<- This flag is set true through FGWinds::StartGust().
 228     double elapsedTime;     ///<- Stores the elapsed time for the ongoing gust.
 229     double startupDuration; ///<- Specifies the time it takes for the gust startup transient.
 230     double steadyDuration;  ///<- Specifies the duration of the steady gust.
 231     double endDuration;     ///<- Specifies the time it takes for the gust to subsude.
 232     OneMinusCosineProfile() ///<- The constructor.
 233     {
 234       elapsedTime = 0.0;
 235       Running = false;
 236       startupDuration = 2;
 237       steadyDuration = 4;
 238       endDuration = 2;
 239     }
 240   };
 241
 242   enum eGustFrame {gfNone=0, gfBody, gfWind, gfLocal};
 243
 244   /// Stores the information about a single one minus cosine gust instance.
 245   struct OneMinusCosineGust {
 246     FGColumnVector3 vWind;                    ///<- The input normalized wind vector.
 247     FGColumnVector3 vWindTransformed;         ///<- The transformed normal vector at the time the gust is started.
 248     double magnitude;                         ///<- The magnitude of the wind vector.
 249     eGustFrame gustFrame;                     ///<- The frame that the wind vector is specified in.
 250     struct OneMinusCosineProfile gustProfile; ///<- The gust shape (profile) data for this gust.
 251     OneMinusCosineGust()                      ///<- Constructor.
 252     {
 253       vWind.InitMatrix(0.0);
 254       gustFrame = gfLocal;
 255       magnitude = 1.0;
 256     };
 257   };
 258
 259   /// Stores information about a specified Up- or Down-burst.
 260   struct UpDownBurst {
 261     double ringLatitude;                           ///<- The latitude of the downburst run (radians)
 262     double ringLongitude;                          ///<- The longitude of the downburst run (radians)
 263     double ringAltitude;                           ///<- The altitude of the ring (feet).
 264     double ringRadius;                             ///<- The radius of the ring (feet).
 265     double ringCoreRadius;                         ///<- The cross-section "core" radius of the ring (feet).
 266     double circulation;                            ///<- The circulation (gamma) (feet-squared per second).
 267     struct OneMinusCosineProfile oneMCosineProfile;///<- A gust profile structure.
 268     UpDownBurst() {                                ///<- Constructor
 269       ringLatitude = ringLongitude = 0.0;
 270       ringAltitude = 1000.0;
 271       ringRadius = 2000.0;
 272       ringCoreRadius = 100.0;
 273       circulation = 100000.0;
 274     }
 275   };
 276
 277   // 1 - Cosine gust setters
 278   /// Initiates the execution of the gust.
 279   virtual void StartGust(bool running) {oneMinusCosineGust.gustProfile.Running = running;}
 280   ///Specifies the duration of the startup portion of the gust.
 281   virtual void StartupGustDuration(double dur) {oneMinusCosineGust.gustProfile.startupDuration = dur;}
 282   ///Specifies the length of time that the gust is at a steady, full strength.
 283   virtual void SteadyGustDuration(double dur) {oneMinusCosineGust.gustProfile.steadyDuration = dur;}
 284   /// Specifies the length of time it takes for the gust to return to zero velocity.
 285   virtual void EndGustDuration(double dur) {oneMinusCosineGust.gustProfile.endDuration = dur;}
 286   /// Specifies the magnitude of the gust in feet/second.
 287   virtual void GustMagnitude(double mag) {oneMinusCosineGust.magnitude = mag;}
 288   /** Specifies the frame that the gust direction vector components are specified in. The
 289       body frame is defined with the X direction forward, and the Y direction positive out
 290       the right wing. The wind frame is defined with the X axis pointing into the velocity
 291       vector, the Z axis perpendicular to the X axis, in the aircraft XZ plane, and the Y
 292       axis completing the system. The local axis is a navigational frame with X pointing north,
 293       Y pointing east, and Z pointing down. This is a locally vertical, locally horizontal
 294       frame, with the XY plane tangent to the geocentric surface. */
 295   virtual void GustFrame(eGustFrame gFrame) {oneMinusCosineGust.gustFrame = gFrame;}
 296   /// Specifies the X component of velocity in the specified gust frame (ft/sec).
 297   virtual void GustXComponent(double x) {oneMinusCosineGust.vWind(eX) = x;}
 298   /// Specifies the Y component of velocity in the specified gust frame (ft/sec).
 299   virtual void GustYComponent(double y) {oneMinusCosineGust.vWind(eY) = y;}
 300   /// Specifies the Z component of velocity in the specified gust frame (ft/sec).
 301   virtual void GustZComponent(double z) {oneMinusCosineGust.vWind(eZ) = z;}
 302
 303   // Up- Down-burst functions
 304   void NumberOfUpDownburstCells(int num);
 305
 306   struct Inputs {
 307     double V;
 308     double wingspan;
 309     double DistanceAGL;
 310     double AltitudeASL;
 311     double longitude;
 312     double latitude;
 313     double planetRadius;
 314     FGMatrix33 Tl2b;
 315     FGMatrix33 Tw2b;
 316     double totalDeltaT;
 317   } in;
 318
 319 private:
 320
 321   double MagnitudedAccelDt, MagnitudeAccel, Magnitude;
 322   double h;
 323   double TurbGain;
 324   double TurbRate;
 325   double Rhythmicity;
 326   double wind_from_clockwise;
 327   double spike, target_time, strength;
 328   FGColumnVector3 vDirectiondAccelDt;
 329   FGColumnVector3 vDirectionAccel;
 330   FGColumnVector3 vDirection;
 331   FGColumnVector3 vTurbulenceGrad;
 332   FGColumnVector3 vBodyTurbGrad;
 333   FGColumnVector3 vTurbPQR;
 334
 335   struct OneMinusCosineGust oneMinusCosineGust;
 336   std::vector <struct UpDownBurst*> UpDownBurstCells;
 337
 338   // Dryden turbulence model
 339   double windspeed_at_20ft; ///< in ft/s
 340   int probability_of_exceedence_index; ///< this is bound as the severity property
 341   FGTable *POE_Table; ///< probability of exceedence table
 342
 343   double psiw;
 344   FGColumnVector3 vTotalWindNED;
 345   FGColumnVector3 vWindNED;
 346   FGColumnVector3 vGustNED;
 347   FGColumnVector3 vCosineGust;
 348   FGColumnVector3 vBurstGust;
 349   FGColumnVector3 vTurbulenceNED;
 350
 351   void Turbulence(double h);
 352   void UpDownBurst();
 353
 354   void CosineGust();
 355   double CosineGustProfile( double startDuration, double steadyDuration,
 356                             double endDuration, double elapsedTime);
 357   double DistanceFromRingCenter(double lat, double lon);
 358
 359   virtual void bind(void);
 360   void Debug(int from);
 361 };
 362
 363 } // namespace JSBSim
 364
 365 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 366 #endif
 367