fft的c语言和matlab对比_应用MatLab对语音信号进行处理

大三时数字信号处理做过一个简单的语音信号处理。复制过来。

课题内容

录制一段语音信号，对录制的信号进行采样；画出采样后语音信号的时域波形和频谱图；给定滤波器的性能指标，采用窗函数法或双线性变换设计滤波器，并画出 滤波器的频率响应；然后用自己设计的滤波器对采集的信号进行滤波，画出滤波后信号的时域波形和频谱，并对滤波前后的信号进行对比，分析信号的变化；最后，实现对原始信号进行处理，包括：加噪声、快放、慢放、变声等。

具体实现
1．语音信号的采集
　　利用手机录音机，录制一段话音（我们所录的是---“现在所录的是男声”）。然后在Matlab软件平台下，利用函数audioread对语音信号进行采样。通过audioread函数的使用，我们很快理解了采样频率、采样位数等概念。
2．原始语音信号的频谱分析
　　我们首先画出语音信号的时域波形；然后对语音号进行快速傅里叶变换，得到信号的频谱特性，从而加深对频谱特性的理解。
程序如下：
[x,fs]=audioread('nansheng.wav'); %打开语音信号
sound(x,fs); %播放语音信号
N=length(x); %长度
n=0:N-1;
w=2*n*pi/N;
y1=fft(x); %对原始信号做FFT变换
subplot(2,1,1);
plot(n,x) %做原始语音信号的时域波形图
title('原始语音信号时域图');
xlabel('时间t');
ylabel('幅值');
subplot(2,1,2); %做原始语音信号的频谱图
plot(w/pi,abs(y1));
title('原始语音信号频谱')
xlabel('频率Hz');
ylabel('幅度');
波形如下：

3.加噪语音信号并对其FFT频谱分析

程序如下：

[x,fs]=audioread('mei.wav');

n=length(x);

x_p=fft(x,n);

f=fs*(0:n/2-1)/n;

figure(1)

subplot(2,1,1);

plot(x);

title('原始语音信号采样后的时域波形');

xlabel('时间轴')

ylabel('幅值A')

subplot(2,1,2);

plot(f,abs(x_p(1:n/2)));

title('原始语音信号采样后的频谱图');

xlabel('频率Hz');

ylabel('频率幅值');

noise=0.005*randn(1,n);

x_z=x+noise';

sound(x_z,fs)

n=length(x);

x_zp=fft(x_z,n);

f=fs*(0:n/2-1)/n;

figure(2)

subplot(2,1,1);

plot(x_z);

title('加噪语音信号时域波形');

xlabel('时间轴')

ylabel('幅值A')

subplot(2,1,2);

plot(f,abs(x_zp(1:n/2)));

title('加噪语音信号频谱图');

xlabel('频率Hz');

ylabel('频率幅值');

波形如下：

4. 滤波器设计

我们用双线性变换法设计了一个数字低通滤波器对原始语音信号进行了滤波处理。

程序如下：

fp=800;fs=1300;rs=35;rp=0.5;Fs=44100;

wp=2*Fs*tan(2*pi*fp/(2*Fs));

ws=2*Fs*tan(2*pi*fs/(2*Fs));

[n,wn]=buttord(wp,ws,rp,rs,'s');

[b,a]=butter(n,wn,'s');

[num,den]=bilinear(b,a,Fs);

[h,w]=freqz(num,den,512,Fs);

figure(1)

plot(w,abs(h));

xlabel('频率/Hz');ylabel('幅值');

title('巴特沃斯低通滤波器幅度特性');

axis([0,5000,0,1.2]); grid on;

[s1,Fs]=audioread('mei.wav');

x1=s1(:,1);

sound(x1,Fs);

N1=length(x1);

Y1=fft(x1,N1);

f1=Fs*(0:N1-1)/N1; t1=(0:N1-1)/Fs;

figure(2)

plot(f1,abs(Y1))

xlabel('频率/Hz');ylabel('幅度');

title('原始信号频谱');

grid on;axis([0 50000 0 200])

y=filter(num,den,x1);

sound(y,Fs);

N2=length(y);

Y2=fft(y,N2);

f2=Fs*(0:N2-1)/N2;

t2=(0:N2-1)/Fs;

figure(3)

plot(f2,abs(Y2))

xlabel('频率/Hz');ylabel('幅度');

title('过滤后信号的频谱'); grid on;

axis([0 50000 0 200])

波形如下：

5.快放：

改变采样频率即可实现。

[x,fs]=audioread('nansheng.wav')

w=2；

M=w*fs; %采样频率加倍

6.慢放：

[x,fs]=audioread('nansheng.wav')

w=0.5；

M=w*fs; %采样频率减半

7.男女变声处理：

[y,fs]=audioread('nansheng.wav'); %读取声音文件

x1=y(:,1); %读入的y矩阵有两列,取第1列

sound(voice(x1,1.5),fs);

N=length(voice(x1,1.5)); %长度

n=0:N-1;

w=2*n*pi/N;

y1=fft(voice(x1,1.5)); %对原始信号做FFT变换

subplot(2,1,1);

plot(n,voice(x1,1.5)) %做原始语音信号的时域波形图

title('变声语音信号时域图');

xlabel('时间t');

ylabel('幅值');

subplot(2,1,2); %做原始语音信号的频谱图

plot(w/pi,abs(y1));

title('变声语音信号频谱')

xlabel('频率');

ylabel('幅度');

※所调用的Voice函数代码：

function Y=voice(x,f) %更改采样率使基频改变 f>1降低;f<1升高

f=round(f*1000);

d=resample(x,f,1000); %时长整合使语音文件恢复原来时长

W=400;

Wov=W/2;

Kmax=W*2;

Wsim=Wov;

xdecim=8;

kdecim=2; X=d';

F=f/1000;

Ss =W-Wov;

xpts = size(X,2);

ypts = round(xpts / F);

Y = zeros(1, ypts);

xfwin = (1:Wov)/(Wov+1);

ovix = (1-Wov):0; newix = 1:(W-Wov);

simix = (1:xdecim:Wsim) - Wsim;

padX = [zeros(1, Wsim), X, zeros(1,Kmax+W-Wov)];

Y(1:Wsim) = X(1:Wsim); lastxpos = 0; km = 0;

for ypos = Wsim:Ss:(ypts-W)

xpos = round(F * ypos);

kmpred = km + (xpos - lastxpos);

lastxpos = xpos;

if (kmpred <= Kmax)

km = kmpred;

else

ysim = Y(ypos + simix);

rxy = zeros(1, Kmax+1);

rxx = zeros(1, Kmax+1);

Kmin = 0;

for k = Kmin:kdecim:Kmax

xsim = padX(Wsim + xpos + k + simix);

rxx(k+1) = norm(xsim);

rxy(k+1) = (ysim * xsim');

end

Rxy = (rxx ~= 0).*rxy./(rxx+(rxx==0));

km = min(find(Rxy == max(Rxy))-1);

end

xabs = xpos+km;

Y(ypos+ovix) = ((1-xfwin).*Y(ypos+ovix)) + (xfwin.*padX(Wsim+xabs+ovix));

Y(ypos+newix) = padX(Wsim+xabs+newix);

end

End

波形如下：

f=0.6(主要实现男声变女声）

f=1.4（实现女声变男声）