吉布斯现象的matlab仿真,吉布斯现象matlab代码

信号与系统实验

一：实验题目

运用matlab验证吉布斯现象 二：实验原理

对于具有不连续点(跳变点)的波形，所取级数项数越多，近似波形的方均误差虽可减小，但在跳变点处的峰起(上冲)值不能减小，此峰起随项数增多向跳变点靠近。(详见课本279页) 三：实验内容

1.计算N取不同值时的方均误差EN的值 2.用matlab画出N取不同值时

S(t)=a0[ancos(n1t)bnsin(n1t)]的图。

n1N

四：实验过程

我们取课本99页的函数f(t)作为实验的目标函数

 -E/2 (T/2tT/4)



f(t)=E/2(T/4tT/4)

-E/2(T/4tT/2)

matlab有专门产生矩形波的函数square(t)，为了显示对称波形，

将它移动四分之一个周期，即square(t+T/2)，然后再乘上我们此题中的振幅,因此得到

f_t=(E/2)*square(t+T/2);

E2

已知周期函数f(t)的傅里叶级数为

f(t)=a0[ancos(n1t)bnsin(n1t)]

n1



因为此题的函数既是偶函数，又是奇谐函数。因此在它的傅里叶级数中只可能含有奇次谐波的余弦项。所以于是

f(t)

2E

an

2Ensin()n2

11

[cos(1t)cos(31t)cos(51t)]35

若取傅里叶级数的前(2N1)项来逼近周期函数f(t)，则有限项傅里叶级数为

S(t)=a0[ancos(n1t)bnsin(n1t)]

n1N

2E2N1sin(nt)

SNTnn1进一步简化为

为了得到对称波形，将它移动四分之一个周期，即

2E

T

2N1sin(nt

SN



n1

n

nT

)

这样用S(t)逼近f(t)引起的误差函数为N(t)f(t)S(t) 方均误差等于

2

ENN(t)

SN就是代码中的S

1t0T12

N(t)dtt0T1

在程序中我们取t=-T:0.0001:T;将t等分成多份，所以不能直接用积分函数int(即不连续)，所以依据定义的方法来求积分(极限法)

lim((t1)t)因为均分所以积分变为t

2

N

t1TT

t0

t1T



T

2N

2

其中N为

t1T

T

2

所有分量t1对应的N值的总和

求EN时应除以总周期即2T数，matlab中为

t

即是将区间分为多少等份的倒2T

1

，其中length(t)为将区间分成的等份数。

length(t)

综上 ，我们得EN

t1T



T

2N

length(t)



sum((f_t-S).^2)

(后面是用在

length(t)

matlab 中的计算方法)

第一问：

运行matlab分别令N1;N3;N5 得n2N1，n为En的计数，出E10.7577;E20.3975;E30.2678…………如下表1

表1 方均误差

第二问：

3

2

相应的函数值f(t)和S(t)

1

-1

-2

-3-4

-3-2-1

0时间t

1234

有限级数 S2

3

2

相应的函数值f(t)和S(t)

1

-1

-2

-3-4

-3-2-1

0时间t

1234

3

2

相应的函数值f(t)和S(t)

1

-1

-2

-3-4

-3-2-1

0时间t

1234

有限级数 S20

相应的函数值f(t)和S(t)

-4

-3-2-1

0时间t

1234

放大跳变点附近的图像得到下图

有限级数 S20

2.352.32.252.22.152.12.051.3

1.35

1.4

1.45

1.51.55时间t

1.6

1.65

1.7

1.75

相应的函数值f(t)和S(t)

五：实验结论与分析

由第一问的表1和第二问的图像可知，随着N的增大，即随着所取级数项数的增多，进似方均误差减小，且峰起随着项数的增多向跳变点靠近，并且从有限级数S20跳变点附近的放大图可知峰起值趋于相同(峰起值趋近于跳变值的9%，这里不再予以验证)。 代码如下：

clear;

E=4; T=pi;

t=-T:0.0001:T;%周期 N=input('N'); i=1;

for i=1:N%设出i，让它变化，从而绘出S1，S2，S3……SN a=0;

for n=1:2:i

b=(E/2)*sin(n*t+n*T/2)/n;%移动四分之一周期以显示对称波形 a=a+b; end

S(i,:)=4*a/T;

f_t=(E/2)*square(t+T/2); plot(t,S,t,f_t) xlabel('时间t')

ylabel('相应的函数值f(t)和S(t)') title('有限级数 S20') end

E_N=sum((f_t-S(N,:)).^2)/length(t)