CN102998690B - Attitude angle direct resolving method based on global position system (GPS) carrier wave double-difference equation - Google Patents

Abstract

The invention provides an attitude angle direct resolving method based on a global position system (GPS) carrier wave double-difference equation. The attitude angle direct resolving method mainly comprises building a carrier wave double-difference equation, guiding carrier attitude information into the carrier wave double-difference equation and solving an attitude angle by utilizing nonlinear least squares. The attitude angle direct resolving method can effectively reduce or eliminate short base line public errors, achieves high-precision resolving, can directly resolve the attitude angle by leading the attitude angle information into the double-difference equation, greatly reduces intermediate evaluated errors, evaluates the attitude angle by utilizing the nonlinear least squares method, improves system resolving speed, finishes real-time resolving of the attitude information and is more suitable for real-time carrier attitude measuring.

Description

A Direct Solution Method of Attitude Angle Based on GPS Carrier Double Difference Equation

technical field

The invention relates to a real-time attitude measurement method of a large carrier, in particular to a method for directly solving the attitude angle based on the GPS carrier double-difference equation.

Background technique

In order to meet the attitude information requirements of different application platforms, a large number of attitude measurement devices have come out, such as star sensors, horizon trackers, sun sensors, and magnetometers for space carrier attitude measurement; Or the magnetic compass and electromagnetic compass that provide heading for underwater vehicles; the inertial devices that provide attitude information for various aircraft, land transportation equipment, ships, submersibles, and space carriers, etc. Among them, the star sensor is easily affected by the weather, terrain or other objective factors, and cannot give the attitude or orientation value in real time; the inertial navigation equipment is complex in structure and expensive, and the error of the inertial device will increase with time when the working time is long. The accumulation causes the decrease of attitude measurement accuracy.

The use of GPS positioning signals for carrier attitude measurement only needs to use low-cost receivers to provide high-precision attitude information, thereby replacing expensive traditional attitude measurement equipment, and can also complete carrier positioning and timing. And it is less affected by the environment, and can perform high-precision attitude measurement tasks for a long time.

GPS carrier attitude measurement is to determine the attitude angle by using the relative position between the receiving antennas installed on the fixed position of the carrier and the coordinate transformation relationship. The specific implementation method is to use the high-precision positioning system installed on the carrier as a reference station, and install three (or more) receiving antennas at a fixed position on the carrier, and use GPS signals to measure the relative position between the antenna and the base station to obtain the baseline in WGS- 84 (World Geodetic System-84) coordinates, the attitude rotation matrix is determined according to the known coordinates of the baseline carrier system, and finally the attitude angle is obtained.

At present, the GPS attitude calculation algorithm is mainly based on the processing of the positioning results to obtain the attitude angle, and the commonly used mathematical estimation method is applied to the attitude calculation. According to the calculation results, it is divided into two types: the calculation of the attitude matrix and the calculation of the Euler angle.

The method of solving the attitude matrix to determine the attitude angle is to use the positioning equation to establish the positioning equations of multiple satellites and multiple epochs, to complete the estimation of the baseline coordinates by solving the overdetermined equation, and to use the coordinate conversion relationship of the baseline in different coordinate systems to obtain The attitude matrix is estimated to determine the attitude angle. This type of method needs to use mathematical estimation means to realize the solution of baseline coordinates and the determination of attitude angle. Commonly used methods include:

(1) According to the coordinate relationship between the baseline in the geographic system and the carrier system, the TRIAD algorithm is used for optimal estimation using the orthogonality of the attitude matrix;

(2) Solving methods based on Wahba problems: such as QUEST (QUaternion ESTimator) algorithm, FOAM (FastOptimal Attitude Matrix) algorithm, Euler-q algorithm, SVD (Singular Value Decomposition) algorithm, etc.;

(3) The least squares direct estimation method using multi-epoch baseline coordinate solution;

(4) Using the pseudorange or carrier observations and the baseline coordinates of the carrier system to establish an observation equation to estimate the attitude matrix.

The direct estimation method of the attitude angle generally uses the known relationship between the antenna position and the rotation attitude angle to solve the attitude angle step by step. Including: set a baseline in the direction of the main axis of the carrier, first obtain the yaw angle and pitch angle, and then use the rotation relationship of the second baseline to obtain the two-antenna attitude measurement method; use the two-antenna attitude measurement formula to determine Yaw angle and pitch angle, the multi-antenna attitude measurement method of roll angle obtained by rotating another antenna twice.

Generally speaking, there are the following problems in GPS attitude measurement:

1. When using the mathematical estimation method to estimate the attitude matrix, the attitude angle result is not directly obtained, and multi-epoch information collection is required to complete the final attitude angle calculation. This introduces further estimation errors on the one hand, and on the other hand It also affects the real-time performance of attitude measurement.

2. Using the direct measurement method to measure the attitude angle does not require the coordinates of the baseline carrier system, nor does it need to calculate the attitude matrix, but the z-axis component of the carrier coordinates of all antenna baselines is zero, and the matrix formed by the coordinates is not of rank, which will lead to attitude results Unreliable and far less accurate than mathematical estimation methods.

Contents of the invention

