WO2009028929A1

WO2009028929A1 - Device and method for calculating position of mobile station

Info

Publication number: WO2009028929A1
Application number: PCT/NL2007/050427
Authority: WO
Inventors: Mitsuru Nakamura; Yasuhiro Tajima; Kazunori Kagawa; Peter Teunissen; Peter Joosten
Original assignee: Toyota Jidosha Kabushiki Kaisha; Delft University
Priority date: 2007-08-29
Filing date: 2007-08-29
Publication date: 2009-03-05

Abstract

A disclosed device for calculating a position of a mobile station using an integer ambiguity comprises a testing unit configured to evaluate reliability of observables on epoch basis, the observables including at least code data and carrier phase obtained at a single epoch. The testing unit includes a first testing unit configured to test residuals as a whole, and a second testing unit configured to test one or more items of the residuals separately.

Description

DESCRIPTION

DEVICE AND METHOD FOR CALCULATING POSITION OF MOBILE STATION

TECHNICAL FIELD

The present invention relates to a device and a method for calculating the position of the mobile station.

BACKGROUND ART

Carrier phase ambiguity resolution is the key to fast and high-precision GNSS (Global Navigation Satellite System) positioning and navigation. It applies to a great variety of current and future models of GPS, modernized GPS and Galileo.

Ambiguity resolution is the process of resolving the unknown cycle ambiguities of double difference carrier phase data (or single difference carrier phase data) as integers . Various methods using various models are known or proposed to estimate the ambiguity. Integer ambiguity estimation using the least-squares method, the Kalman filter, or the like is well-known.

Despite the differences in application of the various methods, their ambiguity resolution problems are intrinsically the same. One of the problems is that they cannot provide the user or analyst with tools to evaluate the quality of the integer solution. It is of importance to be able to evaluate the quality of the integer solution, since unsuccessful ambiguity resolution, when remaining unnoticed, may lead to unacceptable errors in the positioning result.

Unsuccessful ambiguity resolution may be caused by various error factors such as outliers in code data, slips in carrier phase data (referred to as cycle-slip) , or the like. Some of these error factors cannot be removed even when double difference carrier phase data are used. The error factors may also be hidden in the integer ambiguity estimating process. Therefore, it is of importance to be able to identify or narrow the error factors if it is likely that there are errors in the positioning result or lack of reliability in the observation data.

Papers titled "Quality Control in Integrated Navigation Systems" and ^WAN INTEGRITY AND QUALITY CONTROL PROCEDURE FOR USE IN MϋLTI SENSOR INTEGRATION" published in 1990s by P. J.G. Teunissen discloses methods of evaluating the quality of the ambiguity resolution.

The disclosed methods are directed at evaluating the quality of the ambiguity resolution conducted by using the recursive Kalman filter algorithm. Of course, it is of importance to be able to evaluate the quality of the ambiguity resolution in such a multi-epoch application; however, the multi-epoch application has a disadvantage in that the ambiguity resolution at each epoch is affected by observables obtained at the previous epochs. Thus, according to these methods, it is difficult to evaluate the quality of the observables (and thus the quality of the ambiguity resolution conducted by using the observables) on epoch basis.

In particular, in a single epoch application where the ambiguity resolution is conducted by using the observables obtained at a single epoch, it is of importance to be able to evaluate the quality of the integer solution on epoch basis, because in the single epoch application the quality of the integer solution may be significantly affected by observables obtained at each epoch.

DISCLOSURE OF INVENTION It is an object of the present invention to provide a device and a method for calculating the position of a mobile station which can appropriately evaluate reliability of observables on epoch basis.

In order to achieve the above-mentioned objects, according to one aspect of the present invention a method for calculating the position of a mobile station using an integer ambiguity is provided, comprising; a testing unit configured to evaluate reliability of observables on epoch basis, said observables including at least code data and carrier phase obtained at a single epoch, wherein said testing unit includes a first testing unit configured to test residuals as a whole, and a second testing unit configured to test one or more items of the residuals separately.

In this aspect, the first testing unit may be configured to determine whether the sum of squares of the residuals is larger than a reference value, and the second testing unit is configured to perform the^* test . if the sum of squares of the residuals is larger than the reference value.

Preferably, the residuals are derived using the integer ambiguity which is estimated based on the observables . Preferably, the device further includes a third testing unit configured to evaluate reliability of the integer ambiguity estimated based on the observables.

Preferably, the device further a reliability outputting unit configured to output the reliability of the calculated position of the mobile station, said reliability being derived based on the test result.

Preferably, items of the observables to be used for estimating the integer ambiguity are selected according to the test result.

Preferably, the integer ambiguity is estimated using the observables obtained at a single epoch.

Preferably, the integer ambiguity is estimated after modifying a float solution according to the test result.

Preferably, if the test result at a certain epoch does not meet predetermined criteria, the integer ambiguity at the epoch is estimated using an item of observables obtained at another epoch.

According to another aspect of the present invention a method for calculating a position of a mobile station using an integer ambiguity is provided, comprising; a float solution testing unit configured to evaluate reliability of a float solution of the integer ambiguity derived based on observables which are obtained at a single epoch; and an integer solution testing unit configured to evaluate reliability of an integer solution of the integer ambiguity estimated based the float solution.

In this aspect, preferably, if the test result on the float solution doesn't meet predetermined criteria, the estimation of the integer solution based on said float solution is ceased. Preferably, the tests are performed separately on the float solution derived based on double difference observables and on the float solution derived based on single difference observables. Preferably, float solutions to be used for estimating the integer ambiguity are selected according to the test result.

Preferably, the tests are performed separately on the integer solution derived based on double difference observables and on the integer solution derived based on single difference observables.

