DataArray accessor unit_array property returns dimensionless

[DanielAdriaansen]

I have an Xarray Dataset object containing temperature (T) and relative humidity (RH). I am using MetPy's dewpoint_from_relative_humidity() to compute dewpoint temperature (TD). A code snippet looks like this:

import metpy.calc as mpcalc
TD = mpcalc.dewpoint_from_relative_humidity(ds['T'],ds['RH'])

This returns an ndarray object (TD) with pint quantity attached as expected. For convenience, I wanted to add TD to my Dataset. I did that like this:

ds['TD'] = xr.DataArray(TD,dims=dims,coords=ds.coords)

Then I wanted to compute precipitable water (PW) using precipitable_water(). When precipitable_water() is called, it tells me that TD is dimensionless. I traced this to the MetPy DataArray accessor. It's trying to access the units attribute of the DataArray containing TD, however since I created the DataArray myself I did not attach a units attribute since I created it using an ndarray object with pint quantity in it already. I confirm that the pint quantity is still attached by executing:

print ds['TD']

which shows me it has units of degree_Celsius.

In the documentation here, it tells me the unit_array property will obtain a pint.Quantity array of just the data from the DataArray. I think this is only true if the DataArray has an attribute of 'units' and doesn't work even if the DataArray has a pint.Quantity attached to it.

I did a little digging and found @jthielen issue over at the pint repo. So it seems like this is a known issue but perhaps the documentation could be updated to reflect that currently a units attr is required or the unit_array property updated to explicitly check for a units attribute and error out if there isn't one attached to the DataArray?

I also don't know what the difference is between these two items of the DataArrayAccessor:
@property
def unit_array(self):
"""Return the data values of this DataArray as a pint.Quantity."""
return self._data_array.values * self.units

@unit_array.setter
def unit_array(self, values):
"""Set data values from a pint.Quantity."""
self._data_array.values = values.magnitude
self._units = self._data_array.attrs['units'] = str(values.units)

It looks like the setter does explicitly access the units attr, but the property (where I am having trouble) doesn't explicitly access the units attr but it seems implicit based on adding some print statements and observing the behavior.

The trace of the error is:
calc/indices.py where preprocess_xarray() is called prior to precipitable_water() being executed, which is trying to access the unit_array property of the DataArray.

[DanielAdriaansen]

OK, I am using rc1, and just noticed the code in the DataArrayAccessor is different in the "latest" (master) branch. I am going to clone the latest and observe the behavior.

[jthielen]

@DanielAdriaansen Thank you for bringing this up. This is an unfortunate consequence of MetPy's unit functionality on its accessor (at least until #1325) predating xarray's support of Pint Quanities inside of Variables/DataArrays/Datasets (as of v0.13.0). Any currently released version of MetPy (including v0.12.1 and v1.0.0rc1) will have this issue.

After some discussion (starting in #1304 (comment)), it was decided that MetPy should take advantage of Quantity-in-xarray support wherever possible rather than continue with its old behavior (which has the issues you've pointed out) or try to equally support both Quantity-in-xarray and units-on-attribute (which would likely be confusing for users). These changes were made in #1325, and are reflected in the current state of the master branch as you saw. The tutorial, however, has not been brought up to date, since that is waiting on the even more substantial changes of #1353.

If you'd want to learn more about the current state of things, I'd encourage you to read up on these PRs: #1304, #1325, #1353. In short, however, in not too long (i.e., once #1353 is complete and v1.0.0 final is released), your use case will be able to be simplified greatly, as MetPy's functions will return DataArrays when supplied DataArrays:

import metpy.calc as mpcalc
ds['TD'] = mpcalc.dewpoint_from_relative_humidity(ds['T'],ds['RH'])
ds['PW'] = mpcalc.precipitable_water(ds['TD'], ds['P'])

I still think your suggestion of adding a warning in the documentation that pre-v1.0.0 MetPy requires units to be on the xarray attribute and not on a Quantity wrapped by the xarray object is a good idea.

[dopplershift]

Yeah, the xarray+units is a moving target right now (all the more reason to get 1.0 or even maybe an rc2? out). I think a warning in the docs is appropriate.

[dopplershift]

@DanielAdriaansen Would you say this issue was addressed by the warning in the docs?

[DanielAdriaansen]

@dopplershift indeed it does! Thanks to @jthielen for making the update to the docs and for all of the background reading material. It seems I poked the proverbial hornets nest here, but as usual the MetPy developers are on it already. I am glad to hear that MetPy's functions will return DataArrays when a DataArray is supplied in future versions that will be great. Thanks again everyone.