2011/10/20

バイナリダンプから浮動小数点数を概算

今回は浮動小数点数のfloat型とdouble型のバイナリエディタからの見え方、つまりは16進法表示についてです。 16進数の並びから、暗算でどこまで求められるんだろかーと何となく思ったので考えてみました。

このサイトwikipedia浮動小数点数を参考にしています。



この記事はfloat型についてです。double型もほとんど同じですのでリンク先を参考にしてみてください。
具体的なことに入る前に、Intelのアーキテクチャではデータはリトルエンディアンで扱われます。
だから、バイナリエディタから見た場合はまず頭の中でエンディアンを直す必要があるのです。

下のほうにあるスクリプトのfloat_little関数に適当なfloat型の数値を渡すとリトルエンディアンの並びの16進表示が返されます。
float_big関数はその逆でビッグエンディアンで返します。


ここによれば
float の表す値 = (-1)符号部 × 2指数部-127 × 1.仮数部
だそうです。
要するに、
def float_analyze(f):
    print u"符号部", hex((f2i(f) & int("10000000000000000000000000000000", 2)) >> 31)
    print u"指数部", hex(((f2i(f) & int("01111111100000000000000000000000", 2)) >> 23) -127)
    print u"仮数部", hex(f2i(f) & int("00000000011111111111111111111111", 2))
こういう事ですね。
f2i関数は、単純にANDをつかってビットマスクができるように、一度バイナリにパックしたのをアンパックすることで、float型の浮動小数点数をバイナリ表記は同じであるint型の符号なし整数に変換しています
f2i = lambda i:struct.unpack("I", struct.pack("f", i))[0]

まず、2進法で表した最上位ビットが立っていれば符号はマイナスということなので、
すなわち一番上の桁が0x8以上ならマイナス、0x8未満ならプラスとなります。

次は指数部についてです。
16進法で先頭から3桁をみて、3桁目は8未満なら0にして、8以上なら0x8に変えます。
この時、符号がマイナスだったなら一番上の桁から0x8を引いておきます。
そしたらこの3桁を0x8で割るんですが、
0x8 = 0x10 / 0x2
なので0x2をかけて、0x10で割るすなわち桁を右に一つずらす というのと同じです。
指数部はこれでいいのですが、2の指数部-127乗があとで必要になるということなので127を引いた数だけ求めておきます。
0x7F(127) = 0x80 - 0x1
なので、さっき求めた値に0x1を足して0x80を引きます。

なんか簡単に計算できるとか言っておきながら面倒くさいったらありゃしないですねー
そこで、じゃじゃーん!
いまの指数部-127を求めるまでの行程を簡単にしてみました。
まずは上から2桁をみて、最上桁は同じように、0x8以上なら0x8を引いておきます。
その数に0x2をかけます。
今度は上から3桁目の数を見て、0x8未満なら0x1、0x8以上なら0x2を、0x2をかけた数に足します。
最後に0x80を引きます。
これならまださっきよりかはマシになった・・・のか?
1とか2の足し算ならまだ暗算できますけど、0x8を引くってだけでもなかなか難しいですよねぇ


これで符号部と指数部が分かりました。 仮数部の計算は、計算続けるほど面倒だけど精度が上がってく感じなので、仮数部は計算しなくても、元の数の大体の範囲は分かります。
次のプログラムは仮数部を計算していくと精度の上がって行く様子を表示できます。
#!/usr/bin/env python
#-*- coding:utf-8 -*-

import struct

def toBinary(decimal, bit):
    if decimal == 0:
        return '0'*bit    
    bin_str = ""    
    while decimal > 0:
        bin_str += str(decimal % 2)
        decimal >>= 1
        if len(bin_str) < bit:
            return "0"*(bit-len(bin_str)) + bin_str[::-1]
        elif len(bin_str) == bit:
            return bin_str[::-1]
        else:
            raise ValueError, "Parameter should be under %sbit!" % bit
        
float_little = lambda x: "".join([
        hex(ord(struct.pack("f", x)[i]))[1:].replace("x", " ") 
        for i in xrange(4)
        ]).strip()
