matlab振幅谱函数,⽤MATLAB求周期信号的振幅谱和相位谱
函数运⾏出错（两个）...
第⼀个程序：
function[A_sym,B_sym]=gcpfbs
% T      ％输⼊信号周期
%        ％函数输⼊输出都是数值量
% Nf    ％谐波最⼤次数
% Nn    ％输出数据的准确位数
% An    ％cos项系数，第1项是直流项，以后依次为1，2，3…次谐波cos项系数
% Bn    ％sin项系数，第2，3，4…元素依次是1，2，3…次谐波sin项系数
syms y n t k  ％定义符号变量
T=5;
if nargin<2;Nf=input('pleas Input 所需展开的最⾼谐波次数');end
if nargin<3;Nn=32;end
% nargin表⽰输⼊参数个数
y=fb_symfun;
% 调⽤⽅波的符号表达式⼦函数，具体函数见后⾯的⼦函数
A0=2*int(y,t,0,T)/T;
% 根据定义求A0
A=int(2*y*cos(2*pi*n*t/T)/T,t,0,T);
% 求cos项系数的符号积分表达式
B=int(2*y*sin(2*pi*n*t/T)/T,t,0,T);
% 求sin项系数的符号积分表达式
An(1)=double(vpa(A0,Nn));
% 将A0的值以双精度32位形式传给An(1)
for k=1:Nf
% 循环语句，每执⾏⼀次，k加1直到等于Nf
An(k+1)=double(vpa(subs(A,n,k),Nn));
% 根据定义求An(n=2…Nf+1)
Bn(k+1)=double(vpa(subs(B,n,k),Nn));
end
% 根据定义求Bn(n=2,3,…Nf+1)
if nargout==0
% nargout表⽰输出参数，输出参数为零表⽰所有参数都已计算结束
S1=fliplr(An)
% 对An阵左右对称变换，An＝[A2...A(Nf+1)]，所以S1＝[A(Nf+1), (2)
S1(1,k+1)=An(1);
% S1阵列扩展为k＋1项，第k＋1项定为A1，S1＝[A(Nf+1), A(Nf), (1)
S2=fliplr(1/2*S1)
% S1取半然后反折后给S2，S2＝[A1/2,A2/2,…A(Nf+1)/2]
S3=fliplr(1/2*Bn);
% Bn取半然后反折后给S3，S3＝[B(Nf+1)/2,BNf/2,…B2/2]
S3(1,k+1)=0
% S3阵列扩展为k＋1项，第k＋1项定为0，S3＝[B(Nf+1)/2,BNf/2,…B2/2,0]
S4=fliplr(S3)
% 将S3反折后赋给S4， S4＝[0,B2/2,B3/2,…B(Nf)/2,B(N
S5=S2-1i*S4;
% 由三⾓系数得出指数形式的系数，S5=[ A1/2,…A(Nf+1)/2-i B(Nf+1)/2]
S6=fliplr(S5);
% 对得到的系数矩阵反折后传给S6，S6＝[A(Nf+1)/2-i B(Nf+1)/2,…A1/2]
N=Nf*2*pi/T;
% N为最⼤谐波⾓频率，等于所需展开谐波次数乘以，＝
k2=-N:2*pi/T:N;
% 形成－N到N的变量序列，变量之间的间隔为，[－N,－N＋,…N－,N]
S7=[S6,S5(2:end)];
% 形成－N:N的的傅⽴叶指数对称复系数，S7＝[A(Nf+1)/2-iB(Nf+1)/2…    A2/2- iB2/2, A1/2, A2/A(Nf+1)/2-i B(Nf+1)/2]
subplot(2,1,1)
% 创建⼦图1，位于第⼀⾏第⼀列
x=fb_pulse
% 调⽤⽅波脉冲函数画出⽅波波形图
T=5;t=-2*T:0.01:2*T;
% 给定⽅波周期和时间的最⼤最⼩值，X轴区间为[－10，10]，隔0.01取⼀个点
plot(t,x)
% 绘制⽅波波形图
title('周期⽅波脉冲')
% 给⼦图加标题
axis([-10,10,-1,1.2]);
% 修正坐标轴，依次为xmin，xmax，ymin，ymax
line([-10,10],[0,0])
% 添加图形中间的坐标轴
subplot(2,1,2)
% 创建⼦图2，位于第⼆⾏第⼀列
stem(k2,abs(S7));
% 绘制离散频谱图，括号中前者是X轴坐标序列，后者是取绝对值后的系数序列title('周期⽅波脉冲的双边频谱')
% 给⼦图加标题
axis([-30,30,0,0.6])
% 修正坐标轴
end
% 主函数结束
function y=fb_symfun
%⼦函数 由符号变量写成的⽅波脉冲的符号表达式
syms a a1
%定义符号变量
T=input('please input 信号的周期 T=');
M=input('周期与脉冲宽度之⽐ M=');
A=1;
tao=T/M;
a=tao/2;
y1=sym('Heaviside(t+a1)')*A;
y=y1-sym('Heaviside(t-a1)')*A;
y=subs(y,a1,a);
y=simple(y);
function x=fb_pulse
%⼦函数：时域⽅波脉冲函数
%t
%T
%duty：‘占空⽐’,
T=5;t=-2*T:0.01:2*T;duty=50;
%占空⽐为50％即为1:2
x=square(t,duty);
%调⽤系统中的⽅波函数
(第⼀个程序完)
运⾏报错：见图1
程序⼆：
function[A_sym,B_sym]=CTFSigsym
syms t n y
if nargin<3;Nf=input('pleas input 所需展开的最⾼谐波次数：Nf=');end T=input('pleas input 信号的周期 T=');
if nargin<5;Nn=32;end
y=time_fun_s(t);
A0=2*int(y,t,0,T)/T;
As=int(2*y*cos(2*pi*n*t/T)/T,t,0,T);
Bs=int(2*y*sin(2*pi*n*t/T)/T,t,0,T);matlab傅里叶变换的幅度谱和相位谱
A_sym(1)=double(vpa(A0,Nn));
for k=1:Nf
A_sym(k+1)=double(vpa(subs(As,n,k),Nn));
B_sym(k+1)=double(vpa(subs(Bs,n,k),Nn));
end
if nargout==0
S1=fliplr(A_sym)
S1(1,k+1)=A_sym(1)
S2=fliplr(1/2*S1)
S3=fliplr(1/2*B_sym)
S3(1,k+1)=0
S4=fliplr(S3)
S5=S2-1i*S4;
N=Nf*2*pi/T;
k2=0:2*pi/T:N;
x=time_fun_e(t)
subplot 212
stem(k2,abs(S5));
title('连续时间函数周期矩形脉冲的单边幅度谱')
axis([0,80,0,0.12])
line([0,80],[0,0])
line([0,0],[0,0.12])
end
function x=time_fun_e(t)
T=input('pleas input 信号的周期 T='); M=input('周期与脉冲宽度之⽐ M=');
t=-2*T:0.01:2*T;tao=T/M;
x=rectpuls(t,tao);
subplot 211
plot(t,x)
hold on
x=rectpuls(t-T,tao);
plot(t,x)
hold on
x=rectpuls(t+T,tao);
plot(t,x)
title('周期为T，脉宽tao=T/M 的矩形脉冲') axis([-10-T,10+T,0,1.2])
function y=time_fun_s(t)
syms a al
T=input('pleas input 信号的周期 T='); M=input('周期与脉冲宽度之⽐ M=');
A=1;tao=T/M;a=tao/2;
y1=sym('Heaviside(t+al)')*A;
y=y1-sym('Heaviside(t-al)')*A;
y=subs(y,al,a);
y=simple(y);
报错见图⼆
A.jpg
(18.99 KB, 下载次数: 0)
2014-12-1 15:40 上传
第⼀个程序运⾏报错

版权声明：本站内容均来自互联网，仅供演示用，请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198，我们将在24小时内删除。