# Trajectory Estimation of Aircraft in a Double-Satellite Passive Positioning System with the Adjoint Method

Authors: Anzhou Cao; Yanqiu Gao; Jicai Zhang; Xianqing Lv
Journal: Mathematical Problems in Engineering
 (2013)

---

## Abstract

A double-satellite passive positioning system is constructed based on the theory of space geometry, where two observation coordinate systems and a fundamental coordinate system exist. In each observation coordinate system, there exists a ray from the observation satellite to the aircraft. One difficulty lies in that these two rays may not intersect due to the existence of various errors. Under this situation, this work assumes that the middle point of common perpendicular between two rays is the actual position of aircraft. Based on the theory of space geometry, the coordinates of aircraft in the fundamental coordinate system can be determined. A dynamic model with the adjoint method is developed to estimate the trajectory of aircraft during the process of rocket propulsion. By assimilating observations, the trajectory of aircraft can be calculated. Numerical experiments are designed to validate the reasonability and feasibility of this model. Simulated results indicate that even by assimilating a small number of observations, the trajectory of aircraft can be estimated. In addition, the trajectory estimation can become more accurate when more observations are assimilated to the model.


---

## Body

## 1. Introduction

The determination of the position of an aircraft by a passive observer is a typical navigation problem [1]. In this paper, a double-satellite passive positioning system is constructed based on the theory of space geometry. In this system, two observation coordinate systems corresponding to two observation satellites and a fundamental coordinate system exist. In each observation coordinate system, there exists a ray from the observation satellite to the aircraft. In the ideal situation where no error exists in the double-satellite passive positioning system, these two rays intersect and the point of intersection is the aircraft. One difficulty lies in that these two rays may not intersect due to the existence of various errors which contain systemic error of observation satellite and random error of measurement. Under this situation, this work assumes that the middle point of common perpendicular between two rays is the actual position of aircraft. Based on the assumption and theory of space geometry, the coordinates of aircraft in the fundamental coordinate system can be determined.Based on observations of aircraft in the fundamental coordinate system, the trajectory of aircraft could be estimated. Spingarn [1] used an extended Kalman filter to obtain the estimations of optimal filtered position combined with observations. Jamilnia and Naghash [2] used two different approaches, which were separated approach and integrated approach, to optimize trajectories of multistage launch vehicles simultaneously. In this paper, the adjoint method is used to estimate the trajectory of aircraft during the process of rocket propulsion. Based on the theory of inverse problem, the adjoint method is a useful tool for parameter estimation by assimilating observations into numerical models. Lardner [3], Zhang and Lu [4], Guo et al. [5], Cao et al. [6], and Chen et al. [7] studied the inversion of open boundary conditions for ocean dynamic models. As for the inverse problem of other parameters, Lu and Zhang [8] and Zhang et al. [9] inverted spatially varying bottom friction coefficient for the two-dimensional tidal model, and Yu and O’Brien [10] estimated wind drag coefficient and eddy viscosity coefficients in an Ekman layer model. Results of their work demonstrate that the adjoint method is a useful tool for parameter estimation and can improve the accuracy of simulated results.This paper is organized as follows. The double-satellite passive positioning system is constructed in Section2. The dynamic model with the adjoint method is developed to estimate the trajectory of aircraft in Section 3. In Section 4, numerical experiments are designed to validate the reasonability and feasibility of this model. Summaries and conclusions are presented in Section 5.


## 2. Double-Satellite Passive Positioning System

### 2.1. Transformation of Coordinates

