Skip to content

Commit

Permalink
PEP8 clean-up.
Browse files Browse the repository at this point in the history
Ignoring files that will shortly be almost completely replaced.
  • Loading branch information
dopplershift committed Mar 15, 2014
1 parent 838a839 commit d8b5118
Show file tree
Hide file tree
Showing 6 changed files with 2,400 additions and 1,613 deletions.
79 changes: 48 additions & 31 deletions metpy/calc/basic.py
Original file line number Diff line number Diff line change
@@ -1,14 +1,17 @@
'A collection of generic calculation functions.'

__all__ = ['vapor_pressure', 'dewpoint', 'get_speed_dir','get_wind_components',
'mixing_ratio','tke', 'windchill', 'heat_index', 'h_convergence',
'v_vorticity', 'convergence_vorticity', 'advection', 'geostrophic_wind']
__all__ = ['vapor_pressure', 'dewpoint', 'get_speed_dir',
'get_wind_components', 'mixing_ratio', 'tke', 'windchill',
'heat_index', 'h_convergence', 'v_vorticity',
'convergence_vorticity', 'advection', 'geostrophic_wind']

import numpy as np
from numpy.ma import log, exp, cos, sin, masked_array
from scipy.constants import degree, kilo, hour, g

sat_pressure_0c = 6.112 # mb
sat_pressure_0c = 6.112 # mb


def vapor_pressure(temp):
'''
Calculate the saturation water vapor (partial) pressure given
Expand All @@ -26,6 +29,7 @@ def vapor_pressure(temp):
'''
return sat_pressure_0c * exp(17.67 * temp / (temp + 243.5))


def dewpoint(temp, rh):
'''
Calculate the ambient dewpoint given air temperature and relative
Expand All @@ -42,9 +46,10 @@ def dewpoint(temp, rh):
of the result being determined using numpy's broadcasting rules.
'''
es = vapor_pressure(temp)
val = log(rh * es/sat_pressure_0c)
val = log(rh * es / sat_pressure_0c)
return 243.5 * val / (17.67 - val)


def mixing_ratio(part_press, tot_press):
'''
Calculates the mixing ratio of gas given its partial pressure
Expand All @@ -64,22 +69,25 @@ def mixing_ratio(part_press, tot_press):
'''
return part_press / (tot_press - part_press)

def get_speed_dir(u,v,w=None):

def get_speed_dir(u, v, w=None):
'''
Compute the wind speed (horizontal and vector is W is supplied) and
wind direction.
Return horizontal wind speed, vector wind speed, and wind direction in a tuple
* if w is not supplied, returns tuple of horizontal wind speed and wind direction.
Return horizontal wind speed, vector wind speed, and wind direction in
a tuple. If w is not supplied, returns tuple of horizontal wind speed
and wind direction.
'''
hws = np.sqrt(u*u+v*v)
wd = np.arctan2(-u,-v)*180./np.pi
wd[wd<0]=360+wd[wd<0]
hws = np.sqrt(u * u + v * v)
wd = np.rad2deg(np.arctan2(-u, -v))
wd[wd < 0] = 360. + wd[wd < 0]
if w is None:
return hws,wd
return hws, wd
else:
vws = np.sqrt(u*u+v*v+w*w)
return hws,vws,wd
vws = np.sqrt(u * u + v * v + w * w)
return hws, vws, wd


def get_wind_components(speed, wdir):
'''
Expand All @@ -99,9 +107,10 @@ def get_wind_components(speed, wdir):
wdir = wdir * degree
u = -speed * sin(wdir)
v = -speed * cos(wdir)
return u,v
return u, v

def tke(u,v,w):

def tke(u, v, w):
'''
Compute the turbulence kinetic energy (tke) from the time series of the
velocity components u, v, and w.
Expand All @@ -123,13 +132,14 @@ def tke(u,v,w):
wp = w - w.mean()

tke = np.power(np.average(np.power(up, 2)) +
np.average(np.power(vp, 2)) +
np.average(np.power(wp, 2)), 0.5)
np.average(np.power(vp, 2)) +
np.average(np.power(wp, 2)), 0.5)

return tke


