Algorithms with Python / レンジコーダ (range coder) [1]

はじめに

拙作のページシャノン・ファノ符号とハフマン符号で説明したハフマン符号は、各記号の出現確率を調べて頻繁に現れる記号は短いビットで表し、あまり現れない記号は長いビットで表すことで、データを圧縮する古典的な方法です。古典的とはいっても、ほかのアルゴリズムと簡単に組み合わせることができるため、ハフマン符号は今でも現役のアルゴリズムです。

記号の出現確率だけを利用してデータを圧縮する方法は、ハフマン符号のほかに「算術符号」という方法があります。一般に、算術符号はハフマン符号よりも性能が良いのですが、実装方法が難しくて実行速度がハフマン符号よりも遅く、なおかつ特許の問題もあって、実際に使われることはあまりありませんでした。

ところが、近年になって「レンジコーダ (Range Coder)」という方法が注目されるようになりました。レンジコーダは原理的には算術符号と同じ方法ですが、参考文献, URL 2 によると『(おそらく) 特許フリー』とのことで、性能は算術符号に比べるとわずかに劣化しますが、実装方法はとても簡単で実行速度も算術符号より高速です。もちろん、ハフマン符号よりも高性能です。

今回はレンジコーダについて詳しく説明します。そして、実際にファイルを圧縮してみましょう。

●算術符号の基本的な考え方

最初に算術符号の基本的な考え方について説明します。算術符号は記号 (文字) 列全体を一つの符号語にする方法で、1960 年代に P. Elias 氏によって提案されました。

●算術符号の符号化

算術符号は、記号列を実数 0 と 1 の間の区間を用いて表します。たとえば、記号は {a, b, c} の 3 種類があり、出現確率はそれぞれ 0.2, 0.6, 0.2 とします。算術符号は、区間を記号の出現確率に比例した小区間に分割していくことで符号化を行います。それでは、記号列 abbbc を符号化してみましょう。次の図を見てください。

x 以上 y 未満の区間を [x, y) と表すことにします。区間の初期値は [0, 1) です。記号を読み込んだら区間 [0, 1) を分割します。記号が a ならば区間の 0 から 0.2 までの部分、b ならば 0.2 から 0.8 までの部分、c ならば 0.8 から 1.0 までの部分に分割します。

最初の記号は a なので区間は [0, 0.2) となります。次の記号は b なので、区間 [0, 0.2) の 0.2 から 0.8 までの部分 [0.04, 0.16) が新しい区間となります。このように、記号を読み込むたびに区間を分割していくと、記号列 abbbc を表す区間は [0.11296, 0.1216) となります。そして、実際の算術符号の符号語は、この区間に含まれる一つの実数を指定します。

ここで符号語を 2 進数で表して、区間内で小数点以下のビット数の少ない値を選ぶことにしましょう。たとえば、0.1171875 を 2 進数で表すと次のようになります。

0.1171875 = 1/16 + 1/32 + 1/64 + 1/128 = (0.0001111)₂

0001111 の 7 ビットを符号語として出力すれば、記号列 abbbc の 5 文字を 7 ビットに圧縮することができるわけです。この例では 1 文字が 1.4 ビットに圧縮されましたが、記号の出現確率によっては 1 文字が 1 ビット未満に圧縮できる場合もあります。ちなみに、ハフマン符号では 1 文字が 1 ビット未満になることはありえません。

●算術符号の復号

次は復号を説明します。ここでは説明の都合上、符号語は下限値の 0.11296 とします。0.11296 は [0, 0.2) の間にあるので、最初の記号は a であることがわかります。次に、a を表す区間 [0, 0.2) を [0, 1.0) になるように拡大すると、符号語は次のように変換できます。

新しい符号語 = (符号語 - 記号の下限値) / 記号の区間幅
             = (0.11296 - 0) / 0.2
             = 0.5648

新しい符号語 0.5648 は [0.2, 0.8) の間にあるので、次の記号は b であることがわかります。このような操作を繰り返し行うことで、次のように記号列 abbbc を求めることができます。

表：復号の過程
符号語	記号	区間
0.11296	a	[0, 0.2)
0.5648	b	[0.2, 0.8)
0.608	b	[0.2, 0.8)
0.68	b	[0.2, 0.8)
0.8	c	[0.8, 1)
0

●算術符号の問題点