The purpose of the present invention is to provide a method for directly solving the attitude angle based on the GPS carrier double-difference equation that can improve the accuracy of the estimation process and realize the requirement of real-time attitude measurement.

The purpose of the present invention is achieved like this:

A method for directly solving the attitude angle based on the GPS carrier double-difference equation comprises the following steps:

以直角处0号天线的载波观测信息作为基站信号；(1) Use three GPS receiver measurement antennas to collect GPS satellite signals to obtain carrier observation information

The carrier observation information of antenna 0 at the right angle is used as the base station signal;

n=0,1,2，然后对应不同天线，将各卫星单差载波信息向参考卫星的单差载波观测量作差得到双差载波观测量 (2) Select the public visible satellite i, i=1,2,...,M, and take the satellite with the largest elevation angle among the visible satellites as the reference satellite (number it 1), and the receiver obtains the three-antenna carrier observation information Finally, based on the observation information of the base station antenna, the difference between the carrier signal of each antenna relative to the same visible satellite to the base station signal is obtained to obtain the single-difference carrier phase observation

n=0, 1, 2, and corresponding to different antennas, the difference between the single-difference carrier information of each satellite and the single-difference carrier observation of the reference satellite is obtained to obtain the double-difference carrier observation

(3) Using the relationship between carrier observations and attitude angle shown in formula (1), establish the attitude angle solution equation:

表示n号天线的双差载波观测量；λ表示载波波长，单位：m；

表示双差整周模糊度；

表示卫星i到参考的视线向量；R表示地理系到载体系的姿态矩阵；

表示地理系到地心系的姿态转换矩阵；y，r，p表示载体的三个姿态角，分别为航向角、横滚角、俯仰角，即方程未知数；b_n表示载体坐标系下基线向量，为已知量，单位：m。in,

Indicates the double-difference carrier observation of antenna number n; λ indicates the carrier wavelength, unit: m;

Indicates the double-differenced integer ambiguity;

Represents the line-of-sight vector from satellite i to the reference; R represents the attitude matrix from the geographic system to the carrier system;

Represents the attitude transformation matrix from the geographic system to the earth-centered system; y, r, p represent the three attitude angles of the carrier, which are the heading angle, roll angle, and pitch angle, which are the unknowns of the equation; b _n represents the baseline vector in the carrier coordinate system , is a known quantity, unit: m.

(4) The attitude measurement system acquires the single-epoch carrier signal in real time for double-difference processing, and calculates the attitude information of the carrier in real time according to steps (1)-(4).

The present invention may also include:

The method for directly solving the attitude angle based on the GPS carrier double-difference equation is characterized in that: the antenna adopts a right-angled orthogonal layout, that is, a right-angled triangle plane is formed, and the installation distance is more than 10m.

The direct solution method of the attitude angle based on the GPS carrier double-difference equation is characterized in that: the calculation process of the attitude angle solution equation is:

1) Since the ambiguity on the right side of the equation (1) is known, the equation (1) is expressed as:

y _t =f(x)

表示整周模糊度求解后历元时刻t的双差载波观测量； $f\left(x\right)={\mathrm{\λ}}^{-1}\left[\begin{array}{c}-{(I_0^2-I_0^1)}^T\\ -{(I_0^3-I_0^1)}^T\\ ...\\ -{(I_0^M-I_0^1)}^T\end{array}\right]A_R^E{R}^T\left(x{b}_n\right)$ 表示以x为自变量的函数，x＝[y r p]^T。in,

Indicates the double-difference carrier observations at epoch time t after the whole-week ambiguity is solved; $f\left(x\right)={\mathrm{\&lambda;}}^{-1}\left[\begin{array}{c}-{({\mathrm{<figure-callout\; id="0"\; label="I"\; filenames="HDA00002467225800022.png,HDA00002467225800032.png"\; state="\{\{state\}\}">I</figure-callout>}}_0^2-I_0^1)}^T\\ -{({\mathrm{<figure-callout\; id="0"\; label="I"\; filenames="HDA00002467225800022.png,HDA00002467225800032.png"\; state="\{\{state\}\}">I</figure-callout>}}_0^3-I_0^1)}^T\\ ...\\ -{(I_0^m-I_0^1)}^T\end{array}\right]A_R^{\mathrm{E.}}{R}^T\left(x{b}_{\mathrm{no}}\right)$ Represents a function with x as an independent variable, x=[y r p] ^T .

在初始姿态角处对方程f(x)进行一阶泰勒级数展开，令 $\mathrm{\&Delta;x}={\left[\begin{array}{ccc}\mathrm{\&Delta;y}& \mathrm{\&Delta;r}& \mathrm{\&Delta;p}\end{array}\right]}^T=\begin{array}{}\end{array}{\left[\begin{array}{ccc}y-\hat{y}& r-\hat{r}& p-\hat{p}\end{array}\right]}^T$ 为状态变量，方程简化为如下形式：2) Give the initial attitude angle

