Fourier级数变换小程序_Matlab

Fourier级数变换

这里讨论周期函数 f(t) 在对称区间(-L，L)上进行Fourier级数展开.
(1)Fourier级数公式

其中

（2）Euler公式

(3) 利用(2)式整理可得ck如下

当K=0时

K不等于0时

将ak,bk带回Fourier公式，从而可以得到Fourier级数的展开式如下

由余弦函数是偶函数，正弦函数是奇函数的特点，可进一步将上式进行改写成如下的式子.

（4）选取如下的周期函数进行测试

（5）代码部分（包含了对Parseval恒等式的验证）
Parseval恒等式（这里我只考虑了一个函数做内积！也可以改成两个不同函数做内积的情况，(f,g）.

mian.m

% Fourier 级数展开
% 讨论周期函数的傅里叶展开式
clc;clear;close;
syms t k;
f1=2*t+1;
f2=-2*t+1;
L=1/2;       %讨论[-1/2,1/2]
m=10;        %展开的级数项数

% 计算 f的内积
I=(int(f1*f1,t,-L,0)+int(f2*f2,t,0,L))/(2*L);

% 进行傅里叶转化
[a0,ak,bk]=myfourier(f1,f2,L);  % 得到 ak，bk的通项表达式
[ck,s]=sfourier(a0,ak,bk,m,L);  % 级数求和
c0=a0;

%代入数值
dt=0.01;
t=-L:dt:L;
ss=vpa(subs(s)); % vpa将结果转发为小数
ff=ff(t,L);

%验证 Parseval恒等式
a=real(ck);
b=imag(ck);
ck_=a-1i*b;
%sum_ck=sum(abs(ck).^2 )+(c0^2)/2;
sum_ck=2*sum(ck.*ck_)+(c0^2);
eps=vpa(abs(I-sum_ck)) %f的内积与级数系数乘积加和的误差
ck=[c0,ck];

%画图比较真实值与级数展开的结果
subplot(211)
plot(t,ff,'r*',t,ss,'bo');
legend('周期函数图像','级数展开图像')
xlabel('t');
ylabel('f');
subplot(212);
K=0:1:m;  % k 的取值范围
stem(K,abs(ck),'filled');
xlabel('k');
ylabel('频率值 Ck');
legend('频谱图');

ff函数代码如下 ff.m

function y=ff(t,L)
y=(2*t+1).*( t>=-L & t<0 )+(-2*t+1).*(t>=0 & t<L);

计算Fourier计算ak，bk. ( ck=ak-i bk ).ck是Fourier系数
myfourier.m

%计算系数积分表达式,区间为（-L,L)
function [a0,ak,bk]=myfourier(f1,f2,L)
syms t k;
a0=(int(f1,t,-L,0)+int(f2,t,0,L) )/(2*L);
ak=( int(f1*cos(k*pi*t/L),t,-L,0)+int(f2*cos(k*pi*t/L),t,0,L))/(2*L);
bk=( int(f1*sin(k*pi*t/L),t,-L,0)+int(f2*sin(k*pi*t/L),t,0,L))/(2*L);
return

级数求和代码 sfourier.m

% 求信号f在(-L,L)上的傅里叶级数前m项
function [ck,s]=sfourier(a0,ak,bk,m,L)
syms t k;
sk=ak*cos(k*pi*t/L)+bk*sin(k*pi*t/L); % k 取正值
ssk=symsum(sk,k,1,m);
s=a0+2*ssk;  % a0是数值,ssk是含t的表达式
k=1:m;
ak=subs(ak); %k遍历，得到数值ak，不再是符号表达式
bk=subs(bk);
ck=ak-1i*bk;         % 数值
return

%实验结果图

我刚学完Fourier级数这里，看了很多博客，很多文章的Fourier级数的表达式和我写的这个不太一样，可能差一个常数！相应地，Fourier系数表达式做相应的变化就行（和我给的表达式差一个常数）。但其本质和推导过程是相似的！第一次写，不准确的地方欢迎指出来！哈哈哈哈~~~~~~！（接下来有时间就把Fourier变换推导写写！）