ところで、算術符号には 2 つの問題点があります。一つは記号列の最後を判定できないことです。さきほどの復号の例では最後に符号語が 0 になりましたが、このあとも記号 a を復号することができます。つまり、符号語 0.11296 は記号列 abbbc だけではなく、abbbca, abbbcaa, abbcaaa... などにも復号することができるのです。この問題は、終端を表す記号を用意して終端記号を復号したら終了する、または、記号の総数をファイルの先頭に書き込んでおく、などといった方法で解決することができます。

もう一つは、入力する記号列が長くなるほど、より多くの桁数が必要になることです。また、浮動小数点演算の誤差も考慮しなければいけません。これはとても大きな問題点で、解決するまでに 10 年以上の時間がかかりました。

1981 年に C. B. Jones 氏によって発表された算術符号 (Jones 符号) は、実数のかわりに整数で演算するように工夫されています。Jones 符号に興味のある方は参考文献, URL 1 をお読みください。また、『Ｃ言語による最新アルゴリズム事典』 (参考文献 2) にも算術圧縮のプログラム (Ｃ言語) があります。

●レンジコーダの基本的な考え方

レンジコーダは 1998 年に Michael Schindler 氏が発表し、「高性能、高速、特許フリー」の方法として注目を集めるようになりました。Michael Schindler 氏のレンジコーダは計算の途中で「桁上がり」が発生しますが、ロシアの Dmitry Subbotin 氏が発表した「桁上げのないレンジコーダ」は、その名のごとく桁上がりが発生しません。現在、レンジコーダは主にこの 2 種類の形式が存在するようです。

それでは、レンジコーダの基本的な考え方について説明しましょう。なお、ここで説明するレンジコーダは「桁上がり」が発生するバージョンです。

算術符号は区間 [0, 1) を分割していきますが、レンジコーダは [0, 1) を分割するのではなく、最初に大きな区間たとえば [0, 1000) を設定して、それを小さな区間に分割していくことで符号化を行います。レンジコーダは整数で演算するので、記号列が長くなると当然ですが区間が狭くなって分割できなくなります。そのときは区間を引き伸ばすことで対応します。

たとえば、[0, 1000) を分割していくと [123, 124) になりました。もうこれ以上分割できませんね。そこで、区間をたとえば 100 倍して [12300, 12400) を分割することにします。このとき、区間全体の大きさは [0, 1000) ではなく、それを 1000 倍した [0, 100000) と考えるわけです。

単純に考えると、区間を表すために多倍長整数が必要になりますが、区間を引き伸ばすタイミングを定めることにより、通常の整数演算でレンジコーダをプログラムすることができます。また、区間全体の大きさも覚えておく必要はありません。レンジコーダは分割した区間の幅 (range) と下限値 (low) だけで符号化することができます。復号の処理でも、符号化と同じタイミングで range を引き伸ばしていくことで、符号語を記号列に復号することができます。レンジコーダは区間の下限値を符号語として出力します。

●レンジコーダの符号化

それでは具体的に説明しましょう。ここでは説明の都合上、区間の幅 range を 0x1000000 に設定します。実際には range を 0xffffffff (32 ビット) に設定する場合が多いでしょう。下限値 low は 0 に初期化します。区間は [low, low + range) と表すことができるので、最初の区間は [0, 0x1000000) となります。range の初期値が 0x1000000 なので、low の値は 0 から 0xffffff までの範囲 (24 ビット) になります。

記号の出現確率により区間を分割するところは算術符号と同じです。レンジコーダは range が一定の値より小さくなった時点で、range を引き伸ばすところがポイントです。レンジコーダでは、range が初期値の 1/256 (0x10000) より小さくなったら 256 倍します。これは処理をバイト単位で行うための工夫です。次の例を見てください。

[0x123456, 0x123456 + 0xabcd) = 256 倍 => [0x12345600, 0x12345600 + 0xabcd00)

いま low の値が 0x123456 で range の値が 0xabcd だとします。0xabcd < 0x10000 なので range を 256 倍します。このとき、low の値もいっしょに 256 倍すればいいわけです。これで区間を正しく表すことができますが、このままでは low の値がどんどん大きくなる一方ですね。そこで、low の値を一定の範囲内 (24 ビット) に収めることを考えます。

range の値は 24 ビットの範囲内に収まるので、low の計算は 24 ビットの足し算になります。桁上がりの処理を工夫すれば、low を 24 ビットで保持することが可能です。たとえば、次のように low の上位 8 ビット (0x12) をバッファへ出力します。

[0x12345600, 0x12345600 + 0xabcd00) => [0x345600, 0x345600 + 0xabcd00)
low (0x12345600) の上位 8 ビット (0x12) をバッファへ出力 => (0x12)

