]> git.mxchange.org Git - flightgear.git/blobdiff - src/FDM/JSBSim/models/FGAtmosphere.cpp
Sync. withn JSBSim CVS
[flightgear.git] / src / FDM / JSBSim / models / FGAtmosphere.cpp
index 215b3b8a19a25dbcea2cfb8c3a1ee0ed2302f386..8acb17c1b42056c7c112169b7023c5294727a8f8 100644 (file)
@@ -7,7 +7,7 @@
  Purpose:      Models the atmosphere
  Called by:    FGSimExec
 
- ------------- Copyright (C) 1999  Jon S. Berndt (jsb@hal-pc.org) -------------
+ ------------- Copyright (C) 1999  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
@@ -48,12 +48,16 @@ INCLUDES
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
 
 #include "FGAtmosphere.h"
-#include <FGState.h>
-#include <FGFDMExec.h>
+#include "FGState.h"
+#include "FGFDMExec.h"
 #include "FGAircraft.h"
 #include "FGPropagate.h"
 #include "FGInertial.h"
-#include <input_output/FGPropertyManager.h>
+#include "input_output/FGPropertyManager.h"
+#include <iostream>
+#include <cstdlib>
+
+using namespace std;
 
 namespace JSBSim {
 
@@ -80,9 +84,10 @@ FGAtmosphere::FGAtmosphere(FGFDMExec* fdmex) : FGModel(fdmex)
   htab[7]=278385.0; //ft.
 
   MagnitudedAccelDt = MagnitudeAccel = Magnitude = 0.0;
-  SetTurbType( ttCulp );
-  TurbGain = 0.0;
-  TurbRate = 1.7;
+//  SetTurbType( ttCulp );
+  SetTurbType( ttNone );
+  TurbGain = 1.0;
+  TurbRate = 10.0;
   Rhythmicity = 0.1;
   spike = target_time = strength = 0.0;
   wind_from_clockwise = 0.0;
@@ -134,8 +139,10 @@ bool FGAtmosphere::Run(void)
   if (FGModel::Run()) return true;
   if (FDMExec->Holding()) return false;
 
+  RunPreFunctions();
+
   T_dev = 0.0;
-  h = Propagate->Geth();
+  h = Propagate->GetAltitudeASL();
 
   if (!useExternal) {
     Calculate(h);
@@ -144,6 +151,8 @@ bool FGAtmosphere::Run(void)
     CalculateDerived();
   }
 
+  RunPostFunctions();
+
   Debug(2);
   return false;
 }
@@ -274,6 +283,7 @@ void FGAtmosphere::CalculateDerived(void)
 
   vTotalWindNED = vWindNED + vGustNED + vTurbulenceNED;
 
+   // psiw (Wind heading) is the direction the wind is blowing towards
   if (vWindNED(eX) != 0.0) psiw = atan2( vWindNED(eY), vWindNED(eX) );
   if (psiw < 0) psiw += 2*M_PI;
 
@@ -336,6 +346,8 @@ static inline double square_signed (double value)
 }
 
 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+//
+// psi is the angle that the wind is blowing *towards*
 
 void FGAtmosphere::SetWindspeed(double speed)
 {
@@ -357,6 +369,8 @@ double FGAtmosphere::GetWindspeed(void) const
 }
 
 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+//
