High Performance , Low Cost Loran-C Cycle Identification and ECD Estimation

A Loran-C system is a hyperbolic navigation system which works based on time difference of arrival (TDOA). Cycle identification is the task of finding time of arrival of the incoming signal. Finding time of arrival of the signal needs choosing a reference point which is defined to be the third zero crossing of the signal. ECD is the varied Loran-C signal’s envelope from the original pulse. Cycle identification and ECD estimation accuracies have considerable effects on the Loran-C system receivers’ localization accuracy and the error in cycle identification or ECD estimation has direct effect on each other so it is important to estimate the reference point and ECD as precise as possible. In this paper, algorithms for cycle identification and ECD estimation are proposed.Furthermore, this paper addresses the problem of the reference points existence between two samples and proposes two algorithms to estimate the reference points time of arrival between two samples which leads to reach high accuracy using low sampling frequency. The simulation results show that the proposed methods for cycle identification and ECD estimation are robust in noisy conditions and intersample cycle identification algorithms give accurate estimate of the reference points between two samples.


Introduction
A Loran-C system is a hyperbolic navigation system based on the TDOA methods.A Loran-C system consists of one Master station and at least two secondary stations which are organized into a chain.The master station transmits 9 pulses.The first eight pulses are 1 millisecond apart from together and the ninth pulse spaced 2 milliseconds after the eighth pulse.The secondary stations transmit a series of eight pulses, each spaced 1 millisecond apart from together too [1] [2] .Once the last secondary has transmitted, the master transmits again, and the cycle is repeated.The period consists of two series of pulses for each station and it is called the group repetition interval (GRI) which is an important characteristic of the chain.The location of a receiver is computed based on the master and secondary stations time of arrival differences.To obtain the time of arrival, the reference point of each pulse must be detected.The detection accuracy of reference points has substantial effect on the estimation of the receiver's location.The process of obtaining of the reference point is called cycle identification and usually the third zero crossing point is chosen to be the reference point [1].Another important characteristic is the envelope to cycle difference (ECD) which indicates that the envelope of the pulse shifts in time.ECD can cause error in cycle identification process.The main focus of this paper is to develop algorithms to obtain nearest sample to the reference point, estimate the reference point time of arrival between two samples and estimate ECD.This paper modifies and improves our work in [7] and most of topics in [7] are covered.Also, the algorithm in [7] for cycle identification is modified to improve the performance in noisy environments.Alongside with that, two algorithms for estimation of the reference point between two samples are proposed which enables inter sample cycle identification so that, high sampling frequency is not required for high accuracy.The reminder of this paper is organized in the following manner, section 2 is dedicated to cycle identification, section 3 describes two proposed algorithms for intersample cycle identification.section 4 explains methods of ECD estimation and in the section 5 the simulations results are given for different SNRs, eventually the section 6 is the conclusion of the paper.

Cycle Identification
A typical Loran-C pulse is illustrated in Figure 1.The "*" in Figure 1 represents the location of the reference point.The equation of a Loran-C pulse is [1]- [5] [11]: Where t is time in µs, PC is the value of the phase coding which can be zero or π and it is used to distinguish between different slaves and A is the amplitude.Cycle identification is the process of detecting third zero crossing point of a Loran-C signal.This process has considerable effect on the accuracy of a receiver localization.In [3] the delay and sum method is used which uses addition of the signal and delayed version of it to obtain a suitable signal to find the third zero crossing point.This method is described in section 2.1.

Delay and Sum Method
The delay and sum method is also known as the Half Cycle Peak Ratio (HCPR) method.If y(n) is the received signal, the synthetic signal z(n) is created as [3]: Figure 2 shows z(n) and the zero crossing points which are marked with '+'.As it is illustrated in Figure 2, there is a phase reversal at the third zero crossing point which can easily be detected by sampling z(n-5) and z(n+5) of supposed zero crossing point until z(n-5) and z(n+5) have opposite signs [3].For more information about the delay and sum method see [3][1].To investigate another method which uses delayed lock loop (DLL) see [5].

