ECI Cartesian Coordinates to Kepler Orbit Elements

Uriah

New member
Joined
Apr 21, 2015
Messages
5
Reaction score
0
Points
0
Hi,

I'm using JSBSim as the physics engine for simulating space flight in Outerra Anteworld.

In one of the XML files pertaining to the Jupiter 246 / Ares IV launch vehicle, there are a few equations for calculating the semi-major axis, eccentricity, apoapsis, periapsis and orbital period. I would like to complete the calculations for the rest of the orbital elements, but noticed a serious flaw in the method used for eccentricity, which brings into question the rest of the equations used.

At the start of the simulation, where the rocket is on the launch pad, the eccentricity value is initially 1.0, and during ascent decreases slowly towards 0.0. Immediately I knew something was wrong, as the eccentricity should begin at 0.0 and increase during ascent, until you perform a burn to circularize the orbit, then it should return to near 0.0. I've included the XML equations, as well as what I have re-written into JavaScript functions. They are deriving the Keplerian elements from ECI Cartesian coordinates.

Can anyone determine where they went wrong?

Best regards,
Uriah George

XML:
Code:
<property value="14.07644180e15"> guidance/executive/gm </property>

  <channel name="Orbital Parameters">

    <fcs_function name="guidance/executive/semimajor-axis-length-ft">
      <description>
        This equation calculates the semimajor-axis length in feet. a = mu/(2*mu/r - v^2) 
      </description>
      <function>
        <quotient>
          <p> guidance/executive/gm </p>
          <difference>
            <quotient>
              <product>
                <v> 2.0 </v>
                <p> guidance/executive/gm </p>
              </product>
              <p> position/radius-to-vehicle-ft </p>
            </quotient>
            <pow>
              <p> velocities/eci-velocity-mag-fps </p>
              <v> 2.0 </v>
            </pow>
          </difference>
        </quotient>
      </function>
    </fcs_function>

    <fcs_function name="guidance/executive/eccentricity">
      <description>
        This equation calculates eccentricity. e = sqrt (1 - [ r v cos(gamma) ]^2 / {a * mu})
      </description>
      <function>
        <pow>
          <difference>
            <v> 1.0 </v>
            <quotient>
              <pow>
                <product>
                  <property> position/radius-to-vehicle-ft </property>
                  <property> velocities/eci-velocity-mag-fps </property>
                  <cos>
                    <property> flight-path/gamma-rad </property>
                  </cos>
                </product>
                <v> 2.0 </v>
              </pow>
              <product>
                <p> guidance/executive/semimajor-axis-length-ft </p>
                <p> guidance/executive/gm </p>
              </product>
            </quotient>
          </difference>
          <v> 0.5 </v>
        </pow>
      </function>
    </fcs_function>

    <fcs_function name="guidance/executive/apoapsis">
      <function>
        <quotient>
          <difference>
            <product>
              <property>guidance/executive/semimajor-axis-length-ft</property>
              <sum>
                <value>1</value>
                <property>guidance/executive/eccentricity</property>
              </sum>
            </product>
            <property>inertial/sea-level-radius_ft</property>
          </difference>
          <value>5280</value>
        </quotient>
      </function>
    </fcs_function>

    <fcs_function name="guidance/executive/periapsis">
      <function>
        <quotient>
          <difference>
            <product>
              <property>guidance/executive/semimajor-axis-length-ft</property>
              <difference>
                <value>1</value>
                <property>guidance/executive/eccentricity</property>
              </difference>
            </product>
            <property>inertial/sea-level-radius_ft</property>
          </difference>
          <value>5280</value>
        </quotient>
      </function>
    </fcs_function>

    <fcs_function name="guidance/executive/orbital_period">
      <description>
        This equation calculates orbital period in seconds. T = (2 * pi) sqrt ({a^3} / mu)
      </description>)
      <function>
        <product>
          <product>
            <v> 2.0 </v>
            <property> guidance/executive/pi </property>
          </product>
          <pow>
            <quotient>
              <pow>
                <property>guidance/executive/semimajor-axis-length-ft</property>
                <v> 3.0 </v>
              </pow>
              <property> guidance/executive/gm </property>
            </quotient>
            <v>0.5</v>
          </pow>
        </product>
      </function>
    </fcs_function>

JavaScript:

Code:
//This equation calculates the semimajor-axis length in feet. a = GM/(2*GM/r - v^2)
// GM = standard gravitational parameter of the central body, r = orbit radius (ft), v = orbit velocity (fps)
function semi_major_axis(r, v, GM){
	return GM/(2.0 * GM/r - Math.pow(v,2) );
}

//This equation calculates eccentricity. e = sqrt (1 - [ r v cos(gamma) ]^2 / {a * GM})
// GM = standard gravitational parameter of the central body, r = orbit radius (ft), v = orbit velocity (fps), gamma (rad) = 
function eccentricity(a, r, v, gamma, GM){

	return Math.sqrt(1.0 - Math.pow( ( r * v * Math.cos(gamma) ), 2.0) / (a * GM ) );
}

//This equation calculates apoapsis in miles. ap = ((a * {1+e}) - R) / 5280
function apoapsis(a, e, R){
	return ((a * (1+e))-R)/5280.0;
}

//This equation calculates periapsis in miles. pe = ((a * {1-e}) - R) / 5280
function periapsis(a, e, R){
	return ((a * (1-e))-R)/5280.0;
}

//This equation calculates orbital period in seconds. T = (2 * pi) sqrt ({a^3} / GM)
function orbital_period(a, pi, GM){
	return (2.0 * pi) * Math.sqrt(Math.pow(a,3.0)/GM);
}
 

Keithth G

New member
Joined
Nov 20, 2014
Messages
272
Reaction score
0
Points
0
I'm not sure they are wrong.

If you are stationary on the Earth's surface (and ignoring the Earth's rotation for a moment), the eccentricity of your orbit is 1 - i.e., you are at the apogee of a highly eccentric ellipse. The only reason that one doesn't move in this elliptical path is because the you are subject to a resistive force from the ground.

If you throw in the effect of the Earth's rate of rotation while sitting on the Earth's surface, the eccentricity is going to be slightly less than 1 - but not by much.

As you approach a circular orbit, your eccentricity will fall to zero.
 
Last edited:

Uriah

New member
Joined
Apr 21, 2015
Messages
5
Reaction score
0
Points
0
Oh I see now, that makes perfect sense. Because it is angular momentum relative to the earth, and if your velocity is zero and the ground wasn't there, you would fall towards the (point) mass center, so your orbital trajectory would therefore be parabolic, or e = 1.

I now see that at launch in Orbiter that e = 1.

Thanks for the help!

Next I need to determine the inclination from the orbital momentum vector. :)

Regards,
Uriah George
 
Top