FGWinds.h

/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
 Header:       FGWinds.h
 Author:       Jon Berndt, Andreas Gaeb, David Culp
 Date started: 5/2011
 
 ------------- Copyright (C) 2011  Jon S. Berndt (jon@jsbsim.org) -------------
 
 This program is free software; you can redistribute it and/or modify it under
 the terms of the GNU Lesser General Public License as published by the Free
 Software Foundation; either version 2 of the License, or (at your option) any
 later version.
 
 This program is distributed in the hope that it will be useful, but WITHOUT
 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
 FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public License for more
 details.
 
 You should have received a copy of the GNU Lesser General Public License along
 with this program; if not, write to the Free Software Foundation, Inc., 59
 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
 
 Further information about the GNU Lesser General Public License can also be
 found on the world wide web at http://www.gnu.org.
 
HISTORY
--------------------------------------------------------------------------------
5/2011   JSB   Created
 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
SENTRY
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
 
#ifndef FGWINDS_H
#define FGWINDS_H
 
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
INCLUDES
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
 
#include "models/FGModel.h"
#include "math/FGMatrix33.h"
#include <optional>
 
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
 
namespace JSBSim {
 
class FGTable;
 
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
CLASS DOCUMENTATION
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
 
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
CLASS DECLARATION
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
 
class FGWinds : public FGModel {
public:
 
  explicit FGWinds(FGFDMExec*);
  ~FGWinds();
  bool Run(bool Holding) override;
  bool InitModel(void) override;
  enum tType {ttNone, ttStandard, ttCulp, ttMilspec, ttTustin} turbType;
 
  // TOTAL WIND access functions (wind + gust + turbulence)
 
  virtual const FGColumnVector3& GetTotalWindNED(void) const { return vTotalWindNED; }
 
  virtual double GetTotalWindNED(int idx) const {return vTotalWindNED(idx);}
 
  // WIND access functions
 
  virtual void SetWindNED(double wN, double wE, double wD) { vWindNED(1)=wN; vWindNED(2)=wE; vWindNED(3)=wD;}
 
  virtual void SetWindNED(int idx, double wind) { vWindNED(idx)=wind;}
 
  virtual void SetWindNED(const FGColumnVector3& wind) { vWindNED=wind; }
 
  virtual const FGColumnVector3& GetWindNED(void) const { return vWindNED; }
 
  virtual double GetWindNED(int idx) const {return vWindNED(idx);}
 
  virtual double GetWindPsi(void) const { return psiw; }
 
  virtual void SetWindPsi(double dir);
 
  virtual void SetWindspeed(double speed);
 
  virtual double GetWindspeed(void) const;
 
  // GUST access functions
 
  virtual void SetGustNED(int idx, double gust) { vGustNED(idx)=gust;}
 
  virtual void SetTurbNED(int idx, double turb) { vTurbulenceNED(idx)=turb;}
 
  virtual void SetGustNED(double gN, double gE, double gD) { vGustNED(eNorth)=gN; vGustNED(eEast)=gE; vGustNED(eDown)=gD;}
 
  virtual double GetGustNED(int idx) const {return vGustNED(idx);}
 
  virtual double GetTurbNED(int idx) const {return vTurbulenceNED(idx);}
 
  virtual const FGColumnVector3& GetGustNED(void) const {return vGustNED;}
 
  virtual void   SetTurbType(tType tt) {turbType = tt;}
  virtual tType  GetTurbType() const {return turbType;}
 
  virtual void   SetTurbGain(double tg) {TurbGain = tg;}
  virtual double GetTurbGain() const {return TurbGain;}
 
  virtual void   SetTurbRate(double tr) {TurbRate = tr;}
  virtual double GetTurbRate() const {return TurbRate;}
 
  virtual void   SetRhythmicity(double r) {Rhythmicity=r;}
  virtual double GetRhythmicity() const {return Rhythmicity;}
 
  virtual double GetTurbPQR(int idx) const {return vTurbPQR(idx);}
  virtual double GetTurbMagnitude(void) const {return vTurbulenceNED.Magnitude();}
  virtual double GetTurbDirection(void) const {return TurbDirection;}
  virtual const FGColumnVector3& GetTurbPQR(void) const {return vTurbPQR;}
 
  virtual void   SetWindspeed20ft(double ws) { windspeed_at_20ft = ws;}
  virtual double GetWindspeed20ft() const { return windspeed_at_20ft;}
 
  virtual void   SetProbabilityOfExceedence( int idx) {probability_of_exceedence_index = idx;}
  virtual int    GetProbabilityOfExceedence() const { return probability_of_exceedence_index;}
 
