diff --git a/src/metpy/xarray.py b/src/metpy/xarray.py index 9b0b4e0e8fa..d9d0cb7edf6 100644 --- a/src/metpy/xarray.py +++ b/src/metpy/xarray.py @@ -174,6 +174,22 @@ def convert_units(self, units): """ return self.quantify().copy(data=self.unit_array.to(units)) + def convert_to_base_units(self): + """Return new DataArray with values converted to base units. + + See Also + -------- + convert_units + + Notes + ----- + Any cached/lazy-loaded data (except that in a Dask array) will be loaded into memory + by this operation. Do not utilize on moderate- to large-sized remote datasets before + subsetting! + + """ + return self.quantify().copy(data=self.unit_array.to_base_units()) + def convert_coordinate_units(self, coord, units): """Return new DataArray with specified coordinate converted to different units. diff --git a/tests/test_xarray.py b/tests/test_xarray.py index 09d44ec0225..08055017399 100644 --- a/tests/test_xarray.py +++ b/tests/test_xarray.py @@ -153,6 +153,19 @@ def test_convert_units(test_var): assert_almost_equal(result[0, 0, 0, 0], 18.44 * units.degC, 2) +def test_convert_to_base_units(test_ds): + """Test conversion of units.""" + uwnd = test_ds.u_wind.metpy.quantify() + result = (uwnd * (500 * units.hPa)).metpy.convert_to_base_units() + + # Check that units are updated without modifying original + assert result.metpy.units == units('kg s**-3') + assert test_ds.u_wind.metpy.units == units('m/s') + + # Make sure we now get an array back with properly converted values + assert_almost_equal(result[0, 0, 0, 0], -448416.12 * units('kg s**-3'), 2) + + def test_convert_coordinate_units(test_ds_generic): """Test conversion of coordinate units.""" result = test_ds_generic['test'].metpy.convert_coordinate_units('b', 'percent') diff --git a/tutorials/xarray_tutorial.py b/tutorials/xarray_tutorial.py index 784f367de55..8e676ccf20a 100644 --- a/tutorials/xarray_tutorial.py +++ b/tutorials/xarray_tutorial.py @@ -227,6 +227,12 @@ temperature_degc = temperature[0].metpy.convert_units('degC') temperature_degc +######################################################################### +# To base unit conversion: + +temperature_degk = temperature_degc.metpy.convert_to_base_units() +temperature_degk + ######################################################################### # Unit conversion for coordinates: heights_on_hpa_levels = heights.metpy.convert_coordinate_units('isobaric3', 'hPa')