2013年11月1日 星期五

Huffman Coding for DC & AC Coefficients of Grayscale JPEG Image

Joint Photographic Expert Group 是一種利用餘弦轉換 (Cosine Transform)的國際標準壓縮方法,詳細的方法可以參考 [1]www.w3.org/Graphics/JPEG/itu-t81.pdf


比較困難的部份並不是餘弦轉換,而是轉換後的AC值與DC值的壓縮編碼。今天
我們利用Matlab來實作一遍吧!


◎影像切塊成8乘8(partition into 8x8 blocks)
    聽起來這一步應該很簡單,唯一的問題是影像的尺寸如果不是8的倍數時,要如何處理?最簡單的方法就是在影像的最右邊及最下邊補灰度值為128的像素,讓寬度與高度都是8的倍數。


↓利用padarray將區塊補足至8X8

◎平移128 (shift by 128)
    將區塊的每一個灰度值減去128,由於灰度影像的class是uint8,像素值介於0與255之間,在減去128之前,需要影像格式轉為double。

↓轉換格式後再減去128


◎餘弦轉換(cosine transform)
由於Matlab 已有餘弦轉換的函式,直接套用就可以了

↓利用函式dct2作餘弦轉換


◎量化餘弦轉換後的值 (quantization)
    餘弦轉換後的值為實數,不利壓縮。所以將其量化,其實就是除以一個整數,並四捨五入成為一個整數。在JPG標準文件中,量化矩陣Q(其實就是除數)定義為

↓量化矩陣Q

↓餘弦轉換後的8X8區塊除以矩陣Q,再四捨五入

◎儲存量化區塊
    量化區塊的值分成兩類,一個叫DC值,也就是左上角第一個值,其它的數值稱之為AC值。
    DC值只有一個,所以將每個區塊的DC值依掃瞄線順序排列儲存。每個區塊有63個AC值,以Z字形(zigzag)順序儲存。
Z字形順序就是將區塊順時鐘轉45度,再以平常的掃瞄線方式循環排列。

↓zigzag在8X8區塊的順序

在matlab中,我們先以人工的方式建立8x8 zigzag序列的一維座標,再套入二維區塊中。以上圖為例,序號1在一維座標是9,序號2的一維座標是2,序號3的一維座標是3,序號4的一維座標是10,依此類推。matlab的二維座標轉一維座標是先算y軸,再算x軸,也就是 A(y,x) = A((x-1)*8+y)。 EX. A(1,2)=A((2-1)*8+1)=A(9)。


↓建立ZIGZAG的一維座標


↓將量化後的矩陣以zigzag的順序列出


DC值編碼
    由於影像被切成許多8X8區塊,每一個區塊有一個DC值,這個值在餘弦轉換佔有重要地位,解壓縮影像的畫值與其有密切關係。
    將每個DC值排成一維陣列,假設是DC1, DC2, ... DCn。第一個DC值不動,計算第2個DC值與第1個DC值的差,再計算第3個DC值與第2值DC值的差,...,依此類推。求得 DC1, (DC2-D1), (DC3-DC2), ..., (DCn-DCn-1)。
    由於JPG是個國際標準的壓縮方式,對DC值的Huffman Code已有現成表格可供使用。對於每一個DC差值,看它的範圍落在哪個Category,找出對應的Code word來使用。
    另一個小細節:第一個DC值(不是差值),仍然使用與差值相同的方法來編碼!



↓DC差值的Huffman編碼(取自[1]第149頁)

為了容易理解,把DC的差值加入上表中,可得下表

CategoryDC差值範圍Huffman Code編碼長度
00002
1-1, 10103
2-3, -2, 2, 30113
3-7~-4, 4~71003
4-15~-8, 8~151013
5-31~-16, 16~311103
6-63~-32, 32~6311104
7-127~-64, 64~127111105
8-255~-128, 128~2551111106
9-511~-256, 256~51111111107
10-1023~-512, 512~1023111111108
11-2047~-1023, 1023~20471111111109

因此DC差值(假設為DIFF)的編碼流程如下:
如果 DIFF=0 則輸出為 00
如果 DIFF>0 則
    依上表找出對應的Huffman Code(假設為HCode),並求出DIFF的二進位碼(假設為BDIFF),輸出 [HCode BDIFF];
如果DIFF<0 則
     依上表找出對應的Huffman Code(假設為HCode),並求出DIFF的絕對值的二進位碼(假設為BDIFF),輸出 [HCode ~BDIFF] (~代表0變1, 1變0, 也就是補數)。

舉例:

  • DIFF = 2 → Category=2→HCode = 011, BDIFF = 10, 所以輸出為[011 10] = 01110
  • DIFF = 12 →Category=4→HCODE = 101, BDIFF = 1100, 所以輸出為[101 1100] = 1011100
  • DIFF = -4 →Category=3→HCODE = 100, BDIFF = 100, ~BDIFF = 011, 所以輸出為[100 001]= 100011
  • DIFF = -100→Category=7→HCODE = 11110, BDIFF = 1100100, ~BDIFF = 0011011, 所以輸出為[100 001]= 111100011011

↓matlab 的程式碼如下


↓執行看看....結果與上述相同


AC值編號
在一個8X8的區塊中,DC值編碼與另外63個AC值編碼略有不同。利用AC值本身並不是AC值的差值來編號,另外,由於AC值有許多零,所以將零的個數納入編碼的因素,採用 run-length及Huffman Code的方式來編碼。