  // Stores data defining a 1 - cosine gust profile that builds up, holds steady
  // and fades out over specified durations.
  struct OneMinusCosineProfile {
    bool Running;           
    double elapsedTime;     
    double startupDuration; 
    double steadyDuration;  
    double endDuration;     
    OneMinusCosineProfile() 
    {
      elapsedTime = 0.0;
      Running = false;
      startupDuration = 2;
      steadyDuration = 4;
      endDuration = 2;
    }
  };
 
  enum eGustFrame {gfNone=0, gfBody, gfWind, gfLocal};
 
  struct OneMinusCosineGust {
    FGColumnVector3 vWind;                    
    FGColumnVector3 vWindTransformed;         
    double magnitude;                         
    eGustFrame gustFrame;                     
    struct OneMinusCosineProfile gustProfile; 
    OneMinusCosineGust()                      
    {
      vWind.InitMatrix(0.0);
      gustFrame = gfLocal;
      magnitude = 1.0;
    };
  };
 
  struct UpDownBurst {
    double ringLatitude;                           
    double ringLongitude;                          
    double ringAltitude;                           
    double ringRadius;                             
    double ringCoreRadius;                         
    double circulation;                            
    struct OneMinusCosineProfile oneMCosineProfile;
    UpDownBurst() {                                
      ringLatitude = ringLongitude = 0.0;
      ringAltitude = 1000.0;
      ringRadius = 2000.0;
      ringCoreRadius = 100.0;
      circulation = 100000.0;
    }
  };
 
  // 1 - Cosine gust setters
  virtual void StartGust(bool running) {oneMinusCosineGust.gustProfile.Running = running;}
  virtual void StartupGustDuration(double dur) {oneMinusCosineGust.gustProfile.startupDuration = dur;}
  virtual void SteadyGustDuration(double dur) {oneMinusCosineGust.gustProfile.steadyDuration = dur;}
  virtual void EndGustDuration(double dur) {oneMinusCosineGust.gustProfile.endDuration = dur;}
  virtual void GustMagnitude(double mag) {oneMinusCosineGust.magnitude = mag;}
  virtual void GustFrame(eGustFrame gFrame) {oneMinusCosineGust.gustFrame = gFrame;}
  virtual void GustXComponent(double x) {oneMinusCosineGust.vWind(eX) = x;}
  virtual void GustYComponent(double y) {oneMinusCosineGust.vWind(eY) = y;}
  virtual void GustZComponent(double z) {oneMinusCosineGust.vWind(eZ) = z;}
 
  // Up- Down-burst functions
  void NumberOfUpDownburstCells(int num);
 
  struct Inputs {
    double V;
    double wingspan;
    double DistanceAGL;
    double AltitudeASL;
    double longitude;
    double latitude;
    double planetRadius;
    FGMatrix33 Tl2b;
    FGMatrix33 Tw2b;
    double totalDeltaT;
  } in;
 
private:
 
  double MagnitudedAccelDt, MagnitudeAccel, Magnitude, TurbDirection;
  //double h;
  double TurbGain;
  double TurbRate;
  double Rhythmicity;
  double wind_from_clockwise;
  double spike, target_time, strength;
  FGColumnVector3 vTurbulenceGrad;
  FGColumnVector3 vBodyTurbGrad;
  FGColumnVector3 vTurbPQR;
 
  struct OneMinusCosineGust oneMinusCosineGust;
  std::vector <struct UpDownBurst*> UpDownBurstCells;
 
  // Dryden turbulence model
  double windspeed_at_20ft; 
  int probability_of_exceedence_index; 
  FGTable *POE_Table; 
 
  // keep values from last timesteps
  // TODO maybe use deque?
  double xi_u_km1, nu_u_km1;
  double xi_v_km1, xi_v_km2, nu_v_km1, nu_v_km2;
  double xi_w_km1, xi_w_km2, nu_w_km1, nu_w_km2;
  double xi_p_km1, nu_p_km1;
  double xi_q_km1, xi_r_km1;
 
  double psiw;
  FGColumnVector3 vTotalWindNED;
  FGColumnVector3 vWindNED;
  FGColumnVector3 vGustNED;
  FGColumnVector3 vCosineGust;
  FGColumnVector3 vBurstGust;
  FGColumnVector3 vTurbulenceNED;
 
  std::optional<unsigned int> RandomSeed;
  std::shared_ptr<RandomNumberGenerator> generator;
 
  void SetRandomSeed(int sr);
  int  GetRandomSeed(void) const;
 
  void Turbulence(double h);
  void UpDownBurst();
 
  void CosineGust();
  double CosineGustProfile( double startDuration, double steadyDuration,
                            double endDuration, double elapsedTime);
  double DistanceFromRingCenter(double lat, double lon);
 
  virtual void bind(void);
  void Debug(int from) override;
};
 
} // namespace JSBSim
 
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
#endif