When an incoming Global Positioning System (GPS) signal is acquired, pullin search performs a finer search of the Doppler frequency of the incoming signal so that phase lock loop can be quickly stabilized and the receiver can produce an accurate pseudorange measurement. However, increasing the accuracy of the Doppler frequency estimation often involves a higher computational cost for weaker GPS signals, which delays the position fix. In this paper, we show that the Doppler frequency detectable by a long coherent autocorrelation can be accurately estimated using a complexweighted sum of consecutive short coherent autocorrelation outputs with a different Doppler frequency hypothesis, and by exploiting this we propose a noise resistant, lowcost and highly accurate Doppler frequency and phase estimation technique based on a reverse directional application of the finite rate of innovation (FRI) technique. We provide a performance and computational complexity analysis to show the feasibility of the proposed technique and compare the performance to conventional techniques using numerous Monte Carlo simulations.
1. INTRODUCTION
Successful positioning in harsh Global Positioning System (GPS) environments has been one of the most challenging problems in GPS, and the demand for locationbased service (LBS) in harsh environments for GPS, such as dense urban and indoor environments is continuously increasing. To meet the needs for GPS in harsh environments, GPS receivers may need a faster and higher sensitivity acquisition and tracking functions that can produce accurate carrier phase and code phase estimates in various GPS signal conditions. In the literature, a number of acquisition techniques for weak GPS signals have been introduced. For example, a GPS receiver connected to a synchronized cellular communication network (
Kransner 1998
,
van Diggelen 2009
) can make use of the downlink signal measurements and downlink information from the cellular base station to quickly acquire GPS signals with high sensitivity. On the other hand, conventional techniques to enhance the acquisition sensitivity and to reduce the mean acquisition time (MAT) such as long coherent integration and noncoherent accumulation techniques (
Parkinson et al. 1996
) have been continuously studied, and, recently, techniques such as parallel signal search, FFTbased signal search (
Borre et al. 2007
), and twodimensional compressed correlator (
Kong & Kim 2013
) techniques have been developed.
To produce accurate measurements of code phase and carrier phase of incoming GPS signals, GPS receivers employ tracking functions that require initial code phase and Doppler frequency estimates with sufficient accuracy. In practice, the accuracy required for the initial resolutions of the code phase and Doppler frequency are less than a half chip and about ten Hertz (
Tsui 2005
), respectively. Since most of the acquisition techniques perform Doppler frequency search for incoming GPS signals with a search step size about 1/(2
T
) Hz, where the coherent integration time
T
is usually 1~5 ms, the Doppler frequency estimate may have an error of up to 250 Hz, which can be too coarse for a carrier phase tracking loop. However, an initial code phase error of up to a half chip can be tolerable for most of the code phase tracking functions (
Irsigler & Eissfeller 2003
). Therefore, a GPS receiver has a pullin state between the acquisition and tracking states to increase the accuracy of the Doppler frequency estimate to a high enough level to accurately run carrier phase tracking. In the pullin state, a GPS receiver tries to obtain a fine Doppler frequency, carrier phase, and
C
/
N
_{0}
estimates of an incoming GPS signal, which can take about 1.5s (
Kelley et al. 2002
). To improve the pullin search performance, conventional techniques, such as
Tsui (2005)
and
Goiser & Berger (1996)
, perform a fast Fourier transform (FFT) to find the maximum frequency component from code wipedoff sample data and then utilize the amplitude of the component to obtain a fine Doppler frequency estimate. Consecutive fine Doppler frequency estimates can be averaged to increase accuracy in the presence of noise, but, in general, the conventional techniques work well for high
C
/
N
_{0}
. To lessen the noise effect on the Doppler frequency estimation,
Sagiraju et al. (2006)
suggest increasing
T
up to 5 ms, however, the scheme requires a smaller initial Doppler frequency estimation error than 50 Hz, which, in effect increases the complexity of the acquisition state. On the other hand,
Zeng & Li (2010)
suggest a serial Doppler frequency search using correlation in the time domain, and employ a curve fitting algorithm to find a more accurate Doppler frequency estimate. The curve fitting may be able to provide good Doppler frequency and carrier phase estimates at the cost of increased computational complexity. In general, since a GPS receiver does not yet have the knowledge on the received data bit in the pullin state, increasing the accuracy of Doppler frequency and carrier phase estimates of weak GPS signals can be a complex problem due to the computational complexity and algorithmic complexity required to mitigate the effect of bit transition. In practice, the carrier phase estimation is performed by the phase locked loop in the tracking state.
Despite the complexity and importance, there has been little attention paid to the pullin search of GPS receivers in the literature. To reduce the computational complexity and improve the accuracy of the Doppler frequency and carrier phase estimates in the pullin state, we propose a reversedirectional finite rate of innovation (FRI) technique with iterative Cadzow denoising applied to consecutive coherent autocorrelation function (ACF) outputs to reduce the computational complexity in estimating the Doppler frequency and carrier phase of an incoming GPS signal with high accuracy and to increase noise robustness. The complexity of the proposed technique is analyzed and the performance of the proposed technique is compared to the conventional techniques using numerous Monte Carlo simulations.
The rest of this paper is organized as follows. In Section 2, we introduce a general continuous form expression for a coherent ACF output, and in Section 3, we verify that a long coherent ACF output for a Doppler frequency can be estimated with consecutive short coherent ACF outputs at a close Doppler frequency. Section 4 provides the details and complexity of the proposed technique that exploits the findings in Section 3, and Section 5 shows the analysis of the performance of the proposed technique and a performance comparison to conventional techniques using numerous Monte Carlo simulations. And a conclusion is drawn in Section 6.
Throughout this paper, the following conventions are used for notation. Vectors or matrices are denoted by boldface symbols. Small letters are used for scalars and vectors, and capital letters are used for matrices.
2. COHERENT ACF OUTPUT OF GPS SIGNAL
Let
r
(
t
) be the frequency downconverted incoming
L
1 GPS coarse acquisition (C/A) code signal to an IF frequency
f_{IF}
, then the complex IF signal can be expressed as
where Q is the number of GPS satellite signals being received,
A_{q}
,
τ_{q}
,
, and
ø_{q}
are the slowly varying amplitude, code phase, Doppler frequency, and unknown carrier phase of the qth satellite signal, respectively,
D_{q}
(
t
) and
C_{q}
(
t
) are the navigation data bit signal and PseudoRandom Noise code signal of the qth satellite at time t, respectively, and n (t) represents an additive whiteGaussian noise with twosided power spectral density (PSD)
N_{0}
/2. Note that the data bit and chip rates are
R_{b}
(=50 bps) and
R_{c}
(= 1.023 Mcps), respectively, and thet the C/A code has code length
L_{c}
=1023 chips so that the code period is only
T
_{1}
(= 1 ms). When the coherent integration (correlation) length used in the acquisition state is
T
=
L_{T}T
_{1}
, where
L_{T}
is an integer much larger than 1, and code phase search resolution is
T_{c}
/2, the detected code phase
τ_{q}
and Doppler frequency
at the acquisition state can have errors upto
T_{c}
/4 and 1/(4
T
), respectively. Note that the navigation data remains unknown in the pullin state. To detect and track the code phase
τ_{q}
and Doppler frequency
of the qth satellite signal using correlation, the receiver generates a prompt replica signal
where
and
are the code phase and carrier frequency estimates made by the receiver’s acquisition function, respectively. The receiver integrates the product of
r
(
t
) and
(
t
)(complex conjugate of
r_{q}
(
t
)) over the coherent integration interval
T
to see the correlation between the two functions. Defining
and
, the mathematical expression of the ACF output over the twodimensional space is readily available from
Kong & Kim (2013)
.
where
D_{q}
(
t

τ_{q}
) is assumed as constant
D_{q}
for the integration interval
T
, and
is a complex Gaussian noise with twosided power spectral density
N_{0}
/2, since signal powers are negligible when compared to the noise power. When
T
=
T_{b}
(=1/
R_{b}
) and when
ø
,
δτ
, and
δf
are small enough,
R
_{2}
(
δτ
,
δf
)≥0 indicates that
D_{q}
=+1, and
R
_{2}
(
δτ
,
δf
)<0 indicates that
D_{q}
=1. Note that there are several techniques to remove the effect of the navigation data signal
D_{q}
(
t
) from
R
_{2}
(
δτ
,
δf
) ; a number of bit sequences are estimated using Kalman filter (
Psiaki & Jung 2002
), the binary value of
D_{q}
(
t
) is estimated by testing different bit combinations (
Petovello et al. 2008
), or with the assistance information from the cellular base station in the AssistedGPS system (
Kransner 1998
,
van Diggelen 2009
). Since figuring out the data
D_{q}
(
t
) is one of the primary tasks in the tracking state, we neglect the effect of the unknown data in this paper.
3. FREQUENCY ESTIMATION WITH CONSECUTIVE SHORT COHERENT ACF OUTPUTS
In this section, we break
R
_{2}
(
δτ
,
δf
) obtained with a long coherent interval T>>
T
_{1}
into
L_{T}
(=
T/T
_{1}
) consecutive
R
_{2}
(
δτ
,
δf
)s obtained with short coherent interval
T
_{1}
and investigate the possibility of precisely estimating the Doppler frequency of
r
(
t
) with the
L_{T}
R
_{2}
(
δτ
,
δf
) s obtained. In the following analysis, we drop the subscript (·)
_{q}
(unless necessary) for simplicity in algebraic expressions.
Let the input to the pullin state be denoted by
where
is alreadly available from
Kong & Kim (2013)
and (·)
^{l}
represents that the coherent integration from
t
=
lT
_{1}
=
t_{l}
to
t
=(
l
+1)
T
_{1}
=
t_{l}
+
T
_{1}
, and
is a complex Gaussian process with twosided PSD
N_{0}
/2. Note that in Eq. (5) the expected ACF output for 
δτ
≥
T_{c}
is assumed as zero. When the code phase estimation is successful, i.e.,
δτ
=
T_{c}
/2, the expression Eq. (5) can be simplified to a singledimensional function
and let
when the carrier frequency of
r_{q}
(t) is
, it is found from Eq. (7) that
where
which is a complex Gaussian process with the same PSD to Eq. (6). Note that Eq. (10) holds only for a very small Δ
_{f}

