Skip to content

Commit

Permalink
Merge pull request #1253 from dopplershift/fix-smooth-units
Browse files Browse the repository at this point in the history
Fix smoothing unit-handling
  • Loading branch information
dopplershift committed Dec 23, 2019
2 parents c2b87ee + 33a030e commit ff1bb66
Show file tree
Hide file tree
Showing 4 changed files with 91 additions and 22 deletions.
14 changes: 6 additions & 8 deletions src/metpy/calc/basic.py
Expand Up @@ -656,6 +656,7 @@ def sigma_to_pressure(sigma, psfc, ptop):

@exporter.export
@preprocess_xarray
@units.wraps('=A', ('=A',))
def smooth_gaussian(scalar_grid, n):
"""Filter with normal distribution of weights.
Expand Down Expand Up @@ -742,15 +743,13 @@ def smooth_gaussian(scalar_grid, n):
# Assume the last two axes represent the horizontal directions
sgma_seq = [sgma if i > nax - 3 else 0 for i in range(nax)]

# Compute smoothed field and reattach units
res = gaussian_filter(scalar_grid, sgma_seq, truncate=2 * np.sqrt(2))
if hasattr(scalar_grid, 'units'):
res = res * scalar_grid.units
return res
# Compute smoothed field
return gaussian_filter(scalar_grid, sgma_seq, truncate=2 * np.sqrt(2))


@exporter.export
@preprocess_xarray
@units.wraps('=A', ('=A',))
def smooth_n_point(scalar_grid, n=5, passes=1):
"""Filter with normal distribution of weights.
Expand All @@ -764,8 +763,7 @@ def smooth_n_point(scalar_grid, n=5, passes=1):
are 5 and 9. Defaults to 5.
passes : int
The number of times to apply the filter to the grid. Defaults
to 1.
The number of times to apply the filter to the grid. Defaults to 1.
Returns
-------
Expand All @@ -782,7 +780,7 @@ def smooth_n_point(scalar_grid, n=5, passes=1):
masked value or NaN values exists in the array, it will propagate
to any point that uses that particular grid point in the smoothing
calculation. Applying the smoothing function multiple times will
propogate NaNs further throughout the domain.
propagate NaNs further throughout the domain.
"""
if n == 9:
Expand Down
13 changes: 13 additions & 0 deletions src/metpy/testing.py
Expand Up @@ -131,6 +131,17 @@ def check_and_drop_units(actual, desired):
return actual, desired


def check_mask(actual, desired):
"""Check that two arrays have the same mask.
This handles the fact that `~numpy.testing.assert_array_equal` ignores masked values
in either of the arrays. This ensures that the masks are identical.
"""
actual_mask = getattr(actual, 'mask', np.full(np.asarray(actual).shape, False))
desired_mask = getattr(desired, 'mask', np.full(np.asarray(desired).shape, False))
np.testing.assert_array_equal(actual_mask, desired_mask)


def assert_nan(value, units):
"""Check for nan with proper units."""
value, _ = check_and_drop_units(value, np.nan * units)
Expand All @@ -152,6 +163,7 @@ def assert_array_almost_equal(actual, desired, decimal=7):
Wrapper around :func:`numpy.testing.assert_array_almost_equal`
"""
actual, desired = check_and_drop_units(actual, desired)
check_mask(actual, desired)
numpy.testing.assert_array_almost_equal(actual, desired, decimal)


