JP7595435B2 - Positioning device - Google Patents

Description

This disclosure relates to a positioning device that is mounted on a mobile body and determines the position of a radio wave source that is a communication device that performs satellite communication.

To determine the position (location) of a communication device (radio wave source) that communicates via satellite, a method is known that uses the time difference of arrival (TDOA) and frequency difference of arrival (FDOA) of two received signals that are relayed by two different satellites from the communication device to be located and received at a ground station (see, for example, non-patent document 1).

Satellites used for communications include geostationary satellites, quasi-zenith orbit satellites, and low-earth orbit satellites. When two satellites are both geostationary, the frequency shift due to the Doppler effect caused by the movement of the radio wave source, among the frequency difference between the two received signals relayed by the two satellites, may be greater than the frequency shift due to the movement of the satellites. In such cases, conventional positioning devices may not be able to determine iso-FDOA lines that represent the position of the radio wave source that will give the detected frequency difference, and may not be able to locate the position. Even if iso-FDOA lines can be determined, there may be a large error in the iso-FDOA lines, resulting in a significant deterioration in the accuracy of the location position.

Conventional positioning devices do not display the location of the radio wave source if they are unable to locate the location. If they are able to locate the location, they will display it as is even if the accuracy has deteriorated.

D. P. Haworth, et. al, “Interference Localization For Eutelsat Satellites -The First European Transmitter Location System,” International Journal of Satellite Communications, Vol. 15, pp. 155-183, 1997

Conventional positioning devices have a problem in that it is difficult to properly present information regarding the position of a radio wave source when both relay satellites are geostationary satellites.

The position location device according to the present disclosure includes a correlation function calculation unit that receives a first signal and a second signal expressed by complex numbers, and a time difference and a frequency difference between the first signal and the second signal, takes a complex conjugate of the second signal having a time difference with respect to the first signal, multiplies the first signal by a complex number having an absolute value of 1 that rotates by the frequency difference, and integrates the multiplied value by the first signal for a determined time to calculate a two-dimensional correlation value between the first signal and the second signal, and a correlation function calculation unit that calculates the two-dimensional correlation value by changing the time difference and the frequency difference, and calculates a peak time difference and a frequency difference that are the time difference at which the two-dimensional correlation value reaches its maximum peak. a peak detection unit that detects a peak frequency difference, which is the frequency difference at which the two-dimensional correlation value has a maximum peak, and an arrival time, which is the difference between the arrival time of the radio wave from the radio wave source to the first satellite and the arrival time of the radio wave from the radio wave source to the second satellite, from a peak time difference of the two-dimensional correlation value calculated using a first reception signal, in which a radio wave radiated from a radio wave source whose position is unknown is relayed by a first satellite and received at a first observation station, as a first signal, and a second reception signal, in which a radio wave radiated from the radio wave source is relayed by a second satellite different from the first satellite and received at a second observation station, as a second signal, using the positions of the first and second satellites. an arrival frequency difference calculation unit that calculates an arrival time difference, using the positions and velocity vectors of the first and second satellites to calculate an arrival frequency difference, which is a difference between a frequency shift occurring from the radio wave source to the first satellite and a frequency shift occurring from the radio wave source to the second satellite, from a peak frequency difference of a two-dimensional correlation value calculated with the first reception signal as the first signal and the second reception signal as the second signal; an equal arrival time difference line calculation unit that calculates an equal arrival time difference line, which is a line representing the position of the radio wave source from which the arrival time difference is obtained, using the positions of the first and second satellites; a geostationary satellite discrimination unit which discriminates whether the first satellite and the second satellite are geostationary satellites ; a radio wave source speed estimation unit which estimates a radio wave source speed, which is a speed vector representing the movement of the radio wave source, from the arrival frequency difference when the geostationary satellite discrimination unit discriminates that both the first satellite and the second satellite are geostationary satellites; and a second positioning unit which generates radio wave source position estimation information relating to the position of the radio wave source using the radio wave source speed, the equal arrival time difference line, and the positions and speed vectors of the first satellite and the second satellite when the geostationary satellite discrimination unit discriminates that both the first satellite and the second satellite are geostationary satellites.

According to the present disclosure, when both relay satellites are geostationary satellites, information regarding the location of the radio wave source can be presented more appropriately than in the past.

1 is a conceptual diagram illustrating a process of locating the position of a radio wave source whose position is unknown using a conventional position locating device. 3 is a conceptual diagram illustrating a process of locating the position of a radio wave source whose position is unknown, in the position locating device according to the first embodiment. FIG. 1 is a block diagram showing a configuration of a position locating device according to a first embodiment; 2 is a schematic diagram illustrating an XYZ coordinate system used in a position locating device according to the present disclosure. FIG. 4 is a flowchart illustrating an operation of the position locating device according to the first embodiment. 3A to 3C are diagrams illustrating an example of radio wave source position estimation information output by the position locating device according to the first embodiment. 11 is a conceptual diagram illustrating a process of locating the position of a radio wave source whose position is unknown, in the position locating device according to the second embodiment. FIG. FIG. 11 is a block diagram showing a configuration of a position locating device according to a second embodiment. 10 is a flowchart illustrating an operation of the position locating device according to the second embodiment. FIG. 11 is a block diagram showing a configuration of a position locating device according to a third embodiment. 13 is a flowchart illustrating an operation of the position locating device according to the third embodiment. 13 is a diagram for explaining an example of radio wave source position estimation information output by a position locating device according to embodiment 3. FIG. 13 is a diagram for explaining another example of the radio wave source position estimation information output by the position locating device according to the third embodiment. FIG. 13 is a graph showing the presence probability in two examples of radio wave source position estimation information in a position locating device according to embodiment 3. FIG. 13 is a block diagram showing a configuration of a position locating device according to a fourth embodiment. 13 is a diagram for explaining, by way of example, the relationship between the possible direction range, the assumed speed range, and the allowable direction range used by the position locating device according to the fourth embodiment. FIG. 13 is a flowchart illustrating an operation of the position locating device according to the fourth embodiment.