桁上がりはバッファを使うことで簡単に対応することができます。この処理はプログラムを作るときに詳しく説明します。また、桁上がりが発生しないように工夫することができれば、上位 8 ビット (0x12) をそのまま符号語としてファイルなどへ出力することができます。あとは、記号を読み込んで区間の分割と引き伸ばしを繰り返して、最後に low の値 (24 ビット) を出力します。

簡単な例を示しましょう。記号列 "dcbbaaaa" を符号化します。記号の出現確率は次のようになります。

表：記号の出現確率
	a	b	c	d
出現確率	1/2	1/4	1/8	1/8
下限値	0	4/8	6/8	7/8
上限値	4/8	6/8	7/8	8/8

符号化の過程は次のようになります。

記号を読み込むたびに、range の値は小さくなり low の値は増えていきます。d, c, b, b まで記号を読み込むと、range は 0x4000 になり 0x10000 より小さくなります。ここで range と low を 256 倍して、low の上位 8 ビット (0xfa) を出力します。次に記号 a を読み込みます。range の値は小さくなりますが、a の下限値が 0 なので low の値は増えません。最後に low の値を出力して終了です。符号語は [0xfa, 0x80, 0, 0] になります。

●レンジコーダの復号

次は復号について説明します。下限値 low と幅 range は符号化と同様に 0 と 0x1000000 に初期化します。符号語を code とすると、最初 low は 0 なので [0, range) の範囲で code に対応する記号を探すことになります。見つけた記号を c1 とすると、low と range の値は符号化と同様に次式で更新します。

low1   = low (0) + (range * 記号 c1 の下限値)
range1 = range * 記号 c1 の出現確率

こんどは [low1, low1 + range1) の範囲で code に対応する記号を探します。ここで code から下限値の増分を引き算した値 code1 を求めてみます。すると、次の図に示すように code1 は区間 [0, range1) の符号語に対応していることがわかります。つまり、次は [0, range1) の範囲で code1 に対応する記号を探せばよいのです。

このように、符号語 code から下限値の増分を引き算することで、区間を [low, low + range) から [0, range) に変換することができるわけです。したがって、復号処理では下限値 low の値を覚えておく必要はありません。

range が 0x10000 より小さくなったら range を 256 倍するのは符号化と同じです。このとき符号語 code も 256 倍して、新しい符号語を 1 バイト読み込んで code に加算します。これで符号語を復号することができます。

それでは、復号の過程を具体的に説明しましょう。次の図を見てください。

最初に range と code を初期化します。code の範囲は 24 ビットなので、3 バイト読み込んで 0xfa8000 に初期化します。次に、「記号の下限値 <= code / range < 記号の上限値」を満たす記号を探します。この場合、記号は d になります。そして、range を記号 d の出現確率で縮小して、code から (range * d の下限値) を引き算します。今度は 0x200000 の幅の中で 0x1a8000 に相当する記号を探すわけです。

d, c, b, b まで復号すると、range は 0x4000 になり 0x10000 より小さくなります。ここで range と code を 256 倍して、新しい符号語を 1 バイト読み込んで code に足し算します。この場合、符号語は 0 なので code の値は増えません。あとは、同じ処理を繰り返して記号列 "dcbbaaaa" を求めることができます。

●レンジコーダの実装

それではプログラムを作りましょう。最初に、レンジコーダを表すクラスを定義します。

リスト : レンジコーダ

from byteio import *

# 定数
ENCODE = "encode"
DECODE = "decode"
MAX_RANGE = 0xffffffff
MIN_RANGE = 0x1000000
MASK      = 0xffffffff
SHIFT     = 24

# レンジコーダー
class RangeCoder:
    def __init__(self, s, mode):
        self.stream = s
        self.range_ = MAX_RANGE
        self.buff = 0
        self.cnt = 0
        if mode == ENCODE:
            self.low = 0
        elif mode == DECODE:
            # buff の初期値 (0) を読み捨てる
            self.stream.getc()
            # 4 byte read
            self.low = self.stream.getc()
            self.low = (self.low << 8) + self.stream.getc()
            self.low = (self.low << 8) + self.stream.getc()
            self.low = (self.low << 8) + self.stream.getc()
        else:
            raise "RangeCoder mode error"

クラス名は RangeCoder としました。RangeCoder() の引数 s が入出力ファイル (ByteIO のオブジェクト)、mode はモードを表していて、ENCODE が符号化で DECODE が復号です。インスタンス変数 range_ は幅、buff と cnt が桁上がり用のバッファ、low は符号化のときは下限値で復号のときは符号語 code を表します。

