matlab信号处理基础（音频和图像）

提示：文章写完后，目录可以自动生成，如何生成可参考右边的帮助文档

---

前言

【实验目的】
了解音频和图像数据系数特点， 掌握音频和图像文件的离散傅里叶、离散余弦和离散小波变换等基本操作。
【实验环境】
(1) Window 11 操作系统;
(2) Matlab 2020版本软件;
(3) BMP 格式图像文件;
(4) Wav 格式音频文件

一、用离散傅里叶变换分析合成音频和图像

1.分析合成音频文件

①读取音频文件数据

clc;
clear;
len=40000;
[fn,pn]=uigetfile('*.wav','请选择音频文件');
filename=strcat(pn,fn);
[x,fs]=audioread(filename,[1,len]);
 
%新版本的matlab中已经不再支持wavread()函数，替代函数为audioread(filename,N),
%其中N必须为[m,n]格式，如[2,100],且m,n均为正数
%作用：读取微软音频格式 (wav) 文件内容
%输入参数: filename表示音频文件名，字符
%返回值:音频样点，浮点型

②一维离散傅里叶变换

xf=fft(x);%fft()对输入参数进行一维离散傅里叶变换并返回其系数，对应频率从0到fs（采样频率) ，
f1=[0:len-1]*fs/len;
xff=fftshift(xf);%fftshift()将零频对应系数移至中央
h1=floor(len/2);
f2=[-h1:h1]*fs/len;%计算离散样点对应的频率值，以便更好地观察频谱

③一维离散傅里叶逆变换

xsync=ifft(xf);%ifft()对输入参数进行一维离散傅里叶逆变换并返回其系数

④观察结果

figure;%figure(n)表示建第n个图形窗。
subplot(2,2,1);plot(x);title('original audio');
subplot(2,2,2);plot(xsync);title('synthesize audio');
subplot(2,2,3);plot(f1,abs(xf));title('fft coef. of audio');
subplot(2,2,4);plot(f2(1:len),abs(xff));title('fftshift coef. of auio');

结果展示：

2.分析合成图片文件

①读取图像文件数据

[fn,pn]=uigetfile('*.bmp','请选择图像文件');
[x,map]=imread(strcat(pn,fn),'bmp');
I=rgb2gray(x);%rgb2ray()将RGB图像转换为灰度图
 
%函数原型 :A=imread(filename,fmt)
%功能:读取 fmt 指定格式的图像文件内容。
%输入参数 :filename 表示图像文件名，字符串。
%fmt表示图像文件格式名，字符串、函数支持的图像格式包括:JPEG,TIFF,GIF,BMP等，
%当参数不包括文件格式名时，函数尝试推断出文件格式
%返回参数:A表示图像数据内容,整型。

②二维离散傅里叶变换

xf=fft2(I);%fft2()对输入参数进行二维离散傅里叶变换并返回其系数
xff=fftshift(xf);%fftshift()将零频对应系数移至中央

③二维离散傅里叶逆变换

xsync=ifft2(xf);%ifft()对输入参数进行二维离散傅里叶逆变换并返回其系数

④观察结果

figure;
subplot(2,2,1);imshow(x);title('original image');
subplot(2,2,2);imshow(uint8(abs(xsync)));title('synthesize image');
subplot(2,2,3);mesh(abs(xf));title('fft coef. of image');
subplot(2,2,4);mesh(abs(xff));title('fftshift coef. of image');
 
%imshow是二维数据绘图函数，mesh通过三维平面显示数据

结果展示：

二、用离散余弦变换分析合成音频和图像

1.分析合成音频文件数据：

①读取音频文件数据

clc;
clear;
len=40000;
[fn,pn]=uigetfile('*.wav','请选择音频文件');
filename=strcat(pn,fn);
[x,fs]=audioread(filename,[1,len]);

②:一维离散余弦变换

xf=dct(x);

③一维离散余弦逆变换

xsync=idct(xf);
[row,col]=size(x);
xff=zeros(row,col);
xff(1:row,1:col)=xf(1:row,1:col);
y=idct(xff);

④观察结果

figure;
 
%xf=fft(x);
f1=[0:len-1]*fs/len;
%xff=fftshift(xf);
h1=floor(len/2);
f2=[-h1:h1]*fs/len;
 
 
subplot(2,2,1);plot(x);title('original audio');
subplot(2,2,2);plot(xsync);title('synthesize audio');
subplot(2,2,3);plot(f1,abs(xf));title('fft coef. of audio');
subplot(2,2,4);plot(f2(1:len),abs(xff));title('fftshift coef. of auio');

结果展示：

2.分析合成图像文件数据：

①读取图像文件数据

[fn,pn]=uigetfile('*.bmp','请选择图像文件');
[x,map]=imread(strcat(pn,fn),'bmp');
I=rgb2gray(x);

②二维离散余弦变换

xf=dct2(I);

③二维离散余弦逆变换

xsync=uint8(idct(xf));
[row,col]=size(I);
lenr=round(row*4/5);
lenc=round(col*4/5);
xff=zeros(row,col);
xff=zeros(row,col);
xff(1:lenr,1:lenc)=xf(1:lenr,1:lenc);
y=uint8(idct2(xff));

④观察结果