Envelope Based Method for Cycle Identification Using Matched Filter
If s(t) is a Loran-C single pulse, the matched filter of the pulse is defined as [10], [11]: where T is the length of a typical Loran-C pulse.Input signal is defined in (4) and the output of the matched filter can be obtained by: (4) where y(t) is the input signal contaminated by additive Wight Gaussian noise.Substituting (3) into (4) results in: Equation ( 5) is the cross correlation (CC) of the input signal with a typical Loran-C signal and also, it can be described as convolution of the input signal with the matched filter.Cross correlation is used to find the time of arrival of signals in many applications [8][9].According to [10] when a signal is filtered by matched filter, the noise can be efficiently removed and higher SNR is achieved.To see another application of matched filter in Loran-C see [11].In [11] matched filter is used in acquisition stage to improve SNR.However, here matched filter is used for cycle identification.The block diagram of the proposed method is shown in Figure 3 and Q(n) signal is shown in Figure 4 (b).The Loran-C received signal is filtered with g(n) and Hilbert transform is used to obtain the envelope of the signal.The envelope of the signal is used because it does not have local peaks and also it does not depend on phase coding so that the output shape is identical for all Loran-C input signals.To remove spurious peaks caused by noise, a Hamming filter is used and the peaks of the obtained signal show the beginning of each pulse.The third zero crossing point is 30 µs after the beginning of a pulse.

Modification of the Envelope Based Method
There are two problems that limit the performance of envelope based method.The first problem of envelope based method is that although envelope detection removes carrier frequency and makes peak detection straightforward, samples of envelope near the peak sample have amplitudes near maximum value which increases probability of choosing wrong sample as the peak sample.Another disadvantage of the envelope based method is that nonzero ECD can change the estimated sample location so that it needs a post processing stage to reduce the effect of ECD.For more information about this problem see [7].These problems motivate us to find a way to find an alternative for the envelope detection.To improve the previous method, envelope detection and hamming filter are replaced by absolute function and GPD (Global Peak Detection).Block diagram of the modified cycle identification method is shown in Figure 5.The output of the algorithm (after absolute function) is shown in Figure 6 with a global peak at 2000.The output of matched filter cannot be used without absolute function because phase coding changes sample index belongs to maximum of the signal (In case of using envelope like the method described in section 2.2 the output of the matched filter is insensitive to the phase coding).Absolute function makes the output of the matched filter insensitive to the phase coding.There are many peaks in the output of the absolute function and only index of the global peak of each pulse matters so that GPD algorithm is used to find global peak of each pulse.GPD algorithm is an iterative algorithm to find global peak of each pulse and removing local peaks.GPD algorithm pseudo code is shown in Algorithm 1.The GPD algorithm finds the maximum of the signal and remove all local peaks (by setting to zero) near the maximum (local peaks that are closer to the maximum than L).The GPD algorithm repeats this action until all global peaks higher than a predefined threshold are detected.This method mitigates both problems (it is more robust to noise and ECD) and has low complexity so that it can be a suitable candidate for low cost, high performance receivers.

Intersample Cycle Identification
Algorithms described in section 2 only find the nearest sample to the reference location.The task of this section is to introduce two methods to find the exact location of the reference points.Suppose that the sampling frequency is 1 MHz therefore the time duration between two successive samples is 1 µs.Exact location of the third zero crossing can be anywhere between the two nearest sample therefore it can be considered to be a Gaussian random variable between two successive samples.By considering the nearest sample to the exact third zero crossing as the reference location, there will be 0.25 µs error on average which causes great amount of error in estimation of the position of the target.One naïve way to overcome this problem is to increase the sampling frequency.For instance, to reach 20 ns average error 12.5 MHz sampling frequency is required which increases complexity and cost of the system considerably.In this section two methods are proposed to find location of the third zero crossing between two samples without any need for high sampling frequency.In both algorithms it is assumed that the phase coding of the pulse is zero otherwise the input signal must be multiplied by -1.

Linear Model
Sinus function can be estimated by a line near the zero.The important point is that the two samples near zero value must be close enough to be able to use a line as an estimator consequently, the sampling frequency must be high enough to guarantee good performance of this method.At first, location of the nearest sampling point is found by cycle identification methods described in section 2. In the next step, two nearest sampling point to the zero crossing are found.In the end a line crossing these two sampling point is found using (6).
input : output signal of absolute function (y) output : index of global peak sample (GPD_index) variables: threshold : threshold to remove noise i : counter for the loop L : length of a typical Loran-C pulse i=0 while max(y)>threshold i = i+1; GPD_index(i)=index of max(y) set to zero y samples between GPD_index(i)-L, GPD_index(i)+L end Where x 0 and x 1 are the time of two successive samples near the third zero crossing and y 0 , y 1 are their corresponding amplitudes.The zero crossing is the location where y=0 so that reference time can be obtained by (8).
Second order polynomial function can be used instead of the linear function.However, according to our observation linear function suffice.

