10-1.m

《DSP嵌入式系统开发典型案例》此书的源代码

clear all;										% 清所有的变量
clf;											% 清屏
N=1224;										% 设置数据点数
f1=200;										% 设置信号1频率
f2=700;										% 设置信号2频率
f3=1100;										% 设置信号3频率
dt=1/1024;										% 设置采样时间
sample=1/dt;									% 设置采样频率
randn('state',sum(100*clock));						% 产生随机噪声
dis=2*rand(1,N);
dis=dis-mean(dis);								% 将随机噪声的均值设置成0
dis=2*dis;										% 扩大随机噪声的幅度
for k=1:N;										% 产生待检测的信号
    w(k)=1*sin(2*pi*f1*k*dt)+1*sin(2*pi*f2*k*dt)+1*sin(2*pi*f3*k*dt)+dis(k);
end

figure(1);										% 计算待检测信号的均匀采样频谱
fw=fft(w);
s=abs(fw).^2;
f=(1:N-1)*sample/N;
plot(f,s(2:N));
xlabel('Frequency /Hz');							% 横坐标说明
ylabel('Amplitude');								% 纵坐标说明
title('DFT of uniformly sampling signal');				% 标题说明

t=tnonunif(N,1,1);								% 调用非均匀采样时间函数
for k=1:N;										% 产生非均匀采样数据
    w2(k)=1*sin(2*pi*f1*t(k)*dt)+1*sin(2*pi*f2*t(k)*dt)+1*sin(2*pi*f3*t(k)*dt)+dis(k);
end

for k=1:N;										% 计算非均匀采样数据的频谱
    s2(k)=0;
    for n=1:N;									% 计算非均匀采样数据的傅立叶变换
        s2(k)=s2(k)+w2(n)*(cos(2*pi*t(n)*(k-1)/N)-j*sin(2*pi*t(n)*(k-1)/N));
    end
    s2(k)=abs(s2(k)).^2;
end

figure(2);										% 画出非均匀采样数据的频谱
f=(1:N-1)*sample/N;
plot(f,s2(2:N));
xlabel('Frequency /Hz');
ylabel('Amplitude');
title('DFT of additive random sampling signal');

f=(1:N-1)*sample/N;
plot(f,s2(2:N));
xlabel('Frequency /Hz');
ylabel('Amplitude');
title('DFT of additive random sampling signal');