Two coordinate systems, namely, the fundamental and observation coordinate system, are used in a passive positioning system. Figure1 is a snapshot of these coordinate systems. In the fundamental coordinate system (Oc-XcYcZc), the geocenter is taken as the origin. The earth and observation satellite are simplified as two particles and constitute a two-body problem under the effect of gravity. Considering that the time scale of studied problem is much shorter than that of earth rotation, the fundamental coordinate system is regarded as an inertial system. The aircraft is measured by observation satellite in the observation coordinate system (Os-XsYsZs). As shown in the observation coordinate system of Figure 1, the position of observation satellite (Os) is taken as the origin, Xs-axis is along the vector OcOs→, and Zs-axis is perpendicular to Xs-axis and pointing to the north exactly. Ys-axis is then ascertained following right hand system. Generally speaking, observation satellite would not pass the south and north poles, so the observation coordinate is definite and reasonable.Figure 1
Outline of the fundamental(Oc-XcYcZc) and observed (Os-XsYsZs) coordinate system.Considering that observation satellite is passive detective satellite which only measures direction, the aircraft is observed and represented by two dimensionless ratios (α and β) in the observation coordinate:
(1)α=ysxs,β=zsxs,
where xs, ys, and zs are the coordinates of aircraft in the observation coordinate system. It should be noted that xs, ys, and zs are just used to represent the position of aircraft and cannot be measured by one observation satellite.The key to achieve double-satellite passive position is to transform coordinates from the observation coordinate system to the fundamental one. The basic idea of transformation is described as follows. Assume thatQ is an arbitrary point in the observation coordinate system. Three components of vector OsQ→ along OsXs→, OsYs→, and OsZs→ are set to Qx→, Qy→ and Qz→, respectively. Considering that OsXs→, OsYs→, and Oszs→ are well-determined and can be described by the unit vectors i→, j→ and k→ in the fundamental coordinate system, Qx→, Qy→, and Qz→ have the following relation with i→, j→, and k→:
(2)Qx→=A1OsXs→+B1OsYs→+C1OsZs→=a1i→+b1j→+c1k→,Qy→=A2OsXs→+B2OsYs→+C2OsZs→=a2i→+b2j→+c2k→,Qz→=A3OsXs→+B3OsYs→+C3OsZs→=a3i→+b3j→+c3k→,
where A1,A2, A3, B1, B2, B3, C1, C2, C3, a1, a2, a3, b1, b2, b3, c1, c2, and c3 are parameters and will be derived in the following part. Then OsQ→ can be rewritten as
(3)OsQ→=(a1+a2+a3)i→+(b1+b2+b3)j→+(c1+c2+c3)k→,
which is the description of Q in the fundamental coordinate system. Detailed derivations are displayed as follows.Assume thatOs is (x0,y0,z0) in the fundamental coordinate system, Q is (x1,y1,z1) in the observation coordinate system, and (x2,y2,z2) in the fundamental one. Because OsXs→ is along vector OcOs→, the unit vector of OsXs→ in the fundamental coordinate system is
(4)v1→=(x0x02+y02+z02,y0x02+y02+z02,z0x02+y02+z02).
Considering that OsZs→ points to the north exactly, OsZs→ lies in the plane OsOcZc. As a result, direction vector of OsYs→ is also the normal vector of OsOcZc. Then unit vector of OsYs→ in the fundamental coordinate system can be described as
(5)v2→=(-y0x02+y02,x0x02+y02,0).
As unit vectors of OsXs→ and OsYs→ have been determined, that of Oszs→ can be calculated by vector multiplication cross
(6)v3→=(-x0z0x02z02+y02z02+(x02+y02)2,-y0z0x02z02+y02z02+(x02+y02)2,x02+y02x02z02+y02z02+(x02+y02)2).
Consequently, Q in the fundamental coordinate system is described as
(7)(x2,y2,z2)=x1v1→+y1v2→+z1v3→+(x0,y0,z0).


### 2.2. Double-Satellite Passive Positioning System