Mapping Input Signal to a Template Signal Using Least Square Optimization
This algorithm tries to map samples of the received signal to a template signal which has high sampling frequency and it is stored in receiver's memory.By mapping received signal samples in the template signal, accurate location of the third zero crossing can be identified.Suppose by using one of cycle identification methods, 4 nearest sample point to the third zero crossing are chosen so that the following equations can be written: where S is received signal, n 1 is the earliest sample, f s is the received signal sampling frequency and f t is the template signal frequency which is much higher than f s .According to (11), the distance of two successive samples in the input signal equals to M samples in the template signal.If two corresponding samples in the input signal and the template signal are found, finding other corresponding samples will be a trivial task.Suppose n i and nt j points to the same time in the Loran-C pulse (for example third zero crossing) other corresponding samples can be found by: where n q is a sample index of the input signal, nt q is a sample index of the template signal which corresponds to n q, n i and nt i are samples index of the input signal and the template signal respectively that correspond to each other's.To find two corresponding samples, amplitude information cannot be used because maximum amplitude is susceptible to noise and also, it is not highly probable that the input signal is sampled exactly at maximum value.The only available information is (10).( 10) is used to find two corresponding samples by minimizing the cost function below. ) In (13) P is the template (pattern) signal, Ci is obtained from four nearest samples of the input signal to the third zero crossing using (10), M is obtained by (11) and must be an integer.nt1 is the sample index in the template signal which corresponds to n1 (earliest sample in 4 nearest samples to the zero crossing).By minimizing (13), nt1 can be found.After finding nt1, exact third zero crossing time can be calculated using (14).
In (14), n r is the nearest sample index to the third zero crossing in the input signal, nt z is the nearest sample index to the reference point (third zero crossing) in the template signal and nt r is the sample index of the template signal which corresponds to the nearest sample index to the reference point of the input signal (n r ).The prominent advantage of this algorithm is that the accuracy of cycle identification does not depend very much on the received signal sampling frequency and it depends more on the template signal sampling frequency so that this algorithm has advantages over linear model when the received signal has low sampling frequency.Another point is that, formula of the Loran-C signal can be used to obtain the template signal (saved signal with high sampling frequency).However, if the received signal does not obey ideal Loran-C formula, template signal must be obtained from the received signal using interpolation.As an example to show the average of error, suppose the template signal sampling frequency is 10 MHz.The distance between two successive sample is 0.1 µs and maximum of error would be 50 ns thus the average of error is 25 ns.

ECD Estimation
ECD is the varied Loran-C signal's envelope.Loran-C signal in presence of ECD is obtained from (15) [6] where T is the ECD value.In normal conditions, the value of ECD is zero.When the signal travels, the velocity of the propagation changes slightly with frequency and some components of the Loran-C signal travels faster than the others which result in ECD [3]. Figure 7 shows two Loran-C pulse with different values of ECD.Many researchers have done research about ECD.In [6] a method to estimate ECD is introduced which uses optimization to find ECD.In [4] Peterson and Dewalt used the collected data of some flights to estimate the relations between conductivity and ECD.They obtained a curve which shows that increasing in conductivity decreases ECD in nonlinear form.This paper focuses on the obtaining value of ECD directly from the signal.This paper does not intend to model ECD, for ECD modeling see [4] [2].

ECD Estimation Based on Least Square Error (LSE)
This method is described in [6].At first, envelope of the signal is obtained using Hilbert transform (this method can be used without envelope too, however it limits the location of sampling points) after that minimum of a cost function is found by optimization algorithms like steepest descend or Newton-Raphson algorithm.The minimum of the cost function ( 16) is the estimated ECD [6].
In (16), y(s i ) is a sample from the given signal and N is the number of samples which are used for obtaining the cost function.The cost function can be minimized by only one sample (N=1).However, in applications several samples are used to obtain more accurate result.A and T are parameters to be optimized therefore it is a 2D optimization.The cost function ( 16) is a convex function which means that it has only one minimum point.