今回は 32 bit 無符号整数の範囲で low と range_ を表すことにします。range_ は MAX_RANGE (0xffffffff) で初期化します。low は符号化するときは 0 で初期化します。復号するときはファイルから 4 バイト分読み込んで初期化します。符号化のとき、最初に出力される 1 バイトは buff の初期値 0 なので、復号のときは 1 バイト読み捨てることに注意してください。

なお、符号化の処理で一番最初に符号語を出力するとき、buff の値を出力しないようにすれば、復号の初期化で 1 バイト読み捨てる処理は不要になります。また、圧縮後のファイルサイズも 1 バイト少なくなります。興味のある方はプログラムを改造してみてください。

RangeCoder では、次に示すメソッドを用意します。

encode_normalize()
符号化で正規化を行う
decode_normalize()
復号で正規化を行う
finish()
符号化の終了

記号の符号化・復号を行うメソッドは、RangeCoder ではなく出現頻度表を表すクラス Freq のメソッドとして定義します。レンジコーダの場合、出現頻度表の作り方によって、いろいろなモデルを実現することができます。今回はレンジコーダの核となる処理だけをクラス RangeCoder で定義します。

次は記号の出現頻度表を表すクラス Freq を作成します。

リスト : 記号の出現頻度表

class Freq:
    def __init__(self, count):
        size = len(count)
        self.size = size
        m = max(count)
        # 2 バイトに収める
        if m > 0xffff:
            self.count = [0] * size
            n = 0
            while m > 0xffff:
                m >>= 1
                n += 1
            for x in range(size):
                if count[x] != 0:
                    self.count[x] = (count[x] >> n) | 1
        else:
            self.count = count[:]
        self.count_sum = [0] * (size + 1)
        # 累積度数表
        for x in range(size):
            self.count_sum[x + 1] = self.count_sum[x] + self.count[x]

Freq() の引数 count は記号の出現頻度を表す配列です。レンジコーダでファイルを圧縮する場合、ハフマン符号と同様に記号の出現頻度表をファイルに付加する必要があります。記号の出現頻度を 1 バイトで表すと、ファイルに付加するデータは 256 バイトですみますが、これではレンジコーダの力を十分に発揮させることはできません。そこで、今回は記号の出現頻度を 2 バイト (0 - 0xffff) で表すことにします。

出現頻度を 2 バイトに丸める処理は簡単です。count の最大値を求め変数 m にセットします。m が 0xffff 以下になるまで m を右へシフトしていき、その回数を変数 n にセットします。あとは count の各要素を n ビット右へシフトするだけです。このとき、出現頻度が 0 にならないように最下位ビットを 1 にしています。

最後に累積度数表 count_sum を作ります。レンジコーダは整数で演算するので、出現確率は各記号の個数と記号の総数から求めます。また、記号の下限値と上限値は分数ではなく累積度数で表します。この場合、累積度数の合計値は count_sum[size] になります。記号 c の出現確率は、count[c] / count_sum[size] で求めることができます。

●符号化のプログラム

次は符号化を行うメソッド encode() を作ります。

リスト : 符号化

    def encode(self, rc, c):
        temp = rc.range_ // self.count_sum[self.size]
        rc.low += self.count_sum[c] * temp
        rc.range_ = self.count[c] * temp
        rc.encode_normalize()

encode() の引数 rc は RangeCoder のオブジェクト、c が符号化する記号です。区間の幅 range_ を狭めて下限値 low の値を計算します。range_ は記号の出現確率で縮小すればいいので、range_ * count[c] / count_sum[size] となります。low の増分は区間の下限値なので、range_ * count_sum[c] / count_sum[size] となります。

プログラムでは、あらかじめ range_ / count_sum[size] を計算して変数 temp にセットし、その値を使って range_ と low の値を計算しています。なお、記号の合計値 (count_sum[size]) が range_ より大きくなると、temp が 0 になってしまうので、レンジコーダは正常に動作しません。合計値が MIN_RANGE 以上にならないように注意してください。

今回のプログラムでは、記号の出現頻度が 0xffff 以下になるように丸められているので、size が 256 以下であれば合計値が MIN_RANGE をオーバーすることありません。size を 256 より大きくする場合は注意してください。

次は range_ と low の値を正規化する RangeCoder のメソッド encode_normalize() を作ります。