As described before, observation satellite is passive detective satellite and can measure direction only. It means that the position of aircraft can be determined only when two or more observation satellites are available. In this paper, a double-satellite passive positioning system is constructed to study the positioning problem. This means that there exist two observation coordinate systems, of which the origins are the two observation satellites, respectively. In each observation coordinate system, there exists a ray from the observation satellite to the aircraft. In the ideal situation where no error exists in the double-satellite passive positioning system, these two rays intersect and the intersection point is the aircraft. One difficulty lies in that these two rays may not intersect due to the existence of various errors. Under this situation, this work assumes that the middle point of common perpendicular between two rays is the actual position of aircraft.Asα, β, and the trajectory information of observation satellites are known, the equation of the ray from the observation satellite to the aircraft in the fundamental coordinate system can be obtained as long as (1,α,β) is transformed from the observation coordinate system to the fundamental one. According to the assumption, the position of aircraft can be determined based on the theory of space geometry. Detailed derivation is displayed as follows. (It should be noted that the unit vector rather than the equation of the ray is used in the derivation.)Figure2 is the schematic diagram of the double-satellite passive positioning system. M and N represent two satellites, and A and B represent the observations of the aircraft by satellites M and N, respectively. Based on the assumption, the middle point of AB is the actual position of the aircraft, which is represented by C in Figure 2. In the fundamental coordinate system, assume that the coordinates of M and N are (xm,ym,zm) and (xn,yn,zn), and unit vectors of MA→ and NB→ are (φ1,φ2,φ3) and (ϕ1,ϕ2,ϕ3). Consider
(8)θ1=arccos(φ1(xn-xm)+φ2(yn-ym)+φ3(zn-zm)(xn-xm)2+(yn-ym)2+(zn-zm)2),θ2=arccos(ϕ1(xn-xm)+ϕ2(yn-ym)+ϕ3(zn-zm)(xn-xm)2+(yn-ym)2+(zn-zm)2),θ=arccos(φ1ϕ1+φ2ϕ2+φ3ϕ3).
Then,
(9)|MA→|=cosθ1+cosθ2·cosθsin2θ·(xn-xm)2+(yn-ym)2+(zn-zm)2,|NB→|=cosθ2+cosθ1·cosθsin2θ·(xn-xm)2+(yn-ym)2+(zn-zm)2.
Finally, the coordinate of aircraft (C) in the fundamental coordinate system can be calculated as
(10)(xa,ya,za)=(xm,ym,zm)+|MA|→(φ1,φ2,φ3),(xb,yb,zb)=(xn,yn,zn)+|NB|→(ϕ1,ϕ2,ϕ3),(xc,yc,zc)=12[|MA→|(φ1,φ2,φ3)(xm,ym,zm)+(xn,yn,zn)+|MA→|(φ1,φ2,φ3)+|NB→|(ϕ1,ϕ2,ϕ3)].Figure 2
Outline of the double-satellite passive positioning system.


## 2.1. Transformation of Coordinates

