HIVE
/
TriOS-alt
1
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity
TriOS-alt/zubehör/sphinx/spinx100225-ori/sphinx2/fltmath.spn

219 lines
7.4 KiB
Plaintext
Raw Normal View History

2010-11-26 23:58:06 +01:00
''***************************************
''* Floating-Point Math *
''* Single-precision IEEE-754 *
''* Author: Chip Gracey *
''* Copyright (c) 2006 Parallax, Inc. *
''* See end of file for terms of use. *
''***************************************
PUB FFloat(integer) : single | s, x, m
''Convert integer to float
if m := ||integer 'absolutize mantissa, if 0, result 0
s := integer >> 31 'get sign
x := >|m - 1 'get exponent
m <<= 31 - x 'msb-justify mantissa
m >>= 2 'bit29-justify mantissa
return Pack(@s) 'pack result
PUB FRound(single) : integer
''Convert float to rounded integer
return FInteger(single, 1) 'use 1/2 to round
PUB FTrunc(single) : integer
''Convert float to truncated integer
return FInteger(single, 0) 'use 0 to round
PUB FNeg(singleA) : single
''Negate singleA
return singleA ^ $8000_0000 'toggle sign bit
PUB FAbs(singleA) : single
''Absolute singleA
return singleA & $7FFF_FFFF 'clear sign bit
PUB FSqr(singleA) : single | s, x, m, root
''Compute square root of singleA
if singleA > 0 'if a =< 0, result 0
Unpack(@s, singleA) 'unpack input
m >>= !x & 1 'if exponent even, shift mantissa down
x ~>= 1 'get root exponent
root := $4000_0000 'compute square root of mantissa
repeat 31
result |= root
if result ** result > m
result ^= root
root >>= 1
m := result >> 1
return Pack(@s) 'pack result
PUB FAdd(singleA, singleB) : single | sa, xa, ma, sb, xb, mb
''Add singleA and singleB
Unpack(@sa, singleA) 'unpack inputs
Unpack(@sb, singleB)
if sa 'handle mantissa negation
-ma
if sb
-mb
result := ||(xa - xb) <# 31 'get exponent difference
if xa > xb 'shift lower-exponent mantissa down
mb ~>= result
else
ma ~>= result
xa := xb
ma += mb 'add mantissas
sa := ma < 0 'get sign
||ma 'absolutize result
return Pack(@sa) 'pack result
PUB FSub(singleA, singleB) : single
''Subtract singleB from singleA
return FAdd(singleA, FNeg(singleB))
PUB FMul(singleA, singleB) : single | sa, xa, ma, sb, xb, mb
''Multiply singleA by singleB
Unpack(@sa, singleA) 'unpack inputs
Unpack(@sb, singleB)
sa ^= sb 'xor signs
xa += xb 'add exponents
ma := (ma ** mb) << 3 'multiply mantissas and justify
return Pack(@sa) 'pack result
PUB FDiv(singleA, singleB) : single | sa, xa, ma, sb, xb, mb
''Divide singleA by singleB
Unpack(@sa, singleA) 'unpack inputs
Unpack(@sb, singleB)
sa ^= sb 'xor signs
xa -= xb 'subtract exponents
repeat 30 'divide mantissas
result <<= 1
if ma => mb
ma -= mb
result++
ma <<= 1
ma := result
return Pack(@sa) 'pack result
PRI FInteger(a, r) : integer | s, x, m
'Convert float to rounded/truncated integer
Unpack(@s, a) 'unpack input
if x => -1 and x =< 30 'if exponent not -1..30, result 0
m <<= 2 'msb-justify mantissa
m >>= 30 - x 'shift down to 1/2-lsb
m += r 'round (1) or truncate (0)
m >>= 1 'shift down to lsb
if s 'handle negation
-m
return m 'return integer
PRI Unpack(pointer, single) | s, x, m
'Unpack floating-point into (sign, exponent, mantissa) at pointer
s := single >> 31 'unpack sign
x := single << 1 >> 24 'unpack exponent
m := single & $007F_FFFF 'unpack mantissa
if x 'if exponent > 0,
m := m << 6 | $2000_0000 '..bit29-justify mantissa with leading 1
else
result := >|m - 23 'else, determine first 1 in mantissa
x := result '..adjust exponent
m <<= 7 - result '..bit29-justify mantissa
x -= 127 'unbias exponent
longmove(pointer, @s, 3) 'write (s,x,m) structure from locals
PRI Pack(pointer) : single | s, x, m
'Pack floating-point from (sign, exponent, mantissa) at pointer
longmove(@s, pointer, 3) 'get (s,x,m) structure into locals
if m 'if mantissa 0, result 0
result := 33 - >|m 'determine magnitude of mantissa
m <<= result 'msb-justify mantissa without leading 1
x += 3 - result 'adjust exponent
m += $00000100 'round up mantissa by 1/2 lsb
if not m & $FFFFFF00 'if rounding overflow,
x++ '..increment exponent
x := x + 127 #> -23 <# 255 'bias and limit exponent
if x < 1 'if exponent < 1,
m := $8000_0000 + m >> 1 '..replace leading 1
m >>= -x '..shift mantissa down by exponent
x~ '..exponent is now 0
return s << 31 | x << 23 | m >> 9 'pack result
{{
+------------------------------------------------------------------------------------------------------------------------------+
<EFBFBD> TERMS OF USE: MIT License <20>
+------------------------------------------------------------------------------------------------------------------------------<2D>
<EFBFBD>Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation <20>
<EFBFBD>files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, <20>
<EFBFBD>modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software<72>
<EFBFBD>is furnished to do so, subject to the following conditions: <20>
<EFBFBD> <20>
<EFBFBD>The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.<2E>
<EFBFBD> <20>
<EFBFBD>THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE <20>
<EFBFBD>WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR <20>
<EFBFBD>COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, <20>
<EFBFBD>ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. <20>
+------------------------------------------------------------------------------------------------------------------------------+
}}