Spacecraft, Thruster, Tank, Propagator, CoordinateSystem, etc.

Is it possible to access the orbital parameters (cartesian or keplerian) of solar system bodies (Sun, Earth, Luna, etc) in the GMAT script just like we can do for a spacecraft (like DefaultSC.SMA)?
prince

Posts: 16
Joined: Mon Jul 02, 2012 8:02 pm

I don't believe this is possible. You can report parameters, and these are all listed when you click the Edit button on the ReportFile properties window (after adding a ReportFile to the Output folder in the GUI).

One option here is to write a MATLAB script to use SPICE to read the same DE file that GMAT is using, then pass back the state to GMAT. This forum has some extensive topics describing how to use the MATLAB interface (just search for them).
Joel J. K. Parker
Flight dynamics engineer, GMAT team
http://gmatcentral.org
jjkparker

Posts: 617
Joined: Thu Jan 07, 2010 9:48 pm
Location: NASA Goddard Space Flight Center, Greenbelt, MD

There is a way to do this, although it is not as clean as I would like. We're working on defining features for R2013b and the ability to get access to lower level data, in a more elegant manner, is currently on the list.

Say you want body1 w/r/t to body2. Put a spacecraft at the center of a coordinate system that has body1 as the origin (i.e. coordinate [0,0,0]), and ask for the spacecraft state in a coordinate system centered at body2. Assume you need earth w/r/t sun:

Code: Select all
`Create Spacecraft earthSatearthSat.Epoch.UTCGregorian = '04 Jul 2012 12:00:00.000'earthSat.CoordinateSystem   = EarthMJ2000Eq;earthSat.X = 0;earthSat.Y = 0;earthSat.Z = 0;Create CoordinateSystem SunMJ2000EqSunMJ2000Eq.Axes    = MJ2000EqSunMJ2000Eq.Origin  = Sun;Create ReportFile rfBeginMissionSequenceReport rf earthSat.SunMJ2000Eq.X earthSat.SunMJ2000Eq.Y earthSat.SunMJ2000Eq.Z`
shughes

Posts: 443
Joined: Mon Jun 09, 2008 6:27 pm

thanks that's a nice trick!
prince

Posts: 16
Joined: Mon Jul 02, 2012 8:02 pm

shughes wrote:There is a way to do this, although it is not as clean as I would like. We're working on defining features for R2013b and the ability to get access to lower level data, in a more elegant manner, is currently on the list.

Say you want body1 w/r/t to body2. Put a spacecraft at the center of a coordinate system that has body1 as the origin (i.e. coordinate [0,0,0]), and ask for the spacecraft state in a coordinate system centered at body2. Assume you need earth w/r/t sun:

Code: Select all
`Create Spacecraft earthSatearthSat.Epoch.UTCGregorian = '04 Jul 2012 12:00:00.000'earthSat.CoordinateSystem   = EarthMJ2000Eq;earthSat.X = 0;earthSat.Y = 0;earthSat.Z = 0;Create CoordinateSystem SunMJ2000EqSunMJ2000Eq.Axes    = MJ2000EqSunMJ2000Eq.Origin  = Sun;Create ReportFile rfBeginMissionSequenceReport rf earthSat.SunMJ2000Eq.X earthSat.SunMJ2000Eq.Y earthSat.SunMJ2000Eq.Z`

