关于DSP的实验程序源代码

% dm09201 % 绘制周期信号的频谱特性 function CTFS_RP % 以周期矩形脉冲信号为例，计算其频谱特性 % Nf：级数分解的谐波次数，由键盘输入 % Nn：输出数据的准确位数 % a0：直流项系数 % an：第1,2,3,...次谐波余弦项展开系数 % bn：第1,2,3,...次谐波正弦项展开系数 % tao：周期矩形脉冲信号脉宽，由键盘输入 % T：周期矩形脉冲信号周期，由键盘输入 display('Please input the value of T, tao and Nf'); %命令窗口提示用户输入参数 T = input('T = '); tao = input('tao = '); Nf = input('Nf = '); syms t n k x %定义符号变量 Nn = 32; %输出数据的位数为32位 an = zeros(Nf+1,1); %分配an系数数组 bn = zeros(Nf+1,1); %分配bn系数数组 phase = zeros(Nf+1,1); %分配相位数组 s1 = strcat('Heaviside(t+',num2str(tao/2),')'); %构造u(t+tao/2) s2 = strcat('Heaviside(t-',num2str(tao/2),')'); %构造u(t-tao/2) x = sym(s1) - sym(s2); %构造一个周期的脉冲信号u(t+tao/2)-u(t-tao/2) A0 =2*int(x,t,-T/2,T/2)/T; %求出直流项a0 As=2*int(x*cos(2*pi*n*t/T),t,-T/2,T/2)/T; %求出余弦项系数an Bs=2*int(x*sin(2*pi*n*t/T),t,-T/2,T/2)/T; %求出正弦项系数bn an(1) = double(vpa(A0,Nn)); %获取参数组A0所对应的ASCII码数值数组 for k=1:Nf an(k+1)=double(vpa(subs(As,n,k),Nn)); %获取参数组As所对应的ASCII码数值数组 bn(k+1)=double(vpa(subs(Bs,n,k),Nn)); %获取参数组Bs所对应的ASCII码数值数组 end cn = sqrt(an.*an+bn.*bn); %计算幅度谱 for i = 1:Nf %由于bn=0，相位由an的符号决定为0或pi if an(i) >= 0 phase(i) = 0; else phase(i) = pi; end end t = -T*5:0.001:T*5; d = -T*5:T:T*5; xx = pulstran(t,d,'rectpuls',tao); %用pulstran函数生成矩形脉冲信号 subplot(311); %将显示窗口分为3个子窗口，并指向第1个子窗口 plot(t,xx); %绘制周期矩形脉冲信号 axis([-T*5 T*5 0 1.1]); %指定坐标系范围 %title('周期矩形脉冲信号','Fontsize',8); %标注标题 s1 = strcat('周期矩形脉冲信号 T=',num2str(T),' Tao=',num2str(tao),'t'); xlabel(s1,'Fontsize',8); %x轴标签 subplot(312); %指向第2个子窗口 k = 0:Nf; stem(k,cn); %绘制幅度谱 hold on; plot(k,cn); %绘制幅度谱包络线 xlabel('幅度谱 \omega','Fontsize',8); subplot(313); %指向第3个子窗口 stem(k,phase); %绘制相位谱 xlabel('相位谱 \omega','Fontsize',8); % End