clear

%% parameter values

% system parameters
c1=0.5;
c3cr=40*c1/9;

% range of c3 prarameter for globally stable, close to bif and after bif
c3_stab_range=[0.1, 0.5]*c3cr;
% c3_close_range=[0.985, 0.995]*c3cr;
c3_close_range=[0.9, 0.995]*c3cr;
c3_after_range=[1.01, 1.5]*c3cr;

% range of initial conditions for x0 and x0d
x0_range=[0.7, 2];
xd0_range=[0.7, 2];

% number of samples for each range
num_sampl=200;

% time span and time interval
tfinal=700;
t_steps=2000;
tspan=linspace(0,tfinal,t_steps);

%% generate random parameters

c3_stab=rand(num_sampl,1)*(c3_stab_range(2)-c3_stab_range(1))+c3_stab_range(1);
c3_close=rand(num_sampl,1)*(c3_close_range(2)-c3_close_range(1))+c3_close_range(1);
c3_after=rand(num_sampl,1)*(c3_after_range(2)-c3_after_range(1))+c3_after_range(1);
x0(:,1)=rand(num_sampl,1)*(x0_range(2)-x0_range(1))+x0_range(1);
x0(:,2)=rand(num_sampl,1)*(xd0_range(2)-xd0_range(1))+xd0_range(1);

c3=[c3_stab;c3_close;c3_after];
x0=[x0;x0;x0];

%% generate classification target
target(2*num_sampl+1:3*num_sampl,1)=3;
target(1*num_sampl+1:2*num_sampl,1)=2;
target(1:num_sampl,1)=1;
target=categorical(target);

%% perform simulations
options = odeset('RelTol', 1e-6, 'AbsTol', 1e-8);
x=cell(length(c3),1);
for i=1:num_sampl*3
    [~, x{i}] = ode45 (@(t, xt) system_NLD(t, xt, [c1, c3(i)]), tspan, x0(i,:), options);
end

%% save data
clear("i");
% save NLD_4training_raw_simulation_data_c1_05.mat
save NLD_4testing_raw_simulation_data_c1_01.mat

%% normalize data polar_01_log_mvm

clear

load NLD_4training_raw_simulation_data_c1_05.mat
% load NLD_4testing_raw_simulation_data_c1_01.mat

rho=cell(length(x),1);
for i=1:length(x)
    polar_pure=sqrt(x{i}(:,1).^2+x{i}(:,2).^2);
    max_pol=max(polar_pure);
    polar_01=polar_pure/max_pol;
    polar_01_log=log(polar_01);
    for j=20:length(polar_01_log)
        if max([polar_01_log(j-19:j)])<-8
            polar_01_log=polar_01_log(1:j);
            break
        end
    end
    polar_01_log_resumpled=spline(linspace(1,length(polar_01_log),length(polar_01_log)), polar_01_log, linspace(1,length(polar_01_log),1000));
    rho{i}=movmean(polar_01_log_resumpled,round(length(polar_01_log_resumpled)/20));
%     rho{i}=transpose(movmean(polar_01_log,round(length(polar_01_log)/10)));
end

% save normalized data
clear("x");
clear("i");
save NLD_4training_polar_01_log_mvm_c1_05.mat
% save NLD_4testing_polar_01_log_mvm_c1_01.mat

%% normalize data min/max [-1 1]

clear

% load NLD_4training_raw_simulation_data_c1_05.mat
load NLD_4testing_raw_simulation_data_c1_01.mat

rho=cell(length(x),1);
for i=1:length(x)
    massimo=max(x{i}(:,1));
    minimo=min(x{i}(:,1));
    minmax=zeros(length(x{i}(:,1)),1);
    for j=1:length(x{i})
        if x{i}(j,1)>0
            minmax(j)=x{i}(j,1)/massimo;
        else
            minmax(j)=x{i}(j,1)/abs(minimo);
        end
    end
    for j=50:length(minmax)
        if max([minmax(j-49:j)])<1e-8
            minmax=minmax(1:j);
            break
        end
    end
    minmax_resumpled=spline(linspace(1,length(minmax),length(minmax)), minmax, linspace(1,length(minmax),1000));
    rho{i}=minmax_resumpled;
end

% save normalized data
clear("x");
clear("i");
% save NLD_4training_minmax_c1_05.mat
save NLD_4testing_minmax_c1_01.mat

%% normalize data polar [0 1]

clear

% load NLD_4training_raw_simulation_data_c1_05.mat
load NLD_4testing_raw_simulation_data_c1_01.mat

rho=cell(length(x),1);
for i=1:length(x)
    polar_pure=sqrt(x{i}(:,1).^2+x{i}(:,2).^2);
    max_pol=max(polar_pure);
    polar_01=polar_pure/max_pol;
    for j=20:length(polar_01)
        if max([polar_01(j-19:j)])<1e-8
            polar_01=polar_01(1:j);
            break
        end
    end
    polar_01_resumpled=spline(linspace(1,length(polar_01),length(polar_01)), polar_01, linspace(1,length(polar_01),1000));
    rho{i}=polar_01_resumpled;
end

% save normalized data
clear("x");
clear("i");
% save NLD_4training_polar_01_c1_05.mat
save NLD_4testing_polar_01_c1_01.mat

%% normalize data polar [0 1] mvm

clear

% load NLD_4training_raw_simulation_data_c1_05.mat
load NLD_4testing_raw_simulation_data_c1_01.mat

rho=cell(length(x),1);
for i=1:length(x)
    polar_pure=sqrt(x{i}(:,1).^2+x{i}(:,2).^2);
    max_pol=max(polar_pure);
    polar_01=polar_pure/max_pol;
    for j=20:length(polar_01)
        if max([polar_01(j-19:j)])<1e-8
            polar_01=polar_01(1:j);
            break
        end
    end
    polar_01_resumpled=spline(linspace(1,length(polar_01),length(polar_01)), polar_01, linspace(1,length(polar_01),1000));
    rho{i}=movmean(polar_01_resumpled,round(length(polar_01_resumpled)/20));
%     rho{i}=transpose(movmean(polar_01_log,round(length(polar_01_log)/10)));
end

% save normalized data
clear("x");
clear("i");
% save NLD_4training_polar_01_mvm_c1_05.mat
save NLD_4testing_polar_01_mvm_c1_01.mat


%% normalize data polar_01_log

clear

load NLD_4training_raw_simulation_data_c1_05.mat
% load NLD_4testing_raw_simulation_data_c1_05.mat

rho=cell(length(x),1);
for i=1:length(x)
    polar_pure=sqrt(x{i}(:,1).^2+x{i}(:,2).^2);
    max_pol=max(polar_pure);
    polar_01=polar_pure/max_pol;
    polar_01_log=log(polar_01);
    for j=20:length(polar_01_log)
        if max([polar_01_log(j-19:j)])<-8
            polar_01_log=polar_01_log(1:j);
            break
        end
    end
    polar_01_log_resumpled=spline(linspace(1,length(polar_01_log),length(polar_01_log)), polar_01_log, linspace(1,length(polar_01_log),1000));
    rho{i}=polar_01_log_resumpled;
%     rho{i}=transpose(movmean(polar_01_log,round(length(polar_01_log)/10)));
end

% save normalized data
clear("x");
clear("i");
save NLD_4training_polar_01_log_c1_05.mat
% save NLD_4testing_polar_01_log_c1_05.mat






















