傅里叶级数理论详讲&实例应用

巴山学长

发布于 2020-08-06 01:44:22

3K00

代码可运行

文章被收录于专栏：巴山学长巴山学长

运行总次数：0

代码可运行

久久不更，更则必是大招，本期推文有理论有案例，有兴趣的读者仔细阅读。

过冷水之前有和大家讲傅里叶级数，并给出以一个函数用傅里叶级数近似的案例。本期就进一步详讲傅里叶级数。傅里叶级数展开时基底函数取1,cosx、sinx,cos2x、sin2x.....cosnx、sinnx,傅里叶级数一般情况下表示为：

a0、an、bn是展开系数。假定一个周期为2π的函数f(x+2π)=f(x)现在计算其系数。这就需要一点灵活的数学思想来解决问题。我们对上式两边在[0,2π]范围内积分。可得：

我们继续采用两边积分的方法求系数an、bn，对泰勒级数两边同时乘以cos(mx),然后在[0,2π]范围积分，如下：

则当m=n时：同理推bn：

许多实际问题中，函数f(x)是一个定义在有限区间0<x<X上的任意函数。对这样的函数如何展开傅里叶级数？大师告诉我们依旧可以用正弦函数、余弦函数展开：

现在求展开系数的表达式：

讲完傅里叶级数的理论，我们来看看实际傅里叶级数实际应用。过冷水在学习过程中遇到如下一幅图，需要对该函数进行积分求值：

数值图像求积分的问题过冷水的往期推文数值计算——MATLAB数值积分原理详讲、数值优化—三种复杂函数数值积分方法实例演示都有讲。过冷水本期的重点不是求数值积分的方法讨论，就想弄清楚如何很好的给出很好拟合该线段的函数表达式。这个函数图像算是比较复杂的，很少接触到所示的图像函数式，在此只能给出一个较为接近的拟合公式：

