sine curve fitting by recursive method

Find sine regression of periodic signal.

I have long period so I decomposed it in to small signals

Frequency = 50;                     % hertz
StopTime = 1/Frequency;             % seconds
FittingTime = (0:dt:StopTime-dt)';     % seconds
%%Sine wave:
FittingVoltage = sin(2*pi*Fc*FittingTime);

Then I want to do curve fitting by recursive method. whereas

RangeOfVector= 5

for  i = 0:1/Frequency:RangeOfVector
iter(i)=min(TrainTime):(1/Frequency):max(TrainTime)
end

This should process [x 0 0 0 0] then [0 x 0 0 0] then [0 0 x 0 0] and so on

Answer 1

[Mathieu NOE]

hello

this is a little demo you can adapt to your own needs ...

clc
clearvars
close all

% dummy signal (sinus + noise)
dt = 1e-4;
samples = 1000;
f = 50;

t = (0:samples-1)*dt;
s = 0.75*sin(2*pi*f*t) + 0.1 *rand(1,samples);

%%%%%%%%%%%%% main code %%%%%%%%%%%%%%%%%
ym = mean(s);                                                                   % Estimate offset
yu = max(s);
yl = min(s);
yr = (yu-yl);                                                                   % Range of ‘y’
yz = s-ym;
yzs = smoothdata(yz,'gaussian',25);                                             % smooth data to remove noise artifacts (adjust factors)
zt = t(yzs(:) .* circshift(yzs(:),[1 0]) <= 0);                                 % Find zero-crossings
per = 2*mean(diff(zt));                                                         % Estimate period
fre = 1/per;                                                                    % Estimate FREQUENCY

% stationnary sinus fit
fit = @(b,x)  b(1) .* (sin(2*pi*x*b(2) + b(3))) + b(4);   % Objective Function to fit
fcn = @(b) norm(fit(b,t) - s);                                                  % Least-Squares cost function
B = fminsearch(fcn, [yr/2;  fre;  0;  0;]);                                   % Minimise Least-Squares
amplitude = B(1)
frequency_Hz = B(2)
phase_rad = B(3)
DC_offset = B(4)

xp = linspace(min(t),max(t),samples);
yp = fit(B,xp);

figure(1),
plot(t, s, 'db',xp, yp, '-r')
legend('data + noise','model fit');

[MUHAMMAD SULEMAN]

Thank you very much, Mathieu. Its really helpful. May I ask for alogorithm for evaluating the dummy signal through sliding window because my time period is continous and I want to process signal for each cycle of sine wave.

Thanking you in anticipation.

[Mathieu NOE]

hello again

so the code has been modified to split into one period segments and do the fit on those segments;

as a consequence, all parameters of the fitted sine are stored in an array , so if you wish you can export them in a table

this is the updated code :

clc
clearvars
close all
% dummy signal (sinus + noise)
dt = 1e-4;
samples = 1000;
f = 50;
t = (0:samples-1)*dt;
s = 0.75*sin(2*pi*f*t) + 0.1 *rand(1,samples);
%%%%%%%%%%%%% main code %%%%%%%%%%%%%%%%%
ym = mean(s);                                                                   % Estimate offset
yu = max(s);
yl = min(s);
yr = (yu-yl);                                                                   % Range of ‘y’
yz = s-ym;
yzs = smoothdata(yz,'gaussian',25);                                             % smooth data to remove noise artifacts (adjust factors)

ind_zt = find(yzs(:) .* circshift(yzs(:),[1 0]) <= 0);                                 % Find zero-crossings indexes
% we keep only 1 of 2 indexes to have one period signal
ind_zt2 = ind_zt(1:2:length(ind_zt));
zt = t(ind_zt2);                                                                 % Find zero-crossings

for ci = 1:length(ind_zt2)-1
    % extract one period of signal 
    ind_extract = ind_zt2(ci):ind_zt2(ci+1);
    t_extract = t(ind_extract);
    s_extract = s(ind_extract);
    
    per = t(ind_zt2(ci+1)) - t(ind_zt2(ci));                                      % Estimate period
    fre = 1/per;                                                                   % Estimate FREQUENCY
    % stationnary sinus fit
    fit = @(b,x)  b(1) .* (sin(2*pi*x*b(2) + b(3))) + b(4);                     % Objective Function to fit
    fcn = @(b) norm(fit(b,t_extract) - s_extract);                               % Least-Squares cost function
    B = fminsearch(fcn, [yr/2;  fre;  0;  0;]);                                   % Minimise Least-Squares
    amplitude(ci) = B(1)
    frequency_Hz(ci) = B(2)
    phase_rad(ci) = B(3)
    DC_offset(ci) = B(4)
    yp = fit(B,t_extract);
    figure(ci),
    plot(t_extract, s_extract, 'db',t_extract, yp, '-r')
    legend('data + noise','model fit');
end

[MUHAMMAD SULEMAN]

Thank you very much. You are a life saver.

Best Wishes and Regards,

[Mathieu NOE]

You're welcome !