Two coordinate systems, namely, the fundamental and observation coordinate system, are used in a passive positioning system. Figure1 is a snapshot of these coordinate systems. In the fundamental coordinate system (Oc-XcYcZc), the geocenter is taken as the origin. The earth and observation satellite are simplified as two particles and constitute a two-body problem under the effect of gravity. Considering that the time scale of studied problem is much shorter than that of earth rotation, the fundamental coordinate system is regarded as an inertial system. The aircraft is measured by observation satellite in the observation coordinate system (Os-XsYsZs). As shown in the observation coordinate system of Figure 1, the position of observation satellite (Os) is taken as the origin, Xs-axis is along the vector OcOs→, and Zs-axis is perpendicular to Xs-axis and pointing to the north exactly. Ys-axis is then ascertained following right hand system. Generally speaking, observation satellite would not pass the south and north poles, so the observation coordinate is definite and reasonable.Figure 1
Outline of the fundamental(Oc-XcYcZc) and observed (Os-XsYsZs) coordinate system.Considering that observation satellite is passive detective satellite which only measures direction, the aircraft is observed and represented by two dimensionless ratios (α and β) in the observation coordinate:
(1)α=ysxs,β=zsxs,
where xs, ys, and zs are the coordinates of aircraft in the observation coordinate system. It should be noted that xs, ys, and zs are just used to represent the position of aircraft and cannot be measured by one observation satellite.The key to achieve double-satellite passive position is to transform coordinates from the observation coordinate system to the fundamental one. The basic idea of transformation is described as follows. Assume thatQ is an arbitrary point in the observation coordinate system. Three components of vector OsQ→ along OsXs→, OsYs→, and OsZs→ are set to Qx→, Qy→ and Qz→, respectively. Considering that OsXs→, OsYs→, and Oszs→ are well-determined and can be described by the unit vectors i→, j→ and k→ in the fundamental coordinate system, Qx→, Qy→, and Qz→ have the following relation with i→, j→, and k→:
(2)Qx→=A1OsXs→+B1OsYs→+C1OsZs→=a1i→+b1j→+c1k→,Qy→=A2OsXs→+B2OsYs→+C2OsZs→=a2i→+b2j→+c2k→,Qz→=A3OsXs→+B3OsYs→+C3OsZs→=a3i→+b3j→+c3k→,
where A1,A2, A3, B1, B2, B3, C1, C2, C3, a1, a2, a3, b1, b2, b3, c1, c2, and c3 are parameters and will be derived in the following part. Then OsQ→ can be rewritten as
(3)OsQ→=(a1+a2+a3)i→+(b1+b2+b3)j→+(c1+c2+c3)k→,
which is the description of Q in the fundamental coordinate system. Detailed derivations are displayed as follows.Assume thatOs is (x0,y0,z0) in the fundamental coordinate system, Q is (x1,y1,z1) in the observation coordinate system, and (x2,y2,z2) in the fundamental one. Because OsXs→ is along vector OcOs→, the unit vector of OsXs→ in the fundamental coordinate system is
(4)v1→=(x0x02+y02+z02,y0x02+y02+z02,z0x02+y02+z02).
Considering that OsZs→ points to the north exactly, OsZs→ lies in the plane OsOcZc. As a result, direction vector of OsYs→ is also the normal vector of OsOcZc. Then unit vector of OsYs→ in the fundamental coordinate system can be described as
(5)v2→=(-y0x02+y02,x0x02+y02,0).
As unit vectors of OsXs→ and OsYs→ have been determined, that of Oszs→ can be calculated by vector multiplication cross
(6)v3→=(-x0z0x02z02+y02z02+(x02+y02)2,-y0z0x02z02+y02z02+(x02+y02)2,x02+y02x02z02+y02z02+(x02+y02)2).
Consequently, Q in the fundamental coordinate system is described as
(7)(x2,y2,z2)=x1v1→+y1v2→+z1v3→+(x0,y0,z0).


## 2.2. Double-Satellite Passive Positioning System