リスト : 符号化のときの正規化

    def encode_normalize(self):
        if self.low > MAX_RANGE:
            # 桁上がり
            self.buff += 1
            self.low &= MASK
            if self.cnt > 0:
                self.stream.putc(self.buff)
                for _ in range(self.cnt - 1): self.stream.putc(0)
                self.buff = 0
                self.cnt = 0
        while self.range_ < MIN_RANGE:
            if self.low < (0xff << SHIFT):
                self.stream.putc(self.buff)
                for _ in range(self.cnt): self.stream.putc(0xff)
                self.buff = (self.low >> SHIFT) & 0xff
                self.cnt = 0
            else:
                self.cnt += 1
            self.low = (self.low << 8) & MASK
            self.range_ <<= 8

基本的な処理は簡単で、range_ の値が MIN_RANGE 未満の場合は、range_ と low の値を 256 倍して符号語をファイルへ出力します。このときポイントになるのが桁上がりの処理です。

●桁上がりの処理

low を 256 倍するとき、最上位の 8 ビットを符号語としてすぐにファイルへ出力すると、桁上がりに対応することができません。そこで、最上位 8 ビットを変数 buff に格納することにします。low を 256 倍するときは、先に buff の値を符号語として出力し、low の最上位 8 ビットを buff に格納します。そして、low の値は MAX_RANGE 以下 (32 ビット) で保持します。桁上がりが発生したら、buff の値を +1 すればいいわけです。

ここで問題点が一つあります。それは buff の値が 255 (0xff) のとき、桁上がりが発生すると buff が 0x100 になることです。この場合、buff の値は 0 になり、先に出力した符号語を +1 しないといけません。そこで、low の最上位 8 ビットが 0xff の場合は buff を出力しないで、0xff の個数を変数 cnt でカウントすることにします。次の図を見てください。

たとえば、buff が 0x12 で low の最上位 8 ビットが 0xff の場合、buff を出力しないで cnt の値を +1 します。そして、low を 8 ビット左へシフト (256 倍) します。つまり、buff と cnt で [12, ff] を表していることになります。また、最上位 8 ビットの値が 0xff で続く場合もありえます。この場合は、cnt の値を +1 していきます。たとえば cnt が 3 であれば、buff と cnt で [12, ff, ff, ff] を表します。

桁上がりが発生したときは buff を +1 します。このとき、cnt が 0 よりも大きい場合はバッファを出力します。たとえば、buff が 0x12 で cnt が 3 の場合、バッファは [12, ff, ff, ff] を表しています。これに 1 を加えると、[13, 00, 00, 00] になります。つまり、buff に 1 を加えてから出力し、そのあと 0 を出力すればいいわけです。このとき、最後の 0 を buff にセットするので、出力する 0 の個数は cnt - 1 になります。

プログラムの説明に戻ります。encode_normalize() の前半部分が桁上がりの処理です。low が MAX_RANGE より大きい場合は桁上がりが発生しています。buff の値を +1 して、low の値を 32 ビットの範囲に収めます。バッファを出力する場合は、buff を出力したあとで cnt - 1 個の 0 を出力します。そのあと、buff と cnt に 0 をセットします。

次の while ループで low と range_ の値を 256 倍していきます。このとき、low の上位 8 ビットの値をチェックします。low が 0xff000000 未満の場合、low の上位 8 ビットは 0xff ではありません。この場合、buff を出力したあと、0xff を cnt 個出力するだけです。そのあとで、low の上位 8 ビットを buff にセットし、cnt を 0 にします。low の上位 8 ビットが 0xff の場合は cnt を +1 するだけです。最後に low と range_ を 256 倍 (左へ 8 ビットシフト) します。

次はレンジコーダの符号化を終了するメソッド finish() を作ります。

リスト : 符号化の終了

    def finish(self):
        c = 0xff
        if self.low > MAX_RANGE:
            # 桁上がり
            self.buff += 1
            c = 0
        self.stream.putc(self.buff)
        for _ in range(self.cnt): self.stream.putc(c)
        #
        self.stream.putc((self.low >> 24) & 0xff)
        self.stream.putc((self.low >> 16) & 0xff)
        self.stream.putc((self.low >> 8) & 0xff)
        self.stream.putc(self.low & 0xff)

最初に桁上がりをチェックし、そうであれば buff を +1 して変数 c を 0 にします。次に、バッファの内容を出力します。buff を出力したあと、0 または 0xff を cnt 個出力します。最後に low の値 (4 バイト) を出力します。

●レンジコーダによる符号化

最後に、符号化を行う関数 encode() を作ります。

リスト : レンジコーダによる符号化

# ファイルの読み込み
def read_file(fin):
    while True:
        c = fin.getc()
        if c is None: break
        yield c

