How can i save the data of "Long" and "Lat" as a matrix?
clc
clear
format long g
% PROBLEM = GroundTrack Always same
ge = 398600.8; % Earth gravitational constant
% ge = 6.67191*(10^-11) ;
TWOPI = 2*pi;
MINUTES_PER_DAY = 1440.;
MINUTES_PER_DAY_SQUARED = (MINUTES_PER_DAY * MINUTES_PER_DAY);
MINUTES_PER_DAY_CUBED = (MINUTES_PER_DAY * MINUTES_PER_DAY_SQUARED);
fname = 'tle.txt';
fid = fopen(fname, 'r');
while (1)
% read first line
tline = fgetl(fid);
if ~ischar(tline)
break
end
Cnum = tline(3:7); % Catalog Number (NORAD)
SC = tline(8); % Security Classification
ID = tline(10:17); % Identification Number
epoch = str2double(tline(19:32)); % Epoch
TD1 = str2double(tline(34:43)); % first time derivative
TD2 = str2double(tline(45:50)); % 2nd Time Derivative
ExTD2 = tline(51:52); % Exponent of 2nd Time Derivative
BStar = str2double(tline(54:59)); % Bstar/drag Term
ExBStar = str2double(tline(60:61)); % Exponent of Bstar/drag Term
BStar = BStar*1e-5*10^ExBStar;
Etype = tline(63); % Ephemeris Type
Enum = str2double(tline(65:end)); % Element Number
TLEYear = tline(19:20) ;
TLEEpochDate = str2double(tline(21:32));
% read second line
tline = fgetl(fid);
if ~ischar(tline)
break
end
i = str2double(tline(9:16)); % Orbit Inclination (degrees)
raan = str2double(tline(18:25)); % Right Ascension of Ascending Node (degrees)
e = str2double(strcat('0.',tline(27:33))); % Eccentricity
omega = str2double(tline(35:42)); % Argument of Perigee (degrees)
M = str2double(tline(44:51)); % Mean Anomaly (degrees)
no = str2double(tline(53:63)); % Mean Motion
a = ( ge/(no*2*pi/86400)^2 )^(1/3); % semi major axis (m)
rNo = str2double(tline(64:68)); % Revolution Number at Epoch
end
fclose(fid);
satdata.epoch = epoch;
satdata.norad_number = Cnum;
satdata.bulletin_number = ID;
satdata.classification = SC; % almost always 'U'
satdata.revolution_number = rNo;
satdata.ephemeris_type = Etype;
satdata.xmo = M * (pi/180);
satdata.xnodeo = raan * (pi/180);
satdata.omegao = omega * (pi/180);
satdata.xincl = i * (pi/180);
satdata.eo = e;
satdata.xno = no * TWOPI / MINUTES_PER_DAY;
satdata.xndt2o = TD1 * 1e-8 * TWOPI / MINUTES_PER_DAY_SQUARED;
satdata.xndd6o = TD2 * TWOPI / MINUTES_PER_DAY_CUBED;
satdata.bstar = BStar;
GMTDelay = 0 ;
ERS = 360/(24*3600) ; % Earth Rotation Speed in deg/s
mapll = load('coast') ;
mlat = mapll.lat ;
mlong = mapll.long ;
plot (mlong,mlat,'b')
% xlim([-180 180])
% ylim([-90 90])
%
hold on
rgb = imread('Earthmap.jpg'); % read world equidistant projection image
map_image=image(-180:180,90:-1:-90,rgb); %
hold on
day_night_terminator=area(1,1,1,'FaceColor',[0.0,0.0,0.0],'FaceAlpha',0.5,'EdgeColor','None');
tnow = datetime('now') ;
[baseline,coords] = plotdaynightterminator(tnow-(GMTDelay/24)) ;
set(day_night_terminator,'BaseValue',baseline,'XData',coords(:,1),'YData',coords(:,2));
clear mapll ;
rad = 57.2957795130823 ;
% tnow = datetime('now') ;
FirstYear = ['1-Jan-20' TLEYear ' 00:00:00'];
EDay = days(tnow - FirstYear) ;
DeltaDay = (EDay - TLEEpochDate) ;
DeltaMinute = DeltaDay * 24 * 60 ;
step = 30 ; % in second
InitStart = 1 ;
Long=zeros(2,1);
Lat=zeros(2,1);
while (InitStart ~= 0 )
[pos, vel] = sgp4(DeltaMinute, satdata);
[Lat,Long,Alt] = ecef2lla(pos(1)*1000,pos(2)*1000,pos(3)*1000) ;
Lat(InitStart) = Lat*rad;
Long(InitStart) = (Long*rad) - 180 ;
Long = Long - (DeltaMinute*60*ERS) ;
if Long < -180
Long = Long + 360 ;
end
if Long > 180
Long = Long - 360 ;
end
Alt = Alt/1000 ; % convert m to km
DeltaMinute = DeltaMinute + (step/60) ;
hold on
h1 = plot(Long,Lat, 'or','MarkerFaceColor','r');
axis equal
drawnow
title([datestr(tnow) ' - [Time step : ' num2str(step) ' sec]'],'Fontsize', 14);
tnow = tnow+(step/(3600*24)) ;
%pause(0.5)
%xlim([-180 180])
%ylim([-90 90])
%set(h1,'Visible','off')
%InitStart = 0;
end