As described before, observation satellite is passive detective satellite and can measure direction only. It means that the position of aircraft can be determined only when two or more observation satellites are available. In this paper, a double-satellite passive positioning system is constructed to study the positioning problem. This means that there exist two observation coordinate systems, of which the origins are the two observation satellites, respectively. In each observation coordinate system, there exists a ray from the observation satellite to the aircraft. In the ideal situation where no error exists in the double-satellite passive positioning system, these two rays intersect and the intersection point is the aircraft. One difficulty lies in that these two rays may not intersect due to the existence of various errors. Under this situation, this work assumes that the middle point of common perpendicular between two rays is the actual position of aircraft.Asα, β, and the trajectory information of observation satellites are known, the equation of the ray from the observation satellite to the aircraft in the fundamental coordinate system can be obtained as long as (1,α,β) is transformed from the observation coordinate system to the fundamental one. According to the assumption, the position of aircraft can be determined based on the theory of space geometry. Detailed derivation is displayed as follows. (It should be noted that the unit vector rather than the equation of the ray is used in the derivation.)Figure2 is the schematic diagram of the double-satellite passive positioning system. M and N represent two satellites, and A and B represent the observations of the aircraft by satellites M and N, respectively. Based on the assumption, the middle point of AB is the actual position of the aircraft, which is represented by C in Figure 2. In the fundamental coordinate system, assume that the coordinates of M and N are (xm,ym,zm) and (xn,yn,zn), and unit vectors of MA→ and NB→ are (φ1,φ2,φ3) and (ϕ1,ϕ2,ϕ3). Consider
(8)θ1=arccos(φ1(xn-xm)+φ2(yn-ym)+φ3(zn-zm)(xn-xm)2+(yn-ym)2+(zn-zm)2),θ2=arccos(ϕ1(xn-xm)+ϕ2(yn-ym)+ϕ3(zn-zm)(xn-xm)2+(yn-ym)2+(zn-zm)2),θ=arccos(φ1ϕ1+φ2ϕ2+φ3ϕ3).
Then,
(9)|MA→|=cosθ1+cosθ2·cosθsin2θ·(xn-xm)2+(yn-ym)2+(zn-zm)2,|NB→|=cosθ2+cosθ1·cosθsin2θ·(xn-xm)2+(yn-ym)2+(zn-zm)2.
Finally, the coordinate of aircraft (C) in the fundamental coordinate system can be calculated as
(10)(xa,ya,za)=(xm,ym,zm)+|MA|→(φ1,φ2,φ3),(xb,yb,zb)=(xn,yn,zn)+|NB|→(ϕ1,ϕ2,ϕ3),(xc,yc,zc)=12[|MA→|(φ1,φ2,φ3)(xm,ym,zm)+(xn,yn,zn)+|MA→|(φ1,φ2,φ3)+|NB→|(ϕ1,ϕ2,ϕ3)].Figure 2
Outline of the double-satellite passive positioning system.


## 3. Dynamic model

As mentioned in Section2, the coordinates of aircraft in the fundamental coordinate system are determined by the double-satellite passive positioning system, and they will be used as observations in the dynamic model. If adequate observations are obtained, the trajectory of the aircraft can be estimated by using the dynamic model with the adjoint method.


### 3.1. Forward Model (Governing Equation)

The governing equation of the aircraft during the process of rocket propulsion in the fundamental coordinate system can be described as follows:(11)r→¨(t)=F→e+F→T,
where r→(t) is the vector of position, t is time, r→¨(t) is the double derivative in time of r→(t), that is, acceleration, Gm is the constant of earth attraction, F→e=-(Gm/|r→(t)|3)r→(t) is the acceleration caused by universal gravitation, and F→T is the acceleration caused by rocket propulsion, which is taken as the unknown parameter and will be optimized by the adjoint method. Assume that Fx, Fy and Fz, are the x-, y-, and z-components of F→T, respectively, (11) can be changed into
(12)d2xdt2=-Gmr3x+Fx,d2ydt2=-Gmr3y+Fy,d2zdt2=-Gmr3z+Fz,r=x2+y2+z2,
where x, y, and z are three components of r→(t).


### 3.2. Adjoint Model and Corrections

