基于Matlab的帧内DPCM的编解码方法
DIP实验5:DPCM
实验目的
掌握帧内DPCM的编解码方法(预测、量化)及其压缩性能。
实验内容
1)编写一维DPCM编解码程序;
2)编写二维DPCM编解码程序(预测系数[1/2,1/8,1/4,1/8]);
3)分别计算上两种情况中重建图像的PSNR。其中,重建图像和原图像的差值用绝对误差表示,可设置放大因子n以便观察误差分布情况。量化器采用15个量化分层,量化电平分别为(0, 5, 10, 17, 28, 39, 52, 67)。实验图像为cla0。
4)当人为加入单个误码,可观察到怎样的误码图案?
输出图像排列格式如下:
原图像 | 一维DPCM重建图像 | 二维DPCM重建图像 |
---|---|---|
一维DPCM误码图案 | 二维DPCM误码图案 |
需要编写的函数:
1)发送端:线性预测编解码器函数
function [newbuf,Qerror]=Dpcm_code(oldbuf,M,N,dim)
2)求预测值函数
function Pvalue=Predict_Value(row,col,N,dim)
3)求量化值函数
function Qvalue=Quant_Value(err)
4)求重建值函数
function Rvalue=Restor_Value(quan_err,pre_val)
5)限幅函数
function Cvalue=Clip_Value(res_val)
6)接收端:传输误码解码器函数
function newbuf=Error_Code(M,N,dim)
7)峰值信噪比
附:计算预测值,预测公式如下:
- 1_D DPCM:
if j=1(若像素在图像左边第1列,则预测值设为128) ^x[i,j]=128
if j>1 ^x[i,j]=x'[i,j-1]
- 2_D DPCM:
if i=1,j=1
^x[i,j]=128
if i=1,j>1
^x[i,j]=x'[ i,j-1]
if i>1,j=1 or 255
^x[i,j]=x'[i-1,j ]
if i>1,j>1 (Pirsch's predictor)
^x[i,j]=1/2x'[ i,j-1]+1/8x'[i-1,j-1]+1/4x'[i-1,j ]+1/8x'[i-1,j+1]
参考代码
% 预测编码
OI = imread('CLA1.bmp');
subplot(2,3,1);
imshow(OI);
title("原图像");
OI = double(OI);
[M,N] = size(OI);
[newImg1,Qerror1] = Dpcm_code(OI,M,N,1);
newImg1 = uint8(newImg1);
subplot(2,3,2);
imshow(newImg1);
title("一维DPCM重建图像")
% 计算峰值信噪比
psnr1 = calPsnr(OI,newImg1)
[newImg2,Qerror2] = Dpcm_code(OI,M,N,2);
newImg2 = uint8(newImg2);
subplot(2,3,3);
imshow(newImg2);
title("二维DPCM重建图像")
% 计算峰值信噪比
psnr2 = calPsnr(OI,newImg2)
newImg=Error_Code(Qerror1,M,N,1);
newImg = uint8(newImg);
subplot(2,3,5);
imshow(newImg);
title("一维DPCM误码图案")
newImg=Error_Code(Qerror2,M,N,2);
newImg = uint8(newImg);
subplot(2,3,6);
imshow(newImg);
title("二维DPCM误码图案")
%1)发送端:线性预测编解码器函数
function [newImg,Qerror]=Dpcm_code(oldbuf,M,N,dim)
global newbuf;
newbuf = zeros(M,N); % 预测值矩阵,定义为全局变量,因为其他函数中需要用到
Qerror = zeros(M,N);
for i = 1:M
for j = 1:N
preVal = Predict_Value(i,j,N,dim); % 得到该像素点的预测值
error = oldbuf(i,j) - preVal; % 与原来的像素值相减得到预测误差
Qerror(i,j) = Quant_Value(error); % 将预测误差量化
rValue = Restor_Value(Qerror(i,j),preVal); % 用预测值加上量化后的误差
newbuf(i,j) = Clip_Value(rValue); % 限幅函数
end
end
newImg = newbuf;
end
% 2)求预测值函数
function Pvalue=Predict_Value(row,col,N,dim)
global newbuf;
% dim 为1时是一维预测,dim 为2时是二维预测
switch dim
case 1
if col == 1 % 图像的第一列没有前一个预测像素值
Pvalue = 128;
else
Pvalue = newbuf(row,col-1);
end
case 2
if row == 1 && col == 1 % (1,1)像素
Pvalue = 128;
end
if row == 1 && col > 1 % 第一行,列数大于1的像素,它的预测值只有前一个像素
Pvalue = newbuf(row,col-1);
end
if (row > 1 && col == 1) || (row > 1 && col == N) % 第一列或最后一列,且行数大于1
Pvalue = newbuf(row-1,col);
end
if row > 1 && col > 1 && col < N % 其余情况,四个像素值加权
Pvalue = (1/2)*newbuf(row,col-1) + (1/8)*newbuf(row-1,col-1) + (1/4)*newbuf(row-1,col) + (1/8)*newbuf(row-1,col+1);
end
end
end
% 3)求量化值函数
% 量化电平分别为(0,5,10,17,28,39,52,67) 范围对应(3,11,21,32,44,58,66)
function Qvalue=Quant_Value(err)
% 判断量化电平
if abs(err) <= 3
level = 0;
elseif abs(err) <= 11
level = 5;
elseif abs(err) <= 21
level = 10;
elseif abs(err) <= 32
level = 17;
elseif abs(err) <= 44
level = 28;
elseif abs(err) <= 58
level = 39;
elseif abs(err) <= 66
level = 52;
else
level = 67;
end
if err ~= 0
Qvalue = level * (err/abs(err)); % 量化值的正负
else
Qvalue = 0;
end
end
% 4)求重建值函数
function Rvalue=Restor_Value(quan_err,pre_val)
Rvalue = quan_err + pre_val;
end
% 5)限幅函数
function Cvalue=Clip_Value(res_val)
if res_val < 0
Cvalue = 0;
elseif res_val > 255
Cvalue = 255;
else
Cvalue = res_val;
end
end
% 6)接收端:传输误码解码器函数
function newImg=Error_Code(qerror,M,N,dim)
global newbuf;
newbuf = zeros(M,N);
qerror(50,50) = 100; % 将某一量化误差更改
for i = 1:M
for j = 1:N
preVal = Predict_Value(i,j,N,dim);
rValue = Restor_Value(qerror(i,j),preVal);
newbuf(i,j) = Clip_Value(rValue);
end
end
newImg = newbuf;
end
% 7)峰值信噪比
function PSNR = calPsnr(img1,img2)
img1 = double(img1);
img2 = double(img2);
[M,N] = size(img1);
% 峰值信噪比公式
MES = sum(sum((img1-img2).^2))/(M*N);
PSNR = 20 * log10(255/sqrt(MES));
end