float_big = lambda x: "".join([
        hex(ord(struct.pack(">f", x)[i]))[1:].replace("x", " ") 
        for i in xrange(4)
        ]).strip()
f2i = lambda i:struct.unpack("I", struct.pack("f", i))[0]

def HexToFloat(hex_int):    
    print "[-- HexToFloat --]"    
    mask_sign = 0x80000000 # 0b10000000000000000000000000000000    
    mask_exponent = 0x7f800000 # 0b0111/1111/1000/0000/0000/0000/0000/0000    
    mask_fraction = 0x7fffff # 0b0000/0000/0111/1111/1111/1111/1111/1111

    print "Dump(little endian):", float_little(hex_int)
    print "Dump(big endian):", float_big(hex_int)

    sign = (hex_int & mask_sign) >> 31    
    print "Sign(bin):", toBinary(sign, 1)
    exponent = int(((hex_int & mask_exponent) >> 23) - 127)
    print "Exponent(bin):", toBinary(exponent, 8), "(%s = %d)" % (hex(exponent), exponent)    
    fraction = toBinary(hex_int & mask_fraction, 23)
    print "Fraction(bin):", fraction, "(%s)" % hex(hex_int & mask_fraction)

    approx_base = (-1) ** sign * 2 ** exponent    
    print "No fraction:", approx_base, "<= f <", approx_base*2

    fraction_list = [1.]    
    for i in xrange(23): 
        if int(fraction[i]):
            fraction_list.append(1. / 2**(1+i))
        else:
            fraction_list.append(0)

        print "Approximation(fraction:1/2^%d):" % (1+i), \
            approx_base * sum(fraction_list), \
            "<= f <", \
            approx_base * (sum(fraction_list) + 1. / 2**(1+i))

    print
    result = approx_base * sum(fraction_list)
    return result

if __name__ == "__main__":
    import sys
    if len(sys.argv) < 2:
        print "Error!\n"
        print "Usage: %s [floating-point number]\n" % sys.argv[0]
    else:
        HexToFloat(f2i(float(sys.argv[1])))
計算を進めるごとに、引数に渡した浮動小数点数に近似していきます。
"No fraction:"は仮数部を計算しないで指数部までだけの結果です。
何回か試せば分かると思いますが、仮数部なしだと2のn乗とn+1乗の間にあるということしかわかりません。
つまり、整数部分が大きいほど範囲が広がってあまり使えません。 
では、仮数部はどうやって求めていくかいうと、
下から23ビットのデータを十進法の小数点数に直していきます。
23ビットということなので下から6桁で、一番上の桁が0x8より大きいならそこから0x8を引いた部分です。
その一番上の桁を見た時、4,5,6,7のいずれかなら1/2 (0.5)を、
 2,3,5,7のいずれかなら1/4 (0.25)を、
 1.3.5.7のいずれかなら1/8 (0.125)を 仮数部に足します。
これはつまり上から何ビット目がたっているかをみていて、一番上のビットから1/2の位、1/4の位である、と考えることで分かると思います。
仮数部の6桁の上から二番目についても
 8,9,A,B,C,D,E,F -> 1/16 (0.0625)
 4,5,6,7,C,D,E,F -> 1/32 (0.03125)
 2,3,6,7,A,B,E,F -> 1/64 (0.015625)
 1,3,5,7,9,B,D,F -> 1/128 (0.0078125)
といった感じで計算できますがさすがにここまでする必要はないでしょう。

では、練習してみますかw
ちなみに全てリトルエンディアンで表記してあります。
1. [C3 F5 48 40]   答え:3.14
2. [5B 31 86 4E]   答え:152.678192
3. [C7 31 41 4F]   答え:-11.1111
答えは反転して確認してください。 
どうでしょうか? 16進法の繰り上がりとかに慣れたら・・・まぁいけるかな、という感じですかね

結論: 16進法と友達になる気がないなら、誰かが作った便利なツールを使いましょう。

0 件のコメント:

コメントを投稿