METHOD FOR SEPARATING LORAN SKY AND GROUND WAVES BASED ON LEVENBERG-MARQUART ALGORITHM

20220317228 · 2022-10-06

    Inventors

    Cpc classification

    International classification

    Abstract

    Provided herein is a method for separating Loran sky and ground waves based on a Levenberg-Marquart algorithm, including: (1) collecting a plurality of Loran sky-ground wave signals followed by normalization to obtain a normalized signal; (2) preprocessing the normalized signal by inverse Fourier transform method to obtain an initialization parameter; (3) establishing a mathematical model for the Loran sky-ground wave signals in time domain; and (4) solving parameters of the mathematical model using the Levenberg-Marquart algorithm to separate the Loran sky and ground waves.

    Claims

    1. A method for separating Loran sky and ground waves based on a Levenberg-Marquart algorithm, comprising: (1) generating a plurality of simulated Loran sky-ground wave signals by simulation, and collecting a plurality of actual Loran sky-ground wave signals; and normalizing the plurality of simulated Loran sky-ground wave signals and the plurality of actual Loran sky-ground wave signals to obtain a normalized signal; (2) preprocessing the normalized signal to obtain an initialization parameter; (3) according to the initialization parameter, establishing a mathematical model for a Loran received signal in time domain; and (4) fitting parameters in the mathematical model using the Levenberg-Marquart algorithm to obtain a Loran sky-wave time delay, a Loran sky-wave amplitude, a Loran ground-wave time delay and a Loran ground-wave amplitude.

    2. The method of claim 1, wherein the step (1) is performed through steps of: (1.1) generating 1000 simulated Loran sky-ground wave signals under a signal-to-noise ratio of 15 dB; wherein each of the 1000 Loran sky-ground wave signals is composed of ground wave, sky wave and noise; a time delay of the ground wave is 77 μs, and a time delay of the sky wave is 133 μs; a signal ratio of the sky wave to the ground wave is 12 dB; a sampling frequency is 1 MHz; and a single pulse repetition period of the Loran signals is 1 ms; for each Loran sky-ground wave signal, 1000 sampling points are selected, and a total sampling time is 1000 μs; (1.2) generating 1000 simulated Loran sky-ground wave signals under a signal-to-noise ratio of 0 dB, 1000 simulated Loran sky-ground wave signals under a signal-to-noise ratio of 5 dB, 1000 simulated Loran sky-ground wave signals under a signal-to-noise ratio of 10 dB and 1000 simulated Loran sky-ground wave signals under a signal-to-noise ratio of 20 dB; (1.3) collecting the plurality of actual Loran sky-ground wave signals with a sampling frequency of 1 MHz; wherein a single pulse repetition period of each of the plurality of actual Loran sky-ground wave signals is 1 ms; for each of the plurality of actual Loran sky-ground wave signals, 1000 sampling points are selected, and a total sampling time is 1000 μs; and (1.4) normalizing a signal amplitude of a total of 5000 simulated Loran sky-ground wave signals and the plurality of actual Loran sky-ground wave signals to [−1,1] to obtain the normalized signal expressed as (t.sub.i, y.sub.i); wherein i=1,2 . . . m; and m is equal to 1000, and represents the number of samples for a single signal.

    3. The method of claim 1, wherein in step (2), the normalized signal is preprocessed by inverse Fourier transform method; and the initialization parameter comprises: the number of signal multipaths, an initial time delay and an initial amplitude of each multipath signal.

    4. The method of claim 1, wherein in step (3), considering that a Loran signal transmission process only experiences amplitude and delay changes, the mathematical model is established in the time domain according to the number of signal multipaths, expressed as: y ( t ) = .Math. n = 1 N a n q ( t - τ n ) ; and ( 1 ) q ( t ) = t 2 sin ( 0.2 π t ) exp ( - 2 t / 65 ) ; ( 2 ) wherein y(t) is a certain signal sample; q(t) is a standard enhanced Loran signal; t is time, in μs; α.sub.n is a normalized amplitude of a n.sup.th multipath signal; τ.sub.n is a time delay of the n.sup.th multipath signal, in μs; N is the number of signal multipaths; for the plurality of simulated Loran sky-ground wave signals, N is equal to 2; and for the plurality of actual Loran sky-ground wave signals, N is equal to 5.

    5. The method of claim 4, wherein in step (4), a nonlinear equation (1) is solved using the Levenberg-Marquart algorithm, wherein there are 2n parameters to be solved α.sub.1, α.sub.2, . . . , α.sub.n; τ.sub.1, τ.sub.2, . . . , τ.sub.n; let x=[α.sub.1, α.sub.2, . . . , α.sub.n; τ.sub.1, τ.sub.2, . . . , τ.sub.n]=[x.sub.1, x.sub.2, . . . , x.sub.2n]; a model function is expressed as φ(x; t); a residual r(x) is expressed as:
    r(x)=[r.sub.1(x),r.sub.2(x), . . . r.sub.m(x)]  (3); and
    r.sub.i(x)=φ(t.sub.i,x)−y.sub.i,i=1, . . . ,m  (4); an objective function of the Levenberg-Marquardt algorithm is expressed as: f ( x ) = 1 2 r ( x ) T r ( x ) = 1 2 .Math. i = 1 m [ r i ( x ) ] 2 ; ( 5 ) according to a principle of least squares, when solving parameters of the objective function to be fitted, it is necessary to make a residual sum of squares be a minimum value:
    min f(x)  (6); through the Levenberg-Marquart algorithm, multiple iterations are performed to make the parameters to be fitted infinitely approach optimal parameters making the min f (x); and an iteration incremental equation is expressed as:
    x.sub.k+1=x.sub.k−(h+λI).sup.−1g  (7); wherein k represents the current number of iterations; λ is a damping coefficient; I is an identity matrix; h(x.sub.k)=J(x.sub.k).sup.TJ(x.sub.k), g(x.sub.k)=J(x.sub.k).sup.Tr(x .sub.k), J(x.sub.k) is a Jacobi matrix of r(x.sub.k), expressed as: J ( x k ) = [ r 1 x 1 .Math. r 1 x 2 n .Math. .Math. r m x 1 .Math. r m x 2 n ] x k ; ( 8 ) the step (4) is performed through steps of: (4.1) according to the initialization parameter obtained in step (2), setting an initial point x.sub.k=x.sub.0, and setting the damping coefficient λ to be greater than 1, and a convergence accuracy ϵ to be greater than 0; (4.2) calculating r(x.sub.k) and f(x.sub.k); (4.3) calculating J(x.sub.k), h(x.sub.k) and g(x.sub.k); (4.4) solving the incremental equation (h(x.sub.k)+λI)Δx.sub.k=g(x.sub.k) to obtain an increment Δx.sub.k; (4.5) letting x.sub.k+1=x.sub.k+Δx.sub.k, and calculating f(x.sub.k+1); (4.6) if ∥Δx.sub.k∥.sub.2<ϵ, ending iterations; otherwise, proceeding to step (4.7); (4.7) if f(x.sub.k+1)<f(x.sub.k)+βg.sup.TΔx.sub.k, letting λ=λ/v, and proceeding to step (4.8); otherwise, letting λ=λv, repeating steps (4.4)-(4.7); wherein β,v are constraint variables; and (4.8) letting k=k+1, repeating steps (4.2)-(4.8); wherein β=0, v=2, λ=2, and ϵ=10.sup.−4.

    Description

    BRIEF DESCRIPTION OF THE DRAWINGS

    [0041] FIG. 1 is a flow chart of a Levenberg-Marquart algorithm according to an embodiment of the present disclosure;

    [0042] FIG. 2a is a histogram illustrating performance comparison of algorithms under simulated signals in step (1.1);

    [0043] FIG. 2b illustrates performance comparison of algorithms under different signal-to-noise ratios in step (1.2); where a statistical average of 1000 groups of data errors is taken as the delay error;

    [0044] FIG. 2c is a histogram showing performance comparison of algorithms under a signal-to-noise ratio of 0 dB in step (1.2); where the Levenberg-Marquart algorithm is represented by L-M; and

    [0045] FIG. 3 shows separation results of actual signals acquired in step (1.3).

    DETAILED DESCRIPTION OF EMBODIMENTS

    [0046] The present disclosure will be described completely and clearly below with reference to the accompanying drawings and embodiments to make the object, technical solutions, and beneficial effects of the present disclosure clearer.

    [0047] Provided herein is a method for separating Loran sky and ground waves based on a Levenberg-Marquart algorithm, which is performed through the following steps. [0048] (1) A plurality of simulated Loran sky-ground wave signals are generated by simulation, and a plurality of actual Loran sky-ground wave signals are collected. The plurality of simulated Loran sky-ground wave signals and the plurality of actual Loran sky-ground wave signals are normalized to obtain a normalized signal. [0049] (1.1) 1000 simulated Loran sky-ground wave signals are generated under a signal-to-noise ratio of 15 dB.

    [0050] Each of the 1000 Loran sky-ground wave signals is composed of ground wave, sky wave and noise. A time delay of the ground wave is 77 μs, and a time delay of the sky wave is 133 μs. A signal ratio of the sky wave to the ground wave is 12 dB. A sampling frequency is 1 MHz. A single pulse repetition period of the Loran signals is 1 ms. For each Loran sky-ground wave signal, 1000 sampling points are selected, and a total sampling time is 1000 μs. [0051] (1.2) 1000 simulated Loran sky-ground wave signals are generated under a signal-to-noise ratio of 0 dB. 1000 simulated Loran sky-ground wave signals are generated under a signal-to-noise ratio of 5 dB. 1000 simulated Loran sky-ground wave signals are generated under a signal-to-noise ratio of 10 dB. 1000 simulated Loran sky-ground wave signals are generated under a signal-to-noise ratio of 20 dB. [0052] (1.3) The plurality of actual Loran sky-ground wave signals are collected with a sampling frequency of 1 MHz. A single pulse repetition period of actual Loran signals is 1 ms. For each of the actual received Loran signals, 1000 sampling points are selected, and a total sampling time is 1000 μs. [0053] (1.4) A signal amplitude of a total of 5000 simulated Loran sky-ground wave signals and the plurality of actual Loran sky-ground wave signals is normalized to [−1,1] to obtain the normalized signal expressed as (t.sub.i, y.sub.i); where i=1, 2 . . . m; and m is equal to 1000, and represents the number of samples for a single signal. [0054] (2) The normalized signal is preprocessed by inverse Fourier transform method. An initialization parameter includes: the number of signal multipaths, an initial time delay and an initial amplitude of each multipath signal. [0055] (3) Considering that a Loran signal transmission process only experiences amplitude and delay changes, a mathematical model for the plurality groups of Loran sky-ground wave signals is established in a time domain according to the number of signal multipaths, expressed as:

    [00004] y ( t ) = .Math. n = 1 N a n q ( t - τ n ) ; and ( 1 ) q ( t ) = t 2 sin ( 0.2 π t ) exp ( - 2 t / 65 ) ; ( 2 ) [0056] where y(t) is a certain signal sample; q(t) is a standard enhanced Loran signal; t is time, in μs; α.sub.n is a normalized amplitude of the n.sup.th multipath signal; τ.sub.n is a time delay of the n.sup.th multipath signal, in μs; N is the number of signal multipaths; for the plurality of simulated Loran sky-ground wave signals, N is equal to 2; and for the plurality of actual Loran sky-ground wave signals, N is equal to 5. [0057] (4) A nonlinear equation (1) is solved using the Levenberg-Marquart algorithm. There are 2n parameters to be solved α.sub.1, α.sub.2, . . . , α.sub.n; τ.sub.1, τ.sub.2, . . . , τ.sub.n. let x=[α.sub.1, α.sub.2, . . . , α.sub.n; τ.sub.1,τ.sub.2, . . . , τ.sub.n]=[x.sub.1, x.sub.2, . . . , x.sub.2n]. A model function is expressed as φ(x; t). A residual r(x) is expressed as:


    r(x)=[r.sub.1(x), r.sub.2(x), . . . r.sub.m(x)]  (3); and


    r.sub.i(x)=φ(t.sub.i,x)−y.sub.i,i=1, . . . , m  (4).

    [0058] An objective function of the Levenberg-Marquardt algorithm is expressed as:

    [00005] f ( x ) = 1 2 r ( x ) T r ( x ) = 1 2 .Math. i = 1 m [ r i ( x ) ] 2 . ( 5 )

    [0059] According to the principle of least squares, when solving parameters of the objective function to be fitted, it is necessary to make a residual sum of squares be a minimum value:


    min f(x)  (6).

    [0060] Through the Levenberg-Marquart algorithm, multiple iterations are performed to make the parameters to be fitted infinitely approach optimal parameters, which make the min f (x). An iteration incremental equation of an iteration is expressed as:


    x.sub.k+1=x.sub.k−(h+λI).sup.−1g  (7);

    [0061] where k represents the current number of iterations; λ is a damping coefficient; I is an identity matrix; h(x.sub.k)=J(x.sub.k).sup.TJ(x.sub.k), g(x.sub.k)=(x.sub.k).sup.Tr(x.sub.k), J(x.sub.k) is a Jacobi matrix of r(x.sub.k), expressed as:

    [00006] J ( x k ) = [ r 1 x 1 .Math. r 1 x 2 n .Math. .Math. r m x 1 .Math. r m x 2 n ] x k . ( 8 )

    [0062] The step (4) is performed as follows. [0063] (4.1) According to the initialization parameter obtained in step (2), an initial point x.sub.k is set as x.sub.0; the damping coefficient A is set to be greater than 1; and a convergence accuracy ϵ is set to be greater than 0. [0064] (4.2) r(x.sub.k) and f(x.sub.k) are calculated. [0065] (4.3) J(x.sub.k), h(x.sub.k) and g(x.sub.k) re calculated. [0066] (4.4) The incremental equation (h(x.sub.k)+λI)Δx.sub.k=g(x.sub.k) is solved to obtain an increment Δx.sub.k. [0067] (4.5) Let x.sub.k+1=x.sub.k+Δx.sub.k, and f(x.sub.k+1) is calculated. [0068] (4.6) If ∥Δx.sub.k∥.sub.2<ϵ, iterations are ended; otherwise, step (4.7) is proceeded. [0069] (4.7) If f(x.sub.k+1)<f(x.sub.k)+βg.sup.TΔ.sub.k, let λ=λ/v, step (4.8) is proceeded; otherwise, let λ=λv, and steps (4.4)-(4.7) are repeated, where β,v are constraint variables. [0070] (4.8) Let k=k+1, and steps (4.2)-(4.8) are repeated, where β=0, v=2, λ=2, and ϵ=10.sup.−4.