Perform a first-order Taylor series expansion on the equation f(x) at the initial attitude angle, so that $\mathrm{\&Delta;x}={\left[\begin{array}{ccc}\mathrm{\&Delta;y}& \mathrm{\&Delta;r}& \mathrm{\&Delta;p}\end{array}\right]}^T=\begin{array}{}\end{array}{\left[\begin{array}{ccc}\mathrm{the\; y}-\widehat{\mathrm{the\; y}}& r-\hat{r}& p-\hat{p}\end{array}\right]}^T$ is the state variable, the equation simplifies to the following form:

H·Δx= _Δyt (2)

的偏差；观测矩阵H为：In the formula, That is, the estimated result x _k of the state vector in the kth iteration is relative to the measured value

The deviation; the observation matrix H is:

$\mathrm{<span\; class="notranslate">\; h</span>}<span\; class="notranslate">\; =</span>\left[\begin{array}{ccc}& <span\; class="notranslate">\; .</span>& \\ & <span\; class="notranslate">\; .</span>& \\ & <span\; class="notranslate">\; .</span>& \\ <span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span>& <span\; class="notranslate">\; -</span>{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; I</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}^{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; <figure-callout\; id="0"\; label="I"\; filenames="HDA00002467225800022.png,HDA00002467225800032.png"\; state="\{\{state\}\}">I</figure-callout></span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}^{\mathrm{<span\; class="notranslate">\; 1</span>}}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; T</span>}}{\mathrm{<span\; class="notranslate">\; A</span>}}_{\mathrm{<span\; class="notranslate">\; R</span>}}^{\mathrm{<span\; class="notranslate">\; E.</span>}}<span\; class="notranslate">\; (</span>\frac{{<span\; class="notranslate">\; \&PartialD;</span>\mathrm{<span\; class="notranslate">\; R</span>}}^{\mathrm{<span\; class="notranslate">\; T</span>}}}{<span\; class="notranslate">\; \&PartialD;</span>\mathrm{<span\; class="notranslate">\; x</span>}}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; b</span>}}_{\mathrm{<span\; class="notranslate">\; no</span>}}& <span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span>\\ & <span\; class="notranslate">\; .</span>& \\ & <span\; class="notranslate">\; .</span>& \\ & <span\; class="notranslate">\; .</span>& \end{array}\right]$

表示姿态矩阵对各姿态角的偏导数。in,

Indicates the partial derivative of the attitude matrix with respect to each attitude angle.

3) Using least square estimation to solve overdetermined equation (2), the solution is:

Δx＝(H ^T H) ^-1 H ^T Δy _t (3)

4) Set the iteration cut-off threshold and estimate the final value x of the attitude angle:

$\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; =</span>\begin{array}{c}\left[\begin{array}{c}\mathrm{<span\; class="notranslate">\; the\; y</span>}\\ \mathrm{<span\; class="notranslate">\; r</span>}\\ \mathrm{<span\; class="notranslate">\; p</span>}\end{array}\right]\end{array}<span\; class="notranslate">\; =</span>\left[\begin{array}{c}\stackrel{<span\; class="notranslate">\; ^</span>}{\mathrm{<span\; class="notranslate">\; the\; y</span>}}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; \&Delta;y</span>}\\ \stackrel{<span\; class="notranslate">\; ^</span>}{\mathrm{<span\; class="notranslate">\; r</span>}}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; \&Delta;r</span>}\\ \stackrel{<span\; class="notranslate">\; ^</span>}{\mathrm{<span\; class="notranslate">\; p</span>}}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; \&Delta;p</span>}\end{array}\right]<span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 4</span>}<span\; class="notranslate">\; )</span>$

The main features of the method of the present invention are as follows:

(1) GPS carrier observations can provide high-precision positioning information, and the use of double-difference to establish carrier equations can effectively eliminate the common errors at both ends of the attitude measurement baseline and greatly improve the attitude measurement accuracy;

(2) The attitude angle is introduced into the double-difference measurement equation as an unknown quantity, which reduces the intermediate process of the traditional attitude measurement method, can effectively reduce the estimation error, and can directly obtain the attitude angle information through the single epoch information solution;

(3) The introduction of the nonlinear least squares estimation algorithm can realize the estimation of the nonlinear overdetermined equation, realize the fast solution of a single epoch, and greatly improve the real-time estimation of the attitude information of the carrier.

Beneficial effect of the present invention can be verified by following simulation:

1. Attitude Angle Calculation Simulation Verification Model Establishment

The design of the attitude simulation verification model uses the known satellite position, and the attitude algorithm can be verified by setting the geodetic coordinates of the main antenna, the coordinates of the baseline carrier system, and the preset attitude angle.

The design of the verification system includes the following parts:

(1) Satellite coordinate acquisition

Using the almanac files actually broadcast by GPS satellites provided by the navcen.uscg.gov website, calculate the real-time coordinates of all satellites in orbit, and obtain the satellite trajectory with coarse precision. Complete the calculation of elevation angles of all satellites according to the set user position, and obtain the position information of visible satellites.

(2) Carrier signal simulation

After setting the positions of the main antenna and the reference antenna, the geocentric coordinates of all antennas can be obtained according to the predetermined attitude angle, that is, the distance from the satellite to the antenna is known. By adding a certain measurement noise, the carrier signal simulation of each antenna can be realized, and the establishment of the carrier observation equation can be completed (the system assumes that the ambiguity of the whole cycle is known).

