Module: FGAtmosphere.cpp
Author: Jon Berndt
+ Implementation of 1959 Standard Atmosphere added by Tony Peden
Date started: 11/24/98
Purpose: Models the atmosphere
Called by: FGSimExec
HISTORY
--------------------------------------------------------------------------------
11/24/98 JSB Created
+07/23/99 TP Added implementation of 1959 Standard Atmosphere
+ Moved calculation of Mach number to FGTranslation
+********************************************************************************
+COMMENTS, REFERENCES, and NOTES
+********************************************************************************
+[1] Anderson, John D. "Introduction to Flight, Third Edition", McGraw-Hill,
+ 1989, ISBN 0-07-001641-0
********************************************************************************
INCLUDES
#include "FGPosition.h"
#include "FGAuxiliary.h"
#include "FGOutput.h"
+#include "FGDefs.h"
/*******************************************************************************
************************************ CODE **************************************
*******************************************************************************/
+
FGAtmosphere::FGAtmosphere(FGFDMExec* fdmex) : FGModel(fdmex)
{
Name = "FGAtmosphere";
+ h=0;
+ Calculate();
}
bool FGAtmosphere::Run(void)
{
if (!FGModel::Run()) { // if false then execute this Run()
- rho = 0.002377 - 7.0E-08*State->Geth()
- + 7.0E-13*State->Geth()*State->Geth()
- - 2.0E-18*State->Geth()*State->Geth()*State->Geth();
-
- State->SetMach(State->GetVt()/State->Geta());
+ h = State->Geth();
+ Calculate();
+ State->Seta(soundspeed);
} else { // skip Run() execution this time
}
return false;
}
+
float FGAtmosphere::CalcRho(float altitude)
{
- return (0.002377 - 7.0E-08*altitude
+ return (0.00237 - 7.0E-08*altitude
+ 7.0E-13*altitude*altitude
- 2.0E-18*altitude*altitude*altitude);
+}
+
+
+void FGAtmosphere::Calculate(void)
+{
+ //see reference [1]
+
+ float slope,reftemp,refpress,refdens;
+ int i=0;
+ float htab[]={0,36089,82020,154198,173882,259183,295272,344484}; //ft.
+
+ if (h <= htab[0]) {
+ h=0;
+ } else if (h >= htab[7]){
+ i = 7;
+ h = htab[7];
+ } else {
+ while (htab[i+1] < h) {
+ i++;
+ }
+ }
+
+ switch(i) {
+ case 0: // sea level
+ slope = -0.0035662; // R/ft.
+ reftemp = 518.688; // R
+ refpress = 2116.17; // psf
+ refdens = 0.0023765; // slugs/cubic ft.
+ break;
+ case 1: // 36089 ft.
+ slope = 0;
+ reftemp = 389.988;
+ refpress = 474.1;
+ refdens = 0.0007078;
+ break;
+ case 2: // 82020 ft.
+ slope = 0.00164594;
+ reftemp = 389.988;
+ refpress = 52.7838;
+ refdens = 7.8849E-5;
+ break;
+ case 3: // 154198 ft.
+ slope = 0;
+ reftemp = 508.788;
+ refpress = 2.62274;
+ refdens = 3.01379E-6;
+ break;
+ case 4: // 173882 ft.
+ slope = -0.00246891;
+ reftemp = 508.788;
+ refpress = 1.28428;
+ refdens = 1.47035e-06;
+ break;
+ case 5: // 259183 ft.
+ slope = 0;
+ reftemp = 298.188;
+ refpress = 0.0222008;
+ refdens = 4.33396e-08;
+ break;
+ case 6: // 295272 ft.
+ slope = 0.00219459;
+ reftemp = 298.188;
+ refpress = 0.00215742;
+ refdens = 4.21368e-09;
+ break;
+ case 7: // 344484 ft.
+ slope = 0;
+ reftemp = 406.188;
+ refpress = 0.000153755;
+ refdens = 2.20384e-10;
+ break;
+ }
+
+
+ if (slope == 0) {
+ temperature = reftemp;
+ pressure = refpress*exp(-GRAVITY/(reftemp*Reng)*(h-htab[i]));
+ density = refdens*exp(-GRAVITY/(reftemp*Reng)*(h-htab[i]));
+ } else {
+ temperature = reftemp+slope*(h-htab[i]);
+ pressure = refpress*pow(temperature/reftemp,-GRAVITY/(slope*Reng));
+ density = refdens*pow(temperature/reftemp,-(GRAVITY/(slope*Reng)+1));
+ }
+
+ soundspeed = sqrt(SHRATIO*Reng*temperature);
+
+}
+
+
+float FGAtmosphere::GetTemperature(float altitude)
+{
+ if (altitude != h) {
+ h = altitude;
+ Calculate();
+ }
+ return temperature;
+}
+
+
+float FGAtmosphere::GetPressure(float altitude)
+{
+ if (altitude != h) {
+ h = altitude;
+ Calculate();
+ }
+
+ return pressure;
+}
+
+float FGAtmosphere::GetDensity(float altitude)
+{
+ if (altitude != h) {
+ h = altitude;
+ Calculate();
+ }
+
+ return density;
+}
+
+
+float FGAtmosphere::GetSoundSpeed(float altitude)
+{
+ if (altitude != h) {
+ h = altitude;
+ Calculate();
+ }
+ return soundspeed;
}
Header: FGAtmosphere.h
Author: Jon Berndt
+ Implementation of 1959 Standard Atmosphere added by Tony Peden
Date started: 11/24/98
------------- Copyright (C) 1999 Jon S. Berndt (jsb@hal-pc.org) -------------
HISTORY
--------------------------------------------------------------------------------
11/24/98 JSB Created
+07/23/99 TP Added implementation of 1959 Standard Atmosphere
+ Moved calculation of Mach number to FGTranslation
+
********************************************************************************
SENTRY
*******************************************************************************/
-#ifndef FGATMOSPHERE_H
-#define FGATMOSPHERE_H
+#ifndef FGAtmosphere_H
+#define FGAtmosphere_H
/*******************************************************************************
INCLUDES
/*******************************************************************************
COMMENTS, REFERENCES, and NOTES
-*******************************************************************************/
-/**
-The equation used in this model was determined by a third order curve fit using
-Excel. The data is from the ICAO atmosphere model.
-@memo Models the atmosphere.
-@author Jon S. Berndt
-*/
+********************************************************************************
+
+[1] Anderson, John D. "Introduction to Flight, Third Edition", McGraw-Hill,
+ 1989, ISBN 0-07-001641-0
+
/*******************************************************************************
CLASS DECLARATION
*******************************************************************************/
+
using namespace std;
class FGAtmosphere : public FGModel
{
public:
+
FGAtmosphere(FGFDMExec*);
~FGAtmosphere(void);
bool Run(void);
inline float Getrho(void) {return rho;}
float CalcRho(float altitude);
+ inline float GetTemperature(void){return temperature;}
+ inline float GetDensity(void) {return density;} // use only after Run() has been called
+ inline float GetPressure(void) {return pressure;}
+ inline float GetSoundSpeed(void) {return soundspeed;}
+
+ float GetTemperature(float altitude); //Rankine, altitude in feet
+ float GetDensity(float altitude); //slugs/ft^3
+ float GetPressure(float altitude); //lbs/ft^2
+ float GetSoundSpeed(float altitude); //ft/s
+
protected:
private:
float rho;
+
+ float h;
+ float temperature;
+ float pressure;
+ float density;
+ float soundspeed;
+ void Calculate(void);
+
+
};
/******************************************************************************/
// $Log$
-// Revision 1.4 1999/07/31 02:55:24 curt
+// Revision 1.5 1999/08/17 19:18:11 curt
// Updates from Jon.
//
// Revision 1.1 1999/02/13 01:12:03 curt
// $Log$
-// Revision 1.4 1999/07/31 02:55:24 curt
+// Revision 1.5 1999/08/17 19:18:11 curt
// Updates from Jon.
//
// Revision 1.1 1999/02/13 01:12:03 curt
#define INVECCENTSQRD 1.0067395
#define INVECCENTSQRDM1 0.0067395
#define EPS 0.081819221
+#define Reng 1716 //Specific Gas Constant,ft^2/(sec^2*R)
+#define SHRATIO 1.4 //Specific Heat Ratio
+#define RADTODEG 57.29578
+#define DEGTORAD 1.745329E-2
+#define KTSTOFPS 1.68781
+#define FPSTOKTS 0.592484
+
/******************************************************************************/
Rotation = new FGRotation(this);
Position = new FGPosition(this);
Auxiliary = new FGAuxiliary(this);
-// Output = new FGOutput(this);
+ Output = new FGOutput(this);
State = new FGState(this);
if (!Rotation->InitModel()) {cerr << "Rotation model init failed"; Error+=16;}
if (!Position->InitModel()) {cerr << "Position model init failed"; Error+=32;}
if (!Auxiliary->InitModel()) {cerr << "Auxiliary model init failed"; Error+=64;}
-// if (!Output->InitModel()) {cerr << "Output model init failed"; Error+=128;}
+ if (!Output->InitModel()) {cerr << "Output model init failed"; Error+=128;}
+
+ // Schedule a model. The second arg (the integer) is the pass number. For
+ // instance, the atmosphere model gets executed every fifth pass it is called
+ // by the executive. Everything else here gets executed each pass.
Schedule(Atmosphere, 5);
Schedule(FCS, 1);
Schedule(Translation, 1);
Schedule(Position, 1);
Schedule(Auxiliary, 1);
-// Schedule(Output, 1);
+ Schedule(Output, 1);
terminate = false;
frozen = false;
return true;
}
+bool FGFDMExec::RunIC(FGInitialCondition *fgic)
+{
+ float save_dt=State->Getdt();
+ State->Setdt(0.0);
+ State->Initialize(fgic);
+ Run();
+ State->Setdt(save_dt);
+ return true;
+}
+
+
HISTORY
--------------------------------------------------------------------------------
11/17/98 JSB Created
+7/31/99 TP Added RunIC function that runs the sim so that every frame
+ begins with the IC values from the given FGInitialCondition
+ object and dt=0.
********************************************************************************
SENTRY
*******************************************************************************/
#include "FGModel.h"
+#include "FGInitialCondition.h"
using namespace std;
class FGPosition;
class FGAuxiliary;
class FGOutput;
+class FGInitialCondition;
class FGFDMExec
{
bool Initialize(void);
int Schedule(FGModel* model, int rate);
bool Run(void);
+ bool RunIC(FGInitialCondition *fgic);
void Freeze(void) {frozen = true;}
void Resume(void) {frozen = false;}
#include <iostream>
#include <ctime>
-int main(int argc, char** argv)
+void main(int argc, char** argv)
{
FGFDMExec* FDMExec;
}
+//***************************************************************************
+//
+// Reset: Assume all angles READ FROM FILE IN DEGREES !!
+//
+
bool FGState::Reset(string path, string fname)
{
string resetDef;
float U, V, W;
float phi, tht, psi;
- float alpha, beta, gamma;
- float Q0, Q1, Q2, Q3;
- float T[4][4];
resetDef = path + "/" + FDMExec->GetAircraft()->GetAircraftName() + "/" + fname;
resetfile >> h;
resetfile.close();
- Initialize(U, V, W, phi, tht, psi, latitude, longitude, h);
+ Initialize(U, V, W, phi*DEGTORAD, tht*DEGTORAD, psi*DEGTORAD,
+ latitude*DEGTORAD, longitude*DEGTORAD, h);
return true;
} else {
}
}
+//***************************************************************************
+//
+// Initialize: Assume all angles GIVEN IN RADIANS !!
+//
void FGState::Initialize(float U, float V, float W,
float phi, float tht, float psi,
longitude = Longitude;
h = H;
-// Change all angular measurements from degrees (as in config file) to radians
-
gamma = 0.0;
if (W != 0.0)
alpha = U*U > 0.0 ? atan2(W, U) : 0.0;
else
beta = 0.0;
- latitude *= M_PI / 180.0;
- longitude *= M_PI / 180.0;
- phi *= M_PI / 180.0;
- tht *= M_PI / 180.0;
- psi *= M_PI / 180.0;
-
FDMExec->GetTranslation()->SetUVW(U, V, W);
FDMExec->GetRotation()->SetEuler(phi, tht, psi);
FDMExec->GetTranslation()->SetABG(alpha, beta, gamma);
}
+void FGState::Initialize(FGInitialCondition *FGIC)
+{
+
+ float tht,psi,phi;
+ float U,V,W;
+
+ latitude = FGIC->GetLatitudeRadIC();
+ longitude = FGIC->GetLongitudeRadIC();
+ h = FGIC->GetAltitudeFtIC();
+ U = FGIC->GetUBodyFpsIC();
+ V = FGIC->GetVBodyFpsIC();
+ W = FGIC->GetWBodyFpsIC();
+ tht = FGIC->GetThetaRadIC();
+ phi = FGIC->GetPhiRadIC();
+ psi = FGIC->GetPsiRadIC();
+
+ Initialize(U, V, W, phi*DEGTORAD, tht*DEGTORAD, psi*DEGTORAD,
+ latitude*DEGTORAD, longitude*DEGTORAD, h);
+}
+
+
bool FGState::StoreData(string fname)
{
ofstream datafile(fname.c_str());
cout << "U: " << FDMExec->GetTranslation()->GetU() << endl;
cout << "V: " << FDMExec->GetTranslation()->GetV() << endl;
cout << "W: " << FDMExec->GetTranslation()->GetW() << endl;
- cout << "P: " << FDMExec->GetRotation()->GetP() << endl;
- cout << "Q: " << FDMExec->GetRotation()->GetQ() << endl;
- cout << "R: " << FDMExec->GetRotation()->GetR() << endl;
+ cout << "P: " << FDMExec->GetRotation()->GetP()*RADTODEG << endl;
+ cout << "Q: " << FDMExec->GetRotation()->GetQ()*RADTODEG << endl;
+ cout << "R: " << FDMExec->GetRotation()->GetR()*RADTODEG << endl;
cout << "L: " << FDMExec->GetAircraft()->GetL() << endl;
cout << "M: " << FDMExec->GetAircraft()->GetM() << endl;
cout << "N: " << FDMExec->GetAircraft()->GetN() << endl;
cout << "Vt: " << Vt << endl;
cout << "latitude: " << latitude << endl;
cout << "longitude: " << longitude << endl;
- cout << "alpha: " << FDMExec->GetTranslation()->Getalpha() << endl;
- cout << "beta: " << FDMExec->GetTranslation()->Getbeta() << endl;
- cout << "gamma: " << FDMExec->GetTranslation()->Getgamma() << endl;
- cout << "phi: " << FDMExec->GetRotation()->Getphi() << endl;
- cout << "tht: " << FDMExec->GetRotation()->Gettht() << endl;
- cout << "psi: " << FDMExec->GetRotation()->Getpsi() << endl;
- cout << "Pdot: " << FDMExec->GetRotation()->GetPdot() << endl;
- cout << "Qdot: " << FDMExec->GetRotation()->GetQdot() << endl;
- cout << "Rdot: " << FDMExec->GetRotation()->GetRdot() << endl;
+ cout << "alpha: " << FDMExec->GetTranslation()->Getalpha()*RADTODEG << endl;
+ cout << "beta: " << FDMExec->GetTranslation()->Getbeta()*RADTODEG << endl;
+ cout << "gamma: " << FDMExec->GetTranslation()->Getgamma()*RADTODEG << endl;
+ cout << "phi: " << FDMExec->GetRotation()->Getphi()*RADTODEG << endl;
+ cout << "tht: " << FDMExec->GetRotation()->Gettht()*RADTODEG << endl;
+ cout << "psi: " << FDMExec->GetRotation()->Getpsi()*RADTODEG << endl;
+ cout << "Pdot: " << FDMExec->GetRotation()->GetPdot()*RADTODEG << endl;
+ cout << "Qdot: " << FDMExec->GetRotation()->GetQdot()*RADTODEG << endl;
+ cout << "Rdot: " << FDMExec->GetRotation()->GetRdot()*RADTODEG << endl;
cout << "h: " << h << endl;
cout << "a: " << a << endl;
cout << "rho: " << FDMExec->GetAtmosphere()->Getrho() << endl;
#endif
#include "FGDefs.h"
+#include "FGInitialCondition.h"
/*******************************************************************************
DEFINES
*******************************************************************************/
class FGFDMExec;
+
class FGState
{
public:
bool Reset(string, string);
void Initialize(float, float, float, float, float, float, float, float, float);
+ void Initialize(FGInitialCondition *FGIC);
bool StoreData(string);
bool DumpData(string);
bool DisplayData(void);
HISTORY
--------------------------------------------------------------------------------
12/02/98 JSB Created
+ 7/23/99 TP Added data member and modified Run and PutState to calcuate
+ Mach number
********************************************************************************
COMMENTS, REFERENCES, and NOTES
qbar = 0.5*rho*Vt*Vt;
+ mach = Vt / State->Geta();
+
PutState();
} else {
}
Fz = Aircraft->GetFz();
Mass = Aircraft->GetMass();
- rho = Atmosphere->Getrho();
+ rho = Atmosphere->GetDensity();
phi = Rotation->Getphi();
tht = Rotation->Gettht();
{
State->SetVt(Vt);
State->Setqbar(qbar);
+ State->SetMach(mach);
}
private:
float U, V, W; // Body frame velocities owned by FGTranslation
float P, Q, R;
- float Vt, qbar;
+ float Vt, qbar, mach;
float Udot, Vdot, Wdot;
float lastUdot, lastVdot, lastWdot;
float phi, tht, psi;
+EXTRA_DIST = JSBSim.cpp Makefile.solo
+
noinst_LIBRARIES = libJSBsim.a
libJSBsim_a_SOURCES = FGAircraft.cpp FGAircraft.h \
FGDefs.h \
FGFCS.cpp FGFCS.h \
FGFDMExec.cpp FGFDMExec.h \
+ FGInitialCondition.cpp FGInitialCondition.h \
FGModel.cpp FGModel.h \
FGOutput.cpp FGOutput.h \
FGPosition.cpp FGPosition.h \
noinst_PROGRAMS = testJSBsim
-testJSBsim_SOURCES = FGMain.cpp
+testJSBsim_SOURCES = JSBSim.cpp
testJSBsim_LDADD = libJSBsim.a