Need help propogating a near-Earth asteroid

Dynamics and environment models, spacecraft model, solver algorithms, etc.

Need help propogating a near-Earth asteroid

Postby timmyd7777 » Sun Jan 25, 2015 7:33 pm

Hello GMAT folks,

I'm trying to use GMAT to model the path of a near-Earth-passing asteroid, and failing. I've set up a Sun-centered Mj2000 ecliptic coordinate system for performing all the calculations and visualizaton.

I've defined a "spacecraft" called "Asteroid357439", using Keplerian orbital elements from JPL HORIZONS. I've set up a propogator with Sun as the central body, planets Mercury thru Neptune as point masses, using RK89 integration and step size of 60 to 3600 secs (default 600 sec).

And then I run the propogator for a year. The problem is that the asteroid's trajectory is completely wrong. I don't think GMAT is broken. Rather I suspect that I'm entering the initial keplerian elements wrong somehow.

To debug this, here are JPL's original elements for asteroid 357439. Note these are heliocentric, not geocentric!

EPOCH= 2456230.5 ! 2012-Oct-30.00 (CT) Residual RMS= .30844
EC= .4031766150000584 QR= .8964725585602467 TP= 2456339.133332551
OM= 126.7275866134088 W= 311.2334402904601 IN= 23.7426510647147
A= 1.502073445999999 MA= 301.8391930397593 ADIST= 2.107674333439752
PER= 1.84096 N= .5353863679999999 ANGMOM= .019293282
DAN= .99381 DDN= 1.71318 L= 80.48398880000001
B= -17.6252182 MOID= .00862029999999999TP= 2013-Feb-15.6333325510

A= 1.502073445999999 is the semimajor axis in AU - easy to translate that to km.
EC= .4031766150000584 is the eccentricity, no translation needed.
IN= 23.7426510647147 is the inclination of the orbit plane to the MJ2000 ecliptic plane in degrees.
W= 311.2334402904601 is the argument of perihelion in degrees, no translation needed.
OM= 126.7275866134088 is the longitude of the ascending node in degrees.
MA= 301.8391930397593 is the mean anomaly in degrees.

Translating this to GMAT, the "spacecraft" starting elements I've entered are as follows:

Create Spacecraft Asteroid357439;
GMAT Asteroid357439.DateFormat = TDBGregorian;
GMAT Asteroid357439.Epoch = '30 Oct 2012 12:00:00.000';
GMAT Asteroid357439.CoordinateSystem = SunMJ2000Ec;
GMAT Asteroid357439.DisplayStateType = Keplerian;
GMAT Asteroid357439.SMA = 224706989.6072336;
GMAT Asteroid357439.ECC = 0.4031766150000588;
GMAT Asteroid357439.INC = 23.7426510647147;
GMAT Asteroid357439.RAAN = 126.7275866134088;
GMAT Asteroid357439.AOP = 311.2334402904602;
GMAT Asteroid357439.TA = 301.8391930397593;

I've translated the semimajor axis (SMA) from AU to km, as noted above.

I'm assuming that GMAT considers RAAN (Right Ascension of Ascending Node) to be the same thing as longitude of ascending node, since the CoordinateSystem is SunMJ200Ec. Is that correct?

I've also entered the mean anomaly (301.8391930397593) as the true anomaly (TA). I know this is wrong. Is there some easy way to specify a spacecraft starting orbit in GMAT using mean anomaly instead of true anomaly? Converting back-and-forth is a PITA.

But regardless of this known error, I would not expect the asteroid's orbital eccentricity or semimajor axis to change much after a month of propogation. But here is what I see in the command summary after running the propogator for 30 days simulated time:

Keplerian State
--------------------------------
SMA = 320781950.93122 km
ECC = 0.6429275692153
INC = 23.742544187677 deg
RAAN = 126.72702374823 deg
AOP = 319.97582752493 deg
TA = 111.30720566592 deg
MA = 34.307288793120 deg
EA = 68.606043083868 deg

This just doesn't make any sense to me. The orbital elements for a heliocentric body just don't change that fast. I've clearly set something up wrong ... but what?

My entire GMAT script follows. Thanks in advance,