%原图像数据
x=[-0.0171386436557126,2.64214113473764,2.75334303194563,2.92014587775763,2.92134159708245,2.97813826501126,3.03453635983180,3.09133302776062,3.14852826879770,3.20592279638893,3.26311803742601,3.31991470535483,3.37691065983778,3.43350804121246,3.49030470914127,3.54490922497459,3.54371350564977,3.59811873492895,3.65232467765400,3.70633133382491,3.76053727654995,3.81514179238327,3.86974630821658,3.92435082404990,3.97935391299149,4.03415771537894,4.14416389326212,4.25436935769944,4.36457482213675,4.58578289722792,4.69598836166524,4.80738954542737,4.86318978058550,4.97419239123936,5.14019809083481,5.25080212838040,5.36160545248012,5.47280734968812,5.69441399788756,5.80581518164969,5.91721636541183,5.91781422507423,5.97401303334064,6.03001255505291,6.08581279021105,6.19581896809422,6.25082205703581,6.36122680802726,6.41642918352299,6.69343749377230,6.85944319336774,7.02584746607147,7.13625221706292,7.30185934355009,7.57866836724526,7.74527192650312,7.80127144821539,7.85707168337352,7.96727714781084,8.07748261224815,8.18768807668547,8.29829211423105,8.57530042448036,8.63110065963850,8.90810896988780,9.01931086709580,9.07531038880807,9.13130991052034,9.24251180772833,9.29831204288647,9.40851750732378,9.51872297176110,9.68433009824827,9.79533270890213,9.90593674644772,10.0719424460432,10.3489507562925,10.4601526535005,10.6819585882541,10.8475657147412,10.9583690388410,11.0695709360489,11.1805735467028,11.2907790111401,11.4013830486857,11.5677873213894,11.7341915940931,11.8997987205803,12.0660037067299,12.2326072659878,12.3436098766416,12.5100141493453,12.6214153331075,12.6774148548197,12.8432212678611,13.0088283943482,13.1750333804978,13.2302357559935,13.5626457282927,13.6736483389466,13.7286514278882,13.8392554654337,13.9496602164252,14.3374718507742,14.4486737479822,15.1134936925806,15.2804958249467,15.5575041351960,15.6135036569083,15.7247055541163,15.8909105402658,16.0013152912573,16.1117200422488,16.1667231311903,16.2773271687359,16.3325295442317,16.4987345303812,16.5539369058770,16.7755435540764,16.9419478267801,17.0523525777716,17.1629566153172,17.3845632635166,17.6065684848243,17.8283744195779,18.0499810677774,18.2167839135894,18.9370055202376,19.2134159708245,19.3796209569740,19.4906235676279,19.7122302158273,19.8226349668188];
y=[5.95683453237410,5.95683453237410,6.02158273381295,6.11870503597122,6.31294964028777,6.53956834532374,6.70143884892086,6.92805755395683,7.21942446043166,7.54316546762590,7.83453237410072,8.06115107913669,8.32014388489209,8.51438848920863,8.74100719424460,8.61151079136691,8.41726618705036,8.25539568345324,8.06115107913669,7.83453237410072,7.64028776978417,7.51079136690648,7.38129496402878,7.25179856115108,7.18705035971223,7.08992805755396,6.96043165467626,6.86330935251799,6.76618705035971,6.70143884892086,6.60431654676259,6.70143884892086,6.76618705035971,6.79856115107914,6.76618705035971,6.73381294964029,6.73381294964029,6.79856115107914,6.79856115107914,6.89568345323741,6.99280575539568,7.08992805755396,7.21942446043166,7.31654676258993,7.38129496402878,7.25179856115108,7.18705035971223,7.12230215827338,7.08992805755396,7.08992805755396,7.05755395683453,7.08992805755396,7.02517985611511,6.92805755395683,6.89568345323741,6.96043165467626,7.05755395683453,7.12230215827338,7.02517985611511,6.92805755395683,6.83093525179856,6.79856115107914,6.79856115107914,6.86330935251799,6.86330935251799,6.92805755395683,7.02517985611511,7.12230215827338,7.18705035971223,7.25179856115108,7.15467625899281,7.05755395683453,6.96043165467626,6.99280575539568,6.96043165467626,6.92805755395683,6.92805755395683,6.99280575539568,7.02517985611511,6.92805755395683,6.92805755395683,6.99280575539568,7.02517985611511,6.92805755395683,6.89568345323741,6.92805755395683,6.96043165467626,6.86330935251799,6.86330935251799,6.92805755395683,6.96043165467626,6.99280575539568,7.08992805755396,7.18705035971223,7.12230215827338,7.02517985611511,7.02517985611511,6.99280575539568,6.99280575539568,7.02517985611511,6.96043165467626,6.92805755395683,6.86330935251799,6.86330935251799,6.92805755395683,6.92805755395683,7.05755395683453,7.05755395683453,7.15467625899281,7.21942446043166,7.21942446043166,7.15467625899281,7.08992805755396,7.02517985611511,6.99280575539568,6.96043165467626,6.96043165467626,6.92805755395683,6.92805755395683,6.96043165467626,6.89568345323741,6.86330935251799,6.86330935251799,6.92805755395683,6.96043165467626,6.96043165467626,7.05755395683453,7.05755395683453,6.96043165467626,6.96043165467626,6.99280575539568,6.99280575539568,6.92805755395683];
%函数拟合
f=fit(x',y','smoothingspline');
x1=linspace(2,19,250)';
y1=f(x1);
%绘图
figure1 = figure;
axes1 = axes('Parent',figure1);
hold(axes1,'on');
plot(x,y,'LineWidth',2,'Parent',axes1);
plot(x1,y1,'LineWidth',2,'Parent',axes1);
legend('原函数图像','拟合函数图像');
xlabel('$x$','Interpreter','latex');
ylabel('$y$','Interpreter','latex');
xlim(axes1,[0 20]);
set(axes1,'FontSize',12,'LineWidth',2);

函数拟合代码很简单，后期会有详细推文讲解函数拟合问题，本期先先给出代码。整体来看还算拟合的好，部分拟合能够看出来有明显的错差，该方法并没有给出具体的函数解析式，类似于差值法，可以用于原图像的光滑处理。过冷水的实际问题是用具体函数表达式表示出第一个谷底区间的函数表达式，不需要给出整个图像区间的函数解析式，样问题就变得容易了，实际拟合一下图像：

使用函数拟合工具箱直接给出多项式拟合、高斯拟合、傅里叶拟合的函数图像，由拟合图像可知，傅里叶拟合和高斯拟合明显较好，使用起来比较简单。再一次告诫大家，在拟合图像的时候只知道多项式拟合就显得太不成熟了。

在此给出傅里叶拟合的代码：

[~,locs]=findpeaks(-y,'npeaks',1)
x3=x(2:locs);y3=y(2:locs);
f= fit(x3', y3','fourier3');
syms x 
y4=  f.a0 + f.a1*cos(x*f.w) + f.b1*sin(x*f.w) +f.a2*cos(2*x*f.w) + f.b2*sin(2*x*f.w) + f.a3*cos(3*x*f.w) + f.b3*sin(3*x*f.w);