Method for parsing and calculating performance of satellite gravity field measurement by low-to-low satellite-to-satellite tracking

Abstract

A method for calculating performance of satellite gravity field measurement by low-to-low satellite-to-satellite tracking, includes: acquiring parameters of gravity satellite system; calculating an effect of satellite loads on the power spectrum of nonspherical perturbation potential, so as to obtain an degree error variance; comparing degree error variance with degree variance given by Kaula Rule, and when degree error variance equals degree variance, considering that the highest degree of gravity field measurement is obtained, calculating geoid degree error and its accumulative error, gravity anomaly degree error and its accumulative error, so as to obtain the performance of satellite gravity field measurement by low-to-low satellite-to-satellite tracking. The method is capable of evaluating gravity field measurement performance quickly and effectively, obtaining a rule of effects of the gravity satellite system parameters on the gravity field measurement performance, so as to avoid shortcoming caused by numerical simulation.

Claims

1. A method for parsing and calculating performance of satellite gravity field measurement by low-to-low satellite-to-satellite tracking, comprising steps of: step (1): acquiring at least one parameter of a gravity satellite system by low-to-low satellite-to-satellite tracking by the gravity satellite system; step (2): according to the parameter of the gravity satellite system by low-to-low satellite-to-satellite tracking, calculating an effect of a measuring error of a gravity satellite load on a power spectrum of a nonspherical perturbation potential of earth gravity, so as to obtain an degree error variance of a potential coefficient in a gravity field recovery model; step (3): comparing the degree error variance of the potential coefficient in the gravity field recovery model with an degree variance of a potential coefficient given by Kaula Rule, wherein with an increase of an order of the gravity field recovery model, the degree error variance gradually increases and the degree variance gradually decreases; when the degree error variance is equal to the degree variance, considering that a maximum valid order of the gravity field measurement is obtained, wherein the degree error variance, the degree variance and the maximum valid order are obtained by the gravity field model; step (4): according to the degree error variance of the gravity field recovery model, calculating a geoid-order error, an accumulative error, an gravity-anomaly order error and an accumulative error of the gravity field recovery model, wherein the geoid-order error, the accumulative error, the gravity-anomaly order error and the accumulative error are obtained by the gravity field model; and step (5): summarizing the valid degree of the gravity field measurement, the geoid-order error, the accumulative error, the gravity-anomaly order error and the accumulative error, so as to obtain the performance of satellite gravity field measurement by low-to-low satellite-to-satellite tracking, wherein the valid degree of the gravity field measurement, the geoid-order error, the accumulative error, the gravity-anomaly order error and the accumulative error are obtained by the gravity field model; wherein the parameter of the gravity satellite system by low-to-low satellite-to-satellite tracking comprises but not limited to at least one orbit parameter of the gravity satellite system and at least one load indicator of the gravity satellite system; wherein the orbit parameter of the gravity satellite system comprises at least one member of: a maximum valid order of the gravity field recovery N.sub.max, a gravity-anomaly order error .sub.n of an nth order, a geoid-order accumulative error of the nth order, a gravity-anomaly order error g.sub.n of the nth order and a gravity-anomaly accumulative error g of the nth order, wherein N.sub.max, .sub.n, , g.sub.n and g are obtained by the gravity field model; the parameter of the gravity satellite system comprises: a gravity satellite orbit height h and an included angle .sub.0 of satellite-to-satellite geocentric vectors, wherein h and .sub.0 are obtained by the gravity satellite system, the load indicator of the gravity satellite system comprises: an inter-satellite range change rate measurement error ({dot over ()}).sub.m, a satellite orbit determining position error (r).sub.m, a non-gravitational interference F, an inter-satellite range rate data sampling interval (t).sub.{dot over ()}, a non-gravitational interference data interval (t).sub.F, satellite orbit position data sampling interval (t).sub.r and a gravity field measurement service life T, wherein ({dot over ()}).sub.m is obtained by an inter-satellite range measurement device, (r).sub.m is obtained by a spaceborne GPS system, and the F is obtained by an accelerometer, wherein the inter-satellite range measurement device, the spaceborne GPS system and the accelerometer are conventional and all provided on the gravity satellite system; wherein the step (2) specifically comprises steps of: establishing an analytic formula of a low-to-low satellite-to-satellite tracking gravity field measurement degree error variance .sub.n.sup.2 of a potential coefficient: ${}_n^2=\frac{1}{2n+1}\underset{k=0}{\overset{n}{.Math.}}\left[{\left({\stackrel{\_}{C}}_{\mathrm{nk}}\right)}^2+{\left({\stackrel{\_}{S}}_{\mathrm{nk}}\right)}^2\right]=\frac{1}{\underset{k=0}{\overset{n}{.Math.}}\left[B_1\left(r_0,n,k,{}_0\right)+B_2\left(r_0,n,k,{}_0\right)+B_3\left(r_0,n,k,{}_0\right)\right]}\frac{1}{2n+1}{\frac{2\left(n+1\right)r_0^{2n}}{{}^2{\mathrm{Ta}}_e^{2n-2}}\left[\begin{array}{c}-D\sqrt{\left(\frac{{\left(F\right)}^2{T}_{\mathrm{ar}c}^4}{36}{\left(t\right)}_F+{\left(r\right)}_m^2{\left(t\right)}_r\right)f_r\left(n\right)}+\\ \sqrt{\frac{}{r_0}}\sqrt{\left(\begin{array}{c}\frac{{\left(F\right)}^2{T}_{\mathrm{ar}c}^2}{4}{\left(t\right)}_F+\\ {\left(\stackrel{.}{}\right)}_m^2{\left(t\right)}_{\stackrel{.}{}}\end{array}\right)f_{\stackrel{.}{}}\left(n\right)}\end{array}\right]}^2$ wherein: $\begin{array}{cc}\{\begin{array}{c}\begin{array}{c}{B}_1\left(r_0,n,k,{}_0\right)=A_n^2\left(n-1,k,{}_0\right)+\left(\frac{a_e}{r_0}\right).Math.A_n\left(n-1,k,{}_0\right)A_n\left(n,k,{}_0\right).Math.+\\ {\left(\frac{a_e}{r_0}\right)}^2.Math.A_n\left(n-1,k,{}_0\right)A_n\left(n+1,k,{}_0\right).Math.\end{array}\\ \begin{array}{c}{B}_2\left(r_0,n,k,{}_0\right)={\left(\frac{a_e}{r_0}\right)}^2{A}_n^2\left(n,k,{}_0\right)+\\ \left(\frac{a_e}{r_0}\right).Math.A_n\left(n-1,k,{}_0\right)A_n\left(n,k,{}_0\right).Math.+{\left(\frac{a_e}{r_0}\right)}^3.Math.A_n\left(n,k,{}_0\right)A_n\left(n+1,k,{}_0\right).Math.\end{array}\\ \begin{array}{c}{B}_3\left(r_0,n,k,{}_0\right)={\left(\frac{a_e}{r_0}\right)}^4{A}_n^2\left(n+1,k,{}_0\right)+\\ {\left(\frac{a_e}{r_0}\right)}^2.Math.A_n\left(n-1,k,{}_0\right)A_n\left(n+1,k,{}_0\right).Math.+{\left(\frac{a_e}{r_0}\right)}^3.Math.A_n\left(n,k,{}_0\right)A_n\left(n+1,k,{}_0\right).Math.\end{array}\end{array}& \\ A_n\left(l,k,{}_0\right)={}_k{}_0^{}\left[{\stackrel{\_}{P}}_{\mathrm{lk}}\left(\cos \left(+{}_0\right)\right)-{\stackrel{\_}{P}}_{\mathrm{lk}}\left(\cos \right)\right]{\stackrel{\_}{P}}_{\mathrm{nk}}\left(\cos \right)\sin {}_k=\{\begin{array}{cc}2,& k=0\\ 1,& k0\end{array}D=K\left(\frac{\cos \frac{{}_0}{2}}{2}-\frac{4}{3}\right)\frac{\sin {}_0}{r_0^4}{f}_r\left(n\right)=\underset{k=0}{\overset{n}{.Math.}}{}_0^{}{\left[{\stackrel{\_}{P}}_{\mathrm{nk}}\left(\cos \right)\right]}^2{\cos }^2\left(+\right){\sin }^2=\frac{3}{8}\left(2n+1\right)f_{\stackrel{.}{}}\left(n\right)=\left(n+\frac{1}{2}\right)r_0=a_e+h& \end{array}$ wherein .sub.n.sup.2 is an degree error variance of a geopotential coefficient of the gravity field recovery model, r.sub.0 is a geocentric range of a satellite, .sub.0 is a geocentric vector included angle, a.sub.e is an earth average radius, is a product of a gravitation constant and an earth mass, h is the gravity satellite orbit height, T.sub.arc is an integral arc length, F is non-gravitational interference, (t).sub.F is non-gravitational interference data interval, (r).sub.m is a satellite orbit determination position error, (t).sub.r is a satellite orbit data sampling interval, ({dot over ()}).sub.m is an inter-satellite range change rate measurement error, (t).sub.{dot over ()} is an inter-satellite range change rate sampling interval, T is the gravity field measurement service life, wherein ({dot over ()}).sub.m and (t).sub.{dot over ()} are obtained by the inter-satellite range measurement device, (r).sub.m and (t).sub.r are obtained by the spaceborne GPS system, F and (t).sub.F are obtained by the accelerometer; P.sub.nm(x) is the fully normalized associated Legendre polynomial, where n is the degree, m is the order, and x is a specific input ${\stackrel{\_}{P}}_{\mathrm{nm}}\left(x\right)=\sqrt{\left(2n+1\right)k\frac{\left(n-m\right)!}{\left(n+m\right)!}}{P}_{\mathrm{nm}}\left(x\right),k=\{\begin{array}{c}1,m=0\\ 2,m0\end{array}$ P.sub.nm(x) is the associated Legendre polynomial, where n is the degree, m is the order, and x is a specific input; $P_{\mathrm{nm}}\left(x\right)=\frac{{\left(-1\right)}^m}{2^{n}n!}{\left(1-x^2\right)}^{m/2}\frac{{}^{n+m}}{x^{n+m}}{\left(x^2-1\right)}^n$ C.sub.nm and S.sub.nm are cosine and sine terms of normalized gravitational coefficients in a gravity field model, respectively, C.sub.nm or S.sub.nm is a measure of difference between the observed value of gravitational coefficient and the theoretical value, n is the degree, and m is the order.

2. The method for parsing and calculating performance of satellite gravity field measurement by low-to-low satellite-to-satellite tracking, as recited in claim 1, wherein a value of a coefficient K and a phase is respectively: $K=1.3476{10}^{11}{m}^2,=-\frac{}{2}.$

3. The method for parsing and calculating performance of satellite gravity field measurement by low-to-low satellite-to-satellite tracking, as recited in claim 1, wherein regarding to the non-gravitational interference data interval (t).sub.F, if the gravity satellite system measures non-gravitational interference by low-to-low satellite-to-satellite tracking, (t).sub.F is a non-gravitational interference measurement interval; if the gravity satellite system suppresses non-gravitational interference by low-to-low satellite-to-satellite tracking, (t).sub.F is a non-gravitational interference suppression interval.

4. The method for parsing and calculating performance of satellite gravity field measurement by low-to-low satellite-to-satellite tracking, as recited in claim 1, wherein the step (3) specifically comprises a step of establishing a maximum valid order N.sub.max of the gravity field measurement, which meets a formula of: ${}_{N_{\mathrm{ma}x}}^2=\frac{1}{2N_{\mathrm{ma}x}+1}\frac{1.6{10}^{-10}}{N_{\mathrm{ma}x}^3}$ wherein N.sub.max is obtained by the gravity field model.

5. The method for parsing and calculating performance of satellite gravity field measurement by low-to-low satellite-to-satellite tracking, as recited in claim 1, wherein the step (4) specifically comprises: establishing the geoid-order error corresponding to an nth degree:
.sub.n=a.sub.e{square root over ((2n+1).sub.n.sup.2)}, and/or establishing the geoid accumulative error corresponding to the nth degree: $=\sqrt{\underset{k=2}{\overset{n}{.Math.}}{\left({}_k\right)}^2},$ and/or establishing the gravity-anomaly order error corresponding to the nth degree: $g_n=\frac{}{a_e^2}\left(n-1\right)\sqrt{\left(2n+1\right){}_n^2};$ and/or establishing the gravity-anomaly accumulative error corresponding to the nth degree: $g=\sqrt{\underset{k=2}{\overset{n}{.Math.}}{\left(g_k\right)}^2}.$

Description

BRIEF DESCRIPTION OF THE DRAWINGS

(1) FIG. 1 is a sketch view of measuring a satellite formation of gravity field by low-to-low satellite-to-satellite tracking.

(2) FIG. 2 is a sketch view showing relationships among an inter-satellite amplitude and r.sub.0 and .sub.0.

(3) FIG. 3 is a geopotential coefficient degree error variance curve measured in a gravity field under a GRACE system parameter.

(4) FIG. 4 is a GRACE gravity field measurement geoid order error and a geoid accumulative error.

(5) FIG. 5 is a GRACE gravity field measurement degree error of a gravity anomaly error and an accumulative error of gravity anomaly error.

(6) Table 1 shows an inter-satellite range change rate amplitude (m/s) under different orbit height and an inter-satellite geocentric vector included angle.

(7) Table 2 shows physical parameters of gravity field measurement by low-to-low satellite-to-satellite tracking.

(8) Table 3 shows parameters of the GRACE satellite system.

DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENT

(9) These and other objectives, features, and advantages of the present invention will become apparent from the following detailed description, the accompanying drawings, and the appended claims.

(10) One skilled in the art will understand that the embodiment of the present invention as shown in the drawings and described above is exemplary only and not intended to be limiting.

(11) As shown in FIGS. 1-5, a method for parsing and calculating performance of satellite gravity field measurement by low-to-low satellite-to-satellite tracking, comprises following steps of:

(12) step (1): acquiring at least one parameter of a gravity satellite system by low-to-low satellite-to-satellite tracking;

(13) step (2): according to the parameter of the gravity satellite system by low-to-low satellite-to-satellite tracking, calculating an effect of a measuring error of a gravity satellite load on a power spectrum of a nonspherical perturbation potential of earth gravity, so as to obtain an degree error variance of a potential coefficient in a gravity field recovery model;

(14) step (3): comparing the degree error variance of the potential coefficient in the gravity field recovery model with an degree variance of a potential coefficient given by Kaula Rule, wherein with an increase of an order of the gravity field recovery model, the degree error variance gradually increases and the degree variance gradually decreases; when the degree error variance is equal to the degree variance, considering that a maximum valid order of the gravity field measurement is obtained;

(15) step (4): according to the degree error variance of the gravity field recovery model, calculating a geoid-order error and an accumulative error of the geoid-order error, an gravity-anomaly order error and an accumulative error of the gravity-anomaly order error of the gravity field recovery model; and

(16) step (5): summarizing the valid order of the gravity field measurement, the geoid-order error and the accumulative error of the geoid-order error, and the gravity-anomaly order error and the accumulative error of the gravity-anomaly order error, so as to obtain the performance of satellite gravity field measurement by low-to-low satellite-to-satellite tracking.

(17) In the step (1), the parameter of the gravity satellite system by low-to-low satellite-to-satellite tracking is selected from but not limited to the group consisting of at least one orbit parameter of the gravity satellite system and at least one load indicator of the gravity satellite system;

(18) wherein the orbit parameter of the gravity satellite system is selected from the group consisting of a maximum valid order of the gravity field recovery N.sub.max, a gravity-anomaly order error .sub.n of an nth order, a geoid-order accumulative error of the nth order, a gravity-anomaly order error g.sub.n of the nth order and a gravity-anomaly accumulative error g of the nth order;

(19) the parameter of the gravity satellite system is selected from the group consisting of a gravity satellite orbit height h and an included angle .sub.0 of satellite-to-satellite geocentric vectors,

(20) the load indicator of the gravity satellite system is selected from the group consisting of an inter-satellite range change rate measurement error ({dot over ()}).sub.m, a satellite orbit determining position error (r).sub.m, a non-gravitational interference F, an inter-satellite range rate data sampling interval (t).sub.{dot over ()}, satellite orbit position data sampling interval (t).sub.r and a gravity field measurement service life T.

(21) The step (2) specifically comprises steps of:

(22) establishing an analytic formula of a low-to-low satellite-to-satellite tracking gravity field measurement degree error variance .sub.n.sup.2 of a potential coefficient:

(23) ${}_n^2=\frac{1}{2n+1}\underset{k=0}{\overset{n}{.Math.}}\left[{\left({\stackrel{\_}{C}}_{\mathrm{nk}}\right)}^2+{\left({\stackrel{\_}{S}}_{\mathrm{nk}}\right)}^2\right]=\frac{1}{\underset{k=0}{\overset{n}{.Math.}}\left[\begin{array}{c}{B}_1\left(r_0,n,k,{}_0\right)+\\ B_2\left(r_0,n,k,{}_0\right)+B_3\left(r_0,n,k,{}_0\right)\end{array}\right]}\frac{1}{2n+1}{\frac{2\left(n+1\right)r_0^{2n}}{{}^2{\mathrm{Ta}}_e^{2n-2}}\left[\begin{array}{c}-D\sqrt{\left(\frac{{\left(F\right)}^2{T}_{\mathrm{arc}}^4}{36}{\left(t\right)}_F+{\left(r\right)}_m^2{\left(t\right)}_r\right)f_r\left(n\right)}+\\ \sqrt{\frac{}{r_0}}\sqrt{\left(\frac{{\left(F\right)}^2{T}_{\mathrm{arc}}^2}{4}{\left(t\right)}_F+{\left(\stackrel{.}{}\right)}_m^2{\left(t\right)}_{\stackrel{.}{}}\right)f_{\stackrel{.}{}}\left(n\right)}\end{array}\right]}^2$
wherein:

(24) $\{\begin{array}{c}{B}_1\left(r_0,n,k,{}_0\right)=A_n^2\left(n-1,k,{}_0\right)+\\ \left(\frac{a_e}{r_0}\right).Math.A_n\left(n-1,k,{}_0\right)A_n\left(n,k,{}_0\right).Math.+\\ {\left(\frac{a_e}{r_0}\right)}^2.Math.A_n\left(n-1,k,{}_0\right)A_n\left(n+1,k,{}_0\right).Math.\\ B_2\left(r_0,n,k,{}_0\right)={\left(\frac{a_e}{r_0}\right)}^2{A}_n^2\left(n,k,{}_0\right)+\\ \left(\frac{a_e}{r_0}\right).Math.A_n\left(n-1,k,{}_0\right)A_n\left(n,k,{}_0\right).Math.+\\ {\left(\frac{a_e}{r_0}\right)}^3.Math.A_n\left(n,k,{}_0\right)A_n\left(n+1,k,{}_0\right).Math.\\ B_3\left(r_0,n,k,{}_0\right)={\left(\frac{a_e}{r_0}\right)}^4{A}_n^2\left(n+1,k,{}_0\right)+\\ {\left(\frac{a_e}{r_0}\right)}^2.Math.A_n\left(n-1,k,{}_0\right)A_n\left(n+1,k,{}_0\right).Math.+\\ {\left(\frac{a_e}{r_0}\right)}^3.Math.A_n\left(n,k,{}_0\right)A_n\left(n+1,k,{}_0\right).Math.\end{array}{A}_n\left(l,k,{}_0\right)={}_k\underset{0}{\overset{}{}}\left[{\stackrel{\_}{P}}_{\mathrm{lk}}\left(\cos \left(+{}_0\right)\right)-{\stackrel{\_}{P}}_{\mathrm{lk}}\left(\cos \right)\right]{\stackrel{\_}{P}}_{\mathrm{nk}}\left(\cos \right)\sin {}_k=\{\begin{array}{cc}2,& k=0\\ 1,& k0\end{array}D=K\left(\frac{\cos \frac{{}_0}{2}}{2}-\frac{4}{3}\right)\frac{\sin {}_0}{r_0^4}{f}_r\left(n\right)=\underset{k=0}{\overset{n}{.Math.}}\underset{0}{\overset{}{}}{\left[{\stackrel{\_}{P}}_{\mathrm{nk}}\left(\cos \right)\right]}^2{\cos }^2\left(+\right){\sin }^2=\frac{3}{8}\left(2n+1\right)f_{\stackrel{.}{}}\left(n\right)=\left(n+\frac{1}{2}\right)r_0=a_e+h$

(25) wherein .sub.n.sup.2 is an degree error variance of a geopotential coefficient of the gravity field recovery model, r.sub.0 is a geocentric range of a satellite, .sub.0 is a geocentric vector included angle, a.sub.e is an earth average radius, is a product of a gravitation constant and an earth mass, h is the gravity satellite orbit height, T.sub.arc is an integral arc length, F is non-gravitational interference, (t).sub.F is non-gravitational interference data interval, (r).sub.m is a satellite orbit determination position error, (t).sub.r is a satellite orbit data sampling interval, ({dot over ()}).sub.m is an inter-satellite range change rate measurement error, (t).sub. is an inter-satellite range change rate sampling interval, T is the gravity field measurement service life.

(26) Value of a coefficient K and a phase is respectively:

(27) $K=1.3476{10}^{11}{m}^2,=-\frac{}{2}.$

(28) The step (3) specifically comprises a step of establishing a maximum valid order N.sub.max of the gravity field measurement, which meets a formula of:

(29) 0${{}_N^2}_{\max }=\frac{1}{2N_{\max }+1}\frac{1.6{10}^{-10}}{N_{\max }^3}.$

(30) The step (4) specifically comprises:

(31) establishing the geoid-order error corresponding to an nth order:
.sub.n=a.sub.e{square root over ((2n+1).sub.n.sup.2)},

(32) and/or

(33) establishing the geoid accumulative error corresponding to the nth order:

(34) $=\sqrt{\underset{k=2}{\overset{n}{.Math.}}{\left({}_k\right)}^2},$

(35) and/or

(36) establishing the gravity-anomaly order error corresponding to the nth order:

(37) $g_n=\frac{}{a_e^2}\left(n-1\right)\sqrt{\left(2n+1\right){}_n^2};$

(38) and/or

(39) establishing the gravity-anomaly accumulative error corresponding to the nth order.

(40) According to a preferred embodiment, if the gravity satellite system measures non-gravitational interference by low-to-low satellite-to-satellite tracking, (t).sub.F is a non-gravitational interference measurement interval; if the gravity satellite system suppresses non-gravitational interference by low-to-low satellite-to-satellite tracking, (t).sub.F is a non-gravitational interference suppression interval.

(41) Further description of the present invention is illustrated combining with the preferred embodiment.

(42) A low-to-low satellite-to-satellite tracking gravity field satellite system comprises two low-orbit satellites to form a following formation along a track. By measuring a range between the two satellites and a change rate thereof, a gravity field of the earth is recovered. In order to make a measurement accuracy of the data as consist as possible, and facilitate controlling satellite attitude and orbit, the orbit of the satellite is generally selected as a circular orbit or a near-circular orbit. And further, in order to meet the requirement of measuring all over the world, the satellite orbit is usually selected as a polar orbit or a near-polar orbit. Thus, in the method for parsing and calculating performance of satellite gravity field measurement by low-to-low satellite-to-satellite tracking, an eccentricity of the satellite orbit is supposed to be 0, and an orbit inclination is 90.

(43) 1. Energy Conservation Equation of the Gravity Field Measurement by Low-to-Low Tracking

(44) As showing in FIG. 1, two satellites which form the formation of the gravity field measurement by low-to-low satellite-to-satellite tracking are respectively denoted as A and B. B is a guiding satellite and A is a tracking satellite. An energy conservation exists between A and B as follows:

(45) $\begin{array}{cc}{V}_A-V_B=\frac{1}{2}\left({\stackrel{.}{r}}_A^2-{\stackrel{.}{r}}_B^2\right)-\left[{\left(r_x{\stackrel{.}{r}}_y-r_y{\stackrel{.}{r}}_x\right)}_A-{\left(r_x{\stackrel{.}{r}}_y-r_y{\stackrel{.}{r}}_x\right)}_B\right]-{}_{t_0}^t\left(a_A.Math.{\stackrel{.}{r}}_A-a_B.Math.{\stackrel{.}{r}}_B\right)t-V_{t,\mathrm{AB}}-E_{0,\mathrm{AB}}& \left(1\right)\end{array}$

(46) In formula (1), a left of an equal sign is gravitational potential difference of the satellite A and the satellite B. A known gravitational potential spheric harmonics expansion is

(47) $\begin{array}{cc}V\left(r,,\right)=\frac{}{r}\left[1+\underset{n=2}{\overset{}{.Math.}}\underset{k=0}{\overset{n}{.Math.}}{\left(\frac{a_e}{r}\right)}^n\left({\stackrel{\_}{C}}_{\mathrm{nk}}\cos k+{\stackrel{\_}{S}}_{\mathrm{nk}}\sin k\right){\stackrel{\_}{P}}_{\mathrm{nk}}\left(\cos \right)\right]& \left(2\right)\end{array}$

(48) wherein (r, , ) is a spherical coordinate of a Earth-Centered Earth-Fixed coordinate system, is a gravitational constant, a.sub.e is an earth diameter, C.sub.nk and S.sub.nk are respectively a cosine value and a sine value of a nth-order k-gravitational coefficient, P.sub.nk(cos ) is a fully normalized Associated Legendre polynomials. Gravitational potential V(r, , ) is divided into two parts of a central gravitational potential /r and a nonspherical perturbation potential.

(49) $\begin{array}{cc}V\left(r,,\right)=\frac{}{r}+R\left(r,,\right)& \left(3\right)\\ R\left(r,,\right)=\frac{}{r}\underset{n=2}{\overset{}{.Math.}}\underset{k=0}{\overset{n}{.Math.}}{\left(\frac{a_e}{r}\right)}^n\left({\stackrel{\_}{C}}_{\mathrm{nk}}\cos k+{\stackrel{\_}{S}}_{\mathrm{nk}}\sin k\right){\stackrel{\_}{P}}_{\mathrm{nk}}\left(\cos \right)& \left(4\right)\end{array}$

(50) Thus, a left hand of the formula (1) is expressed as:

(51) $\begin{array}{cc}\begin{array}{c}{V}_A\left(r,,\right)-V_B\left(r,,\right)=\frac{}{r_A}+R_A\left(r,,\right)-\frac{}{r_B}-R_B\left(r,,\right)\\ =\left(\frac{}{r_A}-\frac{}{r_B}\right)+\left[R_A\left(r,,\right)-R_B\left(r,,\right)\right].\end{array}& \left(5\right)\end{array}$

(52) In the formula (1), a first item is a difference of kinetic energy of the satellite A and the satellite B. Setting e is a unit vector of the satellite A directing the satellite B,

(53) $\begin{array}{cc}e=\frac{r_{\mathrm{AB}}}{.Math.r_{\mathrm{AB}}.Math.}.& \left(6\right)\end{array}$

(54) A difference of kinetic energy of the satellite A and the satellite B is expressed as:

(55) $\begin{array}{cc}\frac{1}{2}\left({\stackrel{.}{r}}_A^2-{\stackrel{.}{r}}_B^2\right)=-\frac{1}{2}\left({\stackrel{.}{r}}_A+{\stackrel{.}{r}}_B\right).Math.{\stackrel{.}{r}}_{\mathrm{AB}}=-\frac{1}{2}\left({\stackrel{.}{r}}_A+{\stackrel{.}{r}}_B\right).Math.\left\{\left({\stackrel{.}{r}}_{\mathrm{AB}}.Math.e\right)e+\left[{\stackrel{.}{r}}_{\mathrm{AB}}-\left({\stackrel{.}{r}}_{\mathrm{AB}}.Math.e\right)e\right]\right\}.& \left(7\right)\end{array}$

(56) Considering a range between the satellite A and the satellite B is very small, an average velocity thereof is along a line segment of the satellite A and the satellite B,

(57) $\begin{array}{cc}\frac{1}{2}\left({\stackrel{.}{r}}_A+{\stackrel{.}{r}}_B\right).Math.\left[\left({\stackrel{.}{r}}_{\mathrm{AB}}.Math.e\right)e\right]\frac{1}{2}\left({\stackrel{.}{r}}_A+{\stackrel{.}{r}}_B\right).Math.\left[{\stackrel{.}{r}}_{\mathrm{AB}}-\left({\stackrel{.}{r}}_{\mathrm{AB}}.Math.e\right)e\right],& \left(8\right)\end{array}$

(58) wherein the expression (8) has an approximate expression of:

(59) 0$\begin{array}{cc}\frac{1}{2}\left({\stackrel{.}{r}}_A^2-{\stackrel{.}{r}}_B^2\right)-\frac{1}{2}\left({\stackrel{.}{r}}_A+{\stackrel{.}{r}}_B\right).Math.\left[\left({\stackrel{.}{r}}_{\mathrm{AB}}.Math.e\right)e\right]=-\frac{1}{2}\left({\stackrel{.}{r}}_A+{\stackrel{.}{r}}_B\right).Math.\left(\stackrel{.}{}e\right),& \left(9\right)\end{array}$

(60) wherein {dot over ()} is a range change rate between the satellite A and the satellite B, which is equal to a projection of {dot over (r)}.sub.AB on a unit vector e.

(61) Since orbit of the satellite A and the satellite B are assumed to be circular orbit, and geocentric distances of the satellite A and the satellite B are basically equal, an average value is set to be r.sub.0, the formula (9) is further expressed as:

(62) $\begin{array}{cc}\frac{1}{2}\left({\stackrel{.}{r}}_A^2-{\stackrel{.}{r}}_B^2\right)-\sqrt{\frac{}{r_0}}\stackrel{.}{}.& \left(10\right)\end{array}$

(63) In the formula (1), a second item of a right hand of the formula is an energy difference of the satellite A and the satellite B caused by earth rotation, wherein calculation of value indicates that the second item thereof is 4 magnitudes less than a kinetic energy difference, and thus can be ignored in parsing and analyzing. A third item of a right hand of the formula (1) is an energy difference between the satellite A and the satellite B caused by non-gravitational disturbance. In the gravity field measurement, time accumulation of the non-gravitational disturbance causes a result that the satellite orbit shifts a pure gravity orbit, and meanwhile an inter-satellite distance range change rate shifts an inter-satellite range change rate under an effect of pure gravity. An offset of both a pure gravity orbit and a measurement error determine acquiring ability of change rate of the range between the pure-gravity satellite, in such a manner that an effect of the non-gravity disturbance can be introduced in analyzing a pure gravity orbit error and a pure gravity inter-satellite range change rate. A fourth item of the right hand of the item 4 is an energy difference of the satellite A and the satellite B caused by a third body gravitation and a tidal perturbation and etc. The third body gravitation and the tidal perturbation have high-precision model, and are capable of meeting measurement requirements of a static gravity field, and an effect of the third body gravitation and the tidal perturbation can be neglected. A fifth item on a right hand of the formula (1) is an integral constant which has no effect on measurement of the gravity field measurement.

(64) When formula (5) and formula (6) are substituted into the formula (1), it can be obtained that

(65) $\begin{array}{cc}{R}_A\left(r,,\right)-R_B\left(r,,\right)=-\sqrt{\frac{}{r_0}}\stackrel{.}{}-\left(\frac{}{r_A}-\frac{}{r_B}\right)-E_{0,\mathrm{AB}}.& \left(11\right)\end{array}$

(66) Performing variation on the formula (11), it is obtained that

(67) $\begin{array}{cc}\begin{array}{c}{R}_A-R_B=\frac{1}{2}\sqrt{\frac{}{r_0^3}}\stackrel{.}{}\left(r\right)-\sqrt{\frac{}{r_0}}\stackrel{.}{}+\left(\frac{}{r_A^2}-\frac{}{r_B^2}\right)\left(r\right)\\ \frac{1}{2}\sqrt{\frac{}{r_0^3}}\stackrel{.}{}\left(r\right)-\sqrt{\frac{}{r_0}}\stackrel{.}{}+\frac{2}{r_0^3}\left(r_B-r_A\right)\left(r\right);\end{array}& \left(12\right)\end{array}$

(68) wherein r is a pure gravity orbit position error, comprising an orbit determination error and a pure gravity orbit offset caused by non-gravitational interference; {dot over ()} is a pure-gravity orbit inter-satellite range change rate, comprising an inter-satellite range measurement error and an offset of an inter-satellite range change rate caused by non-gravitational interference. A relationship is established between the pure gravity orbit position error r and the pure-gravity orbit inter-satellite range change rate {dot over ()} in the following disclosure.

(69) 2. Geometrical Relationship of the Difference of an Inter-Satellite Geocentric Range and an Inter-Satellite Range Change Rate.

(70) As known that an inter-satellite range change rate of the satellite A and the satellite B is:
{dot over ()}={dot over (r)}.sub.AB.Math.e=|{dot over (r)}.sub.AB| cos {dot over (r)}.sub.AB,e(13)

(71) Polar coordinates are adopted for calculating a geometrical relationship between the difference of the inter-satellite geocentric range r.sub.Br.sub.A and the inter-satellite range change rate {dot over ()}. As shown in FIG. 1, in a plane determined by r.sub.A and r.sub.B, L.sub.2 is set to equally divide an included angle of r.sub.A and r.sub.B and a straight line L.sub.1 is perpendicular to L.sub.2. The core of the earth is selected as a pole, and the straight line is a polar axis, an anti-clockwise direction is an angle positive position, so a position vector of the satellite A and the satellite B is:

(72) $\begin{array}{cc}{r}_A=r_A{}^{\frac{-{}_0}{2}}& \left(14\right)\\ r_B=r_B{}^{\frac{+{}_0}{2}};& \left(15\right)\end{array}$

(73) wherein .sub.0 is a geocentric vector included angle of the satellite A and the satellite B.

(74) A velocity vector of the satellite A and the satellite B is expressed as:

(75) $\begin{array}{cc}{\stackrel{.}{r}}_A={\stackrel{.}{r}}_A{}^{\frac{2-{}_0}{2}};& \left(16\right)\\ {\stackrel{.}{r}}_B={\stackrel{.}{r}}_B{}^{\frac{2+{}_0}{2}}.& \left(17\right)\end{array}$

(76) Thus, a position vector difference and a velocity vector difference of the satellite A and the satellite B is obtained:

(77) $\begin{array}{cc}{r}_{\mathrm{AB}}=r_B-r_A=r_B{}^{\frac{+{}_0}{2}}-r_A{}^{\frac{-{}_0}{2}};& \left(18\right)\\ {\stackrel{.}{r}}_{\mathrm{AB}}={\stackrel{.}{r}}_B-{\stackrel{.}{r}}_A={\stackrel{.}{r}}_B{}^{\frac{2+{}_0}{2}}-{\stackrel{.}{r}}_A{}^{\frac{2-{}_0}{2}}.& \left(19\right)\end{array}$

(78) Under the non-spherical gravitational perturbation of the earth, a geocentric distance between the satellite A and the satellite B has periodic change. Setting a geocentric distance of the satellite A at a t moment is r.sub.0, a geocentric distance of the satellite B is r.sub.0+r.sub.0, an argument of a vector e is calculated as follows. The formula (18) is substituted into the formula (6) to obtain:

(79) $\begin{array}{cc}e=\frac{r_{\mathrm{AB}}}{.Math.r_{\mathrm{AB}}.Math.}=\frac{\left(r_0+r_0\right){}^{\frac{+{}_0}{2}}-r_0{}^{\frac{+{}_0}{2}}}{},& \left(20\right)\end{array}$

(80) Wherein is an inter-satellite distance of the satellite A and the satellite B. Rearranging the formula (20) to obtain;

(81) $\begin{array}{cc}e=\frac{-\left(2r_0+r_0\right)\sin \frac{{}_0}{2}+r_0\cos \frac{{}_0}{2}}{},& \left(21\right)\end{array}$

(82) In such a manner that an argument of e is obtained:

(83) $\begin{array}{cc}\arg \left(e\right)=+\arctan \left(-\frac{r_0}{\left(2r_0+r_0\right)\tan \frac{{}_0}{2}}\right).& \left(22\right)\end{array}$

(84) An argument of {dot over (r)}.sub.AB is obtained, under a consumption of a circular orbit, a velocity value of the satellite A and the satellite B is approximately:

(85) 0$\begin{array}{cc}{\stackrel{.}{r}}_A=\sqrt{\frac{}{r_0}}& \left(23\right)\\ {\stackrel{.}{r}}_B=\sqrt{\frac{}{r_0+r_0}}.& \left(24\right)\end{array}$

(86) The formula (23) and the formula (24) are substituted into the formula (19) to obtain:

(87) $\begin{array}{cc}{\stackrel{.}{r}}_{\mathrm{AB}}={\stackrel{.}{r}}_B-{\stackrel{.}{r}}_A=\sqrt{\frac{}{r_0+r_0}}{}^{\frac{2+{}_0}{2}}-\sqrt{\frac{}{r_0}}{}^{\frac{2-{}_0}{2}}.& \left(25\right)\end{array}$

(88) Rearrange the formula (25) and obtains:

(89) $\begin{array}{cc}{\stackrel{.}{r}}_{\mathrm{AB}}=\left(\sqrt{\frac{}{r_0}}-\sqrt{\frac{}{r_0+r_0}}\right)\cos \frac{{}_0}{2}-\left(\sqrt{\frac{}{r_0+r_0}}+\sqrt{\frac{}{r_0}}\right)\sin \frac{{}_0}{2}.& \left(26\right)\end{array}$

(90) Thus, a ratio of an inscriber to a real part of {dot over (r)}.sub.AB is:

(91) $\begin{array}{cc}\frac{\mathrm{Im}\left({\stackrel{.}{r}}_{\mathrm{AB}}\right)}{\mathrm{Re}\left({\stackrel{.}{r}}_{\mathrm{AB}}\right)}=-\frac{\left(\sqrt{\frac{}{r_0+r_0}}+\sqrt{\frac{}{r_0}}\right)\sin \frac{{}_0}{2}}{\left(\sqrt{\frac{}{r_0}}-\sqrt{\frac{}{r_0+r_0}}\right)\cos \frac{{}_0}{2}}-\frac{\left(4r_0+r_0\right)\tan \frac{{}_0}{2}}{r_0},& \left(27\right)\end{array}$

(92) Thus, an argument of {dot over (r)}.sub.AB is:

(93) $\begin{array}{cc}\arg \left({\stackrel{.}{r}}_{\mathrm{AB}}\right)=\arctan \left(\frac{r_0}{\left(4r_0+r_0\right)\tan \frac{{}_0}{2}}\right)-\frac{}{2}& \left(28\right)\end{array}$

(94) From the formula (22) and (28), an included angle between e and {dot over (r)}.sub.AB is:

(95) $\begin{array}{cc}\begin{array}{c}.Math.{\stackrel{.}{r}}_{\mathrm{AB}},e.Math.=\arg \left({\stackrel{.}{r}}_{\mathrm{AB}}\right)-\arg \left(e\right)\\ =\arctan \left(\frac{r_0}{\left(4r_0+r_0\right)\tan \frac{{}_0}{2}}\right)+\\ \arctan \left(\frac{r_0}{\left(2r_0+r_0\right)\tan \frac{{}_0}{2}}\right)-\frac{3}{2};\end{array}& \left(29\right)\\ \mathrm{so}\cos .Math.{\stackrel{.}{r}}_{\mathrm{AB}},e.Math.=\cos \left[\arctan \left(\frac{r_0}{\left(4r_0+r_0\right)\tan \frac{{}_0}{2}}\right)+\arctan \left(\frac{r_0}{\left(2r_0+r_0\right)\tan \frac{{}_0}{2}}\right)-\frac{3}{2}\right].& \left(30\right)\end{array}$

(96) Simplifying the formula (30) to obtain:

(97) $\begin{array}{cc}\cos .Math.{\stackrel{.}{r}}_{\mathrm{AB}},e.Math.\frac{3r_0}{4r_0\tan \frac{{}_0}{2}}.& \left(31\right)\end{array}$

(98) From the formula (13) and the formula (31), it is obtained that:

(99) $\begin{array}{cc}.Math.{\stackrel{.}{r}}_{\mathrm{AB}}.Math.=\frac{\stackrel{.}{}}{\cos .Math.{\stackrel{.}{r}}_{\mathrm{AB}},e.Math.}=\frac{4r_0\left(\tan \frac{{}_0}{2}\right)\stackrel{.}{}}{3r_0}.& \left(32\right)\end{array}$

(100) From the formula (26), it is obtained that:

(101) $\begin{array}{cc}{.Math.{\stackrel{.}{r}}_{\mathrm{AB}}.Math.}^2={\left(\sqrt{\frac{}{r_0}}-\sqrt{\frac{}{r_0+r_0}}\right)}^2{\cos }^2\frac{{}_0}{2}+{\left(\sqrt{\frac{}{r_0+r_0}}+\sqrt{\frac{}{r_0}}\right)}^2{\sin }^2\frac{{}_0}{2}.& \left(33\right)\end{array}$

(102) From the formula (32) and the formula (33), a relationship between an inter-satellite geocentric range difference .sub.0 and an inter-satellite distance change rate {dot over ()} is:

(103) $\begin{array}{cc}{r}_B-r_A=r_0=-\frac{2r_0}{3\cos \frac{{}_0}{2}}\sqrt{\frac{r_0}{}}\stackrel{.}{}.& \left(34\right)\end{array}$

(104) When the formula (34) is substituted into the formula (12), an energy conservation of the gravity field measurement by low-to-low satellite-to-satellite is:

(105) 0$\begin{array}{cc}{R}_A-R_B=\left(\frac{1}{2}-\frac{4}{3\cos \frac{{}_0}{2}}\right)\sqrt{\frac{}{r_0^3}}\stackrel{.}{}r-\sqrt{\frac{}{r_0}}\stackrel{.}{}.& \left(35\right)\end{array}$

(106) In order to calculate a power spectral density of items in a right hand of the formula (35), a change rule of {dot over ()} with time is needed.

(107) 3. Mathematical Expression of the Inter-Satellite Distance Change Rate

(108) A semi-major of the satellite orbit changes periodically due to perturbation affection of the earth gravity, which causes a result that the inter-satellite change rate changes periodically accordingly, wherein J2 perturbation is a main perturbation item, and an effect thereof is classified as a long term item and a short term item. The long term item of the J2 perturbation has no effect on semi-major of the orbit, and an effect of the short term item on the semi-major is:

(109) $\begin{array}{cc}{a}_s\left(t\right)=\frac{3}{2}\frac{J_2}{a}\left\{\frac{2}{3}\left(1-\frac{3}{2}{\sin }^{2}i\right)\left[{\left(\frac{a}{r}\right)}^3-{\left(1-e^2\right)}^{-3/2}\right]+{\sin }^2{i\left(\frac{a}{r}\right)}^3\cos 2\left(f+\right)\right\};& \left(36\right)\end{array}$

(110) wherein a is a semi-major, i is an orbit inclination, r is an geocentric distance, e is an orbit eccentricity, f is a true anomaly, is a perigee argument. Regarding to a low-to-low satellite-to-satellite tracking gravity field satellite,

(111) $\begin{array}{cc}\{\begin{array}{c}i90\\ e0\\ ar\end{array}& \left(37\right)\end{array}$

(112) so the formula (36) is approximately expressed as:

(113) $\begin{array}{cc}{r}_s\left(t\right)=\frac{3}{2}\frac{J_2}{r}\cos \left(2+\right);& \left(38\right)\end{array}$

(114) wherein is colatitude, is an initial phase angle, so a difference of an instantaneous geocentric distance between two satellites of low-to-low satellite-to-satellite tracking is:

(115) $\begin{array}{cc}{r}_0=\frac{3}{2}\frac{J_2}{r}\left[\cos \left(2+2{}_0+\right)-\cos \left(2+\right)\right]& \left(39\right)\end{array}$

(116) So it is known that an amplitude of change of the difference of the geocentric distance of the two satellites is in direct proportion to

(117) $\frac{\sin {}_0}{r_0},$
wherein .sub.0 is geocentric vector included angle of the two satellites. Considering the formula (34), it is known that:

(118) $\begin{array}{cc}\mathrm{Am}\left(\stackrel{.}{}\right)\frac{\cos \frac{{}_0}{2}}{r_0^2}\sqrt{\frac{}{r_0}}\sin {}_0,& \left(40\right)\end{array}$

(119) wherein Am() represents an amplitude, a change period of the inter-satellite distance change rate is identical to an orbit period, so:

(120) $\begin{array}{cc}\stackrel{.}{}\left(t\right)=K\frac{\cos \frac{{}_0}{2}}{r_0^2}\sqrt{\frac{}{r_0}}\sin {}_0\cos \left(nt+\right),& \left(41\right)\end{array}$

(121) wherein K is an undetermined coefficient, n is an angular velocity of the orbit, =/2 is an initial phase. Since an orbit of the satellite is a polar orbit, so:

(122) $\begin{array}{cc}\stackrel{.}{}\left(\right)=K\frac{\cos \frac{{}_0}{2}}{r_0^2}\sqrt{\frac{}{r_0}}\sin {}_0\cos \left(+\right).& \left(42\right)\end{array}$

(123) In order to obtain a coefficient K, integral of pure gravity orbit is calculated under different orbit height and inter-satellite distance. An amplitude of vibration of the inter-satellite distance change rate is obtained according to a integral orbit, as shown in Table. 1.

(124) TABLE-US-00001 TABLE 1 0.2 0.4 0.6 0.8 1.0 200 km 8.52384717 10.sup.2 1.70054302 10.sup.1 2.54878559 10.sup.1 3.39375389 10.sup.1 4.24052520 10.sup.1 250 km 8.39648169 10.sup.2 1.66955167 10.sup.1 2.49980172 10.sup.1 3.32888860 10.sup.1 4.15999957 10.sup.1 300 km 8.22614068 10.sup.2 1.63879328 10.sup.1 2.45575028 10.sup.1 3.26890250 10.sup.1 4.08533929 10.sup.1 350 km 8.07749317 10.sup.2 1.60715522 10.sup.1 2.40566235 10.sup.1 3.20445000 10.sup.1 4.00451832 10.sup.1 400 km 7.90730919 10.sup.2 1.57561129 10.sup.1 2.35965947 10.sup.1 3.14446876 10.sup.1 3.92894479 10.sup.1 1.2 1.4 1.6 1.8 2.0 200 km 5.08769221 10.sup.1 5.93621975 10.sup.1 6.78081254 10.sup.1 7.62879389 10.sup.1 8.47098944 10.sup.1 250 km 4.99688820 10.sup.1 5.82412568 10.sup.1 6.65497074 10.sup.1 7.48424917 10.sup.1 8.31446183 10.sup.1 300 km 4.90153499 10.sup.1 5.71442505 10.sup.1 6.53303618 10.sup.1 7.34316987 10.sup.1 8.15750088 10.sup.1 350 km 4.80516090 10.sup.1 5.60716101 10.sup.1 6.40541874 10.sup.1 7.20322030 10.sup.1 8.00195000 10.sup.1 400 km 4.71103131 10.sup.1 5.49539036 10.sup.1 6.28048672 10.sup.1 7.06380954 10.sup.1 7.84346464 10.sup.1

(125) $\frac{\cos \frac{{}_0}{2}}{r_0^2}\sqrt{\frac{}{r_0}}\sin {}_0$
is a horizontal ordinate, an amplitude of the inter-satellite distance change rate in the Table 1 is an ordinate, so as to obtain a straight line shown in FIG. 2, which improves that the inter-satellite change rate formula in the expression (42) is rational. Fitting to obtain the coefficient K;
K=1.347610.sup.11 m.sup.2(43)

(126) The formula (42) is substituted into the formula (35), so as to obtain an energy change equation by low-to-low satellite-to-satellite tracking:

(127) 0$\begin{array}{cc}{R}_A-R_B=K\left(\frac{\cos \frac{{}_0}{2}}{2}-\frac{4}{3}\right)\frac{\sin {}_0}{r_0^4}\cos \left(+\right)r-\sqrt{\frac{}{r_0}}\stackrel{.}{}.& \left(44\right)\end{array}$

(128) Based on the formula (44), an degree error variance of the gravity field measurement by low-to-low satellite-to-satellite tracking is established, in such a manner that a valid order of the gravity field measurement, a geoid error and the gravity anomaly error is obtained.

(129) 4. Degree Error Variance of the Gravity Field Measurement by Low-to-Low Satellite-to-Satellite Tracking

(130) From the formula (44), an degree error variance of an nth order is:

(131) $\begin{array}{cc}{P}_n^2\left(R_A-R_B\right)=P_n^2\left[K\left(\frac{\cos \frac{{}_0}{2}}{2}-\frac{4}{3}\right)\frac{\sin {}_0}{r_0^4}\cos \left(+\right)r-\sqrt{\frac{}{r_0}}\stackrel{.}{}\right].& \left(45\right)\end{array}$

(132) Regarding to the formula (45), it is obtained from a power spectrum definition that:

(133) $\begin{array}{cc}{P}_n^2\left(R_A-R_B\right)=\underset{k=-n}{\overset{n}{.Math.}}{\left[\frac{1}{4}{}_0^2{}_0^{}\left(R_A-R_B\right){\stackrel{\_}{Y}}_{\mathrm{nk}}\left(,\right)\sin \right]}^2,& \left(46\right)\end{array}$

(134) wherein:

(135) $\begin{array}{cc}{\stackrel{\_}{Y}}_{\mathrm{nk}}\left(,\right)={\stackrel{\_}{P}}_{n.Math.k.Math.}\left(\cos \right)Q_k\left(\right)& \left(47\right)\\ Q_k\left(\right)=\{\begin{array}{cc}\cos k& k0\\ \sin .Math.k.Math.& k<0\end{array}.& \left(48\right)\end{array}$

(136) Since the two satellites of low-to-low satellite-to-satellite tracking are both operated along a polar orbit and are provide in an identical orbit plane, so a nonspherical perturbation potential of the two satellites are respectively expressed:

(137) $\begin{array}{cc}{R}_A\left(r,,\right)=\frac{}{r_0}\underset{n=2}{\overset{}{.Math.}}\underset{k=0}{\overset{n}{.Math.}}{\left(\frac{a_e}{r_0}\right)}^n\left({\stackrel{\_}{C}}_{\mathrm{nk}}\cos k+{\stackrel{\_}{S}}_{\mathrm{nk}}\sin k\right){\stackrel{\_}{P}}_{\mathrm{nk}}\left(\cos \left(+{}_0\right)\right)& \left(49\right)\\ R_B\left(r,,\right)=\frac{}{r_0}\underset{n=2}{\overset{}{.Math.}}\underset{k=0}{\overset{n}{.Math.}}{\left(\frac{a_e}{r_0}\right)}^n\left({\stackrel{\_}{C}}_{\mathrm{nk}}\cos k+{\stackrel{\_}{S}}_{\mathrm{nk}}\sin k\right){\stackrel{\_}{P}}_{\mathrm{nk}}\left(\cos \right).& \left(50\right)\end{array}$

(138) The formula (49) and the formula (50) are substituted into the formula (46) to obtain:

(139) $\begin{array}{cc}{P}_n^2\left(R_A-R_B\right)=\underset{k=-n}{\overset{n}{.Math.}}{\left[\begin{array}{c}\frac{1}{4}{}_0^2{}_0^{}\left(\begin{array}{c}\frac{}{r_0}\underset{l=2}{\overset{}{.Math.}}\underset{m=0}{\overset{l}{.Math.}}{\left(\frac{a_e}{r_0}\right)}^l\left({\stackrel{\_}{C}}_{lm}\cos m+{\stackrel{\_}{S}}_{lm}\sin m\right)\\ \left[{\stackrel{\_}{P}}_{lm}\left(\cos \left(+{}_0\right)\right)-{\stackrel{\_}{P}}_{lm}\left(\cos \right)\right]\end{array}\right)\\ {\stackrel{\_}{Y}}_{\mathrm{nk}}\left(,\right)\sin \end{array}\right]}^2\frac{1}{16}\frac{{}^2}{r_0^2}\underset{k=0}{\overset{n}{.Math.}}{\left[\begin{array}{c}\underset{l=n-1,n,n+1}{.Math.}\left(k\right)\left({\stackrel{\_}{C}}_{\mathrm{lk}}\right){\left(\frac{a_e}{r_0}\right)}^l{}_0^{}\left[\begin{array}{c}{\stackrel{\_}{P}}_{\mathrm{lk}}\left(\cos \left(+{}_0\right)\right)-\\ {\stackrel{\_}{P}}_{\mathrm{lk}}\left(\cos \right)\end{array}\right]\\ {\stackrel{\_}{P}}_{\mathrm{nk}}\left(\cos \right)\sin \end{array}\right]}^2+\frac{1}{16}\frac{{}^2}{r_0^2}\underset{k=1}{\overset{n}{.Math.}}{\left[\begin{array}{c}\stackrel{}{\underset{l=n-1,n,n+1}{.Math.}}\left({\stackrel{\_}{S}}_{\mathrm{lk}}\right){\left(\frac{a_e}{r_0}\right)}^l{}_0^{}\left[\begin{array}{c}{\stackrel{\_}{P}}_{\mathrm{lk}}\left(\cos \left(+{}_0\right)\right)-\\ {\stackrel{\_}{P}}_{\mathrm{lk}}\left(\cos \right)\end{array}\right]\\ {\stackrel{\_}{P}}_{\mathrm{nk}}\left(\cos \right)\sin \end{array}\right]}^2,& \left(51\right)\end{array}$