The cost function is defined as(13)J=κ2∫t1t2[(x-x^)2+(y-y^)2+(z-z^)2]dt,
where κ is the unit matrix [8] x^, y^, z^, x, x, and z are observations and simulated coordinate of the aircraft, respectively. Then the Lagrangian function can be expressed as
(14)L=J+∫t1t2[λ(d2xdt2+Gmr3x-Fx)+μ(d2ydt2+Gmr3y-Fy)+ν(d2zdt2+Gmr3z-Fz)]dt,
where λ,  μ, and ν are the adjoint variables (namely, Lagrangian multipliers) of x, y, and z, respectively.According to the theory of Lagrangian multiplier method, we obtain the first-order derivatives of Lagrangian function with respect to all the variables,(15)∂L∂λ=0,∂L∂μ=0,∂L∂ν=0,(16)∂L∂x=0,∂L∂y=0,∂L∂z=0,(17)∂L∂Fx=0,∂L∂Fy=0,∂L∂Fz=0.
Equation (15) give the forward model described by (12). From (16) the adjoint model can be deduced as
(18)d2λdt2+Gm(y2+z2-2x2)λ-3x(yμ+zν)r5+κ(x-x^)=0,d2μdt2+Gm(x2+z2-2y2)μ-3y(xλ+zν)r5+κ(y-y^)=0,d2νdt2+Gm(x2+y2-2z2)ν-3z(xλ+yμ)r5+κ(z-z^)=0.
From (17), the corrections of Fx, Fy, and Fz are obtained
(19)Fx-F~x=λ,Fy-F~y=μ,Fz-F~z=ν,
where F~x,  F~y,  F~z and Fx, Fy, Fz are prior and optimized values, receptively. Because the values of cost function decrease along the opposite direction of the gradient, by employing typical steepest descent method [11], Fx, Fy, and Fz can be optimized during iterations.Equations (12), (18) and (19) are discretized using the finite difference method and are shown as follows:
(20)xn+1-2xn+xn-1Δt2=-Gm(rn)3xn+Fxn,yn+1-2yn+yn-1Δt2=-Gm(rn)3yn+Fyn,zn+1-2zn+zn-1Δt2=-Gm(rn)3zn+Fznrn=(xn)2+(yn)2+(zn)2,λn-1-2λn+λn+1Δt2+Gm((yn)2+(zn)2-2(xn)2)λn-3xn(ynμn+znνn)(rn)5+κ(xn-x^n)=0,μn-1-2μn+μn+1Δt2+Gm((xn)2+(zn)2-2(yn)2)μn-3yn(xnλn+znνn)(rn)5+κ(yn-y^n)=0,νn-1-2νn+νn+1Δt2+Gm((xn)2+(yn)2-2(zn)2)νn-3zn(xnλn+ynμn)(rn)5+κ(zn-z^n)=0,Fxn-F~xn=λn,Fyn-F~yn=μn,Fzn-F~zn=νn.The computed process is designed as follows. (a) Run the forward model with prescribed initial values ofFx, Fy, and Fz. (b) Calculate the value of cost function J. (c) The difference between simulated results and observations plays as the external driving force of the adjoint model. The adjoint variables can then be obtained through the backward integration of the adjoint model. (d) Based on the corrections and typical steepest descent method, Fx, Fy, and Fz can be optimized. Repeating the process aforementioned, the difference between simulated results and observations is decreased. When the convergence criterion is met, the optimization stops.The convergence criterion could be that the number of iteration steps is prescribed, the last two values of the cost function are sufficiently close, the magnitude of the gradient norm is sufficiently small, the discrepancy between the updated and old parameters is sufficiently small, or a combination of these. In this work, the convergence criterion is that the number of iteration steps is equal to 10000, which leads to satisfied values of cost function and gradient norm and is obtained through a trial process. In addition, this convergence criterion is reasonable and essential to investigate the influence of the number of observations on the accuracy of simulated results in the numerical experiments.


## 3.1. Forward Model (Governing Equation)

The governing equation of the aircraft during the process of rocket propulsion in the fundamental coordinate system can be described as follows:(11)r→¨(t)=F→e+F→T,
where r→(t) is the vector of position, t is time, r→¨(t) is the double derivative in time of r→(t), that is, acceleration, Gm is the constant of earth attraction, F→e=-(Gm/|r→(t)|3)r→(t) is the acceleration caused by universal gravitation, and F→T is the acceleration caused by rocket propulsion, which is taken as the unknown parameter and will be optimized by the adjoint method. Assume that Fx, Fy and Fz, are the x-, y-, and z-components of F→T, respectively, (11) can be changed into
(12)d2xdt2=-Gmr3x+Fx,d2ydt2=-Gmr3y+Fy,d2zdt2=-Gmr3z+Fz,r=x2+y2+z2,
where x, y, and z are three components of r→(t).


## 3.2. Adjoint Model and Corrections

