Small helper function to represent a Float64 exactly

gistfile1.jl

# reprExact is a small helper function that takes a float (currently only Float64)
# and returns an "exact" representation of that float. Invariant:
# x == eval(parse(reprExact(x))), and even
# bits(x) == bits(eval(parse(reprExact(x))))
# Not tested yet: subnormals, signaling NaNs, NaNs with payload
function reprExact(x::Float64)
    w, p, bias = 11, 53, 1023 # 8,24,127 for float32

    b = bits(x)

    Sbits = b[1:1]
    Spm = Sbits == "0" ? "+" : "-"
    Ebits = b[2:2+w-1]
    Tbits = b[2+w:]

    S = parseint(Sbits,2)
    E = parseint(Ebits,2)
    T = parseint(Tbits,2)

    if E == 2^w-1  # all bits set: NaN or Inf
        if T != 0       # a) NaN
            if Tbits[1] == '1'  # quiet NaN
                return string(Spm,"NaN #q ",Tbits)
            else                # silent NaN
                return string(Spm,"NaN #s ",Tbits)
            end
        else            # b) Inf
            return string(Spm,"Inf")
        end
    elseif E == 0   # no bits set: subnormal or zero
        if T != 0       # d) subnormal
        # v = (-1)^S × 2^(emin)×(0+2^(1-p)×T) = (-1)^S x T x 2^(emin+1-p), emin+1-p = 2+bias-p
            return string(Spm,T,"*2.0^(",2+bias-p,")")
        else            # e) zero
            return string(Spm,"0.0")
        end
    else # normal number!
        # c) If E otherwise, then v = (-1)^S × 2^(E-bias)×(1+2^(1-p)×T)
        return string(Spm,"(1+",T,"/2^",p-1,")*2.0^(",E-bias,")")
    end
end