T
_{1}
=1. From Eqs. (7,9,10), it is found that
which indicates that a long coherent ACF output
R
_{1}
(
δf
Δ
_{f}
) with a long coherent integration interval
T
and a Doppler frequency estimate
Δ
_{f}
can be obtained from the linear sum of
L_{T}
consecutive short coherent ACF outputs
from
L_{T}
integration segments [
lT
_{1}
,(
l
+1)
T
_{1}
](
l
∈{0,1,…,
R_{T}
1}) of [0,
T
_{1}
] with a neighboring Doppler frequency
, when Δ
_{f}
 = 1/
T
_{1}
and 
δf
+Δ
_{f}
 = 1/
T
_{1}
are satisfied. In other words, Eq. (11b) states that it is possible to precisely estimate a long coherent ACF output
R
_{1}
(
δf
Δ
_{f}
) for a Doppler frequency hypothesis using consecutive short coherent ACF outputs
with a close Doppler frequency hypothesis. This observation is exploited in the proposed technique for an accurate and fast GPS pullin search. In addition, it should be noted that the carrier frequency
δf
of
can be estimated by compensating the frequency
δf
with Δ
_{f}
which maximizes the linear sum in Eq.(11b).
Fig. 1
shows 10
^{4}
Monte Carlo simulation results demonstrating the estimation accuracy of the expression Eq. (11b) for a random Doppler frequency
δf
of the incoming signal from 0Hz to 250Hz, when
C
/
N
_{0}
=40 dBHz and
L_{T}
=2
^{5}
.
Fig. 1a
shows that the result of absolute normalized amplitude difference, i.e., │
for
δf
varying from 0 to 250Hz.
Fig. 1b
shows that the standard deviation of the Doppler frequency and the estimated Doppler frequency using Eq. (11b) is small, also. As a result, it is demonstrated that
is a set of samples of a noisy complex carrier signal, with a carrier frequency
δf
, amplitude 1/
L_{T}
, and phase
πδfT
_{1}
, samples at 1/
T
_{1}
[Hz], and that estimating the carrier frequency and phase of
R
_{1}
(
δf
) can provide a mismatch of Doppler frequency and phase estimates between those of the incoming signal and the receiver replica. In addition, when the phase of
is denoted by
ø
as in Eq. (7) the accumulation of
R
_{1}
(
δf
Δ
_{f}
) for all (
l
∈{0,1,…,
R_{T}
1}) in Eq. (11a), i.e., the summation of
L_{T}
consecutive short coherent ACF outputs, will have an overall phase
ø

