3
$\begingroup$

I am currently working on a project from this book. In appendix D, they provide a few projects that can be coded as a review of all the material in the book. I finished the first project, Site/Track and my output is reasonably close (thanks to the numerical imprecision of Python) to the solution. I tried to use the Orbital python library to do some plotting of my solution to add a nice visual accent to the project.

I input the state vectors, converting from Distance Units and Time Units defined in the book to where mu=1, to meters and meters per second (as required by orbital's function conversion function. The output was an eccentricity value of .98, and a semimajor axis about equal to the earths radius. This obviously isn't correct so I proceeded the next data point - same result. Okay maybe it's my code? So I looked up the state vector for the ISS, found here, input that into my code and out comes a beautiful graph of the orbit. So maybe it's just the data in the book. I enter a different example from the book from the chapter that I followed to create the program, where they convert state vectors to orbital elements. It gets that calculation perfect also. I then try one more source here, and these vectors do not work in my program. I also put all 3 vectors into this calculator, the ISS and one from my book worked, but not the others. I'm trying to find the inconsistency that I'm obviously missing...

So my question is. Are there different reference frames for these vectors that could be causing these errors. What other errors could I be inadvertently introducing?

Here is my input data from the book (not functioning) defined in the book as the IJK reference frame (inertial)

r = [ 1779987.13023855 -4944211.65294755  4065801.40507205] 
v = [ 2082.84303416 -1179.76617351   410.70269708]

Output:

KeplerianElements:
    Semimajor axis (a)                           =   3493.823 km
    Eccentricity (e)                             =      0.950399
    Inclination (i)                              =     45.1 deg
    Right ascension of the ascending node (raan) =    160.3 deg
    Argument of perigee (arg_pe)                 =    303.1 deg
    Mean anomaly at reference epoch (M0)         =    144.4 deg
    Period (T)                                   = 0:34:15.239378
    Reference epoch (ref_epoch)                  = J1970.000
        Mean anomaly (M)                         =    144.4 deg
        Time (t)                                 = 0:00:00
        Epoch (epoch)                            = J1970.000

The ISS data I used (proper graph shown) defined on the website as m50 cartesian.

r = [4607312.46, -1531324.39, 4749270.39]
v = [4597.800926, 5516.878978, -2671.990580]

Output:

KeplerianElements:
    Semimajor axis (a)                           =   6794.798 km
    Eccentricity (e)                             =      0.000980
    Inclination (i)                              =     51.4 deg
    Right ascension of the ascending node (raan) =    212.9 deg
    Argument of perigee (arg_pe)                 =     53.7 deg
    Mean anomaly at reference epoch (M0)         =     62.7 deg
    Period (T)                                   = 1:32:54.113176
    Reference epoch (ref_epoch)                  = J1970.000
        Mean anomaly (M)                         =     62.7 deg
        Time (t)                                 = 0:00:00
        Epoch (epoch)                            = J1970.000
$\endgroup$
7
  • 2
    $\begingroup$ I'm not qualified to answer this question, but if I were, I'd want to see your input data and your program's outputs for each of the cases you mention, not just words about them. $\endgroup$ Feb 24 at 22:11
  • 2
    $\begingroup$ "The output was an eccentricity value of .98" Why can I not find this number in your tables? $\endgroup$ Feb 24 at 22:45
  • 1
    $\begingroup$ There's obviously something wrong with your first example, because a semimajor axis less than the radius of the earth isn't possible. Bate, Mueller, and White contains lots of material about ballistic missile trajectories, so it is always good to check whether the data they've given you corresponds to an "orbit" that only lasts for 20 minutes between launch and impact, but a < 3500 km is deep inside the earth at all times. $\endgroup$
    – Ryan C
    Feb 24 at 23:00
  • 1
    $\begingroup$ @RyanC it's a good point; I've hypothesized that it's a reusable first stage. $\endgroup$
    – uhoh
    Feb 24 at 23:15
  • 1
    $\begingroup$ SE Adjusted. Thanks, there is no code because I don’t believe my code is the problem, as all the outputs match. I agree with @uhoh and think the velocity is the problem, which is why I looked up this myreaders.info/05_Satellites_Orbit_Elements.pdf dataset. These still have terribly, almost hyperbolic, orbits but I will plug it into the vis-viva equation to check if the data is just bad for a Satellite. Assuming the data is “bad” I will assume the only state vector that produces a “good” orbit is the ISS. Maybe I’ll check SPICE for some more datasets to test against. $\endgroup$ Feb 25 at 14:44
4
$\begingroup$

Here's a partial answer until you add more information as requested in comments

thanks to the numerical imprecision of Python

I don't think you are anywhere near the limit of Python's floats. Instead let's remember that Keplerian orbits are theoretical approximations only. The biggest deviations come from Earth's equatorial oblateness as expressed by $J_2$ which is about a 1 part per thousand effect. If you just use the central $GM/r^2$ force you're going to be off by of order 0.1%. Then there are higher order multipole terms, gravity from the Sun and Moon, etc.

Here's the initial state vectors that you report. (I think you've converted from unit-less numbers in your source? It would be nice to see them also)

r = [1779987.13023855, -4944211.65294755, 4065801.40507205]  
v = [2082.84303416, -1179.76617351, 410.70269708]  

The vis-viva equation is incredibly useful. For central force (Keplerian) orbits, it relates the scalar values $a$, $r$ and $v$ regardless of what direction those last two vectors are pointing, so it's very handy!

$$v^2 = GM\left( \frac{2}{r} - \frac{1}{a} \right)$$

or

$$\frac{1}{a} = \frac{2}{r} - \frac{v^2}{GM}$$

Earth's standard gravitational parameter (the product $GM$) is about 3.986E+14 m^3/2^2.

You have $r=$ 6644.12 km and $v=$ 2.42874 km/s. The equatorial radius of Earth (6378.137 km) is usually used to assign an altitude value, so the object is at an altitude of 266 km moving only 31% of the 7.746 km/s required for a circular orbit at that distance, so ya it's gonna dive hard into the Earth.

The vis-viva equation above gives $a=$ 3493.823 km exactly as you show in the output which is no surprise, so I won't confirm the other items.

It might be a poorly thought out problem sending the spacecraft into a deep dive, or since it's so slow it might even survive reentry if it has engines, a bit like a reusable first stage at the top of its arc!

$\endgroup$
1
  • 1
    $\begingroup$ Thank you. I knew that the velocity values looked off, I was just trying to make some pretty graphic presentations of the data :) $\endgroup$ Feb 25 at 14:47

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.