Masked arrays and constants

[dennissergeev]

I think there is a bug in calc.equivalent_potential_temperature function (and potentially in a number of similar calc functions).
If I understand correctly, it appears if masked arrays are used as input. The constant kappa is defined as
0.00028562982892500527 joule * second ** 2 / gram / meter ** 2
and not "quite" dimensionless (needs to be multiplied by 1000). Then on the line https://github.com/Unidata/MetPy/blob/master/metpy/calc/thermo.py#L579
masked arrays are raised to the power kappa, which in turn uses numpy.ma.getdata, returning the non-scaled magnitude of 0.00028562982892500527. Thus the conversion to potential temperature becomes incorrect.
Minimal example:

import metpy as mp
import numpy as np

a = mp.units.masked_array(270, data_units='kelvin') * (1000 / np.ma.masked_array([700])) ** (mp.constants.kappa)
# Incorrect:
# a = (270.0275081920714)kelvin

b = mp.units.masked_array(270, data_units='kelvin') * (1000 / np.array([700])) ** (mp.constants.kappa)
# Correct:
# b = (298.95676438544837)kelvin

One solution is to use mp.constants.kappa.to('dimensionless'). I guess that's easiest solution (either in the constants or in calc submodules)

[dopplershift]

😱 ^%!@$%#!$@% unit incompatibility. 💩 👿 🔥

Sorry. I’m better now. Thanks a lot for reporting the problem.

Given that we’re pretty much always wanting to use kappa as an exponent, we should probably just convert it to dimensionless by default. So in constants.py it should read:

kappa = poisson_exponent = (Rd / Cp_d).to(‘dimensionless’)

This would ensure that we don’t miss any places where this is a problem, since silent stripping of units would at least give the nominally expected behavior.

Would you be willing to submit a PR with that fix?

[jrleeman]

Great catch @dennissergeev - units are great when they work 😉

[dennissergeev]

@dopplershift yep, now I better understand your SciPy2017 talk :)

I've just submitted a PR (#625) that should fix this issue.

I also added .to(‘dimensionless’) to the functions dealing with mixing ratio and specific humidity to be on the safe side. However, I'm not sure if I spotted all the cases where this could be an issue...

[dopplershift]

Closed by #625.