mission analysis with the orbit propagator block -凯发k8网页登录

define mission parameters and satellite 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(2019, 1, 4, 12, 0, 0);
mission.duration  = hours(6);

specify keplerian orbital elements for the satellite(s) at the mission.startdate.

mission.satellite.semimajoraxis  = 6786233.13; % meters
mission.satellite.eccentricity   = 0.0010537;
mission.satellite.inclination    = 51.7519;    % deg
mission.satellite.raan           = 95.2562;    % deg
mission.satellite.argofperiapsis = 93.4872;    % deg
mission.satellite.trueanomaly    = 202.9234;   % deg

specify the latitude and longitude of a ground station to use in access analysis below.

mission.groundstation.latitude  =42;  % deg
mission.groundstation.longitude =-71; % deg

open and configure the orbit propagation model

open the included simulink model. this model contains an orbit propagator block connected to output ports. the orbit propagator block supports vectorization. this allows you to model multiple satellites in a single block by specifying arrays of initial conditions in the block parameters window or using set_param. the model also includes a "mission analysis and visualization" section that contains a dashboard callback button. when clicked, this button runs the model, creates a new satellitescenario object in the global base workspace containing the satellite or constellation defined in the orbit propagator block, and opens a satellite scenario viewer window for the new scenario. to view the source code for this action, double click the callback button. the "mission analysis and visualization" section is a standalone workflow to create a new satellitescenario object and is not used as part of this example.

mission.mdl = "orbitpropagatorblockexamplemodel";
open_system(mission.mdl);
snapshotmodel(mission.mdl);

define the path to the orbit propagator block in the model.

mission.satellite.blk = mission.mdl   "/orbit propagator";

set satellite initial conditions. to assign the keplerian orbital element set defined in the previous section, use set_param.

set_param(mission.satellite.blk, ...
    "startdate",      num2str(juliandate(mission.startdate)), ...
    "stateformatnum", "orbital elements", ...
    "orbittype",      "keplerian", ...
    "semimajoraxis",  "mission.satellite.semimajoraxis", ...
    "eccentricity",   "mission.satellite.eccentricity", ...
    "inclination",    "mission.satellite.inclination", ...
    "raan",           "mission.satellite.raan", ...
    "argperiapsis",   "mission.satellite.argofperiapsis", ...
    "trueanomaly",    "mission.satellite.trueanomaly");

set the position and velocity output ports of the block to use the earth-centered earth-fixed frame, which is the international terrestrial reference frame (itrf).

set_param(mission.satellite.blk, ...
    "centralbody",  "earth", ...
    "outportframe", "fixed-frame");

