def assertFloatIdentical(self, x, y):
"""Fail unless floats x and y are identical, in the sense that:
(1) both x and y are nans, or
(2) both x and y are infinities, with the same sign, or
(3) both x and y are zeros, with the same sign, or
(4) x and y are both finite and nonzero, and x == y
"""
msg = 'floats {!r} and {!r} are not identical'
if math.isnan(x) or math.isnan(y):
if math.isnan(x) and math.isnan(y):
return
elif x == y:
if x != 0.0:
return
# both zero; check that signs match
elif math.copysign(1.0, x) == math.copysign(1.0, y):
return
else:
msg += ': zeros have different signs'
self.fail(msg.format(x, y))
python类copysign()的实例源码
def ventogeostrofico(lat,GradientePressione,R):
if R==0.0 or GradientePressione==0.0:
tws=0
else:
omega=2*math.pi/(24*3600)#rotazione terrestre
f=2*omega*math.sin(lat)#parametro di coriolis
densitaaria=1.2#kgm/mc vale per T=20gradi
R=R*1853.0#raggio di curvatura isobare espressa in m
GradientePressione=GradientePressione*100.0/1853.0#Gradiente espresso in Pa/m
segno=GradientePressione/math.copysign(GradientePressione,1)
a=segno/R
b=-f
c=math.copysign(GradientePressione,1)/densitaaria
discriminante=(b**2-4*a*c)
if discriminante >0:
tws=(-b-discriminante**0.5)/(2*a)
tws=tws*3600.0/1853.0#converte tws in knts
else:
tws=0.0#caso del centro di alta pressione
return tws
def ventogeostroficoCentroAzione(PCentroAzione,LatCentroAzione,LonCentroAzione,extension,lat,lon):
#campo di pressione: P=0.5*(PCentro-1013)*cos(pi/extension*r)+(PCentro+1013)*0.5, per r<=extension, P=1013 per r >extension
#r distanza dal centro d'azione e raggio curvatura isobara
losso=lossodromica(LatCentroAzione,LonCentroAzione,lat,lon)
r=losso[0]
brg=losso[1]
tws=0
twd=0
segno=0
exp=0
if r<extension:
if PCentroAzione>1013:
segno=1
else:
segno=-1
#GradientePressione=-0.5*(PCentroAzione-1013)*math.pi/extension*math.sin(math.pi/extension*r)
GradientePressione=0.5*(PCentroAzione-1013)*math.pi/extension*math.sin(math.pi/extension*r)
GradientePressione=math.copysign(GradientePressione,1)
tws=ventogeostrofico(LatCentroAzione,segno*GradientePressione,r)
if LatCentroAzione>0:
twd=riduci360(brg-segno*math.radians(100))#emisfero Nord
else:
twd=riduci360(brg+segno*math.radians(100))#emisfero Sud
return tws,twd
def interno(poligono,pto):
#poligono orientato in senso orario
somma=0
for i in range(0,len(poligono)-1):
v1=poligono[i]
v2=poligono[i+1]
rlv1=math.atan2((v1[1]-pto[1]),(v1[0]-pto[0]))
if rlv1<0:rlv1=2*math.pi+rlv1
rlv2=math.atan2((v2[1]-pto[1]),(v2[0]-pto[0]))
if rlv2<0:rlv2=2*math.pi+rlv2
angolo=math.copysign(rlv2-rlv1,1)
if angolo>math.pi:
angolo=2*math.pi-angolo
somma=somma+angolo
if somma<1.99*math.pi:
Interno=False
else:
Interno=True
return Interno
def loxodrome(latA, lonA, latB, lonB): # lossodromica
# Rhumb line navigation
# Takes two points on earth and returns:
# the distance and constant (true) bearing one would need to
# follow to reach it.
# Doesn't function near poles:
# but you shouldn't be sailing near the poles anyways!
# when latB = -90: math domain error log(0)
# when latA = -90 [zero divison error]
# when A==B returns (0.0,0.0)
# if latA == latB:
if math.copysign(latA-latB, 1) <= math.radians(0.1/3600.0):
q = math.cos(latA)
else:
Df = math.log(math.tan(latB/2+math.pi/4)
/ math.tan(latA/2+math.pi/4), math.e)
q = (latB-latA) * 1.0/Df
Distance = EARTH_RADIUS * ((latA-latB)**2+(q*(lonA-lonB))**2)**0.5
Heading = math.atan2(-q*(lonB-lonA), (latB-latA))
if Heading < 0:
Heading = Heading + 2.0 * math.pi # atan2:[-pi,pi]
return Distance, Heading
def print_lat(lat): # def stampalat(lat):
# returns a string in the format xxDegrees,xxMinutes, N/S
lat_decimal = math.copysign(math.degrees(lat), 1) # latdecimali
lat_degree = int(lat_decimal) # latgradi
lat_minute = (lat_decimal - lat_degree) * 60 # latprimi
if lat_minute > 59.51:
lat_degree = lat_degree + 1
lat_minute = 0
else:
if lat_minute - int(lat_minute) > 0.51:
lat_minute = int(lat_minute) + 1
else:
lat_minute = int(lat_minute)
if lat > 0:
hemisphere = "N" # segno
else:
hemisphere = "S"
gradi = "%2d" % lat_degree # gradi
primi = "%2d" % lat_minute # primi
lat = (gradi.replace(" ", "0") + u"°" + primi.replace(" ", "0")
+ "'" + hemisphere)
return lat
def print_lon(lon): # def stampalon(lon):
# returns a string in the format xxDegrees,xxMinutes, N/S
lon_decimal = math.copysign(math.degrees(lon), 1) # londecimali
lon_degree = int(lon_decimal) # longradi
lon_minutes = (lon_decimal - lon_degree) * 60 # lonprimi
if lon_minutes > 59.51:
lon_degree = lon_degree + 1
lon_minutes = 0
else:
if lon_minutes - int(lon_minutes) > 0.51:
lon_minutes = int(lon_minutes) + 1
else:
lon_minutes = int(lon_minutes)
if lon > 0:
hemisphere = "W"
else:
hemisphere = "E"
degrees = "%3d" % lon_degree # gradi
minutes = "%2d" % lon_minutes
lon = (degrees.replace(" ", "0") + u"°" + minutes.replace(" ", "0")
+ "'" + hemisphere)
return lon
def ventogeostrofico(lat,GradientePressione,R):
if R==0.0 or GradientePressione==0.0:
tws=0
else:
omega=2*math.pi/(24*3600)#rotazione terrestre
f=2*omega*math.sin(lat)#parametro di coriolis
densitaaria=1.2#kgm/mc vale per T=20gradi
R=R*1853.0#raggio di curvatura isobare espressa in m
GradientePressione=GradientePressione*100.0/1853.0#Gradiente espresso in Pa/m
segno=GradientePressione/math.copysign(GradientePressione,1)
a=segno/R
b=-f
c=math.copysign(GradientePressione,1)/densitaaria
discriminante=(b**2-4*a*c)
if discriminante >0:
tws=(-b-discriminante**0.5)/(2*a)
tws=tws*3600.0/1853.0#converte tws in knts
else:
tws=0.0#caso del centro di alta pressione
return tws
def ventogeostroficoCentroAzione(PCentroAzione,LatCentroAzione,LonCentroAzione,extension,lat,lon):
#campo di pressione: P=0.5*(PCentro-1013)*cos(pi/extension*r)+(PCentro+1013)*0.5, per r<=extension, P=1013 per r >extension
#r distanza dal centro d'azione e raggio curvatura isobara
losso=lossodromica(LatCentroAzione,LonCentroAzione,lat,lon)
r=losso[0]
brg=losso[1]
tws=0
twd=0
segno=0
exp=0
if r<extension:
if PCentroAzione>1013:
segno=1
else:
segno=-1
#GradientePressione=-0.5*(PCentroAzione-1013)*math.pi/extension*math.sin(math.pi/extension*r)
GradientePressione=0.5*(PCentroAzione-1013)*math.pi/extension*math.sin(math.pi/extension*r)
GradientePressione=math.copysign(GradientePressione,1)
tws=ventogeostrofico(LatCentroAzione,segno*GradientePressione,r)
if LatCentroAzione>0:
twd=riduci360(brg-segno*math.radians(100))#emisfero Nord
else:
twd=riduci360(brg+segno*math.radians(100))#emisfero Sud
return tws,twd
def interno(poligono,pto):
#poligono orientato in senso orario
somma=0
for i in range(0,len(poligono)-1):
v1=poligono[i]
v2=poligono[i+1]
rlv1=math.atan2((v1[1]-pto[1]),(v1[0]-pto[0]))
if rlv1<0:rlv1=2*math.pi+rlv1
rlv2=math.atan2((v2[1]-pto[1]),(v2[0]-pto[0]))
if rlv2<0:rlv2=2*math.pi+rlv2
angolo=math.copysign(rlv2-rlv1,1)
if angolo>math.pi:
angolo=2*math.pi-angolo
somma=somma+angolo
if somma<1.99*math.pi:
Interno=False
else:
Interno=True
return Interno
def calcolaventogeostrofico(lat,GradientePressione,R):
omega=2*math.pi/(24*3600)#rotazione terrestre
f=2*omega*math.sin(lat)#parametro di coriolis
densitaaria=1.2#kgm/mc vale per T=20gradi
R=R*1853.0#curvatura isobare espressa in m
GradientePressione=GradientePressione*100.0/1853.0#Gradiente espresso in Pa/m
segno=GradientePressione/math.copysign(GradientePressione,1)
a=segno/R
b=-f
c=math.copysign(GradientePressione,1)/densitaaria
discriminante=(b**2-4*a*c)
if discriminante >0:
tws=(-b-discriminante**0.5)/(2*a)
tws=tws*3600.0/1853.0#converte tws in knts
else:
tws=0.0#caso del centro di alta pressione
return tws
def interno(poligono,pto):
#poligono orientato in senso orario
somma=0
for i in range(0,len(poligono)-1):
v1=poligono[i]
v2=poligono[i+1]
rlv1=math.atan2((v1[1]-pto[1]),(v1[0]-pto[0]))
if rlv1<0:rlv1=2*math.pi+rlv1
rlv2=math.atan2((v2[1]-pto[1]),(v2[0]-pto[0]))
if rlv2<0:rlv2=2*math.pi+rlv2
angolo=math.copysign(rlv2-rlv1,1)
if angolo>math.pi:
angolo=2*math.pi-angolo
somma=somma+angolo
if somma<2*math.pi:
Interno=False
else:
Interno=True
return Interno
def get_cubic_root(self):
# We have the equation x^2 D^2 + (1-x)^4 * C / h_min^2
# where x = sqrt(mu).
# We substitute x, which is sqrt(mu), with x = y + 1.
# It gives y^3 + py = q
# where p = (D^2 h_min^2)/(2*C) and q = -p.
# We use the Vieta's substution to compute the root.
# There is only one real solution y (which is in [0, 1] ).
# http://mathworld.wolfram.com/VietasSubstitution.html
# eps in the numerator is to prevent momentum = 1 in case of zero gradient
if np.isnan(self._dist_to_opt) or np.isnan(self._h_min) or np.isnan(self._grad_var) \
or np.isinf(self._dist_to_opt) or np.isinf(self._h_min) or np.isinf(self._grad_var):
logging.warning("Input to cubic solver has invalid nan/inf value!")
raise Exception("Input to cubic solver has invalid nan/inf value!")
p = (self._dist_to_opt + eps)**2 * (self._h_min + eps)**2 / 2 / (self._grad_var + eps)
w3 = (-math.sqrt(p**2 + 4.0 / 27.0 * p**3) - p) / 2.0
w = math.copysign(1.0, w3) * math.pow(math.fabs(w3), 1.0/3.0)
y = w - p / 3.0 / (w + eps)
x = y + 1
if self._verbose:
logging.debug("p %f, denominator %f", p, self._grad_var + eps)
logging.debug("w3 %f ", w3)
logging.debug("y %f, denominator %f", y, w + eps)
if np.isnan(x) or np.isinf(x):
logging.warning("Output from cubic is invalid nan/inf value!")
raise Exception("Output from cubic is invalid nan/inf value!")
return x
def transact_stub(slippage, commission, event, open_orders):
"""
This is intended to be wrapped in a partial, so that the
slippage and commission models can be enclosed.
"""
for order, transaction in slippage(event, open_orders):
if transaction and transaction.amount != 0:
direction = math.copysign(1, transaction.amount)
per_share, total_commission = commission.calculate(transaction)
transaction.price += per_share * direction
transaction.commission = total_commission
yield order, transaction
def update(self, txn):
if self.sid != txn.sid:
raise Exception('updating position with txn for a '
'different sid')
total_shares = self.amount + txn.amount
if total_shares == 0:
self.cost_basis = 0.0
else:
prev_direction = copysign(1, self.amount)
txn_direction = copysign(1, txn.amount)
if prev_direction != txn_direction:
# we're covering a short or closing a position
if abs(txn.amount) > abs(self.amount):
# we've closed the position and gone short
# or covered the short position and gone long
self.cost_basis = txn.price
else:
prev_cost = self.cost_basis * self.amount
txn_cost = txn.amount * txn.price
total_cost = prev_cost + txn_cost
self.cost_basis = total_cost / total_shares
# Update the last sale price if txn is
# best data we have so far
if self.last_sale_date is None or txn.dt > self.last_sale_date:
self.last_sale_price = txn.price
self.last_sale_date = txn.dt
self.amount = total_shares
def __init__(self, dt, sid, amount, stop=None, limit=None, filled=0,
commission=None, id=None):
"""
@dt - datetime.datetime that the order was placed
@sid - stock sid of the order
@amount - the number of shares to buy/sell
a positive sign indicates a buy
a negative sign indicates a sell
@filled - how many shares of the order have been filled so far
"""
# get a string representation of the uuid.
self.id = id or self.make_id()
self.dt = dt
self.reason = None
self.created = dt
self.sid = sid
self.amount = amount
self.filled = filled
self.commission = commission
self._status = ORDER_STATUS.OPEN
self.stop = stop
self.limit = limit
self.stop_reached = False
self.limit_reached = False
self.direction = math.copysign(1, self.amount)
self.type = zp.DATASOURCE_TYPE.ORDER
def from_float(cls, f):
"""Converts a float to a decimal number, exactly.
Note that Decimal.from_float(0.1) is not the same as Decimal('0.1').
Since 0.1 is not exactly representable in binary floating point, the
value is stored as the nearest representable value which is
0x1.999999999999ap-4. The exact equivalent of the value in decimal
is 0.1000000000000000055511151231257827021181583404541015625.
>>> Decimal.from_float(0.1)
Decimal('0.1000000000000000055511151231257827021181583404541015625')
>>> Decimal.from_float(float('nan'))
Decimal('NaN')
>>> Decimal.from_float(float('inf'))
Decimal('Infinity')
>>> Decimal.from_float(-float('inf'))
Decimal('-Infinity')
>>> Decimal.from_float(-0.0)
Decimal('-0')
"""
if isinstance(f, (int, long)): # handle integer inputs
return cls(f)
if _math.isinf(f) or _math.isnan(f): # raises TypeError if not a float
return cls(repr(f))
if _math.copysign(1.0, f) == 1.0:
sign = 0
else:
sign = 1
n, d = abs(f).as_integer_ratio()
k = d.bit_length() - 1
result = _dec_from_triple(sign, str(n*5**k), -k)
if cls is Decimal:
return result
else:
return cls(result)
def __update(self):
area = self.area
if not area:
self.__zero()
return
# Center of mass
# https://en.wikipedia.org/wiki/Center_of_mass#A_continuous_volume
self.meanX = meanX = self.momentX / area
self.meanY = meanY = self.momentY / area
# Var(X) = E[X^2] - E[X]^2
self.varianceX = varianceX = self.momentXX / area - meanX**2
self.varianceY = varianceY = self.momentYY / area - meanY**2
self.stddevX = stddevX = math.copysign(abs(varianceX)**.5, varianceX)
self.stddevY = stddevY = math.copysign(abs(varianceY)**.5, varianceY)
# Covariance(X,Y) = ( E[X.Y] - E[X]E[Y] )
self.covariance = covariance = self.momentXY / area - meanX*meanY
# Correlation(X,Y) = Covariance(X,Y) / ( stddev(X) * stddev(Y) )
# https://en.wikipedia.org/wiki/Pearson_product-moment_correlation_coefficient
correlation = covariance / (stddevX * stddevY)
self.correlation = correlation if abs(correlation) > 1e-3 else 0
slant = covariance / varianceY
self.slant = slant if abs(slant) > 1e-3 else 0
def copysign(x, y):
"""return x with the sign of y (missing from python 2.5.2)"""
if sys.version_info > (2, 6):
return math.copysign(x, y)
return math.fabs(x) * (-1 if y < 0 or (y == 0 and 1/y < 0) else 1)
def f0(cls, x):
return math.sqrt(abs(x))*math.copysign(1,x)*math.sin(1/(abs(x)+pow(10, -9)))
