Max Degree / Order in GravityModel

  • Anonymous - 2013-03-14


    I want to use the GravityModel class to compute accelerations affecting satellites due to the gravitational potential. I wonder if it's possible to define the degree and order up to which I want the field to be calculated during the runtime. As far as I understand the code, it is always computing the gravitational field to the highest possible level?

  • Charles Karney

    Charles Karney - 2013-03-14

    At present a full sum is always done. You could construct model
    files with a subset of the coefficients. Alternatively modifying
    the code to allow the max degree and order to be specified in the
    constructor would be possible. Actually, it looks really easy, just
    modify the two calls to the SphericalHarmonic constructor in the
    constructor for GravityModel and substitute smaller numbers for the
    4th and 5th arguments. This isn't perhaps exactly what you want because
    you would need to instantiate a new GravityModel for each set of
    [max-degree, max-order] you want. But it might be enough to get you

  • Freddy Kay

    Freddy Kay - 2013-03-15

    Hi Charles,

    Thank you for your help! Though I actually wanted to avoid modifiyng the code, I tried you suggestion anyway and it works just fine.
    As far as I can tell there are only two calls to SphericalHarmonic, just as you said, or did I miss one?
    And another question: have you ever crosschecked the results of your Code? Before I tried your library, I already used another implementation which can be found here (see GEOGRAV). The results I get differ by aproximately 0.1%. Do you have any suggestions as to why this is the case?

    Method: GravityModel::V();
    ECI-Coords used [m]:
    X = -5522333.38232309840
    Y = 4238164.28094564470
    Z = 1208232.17500043990

    gx = 6.2485000069713497
    gy = -4.7954346178341929
    gz = -1.3708162469637046

    Results by GEOGRAV (though here they use ECEF Coordinates):
    ax = 6.2481691347188901
    ay = -4.7954745597432824
    az = -1.3707634134685798

  • Charles Karney

    Charles Karney - 2013-03-15

    I guess that you are using the egm96 gravity model??

    I carried out various checks on the Gravity class: specifically
    verifying that the spherical harmonic sums are carried out correctly
    (without underflow or overflow) and cross-checking the results it gives
    for the geoid height against the tabulated geoid heights given by NGA.

    As to what can cause the discrepancy of 0.1% between GeographicLib and
    another software package... Possibilities are:

    • accumulated round-off errors (but this is unlikely with egm96). I ran
      the test point you gave comparing a calculation with double precision
      (53 bit fractions) and long doubles (64 bit fractions). Calling
      V(-5522333.3823230984L, 4238164.2809456447L, 1208232.1750004399L,
      GX, GY, GZ);
      gives for GX, GY, GZ
      53-bit 6.2485000069713497 -4.7954346178341929 -1.3708162469637046
      64-bit 6.2485000069713498 -4.7954346178341933 -1.3708162469637048
      So GeographicLib is calculating the sum accurately. You can try
      reducing the degree and order for both packages to check whether
      round-off is a problem with GEOGRAV.

    • your mysterious distinction between ECI coordinates (which I've never
      heard of) and ECEF coordinates (which I would call geocentric).
      You're on your own on this one.

    • different gravity models. I assert that GeographicLib does things the
      NGA way; but let me know if this is not so. GEOGRAV uses a model that
      is subtly different (different values for a, GM, C21, S21. The
      differences don't quite look big enough to explain the 0.1%
      discrepancy. But you can try fiddling with these constants in one or
      other of the two packages to check.

    You're quite right to worry about a discrepancy of this magnitude. In
    your shoes I would want to isolate the reason for the discrepancy.

    • Freddy Kay

      Freddy Kay - 2013-03-20

      Hey Charles,

      I tried to figure out where the discrepency is coming from. I used STK to calculate the accelerations and fact is I get results which are in compliance with the results in GEOGRAV and show more or less the same discrepency which I mentioned above. As you already checked the sum calculation, I suspect there is a problem with the reference frames.

      ECI (Earth centered inertial) and ECEF (Earth centered Earth Fixed) are both geocentric coordinates. THe ECI frame refers to a system where the z-Axis is pointing towards the northpole/earth rotation axis and the x Axis is pointing towards the vernal equinox, most commonly it is ecpressed in terms of cartesian coordinates but you can also use spherical coordinates (Right ascension, declination and distance). The ECEF frame on the other hand is a system where the z-Axis is again directed towards the northpole but the x Axis is pointing towards the Greenwich meridian, it is therefore rotating in coincidence with the Earth's rotation It is most commonly expressed in terms of longitude, latitude and height above ground/distance from origin.

      Whats weird when I use your library is, that when I use coordinates in the ECI frame, the results I get are of the same magnitude. If I use an ECEF frame I cannot reproduce the results I get with your library. As an ECI frame would need some kind of information on the timeframe in order to evaluate the geopotential of a rotating body, it seems weird that the results of GeographicLib are in the same region as the results I get with Geograv and STK (which include a time dependent conversion).

      I already started to go through your code to find where and how you are converting coordinates in GravityModel but figured it would be easier to ask right away. So, where are you converting coordinates in GravityModel and are you converting them at all? And more important which Reference System are you using? Does it account for the true centre of gravity or is it using the geometric centre of earth? Are you maybe using a Coordinate System that is refering to an Epoch?

      I hope you can tell me where I might have misunderstood your library or confused the coordinates or tell me if I am completely off track.

      Thanks a lot in advance!

  • Charles Karney

    Charles Karney - 2013-03-24

    The coordinates systems that GeographicLib uses are all tied to the
    earth. The only conversions done are simple ones, geocentric to
    geographic and the like. I believe I am faithfully implementing the
    EGM96 and EGM2008 gravity models published by NGA. So you should
    consult the NGA documentation for an explanation of the model:

    I think you'll find that the results GeographicLib gives accurately
    matches the results from the NGA programs, e.g., f477.f or

    Incidentally there's no time dependence incorporated into these models.



Cancel  Add attachments

Get latest updates about Open Source Projects, Conferences and News.

Sign up for the SourceForge newsletter:

No, thanks