MI - Fimex
Functions
vertical_coordinate_transformations.h File Reference
#include "fimex/mifi_constants.h"

Go to the source code of this file.

Functions

int mifi_atmosphere_ln_pressure (size_t n, double p0, const double *lev, double *pressure)
 
int mifi_atmosphere_sigma_pressure (size_t n, double ptop, double ps, const double *level, double *pressure)
 
int mifi_atmosphere_pressure_sigma (size_t n, double ptop, double ps, const double *pressure, double *level)
 
int mifi_atmosphere_hybrid_sigma_pressure (size_t n, double p0, double ps, const double *a, const double *b, double *pressure)
 
int mifi_atmosphere_hybrid_sigma_ap_pressure (size_t n, double ps, const double *ap, const double *b, double *pressure)
 
int mifi_barometric_pressure (size_t n, double P_b, const double *h, double T_b, double *pressure)
 
int mifi_barometric_standard_pressure (size_t n, const double *h, double *pressure)
 
int mifi_barometric_height (size_t n, double P_b, const double *p, double T_b, double *height)
 
int mifi_barometric_standard_altitude (size_t n, const double *p, double *altitude)
 
float mifi_virtual_temperature (float spec_humidity, float t)
 
float mifi_relative_to_specific_humidity (float rh, float t, float p)
 
float mifi_specific_to_relative_humidity (float sh, float t, float p)
 
float mifi_dewpoint_to_relative_humidity (float dew, float t)
 
float mifi_barometric_layer_thickness (float p_low_alti, float p_high_alti, float T)
 
int mifi_ocean_s_g1_z (size_t n, double h, double h_c, double zeta, const double *s, const double *C, double *z)
 
int mifi_ocean_s_g2_z (size_t n, double h, double h_c, double zeta, const double *s, const double *C, double *z)
 
int mifi_omega_to_vertical_wind (size_t n, const double *omega, const double *p, const double *t, double *w)
 
int mifi_vertical_wind_to_omega (size_t n, const double *w, const double *p, const double *t, double *omega)
 

Function Documentation

int mifi_atmosphere_hybrid_sigma_ap_pressure ( size_t  n,
double  ps,
const double *  ap,
const double *  b,
double *  pressure 
)

convert a standard_name="atmosphere_hybrid_sigma_pressure_coordinate" to pressure using the formula p(k) = ap(k) + b(k)*ps

This is the same as mifi_atmosphere_hybrid_sigma_pressure(), but with the reference pressure and a joined already. Choice depends on the model, i.e. available input values.

Parameters
nsize of arrays ap, b and pressure
pssurface pressure - usually varying in time,x,y
appressure level values
bdimensionless level values
pressureoutput values in the same unit as p0 and ps and at the same place as ps
Returns
MIFI_OK on success or MIFI_ERROR on failure
int mifi_atmosphere_hybrid_sigma_pressure ( size_t  n,
double  p0,
double  ps,
const double *  a,
const double *  b,
double *  pressure 
)

convert a standard_name="atmosphere_hybrid_sigma_pressure_coordinate" to pressure using the formula p(k) = a(k)*p0 + b(k)*ps

Parameters
nsize of arrays a, b and pressure
p0reference pressure
pssurface pressure - usually varying in time,x,y
adimensionless level values
bdimensionless level values
pressureoutput values in the same unit as p0 and ps and at the same place as ps
Returns
MIFI_OK on success or MIFI_ERROR on failure
int mifi_atmosphere_ln_pressure ( size_t  n,
double  p0,
const double *  lev,
double *  pressure 
)

convert a standard_name="atmosphere_ln_pressure_coordinate" to pressure using the formula p(k) = p0 * exp(-lev(k))

Parameters
nsize of arrays lev and pressure
p0base pressure
levlevel values
pressureoutput values in the same unit as p0
Returns
MIFI_OK on success or MIFI_ERROR on failure
int mifi_atmosphere_pressure_sigma ( size_t  n,
double  ptop,
double  ps,
const double *  pressure,
double *  level 
)

