I=imread('lenna.bmp');
figure(1);
imshow(I);
title('原始图像')
I1=double(I);%图像存储类型由整型转换为双精度型
T=dctmtx(8); %离散余弦变换矩阵
M=blkproc(I1,[8 8],'P1*x*P2',T,T'); %对原图像进行DCT变换
figure(2);
imshow(uint8(M));
title('DCT变换后的图像')
Q=[16 11 10 16 24 40 51 61
12 12 14 19 26 58 60 55
14 13 16 24 40 57 67 56
14 17 22 29 51 87 80 62
18 22 37 56 68 109 103 77
24 35 55 64 81 104 113 92
49 64 78 87 103 121 120 101
72 92 95 98 112 100 103 99 ];%QF=50的标准量化矩阵
N=blkproc(M,[8 8],'round(x./P1)',Q);%DCT变换后的矩阵除量化因子（四舍五入）
figure(3);
imshow(uint8(N));
title('量化压缩后的图像')
for k=1:8:249
for m=1:8:249
n=1;p=1;
while n<=8
for i=k:k+7
while p<=8
for j=m:m+7
f(n,p)=N(i,j);
p=p+1;
if p+1>9
break
end
end
end
p=1;
n=n+1 ;
if n+1>9
break
end
end
end
zig=[ %zig扫描顺序表
0, 1, 8, 16, 9, 2, 3, 10, ...
17, 24, 32, 25, 18, 11, 4, 5, ...%...表示zig是一维数组，而不是矩阵
12, 19, 26, 33, 40, 48, 41, 34, ...
27, 20, 13, 6, 7, 14, 21, 28, ...
35, 42, 49, 56, 57, 50, 43, 36, ...
29, 22, 15, 23, 30, 37, 44, 51, ...
58, 59, 52, 45, 38, 31, 39, 46, ...
53, 60, 61, 54, 47, 55, 62, 63
];
zig = zig + 1; % 下标加1，MATLAB的下标从1开始
a1 = reshape(f,1,64); % 将输入块变成1x64的向量
b= a1(zig); % 对a1按照查表方式取元素，得到zig-zag扫描结果
fzig=[1, 2, 6, 7, 15 ,16, 28, 29, ... % 反Z扫描顺序表
3, 5 ,8, 14, 17, 27, 30, 43,...
4, 9, 13, 18, 26, 31, 42, 44,...
10, 12, 19, 25, 32,41, 45, 54,...
11, 20, 24, 33, 40, 46, 53, 55,...
21, 23, 34, 39, 47, 52, 56, 61,...
22, 35,38, 48, 51, 57, 60,62,...
36, 37, 49, 50, 58, 59, 63, 64,];
fa1=b(fzig); %反Z变换
fa=reshape(fa1,8,8);% 将1x64的向量变成8x8的块输出
for s=1:8;
for q=1:8
B(s+k-1,q+m-1)=fa(s,q);
end
end
end
end
B1=blkproc(B,[8 8],'x.*P1',Q);
figure(4);
imshow(uint8(B1));% 显示反量化处理的图象
title('反量化后的图像')
B2=blkproc(B1,[8 8],'P1*x*P2',T',T);
B3=uint8(B2);%数据类型反变换，由双精度型转换为整型
figure(5);
imshow(B3);% 显示恢复后的图象
title('解压缩后（反DCT变换）的图像')
t1=(I-B3).^2; % I是原始矩阵 B3是反DCT 变换后的矩阵
t2=sum(sum(t1));
q=255;m=255;n=255;
MSE=t2/(m*n); %计算反DCT变换后图形的均方误差
PSNR=10*log10(q*q/MSE); %计算反DCT变换后峰值信噪比
MSE
PSNR
评论0