実験データから、正弦波のパラメータを同定する。

下図のような実験データから、正弦波のパラメータを同定したいです。

Y = x1*sin(2π*x2 + x3) + x4

x1:振幅

x2:周波数

x3:位相

x4:オフセット

データは、正弦波で、計測ノイズが多少乗っている程度です。

lsqcurvefit関数を使って試してみましたが、初期値が悪かったり、データ点数が多いと、

上手くパラメータフィッティングできません。

何か良い方法がありましたら、アドバイスお願いします。

[Takumi]

信号処理詳しくないので参考になるかわかりませんが，フーリエ変換を利用してはいかがでしょうか．

スペクトルの周波数と振幅でsin波を推定できるかと思います．DC成分は平均でざっくり推定できますね．

（ノイズがあると位相角はうまく求まらなさそうです．．．）

高速フーリエ変換

x1 = 2;     % 振幅
x2 = 10;    % 基本周波数
x3 = pi/4;  % 初期位相
x4 = 0.5;   % オフセット

Fs = 1000;	% サンプリング周期                  
Ts = 1/Fs;	% サンプリング周波数     
L = 500;	% 信号長
t = (0:L-1)*Ts; % 時間配列

y = x1*sin(2*pi*x2*t+x3)+x4;    % 信号
plot(t,y,'.-');
xlabel('t[s]');

% オフセットは平均で推定
x4est = mean(y)

% フーリエ変換
y = y-mean(y);  % DC成分を排除
Y = fft(y);
P2 = abs(Y/L);
P1 = P2(1:L/2+1);
P1(2:end-1) = 2*P1(2:end-1);
f = Fs*(0:(L/2))/L;

tol = 1e-6;
Y(abs(Y) < tol) = 0;
phaseY = unwrap(angle(Y));  % 位相


figure;
subplot(2,1,1); plot(f,P1)
xlabel('f[Hz]');
ylabel('|P1|');

subplot(2,1,2); plot(f,rad2deg(phaseY(1:L/2+1)))
xlabel('f[Hz]');
ylabel('phase[deg]');