DCT變換、DCT反變換、分塊DCT變換
歡迎轉載，但請注明出處！
一、引言
DCT變換的全稱是離散余弦變換(Discrete Cosine Transform)，主要用於將數據或圖像的壓縮，能夠將空域的信號轉換到頻域上，具有良好的去相關性的性能。DCT變換本身是無損的，但是在圖像編碼等領域給接下來的量化、哈弗曼編碼等創造了很好的條件，同時，由於DCT變換時對稱的，所以，我們可以在量化編碼后利用DCT反變換，在接收端恢復原始的圖像信息。DCT變換在當前的圖像分析已經壓縮領域有着極為廣大的用途，我們常見的JPEG靜態圖像編碼以及MJPEG、MPEG動態編碼等標准中都使用了DCT變換。
二、一維DCT變換
一維DCT變換時二維DCT變換的基礎，所以我們先來討論下一維DCT變換。一維DCT變換共有8種形式，其中最常用的是第二種形式，由於其運算簡單、適用范圍廣。我們在這里只討論這種形式，其表達式如下：
其中，f(i)為原始的信號，F(u)是DCT變換后的系數，N為原始信號的點數，c(u)可以認為是一個補償系數，可以使DCT變換矩陣為正交矩陣。
三、二維DCT變換
二維DCT變換其實是在一維DCT變換的基礎上在做了一次DCT變換，其公式如下：
由公式我們可以看出，上面只討論了二維圖像數據為方陣的情況，在實際應用中，如果不是方陣的數據一般都是補齊之后再做變換的，重構之后可以去掉補齊的部分，得到原始的圖像信息，這個嘗試一下，應該比較容易理解。
另外，由於DCT變換高度的對稱性，在使用Matlab進行相關的運算時，我們可以使用更簡單的矩陣處理方式：
接下來利用Matlab對這個過程進行仿真處理：
1 clear; 2 clc; 3 X=round(rand(4)*100) %產生隨機矩陣 4 A=zeros(4); 5 for i=0:3
6 for j=0:3
7 if i==0
8 a=sqrt(1/4); 9 else
10 a=sqrt(2/4); 11 end 12 A(i+1,j+1)=a*cos(pi*(j+0.5)*i/4); 13 end 14 end 15 Y=A*X*A' %DCT變換 16 YY=dct2(X) %Matlab自帶的dct變換
運行結果為：
1 X =
2
3 42 66 68 66
4 92 4 76 17
5 79 85 74 71
6 96 93 39 3
7
8
9 Y =
10
11 242.7500 48.4317 -9.7500 23.5052
12 -12.6428 -54.0659 7.4278 22.7950
13 -6.2500 10.7158 -19.7500 -38.8046
14 40.6852 -38.7050 -11.4653 -45.9341
15
16
17 YY =
18
19 242.7500 48.4317 -9.7500 23.5052
20 -12.6428 -54.0659 7.4278 22.7950
21 -6.2500 10.7158 -19.7500 -38.8046
22 40.6852 -38.7050 -11.4653 -45.9341
由上面的結果我們可以看出，我們采用的公式的方法和Matlab自帶的dct變化方法結果是一致的，所以驗證了我們方法的正確性。
如果原始信號是圖像等相關性較大的數據的時候，我們可以發現在變換之后，系數較大的集中在左上角，而右下角的幾乎都是0，其中左上角的是低頻分量，右下角的是高頻分量，低頻系數體現的是圖像中目標的輪廓和灰度分布特性，高頻系數體現的是目標形狀的細節信息。DCT變換之后，能量主要集中在低頻分量處，這也是DCT變換去相關性的一個體現。
之后在量化和編碼階段，我們可以采用“Z”字形編碼，這樣就可以得到大量的連續的0，這大大簡化了編碼的過程。
四、二維DCT反變換
在圖像的接收端，根據DCT變化的可逆性，我們可以通過DCT反變換恢復出原始的圖像信息，其公式如下：
同樣的道理，我們利用之前的矩陣運算公司可以推導出DCT反變換相應的矩陣形式：
下面我們用Matlab對這個過程進行仿真：
1 clear; 2 clc; 3 X=[ 4 61 19 50 20
5 82 26 61 45
6 89 90 82 43
7 93 59 53 97] %原始的數據 8 A=zeros(4); 9 for i=0:3
10 for j=0:3
11 if i==0
12 a=sqrt(1/4); 13 else
14 a=sqrt(2/4); 15 end 16 A(i+1,j+1)=a*cos(pi*(j+0.5)*i/4); %生成變換矩陣 17 end 18 end 19 Y=A*X*A' %DCT變換后的矩陣 20 X1=A'*Y*A %DCT反變換恢復的矩陣
運行結果為：
1 X =
2
3 61 19 50 20
4 82 26 61 45
5 89 90 82 43
6 93 59 53 97
7
8
9 Y =
10
11 242.5000 32.1613 22.5000 33.2212
12 -61.8263 7.9246 -10.7344 30.6881
13 -16.5000 -14.7549 22.5000 -6.8770
14 8.8322 16.6881 -35.0610 -6.9246
15
16
17 X1 =
18
19 61.0000 19.0000 50.0000 20.0000
20 82.0000 26.0000 61.0000 45.0000
21 89.0000 90.0000 82.0000 43.0000
22 93.0000 59.0000 53.0000 97.0000
我們可以看到反變換后無損的恢復了原始信息，所以證明了方法的正確性。但是在實際過程中，需要量化編碼或者直接舍棄高頻分量等處理，所以會出現一定程度的誤差，這個是不可避免的。
五、分塊DCT變換
在實際的圖像處理中，DCT變換的復雜度其實是比較高的，所以通常的做法是，將圖像進行分塊，然后在每一塊中對圖像進行DCT變換和反變換，在合並分塊，從而提升變換的效率。具體的分塊過程中，隨着子塊的變大，算法復雜度急速上升，但是采用較大的分塊會明顯減少圖像分塊效應，所以，這里面需要做一個折中，在通常使用時，大都采用的是8*8的分塊。
Matlab的 blkproc 函數可以幫我們很方便的進行分塊處理，下面給出我們的處理過程：
1 clear; 2 clc; 3
4 X=imread('pepper.bmp'); 5 X=double(X); 6 [a,b]=size(X); 7 Y=blkproc(X,[8 8],'dct2'); 8 X1=blkproc(Y,[8 8],'idct2'); 9
10 figure 11 subplot(1,3,1); 12 imshow(uint8(X)); 13 title('原始圖'); 14
15 subplot(1,3,2); 16 imshow(uint8(Y)); 17 title('分塊DCT變換圖'); 18
19 subplot(1,3,3); 20 imshow(uint8(X1)); 21 title('分塊DCT恢復圖'); 22
23 Y1=dct2(X); 24 X10=idct2(Y1); 25
26 figure 27 subplot(1,3,1); 28 imshow(uint8(X)); 29 title('原始圖'); 30
31 subplot(1,3,2); 32 imshow(uint8(Y1)); 33 title('DCT變換圖'); 34
35 subplot(1,3,3); 36 imshow(uint8(X10)); 37 title('DCT反變換恢復圖');
運行結果為：
從圖中，我們可以明顯看出DCT變換與分塊DCT變換在使用時的區別。
六、小結
DCT、DWT等是圖像處理的基礎知識，之前一直有用到，但是沒怎么好好整理下，今天在做稀疏編碼的時候正好有用到，就順便整了下，希望能夠給后來者一些提示。
我的新浪微博：http://weibo.com/3109428257/profile?rightmod=1&wvr=5&mod=personinfo
Reference:http://wuyuans.com/2012/11/dct2/