- Email: [email protected]
The use of the Doppler effect to deduce an accurate position for an artificial earth satellite
Planet.
Space Sci., 1962, Vol. 9, pp. 607 to 623.
Pergamon
Press Ltd.
Printed
in Northern
Ireland
THE USE OF THE DOPPLER EFFECT TO DEDUCE AN ACCURATE POSITION FOR AN ARTIFICIAL EARTH SATELLITE E. GBLTON Department of Scientific and Industrial Research, Radio Research Station, Ditton Park, Slough (Received 4 June 1962)
Abstract-Methods by which observations of the Doppler shift of radio signals from. earth satellites may be used to determine the satellites position are discussed. Two methods for the calculation of the theoretical DoppIer shift are given, followed by a method of linear dif%rential orbit correction in terms of Keplerian orbital elements. Limitations on the possibility of solution are discussed, together with the use of single and multiple receiving sites, and limitations due to experimental measurement accuracy. Methods by which corrections for the effect of the ionosphere may be made are briefly considered. INTRODUCTION A knowledge of the position of a satellite at a given time is of impo~an~ both to experimenters interested in the reduction of telemetry data and to those deducing ionospheric information from signal propagation studies. Unfortunately past experience has shown that it is usually many months before reasonably accurate ephemeris becomes available, if at all. This paper describes how observations of the well-known Doppler effect may be used to obtain positioned data of sufficient accuracy, soon after the satellite passage. Little appears to have been published hitherto on the accurate reduction of Doppler shift observations. Since time variation of the magnitude of the effect observed from the rotating earth is a unique, though nonlinear, function of the orbital elements of the satellite, it follows that it should be possible to determine the orbital elements from a set of observations of the Doppler effect. The practical possibility of doing this appears to have been first described by Guier. (l) Otherwise it has been customary to divide the plotted Doppler curve by some geometrical method to obtain a “time of nearest approach” and to use the maximum slope to obtain a “slant range”. However, these calculations are based on tremendously simplifying assumptions, and consequently the time obtained is usually in error by a few seconds, and the range in error by an appreciable fraction. Whilst these values may be very useful in predicting future passages of a satellite, they are of little use in studying the detailed orbital motion of the satellite or in fixing its position accurately. One method by which the orbital elements might be corrected by observation of the Doppler shift is to compare the measured values with theoretical values derived from a trial set of elements, and calculate the mean square deviation: this is carried out for many different sets of elements, and the set giving the smallest mean square deviation is nearest to the correct set of elements. In practice this could mean searching for a minims with perhaps ten values of each of seven elements, i.e. 10’ sets, the Doppler shift being calculated perhaps 200 times for each set of elements. Thus although this method is probably the best for solving a nonlinear system of the type under discussion, it is only of practical use where the largest and fastest electronic computers are available, Another method which might have a more widespread application is to derive a system 607
608
E. GOLTON
of linear equations which relate a variation in Doppler shift, or the equivalent range rate S:, to variations in the orbital elements, xi, of the form:
Differences between an observed Doppler shift and values calculated from a trial set of elements might then be used in an equation solution to give corrections to the etements. N A
Ascending FIG. 1. DIAGRAM
node
SHO~NGELEMENTSOFSA~~LITEORB~.
Because the function S(x,) at any instant is in fact nonlinear in x6 the derived corrections will not be exact, and indeed the system may not provide useful corrections at all. Because of this the nonlinearity of the coefficients in a linear equation system must be investigated, together with the effect of errors in the observed Doppler shift values, which will also limit the accuracy of the first method of solution described. In the present paper expressions for the theoretical Doppler shift will be developed in terms of Keplerian orbital elements defined at the ascending node, assuming a constant rate of change of perigee and ascending node positions. From these will be developed a differential correction scheme, and the necessary investigations mentioned in the previous paragraph will then be carried out. It is obvious that in order to attempt useful Doppler shift measurements one must have a stable transmitter frequency at the satellite, and a stable reference frequency at the ground observing site. At the present time crystal controlled oscillators with stabilities of a few parts in lOlo are in use in certain satellites. The observed Doppler shift will not be quite the same as that observed under free space conditions, because of the effect of the ionosphere: however, the first order effect may be removed by the transmission of two phase locked signals from the satellite, and the second order effect by a study of the ionosphere. These processes will be outlined at the end of the paper. DEFINITIONS
A geocentric co-ordinate system is used, the reference directions being the northward axis of the earth and the first point of Aries. Angles measured eastward are considered
THE USE OF THE DOPPLER EFFECT TO DEDUCE AN ACCURATE POSITION
609
positive. A satellite moving from west to east is considered primarily, although the analysis of the complementary case is almost identical. The satellite motion is assumed to be in the form of a Kepler ellipse, which is allowed to precess about the earth’s axis, and the major axis of which is allowed to rotate in the mean orbit plane. The last two secular additions represent the main perturbations of the orbit due to the oblateness of the earth, e.g. King-Hele.(2) The main symbols used are defined below; others will be introduced as required. Some of the symbols are shown in Fig. 1. The words “ascending node” will be abbreviated to “A.N.” hereafter. t, time of the A.N. T absolute time of the satellite at point S. i inclination of orbit To satellite anomalistic period of revolution e eccentricity angle from A.N. to perigee at time t, Right ascension of A.N. at time t, angle from A.N. to perigee at time T Right ascension of A.N. at time T rate of change of CO rate of change of Sz Right ascension of Greenwich at time T = 0 earth’s sidereal rate of rotation Right ascension of observer at time T time from A.N. to perigee \ time from perigee to S I angle from perigee to S observers radial distance geocentric observers latitude values at P observers longitude radial distance of satellite geocentric latitude of satellite values longitude of satellite Eccentric anomaly between node and perigee Eccentric anomaly between perigee and satellite geocentric angle PBS between observer and satellite In all cases time is measured in terms of the mean solar day. to represent differentiation with respect to time. EQUATION
DESCRIBING
The basic equations
required
The dot notation
will be used
THE SATELLITE MOTION IN A KEPLER ELLIPSE
are as follows:
(1)
(2)
610
E. GOLTON
t = 2 (E - e sin E) sin E = (1 - e2)* *
sin Y (1 + e cos Y)
(4)
- tan E/2 27r (1 i e cos v)2 * = z * (1 - ,213
(6)
We have chosen to work in terms of the period TOrather than the semi-major axis length a, these being connected by the equation TO2= 4n2a3 P where ,u = GM, the product To solve these equations The observation is made The perigee angle is thus
of the gravitational constant and the mass of the earth. we must first consider time. at a time (T -- t,) after the A.N. o = co0+ c&(7’- t,)_ Equation (4) will now give Em. sin w sin E, = (1 - e216(1 + e cos w)
(7)
Substitution of E, into equation 3 gives: t, = 2 (E, - e sin Em) We can now find the time of the satellite from perigee. This is: t = T -
t, -
t,.
(9)
Equation 3 must now be solved using the value t found above to give E. This may be done by a process of iteration since the equation is monotonic. Equation 5 is now used to find the angle v. From this we can calculate r, i, and i by using equations 1,2, and 6. The position of the satellite must now be found. The Right Ascension of the A.N. at time T is CI = Q2,+ h(T - I,): that of Greenwich is (0, + 6 . T). By reference to Fig. 2 we see that the satellite has moved in longitude through an angle A, given by:
Also
tan 2 = tan (w + P) *cos i.
(10)
sin #s = sin ((0 + Y)* sin i.
(11)
Thus the satellite is at geocentric latitude 4s and longitude 0, = Q + Iz - (0, + 8.T). The difference in longitudes between satellite and observer is 8,-08,=a-e+n
(12)
The geocentric angle y between satellite and observer is given by: cos y = sin #@- sin A t_ cos & cos & cos (0, - 0,)
(13)
THE USE OF THE DOPPLER EFFECT TO DEDUCE AN ACCURATE POSITION
Kepler
611
ellipse
Equator
FIG. 2. DIAGRAMSOF ORBIT GEOMETRY.
The distance s between satellite and observer is given by: s2 = rS2 + y02- 2r,r, cos y CALCULATION
OF THE
THEORETICAL
DOPPLER
(14) SHIFT
(A) Method 1. We will first calculate the theoretical Doppler shift by finding the velocity of the satellite relative to the earth, and then taking the component in the direction of the observer. The first part of this method is essential for calculating the Doppler shift along a ray path through the ionosphere and for eliminating the second order ionospheric effect. With reference to Fig. 2 the velocity components relative to the earth are: 3, outwards, r,(C + h) cos B northwards, rs - (i’ + h) sin B - rs(8 - 0) cos & eastwards.
and cos B = tan &/tan (CO+ Y), sin 3 = sin @in (w + Y). Thus the satellite speed is given by V in the following equation, V2 = i,2 + r,“[(i, + 61)~+ (0 - fi)2 * cos2 $, - 2(i + 15) * (6 - 0) sin B * cos #,I (1%
and is directed at elevation a and bearing B given by: tan M.= i,/r,[(G + I%)”+ (6 - h)2 * co? y& - 2(2; + c.!I)(~- $) sin B - cos C&I* (16) tan /3 = [(+ + &) sin B - (6 - $) cos &]/(i, + h) cos B
(17)
With reference to Fig. 3 we now require to find the zenith angle e, and the bearing b, of the observer from the satellite. These are given best by: tanc,=siny/(z-cosy) cos b, = [sin $O - sin & +cos yJ/cos q3,* sin y.
(18) (19)
E. GOLTON
612
0 FIG.
3.
DIAGRAM
FOR
CALCULATION
OF SATELLITE TO OBSERVER
These four angles may now be combined cos x = -cos The actual Doppler
LINE
VELOCITY
COMPONENT
ALONG
SATELLITE
S-P.
to give the angle x between
e, * sin tl + sin e, * cos cc * cos (b, -
V and SP. 8)
(20)
shift is given by: II’-’
f
v.cosx
c
where f is the transmitter frequency refractive index (n = 1 at present).
and c is the speed of light, and n is the local phase
(B) Method 2. The free space Doppler shift may be calculated more easily by considering the rate of or by elimichange of distance to the satellite. From a separate geometrical calculation nating #Q from and substituting equation 12 into equation 13 we obtain : cos y = cos (v + Co)* cos $0 * cos (e -
sz)
+ sin (Y + 0) - [sin 9S0- sin i + cos +O * cos i * sin (f3 -
L?)] (21)
At this stage it will be convenient to introduce letters to represent various terms, the first being B, = cos y as in equation 21. Also the subscript ‘s’ attached to r, will be dropped where possible. Differentiation of equation 14 gives : (sS) = r . J -
r,iB, - rOr. ; @A
613
THE USE OF THE DOPPLER EFFECT TO DEDUCE AN ACCURATE POSITION
Differentiation
of B, and substitution (sS) = rL -
into the above yields
ro3Bz -
r,,r(ci, + +)B, + r&6’ -
fi)B4
(22)
where B, = sin $O sin i. cos (co + v) -
cos $,[sin (w + Y) cos (0 - a) - cos i * cos (o + Y) sin (0 -
Q)]
(23)
and B4 = cos +,,[sin (f3 -
Q) * cos (co + Y) -
cos (0 -
Cl) * sin (w + Y) * cos il
Equation 22 in conjunction with equation 14 enables (S) to be calculated, of change of distance to the satellite, and since
(24)
this being the rate
(S) = -V*cos~ the theoretical
Doppler
shift may be calculated THE DIFFERENTIAL
directly.
CORRECTION SCHEME
A linear relation between an infinitesimal variation A(.+) of the Doppler shift and infinitesimal variations Ax of the orbital elements must now be found. This will be done by partial differentiation of expressions found above. We start with the identity: ~3 - A(S) = s2 . A($ From
equation
-
(8s) . 4 - A(s2)
(25)
14 we obtain: + * A(s2) = r * Ar -
And from equation
r, * B2 * Ar -
r * r,,AB2
(26)
22 :
A@) = Ar * [3 + r,,(d -
$) * B, -
r,-,B2] -
+ Ai[r -
+ AB, * r - r,,(8 -
r,(h + ?‘) - B3]
Afi - r - r, * B4 0) -
(27)
AB, . r - r,,(b + 1’>-
AB, . r,,i
(A& + A+)r * r,,B,
It will be noticed that a variation of ci, and 2;has been allowed, although because these terms are generally small a useful correction may not be obtained. Substitution of equations 26 and 27 into equation 25 gives: s3A(S) = Ar * D, -
AB, . s2 * r - ro(ci) + ti)
AB, * D, + A3 * D, -
+ AB4 . s2 * r * r,(6 -
0) -
Ah * s2 * r - r,, * B4 -
(Ah + Ati)s2 . r * r,, * B,
(28)
where : D, = s2 * i + s2 . r,(8 -
0) * C, -
D, = s2 * r, * i -
The variations Appendix.
s2 * r,,(h + ?j)B, -
(si)r,, - r ;
D, = s2 * r -
(sS) * r + (sS) . r, * B, s2 - r,, * B,.
of the terms B,, B,, B,, are given below, the coefficients
I
(29)
C, being given in the
AB, = C, . AQ + C, * Ai + C, * (Aw + AY) AB, = C, . AQ + C,Ai + C, - (Aw + Av)
(30)
ABd = C, * AG + C, * Ai + C, * (Aw + Av) The variations
of r, t and ?; must now be expressed
as variations
of the orbital
elements.
E. GOLTON
614
This may be carried finally : Ar -= Y
out by partial logarithmic differentiation of equations 1, 2, 6, giving
2 AT, he e + cos Y ---e sm v Av e+ + 1 -+ e cos v 1 - ez 1 + e cos v 3 T, A3 AT, Ae IL= -$cotv*Av + i 3T@ e(l - e2) 3e 2e sin v . Av -?+Ae 12+~~Oysv+-_ 1 e2 (1 + e cos Y) i 0
A+ -= 4
(31)’ (32) (33)
Substitution of these expressions into equation 28 together with equations 29 and 30 give: A(S) = A, *AT, + A, *Ae + A, * A$ + A, * Aci, + A, * AY -+ A, *ASZ+ A, *Ai + A, * AU
(34)
where tbe A, are given in the Appendix. AI1 the above variations are of orbital elements except Av, which besides being a function of Ae and AT,, is a function of At, and Ao. The required dependence may be derived from equations 3, 4. Partial differentiation gives : At=g(l-e*cosE)*AE+$-*At-g*sinE*Ae
(35)
0
and (e -t cos v)
cotE*AE=(l
(36)
+ e cos 7.~)
Elimination of AE now gives: Av =f(v)
* At --f(v).
G . AT, + he. 0
where sin v * cot E * (1 + e +cos v) * 2~ f(v) = (e + cos Y)- (1 - e * cos E) * To
(38)
Equation 37 expresses a variation of Y in terms of variations of e, To, and t. Now t is dependent upon t, and t, as given by equation 9, so that At = -AT,
- At,
(39)
We can write down the dependence of At,,, immediately by analogy to equation 37, viz, : AW =f(~)
e At, -f(u)
. k AT0 + Ae - ie2
+ f(w) * 2 * sin Em)
(40)
These equations result from the dependence of the time taken to travel from ascending node to perigee upon eccentricity, period, and perigee position; this time in turn determines the time available to travel from perigee to the point of observation, at fixed time T. Substitution of equations 39, 40 into equation 37 gives: Av = -f(v)
f(v)
- At, - -
f(v)
aAGO- AT, - r o * 0 + CJ
.f(co)
+ f(v) * 2 (Sin-E + sin E,)
I
(41)
615
THE USE OF THE DOPPLER EFFECT TO DEDUCE AN ACCURATE POSITION
Here t + t, is simply the time from node to point of observation, T - t,,. Equation 41 may after this has been carried out, Q and 03 must be now be substituted into equation 34: replaced by their component values : Cl = Q2, + $(T
-
t,)
ro = CII~+ cb(T -
t,)
(42)
are :
Their variations
AL’=ACI,+(T-t&A&_&At, + (T -
AU = ho, Substitution
of equations
A(i)=
A,-Afo*(T5T,, -
[A&)
+
AZ +
t) n
(1
AT, + [A3 + A&T -
_
sin v + f*
$)
t,>l
*Ah + A,’ . Aq,
1
2
(sin
E
+
sin E,)
d!zf!
A,-_A
equation:
sin w
.f
A,f(v)
correction
t,)] . Aa + A, * AQ,
. At, + [A* + A,‘@- -
+ A,‘h]
1
.-.
1.
+
where A,’ =
C%. At,
41 and 43 into 34 give the final differential
+ A,0
A5
(43)
t,) . A& -
1
* Ae + A, - Ai
(44)
5 f(4’ ORBITAL PERTURBATIONS
A constant rate of change of perigee and ascending node positions are the only perturbations which have so far been included in the correction scheme. Other perturbations are generally small, and inclusion of theoretical expressions for them in the correction scheme would complicate and lengthen it considerably. However, it is quite simple to include theoretical expressions for the perturbations in calculating the elements at a given time, to be used in the correction scheme. For example, one perturbation not included previously is that of the Keplerian relation between period and semi-major axis, which is affected by the oblateness of the earth and is mainly a function of inclination. This can be included by simply modifying the value ,u = GM as a function of inclination. INVESTIGATION
The coefficients out to find. many times, values have satellite past
of the variations
OF THE SOLUTION SYSTEM
of elements
are the partial
derivatives
g,which
we set
Although most are long expressions, simple arithmetical operat;ons are used and the coefficients can be found very quickly on a small computer. Typical been plotted at one minute intervals in Figs. 4, 5 for a descending passage of a an observer at 51”N latitude, using the following elements:
T, 101.66
0.2:14
66!77
CI 352.0428
fi -0901775
,9;0,33
-0;005
t 73”0.7
Time is measured in minutes, and distance in thousands of kilometres. The correction coefficients are thus measured in units of lo3 km/min, and, where applicable, radians. They 2
616
E. GOLTON
show the variations in Doppler shift caused by an infinitesimal change in a particular
as
orbital element. For example, the curve of - in effect shows the slope of the Doppler curve, at, reaching the greatest (negative) value at the point of nearest approach. The coefficients were checked by finding the change in Doppler shift caused by small changes in each element. A3 When these hxi values were plotted, the values at Ax, = 0 were within O-1per cent of the
as
theoretical x
values in all cases, thus checking the correctness of the previous theory.
The soluti:n of a system of linear equations of the form:
. AxZj =
Asi
is affected by the accuracy of both the RHS and the coefficients of the LIB, generally to a greater extent the more variables we include. For a given error in a coefficient the effect will be more serious the worse the conditioning of the set of coefficients, Although we cannot alter the conditioning, it is useful to obtain a measure of it. This may be done by making use of the fact that each linear equation represents a hype~lane in n-dimensional space, and the normalized coefficients are the direction cosines of the vector normal to the hyperplane. The angles between these normal vectors are a direct measure of the conditioning. If the angles are small, the equations are badly conditioned, but if the angles are all near 90” the equations are very well conditioned. A given change in a coefficient results in a tilting of the hyperplane: In a perfectly conditioned set the largest change of one variable occurs when the hyperplanes are parallel to the axes, and is given by Z/zbn2 - b,2 4 sin 19where the b, are the normalized RHS and i is the equation with the chanied coefficient, and 8 is the change in direction of its normal. When the hyperplanes are not parallel to the axes, changes in all variables occur when one plane is tilted, but are of smaller magnitude than in the previous case. The degradation occurring when the hyperplanes are not perpendicular to each other is related to the cosecants of the angles between pairs of normals. The equations with coefficients plotted in Figs. 4,5 were investigated for several numbers of the variables. In particular with seven variables, AT,, Ae, Ai, AC?;,,Ao,, At,, AC, the last being the error in assumed transmitter frequency, and seven equations applying at one minute intervals, the smallest angle between any two hyperplanes was 12”, and the largest 90”. Thus a degradation of only about 5: 1 occurs, and a computer solution will not be limited by rounding errors. In the case of four variables, A!&, Ai, At,, AC, the smallest angle was 7 degrees, rather worse than before. We must now consider the linearity of the correction coefficients. For a given finite set of corrections to the elements and perfect RHS, the appropriate coefficients will not be quite 6% the same as the theoretical coefficients d;r. This is because the function $(x2) is nonlinear in X, at any given instant of time, The nonlinearity has been calculated by finding the change in Doppler shift for a finite change in each element, and expressing its deviation from the theoretical correction as a percentage. These values have been plotted below the corresponding correction coefficients in Figs. 4, 5. A constant percentage nonlinearity would be
THE USE OF THE DOPPLER EFFECT TO DEDUCE AN ACCURATE POSITION
al
g”s NON + I
1
617
618
E. GOLTON
THE USE OF THE DOPPLER EFFECT TO DEDUCE AN ACCURATE POSITION
619
of little consequence, since if all the coefficients of one variable in an equation solution were for example 2 per cent too small, the answer for the variable would simply be 2 per cent too large. However, the curves drawn in Figs. 4, 5, show considerable variation in magnitude and sign. The nonlinearity is dependent of course on the magnitude of the variation of the error in the orbital element, and will set a limit to the maximum permissible error in the trial set of elements from which a sensible convergent solution can be obtained. The RHS of the equation system contain observational errors and ultimately limit the accuracy of orbital element corrections once a convergent sequence has been obtained. The usual way of minimizing random effects is to use more equations than variabIes and to normalize the system. This in theory gives a least-squares fit to the RHS, assuming perfect coefficients, but in our case normalizing the system has the added advantage of smoothing out the effects of the nonlinearity of the coefficients, and increasing the range of convergence. This can be simply explained by the fact that errors in the coefficients of an equation are equivalent to a certain error in the RHS. As an example, it has been found possible to obtain convergence with 4 variables, AI&,, Ai, At,, AC using 8 equations. However, a least-squares solution is not reliable in this case, and this is reflected by the fact that a variation of any one RHS by O-1c/s in the measured Doppler shift at 50 MC/S gives rise to changes of the order of 0.1” and 4 set in time in the elements. However, when the solution was repeated using 47 points equally distributed along the Doppler curve, the same error in any one point gave rise to changes of the order of 0.001 degree and 0.005 set, an improvement of the order of 100 times. This accuracy is not of course always significant, since the other elements for which no solution was attempted are not always known to equivalent accuracy. Regular, rather than random errors impose a rather more serious limitation. In the previous example, using 47 points, a linear drift in transmitter frequency of O-1 c/s over 10 min at 50 MC/S results in corrections of about 0*03” and 0.15 sec. Such Doppler errors might also arise from imperfections in correcting ionospheric effects. At the present stage it has not proved possible to solve usefully for all the orbital elements independently using observations of one passage from one station, partly because of the small range of convergence and partly because of the limited accuracy of observations of the Doppler shift. So far only observations of one passage from one station have been considered. The solution becomes very much easier when more observing stations are available. With observations of several passages from one station it becomes possible to include the period correction AT, with AL&,,Ai, At,, AC, and obtain a strong solution. Usually the difficult elements to correct are perigee position and eccentricity. If the single station is situated at a lower latitude than apex, the observation of ascending and descending passages will cover two sections of the orbit, and useful corrections to these elements can be made. The rate of change 0 becomes important in this case, but rather than try to solve for it as an independent variable it appears more useful to obtain Q2,from each passage separately and then use the mean fi in a combined solution. When several observing sites are available, several different sections of the orbit may be observed, and it appears that all the elements can be corrected independently. If four stations are separated by distances of the order of hundreds of miles so that simultaneous observations are made, it is possible from a consideration of the change in distance between two times, i.e. the integrated Doppler shift, to deduce the position of the satellite at any time. The instantaneous Doppler shifts at any time then give the satellite velocity and the transmitter frequency. Another method for solving this case is given by Carrara et of. This
620
E. GOLTON
type of multistation observation is necessary when the orbit has to be found quickly, possibly in a fraction of a revolution. Another improvement in this classification is to use a transponding satellite and a ground transmitter, thus eliminating an unknown and variable satellite transmitter frequency from the solution. This case and a related solution system has been discussed by Patton and Richard.t4) When long term tracking is the main consideration, it would appear to be more profitable to place observing stations around the globe rather than have several relatively close together, so that several sections of the orbit can be observed at frequent intervals. As an example, correction coefficients have been calculated at 1 min. intervals for an orbit passing near to a station at 51”N latitude and subsequently near to an equatorial station, allowing for the seven variables AS&, Ai, At,, Ae, AT,, AU, AC. Taking one equation from each station, the angle between the corresponding hyperplanes was calculated, this being repeated for all possible pairs, The smallest angle was 30”, most being near to 90”. This is a considerable improvement over the one station example mentioned earlier, where the smallest angle was 12”, and several were small. A worldwide set of Doppler observing stations has been in use for tracking the American Transit navigational satellites G) to very high accuracy, often better than 1 km, on a routine basis, proving beyond doubt the feasibility of high accuracy Doppler tracking. The necessary equipment can be of quite limited extent, and simple aerials can be used, the complete system being of moderate cost compared with high accuracy radar or interferometer systems. Another major advantage of Doppler tracking over interferometer systems is that observations can be made at any satellite frequency with relatively minor equipment changes. The major disadvantage is the requirement of a high stability frequency source in the satellite, but suitable sources are now available at low cost. The equipment inusefor Doppler frequency measurements at the Radio Research Station, Slough, England, has been described by Henderson, (Q later developments will be given in a future memorandum. The orbital correction method described in this paper has been applied to passages of certain satellites, to obtain the satellite position to sufficient accuracy to facilitate the analysis of ionospheric propagation effects soon after the event. Although the use of observations of the Doppler effect from various combinations of stations and passages has been described, the single experimenter will probably only be interested in using observations from a single station. It is suggested that it will probably be adequate to use values for the slowly varying elements such as e, i, LL)~ provided by the prediction services of the world tracking networks. The elements O,, t, can be corrected, and the period obtained either as a correction using different passages in the solution, or from the t, values of individual passage solutions. Only the correction coefficients for the elements to be corrected need be computed, of course. In this way it should be possible to obtain the satellite position in the region of the observer within 20 km for satellites having fairly stable orbits. IONOSPHERIC
EFFECTS
In the previous sections the Doppler shift in free space has been considered. We must now consider briefly the practical case, where the earth’s ionosphere affects the radio wave propagated from the satellite. We will correct the observed Doppler shift to a free-space value, this process being effectively part of the data processing, rather than modify the orbital correction scheme to cope with the observations directly. It is the phase refractive index which needs to be considered in a ray treatment of C.W. satellite transmissions. In the ionosphere the phase refractive index is less than unity for a
THE USE OF THE DOPPLER EFFECT TO DEDUCE
AN ACCURATE POSITION
621
transmitted wave, giving a speed greater than that of light, and an increased wavelength. This gives rise to a phase path between satellite and observer which differs from the free space value, the path also deviating from a straight line. The latter effect is usually the smaller, and a good first order correction can be obtained by considering a straight line path. At frequencies appreciably above the critical frequency of the ionosphere, the local refractive 1 aN index ‘II’ is given approximately by n = 1 - - - where N is the electron density, f is the 2 f2
signal frequency and a is a constant. The total number of cycles of thewave between observer and satellite is given by: c, = ’-
4, di
c s0
taken along the slant path, element df. The time rate of change of this is simply the observed Doppler shift, Df, which after having substituted for IZbecomes D
I
z(D)
f0
_f..k.d
where (Dr)o is the free space Doppler shift. Doppler shift is proportional to the electron proportional to frequency. Because of this two related frequencies, e.g. f and pf, to follows :
c 2f2
‘Ndl
dt s o
Thus the deviation of observed from free space content JN d~along the slant path and inversely Iatter effect it is possible to use transmissions on eliminate the first order ionospheric effect as
Here (AD)f is the differentia1 Doppler shift, a quantity which is easily measured experimentally, being the beat frequency between Doppler shifts on the two frequencies after having removed the difference ratio p. Evaluated at the lower frequency fit is given by: (47) In this way a good first order approximation to the free space Doppler shift is obtained by observation of the Doppler shift on two related frequencies. To obtain a second order approximation it is necessary to consider the deviation of the ray path from a straight line, and effects due to the geomagnetic field. A knowledge of the jonospheric critical frequency is of considerable help, enabling a comparison of the observed differential Doppler frequency to be made with a theoretical value computed from a model ionosphere, including the above mentioned effects. Once the ray path at each frequency has been found, the first method of computing the satellite velocity can be used to find the component along each ray path and thus the exact differential Doppler shift. Typical values for the deviation of the observed Doppler shift from the free space value are di%cult to describe, since they vary during a satellite passage, and with the time of day, but as an indication of magnitude involved the effect at 50 MC/S during daytime at low elevation can be several cycles per second, compared with the Doppler shift of ~1 Kc/s. At night the value is smaller by a factor of 3 to 5 due to the decreased electron content. At this signal frequency, together with a higher frequency, the first order correction process gives
622
E. GOLTON
the free-space value within a fraction of one cycIe per second. The improvement obtained using the second order correction is not yet certain, since the validity of the ray path treatment is in doubt; see Guier.t7) CONCLUSION
A method of calculation has been presented which enables DoppIer frequency shift observations to be used in determining a satellite orbit. An electronic computer is necessary to carry out the method, but a small computer can be used. Factors affecting the accuracy of solution are discussed, together with a consideration of the effect of the ionosphere upon the observations. Acknowledgments-The work described was carried out as part of the programme of the Radia Research Board and is published by permission of the Director of Radio Research of the Department of Scientific and lndustriai Research, Slough, England. REFERENCES 1. W. H. GUI~R,~u~~~e, Lo&. 181, 1525 (1958). 2. D. G. KING-HELE, Pm. I&q. Sm. 247A, 49 (1958). 3. N. CARRARA, P. F. CHEOYA~CI and L. RON&HI, Space Research II, Ml, p, 215, North Holland Pub. Co., New York (1962). 4. R. B. PATTON and V. W. RICHARD, $mce Research II, 1961, p. 218, North Holland Pub. Co., New York (1962). 5. W. H. GULERand G, C. W~I~~E~~AC~, Proc. hf. Radh Eqp. N. Y. 48, 507 11960). 6. R. E. HENDERSON,&if. Comm. Electron. 8, 506 Cf961). 7. W. H. &JIER, Proc. Inst. Rad~a Engrs. N. Y. 49, 1680 (1961).
APPENDIX Ba = cos 4, cos (0 -
Sz) cos (Q + Y) -t sin (~0 + Y) x [sin #e sin i + cos #a cos i sin {O -
a)]
B, = sin (li,sin i cos (CO-+ V) - cos &b,[sin(0 + Y)cos (B - Q)
- cos i cos (ut + Y)sin (0 - Q)] B$ = cos gb,[cos(a f P) sin (0 - Q) - sin (co + Y) cos i cos (0 - Q)] C, = BP
C, = sin (w -!- y)[sin #@cos i - cos +@sin i sin (a - ft>f
C, = B3
c, = c, cos cfo
C, = -cos
C, = -Bz
#&in (03 + V) sin (0 + Q) -I- cos (fB + Y}cos (8 - fz) cos if
C, = cos (cu + Y)[sin +a cos i - cos 4, sin i sin (6 - o-2)] c, = --sin C, I
A, =
(I# + Y) cos i sin (6 -
sin (cu i; v) cos (@ -
- D,r
SS(l - e”)
cos (co + v) cos (B -
Q)
Sz} sin i
pbB& %?(I e”) (1 i- e cos V) + (f? -/- cos (I
e-S
Q) -
It(r -
vr$2&
. .
s
3e 2cos Y 1 + e cos Y -I- (I - e2)i
(
THE USE OF THE DOPPLER EFFECT TO DEDUCE AN ACCURATE
POSITION
tir,rB, 2e sin v
A, =
~(1 + e cos v) - Cd& X yr+
-CID, A, = ~ s3
&w y-
(h + 4 + 5y
-c~~(m+s)+c~~(8-n,,0,~,
- GD, A,=-.-.-s3
c+b
+ Y) + c*y
- fl) cos $”
4734
A* = ~
s3
f
Bz
‘o’(h+Y)+c+CI-n, $
D, = s2L + s2ro(d - fi)C, - s2ro(cb + P)B3 - sir + sir,,B, D, = s2r,,i - sSr,r D, = s2r - s2r,,B2.
_ $)
623
Space Sci., 1962, Vol. 9, pp. 607 to 623.
Pergamon
Press Ltd.
Printed
in Northern
Ireland
THE USE OF THE DOPPLER EFFECT TO DEDUCE AN ACCURATE POSITION FOR AN ARTIFICIAL EARTH SATELLITE E. GBLTON Department of Scientific and Industrial Research, Radio Research Station, Ditton Park, Slough (Received 4 June 1962)
Abstract-Methods by which observations of the Doppler shift of radio signals from. earth satellites may be used to determine the satellites position are discussed. Two methods for the calculation of the theoretical DoppIer shift are given, followed by a method of linear dif%rential orbit correction in terms of Keplerian orbital elements. Limitations on the possibility of solution are discussed, together with the use of single and multiple receiving sites, and limitations due to experimental measurement accuracy. Methods by which corrections for the effect of the ionosphere may be made are briefly considered. INTRODUCTION A knowledge of the position of a satellite at a given time is of impo~an~ both to experimenters interested in the reduction of telemetry data and to those deducing ionospheric information from signal propagation studies. Unfortunately past experience has shown that it is usually many months before reasonably accurate ephemeris becomes available, if at all. This paper describes how observations of the well-known Doppler effect may be used to obtain positioned data of sufficient accuracy, soon after the satellite passage. Little appears to have been published hitherto on the accurate reduction of Doppler shift observations. Since time variation of the magnitude of the effect observed from the rotating earth is a unique, though nonlinear, function of the orbital elements of the satellite, it follows that it should be possible to determine the orbital elements from a set of observations of the Doppler effect. The practical possibility of doing this appears to have been first described by Guier. (l) Otherwise it has been customary to divide the plotted Doppler curve by some geometrical method to obtain a “time of nearest approach” and to use the maximum slope to obtain a “slant range”. However, these calculations are based on tremendously simplifying assumptions, and consequently the time obtained is usually in error by a few seconds, and the range in error by an appreciable fraction. Whilst these values may be very useful in predicting future passages of a satellite, they are of little use in studying the detailed orbital motion of the satellite or in fixing its position accurately. One method by which the orbital elements might be corrected by observation of the Doppler shift is to compare the measured values with theoretical values derived from a trial set of elements, and calculate the mean square deviation: this is carried out for many different sets of elements, and the set giving the smallest mean square deviation is nearest to the correct set of elements. In practice this could mean searching for a minims with perhaps ten values of each of seven elements, i.e. 10’ sets, the Doppler shift being calculated perhaps 200 times for each set of elements. Thus although this method is probably the best for solving a nonlinear system of the type under discussion, it is only of practical use where the largest and fastest electronic computers are available, Another method which might have a more widespread application is to derive a system 607
608
E. GOLTON
of linear equations which relate a variation in Doppler shift, or the equivalent range rate S:, to variations in the orbital elements, xi, of the form:
Differences between an observed Doppler shift and values calculated from a trial set of elements might then be used in an equation solution to give corrections to the etements. N A
Ascending FIG. 1. DIAGRAM
node
SHO~NGELEMENTSOFSA~~LITEORB~.
Because the function S(x,) at any instant is in fact nonlinear in x6 the derived corrections will not be exact, and indeed the system may not provide useful corrections at all. Because of this the nonlinearity of the coefficients in a linear equation system must be investigated, together with the effect of errors in the observed Doppler shift values, which will also limit the accuracy of the first method of solution described. In the present paper expressions for the theoretical Doppler shift will be developed in terms of Keplerian orbital elements defined at the ascending node, assuming a constant rate of change of perigee and ascending node positions. From these will be developed a differential correction scheme, and the necessary investigations mentioned in the previous paragraph will then be carried out. It is obvious that in order to attempt useful Doppler shift measurements one must have a stable transmitter frequency at the satellite, and a stable reference frequency at the ground observing site. At the present time crystal controlled oscillators with stabilities of a few parts in lOlo are in use in certain satellites. The observed Doppler shift will not be quite the same as that observed under free space conditions, because of the effect of the ionosphere: however, the first order effect may be removed by the transmission of two phase locked signals from the satellite, and the second order effect by a study of the ionosphere. These processes will be outlined at the end of the paper. DEFINITIONS
A geocentric co-ordinate system is used, the reference directions being the northward axis of the earth and the first point of Aries. Angles measured eastward are considered
THE USE OF THE DOPPLER EFFECT TO DEDUCE AN ACCURATE POSITION
609
positive. A satellite moving from west to east is considered primarily, although the analysis of the complementary case is almost identical. The satellite motion is assumed to be in the form of a Kepler ellipse, which is allowed to precess about the earth’s axis, and the major axis of which is allowed to rotate in the mean orbit plane. The last two secular additions represent the main perturbations of the orbit due to the oblateness of the earth, e.g. King-Hele.(2) The main symbols used are defined below; others will be introduced as required. Some of the symbols are shown in Fig. 1. The words “ascending node” will be abbreviated to “A.N.” hereafter. t, time of the A.N. T absolute time of the satellite at point S. i inclination of orbit To satellite anomalistic period of revolution e eccentricity angle from A.N. to perigee at time t, Right ascension of A.N. at time t, angle from A.N. to perigee at time T Right ascension of A.N. at time T rate of change of CO rate of change of Sz Right ascension of Greenwich at time T = 0 earth’s sidereal rate of rotation Right ascension of observer at time T time from A.N. to perigee \ time from perigee to S I angle from perigee to S observers radial distance geocentric observers latitude values at P observers longitude radial distance of satellite geocentric latitude of satellite values longitude of satellite Eccentric anomaly between node and perigee Eccentric anomaly between perigee and satellite geocentric angle PBS between observer and satellite In all cases time is measured in terms of the mean solar day. to represent differentiation with respect to time. EQUATION
DESCRIBING
The basic equations
required
The dot notation
will be used
THE SATELLITE MOTION IN A KEPLER ELLIPSE
are as follows:
(1)
(2)
610
E. GOLTON
t = 2 (E - e sin E) sin E = (1 - e2)* *
sin Y (1 + e cos Y)
(4)
- tan E/2 27r (1 i e cos v)2 * = z * (1 - ,213
(6)
We have chosen to work in terms of the period TOrather than the semi-major axis length a, these being connected by the equation TO2= 4n2a3 P where ,u = GM, the product To solve these equations The observation is made The perigee angle is thus
of the gravitational constant and the mass of the earth. we must first consider time. at a time (T -- t,) after the A.N. o = co0+ c&(7’- t,)_ Equation (4) will now give Em. sin w sin E, = (1 - e216(1 + e cos w)
(7)
Substitution of E, into equation 3 gives: t, = 2 (E, - e sin Em) We can now find the time of the satellite from perigee. This is: t = T -
t, -
t,.
(9)
Equation 3 must now be solved using the value t found above to give E. This may be done by a process of iteration since the equation is monotonic. Equation 5 is now used to find the angle v. From this we can calculate r, i, and i by using equations 1,2, and 6. The position of the satellite must now be found. The Right Ascension of the A.N. at time T is CI = Q2,+ h(T - I,): that of Greenwich is (0, + 6 . T). By reference to Fig. 2 we see that the satellite has moved in longitude through an angle A, given by:
Also
tan 2 = tan (w + P) *cos i.
(10)
sin #s = sin ((0 + Y)* sin i.
(11)
Thus the satellite is at geocentric latitude 4s and longitude 0, = Q + Iz - (0, + 8.T). The difference in longitudes between satellite and observer is 8,-08,=a-e+n
(12)
The geocentric angle y between satellite and observer is given by: cos y = sin #@- sin A t_ cos & cos & cos (0, - 0,)
(13)
THE USE OF THE DOPPLER EFFECT TO DEDUCE AN ACCURATE POSITION
Kepler
611
ellipse
Equator
FIG. 2. DIAGRAMSOF ORBIT GEOMETRY.
The distance s between satellite and observer is given by: s2 = rS2 + y02- 2r,r, cos y CALCULATION
OF THE
THEORETICAL
DOPPLER
(14) SHIFT
(A) Method 1. We will first calculate the theoretical Doppler shift by finding the velocity of the satellite relative to the earth, and then taking the component in the direction of the observer. The first part of this method is essential for calculating the Doppler shift along a ray path through the ionosphere and for eliminating the second order ionospheric effect. With reference to Fig. 2 the velocity components relative to the earth are: 3, outwards, r,(C + h) cos B northwards, rs - (i’ + h) sin B - rs(8 - 0) cos & eastwards.
and cos B = tan &/tan (CO+ Y), sin 3 = sin @in (w + Y). Thus the satellite speed is given by V in the following equation, V2 = i,2 + r,“[(i, + 61)~+ (0 - fi)2 * cos2 $, - 2(i + 15) * (6 - 0) sin B * cos #,I (1%
and is directed at elevation a and bearing B given by: tan M.= i,/r,[(G + I%)”+ (6 - h)2 * co? y& - 2(2; + c.!I)(~- $) sin B - cos C&I* (16) tan /3 = [(+ + &) sin B - (6 - $) cos &]/(i, + h) cos B
(17)
With reference to Fig. 3 we now require to find the zenith angle e, and the bearing b, of the observer from the satellite. These are given best by: tanc,=siny/(z-cosy) cos b, = [sin $O - sin & +cos yJ/cos q3,* sin y.
(18) (19)
E. GOLTON
612
0 FIG.
3.
DIAGRAM
FOR
CALCULATION
OF SATELLITE TO OBSERVER
These four angles may now be combined cos x = -cos The actual Doppler
LINE
VELOCITY
COMPONENT
ALONG
SATELLITE
S-P.
to give the angle x between
e, * sin tl + sin e, * cos cc * cos (b, -
V and SP. 8)
(20)
shift is given by: II’-’
f
v.cosx
c
where f is the transmitter frequency refractive index (n = 1 at present).
and c is the speed of light, and n is the local phase
(B) Method 2. The free space Doppler shift may be calculated more easily by considering the rate of or by elimichange of distance to the satellite. From a separate geometrical calculation nating #Q from and substituting equation 12 into equation 13 we obtain : cos y = cos (v + Co)* cos $0 * cos (e -
sz)
+ sin (Y + 0) - [sin 9S0- sin i + cos +O * cos i * sin (f3 -
L?)] (21)
At this stage it will be convenient to introduce letters to represent various terms, the first being B, = cos y as in equation 21. Also the subscript ‘s’ attached to r, will be dropped where possible. Differentiation of equation 14 gives : (sS) = r . J -
r,iB, - rOr. ; @A
613
THE USE OF THE DOPPLER EFFECT TO DEDUCE AN ACCURATE POSITION
Differentiation
of B, and substitution (sS) = rL -
into the above yields
ro3Bz -
r,,r(ci, + +)B, + r&6’ -
fi)B4
(22)
where B, = sin $O sin i. cos (co + v) -
cos $,[sin (w + Y) cos (0 - a) - cos i * cos (o + Y) sin (0 -
Q)]
(23)
and B4 = cos +,,[sin (f3 -
Q) * cos (co + Y) -
cos (0 -
Cl) * sin (w + Y) * cos il
Equation 22 in conjunction with equation 14 enables (S) to be calculated, of change of distance to the satellite, and since
(24)
this being the rate
(S) = -V*cos~ the theoretical
Doppler
shift may be calculated THE DIFFERENTIAL
directly.
CORRECTION SCHEME
A linear relation between an infinitesimal variation A(.+) of the Doppler shift and infinitesimal variations Ax of the orbital elements must now be found. This will be done by partial differentiation of expressions found above. We start with the identity: ~3 - A(S) = s2 . A($ From
equation
-
(8s) . 4 - A(s2)
(25)
14 we obtain: + * A(s2) = r * Ar -
And from equation
r, * B2 * Ar -
r * r,,AB2
(26)
22 :
A@) = Ar * [3 + r,,(d -
$) * B, -
r,-,B2] -
+ Ai[r -
+ AB, * r - r,,(8 -
r,(h + ?‘) - B3]
Afi - r - r, * B4 0) -
(27)
AB, . r - r,,(b + 1’>-
AB, . r,,i
(A& + A+)r * r,,B,
It will be noticed that a variation of ci, and 2;has been allowed, although because these terms are generally small a useful correction may not be obtained. Substitution of equations 26 and 27 into equation 25 gives: s3A(S) = Ar * D, -
AB, . s2 * r - ro(ci) + ti)
AB, * D, + A3 * D, -
+ AB4 . s2 * r * r,(6 -
0) -
Ah * s2 * r - r,, * B4 -
(Ah + Ati)s2 . r * r,, * B,
(28)
where : D, = s2 * i + s2 . r,(8 -
0) * C, -
D, = s2 * r, * i -
The variations Appendix.
s2 * r,,(h + ?j)B, -
(si)r,, - r ;
D, = s2 * r -
(sS) * r + (sS) . r, * B, s2 - r,, * B,.
of the terms B,, B,, B,, are given below, the coefficients
I
(29)
C, being given in the
AB, = C, . AQ + C, * Ai + C, * (Aw + AY) AB, = C, . AQ + C,Ai + C, - (Aw + Av)
(30)
ABd = C, * AG + C, * Ai + C, * (Aw + Av) The variations
of r, t and ?; must now be expressed
as variations
of the orbital
elements.
E. GOLTON
614
This may be carried finally : Ar -= Y
out by partial logarithmic differentiation of equations 1, 2, 6, giving
2 AT, he e + cos Y ---e sm v Av e+ + 1 -+ e cos v 1 - ez 1 + e cos v 3 T, A3 AT, Ae IL= -$cotv*Av + i 3T@ e(l - e2) 3e 2e sin v . Av -?+Ae 12+~~Oysv+-_ 1 e2 (1 + e cos Y) i 0
A+ -= 4
(31)’ (32) (33)
Substitution of these expressions into equation 28 together with equations 29 and 30 give: A(S) = A, *AT, + A, *Ae + A, * A$ + A, * Aci, + A, * AY -+ A, *ASZ+ A, *Ai + A, * AU
(34)
where tbe A, are given in the Appendix. AI1 the above variations are of orbital elements except Av, which besides being a function of Ae and AT,, is a function of At, and Ao. The required dependence may be derived from equations 3, 4. Partial differentiation gives : At=g(l-e*cosE)*AE+$-*At-g*sinE*Ae
(35)
0
and (e -t cos v)
cotE*AE=(l
(36)
+ e cos 7.~)
Elimination of AE now gives: Av =f(v)
* At --f(v).
G . AT, + he. 0
where sin v * cot E * (1 + e +cos v) * 2~ f(v) = (e + cos Y)- (1 - e * cos E) * To
(38)
Equation 37 expresses a variation of Y in terms of variations of e, To, and t. Now t is dependent upon t, and t, as given by equation 9, so that At = -AT,
- At,
(39)
We can write down the dependence of At,,, immediately by analogy to equation 37, viz, : AW =f(~)
e At, -f(u)
. k AT0 + Ae - ie2
+ f(w) * 2 * sin Em)
(40)
These equations result from the dependence of the time taken to travel from ascending node to perigee upon eccentricity, period, and perigee position; this time in turn determines the time available to travel from perigee to the point of observation, at fixed time T. Substitution of equations 39, 40 into equation 37 gives: Av = -f(v)
f(v)
- At, - -
f(v)
aAGO- AT, - r o * 0 + CJ
.f(co)
+ f(v) * 2 (Sin-E + sin E,)
I
(41)
615
THE USE OF THE DOPPLER EFFECT TO DEDUCE AN ACCURATE POSITION
Here t + t, is simply the time from node to point of observation, T - t,,. Equation 41 may after this has been carried out, Q and 03 must be now be substituted into equation 34: replaced by their component values : Cl = Q2, + $(T
-
t,)
ro = CII~+ cb(T -
t,)
(42)
are :
Their variations
AL’=ACI,+(T-t&A&_&At, + (T -
AU = ho, Substitution
of equations
A(i)=
A,-Afo*(T5T,, -
[A&)
+
AZ +
t) n
(1
AT, + [A3 + A&T -
_
sin v + f*
$)
t,>l
*Ah + A,’ . Aq,
1
2
(sin
E
+
sin E,)
d!zf!
A,-_A
equation:
sin w
.f
A,f(v)
correction
t,)] . Aa + A, * AQ,
. At, + [A* + A,‘@- -
+ A,‘h]
1
.-.
1.
+
where A,’ =
C%. At,
41 and 43 into 34 give the final differential
+ A,0
A5
(43)
t,) . A& -
1
* Ae + A, - Ai
(44)
5 f(4’ ORBITAL PERTURBATIONS
A constant rate of change of perigee and ascending node positions are the only perturbations which have so far been included in the correction scheme. Other perturbations are generally small, and inclusion of theoretical expressions for them in the correction scheme would complicate and lengthen it considerably. However, it is quite simple to include theoretical expressions for the perturbations in calculating the elements at a given time, to be used in the correction scheme. For example, one perturbation not included previously is that of the Keplerian relation between period and semi-major axis, which is affected by the oblateness of the earth and is mainly a function of inclination. This can be included by simply modifying the value ,u = GM as a function of inclination. INVESTIGATION
The coefficients out to find. many times, values have satellite past
of the variations
OF THE SOLUTION SYSTEM
of elements
are the partial
derivatives
g,which
we set
Although most are long expressions, simple arithmetical operat;ons are used and the coefficients can be found very quickly on a small computer. Typical been plotted at one minute intervals in Figs. 4, 5 for a descending passage of a an observer at 51”N latitude, using the following elements:
T, 101.66
0.2:14
66!77
CI 352.0428
fi -0901775
,9;0,33
-0;005
t 73”0.7
Time is measured in minutes, and distance in thousands of kilometres. The correction coefficients are thus measured in units of lo3 km/min, and, where applicable, radians. They 2
616
E. GOLTON
show the variations in Doppler shift caused by an infinitesimal change in a particular
as
orbital element. For example, the curve of - in effect shows the slope of the Doppler curve, at, reaching the greatest (negative) value at the point of nearest approach. The coefficients were checked by finding the change in Doppler shift caused by small changes in each element. A3 When these hxi values were plotted, the values at Ax, = 0 were within O-1per cent of the
as
theoretical x
values in all cases, thus checking the correctness of the previous theory.
The soluti:n of a system of linear equations of the form:
. AxZj =
Asi
is affected by the accuracy of both the RHS and the coefficients of the LIB, generally to a greater extent the more variables we include. For a given error in a coefficient the effect will be more serious the worse the conditioning of the set of coefficients, Although we cannot alter the conditioning, it is useful to obtain a measure of it. This may be done by making use of the fact that each linear equation represents a hype~lane in n-dimensional space, and the normalized coefficients are the direction cosines of the vector normal to the hyperplane. The angles between these normal vectors are a direct measure of the conditioning. If the angles are small, the equations are badly conditioned, but if the angles are all near 90” the equations are very well conditioned. A given change in a coefficient results in a tilting of the hyperplane: In a perfectly conditioned set the largest change of one variable occurs when the hyperplanes are parallel to the axes, and is given by Z/zbn2 - b,2 4 sin 19where the b, are the normalized RHS and i is the equation with the chanied coefficient, and 8 is the change in direction of its normal. When the hyperplanes are not parallel to the axes, changes in all variables occur when one plane is tilted, but are of smaller magnitude than in the previous case. The degradation occurring when the hyperplanes are not perpendicular to each other is related to the cosecants of the angles between pairs of normals. The equations with coefficients plotted in Figs. 4,5 were investigated for several numbers of the variables. In particular with seven variables, AT,, Ae, Ai, AC?;,,Ao,, At,, AC, the last being the error in assumed transmitter frequency, and seven equations applying at one minute intervals, the smallest angle between any two hyperplanes was 12”, and the largest 90”. Thus a degradation of only about 5: 1 occurs, and a computer solution will not be limited by rounding errors. In the case of four variables, A!&, Ai, At,, AC, the smallest angle was 7 degrees, rather worse than before. We must now consider the linearity of the correction coefficients. For a given finite set of corrections to the elements and perfect RHS, the appropriate coefficients will not be quite 6% the same as the theoretical coefficients d;r. This is because the function $(x2) is nonlinear in X, at any given instant of time, The nonlinearity has been calculated by finding the change in Doppler shift for a finite change in each element, and expressing its deviation from the theoretical correction as a percentage. These values have been plotted below the corresponding correction coefficients in Figs. 4, 5. A constant percentage nonlinearity would be
THE USE OF THE DOPPLER EFFECT TO DEDUCE AN ACCURATE POSITION
al
g”s NON + I
1
617
618
E. GOLTON
THE USE OF THE DOPPLER EFFECT TO DEDUCE AN ACCURATE POSITION
619
of little consequence, since if all the coefficients of one variable in an equation solution were for example 2 per cent too small, the answer for the variable would simply be 2 per cent too large. However, the curves drawn in Figs. 4, 5, show considerable variation in magnitude and sign. The nonlinearity is dependent of course on the magnitude of the variation of the error in the orbital element, and will set a limit to the maximum permissible error in the trial set of elements from which a sensible convergent solution can be obtained. The RHS of the equation system contain observational errors and ultimately limit the accuracy of orbital element corrections once a convergent sequence has been obtained. The usual way of minimizing random effects is to use more equations than variabIes and to normalize the system. This in theory gives a least-squares fit to the RHS, assuming perfect coefficients, but in our case normalizing the system has the added advantage of smoothing out the effects of the nonlinearity of the coefficients, and increasing the range of convergence. This can be simply explained by the fact that errors in the coefficients of an equation are equivalent to a certain error in the RHS. As an example, it has been found possible to obtain convergence with 4 variables, AI&,, Ai, At,, AC using 8 equations. However, a least-squares solution is not reliable in this case, and this is reflected by the fact that a variation of any one RHS by O-1c/s in the measured Doppler shift at 50 MC/S gives rise to changes of the order of 0.1” and 4 set in time in the elements. However, when the solution was repeated using 47 points equally distributed along the Doppler curve, the same error in any one point gave rise to changes of the order of 0.001 degree and 0.005 set, an improvement of the order of 100 times. This accuracy is not of course always significant, since the other elements for which no solution was attempted are not always known to equivalent accuracy. Regular, rather than random errors impose a rather more serious limitation. In the previous example, using 47 points, a linear drift in transmitter frequency of O-1 c/s over 10 min at 50 MC/S results in corrections of about 0*03” and 0.15 sec. Such Doppler errors might also arise from imperfections in correcting ionospheric effects. At the present stage it has not proved possible to solve usefully for all the orbital elements independently using observations of one passage from one station, partly because of the small range of convergence and partly because of the limited accuracy of observations of the Doppler shift. So far only observations of one passage from one station have been considered. The solution becomes very much easier when more observing stations are available. With observations of several passages from one station it becomes possible to include the period correction AT, with AL&,,Ai, At,, AC, and obtain a strong solution. Usually the difficult elements to correct are perigee position and eccentricity. If the single station is situated at a lower latitude than apex, the observation of ascending and descending passages will cover two sections of the orbit, and useful corrections to these elements can be made. The rate of change 0 becomes important in this case, but rather than try to solve for it as an independent variable it appears more useful to obtain Q2,from each passage separately and then use the mean fi in a combined solution. When several observing sites are available, several different sections of the orbit may be observed, and it appears that all the elements can be corrected independently. If four stations are separated by distances of the order of hundreds of miles so that simultaneous observations are made, it is possible from a consideration of the change in distance between two times, i.e. the integrated Doppler shift, to deduce the position of the satellite at any time. The instantaneous Doppler shifts at any time then give the satellite velocity and the transmitter frequency. Another method for solving this case is given by Carrara et of. This
620
E. GOLTON
type of multistation observation is necessary when the orbit has to be found quickly, possibly in a fraction of a revolution. Another improvement in this classification is to use a transponding satellite and a ground transmitter, thus eliminating an unknown and variable satellite transmitter frequency from the solution. This case and a related solution system has been discussed by Patton and Richard.t4) When long term tracking is the main consideration, it would appear to be more profitable to place observing stations around the globe rather than have several relatively close together, so that several sections of the orbit can be observed at frequent intervals. As an example, correction coefficients have been calculated at 1 min. intervals for an orbit passing near to a station at 51”N latitude and subsequently near to an equatorial station, allowing for the seven variables AS&, Ai, At,, Ae, AT,, AU, AC. Taking one equation from each station, the angle between the corresponding hyperplanes was calculated, this being repeated for all possible pairs, The smallest angle was 30”, most being near to 90”. This is a considerable improvement over the one station example mentioned earlier, where the smallest angle was 12”, and several were small. A worldwide set of Doppler observing stations has been in use for tracking the American Transit navigational satellites G) to very high accuracy, often better than 1 km, on a routine basis, proving beyond doubt the feasibility of high accuracy Doppler tracking. The necessary equipment can be of quite limited extent, and simple aerials can be used, the complete system being of moderate cost compared with high accuracy radar or interferometer systems. Another major advantage of Doppler tracking over interferometer systems is that observations can be made at any satellite frequency with relatively minor equipment changes. The major disadvantage is the requirement of a high stability frequency source in the satellite, but suitable sources are now available at low cost. The equipment inusefor Doppler frequency measurements at the Radio Research Station, Slough, England, has been described by Henderson, (Q later developments will be given in a future memorandum. The orbital correction method described in this paper has been applied to passages of certain satellites, to obtain the satellite position to sufficient accuracy to facilitate the analysis of ionospheric propagation effects soon after the event. Although the use of observations of the Doppler effect from various combinations of stations and passages has been described, the single experimenter will probably only be interested in using observations from a single station. It is suggested that it will probably be adequate to use values for the slowly varying elements such as e, i, LL)~ provided by the prediction services of the world tracking networks. The elements O,, t, can be corrected, and the period obtained either as a correction using different passages in the solution, or from the t, values of individual passage solutions. Only the correction coefficients for the elements to be corrected need be computed, of course. In this way it should be possible to obtain the satellite position in the region of the observer within 20 km for satellites having fairly stable orbits. IONOSPHERIC
EFFECTS
In the previous sections the Doppler shift in free space has been considered. We must now consider briefly the practical case, where the earth’s ionosphere affects the radio wave propagated from the satellite. We will correct the observed Doppler shift to a free-space value, this process being effectively part of the data processing, rather than modify the orbital correction scheme to cope with the observations directly. It is the phase refractive index which needs to be considered in a ray treatment of C.W. satellite transmissions. In the ionosphere the phase refractive index is less than unity for a
THE USE OF THE DOPPLER EFFECT TO DEDUCE
AN ACCURATE POSITION
621
transmitted wave, giving a speed greater than that of light, and an increased wavelength. This gives rise to a phase path between satellite and observer which differs from the free space value, the path also deviating from a straight line. The latter effect is usually the smaller, and a good first order correction can be obtained by considering a straight line path. At frequencies appreciably above the critical frequency of the ionosphere, the local refractive 1 aN index ‘II’ is given approximately by n = 1 - - - where N is the electron density, f is the 2 f2
signal frequency and a is a constant. The total number of cycles of thewave between observer and satellite is given by: c, = ’-
4, di
c s0
taken along the slant path, element df. The time rate of change of this is simply the observed Doppler shift, Df, which after having substituted for IZbecomes D
I
z(D)
f0
_f..k.d
where (Dr)o is the free space Doppler shift. Doppler shift is proportional to the electron proportional to frequency. Because of this two related frequencies, e.g. f and pf, to follows :
c 2f2
‘Ndl
dt s o
Thus the deviation of observed from free space content JN d~along the slant path and inversely Iatter effect it is possible to use transmissions on eliminate the first order ionospheric effect as
Here (AD)f is the differentia1 Doppler shift, a quantity which is easily measured experimentally, being the beat frequency between Doppler shifts on the two frequencies after having removed the difference ratio p. Evaluated at the lower frequency fit is given by: (47) In this way a good first order approximation to the free space Doppler shift is obtained by observation of the Doppler shift on two related frequencies. To obtain a second order approximation it is necessary to consider the deviation of the ray path from a straight line, and effects due to the geomagnetic field. A knowledge of the jonospheric critical frequency is of considerable help, enabling a comparison of the observed differential Doppler frequency to be made with a theoretical value computed from a model ionosphere, including the above mentioned effects. Once the ray path at each frequency has been found, the first method of computing the satellite velocity can be used to find the component along each ray path and thus the exact differential Doppler shift. Typical values for the deviation of the observed Doppler shift from the free space value are di%cult to describe, since they vary during a satellite passage, and with the time of day, but as an indication of magnitude involved the effect at 50 MC/S during daytime at low elevation can be several cycles per second, compared with the Doppler shift of ~1 Kc/s. At night the value is smaller by a factor of 3 to 5 due to the decreased electron content. At this signal frequency, together with a higher frequency, the first order correction process gives
622
E. GOLTON
the free-space value within a fraction of one cycIe per second. The improvement obtained using the second order correction is not yet certain, since the validity of the ray path treatment is in doubt; see Guier.t7) CONCLUSION
A method of calculation has been presented which enables DoppIer frequency shift observations to be used in determining a satellite orbit. An electronic computer is necessary to carry out the method, but a small computer can be used. Factors affecting the accuracy of solution are discussed, together with a consideration of the effect of the ionosphere upon the observations. Acknowledgments-The work described was carried out as part of the programme of the Radia Research Board and is published by permission of the Director of Radio Research of the Department of Scientific and lndustriai Research, Slough, England. REFERENCES 1. W. H. GUI~R,~u~~~e, Lo&. 181, 1525 (1958). 2. D. G. KING-HELE, Pm. I&q. Sm. 247A, 49 (1958). 3. N. CARRARA, P. F. CHEOYA~CI and L. RON&HI, Space Research II, Ml, p, 215, North Holland Pub. Co., New York (1962). 4. R. B. PATTON and V. W. RICHARD, $mce Research II, 1961, p. 218, North Holland Pub. Co., New York (1962). 5. W. H. GULERand G, C. W~I~~E~~AC~, Proc. hf. Radh Eqp. N. Y. 48, 507 11960). 6. R. E. HENDERSON,&if. Comm. Electron. 8, 506 Cf961). 7. W. H. &JIER, Proc. Inst. Rad~a Engrs. N. Y. 49, 1680 (1961).
APPENDIX Ba = cos 4, cos (0 -
Sz) cos (Q + Y) -t sin (~0 + Y) x [sin #e sin i + cos #a cos i sin {O -
a)]
B, = sin (li,sin i cos (CO-+ V) - cos &b,[sin(0 + Y)cos (B - Q)
- cos i cos (ut + Y)sin (0 - Q)] B$ = cos gb,[cos(a f P) sin (0 - Q) - sin (co + Y) cos i cos (0 - Q)] C, = BP
C, = sin (w -!- y)[sin #@cos i - cos +@sin i sin (a - ft>f
C, = B3
c, = c, cos cfo
C, = -cos
C, = -Bz
#&in (03 + V) sin (0 + Q) -I- cos (fB + Y}cos (8 - fz) cos if
C, = cos (cu + Y)[sin +a cos i - cos 4, sin i sin (6 - o-2)] c, = --sin C, I
A, =
(I# + Y) cos i sin (6 -
sin (cu i; v) cos (@ -
- D,r
SS(l - e”)
cos (co + v) cos (B -
Q)
Sz} sin i
pbB& %?(I e”) (1 i- e cos V) + (f? -/- cos (I
e-S
Q) -
It(r -
vr$2&
. .
s
3e 2cos Y 1 + e cos Y -I- (I - e2)i
(
THE USE OF THE DOPPLER EFFECT TO DEDUCE AN ACCURATE
POSITION
tir,rB, 2e sin v
A, =
~(1 + e cos v) - Cd& X yr+
-CID, A, = ~ s3
&w y-
(h + 4 + 5y
-c~~(m+s)+c~~(8-n,,0,~,
- GD, A,=-.-.-s3
c+b
+ Y) + c*y
- fl) cos $”
4734
A* = ~
s3
f
Bz
‘o’(h+Y)+c+CI-n, $
D, = s2L + s2ro(d - fi)C, - s2ro(cb + P)B3 - sir + sir,,B, D, = s2r,,i - sSr,r D, = s2r - s2r,,B2.
_ $)
623