Posits.jl

module Posits
export Posit8, Posit16, Posit32, Posit64, Posit128, AbstractPosit

abstract type AbstractPosit{N, ES} <: Real end

for size in [8, 16, 32, 64, 128]
    Posit = Symbol("Posit", size)
    UInt = Symbol("UInt", size)

    @eval begin
        primitive type $Posit{N, ES} <: AbstractPosit{N, ES}
            $size
        end

        machine_type(::Type{$Posit{N, ES}}) where {N, ES} = $UInt
    end
end

machine_type(x::AbstractPosit) = machine_type(typeof(x))

# The total bit size of the given posit type. Named to match the terminology of
# the posit papers.
nbits(::Type{<:AbstractPosit{N, ES}}) where {N, ES} = N
nbits(x::AbstractPosit) = nbits(typeof(x))

# The number of bits used for the exponent in the given posit type. Named to
# match the terminology of the posit papers.
es(::Type{<:AbstractPosit{N, ES}}) where {N, ES} = ES
es(x::AbstractPosit) = es(typeof(x))

useed(::Type{T}) where T <: AbstractPosit = useed(es(T))
useed(es::Integer) = 2 ^ (2 ^ es)

# Create a mask of the lowest `bits` bits for the given type.
# E.g. `mask(UInt8, 4) == 0x0f`
mask(::Type{T}, bits::Int) where T = (one(T) << bits) - one(T)

function to_bits(x::AbstractPosit)
    # Use a mask to keep only the bits used by the posit `x`. This lets us
    # ignore the unused high bits just in case they are set.
    T = machine_type(x)
    bits = reinterpret(T, x)
    bits & mask(T, nbits(x))
end

Base.bitstring(x::AbstractPosit) =
    string(to_bits(x), base = 2, pad = nbits(x))

function components(x::AbstractPosit)
    # Within this function we will use `size` to refer to the number of bits
    # used by types or posits components.
    T = machine_type(x)
    machine_size = sizeof(T) * 8
    size = nbits(x)
    exp_size = es(x)
    bits = to_bits(x)

    # The number of high bits in the machine integer that aren't used by this
    # posit type.
    unused_size = machine_size - size

    # Extract the sign bit.
    sign = Bool(bits >> (size - 1))

    # Extract the first regime bit and shift the regime to the start of the
    # machine integer so we can count leading zeros or ones.
    first_regime_bit = (bits >> (size - 2)) & one(T)
    regime_leading = bits << (unused_size + 1)

    if first_regime_bit == 1
        regime_count = leading_ones(regime_leading)

        # A string of N ones stands for a regime of N - 1.
        regime = regime_count - 1
    else
        # When counting leading zeros we set the bit after the end of the posit
        # to 1 to deal with the case when the string of zeros extends to the
        # end.
        #
        # The opposite is unecessary when counting leading ones because the bit
        # will be zero automatically.
        regime_leading |= one(T) << (unused_size)
        regime_count = leading_zeros(regime_leading)

        # Finally, a string of N zeros stands for a regime of -N.
        regime = -regime_count
    end

    # The number of bits used by the regime is 1 greater than the unary count,
    # since it includes the opposite bit that ends the run.
    regime_size = regime_count + 1

    # Once we know the regime size, since the sign and exponent sizes are
    # constant, we can finally determine the number of fraction bits.
    frac_size = size - 1 - regime_size - exp_size

    # However, note carefully that this fraction size calculation can turn out
    # negative! Because regime bits can extend to fill the entire posit, they
    # can push all the fraction bits out (making this calculation zero), and
    # even push the exponent bits out (making this calculation negative).
    #
    # Example
    # -------
    # A Posit8{4, 2} value (i.e. a 4-bit posit with 2-bit exponent in
    # an 8-bit machine type):
    #
    # Bits:    0101
    # Meaning: srre
    #
    # s = sign, r = regime, e = exponent, f = fraction (none present)
    #
    # Note that there is only 1 exponent bit present even though this is a
    # 2-bit exponent posit. It is taken as the high bit of the exponent, and
    # the missing low bit(s) are taken to be zero.
    #
    # In this case, it turns out to work just fine for us when extracting the
    # exponent bits to shift right by a possibly-negative amount. If 1 exponent
    # bit has been pushed out by the regime, the apparent fraction size will be
    # -1 and `bits >> -1` will position the exponent in the low bits with the
    # pushed-out zero bits made explicit.
    #
    # Thus, we take this value as our shift amount, but correct the fraction
    # size to be non-negative.
    exp_shift = frac_size
    if frac_size < 0
        frac_size = 0
    end

    exponent = bits >> exp_shift & mask(T, exp_size)
    fraction = bits & mask(T, frac_size)

    (sign = sign, regime_type = Int(first_regime_bit), regime_size = regime_size, regime = regime, exponent = Int(exponent), fraction = Int(fraction))
end

end