Skip to content

Moist Physics

clausius_clapeyron_factor(temp, pressure)

This is the factor \(\alpha\), such that \(dq^*/dT = \alpha q^*\).

I explicitly compute alpha from the formula for saturation_vapor_pressure in the function here.

Parameters:

Name Type Description Default
temp ndarray

Temperature at each coordinate considered. Units: Kelvin.

 required
pressure Union[float, ndarray]

Pressure level in Pa, temperature corresponds to. If all temp are at the lowest atmospheric level, then pressurewill be the lowest level pressure i.e. afloat`.

 required

Returns:

Type Description
ndarray

Clausius clapeyron factor, \(\alpha\). Units: Kelvin\(^{-1}\)
Source code in isca_tools/utils/moist_physics.py 
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
def clausius_clapeyron_factor(temp: np.ndarray, pressure: Union[float, np.ndarray]) -> np.ndarray:
    """
    This is the factor $\\alpha$, such that $dq^*/dT = \\alpha q^*$.

    I explicitly compute alpha from the formula for `saturation_vapor_pressure` in the function here.

    Args:
        temp: Temperature at each coordinate considered. Units: *Kelvin*.
        pressure: Pressure level in *Pa*, temperature corresponds to.
            If all `temp` are at the lowest atmospheric level, then pressure` will be the lowest level pressure i.e.
            a `float`.


    Returns:
        Clausius clapeyron factor, $\\alpha$. Units: *Kelvin$^{-1}$*
    """
    lambda_const = 4302.645 / (temp - 29.65) ** 2
    return lambda_const * pressure / epsilon * sphum_sat(temp, pressure) / saturation_vapor_pressure(temp)

get_density(temp, pressure, sphum=None)

Equation for density using ideal gas equation of state: \(\rho = \frac{p}{RT}\). If specific humidity given, will compute density using virtual temperature, \(T_v\).

Parameters:

Name Type Description Default
temp Union[float, ndarray]

float [n_p_levels] Temperature in K to find density at.

 required
pressure Union[float, ndarray]

float [n_p_levels] Pressure in Pa to find density at.

 required
sphum Optional[Union[float, ndarray]]

float [n_p_levels] Specific humidity in kg/kg to find density at.

 None

Returns:

Type Description
Union[float, ndarray]

Density in units of \(kg/m^3\).
Source code in isca_tools/utils/moist_physics.py 
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
def get_density(temp: Union[float, np.ndarray], pressure: Union[float, np.ndarray],
                sphum: Optional[Union[float, np.ndarray]] = None) -> Union[float, np.ndarray]:
    """
    Equation for density using ideal gas equation of state: $\\rho = \\frac{p}{RT}$.
    If specific humidity given, will compute density using virtual temperature, $T_v$.

    Args:
        temp: `float [n_p_levels]`
            Temperature in *K* to find density at.
        pressure: `float [n_p_levels]`
            Pressure in *Pa* to find density at.
        sphum: `float [n_p_levels]`
            Specific humidity in *kg/kg* to find density at.

    Returns:
        Density in units of $kg/m^3$.
    """
    if sphum is None:
        return pressure / (R * temp)
    else:
        return pressure / (R * virtual_temp(temp, sphum))

mixing_ratio_from_partial_pressure(partial_pressure, total_pressure)

Computes the mixing ratio, \(w\), from partial pressure, \(e\), and total atmospheric pressure, \(p\), according to:

\(w = \epsilon \frac{e}{p-e}\)

Where \(\epsilon = R_d/R_v = 0.622\) is the ratio of molecular weight of water to that of dry air.

This is the same equation used by MetPy.

Parameters:

Name Type Description Default
partial_pressure Union[float, ndarray]

float [n_levels]. Partial pressure at each level, \(e\), in Pa.

 required
total_pressure Union[float, ndarray]

float [n_levels]. Atmospheric pressure at each level, \(p\), in Pa.

 required

Returns:

Type Description
Union[float, ndarray]

Mixing ratio, \(w\), in units of \(kg/kg\).
Source code in isca_tools/utils/moist_physics.py 
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
def mixing_ratio_from_partial_pressure(partial_pressure: Union[float, np.ndarray],
                                       total_pressure: Union[float, np.ndarray]) -> Union[float, np.ndarray]:
    """
    Computes the mixing ratio, $w$, from partial pressure, $e$, and total atmospheric pressure, $p$, according to:

    $w = \epsilon \\frac{e}{p-e}$

    Where $\epsilon = R_d/R_v = 0.622$ is the ratio of molecular weight of water to that of dry air.

    This is the same equation used by
    [MetPy](https://unidata.github.io/MetPy/latest/api/generated/metpy.calc.mixing_ratio.html).

    Args:
        partial_pressure: `float [n_levels]`. Partial pressure at each level, $e$, in *Pa*.
        total_pressure: `float [n_levels]`. Atmospheric pressure at each level, $p$, in *Pa*.

    Returns:
        Mixing ratio, $w$, in units of $kg/kg$.
    """
    return epsilon * partial_pressure / (total_pressure - partial_pressure)

mixing_ratio_from_sphum(sphum)

Computes the mixing ratio, \(w\), from specific humidity, \(q\), according to:

\(w = q/(1-q)\)

Parameters:

Name Type Description Default
sphum Union[float, ndarray]

Specific humidity, \(q\), in units of \(kg/kg\).

 required

Returns:

Type Description
Union[float, ndarray]

Mixing ratio, \(w\), in units of \(kg/kg\).
Source code in isca_tools/utils/moist_physics.py 
50
51
52
53
54
55
56
57
58
59
60
61
62
def mixing_ratio_from_sphum(sphum: Union[float, np.ndarray]) -> Union[float, np.ndarray]:
    """
    Computes the mixing ratio, $w$, from specific humidity, $q$, according to:

    $w = q/(1-q)$

    Args:
        sphum: Specific humidity, $q$, in units of $kg/kg$.

    Returns:
        Mixing ratio, $w$, in units of $kg/kg$.
    """
    return sphum / (1 - sphum)

moist_static_energy(temp, sphum, height, c_p_const=c_p)

Returns the moist static energy in units of kJ/kg.

Parameters:

Name Type Description Default
temp Union[ndarray, float]

float [n_lat, n_p_levels]. Temperature at each coordinate considered. Units: Kelvin.

 required
sphum Union[ndarray, float]

float [n_lat, n_p_levels]. Specific humidity at each coordinate considered. Units: kg/kg.

 required
height Union[ndarray, float]

float [n_lat, n_p_levels] or float. Geopotential height of each level considered. Just a float if only one pressure level considered for each latitude e.g. common to use 2m values. Units: m.

 required
c_p_const float

Heat capacity constant in units of J/K/kg. This gives the option to easily modify the moist static energy but almost always be kept at default c_p.

 c_p

Returns:

Type Description
Union[ndarray, float]

Moist static energy at each coordinate given.
Source code in isca_tools/utils/moist_physics.py 
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
def moist_static_energy(temp: Union[np.ndarray, float], sphum: Union[np.ndarray, float],
                        height: Union[np.ndarray, float],
                        c_p_const: float = c_p) -> Union[np.ndarray, float]:
    """
    Returns the moist static energy in units of *kJ/kg*.

    Args:
        temp: `float [n_lat, n_p_levels]`. Temperature at each coordinate considered. Units: *Kelvin*.
        sphum: `float [n_lat, n_p_levels]`. Specific humidity at each coordinate considered. Units: *kg/kg*.
        height: `float [n_lat, n_p_levels]` or `float`. Geopotential height of each level considered.
            Just a `float` if only one pressure level considered for each latitude e.g. common to use 2m values.
            Units: *m*.
        c_p_const: Heat capacity constant in units of J/K/kg.
            This gives the option to easily modify the moist static energy but almost always be kept at default `c_p`.

    Returns:
        Moist static energy at each coordinate given.
    """
    return (L_v * sphum + c_p_const * temp + g * height) / 1000

partial_pressure_from_sphum(sphum, total_pressure)

Computes the partial pressure, \(e\), from specific humidity, \(q\), and total atmospheric pressure, \(p\), according to:

\(e = wp / (\epsilon+w)\)

where \(w\) is the mixing ratio, calculated using mixing_ratio_from_sphum, and \(\epsilon = R_d/R_v = 0.622\) is the ratio of molecular weight of water to that of dry air.

Parameters:

Name Type Description Default
sphum Union[float, ndarray]

Specific humidity, \(q\), in units of \(kg/kg\).

 required
total_pressure Union[float, ndarray]

float [n_levels]. Atmospheric pressure at each level, \(p\), in Pa.

 required

Returns:

Name Type Description
partial_pressure Union[float, ndarray]

float [n_levels]. Partial pressure at each level, \(e\), in Pa.
Source code in isca_tools/utils/moist_physics.py 
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
def partial_pressure_from_sphum(sphum: Union[float, np.ndarray],
                                total_pressure: Union[float, np.ndarray]) -> Union[float, np.ndarray]:
    """
    Computes the partial pressure, $e$, from specific humidity, $q$, and total atmospheric pressure, $p$, according to:

    $e = wp / (\epsilon+w)$

    where $w$ is the mixing ratio, calculated using `mixing_ratio_from_sphum`, and
    $\epsilon = R_d/R_v = 0.622$ is the ratio of molecular weight of water to that of dry air.

    Args:
        sphum: Specific humidity, $q$, in units of $kg/kg$.
        total_pressure: `float [n_levels]`. Atmospheric pressure at each level, $p$, in *Pa*.

    Returns:
        partial_pressure: `float [n_levels]`.
            Partial pressure at each level, $e$, in *Pa*.
    """
    w = mixing_ratio_from_sphum(sphum)
    return w * total_pressure / (epsilon + w)

rh_from_sphum(sphum, temp, total_pressure)

Relative humidity, \(rh\), is computed from specific humidity, \(q\) according to:

\(rh = e / e_s\)

Where, \(e = pw/(\epsilon + w)\), is the partial pressure, \(e_s\) is the saturation partial pressure and \(w\) is the mixing ratio.

Parameters:

Name Type Description Default
sphum Union[float, ndarray]

float [n_levels]. Specific humidity, \(q\), at each level considered, in units of \(kg/kg\).

 required
temp Union[float, ndarray]

float [n_levels]. Temperature at each level considered. Units: Kelvin.

 required
total_pressure Union[float, ndarray]

float [n_levels]. Atmospheric pressure, \(p\), at each level considered in Pa.

 required

Returns:

Type Description
Union[float, ndarray]

Percentage relative humidity (\(0 < rh < 100\)).
Source code in isca_tools/utils/moist_physics.py 
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
def rh_from_sphum(sphum: Union[float, np.ndarray], temp: Union[float, np.ndarray],
                  total_pressure: Union[float, np.ndarray]) -> Union[float, np.ndarray]:
    """
    Relative humidity, $rh$, is computed from specific humidity, $q$ according to:

    $rh = e / e_s$

    Where, $e = pw/(\epsilon + w)$, is the partial pressure, $e_s$ is the saturation partial pressure and
     $w$ is the mixing ratio.

    Args:
        sphum: `float [n_levels]`. Specific humidity, $q$, at each level considered, in units of $kg/kg$.
        temp: `float [n_levels]`. Temperature at each level considered. Units: *Kelvin*.
        total_pressure: `float [n_levels]`. Atmospheric pressure, $p$, at each level considered in *Pa*.

    Returns:
        Percentage relative humidity ($0 < rh < 100$).
    """
    sat_mix_ratio = mixing_ratio_from_partial_pressure(saturation_vapor_pressure(temp), total_pressure)
    mix_ratio = mixing_ratio_from_sphum(sphum)
    # return 100 * mix_ratio / sat_mix_ratio
    return 100 * mix_ratio / (epsilon + mix_ratio) * (epsilon + sat_mix_ratio) / sat_mix_ratio

saturation_vapor_pressure(temp)

Computes the saturation vapor pressure, \(e_s(T)\), corresponding to a given temperature.

Uses Equation 10 in Bolton 1980. Valid for \(-35^\circ C < T < 35^\circ C\).

Parameters:

Name Type Description Default
temp Union[float, ndarray]

Temperature to compute vapor pressure at. Units: Kelvin.

 required

Returns:

Type Description
Union[float, ndarray]

Saturation vapor pressure, \(e_s(T)\), in units of Pa.
Source code in isca_tools/utils/moist_physics.py 
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
def saturation_vapor_pressure(temp: Union[float, np.ndarray]) -> Union[float, np.ndarray]:
    """
    Computes the saturation vapor pressure, $e_s(T)$, corresponding to a given temperature.

    Uses *Equation 10* in *Bolton 1980*. Valid for $-35^\circ C < T < 35^\circ C$.

    Args:
        temp: Temperature to compute vapor pressure at. Units: *Kelvin*.

    Returns:
        Saturation vapor pressure, $e_s(T)$, in units of *Pa*.
    """
    # Alternative equation from MATLAB exercise M9.2 in Holdon 2004
    # return 611 * np.exp(L_v/R_v * (1/temp_kelvin_to_celsius - 1/temp))
    temp = temp - temp_kelvin_to_celsius  # Convert temperature in kelvin to celsius, as celsius used for this formula.
    # if np.abs(np.asarray(temp)).max() > 35:
    #     warnings.warn('This formula is only valid for $-35^\circ C < T < 35^\circ C$\n'
    #                   'At least one temperature given is outside this range.')
    # Multiply by 100 below to convert from hPa to Pa.
    return 611.2 * np.exp(17.67 * temp / (temp + 243.5))

sphum_from_dew(temp_dew, pressure)

Calculates specific humidity from dew temperature at a given pressure. The dew temperature is defined such as the temperature at which the saturation vapour pressure is equal to the partial pressure i.e. \(e_s(T_d) = e(T)\) where \(T_d\) is the dew temperature, and \(T\) is the actual temperature.

Parameters:

Name Type Description Default
temp_dew Union[float, ndarray]

float [n_lat, n_levels]. Dew temperature at each coordinate considered. Units: Kelvin.

 required
pressure Union[float, ndarray]

float [n_levels]. Atmospheric pressure, \(p\), at each level considered in Pa.

 required

Returns:

Name Type Description
specific_humidity Union[float, ndarray]

float [n_lat, n_levels]. Specific humidity, \(q\), in units of \(kg/kg\).
Source code in isca_tools/utils/moist_physics.py 
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
def sphum_from_dew(temp_dew: Union[float, np.ndarray], pressure: Union[float, np.ndarray]) -> Union[float, np.ndarray]:
    """
    Calculates specific humidity from dew temperature at a given `pressure`.
    The dew temperature is defined such as the temperature at which the saturation vapour pressure
    is equal to the partial pressure i.e. $e_s(T_d) = e(T)$ where $T_d$ is the dew temperature, and
    $T$ is the actual temperature.

    Args:
        temp_dew: `float [n_lat, n_levels]`. Dew temperature at each coordinate considered. Units: *Kelvin*.
        pressure: `float [n_levels]`. Atmospheric pressure, $p$, at each level considered in *Pa*.

    Returns:
        specific_humidity: `float [n_lat, n_levels]`.
            Specific humidity, $q$, in units of $kg/kg$.
    """
    # dew temperature is defined as temp such that saturation vapour pressure equals partial pressure
    e = saturation_vapor_pressure(temp_dew)
    return sphum_from_partial_pressure(e, pressure)

sphum_from_partial_pressure(partial_pressure, total_pressure)

Computes the specific humidity, \(q\), from partial pressure, \(e\), and total atmospheric pressure, \(p\), according to:

\(q = w / (1+w)\)

where \(w\) is the mixing ratio, calculated using mixing_ratio_from_partial_pressure.

Parameters:

Name Type Description Default
partial_pressure Union[float, ndarray]

float [n_levels]. Partial pressure at each level, \(e\), in Pa.

 required
total_pressure Union[float, ndarray]

float [n_levels]. Atmospheric pressure at each level, \(p\), in Pa.

 required

Returns:

Type Description
Union[float, ndarray]

Specific humidity, \(q\), in units of \(kg/kg\).
Source code in isca_tools/utils/moist_physics.py 
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
def sphum_from_partial_pressure(partial_pressure: Union[float, np.ndarray],
                                       total_pressure: Union[float, np.ndarray]) -> Union[float, np.ndarray]:
    """
    Computes the specific humidity, $q$, from partial pressure, $e$, and total atmospheric pressure, $p$, according to:

    $q = w / (1+w)$

    where $w$ is the mixing ratio, calculated using `mixing_ratio_from_partial_pressure`.

    Args:
        partial_pressure: `float [n_levels]`. Partial pressure at each level, $e$, in *Pa*.
        total_pressure: `float [n_levels]`. Atmospheric pressure at each level, $p$, in *Pa*.

    Returns:
        Specific humidity, $q$, in units of $kg/kg$.
    """
    w = mixing_ratio_from_partial_pressure(partial_pressure, total_pressure)
    return w / (1+w)

sphum_sat(temp, pressure)

Returns the saturation specific humidity, \(q^*\), in kg/kg.

Parameters:

Name Type Description Default
temp Union[float, ndarray]

Temperature at each coordinate considered. Units: Kelvin.

 required
pressure Union[float, ndarray]

Pressure level in Pa, temperature corresponds to. If all temp are at the lowest atmospheric level, then pressurewill be the lowest level pressure i.e. afloat`.

 required

Returns:

Type Description
Union[float, ndarray]

Saturation specific humidity at each coordinate given.
Source code in isca_tools/utils/moist_physics.py 
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
def sphum_sat(temp: Union[float, np.ndarray], pressure: Union[float, np.ndarray]) -> Union[float, np.ndarray]:
    """
    Returns the saturation specific humidity, $q^*$, in *kg/kg*.

    Args:
        temp: Temperature at each coordinate considered. Units: *Kelvin*.
        pressure: Pressure level in *Pa*, temperature corresponds to.
            If all `temp` are at the lowest atmospheric level, then pressure` will be the lowest level pressure i.e.
            a `float`.

    Returns:
        Saturation specific humidity at each coordinate given.
    """
    # Saturation specific humidity
    w_sat = mixing_ratio_from_partial_pressure(saturation_vapor_pressure(temp), pressure)
    q_sat = w_sat / (1 + w_sat)
    return q_sat

virtual_temp(temp, sphum)

Equation for virtual temperature using Isca code.

The constants d622, d378, d608 are to match the Isca code.

Parameters:

Name Type Description Default
temp Union[float, ndarray]

float [n_p_levels] Temperature in K to find virtual potential temperature at.

 required
sphum Union[float, ndarray]

float [n_p_levels] Specific humidity of parcel at each pressure level in kg/kg.

 required

Returns:

Type Description
Union[float, ndarray]

Virtual temperature at each pressure level in K.
Source code in isca_tools/utils/moist_physics.py 
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
def virtual_temp(temp: Union[float, np.ndarray], sphum: Union[float, np.ndarray]) -> Union[float, np.ndarray]:
    """
    Equation for virtual temperature using [Isca code](https://github.com/jduffield65/Isca/blob/b9249275469583c1723f12ac62333067f9460fea/isca_source/src/coupler/surface_flux.F90#L463).

    The constants `d622`, `d378`, `d608` are to match the
    [Isca code](https://github.com/jduffield65/Isca/blob/b9249275469583c1723f12ac62333067f9460fea/isca_source/src/coupler/surface_flux.F90#L935-L940).

    Args:
        temp: `float [n_p_levels]`
            Temperature in *K* to find virtual potential temperature at.
        sphum: `float [n_p_levels]`
            Specific humidity of parcel at each pressure level in *kg/kg*.

    Returns:
        Virtual temperature at each pressure level in *K*.
    """
    d622 = R / R_v
    d378 = 1 - d622
    d608 = d378 / d622
    return (1 + d608 * sphum) * temp