图像处理.docx
《图像处理.docx》由会员分享，可在线阅读，更多相关《图像处理.docx（39页珍藏版）》请在冰豆网上搜索。
图像处理
作业一：
编写程序，研究中值滤波算法、均值滤波算法、选择式掩模平滑算法对不同噪声图像的滤波效果。
我编写了中值滤波算法、均值滤波算法、超限邻域平均算法以及选择式掩模平滑算法四个算法。
算法处理效果图如下：
结论：
均值滤波可以起到很好地将噪声点消除，但是均值滤波处理是以图像模糊为代价来减小噪声的。
邻域越大，处理图像越模糊。
为了减小模糊效应，需要力求改进的途径。
超限邻域平均法很好地减少了边缘模糊现象。
对抑制椒盐噪声比较有效。
选择式掩模平滑法是一种自适应的局部平滑滤波。
是以尽量不模糊边缘轮廓为目的的。
不过实验效果并不是很好，我觉得可能与图片本身有关，这个算法也需要进一步改进。
中值滤波是一种非线性的信号处理方法。
它可以把干扰去除，而不刻意让图像模糊。
作业二：
编写程序，提取目标边缘，比较各种算子性能。
我主要编写了canny算子，另外与sobel算子做了对比。
（prewitt算子和高斯-拉普拉斯算子与sobel算子的实现方式大同小异。
）
Canny算子是一种考虑比较全面，效果良好的边缘检测算子。
算法处理效果图如下：
结论：
Robert算子是一种最简单的算子，利用局部差分寻找边缘。
它对陡峭的低噪声图像效果很好，但是对边缘的定位不够准确，边缘较粗。
Sobel算子是滤波算子的形式，用两个模板对图像做卷积，求出x、y方向的导数进一步求出合成方向的导数和方向。
对灰度渐变和噪声较多的图像处理效果较好。
Prewitt算子与sobel算子类似，只是卷积模板不同，理解的邻域权值不同。
对灰度渐变和噪声较多的图像处理效果较好。
Log算子是一种二阶微分算子，经常出现双像素边界，对噪声敏感。
Canny算子不易受噪声干扰，能够检测到真正的弱边缘。
使用时需要调节好高低阈值。
从图片中看出，canny算子检测效果更好。
检测效果与检测的图片有关，也与调节参数有关。
作业三：
编写程序，以一副中国地图和各省GDP、人口作为输入，计算中国的地理中心、经济中心、人口中心。
我编写了地理中心以及地势高度中心的程序。
虽然后两个原理是一样的，但调起来还是需要时间的。
程序效果图：
图3·1地理中心
地理中心大约在甘肃省中下部。
图3·2地势高度中心
图3·3地势高度中心
地势高度中心大约在青海省东部。
图3·3是进行开运算后的处理结果，位置往上移了一点。
作业四：
图像配准、运动跟踪
我编写了矩匹配算法，可以实现模板和图片的匹配。
效果如下：
第一幅图为模板，在第二幅图上进行了匹配。
运动匹配学习了ppt上的算法。
程序代码如下：
%主程序：
读取图像，调用不同函数（处理算法）。
clear
clc
f=imread（'F:
\matlab\shipintuxiangchuli\filter\test.jpg'）;
f=f（:
:
2）;%取灰度值，第三个值定值
figure
imshow（f）;
title（'原始图像'）;
%%subplot（3,1,1）,imshow（f）;
f1=imnoise（f,'gaussian',0.1）;
figure
imshow（f1）;
f2=imnoise（f,'salt&pepper',0.1）;
%figure
%imshow（f2）;
%%subplot（3,1,2）,imshow（f2）;
%title（'加入椒盐噪声'）;
f3=imnoise（f1,'salt&pepper',0.1）;
%figure
%imshow（f3）;
%%subplot（2,3,3）,imshow（f2）;
%title（'加入白噪声和椒盐噪声'）;
%%中值滤波情况
ff1=mid_filter（f1,3）;
f_mid_filter_1=uint8（ff1）;
figure
imshow（f_mid_filter_1）;
%subplot（3,1,3）,imshow（f_mid_filter_1）;
title（'3*3中值滤波'）;
%ff2=mid_filter（f2,5）;
%f_mid_filter_2=uint8（ff2）;
%figure
%imshow（f_mid_filter_2）;
%%subplot（2,3,6）,imshow（f_mid_filter_2）;
%title（'5*5中值滤波'）;
%%均值滤波情况
%ff1=ave_filter（f3,3）;
%f_ave_filter_1=uint8（ff1）;
%figure
%imshow（f_ave_filter_1）;
%ff2=ave_filter（f3,5）;
%f_ave_filter_2=uint8（ff2）;
%figure
%imshow（f_ave_filter_2）;
%超限邻域平均法
%ff1=superave_filter（f2,3）;
%f_superave_filter_1=uint8（ff1）;
%figure
%imshow（f_superave_filter_1）;
%title（'3*3超限邻域平均（加入椒盐噪声,阈值50）'）;
%
%
%%选择式掩膜平滑
%ff1=Sel_masking_smoothing（f1）;
%f_superave_filter_1=uint8（ff1）;
%figure
%imshow（f_superave_filter_1）;
%title（'选择式掩膜平滑（白噪声）'）
%中值滤波算法函数
functionx2=mid_filter（f,n）
x1=double（f）;
x2=x1;
[height,width]=size（x1）;
fori=（（n+1）/2）:
height-（（n+1）/2）
forj=（（n+1）/2）:
width-（（n+1）/2）
%取窗口元素
temp1=x1（i-（（n-1）/2）:
i+（（n-1）/2）,j-（（n-1）/2）:
j+（（n-1）/2））;
temp2=temp1（1,:
）;
foru=2:
n
temp2=[temp2,temp1（u,:
）];%将窗口元素转化成行矩阵
end
m=median（temp2）;%取中值
x2（i,j）=m;%取窗口中值作为像素点
end
end
end
%均值滤波算法函数
functionx2=ave_filter（f,n）
x1=double（f）;
x2=x1;
[height,width]=size（x1）;
fori=（（n+1）/2）:
height-（（n+1）/2）
forj=（（n+1）/2）:
width-（（n+1）/2）
%取窗口元素
temp1=x1（i-（（n-1）/2）:
i+（（n-1）/2）,j-（（n-1）/2）:
j+（（n-1）/2））;
s=sum（sum（temp1））;%求和
x2（i,j）=s/n^2;%取窗口均值作为像素点
end
end
end
%超限邻域平均算法函数
functionx2=superave_filter（f,n）
x1=double（f）;
x2=x1;
a=50;
[height,width]=size（x1）;
fori=（（n+1）/2）:
height-（（n+1）/2）
forj=（（n+1）/2）:
width-（（n+1）/2）
%取窗口元素
temp1=x1（i-（（n-1）/2）:
i+（（n-1）/2）,j-（（n-1）/2）:
j+（（n-1）/2））;
%求周围n*n-1个点的和
s=（sum（sum（temp1））-x1（i,j））/（n^2-1）;
%判断像素点与周围像素的均值之差是否超过阈值
ifabs（x1（i,j）-s）>a
x2（i,j）=s;
%else
%x2（i,j）=x1（i,j）;
end
end
end
end
%选择式掩模平滑算法
functionx2=Sel_masking_smoothing（f）
x1=double（f）;
x2=x1;
[height,width]=size（x1）;
k（1,1:
9）=0;
fori=3:
height-3
forj=3:
width-3
%%%%周围9邻近
a1=[x1（i-1,j-1:
j+1）x1（i,j-1:
j+1）x1（i+1,j-1:
j+1）];
s
（1）=sum（a1）/9;
forb=1:
9
k
（1）=k
（1）+（a1（b）^2-s
（1）^2）;
end
%%%%左7邻近
a2=[x1（i-2,j-1:
j+1）x1（i-1,j-1:
j+1）x1（i,j）];
s
（2）=sum（a2）/7;
forb=1:
7
k
（2）=k
（2）+（a2（b）^2-s
（2）^2）;
end
%%%%上7邻近
a3=[x1（i-1:
i+1,j-2）'x1（i-1:
i+1,j-1）'x1（i,j）];
s（3）=sum（a3）/7;
forb=1:
7
k（3）=k（3）+（a3（b）^2-s（3）^2）;
end
%%%%右7邻近
a4=[x1（i+1,j-1:
j+1）x1（i+2,j-1:
j+1）x1（i,j）];
s（4）=sum（a4）/7;
forb=1:
7
k（4）=k（4）+（a4（b）^2-s（4）^2）;
end
%%%%下7邻近
a5=[x1（i-1:
i+1,j+1）'x1（i-1:
i+1,j+2）'x1（i,j）];
s（5）=sum（a5）/7;
forb=1:
7
k（5）=k（5）+（a5（b）^2-s（5）^2）;
end
%%%%左上7邻近
a6=[x1（i-2,j-2:
j-1）x1（i-1,j-2:
j）x1（i,j-1:
j）];
s（6）=sum（a6）/7;
forb=1:
7
k（6）=k（6）+（a6（b）^2-s（6）^2）;
end
%%%%右上7邻近
a7=[x1（i,j-1:
j）x1（i+1,j-2:
j）x1（i+2,j-2:
j）];
s（7）=sum（a7）/7;
forb=1:
7
k（7）=k（7）+（a7（b）^2-s（7）^2）;
end
%%%%右下7邻近
a8=[x1（i,j:
j+1）x1（i+1,j:
j+2）x1（i+2,j:
j+2）];
s（8）=sum（a8）/7;
forb=1:
7
k（8）=k（8）+（a8（b）^2-s（8）^2）;
end
%%%%左下邻近
a9=[x1（i-2,j+1:
j+2）x1（i-1,j:
j+2）x1（i,j:
j+1）];
s（9）=sum（a9）/7;
forb=1:
7
k（9）=k（9）+（a9（b）^2-s（9）^2）;
end
[maxr,index]=min（k）;
x2（i,j）=s（index）;
end
end
end
%%Canny算子
clear
clc
closeall
f=imread（'1.jpg'）;
figure,imshow（f）;
title（'原图'）;
f=f（:
:
2）;
figure,imshow（f）;
title（'灰度图'）;
q1=double（f）;
J=q1;
[height,width]=size（q1）;
%%一般情况下,使用高斯平滑滤波器卷积降噪.K为size=5的高斯内核
K=1/139*[24542;
491294;
51215125;
491294;
24542];
%%对图像实施高斯滤波
fori=1:
height
forj=1:
width
sum=0;
form=1:
5
forn=1:
5
%sum=sum+K（m,n）*q1（i,j）;
if（i-3+m）>0&&（i-3+m）<=height&&（j-3+n）>0&&（j-3+n）sum=sum+K（m,n）*q1（i-3+m,j-3+n）;
end
end
end
J（i,j）=sum;
end
end
J1=uint8（J）;
figure,imshow（J1）;
title（'高斯滤波后的结果'）;
%%用一阶偏导的有限差分来计算梯度的幅值和方向；
dx=zeros（height,width）;%x方向梯度
dy=zeros（height,width）;%y方向梯度
d=zeros（height,width）;%近似梯度
%fori=1:
height-1
%forj=1:
width-1
%dx（i,j）=（J（i,j+1）-J（i,j）+J（i+1,j+1）-J（i+1,j））/2;
%dy（i,j）=（J（i,j）-J（i+1,j）+J（i,j+1）-J（i+1,j+1））/2;
%d（i,j）=sqrt（dx（i,j）^2+dy（i,j）^2）;
%end
%end
%%用sobel算子计算梯度
fori=2:
height-1
forj=2:
width-1
dx（i,j）=J（i+1,j-1）+2*J（i+1,j）+J（i+1,j+1）-J（i-1,j-1）-2*J（i-1,j）-J（i-1,j+1）;
dy（i,j）=J（i-1,j+1）+2*J（i,j+1）+J（i+1,j+1）-J（i-1,j-1）-2*J（i,j-1）-J（i+1,j-1）;
d（i,j）=sqrt（dx（i,j）^2+dy（i,j）^2）;
end
end
figure,imshow（d,[]）;
title（'求梯度后的结果'）;
%%非极大值抑制
D=d;%记录进行抑制后的梯度
%%首先图像四周边缘肯定不是边缘点
forj=1:
width
D（1,j）=0;
end
forj=1:
width
D（height,j）=0;
end
fori=2:
width-1
D（i,1）=0;
end
fori=2:
width-1
D（i,width）=0;
end
fori=2:
height-1
forj=2:
width-1
ifd（i,j）==0%%当前像素点梯度为0,则一定不是边缘点
D（i,j）=0;
else
gradX=dx（i,j）;%当前点x方向导数
gradY=dy（i,j）;%当前点y方向导数
gradTemp=d（i,j）;%当前点梯度
%如果y方向梯度大,则赋值y梯度
ifabs（gradY）>abs（gradX）
weight=abs（gradX）/abs（gradY）;%权重
grad2=d（i,j-1）;
grad4=d（i,j+1）;
%如果x、y方向导数符号相同
%像素点位置关系
%g1g2
%C
%g4g3
ifgradX*gradY>0
grad1=d（i-1,j-1）;
grad3=d（i+1,j+1）;
else
%如果x、y方向导数符号反
%像素点位置关系
%g2g1
%C
%g3g4
grad1=d（i+1,j-1）;
grad3=d（i-1,j+1）;
end
%如果X方向幅度值较大
else
weight=abs（gradY）/abs（gradX）;%权重
grad2=d（i-1,j）;
grad4=d（i+1,j）;
%如果x、y方向导数符号相同
%像素点位置关系
%g3
%g4Cg2
%g1
ifgradX*gradY>0
grad1=d（i-1,j+1）;
grad3=d（i+1,j-1）;
else
%如果x、y方向导数符号反
%像素点位置关系
%g1
%g4Cg2
%g3
grad1=d（i-1,j-1）;
grad3=d（i+1,j+1）;
end
end
%%利用grad1~grad4权重关系对梯度进行插值
gradTemp1=weight*grad1+（1-weight）*grad2;
gradTemp2=weight*grad3+（1-weight）*grad4;
%%当前像素的梯度是局部的最大值，可能是边缘点
ifgradTemp>=gradTemp1&&gradTemp>=gradTemp2
D（i,j）=gradTemp;
else
%%不可能是边缘点
D（i,j）=0;
end
end
end
end
figure,imshow（D,[]）;
title（'非极大值抑制后的结果'）;
%%滞后阈值:
高阈值和低阈值
EP_MIN=30;
EP_MAX=60;
EdgeLarge=zeros（height,width）;%边缘像素
EdgeBetween=zeros（height,width）;%记录可能存在的边缘点
%dd=zeros（height,width）;
%fori=2:
height
%forj=2:
width
%ifD（i,j）>=D_max;%如果超过高阈值,像素被保留为边缘像素
%%EdgeLarge（i,j）=D（i,j）;
%dd=D（i,j）;
%elseifD（i,j）>=D_min&&（D（i,j+1）>=D_max||D（i,j-1）>=D_max...
%||D（i+1,j）>=D_max||D（i-1,j）>=D_max）%在两个阈值之间,暂时保留在此处,接下来进行连接计算
%%EdgeBetween（i,j）=D（i,j）;
%dd=D（i,j）;
%end%低于低阈值的,被排除
%end
%end
%end
%figure,imshow（dd,[]）
%title（'滞后阈值的结果'）;
fori=1:
height
forj=1:
width
ifD（i,j）>=EP_MAX%大于最大阈值，为边缘点
EdgeLarge（i,j）=D（i,j）;
elseifD（i,j）>=EP_MIN%在两个阈值之间,进行下一步判断
EdgeBetween（i,j）=D（i,j）;
end
end
end
end
%把EdgeLarge的边缘连成连续的轮廓
MAXSIZE=999999;
Queue=zeros（MAXSIZE,2）;%用数组模拟队列
front=1;%队头
rear=1;%队尾
edge=zeros（height,width）;
fori=1:
height
forj=1:
width
ifEdgeLarge（i,j）>0
%强点入队
Queue（rear,1）=i;
Queue（rear,2）=j;
rear=rear+1;
edge（i,j）=EdgeLarge（i,j）;
EdgeLarge（i,j）=0;%避免重复计算
end
whilefront~=rear%队不空
%队头出队
temp_i=Queue（front,1）;
temp_j=Queue（front,2）;
front=front+1;
%8-连通域寻找可能的边缘点
%左上方
ifEdgeBetween（temp_i-1,temp_j-1）>0%把在强点周围的弱点变为强点
EdgeLarge（temp_i-1,temp_j-1）=D（temp_i-1,temp_j-1）;
EdgeBetween（temp_i-1,temp_j-1）=0;%避免重复计算
%入队
Queue（rear,1）=temp_i-1;
Queue（rear,2）=temp_j-1;
rear=rear+1;
end
%正上方
ifEdgeBetween（temp_i-1,temp_j）>0%把在强点周围的弱点变为强点
EdgeLarge（temp_i-1,temp_j）=D（temp_i-1,temp_j）;
EdgeBetween（temp_i-1,temp_j）=0;
%入队
Queue（rear,1）=temp_i-1;
Queue（rear,2）=temp_j;
rear=rear+1;
end
%右上方
ifEdgeBetween（temp_i-1,temp_j+1）>0%把在强点周围的弱点变为强点
EdgeLarge（temp_i-1,temp_j+1）=D（temp_i-1,temp_j+1）;
EdgeBetween（temp_i-1,temp_j+1）=0;
%入队
Queue（rear,1）=temp_i-1;
Queue（rear,2）=temp_j+1;
rear=rear+1;
end
%正左方
ifEdgeBetween（temp_i,temp_j-1）>0%把在强点周围的弱点变为强点
EdgeLarge（temp_i,temp_j-1）=D（temp_i,temp_j-1）;
EdgeBetween（temp_i,temp_j-1）=0;
%入队
Queue（rear,1）=temp_i;
Queue（rear,2）=temp_j-1;
rear=rear+1;
end
%正右方
ifEdgeBetween（temp_i,temp_j+1）>0%把在强点周围的弱点变为强点
EdgeLarge（temp_i,temp_j+1）=D（temp_i,temp_j+1）;
EdgeBetween（temp_i,temp_j+1）=0;
%入队
Queue（rear,1）=temp_i;
Queue（rear,2）=temp_j+1;
rear=rear+1;
end
%左下方
ifEdgeBetween（temp_i+1,temp_j-1）>0%把在强点周围的弱点变为强点
EdgeLarge（temp_i+1,temp_j-1）=D（temp_i+1,temp_j-1）;
EdgeBetween（temp_i+1,temp_j-1）=0;
%入队
Queue（rear,1）=temp_i+1;
Queue（rear,2）=temp_j-1;
rear=rear+1;
end
%正下方
ifEdgeBet
升级会员 