-Tim

%General Mission Analysis Tool(GMAT) Script
%Created: 2015-01-25 11:06:03

%----------------------------------------
%---------- Spacecraft
%----------------------------------------

Create Spacecraft Asteroid357439;
GMAT Asteroid357439.DateFormat = TDBGregorian;
GMAT Asteroid357439.Epoch = '30 Oct 2012 12:00:00.000';
GMAT Asteroid357439.CoordinateSystem = SunMJ2000Ec;
GMAT Asteroid357439.DisplayStateType = Keplerian;
GMAT Asteroid357439.SMA = 224706989.6072336;
GMAT Asteroid357439.ECC = 0.4031766150000588;
GMAT Asteroid357439.INC = 23.7426510647147;
GMAT Asteroid357439.RAAN = 126.7275866134088;
GMAT Asteroid357439.AOP = 311.2334402904602;
GMAT Asteroid357439.TA = 105.5723922486322;
GMAT Asteroid357439.DryMass = 850;
GMAT Asteroid357439.Cd = 2.2;
GMAT Asteroid357439.Cr = 1.8;
GMAT Asteroid357439.DragArea = 15;
GMAT Asteroid357439.SRPArea = 1;
GMAT Asteroid357439.NAIFId = -123456789;
GMAT Asteroid357439.NAIFIdReferenceFrame = -123456789;
GMAT Asteroid357439.OrbitColor = Red;
GMAT Asteroid357439.TargetColor = Teal;
GMAT Asteroid357439.Id = 'SatId';
GMAT Asteroid357439.Attitude = CoordinateSystemFixed;
GMAT Asteroid357439.SPADSRPScaleFactor = 1;
GMAT Asteroid357439.ModelFile = '../data/vehicle/models/aura.3ds';
GMAT Asteroid357439.ModelOffsetX = 0;
GMAT Asteroid357439.ModelOffsetY = 0;
GMAT Asteroid357439.ModelOffsetZ = 0;
GMAT Asteroid357439.ModelRotationX = 0;
GMAT Asteroid357439.ModelRotationY = 0;
GMAT Asteroid357439.ModelRotationZ = 0;
GMAT Asteroid357439.ModelScale = 3;
GMAT Asteroid357439.AttitudeDisplayStateType = 'Quaternion';
GMAT Asteroid357439.AttitudeRateDisplayStateType = 'AngularVelocity';
GMAT Asteroid357439.AttitudeCoordinateSystem = SunMJ2000Ec;
GMAT Asteroid357439.EulerAngleSequence = '321';

%----------------------------------------
%---------- ForceModels
%----------------------------------------

Create ForceModel AsteroidPropogator_ForceModel;
GMAT AsteroidPropogator_ForceModel.CentralBody = Sun;
GMAT AsteroidPropogator_ForceModel.PointMasses = {Earth, Jupiter, Luna, Mars, Mercury, Neptune, Saturn, Uranus, Venus};
GMAT AsteroidPropogator_ForceModel.Drag = None;
GMAT AsteroidPropogator_ForceModel.SRP = Off;
GMAT AsteroidPropogator_ForceModel.RelativisticCorrection = On;
GMAT AsteroidPropogator_ForceModel.ErrorControl = RSSStep;

%----------------------------------------
%---------- Propagators
%----------------------------------------

Create Propagator AsteroidPropogator;
GMAT AsteroidPropogator.FM = AsteroidPropogator_ForceModel;
GMAT AsteroidPropogator.Type = RungeKutta89;
GMAT AsteroidPropogator.InitialStepSize = 600;
GMAT AsteroidPropogator.Accuracy = 9.999999999999999e-012;
GMAT AsteroidPropogator.MinStep = 60;
GMAT AsteroidPropogator.MaxStep = 3600;
GMAT AsteroidPropogator.MaxStepAttempts = 50;
GMAT AsteroidPropogator.StopIfAccuracyIsViolated = true;

%----------------------------------------
%---------- Burns
%----------------------------------------

