fzzyquad.m

精通matlab

function [t,y,S,A,B]=fzzyquad(a,T,Nf,K)
% 用数值积分求Fourier级数展开系数A,B和近似波形S 。如调用本指令时，没有任何
% 输出宗量，那么将自动随积分进程动态地画出原波形和各阶近似波形。
%	a		是被展开函数的时间区间的左端。
%	T		是被展开函数的周期。
%	Nf		是所需展开的最高谐波阶次
%	K		绘波形图所采用的数据长度
%	t,y		是被展开的时间函数数据。
%	S		是一个矩阵。它每行的长度与t,y相同。
%				它的第（n+1)行是所有n次及更低次谐波迭加生成的近似波形。
%	A,B		是展开系数向量。A(n+1),B(n+1)分别存放cos，sin项n次谐波的展开系数。
%				A(1)是直流项
if nargin<4;K=100;end
if (nargin<3|isempty(Nf));Nf=15;end
k=1:K;
t=a+k*T/length(k);
y=time_fun(t,T);
S=zeros(Nf+1,length(k));
a0=mean(y); 
S(1,:)=a0;
A=zeros(1,Nf+1);B=zeros(1,Nf+1);
A(1)=a0;
tolzzy=[0.001,0.018];
for n = 1:Nf
   A(n+1) = quad8('cos_y',a,a+T,tolzzy,[],n,T)/T*2;
   B(n+1) = quad8('sin_y',a,a+T,tolzzy,[],n,T)/T*2;
   S(n+1,:)=S(n,:)+A(n+1)*cos(2*n*pi*t/T)+B(n+1)*sin(2*n*pi*t/T);
   if nargout==0
      plot(t,y,'b');hold on
      plot(t,S(n+1,:),':r');
      aa=axis;dx=aa(2)-aa(1);dy=aa(4)-aa(3);
      text(aa(1)+0.1*dx,aa(3)+.75*dy,'n = ')
      text(aa(1)+0.2*dx,aa(3)+.75*dy,num2str(n));
      pause(1);hold off
   end
end 
%------------------------------
function y=time_fun(t,T)
% 这是需要用户自己编写的待展开的时间函数 y(t+k*T)=y(t)
% 编写程序时，一定要注意使表达式适于数组运算
%	t		是时间数组
%	T		是周期
y=zeros(size(t));ii=find(t>=0.5 & t<=1.5);
y(ii)=ones(size(ii)).*(t(ii)-0.5);y(y==1.0)=0.5;