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