π
Δ
_{f}
T
_{1}
. Therefore, when the frequency error is estimated, the phase of the incoming signal can be readily available from
Validity of the approximation in Eq. (11b).
consequently, we can see thet estimating the frequency error
δf
is to search for a solution Δ
_{f}
that maximizes │
R
_{1}
(
δf
Δ
_{f}
)│ such that
which is equival to the frequency estimation of carrier signal sample
.
4. FREQUENCY AND PHASE ESTIMATION USING THE PROPOSED TECHIQUE
Ther are a number of techniques that can provide a solution for the problem in Eq. (13). For example, instead of searching for Δ
_{f}
, trivial techniques can be the Discrete Fourier Transform (DFT) and the FFT that can estimate the frequency of
R
_{1}
(
δf
) by finding the maximum frequency component such that
which may require only
L_{T}
log
_{2}
L_{T}
complex multiplications. The esimate
is asymptotically unbiased for a very large
L_{T}
(
Hayes 1996
), and the maximum error of the FFTbased scheme Eq. (14) is 1/(2
T
) [Hz]. Since the frequency estimation error with the FFT (or DFT) is uniformly distributed in [1/(2
T
), 1/(2
T
)] [Hz], the absolute mean frequency estimation error 1/(4
T
) [Hz] is not sufficiently accurate when
T
is not very large. In this section, however, we propose a reverse directional FRI technique to increase the frequency estimation accuracy when
T
is not very large.
 4.1 Proposed Reverse Directional FRI Technique
