matlab練習程序（傅里葉擬合）

對於一組數據，通常可以用多項式來擬合，當然對於有周期規律的數據，我們也可以用傅里葉級數來擬合。

傅里葉級數公式形式如下：

當我們確定好n之后，關鍵就是求出A0、an、bn和w即可。

由於有待求系數在非線性函數cos和sin中，我們用非線性最優化方法來求解。

matlab代碼如下：

clear all;close all;clc;

n=7;                                                 %傅里葉級數的階數
x = (0:0.1:5)';
y = x - x.^2+ 10*cos(2*x).*sin(x)+5*cos(2*x);        %構造一個模型
% y = [zeros(25,1);ones(26,1)];                      %試試階躍
plot(x,y,'ro')

epsilon = 1e-6;
epsilon_inv = 1/epsilon;
pre=rand(2*n+2,1);
for i=1:500                                         %非線性迭代優化
    
    f0 = func(pre,x,n);
    g = y - f0;
    
    %數值計算雅克比
    for k = 1:length(pre)
        x_ = pre;
        x_(k) = pre(k) + epsilon;
        jac(:, k) = (func(x_,x,n) - f0) .* epsilon_inv;
    end
    J = jac;
    
    delta = inv(J'*J + eye(2*n+2))*J'* g;
    
    pcur = pre+delta;
    if norm(delta) <1e-16
        break;
    end
    pre = pcur;
end

hold on;
x=0:0.01:5;
plot(x,func(pre,x,n),'g.');
legend('待擬合數據','傅里葉擬合結果')

function y=func(a,x,n)          %傅里葉級數函數
y=a(1);
for i=1:n
    y= y + a(i*2)*cos(i*x*a(end)) + a(i*2+1)*sin(i*x*a(end));
end
end

結果如下：

~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

感興趣的朋友還可以和多項式擬合比較一下，下面是和6階做的對比，代碼我就不貼了，明顯傅里葉擬合效果更好些（當然，如果多項式堆階數效果應該也還可以）。