Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add TT, VT and CT calculate function. #2113

Merged
merged 3 commits into from
Sep 24, 2021
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 3 additions & 0 deletions docs/_templates/overrides/metpy.calc.rst
Original file line number Diff line number Diff line change
Expand Up @@ -71,6 +71,7 @@ Soundings
bunkers_storm_motion
cape_cin
critical_angle
cross_totals
el
k_index
lcl
Expand All @@ -89,6 +90,8 @@ Soundings
storm_relative_helicity
supercell_composite
surface_based_cape_cin
total_totals_index
vertical_totals


Dynamic/Kinematic
Expand Down
4 changes: 4 additions & 0 deletions docs/api/references.rst
Original file line number Diff line number Diff line change
Expand Up @@ -134,6 +134,10 @@ References
.. [Markowski2010] Markowski, P. and Y. Richardson, 2010: *Mesoscale Meteorology in the
Midlatitudes*. Wiley, 430 pp.

.. [Miller1972] Miller, R. C, 1972: *Notes on Analysis and Severe-Storm Forecasting Procedures
of the Air Force Global Weather Central*. `AD0744042
<https://apps.dtic.mil/sti/citations/AD0744042>`_.

.. [Moritz2000] Moritz, H., 2000: Geodetic Reference System 1980.
*Journal of Geodesy* **74**, 128-133, doi:`10.1007/s001900050278
<https://doi.org/10.1007/s001900050278>`_.
Expand Down
131 changes: 131 additions & 0 deletions src/metpy/calc/thermo.py
Original file line number Diff line number Diff line change
Expand Up @@ -3436,3 +3436,134 @@ def showalter_index(pressure, temperature, dewpoint):

# Calculate the Showalter index
return tp500 - prof[-1]


@exporter.export
@preprocess_and_wrap()
@check_units('[pressure]', '[temperature]', '[temperature]')
def total_totals_index(pressure, temperature, dewpoint):
"""Calculate Total Totals Index from the pressure temperature and dewpoint.

Total Totals Index formula derived from [Miller1972]_:
TT = (T850 + Td850) - (2 * T500)

where:

T850 is the temperature at 850 hPa
T500 is the temperature at 500 hPa
Td850 is the dewpoint at 850 hPa

Calculation of the Total Totals Index is defined as the temperature at 850 hPa plus
the dewpoint at 850 hPa, minus twice the temperature at 500 hPa. This index consists of
two components, the Vertical Totals (VT) and the Cross Totals (CT).

Parameters
----------
pressure : `pint.Quantity`
Pressure level(s), in order from highest to lowest pressure

temperature : `pint.Quantity`
Temperature corresponding to pressure

dewpoint : `pint.Quantity`
Dewpoint temperature corresponding to pressure

Returns
-------
`pint.Quantity`
Total Totals Index

"""
# Find temperature and dewpoint at 850 and 500 hPa.
(t850, t500), (td850, _) = interpolate_1d(units.Quantity([850, 500], 'hPa'),
pressure, temperature, dewpoint)

# Calculate total totals index.
tt_index = (t850 - t500) + (td850 - t500)

return tt_index


@exporter.export
@preprocess_and_wrap()
@check_units('[pressure]', '[temperature]')
def vertical_totals(pressure, temperature):
"""Calculate Vertical Totals from the pressure and temperature.

Vertical Totals formula derived from [Miller1972]_:
VT = T850 - T500

where:

T850 is the temperature at 850 hPa
T500 is the temperature at 500 hPa

Calculation of the Vertical Totals is defined as the temperature difference between
850 hPa and 500 hPa. This is a part of the Total Totals Index.

Parameters
----------
pressure : `pint.Quantity`
Pressure level(s), in order from highest to lowest pressure

temperature : `pint.Quantity`
Temperature corresponding to pressure

Returns
-------
`pint.Quantity`
Vertical Totals

"""
# Find temperature at 850 and 500 hPa.
(t850, t500) = interpolate_1d(units.Quantity([850, 500], 'hPa'),
pressure, temperature)

# Calculate vertical totals.
vt = t850 - t500

return vt


