function dsdt = eqmoto_crtbp(t,s)  % quantità adimensionalizzate nel sistema sinodico CR3BP

rho = s(1:3);
vu = s(4:6);

% costanti gravitazionali
mu_e = cspice_bodvrd('399','GM',1);
mu_m = cspice_bodvrd('301','GM',1);
mu = mu_m/(mu_e+mu_m);


% distanze
rho1 = sqrt((rho(1)+mu)^2+rho(2)^2+rho(3)^2); % distanza terra
rho2 = sqrt((rho(1)+mu-1)^2+rho(2)^2+rho(3)^2); % distanza luna

% derivate pseudo potenziale
dUdcsi = rho(1)-(1-mu)*(rho(1)+mu)/rho1^3-mu*(rho(1)-1+mu)/rho2^3;
dUdeta = rho(2)-(1-mu)*rho(2)/rho1^3-mu*rho(2)/rho2^3;
dUdzita = -(1-mu)*rho(3)/rho1^3-mu*rho(3)/rho2^3;

% equazioni del moto
dsdt = [vu; 2*vu(2)+dUdcsi; -2*vu(1)+dUdeta; dUdzita];