I tried but it seems like using this way the state of earth is not propagated with time. See the following script
Code: Select all
`Create Spacecraft earthSatearthSat.Epoch.UTCGregorian = '04 Jul 2012 12:00:00.000'earthSat.CoordinateSystem   = EarthMJ2000Eq;earthSat.X = 0;earthSat.Y = 0;earthSat.Z = 0;Create Spacecraft earthSat1earthSat1.Epoch.UTCGregorian = '04 Jul 2012 12:00:00.000'earthSat1.CoordinateSystem   = EarthMJ2000Eq;earthSat1.X = 6678;earthSat1.Y = 0;earthSat1.Z = 0;Create String newrow;GMAT newrow = ' ';Create ForceModel DefaultProp_ForceModel;GMAT DefaultProp_ForceModel.CentralBody = Earth;GMAT DefaultProp_ForceModel.PrimaryBodies = {Earth};GMAT DefaultProp_ForceModel.Drag = None;GMAT DefaultProp_ForceModel.SRP = Off;GMAT DefaultProp_ForceModel.RelativisticCorrection = Off;GMAT DefaultProp_ForceModel.ErrorControl = RSSStep;GMAT DefaultProp_ForceModel.GravityField.Earth.Degree = 4;GMAT DefaultProp_ForceModel.GravityField.Earth.Order = 4;GMAT DefaultProp_ForceModel.GravityField.Earth.PotentialFile = 'JGM2.cof';GMAT DefaultProp_ForceModel.GravityField.Earth.EarthTideModel = 'None';Create Propagator DefaultProp;GMAT DefaultProp.FM = DefaultProp_ForceModel;GMAT DefaultProp.Type = RungeKutta89;GMAT DefaultProp.InitialStepSize = 60;GMAT DefaultProp.Accuracy = 9.999999999999999e-012;GMAT DefaultProp.MinStep = 0.001;GMAT DefaultProp.MaxStep = 2700;GMAT DefaultProp.MaxStepAttempts = 50;GMAT DefaultProp.StopIfAccuracyIsViolated = true;Create CoordinateSystem SunMJ2000EqSunMJ2000Eq.Axes    = MJ2000EqSunMJ2000Eq.Origin  = Sun;Create ReportFile rfBeginMissionSequenceReport rf earthSat.SunMJ2000Eq.X earthSat.SunMJ2000Eq.Y earthSat.SunMJ2000Eq.ZReport rf earthSat1.SunMJ2000Eq.X earthSat1.SunMJ2000Eq.Y earthSat1.SunMJ2000Eq.ZReport rf earthSat1.EarthMJ2000Eq.X earthSat1.EarthMJ2000Eq.Y earthSat1.EarthMJ2000Eq.ZPropagate DefaultProp(earthSat1) {earthSat1.ElapsedSecs = 120000.0};Report rf newrowReport rf earthSat.SunMJ2000Eq.X earthSat.SunMJ2000Eq.Y earthSat.SunMJ2000Eq.ZReport rf earthSat1.SunMJ2000Eq.X earthSat1.SunMJ2000Eq.Y earthSat1.SunMJ2000Eq.ZReport rf earthSat1.EarthMJ2000Eq.X earthSat1.EarthMJ2000Eq.Y earthSat1.EarthMJ2000Eq.Z`

Results
Code: Select all
`earthSat.SunMJ2000Eq.X   earthSat.SunMJ2000Eq.Y   earthSat.SunMJ2000Eq.Z   33526136.91647737        -136111176.19852         -59006945.77960628       earthSat1.SunMJ2000Eq.X   earthSat1.SunMJ2000Eq.Y   earthSat1.SunMJ2000Eq.Z   33532814.91647737         -136111176.19852          -59006945.77960628        earthSat1.EarthMJ2000Eq.X   earthSat1.EarthMJ2000Eq.Y   earthSat1.EarthMJ2000Eq.Z   6678                        0                           0                           newrow                                        earthSat.SunMJ2000Eq.X   earthSat.SunMJ2000Eq.Y   earthSat.SunMJ2000Eq.Z   33526136.91647737        -136111176.19852         -59006945.77960628       earthSat1.SunMJ2000Eq.X   earthSat1.SunMJ2000Eq.Y   earthSat1.SunMJ2000Eq.Z   36947213.61357012         -135369604.859301         -58683818.51481943        earthSat1.EarthMJ2000Eq.X   earthSat1.EarthMJ2000Eq.Y   earthSat1.EarthMJ2000Eq.Z   3270.819287809069           -5335.765853022342          -577.4975532559364          `
prince

Posts: 16
Joined: Mon Jul 02, 2012 8:02 pm

This is because time isn't universal to the script; it's local to each Spacecraft.

At the beginning of the script you have `earthSat` and `earthSat1`, both of which are at 04 Jul 2012 12:00:00.000. Then you're propagating `earthSat1` forward by 120000 seconds, so now it's at 05 Jul 2012 21:20:00.000. But when you run this line the second time:
Code: Select all
`Report rf earthSat.SunMJ2000Eq.X earthSat.SunMJ2000Eq.Y earthSat.SunMJ2000Eq.Z`

`earthSat` is still at 04 Jul 2012 12:00:00.000.

Try something like this (mission sequence only). Notice the new line after `Report rf newrow`:
Code: Select all
`Report rf earthSat.SunMJ2000Eq.X earthSat.SunMJ2000Eq.Y earthSat.SunMJ2000Eq.ZReport rf earthSat1.SunMJ2000Eq.X earthSat1.SunMJ2000Eq.Y earthSat1.SunMJ2000Eq.ZReport rf earthSat1.EarthMJ2000Eq.X earthSat1.EarthMJ2000Eq.Y earthSat1.EarthMJ2000Eq.ZPropagate DefaultProp(earthSat1) {earthSat1.ElapsedSecs = 120000.0};Report rf newrowearthSat.Epoch = earthSat1.UTCGregorianReport rf earthSat.SunMJ2000Eq.X earthSat.SunMJ2000Eq.Y earthSat.SunMJ2000Eq.ZReport rf earthSat1.SunMJ2000Eq.X earthSat1.SunMJ2000Eq.Y earthSat1.SunMJ2000Eq.ZReport rf earthSat1.EarthMJ2000Eq.X earthSat1.EarthMJ2000Eq.Y earthSat1.EarthMJ2000Eq.Z`
Joel J. K. Parker
Flight dynamics engineer, GMAT team
http://gmatcentral.org
jjkparker