Embodiment 1.
FIG. 1 is a conceptual diagram for explaining a process of locating the position of a radio wave source whose position is unknown by a conventional position locating device. FIG. 1 is based on a method disclosed in Non-Patent Document 1. In the example of FIG. 1, radio waves are transmitted from a radio wave source 70 whose position is unknown. A main lobe radio wave 81 emitted by the radio wave source 70 is transmitted toward a satellite (#1) 71. At the same time, the radio wave source 70 radiates a side lobe radio wave 82. The satellite (#1) 71 relays the radio wave 81 and transmits a radio wave 83. The radio wave 83 from the satellite (#1) 71 is received by an antenna (#1) 72 of a monitoring station. The antenna (#1) 72 receives the radio wave 83 and generates a reception signal 91. The side lobe radio wave 82 emitted by the radio wave source 70 is relayed by a satellite (#2) 73 and is received as a radio wave 84 by an antenna (#2) 74 of a monitoring station. The antenna (#2) 74 receives the radio wave 84 and generates a reception signal 92.

Antenna (#1) 72 and antenna (#2) 74 may be installed in different monitoring stations. The monitoring station where antenna (#1) 72 is installed is the first monitoring station, and the monitoring station where antenna (#2) 74 is installed is the second monitoring station. Received signal 91 is the first received signal relayed by satellite (#1) 71, which is the first satellite, and received at the first observation station. Received signal 92 is the second received signal relayed by satellite (#2) 73, which is the second satellite, and received at the second observation station. When antenna (#1) 72 and antenna (#2) 74 are installed in the same monitoring station, that monitoring station is both the first monitoring station and the second monitoring station.

The received signals 91 and 92 are input to the positioning device 50X. The positioning device 50X is a positioning device that performs positioning using a conventional method. The positioning device 50X estimates the time difference of arrival (TDOA) and frequency difference of arrival (FDOA) by correlating the received signals 91 and 92 relayed and received from satellite (#1) 71 and satellite (#2) 73, respectively. The positioning device 50X estimates the intersection of an iso-TDOA line 61, which is the intersection line between a plane with a constant estimated TDOA and the Earth's surface, and an iso-FDOA line 62, which is the intersection line between a plane with a constant estimated FDOA and the Earth's surface, as the position of the radio wave source 70. The estimated position of the radio wave source 70, i.e., the intersection line between the iso-TDOA line 61 and the iso-FDOA line 62, is called the location position 63.

In the path of radio waves 81 emitted from an unknown radio wave source 70 until they are received as radio waves 83 by antenna (#1) 72, and the path of radio waves 82 until they are received as radio waves 84 by antenna (#2) 74, the path after satellite (#1) 71 and the path after satellite (#2) 73 are known. The time difference and frequency difference obtained by correlating the received signals are subtracted from the time difference and frequency difference that occur on the path after satellite (#1) 71 and after satellite (#2) 73 to calculate the time difference of arrival (TDOA) and frequency difference of arrival (FDOA). TDOA is the difference between the time when radio waves 81 arrive at satellite (#1) 71 and the time when radio waves 82 arrive at satellite (#2) 73. FDOA is the difference between the frequency shift that occurs before radio waves 81 are received by satellite (#1) 71 and the frequency shift that occurs before radio waves 82 are received by satellite (#2) 73.

First, we will explain the conventional method for locating the position of a radio source 70 from TDOA and FDOA.

Define the following variables:
r: Position of the radio wave source 70. Expressed in a three-dimensional coordinate system with the center of the earth as the origin.
r ₁ : Position of satellite (#1) 71.
_r2 : Position of satellite (#2) 73.
r _m1 : Position of antenna (#1) 72.
r _m2 : Position of antenna (#2) 74.
c: speed of light.
v ₁ : Velocity vector of satellite (#1) 71.
v ₂ : Velocity vector of satellite (#2) 73.
f: frequency of the carrier wave of the radio wave from the radio wave source 70.
f ₁ : The difference between the uplink frequency and the downlink frequency in the satellite (#1) 71 .
f ₂ : The difference between the uplink frequency and the downlink frequency in satellite (#2) 73.
Δt: time difference between radio waves relayed by satellite (#1) 71 and satellite (#2) 73.
Δf: frequency difference between the radio waves relayed by satellite (#1) 71 and satellite (#2) 73.
Δtm: time difference occurring on the path from satellite (#1) 71 or satellite (#2) 73.
Δfm: frequency difference occurring on the path from satellite (#1) 71 or satellite (#2) 73.
Δt r : TDOA occurring up to satellite (#1) 71 or satellite (#2) 73.
Δfr: FDOA occurring up to satellite (#1) 71 or satellite (#2) 73.
u ₁ : A unit vector from the radio wave source 70 to the satellite (# 1) 71 .
u ₂ : A unit vector from the radio wave source 70 to the satellite (# 2) 73 .
u _1m : A unit vector directed from the antenna (# 1) 72 to the satellite (# 1) 71 .
u _2m : A unit vector from antenna (# 2) 74 to satellite (# 2) 73 .
s ₁ (t): signal 91 received by antenna (# 1 ) 72 .
s ₂ (t): signal 92 received at antenna (#2) 74 .
R: The radius of the Earth. The Earth is assumed to be a perfect sphere.

In position location, a correlation function (Cross Ambiguity Function: CAF) is calculated using equation (1) between signal _s1 (t) received by antenna (#1) 72 and the conjugate _s2 ^* (t) of the signal received by antenna (#2) 74. Signal _s1 (t) is the first signal in calculating the correlation function, and signal _s2 (t) is the second signal. The correlation function CAF is calculated using the first received signal as the first signal and the second received signal as the second signal.

ここで、Ｔは信号長である。τは、信号ｓ_１(t)と信号ｓ_２(t)の間に設定する時間差である。νは、信号ｓ_１(t)と信号ｓ_２(t)の間に設定する周波数差である。ｓ_１(t)が第一の信号である。信号ｓ_２(t)が第２の信号である。式(1)で計算されるＣＡＦ(τ,ν)は、時間差τおよび周波数νが入力されて、ｓ_１(t)に対して時間差τを持たせ複素共役をとったｓ_２ ^*(t+τ)に、周波数νで回転する絶対値が１である複素数ｅ^{－ｊ２πνｔ}を乗算して得られる信号とｓ_１(t)とを乗算して得られる値を決められた時間Ｔだけ積分した、ｓ_１(t)とｓ_２(t)との間の２次元相関値である。

Here, T is the signal length. τ is the time difference set between signal s ₁ (t) and signal s ₂ (t). ν is the frequency difference set between signal s ₁ (t) and signal s ₂ (t). s ₁ (t) is the first signal. Signal s ₂ (t) is the second signal. CAF(τ,ν) calculated by formula (1) is a two-dimensional correlation value between s 1 (t) and s 2 (t) obtained by inputting a time difference τ and a frequency ν, multiplying s ₂ ^* (t+τ) obtained by giving a time difference τ to s ₁ (t) and taking a complex conjugate, by a complex number e ^−j2πνt having an absolute value of 1 that rotates at a frequency ν, and multiplying the signal obtained by this multiplication by s ₁ (t), and integrating the value obtained by this _{multiplication} by s ₁ (t) for a determined time T.

The time difference Δt is the time when the correlation function CAF calculated by formula (1) reaches a peak value. The frequency difference Δf is the time when the correlation function CAF reaches a peak value. Δt is generated by the difference in the path length between the signals s ₁ (t) and s ₂ (t), and can be expressed by the formula shown below.

The frequency difference Δf is the difference between the Doppler frequency shift generated at antenna (#1) 72 on satellite (#1) 71 and the Doppler frequency shift generated at satellite (#2) 73 and antenna (#2) 74, so it is calculated using the following formula. Conventionally, the frequency difference is calculated without taking into account the speed of radio wave source 70.

The unit vectors used in equation (3) are calculated using the formula below.

Of Δt, Δtm, which is the component due to the difference between the path taken by radio wave 81 after it reaches satellite (#1) 71 and the path taken by radio wave 82 after it reaches satellite (#2) 73, is calculated using formula (8) shown below. As shown in formula (9), Δtm is subtracted from Δt to calculate Δtr. The condition formula for the position r of the radio wave source 70 regarding Δtr is formula (10). The position locating device 50X calculates the position r that satisfies formula (10) for Δtr.

The frequency difference Δfm that occurs regardless of the position r and velocity vector v of the radio wave source 70 is calculated by the following formula (11). As shown in formula (12), Δf
The conditional expression for Δfr with respect to the position r of the radio wave source is given by expression (13). The position locating device 50X calculates the position r that satisfies expression (13) for Δfr.

Assuming that radio wave source 70 exists on the surface of the Earth, the following conditions are satisfied with respect to position r of radio wave source 70.

The conventional position locating device 50X finds an iso-TDOA line 61 that represents a position r that satisfies equations (10) and (14). It also finds an iso-FDOA line 62 that represents a position r that satisfies equations (13) and (14). The intersection of the iso-TDOA line 61 and the iso-FDOA line 62 is taken as the location position 63 of the radio wave source 70.

If the radio wave source 70 is not moving, or if both or either of the satellites (#1) 71 and (#2) 73 are not geostationary satellites, the position of the radio wave source 70 can be located with appropriate accuracy even with the conventional positioning device 50X. If the radio wave source 70 is not moving, the frequency difference can be calculated correctly using equation (13). If both or either of the satellites (#1) 71 and (#2) 73 are not geostationary satellites, the speed of the radio wave source 70 is small compared to the speed of the satellites, so the position r of the radio wave source 70 from which the frequency difference is obtained can be calculated with appropriate accuracy even using equation (13) that does not take the speed of the radio wave source 70 into account.

The speed at which a satellite that is not a geostationary satellite moves is thought to be more than 50 times faster than the speed at which radio wave source 70 moves. This is because a satellite that is not a geostationary satellite orbits the Earth in a period of a few hours or less, and moves at a speed of about 1 km per second or more relative to the Earth's surface. The upper limit of speed for radio wave source 70 is thought to be about 30 m per second. In the case of a ship, the upper limit of speed is said to be 30 knots per hour (55.6 km per hour, 15.4 m per second). In the case of an automobile, for example, the upper limit of speed is thought to be about 100 km per hour (27.8 m per second).

The speed at which a geostationary satellite moves relative to the Earth's surface is about 1 m per second. When the radio wave source 70 moves relatively fast, the speed at which the radio wave source 70 moves becomes greater than the speed at which the geostationary satellite moves. As a result, when both satellite (#1) 71 and satellite (#2) 73 are geostationary satellites, the frequency difference occurs mainly due to the movement of the radio wave source 70. Therefore, when both satellite (#1) 71 and satellite (#2) 73 are geostationary satellites and the radio wave source 70 is moving, it may not be possible to obtain the position of the radio wave source 70 that satisfies equation (13), which does not take into account the frequency difference that occurs due to the movement of the radio wave source 70. Even if the position of the radio wave source 70 that satisfies equation (13) is obtained, it may be a position that is far away from the actual position of the radio wave source 70.

The positioning device according to the present disclosure estimates the velocity vector along which the radio wave source 70 moves, and when both satellite (#1) 71 and satellite (#2) 73 are geostationary satellites and the radio wave source 70 is moving, presents the user with more appropriate positioning information than in the past. Figure 2 is a conceptual diagram illustrating the process of locating the position of a radio wave source whose position is unknown, using the positioning device according to the first embodiment. The biggest difference from the past is that the velocity vector v of the moving radio wave source 70 is estimated.

In the position locating device 50 according to the present disclosure, the frequency difference Δf also takes into account a frequency difference occurring due to the velocity vector v of the radio wave source 70. In order to explain the processing of the position locating device 50, the following variables are defined.
v: the velocity vector of the radio source 70. Also called the radio source velocity.
_FDOA (v,r): The frequency difference caused by the velocity vector v of the radio source 70. This is called the radio source frequency difference.
Δfc(r): the frequency difference occurring between satellite (#1) 71 and satellite (#2) 73.
u _d : A vector representing the difference between u ₁ and u _2. u _d = u ₁ - u ₂ .

The position locator 50 uses the equations shown below for the frequency differences Δf and Δfr. The frequency difference Δf has the components shown in equation (15). For Δf according to equation (15), the relational equation between the position r and velocity vector v of the radio wave source 70 and Δfr calculated using equation (12) is given by equation (16). The third term in equation (16) is the frequency difference that occurs when the movement of the radio wave source 70 is taken into consideration.

Once the velocity vector v of the radio wave source 70 can be estimated, the radio wave source frequency difference _еFDOA (v,r) is calculated using the following formula (17). _u1 (r) and _u2 (r) use the vectors at the position r of the radio wave source 70 used when estimating v. Furthermore, as shown in formula (18), Δfc(r) is calculated by subtracting _еFDOA (v,r) from Δfr. Δfc(r) is the frequency difference obtained by subtracting the radio wave source frequency difference _еFDOA (v,r) from the FDOA Δfr. Δfc(r) is an estimated value of the frequency difference that occurs due to the movement of the satellite. Δfc(r) is called the estimated FDOA. The line that connects the position r of the radio wave source 70 calculated using formula (19) with respect to Δfc(r) is called the estimated iso-FDOA line 64.

The configuration of the positioning device according to the first embodiment will be described with reference to FIG. 3. The positioning device 50 includes frequency conversion units 1A and 1B, AD conversion units 2A and 2B, quadrature detection units 3A and 3B, correlation function calculation unit 4, peak detection unit 5, satellite data acquisition unit 6, TDOA calculation unit 7, FDOA calculation unit 8, iso-TDOA line calculation unit 9, stationary satellite discrimination unit 10, first positioning unit 11, radio wave source speed estimation unit 12, second positioning unit 13, and memory unit 14. The first positioning unit 11 includes an iso-FDOA line calculation unit 15. The second positioning unit 13 includes an FDOA correction unit 16 and an estimated iso-FDOA line calculation unit 17.

Received signals 91 and 92 are input to the position locating device 50. The frequency conversion unit 1A converts the received signal 91 to an intermediate frequency. The AD conversion unit 2A converts the signal output by the frequency conversion unit 1A into a digital value with a determined number of bits at a determined time interval. The quadrature detection unit 3A performs quadrature detection on the signal digitized by the AD conversion unit 2A and outputs a complex baseband signal s ₁ (t). The frequency conversion unit 1B, the AD conversion unit 2B, and the quadrature detection unit 3B perform similar processing on the received signal 92. The quadrature detection unit 3B outputs a complex baseband signal s ₂ (t). The complex baseband signal is a signal expressed by a complex number.

The correlation function calculation unit 4 calculates the correlation function CAF by formula (1) for s ₁ (t) and s ₂ (t) in a determined period T. The peak detection unit 5 calculates the correlation function CAF by changing τ and ν, and detects Δt and Δf, which are τ and ν at which the correlation function CAF has the maximum peak. Δt is the peak time difference at which the correlation function CAF has the maximum peak. Δf is the peak frequency difference at which the correlation function CAF has the maximum peak. From the signals 91 and 92, the relay satellite (#1) 71 and satellite (#2) 73 are determined. The satellite data acquisition unit 6 acquires satellite data including the positions and velocity vectors of the relay satellite (#1) 71 and satellite (#2) 73 from the outside.

The TDOA calculation unit 7 uses the positions _r1 , _r2 , _rm1 , and _rm2 to calculate Δt from Δt using equation (9). Δtr is the time difference of arrival (TDOA). The FDOA calculation unit 8 uses the positions _r1 , _r2 , _rm1 , and _rm2 and the velocity vectors _v1 and _v2 to calculate Δfr from Δf using equation (12). Δfr is the frequency difference of arrival (FDOA). Note that although the positions _r1 , _r2 , _rm1 , and _rm2 are not included in equation (12), the unit vector _u1m is calculated from the positions _r1 and _rm1 using equation (6), and the unit vector _u2m is calculated from the positions _r2 and _rm2 using equation (7). Δfr is calculated indirectly using positions r ₁ , r ₂ , r _m1 , and r _m2 .

The iso-TDOA line calculation unit 9 uses the positions _r1 and _r2 to calculate an iso-TDOA line 61 that represents the position r that satisfies the formula (10) and the formula (14). The iso-FDOA line calculation unit 15 uses the positions _r1 and _r2 and the velocity vectors _v1 and _v2 to calculate an iso-FDOA line 62 that represents the position r that satisfies the formula (13) and the formula (14). The geostationary satellite determination unit 10 determines whether the satellite (#1) 71 is a geostationary satellite, and whether the satellite (#2) 73 is a geostationary satellite. If at least one of the satellite (#1) 71 and the satellite (#2) 73 is not a geostationary satellite, the first positioning unit 11 determines the intersection of the iso-TDOA line 61 and the iso-FDOA line 62 as the location position 63 (position r) of the radio wave source 70.

The iso-TDOA lines 61 are iso-time difference lines that represent the position r of the radio wave source 70 where the time difference of arrival (TDOA) is obtained. The iso-TDOA line calculation unit 9 is an iso-time difference line calculation unit that calculates the iso-time difference of arrival lines. The iso-FDOA lines 62 are iso-frequency difference lines that represent the position r of the radio wave source 70 where the frequency difference of arrival (FDOA) is obtained. The iso-FDOA line calculation unit 15 is an iso-frequency difference line calculation unit that calculates the frequency difference of arrival.

When both satellite (#1) 71 and satellite (#2) 73 are geostationary satellites, the radio wave source velocity estimation unit 12 estimates the velocity vector v of the radio wave source 70. The second positioning unit 13 generates an estimated iso-FDOA line 64 calculated using the iso-TDOA line 61 and the velocity vector v as radio wave source position estimation information 65.

The second positioning unit 13 has an FDOA correction unit 16 and an estimated isoFDOA line calculation unit 17. The FDOA correction unit 16 corrects the FDOA (= Δfr) from the estimated velocity vector v to estimate Δfc(r). A determined number of Δfc(r) are extracted from the estimated Δfc(r) so that they are approximately equally spaced. For the extracted Δfc(r), the estimated isoFDOA line calculation unit 17 calculates an estimated isoFDOA line 64, which is the position r of the radio wave source 70 that satisfies equations (14) and (19).

Δfc(r) is a satellite frequency difference, which is the difference between the frequency shift caused by the movement of satellite (#1) 71 and the frequency shift caused by the movement of satellite (#2) 73. The satellite frequency difference is also a frequency difference obtained by subtracting the radio wave source frequency difference _еFDOA (v,r) from Δfr. The FDOA correction unit 16 is a satellite frequency difference calculation unit that calculates the satellite frequency difference. The estimated iso-FDOA line 64 is an estimated iso-FDOA line that is the position r of the radio wave source 70 from which Δfc(r) is obtained. The estimated iso-FDOA line calculation unit 17 is an estimated iso-FDOA line calculation unit that calculates an estimated iso-FDOA line.

The storage unit 14 stores data necessary for processing such as positions _r1 , _r2 , _rm1 , _rm2 , velocity vectors _v1 , _v2 , and frequencies f, _f1 , _f2 . The storage unit 14 also stores data used during processing such as _u1 , _u2 , _u1m , _u2m , Δtm, and Δfm. The storage unit 14 also stores an assumed location range 31 and an assumed speed range 32. The assumed location range 31 is a range in which the position r of the radio wave source 70 is assumed to exist. The assumed speed range 32 is a range in which the speed vector v of the radio wave source 70 is assumed to exist. The radio wave source speed estimation unit 12 estimates the speed vector v of the radio wave source 70 using the assumed location range 31 and the assumed speed range 32. The radio wave source speed estimation unit 12 estimates v existing in the assumed speed range 32 at a position on the iso-TDOA line 61 included in the assumed location range 31.

A method for the radio wave source velocity estimation unit 12 to estimate the velocity vector v of the radio wave source 70 will be described. The following variables are defined.
r _TDOAk : Positions distributed on the iso-TDOA line 61 included in the assumed position range 31. The number of positions where n is arranged is k=1, 2, ..., n. It is desirable to arrange them evenly.
Q _k (r _TDOAk , v): the error function used in calculating the velocity vector v.
ε: A threshold value for the error function Q _k (r _TDOAk , v).
v _max : the upper limit of |v| defined in the assumed speed range 32.
v _min : The lower limit value of |v| defined in the assumed speed range 32.
(v _x , v _y , v _z ): A vector expressing the velocity vector v in the Cartesian coordinate system.
Δv: The step size for changing each component of the Cartesian coordinate system when searching for the velocity vector v.

The assumed location range 31 is set, for example, such that the latitude is 130 degrees to 140 degrees east and the longitude is 25 degrees to 45 degrees east. The assumed location range 31 may have a boundary where the latitude and longitude change simultaneously. The assumed speed range 32 is set, for example, by an upper limit v _max and a lower limit v _min for a scalar value |v|. The upper limit v _max is set, for example, to 15 meters per second. The lower limit v _min is set, for example, to 0 meters per second. The assumed speed range 32 may also restrict the direction of the speed vector v.

The XYZ coordinate system is used, which is an orthogonal coordinate system that rotates with the earth's rotation with the center of the earth as the origin. The XYZ coordinate system will be described with reference to FIG. 4. FIG. 4 is a schematic diagram illustrating the XYZ coordinate system used in the positioning device according to the present disclosure. Here, Rc is the radius of the geostationary orbit (Rc ≒ 35,786 km). In the XYZ coordinate system, the Z axis is set on the axis on which the earth 60 rotates. The north direction is the positive direction of the Z axis. The plane passing through the equator 77 is the XY plane. The Y axis is set so that satellite (#1) 71 and satellite (#2) 73 are arranged symmetrically with respect to the YZ plane. The side on which satellite (#1) 71 and satellite (#2) 73 exist is the positive direction of the Y axis. In other words, the longitude midway between satellite (#1) 71 and satellite (#2) 73 is zero degrees, and the meridian 78 passing midway between satellite (#1) 71 and satellite (#2) 73 exists on the YZ plane.

The X-axis is set so that it is perpendicular to the Y-axis and Z-axis. The positive direction of the X-axis is 90 degrees east of the positive direction of the Y-axis. Note that this calculation assumes that satellite (#1) 71 and satellite (#2) 73, which are geostationary satellites, are located on geostationary orbit 76. There are cases where geostationary satellites are located in positions that are slightly offset from the geostationary orbit, but the deviation from the geostationary orbit can be accommodated by rotating the XYZ coordinate system to correct for this.

The r _TDOAk is set, for example, at intervals of 0.5 degrees of latitude on the iso-TDOA lines 61 included in the assumed location range 31. The r _TDOAk may be set so that the decimal point of the latitude is 0 or 0.5, or may be set so that it is a value other than 0 or 0.5. An appropriate number of r _TDOAk are set at appropriate intervals according to the size of the assumed location range 31. The r _TDOAk is an assumed position set on the iso-TDOA lines 61 included in the assumed location range 31.

The error function Q _k (r _TDOAk , v) used when determining the velocity vector v of the radio wave source 70 at each r _TDOAk is defined by the following equation.

The radio wave source velocity estimation unit 12 determines the velocity vector v so that Q _k (r _TDOAk , v) is sufficiently small, that is, so that the following formula is satisfied: ε is determined so that the velocity vector v can be obtained with the required accuracy.

Determining the velocity vector v so that the error function _Qk ( _rTDOAk ,v) calculated by equation (20) satisfies equation (21) means determining v so that _еFDOA (v, _rTDOAk ), which is the frequency difference caused by v, satisfies equation (22) shown below. _еFDOA (v,r) is the radio wave source frequency difference, which is the difference between the frequency shift caused by the velocity vector v of the radio wave source 70 moving up to satellite (#1) 71 and the frequency shift caused by the velocity vector v up to satellite (#2) 73. The radio wave source velocity estimation unit 12 estimates the velocity vector v according to the assumed position _rTDOAk . The radio wave source velocity estimation unit 12 estimates the velocity vector v so that _еFDOA (v, _rTDOAk ) becomes the frequency difference obtained by subtracting Δfc( _rTDOAk ) from Δfr.

There may be cases where sufficient calculation time, memory, etc. cannot be secured due to hardware constraints of the device. As shown in equations (23) and (24), V may be calculated so that the average value of _Qk ( _rTDOAk , v) along the TDOA line is sufficiently small.

By using equations (23) and (24), the number of calculated velocity vectors is reduced, which in turn reduces the number of calculated FDOA lines and control points, making it possible to reduce the computational load and the load on the device's memory, etc.

The radio wave source velocity estimation unit 12 calculates the magnitude |v| of the velocity vector v using the following formula:

The expected velocity range 32 limits the magnitude |v| of the velocity vector v, as shown in equation (26) below.

The radio wave source velocity estimation unit 12 changes each of _vx , _vy , and _vz by Δv to calculate all v that satisfy equations (21) and (26). The magnitude of Δv is determined so that v can be calculated with the required accuracy and the calculation time is not excessively long.

The FDOA correction unit 16 calculates _еFDOA (v,r) using equation (17) from _u1 ( _rTDOAk ) and _u2 ( _rTDOAk ) at _rTDOAk used to estimate the velocity vector v and the estimated velocity vector v. Furthermore, it calculates Δfc(r) using equation (18). Δfc(r) is the corrected FDOA. The geostationary satellite is in a geostationary orbit at a distance of about 35,786 km from the center of the earth. Even if the position r of the radio wave source 70 changes by about several hundred km, the changes in _u1 (r) and _u2 (r) are small. The change in the difference vector _u2 (r) _-u1 (r) between _u1 (r) and _u2 (r) when r changes is smaller than the changes in _u1 (r) and _u2 (r). Therefore, _еFDOA (v,r) can be estimated even if the position vector r is not accurately determined.

From the estimated Δfc(r), a number of Δfc(r) according to the size of the assumed location range 31 is extracted. The number of Δfc(r) is determined so that iso-FDOA lines 64 are arranged at appropriate intervals within the assumed location range 31. Δfc(r) are extracted so that the positions r are close to equal intervals. For each extracted Δfc(r), the estimated iso-FDOA line calculation unit 17 calculates an estimated iso-FDOA line 64 which is the position r of the radio wave source 70 that satisfies equations (14) and (19). The second positioning unit 13 outputs the iso-TDOA lines 61 and the estimated iso-FDOA lines 64 within the assumed location range 31 as radio wave source position estimation information.

The operation will be described. Fig. 5 is a flow chart for explaining the operation of the position locating device according to the first embodiment. In step S01, the antenna (#1) 72 receives radio waves 83 relayed by the satellite (#1) 71 to generate a first reception signal _s1 (t). At the same time, the antenna (#2) 74 receives radio waves 84 relayed by the satellite (#2) 73 to generate a second reception signal _s2 (t). In step S02, the peak detection unit 5 detects the peak time difference Δt and peak frequency difference Δf at which the correlation function between _s1 (t) and _s2 (t) has a maximum peak. In step S03, the TDOA calculation unit 7 calculates Δtr, which is the TDOA. The FDOA calculation unit 8 calculates Δfr, which is the FDOA. In step S04, the iso-TDOA line calculation unit 9 calculates iso-TDOA lines 61 for Δtr.

In step S05, the geostationary satellite determination unit 10 determines whether satellite (#1) 71 and satellite (#2) 73 are geostationary satellites. In step S06, it is checked whether both satellite (#1) 71 and satellite (#2) 73 are geostationary satellites. In step S06, if at least one of satellite (#1) 71 and satellite (#2) 73 is not a geostationary satellite (NO in step S06), in step S07, the iso-FDOA line calculation unit 15 calculates iso-FDOA lines 62 for Δfr. In step S08, the first positioning unit 11 determines the intersection of iso-TDOA lines 61 and iso-FDOA lines 62 as the location position 63 (position r) of the radio wave source 70.

If, in S06, both satellite (#1) 71 and satellite (#2) 73 are geostationary satellites (YES in S06), in step S09, the radio wave source velocity estimation unit 12 estimates the velocity vector v of the radio wave source 70. In step S10, the FDOA correction unit 16 corrects Δfr taking into account the velocity vector v and calculates Δfc(r). In step S11, the estimated iso-FDOA line calculation unit 17 calculates the estimated iso-FDOA line 64 for Δfc(r). In step S12, the second positioning unit 13 outputs the iso-TDOA line 61 and the estimated iso-FDOA line 64 included within the assumed position range 31 as radio wave source position estimation information.

FIG. 6 shows an example of the result of position location performed by the position location device 50. In FIG.
The conditions for positioning are as follows:
Position of satellite (#1) 71: Geostationary orbit with longitude of 2 degrees.
Position of Satellite (#2) 73: Geostationary orbit at longitude -2 degrees.
Velocity vector v ₁ of satellite (#1) 71:
v ₁ = (v _x , v _y , v _z ) = (0.5, -0.3, -0.9) [m/sec].
Velocity vector _v2 of satellite (#2) 73:
v ₂ =(v _x , v _y , v _z )=(−0.3, −0.2, 0.9) [m/sec].
Δt r = -0.2 [milliseconds].
Δfr = 24 [Hz].
Expected location range 31: latitude 20 to 35 degrees, longitude -14 to -1 degree.
vmax: 15 [m/sec]

The estimated iso-FDOA lines 64 and iso-TDOA lines 61 within the assumed location range 31 are displayed as radio source position estimation information 65. The estimated iso-FDOA lines 64 are displayed so that the value of Δf c (r), which is the FDOA generated by the movement of satellite (#1) 71 and satellite (#2) 73, can be seen. Only four estimated iso-FDOA lines 64 are displayed to make the diagram easier to see. In reality, an appropriate number of estimated iso-FDOA lines 64 are displayed at appropriate intervals. The iso-TDOA lines 61 may be displayed only within the assumed location range 31. The assumed location range 31 may be displayed with transparent hatching instead of a frame line. When the estimated iso-FDOA lines 64 are displayed only within the assumed location range 31, the assumed location range 31 does not need to be displayed. The estimated iso-FDOA lines 64 may be displayed inside and outside the assumed location range 31, and the assumed location range 31 may also be displayed.

Even if both relaying satellites are geostationary satellites and position location is not possible, the positioning device 50 can estimate the velocity vector v of the radio wave source 70 and present radio wave source position estimation information 65 based on the velocity vector v. When both relaying satellites are geostationary satellites, information regarding the position of the radio wave source 70 can be presented more appropriately than in the past.

In actual positioning, it is unavoidable that there is a frequency shift in the receiver device of the monitoring station or the oscillator of the communication satellite, and that there is an offset (steady deviation) in the orbit (position, velocity) of the communication satellite. Therefore, there is a time or frequency offset in the received signals _s1 (t) and _s2 (t). If there is a time or frequency offset, the accuracy of positioning deteriorates.

The position locating device may not include the iso-FDOA line calculation unit 15, and the first positioning unit 11 may use the positions _r1 , _r2 and the velocity vectors _v1 , _v2 to determine a position on the iso-TODA line 61 that satisfies equation (13), and set it as the located position 63. A position that satisfies equation (13) is a position that satisfies FDOA. Regardless of whether the position locating device includes the iso-FDOA line calculation unit 15, the first positioning unit 11 determines a position on the iso-TODA line 61 that satisfies FDOA, and sets it as the located position 63.
The above also applies to other embodiments.

Embodiment 2.
A position locating device 50A according to the second embodiment is a modification of the first embodiment in which the position of a radio wave source is located using a third reception signal and a fourth reception signal that are received by relaying radio waves emitted from a reference station, which is a reference station whose position is known, by two satellites. The position locating device according to the second embodiment removes the time or frequency offset deviation present in the first reception signal and the second reception signal by subtracting Δt _ref and Δf _ref detected in the third reception signal and the fourth reception signal from the peak time difference Δt and the peak frequency difference Δf detected in the first reception signal and the second reception signal.

Figure 7 is a conceptual diagram for explaining the process of locating the position of a radio wave source whose position is unknown, using a position locating device according to the second embodiment. In Figure 7, radio waves are also emitted from a reference station 75, whose position is known. A main lobe radio wave 85 emitted by the reference station 75 is transmitted toward a satellite (#1) 71. At the same time, the reference station 75 radiates a side lobe radio wave 86. The satellite (#1) 71 relays the radio wave 85 and transmits radio wave 87. The radio wave 87 from the satellite (#1) 71 is received by an antenna (#1) 72 of the monitoring station. The antenna (#1) 72 receives the radio wave 87 and generates a reception signal 93. The reception signal 93 is a third reception signal. The side lobe radio wave 86 emitted by the reference station 75 is relayed by a satellite (#2) 73 and is received as radio wave 88 by an antenna (#2) 74 of the monitoring station. The antenna (#2) 74 receives the radio wave 88 and generates a reception signal 94. Received signal 94 is the fourth received signal.

When using a reference station 75, we will explain how to determine the position of the radio source 70 from the TDOA and frequency difference (FDOA).

Define the following variables:
r _c : the position of the reference station 75 .
Δt _ref : the time difference that occurs when the radio waves emitted by the reference station 75 are relayed by the satellite (# 1 ) 71 and the satellite (# 2 ) 73 .
Δf _ref : A frequency difference that occurs when the radio wave emitted by the reference station 75 is relayed by the satellite (# 1 ) 71 and the satellite (# 2 ) 73 .
Δt _tar-ref : The difference between Δt and Δt _ref .
Δf _tar-ref : The difference between Δf and Δf _ref .
Δts: the time difference that occurs between the radio waves emitted by the reference station 75 and the satellite (#1) 71 and the satellite (#2) 73.
Δfs: the frequency difference that occurs between the satellite (#1) 71 and the satellite (#2) 73 in the radio waves emitted by the reference station 75.
u _1c : A unit vector from the reference station 75 to the satellite (# 1) 71 .
u _2c : A unit vector from the reference station 75 to the satellite (# 2) 73 .
s ₃ (t): signal 93 received by antenna (# 1 ) 72 .
s ₄ (t): signal 94 received by antenna (#2) 74.

A correlation function CAF calculated in the same manner as in equation (1) is also calculated between a signal _s3 (t) radiated from a reference station 75 and received by antenna (#1) 72, and a conjugate _s4 ^* (t) of a signal received by antenna (#2) 74. The correlation function CAF is calculated by taking the third received signal _s3 (t) as the first signal and the fourth received signal _s4 (t) as the second signal.

Regarding the received signal of the radio wave emitted from the reference station 75, τ at which the correlation function calculated by equation (1) reaches a peak value is the time difference Δt _ref . Also, ν at which the correlation function reaches a peak is the frequency difference Δf _ref . Δt _ref occurs due to the difference in the path length between the signal s ₃ (t) and the signal s ₄ (t), and can be expressed by the equation shown below. Also, the time difference Δts that occurs from the reference station 75 to the satellite can be expressed by the equation shown below.

Δts is the reference station time difference, which is the difference between the arrival time of radio waves from the reference station 75 to the satellite (#1) 71 and the arrival time of radio waves from the reference station 75 to the satellite (#3) 73.

The frequency difference Δf _ref and the frequency difference Δf occurring from the reference station 75 to the satellite are calculated by the following formulas.

Δfs is the reference station frequency difference, which is the difference between the frequency shift that occurs from the reference station 75 to satellite (#1) 71 and the frequency shift that occurs from the reference station 75 to satellite (#2) 73.

The unit vectors u _1c (r) and u _2c (r) used in the equations (29) and (30) are calculated using the equations shown below.

In the second embodiment, Δt tar-ref is calculated as the difference between the peak time difference Δt between the received signals s ₁ (t) and s ₂ (t) generated by receiving the radio waves radiated from the radio wave source 70 and the peak time difference Δt ref between the received signals s ₃ (t) and s ₄ (t) generated by receiving the radio waves radiated from the reference station 75. Also, Δf _{tar-ref is} calculated as the difference between the peak frequency difference Δf between the received signals s ₁ (t) and s ₂ (t) and the peak _frequency difference Δt _ref between the received signals s ₃ (t) and s ₄ (t). The peak time difference Δt is called the radio wave source peak time difference, and the peak frequency difference Δf is called the radio wave source peak frequency difference. The peak time difference Δt _ref is called the reference station peak time difference, and the peak frequency difference Δf _ref is called the reference station peak frequency difference.

Δt _tar-ref has components as shown in equation (33) below.

Δf _tar-ref has components as shown in equation (34) below.

Using the positions r _c , r ₁ , and r ₂ of the reference station 75, satellite (#1) 71, and satellite (#2) 73, whose positions are known, Δts is calculated according to equation (28). Δts is the time difference from the reference station 75 to the satellite (#1) 71 and satellite (#2) 73. Furthermore, as shown in equation (35), Δts is subtracted from Δt _tar-ref to obtain Δtr.

Using the positions r _c , r ₁ , and r ₂ of the reference station 75, satellite (#1) 71, and satellite (#2) 73, whose positions are known, and the velocity vectors v ₁ and v ₂ of the satellite (#1) 71 and satellite (#2) 73, Δfs is calculated according to equation (30). Δfs is the frequency difference that occurs from the reference station 75 to the satellite (#1) 71 and satellite (#2) 73. Furthermore, as shown in equation (36), the frequency difference obtained by subtracting Δfs from Δf _tar-ref is defined as Δfr.

In the second embodiment, which uses Δtr calculated using equation (35) and Δfr calculated using equation (36), the position r of the radio wave source 70 can be estimated by eliminating the effects of factors that cause time or frequency differences between the path from satellite (#1) 71 to antenna (#1) 72 and the path from satellite (#2) 73 to antenna (#2) 74, and other factors that cause deviations.

The configuration of the position locating device according to the second embodiment will be described with reference to FIG. 8. In addition to receiving signals 91 and 92, receiving signals 93 and 94 are input to position locating device 50A. For this reason, position locating device 50A has frequency conversion units 1A, 1B, 1C, and 1D, AD conversion units 2A, 2B, 2C, and 2D, and quadrature detection units 3A, 3B, 3C, and 3D. The correlation function calculation unit 4A, peak detection unit 5A, TDOA calculation unit 7A, FDOA calculation unit 8A, and memory unit 14A have been modified.

The frequency conversion unit 1C converts the received signal 93 to an intermediate frequency. The AD conversion unit 3C converts the signal output by the frequency conversion unit 1C into a digital value. The quadrature detection unit 3C performs quadrature detection on the signal digitized by the AD conversion unit 2C and outputs a complex baseband signal _s3 (t). The frequency conversion unit 1D, the AD conversion unit 2D, and the quadrature detection unit 3D perform similar processing on the received signal 94. The quadrature detection unit 3D outputs a complex baseband signal _s4 (t).

The correlation function calculation unit 4A calculates the correlation function CAF by equation (1) for _s1 (t) and _s2 (t) in a determined period T. The peak detection unit 5A detects Δt and Δf at which the correlation function CAF has the maximum peak. Δt is the peak time difference at which the correlation function CAF between _s1 (t) and _s2 (t) has a peak. Δf is the peak frequency difference at which the correlation function CAF between _s1 (t) and _s2 (t) has a peak.

The correlation function calculation unit 4A and the peak detection unit 5A operate similarly for _s3 (t) and _s4 (t) of a determined period T. The correlation function calculation unit 4A calculates the correlation function CAF by equation (1) for _s3 (t) and _s4 (t) of a determined period T. The peak detection unit 5A detects _Δtref and _Δfref at which the correlation function CAF has the maximum peak. _Δtref is the peak time difference between _s3 (t) and _s4 (t) at which the correlation function CAF has a peak. _Δfref is the peak frequency difference between _s3 (t) and _s4 (t) at which the correlation function CAF has a peak.

The TDOA calculation unit 7A calculates Δt _tar-ref by subtracting Δt _ref from Δt as shown in equation (33). Using positions r ₁ , r ₂ , and _rc , it calculates Δts as shown in equation (28). Furthermore, it calculates Δtr by subtracting Δts from Δt _tar-ref as shown in equation (35). The TDOA calculation unit 7A calculates Δtr based on Δt, Δt _ref , and Δts.

The FDOA calculation unit 8A calculates Δf _tar-ref by subtracting Δf _ref from Δf as shown in equation (34). Using the positions r ₁ , r ₂ , _rc and the velocity vectors v ₁ , v ₂ , it calculates Δfs as shown in equation (30). Furthermore, it calculates Δfr by subtracting Δfs from Δf _tar-ref as shown in equation (36). The FDOA calculation unit 8A calculates Δfr based on Δf, Δf _ref , and Δfs.

The storage unit 14A changes the data stored in relation to data necessary for processing and data used during processing, as compared with the storage unit 14. The storage unit 14A does not store positions r _m1 , r _m2 , Δtm, Δfm, vectors u _1m , u _2m, etc. The storage unit 14A stores positions r _c , Δtc, Δfc, vectors u _1c , u _2c, etc.

The operation will be described. FIG. 9 is a flow chart for explaining the operation of the position locating device according to the second embodiment. The difference between FIG. 9 and FIG. 5 for the first embodiment will be described. In step S01A, the antenna (#1) 72 receives the radio wave 83 to generate the first reception signal s ₁ (t), and receives the radio wave 87 to generate the third reception signal s ₃ (t). At the same time, the antenna (#2) 74 receives the radio wave 84 to generate the second reception signal s ₂ (t), and receives the radio wave 88 to generate the fourth reception signal s ₄ (t). In step S02A, the peak detection unit 5 detects the peak time difference Δt and the peak frequency difference Δf between s ₁ (t) and s ₂ (t), and detects the peak time difference Δt _ref and the peak frequency difference Δf _ref between s ₃ (t) and s ₄ (t). In step S03A, the TDOA calculation unit 7A calculates Δtr, which is the TDOA, from Δt and Δt _ref . The FDOA calculation unit 8A calculates Δfr, which is the FDOA, from Δf and Δf _ref .

The operations from step S04 onwards are the same as those of the position locating device 50.

By also using the third and fourth received signals obtained by receiving radio waves transmitted by the reference station, the time or frequency offset present in the first and second received signals is removed. By using the reference station, the accuracy of locating the radio wave source is increased compared to when the reference station is not used.

Embodiment 3.
The third embodiment is based on the first embodiment and is configured to output a position on an iso-TDOA line as radio source position estimation information. Fig. 10 is a block diagram showing the configuration of a position locating device according to the third embodiment. Regarding Fig. 10, the points different from Fig. 3 in the first embodiment will be explained. Similar changes may be made based on the second embodiment or other embodiments.

The position locating device 50B further includes an assumed information input unit 18. The assumed information input unit 18 is used by the user to input the assumed position range 31 and the assumed speed range 32B. The assumed speed range 32B is a range that limits not only the range of the magnitude |v| of the speed vector v, but also the direction of v. The assumed information input unit 18 is an input means by which the user inputs the assumed position range 31 and the assumed speed range 32B.

In the positioning device 50B, the radio wave source velocity estimation unit 12B, the second positioning unit 13B, and the memory unit 14B have been modified. The second positioning unit 13B has a modified FDOA correction unit 16B. The radio wave source velocity estimation unit 12B searches for a velocity vector v within the assumed velocity range 32B. The radio wave source velocity estimation unit 12B reads the assumed velocity range 32B immediately upon startup.

In the second positioning unit 13B, the FDOA correction unit 16B is changed. The FDOA correction unit 16B calculates Δfc(r TDOAk ), which is an estimated FDOA, for each r _TDOAk , which is a position arranged on an iso-TDOA line included in the assumed position range 31. The FDOA correction unit 16B calculates Δfc(r _TDOAk ) for each r _TDOAk by equation ( ₃₇ ).

The second positioning unit 13B also has a presence accuracy calculation unit 19 and a maximum likelihood position estimation unit 20. The presence accuracy calculation unit 19 obtains a cone including all the velocity vectors v calculated for each r _TDOAk . The cone defines a range of directions including all the velocity vectors v. The range of directions of the velocity vectors v defined by the cone is called an allowable direction range 66. Since the radio wave source velocity estimation unit 12B obtains the velocity vectors v within the estimated velocity range 32B, the allowable direction range 66 is included in the estimated velocity range 32B. The presence accuracy calculation unit 19 calculates the proportion (presence accuracy) of the allowable direction range 66 in the assumed velocity range 32B. Furthermore, if the presence accuracy is equal to or greater than a determined lower limit (presence judgment lower limit 33), it is determined that there is a possibility that a radio wave source 70 exists at that point. The r _TDOAk in which there is a possibility that a radio wave source 70 exists is called a possible position 67. The position on the iso-TDOA line 61 closest to the center of gravity obtained for all possible positions 67 using the existence probability as a weighting coefficient is called the most likely position 68. The maximum likelihood position estimation unit 20 obtains the most likely position 68. The most likely position 68 is a position that is considered to be more likely to have a radio wave source 70 present than the possible positions 67. The possible positions 67 and the most likely position 68 are output as radio wave source position estimation information.

The presence determination lower limit value 33 is set to, for example, 10%. The presence determination lower limit value 33 may be defined as a ratio to the maximum presence probability within the assumed position range 31. The presence determination lower limit value 33 may be a fixed value and a ratio to the maximum presence probability, and may be multiple values used in an AND condition or an OR condition. A means for a user to input the presence determination lower limit value 33 may be provided.

For a position r _TDOAk determined as a possible position 67, the estimated iso-FDOA line calculation unit 17 calculates an estimated iso-FDOA line 64 which is the position r at which Δfc(r _TDOAk ), which is the estimated FDOA calculated by equation (37), is obtained.

Memory unit 14B stores assumed position range 31, assumed speed range 32B, and presence determination lower limit value 33. Memory unit 14B stores data different from memory unit 14 as data to be used during processing.

In the assumed velocity range 32B, the elevation angle and azimuth angle in the horizontal coordinate system are used to specify the direction of the velocity vector v. A method of transforming the horizontal coordinate system into the XYZ coordinate system will be described.
Define the following variables:
φ: latitude of position r of the radio source 70.
θ: longitude of position r of the radio wave source 70.
ξ: Elevation angle of the velocity vector v of the radio wave source 70 in the horizontal coordinate system.
ψ: Azimuth angle in the horizontal coordinate system of the velocity vector v of the radio wave source 70. ψ = 0 degrees in the south direction. ψ = 90 degrees in the east direction.
ξmax: The upper limit of the elevation angle of the velocity vector v of the radio wave source 70 in the horizontal coordinate system.
ξmin: The lower limit of the elevation angle of the velocity vector v of the radio wave source 70 in the horizontal coordinate system.
ψmax: the upper limit of the azimuth angle of the velocity vector v of the radio wave source 70 in the horizontal coordinate system.
ψmin: the lower limit of the azimuth angle of the velocity vector v of the radio wave source 70 in the horizontal coordinate system.
[M]: Transformation matrix from horizontal coordinate system to XYZ coordinate system.

In the expected velocity range 32B, the range of the velocity vector v is limited by the formula shown below.

|ｖ|に対する上限値ｖ_maxは、例えば秒速１５ｍに設定する。下限値ｖ_minは、例えば秒速０ｍに設定する。仰角ξに対する上限値ξ_maxは、例えば１２[度]に設定する。下限値ξ_minは、例えば－１２[度]に設定する。これは、自動車が通る道路の勾配が最大で１０度程度であることに合わせている。方位角ψは、範囲を指定しない場合は、ψmin＝－１８０[度]、ψmax＝１８０[度]に設定する。北を含む範囲で方位角ψを制限する場合は、例えばψmin＜１８０[度]＜ψmaxとなるように、ψminおよびψmaxを設定する。

The upper limit v _max for |v| is set to, for example, 15 meters per second. The lower limit v _min is set to, for example, 0 meters per second. The upper limit ξ _max for the elevation angle ξ is set to, for example, 12 degrees. The lower limit ξ _min is set to, for example, -12 degrees. This is in accordance with the fact that the gradient of the road on which cars travel is a maximum of about 10 degrees. If no range is specified for the azimuth angle ψ, then ψmin = -180 degrees, ψmax = 180 degrees. If the azimuth angle ψ is limited to a range that includes north, then ψmin and ψmax are set so that, for example, ψmin < 180 degrees < ψmax.

The velocity vector v expressed in polar coordinates in the horizontal coordinate system can be converted to a value in the XYZ Cartesian coordinate system using the transformation matrix [M] as shown in equation (41) below.

The transformation matrix [M] is expressed by the following equation (42).

When the radio wave source velocity estimation unit 12B is started, it immediately reads the assumed velocity range 32B and converts the assumed velocity range 32B into a range in the XYZ coordinate system using equation (41). When searching for the velocity vector v, the radio wave source velocity estimation unit 12B searches within the assumed velocity range 32B.

The velocity vector V in the Cartesian coordinate system can be transformed into a velocity vector expressed in polar coordinates in the horizontal coordinate system by the transformation matrix [M ⁻¹ ], as shown in equation (43) below.

The transformation matrix [M ⁻¹ ] is expressed by the following equation (44).

In order to explain the processing of the existence probability calculation unit 19, the following variables are defined.
G(r _TDOAk ): the area of the allowable direction range 66 at the position r _TDOAk .
Gdet: Area of assumed velocity range 32B.
γ(r _TDOAk ): The ratio of G(r _TDOAk ) to Gdet. This is called the presence probability.
γmin: a lower limit value for the presence probability γ for determining that the radio wave source 70 exists in r _TDOAk . This is stored in the storage unit 14B as the presence determination lower limit value 33.
r _cent : the weighted centroid of r _TDOAk, which is a possible position 67 .
r _most : represents the most likely position 68. The position on the TODA line 61 closest to _{r cent} .

The presence accuracy calculation unit 19 converts all the velocity vectors v calculated for each r _TDOAk into velocity vectors in the horizontal coordinate system using equation (43). A range that includes all the velocity _vectors v in the horizontal coordinate system is obtained. The obtained range is the allowable direction range 66. As shown in equation (45), the area of the allowable direction range 66, G(r _TDOAk ), is divided by the area of the estimated speed range 32B to calculate the presence accuracy γ(r _TDOAk ).

Here, the area of the estimated velocity range 32B is the area of the range that the estimated velocity range 32B occupies in the ξψ plane that is set so that the elevation angle ξ and the azimuth angle ψ are orthogonal. The allowable direction range 66 is also an area in the ξψ plane.

The possible positions 67 are r _TDOAk whose presence probability γ(r _TDOAk ) is equal to or greater than γmin (γ(r _TDOAk )≧γmin). The maximum likelihood position estimator 20 weights all the possible positions 67, r _TDOAk, to find the center of gravity r _cent .

式(46)において、Σは、γ(ｒ_TDOAk)≧γminであるｒ_TDOAkすなわち各可能位置６７ついて和を取ることを意味する。ｒ_centに至近のＴＯＤＡ線６１上の位置を求めて、最尤位置６８であるｒ_mostとする。

In equation (46), Σ means to take the sum over r _TDOAk where γ(r _TDOAk )≧γmin, that is, over each possible position 67. The position on the TODA line 61 closest to r _cent is found and set as r _most , which is the most likely position 68.

The presence accuracy γ(r _TDOAk ) is a possible speed ratio, which is the ratio of the range of the speed vector v estimated by the radio wave source speed estimation unit 12B to the assumed direction range 32B. The denominator on the right side of equation (46) is the total possible speed ratio, which is the sum of the presence accuracy γ(r _TDOAk ) of the speed vector v estimated at all possible positions 67. The presence accuracy γ(r _TDOAk ) of each possible position 67 divided by the total possible speed ratio is used as a weighting coefficient to determine a position obtained by weighting and averaging each r _TDOAk , i.e., the weighted center of gravity r _cent . The most likely position 68 is determined based on the weighted center of gravity r _cent .

The operation will be described. FIG. 11 is a flow chart for explaining the operation of the position locating device according to the third embodiment. The difference between FIG. 11 and FIG. 5 in the first embodiment will be described. In step S09B, the radio wave source speed estimating unit 12B estimates the speed vector v of the radio wave source 70 in consideration of the assumed speed range 32B that also restricts the direction of the speed vector v. After S09B, in step S13, an allowable direction range 66 including all the speed vectors v calculated for each r _TDOAk is obtained. In step S14, the existence probability γ(r _TDOAk ) of each r _TDOAk is calculated from the allowable direction range 66. In step S15, the r _TDOAk for which γ(r _TDOAk )≧γmin is set as the possible position 67. In step S16, the weighted center of gravity r _cent of the possible position 67 is calculated. In step S17, the position on the TODA line 61 closest to r _cent is obtained and set as r _most , which is the most likely position 68. In step S11B, the estimated iso-FDOA line calculation unit 17 calculates an estimated iso-FDOA line 64 for the possible position 67. In step S12B, the second positioning unit 13B outputs the iso-TDOA line 61, the estimated iso-FDOA line 64, the possible position 67, and the most likely position 68 contained within the assumed position range 31 as radio wave source position estimation information. The estimated iso-FDOA line 64 passing through the possible position 67 or the possible position 67 does not have to be output as the radio wave source position estimation information.

In step S18, a check is made to see if the assumed position range 31 or assumed speed range 32B has been changed. If the assumed position range 31 or assumed speed range 32B has been changed (YES in S18), the process returns to S09B. If the assumed position range 31 or assumed speed range 32B has not been changed (NO in S18), the process ends.

Figures 12 and 13 show examples of the results of positioning performed by the positioning device 50B. Figure 12 shows the case where constraints are set on the magnitude |v| of the velocity vector v and the elevation angle ξ of v, but the azimuth angle ψ of v is not constrained. Figure 13 shows the case where the magnitude |v|, elevation angle ξ, and azimuth angle ψ of the velocity vector v are constrained.

The position locating device 50B also generates possible positions 67 and most likely positions 68 as radio wave source position estimation information 65B. In the figure, the possible positions 67 are displayed as white circles, and the most likely positions 68 are displayed as crosses.

The constraint conditions that were set are as follows: The case shown in FIG.
γ min = 0.1
(Case 1)
|v|≦15 [m/sec]
-12[degrees]≦ξ≦12[degrees]
-180[degrees]≦ψ≦180[degrees]
(Case 2)
13[m/sec]≦|v|≦15[m/sec]
-12[degrees]≦ξ≦12[degrees]
135[degrees]≦ψ≦145[degrees]

The presence probability γ at each position r _TDOAk on the iso-TDOA line 61 included in the assumed location range 31 is as shown in Fig. 14. In Fig. 14, the horizontal axis represents the latitude φ of the position r _TDOAk . In both Case 1 and Case 2, the presence probability γ>γmin at each position r _TDOAk included in the assumed location range 31. Therefore, each position r _TDOAk on the iso-TDOA line 61 included in the assumed location range 31 is a possible position 67. The user can refer to the presence probability γ for each possible position 67.
For ease of viewing, four possible positions 67 are shown. Any suitable number of possible positions 67 may be displayed at suitable intervals.

Between Case 1 and Case 2, the assumed direction range 32B is set narrower in Case 2. Therefore, the existence probability γ of the position r _TDOAk with a large latitude φ is large. As a result, the most likely position 68 in Case 2 is a position with a larger latitude than the most likely position 68 in Case 1. In a situation where the velocity vectors v ₁ and v ₂ of the satellite (#1) 71 and the satellite (#2) 73 are different, if the constraints on the magnitude |v| and direction of the velocity vector v are set narrower, there may be cases where the positions r _TDOAk on the iso-TDOA line 61 included in the assumed position range 31 are not feasible positions 67. In addition, the most likely position 68 changes depending on the set assumed direction range 32B.

In the positioning device 50B, even if both relaying satellites are geostationary satellites and positioning is not possible, the velocity vector v of the radio wave source 70 can be estimated and radio wave source position estimation information 65B based on the velocity vector v can be presented. When both relaying satellites are geostationary satellites, information regarding the position of the radio wave source 70 can be presented more appropriately than in the past.

The estimated iso-FDOA line 64 corresponding to the most likely position 68 may be displayed. Although the estimated iso-FDOA line 64 is displayed corresponding to the possible position 67, the estimated iso-FDOA line 64 does not have to be displayed.

The possible position 67 with the maximum existence probability γ among all possible positions 67 may be determined as the most likely position 68. The possible position 67 with the maximum existence probability γ may be found, for example, by a binary search method.

The presence probability γ(r _TDOAk ) may be calculated as a ratio of the solid angle of the allowable direction range 66 and the assumed speed range 32B in the Cartesian coordinate system. Alternatively, the presence probability γ(r _TDOAk ) may be calculated using a weighted sum. The weighted sum is a function that calculates G(r _TDOAk ) using a weighting function g(ξ, ψ). The weighting function g(ξ, ψ) is a function that returns 0 when a certain direction (ξ, ψ) is not included in the allowable direction range 66, and returns a non-zero value when the direction is included. When the weighting function g(ξ, ψ) returns a non-zero value, it changes the value returned depending on at least one of the values of ξ or ψ. In this way, it is possible to calculate G(r _TDOAk ) by reducing the weight g(ξ, ψ) in the range close to the upper or lower limit of ξ or ψ. In addition, calculating G(r _TDOAk ), which is the area on the ξψ plane of the allowable direction range 66, is equivalent to integrating a function (assuming it is h(ξ, ψ)) that returns 0 if the direction (ξ, ψ) is not included in the allowable direction range 66 and returns 1 if the direction (ξ, ψ) is included in the allowable direction range 66, over ξ and ψ that are included in the assumed speed range 32B.
The above also applies to other embodiments.

Embodiment 4.
The fourth embodiment is based on the third embodiment, but changes the method of estimating the velocity vector v representing the movement of the radio wave source. A position locating device 50C according to the fourth embodiment calculates a radio wave source frequency difference ε _FDOA (v, r _TDOAk ), which is equal to the difference between the FDOA Δfr and the corrected FDOA Δfc(r _TDOAk ), for each r _TDOAk .

Figure 15 is a block diagram showing the configuration of a position locating device 50C according to embodiment 4. Regarding Figure 15, differences from Figure 10 in the case of embodiment 3 will be described. In the position locating device 50C, the radio wave source speed estimating unit 12C and the presence probability calculating unit 19C have been modified.

The radio wave source velocity estimation unit 12C calculates the velocity vector v of the radio wave source 70 that satisfies equation (48) for _eFDOA (v, _rTDOAk ) calculated for each _rTDOAk by equation (47). Note that Δfc( _rTDOAk ) is calculated by equation (37).

To explain how to compute the velocity vector v that satisfies equation (48), we define the following variables:
u _d (r _TDOAk ): the difference vector between u ₁ (r _TDOAk ) and u ₂ (r _TDOAk ).
(u _x , u _y , u _z ): A vector expressing u _d (r _TDOAk ) in Cartesian coordinates.
φd: The angle between u _d (r _TDOAk ) and the XY plane when expressed in polar coordinates.
θd: The angle between u _d (r _TDOAk ) and the YZ plane when expressed in polar coordinates.
ζ: the angle between the velocity vector v and u _d (r _TDOAk ).
ρ: the rotation angle of the velocity vector v around u _d (r _TDOAk ).
ζmax: The upper limit of the angle ζ at which _еFDOA (v, _rTDOAk ) is obtained.
ζmax: The lower limit of the angle ζ at which _еFDOA (v, _rTDOAk ) is obtained.

u _d (r _TDOAk ) is calculated by the following equation (49): u ₁ (r _TDOAk ) and u ₂ (r _TDOAk ) are calculated by using r=r _TDOAk in equations (4) and (5).

The angles φd and θd when expressing u _d (r _TDOAk ) in polar coordinates are calculated by the following formulas.

By substituting equation (49) into equation (48), we obtain equation (53) shown below. Transforming equation (53) gives equation (54).

Solving equation (54) for cosζ gives equation (55).

The angle ζ that satisfies equation (55) is found by considering equation (38), which is the constraint on the magnitude |v| of the velocity vector v. Since the cosine function can only take values less than 1, in order for angle ζ that satisfies equation (55) to exist, vmax, which is the upper limit of |v|, must satisfy equation (56) shown below. If vmax does not satisfy equation (56), then there is no velocity vector v that satisfies equation (55).

When vmax satisfies equation (56), the upper limit ζmax of angle ζ is calculated using equation (57).

When vmax satisfies equation (56), the lower limit ζmin of the angle ζ is calculated using the equation below, taking into account the case-specific value of vmin.

If equation (56) does not hold for each r _TDOAk , the radio wave source velocity estimation unit 12C determines that the velocity vector v does not exist. If equation (56) holds, the upper limit value ζmax of the angle ζ is calculated using equation (57), and the lower limit value ζmin of the angle ζ is calculated using equation (58-1) or equation (58-2). Using ζmax and ζmin, the following ranges are obtained for the velocity vector v that satisfies equations (49) and (38).

The velocity vector v is a vector whose angle ζ with vector ud(r _TDOAk ) is less than ζmax and greater than ζmin, and whose component parallel to ud(r _TDOAk ) is (c/f)(е _FDOA (v,r)/|ud(r _TDOAk )|). If the starting points of all velocity vectors v that satisfy equations (59) to (61) are made to coincide with the starting point of ud(r _TDOAk ), the tip of v will exist in the shape of a circle or donut on a plane perpendicular to ud(r _TDOAk ).

The velocity vector v can be converted to a value in the XYZ coordinate system using the transformation matrix [T] as shown below.

The transformation matrix [T] is expressed as shown below.

The direction of the velocity vector v defined by equations (57) to (60) is the direction of the velocity vector v that satisfies equations (49) and (38). The range of directions of the velocity vector v that satisfy equations (49) and (38) is called the range of possible directions 69.

The presence probability calculation unit 19C converts the possible direction range 69 expressed in the XYZ coordinate system into a range of directions in the horizontal coordinate system using equation (43). The directions expressed by the angles ζ and ρ may be directly converted into directions in the horizontal coordinate system. In that case, the conversion can be performed using the equation shown below.

As shown in FIG. 16, the presence probability calculation unit 19C obtains the overlapping portion of the possible direction range 69 and the assumed speed range 32B as the allowable direction range 66. FIG. 16 is a diagram illustrating the relationship between the possible direction range 69, the assumed speed range 32B, and the allowable direction range 66 by way of an example. In FIG. 16, only the positive range of the azimuth angle ψ is displayed. The azimuth angle ψ changes within a range of 360 degrees (2π rad). The allowable direction range 66 is a direction that is included in the possible direction range 69 and also in the assumed speed range 32B. In FIG. 16, the assumed speed range 32B is set to |ξ|≦12 degrees, 75 degrees≦ψ≦150 degrees. The range of directions that is included in both the possible direction range 69, which is a ring-shaped region, and the assumed direction range 32 is the allowable direction range 66. In FIG. 16, the allowable direction range 66 is displayed with hatching.

The operation of the position locating device according to the fourth embodiment will be described. FIG. 17 is a flow chart for explaining the operation of the position locating device according to the fourth embodiment. The differences between FIG. 17 and FIG. 11 for the third embodiment will be explained. In step S09C, the radio wave source velocity estimation unit 12C estimates the velocity vector v of the radio wave source 70, taking into account the magnitude |v| of the velocity vector v, without taking into account the direction of v. The range including the directions of all the estimated velocity vectors v is set as the possible direction range 69. In step S13C, the existence probability calculation unit 19C determines the overlapping portion of the possible direction range 69 and the assumed velocity range 32B as the allowable direction range 66.
In other respects it operates similarly to position locator 50B.

In the positioning device 50C, even if both relaying satellites are geostationary satellites and positioning is not possible, the velocity vector v of the radio wave source 70 can be estimated and radio wave source position estimation information 65B based on the velocity vector v can be presented. When both relaying satellites are geostationary satellites, information regarding the position of the radio wave source 70 can be presented more appropriately than in the past.

The position locating device 50C analytically determines the possible direction range 69 of the velocity vector v of the radio wave source 70. Therefore, the position locating device 50C can determine the possible direction range 69 more efficiently than the position locating device 50B.

It is possible to freely combine the various embodiments, modify the various embodiments, omit some of the components, or freely combine the various embodiments with some of the components omitted or modified.

50, 50A, 50B, 50C, 50X position location device,

1A, 1B, 1C, 1D Frequency conversion unit,
2A, 2B, 2C, 3D AD conversion section,
3A, 3B, 3C, 3D quadrature detection section,
4, 4A correlation function calculation unit,
5, 5A Peak detection unit,
6 satellite data acquisition unit,
7, 7A TDOA calculation section,
8, 8A FDOA calculation section,
9 Equal TDOA line calculation section,
10 Geostationary satellite identification unit,
11 first positioning unit,
12, 12B, 12C radio wave source speed estimation section,
13, 13B, 13C second positioning unit,
14, 14A, 14B, 14C storage section,
15 IsoFDOA line calculation section,
16, 16B FDOA Amendment Department,
17 Estimated etc. FDOA line calculation unit,
18 assumption information input unit,
19, 19C existence probability calculation unit,
20 Maximum likelihood position estimation unit,

31 expected location range,
32, 32B Expected speed range,
33 Presence judgment lower limit,

60 Earth,
61 Equal TDOA line (equal arrival time difference line),
62 Equal FDOA lines (equal arrival frequency difference lines),
63 location position,
64 Estimated equal FDOA line (estimated equal frequency difference line),
65, 65B Radio source position estimation information,
66 allowable direction range,
67 possible positions,
68 most likely position,
69 Possible direction range,

70 Radio Sources,
71 Satellite (#1) (first satellite),
72 Antenna (#1),
73 Satellite (#2) (Second Satellite),
74 Antenna (#2),
75 Reference Station (Reference Station),
76 geostationary orbit,
77 Equator,
78 a meridian passing through the middle of the satellite (#1) 71 and the satellite (#2) 73;

81 Radio waves,
82 Radio waves,
83 Radio waves,
84 Radio waves,
85 Radio waves,
86 Radio waves,
87 Radio waves,
88 Radio waves,
91 received signal (first received signal),
92 received signal (second received signal),
93 received signal (third received signal);
94 received signal (fourth received signal).

Claims (13)

a correlation function calculation unit which receives a first signal and a second signal expressed by complex numbers, and a time difference and a frequency difference between the first signal and the second signal, calculates a two-dimensional correlation value between the first signal and the second signal by integrating a value obtained by multiplying the first signal by a complex conjugate of the second signal having the time difference with respect to the first signal and multiplying the signal obtained by multiplying the signal by a complex number having an absolute value of 1 rotated by the frequency difference, for a determined time;
a peak detection unit that calculates the two-dimensional correlation value by changing the time difference and the frequency difference, and detects a peak time difference, which is the time difference at which the two-dimensional correlation value has a maximum peak, and a peak frequency difference, which is the frequency difference at which the two-dimensional correlation value has a maximum peak;
an arrival time difference calculation unit that calculates an arrival time difference, which is the difference between the arrival time of the radio wave from the radio wave source to the first satellite and the arrival time of the radio wave from the radio wave source to the second satellite, from the peak time difference of the two-dimensional correlation value calculated using a first reception signal, in which a radio wave radiated from a radio wave source whose position is unknown is relayed by a first satellite and received at a first observation station, as the first signal, and a second reception signal, in which the radio wave radiated from the radio wave source is relayed by a second satellite different from the first satellite and received at a second observation station, as the second signal, using the positions of the first satellite and the second satellite;
an arrival frequency difference calculation unit that calculates an arrival frequency difference, which is a difference between a frequency shift occurring from the radio wave source to the first satellite and a frequency shift occurring from the radio wave source to the second satellite, from the peak frequency difference of the two-dimensional correlation value calculated with the first reception signal as the first signal and the second reception signal as the second signal, using positions and velocity vectors of the first satellite and the second satellite;
an equal arrival time difference line calculation unit that calculates an equal arrival time difference line, which is a line representing the position of the radio wave source from which the arrival time difference is obtained, using the positions of the first satellite and the second satellite;
a geostationary satellite determination unit that determines whether the first satellite and the second satellite are geostationary satellites;
a radio wave source velocity estimation unit that estimates a radio wave source velocity, which is a velocity vector that represents the movement of the radio wave source, from the arrival frequency difference when the geostationary satellite determination unit determines that both the first satellite and the second satellite are geostationary satellites;
a second positioning unit that, when it is determined in the geostationary satellite determination unit that both the first satellite and the second satellite are geostationary satellites, generates radio wave source position estimation information regarding the position of the radio wave source using the radio wave source speed, the equal arrival time difference line, and the positions and velocity vectors of the first satellite and the second satellite. The positioning device according to claim 1, wherein the radio wave source speed estimation unit determines the radio wave source speed that is present in an assumed speed range, which is a range in which the radio wave source speed is assumed to exist, at a position on the equal arrival time difference line that is included in an assumed location range, which is a range in which the position of the radio wave source is assumed to exist. an assumed position range storage unit that stores the assumed position range;
an assumed speed range storage unit that stores the assumed speed range;
3. The position locating device according to claim 2, further comprising an input means for a user to input the assumed position range and the assumed velocity range. the radio wave source velocity estimation unit uses the positions and velocity vectors of the first satellite and the second satellite to estimate the radio wave source velocity so that a difference between a sum of a satellite frequency difference, which is the difference between a frequency shift generated due to the movement of the first satellite and a frequency shift generated due to the movement of the second satellite, at an assumed position set on the equal arrival time difference line, and a radio wave source velocity frequency difference, which is the difference between a frequency shift generated due to the radio wave source velocity up to the first satellite and a frequency shift generated due to the radio wave source velocity up to the second satellite, and the arrival frequency difference falls within a determined range;
The second positioning unit includes an estimated iso-frequency difference line calculation unit that uses the radio wave source speed and the positions and speed vectors of the first satellite and the second satellite to determine the satellite frequency difference at the assumed position as a value obtained by subtracting the radio wave source speed frequency difference from the arrival frequency difference, and calculates an estimated iso-frequency difference line that is the position of the radio wave source from which the satellite frequency difference is obtained;
4. The position locating device according to claim 1, wherein the radio wave source position estimation information is the equal arrival time difference line and the estimated equal frequency difference line. a satellite frequency difference calculation unit that calculates a satellite frequency difference, which is a difference between a frequency shift caused by the movement of the first satellite and a frequency shift caused by the movement of the second satellite, by using positions and velocity vectors of the first satellite and the second satellite with respect to an assumed position where the radio wave source is assumed to exist;
the radio wave source velocity estimation unit estimates the radio wave source velocity included in the assumed velocity range according to the assumed position using positions and velocity vectors of the first satellite and the second satellite so that a radio wave source frequency difference, which is a difference between a frequency shift caused by the radio wave source velocity occurring up to the first satellite and a frequency shift caused by the radio wave source velocity occurring up to the second satellite, is a frequency difference obtained by subtracting the satellite frequency difference for the assumed position set on the equal arrival time difference line included in the assumed position range from the arrival frequency difference;
4. The positioning device according to claim 2, wherein the second positioning unit uses an estimated equal frequency difference line calculated using the equal arrival time difference line and the estimated radio wave source speed as the radio wave source position estimation information. The positioning device according to claim 5, wherein the radio wave source position estimation information is a possible position, which is the assumed position where a possible speed ratio, which is the ratio of the range of the radio wave source speed estimated by the radio wave source speed estimation unit to the assumed speed range, is equal to or greater than a determined lower limit. the second positioning unit includes a maximum likelihood position estimation unit that determines a maximum likelihood position that is estimated to be more likely to include the radio wave source than the possible positions based on the possible speed ratio of the radio wave source speeds estimated at each of the assumed positions set on the equal arrival time difference line included in the assumed position range,
7. The position locating device according to claim 6, wherein the radio wave source position estimation information is the most likely position, or the most likely position and the possible positions. The positioning device according to claim 7, wherein the maximum likelihood position estimation unit determines the maximum likelihood position based on a weighted average position of each of the possible positions, using a weighting coefficient obtained by dividing the possible speed ratio of the radio wave source speed estimated at each of the possible positions by a total possible speed ratio, which is the sum of the possible speed ratios of the radio wave source speeds estimated at all of the possible positions. The positioning device according to claim 7, wherein the maximum likelihood position estimator determines the possible position at which the possible speed ratio is maximized as the maximum likelihood position. the second positioning unit includes an estimated iso-frequency difference line calculation unit that calculates an estimated iso-frequency difference line, which is the position of the radio wave source from which the satellite frequency difference is obtained, by using positions and velocity vectors of the first satellite and the second satellite;
The position locating device according to claim 5 , wherein the radio wave source position estimation information includes the estimated iso-frequency difference line. the arrival time difference calculation unit calculates the arrival time difference based on a reference station peak time difference which is the peak time difference of the two-dimensional correlation value calculated by setting as the first signal a third reception signal in which a radio wave radiated by a reference station whose position is known is relayed by the first satellite and received by the first observation station, and as the second signal a fourth reception signal in which a radio wave radiated by the reference station is relayed by the second satellite and received by the second observation station, the first reception signal being the first signal, and a radio wave source peak time difference which is the peak time difference of the two-dimensional correlation value calculated by setting as the first reception signal the second reception signal being the second signal, and a reference station time difference which is the difference between the arrival time of the radio wave from the reference station to the first satellite and the arrival time of the radio wave from the reference station to the second satellite,
11. The position locating device according to claim 1, wherein the arrival frequency difference calculation unit calculates the arrival frequency difference based on a reference station peak frequency difference which is the peak frequency difference of the two-dimensional correlation value calculated using the third received signal as the first signal and the fourth received signal as the second signal, a radio wave source peak frequency difference which is the peak frequency difference of the two-dimensional correlation value calculated using the first received signal as the first signal and the second received signal as the second signal, and a reference station frequency difference which is the difference between a frequency shift occurring from the reference station to the first satellite and a frequency shift occurring from the reference station to the second satellite. the arrival time difference calculation unit calculates the arrival time difference using positions of the first observation station, the second observation station, the first satellite, and the second satellite;
11. The positioning device according to claim 1, wherein the arrival frequency difference calculation unit calculates the arrival frequency difference using positions of the first observation station and the second observation station, and positions and velocity vectors of the first satellite and the second satellite. 13. The positioning device according to claim 1, further comprising a first positioning unit that, when the geostationary satellite discrimination unit determines that at least one of the first satellite and the second satellite is not a geostationary satellite, uses positions and velocity vectors of the first satellite and the second satellite to determine, as the position of the radio wave source, a position on the equal arrival time difference line where the arrival frequency difference is obtained.

Images