# 符号化
def encode(fin, fout):
    count = [0] * 256
    for x in read_file(fin):
        count[x] += 1
    rc = RangeCoder(fout, ENCODE)
    freq = Freq(count)
    freq.write_count_table(fout)
    fin.stream.seek(0)
    for x in read_file(fin):
        freq.encode(rc, x)
    rc.finish()

引数 fin が入力ファイル、fout が出力ファイルです。どちらも ByteIO のオブジェクトです。最初に記号の出現頻度を求めます。read_file() で fin から記号を読み込み、その個数を配列 count でカウントします。

次に、RangeCoder のオブジェクトを生成して変数 rc にセットします。このとき、出力ファイル fout と ENCODE を指定します。Freq のオブジェクトを生成するときは count を渡します。そして、メソッド write_count_table() で出現頻度表を出力ファイルに付加します。

あとは seek(0) で入力ファイルを巻き戻し、read_file() でファイルから記号を読み込みます。そして、freq.encode() で記号 x を符号化します。最後に rc.finish() でレンジコーダを終了することをお忘れなく。

●復号のプログラム

次は復号を行う Freq のメソッド decode() を作ります。

リスト : 復号

    def decode(self, rc):
        # 記号の探索
        def search_code(value):
            i = 0
            j = self.size - 1
            while i < j:
                k = (i + j) // 2
                if self.count_sum[k + 1] <= value:
                    i = k + 1
                else:
                    j = k
            return i
        #
        temp = rc.range_ // self.count_sum[self.size]
        c = search_code(rc.low // temp)
        rc.low -= temp * self.count_sum[c]
        rc.range_ = temp * self.count[c]
        rc.decode_normalize()
        return c

レンジコーダで記号を復号する場合、累積度数表 count_sum から次式の条件を満たす記号 c を探します。

   count_sum[c] / count_sum[size] <= low / range_ < count_sum[c + 1] / count_sum[size]
=> count_sum[c] <= low * count_sum[size] / range_ < count_sum[c + 1]
   temp = range_ / count_sum[size] とすると
=> count_sum[c] <= low / temp < count_sum[c + 1]

プログラムでは、range_ / count_sum[size] の値を temp にセットし、low / temp を内部関数 search_code() で二分探索しています。レンジコーダは整数で計算するので、割り算の結果が 0 にならないよう計算の順番に注意してください。

記号 c を求めたあと、low と range_ の値を更新します。range_ は記号の出現確率で縮小します。

   range_= range_ * count[c] / count_sum[size] 
=> range_= temp * count[c]

復号の場合は low の値を減らします。

   low -= range_ * count_sum[c] / count_sum[size] 
=> low -= temp * count_sum[c]

最後に、メソッド decode_normalize() を呼び出して、求めた記号 c を返します。

次は low と range_ の値を更新 (正規化) するメソッド decode_normalize() を作ります。

リスト : 復号の正規化

    def decode_normalize(self):
        while self.range_ < MIN_RANGE:
            self.range_ <<= 8
            self.low = ((self.low << 8) + self.stream.getc()) & MASK

decode_normalize() はクラス RangeCoder のメソッドとして定義します。符号化と違って復号の正規化はとても簡単です。range_ が MIN_RANGE よりも小さい場合は range_ と low を 256 倍し、ファイルから 1 記号読み込んで low に加算するだけです。

●レンジコーダによる復号

最後に復号を行う関数 decode() を作ります。

リスト : レンジコーダによる復号

# 出現頻度表の読み込み
def read_count_table(fin):
    count = [0] * 256
    for x in range(256):
        count[x] = (fin.getc() << 8) + fin.getc()
    return count

# 復号
def decode(fin, fout, size):
    freq = Freq(read_count_table(fin))
    rc = RangeCoder(fin, DECODE)
    for _ in range(size):
        fout.putc(freq.decode(rc))

復号は簡単です。関数 read_count_table() で fin から出現頻度表を読み込み、Freq のオブジェクトを生成します。そのあとで、RangeCoder のオブジェクトを生成します。順番を逆にすると正常に動作しません。ご注意ください。あとは、size 個の記号を freq.decode() で復号するだけです。復号した記号は putc() で fout に出力します。

プログラムリスト１ (RangeCoder)
プログラムリスト２ (Freq, main)

●評価結果

それでは、実際に Canterbury Corpus で配布されているテストデータ The Canterbury Corpus を圧縮してみましょう。結果は次にようになりました。

                   表 : レンジコーダの結果
           () はファイルサイズと出現頻度表を引いた値

  ファイル名      サイズ        RangeCoder        下限値   符号化    復号
  ------------------------------------------------------------------------
  alice29.txt    152,089     87,380  ( 86,864)     86,837   0.301   0.453
  asyoulik.txt   125,179     75,770  ( 75,254)     75,235   0.303   0.375
  cp.html         24,603     16,603  ( 16,087)     16,082   0.081   0.082
  fields.c        11,150      7,500  (  6,984)      6,980   0.127   0.041
  grammar.lsp      3,721      2,675  (  2,159)      2,155   0.039   0.045
  kennedy.xls  1,029,744    460,622  (460,106)    459,971   1.937   2.772
  lcet10.txt     426,754    249,679  (249,163)    249,071   0.851   1.215
  plrabn12.txt   481,861    273,569  (273,053)    272,936   0.972   1.384
  ptt5           513,216     78,226  ( 77,710)     77,636   0.791   1.276
  sum             38,240     25,994  ( 25,478)     25,473   0.153   0.141
  xargs.1          4,227      3,109  (  2,593)      2,589   0.034   0.019
  ------------------------------------------------------------------------
  合計         2,810,784  1,281,127 (1,275,451) 1,274,965   5.589   7.803

符号化と復号の単位 : 秒
実行環境 : Python 3.8.10, Ubuntu 20.04 LTS (WSL1), Intel Core i5-6200U 2.30GHz

レンジコーダの圧縮率は圧縮の限界に近い値となりました。出現頻度表が 512 バイトあるので、小さなファイルの圧縮率はハフマン符号よりも悪くなる場合がありますが、大きなファイルではハフマン符号よりも高い圧縮率になりました。特に ptt5 の圧縮率は、レンジコーダの方がとても高くなります。記号に 1 ビット未満の符号語を割り当てることができるレンジコーダ (算術符号) の特徴が結果に出ていると思います。

実行時間ですが、符号化はハフマン符号よりも速くなりましたが、復号はハフマン符号よりも遅くなりました。レンジコーダの場合、記号の探索処理にどうしても時間がかります。復号に時間がかかるアルゴリズムといえるでしょう。

なお、実行時間の結果は M.Hiroi のコーディング、実行したマシン、プログラミング言語などの環境に大きく依存しています。また、これらの環境だけではなく、データの種類によっても実行時間はかなり左右されます。興味のある方は、いろいろなデータをご自分の環境で試してみてください。

●参考文献, URL

植松友彦, 『文書データ圧縮アルゴリズム入門』, ＣＱ出版社, 1994
奥村晴彦, 『データ圧縮の基礎から応用まで』, C MAGAZINE 2002 年 7 月号, ソフトバンク
岡野原大輔, 『高速な算術圧縮を実現する「Range Coder」』, CodeZine
広井誠, 『高性能圧縮ツール bsrc の理論と実装 (後編)』, Interface 2004 年 1 月号, ＣＱ出版社

●プログラムリスト１

#
# rangecoder.py : レンジコーダ (Range Coder)
#
#                 Copyright (C) 2007-2022 Makoto Hiroi
#
from byteio import *

# 定数
ENCODE = "encode"
DECODE = "decode"
MAX_RANGE = 0xffffffff
MIN_RANGE = 0x1000000
MASK      = 0xffffffff
SHIFT     = 24

# レンジコーダー
class RangeCoder:
    def __init__(self, s, mode):
        self.stream = s
        self.range_ = MAX_RANGE
        self.buff = 0
        self.cnt = 0
        if mode == ENCODE:
            self.low = 0
        elif mode == DECODE:
            # buff の初期値 (0) を読み捨てる
            self.stream.getc()
            # 4 byte read
            self.low = self.stream.getc()
            self.low = (self.low << 8) + self.stream.getc()
            self.low = (self.low << 8) + self.stream.getc()
            self.low = (self.low << 8) + self.stream.getc()
        else:
            raise "RangeCoder mode error"

    # 符号化の正規化
    def encode_normalize(self):
        if self.low > MAX_RANGE:
            # 桁上がり
            self.buff += 1
            self.low &= MASK
            if self.cnt > 0:
                self.stream.putc(self.buff)
                for _ in range(self.cnt - 1): self.stream.putc(0)
                self.buff = 0
                self.cnt = 0
        while self.range_ < MIN_RANGE:
            if self.low < (0xff << SHIFT):
                self.stream.putc(self.buff)
                for _ in range(self.cnt): self.stream.putc(0xff)
                self.buff = (self.low >> SHIFT) & 0xff
                self.cnt = 0
            else:
                self.cnt += 1
            self.low = (self.low << 8) & MASK
            self.range_ <<= 8

    # 復号の正規化
    def decode_normalize(self):
        while self.range_ < MIN_RANGE:
            self.range_ <<= 8
            self.low = ((self.low << 8) + self.stream.getc()) & MASK

    # 終了
    def finish(self):
        c = 0xff
        if self.low > MAX_RANGE:
            # 桁上がり
            self.buff += 1
            c = 0
        self.stream.putc(self.buff)
        for _ in range(self.cnt): self.stream.putc(c)
        #
        self.stream.putc((self.low >> 24) & 0xff)
        self.stream.putc((self.low >> 16) & 0xff)
        self.stream.putc((self.low >> 8) & 0xff)
        self.stream.putc(self.low & 0xff)

●プログラムリスト２

#
# rc0.py : レンジコーダ
#
#          Copyright (C) 2007-2022 Makoto Hiroi
#
import argparse, os.path, time
from rangecoder import *

# 出現頻度表
class Freq:
    def __init__(self, count):
        size = len(count)
        self.size = size
        m = max(count)
        # 2 バイトに収める
        if m > 0xffff:
            self.count = [0] * size
            n = 0
            while m > 0xffff:
                m >>= 1
                n += 1
            for x in range(size):
                if count[x] != 0:
                    self.count[x] = (count[x] >> n) | 1
        else:
            self.count = count[:]
        self.count_sum = [0] * (size + 1)
        # 累積度数表
        for x in range(size):
            self.count_sum[x + 1] = self.count_sum[x] + self.count[x]

    #
    def write_count_table(self, fout):
        for x in self.count:
            fout.putc(x >> 8)
            fout.putc(x & 0xff)

    # 符号化
    def encode(self, rc, c):
        temp = rc.range_ // self.count_sum[self.size]
        rc.low += self.count_sum[c] * temp
        rc.range_ = self.count[c] * temp
        rc.encode_normalize()

    # 復号
    def decode(self, rc):
        # 記号の探索
        def search_code(value):
            i = 0
            j = self.size - 1
            while i < j:
                k = (i + j) // 2
                if self.count_sum[k + 1] <= value:
                    i = k + 1
                else:
                    j = k
            return i
        #
        temp = rc.range_ // self.count_sum[self.size]
        c = search_code(rc.low // temp)
        rc.low -= temp * self.count_sum[c]
        rc.range_ = temp * self.count[c]
        rc.decode_normalize()
        return c

# ファイルの読み込み
def read_file(fin):
    while True:
        c = fin.getc()
        if c is None: break
        yield c

# レンジコーダによる符号化
def encode(fin, fout):
    count = [0] * 256
    for x in read_file(fin):
        count[x] += 1
    rc = RangeCoder(fout, ENCODE)
    freq = Freq(count)
    freq.write_count_table(fout)
    fin.stream.seek(0)
    for x in read_file(fin):
        freq.encode(rc, x)
    rc.finish()

# 出現頻度表の読み込み
def read_count_table(fin):
    count = [0] * 256
    for x in range(256):
        count[x] = (fin.getc() << 8) + fin.getc()
    return count

# レンジコーダによる復号
def decode(fin, fout, size):
    freq = Freq(read_count_table(fin))
    rc = RangeCoder(fin, DECODE)
    for _ in range(size):
        fout.putc(freq.decode(rc))

# 符号化
def encode_file(name1, name2):
    size = os.path.getsize(name1)
    with ByteIO(name1, ROPEN) as fin, ByteIO(name2, WOPEN) as fout:
        fout.putc((size >> 24) & 0xff)
        fout.putc((size >> 16) & 0xff)
        fout.putc((size >> 8) & 0xff)
        fout.putc(size & 0xff)
        if size > 0: encode(fin, fout)

# 復号
def decode_file(name1, name2):
    with ByteIO(name1, ROPEN) as fin, ByteIO(name2, WOPEN) as fout:
        size = 0
        for _ in range(4):
            size = (size << 8) + fin.getc()
        if size > 0: decode(fin, fout, size)

# オプション解析
parser = argparse.ArgumentParser(description='RangeCoder符号')
parser.add_argument('name1', help='入力ファイル')
parser.add_argument('name2', help='出力ファイル')
parser.add_argument('-e', '--encode', action='store_true', help='符号化')
parser.add_argument('-d', '--decode', action='store_true', help='復号')
args = parser.parse_args()

# 実行
s = time.time()
if args.encode and not args.decode:
    encode_file(args.name1, args.name2)
elif args.decode:
    decode_file(args.name1, args.name2)
else:
    print('option error')
e = time.time()
print('{:.3f}'.format(e - s))

初版 2007 年 6 月 3 日
改訂 2022 年 9 月 23 日