% parallel coordinate lorenz system as a surface

clf
clear all

dt=.01;
tmax=4;
sigma = 10.0; rho = 28.0; beta = 2.6667;
x=1; y=0; z=0;

index=1;

for t=0:dt:tmax
    %======================================
    % The parallel plot
    
    dx_dt = sigma *(y-x); 
    dy_dt = rho * x - y - x*z;
    dz_dt = x*y - beta * z;
    
    x= x + dx_dt * dt;
    y= y + dy_dt * dt;
    z= z + dz_dt * dt;
    c= [1 0 t/tmax]';
    
    pos(1,index)= x;  
    pos(2,index)= y ;  
    pos(3,index)= z ;
    col(:,index)=c;
    index=index+1;

end

%=================================

subplot(3,1,1)
hold on
for i = 1:length(pos)
    plot(pos(:,i), 'color', col(:,i) ); 
end

%draw the parallel axis
axisx = [1,3; 1,1; 2,2; 3,3; ]';
axisy = [0,0; 0,1; 0,1; 0,1; ]';
line(axisx, axisy, 'color', 'k');
set(gca, 'xticklabel', {'x1', ' ', 'x2',' ','x3'})

%=================================

subplot(3,1,2)

set(gca,'color',[1 1 1])
box on
view(22,63)
axis vis3d
xlabel('x1');ylabel('x2');zlabel('x3');
hold on
for i = 1:length(pos)
    plot3(pos(1,i),pos(2,i),pos(3,i),...
        'markerfacecolor', col(:,i),...
        'markersize',1 ,'color',col(:,i),...
        'marker','o',...
        'linewidth',3)
end


rotate3d on

%=================================

subplot(3,1,3)

set(gca,'color',[1 1 1])
box on
view(22,63)
axis vis3d

xs=1:3;
ys=0:dt:tmax;
zs=pos';
sh=surf(xs,ys,zs)
set(sh, 'edgecolor','none')
colormap(hot)
set(gca, 'zlim', [-20 50])
rotate3d on