Proposed Method for ECD Estimation
At first, exact time of reference point is found so that time of any sample is known.In the next step, two different points of the signal (from single envelope pulse) are sampled according to equation (17) and equation ( 18) as [7]: In ( 17) and (18), only T and A are unknown (the phase coding is known).If equation ( 17) is divided by (18), in the obtained equation only T is unknown and after few mathematical manipulations, equation (19) can be obtained as [7]: )) ( 65 Solving equation ( 19) results in two answers, the nearest result to zero is the value of ECD.This method has low computational complexity and does not need the use of optimization algorithms.The values of C1 and C2 should be enough high to be robust in presence of noise.

Effect of ECD on Cycle Identification
If envelope of the output of the matched filter is used, nonzero ECD can change the estimated sample location.Another point is that cycle identification accuracy has considerable effect on ECD estimation.Envelope based method described in section 2.2 is not immune to ECD so that a post processing algorithm is used.Modified algorithm described in section 2.3 is robust to ECD provided that it is in standard range.Standard range of ECD is between 2.5 µs and -2.5 µs.signal with ECD in this range is usable because ECD effect can be corrected.To reduce effect of ECD on the reference point detection in envelope based method, the nearest zero crossing point to the result of the cycle identification will be the exact point of the reference point.For finding zero crossing point, we searched for a sample which is less than absolute value of previous and the next sample.By using this method, adverse effect of ECD on envelope based cycle identification is diminished.

Simulations and Results
In this section, simulation results of sections 2, 3 and 4 are given.To evaluate the performance of cycle identification methods in presence of white Gaussian noise, the results are simulated for different SNRs.Because of random nature of the noise, the results are evaluated 500 times and the root mean square error of the results are presented.Table 1 shows the values of the different parameters of envelope based algorithm and modified algorithm.The results of the cycle identification methods are given in Table 2 and Figure 8 (a).To obtain results of Table 2 it is assumed that the reference point distance to the nearest sample is zero (only sample error is considered) and ECD is zero as well.The results of Table 2 depict that new proposed method outperform the other methods and also it works well in low SNR condition.Envelope based method of [7] has acceptable results in 10 db and more.Delay and sum method fails to find the true reference points when SNR is 10 db or less.In 30 db and more all methods perform well therefore Delay and sum which has the lowest complexity is preferred.For evaluation of ECD, same procedure has been done and the results are shown in Table 3 and Figure 8 (b) with the known reference points.Results of Table 3 show that proposed algorithm works slightly better than LSE algorithm.Also, the proposed algorithm does not need optimization so that it has lower complexity.4 show that linear model works well for 1 MHz sampling frequency and it is comparable with template mapping 2.
Table 5 shows results of intersample cycle identification for 500 kHz sampling frequency.From results of Table 5 it is evident that linear model fails to find the exact time of the reference point and template mapping is preferred to the linear model.Table 4 and Table 5 illustrate that template mapping algorithm result approximately similar answers for both sampling frequencies.To obtain results of Table 4 and Table 5 the input signal is low noise and ECD is zero.

Conclusion
The cycle identification and ECD estimation accuracies have considerable effects on the Loran-C system receiver's localization accuracy so it is important to estimate the reference point and ECD as precise as possible.White Gaussian noise degrades the performance of cycle identification and ECD estimation.Also, the error in cycle identification or ECD estimation degrades the performance of another one.In this paper proposed methods for cycle identification, ECD estimation and accurate estimation of the reference points between two samples are evaluated and the relation between the cycle identification and ECD estimation is explained.The simulations results showed that the proposed method for cycle identification and ECD estimation are robust in presence of white Gaussian noise.The result of Table 2 and Table 3 showed that the proposed methods outperform the compared methods.Another point is that the proposed methods have acceptable results in presence of noise which makes them suitable to use in very low SNR conditions.Furthermore, two algorithms are proposed to estimate the exact location of the reference points which leads to high accuracy by using low sampling frequency.In this paper the sky waves and inferences issues were not taken into account.The proposed methods are not robust in presence of the sky waves.Also, the proposed algorithms for intersample cycle identification are vulnerable to noise.

Figure 2 :
Figure 2: The Synthetic signal z(n) and four zero crossing point.

Figure 3 :
Figure 3: The block diagram of the envelope based method.

Figure 5 :
Figure 5: The block diagram of the modified method.

Figure 6 :
Figure 6: Output of the absolute function.

Figure 7 :
Figure 7: Loran-C signal with three different ECD value.

Table 4
investigates results of intersample cycle identification for 1 MHz sampling frequency and compares linear model and template mapping method.Results are presented for different distances of the third zero crossing to the nearest sample in µs.Template mapping 1 uses a 20 MHz template signal and template mapping 2 uses a 40 MHz template signal.Results of Table

Table 1 :
The values of parameters.

Table 2 :
Sample based cycle identification evaluation.