In fact, the observation
R
_{1}
(
δf
) is the vector of samples of a signal with FRI (
Vetterli et al. 2002
), i.e., the samples contain only a finite number
K
(<∞) of innovation of a periodic signal. For an example, a low pass filtered input signal
y
(
t
) composed of
T_{P}
periodic distinctive Dirac delta functions with complex weights
x_{k}
(
k
=0,1,...,
k
1) as shown in
Fig. 2a
, the FRI technique (
Vetterli et al. 2002
), is to estimate the delays
t_{k}
(
k
=0,1,...,
k
1) of the Dirac delta functions in the input using the minimum number of samples of
y
(
t
). The estimation of
t_{k}
(
k
=0,1,...,
k
1) is performed in the frequency domain, since the time domain samples of
y
(
t
) are effectively the samples of the linear sum of complex sinusoids. The delays can be accurately estimated from the zeros of the annihilating filter
H
[z] that annihilates the samples of the linear sum of complex sinusoids. And the amplitude
x_{k}
is estimated by solving the Vandermonde system of linear equations (
Vetterli et al. 2002
). On the contrary, the problem in Eq. (13) is to estimate the frequency of a single (
K
=1) complex sinusoid that is equivalently a single Dirac delta function in the frequency domain, which means that the problem in Eq. (13) requires reverse directional application of the FRI technique as shown in
Fig. 2b
and as investigated in
Kim et al. (2011)
, for synthetic aperture radar (SAR) signal processing. In other words,
R
_{1}
(
δf
) can be regarded as a set of noisy time domain samples of a single (i.e.,
K
=1) Dirac delta function in the frequency domain as
Reverse directional application of FRI for GPS.
where the sampling is performed at every
t
=
lT
_{1}
. Therefore, we can expect that the single Dirac delta function in the frequency domain is
periodic, and that the location and amplitude of the corresponding delta function in the frequency domain can indicate the carrier frequency and phase of the
R
_{1}
(
δf
), respectively.
Exploiting
Blu et al. (2008)
for the case of noisy sample input,
R
_{1}
(
δf
), when the input is arranged into a Toeplitz matrix
where
L
≥
K
, an annihilating filter
that minimizez the total least squares
under a constraint 
H

^{2}
=1 provides the frequency estimate
, i.e.,
f
_{1}
=
. The coefficients of the annihilating filter [
h
_{0}
,
h
_{1}
,…,
h_{k}
] can be found from the eigenvector of
A
^{T}
A
corresponding to the smallest eigenvalue, which is related to Pisarenko’s method (
Pisarenko 1973
,
Tufts & Kumaresan 1982
).
To increase the accuracy of the frequency estimate
in the presence of heavy noise in the signal sample input, there can be two ways to increase the noise robustness of the algorithm in Eq. (18). One way is to increase
L_{T}
, which is a trivial solution, and the other way is to apply the iterative Cadzow denoising technique used in
Blu et al. (2008)
. The first step of the technique is to perform a singular value decomposition (SVD) of
A
such that
where
U
,
S
, and
V
are unitary, diagonal with elements in descending order, and unitary matrices of size [(
L_{T}

L
)×(
L
+1)], [(
L
+1)×(
L
+1)], and [(
L
+1)×(
L
+1)], respectively. When the (
K
+1) th largest diagonal element of
S
is not sufficiently smaller than the
K
th largest diagonal element, Cadzow’s iterative denoising technique sets the
L

