読者です 読者をやめる 読者になる 読者になる

ハフマン符号化 : 整理

common lisp algorithm

長さ制限付きハフマン符号化に合わせて、普通の*1ハフマン符号化のソースコードも整理。

###

最初にハフマン符号化本体のソースを載せて、次にその中で利用しているヒープ(順序キュー)の実装を載せる。


ハフマン符号ソース。
参照: nlet, ???-obj-???

;;;;;;;;;;;;;;;;;
;;;; 処理速度優先
(declaim (optimize (speed 3) (safety 0) (debug 0) (compilation-speed 0)))

;; コードの出現頻度表をもとに、各コードに適切な符号ビット長を計算する
;;  heap-pop2: ヒープの先頭から二つの要素を取り出す
;;  heap-push-and-pop2: ヒープに要素を追加した後、先頭から二つの要素を取り出す
;;  ※ ヒープには、出現頻度が低い順に要素が並んでいる
(defun calc-code->bitlength-table (#1=code-frequency-table)
  (declare ((simple-array fixnum) #1#))
  (let ((heap (make-heap #1#)))
    (labels ((iter (1st 2nd end?) ; ※1
               (if end? 
                    ;; ヒープが空になったら、各コードに対応するビット長を数えて、返す
                    (count-code-bitlength (package-obj 1st 2nd) (length #1#))
                  (multiple-value-call #'iter (heap-push-and-pop2 (package-obj 1st 2nd) heap)))))
      (multiple-value-call #'iter (heap-pop2 heap)))))
;; ※1 iter関数呼び出し(ループ)は、(count-if #'plusp #1#)回行われる。
;;      --> これは一回でも出現したコードの個数。この数をNとする。
;;     各iter関数呼び出しで、一回ずつheap-push-and-pop2関数が呼ばれ、
;;     この操作には(実装上)約logN要するので、全体の処理オーダーは、だいたいO(NlogN)となる。


;; ハフマン木(packaged-obj型)を受け取り、各コードに対応する符号ビット長を保持した配列を返す
(defun count-code-bitlength (huffman-tree code-limit)
  (declare (fixnum code-limit))
  (let ((table (make-array code-limit :initial-element 0 :element-type 'fixnum)))
    (nlet self ((tree huffman-tree) (depth 0))
      (declare (fixnum depth))
      (when tree
        (if (not (packaged-obj-p tree))
            (setf (aref table (code-obj-code tree)) depth)
          (destructuring-bind (left . right) (packaged-obj-pair tree)
            (self left  (1+ depth))
            (self right (1+ depth))))))
    table)) 


ヒープの実装。

(declaim (ftype (function (array-index) (simple-array obj)) init-heap-buf)
         (inline upheap downheap heap-push heap-pop heap-pop2 heap-push-and-pop2))

;; 固定サイズの配列によるヒープ実装
;; - buf[i]の子要素は、buf[i+i]及びbuf[i+i+1]
;; - 親要素は子要素よりも常に小さい (obj<関数で)
(defstruct (heap (:constructor heap (size-limit &aux (buf (init-heap-buf size-limit)))))
  (size 0   :type array-index)
  (buf  #() :type (simple-array obj))) 

(defun init-heap-buf (size-limit)
  (make-array (1+ (the fixnum size-limit))
              :initial-element (load-time-value (make-obj :cost 0) t)
              :element-type 'obj))

;; buf[i]にある要素が、その親要素よりも小さい場合、
;; 親要素との交換を繰り返すことで、ヒープ内の適切な位置に再配置する
(defun upheap (i buf)
  (declare ((simple-array obj) buf)
           (fixnum i))
  (loop WITH obj = (aref buf i)
        FOR  i/2 = (ash i -1)
        WHILE (obj< obj (aref buf i/2))
    DO
      (setf (aref buf i) (aref buf i/2)
                      i            i/2)
    FINALLY
      (setf (aref buf i) obj)))

(defun heap-push (obj heap)
  (with-slots (size buf) (the heap heap)
    (setf (aref buf (incf size)) obj)
    (upheap size buf)))

;; buf[i]にある要素が、その子要素よりも大きい場合
;; 子要素との交換を繰り替えすことで、ヒープ内の適切な位置に再配置する
(defun downheap (i size buf &aux (obj (aref buf i)))
  (declare ((simple-array obj) buf)
           (fixnum i size))
  (loop WITH size/2 = (ash size -1)
        WHILE (<= i size/2) DO
    (let ((child (+ i i)))
      (when (and (< child size) 
                 (obj< (aref buf (1+ child)) (aref buf child)))
        (incf child))
      (when (obj< obj (aref buf child))
        (return))
      (setf (aref buf i) (aref buf child)
                      i            child)))
  (setf (aref buf i) obj))

(defun heap-pop (heap)
  (with-slots (size buf) (the heap heap)
    (prog1 (aref buf 1)
      (setf (aref buf 1) (aref buf size))
      (decf size)
      (downheap 1 size buf))))

;; コード出現頻度表を受け取り、ヒープを作成する
(defun make-heap (#1=code-frequency-table)
  (declare ((simple-array fixnum) #1#))
  (let ((heap (heap (count-if #'plusp #1#))))
    (with-slots (buf size) heap
      (loop FOR code FROM 0 BELOW (length #1#)
            FOR freq OF-TYPE fixnum = (aref #1# code)
            WHEN (plusp freq) 
            DO (setf (aref buf (incf size)) (make-code-obj :cost freq :code code)))
      
      (loop FOR k FROM (ash size -1) DOWNTO 1 
            DO (downheap k size buf)))
    heap))

;; ヒープの先頭から二つの要素を取り出す
(defun heap-pop2 (heap)
  (values (heap-pop heap) (heap-pop heap) (zerop (heap-size heap))))

;; ヒープに要素を追加した後、先頭から二つの要素を取り出す
(defun heap-push-and-pop2 (obj heap)
  (heap-push obj heap)
  (heap-pop2 heap))


使用例。
参照: read-binary-file

;; コードの出現頻度カウント関数
(defun count-frequency (file-path)
  (loop WITH freq = (make-array #x100 :element-type 'fixnum :initial-element 0)
        FOR c ACROSS (read-binary-file file-path) 
    DO (incf (aref freq c))
    FINALLY (return freq)))

;; 夏目漱石の『こころ』(UTF-8)のハフマン符号
> (calc-code->bitlength-table (count-frequency "/path/to/kokoro")) 
--> #(0 0 0 0 0 0 0 0 0 0 8 0 0 8 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 18 0 0 0 0 0
      0 0 0 0 0 0 0 11 17 17 16 15 16 16 16 16 17 16 16 16 19 0 0 0 0 0 0 0 0 0 0 0
      0 0 0 19 19 0 0 0 0 0 0 0 0 19 0 0 0 0 19 0 0 0 0 0 0 0 0 17 0 0 0 0 0 19 19
      19 18 0 0 19 0 18 18 0 17 0 17 19 0 17 0 19 19 0 0 0 0 0 4 3 4 8 6 8 7 8 7 7
      6 5 6 7 9 7 8 8 7 7 8 8 8 6 8 7 9 8 9 7 10 6 8 8 9 7 8 8 6 6 7 8 7 6 9 8 6 7
      8 9 9 9 10 9 9 9 7 9 7 8 8 8 7 8 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
      0 0 0 0 0 0 0 0 0 0 0 0 12 2 7 5 6 6 7 7 0 0 0 0 0 10 0 0 0 0 0 0 0 0 0 0 0 0
      0 0 0 0)


長さ制限付きのハフマン符号化との比較。
※ 名前が重複しているので、以下では制限付きの方の関数名を calc-code->bitlength-table-under-limitation に変更している
参照: 制限付きハフマン符号化ソース

;; 『こころ』の各コード出現頻度
> (progn (defvar *kokoro* (count-frequency "/path/to/kokoro")) *kokoro*)
--> #(0 0 0 0 0 0 0 0 0 0 1595 0 0 1595 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 2 0 0 0
      0 0 0 0 0 0 0 0 0 121 4 3 5 18 6 9 8 5 4 7 9 9 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0
      1 1 0 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 4 0 0 0 0 0 1 1 1 2 0 0 1 0 2
      2 0 2 0 3 1 0 3 0 1 1 0 0 0 0 0 24084 103681 33007 1685 8530 2272 4555 2104
      4775 4332 8394 14319 7112 3731 837 3780 1566 2207 3764 4618 2163 2699 1917
      7045 1381 3610 1362 2174 1302 3714 635 8924 2283 2015 722 4773 1841 2708 5828
      7531 4686 1812 5511 5978 841 1890 8909 5330 2120 943 779 1144 672 830 1250
      843 3336 1212 3870 1895 1898 1908 4180 1711 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
      0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 115 138338 5422 13041 8671 8378 5646 3492 0
      0 0 0 0 661 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0)

;; lsコマンドの各コード出現頻度
> (progn (defvar *ls* (count-frequency "/bin/ls")) *ls*)
--> #(11818 1221 510 310 2182 2347 195 236 3460 169 254 63 572 118 54 1193 339 72
      135 41 332 109 42 45 211 105 40 44 279 87 43 105 3023 56 28 40 2020 291 189
      126 113 171 47 81 247 341 162 70 173 274 50 67 233 46 38 27 106 315 73 157
      198 191 32 52 129 143 212 249 817 896 136 100 95 83 39 60 170 234 78 90 223
      38 87 202 217 553 130 113 138 62 45 164 293 533 118 250 144 400 141 290 325
      731 255 157 348 594 26 81 455 222 514 525 263 42 541 637 1279 780 208 196 132
      173 42 44 239 320 78 186 328 104 42 1173 611 1394 81 79 158 3527 19 2846 69
      929 46 28 792 36 25 36 108 361 44 36 83 48 37 26 56 180 50 42 83 198 28 65 57
      168 67 22 87 68 27 30 76 29 32 29 77 20 20 32 185 109 306 48 324 66 138 95
      207 260 87 76 892 302 181 435 215 24 325 1071 154 142 76 31 119 34 58 52 264
      74 194 134 105 32 46 61 304 119 65 210 134 71 63 160 344 157 176 271 463 308
      85 37 1160 717 190 531 490 46 45 57 317 81 107 89 278 55 172 164 474 102 201
      664 838 234 214 7819)

;;;;;;;;;;;;;;;;;
;; 制限なしハフマン
;;  今回の実装では、処理オーダーがO(NlogN)、使用スペースがO(N)
;; 『こころ』
> (time (dotimes (i 100) (huffman::calc-code->bitlength-table *kokoro*)))
Evaluation took:
  0.004 seconds of real time  ; 0.004s
  0.004000 seconds of total run time (0.004000 user, 0.000000 system)
  100.00% CPU
  10,540,573 processor cycles ; 1000万
  579,904 bytes consed        ; 580KB 

;; ls
> (time (dotimes (i 100) (huffman::calc-code->bitlength-table *ls*)))
Evaluation took:
  0.008 seconds of real time  ; 0.008s
  0.008001 seconds of total run time (0.008001 user, 0.000000 system)
  100.00% CPU
  26,334,389 processor cycles ; 2600万
  1,232,120 bytes consed      ; 1.2MB

;;;;;;;;;;;;;;;;;
;; 制限付きハフマン  ※ 最大で15ビット
;;  今回の実装では、処理オーダーがO(NL)、使用スペースが(NL) ※ Lは符号ビット長の最大値
;; 『こころ』
> (time (dotimes (i 100) (huffman::calc-code->bitlength-table-under-limitation *kokoro* 15)))
Evaluation took:
  0.027 seconds of real time  ; 0.027s
  0.024001 seconds of total run time (0.024001 user, 0.000000 system)
  [ Run times consist of 0.004 seconds GC time, and 0.021 seconds non-GC time. ]
  88.89% CPU
  84,188,382 processor cycles ; 8400万
  5,280,032 bytes consed      ; 5.3MB

;; ls
> (time (dotimes (i 100) (huffman::calc-code->bitlength-table-under-limitation *ls* 15)))
Evaluation took:
  0.059 seconds of real time   ; 0.059s
  0.060004 seconds of total run time (0.060004 user, 0.000000 system)
  [ Run times consist of 0.004 seconds GC time, and 0.057 seconds non-GC time. ]
  101.69% CPU
  189,369,133 processor cycles ; 1.9億
  12,390,920 bytes consed      ; 12.4MB

制限の有無で(速度、メモリ共に)一桁程度の差があるが、長さ制限付きハフマン符号化の方の性能も、実用上は十分。

*1:実際には、符号対象となる各コードの符号ビット長だけを保持するカノニカルハフマン符号化の実装となっている。ただし、ソースコード的にはどちらもほとんど変わらないので、特に区別をせずに扱うことにする。