基于偏微分方程的图像处理.docx

1、基于偏微分方程的图像处理基于偏微分方程的图像处理1 图像的基本操作(1) 把一幅彩色图像分解为R、G、B三副单色图像； clear;image_I=imread(rgbtest2.bmp);subplot(2,2,1);imshow(image_I);matrix_R(:,:,1)=image_I(:,:,1);matrix_R(:,:,2)=0;matrix_R(:,:,3)=0;subplot(2,2,2);imshow(matrix_R);title(R分量);matrix_G(:,:,2)=image_I(:,:,2);matrix_G(:,:,1)=0;matrix_G(:,:,3)

2、=0;subplot(2,2,3);imshow(matrix_G);title(G分量);matrix_B(:,:,3)=image_I(:,:,3);matrix_B(:,:,1)=0;matrix_B(:,:,2)=0;subplot(2,2,4);imshow(matrix_B);title(B分量); (2) 把一幅灰度图像分别沿x轴和y轴做反射，扩展为四倍大小； clear;Image=imread(graytest2.bmp);imshow(Image);m,n=size(Image);image11=Image;for i=1:m for j=1:n image12(i,j)=

3、image11(i,n-j+1); endendfor j=1:n for i=1:m image21(i,j)=image11(m-i+1,j); image22(i,j)=image12(m-i+1,j); endendimage1=image11,image12;image2=image21,image22;image=image1;image2;figureimshow(image)原图像：扩展图像： 2 把一幅灰度图像的像素值都变换为原来的1/2使图像质量变差，然后利用Matlab图像处理工具箱中的直方图均衡化函数histeq对图像进行增强，输出原图像、质量变差图像和增强后的图像以及