The cost function is defined as(13)J=κ2∫t1t2[(x-x^)2+(y-y^)2+(z-z^)2]dt,
where κ is the unit matrix [8] x^, y^, z^, x, x, and z are observations and simulated coordinate of the aircraft, respectively. Then the Lagrangian function can be expressed as
(14)L=J+∫t1t2[λ(d2xdt2+Gmr3x-Fx)+μ(d2ydt2+Gmr3y-Fy)+ν(d2zdt2+Gmr3z-Fz)]dt,
where λ,  μ, and ν are the adjoint variables (namely, Lagrangian multipliers) of x, y, and z, respectively.According to the theory of Lagrangian multiplier method, we obtain the first-order derivatives of Lagrangian function with respect to all the variables,(15)∂L∂λ=0,∂L∂μ=0,∂L∂ν=0,(16)∂L∂x=0,∂L∂y=0,∂L∂z=0,(17)∂L∂Fx=0,∂L∂Fy=0,∂L∂Fz=0.
Equation (15) give the forward model described by (12). From (16) the adjoint model can be deduced as
(18)d2λdt2+Gm(y2+z2-2x2)λ-3x(yμ+zν)r5+κ(x-x^)=0,d2μdt2+Gm(x2+z2-2y2)μ-3y(xλ+zν)r5+κ(y-y^)=0,d2νdt2+Gm(x2+y2-2z2)ν-3z(xλ+yμ)r5+κ(z-z^)=0.
From (17), the corrections of Fx, Fy, and Fz are obtained
(19)Fx-F~x=λ,Fy-F~y=μ,Fz-F~z=ν,
where F~x,  F~y,  F~z and Fx, Fy, Fz are prior and optimized values, receptively. Because the values of cost function decrease along the opposite direction of the gradient, by employing typical steepest descent method [11], Fx, Fy, and Fz can be optimized during iterations.Equations (12), (18) and (19) are discretized using the finite difference method and are shown as follows:
(20)xn+1-2xn+xn-1Δt2=-Gm(rn)3xn+Fxn,yn+1-2yn+yn-1Δt2=-Gm(rn)3yn+Fyn,zn+1-2zn+zn-1Δt2=-Gm(rn)3zn+Fznrn=(xn)2+(yn)2+(zn)2,λn-1-2λn+λn+1Δt2+Gm((yn)2+(zn)2-2(xn)2)λn-3xn(ynμn+znνn)(rn)5+κ(xn-x^n)=0,μn-1-2μn+μn+1Δt2+Gm((xn)2+(zn)2-2(yn)2)μn-3yn(xnλn+znνn)(rn)5+κ(yn-y^n)=0,νn-1-2νn+νn+1Δt2+Gm((xn)2+(yn)2-2(zn)2)νn-3zn(xnλn+ynμn)(rn)5+κ(zn-z^n)=0,Fxn-F~xn=λn,Fyn-F~yn=μn,Fzn-F~zn=νn.The computed process is designed as follows. (a) Run the forward model with prescribed initial values ofFx, Fy, and Fz. (b) Calculate the value of cost function J. (c) The difference between simulated results and observations plays as the external driving force of the adjoint model. The adjoint variables can then be obtained through the backward integration of the adjoint model. (d) Based on the corrections and typical steepest descent method, Fx, Fy, and Fz can be optimized. Repeating the process aforementioned, the difference between simulated results and observations is decreased. When the convergence criterion is met, the optimization stops.The convergence criterion could be that the number of iteration steps is prescribed, the last two values of the cost function are sufficiently close, the magnitude of the gradient norm is sufficiently small, the discrepancy between the updated and old parameters is sufficiently small, or a combination of these. In this work, the convergence criterion is that the number of iteration steps is equal to 10000, which leads to satisfied values of cost function and gradient norm and is obtained through a trial process. In addition, this convergence criterion is reasonable and essential to investigate the influence of the number of observations on the accuracy of simulated results in the numerical experiments.


## 4. Numerical Experiments

