How to find minimum value for a loop iteration if we have two variable at some min value i need to save time and velocity ?

1 view (last 30 days)
%file first
function [r_b,v_b,zenit,Elevation,Azi,close,minimum_time]=main(a,e,eta,omega_1,omega_2)
GM= 398600.44;
n = sqrt(GM/(a*a*a)) %mean angular velocity/mean motion
T= 2*pi/(n);
k=1;
for i=0:10:2*T
M= n *(i-0); %mean anomly
E0=M;
err=1;
while(err>10^-12)
E= M+e*sin(E0);
err=abs(E-E0); % I want to understand this step
E0=E; % eccentric anomly
end
v=2*atan(sqrt(1+e)*tan(E/2)/sqrt(1-e)); %atan2 is not working % true anomly
r = a*(1-e*cos(E)); %distance of the satellite
r_b= [r*cos(v);r*sin(v);0] ; %position vector in the orbit system
r_b1(k,:)=[a*(cos(E)-e);a*(sqrt(1-(e*e)))*sin(E);0]; %position vector in the orbit system
v_b(k,:)=[-a*n*sin(E)*a/r;a*n*sqrt(1-(e*e))*cos(E)*a/r]; %velocity vector in the orbit system
r_3_omega_1= [cos(-omega_1) sin(-omega_1) 0;-sin(-omega_1) cos(-omega_1) 0;0 0 1];
r_3_omega_2= [cos(-omega_2) sin(-omega_2) 0;-sin(-omega_2) cos(-omega_2) 0;0 0 1];
r_1=[1 0 0;0 cos(-eta) sin(-eta);0 -sin(-eta) cos(-eta)];
r_i= r_3_omega_1*r_1*r_3_omega_2*r_b;
r_i1(k,:)=transpose(r_i);
w_a=((2*pi)/86164);
phi=w_a*(i-0);
r_phi= [cos(phi) sin(phi) 0;-sin(phi) cos(phi) 0;0 0 1];
r_e=r_phi*r_i ;
r_e1(k,:)=r_phi*r_i ;
z_lat(k,:)= r_e1(k,3);
y_lat(k,:)=r_e1(k,2);
x_lat(k,:)=r_e1(k,1);
lat(k,:)= atan2((z_lat(k,:)),(sqrt((x_lat(k,:)*x_lat(k,:))+(y_lat(k,:)*y_lat(k,:)))));
long(k,:)= atan2(y_lat(k,:),x_lat(k,:)) ;
r_w=[4075.53022;931.78130;4801.61819];
lat_w(k,:)=atan2(r_w(3,1),(sqrt(r_w(1,1)*r_w(1,1)+r_w(2,1)*r_w(2,1))));
long_w(k,:)=atan2(r_w(2,1),r_w(1,1));
A = [-sin(lat_w(k,:))*cos(long_w(k,:)) -sin(long_w(k,:)) cos(lat_w(k,:))*cos(long_w(k,:));-sin(lat_w(k,:))*sin(long_w(k,:)) cos(long_w(k,:)) cos(lat_w(k,:))*sin(long_w(k,:));cos(lat_w(k,:)) 0 sin(lat_w(k,:))];
r_top = transpose(A)*(r_e-r_w);
r_top1(k,:)=transpose(A)*(r_e-r_w);
close=min(r_e-r_w,T)
S=sqrt(r_top(1,1).^2+r_top(2,1).^2+r_top(3,1).^2);
S_1(k,:)=sqrt(r_top1(k,1).^2+r_top1(k,2).^2+r_top1(k,3).^2);
zenit(k,:)=acos(r_top(3,1)/S);
Elevation(k,:)=(pi/2) -zenit(k,:);
Azi(k,:)= atan2(-r_top1(k,2),r_top1(k,1));
k=k+1;
end
% subplot(1,1,1);
figure()
plot3(r_e1(:,1),r_e1(:,2),r_e1(:,3));
save('fig(1)');
% subplot(1,1,2);
figure()
plot3(r_i1(:,1),r_i1(:,2),r_i1(:,3));
save('fig(2)');
%subplot(1,1,3)
figure()
pax=polaraxes;
pax.ThetaDir = 'clockwise';
polarplot(Azi,rad2deg(Elevation))
pax.RDir = 'reverse';
%fig4=polarplot(zenit,rad2deg(Elevation))
% subplot(1,1,4)
figure()
plot(rad2deg(long),rad2deg(lat))
hold on
load coastlines
geoshow(coastlat,coastlon)
hold off
save('fig(5)')
minimum_closest_distance=min(r_e-r_w) ;
minimum_time=(T),
end
%file second, run this file
a_gps=26500;
e_gps=0.01;
eta_gps=deg2rad(55);
omega_1_gps=deg2rad(0);
omega_2_gps=deg2rad(0);
[r_b,v_b,zenit,Elevation,Azi,close]= GPSex1original(a_gps,e_gps,eta_gps,omega_1_gps,omega_2_gps);
a_gnss=42164.140;
e_gnss=0;
eta_gnss=deg2rad(63);
omega_1_gnss=deg2rad(0);
omega_2_gnss=deg2rad(0);
[gnss_r_b,gnss_v_b,gnss_zenit,gnss_Elevation,gnss_Azi]= GPSex1original(a_gnss,e_gnss,eta_gnss,omega_1_gnss,omega_2_gnss);
a=29994;
e=0;
eta=deg2rad(56);
omega_1=deg2rad(0);
omega_2=deg2rad(0);
[galileo_r_b,galileo_v_b,galileo_zenit,galileo_Elevation,galileo_Azi]= GPSex1original(a,e,eta,omega_1,omega_2);
a=6836;
e=0.004;
eta=deg2rad(87);
omega_1=deg2rad(0);
omega_2=deg2rad(0);
[champ_r_b,champ_v_b,champ_zenit,champ_Elevation,champ_Azi]= GPSex1original(a,e,eta,omega_1,omega_2);
a=26554;
e=0.7;
eta=deg2rad(65);
omega_1=deg2rad(245);
omega_2=deg2rad(270);
[molniya_r_b,molniya_v_b,molniya_zenit,molniya_Elevation,molniya_Azi]= GPSex1original(a,e,eta,omega_1,omega_2);
% I have attached files in the attachment please look into it and help me.

Answers (0)

Categories

Find more on Automotive in Help Center and File Exchange

Products


Release

R2019a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!