convert a pressure to standard_name="atmosphere_sigma_coordinate" using the formula p(k) = ptop + sigma(k)*(ps-ptop) -> sigma(k) = (p(k) - ptop)/(ps-ptop)

Parameters
nsize of arrays sigma and pressure
ptoppressure on model top layer (constant for a model)
pssurface pressure - usually varying in time,x,y
pressurepressure level values in the same unit as ptop and ps
levelsigma output values at the same place as ps
Returns
MIFI_OK on success or MIFI_ERROR on failure
int mifi_atmosphere_sigma_pressure ( size_t  n,
double  ptop,
double  ps,
const double *  level,
double *  pressure 
)

convert a standard_name="atmosphere_sigma_coordinate" to pressure using the formula p(k) = ptop + sigma(k)*(ps-ptop)

Parameters
nsize of arrays sigma and pressure
ptoppressure on model top layer (constant for a model)
pssurface pressure - usually varying in time,x,y
levelsigma level values
pressureoutput values in the same unit as ptop and ps and at the same place as ps
Returns
MIFI_OK on success or MIFI_ERROR on failure
int mifi_barometric_height ( size_t  n,
double  P_b,
const double *  p,
double  T_b,
double *  height 
)

convert pressure to height using the inverse formula http://en.wikipedia.org/wiki/Barometric_formula

h(k) = -R*T_b/g*M * log(p(k)/P_b);

with P_b and T_b pressure and temperature at the layer b (i.e. surface)

g = 9.80665 m/s2 M = Molar mass of Earth's air (0.0289644 kg/mol) R = Universal gas constant (8.31432 N·m /(mol·K) )

Parameters
nsize of array h and pressure
P_bpressure at base-layer (i.e. surface, or means-sea-level) - usually varying in time,x,y
ppressure at level
T_btemperature at base layer in K - usually varying in time,x,y
heightoutput values, height above base_layer in m
Warning
This function has not been tested against possibly existing implementations
float mifi_barometric_layer_thickness ( float  p_low_alti,
float  p_high_alti,
float  T 
)

convert pressure difference to layer thickness using the hyspometric equation http://en.wikipedia.org/wiki/Hypsometric_equation

Parameters
p_low_altipressure at lower altitude in hPa
p_high_altipressure at higher altitude in hPa
Ttemperature in K
Returns
layer thickness in m
int mifi_barometric_pressure ( size_t  n,
double  P_b,
const double *  h,
double  T_b,
double *  pressure 
)

convert height above base-layer (usually altitude, i.e. height above MSL) to pressure using the formula http://en.wikipedia.org/wiki/Barometric_formula

P(h) = P_b exp[ -gM/R * h/T_b ]

with P_b and T_b pressure and temperature at the layer b (i.e. surface) and h_b the height above the layer b

g = 9.80665 m/s2 M = Molar mass of Earth's air (0.0289644 kg/mol) R = Universal gas constant (8.31432 N·m /(mol·K) )

Parameters
nsize of array h and pressure
P_bpressure at base-layer (e.g. means-sea-level) - usually varying in time,x,y
hheight in m above base-layer
T_btemperature at base layer in K - usually varying in time,x,y
pressureoutput values in the same unit as p_b and at the same place as ps
Warning
This function has not been tested against possibly existing implementations
int mifi_barometric_standard_altitude ( size_t  n,
const double *  p,
double *  altitude 
)

convert pressure to height using the formula http://en.wikipedia.org/wiki/Barometric_formula and using the international standard atmosphere http://en.wikimedia.org/wiki/International_Standard_Atmosphere

Parameters
altitudeoutput values, i.e. height above mean sea level in m
int mifi_barometric_standard_pressure ( size_t  n,
const double *  h,
double *  pressure 
)

convert altitude to pressure using the formula http://en.wikipedia.org/wiki/Barometric_formula and using the international standard atmosphere http://en.wikipedia.org/wiki/International_Standard_Atmosphere

