求数字滤波的FFT算法实现的MATLAB代码
数字滤波的FFT算法实现的MATLAB代码,谢谢各位了.[ 本帖最后由 eight 于 2007-5-29 12:19 编辑 ] 我也想要的 啊！！！ 一些Matalab相关程序！
现发一些自己做毕业设计做的和一些收集到的相关Matalab程序，共大家参考！
Matlab中的fft分析含噪声的信号，其中有直接FFT，时域同步平均后FFT，频域同步平均后FFT
%%This procedure is to investigate the FFT of sine signal, sine signal with
%%noise.The sine signal with noise is treated by direct FFT,
%%FFT after piece-wise averaging in time domain, and FFT after
%%piece-wise averagingfrequency domain
%%
clear all
pack
w=10;
t=0:(2*pi/w/500):(300*2*pi/w);
x=sin(w*t);
n=length(t);
ss=888;
randn('seed',ss);
z=randn(1,n);
z1=z-mean(z);
z2=z1/max(abs(z1));
z2=z2*500;
figure(1)
plot(t,x)% The sine signal without noise
fre=(1/(2*pi/w/500))*(1:n)/n;
nn=fix(n/2);
y1=fft(x);% The FFT of sine signal without noise
y1=(abs(y1)).^2/n;
y2=y1(2:(nn+1));
fre1=2*pi*fre(1:nn);
figure(2)
plot(fre1,y2)
x=x+z2;
y1=fft(x);% The FFT of sine signal with noise
y1=(abs(y1)).^2/n;
y2=y1(2:(nn+1));
figure(3)
plot(fre1,y2) %The sine signal with noise
%% Averaging in time domain
M=10000;
xx=zeros(1,M);
m1=1;m2=M;
nnnn=fix(300*500/M);
for nnn=1:nnnn
xxx=x(m1:m2);
xx=xx+xxx;
m1=m1+M;
m2=m2+M;
end
fre=(1/(2*pi/w/500))*(1:M)/M;
nn=fix(M/2);
xx=xx/nnnn;
y1=fft(xx);
y1=(abs(y1)).^2/M;
y2=y1(2:(nn+1));
fre1=2*pi*fre(1:nn);
figure(4)
plot(fre1,y2)
% Next is averaging in frequency domain
% In fact, the averaging in time domain is identical with that in frequency
% domain
MM=10000;
xx=zeros(1,MM);
m1=1;m2=MM;
nnnn=fix(300*500/MM);
for nnn=1:nnnn
xxx=x(m1:m2);
xxx=fft(xxx);
%xxx=sqrt((abs(xxx)).^2);
xx=xx+xxx;
m1=m1+MM;
m2=m2+MM;
end
y1=xx/nnnn;
y1=(abs(y1)).^2/MM;
fre=(1/(2*pi/w/500))*(1:MM)/MM;
nn=fix(MM/2);
y2=y1(2:(nn+1));
fre1=2*pi*fre(1:nn);
figure(5)
plot(fre1,y2)
对余弦信号的FFT分析程序:
%f—余弦信号的频率
% M—基2 FFT幂次数 N=2^M为采样点数，这样取值是为了便于作基2的FFT分析
%2. 采样频率Fs
%*******************************************************************%
function samples(f,Fs,M)
N=2^M; % fft点数=取样总点数
Ts=1/Fs; % 取样时间间隔
T=N*Ts; % 取样总时间=取样总点数*取样时间间隔
n=0:N-1;
t=n*Ts;
Xn=cos(2*f*pi*t);
subplot(2,1,1);
stem(t,Xn);
axis();
xlabel('t -->');
ylabel('Xn');
Xk=abs(fft(Xn,N));
subplot(2,1,2);
stem(n,Xk);
axis();
xlabel('frequency -->');
ylabel('!Xk!'); 可以参考http://forum.vibunion.com/forum/thread-44622-1-1.html 很好的学习资料，楼主有做傅里叶逆变换然后再FFT变换的例子没有啊？传上来学习下
页:
[1]