The data used in the numerical experiments, containingα, β, and trajectory information of two observation satellites (satellites 06 and 09), are from problem B of the ninth National Graduate Mathematical Contest in Modeling (http://www.njnet.edu.cn/news/16/). These data are from 50 to 170 s with an interval of 0.2 s, corresponding to the process of rocket propulsion. The trajectories of observation satellites are plotted in Figure 3. α and β are plotted in Figure 4. Then the coordinates of aircraft in the fundamental coordinate system, which are used as observations, are determined according to formulas in Section 2 and are shown in Figure 5.Coordinates of observation satellites in the fundamental coordinate system.
(a)
(b)
(c)α and β of observation satellites from 50 s to 170 s.
(a)
(b)Figure 5
Trajectory (left) and instantaneous coordinates (right) of the aircraft.Four numerical experiments are designed to validate the sensitivity to the number of observations in the model. In experiment one (E1), the trajectory of aircraft is estimated by assimilating all observations. The number of observations is 601 in E1 with internal of 0.2 s. In E2–E4, the number of observations is decreased successively, which is 61, 21, and 13, respectively. In numerical experiments, the initial values ofFx, Fy, and Fz are all equal to 10 m2/s, and the time step is equal to 0.2 s. Simulated results are illustrated as follows.Simulated results show good agreement with the observed ones in all four experiments. However, they are not plotted in a figure for very small differences with each other and the observations. Table1 lists the mean absolute differences (MADs) and mean relative differences (MRDs) between simulated results and observations. In all four experiments, the MADs and MRDs are all very small, indicating the reasonability and feasibility of the model. Log of cost function and gradient norm are plotted in Figure 6. Cost function and gradient norm decrease fast during the iteration, implying the efficiency of the typical steepest method. Careful inspection shows that the MAD and MRD, as well as the value of cost function at the end of iterations, increase from E1 to E4. This is to say, the trajectory estimation of the aircraft becomes worse as the number of observations is decreased. Combining all these results, we can draw two conclusions. One is that even by assimilating a small number of observations, the trajectory of aircraft can be estimated using the dynamic model with the adjoint method. The other is that the trajectory estimation can become more accurate by assimilating more observations.Table 1
Difference between simulated coordinates and observed ones.
Mean absolute difference (m)
Mean relative difference (%)
x
y
z
x
y
z
TE1
116
38
71
0.0097
0.00063
0.0057
TE2
124
40
77
0.0010
0.00065
0.0062
TE3
160
40
103
0.0013
0.00065
0.0083
TE4
244
43
162
0.0020
0.00069
0.0130Log of cost function and gradient norm during iteration.
(a)
(b)


## 5. Summary and Conclusions

In this paper, a double-satellite passive positioning system is constructed based on the theory of space geometry. In this system, two observation coordinate systems corresponding to two observation satellites and a fundamental coordinate system exist. In each observation coordinate system, there exists a ray from the observation satellite to the aircraft. In the ideal situation where no error exists in the double-satellite passive positioning system, these two rays intersect and the point of intersection is the aircraft. One difficulty lies in that these two rays may not intersect due to the existence of various errors. Under this situation, this work assumes that the middle point of common perpendicular between two rays is the actual position of aircraft. Based on the assumption and theory of space geometry, the coordinates of aircraft in the fundamental coordinate system can be determined.In order to estimate the trajectory of aircraft during the process of rocket propulsion, a dynamic model with the adjoint method is constructed. The forward model is used to simulate the motion of aircraft, and the adjoint model is used to optimize parameters. By assimilating observations, which are the coordinates of aircraft in the fundamental coordinate system, the trajectory of aircraft can be estimated. Numerical experiments are designed to validate the reasonability and feasibility of this model. Results indicate that even by assimilating a small number of observations, the trajectory of aircraft can be estimated using the dynamic model with the adjoint method. In addition, the trajectory estimation can become more accurate by assimilating more observations.



---

*Source: 502610-2013-06-15.xml*