K
+1 smallest diagonal elements of
S
to zero to build a new diagonal matrix
S
^{(1)}
in the first step, where (·)
^{(i)}
denotes the
i
th iteration of the Cadzow denoising technique. In the second step, the technique constructs
A
^{(1)}
using
S
^{(1)}
as
The resulting
A
^{(1)}
is not Toeplitz anymore, so the technique obtains the best Toeplitz approximation of
A
^{(1)}
by averaging the diagonals of
A
^{(1)}
in the step. The above three steps are usually iterated
N_{i}
(>1) times until the (
K
+1) th largest diagonal element of
S
^{(Ni)}
becomes
R_{TH}
(>>1) times smaller than the
K
th latgest diagonal element. From experimental observation,
R_{TH}
≥50 is sufficiently high. When it is found that the Cadzow denoising is completed after
N_{i}
th iteration,
A
^{(Ni)}
can be rearranged to form a matrix
B
of size [
L_{T}
×2] as
where
’ is the denoised
by the iterative Cadzow denoising technique, and the eigenvector of
B
^{T}
B
corresponding to the smallest eigenvalue provides the coefficients of the annihilating filter [
h
_{0}
,
h
_{1}
]. From Eq. (17) and for
K
=1,
is the final estimate of the frequency
δf
buried in
R
_{1}
(
δf
) found by the proposed reverse directional FRI technique.
Once the frequency
δf
is estimated, from Eq. (12), the relative phase of the input signal to the receiver generated carrier signal can be easily obtained as
which is the arctangent of the frequency compensated and averaged Inphase and Quadraturephase consecutive short coherent ACF outputs. Note that the proposed technique can estimate both the frequency and phase difference between the incoming signal
r
(
t
) and the receiver replica signal, and that when there is no frequency, i.e.,
δf
=0, the phase estimate of the proposed technique is similar to that of the conventional phase estimation technique.
 4.2 Realization and Receiver Complexity
