第一个程序:
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), …A2]
S1(1