1from __future__ import absolute_import, division, print_function
2"""Math utilities by Russell Owen
3
4History:
52001-01-10 ROwen    Removed inf and nan stuff because it fails on Mac OS X
6                    (and it was known to not work on all platforms).
7                    Maybe someday Python will handle IEEE floating point!
82002-03-04 ROwen    Added checkRange.
92002-07-09 ROwen    Added compareFloats, matchLists and flattenList.
102002-08-01 ROwen    Renamed logEQ to logEq for consistency.
112002-08-08 ROwen    Moved to RO and renamed from ROMath.
122002-08-13 ROwen    Bug fix: wrapCtr gave the wrong result.
132003-03-12 ROwen    Added isSequence; rewrote flattenList to use it.
142003-04-01 ROwen    Added rot2D.
152003-05-29 ROwen    Added asList.
162003-06-27 ROwen    Added removeDups.
172003-07-18 ROwen    Added asSequence.
182003-11-18 ROwen    Moved asList, asSequence, flattenLists, isSequence,
19                    removeDups and matchLists to SeqUtil.
202007-04-24 ROwen    Changed Numeric to numpy in a doc string.
212010-07-30 ROwen    Bug fix: wrapPos returned 360 if input was a very small negative value.
222011-07-21 ROwen    API change and bug fix: rThetaFromXY was documented to raise an exception
23                    if too near the pole, but it did not do this reliably, if at all.
24                    Changed to return theta = NaN (numpy.nan) if too near the pole.
252015-09-24 ROwen    Replace "== None" with "is None" to modernize the code.
262015-11-03 ROwen    Replace "!= None" with "is not None" to modernize the code.
27"""
28__all__ = ["sind", "cosd", "tand", "asind", "acosd", "atand", "atan2d", "compareFloats", "checkRange",
29    "nint", "sign", "logEq", "logNE", "rot2D", "rThetaFromXY", "xyFromRTheta", "vecMag", "wrapCtr", "wrapPos"]
30
31import math
32import numpy
33from RO.PhysConst import RadPerDeg
34import RO.SysConst
35
36DegPerRad = 1.0 / RadPerDeg
37_TinyFloat = numpy.finfo(float).tiny
38
39# The following were commented out 2001-01-10 because inf and nan
40# are not handled on Mac OS X (Python 2.2 from fink).
41# the following probably works on all platforms that support IEEE floating point
42# it has been tested on Solaris and MacOS classic.
43# Probably it's simply an issue of whether Python was built with the -ieee option.
44# Inf = inf = float("inf")
45# NaN = nan = float("nan")
46
47def sind(angDeg):
48    """sine of angle, in degrees"""
49    return math.sin(RadPerDeg * angDeg)
50
51def cosd(angDeg):
52    """cosine of angle, in degrees"""
53    return math.cos(RadPerDeg * angDeg)
54
55def tand(angDeg):
56    """tangent of angle, in degrees"""
57    return math.tan(RadPerDeg * angDeg)
58
59def asind(x):
60    """arcsine of x, in degrees"""
61    return DegPerRad * math.asin(x)
62
63def acosd(x):
64    """arccosine of x, in degrees"""
65    return DegPerRad * math.acos(x)
66
67def atand(x):
68    """arctangent of x, in degrees"""
69    return DegPerRad * math.atan(x)
70
71def atan2d(x, y):
72    """arctangent of y/x, in degrees"""
73    return DegPerRad * math.atan2(x, y)
74
75def compareFloats(a, b, rtol=1.0e-5, atol=RO.SysConst.FAccuracy):
76    """Compares values a and b
77    Returns 0 if the values are approximately equals, i.e.:
78    - |a - b| < atol + (rtol * |a + b|)
79    Else 1 if a > b, -1 if a < b
80
81    Inputs:
82    - a, b: scalars to be compared (int or float)
83    - atol: absolute tolerance
84    - rtol: relative tolerance
85
86    The algorithm used is the same one used by numpy.allclose.
87    """
88    if abs(a - b) < (atol + (rtol * abs(float(a + b)))):
89        return 0
90    return cmp(a, b)
91
92def checkRange(value, minValue, maxValue, valDescr="value"):
93    """Checks that value is in range [minValue, maxValue] and raises a ValueError if not.
94    If minValue or maxValue is None, that limit is not checked.
95    If value is None, nothing is checked.
96    """
97    if value is None:
98        return
99    if maxValue is not None and value > maxValue:
100        raise ValueError("%s too large: %r > %r" % (valDescr, value, maxValue))
101    if minValue is not None and value < minValue:
102        raise ValueError("%s too small: %r < %r" % (valDescr, value, minValue))
103
104def nint(x, n=0):
105    """Returns x rounded to the nearast multiple of 10**-n.
106    Values of n > 0 are treated as 0, so that the result is an integer.
107
108    In other words, just like the built in function round,
109    but returns an integer and treats n > 0 as 0.
110
111    Inputs:
112    - x: the value to round
113    - n: negative of power of 10 to which to round (e.g. -2 => round to nearest 100)
114
115    Error Conditions:
116    - raises OverflowError if x is too large to express as an integer (after rounding)
117    """
118    return int (round (x, min(n, 0)) + 0.5)
119
120def sign(x):
121    """Returns -1 if x < 0, 1 otherwise
122
123    Error Conditions:
124    - raises TypeError if x is not a number
125    """
126    abs(x)  # checks type of argument
127    if x < 0:
128        return -1
129    else:
130        return 1
131
132def logEq(a, b):
133    """Returns 1 if the logical value of a equals the logical value of b, 0 otherwise"""
134    return (not a) == (not b)
135
136def logNE(a, b):
137    """Returns 1 if the logical value of does not equal the logical value of b, 0 otherwise"""
138    return (not a) != (not b)
139
140def rot2D(xyVec, angDeg):
141    """Rotates a 2-dimensional vector by a given angle.
142
143    Inputs:
144    - xyVec     x,y vector to be rotated
145    - angDeg    angle (degrees, 0 along x, 90 along y)
146
147    Outputs:
148    rotVec  x,y rotated vector
149
150    Error Conditions:
151    Raises ValueError if:
152    - xyVec is not two numbers
153    - angDeg is not a number
154
155    Details:
156    Changing coordinate systems:
157    Given a point P whose position in coordinate system A is P_A_xy
158    and another coordinate system B whose angle with respect to A is B_A_ang
159    and whose position with respect to A is B_A_xy,
160    then P_B_xy, the position of P in coordinate system B is:
161    P_B_xy = (P_A_xy - B_A_xy) rotated by -B_A_ang
162
163    History:
164    2003-04-01 Converted to Python from TCC cnv_Rot2D
165    """
166    x, y = xyVec
167    sinAng = sind (angDeg)
168    cosAng = cosd (angDeg)
169
170    return (
171        cosAng * x - sinAng * y,
172        sinAng * x + cosAng * y,
173    )
174
175def rThetaFromXY(xy):
176    """Returns the magnitude and angle of a 2-dim vector.
177
178    Inputs:
179    - xy: cartesian coordinates
180
181    Returns:
182    - r: radius
183    - theta: theta (deg); NaN if cannot be reliably computed
184    """
185    r = math.hypot(*xy)
186    if r < _TinyFloat:
187        theta = numpy.nan
188    else:
189        theta = math.atan2(xy[1], xy[0]) / RO.PhysConst.RadPerDeg
190    return (r, theta)
191
192def xyFromRTheta(rTheta):
193    """Returns the x and y components of a polar vector"""
194    r, theta = rTheta
195    return (
196        r * cosd(theta),
197        r * sind(theta),
198    )
199
200def vecMag(a):
201    """Returns the magnitude of vector a"""
202    sumSq = 0
203    for ai in a:
204        sumSq += ai * ai
205    return math.sqrt(sumSq)
206
207def wrapCtr(angDeg):
208    """Returns the angle (in degrees) wrapped into the range (-180, 180]"""
209    # First wrap into [0, 360]; result is 360 if ctrAng < 0 but so near 0 that adding 360 rounds it
210    ctrAng = angDeg % 360.0
211    if ctrAng > 180.0:
212        ctrAng -= 360.0
213    return ctrAng
214
215def wrapPos(angDeg):
216    """Returns the angle (in degrees) wrapped into the range [0, 360)"""
217    res = angDeg % 360.0
218    # First wrap into [0, 360]; result is 360 if ctrAng < 0 but so near 0 that adding 360 rounds it
219    if res == 360.0:
220        return 0.0
221    return res
222