Matlab数字图像
（数字图像处理课程实验整理）
1、Matlab基本操作
1.1 清除命令窗口的显示
1.2 保存矩阵到指定文件
1.3 清除变量及其存储
1.4 读入文件
1.5 关于矩阵的一些命令
1.6 Matlab的一些常用命令
1.7 一些注意事项
2.图像处理基本操作
2.1 图像处理工具箱
2.2 图像文件的基本操作
3.图像几何变换基本操作
3.1 图像变换之圆柱形
3.2 图像变换之球形
3.3 图像变形切片数
4.图像的简单运算
4.1 图像的逻辑运算
4.2 图像平移
4.3 图像翻转
4.4 图像旋转
4.5 图像缩放
5.图像频率域操作
5.1 傅里叶变换
5.2 均值滤波
5.3 中值滤波
6.图像空间域操作
6.1 均值滤波器
6.2 中值滤波器
7.空间域滤波器
7.1 平滑滤波器
7.2 锐化滤波器
7.3 均值滤波、中值滤波、拉普拉斯子锐化、Robert交叉梯度算子
8.图像复原
8.1 逆滤波
8.2 均值滤波
8.3 统计排序滤波
9.图像分割
一、Matlab基本操作
1.清除命令窗口的显示
clc //清除显示窗口的内容
2.保存矩阵到指定文件中
save 文件 矩阵结果
eg.save A.mat A //将矩阵A中的元素保存在文件A.mat中
图1.1 将矩阵保存文件中
3.清除变量及其存储
clear //清除所有已存储的变量（即显示在工作区中的变量）
图1.2 清除内存
4.读入文件
load A.mat //将文件A.mat的元素读入
5.关于矩阵的一些命令
（1）创建特殊矩阵
A=zeros(n); //生成n*n的全0矩阵
A=zeros(m,n); //生成M*N的全0矩阵
A=zeros(a1,a2,a3……); //生成a1*a2*a3……的全0矩阵
A=zero(size(B)); //生成与矩阵B大小相同的全0矩阵
A=ones(n); //生成n*n的全1矩阵
A=ones(m,n); //生成m*n的全1矩阵
A=ones([m,n]); //生成m*n的全1矩阵
A=ones(a1,a2,a3……); //生成a1*a2*a3……的全1矩阵
A=ones(sinze(B)); //生成与矩阵B大小相同的全1矩阵
A=eye(n); //生成n*n的单位矩阵
…………（同上使用方法）
A=rand(n); //生成n*n的随机矩阵
…………（同上使用方法）
图1.3 矩阵的一些操作
（2）创建特定矩阵
A=[1,2,3;4,5,6];
A=[1 2 3;4 5 6];
//生成2*3的矩阵，相同列之间的数字可以用数字隔开也可以用逗号隔开，不同行之间必须用分号隔开
图1.4 创建特定矩阵
（3）取出矩阵中的特定行特定列元素
eg.A(2,2) //取出矩阵A中第二行第二列的元素
图1.5 取出特定元素
（4）利用矩阵创建矩阵
diag(A); //返回矩阵A对角元素成列向量 diag(v); //以向量v作为对角元素创建矩阵 flipud(A); //矩阵上下翻转 fliplr(A); //矩阵左右翻转 rot90(A); //矩阵逆时针翻转90° tril(A); //提取矩阵A的下三角矩阵 triu(A); //提取矩阵A的上三角矩阵
图1.6 利用矩阵创建矩阵
（5）寻访矩阵中的数据
x=[1 2 3 4 5]; //产生1*5的数组 b=x(3); //寻访数组x的第三个元素 b=x([1 2 5]); //寻访x的第1、2、5元素 b=x(1:3); //寻访前三个元素 b=x(3:end); //寻访第三个元素到最后一个 b=x(find(x>3)); //由大于3的元素构成数组 b=x(3:-1:1); //由前三个数倒排成数组 b=x([1 2 3 4 4 3 2 1]); //对元素的重复访问
图1.7 寻访矩阵中的元素
（6）修改数据
x=[1,2,3;4,5,6;7,8,9]; //产生3×3的矩阵 b=x(find(x>4)) ; //大于4的元素构成数组 x(2,2)=10; //第2行第2列值改为10 x(3,:)=20 ; //第3行的值改为20 x=ones(5) ; //产生5×5全为1的矩阵 x(2:4,2:4)=0 ; //第2、3、4行的第2、3、4列元素改0
图1.8 修改数据
（7）插入、提取、拉长及置空
插入：x=4:6; A=[x-3;x;x+3] %插入创建新矩阵
提取：C=A(1:2,2:3) %提取A的部分，1:2为行数，2:3为列数
拉长：D=A(:) %拉长成列
置空：A(:,2)=[] %删除A的第二列
图1.9 矩阵拉长等操作
（8）矩阵的规模
A=[3:6;1:4] %产生2×4的矩阵
s=size(A) %返回s为[行数 列数]
[r,c]=size(A) %返回r为行数，c为列数
r=size(A,1) %只返回行数
c=size(A,2) %只返回列数
n=length(A) %返回行数、列数最大者
图1.10 矩阵的规模
6.Matlab的一些常用命令
clear； //工作空间中清除所有的变量 clc； //清除命令窗口中显示内容 clf； //清除图形窗口内容 who； //列出当前工作空间中的变量 whos； //列出当前工作空间中的变量及信息 help<命令名> //获得在线帮助
图1.11 常用命令
7.一些注意事项
（1）M文件命名时不要使用纯数字，这样会导致错误的结果
（2）函数文件保存一定要和函数名同名
（3）命令行窗口中，不加分号执行完以后会显示结果，加分号不会显示出最终的结果
二、图像处理基本操作
Matlab之图像处理（写入/保存：imwrite/save test.mat）
1.图像处理工具箱
（1）图像处理工具箱支持4种图像类型：
- 真彩色图像
- 索引色图像
- 灰度图像
- 二值图像
（2）图像序列
A=cat(1,A1,A2,A3,A4,A5) //一维，纵向排列 A=cat(2,A1,A2,A3,A4,A5) //一维，横向排列 A=cat(3,A1,A2,A3,A4,A5) //一维，竖向排列
2.图像文件的基本操作
imread(filename,fmt) //读取文件格式的图像
imwrite(A,map,flename,fmt) //写入图形文件格式的图像
imfinfo('autumn.tif') //获取图像的信息
load\save //以Mat文件加载保存矩阵数据
imshow(I) //显示加载到Matlab中的图像
imshow(I,n)
imshow(I,[low high]) //[low high]为图像数据的值域
说明：
- 用imread读入的图像会以256*256*3的矩阵大小进行存储
- imfinfo函数所获取的信息至少包含：文件名、文件格式、文件格式的版本号、文件修改时间、文件的字节大小、图像的宽度、图像的长度(像素）、每个像素的位数、图像类型
- 对于double型数组为[0,1]，对于uint8型数组为[0,255]
- 直接从磁盘显示图像：imshow('rice.png')
- 可以在显示图像之前用输入figure;这样就可以在matlab显示图像时显示多张
subplot(2,2,1);imshow(a); //创建窗口为2*2的矩阵，在第一个位置显示图像a subplot(2,2,2);imshow(b); //在第二个位置显示图像b
效果：
图2.1 结果
I=imread('lena.jpg');
J1=imrotate(I,45); J2=imrotate(I,45,'crop'); //将图像旋转45°
subplot(131);imshow(I);
subplot(132);imshow(J1);
subplot(133);imshow(J2);
效果：
图2.2 结果
《截止到幻灯片91页》
三、图像几何变换基本操作
1.图像的基本操作
imfinfo('文件名称'); %查询图像信息
结果显示：
图3.1 结果
imshow(I) //灰度图像
imshow(I,n)
imshow(I,[low high])
imshow(BW) //二值图像
imshow(RGB) //RGB图像
结果显示：
图3.2 结果显示
[x,y,z]=cylinder;%创建一个圆柱面
I=imread('lena.jpg');%读取一个真彩色图像
warp(x,y,z,I);%把图像映射成圆柱体的表面纹理
结果显示：
图3.3 结果显示
[x,y,z]=sphere(50);%创建一个圆面
I=imread('lena.jpg');%读取一个真彩色图像
warp(x,y,-z,I);%把图像映射成圆柱体的表面纹理
结果显示：
图3.4 图像变球形
sphere(n); //表示生成的切片数
结果显示：
图3.5 切片显示
四、图像的简单运算
1.对简单图像进行逻辑运算
常见的图像逻辑运算有与(AND)、或(OR)、非(NOT)
样例：
A=zeros(128);A(40:67,60:100)=1;
figure
figure(1);imshow(A);
B=zeros(128);B(50:80,40:70)=1;
figure(2);imshow(B);
C=and(A,B);figure(3);imshow(C);
D=or(A,B);figure(4);imshow(D);
E=not(A);figure(5);imshow(E);
图4.1 结果显示1
图4.2 结果显示2
2.图像的平移
样例：
i=imread('moon.tif');
subplot(121);
imshow(i);
title('原始图像');
[M,N]=size(i);
g=zeros(M,N);
a=70;b=70;
for is=1:M
for j=1:N
if((is-a>0)&(is-a<M)&(j-b>0)&(j-b<N))
g(is,j)=i(is-a,j-b);
else
g(is,j)=0;
end
end
end
subplot(122);imshow(uint8(g));title('平移后的图像');
图4.3 图像的平移
3.图像的翻转
样例
i=imread('moon.tif');
subplot(121);
imshow(i);
title('原始图像');
[M,N]=size(i);
g=zeros(M,N);
a=70;b=70;
for is=1:M
for j=1:N
g(is,j)=i(M-is+1,j);
end
end
subplot(122);imshow(uint8(g));title('翻转后的图像');
图4.4 图像的翻转
4.图像的旋转
样例：
i=imread('moon.tif');j=imrotate(i,60,'bilinear');
k=imrotate(i,60,'bilinear','crop');
subplot(1,3,1),imshow(i);subplot(1,3,2),imshow(j);subplot(1,3,3),imshow(k);
图4.5 图像的旋转
5.图像的缩放
样例:
i=imread('moon.tif');
j1=imresize(i,10,'nearest');
j2=imresize(i,10,'bilinear');
j3=imresize(i,10,'bicubic');
subplot(2,2,1);
imshow(i);
title('原始图像');
subplot(2,2,2);
imshow(j1);
title('最近邻法');
subplot(2,2,3);
imshow(j2);
title('双线性插值法');
subplot(2,2,4);
imshow(j3);
title('三次内插法');
图4.6 图像的缩放
五、图像频率域操作
1.傅里叶变换
clear
clc
img=imread('peppers.png');
subplot(2,2,1);imshow(img);title('原图');
f=rgb2gray(img); %对于RGB图像必须做的一步，也可以用im2double函数
F=fft2(f); %傅里叶变换
F1=log(abs(F)+1); %取模并进行缩放
subplot(2,2,2);imshow(F1,[]);title('傅里叶变换频谱图');
Fs=fftshift(F); %将频谱图中零频率成分移动至频谱图中心
S=log(abs(Fs)+1); %取模并进行缩放
subplot(2,2,3);imshow(S,[]);title('频移后的频谱图');
fr=real(ifft2(ifftshift(Fs))); %频率域反变换到空间域，并取实部
ret=im2uint8(mat2gray(fr)); %更改图像类型
subplot(2,2,4);imshow(ret),title('逆傅里叶变换');
结果显示：
图5.1 结果显示
相关见：https://blog.csdn.net/weixin_43637490/article/details/89196212
说明：对图片进行im2double()函数或者rgb2gray()很重要的一步，因为存入的时候系统默认是三维真彩图像
2.通过对高频进行简单的操纵，实现简单的去噪
II=imread('cameraman.tif');
I=imnoise(II,'gaussian',0,0.01);%加噪
J=fft2(I);%对加噪的照片进行傅里叶变换
J1=fftshift(J);%对加噪的图片进行频移
figure; subplot(1,3,1); imshow(I);
subplot(1,3,2); imshow(log(abs(J)),[8 10]);
subplot(1,3,3); imshow(log(abs(J1)),[8,10]);
J2=J;
J2(10:220,10:220)=0;%对未频移的傅里叶图像进行简单处理
I2=ifft2(J2);%逆变换
I3=abs(I2);
figure; subplot(2,3,1);imshow(I2,[]);
subplot(2,3,2);imshow(I3,[0,255]);
%采用MATLAB中的函数filter2对受噪声干扰的图像进行均值滤波
K1= medfilt2(I);%采用二维中值滤波函数medfilt2对受椒盐噪声干扰的图像滤波
%采用MATLAB中的函数filter2对受噪声干扰的图像进行中值滤波
K2=filter2(fspecial('average',3),I)/255; %模板尺寸为3
K3=filter2(fspecial('average',5),I)/255;% 模板尺寸为5
K4= filter2(fspecial('average',9),I)/255; %模板尺寸为9
subplot(2,3,3);imshow(K1);
title('中值滤波');
subplot(2,3,4); imshow(K2);
title('均值滤波，尺寸3');
subplot(2,3,5);imshow(K3);
title('均值滤波，尺寸5');
subplot(2,3,6);imshow(K4);
title('均值滤波，尺寸9');
图5.2 均值滤波结果
图5.3 结果显示
六、图像增强之空间滤波器基础
（1)简单均值滤波器
clc;
clear;
close all;
img=imread('pears.png');
f=rgb2gray(img);
f=im2double(f);
for i=2:485
for j=2:731
f(i,j)=(f(i-1,j-1)+f(i-1,j)+f(i-1,j+1)+f(i,j-1)+f(i,j)+f(i,j+1)+f(i+1,j-1)+f(i+1,j)+f(i+1,j+1))/9;
end;
end;
subplot(1,2,1);imshow(f);
subplot(1,2,2);imshow(img);
图6.1 均值滤波器
(2)中值滤波器
t=imread('a1.jpg'); t=rgb2gray(t); t1=imnoise(t,'salt & pepper',0.3); subplot(1,2,1),imshow(t1),title('加入椒盐噪声后') t2=medfilt2(t1,[3 3]); subplot(1,2,2),imshow(t2),title('中值滤波后')
图6.2 中值滤波器
七、空间域滤波器
一.知识背景
1.空间平滑滤波器
（1）平滑均值滤波
（2）平滑中值滤波
主要功能：去除椒盐噪声
2.空间锐化滤波器
（1）拉普拉斯算子
（2）反锐化掩蔽与高频提升滤波
（3）Robert交叉梯度算子
（4）Sobel梯度算子
二、样例
clc;
clear;
I=imread('cameraman.tif');
subplot(2,3,1);imshow(I);title('原图');
I1=imnoise(I,'salt & pepper');
subplot(2,3,2);imshow(I1);title('加椒盐噪声');
[w h]=size(I1);
I=double(I);
I1=double(I1);
I_out=I1;
I_out2=I1;
I_out3=I1;
I_out4=I1;
for i = 2:w-1
for j = 2:h-1
%均值滤波
I_out(i,j)=[I1(i-1,j-1)+I1(i-1,j)+I1(i-1,j+1)+...
+I1(i,j-1)+I1(i,j)+I1(i,j+1)+...
+I1(i+1,j-1)+I1(i+1,j)+I1(i+1,j+1)]/9;
%中值滤波
x=[I1(i-1,j-1) I1(i-1,j) I1(i-1,j+1) I1(i,j-1) I1(i,j) I1(i,j+1) I1(i+1,j-1) I1(i+1,j) I1(i+1,j+1)];
I_out2(i,j)=median(x);
%拉普拉斯算子
I_out3(i,j)=[I(i-1,j-1)+I(i-1,j)+I(i-1,j+1)+...
+I(i,j-1)-I(i,j)*8+I(i,j+1)+...
+I(i+1,j-1)+I(i+1,j)+I(i+1,j+1)];
%Robert梯度算子
I_out4(i,j)=max(abs(I(i+1,j+1)-I(j,j)),abs(I(j+1,j)-I(j,j+1)));
end
end
%要加范围[0 255]
subplot(2,3,3);imshow(I_out,[0 255]);title('均值滤波');
subplot(2,3,4);imshow(I_out2,[0 255]);title('中值滤波');
subplot(2,3,5);imshow(I_out3,[0 255]);title('拉普拉斯算子锐化');
subplot(2,3,6);imshow(I_out4,[0 255]);title('Robert交叉梯度算子');
图7.1 结果显示
八、图像复原
逆滤波：就是我们知道了原图受到了什么算法(叫做退化函数 H )带来的噪声，那么我们从噪声图反着用这个算法就可以推算出原图。不过即使知道了退化函数，也还是不能完全恢复到原图。
I = imread('cameraman.tif');
I = im2double(I);
[m, n] = size(I);
M = 2 * m;
N = 2 * n;
u = -m / 2:(m / 2 - 1);
v = -n / 2:(n / 2 - 1);
[U, V] = meshgrid(u, v); %基于向量 u 和 v 中包含的坐标返回二维网格坐标
D = sqrt(U.^2 + V.^2);
D0 = 130; %滤波器截止频率为20
H = exp(-(D.^2)./(2 * (D0^2)));
N = 0.01 * ones(size(I, 1), size(I, 2));
N = imnoise(N, 'gaussian', 0, 0.001);
J = fftfilter(I, H) + N;
HC = zeros(m, n);
M1 = H > 0.1;
HC(M1) = 1 ./ H(M1);
K = fftfilter(J, HC);
HC = zeros(m, n);
M2 = H > 0.01;
HC(M2) = 1 ./ H(M2);
L = fftfilter(J, HC);
figure;
subplot(221);imshow(J);
subplot(222);imshow(J, []);
subplot(223);imshow(K, []);
subplot(224);imshow(L, []);
filter函数：
function Z = fftfilter(X, H)
F = fft2(X, size(H, 1), size(H, 2));
Z = H .* F;
Z = ifftshift(Z);
Z = abs(ifft2(Z));
Z = Z(1:size(X, 1), 1:size(X, 2));
end
图8.1 结果显示
（2）均值滤波
I = imread('火影9.jpg'); I = im2double(I); I = imnoise(I, 'gaussian', 0.05); %加入高斯噪声 %算术均值滤波 PSF = fspecial('average', 3); %创建预定义的二维滤波器 J = imfilter(I, PSF); %几何均值滤波 K = exp(imfilter(log(I), PSF)); figure; subplot(131);imshow(I);title('加入高斯噪声'); subplot(132);imshow(J);title('算术均值滤波'); subplot(133);imshow(K);title('几何均值滤波');
图8.2 结果显示
（3）统计排序滤波
I = imread('火影9.jpg'); I = im2double(I); I = imnoise(I, 'salt & pepper', 0.05); %加入椒盐噪声 J = medfilt2(I, [3, 3]); %中值滤波 figure; subplot(121);imshow(I);title('加入椒盐噪声'); subplot(122);imshow(J);title('中值滤波');
图8.3 统计排序滤波
I = imread('火影9.jpg'); I = rgb2gray(I); I = im2double(I); I = imnoise(I, 'salt & pepper', 0.1); %加入椒盐噪声 domain = [1 1 1 1; 1 1 1 1; 1 1 1 1; 1 1 1 1]; J1 = ordfilt2(I, 1, domain); J2 = ordfilt2(I, 8, domain); J3 = ordfilt2(I, 16, domain); figure; subplot(221);imshow(I);title('加入椒盐噪声'); subplot(222);imshow(J1);title('最小值滤波'); subplot(223);imshow(J2);title('中值滤波'); subplot(224);imshow(J3);title('最大值滤波');
图8.4 结果显示
九、图像分割
图像分割是根据图像的直方图和结构特性或者一些具体的应用需求将图像划分成两个或多个互不相交的子区域的过程，这些子区域是在特定意义下的具有相同属性的像素的连通集合。例如，一幅图像中不同目标物体所占的图像区域、背景所占的背景区域等都属于这样的连通集合概念。因为与上边缘检测等相似，见博客
见：(25条消息) 基于matlab的数字图像处理之彩色图像处理_qq_45840242的博客-CSDN博客
