第一篇:实验8 医学图像频域滤波
实验8 医学图像频域滤波与图像复原
实验目的:
1.熟悉医学图像离散傅里叶变化的原理和方法;
2.掌握医学图像频域滤波的原理;
3.掌握使用Matlab中的函数实现医学图像进行频域滤波的方法;
4.掌握使用Matlab中的图像退化与复原的方法;
实验内容:
一、医学图像频域滤波方法与实现
使用imnoise给图像BMRI1_24bit.bmp添加概率为0.2的椒盐噪声,对原图像和加噪声后的图像进行傅氏变换并显示变换后的移中频谱图,然后分别使用Butterworth低、高通滤波器对噪声图像进行低通和高通滤波,显示D0为5,10,20,40时的滤波效果图,并说明存两种滤波效果中所存在的差异及原因。
二、医学图像退化及复原
1、产生带运动和噪声的退化图像
首先使用fspecial 产生一个了len=7,theta=45的移动退化滤波器,然后使用产生的移动退化滤波器作为imfilter的滤波模板,边界填充选项(boundary_options)选择‘circular’,为图像BMRI1_24bit.bmp添加一个移动退化,再使用imnoise产生一个均值为0,方差为0.001的高斯噪声,并将高斯噪声叠加到已经产生移动退化的图像MRIBrain_10.bmp上。
2、对退化图像进行复原
使用DECONVWNR函数对产生退化的MRIBrain_10.bmprp按下面的四种方法进行复原。
(1)省略DECONVWNR中参数NSR;
(2)NSR取一个大于0的数值;
(3)NSR取一个等于1的数值;
(4)NSR取一个大于1的数值;
通过比较上面四种方法的结果,指出哪一种方法的效果最好。
提示:当省略NSR时表示进行直接逆滤波,不省略时表示进行维纳滤波
第二篇:Matlab频域滤波练习
%% fftshift 对数变换,所应用的图片本身很简单,就只有黑白2种颜色 clc clear
f = imread('.imagesdipum_images_ch04Fig0403(a)(image).tif');imshow(f)
title('原始图像')
imfinfo('.imagesdipum_images_ch04Fig0403(a)(image).tif');%此处如果用Imfinfo(f)就会报错fft
%没有居中的傅里叶频谱
F=fft2(f);%进行二维快速傅里叶变换,其结果和DFT的一样,只是计算机的计算速度变快了而已,因而叫fft S=abs(F);%求傅里叶变换后的幅值 figure,subplot(121),imshow(S,[])
title('傅里叶频谱图像1');%title函数一定要放在坐标显示的下一句才有效。
subplot(122),imshow(S),title('傅里叶频谱图像2')
%居中的傅里叶频谱
Fc=fftshift(F);%将频谱图像原点移至图像矩形中间 S1=abs(Fc);
figure,subplot(121),imshow(S1,[]);%加了第二个参数后显示的图像正常
subplot(122),imshow(S1);%当没有第二个参数时,显示的图像为竖线加一些孤立的黑点,why?
%使用对数后视觉增强后的傅里叶频谱 S2=log(1+S1);imshow(S2,[]);
%%用fftshift对数变换显示稍复杂的图片
clc clear
f=imread('.imagesdipum_images_ch04Fig0409(a)(bld).tif');imshow(f);
f=double(f);%其实这句不要试过对后面的变换结果也没有影响 F=fft2(f);
Fc=fftshift(F);S=abs(Fc);
S2=log(1+S);%如果没有这句的话,那么根本看不到细节的图,所以一定要用对数压缩增加对比度
figure,imshow(S2,[])
%%试一下ifft功能
clc clear
f=imread('.imagesdipum_images_ch04Fig0409(a)(bld).tif');imshow(f);
f=double(f);f(1:8,1:8)F=fft2(f);
f=ifft2(F);%取傅里叶变换后的反傅里叶变换,直接是整数,并不需要像某些代码一样取real部分 f(1:8,1:8)
%%先0扩充再滤波 clc clear f=imread('.imagesdipum_images_ch04Fig0405(a)(square_original).tif');imshow(f);
[M N]=size(f);
F=fft2(f);%没有经过0扩充直接计算fft sig=10;%高斯滤波参数
H=lpfilter('gaussian',M,N,sig);
G=H.*F;%加了.号的乘法表示对应每个元素都相乘,没加.号的乘法说明是真正的矩阵乘法
g=real(ifft2(G));
figure,imshow(g,[]);
PQ=paddedsize(size(f));%PQ为0扩展的面积,这里设置为与图像的大小相同?这里size(f)的大小为256*256
Fp=fft2(f,PQ(1),PQ(2));%如果fft2有3个参数的话,则是进行了0扩展了
Hp=lpfilter('gaussian',PQ(1),PQ(2),2*sig);%构造高斯滤波器 Gp=Hp.*Fp;%对0扩展后的图像进行高斯滤波 gp=real(ifft2(Gp));figure,imshow(gp,[])'
gpc=gp(1:size(f,1),1:size(f,2));%取gp图像的1到f的第一维的行,1到f第二维的列部分,即取与图像大小相同的部分
figure,imshow(gpc,[]);%此时显示的结果应该与g图像一样。
%直接使用空间域滤波
h=fspecial('gaussian',15,7);gs=imfilter(f,h);
figure,imshow(gs,[])
%% clc clear
f=imread('.imagesdipum_images_ch04Fig0405(a)(square_original).tif');f=double(f);imshow(f);
PQ=paddedsize(size(f));%size(f)的大小为256*256
sig=10;
H = lpfilter('gaussian',PQ(1),PQ(2),2*sig);% PQ=[512 512] figure,mesh(abs(H(1:10:512,1:10:512)));%画出滤波器的三维空间图 g=dftfilt(f,H);
figure,imshow(g,[])
%%空间滤波和频域滤波的比较
clc clear
f=imread('.imagesdipum_images_ch04Fig0409(a)(bld).tif');imshow(f);
F=fft2(f);
S=fftshift(log(1+abs(F)));
S=gscale(S);%用gscale来进行将数据缩放到默认值0~255 figure,imshow(S);
h=fspecial('sobel');%构造sobel空间滤波器 figure,freqz2(h);%查看滤波器的图形
PQ=paddedsize(size(f));
H=freqz2(h,PQ(1),PQ(2));%将空间滤波器h转换成频域滤波器H,但是滤波器的原点在矩阵的中心处
H1=ifftshift(H);%原点迁移到左上角
figure,mesh(abs(H1(1:20:1200,1:20:1200)));figure,subplot(121),imshow(abs(H),[]);subplot(122),imshow(abs(H1),[]);
gs=imfilter(double(f),h);%空间滤波,其实就是边缘检测 gf=dftfilt(f,H);%频域滤波
subplot(121),imshow(gs,[]),subplot(122),imshow(gf,[]);%将2种滤波结果显示出来
figure,imshow(abs(gs)>0.2*abs(max(gs(:))));%只显示最大值的20%的边缘
figure,imshow(abs(gf)>0.2*abs(max(gf(:))));
d=abs(gs-gf);a=max(d(:))b=min(d(:))
%%构造低通滤波器
clc clear
f=imread('.imagesdipum_images_ch04Fig0413(a)(original_test_pattern).tif');
imshow(f)
F1=fft2(f);
imshow(log(1+abs(fftshift(F1))),[]);%直接显示取对数后的傅里叶变换
PQ=paddedsize(size(f));
[U V]=dftuv(PQ(1),PQ(2));%这个函数还真不明白什么意思,help中解释的是计算网格频率矩阵U和V,这2个矩阵是用来频率滤波的 D0=0.05*PQ(2);
F=fft2(f,PQ(1),PQ(2));%用0扩充大小为PQ(1)*PQ(2),然后fft变换 figure,imshow(log(1+abs(fftshift(F))),[]);%显示0扩充后的fft变换图
H=exp(-(U.^2+V.^2)/(2*(D0^2)));%构造频域滤波算子 figure,imshow(fftshift(H),[]);%显示平移后滤波器算子 g=dftfilt(f,H)
figure,imshow(g,[]);
%%绘制一个低通滤波器的mesh图 clc clear
H1=lpfilter('gaussian',500,500,50)%构造一个中心点在左上角的高斯低通滤波器,size为500*500
H2=fftshift(H1);%将中心点平移至矩阵中心
mesh(H1(1:10:500,1:10:500)),figure,mesh(H2(1:10:500,1:10:500));%绘出2者的mesh图
axis([0 50 0 50 0 1]);%xy轴范围0~50,z轴范围0~1
figure,subplot(121),imshow(H1,[]),subplot(122),imshow(H2,[]);%绘出2者的灰度图
%%绘制一个低通滤波器的mesh图 clc clear
H=fftshift(hpfilter('ideal',500,500,100));%构造理想的高通滤波器 mesh(H(1:10:500,1:10:500));axis([1 50 1 50 0 1]);colormap([0 0 0]);%mesh网格全部用黑色画 axis off%关掉坐标轴显示 grid off%关掉网格显示 figure,imshow(H,[])
%%使用高通滤波器
clc clear
f=imread('.imagesdipum_images_ch04Fig0413(a)(original_test_pattern).tif');
imshow(f)
title('原始图像')
PQ=paddedsize(size(f));D0=0.05*PQ(1);
H=hpfilter('gaussian',PQ(1),PQ(2),D0);g=dftfilt(f,H);
figure,imshow(g,[])
%%将高通滤波和直方图均衡结合起来
clc clear
f=imread('.imagesdipum_images_ch04Fig0419(a)(chestXray_original).tif');imshow(f)
title('原始图像')
PQ = paddedsize(size(f));D0 = 0.05*PQ(1);
HBW=hpfilter('btw',PQ(1),PQ(2),D0,2);%构造高通Butterworth滤波器,截断频率为D0 H=0.5+2*HBW;
gdw=dftfilt(f,HBW);
gbw=gscale(gdw);%将灰度值自动缩放到0~255之间 figure,subplot(121),imshow(gdw,[]);title('高通滤波后图像');
gbe=histeq(gbw,256);%直方图均衡化 subplot(122),imshow(gbe,[]);title('高通滤波且直方图均衡化后图像');
ghf=dftfilt(f,H);%H是HBW的线性变换 ghf=gscale(ghf);
figure,subplot(121),imshow(ghf,[]);title('高频强调滤波后图像');
ghe=histeq(ghf,256);
subplot(122),imshow(ghe,[]);title('高频强调且均衡化后图像');
%%fftshift和ifftshift的加深理解 clc clear
A=[2 0 0 1
0 0 0 0
0 0 0 00 0 4]
B=fftshift(A)%中心点平移
C=fftshift(fftshift(A))%还原成A
D=ifftshift(fftshift(A))%也同样还原成A E=ifftshift(A)%平移结果和B一样
第三篇:图像滤波总结
数字图像处理:各种变换滤波和噪声的类型和用途总结
一、基本的灰度变换函数 1.1.图像反转
适用场景:增强嵌入在一幅图像的暗区域中的白色或灰色细节,特别是当黑色的面积在尺寸上占主导地位的时候。
1.2.对数变换(反对数变换与其相反)
过程:将输入中范围较窄的低灰度值映射为输出中较宽范围的灰度值。用处:用来扩展图像中暗像素的值,同时压缩更高灰度级的值。特征:压缩像素值变化较大的图像的动态范围。
举例:处理傅里叶频谱,频谱中的低值往往观察不到,对数变换之后细节更加丰富。
1.3.幂律变换(又名:伽马变换)
过程:将窄范围的暗色输入值映射为较宽范围的输出值。
用处:伽马校正可以校正幂律响应现象,常用于在计算机屏幕上精确地显示图像,可进行对比度和可辨细节的加强。
1.4.分段线性变换函数
缺点:技术说明需要用户输入。优点:形式可以是任意复杂的。
1.4.1.对比度拉伸:扩展图像的动态范围。
1.4.2.灰度级分层:可以产生二值图像,研究造影剂的流动。1.4.3.比特平面分层:原图像中任意一个像素的值,都可以类似的由这些比特平面对应的二进制像素值来重建,可用于压缩图片。
1.5.直方图处理
1.5.1直方图均衡:增强对比度,补偿图像在视觉上难以区分灰度级的差别。作为自适应对比度增强工具,功能强大。
1.5.2直方图匹配(直方图规定化):希望处理后的图像具有规定的直方图形状。在直方图均衡的基础上规定化,有利于解决像素集中于灰度级暗端的图像。
1.5.3局部直方图处理:用于增强小区域的细节,方法是以图像中的每个像素邻域中的灰度分布为基础设计变换函数,可用于显示全局直方图均衡化不足以影响的细节的显示。1.5.4直方图统计:可用于图像增强,能够增强暗色区域同时尽可能的保留明亮区域不变,灵活性好。
二、基本的空间滤波器 2.1.平滑空间滤波器
2.1.1平滑线性滤波器(均值滤波器)
输出:包含在滤波器模板邻域内的像素的简单平均值,用邻域内的平均灰度替代了图像中每个像素的值,是一种低通滤波器。结果:降低图像灰度的尖锐变化。
应用:降低噪声,去除图像中的不相关细节。负面效应:边缘模糊。
2.1.2统计排序滤波器(非线性滤波器)举例:中值滤波器。过程:以滤波器包围的图像区域中所包含图像的排序为基础,然后使用统计排序结果决定的值取代中心区域的值。
用处:中值滤波器可以很好的解决椒盐噪声,也就是脉冲噪声。
2.2.锐化空间滤波器
2.2.1拉普拉斯算子(二阶微分)
作用:强调灰度的突变,可以增强图像的细节。
2.2.2非锐化掩蔽和高提升滤波
原理:原图像中减去一幅非锐化(平滑处理)的版本。背景:印刷和出版界使用多年的图像锐化处理。
高提升滤波:原图减去模糊图的结果为模板,输出图像等于原图加上加权后的模板,当权重为1得到非锐化掩蔽,当权重大于1成为高提升滤波。
2.2.3梯度锐化(一阶微分对)
含义:梯度指出了在该位置的最大变化率的方向。
用处:工业检测,辅助人工检测产品的缺陷,自动检测的预处理。
三、基本的频率滤波器 3.1.1理想低(高)通滤波器 特性:振铃现象,实际无法实现。
用处:并不实用,但是研究滤波器的特性很有用。
3.1.2布特沃斯低(高)通滤波器
特点:没有振铃现象,归功于在低频和高频之间的平滑过渡,二阶的布特沃斯低通滤波器是很好的选择。
效果:比理想低(高)通滤波器更平滑,边缘失真小。截止频率越大,失真越平滑。
3.1.3高斯低(高)通滤波器 特点:没有振铃。
用处:任何类型的人工缺陷都不可接受的情况(医学成像)。
3.1.4钝化模板,高提升滤波,高频强调滤波 用处:X射线,先高频强调,然后直方图均衡。
3.1.5同态滤波
原理:图像分为照射分量和反射分量的乘积。
用处:增强图像,锐化图像的反射分量(边缘信息),例如PET扫描。
3.1.6选择性滤波
3.1.6.1带阻滤波器和带通滤波器。作用:处理制定频段和矩形区域的小区域。
3.1.6.2陷阱滤波器
原理:拒绝或通过事先定义的关于频率矩形中心的一邻域。应用:选择性的修改离散傅里叶变换的局部区域。
优点:直接对DFT处理,而不需要填充。交互式的处理,不会导致缠绕错误。用途:解决莫尔波纹。
四、重要的噪声概率密度函数 4.1.高斯噪声
特点:在数学上的易处理性。
4.2瑞利噪声
特点:基本形状向右变形,适用于近似歪斜的直方图。
4.3爱尔兰(伽马)噪声
特点:密度分布函数的分母为伽马函数。
4.4指数噪声
特点:密度分布遵循指数函数。
4.5均匀噪声 特点:密度均匀。
4.6脉冲噪声(双极脉冲噪声又名椒盐噪声)
特点:唯一一种引起退化,视觉上可以区分的噪声类型。
五、空间滤波器还原噪声 5.1均值滤波器 5.1.1算术均值滤波器
结果:模糊了结果,降低了噪声。适用:高斯或均匀随机噪声。5.1.2几何均值滤波器
结果:和算术均值滤波器相比,丢失的图像细节更少。适用:更适用高斯或均匀随机噪声。
5.1.3谐波均值滤波器
结果:对于盐粒噪声(白色)效果较好,但不适用于胡椒噪声(黑色),善于处理高斯噪声那样的其他噪声。
5.1.4逆谐波均值滤波器
结果:适合减少或在实际中消除椒盐噪声的影响,当Q值为正的时候消除胡椒噪声,当Q值为负的时候该滤波器消除盐粒噪声。但不能同时消除这两种噪声。适用:脉冲噪声。
缺点:必须知道噪声是明噪声还是暗噪声。
5.2统计排序滤波器 5.2.1中值滤波器
适用:存在单极或双极脉冲噪声的情况。
5.2.2最大值滤波器
作用:发现图像中的最亮点,可以降低胡椒噪声。
5.2.2最小值滤波器
作用:对最暗点有用,可以降低盐粒噪声。
5.2.3中点滤波器
作用:结合统计排序和求平均,对于随机分布噪声工作的很好,如高斯噪声或均匀噪声。5.2.4修正的阿尔法均值滤波器
作用:在包括多种噪声的情况下很有用,例如高斯噪声和椒盐噪声混合。
5.3自适应滤波器
5.3.1自适应局部降低噪声滤波器
作用:防止由于缺乏图像噪声方差知识而产生的无意义结果,适用均值和方差确定的加性高斯噪声。
5.3.1自适应中值滤波器
作用:处理更大概率的脉冲噪声,同时平滑非脉冲噪声时保留细节,减少诸如物体边界粗化或细化等失真。
5.4频率域滤波器消除周期噪声 5.4.1带阻滤波器
应用:在频率域噪声分量的一般位置近似已知的应用中消除噪声
5.4.2带通滤波器
注意:不能直接在一张图片上使用带通滤波器,那样会消除太多的图像细节。用处:屏蔽选中频段导致的结果,帮助屏蔽噪声模式。
5.4.3陷阱滤波器
原理:阻止事先定义的中心频率的邻域内的频率。作用:消除周期性噪声。
5.4.4最佳陷阱滤波
作用:解决存在多种干扰分量的情况。
第四篇:实验四频域分析
实验四连续信号与系统的频域分析
一、实验目的:
1、绘制非周期信号的频谱。
2、绘制系统的幅频及相频响应曲线。
二、实验内容
1、非周期信号的频谱
调出下列程序,并观察信号的频谱。
例:求单边指数信号f(t)e2tu(t)的傅里叶变换,并画出f(t)及其幅度谱和相位谱图。syms t w phase im re;%定义符号变量
f=exp(-2*t)*sym('heaviside(t)');
F=fourier(f);
subplot(311)
ezplot(f);
axis([-1 2.5 0 1.1]);
subplot(312)
ezplot(abs(F));%绘制幅度谱
im=imag(F);%计算F(jw)虚部
re=real(F);%计算F(jw)实部
phase=atan(im/re)%计算相位谱
subplot(313)
ezplot(phase);
1作业1:试画出矩形信号g(t)0t12的幅度频谱,观察其频率特性。1t22、MATLAB提供了函数freqs来实现连续系统频率响应H(j)的分析。该函数可以求出系统频率响应的数值解,并可绘出系统的幅频及相频响应曲线。调用格式如下:
(1)H=freqs(B,A,W);B为系统频率响应分子多项式系数,或者微分方程的右端系数,A为系统频率响应分母多项式系数,或微分方程左端系数。W为形如W1:P:W2的频率范围,P为频率采样间隔。输出参量H为返回在W所定义的频率点上,系统频率响应的样值。
abs(H):求H的幅度响应;angle(H):求H的相位响应。
作业2:某连续时间系统的频率响应为H(j)j3,求系统频率响应的样值,并绘出幅度2j3j2
响应曲线和相位响应曲线。
(2)freqs(B,A);该调用格式并不返回系统频率响应的样值,而是以伯特图的方式绘出系统的幅度响应和相位响应曲线。
j222500作业3:已知某系统的频率响应H(j),用以上命令绘制幅度响应和相位响应j2200j20000
曲线,并分析该系统的频率特性。
第五篇:0904057医学图像处理教学大纲
《医学图像处理》课程教学大纲
一、课程基本信息
课程编号:0904057 课程中文名称:医学图像处理
课程英文名称:Medical Image Processing 课程性质:专业主干课程 考核方式:考试 开课专业:生物医学工程 开课学期:7 总学时: 48(其中理论38学时,上机10学时)总学分:3
二、课程目的
本课程是一门专业基础课,目的是为了加强生物医学工程专业学生用计算机进行医学图像处理能力的培养。通过讲述bmp图像、锐化增强、图像分割、图像变换、图像的识别等内容,使学生了解医学图像处理所必需的基础知识,掌握医学图像处理的基本技能和实际应用的方法,为今后结合本专业开展相应的研究打下良好的基础。
三、教学基本要求(含素质教育与创新能力培养的要求)
本门课程分为理论学习和上机实践两个部分,要求学生掌握医学图像处理所必需的基础知识的同时,重点培养学生的实际应用能力。
1、对用计算机进行医学图像处理的重要性和特殊性有明确的认识。
2、熟悉bmp图像显示的方法和相应的程序设计。
3、熟悉图像的点运算、图像的代数运算、图像的几何运算的算法和相应的程序设计。
4、熟悉直方图增强、锐化增强、局部增强、伪彩色增强、平滑的算法和相应的程序设计。
5、熟悉基于边界的图像分割、阈值分割、基于区域增长或分裂的分割算法和相应的程序设计。
四、教学内容与学时分配
第一章
数字图像的形成和图像处理系统(4学时)数字图像的形成、图像处理系统的构成,不同成像技术产生的医学图像 第二章
图像的运算(6学时)
图像的点运算、图像的代数运算,图像的几何运算 第三章
图像增强(6学时,课程重点)
直方图增强,图像平滑,锐化增强、局部增强、伪彩色增强 第四章
图像的分割(6学时,课程重点、课程难点)
基于边界的图像分割、阈值分割,基于区域增长或分裂的分割、分割效果的评价 第五章
图像的表达与描述(4学时,课程难点)
目标外特性的表达与描述,目标内特性的表达与描述、目标特性描述的标定 第六章
图像变换(4学时)
傅立叶变换,傅立叶变换在图像处理中的应用 第七章
图像的识别(4学时,课程难点)图像相似性的测量,图像的特征,图像的分类 第八章
医学虚拟现实系统概述(2学时)医学虚拟现实系统的组成及其主要应用等 考试(2学时)
五、教学方法及手段(含现代化教学手段及研究性教学方法)
所运用的教学方法和教学手段有理论教学(课堂讲授)、使用多媒体教学、课堂讨论、上机实践等。
六、实验(或)上机内容
1、编写程序,对医学图像进行平滑处理
2、编写程序,对医学图像进行滤波处理
3、编写程序,对医学图像进行基于阈值的图像分割
4、显微医学图像分析
5、医学虚拟现实演示系统
七、先修课程
先修课程:程序设计基础(C语言)、大学计算机基础、计算机软件基础。
八、教材及主要参考资料
[1] 章鲁,顾顺德,陈瑛.医学图像处理[M].上海:上海科学技术出版社,2002.[2] 谷口庆治.数字图像处理[M].北京:科学出版社&共立出版社,2002.九、课程考核方式
平时成绩占20%,闭卷笔试占80%。
撰写人签字:
院(系)教学院长(主任)签字: