Abstract:
A deconvolution approach to adaptive signal processing has been applied to the elimination of signal multipath errors as embodied in one preferred embodiment in a global positioning system receiver. The method and receiver of the present invention estimates then compensates for multipath effects in a comprehensive manner. Application of deconvolution, along with other adaptive identification and estimation techniques, results in completely novel GPS receiver architecture.

Description:
ORIGIN OF THE INVENTION 
     The invention described herein was made in the performance of work under a NASA contract, and is subject to the provisions of Public Law 96-517 (35 U.S.C. § 202) in which the Contractor has elected not to retain title. 
    
    
     TECHNICAL FIELD 
     The present invention relates to signal error reduction, and more specifically to a method and apparatus for reducing Global Positioning Satellites (GPS) multipath signal error using deconvolution to determine the attitude and position of kinetic platforms. 
     BACKGROUND ART 
     The superiority of GPS receiver technology over other competing technologies, in terms of the cost, reliability, mass, size, and power considerations, as proven in its various applications such as navigation, spacecraft orbit determination, and surveying, has constantly prompted the extension of GPS to other important areas that require increasingly more exacting performance from the GPS system. One such area, for example, is the application of GPS to the attitude determination of aircraft and earth orbiting satellites. 
     For the GPS receivers to provide the requisite precision, the existing error sources in the GPS system must be eliminated or greatly reduced. At the current state of GPS technology, the most significant error source is the signal multipath propagation. For example, in spacecraft operational environments, the GPS signal is reflected from various structural components of the spacecraft and these reflected components are received by the GPS receiver along with the desired direct line-of-sight (LOS) path. The reflected signals differ from the desired LOS path signal in terms of their delays, amplitudes, and phases. The carrier phase tracking loop provides no inherent discrimination against the multipath signals and thus tracks the phase of the composite signal corrupted by multipath components. The resulting differential carrier phase estimation error can be orders of magnitude higher compared to the case of no multipath propagation in many GPS applications. For example,, the measurements obtained by the RADCAL satellite GPS-ADS (Attitude Determination System) experiments have shown that the differential range error in such environment is of the order of 1 cm corresponding to an attitude determination error of about 0.5 degree. Thus for the GPS receivers to provide precision pointing knowledge (order of 1 arcmin or better with 1 meter antenna baseline) or a differential range accuracy of about 0.3 mm or better, the multipath effects must be suppressed by orders of magnitudes. Similar accuracy may also be desirable in other GPS precision applications in the presence of multipath signals such as GPS based geophysical measurements and precision surveying. 
     Among the past approaches to deal with the multipath problem, one approach involves reducing the early-late delay spacing among the correlators in the GPS receiver code lock loop. However while this approach reduces the code range errors to some extent, it does not aid in the carrier phase measurements accuracy that is the basis of most GPS precision applications. Moreover even the reduction in the code range error is limited and if the early-late spacing is smaller than the initial delay error due to multipath signals (easily the case with many multipath situations), then the loop error can be very high and the loop may not even track. 
     In the work of VanNee et al., a set of nonlinear implicit equations are derived on the basis of maximum likelihood estimation theory. These equations contain the amplitude,, phase and delay of multipath signals. The paper proposes to solve these highly nonlinear equations in a recursive manner, but there is no explicit solution presented. The method has been implemented in a Nov Atel receiver using multiple correlators. Although this approach is useful for code phase measurement and a benign multipath environment, it may not be suitable for severe multipath environment and for precision GPS applications based on carrier phase measurements. 
     STATEMENT OF THE INVENTION 
     The present invention provides a method and apparatus for simultaneous estimation and compensation of the multipath estimation errors in both the carrier phase lock and the code delay lock loops in GPS receivers. The preceding is accomplished by disclosing new techniques for dealing with the multipath problem and includes these techniques in an improved GPS receiver. These techniques are based on the application of the optimal deconvolution approach in a somewhat unconventional manner, as compared to its application in other fields such as seismology and telecommunications. The proposed method consists of first estimating the impulse response of the effective multipath channel by a least squares algorithm. This step is followed by obtaining an inverse filter which equalizes the multipath channel response to the desired ideal multipath free response to the maximum extent possible within the specified constraints of the implementation complexity. From the equalized response, estimates of the true carrier phase and code delay can be made. The simulation results demonstrate that the proposed method is capable of reducing the multipath distortion by several orders of magnitude in an environment that is more severe than may possibly take place in any realistic precision GPS applications environment. 
     The invention also is embodied in a receiver comprising a series of correlators which generate demodulated signals from known reference signals, and the demodulated reference signals are then used to develop a discriminator function from which multipath error can be approximately determined. The present invention then compensates for the multipath error, and produces a final differential phase and code range estimate after compensating for the multipath error. 
     The present invention presents an optimum solution to the inherent tradeoffs among hardware/software implementation complexity, the extent of the multipath expected in the application, and the degree of multipath cancellation. 
    
    
     BRIEF DESCRIPTION OF THE DRAWINGS 
     The objects and features of the present invention, which are believed to be novel, are set forth with particularity in the appended claims. The present invention, both as to its organization and manner of operation, together with further objects and advantages, may best be understood by reference to the following description, taken in connection with the accompanying drawings. 
     FIG. 1 is a block diagram of the elements of a receiver in a preferred embodiment of the present invention; 
     FIG. 2 is a flow diagram of the channel impulse response estimator and equalizer of FIG. 1; 
     FIG. 3 is a graph illustrating the tracking of the Discriminator Function amplitude vs. delay in simulation example 1; 
     FIG. 4 is a graph illustrating the Discriminator response vs. delay in simulation example 1; 
     FIG. 5 is a graph of the normalized delay error magnitude vs. iterations for simulation example 1; 
     FIG. 6 is a graph of the pointing error vs. iterations for simulation example 1; 
     FIG. 7 is a graph of the pointing error vs. iterations for simulation example 1 with Delay&gt;0.9; 
     FIG. 8 is a graph of the normalized delay error magnitude vs. iterations for simulation example 2; 
     FIG. 9 is a graph of the pointing error vs. iterations for simulation example 2; 
     FIG. 10 is a graph of the pointing error vs. iterations for simulation example 3; and 
     FIG. 11 is a graph of the normalized delay error vs. iteration number for simulation example 3. 
    
    
     DETAILED DESCRIPTION OF THE INVENTION 
     The following description is provided to enable any person skilled in the art to make and use the invention and sets forth the best modes contemplated by the inventors of carrying out their invention. Various modifications, however, will remain readily apparent to those skilled in the art, since the generic principles of the present invention have been defined herein specifically to provide a method and apparatus for reducing multipath error especially in GPS applications. 
     In the more conventional telecommunication applications of the deconvolution approach or the equalizer theory, the multipath propagation channel is modeled as, ##EQU1## where {u j  } represents the transmitted symbol sequence, {y j  } is the channel output sequence, and {h -q .sbsb.1, h -q .sbsb.1 +1 , . . . , h 0 ,h 1 , . . . h q .sbsb.2 } represents the discrete channel impulse response. The additive noise sequence {v j  } is usually assumed to be a zero-mean white Gaussian. The basic process of deconvolution involves the estimation of the transmitted input sequence {u k  } on the basis of noisy observations {y k  } assuming that the discrete channel impulse response {h k  } is known to the receiver. In the case of unknown channel response, adaptive equalization techniques may be used wherein first an approximate estimate of {h k  } is obtained on the basis of a training sequence known in advance to the receiver and the channel output {y k  } and subsequently the estimate of {h k  } is refined adaptively with u k  replaced by its estimated/detected version in the adaptive algorithm. From the real time estimate of {h k  }, a time-varying inverse filter is derived which then filters {y k  } to obtain the estimate of {u k  }. In practice the two steps of estimating {h k  } and then finding the corresponding inverse filter are combined into a single step of finding directly the adaptive equalizer coefficients. The problem wherein the adaptive equalization is achieved without any training sequence is relatively more difficult and has also received considerable attention in the literature. 
     In the following, the deconvolution approach is extended to the problem of multipath elimination in the GPS receiver code tracking and carrier phase lock loops. Referring to the drawing of GPS receiver multipath compensator 100 shown in FIG. 1, the antenna 1 receives the GPS signal which is filtered and amplified by the RF front end 2. In the absence of the multipath signals, the output 3 of the RF front end 2 is given by: ##EQU2## where the receiver noise is not considered in the first instance, A c  is the received signal amplitude, ω c  is the GPS signal frequency, and a(t)=±1 is the pseudo-random code waveform that phase modulates the carrier. It is assumed that the data modulation is removed in a decision-directed manner. In the presence of N multipaths in addition to the direct line-of-sight path, the input signal may be characterized as 
     
         s.sub.m (t)=s(t)+α.sub.1 s(t-τ.sub.1)+. . . +α.sub.N s(t-τ.sub.N)                                          (3) 
    
     where α i  and τ i  denote respectively the amplitude and delay of the i th multipath for i=1,2, . . . ,N, and s(t) is the direct path or line-of-sight component of s m  (t). After substituting (2) in (3), the composite received signal may be expressed as follows: ##EQU3## 
     This signal s m  (t) is translated from the carrier frequency to an intermediate frequency (IF) by converter 4 to yield IF signal s IF  (t) at an intermediate frequency ω 1 . Thus, equation (4) can then be represented as ##EQU4## In the conventional delay-lock discriminator the signal in (5) is correlated with the reference signals ##EQU5## where t is the delay tracking error, τ d  is the offset delay referred to as early and late correlator delay and θ c  is the reference oscillator phase. The resulting correlation functions are given by 
     
         s.sub.IF (t).s.sub.L (t)=a.sub.c α(t).α(t-τ-τ.sub.d)cos(θ.sub.m.sbsb.0 -θ.sub.C)+. . . +A.sub.c α.sub.N α(t-τ.sub.N).α(t-τ-τ.sub.d).cos(θ.sub.m.sbsb.N -θ.sub.C).                                        (7) 
    
     The resulting signal denoted by R 1i  (t) may be written in the form below. 
     
         R.sub.1i (τ)=h.sub.0i R.sub.c (τ+τ.sub.d)+h.sub.1i R.sub.c (τ+τ.sub.d -τ.sub.1)+. . . +h.sub.Ni R.sub.c (τ+τ.sub.d -τ.sub.N)h.sub.ki =α.sub.k A.sub.c cos(θ.sub.m.sbsb.k -θ.sub.c);k=0,1, . . . ,N  (8) 
    
     In equation (8) R c  (τ) represents the code autocorrelation function and α 0  =1. Similarly the corresponding signal obtained with s L  replaced by S E  in (7) and denoted by R 2i  (t) may be written in the following form. 
     
         R.sub.2i (τ)=h.sub.0i R.sub.c (τ-τ.sub.d)+h.sub.1i R.sub.c (τ-τ.sub.d -τ.sub.1)+. . . +h.sub.Ni R.sub.c (τ-τ.sub.d -τ.sub.N)                          (9) 
    
     The difference between R 1i  and R 2i , the so called discriminator function D i  (τ) is now given by 
     
         D.sub.i (τ)=h.sub.0i g.sub.c (τ)+h.sub.1i g.sub.c (τ-τ.sub.1)+. . . +h.sub.Ni g.sub.c (τ-τ.sub.N) (10) 
    
     where 
     
         g.sub.c (τ)=R.sub.c (τ+τ.sub.)-R.sub.c (τ-τ.sub.d 
    
     is the discriminator function in the ideal case. Clearly if α k  =0 for all k≠0 then D i  (τ)=g c  (τ) corresponding to the case of no multipath propagation. The standard delay lock loop converges to the solution τ m  of the equation D i  (τ)=0 instead of converging to 0 which is the solution of g c  (τ)=0. The multipath error τ m  can be excessive depending upon the actual multipath environment encountered as will be illustrated by several simulation examples. 
     In the approach of the present invention the discriminator function is measured first and then with equation (10) is used to estimate all the unknown variables including the multipath delays, amplitudes, and phases and the estimates of the multipath errors both in the phase lock and code lock loops via deconvolution. The multipath errors are then compensated for in arriving at the final differential phase and code range estimates. 
     To accomplish these objectives, multiple correlators corresponding to multiple values of τ d  are used in the present invention instead of just two correlators. Returning to FIG. 1, signal 5 is correlated with carrier signal 43 C i  (t) in correlator 6, and the output 7 is correlated in correlators 8,12,16 with signals 54, 53, and 55, respectively. The outputs 9,13, and 17 are received by integrate/dump units 10, 14, and 18 and each output is integrated over T seconds to generate the correlator functions represented in equation (8). In the same manner, all the correlator functions R i  (τ+iΔ) for all integers i in the range of -Q 1   to Q 2   are generated using additional correlators and integrate/dump units. Output 11, R i  (τ+Q 2   Δ) is the same as R 11  (τ) of equation (8) with τ d  in (8) replaced by Q 2   and the carrier signal 43 C i  (t) is given by 
     
         C.sub.i (t)=2cos (ω.sub.1 t+θ.sub.c).          (11) 
    
     The multipath resolution parameter Δ is selected sufficiently small to provide the required multipath resolution and the number of correlators (Q 1  +Q 2  +1) is dependent on the extent of multipath delay spread and Δ. The linear combiner 20 receives output 11, 15, and 19 from the respective integrate/dump units 10,14, and 18, as well as any other integrate/dump units (not shown) and generates the sampled values of the function D i  (τ) given in (10) by differencing pairs of its input, that is, 
     
         D.sub.i (τ+kΔ)=R.sub.i (τ+(k+I.sub.d)Δ)-R.sub.i (τ+(k-I.sub.d)Δ);k=-M.sub.i, . . . , 0, . . . , M.sub.2. (12) 
    
     where τ d  is selected equal to I d  Δ for some integer I d  and M 1  and M 2  values are dictated by the extent of multipath. Thus, in terms of M 1  and M 2 , the Q 1   and Q 2   are given by Q 1  =M 1  +I d , Q 2  =M 2  +I d . The output 36 from the linear combiner 20 D i  is the (M 1  +M 2  +1) vector whose k th component is given by D i  (τ+kΔ) in (12). 
     The signal 5, S IF  (t) is correlated in 21 with the carrier signal C q  (t) 45 which is the quadrature phase version of C i  (t) and is given by: 
     
         C.sub.q (t)=-1 sin (ω.sub.1 t+θ.sub.c)         (13) 
    
     Correlating output 22 of the correlator 21 with signal a(t-τ) in correlator 27 and integrating the output 28 over an interval of T seconds in the integrate/dump unit 29 generates R q  (τ) which is related to the parameters of interest by equation (14) as 
     
         R.sub.q (τ)=h.sub.0q R.sub.c (τ)+h.sub.1q R.sub.c (τ-τ.sub.1)+. . . +h.sub.Nq R.sub.c (τ-τ.sub.N) h.sub.kq =α.sub.k A.sub.c sin(θ.sub.m.sbsb.k -θ.sub.c);k=0,1, . . . ,N                                                      (14) 
    
     In the same manner all of the correlator functions R q  (τ+iΔ) for all integers i in the range of -Q 1   to Q 2   are generated using additional correlators exemplified by 23,31 and whose outputs 24,32 are received by integrate/dump units 25,33 as shown in FIG. 1. Outputs 26,30,34 are transferred to the linear combiner 35 which generates (M 1  +M 2  +1) vector D q  whose k th component D q  (τ+kΔ) is given by: 
     
         D.sub.q (τ+kΔ)=R.sub.q (τ+(k+I.sub.d)Δ)-R.sub.q (τ+(k-I.sub.d)Δ; k=-M.sub.1, . . . , 0, . . . , M.sub.2. (15) 
    
     From equation (15) it follows that D q  (τ) may be expressed as: 
     
         D.sub.q (τ)=h.sub.0q g.sub.c (τ)+h.sub.1q g.sub.c (τ-τ.sub.i)+. . . +h.sub.Nq g.sub.c (τ-τ.sub.N). (16) 
    
     For the purpose of estimating and filtering, the outputs 36 and 37 are combined into the vector D given by 
     
         D=D.sub.i +jD.sub.q                                        (17) 
    
     where j=√-1 and the k th component of D denoted by D k  is given by 
     
         D.sup.k =D(τ+kΔ)=D.sub.i (τ+kΔ)+jD.sub.q (τ+kΔ)                                          (18) 
    
     From equations (10), (16), and (18) it follows that 
     
         D(τ)=h.sub.0 g.sub.c (τ)+h.sub.1 g.sub.c (τ-τ.sub.1)+. . . +h.sub.N g.sub.c (τ-τ.sub.N)+n(τ) h.sub.k =h.sub.ki +jh.sub.kq ; k=1, . . . , N.                                         (19) 
    
     Here, D(τk) is the complex-valued discriminator function, and n(t) represents the noise at the correlator output corresponding to the noise at the receiver input. 
     FIG. 2 illustrates a flow chart of the steps executed by the channel impulse response estimator and equalizer 38. Step 200 involves obtaining the two discriminator functions D i  and D q  from the linear combiners 20,35 and solves for possible channel impulse responses h k  based on delay kΔ. The output from step 200 is a series of channel impulse responses h -N  through h N , which are then processed in step 210. Step 210 selects the impulse response which best fits the measurement for the inputs D i  and D q , denoted as h k .sbsp.0. In step 220, the corresponding deconvolution filter coefficients ƒ are computed from the results of step 210, and the coefficients are used in step 230 to remove the multipath contribution by convolving measurements {D j  } with {ƒ j  }. The resultant discriminator function D eq  is then solved in step 240 for Re{D eq  (τ)}=0 by linear interpolation of the sampled version of D eq . The delay obtained, τ is used in step 250 to solve for the channel impulse response h.sup.ƒ and to compute the argument (phase) of the central coefficient h 0 .sup.ƒ in the final step 260. 
     Turning to the preferred solution, in FIG. 1 estimator 38 processes the input vectors D i  and D q  to estimate the multipath channel impulse response, obtain the coefficients of a deconvolution filter, and from these derive an estimate of the delay τ and the phase difference θ=(θ m0  -θ C ) such that these estimates are not corrupted by the multipath signals. For an understanding of the estimator 38, Equation (19) is rewritten in its discrete form (20), that is, substitute τ=jΔ where j is an integer into equation (19): ##EQU6## Here y j  =D(jΔ); g j  =g c  (jΔ); n j  =n(jΔ); and (-M 1 ,M 2 ) represents the interval over which the measured discriminator function is significant. Note that equation (20) is of somewhat more general form than (19) in that this form may also include negative values of the multipath delays. The delays τ 1 , . . . , τ N  are approximated by integer multiples of Δ which can be achieved by selecting Δ sufficiently small. A least squares estimation technique is now applied to estimate the channel impulse response coefficients {h k  } and the corresponding deconvolution or equalizer filter coefficients from {y j  }. Letting 
     
         h.sup.T = h.sub.-q.sbsb.1, . . . , h.sub.-1, h.sub.0, h.sub.1, . . . , h.sub.q.sbsb.2 ! 
    
     
         g.sub.j.sup.H = g.sub.j+q.sbsb.1, . . . , g.sub.j+1, g.sub.j, g.sub.j-1, . . . , g.sub.j-q.sbsb.2 !; j=-M.sub.1, . . . ,0, . . . , M.sub.2 
    
     the least squares estimate of the unknown channel impulse response vector is given by ##EQU7## where σ n   2  =E n j   2  ! and I is the identity matrix. Now denoting by χ 0  the truncated (or zero padded depending upon whether q i  &gt;K i  or vice versa; i=1,2) version of h, i.e., with 
     χ 0  = h K .sbsb.1, . . . , h 0 , . . . h -K .sbsb.2 ! 
     where K 1  and K 2  are the selected numbers of the deconvolution filter coefficients, and by χ j  the j times shifted version of χ 0 , 
     X j  = h K .sbsb.1 +j , . . . , h j , . . . , h -K .sbsb.2 +j  ! 
     the optimum parameter vector f is given by ##EQU8## Note that in the definition of χ j  the entries for h j  are set equal to zero if j falls outside the interval  -q 1 ,q 2  !. Also the constant κ is given by ##EQU9## The equalized channel response z is obtained by convolving the two sequences {h -q .sbsb.1, . . . ,h 0 , . . . ,h q .sbsb.2 } and {ƒ -K .sbsb.1, . . . ,ƒ 0 , . . . ,ƒ K .sbsb.2 }. The elements of z are denoted by {Z - (q.sbsb.1 +K .sbsb.1.sub.), . . . ,z 0 , . . . ,z.sub.(q.sbsb.2 +K .sbsb.2.sub.) }. Therefore the equalized discriminator response obtained by convolving {y j  } with {ƒ j  } sequence is given by ##EQU10## For the case of perfect equalization z 0  =1, z i  =0 for i≠0 and D eq  (τ)=g c  (τ). In practice Re{D eq  (τ)}≈g c  (τ) and the solution of the equation 
     
         Re{D.sub.eq (τ)}=0                                     (24) 
    
     is the estimate of the delay τ. When the initial delay uncertainty τ is smaller compared to Δ/2 in magnitude, the j th components of the vector D i  and D q  at the input of 38 are approximately equal to the required real and imaginary parts of y j  in equation (20), and thus the deconvolution of the multipath channel is achieved by replacing y j  by the j th component of the vector D.tbd.D i  +√-1D q . 
     Delay Estimation 
     When the initial delay uncertainty τ is not small compared with Δ, it is explicitly taken into account by the estimator 38. Denoting the initial uncertainty by τ.sub.ε -τ.sub.ρ, the equation (19) is modified by replacing τ by (τ.sub.ε -τ.sub.ρ), yielding 
     
         D(τ)=h.sub.0 g.sub.c (τ.sub.ε -τ.sub.ρ)+h.sub.1 g.sub.c (τ.sub.ε -τ.sub.1 -τ.sub.ρ)+. . . +h.sub.N g.sub.c (τ.sub.ε -τ.sub.N -τ.sub.ρ)+n(τ). (25) 
    
     In (25) τ.sub.ρ represents the negative of the signed integer part of (τ/Δ) and τ.sub.ε is the fractional part of (τ/Δ) such that 
     
         τ=τ.sub.ε -τ.sub.ρ =τ.sub.ε -k.sub.0 Δ; |τ.sub.ε |&lt;Δ/2 
    
     for some signed integer k 0  and with D selected above. The equation (25) may be rewritten as 
     
         D(τ)=h.sub.0 g.sub.k.sbsb.0 (τ.sub.ε)+h.sub.1 g.sub.k.sbsb.0 (τ.sub.ε -τ.sub.1)+. . . +h.sub.N g.sub.k.sbsb.0 ·(τ.sub.ε -τ.sub.N)+n(τ) g.sub.k.sbsb.0 (τ.sub.ε)= g.sub.c (τ.sub.ε -k.sub.0 Δ)                                         (26) 
    
     Assuming that |k 0  |&lt;N for some integer N, then k 0  is the integer k that minimizes the following index ##EQU11## where h k  is the channel response obtained on the basis of replacing k 0  by k in (26), 
     
         g.sup.k (τ)= g.sub.k (τ.sub.ε),g.sub.k (τ.sub.ε -τ.sub.1), . . . ,g.sub.k (τ.sub.ε -τ.sub.N)!.sup.T (28) 
    
     and g k  (τ.sub.ε) is approximated by g c  (-kΔ). In the estimation of h k  for any integer k, the discretization and estimation procedure of (20)-(22) is followed. Thus the estimator 38 minimizes the following index ##EQU12## where g 0   k  is the vector of dimension (q 1  +q 2  +1) whose lth component is equal to g c  (-(l+k)Δ), g j   k  is equal to the j time shifted version of g k , and h k  is the solution of equation (21) with g j  replaced g j   k . 
     If k 0  is the solution of the minimization in (29), then the overall channel estimate is given by h=h k .sbsp.0. Substitution of h in (22) provides the coefficients of the inverse filter and the solution of (24) yields τ, signal 46 in FIG. 1. Note that if it is known that the multipath error τ is smaller than Δ/2 in magnitude, then this additional procedure is not necessary. However for k 0  &gt;0, considerable error can otherwise ensue in the estimate of h. 
     Multipath Phase Estimation 
     With τ denoting the estimate of the multipath delay error, one now solves equation (20) for h but with g j  equal to g c  (jΔ-τ) by least squares estimation procedure. Denoting the resulting least squares solution by h.sup.ƒ, then the multipath phase error estimate is given by the argument of the zeroth component of the estimated impulse response, i.e., θ=+arg (h 0 .sup.ƒ) which corresponds to signal 39 in FIG. 1. This error is compensated from the carrier phase-lock loop phase estimate in the following estimation cycle. 
     Output 39 of estimator 38 is input to a loop filter 40, which may be a digital filter similar to the ones used in phase-locked loops. In its simple form it is a multiplier by a gain (1/G c  T) where is G c  the sensitivity of the numerically controlled oscillator (NCO) 42. Output 41 of the carrier loop filter 40 controls the phase of the carrier signal 43 (C i  (t)) at the output of NCO 42. Signal 43 is inputted into phase shifter 44 which generates a π/2 phase shifted signal 45 (C q  (t)). Output 46 of estimator 38 is input to a code loop filter 47, which may be a digital filter similar to the ones used in delay lock loops or phase lock loops. In its simple form it is a multiplier by a gain (1/G d  T) where G d  is the sensitivity of the NCO 49. Output 48 of the code loop filter 47 controls the phase (delay) of the output 50 of the code NCO 49. The output of code NCO is inputted to a code generator 51 whose output is the reference code a(t-τ+Q 1  Δ) which is the binary pseudo random code and is the delayed version of the code a(t) appearing in the received signal 3 (S m  (t). 
     The delay τ of this code is controlled by signal 48 at the input of code NCO 49. The code delay signal 52 is inputted to a tapped delay line 59 implemented by a shift register (not shown). The (Q 1  +Q 2  +1) outputs of this tapped delay exemplified by 52,53, and 54 are input to various correlators such as 8, 12, and 16, and to the correlators 23, 27, and 31. 
     In the close-loop estimator, τ is the estimate of the delay error τ between the code phase of the direct component of the incoming signal 5 and the phase of locally generated signal 53 at one sampling time instance earlier. This error is corrected via the signal 48 at the input to the code NCO 49. With the recursive procedure in the limit the error τ will approach zero and the code signal 53 will achieve lock with the direct path component of the incoming signal, thus eliminating the effect of multipath signal components present in signal 5. 
     Similarly, θ is the estimate of the phase error between the phase of the direct path component of the incoming signal 5 and the local reference signal 43 (C i  (t)), and by generating signal 41 via lock loop filter 40, the control of the phase of carrier NCO 42 is such as to drive this error to zero. With the recursive procedure the carrier NCO output signal 43 will achieve phase lock with the direct path component of the incoming signal, thus eliminating the effect of the multipath signal components on the carrier phase tracking. The accumulated value of signal 41 multiplied by G c  T or the accumulated value of signal 39 when the carrier loop filter 40 is a multiplier with gain (1/G c  T), is the estimate of the difference between the phase of the incoming signal direct path component and the initial carrier NCO signal phase. The value of signal 46 when the code loop filter is a multiplier with gain (1/G d  T), is the estimate of the delay between the code phase of the incoming signal direct path component and the initial code phase of the code NCO output 50. 
     For the case of GPS signals, this delay difference multiplied by the velocity of light is equal to the multipath error free estimate of the code phase range. Similarly the carrier phase difference multiplied by (λ/2π) where λ is the carrier signal wavelength, is equal to the multipath error free estimate of the carrier phase range, also known as the accumulated delta range. 
     SIMULATION RESULTS 
     In this section some simulation examples are presented depicting the performance of the proposed multipath cancellation algorithm. First some additional notations are introduced to present these results. Let T and T c  denote the sampling period for the signal processing system and the code chip period respectively. Also let τ d  denote the normalized delay τ d  /T c  and g be the sample signal-to-noise power ratio given by γ=(P c  T/N 0 ) where P c  =A c   2  /2 is the signal power received by the direct or line-of-sight path and N 0  is the one-sided power spectral density of the receiver thermal noise. Further denote by a and θ m  the vectors consisting of the amplitudes and phases of the discrete channel response h in equation (26). One motivation behind this invention has been the GPS application to precision attitude determination of LEO (low earth orbit) spacecrafts. Therefore the achievable carrier phase estimation errors are interpreted in terms of the equivalent attitude pointing errors. It may be easily seen that for 1 meter antenna baseline the conversion factor between the carrier phase error in radians and the corresponding pointing error in arcmin is approximately equal to 109.4. Several simulation examples are presented below in terms of the notations introduced in this section. Recall that (q 1  +q 2 ) is equal to the number of discrete multipaths , (K 1  +K 2  +1) is the number of filter taps and g is the sample SNR. 
     Example 1 
     q 1  =0; q 2  =8; K 1  =12; K 2  =4; γ=10 3  ;τ d  =0.5 
     a= 1 0.2 0.5 0.8 0.9 0.7 1 0.2 0.9! 
     θ m  = 8 0 -0.5 0.7 1.2 0.8 -0.6 1.3 0.2 -1.1! 
     FIG. 3 plots the ideal discriminator response g c  (τ), distorted response D(t) of (19), and the discrete version of the equalized response D eq  (τ). Clearly while the zero-crossing of D(t) is about 0.4 T c , the zero-crossings of both g c  (τ) and D eq  (τ) are equal to zero. Reducing the tap spacing to τ d  =0.1 does not make any significant difference in the zero-crossing of the function D(t) as shown in FIG. 4, which plots the various discriminator functions for the case of τ d  =0.1. FIG. 5 shows the convergence of the delay error as a function of the number of samples processed with an initial normalized delay of -0.498 chips. 
     It can be seen that the steady state error without multipath correction will remain approximately equal to 0.4 T c  while with the correction algorithm the error is approximately equal to 0.001 T c . Note also that with a 50 dB-Hz CNR at the receiver input this period of 100 iterations corresponds to only 1 sec of real time as per the definition of g. FIG. 6 shows the carrier phase error expressed in arcmin for the GPS attitude determination application. Note that the initial error of more than 30 arcmin is reduced to about 0.1 arcmin. 
     Example 2 
     q 1  =0; q 2  =25; K 1  =35; K 2  =4; γ=10 4  ; τ d  =0.1 
     a= 1 0.8 0.5 0.3 0.9 0.2 0.5 0.3 0.7 0.1 0.05 0 0.08 0.1 0.03 0.02 0.1 0.07 0.05 0.04 0.1 0 0.05 0.08 0.01! 
     θ m  = 0 0.5 0.7 -0.7 1 0.3 0.5 -0.6 1.5 1 0.8 0.2 0.1 0 1 -0.5 -0.3 0.6 0.1 0.5 -0.8 -0.5 0.9 1.2 -0.6 1! 
     FIGS. 7 and 8 depict respectively the convergence of carrier phase and the code delay errors versus the number of samples. In this example the sample SNR is equal to 40 dB and thus corresponds to an update period of 0.1 sec for a 50 dB-Hz CNR at the GPS receiver input. Note that the multipaths with delay greater than 9 D are weak in that their relative amplitudes are less than 0.1 of the direct path which makes the identification of these multipaths more difficult. As may be inferred from these figures a residual pointing error of less than 1 arcmin and a normalized delay error of 0.005 chips is obtained in less than 50 iterations. Since the amplitudes of multipaths for delays exceeding 9 D are relatively small, one may consider ignoring these and thus select a filter order that is much smaller than q. FIG. 9 shows the resulting error performance when the filter order K 1  is selected to be 10. As is apparent from the figure the residual error, while much smaller compared to the initial error of about 30 arcmin, is much higher than that achieved with K 1  =35 filter case. 
     Example 3 
     q 1  0; q 2  =35; K 1  =35; K 2  =4; γ=10 3  ; τ d  =0.1 
     a= 1 0.8 0.5 0.3 0.9 0 0.2 0.5 0.3 0.7 0.1 0.5 0.6 0.2 0 1 0.7 0.5 0.2 0.3 0.5 0.7 1 0.4 0.8 0.3 0.7 1 0.7 0.8 0.6 1 0.5 0.3 0 0.7! 
     θ m  = 0 0.5 0.7 -0.7 1 0.3 0.5 -0.6 1.5 1 0.8 0.2 0.1 0 1 -0.5 -0.3 0.6 0.1 0.5 -0.8 -0.5 0.9 1.2 -0.6 1 0.7 -1 -0.7 0.8 0.6 0.5 -0.8 1 0.4 -0.9! 
     FIGS. 10 and 11 plot the residual carrier phase and code delay errors respectively, for  g!=30 dB. As may be inferred from these figures, for this case of very severe multipath propagation, the pointing error is reduced from about 32 arcmin to only 4.5 arcmin while the residual delay error is about 0.01 chips. 
     It will be understood that the embodiment described herein are merely exemplary and that a person skilled in the art may make many variations and modifications without departing from the spirit and scope of the invention. For example, one such modification comprises sequential time sharing among a relatively small number of correlators to obtain correlation signals over a much larger range of delay offsets. This is achieved by sequentially offsetting the code NCO phase over a range of delay offsets and computing correlation signals for each such offsets instead of computing all the correlation signals in parallel as described above in the preferred embodiment. Such sequential method may be attractive if the phase and delay of the input signal vary relatively slowly with time. All such variations and modifications are intended to be included within the scope of the invention as defined in the appended claims.