大家好!欢迎来到本次的总结性的一篇文章,因为咸鱼哥这几个月是真的有点小忙(参加了点小比赛,准备考试等等)所以,在数字图像学习后,我来写一个总结性的文章,同时帮助大家学习,相关细节的东西我会后续写对应内容的详解,写文章的时间有点靠后,所以更新有点慢,多说无益,所以那我们开始吧!
注:咸鱼哥这里的学习是根据(数字图像处理 冈萨雷斯 第四版)来进行的,其实学习中用那本书都可以,笔者也看过其他的基本常见教材,内容基本一致。
一、生成灰度图像
% 读取图像
A = imread('1.1.png');
B = imread('1.1.png');
% 显示原始图像
%figure,imshow(A),title('Original Image A');
% 提取图像的红色、绿色和蓝色通道
Fxy_r = A(:, :, 1); % 红色通道
Fxy_g = A(:, :, 2); % 绿色通道
Fxy_b = A(:, :, 3); % 蓝色通道(注意这里修正了原代码中的语法错误)
% 显示红色、绿色和蓝色通道
% figure;
% subplot(1, 3, 1), imshow(Fxy_r), title('Red Channel');
% subplot(1, 3, 2), imshow(Fxy_g), title('Green Channel');
% subplot(1, 3, 3), imshow(Fxy_b), title('Blue Channel');
% 将图像转换为灰度图
A_gray = rgb2gray(A);
B_gray = rgb2gray(B);
% 显示灰度图像
%figure,imshow(A_gray),title('Grayscale Image');
这里基本上在数字图像处理中都会用到转到灰度这一块,所以要记得灰度在数字图像处理这一节中的重要性(这里难度不是很大的)。
二、建立图像的坐标系
[Rw, Cl] = size(A_gray); % 获取灰度图像的尺寸
[X, Y] = meshgrid(1:Cl, 1:Rw);% 生成网格坐标矩阵
Z = double(A_gray);
% figure, surf(X, Y, Z), shading interp, title('3D');
这里的理解可以简单的论述为:一张图片是由很多的像素点构成的,那么,我们要处理图片局部或者分割图片的操作时,就要先知道区域的位置,而这里的操作也是一样。
三、Sobel边缘检测
预操作(注意其中的double类型):
% 6) Sobel边缘检测
imgGrayDouble = double(rgb2gray(A)); % 将RGB图像转换为灰度图像,并转换为double类型
注意:double类型的灰度,在数值上很精确,在相关处理上效果是较好,在后续的增强明暗对比度的二值化时同样适用,但是用double类型一定要注意数值维度相符合,否则程序会报错
1)使用自定义的算子蒙板,进行垂直边缘的检测
% 垂直Sobel边缘检测
%sobelVertical = imfilter(imgGrayDouble, [-1 0 1; -2 0 2; -1 0 1], 'replicate'); % 应用垂直Sobel滤波器
2)使用自定义的算子蒙板,进行水平边缘的检测
% 水平Sobel边缘检测
%sobelHorizontal = imfilter(imgGrayDouble, [-1 -2 -1; 0 0 0; 1 2 1], 'replicate'); % 应用水平Sobel滤波器
这里推荐下面的方式来进行边缘检测(用matlab内置函数进行):
%sobel边缘检测
edges1 = edge(A_gray,'Sobel');
%canny边缘检测
edges2 = edge(A_gray,'Canny');
但是注意一点的是,在函数中的参数一为灰度图像。
这里是显示图像的操作,uu们可以跳过相关显示的代码,因为在这里只需要知道imshow、figure、title函数的作用和效果就行,用法简单,所以简单描述一下。
% 显示原始灰度图和边缘检测结果
% figure,imshow(uint8(abs(sobelVertical)), []), title('垂直Sobel边缘检测'); % 显示垂直边缘检测结果
% figure,imshow(uint8(abs(sobelHorizontal)), []), title('水平Sobel边缘检测'); % 显示水平边缘检测结果
% figure,imshow(imgGrayDouble, []), title('原始灰度图像'); % 可选:显示原始灰度图像
% figure,imshow(edges1),title('sobel');
% figure,imshow(edges2),title('canny');
四、计算灰度图像
%计算图像直方图
[m,n] = imhist(A);
%计算CDF值
cdf_A = cumsum(m)/numel(A);
%计算CDF增强图
cdf_ip = (1:length(I)/length(m));
% figure,imshow(I),title('gray _ histeq');
%(error)figure,imshow(I1),title('hist');
其CDF是分段 函数,一般为分段函数比如举例中的掷硬币随机变量,求它的CDF:
相关理论概念我会在其他的章节内容中讲述,在这里只需要知道CDF是什么东西即可。
实验四,直方图部分代码分析
A_color = imread('1.jfif'); % 读取彩色图像
A_eq_gray = my_color_to_gray_histeq(A_color); % 转换为灰度图像并进行直方图均衡化
% 显示结果
figure;
subplot(1, 2, 1), imshow(rgb2gray(A_color)), title('转换后的灰度原图'); % 显示转换后的灰度图像作为对比
subplot(1, 2, 2), imshow(A_eq_gray), title('均衡化后的灰度图像');
function A_eq_gray = my_color_to_gray_histeq(A_color)
% 将彩色图像转换为灰度图像
A_gray = rgb2gray(A_color);
% 获取灰度图像的尺寸
[rows, cols] = size(A_gray);
% 计算灰度图像的直方图
[counts, binLocations] = imhist(A_gray);
% 计算累积分布函数(CDF)
cdf = cumsum(counts) / (rows * cols); % 归一化到[0, 1]范围
% 创建一个查找表(LUT),将原始灰度级映射到新的灰度级
% 注意:MATLAB的索引从1开始,但灰度级从0开始,需要相应调整
lut = zeros(1, 256); % 初始化查找表
for i = 0:255
% 由于CDF是单调递增的,找到第一个大于或等于当前灰度级/255的值所对应的索引
% 使用find的第一个输出,并处理边界情况(即当所有值都小于当前灰度级/255时)
idx = find(cdf >= i/255, 1, 'first');
if isempty(idx)
% 如果找不到,则映射到最大值255(索引256,但减1得到灰度级255)
lut(i+1) = 255;
else
% 否则,映射到对应的灰度级
lut(i+1) = idx - 1; % -1是因为lut的索引从1开始,而灰度级从0开始
end
end
% 应用查找表来均衡化灰度图像
% 注意:这里使用uint8转换来确保输出是8位灰度图像
A_eq_gray = uint8(lut(double(A_gray) + 1)); % +1是因为A_gray的灰度级是从0开始的
end
% 使用示例
五、高斯平滑处理
sigma = 0.1
img_guassion = imgaussfilt(A,sigma);
% figure,imshow(img_guassion),title('guassion_smooth');
六、改变图像大小
img_resize = imresize(A,[64 64]);
% figure,imshow(img_resize);
七、为图像添加不同噪声
D=0.05;
img_noise = imnoise(A_gray,'gaussian',0,0.1);
img_noise1 = imnoise(A_gray,'salt & pepper',D);
img_noise2 = imnoise(A_gray,'poisson' );
% figure,subplot(3,2,1),imshow(img_noise);
% subplot(3,2,2),imshow(img_noise1);
% subplot(3,2,3),imshow(img_noise2);
注:这里是单独的调用函数,实际的学习中个人建议还是学习一下底层的关系。
八、进行图像的滤波处理
1)均值滤波
I_avr = filter2(fspecial('average', [3 3]), double(img_noise));
2)加权均值滤波
I_add = filter2(fspecial('gaussian', [3 3], 0.5), double(img_noise));
h_1=[0.1,0.2,0.1;0.2 0.4 0.2;0.1 0.2 0.1];%创建加权滤波模板
I_add2 = imfilter(A,h_1,'same');
% figure;
% imshow(I_add2);
3)中值滤波
I_mid = medfilt2(double(img_noise), [3 3]);
4)最大值滤波
I_max = ordfilt2(double(img_noise), 9, ones(3, 3));
5)最小值滤波
I_min = ordfilt2(double(img_noise), 1, ones(3, 3));
这里是显示出滤波之后的效果图的,和第三部分生成窗口作用一样
% figure;
% subplot(3,2,1),imshow(I_avr);
% subplot(3,2,2),imshow(I_add);
% subplot(3,2,3),imshow(I_max);
% subplot(3,2,4),imshow(I_min);
九、确保A和B具有相同的尺寸
if size(A) ~= size(B)
% 选择一个图像作为参考尺寸,这里以B为例
A_resized = imresize(A, size(B));
else
A_resized = A;
end
%裁剪图像
A_cropped = A(1:size(B,1), 1:size(B,2));
十、图像二值化
brinary_img = imbinarize(A_resized);
十一、相关基本操作
brinary_img = imbinarize(A_cropped);%裁剪完成的二值化操作
brinary_img2= imbinarize(B);
%自适应阈值
brinary_img_1 = imbinarize(A, 'adaptive');
brinary_img_2 = imbinarize(B, 'adaptive');
%求补集
cpl_img = imcomplement(brinary_img);
%求并集
union_img = brinary_img | brinary_img2;
union_img2 = brinary_img_1| brinary_img_2;%自适应的阈值操作
%
% figure,imshow(cpl_img);
%(error)subplot(1,3,2),imshow(union_img);
%图像反转、对数变换、幂次变换和二值化
I_1=255-A_gray;
a=1;
I_2=a*log(1+double(A_gray));
b=0.5;
I_3=a*(double(A_gray) .^b);%b<1,图像变亮,b>1,图像变暗
I_3=unit8(I_3);
% 二值化
k = 2; % 阈值倍数
I_4 = double(A_gray > mean(A_gray) * k);
% figure,subplot(2,2,1),imshow(I_1);
% subplot(2,2,2),imshow(I_2);
% subplot(2,2,3),imshow(I_3);
% subplot(2,2,4),imshow(I_4);
十二、维纳滤波(Wiener Filter)
LEN = 25;
THETA = 20;
PSF = fspecial('motion', LEN, THETA);
J = imfilter(A_gray, PSF, 'conv', 'circular');
NSR = 0;
K = deconvwnr(J, PSF, NSR);
% figure;
% subplot(131); imshow(I);
% subplot(132); imshow(J);
% subplot(133); imshow(K);
十三、锐化
1)(函数)使用 imsharpen 函数进行锐化
sharpened_img = imsharpen(B);
% 显示原始图像和锐化后的图像
% figure;
% subplot(1,2,1), imshow(B), title('原始图像');
% subplot(1,2,2), imshow(sharpened_img), title('锐化后的图像');
2)锐化(laplacian)拉普拉斯滤波器
% 定义拉普拉斯滤波器
laplacian_filter = [0 -1 0; -1 4 -1; 0 -1 0];
% 应用拉普拉斯滤波器
laplacian_img = filter2(laplacian_filter, double(B_gray), 'same');
% 将拉普拉斯结果加回到原始图像以锐化
sharpened_img = B_gray + uint8(laplacian_img);
% 显示原始图像和锐化后的图像
% figure;
% subplot(1,2,1), imshow(B), title('原始图像');
% subplot(1,2,2), imshow(sharpened_img), title('锐化后的图像');
十四、使用高斯混合模型(GMM)进行图像分割
B = im2double(B);% 将图像转换为双精度型,以便进行计算
% 使用高斯混合模型进行图像分割
num_clusters = 2;
gmm_model = fitgmdist(double(B_gray(:)), num_clusters);
% 预测每个像素的类别
pixel_labels = cluster(gmm_model, double(B_gray(:)));
pixel_labels = reshape(pixel_labels,size(B_gray));
% 提取前景区域和背景区域
foreground =pixel_labels == 1;
background =pixel_labels == 2;
% 计算前景区域的距离变换
distance_transform = bwdist(~foreground);
% 计算前景边缘
edges = edge(distance_transform, 'Canny');
%显示结果
figure
% subplot(235)
imshow(distance_transform);
% title('距离变换结果');
colormap(jet);
colorbar;
subplot(2,3,[1:2])
imshow(B);
title('原始图像')
subplot(233)
imshow(foreground)
title('前景图')
subplot(234)
imshow(background)
title('背景图')
subplot(236)
imshow(edges); % 显示边缘图
title('前景边缘显示结果');
十五、滤波器函数的设定
预处理:
[ra,ca] = size(A_gray);
P = 2 * ra,Q =2 * ca;%设置填充后的图像大小
A_2 = zeros(P,Q);%用零填充图像;
A_2(1:ra,1:ca)=A_gray;
[Ra,Ca]=size(A_2);
[y,x]=meshgrid(0:Ca-1,0:Ra-1);%建立空域二维坐标系
F_B=fft2(A_2.*(-1).^(x+y));%傅里叶变化乘以(-1)^(x+y)
1)理想低通滤波器
D_1=12;
H_1=zeros(P,Q);
for n=1:P
for m = 1:Q
if((n-P/2-1)^2+(m-Q/2-1)^2 <=(D_1)^2);
H_1(n,m)=1;%赋值给1
else
H_1(n,m)=0;%赋值给0
end
end
end
figure,mesh(H_1),title('滤波器透视图,理想低通D0=10');
G=F_B.*H_1;
fa=ifft2(G);
fa=real(fa).*(-1).^(x+y);
fA=fa(1:ra,1:ca);%去除补零的部分;
figure,subplot(1,3,2),imshow(fA,[]);
2)理想高通滤波器
figure
subplot(4,2,1),imshow(A),title('原图');
H_2= 1-H_1;
subplot(4,2,3),mesh(H_2),title('滤波器透视图,理想高通D0=10');
G1=F_B.*H_2;
fa1=ifft2(G1);
fa1=real(fa1).*(-1).^(x+y);
fA1=fa1(1:ra,1:ca);
subplot(1,3,3),imshow(fA1,[]);
3)高斯低通滤波器
figure
subplot(4,2,1),imshow(A),title('原图');
U = 0:P-1;
V = 0:Q-1;%赋值傅里叶频谱频率值
[V,U] = meshgrid(0:Q-1,0:P-1);%生成频域二维坐标
D_mid = sqrt((U-P/2).^2+(V-Q/2).^2);%赋值频域点到中心距离
H_uv = exp(-D_mid.^2/2/D_0/D_0);%生成滤波器
subplot(4,2,3),mesh(H_uv),title('滤波器透视图,高斯低D0=10)');
G_pxy=F_B.*H_uv;%用滤波器滤波
g_pxy=ifft2(G_pxy);%做逆傅里叶变换
g_pxy=real(g_pxy).*(-1).^(x+y);%取实部后乘以(-1)的x+y次方
g_xy=g_pxy(1:ra,1:ca);%除去补零部分
subplot(4,2,4),imshow(g_xy,[]),title('滤波结果图,高斯低D0=10)');
4)高斯高通滤波器
figure
subplot(4,2,1),imshow(A),title('原图');
H_uv_3=1-H_uv;%用1减去低通滤波器得到高通滤波器
subplot(4,2,3),mesh(H_uv_3),title('滤波器透视图,高斯高0=10)');
G_pxy_3=F_B.*H_uv_3;%用滤波器滤波
g_pxy_3=ifft2(G_pxy_3);%做逆傅里叶变换
g_pxy_3=real(g_pxy_3).*(-1).^(x+y);%取实部后乘以(-1)的x+y次方
g_xy_3=g_pxy_3(1:ra,1:ca);%除去补零部分
subplot(4,2,4),imshow(g_xy_3,[]),title('滤波结果图,高斯高D0=10)');
5)巴特沃斯低通滤波器
U=0:P-1;V=0:Q-1;%赋值傅里叶频谱频率值
[V,U]=meshgrid(0:Q-1,0:P-1);%生成频域二维坐标
D_mid=sqrt((U-P/2).^2+(V-Q/2).^2);%赋值频域点到中心距离
figure;
subplot(4,2,1),imshow(A),title('原图');
D_3=20;%赋值截止频率
n=1;%赋值巴特沃斯滤波器阶数
H_UV=1./(1+(D_mid/D_3).^(2*n));%生成滤波器
subplot(4,2,3),mesh(H_UV),title('滤波器透视图,巴低通D0=20,n=1)');
G_PXY=F_B.*H_UV;%用滤波器滤波
f_PXY=ifft2(G_PXY);%做逆傅里叶变换
f_PXY=real(f_PXY).*(-1).^(x+y);%取实部后乘以(-1)的x+y次方
f_XY=f_PXY(1:ra,1:ca);%除去补零部分
subplot(4,2,4),imshow(f_XY,[]),title('滤波结果图,巴低通D0=20,n=1)');
6)巴特沃斯高通滤波器
figure
subplot(4,2,1),imshow(A),title('原图');
n=1;%赋值巴特沃斯滤波器阶数
H_UV=1./(1+(D_3./D_mid).^(2*n));%生成滤波器
subplot(4,2,3),mesh(H_UV),title('滤波器透视图,巴高通D0=20,n=1)');
G_PXY=F_B.*H_UV;%用滤波器滤波
f_PXY=ifft2(G_PXY);%做逆傅里叶变换
f_PXY=real(f_PXY).*(-1).^(x+y);%取实部后乘以(-1)的x+y次方
f_XY=f_PXY(1:ra,1:ca);%除去补零部分
subplot(4,2,4),imshow(f_XY,[]),title('滤波结果图,巴高通D0=20,n=1)');
此外还有陷波滤波器等
其中陷波滤波器、巴特沃斯滤波器、高斯滤波器的差异都很小,只是在于滤波器使用的公式不一样,所以我们在学习的时候可以举一反三的学习,这也是学习的一种良好的方式。
十六、形态学相关操作
1)腐蚀图像
se1=[1;1;1]; %线型结构元素
FuShi=imerode(BianYuan,se1);
2)图像聚类、填充图像
se2=strel('rectangle',[25,25]); %矩形结构元素
TianChong=imclose(FuShi,se2);%图像聚类、填充图像
3)移除指定面积图像
YuanShiLvBo=bwareaopen(TianChong,2000);%从对象中移除面积小于2000的小对象
十七、RGB to HSI
rgbImage = imread('image.jpg');
rgbImage = im2double(rgbImage); % Normalize to [0,1]
%Separate color channels
R=rgbImage(:,:,1);
G=rgbImage(:,:,2);
B=rgbImage(:,:,3);
%Calculate intensity
I = (R +G+B)/3;
%Calculate saturation
minRGB =min(min(R, G), B);
S =1-(minRGB ./ I);
S(I == 0) = 0; % Handle division by zero
%Calculate hue
theta = acos(0.5 * ((R- G) + (R- B)) ./ sqrt((R- G).^2 + (R- B).*(G- B)));
H=theta;
H(B >G) =2*pi- H;
%Normalize H to [0,1]
H=H/(2*pi);
%Combine channels into HSI image
hsiImage = cat(3, H, S, I);
%Display the result
imshow(hsiImage);
title('HSI Image')
十八、RGB to HSV
%Read the image
rgbImage = imread('image.jpg');
%Convert to HSV
hsvImage = rgb2hsv(rgbImage);
%Display the result
imshow(hsvImage);
title('HSV Image');
%Assuming `hsiImage` is the HSI image you already have
%Extract individual channels
H=hsiImage(:,:,1); % Hue
S =hsiImage(:,:,2); % Saturation
I = hsiImage(:,:,3); % Intensity
%Display the individual components
figure;
subplot(1, 3, 1);
imshow(H);
title('Hue');
subplot(1, 3, 2);
imshow(S);
title('Saturation');
subplot(1, 3, 3);
imshow(I);
title('Intensity')
十九、彩图的单色图片展示
% 确保图像是RGB格式
if size(img, 3) == 3
% 创建单色图像
redImg = img;
greenImg = img;
blueImg = img;
% 将其他颜色通道设置为0
redImg(:,:,2:3) = 0; % 红色通道图像,G和B通道为0
greenImg(:,:,3:3) = 0; % 绿色通道图像,R和B通道为0
blueImg(:,:,1:2) = 0; % 蓝色通道图像,R和G通道为0
% 显示单色图像
figure;
subplot(1,3,1);
imshow(redImg);
title('Red Channel Image');
subplot(1,3,2);
imshow(greenImg);
title('Green Channel Image');
subplot(1,3,3);
imshow(blueImg);
title('Blue Channel Image');
else
disp('The image is not in RGB format.');
end
二十、绘制区域图像
p=ones(256,256);
p(96:105,:)=0;
p(:,96:105)=0;
figure;
subplot(1,3,1),imshow(p);
p2=ones(256,256);
p2(96:105,:)=0;
p2(:,96:105)=0;
subplot(1,3,2),imshow(p2);
j=union(p,p2);
subplot(1,3,3),imshow(j);
好了,到这里,我们的matlab版对于数字图像处理课程学习的总结就结束了,我们可以发现基本上的所有操作都需要灰度的处理(除需要彩色图像时候,我们采用分图层或者显示同色图像),其他的过程中会有底层编写或者函数的调用,这里咸鱼哥建议大家有时间的话也可以使用matlab工具箱(2022版以上的都有)来进行其他好玩的操作,比如:无人机算法路径规划实现,算法验证等,其实咸鱼哥个人感觉还是很有意思的。
关于理论的分析,咸鱼哥会在后续一个个总结写出,相当于一种总结和学习的经历,望诸君共勉,我们共同学习!
uu们!我们下次见!
另记:以上代码均测试过,可以在matlab中运行,所以初学者也可以放心学习,如有不会之处可以在评论区发表,我会认真解答的。(此外,以上代码结果均太多,而且大部分相同,所以希望大家抽时间在matlab中运行,这样起到学习和理解,也方便我们看清数字图像处理技术有什么作用)。