Question: Plotting gait diagrams

Hello,

I would like to plot gait diagrams (the lines you can see on the picture belowà from the solutions obtained with a NL oscillator (composed with 8 coupled odes). Here the result that I would like to obtain.

Initial plot:

 

Desired plot

 

 

I would like to obtain 4 lines corresponding to the 4 elliptic trajectories obtained with the NL oscillator. The four lines should be done like this. When the trajectory is above 0, the line should be colored in green. When the trajectory is below 0, the line should be colored in black. 

May you help me to define this kind of graph called gait diagrams from the solution of the NL oscillator ?

Here you can find my maple code:

K:=Matrix([<0, -1, 1, -1>,<-1, 0, -1, 1>,<-1, 1, 0,-1>,<1, -1, -1,0>]);

for i to 4
do
r[i]:=sqrt((u[i](t))^2+(v[i](t))^2):
omega[i]:=omega[sw]/(1+exp(b*v[i](t)))+omega[st]/(1+exp(-b*v[i](t))):
Equ[i]:=diff(u[i](t),t)=Au*(1-r[i]^2)*u[i](t)-omega[i]*v[i](t):
Eqv[i]:=diff(v[i](t),t)=Av*(1-r[i]^2)*v[i](t)+omega[i]*u[i](t)+MatrixVectorMultiply(K,<seq(v[i](t),i=1..4)>)[i]:
EqSys[i]:=[Equ[i],Eqv[i]]:
end do:

paramsCycle:=omega[st]=4*2*Pi,omega[sw]=2*Pi,Au=5,Av=50,b=100;
params:=paramsCycle;

Differential system 
sys:=map(op,eval([seq(EqSys[i],i=1..4)],[params]));
ic:=[u[1](0)=0, v[1](0)=0,u[2](0)=0, v[2](0)=-0.1,u[3](0)=0, v[3](0)=0.1,u[4](0)=0, v[4](0)=0.1];
Résolution1
res:=dsolve([sys[],ic[]],numeric):
Initial boundaries
tcalc:=4;
ic2:=[seq(u[i](0)=eval(u[i](t), res(tcalc)),i=1..4),seq(v[i](0)=eval(v[i](t), res(tcalc)),i=1..4)];
Résolution2
res:=dsolve([sys[],ic2[]],numeric):

tmax:= 40:
numpts:=100*tmax:
plots:-odeplot(res,[t,v[1](t)],0..tmax,thickness=2, view=[0..5, -1.5..1.5],numpoints = numpts);
plots:-odeplot(res,[t,v[2](t)],0..tmax,thickness=2, view=[0..5, -1.5..1.5],numpoints = numpts);
plots:-odeplot(res,[t,v[3](t)],0..tmax,thickness=2, view=[0..5, -1.5..1.5],numpoints = numpts);
plots:-odeplot(res,[t,v[4](t)],0..tmax,thickness=2, view=[0..5, -1.5..1.5],numpoints = numpts);
plots:-odeplot(res,[seq([t,v[i](t)+i*5],i=1..4)],0..tmax,thickness=2,view=[0..5,0..25], numpoints = numpts);

Thanks a lot for your help

Please Wait...