% Calcul des trajectoires
% XX et YY sont les dérivées de X et Y
DAA1=['2*(' XX ').*(' YY ')'];
DAA2=['-(' X ').^2+(' Y ').^2+Qy*Py-Qy*(' Y ')+Px*(' X ')+Qx*(' X ')'];
DAA21=['-Px*Qx-Py*(' Y ')'];
DAA3=['((' YY ')-(' XX ')).*((' YY ')+(' XX '))'];
DAA4=['Qy*Px-Qy*(' X ')-Py*(' X ')+2*(' X ').*(' Y ')+Py*Qx-Px*(' Y ')'];
DAA41=['-Qx*(' Y ')'];
DAA=['((' DAA1 ').*(' DAA2 '' DAA21 '))-((' DAA3 ').*(' DAA4 '' DAA41 '))'];
% clf
% % Graphique des solutions de T
% axis('square')
 if choix == 3
   u=4*pi:pi/300:6*pi;
   DA=eval(DAA);
%   plot(u,DA,[4*pi 6*pi],[0 0])
 else
   u=0:pi/300:2*pi;
   DA=eval(DAA);
%   plot(u,DA,[0 2*pi],[0 0])
 end
% title('Les zéros de cette fonction sont le(s) T possible(s) (rebond)');
% xlabel('Appuyez sur une touche pour continuer')
% pause


% Recherche des zeros
mm=length(DA);
S=[];
for i=1:mm-1
  if sign(DA(i))*sign(DA(i+1)) == -1 | 0
    S=[S (u(i)+u(i+1))/2];
  end
end
u=S;
% Définition de T avec le u trouvé plus haut
Tx=eval(X);
Ty=eval(Y);
T=[Tx;Ty];
