The figures 26a and 27 are calculated by using all 3 components of the Roessler system. The data should be in ASCII format, 3 columns, when using the commandline RP software. You can check whether the software has read the data correctly by the output:
Code: Select all
./rp_x86_64i -n MAX -w 2 -i roessler -c -p roessler_hist -e 0.862
Output:
Code: Select all
Used file: roessler (30003 data points in 3 columns read).
BTW: The histogram will include all lines, thus, the option
-l will not work here. This option is only for the RQA measures.
Here is some Matlab code to reproduce the results of the figures. It is also using the commandline RP software.
Code: Select all
%% Line length distribution for Roessler
% Fig. 26A
clear
options=odeset('initialstep',.01,'MaxStep',.1);
% sampling time difference is 0.2
time=0:.2:2050;
% integrate Roessler system
% roessler_abc = [0.2 0.2 5.7];
[t, x]=ode45('roessler',time,[0.1 0.4 0.3],options);
% save data as 3-column vector, ASCII file
x1 = x(251:end,:); save roessler x1 -ascii -tabs
% calculate the recurrence plot and cumulatove line length distribution
% using commandline RP tool
clear, e = 0.862
unix(['./rp_x86_64i -n MAX -w 2 -i roessler -c -p roessler_hist -e ',num2str(e)]);
% load histogram into Matlab
h=load(['roessler_hist']); p = h(:,2);
% plot the histogram
clf
semilogy(1:200,p(1:200,1:4:end),'k'), xlim([1 200])
xlabel('Length l (units)'), ylabel('Total number')
set(gca,'xtick',[0:50:200],'xtickl',[0:50:200]), ylim([1e1 1e8])
Output:
Code: Select all
./rp_x86_64i -n MAX -w 2 -i roessler -c -p roessler_hist -e 0.862
Used file: roessler (30003 data points in 3 columns read).
Used parameters: embedding dimension m = 1
embedding delay t = 1
recurrence threshold e = 0.862
Theiler window w = 2
minimal diagonal line l_min = 2
minimal vertical line v_min = 2
maximum norm used
Calculate recurrence points% a=.2; % periodic
% b=.7;
% c=3.6;
|******************************| 100%
Computation time: 0 min 2 sec
Recurrence quantification analysis:
RR: 0.009573
DET: 0.9717 LAM: 0.522
DET/RR: 101.5 LAM/DET: 0.5372 W_prob: 87
L_max: 404 V_max: 4 W_max: 5241
L_mean: 14.02 TT: 2.255 W_mean: 140.9
L_entr: 3.048 V_entr: 0.6153 W_entr: 3.26
DIV: 0.002475 T1: 101.9 F_min: 0.0001908
T2: 144.2
Some code for Fig. 27:
Code: Select all
%% Line length distribution for Roessler
% Fig. 27
% integrate Roessler system
clear,options=odeset('initialstep',.01,'MaxStep',.1);
% sampling time difference is 0.2
time=0:.2:2050;
% Roessler with C = 9
roessler_abc = [0.1 0.1 9];
[t, x]=ode45('roessler',time,[0.1 0.4 0.3],options,roessler_abc);
% save data as 3-column vector, ASCII file
x1 = x(251:end,:); save roessler9 x1 -ascii -tabs
% Roessler with C = 30
roessler_abc = [0.1 0.1 30];
[t, x]=ode45('roessler',time,[0.1 0.4 0.3],options,roessler_abc);
% save data as 3-column vector, ASCII file
x1 = x(251:end,:); save roessler30 x1 -ascii -tabs
clear, e = 2;
unix(['./rp_x86_64i -s -n MAX -w 2 -i roessler30 -c -p roessler_hist -e ',num2str(e)]);
h=load(['roessler_hist']);
p(:,1) = h(:,2);
unix(['./rp_x86_64i -s -n MAX -w 2 -i roessler9 -c -p roessler_hist -e ',num2str(e)]);
h=load(['roessler_hist']);
p(:,2) = h(:,2);
clf
semilogy(.2*(6:200),p(6:200,1),'k--'), hold on
semilogy(.2*(6:200),p(6:200,2),'k-')
xlabel('Length l (s)'), ylabel('Total number')
yl = ylim;line([12 12],yl,'linestyle',':','color',[0 0 0])
line([24.8 24.8],yl,'linestyle',':','color',[0 0 0])
You will need the file
roessler.m for integrating the Roessler system:
Code: Select all
function dy=roessler(t,y, dummy, abc)
dy=zeros(3,1);
if nargin < 4
a=.2;
b=.2;
c=5.7;
else
a = abc(1);
b = abc(2);
c = abc(3);
end
dy(1)=-y(2)-y(3);
dy(2)=y(1)+a*y(2);
dy(3)=b+y(3)*(y(1)-c);
However, I have just realised that in Fig 27 the two curves are exchanged: c=30 is the solid one and c=9 the dashed one.