Preferably, integer solutions to be used for calculating the position of a mobile station are selected according to the test result.

According to another aspect of the present invention, a method for calculating a position of a mobile station using an integer ambiguity is provided, comprising; a float solution testing step of evaluating reliability of a float solution of the integer ambiguity derived based on observables which are obtained at a single epoch; and an integer solution testing step of evaluating reliability of an integer solution of the integer ambiguity estimated based the float solution.

According to another aspect of the present invention, a computer program is provided for calculating a position of a mobile station using an integer ambiguity, said computer program causing a computer to execute the following step; a testing step for evaluating reliability of observables including at least code data and carrier phase data, said observables being obtained at a single epoch, wherein said testing step includes a first testing step in which residuals of variables with respect to the observables as a whole are evaluated, and a second testing step in which a statistical test for detecting at least one of error factors is performed.

According to another aspect of the present invention a computer program is provided for calculating the position of a mobile station using an integer ambiguity, said computer program causing a computer to execute the following step: a float solution testing step of evaluating reliability of a float solution of the integer ambiguity derived based on observables obtained at a single epoch; and an integer solution testing step of evaluating reliability of an integer solution of the integer ambiguity estimated based the float solution.

BRIEF DESCRIPTION OF THE DRAWINGS

These and other objects, features, and advantages of the present invention will become more apparent from the following detailed description of preferred embodiments given with reference to the accompanying drawings, in which:

Fig. 1 is a schematic diagram of a carrier phase GPS positioning device according to an embodiment of the present invention; Fig. 2 is a diagram showing a configuration of the carrier phase GPS positioning device in FIG. 1;

Fig. 3 is a block diagram showing an embodiment of a carrier phase GPS positioning device 34 installed in the mobile station 30 according to the present invention; Fig. 4 is a diagram illustrating the definitions of coordinate systems used in descriptions;

Fig. 5 is a flowchart illustrating the input test according to an embodiment of the present embodiment; Fig. 6 is a flowchart illustrating a method of giving reliability to the calculated position of the mobile station 30;

Fig. 7 is a flowchart illustrating a method of calculating the position of the mobile station 30 for the epoch when the flag 2 and/or the flag 1 is set to "1";

Fig. 8A shows the curve of the positions calculated according to the conventional single epoch application; Fig. 8B shows the curve of the positions calculated according to an embodiment of the present invention;

Fig. 9 is a flowchart illustrating a method of evaluating the quality of the calculated position of the mobile station 30 according to another embodiment of the present invention; and

Fig. 10 is a flowchart illustrating a method of improving the accuracy in the calculated position of the mobile station 30 based on the test results according to an embodiment of the present invention.

BEST MODE FOR CARRYING OUT THE INVENITON

Hereafter, the preferred embodiments according to the present invention are explained with reference to the drawings.

FIG. 1 is a schematic diagram of a carrier phase GPS positioning system according to an embodiment of the present invention.

As illustrated in FIG. 1, the carrier phase GPS positioning system includes GPS satellites 10 orbiting around the earth, a reference station 20 located at a fixed position (known position) , and a mobile station 30 that is on the earth and is able to move on the earth. Each of the GPS satellites 10 broadcasts navigation messages toward the earth continuously. The navigation messages include orbital information of the corresponding GPS satellite 10, clock correction values, and correction coefficients of the ionospheric layer. The navigation messages are spread using a PRN code, such as C/A code or P code. At present, the C/A code is carried on a Ll carrier (frequency: 1575.42 MHz) and the P codes are carried on the Ll carrier and a L2 carrier (frequency: 1227,6 MHz), and are broadcast toward the earth.

Presently, there are 24 GPS satellites orbiting the earth at an altitude of 20,000 km. Every four GPS satellites are equally arranged on one of six orbital planes of the earth, which are inclined 55 degrees relative to each other. Therefore, at least five satellites are always observable from a position as long as the position is open to the sky, no matter where the position is on the earth.

In FIG. 2, the mobile station 30 has a GPS receiver 32. In the GPS receiver 32, there is an oscillator (not illustrated) having an oscillating frequency equal to the carrier frequency of the GPS satellite 10. The GPS receiver 32 converts an electromagnetic wave, which is emitted from the GPS satellite 10 and is received by the GPS receiver 32 via a GPS antenna 32a, performs C/A code synchronization using the C/A codes generated in the GPS receiver 32, and extracts the navigation messages .

The GPS receiver 32 calculates a carrier phase accumulation value Φi_U of the carrier waves from the GPS satellites 1Oi. Here, in the phase accumulation value Φi_u, the subscript i(=l, 2, ...) represents the numbers assigned to the GPS satellite 1Oi, and the subscript u represents _Q—

that the accumulation value is calculated on the side of the mobile station 30.

The phase accumulation value Φχ_u can be described as the difference between a phase Θi_U(t) of the oscillator at the time t of receiving the carrier wave and a phase Θj_.u(t-τ) of the carrier wave when the satellite signal from the GPS satellite 1Oi is generated, as shown by the following formula (1) .

Here, τ_u represents travel time from the GPS satellite 10 to the GPS receiver 32, and ε i_U represents noise (uncertainty) . Further, at the time when starting observing the phase difference, the GPS receiver 32 can accurately determine the carrier phase within one wavelength of the carrier wave, but cannot determine what number of the wavelengths the present wavelength is. For this reason, in the phase accumulation value Φj_.u(t), as shown in the formula (1), there is an uncertainty factor Ni_U, known as "integer ambiguity".

The phase accumulation value Φi_U can be calculated for both of .the Ll wave and the L2 wave. In this case, two phase accumulation values ~Φ-_LU can be obtained at each epoch.