@exporter.export
@preprocess_and_wrap()
@check_units('[pressure]', '[temperature]', '[temperature]')
def cross_totals(pressure, temperature, dewpoint):
"""Calculate Cross Totals from the pressure temperature and dewpoint.

Cross Totals formula derived from [Miller1972]_:
CT = Td850 - T500

where:

Td850 is the dewpoint at 850 hPa
T500 is the temperature at 500 hPa

Calculation of the Cross Totals is defined as the difference between dewpoint
at 850 hPa and temperature at 500 hPa. This is a part of the Total Totals Index.

Parameters
----------
pressure : `pint.Quantity`
Pressure level(s), in order from highest to lowest pressure

temperature : `pint.Quantity`
Temperature corresponding to pressure

dewpoint : `pint.Quantity`
Dewpoint temperature corresponding to pressure

Returns
-------
`pint.Quantity`
Cross Totals

"""
# Find temperature and dewpoint at 850 and 500 hPa
(_, t500), (td850, _) = interpolate_1d(units.Quantity([850, 500], 'hPa'),
pressure, temperature, dewpoint)

# Calculate vertical totals.
ct = td850 - t500

return ct
79 changes: 74 additions & 5 deletions tests/calc/test_thermo.py
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,7 @@
import xarray as xr

from metpy.calc import (brunt_vaisala_frequency, brunt_vaisala_frequency_squared,
brunt_vaisala_period, cape_cin, density, dewpoint,
brunt_vaisala_period, cape_cin, cross_totals, density, dewpoint,
dewpoint_from_relative_humidity, dewpoint_from_specific_humidity,
dry_lapse, dry_static_energy, el, equivalent_potential_temperature,
exner_function, gradient_richardson_number, InvalidSoundingError,
Expand All @@ -27,10 +27,10 @@
specific_humidity_from_dewpoint, specific_humidity_from_mixing_ratio,
static_stability, surface_based_cape_cin,
temperature_from_potential_temperature, thickness_hydrostatic,
thickness_hydrostatic_from_relative_humidity, vapor_pressure,
vertical_velocity, vertical_velocity_pressure,
virtual_potential_temperature, virtual_temperature,
wet_bulb_temperature)
thickness_hydrostatic_from_relative_humidity, total_totals_index,
vapor_pressure, vertical_totals, vertical_velocity,
vertical_velocity_pressure, virtual_potential_temperature,
virtual_temperature, wet_bulb_temperature)
from metpy.calc.thermo import _find_append_zero_crossings
from metpy.testing import assert_almost_equal, assert_array_almost_equal, assert_nan
from metpy.units import masked_array, units
Expand Down Expand Up @@ -1950,3 +1950,72 @@ def test_showalter_index():

result = showalter_index(pressure, temps, dewp)
assert_almost_equal(result, units.Quantity(7.6024, 'delta_degC'), 4)


def test_total_totals_index():
"""Test the Total Totals Index calculation."""
pressure = np.array([1008., 1000., 947., 925., 921., 896., 891., 889., 866.,
858., 850., 835., 820., 803., 733., 730., 700., 645.,
579., 500., 494., 466., 455., 441., 433., 410., 409.,
402., 400., 390., 388., 384., 381., 349., 330., 320.,
306., 300., 278., 273., 250., 243., 208., 200., 196.,
190., 179., 159., 151., 150., 139.]) * units.hPa
temperature = np.array([27.4, 26.4, 22.9, 21.4, 21.2, 20.7, 20.6, 21.2, 19.4,
19.1, 18.8, 17.8, 17.4, 16.3, 11.4, 11.2, 10.2, 6.1,
0.6, -4.9, -5.5, -8.5, -9.9, -11.7, -12.3, -13.7, -13.8,
-14.9, -14.9, -16.1, -16.1, -16.9, -17.3, -21.7, -24.5, -26.1,
-28.3, -29.5, -33.1, -34.2, -39.3, -41., -50.2, -52.5, -53.5,
-55.2, -58.6, -65.2, -68.1, -68.5, -72.5]) * units.degC
dewpoint = np.array([24.9, 24.6, 22., 20.9, 20.7, 14.8, 13.6, 12.2, 16.8,
16.6, 16.5, 15.9, 13.6, 13.2, 11.3, 11.2, 8.6, 4.5,
-0.8, -8.1, -9.5, -12.7, -12.7, -12.8, -13.1, -24.7, -24.4,
-21.9, -24.9, -36.1, -31.1, -26.9, -27.4, -33., -36.5, -47.1,
-31.4, -33.5, -40.1, -40.8, -44.1, -45.6, -54., -56.1, -56.9,
-58.6, -61.9, -68.4, -71.2, -71.6, -77.2]) * units.degC

