SNAPSHOT GNSS RECEIVER AND METHOD USING SUPER-LONG COHERENT INTEGRATION AND FRACTIONAL FOURIER TRANSFORM

Abstract

Snapshot receiver that comprises a correlator module for correlating incoming GNSS signals and a transform module that transforms resulting correlated outputs using a Fractional Fourier Transform (FrFT) process, to thereby compensate for weak and dynamic signals. The output of the transform module is an estimated Doppler rate with high accuracy. The estimated Doppler rate is passed to a super-resolution-measurement (SRM) module, which outputs error values of the estimated Doppler rate and a pseudorange (measured distance-to-satellite) via a phase dead reckoning (DR) calculation for the snapshot receiver. The pseudorange and error values are passed to a navigator module that determines position information based on those inputs. In some embodiments, the SRM module comprises a maximum likelihood estimator (MLE).

Claims

1. A snapshot receiver for determining position information, said snapshot receiver comprising: a receiving module for receiving a GNSS signal from a satellite; a processing module comprising: a correlator module for correlating said GNSS signal using time-domain integration of a plurality of samples of said GNSS signal, to thereby produce a correlated signal; a transform module for applying a Fractional Fourier Transform (FrFT) process to said correlated signal to thereby produce an estimated Doppler rate; and a super-resolution measurement (SRM) module for receiving the estimated Doppler rate and for determining error values of said estimated Doppler rate and for determining a measured distance between said snapshot receiver and said satellite; and a navigator module for determining said position information based on said error values and said measured distance.

2. The snapshot receiver according to claim 1, wherein said SRM module comprises a phase dead reckoning module for determining said error values and said measured distance.

3. The snapshot receiver according to claim 1, wherein said SRM module comprises a maximum likelihood estimator (MLE) for determining said error values and said measured distance, wherein said error values and said measured distance are values that maximize a probability of convergence between said estimated Doppler rate and a modelled reference Doppler rate.

4. The snapshot receiver according to claim 1, wherein said error values comprise a Doppler rate error, a Doppler frequency error, a carrier phase error, and a code phase error of said estimated Doppler rate.

5. The snapshot receiver according to claim 1, wherein said correlator comprises a plurality of correlation channels.

6. The snapshot receiver according to claim 1, wherein said correlator applies super-long coherent integration (S-LCI).

7. The snapshot receiver according to claim 1, wherein: said transform module comprises a down-sampling submodule; said correlated signal comprises fast-time correlator outputs; said down-sampling submodule down-samples said fast-time correlator outputs to thereby produce slow-time correlator outputs; and said FrFT is applied to said slow-time correlator outputs.

8. The snapshot receiver according to claim 1, wherein said estimated Doppler rate is passed through a pre-processing module before being passed to the SRM module.

9. The snapshot receiver according to claim 3, wherein said MLE applies a gradient descent optimization process when determining said error values and said measured distance.

10. A method for determining position information, said method comprising the steps of: receiving a GNSS signal from a satellite; correlating said GNSS signal using time-domain integration of a plurality of samples of said GNSS signal, to thereby produce a correlated signal; applying a Fractional Fourier Transform (FrFT) process to said correlated signal to thereby produce an estimated Doppler rate, wherein said estimated Doppler rate has an associated probability distribution; and determining error values of said estimated Doppler rate; based on the estimated Doppler rate, determining a measured distance between said snapshot receiver and said satellite; and determining said position information based on said error values and said measured distance.

11. The method according to claim 10, wherein determining said error values and said measured distance uses a dead-reckoning phase based on the super-resolution (SR) Doppler estimation.

12. The method according to claim 10, wherein determining said error values and said measured distance uses maximum likelihood estimation (MLE), such that said error values and said measured distance are values that maximize a probability of convergence between said estimated Doppler rate, an estimated Doppler frequency, an estimated carrier phase, an estimated code phase, an estimated signal amplitude, and a modelled reference Doppler rate, a modelled Doppler frequency, a modelled carrier phase, a modelled code phase, and a modelled signal amplitude, respectively.

13. The method according to claim 10, wherein said error values comprise a Doppler rate error, a Doppler frequency error, a carrier phase error, and a code phase error of said estimated Doppler rate.

14. The method according to claim 10, wherein a plurality of correlation channels are used in said step of correlating.

15. The method according to claim 10, wherein said time-domain integration is super-long coherent integration (S-LCI).

16. The method according to claim 10, wherein said correlated signal comprises fast-time correlator outputs, and wherein said method further comprises down-sampling said fast-time correlator outputs to thereby produce slow-time correlator outputs before applying said FrFT.

17. The method according to claim 10, further comprising pre-processing said estimated Doppler rate before determining said error values and said measured distance.

18. The method according to claim 12, wherein said MLE applies a gradient descent optimization process when determining said error values and said measured distance.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

[0028] The present invention will now be described by reference to the following figures, in which identical reference numerals refer to identical elements and in which:

[0029] FIG. 1 is a schematic diagram of a snapshot receiver according to one aspect of the invention;

[0030] FIG. 2 is a process diagram of a phase-dead-reckoning- and super-resolution (SR)-Doppler-based tracking process for use by a receiver according to an embodiment of the invention;

[0031] FIG. 3 is a process diagram of a maximum likelihood estimator (MLE)-based tracking process for use by a receiver according to an embodiment of the invention;

[0032] FIG. 4 is a diagram showing Sigma correlation according to the prior art;

[0033] FIG. 5 is a schematic diagram showing conventional Fourier transforms and fractional Fourier transforms according to the prior art;

[0034] FIG. 6 is a diagram of a super-long coherent integration process with partially matched filtering-FrFT according to the present invention;

[0035] FIG. 7 is a diagram of processing stages in an MLE based implementation of the invention;

[0036] FIG. 8 is a diagram of a joint Doppler effect/time-of-arrival (TOA) MLE in the baseband, according to an embodiment of the invention;

[0037] FIG. 9 is a diagram of pre-processing processes according to an embodiment of the present invention;

[0038] FIG. 10 is a diagram of matching filtering process up correlation peaks about the code domain (CD), frequency domain (FD), and fractional frequency domain (FrFD); and

[0039] FIG. 11 is a flow chart detailing a method according to one aspect of the invention.

DETAILED DESCRIPTION

[0040] To better understand the present invention, the reader is directed to the listing of citations at the end of this description. For ease of reference, these citations and references have been referred to by their listing number throughout this document. The contents of the citations in the list at the end of this description are hereby incorporated by reference herein in their entirety.

[0041] The present invention discloses a snapshot receiver that correlates incoming GNSS signals and transforms resulting correlated outputs using a Fractional Fourier Transform (FrFT). An estimated Doppler rate is provided as the transformed output, and values associated with the estimated Doppler rate are determined by a super-resolution measurement (SRM) module and used to determine position information of the receiver.

[0042] Referring now to FIG. 1, a snapshot receiver 10 according to an aspect of the invention is shown. (As should be clear, FIG. 1 is not to scale.) A receiving module 20 receives a GNSS signal from a satellite 30. The signal is passed to a processing module 40. The processing module 40 comprises a correlator module 50. The correlator module 50 performs time-domain integration of a plurality of samples of the signal, to thereby produce a correlated signal. The processing module 40 further comprises a transform module 60, which applies a Fractional Fourier Transform (FrFT) to the correlated signal to reduce noise. The output of the transform module 60 is an estimated Doppler rate and Doppler shift of the signal. The processing module 40 further comprises a super-resolution measurement (SRM) module 70 that determines error values of the estimated Doppler rate, Doppler shift, carrier phase, and a pseudorange of the signal (i.e., a distance between the satellite 30 and the receiver 10, as measured at the receiver 10), based on the estimated Doppler rate. Methods of determining the pseudorange of a receiver are well-known in the art. The error values and pseudorange/measured distance are then passed from the processing module 40 to a navigator module 80. The navigator module 80 determines position information based on the error values and pseudorange.

