Orbital Coordinate Systems, Part IIIBy Dr. T.S. Kelso 



Orbital Coordinate Systems, Part IIIBy Dr. T.S. Kelso 
Last time, we worked through the process of calculating the ECI (EarthCentered Inertial) coordinates of an observer's position on the Earth's surface, starting with the observer's latitude and longitude. Then, we used those coordinates to calculate look angles (azimuth and elevation) from the observer's position to an orbiting satellite. The most difficult part of that process was in calculating the sidereal time, a quantity necessary to determines the Earth's orientation in inertial space.
In the process of performing those calculations, however, we made one simplifying assumption: that the Earth is a sphere. Unfortunately, this assumption is not a good one. Ignoring the fact that the Earth's shape can more accurately be described as an oblate spheroid (a flattened sphere) can have a significant effect in certain types of satellite tracking applications. In this column, we will examine the implications of our initial assumption by modifying our calculations to allow for the Earth's flattening at the poles and then tackle the related problem of determining the subpoint of an orbiting satellite. Let's start by looking at a crosssection of the Earth and defining some terms.
Figure 1. CrossSection of Oblate Earth
Figure 1 is an exaggerated view of the crosssection of the Earth. For an observer on the Earth's surface, we can define a couple of terms fairly easily. The first is the local zenith. The local zenith direction is just a fancy way of saying "straight up." It is the direction away from a point on the Earth's surface perpendicular (at a right angle to) the local horizon. On a sphere, this direction is always directly away from the Earth's center. However, on an oblate spheroid, this is not the case since a line from the center of the Earth to the observer's position would not point to the local zenith (except on the equator and at the poles).
Since the local zenith direction depends upon the local horizon, let's take some time to better define it, as well. The local horizon is a plane which is tangent (touching at a point) to the Earth's surface at the observer's position. For our purposes, we will consider the local horizon to be the plane tangent to the reference spheroid. The term reference spheroid is used to define the oblate spheroid which 'best' defines the shape of the Earth. How 'best' is defined is a complicated process and depends upon whether the fit of the reference spheroid is regional or global. We will use the reference spheroid defined in WGS72 (World Geodetic System, 1972) for our standard.
In WGS72, the Earth's equatorial radius, a, is defined to be 6,378.135 km. The Earth's polar radius, b, is related to the equatorial radius by something called the flattening, f, where
b = a(1  f)
The flattening term, as defined in WGS72, is only 1/298.26—a very small deviation from a perfect sphere. Using this value, the Earth's polar radius would be 6,356.751 km—only 22 kilometers difference from the equatorial radius.
The first real significance of using an oblate spheroid instead of a sphere to define the Earth's shape comes in determining the observer's latitude. On a sphere, latitude is defined as the angle between the line going from the center of the Earth to the observer and the Earth's equatorial plane. However, on an oblate spheroid, geodetic latitude is the angle between the local zenith direction and the Earth's equatorial plane. This angle, φ, is the latitude used on maps; the angle formed by the observer's position, the Earth's center, and the equatorial plane is more properly referred to as the geocentric latitude, φ'.
The impact of this change is that in order to calculate the observer's ECI position, we must determine the geocentric latitude from the geodetic latitude. Knowing the geocentric latitude, φ', we can then calculate the geocentric radius, ρ, and from that calculate the z coordinate (ρ sin φ') and the projection in the equatorial plane (ρ cos φ'). Let's start by developing the relationship between φ and φ' since we'll usually be given φ.
From the basic definition of an ellipse,
where
R' = ρ cos(φ')
and
z' = ρ sin(φ').
Now,
tan(φ') =
and
tan(φ) =
(that is, the normal to the tangent of the spheroid). Differentiating the equation of the ellipse,
and rearranging terms,
which can be written as,
tan(φ') = tan(φ) = (1f)^{2} tan(φ).
So, knowing the geodetic latitude and the flattening, we can now determine the geocentric latitude. Now, let's see how much of a difference results from using an oblate spheroid. Figure 2 plots the difference between geodetic and geocentric latitude as a function of geodetic latitude.
Figure 2. Geocentric vs. Geodetic Latitude
That's it? All that work and the maximum error is less than twotenths of a degree? It would hardly seem worth the effort to perform the calculation. But let's explore a little further.
Although the development is too complicated to present here, it can be shown that
ρ sin(φ') = z' = a S sin(φ)
and
ρ cos(φ') = R' = a C cos(φ)
where
.
Our ECI coordinates, are now
x' = a C cos(φ) cos(θ)
y' = a C cos(φ) sin(θ)
z' = a S sin(φ).
Using the example of calculating the ECI coordinates of 40° N (geodetic) latitude, 75° W longitude on 1995 October 01 at 9^{h} UTC,
x' = 1703.295 km, y' = 4586.650 km, z' = 4077.984 km.
Although close to our calculations assuming a spherical Earth, we find this simplification resulted in a position error of 22.8 km.
What we really want to know, however, is just how big an error will result when generating look angles to a satellite from an observer's position on the Earth's surface if we assume a spherical Earth. From Figure 2, we would expect to have the largest errors for observers around 45° N latitude, so let's use a location near Minneapolis at 45° N latitude and 93° W longitude for our example. On a pass of the Mir space station over Minneapolis on 1995 November 18, Mir passed almost directly overhead. At 12^{h} 46^{m} UTC, its ECI position was calculated to be: x = 4400.594 km, y = 1932.870 km, z = 4760.712 km. Calculating the look angles for both a spherical and oblate Earth yields the results shown in Table 1.
Table 1. Look Angles for Spherical vs. Oblate Earth
Spherical Earth  Oblate Earth  
Azimuth  118.80°  100.36° 
Elevation  80.24°  81.52° 
The pointing error produced by assuming a spherical Earth is 3.17 degrees. For most applications, this error might not be significant. However, in applications involving tracking with highgain, typically narrowbeamwidth, antennas, an error of 3 degrees can result in a loss of communications.
So, now that we've completed the calculation of a satellite look angle for an oblate Earth, let's look at how to calculate the subpoint of a satellite in Earth orbit. We'll begin by examining the calculations for a spherical Earth first before looking at the case for an oblate Earth.
First, let's be sure we understand what we're looking for. The satellite subpoint is that point on the Earth's surface directly below the satellite. For the case of a spherical Earth, this point is the intersection of the line from the center of the Earth to the satellite and the Earth's surface, as shown in Figure 3.
Figure 3. Calculating Satellite SubPoint—Spherical Earth
Given the ECI position of the satellite to be [x, y, z], the latitude is
and the (East) longitude is
where θ_{g} is the Greenwich Mean Sidereal Time (GMST). The altitude of the satellite would be
where R_{e} is the Earth's circular radius.
As seen in Figure 4, the calculation for an oblate Earth is somewhat more complicated. The first thing we notice is that our definition of satellite subpoint requires some refinement. The point on the Earth's surface directly below the satellite is not on a line joining the satellite and the center of the Earth. Instead, it is that point on the Earth's surface where the satellite would appear at the zenith.
Figure 4. Calculating Satellite SubPoint—Oblate Earth
Calculating the longitude of the satellite's subpoint doesn't change. However, to calculate the geodetic latitude of the satellite subpoint, we'll want to begin by approximating φ with φ' (as calculated above) and letting (for computational efficiency). Then, we'll want to loop through the following calculations
until is within the desired tolerance. To compute the altitude of the satellite above the subpoint,
.
Using our example of Mir passing over Minneapolis on 1995 November 18 at 12^{h} 46^{m} UTC yields a subpoint at 44.91° N (geodetic) latitude, 92.31° W longitude, and 397.507 km altitude. And while we cannot solve for the subpoint directly, the number of iterations required is typically quite small. For this example, the value of after the first iteration is 0.180537 degrees, after the second iteration it's 0.000574 degrees, and after the third iteration it's 0.000002 degrees.
Admittedly, some of the differences we've found may seem small, but that will depend upon your tracking requirements. And, since they are not that much more difficult to calculate, there is little reason not to use them. As always, if you have questions or comments on this column, feel free to send me email at TS.Kelso@celestrak.com or write care of Satellite Times. Until next time, keep looking up!