Expand All @@ -161,6 +173,7 @@ def assert_array_equal(actual, desired):
Wrapper around :func:`numpy.testing.assert_array_equal`
"""
actual, desired = check_and_drop_units(actual, desired)
check_mask(actual, desired)
numpy.testing.assert_array_equal(actual, desired)


Expand Down
62 changes: 48 additions & 14 deletions tests/calc/test_basic.py
Expand Up @@ -165,10 +165,9 @@ def test_windchill_invalid():

wc = windchill(temp, speed)
# We don't care about the masked values
truth = np.array([2.6230789, np.nan, np.nan, np.nan, np.nan, np.nan]) * units.degF
mask = np.array([False, True, True, True, True, True])
truth = np.ma.array([2.6230789, np.nan, np.nan, np.nan, np.nan, np.nan],
mask=[False, True, True, True, True, True]) * units.degF
assert_array_almost_equal(truth, wc)
assert_array_equal(wc.mask, mask)


def test_windchill_undefined_flag():
Expand Down Expand Up @@ -249,7 +248,7 @@ def test_heat_index_vs_nws():
temp = units.Quantity(np.array([86, 111, 40, 96]), units.degF)
rh = np.ma.array([45, 27, 99, 60]) * units.percent
hi = heat_index(temp, rh)
truth = units.Quantity(np.ma.array([87, 121, 40, 116]), units.degF)
truth = np.ma.array([87, 121, 40, 116], mask=[False, False, True, False]) * units.degF
assert_array_almost_equal(hi, truth, 0)


Expand Down Expand Up @@ -413,8 +412,9 @@ def test_apparent_temperature():
[10, 10, 10]]) * units.percent
wind = np.array([[5, 3, 3],
[10, 1, 10]]) * units.mph
truth = np.array([[99.6777178, 86.3357671, 70],
[8.8140662, 20, 60]]) * units.degF
truth = np.ma.array([[99.6777178, 86.3357671, 70],
[8.8140662, 20, 60]],
mask=[[False, False, True], [False, True, True]]) * units.degF
res = apparent_temperature(temperature, rel_humidity, wind)
assert_array_almost_equal(res, truth, 6)

Expand Down Expand Up @@ -589,24 +589,58 @@ def test_smooth_n_pt_9_repeat():
[5816., 5784., 5744., 5716., 5684.]])
shght = smooth_n_point(hght, 9, 2)
s_true = np.array([[5640., 5640., 5640., 5640., 5640.],
[5684., 5675.4375, 5666.9375, 5658.8125, 5651.],
[5728., 5710.875, 5693.875, 5677.625, 5662.],
[5772., 5746.375, 5720.625, 5696.375, 5673.],
[5816., 5784., 5744., 5716., 5684.]])
[5684., 5675.4375, 5666.9375, 5658.8125, 5651.],
[5728., 5710.875, 5693.875, 5677.625, 5662.],
[5772., 5746.375, 5720.625, 5696.375, 5673.],
[5816., 5784., 5744., 5716., 5684.]])
assert_array_almost_equal(shght, s_true)


def test_smooth_n_pt_wrong_number():
"""Test the smooth_n_pt function using wrong number of points."""
hght = np.array([[5640., 5640., 5640., 5640., 5640.],
[5684., 5676., 5666., 5659., 5651.],
[5728., 5712., 5692., 5678., 5662.],
[5772., 5748., 5718., 5697., 5673.],
[5816., 5784., 5744., 5716., 5684.]])
[5684., 5676., 5666., 5659., 5651.],
[5728., 5712., 5692., 5678., 5662.],
[5772., 5748., 5718., 5697., 5673.],
[5816., 5784., 5744., 5716., 5684.]])
with pytest.raises(ValueError):
smooth_n_point(hght, 7)


def test_smooth_n_pt_temperature():
"""Test the smooth_n_pt function with temperature units."""
t = np.array([[2.73, 3.43, 6.53, 7.13, 4.83],
[3.73, 4.93, 6.13, 6.63, 8.23],
[3.03, 4.83, 6.03, 7.23, 7.63],
[3.33, 4.63, 7.23, 6.73, 6.23],
[3.93, 3.03, 7.43, 9.23, 9.23]]) * units.degC

smooth_t = smooth_n_point(t, 9, 1)
smooth_t_true = np.array([[2.73, 3.43, 6.53, 7.13, 4.83],
[3.73, 4.6425, 5.96125, 6.81124, 8.23],
[3.03, 4.81125, 6.1175, 6.92375, 7.63],
[3.33, 4.73625, 6.43, 7.3175, 6.23],
[3.93, 3.03, 7.43, 9.23, 9.23]]) * units.degC
assert_array_almost_equal(smooth_t, smooth_t_true, 4)


def test_smooth_gaussian_temperature():
"""Test the smooth_gaussian function with temperature units."""
t = np.array([[2.73, 3.43, 6.53, 7.13, 4.83],
[3.73, 4.93, 6.13, 6.63, 8.23],
[3.03, 4.83, 6.03, 7.23, 7.63],
[3.33, 4.63, 7.23, 6.73, 6.23],
[3.93, 3.03, 7.43, 9.23, 9.23]]) * units.degC

smooth_t = smooth_gaussian(t, 3)
smooth_t_true = np.array([[2.8892, 3.7657, 6.2805, 6.8532, 5.3174],
[3.6852, 4.799, 6.0844, 6.7816, 7.7617],
[3.2762, 4.787, 6.117, 7.0792, 7.5181],
[3.4618, 4.6384, 6.886, 6.982, 6.6653],
[3.8115, 3.626, 7.1705, 8.8528, 8.9605]]) * units.degC
assert_array_almost_equal(smooth_t, smooth_t_true, 4)


def test_altimeter_to_station_pressure_inhg():
"""Test the altimeter to station pressure function with inches of mercury."""
altim = 29.8 * units.inHg
Expand Down
24 changes: 24 additions & 0 deletions tests/test_testing.py
@@ -0,0 +1,24 @@
# Copyright (c) 2019 MetPy Developers.
# Distributed under the terms of the BSD 3-Clause License.
# SPDX-License-Identifier: BSD-3-Clause
"""Test MetPy's testing utilities."""

import numpy as np
import pytest

from metpy.testing import assert_array_almost_equal


# Test #1183: numpy.testing.assert_array* ignores any masked value, so work-around
def test_masked_arrays():
"""Test that we catch masked arrays with different masks."""
with pytest.raises(AssertionError):
assert_array_almost_equal(np.array([10, 20]),
np.ma.array([10, np.nan], mask=[False, True]), 2)


def test_masked_and_no_mask():
"""Test that we can compare a masked array with no masked values and a regular array."""
a = np.array([10, 20])
b = np.ma.array([10, 20], mask=[False, False])
assert_array_almost_equal(a, b)

0 comments on commit ff1bb66

Please sign in to comment.