(3) Carrier posture result verification

Using the simulation parameters as the known information, the equations are established according to different attitude calculation methods, the attitude angle information is solved, and then compared with the preset attitude angle, so as to realize the verification of the algorithm.

The attitude measurement array of three antennas is established in the carrier coordinate system, and the distribution is shown in Figure 3.

Assuming that the reference antenna is at the carrier rotation center, the satellite cut-off elevation angle is 5°, and other simulation conditions are set as shown in Table 1.

Table 1 Simulation parameter settings

2. Real-time verification of attitude angle calculation method

Set the initial GPS time as 0, sample 100 epochs, and observe the real-time calculation effect of the algorithm in three states: the carrier is stationary, rotates at a fixed angular rate, and rotates sinusoidally with time. One of the iterations of a single epoch is shown in Figure 4, and it is found that the attitude results are stable after 4 iterations. The verification results of 100 solution epochs are shown in Figure 5, and the solution error situation is shown in Figure 6.

Analyzing the attitude result curve in Figure 6, under static conditions, the accuracy of the direct attitude angle calculation method based on the carrier double-difference equation can reach the order of 10 ^-4 degrees, and can be effectively solved within the simulation time; when the carrier is fixed When the angular rate rotation heading angle changes with time, the dynamic heading angle accuracy also reaches the order of 10 ^-4 degrees; when the carrier pitch angle changes with time according to the sinusoidal law, the solution accuracy reaches the order of 10 ^-2 degrees.

Description of drawings

Fig. 1 is the flow chart of method solution of the present invention;

Fig. 2 is a schematic diagram of carrier signal propagation characteristics;

Fig. 3 is the distribution situation of the simulation antenna of attitude calculation of the present invention;

Fig. 4 is the single epoch least squares iterative estimation process curve of the present invention;

Fig. 5a is a curve diagram of the calculation results of the attitude angle of the carrier in the static state of the present invention;

Fig. 5b is a curve diagram of the solution result of the attitude angle of the carrier of the present invention in the state of rotation at a fixed angular rate;

Fig. 5c is a curve diagram of the calculation results of the attitude angle of the carrier in the sinusoidal rotation state of the present invention;

Fig. 6a is a curve diagram of the error calculation error of the carrier of the present invention in the static state attitude angle;

Fig. 6b is a curve diagram of the attitude angle calculation error curve of the carrier of the present invention in a fixed angular rate rotation state;

Fig. 6c is a curve diagram of the attitude angle calculation error curve of the carrier of the present invention in a sinusoidal rotation state.

Detailed ways

The present invention is described in detail below in conjunction with accompanying drawing example:

1. The implementation process of the attitude angle direct solution algorithm

(1) Establish the relationship between the baseline vector and carrier observations

GPS carrier three-dimensional attitude measurement generally uses three receiving antennas to form two baseline vectors, so the attitude measurement system involved in this patent is composed of three antennas, of which antenna No. 0 is set as the reference antenna, and antenna No. 1 and No. 2 are slave antennas. Assuming that three antennas track M satellites at the same time, for any one of the satellites i, the carrier observation equation of the receiver to the satellite can be established:

表示天线至卫星的实际距离，单位：m；表示电离层延迟,单位：m；

表示对流层延迟,单位：m；δt_n表示接收机钟差，单位：s；δtⁱ表示卫星钟差，单位：s；表示整周模糊度；

表示载波测量噪声。In the formula, Indicates the carrier observation amount of the i-th star received by the antenna n; λ indicates the carrier wavelength, unit: m;

Indicates the actual distance from the antenna to the satellite, unit: m; Indicates ionospheric delay, unit: m;

Indicates tropospheric delay, unit: m; δt _n indicates receiver clock error, unit: s; δt ⁱ indicates satellite clock error, unit: s; Indicates the ambiguity of the whole week;

Indicates the carrier measurement noise.

A single-difference equation is established for each baseline to eliminate error terms such as ionospheric delay, tropospheric delay, and satellite clock bias. Since the length of the baseline is much smaller than that of the satellite, it is considered that the satellite line-of-sight vectors at both ends of the baseline are the same, and the satellite signal propagation is shown in Figure 2.

According to the carrier observation equation of a single receiver to a satellite, the single-difference observation equation of two stations to the same satellite i is established:

Further eliminate the clock error of the receiver, establish the single-difference equation of another satellite j observed by the three antennas at the same time, and combine the two single-difference equations to obtain the double-difference carrier observation equation:

则主从天线到卫星i的单差几何距离

（单位：m），等于基线向量在主天线对卫星i观测方向上投影长度的相反数，即：Relate the baseline vector b _n0 to the double-difference carrier observation equation. The line-of-sight vectors from the satellite to the three antennas are

Then the single-difference geometric distance from the master-slave antenna to satellite i

(unit: m), which is equal to the opposite number of the projection length of the baseline vector in the observation direction of the satellite i by the main antenna, that is:

