8000 Initial UQuantity and tests. · astropy/astropy@113f8d1 · GitHub
[go: up one dir, main page]

Skip to content

Commit 113f8d1

Browse files
committed
Initial UQuantity and tests.
1 parent f5e4856 commit 113f8d1

File tree

5 files changed

+190
-1
lines changed

5 files changed

+190
-1
lines changed

astropy/uncertainty/__init__.py

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,3 +1,5 @@
11
# Licensed under a 3-clause BSD style license - see LICENSE.rst
2+
from __future__ import absolute_import
23

3-
from core import *
4+
from .core import *
5+
from .uquantity import *

astropy/uncertainty/core.py

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,3 +1,4 @@
1+
# -*- coding: utf-8 -*-
12
# Licensed under a 3-clause BSD style license - see LICENSE.rst
23

34
from __future__ import (absolute_import, division, print_function,
@@ -47,12 +48,15 @@
4748
np.sin: np.cos,
4849
np.sinh: np.cosh,
4950
np.sqrt: lambda x: 0.5/np.sqrt(x),
51+
np.square: lambda x: 2.0*x,
5052
np.tan: lambda x: 1+np.tan(x)**2,
5153
np.tanh: lambda x: 1-np.tanh(x)**2}
5254

5355
UFUNC_DERIVATIVES[np.divide] = UFUNC_DERIVATIVES[np.true_divide]
5456
UFUNC_DERIVATIVES[np.abs] = UFUNC_DERIVATIVES[np.fabs]
5557

58+
if hasattr(np, 'cbrt'):
59+
UFUNC_DERIVATIVES[np.cbrt] = lambda x: 1./(3.*np.cbrt(x)**2)
5660

5761

5862
class Uncertainty(object):

astropy/uncertainty/tests/__init__.py

Whitespace-only changes.
Lines changed: 116 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,116 @@
1+
# coding: utf-8
2+
# Licensed under a 3-clause BSD style license - see LICENSE.rst
3+
"""
4+
Test the Quantity class and related.
5+
"""
6+
7+
from __future__ import (absolute_import, unicode_literals, division,
8+
print_function)
9+
10+
import numpy as np
11+
from ...tests.helper import pytest
12+
from ... import units as u
13+
from ... import constants as c
14+
from .. import UQuantity
15+
16+
17+
def test_initialisation():
18+
v1 = UQuantity(5, 2, u.km)
19+
assert v1.value == 5
20+
assert v1.unit == u.km
21+
assert v1.uncertainty == 2
22+
v2 = UQuantity(5 * u.km, 2)
23+
assert v2.value == 5
24+
assert v2.unit == u.km
25+
assert v2.uncertainty == 2
26+
v3 = UQuantity(5 * u.km, 2000 * u.m)
27+
assert v3.value == 5
28+
assert v3.unit == u.km
29+
assert v3.uncertainty == 2
30+
v4 = UQuantity(np.arange(5.), 2., u.km)
31+
assert np.all(v4.value == np.arange(5.))
32+
assert v4.unit == u.km
33+
assert v4.uncertainty == 2
34+
v5 = UQuantity(np.arange(5.), np.array([1., 2., 1., 2., 1.]), u.km)
35+
assert np.all(v5.value == np.arange(5.))
36+
assert v5.unit == u.km
37+
assert np.all(v5.uncertainty == np.array([1., 2., 1., 2., 1.]))
38+
39+
class TestBasics():
40+
def setup(self):
41+
self.v = UQuantity(5., 2., u.km)
42+
self.a = UQuantity(np.arange(1., 5.), 1., u.s)
43+
self.b = UQuantity(np.array([1., 2., 3.]), np.array([0.1, 0.2, 0.1]),
44+
u.m)
45+
def test_addition(self):
46+
c1 = self.v + UQuantity(12, 5, self.v.unit)
47+
assert c1.value == self.v.value + 12
48+
assert c1.unit == u.km
49+
# Uncertainties under addition add in quadrature
50+
assert np.allclose(c1.uncertainty,
51+
np.sqrt(self.v.uncertainty**2 + 5**2))
52+
# now with different units
53+
c2 = self.v + UQuantity(12000., 5000., u.m)
54+
assert c2.value == self.v.value + 12
55+
assert c2.unit == u.km
56+
assert np.allclose(c2.uncertainty,
57+
np.sqrt(self.v.uncertainty**2 + 5**2))
58+
# try array
59+
c3 = self.v + self.b
60+
assert np.all(c3.nominal_value ==
61+
self.v.nominal_value + self.b.nominal_value)
62+
assert np.allclose(c3.uncertainty,
63+
np.sqrt((self.v.uncertainty * self.v.unit)**2 +
64+
(self.b.uncertainty * self.b.unit)**2)
65+
.to(c3.unit).value)
66+
# try adding regular Quantity
67+
c4 = self.v + 10. * self.v.unit
68+
assert c4.value == self.v.value + 10.
69+
assert c4.uncertainty == self.v.uncertainty
70+
71+
def test_subtraction(self):
72+
c1 = self.v - UQuantity(12, 5, self.v.unit)
73+
assert c1.value == self.v.value - 12
74+
assert c1.unit == u.km
75+
# Uncertainties under addition add in quadrature
76+
assert np.allclose(c1.uncertainty,
77+
np.sqrt(self.v.uncertainty**2 + 5**2))
78+
79+
def test_multiplication(self):
80+
c1 = self.v * self.a
81+
assert np.all(c1.nominal_value ==
82+
self.v.nominal_value * self.a.nominal_value)
83+
84+
# Fractional uncertainties under multiplication add in quadrature
85+
assert np.allclose(c1.uncertainty/c1.value,
86+
np.sqrt((self.v.uncertainty/self.v.value)**2 +
87+
(self.a.uncertainty/self.a.value)**2))
88+
# Test multiplication with straight Quantity
89+
c2 = self.a * (10. * u.s)
90+
assert np.all(c2.value == self.a.value * 10.)
91+
assert c2.unit == self.a.unit * u.s
92+
assert np.all(c2.uncertainty == self.a.uncertainty * 10.)
93+
94+
def test_division(self):
95+
c1 = self.v / self.a
96+
assert np.allclose(c1.value, (self.v.nominal_value /
97+
self.a.nominal_value).to(c1.unit).value)
98+
# Fractional uncertainties under division add in quadrature
99+
assert np.allclose(c1.uncertainty/c1.value,
100+
np.sqrt((self.v.uncertainty/self.v.value)**2 +
101+
(self.a.uncertainty/self.a.value)**2))
102+
103+
def test_more_complex():
104+
G = UQuantity(c.G, subok=False)
105+
m1 = UQuantity(1e15, 1e5, u.kg)
106+
m2 = UQuantity(100, 10, u.kg)
107+
r = UQuantity(10000, 500, u.m)
108+
F = G * (m1 * m2) / r**2
109+
assert np.isclose(F.value, 6.67384e-11 * (1e15 * 100< 10000 /span>) / 10000**2)
110+
assert F.unit == u.N
111+
# Uncertainties calculated using partial derivative method
112+
assert np.isclose(F.uncertainty, np.sqrt(
113+
(m1.value*m2.value/(r.value**2)*G.uncertainty)**2 +
114+
(G.value*m2.value/(r.value**2)*m1.uncertainty)**2 +
115+
(G.value*m1.value/(r.value**2)*m2.uncertainty)**2 +
116+
(-2*G.value*m1.value*m2.value/(r.value**3)*r.uncertainty)**2))

