CF units raising temperature by decimal power

38 views
Skip to first unread message

Daniel Brener

unread,
Jul 8, 2019, 9:35:21 AM7/8/19
to SciTools (iris, cartopy, cf_units, etc.) - https://github.com/scitools
I want to raise the temperatures in my cube by the power 3.5.

CF units won't let me (it screens lots of exceptions no matter what I try).

How can I do this or bypass the CF?

Cheers

Daniel

Klaus Zimmermann

unread,
Jul 9, 2019, 10:02:48 AM7/9/19
to scitoo...@googlegroups.com
Hi Daniel,

could you provide us with a minimal working example to see the
exceptions in action?

Also, what are you actually trying to do? Units of (degC)^(3.5) are
rather unusual, so here a little more context might help us to advise
better.

Cheers
Klaus

Daniel Brener

unread,
Jul 10, 2019, 4:26:49 AM7/10/19
to SciTools (iris, cartopy, cf_units, etc.) - https://github.com/scitools
Hi Klaus,

Thank you for replying.

Here's a short example:
import iris
filename = '20110827.oper00.pc000'
vars = ['air_temperature']

temperature = iris.load(filename, vars)

t1 = temperature[0] ** 3.5

print(t1)

Gives the following:

(base) C:\Users\...>python poc_cf_error.py
Traceback (most recent call last):
  File "poc_cf_error.py", line 7, in <module>
    t1 = temperature[0] ** 3.5
  File "C:\Anaconda\lib\site-packages\iris\analysis\maths.py", line 547, in exponentiate
    return _math_op_common(cube, power, cube.units ** exponent, new_dtype,
  File "C:\Anaconda\lib\site-packages\cf_units\__init__.py", line 1701, in __pow__
    raise ValueError(msg)
ValueError: Cannot raise a unit by a decimal (got 3.5).

In my attempts to debug this I tried separating the operations (e.g. cubing then rooting or vice versa). The cubing worked but rooting does not e.g.

import iris
filename = '20110827.oper00.pc000'
vars = ['air_temperature']

temperature = iris.load(filename, vars)

t1 = temperature[0] ** 0.5

print(t1)

Gives:

(base) C:\Users\...>python poc_cf_error.py
productRoot(): It's meaningless to take the 2nd root of "K"
Traceback (most recent call last):
  File "poc_cf_error.py", line 7, in <module>
    t1 = temperature[0] ** 0.5
  File "C:\Anaconda\lib\site-packages\iris\analysis\maths.py", line 547, in exponentiate
    return _math_op_common(cube, power, cube.units ** exponent, new_dtype,
  File "C:\Anaconda\lib\site-packages\cf_units\__init__.py", line 1695, in __pow__
    result = self.root(root)
  File "C:\Anaconda\lib\site-packages\cf_units\__init__.py", line 1419, in root
    six.raise_from(value_error, None)
  File "<string>", line 3, in raise_from
ValueError: [UT_MEANINGLESS] Failed to take the root of Unit('K')

I have tested this on both Windows and UNIX based systems.

Cheers

Daniel

Klaus Zimmermann

unread,
Jul 10, 2019, 4:53:25 AM7/10/19
to scitoo...@googlegroups.com
Hi Daniel,

thanks!

This seems rather consistent. I also cannot think of an interpretation
of roots of kg.

What are you actually in terms of physical meaning trying to do?

Cheers
Klaus
> --
> You received this message because you are subscribed to the Google
> Groups "SciTools (iris, cartopy, cf_units, etc.) -
> https://github.com/scitools" group.
> To unsubscribe from this group and stop receiving emails from it, send
> an email to scitools-iri...@googlegroups.com
> <mailto:scitools-iri...@googlegroups.com>.
> To view this discussion on the web, visit
> https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com
> <https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com?utm_medium=email&utm_source=footer>.
> For more options, visit https://groups.google.com/d/optout.

Daniel Brener

unread,
Jul 10, 2019, 5:05:01 AM7/10/19
to SciTools (iris, cartopy, cf_units, etc.) - https://github.com/scitools
Hi Klaus,

Square root of kelvin should be a reasonable operation to perform, it is required when calculating the dynamic viscosity of a fluid using Sutherlands formula which has T^(3/2).

Is there someway of overriding the exception?

Cheers

Daniel

Klaus Zimmermann

unread,
Jul 10, 2019, 5:13:11 AM7/10/19
to scitoo...@googlegroups.com
Hi Daniel,

correct me if I misunderstand you, but as far as I know Sutherland's
formula, it has

(T/T0)^(3/2)

with reference temperature T0.

Of course the ratio is dimensionless and there is no problem.

Cheers
Klaus
> > an email to scitoo...@googlegroups.com <javascript:>
> > <mailto:scitools-iri...@googlegroups.com <javascript:>>.
> <https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com?utm_medium=email&utm_source=footer
> <https://groups.google.com/d/optout>.
>
> --
> You received this message because you are subscribed to the Google
> Groups "SciTools (iris, cartopy, cf_units, etc.) -
> https://github.com/scitools" group.
> To unsubscribe from this group and stop receiving emails from it, send
> an email to scitools-iri...@googlegroups.com
> <mailto:scitools-iri...@googlegroups.com>.
> To view this discussion on the web, visit
> https://groups.google.com/d/msgid/scitools-iris/60c23a28-8f2d-41b0-9045-88759450ad72%40googlegroups.com
> <https://groups.google.com/d/msgid/scitools-iris/60c23a28-8f2d-41b0-9045-88759450ad72%40googlegroups.com?utm_medium=email&utm_source=footer>.

Daniel Brener

unread,
Jul 10, 2019, 5:22:42 AM7/10/19
to SciTools (iris, cartopy, cf_units, etc.) - https://github.com/scitools
Hi Klaus,

That's correct, that is Sutherlands formula but I'm actually wanting to use a formula form the NASA 1976 standard atmosphere:

dyn viscosity = beta * T^(3/2) / (T + S)

p.g. 35 in this PDF
 
My understanding from this manual is that T here is T(z) and therefore not dimensionless (2nd column p.g. 32)

Unless you know of a better empirical formula for dynamic viscosity in the atmosphere?

Cheers

Daniel
>     > <mailto:scitools-iris+unsub...@googlegroups.com <javascript:>>.

Klaus Zimmermann

unread,
Jul 10, 2019, 5:45:46 AM7/10/19
to scitoo...@googlegroups.com
Hi Daniel,

I see. Ok, in that case the end result contains only integer powers of
units, but intermediate quantities have fractional ones. I would deal
with this by writing a small function like

```
from cf_units import Unit

def dynamic_viscosity(cube):
# Sutherland constant in K
SUTHERLAND_CONSTANT = 110.4
# empirical beta constant in kg/(m s K^.5) see <reference>
BETA = 1.458e-6

new_cube = cube.copy()
# Make sure we have K
new_cube.convert_units('K')
# Get potentially lazy data to profit from dask
temp = new_cube.core_data()
# Perform calculation on data only
new_cube.data = (BETA * temp ** 1.5) / (temp + SUTHERLAND_CONSTANT)
# Fix assign known correct units
new_cube.units = Unit('N s m-2')
# Fix other metadata (standard_name, long_name, ...) here
return new_cube
```

Of course you were quite right to try and have the units calculated all
along the way. That would be cleaner and allow you to easily deal with
data coming in with compatible units, say degC. Unfortunately, I think
in this case this is the best one can do.

Hope that helps.

Cheers
Klaus
> >     > <mailto:scitools-iri...@googlegroups.com
> <javascript:> <javascript:>>.
> >     > To view this discussion on the web, visit
> >     >
> >    
> https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com
> <https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com>
>
> >    
> <https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com
> <https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com>>
>
> >
> >     >
> >    
> <https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com?utm_medium=email&utm_source=footer
> <https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com?utm_medium=email&utm_source=footer>
>
> >    
> <https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com?utm_medium=email&utm_source=footer
> <https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com?utm_medium=email&utm_source=footer>>>.
>
> >
> >     > For more options, visit https://groups.google.com/d/optout
> <https://groups.google.com/d/optout>
> >     <https://groups.google.com/d/optout
> <https://groups.google.com/d/optout>>.
> >
> > --
> > You received this message because you are subscribed to the Google
> > Groups "SciTools (iris, cartopy, cf_units, etc.) -
> > https://github.com/scitools" group.
> > To unsubscribe from this group and stop receiving emails from it,
> send
> > an email to scitoo...@googlegroups.com <javascript:>
> > <mailto:scitools-iri...@googlegroups.com <javascript:>>.
> > To view this discussion on the web, visit
> >
> https://groups.google.com/d/msgid/scitools-iris/60c23a28-8f2d-41b0-9045-88759450ad72%40googlegroups.com
> <https://groups.google.com/d/msgid/scitools-iris/60c23a28-8f2d-41b0-9045-88759450ad72%40googlegroups.com>
>
> >
> <https://groups.google.com/d/msgid/scitools-iris/60c23a28-8f2d-41b0-9045-88759450ad72%40googlegroups.com?utm_medium=email&utm_source=footer
> <https://groups.google.com/d/msgid/scitools-iris/60c23a28-8f2d-41b0-9045-88759450ad72%40googlegroups.com?utm_medium=email&utm_source=footer>>.
>
> > For more options, visit https://groups.google.com/d/optout
> <https://groups.google.com/d/optout>.
>
> --
> You received this message because you are subscribed to the Google
> Groups "SciTools (iris, cartopy, cf_units, etc.) -
> https://github.com/scitools" group.
> To unsubscribe from this group and stop receiving emails from it, send
> an email to scitools-iri...@googlegroups.com
> <mailto:scitools-iri...@googlegroups.com>.
> To view this discussion on the web, visit
> https://groups.google.com/d/msgid/scitools-iris/d92a906c-22d6-4fd4-8554-0bb78fcb23a6%40googlegroups.com
> <https://groups.google.com/d/msgid/scitools-iris/d92a906c-22d6-4fd4-8554-0bb78fcb23a6%40googlegroups.com?utm_medium=email&utm_source=footer>.

Daniel Brener

unread,
Jul 10, 2019, 5:50:45 AM7/10/19
to SciTools (iris, cartopy, cf_units, etc.) - https://github.com/scitools
Hi Klaus,

Wonderful, that's very kind of you to sketch it out for me, thank you.

That's a nice trick to get around the units problem - I'm sure I'll be using that idea in the future for other problems and no doubt this will likely help others!

Thanks again,

Daniel
>     >     > <mailto:scitools-iris+unsub...@googlegroups.com
>     <javascript:> <javascript:>>.
>     >     > To view this discussion on the web, visit
>     >     >
>     >    
>     https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com
>     <https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com>
>
>     >    
>     <https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com
>     <https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com>>
>
>     >
>     >     >
>     >    
>     <https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com?utm_medium=email&utm_source=footer
>     <https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com?utm_medium=email&utm_source=footer>
>
>     >    
>     <https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com?utm_medium=email&utm_source=footer
>     <https://groups.google.com/d/msgid/scitools-iris/7cee1915-96a1-4ae0-a22a-cc58f898c0b2%40googlegroups.com?utm_medium=email&utm_source=footer>>>.
>
>     >
>     >     > For more options, visit https://groups.google.com/d/optout
>     <https://groups.google.com/d/optout>
>     >     <https://groups.google.com/d/optout
>     <https://groups.google.com/d/optout>>.
>     >
>     > --
>     > You received this message because you are subscribed to the Google
>     > Groups "SciTools (iris, cartopy, cf_units, etc.) -
>     > https://github.com/scitools" group.
>     > To unsubscribe from this group and stop receiving emails from it,
>     send
>     > an email to scitoo...@googlegroups.com <javascript:>
>     > <mailto:scitools-iris+unsub...@googlegroups.com <javascript:>>.

Klaus Zimmermann

unread,
Jul 10, 2019, 7:15:20 AM7/10/19
to scitoo...@googlegroups.com
Hi Daniel,

I am happy you find it helpful :)

Cheers
Klaus
> https://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/19770009539.pdf <https://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/19770009539.pdf>
> >     >     > <mailto:scitools-iri...@googlegroups.com
> <javascript:>
> >     > <mailto:scitools-iri...@googlegroups.com
> <javascript:> <javascript:>>.
> >     > To view this discussion on the web, visit
> >     >
> >    
> https://groups.google.com/d/msgid/scitools-iris/60c23a28-8f2d-41b0-9045-88759450ad72%40googlegroups.com
> >     <https://groups.google.com/d/optout
> <https://groups.google.com/d/optout>>.
> >
> > --
> > You received this message because you are subscribed to the Google
> > Groups "SciTools (iris, cartopy, cf_units, etc.) -
> > https://github.com/scitools" group.
> > To unsubscribe from this group and stop receiving emails from it,
> send
> > an email to scitoo...@googlegroups.com <javascript:>
> > <mailto:scitools-iri...@googlegroups.com <javascript:>>.
> > To view this discussion on the web, visit
> >
> https://groups.google.com/d/msgid/scitools-iris/d92a906c-22d6-4fd4-8554-0bb78fcb23a6%40googlegroups.com
> <https://groups.google.com/d/msgid/scitools-iris/d92a906c-22d6-4fd4-8554-0bb78fcb23a6%40googlegroups.com>
>
> >
> <https://groups.google.com/d/msgid/scitools-iris/d92a906c-22d6-4fd4-8554-0bb78fcb23a6%40googlegroups.com?utm_medium=email&utm_source=footer
> <https://groups.google.com/d/msgid/scitools-iris/d92a906c-22d6-4fd4-8554-0bb78fcb23a6%40googlegroups.com?utm_medium=email&utm_source=footer>>.
>
> > For more options, visit https://groups.google.com/d/optout
> <https://groups.google.com/d/optout>.
>
> --
> You received this message because you are subscribed to the Google
> Groups "SciTools (iris, cartopy, cf_units, etc.) -
> https://github.com/scitools" group.
> To unsubscribe from this group and stop receiving emails from it, send
> an email to scitools-iri...@googlegroups.com
> <mailto:scitools-iri...@googlegroups.com>.
> To view this discussion on the web, visit
> https://groups.google.com/d/msgid/scitools-iris/3033803e-78b7-483e-ba5d-4f810584f47b%40googlegroups.com
> <https://groups.google.com/d/msgid/scitools-iris/3033803e-78b7-483e-ba5d-4f810584f47b%40googlegroups.com?utm_medium=email&utm_source=footer>.
Reply all
Reply to author
Forward
0 new messages