${\mathrm{<span\; class="notranslate">\; r</span>}}_{\mathrm{<span\; class="notranslate">\; no</span>}\mathrm{<span\; class="notranslate">\; 0</span>}}^{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; =</span><span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; b</span>}}_{\mathrm{<span\; class="notranslate">\; no</span>}\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; \&CenterDot;</span>{\mathrm{<span\; class="notranslate">\; I</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}^{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 4</span>}<span\; class="notranslate">\; )</span>$

Therefore, the relationship between the double-difference geometric distance and the baseline in the double-difference equation can be expressed as:

${\mathrm{<span\; class="notranslate">\; r</span>}}_{\mathrm{<span\; class="notranslate">\; no</span>}\mathrm{<span\; class="notranslate">\; 0</span>}}^{\mathrm{<span\; class="notranslate">\; ij</span>}}<span\; class="notranslate">\; =</span><span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; b</span>}}_{\mathrm{<span\; class="notranslate">\; no</span>}\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; \&CenterDot;</span>{\mathrm{<span\; class="notranslate">\; I</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}^{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; b</span>}}_{\mathrm{<span\; class="notranslate">\; no</span>}\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; \&CenterDot;</span>{\mathrm{<span\; class="notranslate">\; I</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}^{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; =</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; I</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}^{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; I</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}^{\mathrm{<span\; class="notranslate">\; j</span>}}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; \&CenterDot;</span>{\mathrm{<span\; class="notranslate">\; b</span>}}_{\mathrm{<span\; class="notranslate">\; no</span>}\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 5</span>}<span\; class="notranslate">\; )</span>$

Then the carrier double-difference equation is transformed into:

是双差载波相位测量值，为已知量；b_n0是待求的三维基线向量，单位：m；

是双差整周模糊度，为未知整数。（2）建立姿态角解算模型Equation (6) gives the double difference The relationship with the baseline vector b _n0 . In the formula,

is the double-difference carrier phase measurement value, which is a known quantity; b _n0 is the three-dimensional baseline vector to be obtained, unit: m;

is the double-differenced integer ambiguity, which is an unknown integer. (2) Establish an attitude angle calculation model

In order to simplify the calculation steps of the attitude measurement algorithm, avoid the estimation error introduced in the intermediate process, and improve the calculation accuracy, the relationship between the attitude angle and the double-difference observation is directly established as follows:

Assume that the antenna of the GPS attitude measurement system is fixedly connected with the carrier, that is, the coordinates of the antenna under the carrier system are unchanged and known, and they are r _{0, B} , r _1,B , r _2,B . Two baseline vectors are formed as b _n =r _n,B -r _0,B , where n=1,2. In the geographic coordinate system,

Using the conversion relationship of the coordinate system, the baseline vector is expressed in the carrier system through the attitude matrix, then the formula (6) is transformed into the formula (7):

是地理系到地心系的姿态转换矩阵。In the formula, R is the attitude matrix of the geographic system to the carrier system,

is the attitude transformation matrix from the geographic system to the geocentric system.

When the double-difference integer ambiguity is determined, only three attitude angles are unknown in the equation, and the attitude can be directly estimated by using the relationship between the carrier double-difference observation and the Euler angle in (7) to establish the observation equation matrix horn.

M-1 double-difference equations can be established for each baseline, and the equation matrix form is:

In formula (8), set No. 1 satellite as the reference satellite. According to the relationship between the attitude angle and the equation, the equation is linearized, and then the attitude angle is estimated by the method of least squares, and finally a stable result is obtained.

(3) Least square estimation attitude angle solution process

The three rotation attitude angles are used as unknowns and included in the attitude matrix R, and the attitude matrix has the following form:

$\mathrm{<span\; class="notranslate">\; R</span>}<span\; class="notranslate">\; =</span>\left(\begin{array}{ccc}\mathrm{<span\; class="notranslate">\; CrCy</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; SrSpSy</span>}& \mathrm{<span\; class="notranslate">\; CrSy</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; QUR</span>}& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; SrCp</span>}\\ <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; CpS</span>}& \mathrm{<span\; class="notranslate">\; CpCy</span>}& \mathrm{<span\; class="notranslate">\; Sp</span>}\\ \mathrm{<span\; class="notranslate">\; SrCy</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; CrSpSy</span>}& \mathrm{<span\; class="notranslate">\; SrS</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; CrSpCy</span>}& \mathrm{<span\; class="notranslate">\; CrCp</span>}\end{array}\right)<span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 9</span>}<span\; class="notranslate">\; )</span>$

In the formula, S represents sin; C represents cos; y, r, p are respectively the yaw angle (Yaw) of the carrier around the z-axis of the local horizontal coordinate system, the roll angle (Roll) around the y-axis, and the roll angle (Roll) around the x-axis of the carrier, respectively. The pitch angle of the rotation (Pitch). Since the model has a nonlinear relationship with the attitude angle, nonlinear least squares is used for estimation to determine the attitude angle.

It is considered that the ambiguity on the right side of the equation (8) is known, and the equation is expressed as:

y _t =f(x)

