350 likes | 640 Vues
第 7 章 MATLAB 工具箱简介. 7.1 信号处理工具箱. 信号处理工具箱式算法文件的集合,大部分为 M 文件,覆盖了经典信号处理理论的大多数内容,在语音信号处理、生物医学工程、实时控制理论、雷达信号处理等多个研究领域得到成功的应用。 按功能,信号处理工具箱中的函数可分为: 1 )信号与波形的产生; 2 )信号变换处理; 3 )滤波器的设计、分析与实现; 4 )随机信号处理与谱估计; 5 )时 - 频分析; 6 )到谱分析;等等。. 信号的表示与产生. 信号处理的基本数据可以采用矩阵来描述:每一列代表一个通道,每一行对应一个采样点。
E N D
7.1 信号处理工具箱 • 信号处理工具箱式算法文件的集合,大部分为M文件,覆盖了经典信号处理理论的大多数内容,在语音信号处理、生物医学工程、实时控制理论、雷达信号处理等多个研究领域得到成功的应用。 • 按功能,信号处理工具箱中的函数可分为: 1)信号与波形的产生; 2)信号变换处理; 3)滤波器的设计、分析与实现; 4)随机信号处理与谱估计; 5)时-频分析; 6)到谱分析;等等。
信号的表示与产生 • 信号处理的基本数据可以采用矩阵来描述:每一列代表一个通道,每一行对应一个采样点。 • 例如:产生一个采样频率为1000Hz的离散信号 1)定义时间轴: fs=1000;t=0:1/fs:1; 2) 定义采样得到的信号: y=sin(2*pi*50*t)+2*sin(2*pi*120*t+0.5*randn(size(t)); (注:randn(1,n)产生均值为0,方差为1的白噪声) 3) 绘图:plot(t,y)
经典信号产生函数 • X=sawtooth(t):产生周期为2*pi的锯齿波; • X=square(t,duty):产生周期为2*pi的方波; • Y=sinc(x):产生sinc函数的波形; • 单位冲激信号x=zeros(1,n);x(1)=1; • 单位阶跃信号x=ones(1,n); • 复指数信号n=0:N-1;x=exp((sigma+jw)*n); • ….
信号的基本运算 • 信号加 x=x1+x2; • 信号乘 x=x1.*x2; • 幅度变化 y=alpha*x; • 位移 y=[zeros(1,k),x] 对应于 y(n)=x(n-k); • 折叠 y=fliplr(x) 对应于 y(n)=x(-n); • 采样和 y=sum(x(n1:n2)) 对应于 y(n)=x(n1)+..+x(n2); • 采样积 y=prod(x(n1:n2)) 对应于 y(n)=x(n1)*…*x(n2); • 信号能量 E=sum(abs(x).^2) • 信号功率 p=sum(abs(x).^2)/N • …
离散傅立叶变换DFT • A=dftmtx(n)返回nXn的DFT变换矩阵,y=x*A则返回x的DFT变换至y; • Ai=conj(dftmtx(n))/n返回离散傅立叶反变换矩阵。 t=0:.01:1; x=sinc(2*pi*5*t); A=dftmtx(length(t)); y=x*A; subplot(1,2,1); plot(t,x); subplot(1,2,2); plot(t,y); Ai=conj(A)/length(t); yi=y*Ai; figure; plot(t,yi)
快速离散傅立叶变换FFT y=fft(x,n):计算n点FFT; y=fftshift(x):将零频位置移到频谱中心; y=ifft(x,n):快速离散傅立叶反变换。 clear t=linspace(0,pi,1000); data=exp(-t./2); power=fft(data); mag=abs(power);%振幅 phs=unwrap(angle(power));%相位 subplot(2,1,1); fre=[0,2,4,6,8,10]; stem(fre,mag(1:6),'r.'); subplot(2,1,2); stem(fre,phs(1:6),'b.');
离散余弦变换DCT • 与傅立叶变换相比,离散余弦变换具有更好的压缩性能,在数据压缩和数据通信中得到广泛的应用。 • y=dct(x,n):计算指定的n点DCT; • y=idct(x,n):计算指定的n点DCT反变换。 load mtlb subplot(2,2,1),plot(mtlb); x=dct(mtlb); subplot(2,2,2),plot(x); x(abs(x)<0.5)=0; subplot(2,2,3),plot(x); y=idct(x); subplot(2,2,4),plot(y);
希尔伯特变换HCT • y=hilbert(x):返回x的解析信号y,其实部为原信号,虚部为DCT得到的信号。 t=0:1/1023:1; x=sin(2*pi*60*t); y=Hilbert(x); plot(t(1:50),real(y(1:50))); hold on plot(t(1:50),imag(y(1:50)),':');
7.2 数字图像处理工具箱 • MATLAB的图像处理工具箱由一系列支持图像处理操作的函数组成,按功能可以分为: 图像显示;图像文件输入输出;几何操作;图像分析与增强;图像变换;图像滤波;二值图像处理;等等。 • MATLAB的图像是将数据阵列元素以double(64位)的浮点类型存储,也支持uint8和uint16类型。 • MATLAB支持5种图像类型:索引图象、灰度图像、二值图像、RGB图像和多帧图像阵列。 • 多帧图象阵列是由多帧图像组成,但组成的图像必须为前面中的一种,cat函数可以用来将相同大小的几个独立图像形成多帧文件。
索引图象举例 Image matrix colormap
1 2 Columns 1 2 Rows 图像坐标系统:像素坐标系
图像坐标系统:空间坐标系 将像素看成是一个有区域的矩形块。
改变矩阵图像的坐标原点 >> A=magic(5); >> x=[19.5 23.5]; >> y=[8.0 12.0]; >> image(A,'XData',x,'YData',y); >> axis image >> colormap(jet(25)); >> A A = 17 24 1 8 15 23 5 7 14 16 4 6 13 20 22 10 12 19 21 3 11 18 25 2 9
图像文件的读写 • A=imread(filename,fmt):读文件数据至A。 fmt的可能取值有:‘jpg’或‘jpeg’(联合影像专家组格式JPEG);’tif’或‘tiff’(标志图像文件格式TIFF);’gif’(图形交换格式GIF);’bmp’(位图格式BMP);’png’(可移动网络图形格式);等等。 • imwrite(A,filename,fmt):写图像数据A至文件。 • imfinfo filename:filename要包含完整路径,例如: >> imfinfo c:\MATLAB6P5\toolbox\simulink\blocks\blocksets.bmp
图像类型判断及转换 • isbw:是否为二值图像;isgray:是否为灰度图像;isind:是否为索引图象;isrgb:是否为RGB图像。 • gray2ind:灰度图像转换为索引图象;im2bw:设置阀值转换为二值图像;im2double(uint8, uint16):转换为double(uint8或uint16)类型;mat2gray矩阵转换为灰度图像;其他: ind2gray, ind2rgb, mat2gray, rgb2gray, rgb2ind等 • 还可用cat函数将灰度图像A转换为RGB图像,为: RGB=cat(3,A,A,A);
索引图转换为灰度图 >> load trees >> A=ind2gray(X,map); >> imshow(X,map); >> figure %在同一窗口显示有干扰 %可以采用subimage来在同一窗口显 %示多幅图像。 >> imshow(A)
图像显示 • 一般显示:imshow或image或imagesc,例如 A=imread(‘flowers.tif’);image(A); • 投影显示(纹理成图):warp,例如: I=imread('flowers.tif'); [x,y,z]=cylinder; subplot(121),warp(x,y,z,I); [x,y,z]=sphere(50); subplot(122),warp(x,y,z,I);
多帧图像阵列逐帧显示 • 利用imshow,image或subimage显示一帧图像,例如: >> load mri >> subplot(221); >> colormap(map); >> subimage(D(:,:,:,1)); >> subplot(222); >> subimage(D(:,:,:,2)); >> subplot(223); >> subimage(D(:,:,:,3)); >> subplot(224); >> subimage(D(:,:,:,4));
多帧图像阵列拼接显示 • montage(I):拼接灰度图像;montage(X,map):拼接索引图象;montage(RGB):拼接RGB图像,例如: >> load mri >> montage(D,map)
多帧图像阵列转化成电影显示 • mov=immovie(X,map):将索引图象转换成电影。 >> load mri >> mov=immovie(D,map); >> colormap(map),movie(mov);
图像几何运算 • 改变图像尺寸:imresize,例如:B=imresize(A,[20,20]); • 图像旋转: imrotate,例如:B=imrotate(A,30); • 图像裁剪:imcrop,例如:B=imcrop(A,[10,10,150,150]);其中[10,10,150,150]为用户指定的矩形裁剪区[xmin,ymin,width,height];
图像的数学变换:FFT • B=fft2(A):二维快速傅立叶变换 >> load imdemos saturn2 %提取imdemos中的saturn2 >> imshow(saturn2) >> B=fftshift(fft2(saturn2)); %将变换后的频谱中心从矩阵的始点移到矩阵中心。 >> figure >> imshow(log(abs(B)),[]); >> colormap(jet(64)) >> colorbar 傅立叶振幅谱 原始图像
图像的数学变换:DCT • B=dct2(A):计算A的离散余弦变换。例如: clear; RGB=imread('autumn.tif'); A=rgb2gray(RGB); B=dct2(A); imshow(A); figure; imshow(log(abs(B)),[]); colormap(jeT(64)); colorbar; B(abs(B)<10)=0; %除掉值小于10的数据(压缩) C=idct2(B)/255;%dct反变换 figure imshow(C); A C
像素选取 • pixval:用于交互显示像素数据值 • impixel:返回被选中的若干点的数据值,选点方式可以是提供坐标或用鼠标选择(双击左键结束选择)。 >> imshow canoe.tif >> pixval >> imshow canoe.tif >> vals=impixel vals = 0.5490 0.4824 0.3216 0.1608 0.1608 0.1922 0.2902 0.2588 0.2235 0.2588 0.2902 0.2235
图像轮廓 • imcontour:显示一幅图像的轮廓图。 >> A=imread('blood1.tif'); >> subplot(121),imshow(A); >> subplot(122),imcontour(A);
图像边缘检测 • 边缘检测的基本思想是首先利用边缘增强算子,突出图像中的局部边缘,然后定义像素的“边缘强度”,通过设置门限的方法提取边缘点集。 • 常用的边缘算子有robert、sobel、prewitt、log和canny五个算子。 >> A=imread('blood1.tif'); >> bw1=edge(A,'sobel'); >> bw2=edge(A,'robert'); >> bw3=edge(A,'log'); >> subplot(221),imshow(A); >> subplot(222),imshow(bw1); >> subplot(223),imshow(bw2); >> subplot(224),imshow(bw3);
图像增强 • 图像增强包括:进行灰度级调整和进行噪声消除。 • 灰度级调整是将图像的灰度级映射到一个新的范围,在新的范围内能够容纳图像多数的灰度级,是图像的对比更加鲜明。 • 消除噪声的方法主要有:均值滤波(用像素邻域中所有像素的均值来取代该像素的值)、中值滤波(用像素邻域中所有像素的中间值来取代该像素的值)、自适应滤波(依据局部图像的差异来调整滤波器的参数,对局部差异大的地方进行小的平滑操作,反之进行大的平滑操作)。
灰度级调整 • 1)显示图像及灰度级信息; • 2)依据灰度级信息进行调整。 J=imadjust(I,[low_in high_in],[low_out high_out],gama) 灰度级从[low_in high_in]调整到[low_out high_out],gama指定图像I和J之间的映射关系,为1表示线性映射,大于1表示强映射,小于1表示弱映射。 MATLAB提供了histeq函数用来自动调整直方图。 A=imread('pout.tif'); subplot(221); imshow(A); subplot(222); imhist(A);%显示图像对应直方图 %原图绝大多数灰度级集中在75~170 B=imadjust(A,[75/255 170/255],[]); %将灰度级扩展到[0,255] subplot(223); imshow(B); subplot(224); imhist(B);
噪声消除 • MATLAB工具箱提供imnoise函数,用来引入各种噪声。 • 线性滤波 B=filter2(h,A),其中h为计算模块,可由h=fspecial(type)等函数生成,A为图像矩阵。type的类型有‘gaussian’(Gussian低通滤波)、‘average’(均值滤波)等。 • 中值滤波:B=medfilt2(A). • 自适应滤波: B=wiener2(A,[m,n]):二维维纳滤波,[m,n]指定滤波器的窗口大小为mXn。 A=imread('blood1.tif'); B=imnoise(A,'gaussian',0,0.005); %均值为0,方差为0.005 h=fspecial('average',5); %h为一个5X5的均值滤波器 C=uint8(round(filter2(h,B))); D=medfilt2(B); E=wiener2(B,[5 5]); subplot(221),imshow(A); subplot(222),imshow(B); subplot(223),imshow(C); subplot(224),imshow(D);