The GPS receiver 32 also calculates a pseudo range Pi_U(t) based on the C/A codes carried by the carrier waves from the GPS satellites 1Oi. The pseudo range p i_u(t) calculated here includes errors such as a range error, as shown by the following formula (1-1) .

Here, c denotes speed of light, and b, which is also referred to as a clock bias, corresponds to a range error due to a clock error in GPS receiver 32.

Similarly, the pseudo ranges Pi_U(t) can be measured using the P codes carried by the Ll wave and the L2 wave. In this case, two pseudo ranges Pi_U(t) based on the P codes can be obtained at each epoch.

The mobile station 30 also includes a communication device 33, such as a mobile phone. As described below, the communication device 33 is capable of communicating with a communication facility 23 installed on the reference station 20 side, such as a base station for mobile phones, by bi-directional communication.

A GPS receiver 22 having a GPS antenna 22a is installed in the reference station 20. The GPS receiver 22, the same as the GPS receiver 32 in the mobile station 30, calculates a carrier phase accumulation value Φi_b at time t based on the carrier waves from the GPS satellites 10i, as shown by the following formula (2) .

Here, Nu₃ is an integer ambiguity, and e _H3 represents noise (uncertainty) . In the ptiase accumulation value Φib, the subscript b represents that the accumulation value is calculated on the side of the reference station 20.

Similarly, the phase accumulation value Φ_ϋ can be measured using both of the Ll wave and the L2 wave. In this case, two phase accumulation values Φjjj can be obtained at each epoch.

The GPS receiver 22 also calculates a pseudo range pa_>{t) based on the C/A codes carried by the carrier waves from the GPS satellites 10i, as shown by the following formula (2-1) .

Similarly, the pseudo range p^ (t) can be measured using the P codes carried by the Ll wave and the

L2 wave. In this case, two the pseudo ranges Pϋ_> (t) based on the P codes can be obtained at each epoch. The reference station 20 transmits observation data including the obtained carrier phase accumulation value Φi_b and the pseudo range p a, to the mobile station 30 via the communication facility 23. More than one reference station 20 may be installed in a specified region. As illustrated in FIG. 2, each of the reference stations 20 may be connected to one or more communication facilities 23 through the Internet or other networks, or a communication facility 23 may be installed in each of the reference stations 20. In the former case, as long as the mobile station 30 is able to communicate with the communication facility 23, the mobile station 30 can obtain the information received by each of the reference stations 20.

FIG. 3 is a block diagram showing an embodiment of a carrier phase GPS positioning device 34 installed in the mobile station 30 according to an embodiment of the present invention.

The carrier phase GPS positioning device 34 of the present embodiment includes a calculation unit 40, which is connected to the GPS receiver 32 and the communication device 33, and further, to various sensors 50 provided in the mobile station 30. The calculation unit 40 may also be installed in the GPS receiver 32. When the mobile station is a vehicle, the GPS receiver 32, the calculation unit 40 and/or the communication device 33 may also be mounted in a navigation device.

The calculation unit 40 consists mainly of a 5 microcomputer. A microcomputer consists of a CPU for processing data, memory in which computer programs and data for implementing the process described below or processed data are stored, interfaces, and the like.

The calculation unit 40, based on the orbital

10 information in the navigation messages received by the GPS receiver 32, calculates positions (X±(t), Yi (t) , Zi (t) ) of all observable GPS satellites 1Oi at time t in a global coordinate system.

FIG. 4 is a diagram illustrating the

15 definitions of coordinate systems used in the following description. FIG. 4 shows relationships between the global coordinate system, a local coordinate system, and a body coordinate system.

The body coordinate system is defined on the

20 body of the vehicle.

Because movement of each of the GPS satellites 10 is confined to an orbital plane passing through the center of gravity of the earth, and the orbit of each of the GPS satellites 10 is an ellipse with the center of

25 gravity of the earth as a focus, positions of each of the GPS satellites 10 in the orbital plane can be calculated by successive numerical solutions of Kepler's equation. Because the orbital planes of each of the GPS satellites and the equatorial plane in the global coordinate system

30. satisfy a rotational transformation relationship, positions (Xi(t), Yi(t), Zi(t)) of the GPS satellites 10 at the time t of receiving the carrier waves can be calculated by three dimensional rotational coordinate transformation of the positions of the GPS satellites 10 on the orbital planes.

The calculation unit 40, based on the output signals of the various sensors 50 input periodically, calculates quantities related to movement of the mobile station 30.

For example, if the mobile station 30 is a vehicle, the calculation unit 40 calculates the speed Vx (t) (speed in the forward and backward directions) and Vy (t) (speed in the right and left directions) at the time t of receiving the carrier wave based on output signals from various sensors 50, for example, two wheel speed sensors mounted on the driven wheels of the vehicle, a steering sensor, a yaw rate sensor, a left and right G acceleration sensor, and an azimuth meter.

Because the speed vector (Vx (t) , Vy (t) ) of the vehicle is defined in the body coordinate system whose origin is on the body of the vehicle, it is necessary for the calculation unit 40 to transform the speed vector (Vx (t), Vy (t) ) from the body coordinate system to the global coordinate system via the local coordinate system. Usually, the rotational transformation of coordinates can be performed by using Euler angles. In the present embodiment, the transformation from the body coordinate system to the local coordinate system is performed using only a yaw angle φ (t) since the roll angle and pitch angle are small. Depending on the situation, the roll angle and pitch angle may also be considered, or the yaw angle may also be ignored. The transformation from the local coordinate system to the global coordinate system is performed by using the longitude φ (t) and latitude λ (t) of the position of the vehicle.

