根据傅里叶级数用matlab实现周期信号分解

• 周期方波信号
  $x(t) = \sum_{n=1}^{\infty}b_nsin(n \omega t)$x(t)=n=1∑∞​bn​sin(nωt)
  $b_n = \begin{cases} 0, n=2,4,6,...\\ \frac{4}{n\pi}, n=1,3,5,... \end{cases}$bn​={0,n=2,4,6,...nπ4​,n=1,3,5,...​

% 周期方波信号的傅里叶级数
t = -1:0.001:1;
w0 = 2 * pi;
y = square(2 * pi * t,50);
%周期方波信号
plot(t,y), grid on;
axis([-1 1 -1.5 1.5]);
n_max = [1 3 5 7 31];
N = length(n_max);
for k = 1:N
    n = 1:2:n_max(k);
    b = 4./(pi*n);
    x = b * sin(w0*n'*t);
    figure;
    plot(t,y);
    hold on;
    plot(t,x);
    hold off;
    axis([-1 1 -1.5 1.5]),grid on;
    title(['最大谐波数=',num2str(n_max(k))]);
end

• 周期锯齿脉冲信号
  $b_n = (-1)^{n+1} \frac{1}{n\pi}$bn​=(−1)n+1nπ1​
  $x(t) = \frac{1}{\pi}\sum_{n=1}^{\infty}(-1)^{n+1} \frac{1}{n}sin(n \omega t)$x(t)=π1​n=1∑∞​(−1)n+1n1​sin(nωt)

%周期锯齿脉冲信号的傅里叶级数
%取谐波次数为35
t = -3:0.001:3;
w = pi ;%w的确定
x = sawtooth(pi*(t+1));
figure;
plot(t,x);grid on;
axis([-3 3 -1.2 1.2])
Nf = 35;
bn(1) = 0;
for i = 1:Nf
    bn(i+1) = (-1)^(i+1)*1/(i*pi);
    cn(i+1) = abs(bn(i+1));
end
n_max = [1 2 3 7 21 35];
N = length(n_max);
for k = 1:N
    n = 1:n_max(k);
    b = (-1).^(n+1).*1./(n.*pi); 
    y = 2*b * sin(n'*w*t);%系数的设置，重新计算
    figure
    plot(t,x);hold on
    axis([-3 3 -1.2 1.2])
    plot(t,y);hold off
    title(['最大谐波数=',num2str(n_max(k))]);
end
figure;
k = 0:Nf;
stem(k,cn);
hold on;
plot(k,cn);
title('幅度频谱');

以上内容有参考信号的matlab辅导书，频率周期系数的关系还有一点令人迷惑的地方，下午约了同学一起打羽毛球，三角脉冲信号以后再说吧，大概率会鸽