fly/flightgear
1
0
Fork 0
Code Activity

The WeatherDatabase doesn't need the voronoi code anymore but uses

Browse source 
Dave Eberly's spherical interpolation code (found in the Lib/Math
directory). So it would be great if you could give him also a place
in the thanks file.  Changing the WeatherDatabse made actually a heavy
internal redesign necessary but no code outside the database is
affected (isn't code hiding great?).
This commit is contained in:
curt 1999-12-23 16:54:54 +00:00
parent 1963e006c3
commit 210e87ec3a
14 changed files with 517 additions and 330 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

18
src/WeatherCM/FGAirPressureItem.h
View file

@ -37,6 +37,10 @@ HISTORY
suggestion suggestion
19.10.1999 Christian Mayer change to use PLIB's sg instead of Point[2/3]D 19.10.1999 Christian Mayer change to use PLIB's sg instead of Point[2/3]D
and lots of wee code cleaning and lots of wee code cleaning
15.12.1999 Christian Mayer changed the air pressure calculation to a much
more realistic formula. But as I need for that
the temperature I moved the code to
FGPhysicalProperties
*****************************************************************************/ *****************************************************************************/
/****************************************************************************/ /****************************************************************************/
@ -62,25 +66,23 @@ FGAirPressureItem operator-(const FGAirPressureItem& arg);
/* CLASS DECLARATION */ /* CLASS DECLARATION */
/* NOTE: The value stored in 'value' is the air preasure that we'd have at */ /* NOTE: The value stored in 'value' is the air preasure that we'd have at */
/* an altitude of 0.0 The correct airpreasure at the stored altitude */ /* an altitude of 0.0 The correct airpreasure at the stored altitude */
/* gets calulated when getValue() is called. */ /* gets calulated in FGPhyiscalProperties as I need to know the */
/* temperatures at the different altitudes for that. */
/****************************************************************************/ /****************************************************************************/
class FGAirPressureItem class FGAirPressureItem
{ {
private: private:
WeatherPrecision value; WeatherPrecision value; //that's the airpressure at 0 metres
protected: protected:
public: public:
FGAirPressureItem(const WeatherPrecision v) {value = v;} FGAirPressureItem(const WeatherPrecision v) {value = v; }
FGAirPressureItem() {value = FG_WEATHER_DEFAULT_AIRPRESSURE;} FGAirPressureItem() {value = FG_WEATHER_DEFAULT_AIRPRESSURE;}
WeatherPrecision getValue(const WeatherPrecision& alt) const WeatherPrecision getValue(void) const
{ {
return (WeatherPrecision)((value / 101325.0) * return value;
(
1.01325e5 + alt * (-1.19459535223623e1 + alt * (5.50461110007561e-4 + alt * (-1.13574703113648e-8 + alt * 8.61601726143988e-14)))
));
}; };
FGAirPressureItem& operator*=(const WeatherPrecision arg); FGAirPressureItem& operator*=(const WeatherPrecision arg);

247
src/WeatherCM/FGLocalWeatherDatabase.cpp
View file

@ -38,6 +38,10 @@ HISTORY
suggestion suggestion
19.10.1999 Christian Mayer change to use PLIB's sg instead of Point[2/3]D 19.10.1999 Christian Mayer change to use PLIB's sg instead of Point[2/3]D
and lots of wee code cleaning and lots of wee code cleaning
14.12.1999 Christian Mayer Changed the internal structure to use Dave
Eberly's spherical interpolation code. This
stops our dependancy on the (ugly) voronoi
code and simplyfies the code structure a lot.
*****************************************************************************/ *****************************************************************************/
/****************************************************************************/ /****************************************************************************/
@ -49,7 +53,6 @@ HISTORY
#include <Aircraft/aircraft.hxx> #include <Aircraft/aircraft.hxx>
#include "FGLocalWeatherDatabase.h" #include "FGLocalWeatherDatabase.h"
#include "FGVoronoi.h"
#include "FGWeatherParse.h" #include "FGWeatherParse.h"
@ -57,37 +60,6 @@ HISTORY
/********************************** CODE ************************************/ /********************************** CODE ************************************/
/****************************************************************************/ /****************************************************************************/
/****************************************************************************/
/* return the index (better: ID) of the area with point p */
/****************************************************************************/
unsigned int FGLocalWeatherDatabase::AreaWith(const sgVec2& p) const
{
for (FGMicroWeatherList::size_type i = 0; i != WeatherAreas.size(); i++)
{
if (WeatherAreas[i].hasPoint(p) == true)
return i+1;
}
return 0; //nothing found
}
/****************************************************************************/
/* make tiles out of points on a 2D plane */
/****************************************************************************/
void FGLocalWeatherDatabase::tileLocalWeather(const FGPhysicalProperties2DVector& EntryList)
{
FGVoronoiInputList input;
for (FGPhysicalProperties2DVectorConstIt it1 = EntryList.begin(); it1 != EntryList.end(); it1++)
input.push_back(FGVoronoiInput(it1->p, *it1));
FGVoronoiOutputList output = Voronoiate(input);
for (FGVoronoiOutputList::iterator it2 = output.begin(); it2 != output.end(); it2++)
WeatherAreas.push_back(FGMicroWeather(it2->value, it2->boundary));
}
FGLocalWeatherDatabase* FGLocalWeatherDatabase::theFGLocalWeatherDatabase = 0; FGLocalWeatherDatabase* FGLocalWeatherDatabase::theFGLocalWeatherDatabase = 0;
FGLocalWeatherDatabase *WeatherDatabase; FGLocalWeatherDatabase *WeatherDatabase;
@ -98,7 +70,6 @@ void FGLocalWeatherDatabase::init(const WeatherPrecision visibility, const Datab
if (theFGLocalWeatherDatabase) if (theFGLocalWeatherDatabase)
{ {
//FG_LOG( FG_GENERAL, FG_ALERT, "Error: only one local weather allowed" );
cerr << "Error: only one local weather allowed"; cerr << "Error: only one local weather allowed";
exit(-1); exit(-1);
} }
@ -106,7 +77,7 @@ void FGLocalWeatherDatabase::init(const WeatherPrecision visibility, const Datab
setWeatherVisibility(visibility); setWeatherVisibility(visibility);
DatabaseStatus = type; DatabaseStatus = type;
global = 0; //just get sure... database = 0; //just get sure...
Thunderstorm = false; Thunderstorm = false;
//I don't need to set theThunderstorm as Thunderstorm == false //I don't need to set theThunderstorm as Thunderstorm == false
@ -115,33 +86,42 @@ void FGLocalWeatherDatabase::init(const WeatherPrecision visibility, const Datab
{ {
case use_global: case use_global:
{ {
global = new FGGlobalWeatherDatabase; //initialize GlobalDatabase cerr << "Error: there's no global database anymore!\n";
global->setDatabaseStatus(FGGlobalWeatherDatabase_working); exit(-1);
tileLocalWeather(global->getAll(last_known_position, WeatherVisibility, 3));
} }
break; break;
case use_internet: case use_internet:
{ {
FGWeatherParse parsed_data; FGWeatherParse *parsed_data = new FGWeatherParse();
parsed_data.input( "weather/current.gz" ); parsed_data->input( "weather/current.gz" );
global = new FGGlobalWeatherDatabase; //initialize GlobalDatabase unsigned int n = parsed_data->stored_stations();
global->setDatabaseStatus(FGGlobalWeatherDatabase_working);
sgVec2 *p = new sgVec2 [n];
FGPhysicalProperties *f = new FGPhysicalProperties[n];
// fill the database // fill the database
for (unsigned int i = 0; i != parsed_data.stored_stations(); i++) for (unsigned int i = 0; i < n; i++)
{ {
global->add( parsed_data.getFGPhysicalProperties2D(i) ); f[i] = parsed_data->getFGPhysicalProperties(i);
//cerr << parsed_data.getFGPhysicalProperties2D(i); parsed_data->getPosition(i, p[i]);
if ( (i%100) == 0) if ( (i%100) == 0)
cerr << "."; cerr << ".";
} }
cerr << "\n";
//and finally tile it // free the memory of the parsed data to ease the required memory
tileLocalWeather(global->getAll(last_known_position, WeatherVisibility, 3)); // for the very memory consuming spherical interpolation
delete parsed_data;
//and finally init the interpolation
cerr << "\nInitialiating Interpolation. (2-3 minutes on a PII-350)\n";
database = new SphereInterpolate<FGPhysicalProperties>(n, p, f);
//and free my allocations:
delete[] p;
delete[] f;
cerr << "Finished weather init.\n"; cerr << "Finished weather init.\n";
} }
@ -154,9 +134,11 @@ void FGLocalWeatherDatabase::init(const WeatherPrecision visibility, const Datab
case manual: case manual:
case default_mode: case default_mode:
{ {
double x[2] = {0.0, 0.0}; //make an standard weather that's the same at the whole world
vector<sgVec2Wrap> emptyList; double y[2] = {0.0, 0.0}; //make an standard weather that's the same at the whole world
WeatherAreas.push_back(FGMicroWeather(FGPhysicalProperties2D(), emptyList)); //in these cases I've only got one tile double z[2] = {1.0, -1.0}; //make an standard weather that's the same at the whole world
FGPhysicalProperties f[2]; //make an standard weather that's the same at the whole world
database = new SphereInterpolate<FGPhysicalProperties>(2,x,y,z,f);
} }
break; break;
@ -168,13 +150,7 @@ void FGLocalWeatherDatabase::init(const WeatherPrecision visibility, const Datab
FGLocalWeatherDatabase::~FGLocalWeatherDatabase() FGLocalWeatherDatabase::~FGLocalWeatherDatabase()
{ {
//Tidying up: //Tidying up:
delete database;
//delete every stored area
WeatherAreas.erase(WeatherAreas.begin(), WeatherAreas.end());
//delete global database if necessary
if (DatabaseStatus == use_global)
delete global;
} }
/****************************************************************************/ /****************************************************************************/
@ -182,16 +158,7 @@ FGLocalWeatherDatabase::~FGLocalWeatherDatabase()
/****************************************************************************/ /****************************************************************************/
void FGLocalWeatherDatabase::reset(const DatabaseWorkingType type) void FGLocalWeatherDatabase::reset(const DatabaseWorkingType type)
{ {
//delete global database if necessary cerr << "FGLocalWeatherDatabase::reset isn't supported yet\n";
if ((DatabaseStatus == use_global) && (type != use_global))
delete global;
DatabaseStatus = type;
if (DatabaseStatus == use_global)
tileLocalWeather(global->getAll(last_known_position, WeatherVisibility, 3));
//delete every stored area
WeatherAreas.erase(WeatherAreas.begin(), WeatherAreas.end());
} }
/****************************************************************************/ /****************************************************************************/
@ -199,45 +166,29 @@ void FGLocalWeatherDatabase::reset(const DatabaseWorkingType type)
/****************************************************************************/ /****************************************************************************/
void FGLocalWeatherDatabase::update(const WeatherPrecision dt) void FGLocalWeatherDatabase::update(const WeatherPrecision dt)
{ {
if (DatabaseStatus==use_global) //if (DatabaseStatus==use_global)
global->update(dt); // global->update(dt);
} }
void FGLocalWeatherDatabase::update(const sgVec3& p) //position has changed void FGLocalWeatherDatabase::update(const sgVec3& p) //position has changed
{ {
sgCopyVec3(last_known_position, p); sgCopyVec3(last_known_position, p);
if ( AreaWith(p) == 0 )
{ //I have moved out of my local area...
//now I should erase all my areas and get the new ones
//but that takes too long :-( I have to take care about that soon
}
//uncomment this when you are using the GlobalDatabase //uncomment this when you are using the GlobalDatabase
/* /*
cerr << "****\nupdate(p) inside\n"; cerr << "****\nupdate(p) inside\n";
cerr << "Parameter: " << p[0] << "/" << p[1] << "/" << p[2] << "\n"; cerr << "Parameter: " << p[0] << "/" << p[1] << "/" << p[2] << "\n";
sgVec2 p2d; sgVec2 p2d;
sgSetVec2( p2d, p[0], p[1] ); sgSetVec2( p2d, p[0], p[1] );
cerr << FGPhysicalProperties2D(global->get(p2d), p2d); cerr << FGPhysicalProperties2D(get(p2d), p2d);
cerr << "****\n"; cerr << "****\n";
*/ */
} }
void FGLocalWeatherDatabase::update(const sgVec3& p, const WeatherPrecision dt) //time and/or position has changed void FGLocalWeatherDatabase::update(const sgVec3& p, const WeatherPrecision dt) //time and/or position has changed
{ {
sgCopyVec3(last_known_position, p); sgCopyVec3(last_known_position, p);
if ( AreaWith(p) == 0 )
{ //I have moved out of my local area...
//now I should erase all my areas and get the new ones
//but that takes too long :-( I have to take care about that soon
}
if (DatabaseStatus==use_global)
global->update(dt);
} }
/****************************************************************************/ /****************************************************************************/
@ -245,154 +196,46 @@ void FGLocalWeatherDatabase::update(const sgVec3& p, const WeatherPrecision dt)
/****************************************************************************/ /****************************************************************************/
FGPhysicalProperty FGLocalWeatherDatabase::get(const sgVec3& p) const FGPhysicalProperty FGLocalWeatherDatabase::get(const sgVec3& p) const
{ {
unsigned int a = AreaWith(p); return FGPhysicalProperty(database->Evaluate(p), p[3]);
if (a != 0)
return WeatherAreas[a-1].get(p[3]);
else //point is outside => ask GlobalWeatherDatabase
return global->get(p);
} }
FGPhysicalProperties FGLocalWeatherDatabase::get(const sgVec2& p) const FGPhysicalProperties FGLocalWeatherDatabase::get(const sgVec2& p) const
{ {
sgVec3 temp; return database->Evaluate(p);
sgSetVec3(temp, p[0], p[1], 0.0);
unsigned int a = AreaWith(temp);
if (a != 0)
return WeatherAreas[a-1].get();
else //point is outside => ask GlobalWeatherDatabase
return global->get(p);
} }
WeatherPrecision FGLocalWeatherDatabase::getAirDensity(const sgVec3& p) const WeatherPrecision FGLocalWeatherDatabase::getAirDensity(const sgVec3& p) const
{ {
FGPhysicalProperty dummy; FGPhysicalProperty dummy(database->Evaluate(p), p[3]);
unsigned int a = AreaWith(p);
if (a != 0)
dummy = WeatherAreas[a-1].get(p[3]);
else //point is outside => ask GlobalWeatherDatabase
dummy = global->get(p);
return return
(dummy.AirPressure*FG_WEATHER_DEFAULT_AIRDENSITY*FG_WEATHER_DEFAULT_TEMPERATURE) / (dummy.AirPressure*FG_WEATHER_DEFAULT_AIRDENSITY*FG_WEATHER_DEFAULT_TEMPERATURE) /
(dummy.Temperature*FG_WEATHER_DEFAULT_AIRPRESSURE); (dummy.Temperature*FG_WEATHER_DEFAULT_AIRPRESSURE);
} }
/****************************************************************************/
/* Add a weather feature at the point p and surrounding area */
/****************************************************************************/
void FGLocalWeatherDatabase::addWind(const WeatherPrecision alt, const sgVec3& x, const sgVec2& p)
{
unsigned int a = AreaWith(p);
if (a != 0)
WeatherAreas[a-1].addWind(alt, x);
}
void FGLocalWeatherDatabase::addTurbulence(const WeatherPrecision alt, const sgVec3& x, const sgVec2& p)
{
unsigned int a = AreaWith(p);
if (a != 0)
WeatherAreas[a-1].addTurbulence(alt, x);
}
void FGLocalWeatherDatabase::addTemperature(const WeatherPrecision alt, const WeatherPrecision x, const sgVec2& p)
{
unsigned int a = AreaWith(p);
if (a != 0)
WeatherAreas[a-1].addTemperature(alt, x);
}
void FGLocalWeatherDatabase::addAirPressure(const WeatherPrecision alt, const WeatherPrecision x, const sgVec2& p)
{
unsigned int a = AreaWith(p);
if (a != 0)
WeatherAreas[a-1].addAirPressure(alt, x);
}
void FGLocalWeatherDatabase::addVaporPressure(const WeatherPrecision alt, const WeatherPrecision x, const sgVec2& p)
{
unsigned int a = AreaWith(p);
if (a != 0)
WeatherAreas[a-1].addVaporPressure(alt, x);
}
void FGLocalWeatherDatabase::addCloud(const WeatherPrecision alt, const FGCloudItem& x, const sgVec2& p)
{
unsigned int a = AreaWith(p);
if (a != 0)
WeatherAreas[a-1].addCloud(alt, x);
}
void FGLocalWeatherDatabase::setSnowRainIntensity(const WeatherPrecision x, const sgVec2& p) void FGLocalWeatherDatabase::setSnowRainIntensity(const WeatherPrecision x, const sgVec2& p)
{ {
unsigned int a = AreaWith(p); /* not supported yet */
if (a != 0)
WeatherAreas[a-1].setSnowRainIntensity(x);
} }
void FGLocalWeatherDatabase::setSnowRainType(const SnowRainType x, const sgVec2& p) void FGLocalWeatherDatabase::setSnowRainType(const SnowRainType x, const sgVec2& p)
{ {
unsigned int a = AreaWith(p); /* not supported yet */
if (a != 0)
WeatherAreas[a-1].setSnowRainType(x);
} }
void FGLocalWeatherDatabase::setLightningProbability(const WeatherPrecision x, const sgVec2& p) void FGLocalWeatherDatabase::setLightningProbability(const WeatherPrecision x, const sgVec2& p)
{ {
unsigned int a = AreaWith(p); /* not supported yet */
if (a != 0)
WeatherAreas[a-1].setLightningProbability(x);
}
void FGLocalWeatherDatabase::addProperties(const FGPhysicalProperties2D& x)
{
if (DatabaseStatus==use_global)
{
global->add(x);
//BAD, BAD, BAD thing I'm doing here: I'm adding to the global database a point that
//changes my voronoi diagram but I don't update it! instead I'm changing one local value
//that could be anywhere!!
//This only *might* work when the plane moves so far so fast that the diagram gets newly
//calculated soon...
unsigned int a = AreaWith(x.p);
if (a != 0)
WeatherAreas[a-1].setStoredWeather(x);
}
else
{
unsigned int a = AreaWith(x.p);
if (a != 0)
WeatherAreas[a-1].setStoredWeather(x);
}
} }
void FGLocalWeatherDatabase::setProperties(const FGPhysicalProperties2D& x) void FGLocalWeatherDatabase::setProperties(const FGPhysicalProperties2D& x)
{ {
if (DatabaseStatus==use_global) /* not supported yet */
{
global->change(x);
//BAD, BAD, BAD thing I'm doing here: I'm adding to the global database a point that
//changes my voronoi diagram but I don't update it! Instead I'm changing one local value
//that could be anywhere!!
//This only *might* work when the plane moves so far so fast that the diagram gets newly
//calculated soon...
unsigned int a = AreaWith(x.p);
if (a != 0)
WeatherAreas[a-1].setStoredWeather(x);
}
else
{
unsigned int a = AreaWith(x.p);
if (a != 0)
WeatherAreas[a-1].setStoredWeather(x);
}
} }
void fgUpdateWeatherDatabase(void) void fgUpdateWeatherDatabase(void)
{ {
//cerr << "FGLocalWeatherDatabase::update()\n";
sgVec3 position; sgVec3 position;
sgSetVec3(position, sgSetVec3(position,
current_aircraft.fdm_state->get_Latitude(), current_aircraft.fdm_state->get_Latitude(),

77
src/WeatherCM/FGLocalWeatherDatabase.h
View file

@ -38,6 +38,10 @@ HISTORY
suggestion suggestion
19.10.1999 Christian Mayer change to use PLIB's sg instead of Point[2/3]D 19.10.1999 Christian Mayer change to use PLIB's sg instead of Point[2/3]D
and lots of wee code cleaning and lots of wee code cleaning
14.12.1999 Christian Mayer Changed the internal structure to use Dave
Eberly's spherical interpolation code. This
stops our dependancy on the (ugly) voronoi
code and simplyfies the code structure a lot.
*****************************************************************************/ *****************************************************************************/
/****************************************************************************/ /****************************************************************************/
@ -49,21 +53,16 @@ HISTORY
/****************************************************************************/ /****************************************************************************/
/* INCLUDES */ /* INCLUDES */
/****************************************************************************/ /****************************************************************************/
//This is only here for smoother code change. In the end the WD should be clean
//of *any* OpenGL:
#ifdef HAVE_WINDOWS_H
# include <windows.h>
#endif
#include <GL/glut.h>
#include <XGL/xgl.h>
#include <vector> #include <vector>
#include "sg.h" #include <sg.h>
#include <Math/sphrintp.h>
#include "FGPhysicalProperties.h" #include "FGPhysicalProperties.h"
#include "FGGlobalWeatherDatabase.h" #include "FGPhysicalProperty.h"
#include "FGMicroWeather.h"
#include "FGWeatherFeature.h" #include "FGWeatherFeature.h"
#include "FGWeatherDefs.h" #include "FGWeatherDefs.h"
#include "FGThunderstorm.h" #include "FGThunderstorm.h"
@ -81,10 +80,7 @@ class FGLocalWeatherDatabase
{ {
private: private:
protected: protected:
FGGlobalWeatherDatabase *global; //point to the global database SphereInterpolate<FGPhysicalProperties> *database;
typedef vector<FGMicroWeather> FGMicroWeatherList;
typedef FGMicroWeatherList::iterator FGMicroWeatherListIt;
typedef vector<sgVec2> pointVector; typedef vector<sgVec2> pointVector;
typedef vector<pointVector> tileVector; typedef vector<pointVector> tileVector;
@ -92,10 +88,6 @@ protected:
/************************************************************************/ /************************************************************************/
/* make tiles out of points on a 2D plane */ /* make tiles out of points on a 2D plane */
/************************************************************************/ /************************************************************************/
void tileLocalWeather(const FGPhysicalProperties2DVector& EntryList);
FGMicroWeatherList WeatherAreas;
WeatherPrecision WeatherVisibility; //how far do I need to simulate the WeatherPrecision WeatherVisibility; //how far do I need to simulate the
//local weather? Unit: metres //local weather? Unit: metres
sgVec3 last_known_position; sgVec3 last_known_position;
@ -103,23 +95,11 @@ protected:
bool Thunderstorm; //is there a thunderstorm near by? bool Thunderstorm; //is there a thunderstorm near by?
FGThunderstorm *theThunderstorm; //pointer to the thunderstorm. FGThunderstorm *theThunderstorm; //pointer to the thunderstorm.
/************************************************************************/
/* return the index of the area with point p */
/************************************************************************/
unsigned int AreaWith(const sgVec2& p) const;
unsigned int AreaWith(const sgVec3& p) const
{
sgVec2 temp;
sgSetVec2(temp, p[0], p[1]);
return AreaWith(temp);
}
public: public:
static FGLocalWeatherDatabase *theFGLocalWeatherDatabase; static FGLocalWeatherDatabase *theFGLocalWeatherDatabase;
enum DatabaseWorkingType { enum DatabaseWorkingType {
use_global, //use global database for data use_global, //use global database for data !!obsolete!!
use_internet, //use the weather data that came from the internet use_internet, //use the weather data that came from the internet
manual, //use only user inputs manual, //use only user inputs
distant, //use distant information, e.g. like LAN when used in distant, //use distant information, e.g. like LAN when used in
@ -186,19 +166,12 @@ public:
/************************************************************************/ /************************************************************************/
/* Add a weather feature at the point p and surrounding area */ /* Add a weather feature at the point p and surrounding area */
/************************************************************************/ /************************************************************************/
// !! Adds aren't supported anymore !!
void addWind (const WeatherPrecision alt, const sgVec3& x, const sgVec2& p);
void addTurbulence (const WeatherPrecision alt, const sgVec3& x, const sgVec2& p);
void addTemperature (const WeatherPrecision alt, const WeatherPrecision x, const sgVec2& p);
void addAirPressure (const WeatherPrecision alt, const WeatherPrecision x, const sgVec2& p);
void addVaporPressure(const WeatherPrecision alt, const WeatherPrecision x, const sgVec2& p);
void addCloud (const WeatherPrecision alt, const FGCloudItem& x, const sgVec2& p);
void setSnowRainIntensity (const WeatherPrecision x, const sgVec2& p); void setSnowRainIntensity (const WeatherPrecision x, const sgVec2& p);
void setSnowRainType (const SnowRainType x, const sgVec2& p); void setSnowRainType (const SnowRainType x, const sgVec2& p);
void setLightningProbability(const WeatherPrecision x, const sgVec2& p); void setLightningProbability(const WeatherPrecision x, const sgVec2& p);
void addProperties(const FGPhysicalProperties2D& x); //add a property
void setProperties(const FGPhysicalProperties2D& x); //change a property void setProperties(const FGPhysicalProperties2D& x); //change a property
/************************************************************************/ /************************************************************************/
@ -231,34 +204,10 @@ void inline FGLocalWeatherDatabase::setWeatherVisibility(const WeatherPrecision
WeatherVisibility = visibility; WeatherVisibility = visibility;
else else
WeatherVisibility = MINIMUM_WEATHER_VISIBILITY; WeatherVisibility = MINIMUM_WEATHER_VISIBILITY;
#if 0
//This code doesn't belong here as this is the optical visibility and not
//the visibility of the weather database (that should be bigger...). The
//optical visibility should be calculated from the vapor pressure e.g.
//But for the sake of a smoother change from the old way to the new one...
GLfloat fog_exp_density;
GLfloat fog_exp2_density;
// for GL_FOG_EXP
fog_exp_density = -log(0.01 / WeatherVisibility);
// for GL_FOG_EXP2
fog_exp2_density = sqrt( -log(0.01) ) / WeatherVisibility;
// Set correct opengl fog density
xglFogf (GL_FOG_DENSITY, fog_exp2_density);
// FG_LOG( FG_INPUT, FG_DEBUG, "Fog density = " << w->fog_density );
//cerr << "FGLocalWeatherDatabase::setWeatherVisibility(" << visibility << "):\n";
//cerr << "Fog density = " << fog_exp_density << "\n";
#endif
} }
WeatherPrecision inline FGLocalWeatherDatabase::getWeatherVisibility(void) const WeatherPrecision inline FGLocalWeatherDatabase::getWeatherVisibility(void) const
{ {
//cerr << "FGLocalWeatherDatabase::getWeatherVisibility() = " << WeatherVisibility << "\n";
return WeatherVisibility; return WeatherVisibility;
} }

154
src/WeatherCM/FGPhysicalProperties.cpp
View file

@ -71,7 +71,7 @@ FGPhysicalProperties::FGPhysicalProperties()
AirPressure = FGAirPressureItem(101325.0); AirPressure = FGAirPressureItem(101325.0);
VaporPressure[ 0.0] = FG_WEATHER_DEFAULT_VAPORPRESSURE; //in Pa (I *only* accept SI!) VaporPressure[-1000.0] = FG_WEATHER_DEFAULT_VAPORPRESSURE; //in Pa (I *only* accept SI!)
VaporPressure[10000.0] = FG_WEATHER_DEFAULT_VAPORPRESSURE; //in Pa (I *only* accept SI!) VaporPressure[10000.0] = FG_WEATHER_DEFAULT_VAPORPRESSURE; //in Pa (I *only* accept SI!)
//Clouds.insert(FGCloudItem()) => none //Clouds.insert(FGCloudItem()) => none
@ -126,7 +126,7 @@ ostream& operator<< ( ostream& out, const FGPhysicalProperties2D& p )
out << "\n"; out << "\n";
out << "Stored AirPressure: "; out << "Stored AirPressure: ";
out << p.AirPressure.getValue(0)/100.0 << " hPa at " << 0.0 << "m; "; out << p.AirPressure.getValue()/100.0 << " hPa at " << 0.0 << "m; ";
out << "\n"; out << "\n";
out << "Stored VaporPressure: "; out << "Stored VaporPressure: ";
@ -140,4 +140,154 @@ ostream& operator<< ( ostream& out, const FGPhysicalProperties2D& p )
} }
inline double F(const WeatherPrecision factor, const WeatherPrecision a, const WeatherPrecision b, const WeatherPrecision r, const WeatherPrecision x)
{
const double c = 1.0 / (-b + a * r);
return factor * c * ( 1.0 / (r + x) + a * c * log(abs((r + x) * (b + a * x))) );
}
WeatherPrecision FGPhysicalProperties::AirPressureAt(const WeatherPrecision x) const
{
const double rho0 = (AirPressure.getValue()*FG_WEATHER_DEFAULT_AIRDENSITY*FG_WEATHER_DEFAULT_TEMPERATURE)/(TemperatureAt(0)*FG_WEATHER_DEFAULT_AIRPRESSURE);
const double G = 6.673e-11; //Gravity; in m^3 kg^-1 s^-2
const double m = 5.977e24; //mass of the earth in kg
const double r = 6368e3; //radius of the earth in metres
const double factor = -(rho0 * TemperatureAt(0) * G * m) / AirPressure.getValue();
double a, b, FF = 0.0;
//ok, integrate from 0 to a now.
if (Temperature.size() < 2)
{ //take care of the case that there aren't enough points
//actually this should be impossible...
if (Temperature.size() == 0)
{
cerr << "ERROR in FGPhysicalProperties: Air pressure at " << x << " metres altiude requested,\n";
cerr << " but there isn't enough data stored! No temperature is aviable!\n";
return FG_WEATHER_DEFAULT_AIRPRESSURE;
}
//ok, I've got only one point. So I'm assuming that that temperature is
//the same for all altitudes.
a = 1;
b = TemperatureAt(0);
FF += F(factor, a, b, r, x );
FF -= F(factor, a, b, r, 0.0);
}
else
{ //I've got at least two entries now
//integrate 'backwards' by integrating the strip ]n,x] first, then ]n-1,n] ... to [0,n-m]
if (x>=0.0)
{
map<WeatherPrecision, WeatherPrecision>::const_iterator temp2 = Temperature.upper_bound(x);
map<WeatherPrecision, WeatherPrecision>::const_iterator temp1 = temp2; temp1--;
if (temp1->first == x)
{ //ignore that interval
temp1--; temp2--;
}
bool first_pass = true;
while(true)
{
if (temp2 == Temperature.end())
{
//temp2 doesn't exist. So cheat by assuming that the slope is the
//same as between the two earlier temperatures
temp1--; temp2--;
a = (temp2->second - temp1->second)/(temp2->first - temp1->first);
b = temp1->second - a * temp1->first;
temp1++; temp2++;
}
else
{
a = (temp2->second - temp1->second)/(temp2->first - temp1->first);
b = temp1->second - a * temp1->first;
}
if (first_pass)
{
FF += F(factor, a, b, r, x);
first_pass = false;
}
else
{
FF += F(factor, a, b, r, temp2->first);
}
if (temp1->first>0.0)
{
FF -= F(factor, a, b, r, temp1->first);
temp1--; temp2--;
}
else
{
FF -= F(factor, a, b, r, 0.0);
return AirPressure.getValue() * exp(FF);
}
}
}
else
{ //ok x is smaller than 0.0, so do everything in reverse
map<WeatherPrecision, WeatherPrecision>::const_iterator temp2 = Temperature.upper_bound(x);
map<WeatherPrecision, WeatherPrecision>::const_iterator temp1 = temp2; temp1--;
bool first_pass = true;
while(true)
{
if (temp2 == Temperature.begin())
{
//temp1 doesn't exist. So cheat by assuming that the slope is the
//same as between the two earlier temperatures
temp1 = Temperature.begin(); temp2++;
a = (temp2->second - temp1->second)/(temp2->first - temp1->first);
b = temp1->second - a * temp1->first;
temp2--;
}
else
{
a = (temp2->second - temp1->second)/(temp2->first - temp1->first);
b = temp1->second - a * temp1->first;
}
if (first_pass)
{
FF += F(factor, a, b, r, x);
first_pass = false;
}
else
{
FF += F(factor, a, b, r, temp2->first);
}
if (temp2->first<0.0)
{
FF -= F(factor, a, b, r, temp1->first);
if (temp2 == Temperature.begin())
{
temp1 = Temperature.begin(); temp2++;
}
else
{
temp1++; temp2++;
}
}
else
{
FF -= F(factor, a, b, r, 0.0);
return AirPressure.getValue() * exp(FF);
}
}
}
}
return AirPressure.getValue() * exp(FF);
}

19
src/WeatherCM/FGPhysicalProperties.h
View file

@ -38,6 +38,10 @@ HISTORY
suggestion suggestion
19.10.1999 Christian Mayer change to use PLIB's sg instead of Point[2/3]D 19.10.1999 Christian Mayer change to use PLIB's sg instead of Point[2/3]D
and lots of wee code cleaning and lots of wee code cleaning
15.12.1999 Christian Mayer changed the air pressure calculation to a much
more realistic formula. But as I need for that
the temperature I moved the code to
FGPhysicalProperties
*****************************************************************************/ *****************************************************************************/
/****************************************************************************/ /****************************************************************************/
@ -63,7 +67,7 @@ HISTORY
#include <vector> #include <vector>
#include <map> #include <map>
#include "sg.h" #include <sg.h>
#include "FGWeatherDefs.h" #include "FGWeatherDefs.h"
@ -107,7 +111,7 @@ public:
void WindAt (sgVec3 ret, const WeatherPrecision a) const; void WindAt (sgVec3 ret, const WeatherPrecision a) const;
void TurbulenceAt (sgVec3 ret, const WeatherPrecision a) const; void TurbulenceAt (sgVec3 ret, const WeatherPrecision a) const;
WeatherPrecision TemperatureAt (const WeatherPrecision a) const; WeatherPrecision TemperatureAt (const WeatherPrecision a) const;
WeatherPrecision AirPressureAt (const WeatherPrecision a) const; WeatherPrecision AirPressureAt (const WeatherPrecision x) const; //x is used here instead of a on purpose
WeatherPrecision VaporPressureAt(const WeatherPrecision a) const; WeatherPrecision VaporPressureAt(const WeatherPrecision a) const;
//for easier access to the cloud stuff: //for easier access to the cloud stuff:
@ -214,7 +218,7 @@ inline FGPhysicalProperties& FGPhysicalProperties::operator += (const FGPhysical
if (!Temperature.insert(*TemperatureIt).second) if (!Temperature.insert(*TemperatureIt).second)
Temperature[TemperatureIt->first] += TemperatureIt->second; Temperature[TemperatureIt->first] += TemperatureIt->second;
AirPressure += p.AirPressure.getValue(0.0); AirPressure += p.AirPressure.getValue();
for (scalar_iterator VaporPressureIt = p.VaporPressure.begin(); for (scalar_iterator VaporPressureIt = p.VaporPressure.begin();
VaporPressureIt != p.VaporPressure.end(); VaporPressureIt != p.VaporPressure.end();
@ -249,7 +253,7 @@ inline FGPhysicalProperties& FGPhysicalProperties::operator -= (const FGPhysical
if (!Temperature.insert( make_pair(TemperatureIt->first, -TemperatureIt->second) ).second) if (!Temperature.insert( make_pair(TemperatureIt->first, -TemperatureIt->second) ).second)
Temperature[TemperatureIt->first] -= TemperatureIt->second; Temperature[TemperatureIt->first] -= TemperatureIt->second;
AirPressure -= p.AirPressure.getValue(0.0); AirPressure -= p.AirPressure.getValue();
for (scalar_iterator VaporPressureIt = p.VaporPressure.begin(); for (scalar_iterator VaporPressureIt = p.VaporPressure.begin();
VaporPressureIt != p.VaporPressure.end(); VaporPressureIt != p.VaporPressure.end();
@ -261,7 +265,6 @@ inline FGPhysicalProperties& FGPhysicalProperties::operator -= (const FGPhysical
return *this; return *this;
} }
inline void FGPhysicalProperties::WindAt(sgVec3 ret, const WeatherPrecision a) const inline void FGPhysicalProperties::WindAt(sgVec3 ret, const WeatherPrecision a) const
{ {
typedef map<FGPhysicalProperties::Altitude, FGWindItem>::const_iterator vector_iterator; typedef map<FGPhysicalProperties::Altitude, FGWindItem>::const_iterator vector_iterator;
@ -302,10 +305,8 @@ inline WeatherPrecision FGPhysicalProperties::TemperatureAt(const WeatherPrecisi
return ( (it2->second - it->second)/(it2->first - it->first) ) * (a - it2->first) + it2->second; return ( (it2->second - it->second)/(it2->first - it->first) ) * (a - it2->first) + it2->second;
} }
inline WeatherPrecision FGPhysicalProperties::AirPressureAt(const WeatherPrecision a) const //inline WeatherPrecision FGPhysicalProperties::AirPressureAt(const WeatherPrecision x) const
{ //moved to FGPhysicalProperties.cpp as it got too complex to inline
return AirPressure.getValue(a);
}
inline WeatherPrecision FGPhysicalProperties::VaporPressureAt(const WeatherPrecision a) const inline WeatherPrecision FGPhysicalProperties::VaporPressureAt(const WeatherPrecision a) const
{ {

2
src/WeatherCM/FGPhysicalProperty.h
View file

@ -62,7 +62,7 @@ HISTORY
#include <vector> #include <vector>
#include "sg.h" #include <sg.h>
#include "FGWeatherDefs.h" #include "FGWeatherDefs.h"
#include "FGPhysicalProperties.h" #include "FGPhysicalProperties.h"

2
src/WeatherCM/FGTurbulenceItem.h
View file

@ -58,7 +58,7 @@ HISTORY
# include <windows.h> # include <windows.h>
#endif #endif
#include "sg.h" #include <sg.h>
#include "FGWeatherDefs.h" #include "FGWeatherDefs.h"

2
src/WeatherCM/FGWeatherFeature.h
View file

@ -51,7 +51,7 @@ HISTORY
/****************************************************************************/ /****************************************************************************/
#include <Include/compiler.h> #include <Include/compiler.h>
#include "sg.h" #include <sg.h>
#include "FGWeatherDefs.h" #include "FGWeatherDefs.h"

47
src/WeatherCM/FGWeatherParse.cpp
View file

@ -50,12 +50,13 @@ You can also visit his homepage at http://www.wetterzentrale.de
HISTORY HISTORY
------------------------------------------------------------------------------ ------------------------------------------------------------------------------
18.10.1999 Christian Mayer Created 18.10.1999 Christian Mayer Created
14.12.1999 Christian Mayer minor internal changes
*****************************************************************************/ *****************************************************************************/
/****************************************************************************/ /****************************************************************************/
/* INCLUDES */ /* INCLUDES */
/****************************************************************************/ /****************************************************************************/
#include "Include/fg_constants.h" #include <Include/fg_constants.h>
#include "FGWeatherParse.h" #include "FGWeatherParse.h"
#include "FGWeatherUtils.h" #include "FGWeatherUtils.h"
@ -82,6 +83,14 @@ void FGWeatherParse::input(const char *file)
in.open( file ); in.open( file );
if (! in.is_open() )
{
cerr << "Couldn't find that file!\nExiting...\n";
exit(-1);
}
else
{
bool skip = false;
while ( in ) while ( in )
{ {
entry temp; entry temp;
@ -103,26 +112,45 @@ void FGWeatherParse::input(const char *file)
in >> temp.airpressure_history[2]; in >> temp.airpressure_history[2];
in >> temp.airpressure_history[3]; in >> temp.airpressure_history[3];
for (int i = 0; i < weather_station.size(); i++)
{
if ((weather_station[i].lat == temp.lat) && (weather_station[i].lon == temp.lon))
{
// Two weatherstations are at the same positon
// => averageing both
// just taking care of the stuff that metters for us
weather_station[i].x_wind += temp.x_wind; weather_station[i].x_wind *= 0.5;
weather_station[i].y_wind += temp.y_wind; weather_station[i].y_wind *= 0.5;
weather_station[i].temperature += temp.temperature; weather_station[i].temperature *= 0.5;
weather_station[i].dewpoint += temp.dewpoint; weather_station[i].dewpoint *= 0.5;
weather_station[i].airpressure += temp.airpressure; weather_station[i].airpressure *= 0.5;
skip = true;
}
}
if (skip == false)
weather_station.push_back( temp ); weather_station.push_back( temp );
skip = false;
// output a point to ease the waiting // output a point to ease the waiting
if ( ((nr++)%100) == 0 ) if ( ((nr++)%100) == 0 )
cerr << "."; cerr << ".";
} }
cerr << "\n" << nr << " stations read\n"; cerr << "\n" << nr << " stations read\n";
}
} }
FGPhysicalProperties2D FGWeatherParse::getFGPhysicalProperties2D(const unsigned int nr) const FGPhysicalProperties FGWeatherParse::getFGPhysicalProperties(const unsigned int nr) const
{ {
FGPhysicalProperties2D ret_val; FGPhysicalProperties ret_val;
//chache this entry //chache this entry
entry this_entry = weather_station[nr]; entry this_entry = weather_station[nr];
//set the position of the station
sgSetVec2( ret_val.p, this_entry.lat * DEG_TO_RAD, this_entry.lon * DEG_TO_RAD );
ret_val.Wind[-1000.0] = FGWindItem(this_entry.x_wind, this_entry.y_wind, 0.0); ret_val.Wind[-1000.0] = FGWindItem(this_entry.x_wind, this_entry.y_wind, 0.0);
ret_val.Wind[10000.0] = FGWindItem(this_entry.x_wind, this_entry.y_wind, 0.0); ret_val.Wind[10000.0] = FGWindItem(this_entry.x_wind, this_entry.y_wind, 0.0);
ret_val.Temperature[0.0] = Celsius( this_entry.temperature ); ret_val.Temperature[0.0] = Celsius( this_entry.temperature );
@ -138,6 +166,9 @@ FGPhysicalProperties2D FGWeatherParse::getFGPhysicalProperties2D(const unsigned
return ret_val; return ret_val;
} }
void FGWeatherParse::getPosition(const unsigned int nr, sgVec2 pos) const
{
//set the position of the station
sgSetVec2( pos, weather_station[nr].lat * DEG_TO_RAD, weather_station[nr].lon * DEG_TO_RAD );
}

4
src/WeatherCM/FGWeatherParse.h
View file

@ -49,6 +49,7 @@ You can also visit his homepage at http://www.wetterzentrale.de
HISTORY HISTORY
------------------------------------------------------------------------------ ------------------------------------------------------------------------------
18.10.1999 Christian Mayer Created 18.10.1999 Christian Mayer Created
14.12.1999 Christian Mayer minor internal changes
*****************************************************************************/ *****************************************************************************/
/****************************************************************************/ /****************************************************************************/
@ -121,7 +122,8 @@ public:
return weather_station[nr]; return weather_station[nr];
} }
FGPhysicalProperties2D getFGPhysicalProperties2D(const unsigned int nr) const; FGPhysicalProperties getFGPhysicalProperties(const unsigned int nr) const;
void getPosition(const unsigned int nr, sgVec2 pos) const;
}; };
/****************************************************************************/ /****************************************************************************/

2
src/WeatherCM/FGWeatherVectorWrap.h
View file

@ -43,7 +43,7 @@ HISTORY
/****************************************************************************/ /****************************************************************************/
/* INCLUDES */ /* INCLUDES */
/****************************************************************************/ /****************************************************************************/
#include "sg.h" #include <sg.h>
/****************************************************************************/ /****************************************************************************/
/* DEFINES */ /* DEFINES */

2
src/WeatherCM/FGWindItem.h
View file

@ -58,7 +58,7 @@ HISTORY
# include <windows.h> # include <windows.h>
#endif #endif
#include "sg.h" #include <sg.h>
#include "FGWeatherDefs.h" #include "FGWeatherDefs.h"

3
src/WeatherCM/Makefile.am
View file

@ -3,9 +3,7 @@ noinst_LIBRARIES = libWeatherCM.a
libWeatherCM_a_SOURCES = \ libWeatherCM_a_SOURCES = \
FGAirPressureItem.cpp FGAirPressureItem.h \ FGAirPressureItem.cpp FGAirPressureItem.h \
FGCloud.h FGCloudItem.cpp FGCloudItem.h \ FGCloud.h FGCloudItem.cpp FGCloudItem.h \
FGGlobalWeatherDatabase.cpp FGGlobalWeatherDatabase.h \
FGLocalWeatherDatabase.cpp FGLocalWeatherDatabase.h \ FGLocalWeatherDatabase.cpp FGLocalWeatherDatabase.h \
FGMicroWeather.cpp FGMicroWeather.h \
FGPhysicalProperties.cpp FGPhysicalProperties.h \ FGPhysicalProperties.cpp FGPhysicalProperties.h \
FGPhysicalProperty.cpp FGPhysicalProperty.h \ FGPhysicalProperty.cpp FGPhysicalProperty.h \
FGSnowRain.h \ FGSnowRain.h \
@ -13,7 +11,6 @@ libWeatherCM_a_SOURCES = \
FGThunderstorm.cpp FGThunderstorm.h \ FGThunderstorm.cpp FGThunderstorm.h \
FGTurbulenceItem.cpp FGTurbulenceItem.h \ FGTurbulenceItem.cpp FGTurbulenceItem.h \
FGVaporPressureItem.cpp FGVaporPressureItem.h \ FGVaporPressureItem.cpp FGVaporPressureItem.h \
FGVoronoi.cpp FGVoronoi.h \
FGWeatherDefs.h FGWeatherFeature.h FGWeatherUtils.h \ FGWeatherDefs.h FGWeatherFeature.h FGWeatherUtils.h \
FGWeatherParse.cpp FGWeatherParse.h \ FGWeatherParse.cpp FGWeatherParse.h \
FGWeatherVectorWrap.h \ FGWeatherVectorWrap.h \

212
src/WeatherCM/air-pressure-explanation.txt Normal file
View file

@ -0,0 +1,212 @@
The formula p(x) for calculating the air pressure at a given altitude
---------------------------------------------------------------------
Well known is the baromertic(?) formula
rho0
------ * g * x
p0
p(x) = p0 * e
with p0 being the airpressure and rho0 being the air density at an altitude of
0 metres above sea level and g being the gravity constant of 9.81 m/sq. sec
This formula can easily be derivated when you know, that:
* the pressure difference is
dp = - rho * g * dx
* Boyle-Mariotte says:
p0 : p = rho0 : rho
Combinig the terms and changing them around I get:
dp [ rho0 ]
-- = - rho * g = - [ ------ * p(x) ] * g
dx [ p0 ]
rho0
p'(x) = - ------ * p(x) * g
p0
Solving that differential equation and knowing that p(0) = p0 I get:
rho0
- ------ * g * x
p0
p(x) = p0 * e
q.e.d.
-------------------------------------------------------------------------------
The problem with that equation is that it doesn't take different temperatures
at different altitudes into account. And the inaccuracies due to it are huge.
That's why this formula is only used in very low altitudes.
So to get a usefull formula for FlightGear I need to extend it. And as I'm
already 'recreating' that formula I'm taking the change in g into account, too.
This doesn't make such a dramatic difference to the result as the inclusion of
temperature change does, but it doesn't complicate the final formula too much.
So I get three formulas that I'm combining in the end:
* the change of g with the altitude:
G * m
g(x) = -----------
(x + r)^2
G is the universal gravity constant(?) and is 6.673e-11 m^3 kg^-1 s^-2
m is the mass of the earth and is 5.977e24 kg
r is the radius of the earth and is 6368 km
* The pressure difference stays the same:
dp = - rho * g(x) * dx
* If I combine Boyle-Mariotte with Gay-Lussac I get:
rho0 * T0 p
rho = ----------- * ---
p0 T
Combining the terms again I get this time:
dp [ rho0 * T0 p(x) ]
-- = - rho * g(x) = - [ ----------- * ------ ] * g(x)
dx [ p0 T(x) ]
rho0 * T0 p(x) * g(x)
p'(x) = - ----------- * -------------
p0 T(x)
This DE isn't that easy to solve as the one above, it by looking into the right
books you'll see the general solution for:
y' + f(x)*y = 0
is
x
/\
- | f(x) dx
\/
n
y = m * e
and P(m,n) will be a point on the graph.
For q = n = 0 metres altitude we get y = m. As y is p(x) we know that m has to
be p0.
So our final formuala is
ho0 * T0 g(x)
f1(x) = ----------- * ------
p0 T(x)
x x
/\ /\
- | f1(x) dx | f(x) dx
\/ \/
0 0 F(x) - F(0)
p(x) = p0 * e = p0 * e = p0 * e
The only disturbing thing we've got left is the integral. Luckily there is a
great service at http://integrals.wolfram.com/ that helps me doing it :-)
But the f(x) is still too general so I'm substituting:
rho0 * T0 * G * m
f(x) = - -----------------------
p0 * (x + r)^2 * T(x)
but even that isn't good enough. But as I'm linearily interpolating between
two different temperatures I can say that T(x) = a*x + b for the x inbetween
two different stored temperatures. So I just need to integrate every pice
independandly. But anyway, I get:
rho0 * T0 * G * m
f(x) = - ------------------------------
p0 * (x + r)^2 * (a * x + b)
Integrating that I get:
rho0 * T0 * G * m [ 1
F(x) = - ------------------- * [ ------------------------ -
p0 [ (-b + a * r) * (r + x)
a * log|r + x| a * log|b + a * x| ]
---------------- + -------------------- ]
(b - a * r)^2 (b - a * r)^2 ]
To lower the computional cost I can transfere the equation.
* I'm defining
rho0 * T0 * G * m
factor = - -------------------
p0
1
c = --------------
(-b + a * r)
* now I can write
[ c ]
F(x) = factor * [ --------- - a * c * c * [log|r + x| + log|b + a * x|] ]
[ (r + x) ]
* and simplyfy it to
[ 1 ]
F(x) = factor * c * [ --------- - a * c * log|(r + x) * (b + a * x)| ]
[ (r + x) ]
-------------------------------------------------------------------------------
The following table shows quite nicely how accurate my formula is:
Altitude[m] | Airpressure [hPa] | Error [%]
| Official | My formula |
------------+---------------+---------------+---------------
-200 | 1037.51 | 1037.24 | 0.0260
-100 | 1025.32 | 1025.19 | 0.0127
0 | 1013.25 | 1013.25 | 0.0
500 | 954.59 | 955.224 | 0.0664
1000 | 898.70 | 899.912 | 0.1349
2000 | 794.88 | 797.042 | 0.2720
3000 | 700.99 | 703.885 | 0.4130
4000 | 616.28 | 619.727 | 0.5593
5000 | 540.07 | 543.89 | 0.7073
6000 | 471.67 | 475.731 | 0.8610
7000 | 410.46 | 414.643 | 1.0191
8000 | 355.84 | 360.054 | 1.1842
9000 | 307.27 | 311.422 | 1.3513
10000 | 264.21 | 268.238 | 1.5245
20000 | 54.670/55.3 | 55.7971 | 2.0616/0.8989
30000 | 11.8 | 11.3149 | 1.5441
40000 | 3.0 | 2.74665 | 18.9703
50000 | 0.88 | 0.753043 | 41.9183
60000 | 0.257 | 0.221907 | 57.9802
70000 | 0.0602 | 0.0530785 | 61.9153
80000 | 0.0101 | 0.00905461 | 51.5725
100000 | 2.14e-4 | 2.03984e-4 | 5.5131
The official values are from the CINA atmosphere which assumes a air pressure
of 1013.25 hPa and a temperature of 15 degC at sea level and a temperature
gradient of -6.5 deg/km. The CINA atmosphere gives only values for altiudes
up to 20 km. The values for 20 km and above are from the 1959 ARDC atmosphere.
That's why I've got two values at 20000 metres.
The temperature changes dramtically in the altitudes over 20 km which I didn't
take care of in my calculations which explains the huge errors at that altitude
range. But you can see nicely that the values are at least quite close to the
official values.
Using a better temperature model for the altitudes above 20 km should
dramatically increase the accuracy there.