def windchill(temp, speed, metric=True, face_level_winds=False,
mask_undefined=True):
mask_undefined=True):
'''
Calculate the Wind Chill Temperature Index (WCTI) from the current
temperature and wind speed.
Expand Down Expand Up @@ -178,25 +188,26 @@ def windchill(temp, speed, metric=True, face_level_winds=False,
if metric:
# Formula uses wind speed in km/hr, but passing in m/s makes more
# sense. Convert here.
temp_limit, speed_limit = 10., 4.828 #Temp in C, speed in km/h
temp_limit, speed_limit = 10., 4.828 # Temp in C, speed in km/h
speed = speed * hour / kilo
speed_factor = speed ** 0.16
wcti = (13.12 + 0.6215 * temp - 11.37 * speed_factor
+ 0.3965 * temp * speed_factor)
+ 0.3965 * temp * speed_factor)
else:
temp_limit, speed_limit = 50., 3.
speed_factor = speed ** 0.16
wcti = (35.74 + 0.6215 * temp - 35.75 * speed_factor
+ 0.4275 * temp * speed_factor)
+ 0.4275 * temp * speed_factor)

#See if we need to mask any undefined values
# See if we need to mask any undefined values
if mask_undefined:
mask = np.array((temp > temp_limit) | (speed <= speed_limit))
if mask.any():
wcti = masked_array(wcti, mask=mask)

return wcti


def heat_index(temp, rh, mask_undefined=True):
'''
Calculate the Heat Index from the current temperature and relative
Expand Down Expand Up @@ -225,14 +236,14 @@ def heat_index(temp, rh, mask_undefined=True):
temperature-humidity index based on human physiology and clothing
science. J. Appl. Meteor., 18, 861-873.
'''
rh2 = rh**2
temp2 = temp**2
rh2 = rh ** 2
temp2 = temp ** 2

# Calculate the Heat Index
HI = (-42.379 + 2.04901523 * temp + 10.14333127 * rh
- 0.22475541 * temp * rh - 6.83783e-3 * temp2 - 5.481717e-2 * rh2
+ 1.22874e-3 * temp2 * rh + 8.5282e-4 * temp * rh2
- 1.99e-6 * temp2 * rh2)
- 0.22475541 * temp * rh - 6.83783e-3 * temp2 - 5.481717e-2 * rh2
+ 1.22874e-3 * temp2 * rh + 8.5282e-4 * temp * rh2
- 1.99e-6 * temp2 * rh2)

# See if we need to mask any undefined values
if mask_undefined:
Expand All @@ -242,12 +253,14 @@ def heat_index(temp, rh, mask_undefined=True):

return HI


def _get_gradients(u, v, dx, dy):
#Helper function for getting convergence and vorticity from 2D arrays
# Helper function for getting convergence and vorticity from 2D arrays
dudx, dudy = np.gradient(u, dx, dy)
dvdx, dvdy = np.gradient(v, dx, dy)
return dudx, dudy, dvdx, dvdy


def v_vorticity(u, v, dx, dy):
'''
Calculate the vertical vorticity of the horizontal wind. The grid
Expand All @@ -269,6 +282,7 @@ def v_vorticity(u, v, dx, dy):
dudx, dudy, dvdx, dvdy = _get_gradients(u, v, dx, dy)
return dvdx - dudy


def h_convergence(u, v, dx, dy):
'''
Calculate the horizontal convergence of the horizontal wind. The grid
Expand All @@ -290,6 +304,7 @@ def h_convergence(u, v, dx, dy):
dudx, dudy, dvdx, dvdy = _get_gradients(u, v, dx, dy)
return dudx + dvdy


def convergence_vorticity(u, v, dx, dy):
'''
Calculate the horizontal convergence and vertical vorticity of the
Expand All @@ -313,6 +328,7 @@ def convergence_vorticity(u, v, dx, dy):
dudx, dudy, dvdx, dvdy = _get_gradients(u, v, dx, dy)
return dudx + dvdy, dvdx - dudy


def advection(scalar, wind, deltas):
'''
Calculate the advection of *scalar* by the wind. The order of the
Expand Down Expand Up @@ -345,10 +361,11 @@ def advection(scalar, wind, deltas):

# Make them be at least 2D (handling the 1D case) so that we can do the
# multiply and sum below
grad,wind = np.atleast_2d(grad, wind)
grad, wind = np.atleast_2d(grad, wind)

return (-grad * wind).sum(axis=0)


def geostrophic_wind(heights, f, dx, dy, geopotential=False):
'''
Calculate the geostrophic wind given from the heights. If geopotential
Expand Down Expand Up @@ -389,5 +406,5 @@ def geostrophic_wind(heights, f, dx, dy, geopotential=False):
deltas = deltas + [1.] * (heights.ndim - 2)

grad = np.gradient(heights, *deltas)
dx,dy = grad[0], grad[1] # This throws away unused gradient components
dx, dy = grad[0], grad[1] # This throws away unused gradient components
return -norm_factor * dy, norm_factor * dx
28 changes: 19 additions & 9 deletions metpy/calc/tests/test_calc.py
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,7 @@
from metpy.calc import *
from metpy.constants import g


class TestVaporPressure(TestCase):
def test_basic(self):
temp = np.array([5, 10, 18, 25])
Expand All @@ -13,6 +14,7 @@ def test_scalar(self):
es = vapor_pressure(0)
assert_almost_equal(es, 6.112, 3)


