%《MATLAB及其在理工課程中的應(yīng)用指南》(十一五規(guī)劃版)第五章例5-3-6程序
% 在[-pi,pi]區(qū)間內(nèi)任意函數(shù)的傅立葉級(jí)數(shù)展開(kāi)
% 西安電子科技大學(xué)出版社出版  陳懷琛編著 2007年7月
%
clear,clf,format compact,format short
x=linspace(-pi,pi,1001);dx=2*pi/1000;; 	% [-pi,pi]內(nèi)長(zhǎng)為1001的x數(shù)組步長(zhǎng)
f=input('輸入f=(長(zhǎng)度為1001點(diǎn)的數(shù)組) '); % 用戶輸入長(zhǎng)為1001的f數(shù)組
if isempty(f) error('未給出函數(shù)數(shù)組'), end
n=input('傅立葉系數(shù)的階數(shù)n= ')          % 用戶給出所需傅立葉系數(shù)的階數(shù)
if isempty(n) n=9, end                  % 若不給出階數(shù)，按9階計(jì)算
a0=trapz(f)/pi*dx                       % 計(jì)算傅立葉系數(shù)a0
disp('      k        a(k)       b(k)')  % 給出表頭
for k=1:n 
    a(k)=trapz(f.*cos(k*x))/pi*dx;      % 計(jì)算傅立葉系數(shù)a(k)
    b(k)=trapz(f.*sin(k*x))/pi*dx;      % 計(jì)算傅立葉系數(shù)b(k)
    disp([k,a(k),b(k)])                 % 顯示各階系數(shù)
end
pause,disp('求這些傅立葉系數(shù)構(gòu)成的級(jí)數(shù)的波形')
pause,f1=a0/2*ones(size(x));            % 以a0為基礎(chǔ)，構(gòu)造傅立葉級(jí)數(shù)
for k=1:n
    f1=f1+a(k)*cos(k*x)+b(k)*sin(k*x);  % 累加各項(xiàng)傅立葉級(jí)數(shù) 
end
subplot(1,2,2),plot(x,f1,'linewidth',2),              % 畫出傅立葉級(jí)數(shù)的波形
v=axis;grid on
subplot(1,2,1),plot(x,f,'linewidth',2)                % 畫出原輸入函數(shù)的波形
axis(v),grid on,shg
set(gcf,'color','w')