Such dynamic information indicating dynamic behavior of the mobile station 30 may be used to generate or modify models for calculating the position of the mobile station 30, or may be used to modify the calculated position of the mobile station 30. The calculation unit 40 performs an input test on the observation data obtained at the mobile station 30 and the reference station 20, or other sensor data from the various sensors 50. The observation data include the carrier phase accumulation values Φi_U, Φϋ_>/ and the pseudo ranges p i_U, p ϋ.

The input test is performed for checking reliability of the observation data, as described below. The input test is performed in two stages. The first stage is for diagnosing whether an unspecified model error has occurred. If the error indeed present, then the second stage is initiated for detecting a potential source (error factor) of the model error detected in the first stage .

Fig. 5 is a flowchart illustrating the input test according to an embodiment of the present embodiment. The process routine shown in Fig. 5 is performed at each epoch independently. The process routine shown in Fig. 5 is executed by the calculation unit 40.

In step SlOO, a residual vector v consisting of residuals is derived using an appropriate model. The various models known or proposed are applicable. In this embodiment, the simple model referred to as "geometry-free model" and/or "geometry-based model" is used. The geometry-free model is given as follows.

Or the

follows.

Here, Φ (t) is an observable related to the carrier phase accumulation value, and p (t) is an observable related to the pseudo range. D_p(t), D_φ (t) are each an unknown quantity related to a station-satellite range. N is an unknown quantity related to the integer ambiguity. λ denotes a wavelength of the carrier wave. n_p (t) and n_φ(t) each denotes an observation error, including the clock error in the case of the former equation.

According to this embodiment, since the model is established by using the pseudo ranges pi_u, p_ϋ,in addition to the carrier phase accumulation values Φi_U/ Φ i_b , it becomes possible to have the redundancy that is required to accurately derive the optimal solutions, even if the number of the GPS satellites 10 that is being tracked is relatively small. In other words, it becomes possible to perform the integer resolution on epoch basis, even if the number of the GPS satellites 10 that is being tracked is relatively small.

In one embodiment, the calculation unit 40 calculates^"double differences of the observation data, that is, double differences of the carrier phase accumulation values and the pseudo ranges.

The double difference of the carrier phase accumulation values related to the GPS satellites 10_j and 1O_h (j is not equal to h) at time t (a certain epoch) can be expressed by the following formula (3) .

The double difference of the pseudo ranges related to the GPS satellites 1Oj and 1Oh (j is not equal to h) at time t (a certain epoch) can be expressed by the following formula (4) .

In this embodiment, the phase accumulation value Φ_jMw of the carrier phase accumulation values is used as observables Φ (t) in the model, and the double difference P_jhbuOf pseudo ranges is used as observables p (t) in the model .

In this embodiment, the least-squares method is applied to derive the residual vector v . The least- squar method is iven as follows.

Here, y denotes a vector of observables consisting of Φjhbu and Pjhbu_f x denotes a vector of unknown quantities consisting of D_p (t) , D_φ (t) , and N, n denotes an observation noise including errors which cannot be represented even by approximations, and W denotes an design matrix. The design matrix W is given as follows.

Here, W is a combination of an design matrix WX for the observables Φjhbu and an design matrix Wl for the observables Pjhbu- The design matrix W may be adapted by the solutions of the variables x in the case of the design matrix W being dependent on the variables x. In the case of double difference the same 5 satellites being observable at the mobile station 30 and the reference station 20 and observation data being generated from the Ll wave and the L2 wave, the number of the observables is 12, and the number of the unknown quantities is 11.

As is well known, the optimum solution x of x is given as follows.

As a result, the residual vector v is given as follows;

that in order to derive the residual vector v the single differences of the observation data or the observation data itself can be used as observables. In the case of the double phase difference, it is possible to eliminate influence of the initial phase of oscillators in the GPS receivers 22 and 32, and clock uncertainties. In the case of the single phase difference, it is possible to eliminate influence of the initial phase of oscillators in the GPS satellite 10, and the GPS clock error.

It is also noted that in order to derive the residual vector v , the combination of the single differences of the observation data and the double differences of the observation data can be used as observables.

In single difference application, the single difference of the carrier phase accumulation values related to the GPS satellites 1Oj at time t can be expressed by the following formula (5) .

The single difference of the pseudo ranges related to the GPS satellites 1Oj at time t can be expressed by the following formula (6) .

It is also noted that in order to derive the residual vector v , another model can be used. In another embodiment, a geometry-based model is used. In this embodiment, the distance between the GPS satellite 1Oi and the GPS receiver 22 or 32 equals the wavelength λ of the carrier wave multiplied by the phase accumulation value, and the double difference Φjh_bu of the phase accumulation satisfies the following formula (7) .

In formula (7), [Xj₃ (t) , Y_b(t), Z_b(t)] are coordinates (known) of the reference station 20 at time t in the global coordinate system, [X_u(t), Y_u(t), Z_u(t)] are coordinates (unknown) of the mobile station 30 at time t, and [X,(t)_f Yj (t), Zj (t)] and [X₁, (t), Y_h(t), Z_h(t)] are coordinates of the GPS satellites 10_j and 1O_h at time t calculated by the calculation unit 40. Njht_>u represents the double difference of the integer ambiguity, that is, Nj_h__>u

he double difference Pj_hbuOf the phase accumulation satisfies the following formula (8) .

In formula (8), [X_b(t), Y_b(t), Z_b(t)], [X_u(t), Yu(t), Z_u(t)], [Xj (t), Yj (t), Zj (t)] and [X_h(t), Y_h(t), Z_h (t) ] are the same as mentioned previously. bj_hbU represents the double difference of the clock biases. Similarly, in this embodiment, the least- squares method is applied to derive the residual vector v The least-squares method is given as follows.

Here, the vector of observables y consists of the double difference Φj_ht>u of the phase accumulation