in, Indicates the double-difference carrier observations at epoch time t after the whole-week ambiguity is solved; $f\left(x\right)={\mathrm{\&lambda;}}^{-1}\left[\begin{array}{c}-{({\mathrm{<figure-callout\; id="0"\; label="I"\; filenames="HDA00002467225800022.png,HDA00002467225800032.png"\; state="\{\{state\}\}">I</figure-callout>}}_0^2-{\mathrm{<figure-callout\; id="0"\; label="I"\; filenames="HDA00002467225800022.png,HDA00002467225800032.png"\; state="\{\{state\}\}">I</figure-callout>}}_0^1)}^T\\ -{({\mathrm{<figure-callout\; id="0"\; label="I"\; filenames="HDA00002467225800022.png,HDA00002467225800032.png"\; state="\{\{state\}\}">I</figure-callout>}}_0^3-{\mathrm{<figure-callout\; id="0"\; label="I"\; filenames="HDA00002467225800022.png,HDA00002467225800032.png"\; state="\{\{state\}\}">I</figure-callout>}}_0^1)}^T\\ ...\\ -{(I_0^m-I_0^1)}^T\end{array}\right]A_R^{\mathrm{E.}}{R}^T\left(x\right)b_{\mathrm{no}}$ Represents a function with x as an independent variable, x=[y r p] ^T .

在初始姿态角处对方程f(x)进行一阶泰勒级数展开，忽略高阶项。令 $\mathrm{\&Delta;x}={\left[\begin{array}{ccc}\mathrm{\&Delta;y}& \mathrm{\&Delta;r}& \mathrm{\&Delta;p}\end{array}\right]}^T=\begin{array}{}\end{array}{\left[\begin{array}{ccc}y-\hat{y}& r-\hat{r}& p-\hat{p}\end{array}\right]}^T$ 为状态变量，则非线性方程组可近似转化为以下用矩阵形式表达的线性方程组：Set the initial attitude angle

A first-order Taylor series expansion of the equation f(x) is performed at the initial attitude angle, ignoring higher-order terms. make $\mathrm{\&Delta;x}={\left[\begin{array}{ccc}\mathrm{\&Delta;y}& \mathrm{\&Delta;r}& \mathrm{\&Delta;p}\end{array}\right]}^T=\begin{array}{}\end{array}{\left[\begin{array}{ccc}\mathrm{the\; y}-\widehat{\mathrm{the\; y}}& r-\hat{r}& p-\hat{p}\end{array}\right]}^T$ is the state variable, the nonlinear equations can be approximately transformed into the following linear equations expressed in matrix form:

H·Δx= _Δyt (10)

即第k次迭代中状态向量的估计结果x_k相对测量值

的偏差。观测矩阵H表示为：In the formula,

That is, the estimated result x _k of the state vector in the kth iteration is relative to the measured value

deviation. The observation matrix H is expressed as:

$\mathrm{<span\; class="notranslate">\; h</span>}<span\; class="notranslate">\; =</span>\left[\begin{array}{ccc}& <span\; class="notranslate">\; .</span>& \\ & <span\; class="notranslate">\; .</span>& \\ & <span\; class="notranslate">\; .</span>& \\ <span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span>& <span\; class="notranslate">\; -</span>{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; I</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}^{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; I</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}^{\mathrm{<span\; class="notranslate">\; 1</span>}}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; T</span>}}{\mathrm{<span\; class="notranslate">\; A</span>}}_{\mathrm{<span\; class="notranslate">\; R</span>}}^{\mathrm{<span\; class="notranslate">\; E.</span>}}<span\; class="notranslate">\; (</span>\frac{{<span\; class="notranslate">\; \&PartialD;</span>\mathrm{<span\; class="notranslate">\; R</span>}}^{\mathrm{<span\; class="notranslate">\; T</span>}}}{<span\; class="notranslate">\; \&PartialD;</span>\mathrm{<span\; class="notranslate">\; x</span>}}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; b</span>}}_{\mathrm{<span\; class="notranslate">\; no</span>}}& <span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span><span\; class="notranslate">\; .</span>\\ & <span\; class="notranslate">\; .</span>& \\ & <span\; class="notranslate">\; .</span>& \\ & <span\; class="notranslate">\; .</span>& \end{array}\right]$

表示姿态矩阵对各姿态角的偏导数，形式如下：In the formula,

Indicates the partial derivative of the attitude matrix to each attitude angle, the form is as follows:

$\frac{{<span\; class="notranslate">\; \&PartialD;</span>\mathrm{<span\; class="notranslate">\; R</span>}}^{\mathrm{<span\; class="notranslate">\; T</span>}}}{<span\; class="notranslate">\; \&PartialD;</span>\mathrm{<span\; class="notranslate">\; the\; y</span>}}<span\; class="notranslate">\; =</span>\left(\begin{array}{ccc}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; CrSy</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; CaS</span>}& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; CyCp</span>}& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; SrS</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; CySpCr</span>}\\ \mathrm{<span\; class="notranslate">\; CyCr</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; SySpSr</span>}& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; SyCp</span>}& \mathrm{<span\; class="notranslate">\; <figure-callout\; id="0"\; label="CaS\; +\; SySpCr"\; filenames="HDA00002467225800022.png,HDA00002467225800032.png"\; state="\{\{state\}\}">CaS</figure-callout></span><figure-callout\; id="0"\; label="CaS\; +\; SySpCr"\; filenames="HDA00002467225800022.png,HDA00002467225800032.png"\; state="\{\{state\}\}">\; </figure-callout>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; SySpCr</span>}\\ \mathrm{<span\; class="notranslate">\; 0</span>}& \mathrm{<span\; class="notranslate">\; 0</span>}& \mathrm{<span\; class="notranslate">\; 0</span>}\end{array}\right)$

