周期方波信号分解

pdf里面的代码：

1
2
3
4
5
6
7
8
9
figure
for i = 1:5
    si=tf([1000*i 0],[1 1000*i ((2*i*pi*(10^4))^2)]);
    bode(si);
    hold on
    str{i} = [ num2str(10*i) 'KHz'];
end
legend(str);

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
clc,clear
%设置参数
R=10;
L=10*10^(-3);
C=25.33*10^(-9);
%固定的输入频率、角频率
f=10^4;
w0=2*pi*f;
%设置时间
t = 0:1/(10^4*f):1/f;
%设置分割数
k=1:2:1001;
%sin前的系数
Ak=4./(pi*k);
%sin·分量
Bk=sin(w0*k'*t);
%输入U
U=Ak*Bk;
plot(t,U);
hold on
%滤波器处理函数
Hjw=1./(1+1j*(L*w0.*k./R - 1./(C*R*w0.*k)));
%输出UR,这里加角度是因为不是用atan分母的角，而是直接用angle得到角
Hs=abs(Hjw)'.*sin(w0*k'*t+angle(Hjw)');
UR=Ak*Hs;
plot(t, UR);
legend('输入方波波形','10KHz滤波器')
axis([0 1/f -inf inf])