tt = total_totals_index(pressure, temperature, dewpoint)
assert_almost_equal(tt, 45.10 * units.delta_degC, 2)


def test_vertical_totals():
"""Test the Vertical Totals calculation."""
pressure = np.array([1008., 1000., 947., 925., 921., 896., 891., 889., 866.,
858., 850., 835., 820., 803., 733., 730., 700., 645.,
579., 500., 494., 466., 455., 441., 433., 410., 409.,
402., 400., 390., 388., 384., 381., 349., 330., 320.,
306., 300., 278., 273., 250., 243., 208., 200., 196.,
190., 179., 159., 151., 150., 139.]) * units.hPa
temperature = np.array([27.4, 26.4, 22.9, 21.4, 21.2, 20.7, 20.6, 21.2, 19.4,
19.1, 18.8, 17.8, 17.4, 16.3, 11.4, 11.2, 10.2, 6.1,
0.6, -4.9, -5.5, -8.5, -9.9, -11.7, -12.3, -13.7, -13.8,
-14.9, -14.9, -16.1, -16.1, -16.9, -17.3, -21.7, -24.5, -26.1,
-28.3, -29.5, -33.1, -34.2, -39.3, -41., -50.2, -52.5, -53.5,
-55.2, -58.6, -65.2, -68.1, -68.5, -72.5]) * units.degC

vt = vertical_totals(pressure, temperature)
assert_almost_equal(vt, 23.70 * units.delta_degC, 2)


def test_cross_totals():
"""Test the Cross Totals calculation."""
pressure = np.array([1008., 1000., 947., 925., 921., 896., 891., 889., 866.,
858., 850., 835., 820., 803., 733., 730., 700., 645.,
579., 500., 494., 466., 455., 441., 433., 410., 409.,
402., 400., 390., 388., 384., 381., 349., 330., 320.,
306., 300., 278., 273., 250., 243., 208., 200., 196.,
190., 179., 159., 151., 150., 139.]) * units.hPa
temperature = np.array([27.4, 26.4, 22.9, 21.4, 21.2, 20.7, 20.6, 21.2, 19.4,
19.1, 18.8, 17.8, 17.4, 16.3, 11.4, 11.2, 10.2, 6.1,
0.6, -4.9, -5.5, -8.5, -9.9, -11.7, -12.3, -13.7, -13.8,
-14.9, -14.9, -16.1, -16.1, -16.9, -17.3, -21.7, -24.5, -26.1,
-28.3, -29.5, -33.1, -34.2, -39.3, -41., -50.2, -52.5, -53.5,
-55.2, -58.6, -65.2, -68.1, -68.5, -72.5]) * units.degC
dewpoint = np.array([24.9, 24.6, 22., 20.9, 20.7, 14.8, 13.6, 12.2, 16.8,
16.6, 16.5, 15.9, 13.6, 13.2, 11.3, 11.2, 8.6, 4.5,
-0.8, -8.1, -9.5, -12.7, -12.7, -12.8, -13.1, -24.7, -24.4,
-21.9, -24.9, -36.1, -31.1, -26.9, -27.4, -33., -36.5, -47.1,
-31.4, -33.5, -40.1, -40.8, -44.1, -45.6, -54., -56.1, -56.9,
-58.6, -61.9, -68.4, -71.2, -71.6, -77.2]) * units.degC

ct = cross_totals(pressure, temperature, dewpoint)
assert_almost_equal(ct, 21.40 * units.delta_degC, 2)