$\frac{{<span\; class="notranslate">\; \&PartialD;</span>\mathrm{<span\; class="notranslate">\; R</span>}}^{\mathrm{<span\; class="notranslate">\; T</span>}}}{<span\; class="notranslate">\; \&PartialD;</span>\mathrm{<span\; class="notranslate">\; r</span>}}<span\; class="notranslate">\; =</span>\left(\begin{array}{ccc}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; CaS</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; <figure-callout\; id="0"\; label="SySpCr"\; filenames="HDA00002467225800022.png,HDA00002467225800032.png"\; state="\{\{state\}\}">SySpCr</figure-callout></span>}& \mathrm{<span\; class="notranslate">\; 0</span>}& \mathrm{<span\; class="notranslate">\; CyCr</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; SySpSr</span>}\\ <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; <figure-callout\; id="0"\; label="SySr\; +\; CySpCr"\; filenames="HDA00002467225800022.png,HDA00002467225800032.png"\; state="\{\{state\}\}">SySr</figure-callout></span><figure-callout\; id="0"\; label="SySr\; +\; CySpCr"\; filenames="HDA00002467225800022.png,HDA00002467225800032.png"\; state="\{\{state\}\}">\; </figure-callout>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; CySpCr</span>}& \mathrm{<span\; class="notranslate">\; 0</span>}& \mathrm{<span\; class="notranslate">\; SyCr</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; CaS</span>}\\ <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; CpCr</span>}& \mathrm{<span\; class="notranslate">\; 0</span>}& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; CpSr</span>}\end{array}\right)$

$\frac{{<span\; class="notranslate">\; \&PartialD;</span>\mathrm{<span\; class="notranslate">\; R</span>}}^{\mathrm{<span\; class="notranslate">\; T</span>}}}{<span\; class="notranslate">\; \&PartialD;</span>\mathrm{<span\; class="notranslate">\; p</span>}}<span\; class="notranslate">\; =</span>\left(\begin{array}{ccc}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; SyCpSr</span>}& \mathrm{<span\; class="notranslate">\; SySp</span>}& \mathrm{<span\; class="notranslate">\; SyCpCr</span>}\\ \mathrm{<span\; class="notranslate">\; CyCpSr</span>}& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; CySp</span>}& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; CyCpCr</span>}\\ \mathrm{<span\; class="notranslate">\; SpR</span>}& \mathrm{<span\; class="notranslate">\; Cp</span>}& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; SpCr</span>}\end{array}\right)$

Then the least squares solution is:

Δx＝(H ^T H) ^-1 H ^T Δy _t (11)

By establishing the nonlinear least squares estimation equation above, the estimation of unknown parameters can be completed under the condition of a given initial attitude angle. That is, the direct estimation of the attitude angle by using the carrier double-difference observation equation is realized.

2. Error analysis of the direct solution method of attitude angle

联立卫星载波双差观测方程，通过最小二乘方法进行超定方程求解。超定方程的一般形式为：Assuming that all satellite carrier observation errors are the same as and are independent of each other, the mean value is 0, then the error of the double-difference measurement value is

Simultaneous satellite carrier double-difference observation equations are solved by the least squares method for overdetermined equations. The general form of the overdetermined equation is:

Ax=b (12)

Then, there is a unique least squares solution:

${\mathrm{<span\; class="notranslate">\; A</span>}}^{\mathrm{<span\; class="notranslate">\; h</span>}}\mathrm{<span\; class="notranslate">\; Ax</span>}<span\; class="notranslate">\; =</span>{\mathrm{<span\; class="notranslate">\; A</span>}}^{\mathrm{<span\; class="notranslate">\; h</span>}}\mathrm{<span\; class="notranslate">\; b</span>}<span\; class="notranslate">\; \&DoubleRightArrow;</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; =</span>{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; A</span>}}^{\mathrm{<span\; class="notranslate">\; h</span>}}\mathrm{<span\; class="notranslate">\; A</span>}<span\; class="notranslate">\; )</span>}^{<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}}{\mathrm{<span\; class="notranslate">\; A</span>}}^{\mathrm{<span\; class="notranslate">\; h</span>}}\mathrm{<span\; class="notranslate">\; b</span>}$

In formula (12), the error of b is δb, and the error of A is δA, and their influence on the solution of the equation is proportional to the square of the condition number of A, that is, the condition number of the overdetermined equation will increase in a quadratic relationship:

cond(A ^H A)=[cond(A)] ² (13)

According to the law of error propagation, the error of the least squares solution is

Due to the method of directly solving the attitude angle based on the double-difference equation of the carrier wave, only one-step least squares estimation is performed to complete the calculation of the attitude angle. From the previous analysis, we can see that the estimation error of the vector x is:

cov[Δx]=E[ΔxΔx ^T ]=(H ^T H) ^-1 H ^T E[ΔyΔy ^T ]H(H ^T H) ^-1 (14)

均值为0，且相互独立。因此，其中I为单位阵。The variance of the measurement error Δy is

have a mean of 0 and are independent of each other. therefore, where I is the identity matrix.

define a dimensionless matrix but

Therefore, the attitude angle estimation error is related to the measurement error as well as the measurement matrix and baseline length. Among them, the longer the baseline, the higher the accuracy without considering other error factors; and from formula (15), it can be seen that the observation matrix is only related to the geometric position of the satellite. When the geometric distribution is better, the influence of measurement error on the estimation of attitude angle smaller.

Claims (1)

1. the attitude angle direct solving method based on the two eikonal equations of gps carrier, is characterized in that comprising the following steps:

(1) adopt three GPS receivers to measure antenna collection gps satellite signal and obtain carrier wave observation information

using the carrier wave observation information of No. 0 antenna in place, right angle as base station signal;

(2) choose public satellites in view i, i=1,2, M, and using the satellite of elevation angle maximum in satellites in view as with reference to satellite, receiver obtains after triantennary carrier wave observation information, take antenna for base station observation information as benchmark, respectively the carrier signal of the relatively same visible satellite of each antenna is done to poor single poor carrier phase observed quantity of obtaining to base station signal

n=0,1,2, then corresponding different antennae, does poor two poor carrier wave observed quantities that obtain by the poor carrier information of each satellite list to the poor carrier wave observed quantity of list of reference satellite

(3) utilize formula

shown in the relational expression of carrier wave observed quantity and attitude angle, set up solving of attitude equation, wherein,

represent the two poor carrier wave observed quantity of n antenna; λ represents carrier wavelength, unit: m;

represent two poor integer ambiguities;

represent that satellite i is to the sight line vector of reference antenna; R represents that geography is tied to the attitude matrix of carrier system;

represent the geographical attitude transition matrix that is tied to the earth's core system; Y, r, p represents three attitude angle of carrier, is respectively course angle, roll angle, the angle of pitch, i.e. equation unknown number; b _nrepresenting the baseline vector under carrier coordinate system, is known quantity, unit: m;

(4) attitude measurement system Real-time Obtaining list carrier signal epoch is carried out two poor processing, by step (1)-(4), calculates in real time the attitude information of carrier;

Described antenna adopts right angle orthogonal formula layout, i.e. form right angle triangle projective planum, and mounting distance is more than 10m;

Described solving of attitude equation calculation process is:

1) due to formula

equal sign the right integer ambiguity is known, therefore will

be expressed as:

y _t=f(x)

Wherein,

represent that integer ambiguity solves rear epoch of the two poor carrier wave observed quantity of t constantly; $f\left(x\right)={\mathrm{\&lambda;}}^{-1}\left[\begin{array}{c}{-(I_0^2-I_0^1)}^T\\ {-(I_0^3-I_0^1)}^T\\ ...\\ {-(I_0^M-I_0^1)}^T\end{array}\right]A_R^E{R}^T\left(x\right)b_n$ The function that x is independent variable is take in expression, x=[y r p] ^t;

2) provide initial attitude angle $x_0={\left[\begin{array}{ccc}\hat{y}& \hat{r}& \hat{p}\end{array}\right]}^T,$ At initial attitude angle, place carries out single order Taylor series expansion to Equation f (x), order $\mathrm{\&Delta;x}=\begin{array}{c}{\left[\begin{array}{ccc}\mathrm{\&Delta;y}& \mathrm{\&Delta;r}& \mathrm{\&Delta;p}\end{array}\right]}^T={\left[\begin{array}{ccc}y-\hat{y}& r-\hat{r}& p-\hat{p}\end{array}\right]}^T\end{array}$ For state variable, equation simplification is following form:

H·Δx=Δy _t

In formula,

the i.e. estimated result x of state vector in the k time iteration _krelative measurement value

deviation; Observing matrix H is:

$H=\left(\begin{array}{ccc}& \&CenterDot;& \\ & \&CenterDot;& \\ & \&CenterDot;& \\ ...& -{(I_0^i-I_0^1)}^T{A}_R^E\left(\frac{{\&PartialD;R}^T}{\&PartialD;x}\right)b_n& ...\\ & \&CenterDot;& \\ & \&CenterDot;& \\ & \&CenterDot;& \end{array}\right)$

Wherein,

represent the partial derivative of attitude matrix to each attitude angle;

3) adopt least-squares estimation solution overdetermined equation H Δ x=Δ y _t, its solution is:

Δx=(H ^TH) ^-1H ^TΔy _t

4) set iteration cut-off threshold, estimate attitude angle final value x:

$x=\left[\begin{array}{c}y\\ r\\ p\end{array}\right]=\begin{array}{c}\left[\begin{array}{c}\hat{y}+\mathrm{\&Delta;y}\\ \hat{r}+\mathrm{\&Delta;r}\\ \hat{p}+\mathrm{\&Delta;p}\end{array}\right].\end{array}$

Images