+// psi is the angle that the wind is blowing *towards*
 
 void FGAtmosphere::SetWindPsi(double dir)
 {
@@ -369,9 +383,11 @@ void FGAtmosphere::SetWindPsi(double dir)
 
 void FGAtmosphere::Turbulence(void)
 {
+  double DeltaT = rate*State->Getdt();
+
   switch (turbType) {
   case ttStandard: {
-    TurbGain = TurbGain * TurbGain * 100.0;
+    // TurbGain = TurbGain * TurbGain * 100.0; // what is this!?
 
     vDirectiondAccelDt(eX) = 1 - 2.0*(double(rand())/double(RAND_MAX));
     vDirectiondAccelDt(eY) = 1 - 2.0*(double(rand())/double(RAND_MAX));
@@ -382,8 +398,8 @@ void FGAtmosphere::Turbulence(void)
                                 // away from the peaks
     MagnitudedAccelDt = ((MagnitudedAccelDt - Magnitude) /
                          (1 + fabs(Magnitude)));
-    MagnitudeAccel    += MagnitudedAccelDt*rate*TurbRate*State->Getdt();
-    Magnitude         += MagnitudeAccel*rate*State->Getdt();
+    MagnitudeAccel    += MagnitudedAccelDt*TurbRate*DeltaT;
+    Magnitude         += MagnitudeAccel*DeltaT;
     Magnitude          = fabs(Magnitude);
 
     vDirectiondAccelDt.Normalize();
@@ -392,9 +408,9 @@ void FGAtmosphere::Turbulence(void)
     vDirectiondAccelDt(eX) = square_signed(vDirectiondAccelDt(eX));
     vDirectiondAccelDt(eY) = square_signed(vDirectiondAccelDt(eY));
 
-    vDirectionAccel += vDirectiondAccelDt*rate*TurbRate*State->Getdt();
+    vDirectionAccel += vDirectiondAccelDt*TurbRate*DeltaT;
     vDirectionAccel.Normalize();
-    vDirection      += vDirectionAccel*rate*State->Getdt();
+    vDirection      += vDirectionAccel*DeltaT;
 
     vDirection.Normalize();
 
@@ -443,42 +459,39 @@ void FGAtmosphere::Turbulence(void)
   }
   case ttBerndt: { // This is very experimental and incomplete at the moment.
 
-    TurbGain = TurbGain * TurbGain * 100.0;
-  
-    vDirectiondAccelDt(eX) = 1 - 2.0*(double(rand())/double(RAND_MAX));
-    vDirectiondAccelDt(eY) = 1 - 2.0*(double(rand())/double(RAND_MAX));
-    vDirectiondAccelDt(eZ) = 1 - 2.0*(double(rand())/double(RAND_MAX));
-
-
-    MagnitudedAccelDt = 1 - 2.0*(double(rand())/double(RAND_MAX)) - Magnitude;
-    MagnitudeAccel    += MagnitudedAccelDt*rate*State->Getdt();
-    Magnitude         += MagnitudeAccel*rate*State->Getdt();
+    vDirectiondAccelDt(eX) = GaussianRandomNumber();
+    vDirectiondAccelDt(eY) = GaussianRandomNumber();
+    vDirectiondAccelDt(eZ) = GaussianRandomNumber();
+/*
+    MagnitudedAccelDt = GaussianRandomNumber();
+    MagnitudeAccel    += MagnitudedAccelDt * DeltaT;
+    Magnitude         += MagnitudeAccel * DeltaT;
+*/
+    Magnitude         += GaussianRandomNumber() * DeltaT;
 
     vDirectiondAccelDt.Normalize();
-    vDirectionAccel += vDirectiondAccelDt*rate*State->Getdt();
+    vDirectionAccel += TurbRate * vDirectiondAccelDt * DeltaT;
     vDirectionAccel.Normalize();
-    vDirection      += vDirectionAccel*rate*State->Getdt();
+    vDirection      += vDirectionAccel*DeltaT;
 
-                                // Diminish z-vector within two wingspans
-                                // of the ground
+    // Diminish z-vector within two wingspans of the ground
     double HOverBMAC = Auxiliary->GetHOverBMAC();
-    if (HOverBMAC < 2.0)
-        vDirection(eZ) *= HOverBMAC / 2.0;
+    if (HOverBMAC < 2.0) vDirection(eZ) *= HOverBMAC / 2.0;
 
     vDirection.Normalize();
 
     vTurbulenceNED = TurbGain*Magnitude * vDirection;
     vTurbulenceGrad = TurbGain*MagnitudeAccel * vDirection;
 
-    vBodyTurbGrad = Propagate->GetTl2b()*vTurbulenceGrad;
-    vTurbPQR(eP) = vBodyTurbGrad(eY)/Aircraft->GetWingSpan();
+    vBodyTurbGrad = Propagate->GetTl2b() * vTurbulenceGrad;
+    vTurbPQR(eP) = vBodyTurbGrad(eY) / Aircraft->GetWingSpan();
     if (Aircraft->GetHTailArm() > 0)
-      vTurbPQR(eQ) = vBodyTurbGrad(eZ)/Aircraft->GetHTailArm();
+      vTurbPQR(eQ) = vBodyTurbGrad(eZ) / Aircraft->GetHTailArm();
     else
-      vTurbPQR(eQ) = vBodyTurbGrad(eZ)/10.0;
+      vTurbPQR(eQ) = vBodyTurbGrad(eZ) / 10.0;
 
     if (Aircraft->GetVTailArm() > 0)
-      vTurbPQR(eR) = vBodyTurbGrad(eX)/Aircraft->GetVTailArm();
+      vTurbPQR(eR) = vBodyTurbGrad(eX) / Aircraft->GetVTailArm();
     else
       vTurbPQR(eR) = vBodyTurbGrad(eX)/10.0;
 
@@ -565,7 +578,9 @@ void FGAtmosphere::bind(void)
 {
   typedef double (FGAtmosphere::*PMF)(int) const;
   typedef double (FGAtmosphere::*PMFv)(void) const;
+  typedef int (FGAtmosphere::*PMFt)(void) const;
   typedef void   (FGAtmosphere::*PMFd)(int,double);
+  typedef void   (FGAtmosphere::*PMFi)(int);
   PropertyManager->Tie("atmosphere/T-R", this, (PMFv)&FGAtmosphere::GetTemperature);
   PropertyManager->Tie("atmosphere/rho-slugs_ft3", this, (PMFv)&FGAtmosphere::GetDensity);
   PropertyManager->Tie("atmosphere/P-psf", this, (PMFv)&FGAtmosphere::GetPressure);
@@ -578,7 +593,7 @@ void FGAtmosphere::bind(void)
   PropertyManager->Tie("atmosphere/sigma", this, &FGAtmosphere::GetDensityRatio);
   PropertyManager->Tie("atmosphere/delta", this, &FGAtmosphere::GetPressureRatio);
   PropertyManager->Tie("atmosphere/a-ratio", this, &FGAtmosphere::GetSoundSpeedRatio);
-  PropertyManager->Tie("atmosphere/psiw-rad", this, &FGAtmosphere::GetWindPsi);
+  PropertyManager->Tie("atmosphere/psiw-rad", this, &FGAtmosphere::GetWindPsi, &FGAtmosphere::SetWindPsi);
   PropertyManager->Tie("atmosphere/delta-T", this, &FGAtmosphere::GetDeltaT, &FGAtmosphere::SetDeltaT);
   PropertyManager->Tie("atmosphere/T-sl-dev-F", this, &FGAtmosphere::GetSLTempDev, &FGAtmosphere::SetSLTempDev);
   PropertyManager->Tie("atmosphere/density-altitude", this, &FGAtmosphere::GetDensityAltitude);
@@ -589,8 +604,6 @@ void FGAtmosphere::bind(void)
                                                           (PMFd)&FGAtmosphere::SetWindNED);
   PropertyManager->Tie("atmosphere/wind-down-fps",  this, eDown, (PMF)&FGAtmosphere::GetWindNED,
                                                           (PMFd)&FGAtmosphere::SetWindNED);
-  PropertyManager->Tie("atmosphere/wind-from-cw", this, &FGAtmosphere::GetWindFromClockwise,
-                                                        &FGAtmosphere::SetWindFromClockwise);
   PropertyManager->Tie("atmosphere/wind-mag-fps", this, &FGAtmosphere::GetWindspeed,
                                                         &FGAtmosphere::SetWindspeed);
   PropertyManager->Tie("atmosphere/total-wind-north-fps", this, eNorth, (PMF)&FGAtmosphere::GetTotalWindNED);
@@ -604,9 +617,17 @@ void FGAtmosphere::bind(void)
   PropertyManager->Tie("atmosphere/gust-down-fps",  this, eDown, (PMF)&FGAtmosphere::GetGustNED,
                                                           (PMFd)&FGAtmosphere::SetGustNED);
 
+  PropertyManager->Tie("atmosphere/turb-north-fps", this, eNorth, (PMF)&FGAtmosphere::GetTurbNED,
+                                                          (PMFd)&FGAtmosphere::SetTurbNED);
+  PropertyManager->Tie("atmosphere/turb-east-fps",  this, eEast, (PMF)&FGAtmosphere::GetTurbNED,
+                                                          (PMFd)&FGAtmosphere::SetTurbNED);
+  PropertyManager->Tie("atmosphere/turb-down-fps",  this, eDown, (PMF)&FGAtmosphere::GetTurbNED,
+                                                          (PMFd)&FGAtmosphere::SetTurbNED);
+
   PropertyManager->Tie("atmosphere/p-turb-rad_sec", this,1, (PMF)&FGAtmosphere::GetTurbPQR);
   PropertyManager->Tie("atmosphere/q-turb-rad_sec", this,2, (PMF)&FGAtmosphere::GetTurbPQR);
   PropertyManager->Tie("atmosphere/r-turb-rad_sec", this,3, (PMF)&FGAtmosphere::GetTurbPQR);
+  PropertyManager->Tie("atmosphere/turb-type", this, (PMFt)&FGAtmosphere::GetTurbType, (PMFi)&FGAtmosphere::SetTurbType);
   PropertyManager->Tie("atmosphere/turb-rate", this, &FGAtmosphere::GetTurbRate, &FGAtmosphere::SetTurbRate);
   PropertyManager->Tie("atmosphere/turb-gain", this, &FGAtmosphere::GetTurbGain, &FGAtmosphere::SetTurbGain);
   PropertyManager->Tie("atmosphere/turb-rhythmicity", this, &FGAtmosphere::GetRhythmicity,