A drawback of the proposed technique for the pullin search maybe the high computational cost for multiple SVD’s in the iterative Cadzow denoising. To reduce the computational cost, computationally efficient SVD algorithms such as Riemannian SVD (RSVD) and GolubReinsch SVD (GRSVD), introduced in
Chan (1982)
can be applied as a solution. For an SVD of the matrix
A
of size[(
L_{T}

L
)×(
L
+1)], RSVD and GRSVD require
respectively (
Golub & van Loan 1996
). Assuming
L
=
L_{T}
/4 for Eq. (24a) and
L
=
L_{T}
/2 for Eq. (24b), and
L
>>1 for both cases, we can find that
M
_{1}
≈
M
_{2}
/3. Therefore, 1<
L
<<
L_{T}
/2 is a preferable choice, and the complexity decreases as
L
decreases in RSVD. Letting 1/
L_{T}
<<
α
_{l}
<1/2 and
L
=
α
_{l}
L_{T}
, and assuming
N_{i}
iterations of Cadzow denoising, the total computational complexity for the iterative Cadzow denoising is about
where
Ni
≤5 is usually observed in experimentations (
Blu et al. 2008
). Notice that
M_{CD}
in Eq, (25) is a monotonously increasing function with respect to
α
, so that smaller
α
results in less computational cost.
Fig. 3
shows the results of 10
^{4}
Monte Carlo simulations for a fixed
R_{TH}
=50 and for various
α
to estimate the total number of multiplications
M_{CD}
in the proposed technique employing iterative Cadzow denoising with respect to GPS signal strength.
Fig. 3a
and
b
show the simulation results for
L_{T}
=16 and
L_{T}
32, respectively. In both figures, 1/4≤
α
≤1/2 results in a similar total number of multiplications, but
α
=1/8 results in the largest total number of multiplications, which is because the number of iterations
N_{i}
is increased for smaller
α
. Overall, the computational cost is the lowest or near lowest when
α
=1/3, so
L
≈
L_{T}
/3 is a choice leading to a computationally efficient application of the iterative Cadzow denoising algorithm. Therfore, in the following, we use
L
≈
L_{T}
/3 for the proposed technique.
Complexity comparison.
5. PERFORMANCE ANALYSIS AND NUMERICAL SIMULATIONS
In this section, we compare the performance of the proposed technique and other techniques including the conventional technique (
Tsui 2005
).
Exploiting the analysis on the uncertainty of the estimate with the FRI technique for a single Dirac delta function, the variance of the unbiased frequency estimate
Eq. (22) is found in
Blu et al. (2008)
and is lower bounded by the CramerRao lower bound (CRLB) as
where
is the period of the frequency domain spectrum, and
PSNR
represents the peak signaltonoise ratio such that
Fig. 4a
shows the standard deviation of frequency estimation errors of the proposed technique and other techniques with respect to the carrier to noise denisty ratio (
C
/
N
_{0}
). The performance of the proposed techique is almost the same to the theoretical bound CRLB when
C
/
N
_{0}
is not too weak. Note that when the FFTbased frequency resolution governs the estimation accuracy so that the 64point FFT (i.e.,
L_{T}
=64) has twice the constant resolution error of the 123 point FFT as shown in
Fig. 4a
, where ‘
N_{F}
pt FFT (
N_{F}
ms)’ in the legend denote the
N_{F}
point FFT from
N_{F}
ms of data. Similarly, the performance of the serial Doppler frequency search scheme with coherent correlation length
T_{CO}
=
N_{CO}
T
_{1}
ms if lowerbounded by the Doppler frequency search step size, where
N_{CO}
=20. In the serial Doppler frequency search simulations, we chose the search step size of 10
^{3}
/64 Hz to compare with the 64 point FFTbased technique. When the curve fitting technique (
Zeng & Li 2010
) is applied to the serial Doppler frequency search result obtained with search step size of 1/20 Hz, the performance is improved and shows that it has the same slope to the CRLB in the moderate and high
C
=
N
_{0}
region.
Fig. 4a
shows the performance of the conventional techique (
Tsui2005
) also; the conventional frequency discriminator uses two consecutive ACF outputs, with coherent integration length
T
_{1}
= ms, to compute
Performamce comparison
where
I_{m}
and
Q_{m}
represent the inphase and quadraturephase component of the ACF output at
t
=
mT
_{1}
, and the moving average length of
N
1 is assumed. Note that the performance of the conventional technique is the worst among the techniques shown in
Fig. 4
for
C
=
N
_{0}
not high enough, which is due to the short coherent integration length
T
_{1}
for each ACF output, and that we can expert 3 dB more performance gain wher the coherent integration for each ACF output is doubled. Note also thet the ideal frequency estimate in Eq. (28) does not include any measurement error due to the approximation funstion. After some algebraic manipulations, the frequency error variance can be found as
Notice that when
C
/
N
_{0}
=45 dBHz the frequency errors of the proposed technique for
L_{T}
=16,
L_{T}
=32 and the curve fitting technique are about 1 Hz, 0.3 Hz, and 0.8 Hz, respectively, and the errors increase much slowly than in the conventional technique as the
C
/
N
_{0}
decreases. The resulting phase estimation performnce of the techniques compared in
Fig. 4a
shows some similarity as shown in
Fig. 4b
. Note that the phase estimation is performed with a shorter data length than the data length used for frequency estimation. For example, when the frequency is estimated with
N_{F}
point FFT (
N_{F}
) ms, the maximum frequency error is 10
^{3}
/2
N_{F}
Hz with which the phase error can increase to
π
[rad]. Since a phase estimation with an arctangent function may result in a large arithmetic error in the presence of noise, we limit the maximum possible phase estimation error by taking a shorter data length than that for the frequency estimation. In
Fig. 4b
, a legend, for example, ‘64pt FFT (64 ms)+32 ms’ denotes the phase estimation with 32 ms of data when frequency estimation is performed with
L_{T}
=64, and it can be found that a shorter data length for phase estimation results in a smaller phase estimation with the proposed technique and with the curve fitting technique is much better than with other techniques for moderate and strong signals.
Fig. 4c
shows the computational complexity (i.e., the number of complex multiplications) of the frequency estimation techniques, such as
N_{F}
point FFTbased technique (
N_{F}
=64), serial Doppler frequency search techique with
T_{CO}
=20 ms and search step size 10
^{3}
/
N_{F}
Hz, conventional technique Eq. (28) with
N
=20, and the proposed technique with
L_{T}
=16 ms and 32 ms, discussed in this section. When the sampling rate
f_{s}
=2
R_{c}
=2×1.023 MHz, the computational costs for the frequency estimation techniques can be expressed as
where
M_{FFT}
,
M_{SS}
,
M_{conv}
, and
M_{FRI}
represent the computational complexity of the
N_{F}
point FFT based technique, the serial search technique with coherent integration length
T_{co}
testing
N_{DF}
=
N_{F}
/2=32 Doppler frequency hypotheses tested with a smaller coherent integration length
T
_{1}
(i.e., 500 Hz), the conventional technique Eq. (28) with N=20, and the proposed technique, respectively. Notice that computational complexity for computing the arctangent function is not counted in the express ions in Eq. (30), and that, for high
C
=
N
_{0}
,
M_{FFT}
=64×(2046+
log
_{2}
64)=131328,
M_{SS}
=20×2046×32=1350360,
M_{conv}
=20×2046×
log
_{2}
2046=450062,
M_{FRI}
=2046×16+5×[4×11
^{2}
+22×6
^{2}
]×6+[4×15
^{2}
+22×2
^{2}
]×2=72992 for
L_{T}
=16, where
N_{i}
=5, and
M_{FRI}
=2046×32+2×[4×21
^{2}
+22×12
^{2}
]×12+[4×31
^{2}
+22×2
^{2}
]×2=191704 for
L_{T}
=32, where
N_{i}
=2, which match the simulation results for
C
/
N
_{0}
=45 dBHz shown in
Fig. 4c
.
Obviously, the total complex multiplications of the curve fitting technique are larger than the serial search technique due to the additional curve fitting process. Furthermore, the three figures in
Fig. 4
shows that the proposed technique can produce much more accurate Doppler frequency and phase estimates than other techniques and that the computational complexity of the proposed technique is only around 2×10
^{5}
multiplications for
L_{T}
ms of data, while other techniques have multiple times of estimation error and computational cost in the GPS pullin state.
6. CONCLUSIONS
It has been found that a reverse directional finite rate of innovation technique can be applied to GPS pullin search and improves the accuracy and computational cost of the pullin search in comparison to the conventional techniques. Analytical expressions are derived to obtain the theoretical estimates of performance improvement and computational complexity with the proposed technique. The simulations have demonstrated that the proposed technique has a performance improvement when SNR is moderate or high, and that the improvement becomes significant as the resolution of frequency estimation of the proposed technique increases linearly as SNR increases, which is not possible with the conventional FFTbased frequency estimation techniques. In addition to the improvement of the frequency estimation, the proposed technique contributes to the accuracy in the phase estimation, and, therefore, improves the phase tracking performance of a receiver.
Acknowledgements
This work is supported by the Electronic Warfare Research Center (EWRC) grant funded by the Agency for Defense Development (ADD).
BIO
SeungHyun Kong received a B.S.E.E. from the Sogang University, Korea, in 1992, an M.S.E.E. from the Polytechnic University, New York, in 1994, and a Ph.D. degree in Aeronautics and Astronautics from Stanford University, CA, in 2006. From 1997 to 2004, he worked with Samsung Electronics Inc. and Nexpilot Inc., both in Korea, where his research focus was on 2G CDMA and 3G UMTS PHY and mobile positioning technologies. In 2006, he was involved with the development of hybrid positioning technology using wireless location signature and Assisted GNSS at Polaris Wireles, Inc. and from 2007 to 2009, he was a researcher at Qualcomm Research Center, San Diego, CA, where his R&D focus was on indoor location technologies and advanced GNSS technologies. Since 2010, he as been an Assistant Professor at the CCS Graduate School for Green Transportation in the Korea dvanced Institute of Science and Technology (KAIST). His research interests include superresolution signal processing, detection and estimation for navigation systems, and Assisted GNSS in wireless communication systems.
Kyungwoo Yoo received a B.S.E.E from Inha University, Korea, in 2012, and a M.S. degree from the CCS Graduate School for Green Transportation in the Korea Advanced Institute of Science and Technology (KAIST), Korea, in 2014. He is currently pursuing a Ph.D. degree with the same institution. His research interests include antijamming technique and signal processing for GNSS and navigation system.
Blu T.
,
Dragotti P.L.
,
Vetterli M.
,
Marziliano P.
,
Coulot L.
2008
Sparse Sampling of Signal Innovations: Theory, Algorithms and Performance Bounds
IEEE Signal Process Mag.
http://dx.doi.org/10.1109/MSP.2007.914998
25
31 
40
DOI : 10.1109/MSP.2007.914998
Borre K.
,
Akos D.
,
Bertelsen N.
,
Rinder P.
,
Jensen S.
2007
A SoftwareDefined GPS and Galileo Receiver: A SingleFrequency Approach
Birkhauser
Boston, MA
Chan T. F.
1982
An Improved Algorithm for Computing the Singular Value Decomposition
ACM Trans. Math. Softw.
http://dx.doi.org/10.1145/355984.355990
8
72 
83
DOI : 10.1145/355984.355990
Goiser M. J.
,
Berger G. L.
1996
Synchronizing a digital GPS receiver
MELCON ’96, Proc. of the 8th Mediterranean Electrotechnical Conference
May 1996
http://dx.doi.org/10.1109/MELCON.1996.551387
1043 
1046
DOI : 10.1109/MELCON.1996.551387
Golub G. H.
,
van Loan C. F.
1996
Matrix Computations
3rd ed
The Johns Hopkins Univ. Press
Baltimore
Hayes M. H.
1996
Statistical Digital Signal Processing and Modeling
John Wiley & Sons
New York
Irsigler M.
,
Eissfeller B.
2003
Comparison of Multipath Mitigation Techniques with Consideration of Future Signal Structures
Proc. of ION GPS/GNSS 2003
Portland, OR
Sep 2003
Kelley C.
,
Cheng J.
,
Barnes J.
2002
Open source software for learning about GPS
Proc. of ION GPS/GNSS 2002
Portland, OR
Sep 2002
Kim B.
,
Muchkaev A.
,
Kong S.H.
2011
SAR Image Processing Using SuperResolution Spectral Estimation with Annihilating Filter
2011 3rd International AsiaPacific Conference on Synthetic Aperture Radar, APSAR
Seoul, Korea
14 Sep 2011
Kong S.H.
,
Kim B.
2013
TwoDimensional Compressed Correlator for Fast PN Code Acquisition
IEEE Trans. Wirel. Commun.
http://dx.doi.org/10.1109/TWC.2013.092313.130407
12
5859 
5867
DOI : 10.1109/TWC.2013.092313.130407
Kransner N. F.
1998
GPS receiver and method for processing GPS signals
US Patent, PN 5,781,156
Parkinson B.
,
Spilker J.
,
Axelrad P.
,
Enge P.
1996
Global Positioning System: Theory and Applications
American Institute of Aeronautics and Astronautics
Washington, DC
Petovello M.
,
O&Driscoll C.
,
Lachapelle G.
2008
Weak Signal Carrier Tracking Using Extended Coherent Integration with an UltraTight GNSS/IMU Receiver
Proc. of ENC 2008
Toulouse
Apr 2008
Pisarenko V. F.
1973
The Retrieval of Harmonics from a Covariance Function
Geophys. J. R. Asfr. Soc.
347 
366
Psiaki M. L.
,
Jung H.
2002
Extended Kalman Filter Methods for Tracking Weak GPS Signals
Proc. of ION GPS/ITM 2002
Portland, OR
Sep 2002
2539 
2553
Sagiraju P. K.
,
Akopian D.
,
Valio H.
2006
Fine Frequency Estimation in Weak Signals for GPS Receivers
Proc. of ION GPS/NTM 2006
Monterey, CA
Jan 2006
2524 
2533
Tsui J. B. Y.
2005
Fundamentals of Global Positioning Receivers: A Software Approach
2nd Ed
John Willey & Sons
New York
Tufts D. W.
,
Kumaresan R.
1982
Estimation of frequencies of multiple sinusoids: Making linear prediction perform like maximum likelihood
Proc. IEEE
http://dx.doi.org/10.1109/PROC.1982.12428
975 
989
DOI : 10.1109/PROC.1982.12428
van Diggelen F.
2009
AGPS: Assisted GPS, GNSS, and SBAS
Artech House
Norwood, MA
Vetterli M.
,
Marziliano P.
,
Blu T.
2002
Sampling Signals with Finite Rate of Innovation
IEEE Trans. Signal Process
http://dx.doi.org/10.1109/TSP.2002.1003065
50
1417 
1428
DOI : 10.1109/TSP.2002.1003065
Zeng D.
,
Li J.
2010
GPS Signal Fine Acquisition Algorithm, Information Science and Engineering (ICISE)
2010 2nd International Conference on IEEE
Dec 2010
http://dx.doi.org/10.1109/ICISE.2010.5691791
3729 
3732
DOI : 10.1109/ICISE.2010.5691791