因此決定編碼的索引值就是AC值有幾個零!AC值的編碼與DC差值編碼使用的方法完全一樣,在JPG標準中,以Size來代表AC值的範圍,其範圍如下表:

↓AC值的Size 與DC值的Category意義相同

Size
AC值範圍
1
-1, 1
2
-3, -2, 2, 3
3
-7~-4, 4~7
4
-15~-8, 8~15
5
-31~-16, 16~31
6
-63~-32, 32~63
7
-127~-64, 64~127
8
-255~-128, 128~255
9
-511~-256, 256~511
10=A
-1023~-512, 512~1023

上表與DC差值的對應表幾乎一樣,一個是沒有Size=0,一個是只考慮到第10類,沒有第11類。這是由於AC值為零使用run-length的概念來取代,而且AC值的範圍也沒有DC值這麼廣。

另一個索引值就是零的個數,以Run來代表。在8X8區塊的灰度影像中,JPG標準最多只針對15個連續的零來編碼,因此Size的範圍就是0至15。如果有連續16個零就以ZRL來編碼。

因此(Run/Size)的範圍就是從(0/1)至(15/10) = (F/A),對應表項目很多,只列出一小段。


↓AC值(侷部)的Huffman編碼(取自[1]第150頁),共四頁


在這個很長很長的表格中,有兩個特殊的編碼,一個是0/0,另一個是F/0。0/0代表區塊結束(End Of Block EOB),F/0代表前面有15個零,而且自己也是零(ZRL, 共有16個零)。其它的編碼都沒有Size=0的索引值。

現在開始來寫AC值的編碼流程:

輸入值:AC(1:63) 是一個具有63個AC值的一維陣列,而且已按Zigzag順序排列。

首先找出AC陣列中最後一個不為零的序號值K,即AC(K)不為零,AC(K+1)至AC(63)都是零。
接著依序從i=1至K處理AC陣列的值
B=AC(i)
如果B大於零
    依上表求出其size值,設為size(B);
    求出前面有幾個連續的零,設為zero(B);
    由 zero(B)/size(B)=run/size 查表找出對應的Huffman Code, 設為HC
   輸出 HC及B的二進位碼, 並將零的個數歸零;
如果B小於零
    依上表求出其size值,設為size(B);
     求出前面有幾個連續的零,設為zero(B);
    由 zero(B)/size(B)=run/size 查表找出對應的Huffman Code, 設為HC
   計算絕對值B的二進位碼並計算此二進位的補數,輸出 HC及補數,並將零的個數歸零;
如果B為零
    如果前面零的個數小於15,將零的個數加1
    如果前面零的個數等於15,輸出ZRL的代碼,並將零的個數歸零。
依序完成以上K個AC值,最後輸出EOB的代碼

matlab的程式如下:
But 程式碼太長了,分段顯示


↓程式碼的上半段

上圖中的程式中,Hcode是一個二維的cell,儲存(run/size)對應的Huffman 的編號。利用find尋找最後一個非零的AC值,並將其索引值儲存至K。如果K是空矩陣,代表63個值都是零,就輸出EOB,並結束程式。


↓程式碼的下半段



現在來試試看吧!


↓先SHOW一下AC值的內容

63個數字若排成一行,不易閱讀,利用reshape,用一個9x7的陣列方式顯示。由於matlab儲存二維陣列的順序是先Y軸再X軸,所以在AC一維陣列的順序是7 2 1 -3 -4 -1 1 -2 1 0....最後一個非零的值是-1,位置是25 (AC(25)=-1)。



↓執行結果

解釋一下:

  • AC(1)=7 → index=2, NOZERO=0 →Hcode{0+1,2+1}=100, dec2bin(7)=111→輸出100111
  • AC(2)=2 →index=1, NOZERO=0 →Hcode{1,2}=01, dec2bin(2)=10 →輸出0110
  • AC(4)=-3 → index=1, NOZERO=0 →Hcode{1,2}=01,  dec2bin(abs(-3))=11 補數00→輸出0100
  • AC(11)=2→index=1, NOZERO=1 →Hcode{2,2}=11011, 
  • dec2bin(2)=10 →輸出1101110
  • AC(18)=1index=0, NOZERO=3 →Hcode{4,1}=111010, 
  • dec2bin(1)=1 →輸出1110101
  • AC(24)=-1 →index=0, NOZERO=3 →Hcode{4,1}=111010,  dec2bin(abs(-1))=1 補數0→輸出1110100




↓再雞婆一點,每一個AC(K)及其編碼如下表

上表中,空白部份 (K=10, 15-17, 21-23)代表AC值為零。

大功告成, 累~~~
對了!JPG有幾種變形,以上是簡單的一種,稱之為BASELINE。還有彩色的...漸進壓縮...無損壓縮....更累~~~
--- END ---

4 則留言:

  1. 作者已經移除這則留言。

    回覆刪除
  2. 大大您好:
    請問在DC值編碼的時候
    為什麼不能對負數做一個編碼
    跑出的error是

    Improper index matrix reference.
    Error in getdc (line 19)
    mycode=strcat(Hcode(index+2).bdc);

    不知道是什麼問題
    能請大大做一個解釋嗎!!?? 謝謝你

    回覆刪除
  3. 哦哦哦哦有這邊文章真是太感謝了~

    回覆刪除
  4. 非常感謝,目前正在大學修相關課程,這篇非常有幫助!!

    回覆刪除