(refer to formula (7)) and the double difference p ^_u of pseudo ranges (refer to formula (8)). x denotes a vector of unknown quantities consisting of variables X₁₁, Y_u, Z_u, and the double phase difference Njh_bu of the integer ambiguity.

Since the vector of observables y is non-linear relative to the variables X_u, Y_u, Z_n, the items in the formulas (7) and (8) are partially differentiated (linearized) relative to X_u, Y_u, Z_u. Sirtiilarly, in the embodiment, the optimum solution x of x and the residual vector v are given as follows .

In step SIlO, a test statistic T_1n for testing the overall validity of observables is derived and model errors are checked using the test statistic T_1n . The test statistic T₀₁ , which is the sum of squares of the residuals, is given as follows.

Here, Q_v denotes an error matrix (i.e., a variance-covariance matrix of residuals) . m denotes the number of unknown quantities. c_ip (i =l,....,<gr, ) is a vector for identifying an expected error factor e_t ( /=1,....,#, ) that could contribute to the test statistic T_m .

In the case where the observables consist of the pseudo ranges based on the C/A code, the P code in the Ll wave, the P code in the L2 wave, and the carrier phase accumulation values based on the Ll wave and the L2 wave, and the geometry-free model is used, the variable vector x for a minimum unit (i.e., one pair of the GPS satellites 10j and 1Oh ) can be represented as follows.

Then, residual vector v is given as follows.

he residual vector v used for calculating is modified by omitting the items associated

with integer ambiguity as follows. e

Here, represents a vector

for identifying an error in Φ_u, that is, error in the carrier phase data of the Ll wave, such as Φ (t) for the

Ll wave . epresents a vector for

identifying an error in Φ^ , that is, error in the carrier phase data of the L2 wave, such as Φ (t) for the L2 wave.

represents a vector for identifying

errors that occur simultaneously in Φ_Ll and Φ_L2.

represents a vector for

identifying an error proportional to the ratio between the frequency Ll and the frequency L2, such as an ionospheric

error. represents a vector for

identifying an error proportional to a ratio 9/7.

represents a vector for identifying an

error in code data (pseudo ranges based on the C/A code) .

In step S115, a test is performed on the test statistic T_m . In one embodiment, the test for testing the overall validity of the hypothesis is performed as follows .

Here, χ² _a(b,0) is the upper a probability point of the central % ²-distribution with b degrees of freedom.

If T_m≥χl(b,0) , it is concluded that unspecified model errors have occurred, and then the test process of Fig. 5 proceeds to step 120. Otherwise it is concluded that no unspecified model errors have occurred and therefore the test process finishes.

In step S120, p indicating the number of repetitions is incremented as follows. p = p +l.

In step S130, a test statistic t_lp is derived for each of the expected error factors by applying each corresponding vector C₁ , and the maximum test statistic

If₁-J is selected as follows.

In step S140, a test is performed on the maximum

test statistic W_j A ■ In one embodiment, the test is

performed using a normal distribution with variance 0 and mean 1 as follows.

est

process of Fig. 5 finishes, concluding that all errors have been found. In this case, the vector c_Jp related to

the maximum test statistic |f_Λj_ selected in step 130 is-

considered to indicate the error factor. In other words, it is concluded that the vector c_Jp related to the maximum

test statistic k- J indicates the error factor (error

type) that has contributed to a deterioration of the test statistic T_1n . Otherwise the test process proceeds to step 150.

In step S150, the redundancy is checked, as follows .

If m≤p, the test process of Fig. 5 finishes, concluding that there have been too many errors in the observation data to be identified. Otherwise the test process proceeds to step 160. In step S160, in order to eliminate or reduce

the effect of the maximum test statistic |f,J on the test

statistic T_n , the test statistic T_1n is modified as follows .

In step S170, the test is performed again on the modified test statistic T^* _m . In one embodiment, the test is performed as follows.

Here, F_a{m,∞,ϋ) is the upperα probability point of the central F-distribution with m, ∞ degrees of freedom. If T'm≥T Thr, it is concluded that there may be error factors other than the error factor specified by the vector C_{j p}, and the test process of Fig. 5 proceeds to step 180. Otherwise it becomes clear that the error factor specified by the vector c_Jp is the last error factor to be found, and thus the test process finishes, concluding that all errors have been found.

In step S180, in order to search for the error factors in the different direction, the vector c_ip is orthogonalized to be used for the next routine, as follows;

Then, the test process of Fig. 5 returns to step

120 with the orthogonalized vector c_ip , and continues recursively until it is concluded that all errors have been found or there have been too many errors.

According to the embodiment, it becomes possible to evaluate the quality of the observation data as a whole as well as each of the error factors separately.

It is noted that if the four GPS satellites can be tracked, the aforementioned tests may be performed on the observation data for three pairs of the GPS satellites (if one of GPS satellites is regarded as a reference satellite) . In this case, it is also possible to perform the aforementioned tests simultaneously for 3 pairs of the GPS satellites by integrating the corresponding vectors and the matrixes. In this embodiment, it is also advantageous to modify the float solution and/or the variance according to the test result. The modification is implemented by using matrixes Δ , Q_A , and K . The matrixes Δ , β_Δ , and K are given as follows.

Here, c_J ^T _p is the vector that maximizes the test

statistic W_j A (see step 130 in Fig. 5) . Q_x is a variance-

covariance matrix of the estimated x, and is represented as follows.

If the error factor specified by cj_p is outliers

in the observation data (typically, in the code data) , the optimum solution (float solution) Jc and the variance- covariance matrix of the observables Q_v are modified by being derived as follows.

On the other hand, if the error factor specified

by c^_p is slip in the carrier phase data (i.e. cycle-slip),

x and Q_v are modified by being derived as follows .

Here, [Δ] represents rounded-off components of the matrix Δ .

Next, referring to Fig. 6, a method of giving reliability representative of the quality of the calculated position of the mobile station 30 to the output result is described. The process routine shown in Fig. 6 is performed at each epoch independently. The process routine shown in Fig. 6 is executed by the calculation unit 40.

In step 200, two flags (flag 1 and flag 2) are provided. Both flags have initial values ^M0" indicating that the quality of the calculated position of the mobile station 30 meets predetermined criteria.

In step 210, the data to be used to calculate the position of the mobile station 30 are received (input) The input data may include the observation data, such as carrier phase accumulation values and pseudo ranges, and the dynamic information indicating dynamic behavior of the mobile station 30 derived based on the output signals from the various sensors 50. In step 220, the test is performed on the input data. The test is performed according to the method described in Fig. 5.

In step 220, if the test result doesn't meet predetermined criteria, the flag 1 is set to Λ^Λ1" indicating that there may be an error in the input data that could lead to an unsuccessful ambiguity resolution. For example, if in the step 150 (see Fig. 5) it is concluded that there have been too many errors in the observation data, the flag 1 is set to "1". Otherwise, the process proceeds to step 230.

In step 230, the ambiguity resolution is performed. The ambiguity resolution process includes the process to derive the float solution and the process to estimate the integer solution (also referred to as ^Λfixed solution' ) based on the float solution. In the process to derive the .float solution, a least-squares method or a modified Kalman filter (described later) or the like can be used. In the process to estimate the integer ambiguity, the LA-MBDA, for example, can be used, which de-correlates the integer ambiguities, and narrows the searching space of the integer solutions so as to facilitate finding the integer solution.

In step 240, an output test is performed on the derived float solution and/or the estimated integer ambiguity. The output test may be performed using conditional variance, a ratio test, or the like. The probability of estimating the correct integers is given as follows .

Here, σ_v denotes diagonal components (variance components) of the error matrix Q_v .

The ratio test is given as follows.

Here, Q_a, denotes an appropriate weight matrix, a denotes a subset of the float -solution of the variables x, a_s denotes the integer solution of the variables x with the highest probability and 5_S2 denotes the integer solution with the second highest probability. The optimum solution Jc or the modified solution Jc^* can be used as a. It is noted that other ways of conducting the output test are also applicable. For example, the output test using the innovation or residual vector v and/or the covariance matrix Q_x is applicable.

In step 240, if the test result of the output test doesn't meet the predetermined criteria, the flag 2 is set to "1" indicating that there may be an error in the estimated integer ambiguity that could lead to incorrect position determination of the mobile station 30. For

example, the flag 2 is set to

W η //

J- .

After the integer ambiguity is resolved, the position (coordinates) of the mobile station 30 can be determined using various differential positioning methods, for example.

In step 250, the calculated position of the mobile station 30 together with the flags are output as a result. The position of the mobile station 30 obtained in this way may be used in various controls or be presented as information, for example, it may be output and displayed on a screen of a navigation device, or displayed in a map shown on the screen of a mobile phone. In this way, it becomes possible for a user or analyst to check the quality of the calculated position by checking the value of the flags. For example, if the flag 1 is set to "0" and the flag 2 is set to "1", user or analyst may understand there is a possibility of errors in the calculated position as well as that such error would result from the ambiguity resolution process. It is noted that more flags can be used to indicate the type of error that has been found by the vector c_ip , for example.

Next, referring to Fig. 7, one embodiment of a method of calculating the position of the mobile station 30 at the epoch in which the flag 2 and the flag 1 are set to "1" is described. The process routine shown in Fig. 7 is executed by the calculation unit 40. The process routine shown in Fig. 7 may be performed before the step 250 in Fig. 6 in the case where the flag 2 and the flag 1 are set to "1". It is noted that the process routine shown in Fig. 7 may be performed in the case where either the flag 2 or the flag 1 is set to "1".

In step 300, Ar—Ϊ < 0 is checked. Here, k is the current epoch number and J = O₅I,...., (initial value of i is 0) . If k~i>0, the process proceeds to step 310. Otherwise the process finishes, setting the flag 1 to "k-i(=0)".

This indicates that the process routine has ended without any modification.

In step 310, the data having been obtained at the epoch k-i is read out from the memory. In other words, the previous epoch data are read out from the memory.

In step 320, the ambiguity resolution is performed based on the data obtained from epoch k-i to epoch k. In this case, the ambiguity resolution may be performed by the observation data obtained over more than 2 epochs. Therefore, the recursive Kalman filter can be applied to the ambiguity resolution. If applying the Kalman filter, the following equations can be obtained. For updating of the epoch,

For updating of the observables,

Here, Q and^' R represent the covariance matrix of the external noise and the covariance matrix of the observation noise, respectively. Ki denotes Kalman gain (corrective gain). P(k)^(+> and P(k)^M are the covariance matrix of the expected errors, and the covariance matrix of the estimated errors, respectively. The formulae (9) and (12) are filter equations, and the formulae (10) and (13) are covariance equations. Here, the superscript ^("' and ^!+! indicate time before and after the updating, respectively. In the static model, the item U(k-l) derived from the dynamic information in formula (9) does not exist. In the case of the modified Kalman filter for the single epoch application, the P(k) is initialized on epoch basis and thus the covariance is not took^* over for^{^}the next epoch. It is noted that the ambiguity resolution may be performed by the least-square method using observation data obtained at more than two epochs .

In step 330, the statistical test is performed on the float solution and/or the integer ambiguity which are derived or estimated based on the observation data obtained over more than 2 epochs. The method of the statistical test may be the same as described with reference to step 240' in Fig. 6. The test using the innovation (y(k)- W(k) *x(k) ^l~}) or residual vector v is also applicable. Further, the input test as described with reference to Fig. 5 may be performed on the observation data obtained over more than 2 epochs. If the result of the statistical test doesn't meet the predetermined criteria, the process returns to step 300 and continues recursively until the result of the statistical test meets the predetermined criteria, or k-i≤Q . In step 340, the calculated position of the mobile station 30 together with the flags are output as a result. In this way, it becomes possible for the user or analyst to check the quality of the calculated position by checking the value of the flags. According to the embodiment, while the ambiguity resolution may be performed on epoch basis, the ambiguity resolution on a multi-epoch basis is performed in exceptional cases where the flag 2 and/or the flag 1 is set to "1". Therefore, it is possible to output the calculated position of the mobile station 30 with high accuracy without degrading the advantage of the single epoch application.

It is noted that it is also possible to use more than two GPS antennas 32a, which are preferably provided at different places on the mobile station 30, in order to have more redundancy to ensure reliability.

Fig. 8A graphs the positions of the mobile station 30 at the plural epochs determined according to the conventional single epoch application. Fig. 8B graphs the positions at the plural epochs determined according to the embodiment shown in Fig. 7 of the present invention. The same observation data were used for both applications. In the conventional single epoch application, at 2 epochs, there are hops in the determined position with respect to the neighboring epochs, as shown in Fig. 8A. To the contrary, according to the embodiment shown in Fig, 7 of the present invention, there are no such hops in the determined position, as shown in Fig. 8B. The hop resulting from the epoch kl on the left side in Fig. 8A is avoided by using the observation data obtained over two epochs (epochs kl and kl-1) . The hop resulting from the epoch k2 on the right side in Fig. 8A is avoided by using the observation data obtained over 6 epochs

(epochs k2 and kl-6) . It is noted that in Fig. 8B the positions of the mobile station 30 at other epochs were determined by the observation data obtained at the single corresponding epoch. Next, referring to Fig. 9, a method of evaluating the quality of the calculated position of the mobile station 30 according to another embodiment of the present invention is described. The process routine shown in Fig. 9 is performed at each epoch independently. The process routine shown in Fig. 9 is executed by the calculation unit 40.

In step 400', two flags (flag 3 and flag 4) are provided. The flag 3 indicates reliability of the derived float solutions derived. The flag 4 indicates reliability of the derived fixed solutions. Both flags have initial values "0" indicating that the reliability is relatively high.

In step 410, the float solutions are derived based on the observation data. The float solutions can be derived using the least-squares method or the like.

In step 420, the innovation or residual vector v and/or the covariance matrix Q_x are derived using the derived float solutions. In step 430, the statistical test is performed on the derived float solutions by checking the innovation or residual vector v and/or the covariance matrix Q_x . The statistical test according to the Fig. 5 may be performed by using the residual vector v . The statistical test may be performed by checking the conditional variance.

If the test result doesn't meet the predetermined criteria, the flag 3 is set to "1" indicating that there may be an error in the derived float solutions. In this case, the derived float solutions are not reliable enough to be used for the integer estimating process (step 440) . Thus, the process at this epoch finishes by skipping the integer estimating process. Since the integer estimating process requires a relatively heavy calculation work load, it is useful to check the quality before the integer estimating process. It is noted that it is also possible to re-derive the float solutions using the previous epoch data, as is described referring to Fig. 7. On the other hand, if the test result meets the predetermined criteria, the process proceeds to step 440.

In step 440, the fixed solutions are estimated from the derived float solutions. The fixed solutions can be estimated using the LA-MBDA or the like.

In step 450, the innovation or residual vector v and/or the covariance matrix Q_x is derived using the estimated fixed solutions.

In step 460, the statistical test is performed on the estimated fixed solutions by checking innovation or residual vector v and/or the covariance matrix Q_x . The statistical test according to the Fig. 5 may be performed by using the residual vector v . The statistical test may be performed by checking the conditional variance and/or by using the ratio test as described previously.

If the test result doesn't meet the predetermined criteria, the flag 4 is set to "1" indicating that there may be an error in the estimated fixed solutions, and the process proceeds to step 470. If the test result meets the predetermined criteria, the process proceeds to step 470 without changing the values of the flags. In step 470, the calculated position of the mobile station 30 together with the flags are output as a result. In this way, it becomes possible for a user or analyst to check the quality of the calculated position by checking the value of the flags. Further, if the flag 3 is "1", the user or analyst can understand that the integer estimation is skipped because of the lack of reliability in the solutions.

According to the embodiment, since the reliability of the estimated fixed solutions is checked in addition to the reliability of the derived float solutions, it becomes possible to detect the errors that could occur in the integer estimating process.

Next, referring to Fig. 10, a method of improving the accuracy in the calculated position of the mobile station 30 based on the test results according to an embodiment of the present invention is described. The process routine shown in Fig. 10 is performed at each epoch independently. The process routine shown in Fig. 10 is executed by the calculation unit 40. In step 500, the ambiguity resolution is performed using the single difference (SD) observation data (see the formulas (5) and (6) ) .

In step 510, the innovation or residual vector v and/or the covariance matrix Q_x are derived using the fixed solutions and/or the float solution derived at step 500, as is described referring to Fig. 9.

In step 520, the ambiguity resolution is performed using the double difference (DD) observation data (see the formulas (3) and (4)).

In step 530, the innovation or residual vector v and/or the covariance matrix Q_x are derived using the fixed solutions and/or the float solution derived at step 520, as is described referring to Fig. 9.

In step 540, the reliability of the estimated fixed solutions derived at step 500 is evaluated by checking innovation or residual vector v and/or the covariance matrix Q₁ derived at 510. Similarly, the reliability of the estimated fixed solutions derived at step 520 is evaluated by checking innovation or residual vector v and/or the covariance matrix Q_x derived at step 530. Then, the reliabilities for both of the ambiguity resolutions are compared. In step 550, the more reliable estimated fixed solutions are selected to calculate the position of the mobile station 30. Then, the calculated position of the mobile station 30 is output as a result. In this way, it becomes possible to have more redundancy to ensure reliability and to appropriately select a reliable output by using the test result.

In this embodiment, it is also possible to derive the conditional variance or the ratio used in the ratio test for each of the ambiguity resolutions and then compare them so as to output the more reliable estimated fixed solutions.

In this embodiment, it is also possible to compare the reliability of the float solutions based on single difference observation data with the reliability of the float solutions based on double difference observation data before the integer estimating process. In this case, the more reliable float solutions are used to estimate the fixed solutions. Since the integer estimating process requires a relatively heavy calculation work load, it is useful to compare the reliability levels of the float solutions before the integer estimating process. It is noted that it is also possible to calculate the attitude or the azimuth angle of the mobile station 30 based on the observation data, as is the case with the position of the mobile station 30. The above- described embodiments can be applicable in such an application. The present invention is disclosed with reference to the preferred embodiments. However, it should be understood that the present invention is not limited to the above-described embodiments, and variations and modifications may be made without departing from the scope of the present invention.

For example, in the above embodiments, the influence of the ionospheric layer refraction effect, tropospheric bending effect, and multi-path are not considered^" in the model, but it^* is also possible to modόl these error factors.

In the above embodiment, a vehicle is illustrated as an example of the mobile station 30. The mobile station 30 may also include a folk lift or a robot with the receiver 32 and the calculation unit 40, for example.

Claims

1. A device for calculating a position of a mobile station using an integer ambiguity, comprising; a testing unit configured to evaluate reliability of observables on epoch basis, said observables including at least code data and carrier phase obtained at a single epoch, wherein said testing unit includes a first testing unit configured to test residuals as a whole, and a second testing unit configured to test one or more items of the residuals separately.

2. The device as claimed in claim 1, wherein the first testing unit is configured to test whether the sum of squares of the residuals is larger than a reference value, and the second testing unit is configured to perform the test if the sum of squares of the residuals is larger than the reference value.

3. The device as claimed in claim 1, wherein the residuals are derived using the integer ambiguity which is estimated based on the^* observables.

4. The device as claimed in claim 1, further comprising a third testing unit configured to evaluate reliability of the integer ambiguity estimated based on the observables .

5. The device as claimed in claim 1, further comprising a reliability outputting unit configured to output a reliability of the calculated position of the mobile station, said reliability being derived based on the test result.

6. The device as claimed in claim 1, wherein items of the observables to be used for estimating the integer ambiguity are selected according to the test result .

7. The device as claimed in claim 1, wherein the integer ambiguity is estimated using the observables obtained at a single epoch.

8. The device as claimed in claim 1, wherein the integer ambiguity is estimated after modifying a float solution according to the test result.

9. The device as claimed in claim 7, wherein if the test result at a certain epoch does not meet predetermined criteria, the integer ambiguity at the epoch is estimated using an item of observables obtained at another epoch.

10. A device for calculating a position of a mobile station using an integer^' ambiguity, comprising; a float solution testing unit configured to evaluate reliability of a float solution of the integer ambiguity derived based on observables which are obtained at a single epoch; and an integer solution testing unit configured to evaluate reliability of an integer solution of the integer ambiguity estimated based the float solution.

11. The device as claimed in claim 10, wherein if the test result on the float solution doesn' t^' meet predetermined criteria, the estimation of the integer solution based on said float solution is ceased.

12. The device as claimed in claim 10, wherein the tests are performed separately on the float solution derived based on double difference observables and on the float solution derived based on single difference observables .

13. The device as claimed in claim 12, wherein float solutions to be used for estimating the integer ambiguity are selected according to the test result.

14. The device as claimed in claim 10, wherein the tests are performed separately on the integer solution derived based on double difference observables and on the integer solution derived based on single difference observables .

15. The device as claimed in claim 14, wherein integer solutions to be used for calculating the position of a mobile station are selected according to the test result .

lβ. A method for calculating a position of a mobile station using an integer ambiguity, comprising; a testing step for evaluating reliability of observables including at least code data and carrier phase data, said observables being obtained at a single epoch, wherein said testing step includes a first testing step in which residuals as a whole are evaluated, and a second testing step in which a statistical test for detecting at least one of error factors is performed.

17. A method for calculating a position of a mobile station using an integer ambiguity, comprising; a float solution testing step of evaluating reliability of a float solution of the integer ambiguity derived based on observables which are obtained at a single epoch; and an integer solution testing step of evaluating reliability of an integer solution of the integer ambiguity estimated based the float solution.

18. A computer program for calculating a position of a mobile station using an integer ambiguity, said computer program causing a computer to execute the following step; a testing step for evaluating reliability of observables including at least code data and carrier phase data, said observables being obtained at a single epoch, wherein said testing step includes a first testing step in which residuals of variables with respect to the observables^" as a whole are evaluated, an3 a second testing step in which a statistical test for detecting at least one of error factors is performed.

19. A computer program for calculating a position of a mobile station using an integer ambiguity, said computer program causing a computer to execute the following steps; a float solution testing step of evaluating reliability of a float solution of the integer ambiguity derived based on observables obtained at a single epoch; and an integer solution testing step of evaluating reliability of an integer solution of the integer ambiguity estimated based the float solution.