[0043] As noted above, the correlator module 50 performs time-domain integration using a plurality of samples of the incoming GNSS signal. As previously described, many such time-domain integration approaches are known in the art, including long coherent integration (LCI) and super-long coherent integration (S-LCI). In some embodiments, the correlator module 50 uses LCI. In other embodiments, the correlator module 50 uses S-LCI. In still other embodiments, the correlator module 50 uses a combination of LCI and S-LCI. However, in some embodiments, other time-domain integration methods and approaches may be used.

[0044] Details of the application of the FrFT to the correlated signal by the transform module 60 will also be described in detail below.

[0045] In some embodiments, the SRM module 70 comprises a phase dead reckoning (phase DR) module that determines super-resolution Doppler and phase values. Further description of such a phase DR module is provided below.

[0046] In other embodiments, the SRM module 70 comprises a maximum likelihood estimator (MLE). MLEs are well-known statistical methods for determining values that maximize the likelihood of a specific result in a particular probability distribution. In this case, the MLE determines error values and a pseudorange that maximize the probability of convergence between estimated rates and modelled reference rates. That is, the error values and pseudorange determined by the MLE comprise values that maximize the probability of convergence between the estimated Doppler rate, an estimated Doppler frequency, an estimated carrier phase, an estimated code phase, and an estimated signal amplitude, and the modelled Doppler rate, the modelled Doppler frequency, the modelled carrier phase, the modelled code phase, and the modelled signal amplitude, respectively. Further description of MLE-based implementations of the SRM module 70 is provided below.

[0047] One embodiment of the snapshot receiver 10 is shown in FIG. 2 (i.e., as a block diagram). This embodiment uses a phase-DR-based SRM module 70. The block denoted Front-end in this Figure corresponds to the receiving module 20 of FIG. 1. The block denoted Proposed L5Q baseband signal processing corresponds to the correlator module 50 and transform module 60 of FIG. 1. The block denoted SRM: Doppler frequency and pseudorange corresponds to the SRM module 70 of FIG. 1.

[0048] Another embodiment of the snapshot receiver 10 is shown in FIG. 3. This embodiment uses an MLE-based SRM module 70. As can be seen, the various modules of FIG. 1 can, in some embodiments, be embodied in separate and/or multiple elements of the receiver. For example, in this figure, the receiving module 20 is not explicitly shown. As a further example, as can be clearly seen, in this embodiment, the correlator module 50 comprises the initialization block on the left as well as the down-sampling submodule denoted as PMF (Down-sampling) in block 2.2. In this embodiment, the correlator module 50 first produces fast-time outputs. As will be further described below, the fast-time outputs are down-sampled by the down-sampling submodule to thereby produce slow-time outputs for use in S-LCI.

[0049] Further, as can be seen from FIG. 3, in some embodiments, the receiver 10 comprises a preprocessing module that preprocesses the estimated Doppler rate and Doppler shift after the FrFT is applied and before the SRM module 70 receives the transformed signal. Further, in some embodiments, the receiver 10 also comprises replication and feedback modules and submodules that create code/carrier local replicas based on the output of the SRM module 70. The code/carrier local replicas are, in some embodiments, fed back to the correlator module 50.

[0050] As well, as should be understood, the diagram shown in FIG. 3 details the receiver 10 in a single channel. As is well known, receivers may have a plurality of channels active at any time. As would be clear to the person skilled in the art, in the embodiment shown in FIG. 3, signals in each channel of the receiver 10 would be processed in an identical manner. The outputs of each channel are passed to the navigator module 80, which determines positioning information. Again, the outputs of each channel are the applicable pseudorange, carrier phase, and error values.

[0051] Again, the error values comprise the Doppler rate error, Doppler frequency error, the carrier phase error, the code phase error, and the signal amplitude error.

[0052] Similarly, in some embodiments, the correlator module 50 comprises a plurality of correlation channels. Each of the plurality of correlation channels performs a separate time-domain integration process on samples of the signal. Outputs of each separate channel are then merged to form the correlated signal for further processing within the processing module 40.

[0053] In some embodiments, the MLE applies the well-known gradient descent optimization process. Of course, other optimization processes are possible depending on the embodiment and implementation.

Mathematical Implementations and Approaches

[0054] A dynamic GNSS signal model may be used as a basis for understanding the various processes and description that follow. As should be understood, this model is intended to mimic the behaviour of an incoming GNSS signal. In particular, the model below models the intermediate frequency (IF) GPS L1 C/A signal, as follows:

[00001]$s\left(t\right)=2\left(1+?a\right)d\left(t-{?}_{\mathrm{PRN}}+{?}_{\mathrm{dyn}}\left(t\right)\right)C\left(t-{?}_{\mathrm{PRN}}+{?}_{\mathrm{dyn}}\left(t\right)\right)?\cos \left(2?f_{I}t+\left({?}_I+{?}_{\mathrm{dyn}}\left(t\right)\right)\right)$

where

[00002]${?}_{\mathrm{dyn}}\left(t\right)\stackrel{?}{=}\left(\mathrm{vt}+\frac{?}{2}{t}^2\right)c^{-1}$${?}_{\mathrm{dyn}}\left(t\right)\stackrel{?}{=}2?f_r{?}_{\mathrm{dyn}}\left(t\right)\stackrel{?}{=}2?\left(\mathrm{ft}+\frac{?}{2}{t}^2\right)$$\begin{array}{cc}f\stackrel{?}{=}{f}_r{c}^{-1}v,& ?\stackrel{?}{=}{f}_r{c}^{-1}?\end{array}$${?}_I\stackrel{?}{=}{?}_0-{?}_{\mathrm{LO}}-2?\left(f_r+f\right){?}_{\mathrm{PRN}}$

and where ?? is the amplitude error of the baseband signal; t is the time variable; D(t) is the navigation data; C(t) is the spreading code; f.sub.i is the IF; ?.sub.i is the initial carrier phase; ?.sub.PRN is the propagation delay about the distance between the satellite and the user/receiver at the time of signal emission; ?.sub.dyn(t) and ?.sub.dyn(t) are the respective propagation delay residual and carrier phase residual determined by line-of-sight (LOS) dynamics and signal propagation time; ? is the dynamic term related to the velocity of the changing distance while ? is the counterpart related to acceleration; c is the speed of light; f is the initial Doppler frequency shift; f.sub.r is the radio frequency; ? is the Doppler rate; ?.sub.0 is the initial carrier phase; and ?.sub.Lo is the carrier phase caused by the frequency-mixing process.

[0055] Note as well that, in traditional GNSS baseband signal processing, the higher-order dynamic term ? related to the Doppler rate in ?.sub.dyn(t) and ?.sub.dyn(t) are not usually modelled. As a result, the local signal cannot be adequately replicated as the real signal, as higher-order dynamics cannot be considered. This has a greater impact for longer correlation periods: i.e., the longer the coherent integration interval is, the more apparent the discrepancy between the real signal and the local replica becomes.

[0056] Mathematical implementation details that follow referencing this signal model as needed. As should be understood, the mathematical implementations that follow are simply possible implementations. As described above and as known in the art, multiple approaches and methods may be suitable for each step.

Correlation and Digital Fractional Fourier Transform (FrFT)

[0057] Time-domain integration (correlation) processes are well known. In general, the correlation processes comprise determining average values of a signal over preset time intervals. A diagram of traditional correlation processes in traditional GNSS tracking is shown in FIG. 4. An initial incoming IF signal is locally replicated. A plurality of samples and their corresponding local replicas are correlated together. It should be noted that FIG. 4 depicts three (3) samples; as would be understood by the person skilled in the art, the use of 3 samples in this figure is merely for visual simplicity. In practical cases, the actual sample number may be well into the thousands, if not higher. The averaged samples are summed to provide a single correlated output.

[0058] As is known in the art, LCI and S-LCI involve the use of large numbers of samples. In some embodiments of the present invention, the initial samples are further down-sampled to increase the number available for integration.

[0059] A diagram of FrFT, as compared to standard Fourier transforms (FTs), is shown in FIG. 5. The FrFT is a form of FT that is equivalent to the FT process when p is equal to 1. Assuming that the signal ?(t) satisfies Dirichlet conditions, the FrFT for ?(t) can be defined as an integral transform as follows:

[00003]$?\left(p,u\right)={?}_{-?}^{?}?\left(t\right)K_p\left(t,u\right)\mathrm{dt}$

with

[00004]$K_p\left(t,u\right)=\{\begin{array}{cc}{C}_{?}{e}^{j?\left[\left(t^2+u^2\right)\cot ?-2\mathrm{ut}\csc ?\right]}& ??n?\\ ?\left(t-u\right),& \begin{array}{cc}?=2?n,& n??\end{array}\\ ?\left(t+u\right),& ?=\left(2n+1\right)?\end{array}$

where K.sub.p(t, u) represents the transformation kernel of FRFT;

[00005]$C_{?}=\sqrt{1-j\cot ?};p$

is the FrFT order; ? denotes a rotation angle from the time-frequency plane in terms of the ordinary FT to an extended plane; and where ? satisfies

[00006]$?=\frac{?}{2}p$

FIG. 5 shows that the ordinary FT transfers a dynamic signal from the time to frequency domain, but the signal power cannot be concentrated in this normal time-frequency plane. In contrast, after an FrFT, i.e., a rotation with an angle ? for this plane, the dynamic signal power is concentrated in the FrFD (i.e., the fractional Fourier domain).

[0060] FIG. 6 illustrates the work process of the correlators based on the digital FrFT in the GNSS baseband, in some embodiments of the invention. (That is, in embodiments using S-LCI and multiple correlation channels.) As well, in the illustrated embodiment, a partially matched filter (PMF) technique is used for the digital FrFT to reduce computational burden. This is similar to the PMF-FFT (i.e., partially matched filter-fast Fourier transform) technique described in reference [9]. The input of the digital PMF-FrFT is slow-time correlations (the larger unfilled circles).

[0061] As can be seen, in the embodiment illustrated in FIG. 6, three types of correlator outputs are formed as follows: [0062] First, a sequence of fast-time correlations (the smaller unfilled circles without extending line segments) is formed in the time domain (TD) after an integrating and dumping (I&D) module (a submodule of the correlator module 50) processes the incoming IF samples. [0063] Second, a sequence of slow-time correlations is formed in the TD when the fast-time correlation sequence passes through the N.sub.s-point summation operators. One summation process represents one PMF, and is a down-sampling process (i.e., performed by the down-sampling submodule). As is typical in GNSS tracking application, N.sub.s is much smaller than N.sub.t, and the values satisfy N.sub.t=N.sub.sP, where P is the number of IF samples within one PMF and N.sub.s is the number of PMFs. [0064] Third, a sequence of S-LCI correlations is formed in the FrFD after the PMF-FrFT processes the slow-time correlator outputs. Thus, compared to the original digital FrFT (i.e., without PMF processes), which has a computational load of O(N.sub.t log.sub.2 N.sub.t), the digital PMF-FrFT process has a significantly lighter computational load of O(N.sub.s log.sub.2 N.sub.s).

[0065] Closed-form models for the S-LCI correlations in the FrFD can then be derived as follows.

Closed-Form Models for the S-LCI Correlations in the FrFD

[0066] The expression for the sampled and digitalized s(t) in the dynamic GNSS model can be written as

[00007]$\begin{array}{c}s\left(n_F{T}_F\right)=D\left(n_F{T}_F-{?}_{\mathrm{PRN}}+{?}_{\mathrm{dyn}}\left(n_F{T}_F\right)\right)\\ ?2\left(1+?a\right)C\left(n_F{T}_F-{?}_{\mathrm{PRN}}+{?}_{\mathrm{dyn}}\left(n_F{T}_F\right)\right)\\ ?\cos \left(2?f_I{n}_F{T}_F+\left({?}_I+2?{\mathrm{fn}}_F{T}_F+\mathrm{??}{n}_F^2{T}_F^2\right)\right)\end{array}$

where n.sub.F is the index of incoming IF samples; T.sub.F is the IF sampling interval; s[n.sub.F] indicates the sequence s(n.sub.FT.sub.F); and, letting ?[n.sub.F, m.sub.c, m.sub.f] be the notation of the sequence of local replicas and ?.sup.+[n.sub.S, m.sub.c, m.sub.f] be the sequence of correlator outputs in the time-code-frequency domain (where n.sub.s is the index of the outputs, m.sub.c is the index of local code replica sequence about code chip offsets, and m.sub.f is the index of local carrier replica sequence about Doppler frequency offsets), the correlator outputs may be obtained by summing the outputs of the I&D operator between s [n.sub.F] and ?[n.sub.F, m.sub.c, m.sub.f], the expression of which is known as the cross-ambiguity function (CAF). The outputs of the I&D operator are herein referred to as the fast-time correlation sequence, while ?.sup.+[n.sub.S, m.sub.c, m.sub.f] is herein referred to as the slow-time correlation sequence.

[0067] The digital PMF-FrFT of the slow-time correlation sequence then gives

[00008]$x_b^{+}\left[m_c,m_f,m_{f_p},m_p\right]=\frac{\sqrt{1-j\cot \left(?\left[m_p\right]\right)}}{2\sqrt{N_S}}?\underset{n_s=-N_s}{\overset{N_s}{.Math.}}{?}^{+}\left[\frac{n_S}{2\sqrt{N_S}},m_c,m_f\right]e^{\left(j?\frac{\cot ?\left[m_p\right]m_u^2}{4N_S}-j?\frac{\csc ?\left[m_p\right]m_u{n}_S}{2N_S}+j?\frac{\cot ?\left[m_p\right]n_S^2}{4N_S}\right)}\mathrm{for}{m}_c=-\frac{M_c-1}{2},.Math.,\frac{M_c-1}{2},m_f=-\frac{M_f-1}{2},.Math.,m_{f_p}=-\frac{M_{f_p}-1}{2},.Math.,\frac{M_{f_p}-1}{2},m_p=-\frac{M_p-1}{2},.Math.,\frac{M_p-1}{2},m_c,m_f,m_{f_p},m_p??$

with

[00009]$m_{f_p}\stackrel{?}{=}\frac{m_u}{2\sqrt{N_S}}$

where m.sub.f.sub.p and m.sub.p represent the indexes of the match filters for frequency shift and FrFT order (related to Doppler rate) in the FrFD, respectively; ?p is the discrete search step for the FRFT order; and M.sub.c, M.sub.f, M.sub.f.sub.p and M.sub.p are respectively the numbers of match filters in the CD, FD, FrFD (matching frequency shift), and FrFD (matching frequency rate).

[0068] Then, the matched frequency shift and rate in the FrFD is computed as

[00010]$f_{n_c}\left[m_{f_p},m_p\right]=\frac{m_{f_p}}{N_S{T}_S\sin \left(?\left[m_p\right]\right)}$$?\left[m_p\right]=-\frac{\cot \left(?\left[m_p\right]\right)}{N_S{T}_S^2}$

where

[00011]$?\left[m_p\right]=\frac{?}{2}\left(1+m_p?p\right)$

and is the matched rotation angle.

[0069] The equation for ?.sub.b.sup.+ can then be expressed as

[00012]$\begin{array}{c}{\left[{?}_b^{+}\right]}_{m_c,m_f,m_{f_p},m_p}=\sqrt{1-j\cot \left(?\left[m_p\right]\right)}{e}^{j?\frac{f_{n_C}^2\left(?f\left[m_f\right],??\right)\cot \left(?\left[m_p\right]\right)}{N_S}}\\ ?e^{j{?}_{n_C}\left(??,?f\left[m_f\right],?\right)}.Math.{\left[{?}_b^{+}\right]}_{m_c,m_f,m_{f_p},m_p}.Math.\end{array}$

with

[00013]${?}_{n_C}\left(\mathrm{??},?f\left[m_f\right],?\right)=\mathrm{??}+\mathrm{??}f\left[m_f\right]N_S{T}_S+\frac{?}{4}?N_S^2{T}_S^2$$f_{n_C}\left(?f\left[m_f\right],?\right)=?f\left[m_f\right]+\frac{?}{2}{N}_S{T}_S$$\begin{array}{c}.Math.{\left[{?}_b^{+}\right]}_{m_c,m_f,m_{f_p},m_p}.Math.=R_{.Math.}\left({{\mathrm{??}}_{b,m}}_c\left({\mathrm{??}}_b\right),?f\left[m_f\right],?\right)\\ ?\left(1+?a\right){?}_1^f\left(?f\left[m_f\right],?\right)\\ ?{?}_2^f\left({\mathrm{?u}}_{m_{f_p}},m_p,?f\left[m_f\right],?\right)\\ ?{?}_2^{?}\left({\mathrm{??}}_{m_p}\left(?\right)\right)\end{array}$${\mathrm{??}}_{b,m_c}\left({\mathrm{??}}_b\right)={\mathrm{??}}_b-\frac{m_c{d}_c}{2f_c}$$?u_{m_{f_p},m_p}\left(?f\left[m_f\right]\right)=?f\left[m_f\right]-f_{n_C}\left[m_{f_p},m_p\right]$${\mathrm{??}}_{m_p}\left(?\right)=?-?\left[m_p\right]$$\mathrm{??}={?}_I-{\widehat{?}}_I,{\mathrm{??}}_b={?}_{\mathrm{PRN}}-{\widehat{?}}_{\mathrm{PRN}},?f\left[m_f\right]=f-\left(\hat{f}+m_f?f_{\mathrm{nco}}\right)$

where {circumflex over (?)}.sub.PRN, {circumflex over (f)} and {circumflex over (?)}.sub.1 are, respectively, the estimated code delay in seconds, Doppler frequency in Hz, and carrier phase in radians from the NCO processes; ??.sub.b and ?? are the expressions of initial code phase error and initial carrier phase, respectively, ?f[m.sub.f] is the expression of initial frequency error about the NCO frequency; f.sub.c is the code frequency; f.sub.n.sub.c(?f[m.sub.f], ?) and ?.sub.n.sub.c(??, ?f [m.sub.f], ?) are the closed-form expressions for the averaging frequency error and carrier phase over T.sub.t; and |?.sub.b.sup.+| is the closed-form correlation amplitude.

[0070] In summary, implementing digital PMF-FrFT to the slow-time correlation sequence gives the S-LCI correlation sequence varied with both matched frequency shift and rate in the FrFD. Thus, a 2D set of the S-LCI correlations in the FrFD is formed and can be written in a matrix as:

[00014]$\left[\begin{array}{ccc}{?}_b^{+}\left[m_c,m_f,m_{f_p}-\frac{M_{f_p}-1}{2},m_p-\frac{M_p-1}{2}\right]& .Math.& {?}_b^{+}\left[m_c,m_f,m_{f_p}+\frac{M_{f_p}-1}{2},m_p-\frac{M_p-1}{2}\right]\\ .Math.& ?& .Math.\\ {?}_b^{+}\left[m_c,m_f,m_{f_p}-\frac{M_{f_p}-1}{2},m_p+\frac{M_p-1}{2}\right]& .Math.& {?}_b^{+}\left[m_c,m_f,m_{f_p}+\frac{M_{f_p}-1}{2},m_p+\frac{M_p-1}{2}\right]\end{array}\right]$

[0071] The S-LCI correlation sequences in the FrFD, i.e., after FrFT, can then be used by the SRM module 70 of the receiver 10.

SRM Module 70

[0072] The SRM module returns high-resolution error values based on the output(s) of the transform module. As described above, the SRM module 70, in some embodiments, comprises a phase dead reckoning (phase DR) module. In other embodiments, the SRM module 70 comprises an MLE.

Phase DR Implementation

[0073] In implementations using super-resolution measurements, the closed-form models of the matched frequency shift error and rate can be written as

[00015]$?f_{n_C}\left[{\hat{m}}_u^{},{\hat{m}}_p\right]=\frac{{\hat{m}}_u^{}}{N_S{T}_S\sin \left(?\left[{\hat{m}}_p\right]\right)}$${?}_{n_C}\left[{\hat{m}}_p\right]=-\frac{\cot \left(?\left[{\hat{m}}_p\right]\right)}{N_S{T}_S^2}$

where

[00016]${\stackrel{?}{m}}_u^{}\stackrel{?}{=}{\stackrel{?}{m}}_u/2/\sqrt{N_S};?\left[{\stackrel{?}{m}}_p\right]=\frac{?}{2}\left(1+{\stackrel{?}{m}}_p\mathrm{?p}\right);$

? is the matched rotation angle; and the subscript n.sub.C denotes the epoch index for the S-LCI correlations in the FrFD.

[0074] A frequency discriminator is used to reduce the estimated frequency error further, as

[00017]$?{\tilde{f}}_{\mathrm{nC},\mathrm{disc}}=\frac{\left(?f_{n_C}\left[{\hat{m}}_u^{}-1,{\hat{m}}_p\right]P_l+?f_{n_C}\left[{\hat{m}}_u^{}.Math.{\hat{m}}_p\right]p_n+?f_{n_C}\left[{\hat{m}}_u^{}+1,{\hat{m}}_p\right]P_h\right)}{P_l+P_n+P_h}$$\mathrm{with}$$P_l=.Math.{\widetilde{?}}_b^{+}\left[m_u^{}-1,{\hat{m}}_p\right]{.Math.}^2,P_n=.Math.{\widetilde{?}}_b^{+}\left[m_u^{},{\hat{m}}_p\right]{.Math.}^2,$$P_h=.Math.{\widetilde{?}}_b^{+}\left[m_u^{}+l,{\hat{m}}_p\right]{.Math.}^2$

[0075] The matched code phase offset is modelled as

[00018]$?{\widehat{?}}_{n_s}\left[{\hat{m}}_c\right]\stackrel{?}{=}{?}_{n_s}-\left({\widehat{?}}_{n_S,0}+\frac{{\hat{m}}_c{d}_c}{2f_c}\right)$

and the carrier phase is discriminated by

[00019]$?{\widehat{?}}_{\mathrm{disc},n_C}=a\tan 2\left[\frac{?\left({\widehat{?}}_b^{+}\left[{\hat{m}}_c,{\hat{m}}_u,{\hat{m}}_p\right]\right)}{?\left({\widehat{?}}_b^{+}\left[{\hat{m}}_c,{\hat{m}}_u,{\hat{m}}_p\right]\right)}\right]-\arg \left[\sqrt{1-j\cot \left(?\left[{\hat{m}}_p\right]\right)}{e}^{j?\frac{?f_{n_C}^2\left(?{\hat{f}}_0,{?}_0\right)\cot \left(?\left[{\hat{m}}_p\right]\right)}{N_S}}\right]$

where ?{circumflex over (f)}.sub.0 and {circumflex over (?)}.sub.0 represent the initial frequency offset and Doppler rate estimates, respectively.

[0076] The carrier NCO frequency, in such an implementation, is then computed as

[00020]${\hat{f}}_{?,\mathrm{NCO},n_C}=f_1+{\hat{f}}_{n_C}$$\mathrm{with}$${\hat{f}}_{n_C}={\hat{f}}_{n_C-1}+{?}_{n_C}\left[{\hat{m}}_p\right]T_t+?f_{n_C}\left[{\hat{m}}_u^{}.Math.{\hat{m}}_p\right]+?{\hat{f}}_{\mathrm{disc},n_C}+?{\widehat{?}}_{\mathrm{disc},n_C}$

where {circumflex over (f)}.sub.n.sub.c and {circumflex over (f)}.sub.n.sub.c.sub.?1 are the estimated Doppler frequency in the current and previous S-LCI epochs, respectively.

[0077] The code NCO frequency is solely assisted by the carrier NCO frequency, to alleviate the multipath interference on the code phase estimates. As such, the code NCO frequency is given by

[00021]${\hat{f}}_{?,?\mathrm{CO},n_C}=f_C+f_C?{\widehat{?}}_{n_C}-\frac{f_c}{f_r}{\hat{f}}_{n_C}$

where ?{circumflex over (?)}.sub.n.sub.c is the estimated code delay offset in the previous epoch and f.sub.r is the radio frequency (RF). Then, the code phase (which provides the pseudorange) is extracted via the code NCO frequency (which, as seen, is solely determined by the super-resolution carrier frequency). This is a code phase DR implementation through the integration of the super-resolution carrier frequency.

Doppler Frequency

[0078] The Doppler frequency measurement in an SRM-based embodiment is the averaging value over the current S-LCI interval, and can be obtained from

[00022]${\tilde{f}}_{\mathrm{dop},n_c}={\hat{f}}_{n_C}$

Then, the matched Doppler rate can be used to estimate the instantaneous Doppler in the FrFD and the current Doppler measurement value. Hence, the instantaneous Doppler at the beginning of the n.sub.Cth epoch is computed as

[00023]$\hat{f}\mathrm{dop},n_{C,}0=\tilde{f}\mathrm{dop},n_C-\frac{1}{2}{?}_{n_C}-1\left[{\hat{m}}_p\right]T_t$

Error values can be determined from the averaging Doppler rate, the instantaneous Doppler shift in the FrFD, and the current Doppler measurement value.

Pseudorange

[0079] For n.sub.s=0, . . . , N.sub.s?1, n.sub.s?, the pseudorange measurements are computed through

[00024]${\widetilde{?}}_{n_S}^i=c\left({\hat{i}}_{n_S}^i-{\widehat{?}}_{\mathrm{rem},n_S}^i-\frac{{\hat{f}}_{?,\mathrm{NCO},n_C}^i}{f_C}{T}_S\right)$

where n.sub.s is the index of the CAF samples since the S-LCI epoch n.sub.C+1; the superscript i is the satellite number; {circumflex over (?)}.sub.rem,n.sub.s.sup.i is the remained code delay in seconds at the previous epoch; {circumflex over (?)}.sub.n.sub.s.sup.i is the local time count for the time-of-arrival (TOA) of satellite i, measured with a local oscillator.

[0080] The code NCO frequency value is unchangeable over the period between epoch n.sub.C+1 and n.sub.C+2; similarly, the carrier frequency is also unchanged over this period.

MLE Implementation

[0081] The following discusses a potential joint Doppler-effect/time-of-arrival (TOA) MLE for embodiments where the SRM module 70 is based on an MLE. A diagram of the process is shown in FIG. 7. As can be seen, 1000 slow-time correlation samples are generated (as described above) with the unchanged carrier and code NCO frequencies. Further, there are five correlator channels to process one satellite signal: besides the conventional early-, prompt-, and late-code delay channels, two extra channels (namely, low- and high-frequency channels) are also provided. Further, in the embodiment illustrated, a pre-processing process takes the place of the conventional code and carrier discriminators to process the correlator outputs from the five channels. Then, the joint Doppler effect/TOA MLE produces the final Doppler frequency error, carrier phase error, and the code phase error, to form the sources of GNSS measurements. As well, in the implementation depicted in FIG. 7, the slow-time interval was taken as 1 ms and a 1s signal sequence was taken as an exemplary value. As should be understood, each of these parameters may be varied as considered suitable. For example, more or fewer than 1000 slow-time correlation samples may be generated, and more or fewer than 5 correlation channels may be used. The person skilled in the art would be able to identify suitable numbers of each. For clarity, it should be noted that the MLE implementation requires at least three channels, as follows: a frequency change channel for the frequency error estimation; a code phase/delay change channel for estimation of the code phase error; and a prompt channel, without frequency or code phase offsets, for estimating the Doppler rate with signal sequences based on the digital FRFT. The number of slow-time correlation samples may be any suitable number, from 1 up to the total number of IF sample points, which may be thousands or more. Of course, the person skilled in the art would understand that using very few slow-time correlation samples might not provide optimal results, while using very many slow-time correlation samples might create unnecessary computational burden.

[0082] Another diagram is of the joint Doppler effect/TOA MLE is shown in FIG. 8. A linearized system model of the MLE can be created, corresponding to this figure. The state vector is given by

[00025]$?\stackrel{?}{=}{\left[{\mathrm{??}}_C,{\mathrm{??}}_{b,c},?f,?,??\right]}^{?}$

with

[00026]${\mathrm{??}}_C\stackrel{?}{=}\frac{\mathrm{??}}{2?},{\mathrm{??}}_{b,c}\stackrel{?}{=}{f}_c{\mathrm{??}}_b$

[0083] The first three states represent initial carrier phase error in cycle, code phase error in chip, and Doppler frequency error in Hz, respectively. ? is the Doppler rate state in Hz/s, and ?? (gain) is the state of the normalized amplitude error.

[0084] The closed-form vector about the measurements of the MLE is defined as

[00027]$f\left(?\right)\stackrel{?}{=}\left[\begin{array}{c}{f}_{A0}({\mathrm{??}}_{b,c},?f,?\\ f_{A1}({\mathrm{??}}_{b,c},?f,?\\ f_{A2}({\mathrm{??}}_{b,c},?f,?\\ f_{A3}({\mathrm{??}}_{b,c},?f,?\\ f_B\left({\mathrm{??}}_c\right)\\ f_C\left(?a\right)\end{array}\right]$$\mathrm{with}$$f_{A0}\left({\mathrm{??}}_{b,c},?f,?\right)\stackrel{?}{=}.Math.{\left[{?}_b^{+}\right]}_{{\hat{m}}_c{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p}.Math.$$f_{A1}\left({\mathrm{??}}_{b,c},?f,?\right)\stackrel{?}{=}{\left[.Math.{\left[{?}_b^{+}\right]}_{{\hat{m}}_c{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_{p-1}}.Math.,.Math.{\left[{?}_b^{+}\right]}_{{\hat{m}}_c{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_{p+1}}.Math.\right]}^T$$f_{A2}\left({\mathrm{??}}_{b,c},?f,?\right)\stackrel{?}{=}{\left[.Math.{\left[{?}_b^{+}\right]}_{{\hat{m}}_c\hat{m}f-\frac{1}{2},{\hat{m}}_{f_p},{\hat{m}}_p}.Math.,.Math.{\left[{?}_b^{+}\right]}_{{\hat{m}}_c\hat{m}f+\frac{1}{2},{\hat{m}}_{f_p},{\hat{m}}_p}.Math.\right]}^T$$f_{A3}\left({\mathrm{??}}_{b,c},?f,?\right)\stackrel{?}{=}{\left[.Math.{\left[{?}_b^{+}\right]}_{{\hat{m}}_c-\frac{1}{2},\hat{m}f{\hat{m}}_{f_p},{\hat{m}}_p}.Math.,.Math.{\left[{?}_b^{+}\right]}_{{\hat{m}}_{c+}\frac{1}{2},{\hat{m}}_f{\hat{m}}_{f_p},{\hat{m}}_p}.Math.\right]}^T$$f_B\left({\mathrm{??}}_C\right)\stackrel{?}{=}[{?}_{n_C}\left(\mathrm{??},?f\left[{\hat{m}}_f\right],?\right),{?}_{n_C}\left(\mathrm{??},?f\left[{\hat{m}}_f-\frac{1}{2}\right],?\right),$${{?}_{n_C}\left(\mathrm{??},?f\left[{\hat{m}}_f+\frac{1}{2}\right],?\right)]}^T$$f_C\left(?a\right)\stackrel{?}{=}\sqrt{{\left[C/N_0\right]}_{\mathrm{Hz}}\left(?a\right)}=\sqrt{\frac{{\left(1+?a\right)}^2}{2{?}_{N_S}^2{N}_S{T}_S}}$

where ?.sub.N.sub.s.sup.2 is the noise variance of the in-phase/quadrature summations of the measured S-LCI correlation in a noise channel. In a typical baseband processing manner of radio signals, when the signal is tracked, the baseband signals consist of an in-phase component (primarily containing useful signal power) and a quadrature component (primarily containing noise power). So, in some cases, the noise variance can also be estimated from the quadrature component in the signal channel.

[0085] As the S-LCI correlation sequence is known to follow a Rayleigh distribution when there is no signal present in the signal processing channel, a noise channel or the quadrature component of a signal channel can be used to estimate the noise variance.

[0086] Then, as the unknown state vector is nonlinear with f(?), the linearized form should be derived to allow the estimates to be updated in an iterative way, and is given by

[00028]$f\left(?\right)?f\left({\widehat{?}}_0\right)+\frac{?f\left(?\right)}{?{?}^{?}}{.Math.}_{?={\widehat{?}}_0}-??$$\mathrm{with}$$??=?-{\widehat{?}}_0$

where {circumflex over (?)}.sub.0 denotes the initial guess (or the previous-epoch estimation) of the state vector. Then, the design matrix satisfies

[00029]$H\stackrel{?}{=}\frac{?f\left(?\right)}{?{?}^{?}}{.Math.}_{?={\widehat{?}}_0}$

[0087] The measurement vector is given by

[00030]$\tilde{f}\stackrel{?}{=}{\left[{\tilde{f}}_{A0},{\tilde{f}}_{A1}^T,{\tilde{f}}_{A2}^T,{\tilde{f}}_{A3}^T,{\tilde{f}}_B^T,{\tilde{f}}_C^T\right]}^T$$\mathrm{with}$${\tilde{f}}_{A0}\stackrel{?}{=}.Math.{\widetilde{?}}_{b,\mathrm{FrFD}}^{+,\left(P,N,M\right)}.Math.$${\tilde{f}}_{A1}\stackrel{?}{=}{\left[.Math.{\widetilde{?}}_{b,\mathrm{FrFD}}^{+,\left(P,N,\mathrm{La}\right)}.Math.,.Math.{\widetilde{?}}_{b,\mathrm{FrFD}}^{+,\left(P,N,\mathrm{Sm}\right)}.Math.\right]}^T$${\tilde{f}}_{A2}\stackrel{?}{=}{\left[.Math.{\widetilde{?}}_{b,\mathrm{FrFD}}^{+,\left(P,\mathrm{Lo},M\right)}.Math.,.Math.{\widetilde{?}}_{b,\mathrm{FrFD}}^{+,\left(P,\mathrm{Hi},M\right)}.Math.\right]}^T$${\tilde{f}}_{A3}\stackrel{?}{=}{\left[.Math.{\widetilde{?}}_{b,\mathrm{FrFD}}^{+,\left(E,N,M\right)}.Math.,.Math.{\widetilde{?}}_{b,\mathrm{FrFD}}^{+,\left(L,N,M\right)}.Math.\right]}^T$${\tilde{f}}_B\stackrel{?}{=}{\left[{\widetilde{?}}_{b,\mathrm{FrFD}}^{\left(P,N,M\right)},{\widetilde{?}}_{b,\mathrm{FrFD}}^{\left(P,\mathrm{Lo},M\right)},{\widetilde{?}}_{b,\mathrm{FrFD}}^{\left(P,\mathrm{Hi},M\right)}\right]}^T$${\tilde{f}}_C\stackrel{?}{=}\sqrt{}$

Gradient Descent Optimization

[0088] One embodiment of the MLE, as mentioned above, applies gradient descent optimization to determine the error values, carrier phase and pseudorange values. Gradient descent optimization is well-known and can be applied to the present context as follows.

[0089] Let {circumflex over (?)} represent an estimate of the state vector ? with the measurement vector {tilde over (f)}. A loss function in terms of the S-LCI correlation amplitude (related to ??.sub.b,c, ?f, ?), the carrier phase (related to ??.sub.c) and the carrier-to-noise power density C/N.sub.0 measurements (related to ??) can be written as follows:

[00031]$\mathrm{Minimize}:J\left(?\right)={\left[f\left(\widehat{?}\right)-\tilde{f}\right]}^T{C}_s^{-1}\left[f\left(\widehat{?}\right)-\tilde{f}\right]$

where C.sub.s is the priori covariance matrix.

[0090] Then, that the optimization problem is unbiased and measurement noise is additive white Gaussian noise (AWGN), the minimization of the loss function is equivalent to maximizing the linearized joint likelihood function, as follows:

Maximize:

[0091] [00032]$p_x\left(x;?\right)=\frac{1}{{\left(2?\right)}^{\frac{11}{2}}{\det }^{\frac{1}{2}}\left[C_5\right]}?\exp \left[-\frac{1}{2}{\left(x-s-H??\right)}^T{C}_s^{-1}\left(x-s-H??\right)\right]$

where det[.Math.] is the determinant function; x?{tilde over (f)}(?), s?f({circumflex over (?)}.sub.0), and H is defined as above.

[0092] Maximizing with

[00033]$\frac{?\ln p_x\left(x;?\right)}{\mathrm{??}}=0$

gives an estimation for the MLE of

[00034]$?\widehat{?}={\left(H^T{C}_s^{-1}H\right)}^{-1}{H}^T{C}_s^{-1}\left(x-s\right)$

with ?? as above.

[0093] The iterative estimation of the state vector is given by

[00035]$\widehat{?}={\widehat{?}}_0+?\widehat{?}$

Synthesized Carrier/Code NCO Frequencies and Carrier Phase/TOA Updates

[0094] As mentioned above, in some embodiments, local carrier/code NCO replicas are generated in by the processing module 40 and fed back to the correlator module 50. In some embodiments, these are created as follows:

[00036]${\hat{f}}_{?.\mathrm{NCO},n_C+1}=f_I+{\hat{f}}_{n_C+1}$${\hat{f}}_{?,\mathrm{NCO},n_C+1}=f_C-\frac{f_c}{f_r}{\hat{f}}_{n_C+1}$$\mathrm{with}$${\hat{f}}_{n_C+1}={\hat{f}}_{n_C}+?{\hat{f}}_{n_C}+?{\widehat{?}}_{n_C}{T}_t$

where subscript n.sub.C is the epoch index for a working snapshot receiver; {circumflex over (f)}.sub.n.sub.c, ?{circumflex over (f)}.sub.n.sub.c, and {circumflex over (?)}.sub.n.sub.c are the estimated Doppler frequency, the estimated frequency error and the estimated Doppler rate respectively; and ? is the coefficient for the Doppler rate estimation which may be determined based on empirical data.

[0095] Both the carrier and code phases of the local replicas are also directly adjusted by the estimated carrier and code phase errors as follows,

[00037]${\widehat{?}}_{0,n_C+1}={\widehat{?}}_{0,n_C}+2?\left(?{\widehat{?}}_{c,n_C}+?{\hat{f}}_{n_C}{T}_t+\frac{1}{2}{\widehat{?}}_{n_C}{T}_t^2\right)$${\widehat{?}}_{\mathrm{PRN},n_C+1}={\widehat{?}}_{\mathrm{PRN},n_C}+\frac{1}{f_c}?{\widehat{?}}_{b,c,n_C}=\frac{1}{f_r}\left(?{\hat{f}}_{n_C}{T}_t+\frac{1}{2}{\widehat{?}}_{n_C}{T}_t^2\right)$

where ?{circumflex over (?)}.sub.c,n.sub.C is the estimated carrier phase error and {circumflex over (?)}.sub.0,n.sub.C.sub.+1 and {circumflex over (?)}.sub.PRN,n.sub.C.sub.+1 are the updated carrier phase and TOA estimates.

Preprocessing

[0096] As mentioned above, in some embodiments, preprocessing is performed between the transform module 60 and the SRM module 70. One exemplary preprocessing process is illustrated in FIG. 9. As can be seen, the illustrated preprocessing process includes a Decision-making step. Before the Decision-making step, however, search steps for matching the code phase, frequency shift and rate should first be determined, to help produce non-singular and more accurate estimations ([10]). Rule-of-thumb early-late code and low-high frequency spacings are adopted: 0.5 chips and 2(3N.sub.ST.sub.S).sup.?1 Hz, respectively, which provide the lowest tolerable correlating power peak. Similarly, the small-large rotation angle related to Doppler rate estimation is ??p.sub.0 where ?p.sub.0 is an optimum search step of FrFT order.

Then, for

[0097] [00038]${?}_{m_c}=\left\{x.Math.-\frac{1}{d_c}?x?\frac{1}{d_c},x??\right\}$${?}_{m_f}=\left\{x.Math.-\frac{1}{T_t?f_{nco}}?x?\frac{1}{T_t?f_{nco}},x??\right\}$${?}_{m_{f_P}}=\left\{x.Math.-\frac{N_s}{2}?x?\frac{N_s}{2}-1,x??\right\}$$\mathrm{and}$${?}_{m_p}=\left\{x.Math.-\frac{1}{?p_0}?x?\frac{1}{?p_0},x??\right\}$

the decision-making process in the CD, FD, and FrFD are given by

[00039]$\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]=\underset{\begin{array}{c}{m}_c?{?}_{m_c},m_f?{?}_{m_f},\\ m_{f_p}?{?}_{m_{f_o},m_p?{?}_{m_p}}\end{array}}{\mathrm{argmax}}\left\{.Math.{\mathrm{DF}}^p\left\{{\widetilde{?}}^{+}\left[n_s,m_c,m_f\right]\right\}.Math.\right\}$

where DF.sup.p{.Math.} is the digital operator of the FrFT process with respect to a FrFT order p and {circumflex over (?)}.sup.+[n.sub.S, m.sub.c, m.sub.f] is the measured noisy ?.sup.+[n.sub.S, m.sub.c, m.sub.f] in the baseband.

[0098] The matching process in the various domains is shown in FIG. 10.

[0099] The values of {circumflex over (m)}.sub.c, {circumflex over (m)}.sub.f, {circumflex over (m)}.sub.f, and my are then conveyed to the early, late, low, and high channels. The measured S-LCI correlations from the prompt, early, late, low, and high channels are denoted as follows:

[00040]$\mathrm{Prompt}.fwdarw.{\widetilde{?}}_b^{+,P}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]$$\mathrm{Early}.fwdarw.{\widetilde{?}}_b^{+,E}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]$$\mathrm{Late}.fwdarw.{\widetilde{?}}_b^{+,L}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]$$\mathrm{Low}.fwdarw.{\widetilde{?}}_b^{+,\mathrm{Lo}}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]$$\mathrm{High}.fwdarw.{\widetilde{?}}_b^{+,\mathrm{Hi}}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]$

[0100] With N.sub.F as the IF sample number over T.sub.S, the normalized amplitudes of each set of correlations are computable as:

[00041]$.Math.{\widetilde{?}}_{b,\mathrm{FrFD}}^{+\left(P,N,M\right)}.Math.=\frac{1}{N_F\sqrt{N_S}}\sqrt{{?}^2\left[{\widetilde{?}}_b^{+,P}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]\right]+{?}^2\left[{\widetilde{?}}_b^{+,P}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]\right]}$$.Math.{\widetilde{?}}_{b,\mathrm{FrFD}}^{+\left(E,N,M\right)}.Math.=\frac{1}{N_F\sqrt{N_S}}\sqrt{{?}^2\left[{\widetilde{?}}_b^{+,E}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]\right]+{?}^2\left[{\widetilde{?}}_b^{+,E}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]\right]}$$.Math.{\widetilde{?}}_{b,\mathrm{FrFD}}^{+\left(L,N,M\right)}.Math.=\frac{1}{N_F\sqrt{N_S}}\sqrt{{?}^2\left[{\widetilde{?}}_b^{+,L}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]\right]+{?}^2\left[{\widetilde{?}}_b^{+,L}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]\right]}$$.Math.{\widetilde{?}}_{b,\mathrm{FrFD}}^{+\left(P,\mathrm{Lo},M\right)}.Math.=\frac{1}{N_F\sqrt{N_S}}\sqrt{{?}^2\left[{\widetilde{?}}_b^{+,\mathrm{Lo}}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]\right]+{?}^2\left[{\widetilde{?}}_b^{+\mathrm{Lo}}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]\right]}$$.Math.{\widetilde{?}}_{b,\mathrm{FrFD}}^{+\left(P,\mathrm{Hi},M\right)}.Math.=\frac{1}{N_F\sqrt{N_S}}\sqrt{{?}^2\left[{\widetilde{?}}_b^{+,\mathrm{Hi}}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]\right]+{?}^2\left[{\widetilde{?}}_b^{+,\mathrm{Hi}}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]\right]}$

