function [eo]=eloSF(y0,t) %Función para obtener los elementos orbitales a partir de la posición y velocidad
%Dadas la posicion y velocidad en el sistema de referencia movil, obtiene
%los elementos orbitales en el sistema de referencia fijo
global wt mu
Mtheta=[cos(-wt*t) -sin(-wt*t) 0;sin(-wt*t) cos(-wt*t) 0;0 0 1];
Vm=y0(4:6);
rp=sqrt(y0(1)^2+y0(2)^2);vpm=rp*wt;
rpx=y0(1);rpy=y0(2);
if rpy>=0
theta=acos(rpx/rp);
else
theta=2*pi-acos(rpx/rp);
end
etheta=[-sin(theta) cos(theta) 0];
Vpm=vpm*etheta;
Vf=Vm+Vpm;
% Vm=Vf-Vpm;
y0(4:6)=Vf;
r=(Mtheta^-1*y0(1:3)')';
v=(Mtheta^-1*y0(4:6)')';
% r=y0(1:3);
% v=y0(4:6);
eo=zeros(1,6);%Inicializamos el vector de los elementos orbitales
r;
v;
h=cross(r,v);
%EC
EC=cross(v,h)/mu-r/norm(r);
eo(1)=norm(EC);
%SMA(m)
epsilon=norm(v)^2/2-mu/norm(r);
SMA=-mu/2/epsilon;
eo(2)=SMA;
%TA(grados)
if r*v'>=0
TA=acos(r*EC'/norm(r)/norm(EC))*180/pi;
else
TA=(2*pi-acos(r*EC'/norm(r)/norm(EC)))*180/pi;
end
eo(6)=TA;
%INC(grados)
INC=acos(h(3)/norm(h))*180/pi;
eo(3)=INC;
%%OM(grados)
K=[0 0 1]';
n=cross(K,h)/norm(cross(K,h));
if n(2)>0
OM=acos(n(1))*180/pi;
else
OM=(2*pi-acos(n(1)))*180/pi;
end
eo(4)=OM;
%%W(grados)
if EC(3)>0
W=acos(EC*n'/norm(EC))*180/pi;
else
W=(2*pi-acos(EC*n'/norm(EC)))*180/pi;
end
eo(5)=W;
% % eo=[0.15 a+2e3 147.22*pi/180 90*pi/180 0*pi/180 180*pi/180]
% % [coord]=rv(eo)
% % y0=coord
% % Vf=y0(4:6);
% % rp=sqrt(y0(1)^2+y0(2)^2);vpm=rp*wt;
% % rpx=y0(1)^2;rpy=y0(2);
% % if rpx>0
% % theta=atan(rpy/rpx);
% % else
% % theta=pi+atan(rpy/rpx);
% % end
% % etheta=[-sin(theta);cos(theta);0];
% % Vpm=vpm*etheta;
% %
% % Vm=Vf-Vpm;
% %
% % y0(4:6)=Vm
e-REdING. Biblioteca de la Escuela Superior de Ingenieros de Sevilla.