float mifi_dewpoint_to_relative_humidity ( float  dew,
float  t 
)

Calculate relative humidity from dewpoint temperature and temperature. Derived from Bolton, 1980: http://www.eol.ucar.edu/projects/ceop/dm/documents/refdata_report/eqns.html

Parameters
dewdewpoint temperature in K
ttemperature in K
Returns
relative humidity in %
int mifi_ocean_s_g1_z ( size_t  n,
double  h,
double  h_c,
double  zeta,
const double *  s,
const double *  C,
double *  z 
)

convert a standard_name="ocean_s_coordinate_g1" to z using the formula z(k) = h_c*sigma(k) + C(k)*(h - h_c) + zeta*(1+(h_c*sigma(k)+C(k)*(h-h_c)/h))

Parameters
nsize of arrays s and C pressure
hunperturbed water column thickness - varying in x,y (might be varying in time for sediment applications)
h_ccritical depth, usually min(h(x,y))
zetatime-varying free surface - varying in time,x,y
ss(k) dimensionless s coordinate
Czeta(k) streching function
zoutput values in the same unit as h and h_c and at the same place as h
Returns
MIFI_OK on success or MIFI_ERROR on failure
See also
https://www.myroms.org/wiki/index.php?title=Vertical_S-coordinate
int mifi_ocean_s_g2_z ( size_t  n,
double  h,
double  h_c,
double  zeta,
const double *  s,
const double *  C,
double *  z 
)

convert a standard_name="ocean_s_coordinate_g1" to z using the formula z(k) = zeta + (zeta +h)*(h_c*sigma + h*C(k))/(h_c+h)

Parameters
nsize of arrays s and C pressure
hunperturbed water column thickness - varying in x,y (might be varying in time for sediment applications)
h_ccritical depth, usually min(h(x,y))
zetatime-varying free surface - varying in time,x,y
ss(k) dimensionless s coordinate
Czeta(k) streching function
zoutput values in the same unit as h and h_c and at the same place as h
Returns
MIFI_OK on success or MIFI_ERROR on failure
See also
https://www.myroms.org/wiki/index.php?title=Vertical_S-coordinate
int mifi_omega_to_vertical_wind ( size_t  n,
const double *  omega,
const double *  p,
const double *  t,
double *  w 
)

convert the vertical pressure change omega to vertical wind-speed using omega = - rho * g * w

Parameters
nsize of the array omega, p, t, w
omegavertical flow in Pa/s
ppressure in Pa
ttemperature in K
woutput-array for vertical wind speed in m/s
Returns
MIFI_OK status
float mifi_relative_to_specific_humidity ( float  rh,
float  t,
float  p 
)

Calculate specific humidity from relative humidity and temperature. Derived from HIRLAM version 2.6.0, by G. Cats

Parameters
rhrelative humidity in %
temperaturein K
ppressure in hPa
Returns
specific humidity in kg/kg
float mifi_specific_to_relative_humidity ( float  sh,
float  t,
float  p 
)

Calculate specific humidity from relative humidity and temperature. Derived from HIRLAM version 2.6.0, by G. Cats

Parameters
shspecific humidity in kg/kg
temperaturein K
ppressure in hPa
Returns
relative humidity in %
int mifi_vertical_wind_to_omega ( size_t  n,
const double *  w,
const double *  p,
const double *  t,
double *  omega 
)

convert the vertical wind-speed in m/s to vertical pressure change omega using omega = - rho * g * w

Parameters
nsize of the array omega, p, t, w
wvertical wind speed in m/s
ppressure in Pa
ttemperature in K
omegaoutput-array of vertical flow in Pa/s
Returns
MIFI_OK status
float mifi_virtual_temperature ( float  spec_humidity,
float  t 
)

calculate virtual temperature from specific humidity and temperature

Parameters
spec_humidityspecific humidity in kg/kg
temperaturein K
Returns
virtual temperature in K