Create ImpulsiveBurn DefaultIB;
GMAT DefaultIB.CoordinateSystem = Local;
GMAT DefaultIB.Origin = Earth;
GMAT DefaultIB.Axes = VNB;
GMAT DefaultIB.Element1 = 0;
GMAT DefaultIB.Element2 = 0;
GMAT DefaultIB.Element3 = 0;
GMAT DefaultIB.DecrementMass = false;
GMAT DefaultIB.Isp = 300;
GMAT DefaultIB.GravitationalAccel = 9.810000000000001;

%----------------------------------------
%---------- Coordinate Systems
%----------------------------------------

Create CoordinateSystem SunMJ2000Ec;
GMAT SunMJ2000Ec.Origin = Sun;
GMAT SunMJ2000Ec.Axes = MJ2000Ec;

%----------------------------------------
%---------- Subscribers
%----------------------------------------

Create OrbitView DefaultOrbitView;
GMAT DefaultOrbitView.SolverIterations = Current;
GMAT DefaultOrbitView.UpperLeft = [ 0.002772002772002772 0.007623888182973317 ];
GMAT DefaultOrbitView.Size = [ 0.4968814968814969 0.8081321473951716 ];
GMAT DefaultOrbitView.RelativeZOrder = 98;
GMAT DefaultOrbitView.Maximized = false;
GMAT DefaultOrbitView.Add = {Asteroid357439, Earth, Sun, Mercury, Venus, Mars, Jupiter, Saturn, Uranus, Neptune};
GMAT DefaultOrbitView.CoordinateSystem = SunMJ2000Ec;
GMAT DefaultOrbitView.DrawObject = [ true true true true true true true true true true ];
GMAT DefaultOrbitView.DataCollectFrequency = 1;
GMAT DefaultOrbitView.UpdatePlotFrequency = 50;
GMAT DefaultOrbitView.NumPointsToRedraw = 0;
GMAT DefaultOrbitView.ShowPlot = true;
GMAT DefaultOrbitView.ShowLabels = true;
GMAT DefaultOrbitView.ViewPointReference = Sun;
GMAT DefaultOrbitView.ViewPointVector = [ 0 0 300000000 ];
GMAT DefaultOrbitView.ViewDirection = Sun;
GMAT DefaultOrbitView.ViewScaleFactor = 1;
GMAT DefaultOrbitView.ViewUpCoordinateSystem = SunMJ2000Ec;
GMAT DefaultOrbitView.ViewUpAxis = Z;
GMAT DefaultOrbitView.EclipticPlane = Off;
GMAT DefaultOrbitView.XYPlane = On;
GMAT DefaultOrbitView.WireFrame = Off;
GMAT DefaultOrbitView.Axes = On;
GMAT DefaultOrbitView.Grid = Off;
GMAT DefaultOrbitView.SunLine = Off;
GMAT DefaultOrbitView.UseInitialView = On;
GMAT DefaultOrbitView.StarCount = 7000;
GMAT DefaultOrbitView.EnableStars = On;
GMAT DefaultOrbitView.EnableConstellations = On;


%----------------------------------------
%---------- Mission Sequence
%----------------------------------------

BeginMissionSequence;
Propagate AsteroidPropogator(Asteroid357439) {Asteroid357439.ElapsedDays = 30};
timmyd7777
 
Posts: 3
Joined: Wed Jan 14, 2015 2:18 am

Re: Need help propogating a near-Earth asteroid

Postby timmyd7777 » Mon Jan 26, 2015 7:19 am

OK, problem (mostly) solved - my propogator was missing the Sun from the list of point masses. Somehow I assumed that since the Sun was the central body, it wouldn't have to also be listed as a point mass. Nevertheless, it does make sense that you need to include the Sun, if you're trying to model the trajectory of an object orbiting the Sun.

I found an online Kepler's Equation calculator for quick-and-dirty mean anomaly -> true anomaly conversions here:
http://www.jgiesen.de/kepler/kepler.html

However, it would still be fantastic if the mean anomaly could be specified when originally inputting the asteroid's (er, spacecraft's) orbit. Where do I submit this request?

-Tim
timmyd7777
 
Posts: 3
Joined: Wed Jan 14, 2015 2:18 am


Return to Mathematics and Modeling

Who is online

Users browsing this forum: No registered users and 2 guests