astropy/uncertainty/uquantity.py

Lines changed: 67 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,67 @@
1+
# -*- coding: utf-8 -*-
2+
# Licensed under a 3-clause BSD style license - see LICENSE.rst
3+
4+
from __future__ import (absolute_import, division, print_function,
5+
unicode_literals)
6+
7+
import numpy as np
8+
9+
from ..units import Quantity
10+
from .core import Variable
11+
12+
__all__ = ["UQuantity"]
13+
14+
15+
class UQuantity(Variable, Quantity):
16+
17+
def __new__(cls, value, uncertainty=None, unit=None, **kwargs):
18+
if uncertainty is None:
19+
uncertainty = getattr(value, 'uncertainty', None)
20+
value = getattr(value, 'nominal_value', value)
21+
subok = kwargs.pop('subok', True)
22+
if subok and isinstance(value, Quantity):
23+
q_cls = type(value)
24+
cls = _subclasses[q_cls]
25+
else:
26+
q_cls = Quantity
27+
value = q_cls(value, unit, **kwargs)
28+
if uncertainty is not None:
29+
uncertainty = q_cls(uncertainty, value.unit)
30+
return super(UQuantity, cls).__new__(cls, value, uncertainty, **kwargs)
31+
32+
@property
33+
def uncertainty(self):
34+
return (0 if self._uncertainty is None
35+
else self._uncertainty().to(self.unit).value)
36+
37+
def __quantity_subclass__(self, unit):
38+
q_cls = self._nominal_value.__quantity_subclass__(unit)[0]
39+
return _subclasses[q_cls], True
40+
41+
def __str__(self):
42+
return '{0}±{1}{2:s}'.format(self.value, self.uncertainty,
43+
self._unitstr)
44+
45+
def __repr__(self):
46+
prefix1 = '<' + self.__class__.__name__ + ' '
47+
if self.isscalar:
48+
return '{0}{1}±{2}{3:s}>'.format(prefix1, self.value,
49+
self.uncertainty, self._unitstr)
50+
51+
prefix2 = ' ' * (len(prefix1) - 1) + '±'
52+
valstr = np.array2string(self.value, separator=',', prefix=prefix1)
53+
errstr = np.array2string(self.uncertainty, separator=',',
54+
prefix=prefix2)
55+
return '{0}{1}\n{2}{3}{4:s}>'.format(prefix1, valstr,
56+
prefix2, errstr, self._unitstr)
57+
58+
class _SubclassDict(dict):
59+
def __getitem__(self, q_cls):
60+
if q_cls not in self:
61+
# Use str('U') to avoid unicode for class name in python2.
62+
self[q_cls] = type(str('U') + q_cls.__name__,
63+
(UQuantity, Variable, q_cls), {})
64+
return super(_SubclassDict, self).__getitem__(q_cls)
65+
66+
_subclasses = _SubclassDict()
67+
_subclasses[Quantity] = UQuantity

0 commit comments

Comments
 (0)
0