[0101] The S-LCI correlations related to small and large (rotation angles) channels come from the outputs of the digital PMF-FrFT in the prompt channel, and the related amplitudes are calculated as

[00042]$.Math.{\widetilde{?}}_{b,\mathrm{FrFD}}^{+\left(P,N,\mathrm{Sm}\right)}.Math.=\frac{1}{N_F\sqrt{N_S}}?\sqrt{{?}^2\left[{\widetilde{?}}_b^{+,P}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p-1\right]\right]+{?}^2\left[{\widetilde{?}}_b^{+,P}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p-1\right]\right]}$

[00043]$.Math.{\widetilde{?}}_{b,\mathrm{FrFD}}^{+\left(P,N,\mathrm{La}\right)}.Math.=\frac{1}{N_F\sqrt{N_S}}?\sqrt{{?}^2\left[{\widetilde{?}}_b^{+,P}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p+1\right]\right]+{?}^2\left[{\widetilde{?}}_b^{+,P}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p+1\right]\right]}$

[0102] The carrier phase measurements for the respective normal, low, and high (frequency) channels are calculated through a four-quadrant arctangent/pure PLL discriminator as follows:

[00044]${\widetilde{?}}_{b,\mathrm{FrFD}}^{P,N,M)}=a\tan 2\left[?\left[{\widetilde{?}}_b^{+,P}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]\right],?\left[{?}_b^{+,P}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]\right]\right]-\arg \left[\sqrt{1-j\cot \left(?\left[{\tilde{m}}_p\right]\right)}{e}^{j?\frac{f_{n_C}^2\left(?{\hat{f}}_o,{\hat{u}}_0\right)\cot \left(?\left[{\tilde{m}}_p\right]\right)}{N_S}}\right]$${\widetilde{?}}_{b,\mathrm{FrFD}}^{P,\mathrm{Lo},M)}=a\tan 2\left[?\left[{\widetilde{?}}_b^{+,\mathrm{Lo}}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]\right],?\left[{?}_b^{+,\mathrm{Lo}}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]\right]\right]-\arg \left[\sqrt{1-j\cot \left(?\left[{\tilde{m}}_p\right]\right)}{e}^{j?\frac{f_{n_C}^2\left(?\widehat{f_0}-\frac{1}{2}?f_{\mathrm{nco}},{\widehat{?}}_0\right)\cot \left(?\left[{\tilde{m}}_p\right]\right)}{N_S}}\right]$${\widetilde{?}}_{b,\mathrm{FrFD}}^{P,\mathrm{Hi},M)}=a\tan 2\left[?\left[{\widetilde{?}}_b^{+,\mathrm{Hi}}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]\right],?\left[{?}_b^{+,\mathrm{Hi}}\left[{\hat{m}}_c,{\hat{m}}_f,{\hat{m}}_{f_p},{\hat{m}}_p\right]\right]\right]-\arg \left[\sqrt{1-j\cot \left(?\left[{\tilde{m}}_p\right]\right)}{e}^{j?\frac{f_{n_C}^2\left(?\widehat{f_0}+\frac{1}{2}?f_{\mathrm{nco}},{\widehat{?}}_0\right)\cot \left(?\left[{\tilde{m}}_p\right]\right)}{N_S}}\right]$