4、它们的直方图。 clear;Image=imread(graytest2.bmp);Image1=0.5*double(Image);subplot(3,2,1);imshow(Image);subplot(3,2,2);imhist(Image);subplot(3,2,3);imshow(uint8(Image1);image1=uint8(Image1);imwrite(image1,graytest2_1.bmp)subplot(3,2,4);imhist(image1);Image_hist=histeq(image1);subplot(3,2,5);imshow(Image_his

5、t);imwrite(Image_hist,graytest2_2.bmp)subplot(3,2,6);imhist(Image_hist);3 图像的水平集、水平线和等高线的显示和处理(1) 显示一幅灰度图像在某一水平的水平集、水平线和等高线；function level(image_name,N)%求一图像的水平集、线%输入参数：image_name图像文件名% ：N阀值（0255）if nargin1 disp(请输入文件名!); return;elseif nargin N matrix_temp(row,col)=1; end endendsubplot(2,2,2);imshow

6、(matrix_temp,);title(图像的水平集);imwrite(matrix_temp,level_setzhan.bmp);% 图像矩阵扩展 赋值 便于处理边界matrix_ex=zeros(size_m+2,size_n+2);for row=1:size_m for col=1:size_n matrix_ex(row+1,col+1)=matrix_temp(row,col); endend% 四邻域反填充 得水平线matrix_new=matrix_temp;for row=2:size_m+1 for col=2:size_n+1 if matrix_ex(row+1,c

7、ol)=0 & matrix_ex(row-1,col)=0 & matrix_ex(row,col+1)=0 & matrix_ex(row,col-1)=0 matrix_new(row-1,col-1)=1; end endendsubplot(2,2,3);imshow(matrix_new,);title(图像的水平线);%imwrite(matrix_new,level_line.bmp);% 求图像的等高线contour=zeros(size_m,size_n);for row=1:size_m for col=1:size_n if image_input(row,col)=N

8、 contour(row,col)=1; end contour(row,col)=1-contour(row,col); endendsubplot(2,2,4);imshow(contour,);title(图像的等高线);%imwrite(contour,contour.bmp);function level_line_set(image_name)% 求一图像的阀值N取各值所对应的水平线的叠加 % 输入参数：image_name图像文件名 image_input=imread(image_name); %不变量size_m,size_n=size(image_input); %不变量m

9、atrix_out=zeros(size_m,size_n); %循环中不需重新初始化的变量% 输出矩阵初始化（白）for row=1:size_m for col=1:size_n matrix_out(row,col)=1; endendmatrix_ex=zeros(size_m+2,size_n+2); % 扩展矩阵subplot(1,2,1);imshow(image_input);title(原图像);for N=1:32:255 % N为相应阀值 matrix_temp=zeros(size_m,size_n); % 重新初始化 （黑） for row=1:size_m for

10、col=1:size_n if image_input(row,col) =N matrix_temp(row,col)=1; % 求水平集 end end end % 图像矩阵扩展 赋值 便于处理边界 for row=1:size_m for col=1:size_n matrix_ex(row+1,col+1)=matrix_temp(row,col); end end % 四邻域（黑区域填白）反填充 得水平线 matrix_new=matrix_temp; for row=2:size_m+1 for col=2:size_n+1 if matrix_ex(row+1,col)=0 &

11、matrix_ex(row-1,col)=0 & matrix_ex(row,col+1)=0 & matrix_ex(row,col-1)=0 matrix_new(row-1,col-1)=1; % 填白 end end end matrix_out=min(matrix_out,matrix_new); %各级水平线叠加endmatrix_out(1,1)=0;matrix_out(2,1)=1;subplot(1,2,2);imshow(matrix_out,);title(图像部分水平集的叠加);level(graytest2.bmp,100);(2) 显示图像并显示它的一族水平线（

12、要求阀值为32, 64, 96, 128, 160, 192, 224）。function level_line_set(image_name)% 求一图像的阀值N取各值所对应的水平线的叠加 % 输入参数：image_name图像文件名 image_input=imread(image_name); %不变量size_m,size_n=size(image_input); %不变量matrix_out=zeros(size_m,size_n); %循环中不需重新初始化的变量% 输出矩阵初始化（白）for row=1:size_m for col=1:size_n matrix_out(row,

13、col)=1; endendmatrix_ex=zeros(size_m+2,size_n+2); % 扩展矩阵subplot(1,2,1);imshow(image_input);title(原图像);for N=1:32:255 % N为相应阀值（要求阀值为32, 64, 96, 128, 160, 192, 224）。 matrix_temp=zeros(size_m,size_n); % 重新初始化 （黑） for row=1:size_m for col=1:size_n if image_input(row,col) =N matrix_temp(row,col)=1; % 求水平

14、集 end end end % 图像矩阵扩展 赋值 便于处理边界 for row=1:size_m for col=1:size_n matrix_ex(row+1,col+1)=matrix_temp(row,col); end end % 四邻域（黑区域填白）反填充 得水平线 matrix_new=matrix_temp; for row=2:size_m+1 for col=2:size_n+1 if matrix_ex(row+1,col)=0 & matrix_ex(row-1,col)=0 & matrix_ex(row,col+1)=0 & matrix_ex(row,col-1

15、)=0 matrix_new(row-1,col-1)=1; % 填白 end end end matrix_out=min(matrix_out,matrix_new); %各级水平线叠加endmatrix_out(1,1)=0;matrix_out(2,1)=1;subplot(1,2,2);imshow(matrix_out,);title(图像部分水平集的叠加);4 显示一幅灰度图像及其加高斯白噪声（）和椒盐噪声（噪声密度）后的图像及直方图。function noise(filename)Image=imread(filename);subplot(2,2,1);imshow(Imag

16、e,); % 原图像title(原图像);noise_g=imnoise(Image,gaussian); % 增加高斯白噪声subplot(2,2,2);imshow(noise_g,); %imwrite(noise_g,noise_g.bmp)title(高斯白噪声);noise_s=imnoise(Image,salt & pepper);subplot(2,2,3);imshow(noise_s,); title(黑白象素点噪声(椒盐);%imwrite(noise_s,noise_s.bmp)noise_p=imnoise(Image,speckle); subplot(2,2,4

17、);imshow(noise_p,); %imwrite(noise_p,noise_p.bmp)title(乘积性噪声);Image=imread(graytest2.bmp);subplot(3,2,1);imshow(Image,); % 原图像title(原图像);subplot(3,2,2);imhist(Image);noise_g=imnoise(Image,gaussian,0,0.01); % 增加高斯白噪声subplot(3,2,3);imshow(noise_g,); title(高斯白噪声);imwrite(noise_g,graytest2_g.bmp);subplo

18、t(3,2,4);imhist(noise_g);noise_s=imnoise(Image,salt & pepper,0.02);subplot(3,2,5);imshow(noise_s,); title(黑白象素点噪声(椒盐);imwrite(noise_s,graytest2_s.bmp)subplot(3,2,6);imhist(noise_s);5 分别利用算术平均滤波、加权的线性光滑和热传导方程滤波方法对4题中得到的两幅噪声图像进行去噪处理。算术平均滤波： function aver_filter( Image_name,N,out_filename )%非线性平滑滤波器G1

19、平均%输入参数 文件名 N为模板大小(正奇数)%读取图像if nargin3 disp(输入参数错误！); return;endimage=imread(Image_name);if isrgb(image) image=rgb2gray(image);end% 图像扩展size_m,size_n=size(image);image_ex=zeros(size_m+N-1,size_n+N-1);% 赋值for row=1:size_m for col=1:size_n image_ex(row+(N-1)/2,col+(N-1)/2)=image(row,col); endend% 对每个象

20、素pattern_plate=zeros(N,N); %模板image_new=zeros(size_m,size_n); %滤后图像矩阵 for row=1:size_m for col=1:size_n % 取模板的元素 for row_p=1:N for col_p=1:N pattern_plate(row_p,col_p)=image_ex(row+row_p-1,col+col_p-1); end end % 求模板的平均值得滤后图像矩阵 image_new(row,col)=mean2(pattern_plate); end end% 输出imshow(image_new,);t

21、itle(滤后图像);imwrite( uint8(image_new),out_filename )加权线性光滑：function weight_filter( Image_name,flag,out_filename )% 用各模板进行去噪和锐化 %输入参数 Image_name 文件名 % flag 模板类型标志14if nargin3 disp(请按参数格式输入参数！); return;endif flag=1 & flag=2 & flag=3 & flag=4 disp(参数flag输入错误！); return;endN=3;pattern_plate=zeros(N,N); %模

22、板% 各类模板if flag=1 pattern_plate=1/9 1/9 1/9; 1/9 1/9 1/9; 1/9 1/9 1/9;elseif flag=2 pattern_plate=1/16 1/8 1/16; 1/8 1/4 1/8; 1/16 1/8 1/16; % 平滑滤波器elseif flag=3 pattern_plate=0 -1 0; -1 5 -1; 0 -1 0; % 反锐化掩蔽滤波器elseif flag=4 A=1.7; pattern_plate=0 -1 0; -1 A+5 -1; 0 -1 0; % 高提升滤波器end% 读取图像 image_I=im

23、read(Image_name);if isrgb(image_I) image_I=rgb2gray(image_I);endsubplot(1,2,1);imshow(image_I,);title(原图像);%图像扩展size_m,size_n=size(image_I);image_ex=zeros(size_m+N-1,size_n+N-1);%赋值for row=1:size_m for col=1:size_n image_ex(row+(N-1)/2,col+(N-1)/2)=image_I(row,col); endend%对每个象素temp_plate=zeros(N,N)

24、; %模板image_new=zeros(size_m,size_n); %滤后图像矩阵for row=1:size_m for col=1:size_n %取模板的元素 for row_p=1:N for col_p=1:N temp_plate(row_p,col_p)=image_ex(row+row_p-1,col+col_p-1); end end image_new(row,col)=sum(sum(temp_plate.*pattern_plate); endend%输出subplot(1,2,2);imshow(uint8(image_new);title(滤波后图像);imw

25、rite(uint8(image_new),out_filename);热传导方程滤波：function Heat_Equation( image_name,t,out_filename )% 热传导方程 Ut=U0+t*Delta(U0) % image_name : 原图像文件名% t : 迭代次数 % out_filename: 输出文件名if nargin0) N=t/0.1;else N=(-t)/0.1% 每次叠代的步长为0.1end for index=1:N image_ex=matrix_expand(image_I); % 边界拓展 diff_two_order=zeros

26、(size_r,size_c); for row=2:size_r+1 for col=2:size_c+1 diff_two_order(row-1,col-1)=-4*image_ex(row,col)+image_ex(row-1,col)+image_ex(row+1,col)+image_ex(row,col-1)+image_ex(row-1,col+1); end end if(t0) image_I=image_I+0.1*diff_two_order; else image_I=image_I-0.1*diff_two_order; end imwrite( uint8(im

27、age_I),strcat(Heat_Eq_,num2str(index),.bmp);end cd(.)subplot(1,2,2);imshow(uint8(image_I);if t0 title(热传导方程光滑后图像);else title(Gabor锐化后的图像);endimwrite( uint8(image_I),out_filename );% end%function matrix_ex=matrix_expand(matrix)% 对矩阵进行延拓 size_r,size_c=size(matrix); matrix_ex=zeros(size_r+2,size_c+2); % 中 for row=1:size_r for col=1:size_c matrix_ex(row+1,col+1)=matrix(row,col); end end % 左，右 for