Add USSA 1976 model

This commit is contained in:
Milo Priegnitz 2024-06-21 18:20:06 +02:00
commit f3f7677108
5 changed files with 238 additions and 0 deletions

4
README.md Normal file
View File

@ -0,0 +1,4 @@
# Models
Collecion of Altitude models.
- U.S. Standard Atmosphere 1976

View File

@ -0,0 +1,58 @@
#ifndef INC_MODELS_ALTITUDE_
#define INC_MODELS_ALTITUDE_
#include <sta/math/probability/randomVariable.hpp>
namespace sta
{
namespace models {
class BarometricFormula
{
/**
* @brief Get the altitude from the pressure using the barometric formula. The model is valid for altitudes up to leat 9 km.
*
* @param p Pressure in hPa
* @return Altitude in meters
*/
public:
/*
* @brief Constructor for the Barometric Formula
* @param pRef Reference pressure in hPa
* @param hRef Reference altitude in meters
*/
BarometricFormula(float pRef = 1013.25, float hRef = 0);
/*
* @brief Get the altitude from the pressure using the barometric formula
* @param p Pressure in hPa
* @return Altitude in meters
*/
float getAltitude(float p);
private:
float p0_;
};
namespace USSA1976
{
/**
* @brief Get the altitude from the pressure using the U.S. Standard Atmosphere 1976 model (https://ntrs.nasa.gov/search.jsp?R=19770009539)
* and the barometric formula for multiple layers. The model is valid for altitudes up to at leat 71 km. See https://en.wikipedia.org/wiki/Barometric_formula for more information.
*
* @param p Pressure in hPa
* @return Altitude in meters
*/
float getAltitude(float p);
/**
* @brief Get the altitude as Gaussian Random Variable from the pressure as Gaussian Random Variable using the U.S. Standard Atmosphere 1976 model (https://ntrs.nasa.gov/search.jsp?R=19770009539)
* and the barometric formula for multiple layers. The model is valid for altitudes up to at leat 71 km. See https://en.wikipedia.org/wiki/Barometric_formula for more information.
*
* @param p Pressure in hPa as a Gaussian Random Variable
* @return Altitude in meters as a Gaussian Random Variable
*/
sta::math::probability::grv getAltitude(sta::math::probability::grv p);
} // namespace USSA1976
} // namespace models
} // namespace sta
#endif //INC_MODELS_ALTITUDE_

View File

@ -0,0 +1,12 @@
#ifndef STA_MODELS_CONSTANTS_HPP
#define STA_MODELS_CONSTANTS_HPP
//Altitude Constants
#define STA_CONSTANT_R 8.3144598 //J/(mol*K)
#define STA_CONSTANT_G_0 9.80665 //m/s^2
#define STA_CONSTANT_M 0.028964425278793993 //kg/mol
#endif //STA_MODELS_CONSTANTS_HPP

11
library.json Normal file
View File

@ -0,0 +1,11 @@
{
"owner" : "sta",
"name": "sta-models",
"version": "0.1.0",
"dependencies": [
{
"url": "git@gitlab.com:sta-git/alpaka/sta-peak.git",
"ref": "main"
}
]
}

153
src/models/altitude.cpp Normal file
View File

@ -0,0 +1,153 @@
#include <sta/models/altitude.hpp>
#include <cmath>
#include <sta/models/constants.hpp>
namespace sta
{
namespace models {
BarometricFormula::BarometricFormula(float pRef, float hRef)
{
p0_ = pRef / std::pow(1 - hRef / 44330, 5.255);
}
float BarometricFormula::getAltitude(float p){
return 44330 * (1 - std::pow(p / p0_, 1/5.255));
}
float USSA1976::getAltitude(float p){
float T_b = 288.15;
float L_b = 0.0065;
float P_b = 1013.25;
float exp_b = 5.2561;
float H_b = 0;
if(p > 226.321){
T_b = 288.15;
L_b = 0.0065;
P_b = 1013.25;
exp_b = 5.2561;
H_b = 0;
}
else if(p > 54.7489){
T_b = 216.65;
L_b = 0;
P_b = 226.321;
exp_b = 0;
H_b = 11000;
}
else if(p > 8.6802){
T_b = 216.65;
L_b = -0.001;
P_b = 54.7489;
exp_b = -34.1626;
H_b = 20000;
}
else if(p > 1.1091){
T_b = 228.65;
L_b = -0.0028;
P_b = 8.6802;
exp_b = -12.2009;
H_b = 32000;
}
else if(p > 0.6699){
T_b = 270.65;
L_b = 0;
P_b = 1.1091;
exp_b = 0;
H_b = 47000;
}
else if(p > 0.0399){
T_b = 270.65;
L_b = 0.0028;
P_b = 0.6699;
exp_b = 12.2009;
H_b = 51000;
}
else{
T_b = 214.65;
L_b = 0.002;
P_b = 0.0399;
exp_b = 17.0816;
H_b = 71000;
}
if(L_b == 0){
return H_b - STA_CONSTANT_R*T_b*std::log(p/P_b) / (STA_CONSTANT_G_0* STA_CONSTANT_M);
}
return T_b/L_b * (-std::pow((p/P_b),(1/exp_b)) + 1) + H_b;
}
sta::math::probability::grv USSA1976::getAltitude(sta::math::probability::grv p){
float T_b = 288.15;
float L_b = 0.0065;
float P_b = 1013.25;
float exp_b = 5.2561;
float H_b = 0;
if(p.mean > 226.321){
T_b = 288.15;
L_b = 0.0065;
P_b = 1013.25;
exp_b = 5.2561;
H_b = 0;
}
else if(p.mean > 54.7489){
T_b = 216.65;
L_b = 0;
P_b = 226.321;
exp_b = 0;
H_b = 11000;
}
else if(p.mean > 8.6802){
T_b = 216.65;
L_b = -0.001;
P_b = 54.7489;
exp_b = -34.1626;
H_b = 20000;
}
else if(p.mean > 1.1091){
T_b = 228.65;
L_b = -0.0028;
P_b = 8.6802;
exp_b = -12.2009;
H_b = 32000;
}
else if(p.mean > 0.6699){
T_b = 270.65;
L_b = 0;
P_b = 1.1091;
exp_b = 0;
H_b = 47000;
}
else if(p.mean > 0.0399){
T_b = 270.65;
L_b = 0.0028;
P_b = 0.6699;
exp_b = 12.2009;
H_b = 51000;
}
else{
T_b = 214.65;
L_b = 0.002;
P_b = 0.0399;
exp_b = 17.0816;
H_b = 71000;
}
float h;
float variance_h;
if(L_b == 0){
h = H_b - STA_CONSTANT_R*T_b*std::log(p.mean/P_b) / (STA_CONSTANT_G_0* STA_CONSTANT_M);
variance_h = std::pow(( STA_CONSTANT_R * T_b / (STA_CONSTANT_G_0* STA_CONSTANT_M*p.mean)),2) * p.variance;
}else{
float inner = -std::pow((p.mean/P_b),(1/exp_b));
h = T_b/L_b * (inner + 1) + H_b;
variance_h = std::pow(inner*(T_b/L_b)*(1/exp_b) / p.mean ,2) * p.variance;
}
return {h, variance_h};
}
} // namespace models
} // namespace sta