where ?{circumflex over (f)}.sub.0 and {circumflex over (?)}.sub.0, respectively, represent the initial frequency error and Doppler rate estimates.

[0103] A measured carrier-to-noise ratio density related to the incoming signal amplitude is also used in the MLE. In practice, this value can be obtained from previous-epoch measurements, empirical models, etc.

Flowcharts

[0104] FIG. 11 is a flowchart detailing a method according to one aspect of the invention. At step 1100, a GNSS signal is received. At step 1110, the signal is correlated as detailed above. At step 1120, a fractional Fourier transform (FrFT) is applied. At steps 1130A and 1130B, error values of an estimated Doppler rate and a measured distance (pseudorange) between the relevant receiver and a satellite transmitting the GNSS signal are determined. Once these values are determined (e.g., via Phase DR computations or a joint Doppler/TOA MLE, etc), the position information is determined at step 1140.

[0105] As noted above, for a better understanding of the present invention, the following references may be consulted. Each of these references is hereby incorporated by reference in its entirety: [0106] [1] Dierendonck, A. J. Van. (1996). GPS Receivers. In B. W. Parkinson, J. J. Spilker Jr, P. Axelrad, & P. Enge (Eds.), Global Positioning System: Theory and Applications, Volume 1. American Institute of Aeronautics and Astronautics, Inc. [0107] [2] Lashley, M., Bevly, D. M., & Hung, J. Y. (2009). Performance analysis of vector tracking algorithms for weak GPS signals in high dynamics. IEEE Journal on Selected Topics in Signal Processing, 3(4), 661-673. https://doi.org/10.1109/JSTSP.2009.2023341 [0108] [3] Pany, T., & Eissfeller, B. (2006). Use of a Vector Delay Lock Loop Receiver for GNSS Signal Power Analysis in Bad Signal Conditions. 2006 IEEE/ION Position, Location, And Navigation Symposium, 2006, 893-903. https://doi.org/10.1109/PLANS.2006.1650689 [0109] [4] Maier, D. S., & Pany, T. (2021). Multipath and Attitude Estimation Phase Lock Loop for Antenna Array Signal Processing. In proceedings of ION GNSS+2021, pp. 3402-3421. https://doi.org/10.33012/2021.17979 [0110] [5] Gowdayyanadoddi, N. S., Broumandan, A., Curran, J. T., and Lachapelle, G., Benefits of an ultra stable oscillator for long coherent integration, in Proceedings of the 27th International Technical Meeting of The Satellite Division of the Institute of Navigation (ION GNSS+2014), 2014, pp. 1578-1594. [0111] [6] Soloviev, A. and Dickman, J. Extending GPS carrier phase availability indoors with a deeply integrated receiver architecture, IEEE Wireless Communications, Vol. 18, No. 2, pp. 36-44, April 2011. [Online]. Available: http://ieeexplore.ieee.org/document/5751294/[7] [0112] [7] Faragher, R., Powe, M., Esteves, P., Couronneau, N., Crockett, M., Martin, H., Ziglioli, E., and Higgins, C., Supercorrelation as a Service: S-GNSS Upgrades for Smartdevices, in Proceedings of the 32nd International Technical Meeting of the Satellite Division of The Institute of Navigation (ION GNSS+2019), September 16-20, Miami, Florida, USA, 2019. [0113] [8] Groves, P. D., Zhong, Q., Faragher, R., and Esteves, P., Combining Inertially-aided Extended Coherent Integration (Supercorrelation) with 3D-Mapping-Aided GNSS, in Proceedings of the 33rd International Technical Meeting of the Satellite Division of The Institute of Navigation (ION GNSS+2020), September 2020, pp. 2327-2346. [0114] [9] Stirling-Gallacher, R. A., Hulbert, A. P., & Povey, G. J. R. (1996). A fast acquisition technique for a direct sequence spread spectrum signal in the presence of a large Doppler shift. Proceedings of ISSSTA '95 International Symposium on Spread Spectrum Techniques and Applications, 1, 156-160. https://doi.org/10.1109/ISSSTA.1996.563761 [0115] [10] Li, X., & Pahlavan, K. (2004). Super-Resolution TOA Estimation With Diversity for Indoor Geolocation. IEEE Transactions on Wireless Communications, 3(1), 224-234. https://doi.org/10.1109/TWC.2003.819035 [0116] [11] Ozaktas, H. M., Arikan, O., Kutay, M. A., & Bozdagt, G. (1996). Digital computation of the fractional Fourier transform. IEEE Transactions on Signal Processing, 44(9), 2141-2150. [0117] [12] Luo, Y., Li, J., Yu, C., Lyu, Z., Yue, Z., & El-Sheimy, N. (2019). A GNSS software-defined receiver with vector tracking techniques for land vehicle navigation. Proc. ION 2019 Pacific PNT Meeting, Honolulu, Hawaii, USA, Apr. 8-11, 2019-April, 713-727. https://doi.org/10.33012/2019.16834 [0118] [13] Luo, Y., Yu, C., Chen, S., Li, J., Ruan, H., & El-Sheimy, N. (2019). A Novel Doppler Rate Estimator Based on Fractional Fourier Transform for High-Dynamic GNSS Signal. IEEE Access, 7, 29575-29596. https://doi.org/10.1109/ACCESS.2019.2903185 [0119] [14] Luo, Y., Yu, C., Xu, B., Li, J., Tsai, G.-J., Li, Y., & El-Sheimy, N. (2019). Assessment of Ultra-Tightly Coupled GNSS/INS Integration System towards Autonomous Ground Vehicle Navigation Using Smartphone IMU. 2019 IEEE International Conference on Signal, Information and Data Processing (ICSIDP), 1-6. https://doi.org/10.1109/ICSIDP47821.2019.9173292 [0120] [15] Luo, Y., Zhang, L., & El-Sheimy, N. (2018). An improved DE-KFL for BOC signal tracking assisted by FRFT in a highly dynamic environment. 2018 IEEE/ION Position, Location and Navigation Symposium, PLANS 2018-Proceedings, 1525-1534. https://doi.org/10.1109/PLANS.2018.8373547 [0121] [16] Luo, Y., Zhang, L., & Ruan, H. (2018). An Acquisition Algorithm Based on FRFT for Weak GNSS Signals in A Dynamic Environment. IEEE Communications Letters, 22(6), 1212-1215. https://doi.org/10.1109/LCOMM.2018.2828834

