データ圧縮のお話です。拙作のページ シャノン・ファノ符号とハフマン符号 で説明したハフマン符号は各記号の出現確率を調べ、頻繁に現れる記号は短いビットで表し、あまり現れない記号は長いビットで表すことで、データを圧縮する古典的な方法です。古典的とはいっても、ほかのアルゴリズムと簡単に組み合わせることができるため、ハフマン符号は今でも現役のアルゴリズムです。
記号の出現確率だけを利用してデータを圧縮する方法は、ハフマン符号のほかに「算術符号」という方法があります。一般に、算術符号はハフマン符号よりも性能が良いのですが、実装方法が難しくて実行速度がハフマン符号よりも遅く、なおかつ特許の問題もあって、実際に使われることはあまりありませんでした。
ところが最近になって、「レンジコーダ (Range Coder) 」という方法が注目されています。レンジコーダは原理的には算術符号と同じ方法ですが、参考文献, URL 2 によると 『（おそらく）特許フリー』 とのことで、性能は算術符号に比べるとわずかに劣化しますが、実装方法はとても簡単で実行速度も算術符号より高速です。もちろん、ハフマン符号よりも高性能です。
今回はレンジコーダについて詳しく説明します。そして、実際にファイルを圧縮してみましょう。
最初に算術符号の基本的な考え方について説明します。算術符号は記号（文字）列全体を一つの符号語にする方法で、1960 年代に P. Elias によって提案されました。
なお、算術符号については拙作のページ Common Lisp 入門番外編 Lisp で算術符号 でも説明しています。内容は重複していますが、ご了承くださいませ。
算術符号は、記号列を実数 0 と 1 の間の区間を用いて表します。たとえば、記号は {a, b, c} の 3 種類があり、出現確率はそれぞれ 0.2, 0.6, 0.2 とします。算術符号は、区間を記号の出現確率に比例した小区間に分割していくことで符号化を行います。それでは、記号列 abbbc を符号化してみましょう。次の図を見てください。
1.0 ┬ ┌→ 0.2 ┬ ┌→ 0.16 ┬ ┌→ 0.136 ┬ ┌→ 0.1216 ┬ ─→ 0.1216 ┬
c │ │ │ │ │ │ │ │ │ │
0.8 ┼ │ 0.16 ┼ ┘ 0.136 ┼ ┘ 0.1216 ┼ ┘ 0.11296 ┼ ┐ 0.1198728 ┼
│ │ │ │ │ │ │ │
b │ │ │ │ │ │ │ │
│ │ │ │ │ │ │ │
0.2 ┼ ┘ 0.04 ┼ ┐ 0.064 ┼ ┐ 0.0784 ┼ ┐ 0.08704 ┼ │ 0.114688 ┼
a │ │ │ │ │ │ │ │ │ │
0.0 ┴ ─→ 0.0 ┴ └→ 0.04 ┴ └→ 0.64 ┴ └→ 0.0784 ┴ └→ 0.11296 ┴
記号 a ab abb abbb 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)2
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) にも算術圧縮のプログラム (Ｃ言語) があります。
このほかに、拙作のページ Common Lisp 入門番外編 Lisp で算術符号 では、Common Lisp で算術符号のプログラムを作成しています。学習用のプログラムなので実用性はありませんが、Lisp に興味のある方は読んでみてください。
次はレンジコーダについて説明します。なお、拙作のページ Common Lisp 入門番外編 Lisp でレンジコーダ と内容が重複していますが、ご了承くださいませ。
レンジコーダは 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 |
符号化の過程は次のようになります。
low = low + (range * 記号の下限値)
range = range * 記号の出現確率
[ low, range] [ low, range] (数値は 16 進数)
[ 0,1000000] - d -> [e00000, 200000]
[e00000, 200000] - c -> [f80000, 40000]
[f80000, 40000] - b -> [fa0000, 10000]
[fa0000, 10000] - b -> [fa8000, 4000] 256 倍して fa を出力
[800000, 400000] - a -> [800000, 200000]
[800000, 200000] - a -> [800000, 100000]
[800000, 100000] - a -> [800000, 80000]
[800000, 80000] - a -> [800000, 40000] low [80, 00, 00] を出力
符号語 => [fa, 80, 00, 00]
図 : 符号化の過程（レンジコーダ）
記号を読み込むたびに、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 に対応する記号を探せばよいのです。
0 low1 code low1+range1 range
├──────┼────┼───────┼───────┤
┌──────┘ │ │
│ ┌──────┘ │
│ │ ┌──────┘
↓ ↓ ↓
├────┼───────┤
0 code1 range1
図 : 区間の更新
このように、符号語 code から下限値の増分を引き算することで、区間を [low, low + range) から [0, range) に変換することができるわけです。したがって、復号処理では下限値 low の値を覚えておく必要はありません。
range が 0x10000 より小さくなったら range を 256 倍するのは符号化と同じです。このとき符号語 code も 256 倍して、新しい符号語を 1 バイト読み込んで code に加算します。これで符号語を復号することができます。
それでは、復号の過程を具体的に説明しましょう。次の図を見てください。
code = code - (range * 記号の下限値)
range = range * 記号の出現確率
符号語を 3 バイト [fa, 80, 00] 読み込み code を初期化
[ code, range] [ code, range] (数値は 16 進数)
[fa8000, 1000000] - d -> [1a8000, 200000]
[1a8000, 200000] - c -> [ 28000, 40000]
[ 28000, 40000] - b -> [ 8000, 10000]
[ 8000, 10000] - b -> [ 0, 4000] 256 倍する
符号語を 1 バイト (00) code に加算
[ 0, 400000] - a -> [ 0, 200000]
[ 0, 200000] - a -> [ 0, 100000]
[ 0, 100000] - a -> [ 0, 80000]
[ 0, 80000] - a -> [ 0, 40000]
記号列 => "dcbbaaaa"
図 : 復号の過程（レンジコーダ）
最初に 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" を求めることができます。
それではプログラムを作りましょう。最初に、レンジコーダを表すクラスを定義します。
リスト : レンジコーダ
# 定数
ENCODE = "encode"
DECODE = "decode"
MAX_RANGE = 0x100000000
MIN_RANGE = 0x1000000
MASK = 0xffffffff
SHIFT = 24
class RangeCoder:
def __init__(self, file, mode):
self.file = file
self.range = MAX_RANGE
self.buff = 0
self.cnt = 0
if mode == ENCODE:
self.low = 0
elif mode == DECODE:
# buff の初期値 (0) を読み捨てる
getc(self.file)
# 4 byte read
self.low = getc(self.file)
self.low = (self.low << 8) + getc(self.file)
self.low = (self.low << 8) + getc(self.file)
self.low = (self.low << 8) + getc(self.file)
else:
raise "RangeCoder mode error"
クラス名は RangeCoder としました。RangeCoder() の引数 file が入出力ファイル、mode はモードを表していて、ENCODE が符号化で DECODE が復号です。当然ですがファイルは、符号化のときはライトモード ("wb") で、復号のときはリードモード ("rb") でオープンしてください。インスタンス変数 range は幅、buff と cnt が桁上がり用のバッファ、low は符号化のときは下限値で復号のときは符号語 code を表します。
Python は多倍長整数をサポートしているので、range は MAX_RANGE (0x100000000) で初期化します。符号化の場合、low は 0 で初期化します。復号の場合、low の範囲は 32 ビットなので、ファイルから 4 バイト分読み込んで初期化します。符号化のとき、最初に出力される 1 バイトは buff の初期値 0 なので、復号のときは 1 バイト読み捨てることに注意してください。
なお、符号化の処理で一番最初に符号語を出力するとき、buff の値を出力しないようにすれば、復号の初期化で 1 バイト読み捨てる処理は不要になります。また、圧縮後のファイルサイズも 1 バイト少なくなります。興味のある方はプログラムを改造してみてください。
RangeCoder では、次に示すメソッドを用意します。
記号の符号化・復号を行うメソッドは、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 xrange(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 xrange(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 になってしまうので、レンジコーダは正常に動作しません。count_sum[size] が MIN_RANGE 以上にならないように注意してください。今回のプログラムは記号の出現頻度が 0xffff 以下になるように丸められているので、size が 256 以下であれば count_sum[size] が 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:
putc(self.file, self.buff)
for _ in xrange(self.cnt - 1): putc(self.file, 0)
self.buff = 0
self.cnt = 0
while self.range < MIN_RANGE:
if self.low < (0xff << SHIFT):
putc(self.file, self.buff)
for _ in xrange(self.cnt): putc(self.file, 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, cnt low
------------------------------
12, 0 <= [ff, 34, 56, 78]
------------------------------
12, 1 <= [34, 56, 78, 00]
図 : バッファの動作
たとえば、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
putc(self.file, self.buff)
for _ in xrange(self.cnt): putc(self.file, c)
# low を出力
putc(self.file, (self.low >> 24) & 0xff)
putc(self.file, (self.low >> 16) & 0xff)
putc(self.file, (self.low >> 8) & 0xff)
putc(self.file, self.low & 0xff)
最初に桁上がりをチェックし、そうであれば buff を +1 して変数 c を 0 にします。次に、バッファの内容を出力します。buff を出力したあと、0 または 0xff を cnt 個出力します。最後に low の値 (4 バイト) を出力します。
最後に、符号化を行う関数 encode を作ります。
リスト : レンジコーダによる符号化
# 記号の読み込み
def read_file(fin):
while True:
c = getc(fin)
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.seek(0)
for x in read_file(fin):
freq.encode(rc, x)
rc.finish()
引数 fin が入力ファイル、fout が出力ファイルです。どちらのファイルもバイナリモードでオープンしてください。最初に記号の出現頻度を求めます。関数 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]
レンジコーダは整数で計算するので、割り算の結果が 0 にならないよう計算の順番に注意してください。プログラムでは、range / count_sum[size] の値を temp にセットし、low / temp を内部関数 search_code で二分探索しています。
記号 c を求めたあと、low と range の値を更新します。range は記号の出現確率で縮小すればいいので、range * count[c] / count_sum[size] = temp * count[c] となります。復号の場合は low の値を減らします。その値は range * count_sum[c] / count_sum[size] = 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) + getc(self.file)) & MASK
decode_normalize はクラス RangeCoder のメソッドとして定義します。符号化と違って復号の正規化はとても簡単です。range が MIN_RANGE よりも小さい場合は range と low を 256 倍し、ファイルから 1 記号読み込んで low に加算するだけです。
最後に復号を行う関数 decode を作ります。
リスト : レンジコーダによる復号
# 出現頻度表の読み込み
def read_count_table(fin):
count = [0] * 256
for x in xrange(256):
count[x] = (getc(fin) << 8) + getc(fin)
return count
# 復号
def decode(fin, fout, size):
freq = Freq(read_count_table(fin))
rc = RangeCoder(fin, DECODE)
for _ in xrange(size):
putc(fout, freq.decode(rc))
復号は簡単です。関数 read_count_table で fin から出現頻度表を読み込み、Freq のオブジェクトを生成します。そのあとで、RangeCoder のオブジェクトを生成します。順番を逆にすると正常に動作しません。ご注意ください。あとは、size 個の記号を freq.decode で復号するだけです。復号した記号は putc で fout に出力します。
それでは、実際に Canterbury Corpus で配布されているテストデータ The Canterbury Corpus を圧縮してみましょう。結果は次にようになりました。
表 : レンジコーダの結果
() はファイルサイズと出現頻度表を引いた値
ファイル名 サイズ RangeCoder 符号化 復号 下限値
---------------------------------------------------------------------
alice29.txt 152,089 87,380 ( 86,864) 2.16 2.95 86,837
asyoulik.txt 125,179 75,770 ( 75,254) 1.80 2.44 75,235
cp.html 24,603 16,603 ( 16,087) 0.36 0.49 16,082
fields.c 11,150 7,500 ( 6,984) 0.16 0.22 6,980
grammar.lsp 3,721 2,675 ( 2,159) 0.06 0.08 2,155
kennedy.xls 1,029,744 460,622 (460,106) 13.36 18.55 459,971
lcet10.txt 426,754 249,679 (249,163) 6.07 8.29 249,071
plrabn12.txt 481,861 273,569 (273,053) 6.78 9.32 272,936
ptt5 513,216 78,226 ( 77,710) 5.39 8.36 77,636
sum 38,240 25,994 ( 25,478) 0.57 0.75 25,473
xargs.1 4,227 3,109 ( 2,593) 0.06 0.10 2,589
---------------------------------------------------------------------
合計 2,810,784 1,281,127 (1,275,451) 36.77 51.55 1,274,965
# 符号化と復号の単位 : 秒
実行環境 : Windows XP, celeron 1.40 GHz, Python 2.4.2
レンジコーダの圧縮率は圧縮の限界に近い値となりました。出現頻度表が 512 バイトあるので、小さなファイルの圧縮率はハフマン符号よりも悪くなる場合がありますが、大きなファイルではハフマン符号よりも高い圧縮率になりました。特に ptt5 の圧縮率は、レンジコーダの方がとても高くなります。記号に 1 ビット未満の符号語を割り当てることができるレンジコーダ (算術符号) の特徴が結果に出ていると思います。
実行時間ですが、符号化はハフマン符号よりも少し速くなりましたが、復号はハフマン符号よりもかなり遅くなりました。レンジコーダの場合、記号の探索処理にどうしても時間がかります。LZT 符号よりも遅いので、復号に時間がかかるアルゴリズムといえるでしょう。
なお、実行時間の結果は M.Hiroi のコーディング、実行したマシン、プログラミング言語などの環境に大きく依存しています。また、これらの環境だけではなく、データの種類によっても実行時間はかなり左右されます。興味のある方は、いろいろなデータをご自分の環境で試してみてください。
# coding: utf-8
#
# rangecoder.py : レンジコーダ (Range Coder)
#
# Copyright (C) 2007 Makoto Hiroi
#
# 定数
ENCODE = "encode"
DECODE = "decode"
MAX_RANGE = 0x100000000
MIN_RANGE = 0x1000000
MASK = 0xffffffff
SHIFT = 24
# バイト単位の入出力
def getc(f):
c = f.read(1)
if c == '': return None
return ord(c)
def putc(f, x):
f.write(chr(x & 0xff))
#
class RangeCoder:
def __init__(self, file, mode):
self.file = file
self.range = MAX_RANGE
self.buff = 0
self.cnt = 0
if mode == ENCODE:
self.low = 0
elif mode == DECODE:
# buff の初期値 (0) を読み捨てる
getc(self.file)
# 4 byte read
self.low = getc(self.file)
self.low = (self.low << 8) + getc(self.file)
self.low = (self.low << 8) + getc(self.file)
self.low = (self.low << 8) + getc(self.file)
else:
raise "RangeCoder mode error"
# 符号化の正規化
def encode_normalize(self):
if self.low >= MAX_RANGE:
# 桁上がり
self.buff += 1
self.low &= MASK
if self.cnt > 0:
putc(self.file, self.buff)
for _ in xrange(self.cnt - 1): putc(self.file, 0)
self.buff = 0
self.cnt = 0
while self.range < MIN_RANGE:
if self.low < (0xff << SHIFT):
putc(self.file, self.buff)
for _ in xrange(self.cnt): putc(self.file, 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) + getc(self.file)) & MASK
# 終了
def finish(self):
c = 0xff
if self.low >= MAX_RANGE:
# 桁上がり
self.buff += 1
c = 0
putc(self.file, self.buff)
for _ in xrange(self.cnt): putc(self.file, c)
#
putc(self.file, (self.low >> 24) & 0xff)
putc(self.file, (self.low >> 16) & 0xff)
putc(self.file, (self.low >> 8) & 0xff)
putc(self.file, self.low & 0xff)
# coding: utf-8
#
# rc0.py : レンジコーダ
#
# Copyright (C) 2007 Makoto Hiroi
#
import time, sys, getopt, os.path
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 xrange(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 xrange(size):
self.count_sum[x + 1] = self.count_sum[x] + self.count[x]
#
def write_count_table(self, fout):
for x in self.count:
putc(fout, x >> 8)
putc(fout, 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 = getc(fin)
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.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 xrange(256):
count[x] = (getc(fin) << 8) + getc(fin)
return count
# レンジコーダによる復号
def decode(fin, fout, size):
freq = Freq(read_count_table(fin))
rc = RangeCoder(fin, DECODE)
for _ in xrange(size):
putc(fout, freq.decode(rc))
# 符号化
def encode_file(name1, name2):
size = os.path.getsize(name1)
infile = open(name1, "rb")
outfile = open(name2, "wb")
putc(outfile, (size >> 24) & 0xff)
putc(outfile, (size >> 16) & 0xff)
putc(outfile, (size >> 8) & 0xff)
putc(outfile, size & 0xff)
if size > 0: encode(infile, outfile)
infile.close()
outfile.close()
# 復号
def decode_file(name1, name2):
infile = open(name1, "rb")
outfile = open(name2, "wb")
size = 0
for _ in xrange(4):
size = (size << 8) + getc(infile)
if size > 0: decode(infile, outfile, size)
infile.close()
outfile.close()
#
def main():
eflag = False
dflag = False
opts, args = getopt.getopt(sys.argv[1:], 'ed')
for x, y in opts:
if x == '-e' or x == '-E':
eflag = True
elif x == '-d' or x == '-D':
dflag = True
if eflag and dflag:
print 'option error'
elif eflag:
encode_file(args[0], args[1])
elif dflag:
decode_file(args[0], args[1])
else:
print 'option error'
#
s = time.clock()
main()
e = time.clock()
print "%.3f" % (e - s)