ФИЗИКО-МАТЕМАТИЧЕСКИЕ НАУКИ
Mathematical models for systems of high-precision
attitude determination
1 2 Akhmedov D. , Raskaliyev A.
Математические модели в системах определения высокоточной пространственной ориентации Ахмедов Д. Ш.1, Раскалиев А. С.2
1Ахмедов Даулет Шафигуллович /Akhmedov Daulet Shafigullovich - доктор технических наук,
директор;
2Раскалиев Алмат Серикович /Raskaliyev Almat Serikovich - научный сотрудник, лаборатория спутниковых навигационных технологий, Институт космической техники и технологий, г. Алматы, Республика Казахстан
Abstract: this article describes basic mathematical models of systems that have a significant impact on the solution of the problem of attitude determination. The satellite navigation technique employed is based on determining attitude angles of the object under study based on coordinate information of three GNSS antenna phase centers attached to the object. Аннотация: в данной статье описываются основные математические модели систем, которые оказывают существенное влияние на решение задачи определения пространственной ориентации. Предложенный авторами метод спутниковой навигации основан на определении углов ориентации изучаемого объекта на основе информации о координатах фазовых центров трех ГНСС антенн, закрепленных на искомом объекте.
Keywords: GNSS, attitude determination, least squares method.
Ключевые слова: ГНСС, определение пространственной ориентации, метод наименьших квадратов.
Determination of the coordinates of GNSS navigation satellites has significant impact on accuracy of the solution of this problem. This issue is well documented in the literature [1, 2].
The coordinates are calculated using a variety of satellite ephemeris that is transmitted as part of the navigational message from each of the satellites involved in the GNSS constellation. Ephemeris consists of Keplerian orbital elements and correction factors that take into account the rate of change of inclination and longitude of the ascending node. There are in total 16 parameters per satellite. They can be used to compute within margin of error approximately 1-2 meters of satellite coordinates in the geocentric coordinate system associated at any given time. In addition, each satellite navigational message contains 5 coefficients of an onboard clock model that are used to calculate absolute error of onboard clock and its drift occurring over time. Due to the relativity effect, which occurs because of non-zero eccentricity, one also requires the relativity correction, which is calculated in the middle of the process of deriving coordinates of the satellite from its ephemeris [2].
After calculation of satellite coordinates within the ECEF system, a navigational problem is solved as described in [1, 2]. For example, in order to complete preliminary calculations at least four known coordinates of satellites and pseudo ranges are necessary: coordinates of phase centers of the navigational antenna with respect to common terrestrial ellipsoid with an accuracy of about 15 m; speed with an accuracy of 0.72 km/h; and UTC time with an accuracy of 10-6 s at any time, regardless of the weather conditions.
Before using pseudo-range measurements to determine antenna phase center coordinates, one should take into account the effect of disruptive factors in order to measure the distances to a particular satellite. The measurement result is random, and the measurement error of the
carrier phase and the code is also a random variable. Furthermore, these random errors are normally distributed. In order to solve a basic navigational problem, it is necessary to resolve a number of carrier phase ambiguities. The least squares method (LSM) should be applied to resolve any of these ambiguities. LSM is based on the assumption that the «error» caused by the replacement of an accurate (unknown) value of the physical parameter on its approximate value X, calculated from observations, is proportional to the square of the error (X - a)2. Under these conditions, the optimal estimate is virtually recognized as the random value X, for which the average value of the «error» is minimal.
The relationship between the pseudo-range measurements, satellite coordinates, and unknown consumer coordinates can be written as equation (1):
P =
XC -Xr)2 +(rC -YR)2 + (ZÇ -ZR)2 +AP.
(1)
LR/ 1 y-i y-'i R,
c c c
Where X ., Y ,, Z . (satellite coordinates at the time of receiving navigation signal by the consumer);
X , Y , Z (coordinates of the consumer at the time of receiving navigation signal from the satellite);
AP. (pseudo range error caused by shift of satellite time scale from the customer);
i (sequence number of the satellite).
In order to uniquely find solutions, it is necessary to use approximate coordinates of the consumer ^p, YR ^p, ZRApp. Once we determine the coordinates, we can get the exact solution through correction computations for the approximate coordinates by using an iterative process of LSM. The iteration converges as faster and as more accurate, once the approximate consumer coordinates is identified. A satellite is considered visible if it is operational and its signal is received by the consumer and situated above line of horizon in reference to the consumer. Then we calculate and form a vector-gradient of range by defined
Yr
satellites of involved GNSS constellation as shown below:
parameters Xfl, , according to formulas (2) - (6) by using all visible navigational
dPi [Xr_ _ APP X )
dXR, P,
dP;. Y _ app — y,c )
dYri, P
dP, [ZR _ app — zC )
dzri P
dP = 1 (5)
(2)
(3)
(4)
= [
dAp
dP dP;
dPt
dXRi dYRi dZRi
1]?
(6)
Based on vector-gradient (6), we form a gradient matrix according to the formula (7):
gp=fe g2 ) (7)
Then, we form an error matrix L (column-matrix) as the difference between computed values of pseudo-ranges P. by formula (1) and measured by the receiver pseudo-ranges
P
in accordance with the formula (8):
L = P - P
(8)
Then, we implement an iterative process for determining consumer coordinates by Least Squares Method according to the formula (9):
Rk+1 = Rk + (Gp • K-; • Gpp )1 • Gp • K-; • L, (9)
where Kp (diagonal covariance matrix of the observations with priori variance measurements of pseudo-range by code). We assume that the measurements of pseudorange are equally accurate and uncorrelated;
k (sequence number of the iteration);
R (column matrix of the objective user coordinates on the k-th iteration in accordance with the formula (10)):
rk = (xr, Yr , ZR j (10)
The initial value is found according to the formula (11):
R0 = (X r _ app , yr _ app , zr _ app y (11)
The iterative process terminates under the condition (12):
Rk+i" Rk| <S, (12)
Where the value of s is selected small enough to find a solution to the problem with a good accuracy, but large enough that the iterative process converges to a certain value and does not require huge computations.
To assess the accuracy of positioning through inner convergence, one calculates the error covariance matrix of estimated parameters in accordance with formula (13):
Kr = (Gp • KP 1 • GP y1 (13)
Where the error matrix can be drawn in the formula (14):
k ar =
aa x y aa x z axaAP
a a y x ay2 aa yz ayaAP
azay az2 azaAP
aAPax aAPay aAPaz _2 aAP
(14)
An initial approximate solution for this algorithm is based on the coordinates of Earth's center (XRApp, YRApp, ZRApp ) = (0, 0, 0) and the value of s is 10 meters. A diagonal covariance matrix of the observations with priori variance code measurements is assumed as a unit, when «weighting» factors for all operational satellites are equal to each other.
Since mathematical models of differential correction were described in detail in [1-4], they are not included here. One of these models used in the high-precision attitude determination system is a mathematical model for determining attitude angles of the object. The coordinate points of A, B and C serve as input data for solving the problem of attitude determination (Figure 1). The high-precision coordinates of points A, B and C are obtained by correction variables from the local system of differential correction.
Fig. 1. Attitude determination of the object by means of coordinates of three points
Let us consider the coordinate system OXYZ. We chose a related coordinate system (CS) 0 'x'y z ' so that points A, B, and C, coordinates of which are known in mobile CS 0XYZ,
proved to align correspondingly on axes 0 'x' , 0 'y', and 0'z' (Figure 1).
These distances are known as AB = a, BC = b and AC = c.
The problem objective is to determine the following values: 0 ' A = d, 0 ' B = e and
0 ' C = f.
Triangles A0 'B, A0 C, B0'C are rectangular. Then, using Picadors' theorem, we
obtain the system of equations (15):
( /-2 , i2 2 j + d = c
d2 + e2 = a2 (15)
e2 + f2 = b2
From this system of equations, we derive another system (16):
d2 = 1 /2 * (a2 + c2 - b2) e2 = 1 /2 * (a2 + b2 - c2) f2 = 1 /2 * (b2 + c2 - a2)
(16)
Then, coordinates of points A, B and C in CS 0 ' x 'y ' z ' will be computed according to the formula (17):
A ' (d ,0,0), B ' (0, e,0), C ' (0,0, f ) (17) Taking into account coordinates of points A, B, and C in CS 0XYZ, one can determine
coordinates of the point 0 ' (an origin of related CS) in the given CS, based on the condition of constant distances d, e, and f. Therefore, one can assume the constant position of CS 0 ' x 'y 'z ' regarding to CS 0XYZ. Let's write down the system (18):
(X0 -XA )2 + O0 ~yA )2 + (z0 Z A )2 = d (x0 -xb )2 + (yo - yB )2 + (z0 -zB )2 = e 2,
(x0 -xc )2 + (y0-yc )2 + (z 0 -zc )2 = f2,
Where X0, y0, z0 represent the coordinates of point 0' in CS OXYZ; and
x^, yA, z^ represent the coordinates of point A and etc. Let us indicate the distance
between the centers of the systems 00' = D which is determined according to the formulas (19), (20):
D2 = X02 + y0 + Z02 (19)
2 2 2 2 rA = X A + yA + ZA 2 2,2,2 rB = xB + yB + ZB (20) 2 2 2 2 rc = xc + yc + zc
Let us express D through the known terms of , rB, rc, d, e, f. In order to find these known terms, we put the given formulas (18) and (19) into (20). Then, we obtain the system (21):
d 2 = ra + d 2 - 2*(x0 xa + y0 yA + z0 za )
e2 = rb + d2 - 2*(x0xb + y0yb + z0zb)
(21)
(22)
f = rc + D -2*(x0 xc + y0 yc + z0 zc ) Let us write down the differences as a system (22):
+ {yA-y
+ (yA -yc K +(zA
+ (yB-yc )y0 +(zB-zc )z0 ] where terms (22) are incorporated into the system of linear equations in relation to X0, y0, z0 . Let us write down them in the form (23):
Where F =
d2 -e2 = r 2 rA -r 2 rB -2* | (XA - XB
d2 -f = r 2 rA -r 2 'c -2* Kxa - xc
e2 -f2 = r 2 ' B -r 2 'c -2* | (xb - xc
^0 + (zA
FU= G, (23)
X A - XB yA -yB z A zB au ai2 ai3
X A - xc yA -yc z A zc = a21 a22 a23
XB - xc yA yc z A zc a3i a32 a33
U =
"v0 y0
z„
and
c =
e2 -d2 + r 2 + 'a -r 2 ' B
f -d2 2 + r 2 + 'a -r 2 ' c
f2 e2 2 + r 2 + 'b -r 2 ' c
2
c
c
2
c
3
aii ai3 x0 C1 - ai2y0
1 a33 z 0 C3 - a32y0
The obtained system of equations is singular. Therefore, the first and third equations of (23) can be presented as a system (24):
(24)
One solves this system with respect to x0, z0, and then the values x0 and y0 are determined by the formula (25):
f x0 = Jo + b
lz0 = k 2 Jo + b2
(25)
Where k1 =-; and
q
bl = -33J-— and
q
_ ai2a31 - aiia32
k 2 =- and
q
_ aiic3 - a31C1
b2 =-and
q
q = a33aii - ai3a3i.
One inserts the formula (25) into the first equation of the system (24) and obtains the formula (26):
Ay o2 - 2By o + c = 0
(26)
Where A = 1 + kf + £2, B = ki(Xa - ¿i) + k2(ZA - b2) + yA . , and
C = (¿i - Xa )2 + (¿2 - Za )2 + yA - d2.
The equation (26) has two solutions. As one has chosen to use the coordinate system 0 'x' y z ', it is necessary to pick the solution that satisfies the condition y0 > 0.
After determining the coordinates (x0, y0, z0 ) of point 0', coordinates of any point G, given in CS 0 ' x 'y z ', can be computed in system 0XYZ according to the formula (27):
g = d + m * r
(27)
Where G =
G G G.
D = 1x0,y0,Z0I , and
rb = \xg , yg , zg
M =
m.
m.
, and m
13
m
21
m
22
m
23
m
m
m
the rotation matrix, elements of which is determined
3i 32 33
according to the formula (28):
t
X fi Хл
шъх - C 0
■'•11 , ' ""21 ' 31 Г
d e J
m У а - Уо _ Уд ~ Уо ш - yc ~ y0
Ш12 - -1-> Ш22 - -> Ш32 - . , (28)
d e J
ZA Z0 ,„„ ZB Z 0 m — ^C Zo Ш13 - -1-> Ш23 - -> Ш33 - r
d e J
Another possible alternative method of attitude determination of the examined object is to use the theory of quaternion. This method of quaternion requires more processing time, but it solves potential singularity issues that can emerge during the algorithmic computations while converting angles to radians.
References
XA X0
XB X0
1. Soloviyov U. A. Systems of satellite navigation. Moscow: Eco trends, 2000. - 270 p.
2. Genike А. А., Pobedinskiy G. G. Global satellite systems of positioning and their application in geodesy. Moscow: Mapgeocenter, 2004. - 355 p.
3. Аkhmedov D. Sh., Moldabekov М. M., Raskaliyev А. S., Shabelnikov Ye. А. A mathematical model of code-phase differential correction for signals of global navigation satellite system. Almaty: Herald of the National Engineering Academy of the Republic of Kazakhstan, 2012. - Volume 2, 44 p.
4. Аkhmedov D. Sh., Moldabekov М. M., Raskaliyev А. S., Shabelnikov Ye. А. Mathematical models of differential correction in high-precision navigation systems. Almaty: Applied space research in Kazakhstan, LLP «Dyke-Press», 2010. - Volume 6, pp. 227-231.
Methods and algorithms used for calibration of INS Raskaliyev A.
Методы и алгоритмы, используемые при калибровке ИНС
Раскалиев А. С.
Раскалиев Алмат Серикович /Raskaliyev Almat Serikovich - научный сотрудник, лаборатория спутниковых навигационных технологий, Институт космической техники и технологий, г. Алматы, Республика Казахстан
Abstract: this article describes some practical and theoretical aspects of mathematical modeling of Inertial Navigation System (INS). It provides overview of methods and algorithms used for calibration in micro-electro-mechanical systems (MEMS) sensors. The author of this article also set forth implementation issues of the complementary filter and exponential noise smoothing algorithm that were encountered when calibrating MEMS sensors by six-position method (SPM) method.
Аннотация: в данной статье описываются некоторые практические и теоретические аспекты математического моделирования инерциальной навигационной системы (ИНС). Статья проводит обзор методов и алгоритмов, используемых для калибровки датчиков на основе микро-электромеханических систем (МЭМС). Автор данной статьи также обсуждает вопросы реализации комплиментарного фильтра и экспоненциального алгоритма сглаживания шумов, которые возникли при калибровке МЭМС датчиков шестипозиционным методом (SPM).