图像信号处理一

实验一 图像变换

实验目的

学会对图像进行傅立叶等变换，在频谱上对图像进行分析，增进对图像频域上的感性认识，并用图像变换进行压缩。

实验内容

对Lena或cameraman图像进行傅立叶、离散余弦、哈达玛变换。在频域，对比他们的变换后系数矩阵的频谱情况，进一步，通过逆变换观察不同变换下的图像重建质量情况。

实验要求

实验采用获取的图像，为灰度图像，该图像每象素由8比特表示。具体要求如下：

（1）对图像进行傅立叶变换、获得变换后的系数矩阵；

（2）将傅立叶变换后系数矩阵的频谱用图像输出，观察频谱；

（3）通过设定门限，将系数矩阵中95%的（小值）系数置为0，对图像进行反变换，获得逆变换后图像；

（4）观察逆变换后图像质量，并比较原始图像与逆变后的峰值信噪比（PSNR）。

（5）对输入图像进行离散余弦、哈达玛变换，重复步骤1-5；

（6）比较三种变换的频谱情况、以及逆变换后图像的质量（PSNR）。

比较相同图像质量（PSNR）下，步骤（3）中三种变换的系数矩阵，最少需要保留的非零元素个数。一般而言，需要保留的非零元素越少，则将这种变换用于图像编码，越易于节省码率。

实验原理

傅里叶变换、离散余弦和哈达玛变换都可以减小图片像素之间的相关性，并把更多的能量压缩在低频段，这使图像压缩成为可能，在把变换后小于95%的值置为0后，仍可以恢复图像，但会有不同程度的失真。

实验结果与分析

实验结果：

图 1 傅里叶变换与逆变换

图 2 离散余弦变换与逆变换

图 3 哈达玛变换与逆变换

图像的质量（PSNR）:

	傅里叶逆变换后图像的质量（PSNR）	离散余弦逆变换后图像的质量（PSNR）	哈达玛逆变换后图像的质量（PSNR）
LENA512.bmp	33.0162dB	34.4000dB	28.5244dB
cameraman.bmp	34.7764dB	36.4848dB	27.7271dB

取最终PSNR = 18 dB，计算不同变换需要保留的点数：

	傅里叶变换需要保留的点数	离散余弦变换需要保留的点数	哈达玛变换需要保留的点数
LENA512.bmp	53个	40个	63个
cameraman.bmp	92个	63个	116个

结果分析：

在输入图像信噪比相同的情况下，经过离散余弦处理后的能量比傅里叶变换和哈达玛变换更加集中，这也是在把95%的小能量置零并经过离散余弦逆变换后，其图像的PSNR最高的原因。通过对比三种变换的PSNR可知，这三种变换方式中，离散余弦变换的效果最好，傅里叶变换效果次之，哈达玛变换的效果最差。同样也可又在保证相同PSNR的情况下，各种变换需要保留的点的个数看出，离散余弦变换需要保留的点个数最少，而哈达玛变换需要保留的点最多。

主要代码

主程序：

clc
clear all

gs = 0.95;
% 取最终PSNR = 27dB，计算需要保留的点数

LINA = imread('LENA512.bmp');
LINA = rgb2gray(LINA);
% 傅里叶变换LINAF
F=fft2(LINA);
Fs=fftshift(F);
LINAF=log10(log10(abs(Fs)+1));
imwrite(LINAF,'LINAF.bmp')

% 离散余弦LINAC
DCTI=dct2(double(LINA));%计算余弦变换并移位
DCTIA=abs(DCTI);
LINAC=log10(DCTIA+1);
imwrite(LINAC,'LINAC.bmp')

% 哈达玛变换LINAH
I = LINA;
I=im2double(I);
h1=size(I,1);      %图像的行
h2=size(I,2);      %图像的列
H1=hadamard(h1);       %Hadamard变换矩阵
H2=hadamard(h2);       %Hadamard变换矩阵
J=H1*I*H2/sqrt(h1*h2); %Hadamard变换
LINAH = J;
imwrite(LINAH,'LINAH.bmp')

% 通过设定门限，将系数矩阵中95%的（小值）系数置为0
Y1 = abs(F(:));
[Y,I]=sort(Y1);
jiezhi = Y(round(gs*length(Y1)));
LINAF0 = F;
LINAF0(abs(F)<jiezhi) = 0;
LINAf0 = uint8(ifft2((LINAF0)));
Fs=fftshift(LINAF0);
LINAF0=log10(log10(abs(Fs)+1));
imwrite(LINAF0,'LINAF0.bmp')
imwrite(LINAf0,'iLINAf0.bmp')
Diff = double(LINA)-double(LINAf0);
MSE=sum(Diff(:).^2)/(512*512);
PSNR_LINAf = 10*log10(255^2/MSE(1))

Y1 = DCTIA(:);
[Y,I]=sort(Y1);
jiezhi = Y(round(gs*length(Y)));
LINAC0 = DCTI;
LINAC0(DCTIA<jiezhi) = 0;
imwrite(LINAC0,'LINAC0.bmp')
LINAc0 = uint8(idct2(LINAC0));
imwrite(LINAc0,'iLINAc0.bmp')
Diff = double(LINA)-double(LINAc0);
MSE=sum(Diff(:).^2)/(512*512);
PSNR_LINAc = 10*log10(255^2/MSE(1))