class TestDewpoint(TestCase):
def test_basic(self):
temp = np.array([30, 25, 10, 20, 25])
Expand All @@ -25,14 +27,15 @@ def test_scalar(self):
td = dewpoint(10.6, .37) * 1.8 + 32.
assert_almost_equal(td, 26, 0)


class TestWindComps(TestCase):
def test_basic(self):
'Test the basic calculation.'
speed = np.array([4, 4, 4, 4, 25, 25, 25, 25, 10.])
dirs = np.array([0, 45, 90, 135, 180, 225, 270, 315, 360])
s2 = np.sqrt(2)

u,v = get_wind_components(speed, dirs)
u, v = get_wind_components(speed, dirs)

true_u = np.array([0, -4/s2, -4, -4/s2, 0, 25/s2, 25, 25/s2, 0])
true_v = np.array([-4, -4/s2, 0, 4/s2, 25, 25/s2, 0, -25/s2, -10])
Expand All @@ -44,6 +47,7 @@ def test_scalar(self):
comps = np.array(get_wind_components(8, 150))
assert_array_almost_equal(comps, np.array([-4, 6.9282]), 3)


class TestWindChill(TestCase):
def test_basic(self):
'Test the basic wind chill calculation.'
Expand Down Expand Up @@ -85,6 +89,7 @@ def test_undefined_flag(self):
mask = np.array([False]*6)
assert_array_equal(wc.mask, mask)


class TestHeatIndex(TestCase):
def test_basic(self):
'Test the basic heat index calculation.'
Expand Down Expand Up @@ -117,6 +122,7 @@ def test_undefined_flag(self):
mask = np.array([False]*6)
assert_array_equal(hi.mask, mask)


#class TestIrrad(TestCase):
# def test_basic(self):
# 'Test the basic solar irradiance calculation.'
Expand Down Expand Up @@ -151,11 +157,12 @@ def test_undefined_flag(self):
# False, False, True, True, True])
# assert_array_equal(s.mask, mask)


class TestGradients(TestCase):
def test_basic(self):
'Basic braindead test of vorticity and divergence calculation'
u = np.ones((3,3))
c,v = convergence_vorticity(u, u, 1, 1)
u = np.ones((3, 3))
c, v = convergence_vorticity(u, u, 1, 1)
truth = np.zeros_like(u)
assert_array_equal(c, truth)
assert_array_equal(v, truth)
Expand All @@ -164,7 +171,7 @@ def test_basic2(self):
'Basic test of vorticity and divergence calculation'
a = np.arange(3)
u = np.c_[a, a, a]
c,v = convergence_vorticity(u, u.T, 1, 1)
c, v = convergence_vorticity(u, u.T, 1, 1)
true_c = 2. * np.ones_like(u)
true_v = np.zeros_like(u)
assert_array_equal(c, true_c)
Expand All @@ -174,12 +181,13 @@ def test_basic3(self):
'Basic test of vorticity and divergence calculation'
a = np.arange(3)
u = np.c_[a, a, a]
c,v = convergence_vorticity(u, u, 1, 1)
c, v = convergence_vorticity(u, u, 1, 1)
true_c = np.ones_like(u)
true_v = np.ones_like(u)
assert_array_equal(c, true_c)
assert_array_equal(v, true_v)


class TestAdvection(TestCase):
def test_basic(self):
'Basic braindead test of advection'
Expand Down Expand Up @@ -207,22 +215,23 @@ def test_basic3(self):

def test_2dbasic(self):
'Basic 2D braindead test of advection'
u = np.ones((3,3))
v = np.ones((3,3))
u = np.ones((3, 3))
v = np.ones((3, 3))
s = np.ones_like(u)
a = advection(s, u, (1,))
truth = np.zeros_like(u)
assert_array_equal(a, truth)

def test_2dbasic2(self):
'Basic 2D test of advection'
u = np.ones((3,3))
v = 2 * np.ones((3,3))
u = np.ones((3, 3))
v = 2 * np.ones((3, 3))
s = np.array([[1, 2, 1], [2, 4, 2], [1, 2, 1]])
a = advection(s, [u, v], (1, 1))
truth = np.array([[-3, -2, 1], [-4, 0, 4], [-1, 2, 3]])
assert_array_equal(a, truth)


class TestGeos(TestCase):
def test_basic(self):
'Basic test of geostrophic wind calculation'
Expand All @@ -234,5 +243,6 @@ def test_basic(self):
assert_array_equal(ug, true_u)
assert_array_equal(vg, true_v)


if __name__ == '__main__':
run_module_suite()
Loading

0 comments on commit d8b5118

Please sign in to comment.