File: Thickness_Hydrostatic.py

package info (click to toggle)
metpy 1.7.1%2Bds-1
  • links: PTS, VCS
  • area: main
  • in suites: forky, sid
  • size: 18,584 kB
  • sloc: python: 41,853; makefile: 111; javascript: 57
file content (51 lines) | stat: -rw-r--r-- 2,069 bytes parent folder | download | duplicates (2)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
# Copyright (c) 2022 MetPy Developers.
# Distributed under the terms of the BSD 3-Clause License.
# SPDX-License-Identifier: BSD-3-Clause
"""
=====================
Hydrostatic Thickness
=====================

Use functions from `metpy.calc` as well as pint's unit support to perform calculations.

The code below uses example data from our test suite to calculate the hydrostatic thickness
between the surface and 500-hPa level for the provided sounding data.
"""
import pandas as pd

from metpy.calc import (mixing_ratio_from_relative_humidity, relative_humidity_from_dewpoint,
                        thickness_hydrostatic)
from metpy.cbook import get_test_data
from metpy.units import units

###########################################
# Upper air data can be obtained using the siphon package, but for this example we will use
# some of MetPy's sample data.

# Set column names
col_names = ['pressure', 'height', 'temperature', 'dewpoint', 'direction', 'speed']

# Read in test data using col_names
df = pd.read_fwf(get_test_data('jan20_sounding.txt', as_file_obj=False),
                 skiprows=5, usecols=[0, 1, 2, 3, 6, 7], names=col_names)

###########################################
# Drop any rows with all NaN values for T, Td, winds
df = df.dropna(subset=('temperature', 'dewpoint', 'direction', 'speed'
                       ), how='all').reset_index(drop=True)

###########################################
# Isolate pressure, temperature, and dewpoint and add units
p = df['pressure'].values * units.hPa
T = df['temperature'].values * units.degC
Td = df['dewpoint'].values * units.degC

###########################################
# Calculate the relative humidity to compute the mixing ratio
rh = relative_humidity_from_dewpoint(T, Td)
mixrat = mixing_ratio_from_relative_humidity(p, T, rh)

###########################################
# Calculate the thickness from the pressure, temperature, and mixing ratio
# for the layer from the surface pressure to 500-hPa
print(thickness_hydrostatic(p, T, mixing_ratio=mixrat, depth=p[0] - 500 * units.hPa))