configure the propagator. this example uses a numerical propagator for higher position accuracy. use numerical propagators to model earth gravitational potential using the equation for universal gravitation ("pt-mass"), a second order zonal harmonic model ("oblate ellipsoid (j2)"), or a spherical harmonic model ("spherical harmonics"). spherical harmonics are the most accurate, but trade accuracy for speed. for increased accuracy, you can also specify whether to use earth orientation parameters (eop's) in the internal transformations between inertial (icrf) and fixed (itrf) coordinate systems.

set_param(mission.satellite.blk, ...
    "propagator",   "numerical (high precision)", ...
    "gravitymodel", "spherical harmonics", ...
    "earthsh",      "egm2008", ... % earth spherical harmonic potential model
    "shdegree",     "120", ... % spherical harmonic model degree and order
    "useeops",      "on", ... % use eop's in eci to ecef transformations
    "eopfile",      "aeroiersdata.mat"); % eop data file

apply model-level solver setting using set_param. for best performance and accuracy when using a numerical propagator, use a variable-step solver.

set_param(mission.mdl, ...
    "solvertype", "variable-step", ...
    "solvername", "variablestepauto", ...
    "reltol",     "1e-6", ...
    "abstol",     "1e-7", ...
    "stoptime",   string(seconds(mission.duration)));

save model output port data as a dataset of time series objects.

set_param(mission.mdl, ...
    "saveoutput", "on", ...
    "outputsavename", "yout", ...
    "saveformat", "dataset");

run the model and collect satellite ephemerides

simulate the model. in this example, the orbit propagator block is set to output position and velocity states in the ecef (itrf) coordinate frame.

mission.simoutput = sim(mission.mdl);

extract position and velocity data from the model output data structure.

mission.satellite.timeseriesposecef = mission.simoutput.yout{1}.values;
mission.satellite.timeseriesvelecef = mission.simoutput.yout{2}.values;

set the start data from the mission in the timeseries object.

mission.satellite.timeseriesposecef.timeinfo.startdate = mission.startdate;
mission.satellite.timeseriesvelecef.timeinfo.startdate = mission.startdate;

load the satellite ephemerides into a satellitescenario object

create a satellite scenario object to use during the analysis portion of this example.

scenario = satellitescenario;

add the satellites to the satellite scenario as ecef position and velocity timeseries using the satellite method.

sat = satellite(scenario, mission.satellite.timeseriesposecef, mission.satellite.timeseriesvelecef, ...
    "coordinateframe", "ecef")

sat = 
  satellite with properties:
               name:  satellite
                 id:  1
     conicalsensors:  [1x0 matlabshared.satellitescenario.conicalsensor]
            gimbals:  [1x0 matlabshared.satellitescenario.gimbal]
       transmitters:  [1x0 satcom.satellitescenario.transmitter]
          receivers:  [1x0 satcom.satellitescenario.receiver]
           accesses:  [1x0 matlabshared.satellitescenario.access]
        groundtrack:  [1x1 matlabshared.satellitescenario.groundtrack]
              orbit:  [1x1 matlabshared.satellitescenario.orbit]
    orbitpropagator:  ephemeris
        markercolor:  [0.059 1 1]
         markersize:  6
          showlabel:  true
     labelfontcolor:  [1 1 1]
      labelfontsize:  15

disp(scenario)

  satellitescenario with properties:
         starttime: 04-jan-2019 12:00:00
          stoptime: 04-jan-2019 18:00:00
        sampletime: 60
      autosimulate: 1
        satellites: [1×1 matlabshared.satellitescenario.satellite]
    groundstations: [1×0 matlabshared.satellitescenario.groundstation]
           viewers: [0×0 matlabshared.satellitescenario.viewer]
          autoshow: 1

preview latitude (deg), longitude (deg), and altitude (m) for each satellite. use the states method to query satellite states at each scenario time step.

for idx = numel(sat):-1:1
    % retrieve states in geographic coordinates
    [lladata, ~, llatimestamps] = states(sat(idx), "coordinateframe","geographic");
    % organize state data for each satellite in a seperate timetable
    mission.satellite.llatable{idx} = timetable(llatimestamps', lladata(1,:)', lladata(2,:)', lladata(3,:)',...
        'variablenames', {'lat_deg','lon_deg', 'alt_m'});
    mission.satellite.llatable{idx}
end

ans=361×3 timetable
            time            lat_deg    lon_deg      alt_m   
    ____________________    _______    _______    __________
    04-jan-2019 12:00:00    -44.804    120.35     4.2526e 05
    04-jan-2019 12:01:00    -42.797    124.73     4.2229e 05
    04-jan-2019 12:02:00    -40.626    128.77     4.2393e 05
    04-jan-2019 12:03:00    -38.322    132.53     4.2005e 05
    04-jan-2019 12:04:00    -35.848    136.07     4.2004e 05
    04-jan-2019 12:05:00    -33.289    139.35      4.203e 05
    04-jan-2019 12:06:00    -30.655    142.41      4.187e 05
    04-jan-2019 12:07:00    -27.884    145.34     4.1982e 05
    04-jan-2019 12:08:00    -25.069    148.09     4.1831e 05
    04-jan-2019 12:09:00    -22.234    150.68     4.1404e 05
    04-jan-2019 12:10:00    -19.297    153.19     4.1829e 05
    04-jan-2019 12:11:00    -16.343    155.58     4.1713e 05
    04-jan-2019 12:12:00    -13.388    157.89       4.07e 05
    04-jan-2019 12:13:00    -10.354    160.15      4.104e 05
    04-jan-2019 12:14:00    -7.3077    162.37     4.1291e 05
    04-jan-2019 12:15:00    -4.2622    164.55     4.0487e 05
      ⋮

clear lladata llatimestamps;

display satellite trajectories over the 3d globe

to display the satellite trajectories over earth (wgs84 ellipsoid), use helper function plot3dtrajectory.

mission.colormap = lines(256); % define colormap for satellite trajectories 
mission.colormap(1:3,:) = [];
plot3dtrajectories(mission.satellite, mission.colormap);

display global and regional 2d ground traces

view the global ground trace as a 2d projection using helper function plot2dtrajectories:

plot2dtrajectories(mission.satellite, mission.groundstation, mission.colormap);

view regional ground trace. select the region of interest from the dropdown menu:

plot2dtrajectories(mission.satellite, mission.groundstation, mission.colormap, "usa");

compute satellite to ground station access (line-of-sight visibility)

add the ground station to the satellitescenario object using the groundstation method.

gs = groundstation(scenario, mission.groundstation.latitude, mission.groundstation.longitude, ...
    "minelevationangle", 10, "name", "ground station")

gs = 
  groundstation with properties:
                 name:  ground station
                   id:  2
             latitude:  42 degrees
            longitude:  -71 degrees
             altitude:  0 meters
    minelevationangle:  10 degrees
       conicalsensors:  [1x0 matlabshared.satellitescenario.conicalsensor]
              gimbals:  [1x0 matlabshared.satellitescenario.gimbal]
         transmitters:  [1x0 satcom.satellitescenario.transmitter]
            receivers:  [1x0 satcom.satellitescenario.receiver]
             accesses:  [1x0 matlabshared.satellitescenario.access]
          markercolor:  [1 0.4118 0.1608]
           markersize:  6
            showlabel:  true
       labelfontcolor:  [1 1 1]
        labelfontsize:  15

attach line-of-sight access analyses between all individual satellites and the ground station using the access method.

ac = access(sat, gs);
ac.linecolor = "green";

display access intervals

display access intervals for each satellite as a timetable. use accessstatus and accessintervals satellite methods to interact with access analysis results.

for idx = numel(ac):-1:1
    mission.satellite.accessstatus{idx} = accessstatus(ac(idx));
    mission.satellite.accesstable{idx} = accessintervals(ac(idx));
    % use local function addllatotimetable to add geographic positions and
    % closest approach range to the access intervals timetable
    mission.satellite.accesstable{idx} = addllatotimetable(...
        mission.satellite.accesstable{idx}, mission.satellite.llatable{idx}, mission.groundstation);
end
clear idx;

display access intervals overlaying 2d ground traces of the satellite trajectories using helper function plotaccessintervals.

plotaccessintervals(mission.satellite, mission.groundstation, mission.colormap);

mission.satellite.accesstable{:}

ans=2×8 table
      source            target         intervalnumber         starttime                endtime           duration    lla (deg, deg, m)    closestapproach (m)
    ___________    ________________    ______________    ____________________    ____________________    ________    _________________    ___________________
    "satellite"    "ground station"          1           04-jan-2019 12:44:00    04-jan-2019 12:50:00      360         {6×3 double}           5.0087e 05     
    "satellite"    "ground station"          2           04-jan-2019 14:21:00    04-jan-2019 14:25:00      240         {4×3 double}            1.102e 06     

further analysis

play the satellitescenario object to open and animate the scenario in a satellitescenarioviewer window.

play(scenario);
disp(scenario.viewers(1))

  viewer with properties:
                       name: 'satellite scenario viewer'
                   position: [560 300 800 600]
                    basemap: 'satellite'
    playbackspeedmultiplier: 50
       camerareferenceframe: 'ecef'
                currenttime: 04-jan-2019 12:00:32
                showdetails: 1
                  dimension: '3d'

show the satellite ground track in the viewer.

groundtrack(sat);

references

[1] wertz, james r, david f. everett, and jeffery j. puschell. space mission engineering: the new smad. hawthorne, ca: microcosm press, 2011. print.