[0122] As used herein, the expression at least one of [x] and [y] means and should be construed as meaning [x], [y], or both [x] and [y].

[0123] It should be clear that various aspects of the present invention may be implemented as software modules in an overall software system. As such, some aspects of the present invention may thus take the form of computer executable instructions that, when executed, implements various software modules with predefined functions.

[0124] Similarly, some embodiments of the invention may be executed by a computer processor or similar device programmed in the manner of method steps, or may be executed by an electronic system which is provided with means for executing these steps. Similarly, an electronic memory means such as computer diskettes, CD-ROMs, Random Access Memory (RAM), Read Only Memory (ROM) or similar computer software storage media known in the art, may be programmed to execute such method steps. As well, electronic signals representing these method steps may also be transmitted via a communication network.

[0125] Such embodiments of the invention may be implemented in any conventional computer programming language. For example, preferred embodiments may be implemented in a procedural programming language (e.g., C or Go) or an object-oriented language (e.g., C++, java, PHP, PYTHON or C#). Alternative embodiments of the invention may be implemented as pre-programmed hardware elements, other related components, or as a combination of hardware and software components.

[0126] Further, some embodiments can be implemented as a computer program product for use with a computer system and GNSS device(s). Such implementations may comprise a series of computer instructions fixed either on a tangible medium, such as a computer readable medium (e.g., a diskette, CD-ROM, ROM, or fixed disk) or transmittable to a computer system, via a modem or other interface device, such as a communications adapter connected to a network over a medium. The medium may be either a tangible medium (e.g., optical or electrical communications lines) or a medium implemented with wireless techniques (e.g., microwave, infrared or other transmission techniques). The series of computer instructions embodies all or part of the functionality previously described herein. Those skilled in the art should appreciate that such computer instructions can be written in a number of programming languages for use with many computer architectures or operating systems. Furthermore, such instructions may be stored in any memory device, such as semiconductor, magnetic, optical or other memory devices, and may be transmitted using any communications technology, such as optical, infrared, microwave, or other transmission technologies. It is expected that such a computer program product may be distributed as a removable medium with accompanying printed or electronic documentation (e.g., shrink-wrapped software), preloaded with a computer system (e.g., on system ROM or fixed disk), or distributed from a server over a network (e.g., the Internet or World Wide Web). Of course, some embodiments of the invention may be implemented as a combination of both software (e.g., a computer program product) and hardware. Still other embodiments of the invention may be implemented as entirely hardware, or entirely software (e.g., a computer program product).

[0127] A person understanding this invention may now conceive of alternative structures and embodiments or variations of the above all of which are intended to fall within the scope of the invention as defined in the claims that follow.