modeling satellite constellations using ephemeris data -凯发k8网页登录
this example demonstrates how to add time-stamped ephemeris data for a constellation of 24 satellites (similar to esa galileo gnss constellation) to a satellite scenario for access analysis. the example uses data generated by the aerospace blockset orbit propagator
block. for more information, see the aerospace blockset example constellation modeling with the orbit propagator block.
the satellitescenario
object supports loading previously generated, time-stamped satellite ephemeris data into a scenario from a timeseries
or timetable
object. an ephemeris is a table containing position (and optionally velocity) state information of a satellite during a given period of time. ephemeris data used to add satellites to the scenario object is interpolated via the makima
interpolation method to align with the scenario time steps. this allows you to incorporate data generated by a simulink model into either a new or existing satellitescenario.
define mission parameters and constellation initial conditions
specify a start date and duration for the mission. this example uses matlab structures to organize mission data. these structures make accessing data later in the example more intuitive. they also help declutter the global base workspace.
mission.startdate = datetime(2020, 11, 30, 22, 23, 24); mission.duration = hours(24);
the constellation in this example is a walker-delta constellation modeled similar to galileo, the european gnss (global navigation satellite system) constellation. the constellation consists of 24 satellites in medium earth orbit (meo). the satellites' keplerian orbital elements at the mission start date epoch are:
mission.constellationdefinition = table( ... 29599.8e3 * ones(24,1), ... % semi-major axis (m) 0.0005 * ones(24,1), ... % eccentricity 56 * ones(24,1), ... % inclination (deg) 350 * ones(24,1), ... % right ascension of the ascending node (deg) sort(repmat([0 120 240], 1,8))', ... % argument of periapsis (deg) [0:45:315, 15:45:330, 30:45:345]', ... % true anomaly (deg) 'variablenames', ["a (m)", "e", "i (deg)", ... "ω (deg)", "ω (deg)", "ν (deg)"]); mission.constellationdefinition
ans=24×6 table
a (m) e i (deg) ω (deg) ω (deg) ν (deg)
________ ______ _______ _______ _______ _______
2.96e 07 0.0005 56 350 0 0
2.96e 07 0.0005 56 350 0 45
2.96e 07 0.0005 56 350 0 90
2.96e 07 0.0005 56 350 0 135
2.96e 07 0.0005 56 350 0 180
2.96e 07 0.0005 56 350 0 225
2.96e 07 0.0005 56 350 0 270
2.96e 07 0.0005 56 350 0 315
2.96e 07 0.0005 56 350 120 15
2.96e 07 0.0005 56 350 120 60
2.96e 07 0.0005 56 350 120 105
2.96e 07 0.0005 56 350 120 150
2.96e 07 0.0005 56 350 120 195
2.96e 07 0.0005 56 350 120 240
2.96e 07 0.0005 56 350 120 285
2.96e 07 0.0005 56 350 120 330
⋮
load ephemeris timeseries data
the timeseries objects contain position and velocity data for all 24 satellites in the constellation. the data is referenced in the international terrestrial reference frame (itrf), which is an earth-centered earth-fixed (ecef) coordinate system. the data was generated using the aerospace blockset orbit propagator
block. for more information, see the aerospace blockset example constellation modeling with the orbit propagator block.
mission.ephemeris = load("satellitescenarioephemerisdata.mat", "timeseriespositrf", "timeseriesvelitrf"); mission.ephemeris.timeseriespositrf
timeseries common properties: name: '' time: [57x1 double] timeinfo: [1x1 tsdata.timemetadata] data: [24x3x57 double] datainfo: [1x1 tsdata.datametadata] more properties, methods
mission.ephemeris.timeseriesvelitrf
timeseries common properties: name: '' time: [57x1 double] timeinfo: [1x1 tsdata.timemetadata] data: [24x3x57 double] datainfo: [1x1 tsdata.datametadata] more properties, methods
load the satellite ephemerides into a satellitescenario o
bject
create a satellite scenario object for the analysis.
scenario = satellitescenario(mission.startdate, mission.startdate hours(24), 60);
use the satellite
method to add all 24 satellites to the satellite scenario from the ecef position and velocity timeseries objects. this example uses position and velocity information; however satellites can also be added from position data only and velocity states are then estimated. available coordinate frames for name-value pair coordinateframe
are "ecef"
, "inertial"
, and "geographic"
. if the timeseries object contains a value for ts.timeinfo.startdate
, the method uses that value as the epoch for the timeseries object. if no startdate
is defined, the method uses the scenario start date by default.
sat = satellite(scenario, mission.ephemeris.timeseriespositrf, mission.ephemeris.timeseriesvelitrf, ... coordinateframe="ecef", name="galileo " (1:24))
sat = 1x24 satellite array with properties: name id conicalsensors gimbals transmitters receivers accesses groundtrack orbit orbitpropagator markercolor markersize showlabel labelfontcolor labelfontsize
disp(scenario)
satellitescenario with properties: starttime: 30-nov-2020 22:23:24 stoptime: 01-dec-2020 22:23:24 sampletime: 60 autosimulate: 1 satellites: [1×24 matlabshared.satellitescenario.satellite] groundstations: [1×0 matlabshared.satellitescenario.groundstation] viewers: [0×0 matlabshared.satellitescenario.viewer] autoshow: 1
alternatively, satellites can also be added as ephemerides to the satellite scenario as a matlab timetable
, table
, or tscollection
. for example, a timetable
containing the first 3 satellites of the position timeseries
object in the previous section, formatted for use with satellitescenario
objects is shown below.
satellites are represented by variables (column headers).
each row contains a position vector associated with the row's
time
property.
timetable(... datetime(getabstime(mission.ephemeris.timeseriespositrf), locale="en_us"), ... squeeze(mission.ephemeris.timeseriespositrf.data(1,:,:))', ... squeeze(mission.ephemeris.timeseriespositrf.data(2,:,:))', ... squeeze(mission.ephemeris.timeseriespositrf.data(3,:,:))',... variablenames=["satellite_1", "satellite_2", "satellite_3"])
ans=57×3 timetable
time satellite_1 satellite_2 satellite_3
____________________ ________________________________________ ________________________________________ ________________________________________
30-nov-2020 22:23:24 1.8249e 07 -2.2904e 07 -4.2009e 06 2.3678e 07 -1.075e 07 1.4119e 07 1.5239e 07 7.7076e 06 2.4177e 07
30-nov-2020 22:23:38 1.8252e 07 -2.2909e 07 -4.1563e 06 2.3662e 07 -1.0735e 07 1.4156e 07 1.5214e 07 7.7334e 06 2.4184e 07
30-nov-2020 22:24:53 1.8268e 07 -2.2937e 07 -3.933e 06 2.3584e 07 -1.0663e 07 1.434e 07 1.5088e 07 7.8627e 06 2.4222e 07
30-nov-2020 22:31:05 1.8326e 07 -2.3055e 07 -2.8121e 06 2.3185e 07 -1.028e 07 1.5243e 07 1.4466e 07 8.5229e 06 2.4378e 07
30-nov-2020 22:48:39 1.8326e 07 -2.3223e 07 3.9182e 05 2.2005e 07 -8.9966e 06 1.7621e 07 1.2798e 07 1.0506e 07 2.4539e 07
30-nov-2020 23:08:30 1.8076e 07 -2.3078e 07 3.9992e 06 2.0643e 07 -7.2057e 06 1.9943e 07 1.1124e 07 1.2894e 07 2.4217e 07
30-nov-2020 23:28:27 1.7624e 07 -2.2538e 07 7.5358e 06 1.9321e 07 -5.0678e 06 2.1838e 07 9.7076e 06 1.5379e 07 2.3362e 07
30-nov-2020 23:50:59 1.6968e 07 -2.1428e 07 1.1328e 07 1.7977e 07 -2.3021e 06 2.34e 07 8.4636e 06 1.8183e 07 2.1782e 07
01-dec-2020 00:14:27 1.6244e 07 -1.9712e 07 1.4937e 07 1.6838e 07 8.7771e 05 2.4329e 07 7.5789e 06 2.0966e 07 1.9489e 07
01-dec-2020 00:38:42 1.5585e 07 -1.7375e 07 1.8189e 07 1.6017e 07 4.355e 06 2.4512e 07 7.0779e 06 2.3551e 07 1.6498e 07
01-dec-2020 01:04:35 1.5124e 07 -1.4345e 07 2.1006e 07 1.5585e 07 8.1065e 06 2.383e 07 6.9314e 06 2.5831e 07 1.2718e 07
01-dec-2020 01:31:17 1.5035e 07 -1.079e 07 2.3096e 07 1.562e 07 1.1816e 07 2.2205e 07 7.0715e 06 2.7527e 07 8.3282e 06
01-dec-2020 01:58:58 1.5443e 07 -6.8501e 06 2.4303e 07 1.6102e 07 1.5274e 07 1.9601e 07 7.348e 06 2.8484e 07 3.4363e 06
01-dec-2020 02:27:08 1.6406e 07 -2.8152e 06 2.4478e 07 1.6925e 07 1.8197e 07 1.6103e 07 7.5521e 06 2.8587e 07 -1.6897e 06
01-dec-2020 02:55:18 1.7869e 07 1.001e 06 2.3582e 07 1.7894e 07 2.0376e 07 1.1901e 07 7.4614e 06 2.7856e 07 -6.7427e 06
01-dec-2020 03:23:29 1.9711e 07 4.381e 06 2.1653e 07 1.8787e 07 2.1739e 07 7.1754e 06 6.8858e 06 2.6405e 07 -1.1504e 07
⋮
set graphical properties on the satellites
set satellite in the same orbital plane to have the same orbit color.
set(sat(1:8), markercolor="#ff6929"); set(sat(9:16), markercolor="#139fff"); set(sat(17:24), markercolor="#64d413"); orbit = [sat(:).orbit]; set(orbit(1:8), linecolor="#ff6929"); set(orbit(9:16), linecolor="#139fff"); set(orbit(17:24), linecolor="#64d413");
add ground stations to scenario
to provide accurate positioning data, a location on earth must have access to at least 4 satellites in the constellation at any given time. in this example, use three locations to compare total constellation access over the 1 day analysis window to different regions of earth:
natick, massachusetts, usa (42.30048°, -71.34908°)
münchen, germany (48.23206°, 11.68445°)
bangalore, india (12.94448°, 77.69256°)
gsus = groundstation(scenario, 42.30048, -71.34908, ... minelevationangle=10, name="natick"); gsus.markercolor = "red"; gsde = groundstation(scenario, 48.23206, 11.68445, ... minelevationangle=10, name="munchen"); gsde.markercolor = "red"; gsin = groundstation(scenario, 12.94448, 77.69256, ... minelevationangle=10, name="bangalore"); gsin.markercolor = "red"; figure geoscatter([gsus.latitude gsde.latitude gsin.latitude], ... [gsus.longitude gsde.longitude gsin.longitude], "red", "filled") geolimits([-75 75], [-180 180]) title("ground stations")
compute ground station to satellite access (line-of-sight visibility)
calculate line-of-sight access between the ground stations and each individual satellite using the access
method.
accessus = access(gsus, sat); accessde = access(gsde, sat); accessin = access(gsin, sat);
set access colors to match orbital plane colors assigned earlier in the example.
set(accessus, linewidth="1"); set(accessus(1:8), linecolor="#ff6929"); set(accessus(9:16), linecolor="#139fff"); set(accessus(17:24), linecolor="#64d413"); set(accessde, linewidth="1"); set(accessde(1:8), linecolor="#ff6929"); set(accessde(9:16), linecolor="#139fff"); set(accessde(17:24), linecolor="#64d413"); set(accessin, linewidth="1"); set(accessin(1:8), linecolor="#ff6929"); set(accessin(9:16), linecolor="#139fff"); set(accessin(17:24), linecolor="#64d413");
view the full access table between each ground station and all satellites in the constellation as tables. sort the access intervals by interval start time. satellites added from ephemeris data do not display values for startorbit
and endorbit
.
intervalsus = accessintervals(accessus); intervalsus = sortrows(intervalsus, "starttime", "ascend")
intervalsus=40×8 table
source target intervalnumber starttime endtime duration startorbit endorbit
________ ____________ ______________ ____________________ ____________________ ________ __________ ________
"natick" "galileo 1" 1 30-nov-2020 22:23:24 01-dec-2020 04:04:24 20460 nan nan
"natick" "galileo 2" 1 30-nov-2020 22:23:24 01-dec-2020 01:24:24 10860 nan nan
"natick" "galileo 3" 1 30-nov-2020 22:23:24 30-nov-2020 22:57:24 2040 nan nan
"natick" "galileo 12" 1 30-nov-2020 22:23:24 01-dec-2020 00:00:24 5820 nan nan
"natick" "galileo 13" 1 30-nov-2020 22:23:24 30-nov-2020 23:05:24 2520 nan nan
"natick" "galileo 18" 1 30-nov-2020 22:23:24 01-dec-2020 04:00:24 20220 nan nan
"natick" "galileo 19" 1 30-nov-2020 22:23:24 01-dec-2020 01:42:24 11940 nan nan
"natick" "galileo 20" 1 30-nov-2020 22:23:24 30-nov-2020 22:46:24 1380 nan nan
"natick" "galileo 11" 1 30-nov-2020 22:25:24 01-dec-2020 00:18:24 6780 nan nan
"natick" "galileo 17" 1 30-nov-2020 22:50:24 01-dec-2020 05:50:24 25200 nan nan
"natick" "galileo 8" 1 30-nov-2020 23:20:24 01-dec-2020 07:09:24 28140 nan nan
"natick" "galileo 7" 1 01-dec-2020 01:26:24 01-dec-2020 10:00:24 30840 nan nan
"natick" "galileo 24" 1 01-dec-2020 01:40:24 01-dec-2020 07:12:24 19920 nan nan
"natick" "galileo 14" 1 01-dec-2020 03:56:24 01-dec-2020 07:15:24 11940 nan nan
"natick" "galileo 6" 1 01-dec-2020 04:05:24 01-dec-2020 12:14:24 29340 nan nan
"natick" "galileo 23" 1 01-dec-2020 04:10:24 01-dec-2020 08:03:24 13980 nan nan
⋮
intervalsde = accessintervals(accessde); intervalsde = sortrows(intervalsde, "starttime", "ascend")
intervalsde=40×8 table
source target intervalnumber starttime endtime duration startorbit endorbit
_________ ____________ ______________ ____________________ ____________________ ________ __________ ________
"munchen" "galileo 2" 1 30-nov-2020 22:23:24 01-dec-2020 04:34:24 22260 nan nan
"munchen" "galileo 3" 1 30-nov-2020 22:23:24 01-dec-2020 01:58:24 12900 nan nan
"munchen" "galileo 4" 1 30-nov-2020 22:23:24 30-nov-2020 23:05:24 2520 nan nan
"munchen" "galileo 10" 1 30-nov-2020 22:23:24 30-nov-2020 23:58:24 5700 nan nan
"munchen" "galileo 19" 1 30-nov-2020 22:23:24 01-dec-2020 01:36:24 11580 nan nan
"munchen" "galileo 20" 1 30-nov-2020 22:23:24 01-dec-2020 00:15:24 6720 nan nan
"munchen" "galileo 21" 1 30-nov-2020 22:23:24 30-nov-2020 22:28:24 300 nan nan
"munchen" "galileo 9" 1 30-nov-2020 22:34:24 01-dec-2020 02:22:24 13680 nan nan
"munchen" "galileo 18" 1 30-nov-2020 22:41:24 01-dec-2020 02:31:24 13800 nan nan
"munchen" "galileo 1" 1 30-nov-2020 23:05:24 01-dec-2020 06:42:24 27420 nan nan
"munchen" "galileo 16" 1 30-nov-2020 23:29:24 01-dec-2020 04:47:24 19080 nan nan
"munchen" "galileo 15" 1 01-dec-2020 00:50:24 01-dec-2020 07:27:24 23820 nan nan
"munchen" "galileo 17" 1 01-dec-2020 01:05:24 01-dec-2020 03:00:24 6900 nan nan
"munchen" "galileo 8" 1 01-dec-2020 01:57:24 01-dec-2020 08:25:24 23280 nan nan
"munchen" "galileo 14" 1 01-dec-2020 02:36:24 01-dec-2020 10:19:24 27780 nan nan
"munchen" "galileo 7" 1 01-dec-2020 04:35:24 01-dec-2020 09:43:24 18480 nan nan
⋮
intervalsin = accessintervals(accessin); intervalsin = sortrows(intervalsin, "starttime", "ascend")
intervalsin=31×8 table
source target intervalnumber starttime endtime duration startorbit endorbit
___________ ____________ ______________ ____________________ ____________________ ________ __________ ________
"bangalore" "galileo 3" 1 30-nov-2020 22:23:24 01-dec-2020 05:12:24 24540 nan nan
"bangalore" "galileo 4" 1 30-nov-2020 22:23:24 01-dec-2020 02:59:24 16560 nan nan
"bangalore" "galileo 5" 1 30-nov-2020 22:23:24 01-dec-2020 00:22:24 7140 nan nan
"bangalore" "galileo 9" 1 30-nov-2020 22:23:24 01-dec-2020 03:37:24 18840 nan nan
"bangalore" "galileo 10" 1 30-nov-2020 22:23:24 01-dec-2020 00:09:24 6360 nan nan
"bangalore" "galileo 16" 1 30-nov-2020 22:23:24 01-dec-2020 08:44:24 37260 nan nan
"bangalore" "galileo 21" 1 30-nov-2020 22:23:24 30-nov-2020 23:25:24 3720 nan nan
"bangalore" "galileo 22" 1 30-nov-2020 22:23:24 30-nov-2020 22:58:24 2100 nan nan
"bangalore" "galileo 15" 1 01-dec-2020 00:17:24 01-dec-2020 11:16:24 39540 nan nan
"bangalore" "galileo 2" 1 01-dec-2020 00:25:24 01-dec-2020 07:10:24 24300 nan nan
"bangalore" "galileo 22" 2 01-dec-2020 00:48:24 01-dec-2020 05:50:24 18120 nan nan
"bangalore" "galileo 21" 2 01-dec-2020 01:32:24 01-dec-2020 08:29:24 25020 nan nan
"bangalore" "galileo 1" 1 01-dec-2020 03:06:24 01-dec-2020 07:17:24 15060 nan nan
"bangalore" "galileo 20" 1 01-dec-2020 03:36:24 01-dec-2020 12:38:24 32520 nan nan
"bangalore" "galileo 14" 1 01-dec-2020 05:48:24 01-dec-2020 13:29:24 27660 nan nan
"bangalore" "galileo 19" 1 01-dec-2020 05:53:24 01-dec-2020 17:06:24 40380 nan nan
⋮
view the satellite scenario
open a 3-d viewer window of the scenario. the viewer window contains all 24 satellites and the three ground stations defined earlier in this example. a line is drawn between each ground station and satellite during their corresponding access intervals. hide the details of the satellites and ground stations by setting the showdetails
name-value pair to false
. show satellite orbits and labels for the ground station locations.
viewer3d = satellitescenarioviewer(scenario, showdetails=false); show(sat.orbit); gsus.showlabel = true; gsus.labelfontsize = 11; gsde.showlabel = true; gsde.labelfontsize = 11; gsin.showlabel = true; gsin.labelfontsize = 11;
compare access between ground stations
calculate access status between each satellite and ground station using the accessstatus method. each row of the output array corresponds with a satellite in the constellation. each column corresponds with time steps in the scenario. a value of true
indicates that the satellite can access the aircraft at that specific time sample. the second output of accessstatus
contains the time steps of the scenario. plot cumulative access for each ground station over the one day analysis window.
[statusus, timesteps] = accessstatus(accessus); statusde = accessstatus(accessde); statusin = accessstatus(accessin); % sum cumulative access at each timestep statusus = sum(statusus, 1); statusde = sum(statusde, 1); statusin = sum(statusin, 1); subplot(3,1,1); stairs(timesteps, statusus); title("natick to galileo") ylabel("# of satellites") subplot(3,1,2); stairs(timesteps, statusde); title("münchen to galileo") ylabel("# of satellites") subplot(3,1,3); stairs(timesteps, statusin); title("bangalore to galileo") ylabel("# of satellites")
collect access interval metrics for each ground station in a table for comparison.
statustable = [table(height(intervalsus), height(intervalsde), height(intervalsin)); ... table(sum(intervalsus.duration)/3600, sum(intervalsde.duration)/3600, sum(intervalsin.duration)/3600); ... table(mean(intervalsus.duration/60), mean(intervalsde.duration/60), mean(intervalsin.duration/60)); ... table(mean(statusus, 2), mean(statusde, 2), mean(statusin, 2)); ... table(min(statusus), min(statusde), min(statusin)); ... table(max(statusus), max(statusde), max(statusin))]; statustable.properties.variablenames = ["natick", "münchen", "bangalore"]; statustable.properties.rownames = ["total # of intervals", "total interval time (hrs)",... "mean interval length (min)", "mean # of satellites in view", ... "min # of satellites in view", "max # of satellites in view"]; statustable
statustable=6×3 table
natick münchen bangalore
______ _______ _________
total # of intervals 40 40 31
total interval time (hrs) 167.88 169.95 180.42
mean interval length (min) 251.82 254.93 349.19
mean # of satellites in view 7.018 7.1041 7.5337
min # of satellites in view 5 5 5
max # of satellites in view 9 10 9
walker-delta constellations like galileo are evenly distributed across longitudes. natick and münchen are located at similar latitudes, and therefore have very similar access characteristics with respect to the constellation. bangalore is at a latitude closer to the equator. despite having a lower number of individual access intervals, it has the highest average number of satellites in view, the highest overall interval time, and the longest average interval duration (by about 95 minutes). all locations always have at least 4 satellites in view, as is required for gnss trilateration.
references
[1] wertz, james r, david f. everett, and jeffery j. puschell. space mission engineering: the new smad. hawthorne, ca: microcosm press, 2011. print.
[2] the european space agency: galileo facts and figures. https://www.esa.int/applications/navigation/galileo/facts_and_figures
see also
objects
satellitescenario
| | | |
functions
- | | |