Posts: 617
Joined: Thu Jan 07, 2010 9:48 pm
Location: NASA Goddard Space Flight Center, Greenbelt, MD

thanks it works.

Also I saw a minor issue with STM propagation that if we propagate the same states with and without STM, the results are not exactly same. I have attached the script. The error is of the order of 1e-4 when I propagated for 120000 secs.
Attachments
test.script
prince

Posts: 16
Joined: Mon Jul 02, 2012 8:02 pm

The difference is because the adaptive step algorithm has 36 additional states that can affect the step size and so the integrator is taking slightly different steps and arriving at a slightly different answer that still meets the integration tolerance of 1e-11. I tried two things. (1) I reduced the Acurracy setting to 1e-14 and get sub millimeter agreement as shown below. I set to fixed step (DefaultProp_ForceModel.ErrorControl = None;) and see identical results.

Code: Select all
`earthSat1.SunMJ2000Eq.X   earthSat1.SunMJ2000Eq.Y   earthSat1.SunMJ2000Eq.Z   earthSat1.SunMJ2000Eq.VX  36947213.61359383         -135369604.8592888        -58683818.51481693        35.44631554785417          10.04496492114358   earthSat1.EarthMJ2000Eq.X   earthSat1.EarthMJ2000Eq.Y   earthSat1.EarthMJ2000Eq.Z   earthSat1.EarthMJ2000Eq.VX      3270.819311521271           -5335.765840896569          -577.4975507577153          7.044415173766448  earthSat2.SunMJ2000Eq.X   earthSat2.SunMJ2000Eq.Y   earthSat2.SunMJ2000Eq.Z   earthSat2.SunMJ2000Eq.VX   36947213.61359354         -135369604.859289         -58683818.51481707        35.44631554807062          10.0449649207858earthSat2.EarthMJ2000Eq.X   earthSat2.EarthMJ2000Eq.Y   earthSat2.EarthMJ2000Eq.Z   earthSat2.EarthMJ2000Eq.VX    3270.819311230655           -5335.765841041648          -577.4975508946158          7.044415173982896                     `
shughes

Posts: 443
Joined: Mon Jun 09, 2008 6:27 pm

that makes sense. However I am seeing one more behavior that I am not able to understand fully. I propagate 2 satellites using exactly same propagators. The satellites initial conditions however are not exactly same. After the propagation, I assign state (X,Y,Z,VX,Vy,VZ) of 1 satellite to another and immediately subtract the states of 2 satellites in a different coordinate system. The difference is not zero. If I subtract in the same coordinate system then difference is 0 as expected. I am not sure what I am doing wrong. I have attached my test script.
Attachments
test.script
prince

Posts: 16
Joined: Mon Jul 02, 2012 8:02 pm

Here are my first thoughts, but I am not convinced I fully understand this issue yet. For starters, if you use a single propagate command like the code snippet below, which ensures the epochs of the spacecraft are precisely synchronized, then this issue goes away:

Code: Select all
`Propagate DefaultProp(DefaultSC,DefaultSC3) {DefaultSC.ElapsedSecs = 60};`

Here are some thoughts on the behavior for separate propagations: GMAT models time using a double precision real number in a modified Julian date format (less than ideal and we want to change that eventually). I can't prove it yet, but since the synchronized propagation avoids this issue, I believe the epochs of the separately propagated spacecraft are very slightly different due to ~5 microsecond time discretization that occurs when modelling time as a modified Julian date using a double precision real. When you take into account that the Earth moves about 30 km/s and error in time is about 5 microseconds you get 0.5e-5 s * 30 km/s = 1.5 cm. The RSS of difference of the actual data shown below is 1.1 cm. The precision of a double may come into play here in more ways though because for these values there are 9 digits to the left of the decimal leaving about 6 useful digits to the right which is what we see in the difference.

Code: Select all
`DefaultSC.TAIModJulian    DefaultSC.SunMJ2000Ec.X    DefaultSC.SunMJ2000Ec.Y    DefaultSC.SunMJ2000Ec.Z   21545.00069444444        -26494694.3116229           144697736.8044981          457.9547021538019                DefaultSC3.TAIModJulian   DefaultSC3.SunMJ2000Ec.X   DefaultSC3.SunMJ2000Ec.Y   DefaultSC3.SunMJ2000Ec.Z  21545.00069444444         -26494694.31161353         144697736.8044998          457.9547021538019                X2(1,1)                   X2(2,1)                    X2(3,1)                    X2(4,1)                              -9.365379810333252e-06   -1.728534698486328e-06      0  `
shughes

Posts: 443
Joined: Mon Jun 09, 2008 6:27 pm

Next