converts 64bit IBM floating point numbers to IEEE754 64bit floats

IBMFloat64.jl

import Base.convert

bitstype 64 IBMFloat64

# ibm float format is hex based:
#  * bit 1 = sign
#  * bits 2-8 = exponent to the power of 16, with bias of 64
#  * bits 9-64 = mantissa: no implied leading 1 (unlike ieee),
#    typically normalised (hex of bits 9-12 is non-zero), though this need
#    not be the case
#  * zero is indicated by all bits being 0.
#  * special values arise when mantissa is zero, and sign|exponent is not
#    (therefore, no negative zero, no Inf's): currrently, all converted to NaN
#    can extract by
#      convert(Uint8, reinterpret(Uint64,x) >> 56)

function convert(::Type{Float64}, x::IBMFloat64)
    ibm = reinterpret(Uint64,x)
    if ibm == 0
        return zero(Float64)
    end

    mantissa = ibm & 0x00ff_ffff_ffff_ffff
    if mantissa == 0
        # special value, convert to NaN
        # in future, special value could be stored in NaN payload
        return nan(Float64)
    end

    # normalise mantissa to base 2.
    offset = 1
    while mantissa & 0x0080_0000_0000_0000 == 0
        offset += 1
        mantissa <<= 1
    end
    # drop leading 1 (8th bit), and shift remainder to 12th bit.
    mantissa = (mantissa $ 0x0080_0000_0000_0000) >> 3

    sign = ibm & 0x8000_0000_0000_0000

    exponent = ibm & 0x7f00_0000_0000_0000
    # convert from base 16 to base 2, adjust for normalisation.
    exponent = ((((exponent >> 56) - 64) << 2) - offset + 1023) << 52

    return reinterpret(Float64, sign | exponent | mantissa)
end