figure;
subplot(2,2,1);imshow(x);title('original image');
subplot(2,2,2);imshow(uint8(abs(xsync)));title('synthesize image');
subplot(2,2,3);imshow(uint8(abs(y)));title('part synthesize image');
subplot(2,2,4);mesh(abs(xff));title('fftshift coef. of image');

结果展示：

三、用离散小波变换分析合成音频和图像

1.分析合成音频文件

①读取音频文件数据

clc;
clear;
len=40000;
[fn,pn]=uigetfile('*.wav','请选择音频文件');
filename=strcat(pn,fn);
[x,fs]=audioread(filename,[1,len]);
②:一维离散小波变换
x=x(:,1);%将多声道转换为单声道
[C,L]=wavedec(x,2,'db4');

③一维离散小波逆变换

xsync=waverec(C,L,'db4');
cA2=appcoef(C,L,'db4',2);
cD2=detcoef(C,L,2);
cD1=detcoef(C,L,1);

④观察结果

figure;
subplot(2,3,1);plot(x);title('original audio');
subplot(2,3,2);plot(xsync);title('synthesize audio');
subplot(2,3,4);plot(cA2);title('app coef. of audio');
subplot(2,3,5);plot(cD2);title('det coef. of auio');
subplot(2,3,6);plot(cD1);title('det coef. of auio');

结果展示：

2.分析合成图像文件

①读取图像文件数据

[fn,pn]=uigetfile('*.bmp','请选择图像文件');
[x,map]=imread(strcat(pn,fn),'bmp');
I=rgb2gray(x);

②二维离散小波变换

sx=size(I);
[cA1,cH1,cV1,cD1]=dwt2(I,'bior3.7');

③二维离散小波逆变换

xsync=uint8(idwt2(cA1,cH1,cV1,cD1,'bior3.7',sx));
A1=uint8(idwt2(cA1,[],[],[],'bior3.7',sx));
H1=uint8(idwt2([],cH1,[],[],'bior3.7',sx));
V1=uint8(idwt2([],[],cV1,[],'bior3.7',sx));
D1=uint8(idwt2([],[],[],cD1,'bior3.7',sx));

④观察结果

figure;
subplot(2,3,1);imshow(x);title('original image');
subplot(2,3,2);imshow(uint8(abs(xsync)));title('synthesize image');
subplot(2,3,3);mesh(A1);title('app coef. of image');
subplot(2,3,4);mesh(H1);title('hor coef. of image');
subplot(2,3,5);mesh(V1);title('ver coef. of image');
subplot(2,3,6);mesh(D1);title('dia coef. of image');
 

结果展示：

总结

拓展

图像傅立叶变换的物理意义

图像的频率是表征图像中灰度变化剧烈程度的指标，是灰度在平面空间上的梯度。如：大面积的沙漠在图像中是一片灰度变化缓慢的区域，对应的频率值很低；而对于地表属性变换剧烈的边缘区域在图像中是一片灰度变化剧烈的区域，对应的频率值较高。傅立叶变换在实际中有非常明显的物理意义，设f是一个能量有限的模拟信号，则其傅立叶变换就表示f的谱。从纯粹的数学意义上看，傅立叶变换是将一个函数转换为一系列周期函数来处理的。从物理效果看，傅立叶变换是将图像从空间域转换到频率域，其逆变换是将图像从频率域转换到空间域。换句话说，傅立叶变换的物理意义是将图像的灰度分布函数变换为图像的频率分布函数，傅立叶逆变换是将图像的频率分布函数变换为灰度分布函数。

傅立叶变换以前，图像（未压缩的位图）是由对在连续空间（现实空间）上的采样得到一系列点的集合，我们习惯用一个二维矩阵表示空间上各点，则图像可由z=f(x,y)来表示。由于空间是三维的，图像是二维的，因此空间中物体在另一个维度上的关系就由梯度来表示，这样我们可以通过观察图像得知物体在三维空间中的对应关系。为什么要提梯度？因为实际上对图像进行二维傅立叶变换得到频谱图，就是图像梯度的分布图，当然频谱图上的各点与图像上各点并不存在一一对应的关系，即使在不移频的情况下也是没有。傅立叶频谱图上我们看到的明暗不一的亮点，实际上图像上某一点与邻域点差异的强弱，即梯度的大小，也即该点的频率的大小（可以这么理解，图像中的低频部分指低梯度的点，高频部分相反）。一般来讲，梯度大则该点的亮度强，否则该点亮度弱。这样通过观察傅立叶变换后的频谱图，也叫功率图，我们首先就可以看出，图像的能量分布，如果频谱图中暗的点数更多，那么实际图像是比较柔和的（因为各点与邻域差异都不大，梯度相对较小），反之，如果频谱图中亮的点数多，那么实际图像一定是尖锐的，边界分明且边界两边像素差异较大的。对频谱移频到原点以后，可以看出图像的频率分布是以原点为圆心，对称分布的。将频谱移频到圆心除了可以清晰地看出图像频率分布以外，还有一个好处，它可以分离出有周期性规律的干扰信号，比如正弦干扰，一副带有正弦干扰，移频到原点的频谱图上可以看出除了中心以外还存在以某一点为中心，对称分布的亮点集合，这个集合就是干扰噪音产生的，这时可以很直观的通过在该位置放置带阻滤波器消除干扰。

音频信号处理的目的：时域–>频域，便于进一步处理。