Y1 = LINAH(:);
[Y,I]=sort(abs(Y1));
jiezhi = Y(round(gs*length(Y)));
LINAH0 = LINAH;
LINAH0(abs(LINAH)<jiezhi) = 0;
imwrite(LINAH0,'LINAH0.bmp')
I = LINAH0;
h1=size(I,1);      %图像的行
h2=size(I,2);      %图像的列
H1=hadamard(h1);       %Hadamard变换矩阵
H2=hadamard(h2);       %Hadamard变换矩阵
J=H1*I*H2/sqrt(h1*h2); %Hadamard变换
LINAh0 = uint8(J.*512/2);
imwrite(LINAh0,'iLINAh0.bmp')
Diff = double(LINA)-double(LINAh0);
MSE=sum(Diff(:).^2)/(512*512);
PSNR_LINAh = 10*log10(255^2/MSE(1))

CAM = imread('cameraman.bmp');
% 傅里叶变换CAMF
F=fft2(CAM);
Fs=fftshift(F);
CAMF=log10(log10(abs(Fs)+1));
imwrite(CAMF,'CAMF.bmp')

% 离散余弦CAMC
DCTI=dct2(double(im2gray(CAM)));%计算余弦变换并移位
DCTIA=abs(DCTI);
CAMC=log10(DCTIA+1);
imwrite(CAMC,'CAMC.bmp')

% 哈达玛变换CAMH
I = CAM;
I=im2gray(I);
I=im2double(I);
h1=size(I,1);      %图像的行
h2=size(I,2);      %图像的列
H1=hadamard(h1);       %Hadamard变换矩阵
H2=hadamard(h2);       %Hadamard变换矩阵
J=H1*I*H2/sqrt(h1*h2); %Hadamard变换
CAMH = J;
imwrite(CAMH,'CAMH.bmp')

% 通过设定门限，将系数矩阵中95%的（小值）系数置为0
Y1 = abs(F(:));
[Y,I]=sort(Y1);
jiezhi = Y(round(gs*length(Y1)));
CAMF0 = F;
CAMF0(abs(F)<jiezhi) = 0;
CAMf0 = uint8(ifft2((CAMF0)));
Fs=fftshift(CAMF0);
CAMF0=log10(log10(abs(Fs)+1));
imwrite(CAMF0,'CAMF0.bmp')
imwrite(CAMf0,'iCAMf0.bmp')
Diff = double(CAM)-double(CAMf0);
MSE=sum(Diff(:).^2)/(512*512);
PSNR_CAMf = 10*log10(255^2/MSE(1))

Y1 = DCTIA(:);
[Y,I]=sort(Y1);
jiezhi = Y(round(gs*length(Y)));
CAMC0 = DCTI;
CAMC0(DCTIA<jiezhi) = 0;
imwrite(CAMC0,'CAMC0.bmp')
CAMc0 = uint8(idct2(CAMC0));
imwrite(CAMc0,'iCAMc0.bmp')
Diff = double(CAM)-double(CAMc0);
MSE=sum(Diff(:).^2)/(512*512);
PSNR_CAMc = 10*log10(255^2/MSE(1))

Y = abs(CAMH(:));
[Y,I]=sort(im2gray(Y));
jiezhi = Y(round(gs*length(Y)));
CAMH0 = CAMH;
CAMH0(abs(CAMH)<jiezhi) = 0;
imwrite(CAMH0,'CAMH0.bmp')
I = CAMH0;
h1=size(I,1);      %图像的行
h2=size(I,2);      %图像的列
H1=hadamard(h1);       %Hadamard变换矩阵
H2=hadamard(h2);       %Hadamard变换矩阵
J=H1*I*H2/sqrt(h1*h2); %Hadamard变换
CAMh0 = uint8(J.*512/2);
imwrite(CAMh0,'iCAMh0.bmp')
Diff = double(CAM)-double(CAMh0);
MSE=sum(Diff(:).^2)/(512*512);
PSNR_CAMh = 10*log10(255^2/MSE(1))


% 取最终PSNR = 18dB，计算需要保留的点数

计算点数

clc
clear all

temPSNR_LINAc =29.5247;
    % 取最终PSNR = 27dB，计算需要保留的点数
for gs = 0.9998:0.00001:0.99999
    
    LINA = imread('LENA512.bmp');
    LINA = rgb2gray(LINA);
    
    DCTI=dct2(double(LINA));%计算余弦变换并移位
    DCTIA=abs(DCTI);
    LINAC=log10(DCTIA+1);
    imwrite(LINAC,'LINAC.bmp')

    % 通过设定门限，将系数矩阵中95%的（小值）系数置为0
    Y1 = DCTIA(:);
    [Y,I]=sort(Y1);
    jiezhi = Y(round(gs*length(Y)));
    LINAC0 = DCTI;
    LINAC0(DCTIA<jiezhi) = 0;
    imwrite(LINAC0,'LINAC0.bmp')
    LINAc0 = uint8(idct2(LINAC0));
    imwrite(LINAc0,'iLINAc0.bmp')
    Diff = double(LINA)-double(LINAc0);
    MSE=sum(Diff(:).^2)/(512*512);
    PSNR_LINAc = 10*log10(255^2/MSE(1))
    
    if (PSNR_LINAc<18)&&(temPSNR_LINAc>18)
    jg = (1-gs) * 512*512
    end
    temPSNR_LINAc = PSNR_LINAc;
end