galvez_davison_index#
- metpy.calc.galvez_davison_index(pressure, temperature, mixing_ratio, surface_pressure, vertical_dim=0)[source]#
Calculate GDI from the pressure, temperature, mixing ratio, and surface pressure.
Calculation of the GDI relies on temperatures and mixing ratios at 950, 850, 700, and 500 hPa. These four levels define three layers: A) Boundary, B) Trade Wind Inversion (TWI), C) Mid-Troposphere.
GDI formula derived from [Galvez2015]:
\[GDI = CBI + MWI + II + TC\]where:
\(CBI\) is the Column Buoyancy Index
\(MWI\) is the Mid-tropospheric Warming Index
\(II\) is the Inversion Index
\(TC\) is the Terrain Correction [optional]
# GDI Value
Expected Convective Regime
>=45
Scattered to widespread thunderstorms likely.
35 to 45
Scattered thunderstorms and/or scattered to widespread rain showers.
25 to 35
Isolated to scattered thunderstorms and/or scattered showers.
15 to 25
Isolated thunderstorms and/or isolated to scattered showers.
5 to 10
Isolated to scattered showers.
<5
Strong TWI likely, light rain possible.
- Parameters:
pressure (
pint.Quantity
) – Pressure level(s)temperature (
pint.Quantity
) – Temperature corresponding to pressuremixing_ratio (
pint.Quantity
) – Mixing ratio values corresponding to pressuresurface_pressure (
pint.Quantity
) – Pressure of the surface.vertical_dim (int, optional) – The axis corresponding to vertical, defaults to 0. Automatically determined from xarray DataArray arguments.
- Returns:
pint.Quantity
– GDI Index
Examples
>>> from metpy.calc import mixing_ratio_from_relative_humidity >>> from metpy.units import units >>> # pressure >>> p = [1008., 1000., 950., 900., 850., 800., 750., 700., 650., 600., ... 550., 500., 450., 400., 350., 300., 250., 200., ... 175., 150., 125., 100., 80., 70., 60., 50., ... 40., 30., 25., 20.] * units.hPa >>> # temperature >>> T = [29.3, 28.1, 23.5, 20.9, 18.4, 15.9, 13.1, 10.1, 6.7, 3.1, ... -0.5, -4.5, -9.0, -14.8, -21.5, -29.7, -40.0, -52.4, ... -59.2, -66.5, -74.1, -78.5, -76.0, -71.6, -66.7, -61.3, ... -56.3, -51.7, -50.7, -47.5] * units.degC >>> # relative humidity >>> rh = [.85, .65, .36, .39, .82, .72, .75, .86, .65, .22, .52, ... .66, .64, .20, .05, .75, .76, .45, .25, .48, .76, .88, ... .56, .88, .39, .67, .15, .04, .94, .35] * units.dimensionless >>> # calculate mixing ratio >>> mixrat = mixing_ratio_from_relative_humidity(p, T, rh) >>> galvez_davison_index(p, T, mixrat, p[0]) <Quantity(-8.78797532, 'dimensionless')>