AIS signal constrained least square angle estimation method, apparatus and device, and storage medium

Abstract

The invention provides an AIS signal constrained least square angle estimation method, device and equipment and a storage medium, and the method comprises the steps: collecting time domain, frequency domain and phase domain TDOA observation data, and evaluating the observation precision of each path based on a signal-to-noise ratio, goodness of fit and phase stability; then, mapping multi-source observation values into a unified nonlinear observation equation, and constructing a weighting matrix according to each observation precision, so that a high-precision channel obtains a larger weight; thirdly, forming an equation and inequality constraint set by combining the physical incidence angle range constraint of the double / three-antenna array, the arrival time difference range constraint determined by the baseline length and the propagation speed and the three-antenna closed TDOA geometric constraint; and finally, constructing an optimization target in a constrained least square form, carrying out iterative solution through a sequential quadratic programming (SQP) algorithm, outputting an optimal angle estimation value on the premise of meeting all physical and geometric constraints, and realizing high-precision estimation of the AIS signal arrival angle in a complex electromagnetic environment.

Description

Technical Field

[0001] The present invention relates to the field of AIS signal angle estimation, and in particular to an AIS signal constrained least squares angle estimation method, device, equipment and storage medium. Background Art

[0002] In modern maritime traffic management and ship navigation systems, the Automatic Identification System (AIS) has become a crucial technology for ensuring maritime safety and traffic regulation. Estimating the angle of AIS signals is a key technology for achieving functions such as ship positioning, track tracking, and collision warning. Its accuracy directly impacts the reliability and safety of maritime traffic management systems.

[0003] Traditional AIS signal angle estimation methods are mainly based on single-point estimation or simple peak detection algorithms. These methods have the following technical limitations: First, they lack statistical optimality guarantees and are not based on rigorous statistical theory design, which cannot guarantee the optimality of the estimation results. Second, the uncertainty quantification capability is insufficient, and reliable confidence intervals and accuracy assessments cannot be provided, making it difficult to meet the needs of safety-critical applications. Third, constraints are not fully utilized, and physical and geometric constraints are not effectively used to improve estimation performance and eliminate ill-conditioned solutions. Finally, the multi-observation information fusion method is simple, often using simple averaging or weighting methods, and lacks a scientific fusion strategy that considers differences in observation quality.

[0004] In view of this, this application is filed. Summary of the Invention

[0005] The present invention discloses an angle estimation method, apparatus, device and storage medium for AIS signal constrained least squares, aiming to achieve high-precision estimation of the arrival angle of AIS signals in a complex electromagnetic environment.

[0006] A first embodiment of the present invention provides an angle estimation method for an AIS signal using constrained least squares, including: Receive AIS signal arrival time difference observation data from multiple observation domains, and perform quality assessment on each observation data to obtain the measurement accuracy of each observation data; Establishing a unified observation equation based on the arrival time difference observation data, and constructing an observation weight matrix based on the measurement accuracy of each observation data, wherein the observation equation is used to express the functional relationship between the observation value and the arrival angle of the signal to be estimated; Establishing constraints based on the physical characteristics and geometric configuration of the antenna array, wherein the constraints include physical constraints that limit the physically feasible range of the estimated angle and geometric constraints for ensuring consistency of multi-antenna measurements; Combining the observation equation, the observation weight matrix, and the constraint conditions, constructing a constrained least squares optimization problem; The constrained optimization algorithm is used to solve the constructed optimization problem, and the angle estimation value that satisfies all constraints and minimizes the objective function is obtained through iterative calculation.

[0007] Preferably, the unified observation equation is established based on the arrival time difference observation data, specifically: The m observation values ​​are combined into an observation vector z=[z1,z2,...,z m ] T , where each observation zᵢ is the arrival time difference τ obtained by time domain cross correlation i _time, or arrival time difference τ obtained by frequency domain fitting i _freq or carrier phase to obtain the arrival time difference τ i _phase; Establish the theoretical model function h(θ)=[h1(θ),h2(θ),...,h m (θ)] T , where θ is the angle of arrival of the signal to be estimated; For a dual-antenna system, the theoretical model function is expressed as h i (θ)=(d / c)×sin(θ), where d is the antenna baseline length and c is the signal propagation speed; For a three-antenna system, the theoretical model function is expressed as h1(θ)=(d 12 / c)×sin(θ),h2(θ)=(d 13 / c)×sin(θ-α 13 )、h3(θ)=(d 23 / c)×sin(θ-α 23 ), where d ij is the baseline length between the i-th antenna and the j-th antenna, α ij is the angle between the corresponding baseline and the reference direction; Establish the observation equation z=h(θ)+ε, where ε is the observation error vector.

[0008] Preferably, the observation weight matrix is ​​constructed according to the measurement accuracy of each observation data, specifically: The quality of the observation data from the time domain, frequency domain and phase domain is evaluated separately: the measurement accuracy σ of the time domain observation data is i The measurement accuracy σ of the frequency domain observation data is determined according to the signal-to-noise ratio SNR and the interpolation quality. i The measurement accuracy σ of the phase domain observation data is determined based on the goodness of fit and coherence. i According to the phase stability and ambiguity confidence, the inverse square of the precision of each observation data is used to construct the initial weight matrix: W=diag([1 / σ1 2 ,1 / σ22 ,...,1 / σ m 2 ]), when there is correlation between observations, the weight matrix is ​​modified to W=R -1 , where R is the observation covariance matrix; During the iteration process, the weight w is adaptively adjusted through the robust weight function according to the size of the residual ᵢ i =w i0 ×robust_weight(residual i ).

[0009] Preferably, the constraints include: Angle range constraints: -90°≤θ≤90° for a two-antenna system, 0°≤θ<360° for a three-antenna system; Time difference of arrival range constraint: The absolute value of the time difference of arrival between any two antennas does not exceed the ratio of the baseline length to the signal propagation speed, which is expressed as |τ ij |≤d ij / c; Three-antenna closure constraint: The time difference from the first antenna to the third antenna is equal to the time difference from the first antenna to the second antenna plus the time difference from the second antenna to the third antenna, which is expressed as τ 13 =τ 12 +τ 23 ; The inequality constraint is expressed as g(θ)≤0, where g(θ)=[θ-θ max ,θ min -θ,|h(θ)|-τ max ] T ; The equality constraint is expressed as c(θ)=0, and for a three-antenna system, c(θ)=h3(θ)-h1(θ)-h2(θ).

[0010] Preferably, constructing a constrained least squares optimization problem includes: Define the weighted least squares objective function f(θ)=(zh(θ)) T W(zh(θ)), represents the weighted sum of squares of the differences between the observed values ​​and the theoretical model values; Construct the Lagrangian function L(θ,λ,μ)=f(θ)+λ T c(θ)+μ T g(θ), where λ is the Lagrange multiplier vector for equality constraints and μ is the Lagrange multiplier vector for inequality constraints; Establish KKT optimality conditions: gradient condition ∇θL=∇f+A T λ+B Tμ=0, equality constraint c(θ)=0, inequality constraint g(θ)≤0, non-negativity constraint μ≥0, complementary slack condition μ T g(θ)=0, where A=∂c / ∂θ is the Jacobian matrix for equality constraints, and B=∂g / ∂θ is the Jacobian matrix for inequality constraints.

[0011] Preferably, the constrained optimization algorithm is used to solve the constructed optimization problem, and the angle estimation value that satisfies all constraints and minimizes the objective function is obtained by iterative calculation, specifically: At the kth iteration, estimate θ based on the current angle k and its multiplier λ k 、μ k Construct a quadratic subproblem: min(1 / 2)d T H k d+∇f k T d, the constraint condition is A k d+c(θ k )=0 and B k d+g(θ k )≤0, where d=θ k+1 -θ k is the search direction, H k is the Lagrangian function L(θ,λ,μ) at θ k The Hessian matrix approximation matrix at ∇f k is the objective function gradient; Update the Hessian matrix using the BFGS formula: H k+1 =H k +(y k y k T ) / (y k T s k )-(H k s k s k T H k ) / (s k T H k s k ), where s k =θ k+1 -θ k ,y k =∇L(θ k+1 ,λ k+1 ,μ k+1 )-∇L(θ k ,λ k ,μ k ) Determine the step size α using the Armijo criterion for line search k , so that merit_function(θ k +α k d k )≤merit_function(θ k )+c1α k ∇merit_function k T d k , where c1 is the line search parameter, merit(·) is the comprehensive evaluation function combining the objective function and the penalty term; When the following convergence conditions are met at the same time, the iteration terminates and the optimal estimate is output: Gradient norm ||∇L||<ε0; Parameter increment ||θ k+1 -θ k ||<ε θ ; Objective function change |f(θ k+1 )-f(θ k )|<ε f ; Constraint satisfaction ||c(θ k )||<ε c And max(g(θ k ))<ε g .

[0012] Among them, ε0, ε θ , ε f , ε c and ε g are the preset convergence thresholds respectively.

[0013] A second embodiment of the present invention provides an AIS signal constrained least squares angle estimation device, comprising: A quality assessment unit is used to receive AIS signal arrival time difference observation data from multiple observation domains and perform quality assessment on each observation data to obtain the measurement accuracy of each observation data; an observation model construction unit, configured to establish a unified observation equation based on the arrival time difference observation data, and to construct an observation weight matrix based on the measurement accuracy of each observation data, wherein the observation equation is used to express the functional relationship between the observation value and the arrival angle of the signal to be estimated; A constraint condition establishing unit, configured to establish constraint conditions based on the physical characteristics and geometric configuration of the antenna array, wherein the constraint conditions include physical constraints that limit the physically feasible range of the estimated angle and geometric constraints for ensuring consistency of multi-antenna measurement; An optimization problem construction unit, configured to construct a constrained least squares optimization problem by combining the observation equation, the observation weight matrix, and the constraint conditions; The angle estimation value calculation unit is used to solve the constructed optimization problem using a constrained optimization algorithm, and iteratively calculate the angle estimation value that satisfies all constraints and minimizes the objective function.

[0014] A third embodiment of the present invention provides an AIS signal constrained least squares angle estimation device, including a memory and a processor, wherein the memory stores a computer program, and the computer program can be executed by the processor to implement an AIS signal constrained least squares angle estimation method as described in any one of the above.

[0015] A fourth embodiment of the present invention provides a computer-readable storage medium, characterized in that it stores a computer program, and the computer program can be executed by a processor of the device where the computer-readable storage medium is located to implement an AIS signal constrained least squares angle estimation method as described in any one of the above items.

[0016] The present invention provides an angle estimation method, apparatus, device, and storage medium for AIS signal constrained least squares. By collecting TDOA observation data in the time, frequency, and phase domains, the accuracy of each observation channel is evaluated based on signal-to-noise ratio, goodness-of-fit, and phase stability. The multi-source observations are then mapped into a unified nonlinear observation equation, and a weighting matrix is ​​constructed based on the accuracy of each observation, giving higher weight to the high-precision channel. A set of equality and inequality constraints is then formed by combining the physical incident angle range constraints of the dual / triple antenna array, the arrival time difference range constraints determined by the baseline length and propagation velocity, and the three-antenna closed TDOA geometric constraints. Finally, an optimization objective is constructed in the form of constrained least squares and iteratively solved using a sequential quadratic programming (SQP) algorithm. The optimal angle estimate is output while satisfying all physical and geometric constraints, achieving high-precision estimation of the AIS signal arrival angle in complex electromagnetic environments. BRIEF DESCRIPTION OF THE DRAWINGS

[0017] Figure 1 1 is a flow chart of an angle estimation method for an AIS signal using constrained least squares method provided by the first embodiment of the present invention; Figure 2 1 is a schematic diagram of a module of an AIS signal constrained least squares angle estimation device provided by a second embodiment of the present invention. DETAILED DESCRIPTION

[0018] The following will clearly and completely describe the technical solutions in the embodiments of the present invention in conjunction with the accompanying drawings. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making creative efforts are within the scope of protection of the present invention.

[0019] In order to better understand the technical solution of the present invention, the embodiments of the present invention are described in detail below with reference to the accompanying drawings.

[0020] It should be understood that the embodiments described are only a portion of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by persons of ordinary skill in the art without creative work are within the scope of protection of the present invention.

[0021] The terms used in the embodiments of the present invention are only for the purpose of describing specific embodiments and are not intended to limit the present invention. The singular forms "a", "an", "the" and "the" used in the embodiments of the present invention and the appended claims are also intended to include plural forms unless the context clearly indicates otherwise.

[0022] It should be understood that the term "and / or" as used herein is merely a description of the relationship between associated objects, indicating that three possible relationships exist. For example, "A and / or B" can represent: A exists alone, A and B exist simultaneously, or B exists alone. Furthermore, the character " / " in this document generally indicates that the associated objects are in an "or" relationship.

[0023] The word "if," as used herein, may be interpreted as "at the time of" or "when" or "in response to determining" or "in response to detecting," depending on the context. Similarly, the phrases "if it is determined" or "if (stated condition or event) is detected" may be interpreted as "when it is determined" or "in response to the determination" or "when detecting (stated condition or event)" or "in response to detecting (stated condition or event)," depending on the context.

[0024] The "first" and "second" mentioned in the embodiments are merely used to distinguish similar objects and do not represent a specific ordering of the objects. It is understood that the specific order or precedence of "first" and "second" can be interchanged where appropriate. It should be understood that the objects distinguished by "first" and "second" can be interchanged where appropriate, so that the embodiments described herein can be implemented in an order other than that illustrated or described herein.

[0025] The specific embodiments of the present invention are described in detail below with reference to the accompanying drawings.

[0026] The present invention discloses an angle estimation method, apparatus, device and storage medium for AIS signal constrained least squares, aiming to achieve high-precision estimation of the arrival angle of AIS signals in a complex electromagnetic environment.

[0027] A first embodiment of the present invention provides an angle estimation method for an AIS signal using constrained least squares, which can be performed by an estimation device, and in particular, by one or more processors within the estimation device, to implement at least the following steps: S101, receiving AIS signal arrival time difference observation data from multiple observation domains, and performing quality assessment on each observation data to obtain the measurement accuracy of each observation data; In this embodiment, the estimation device (or system) may be a terminal with data processing capabilities, such as a workstation, desktop computer, laptop computer, server, etc. The estimation device may be installed with a corresponding operating system and application software, and the functions required by this embodiment may be implemented through the combination of the operating system and application software.

[0028] In this embodiment, the estimation device is equipped with three receiving antennas, forming three baselines for measuring the time difference of arrival (TDOA) of AIS signals. When receiving an AIS signal from a ship, the estimation device simultaneously extracts TDOA observations from three different processing domains: time, frequency, and phase.

[0029] Specifically, the system first receives AIS signals through the antenna array. The signals received by each antenna are processed by the RF front end and then digitally sampled. For time domain processing, the system uses a cross-correlation algorithm to calculate the time delay between the signals received by different antennas, and determines the arrival time difference τ by finding the peak position of the cross-correlation function. i During this process, the system records the sharpness of the cross-correlation peak and the signal-to-noise ratio (SNR). When the SNR is high and the cross-correlation peak is sharp, it indicates that the quality of the time domain observation is good. The system comprehensively evaluates the measurement accuracy of the time domain observation based on the SNR value and the accuracy of the interpolation algorithm. i For example, when the SNR is 20 dB, combined with the accuracy of parabolic interpolation, the measurement accuracy can be obtained to be approximately 0.1 sampling period.

[0030] For frequency domain processing, the system converts the received signal into the frequency domain through fast Fourier transform and uses the phase spectrum information of the signal to estimate the delay. In the frequency domain, the phase of the cross-power spectrum of the signals received by the two antennas shows a linear change relationship, and its slope corresponds to the arrival time difference τ i _freq. The system extracts this slope through least squares linear fitting and calculates the goodness of fit R 2 and spectral coherence index. When the goodness of fit R 2When it is close to 1 and the coherence is high, it means that the frequency domain observation quality is good. The system determines the measurement accuracy σ of the frequency domain observation based on the fitting residual and the coherence coefficient. i _freq, which can usually achieve an accuracy of 0.05 sampling periods under high coherence conditions.

[0031] For phase domain processing, the system extracts the carrier phase information of the AIS signal and obtains high-precision arrival time difference τ through phase difference measurement. i _phase. Due to the 2π ambiguity problem in phase measurement, the system needs to first use a rough estimate in the time domain or frequency domain to resolve the ambiguity. When the phase is continuous and stable, the phase domain can provide the highest measurement accuracy. The system evaluates the measurement accuracy of phase domain observations by monitoring phase stability indicators such as phase change rate and Allan variance, as well as the confidence level of ambiguity resolution. i Under ideal conditions, the phase domain measurement accuracy can reach 0.01 sampling period.

[0032] In actual implementation, the system organizes all observation data into observation vector z=[τ 12 _time,τ 13 _time,τ 23 _time,τ 12 _freq,τ 13 _freq,τ 23 _freq,τ 12 _phase,τ 13 _phase,τ 23 _phase] T , where the subscripts represent antenna pairs, and a total of nine observations are included. Each observation has a corresponding measurement accuracy assessment result. Furthermore, the system implements a real-time monitoring mechanism for observation quality. When the quality indicator of an observation domain falls below a preset threshold, the system automatically reduces the weight of the observation or removes it, ensuring that only reliable observations contribute to the final angle estimation. For example, when the time-domain cross-correlation SNR falls below 10dB, the system deems the observation unreliable and sets its weight to 0 or uses robust estimation methods to reduce its impact.

[0033] S102: establishing a unified observation equation based on the arrival time difference observation data, and constructing an observation weight matrix based on the measurement accuracy of each observation data, wherein the observation equation is used to express the functional relationship between the observation value and the arrival angle of the signal to be estimated; It should be noted that the system organizes all time difference of arrival measurements from different processing domains and different antenna pairs into a comprehensive observation vector. In specific implementation, assume that the system is configured with three receiving antennas, forming three baselines (1-2, 1-3, 2-3). Each baseline can provide time difference of arrival measurements in all three processing domains (time domain, frequency domain, phase domain). The observation vector contains 9 elements, represented as z = [z1,z2,...,z9] T Among them, z1 to z3 correspond to the time domain arrival time difference τ of the three baselines respectively 12 _time,τ 13 _time,τ 23 _time, z4 to z6 correspond to the frequency domain measurement value τ 12 _freq,τ 13 _freq,τ 23 _freq, z7 to z9 correspond to the phase domain measurement value τ 12 _phase,τ 13 _phase,τ 23 _phase. Vectorized representation enables observation information from different domains to be fused and processed under the same mathematical framework.

[0034] In order to establish a mathematical relationship between the observed value and the angle to be estimated, the system constructs a theoretical model function h(θ). This function is based on the physical law of electromagnetic wave propagation and describes the theoretical arrival time difference that should be observed for each antenna when the signal is incident from the angle θ. In the three-antenna system of this embodiment, antenna 1 is set as the reference antenna and is located at the coordinate origin. Antennas 2 and 3 are located at different positions, forming a baseline length d 12 =10 m, d 13 =10 m, d 23 =10 meters equilateral triangle array. Baseline 1-2 is along the positive direction of the x-axis, and the angle α between baseline 1-3 and the x-axis 13 =60°, the angle between baseline 2-3 and x-axis is α 23 =120°.

[0035] The specific form of the theoretical model function is: for the theoretical time difference h of baseline 1-2 i (θ)=(d / c)×sin(θ), where c=3×10 8 Meters per second is the propagation speed of electromagnetic waves. For baselines 1-3, since there is a 60° angle between them and the reference direction, the theoretical time difference h2(θ)=(d 13 / c)×sin(θ-α 13 ); For baseline 2-3, the theoretical time difference h3(θ)=(d 23 / c)×sin(θ-α 23Since the theoretical value of each baseline in the three processing domains is the same, the complete theoretical model vector is: h(θ)=[h1(θ),h2(θ),h3(θ),h1(θ),h2(θ),h3(θ),h1(θ),h2(θ),h3(θ)] T .

[0036] The observation equation z = h(θ) + ε establishes the relationship between actual observations and the theoretical model, where ε represents the observation bias caused by various error factors, including thermal noise, multipath effects, and systematic errors. Ideally, if the true incident angle is known and there are no errors, the observed value should be exactly the same as the theoretical value. However, in reality, deviations always exist, and it is precisely by minimizing these deviations that the most likely incident angle is estimated.

[0037] When constructing the observation weight matrix, the system makes full use of the above-mentioned quality assessment results. Taking a specific measurement as an example, when the time domain measurement of baseline 1-2 shows an SNR of 25dB and good interpolation quality, the system evaluates its measurement accuracy σ1=0.08 sampling periods; while the frequency domain measurement of the same baseline has a goodness of fit R 2 =0.98, the coherence coefficient reaches 0.95, and a higher accuracy evaluation of σ4=0.05 sampling period is obtained; when the phase domain measurement is stable and the ambiguity is correctly resolved, the accuracy can reach σ7=0.01 sampling period.

[0038] Based on these accuracy assessments, the system constructs a diagonal weight matrix W = diag([1 / 0.08 2 ,1 / 0.10 2 ,1 / 0.09 2 ,1 / 0.05 2 ,1 / 0.06 2 ,1 / 0.05 2 ,1 / 0.01 2 ,1 / 0.015 2 ,1 / 0.012 2 This weight distribution ensures that the high-precision phase domain observations dominate the estimation, while the less accurate time domain observations still provide auxiliary information, achieving the optimal fusion of multi-domain information.

[0039] In actual processing, the system also considers the possible correlation between observations. For example, when the same signal segment is used for time domain and frequency domain processing, the errors of the two may be correlated. In this case, the system obtains the observation covariance matrix R by analyzing historical data or theoretical derivation, and corrects the weight matrix to W=R -1 , ensuring that dependencies are handled correctly.

[0040] More importantly, the system implements an adaptive weight adjustment mechanism. During the iterative solution process, each time a new angle estimate is obtained, the system calculates the residual of each observation: residualᵢ=zᵢ-hᵢ(θ). When the residual of an observation is abnormally large, it may indicate that the observation is interfered with or has a gross error. The system uses a Huber-type robust weight function. When |residual i |≤1.345σ i When | residual i |>1.345σ i When the weight is w i =w i0 ×1.345σ i / |residual i This adaptive mechanism enables the algorithm to automatically identify and suppress abnormal observations.

[0041] S103, establishing constraints based on the physical characteristics and geometric configuration of the antenna array, wherein the constraints include physical constraints that limit the physically feasible range of the estimated angle and geometric constraints for ensuring consistency of multi-antenna measurement; In this embodiment, the system establishes a constraint system based on the actual antenna array configuration and the propagation characteristics of the AIS signal. For a three-antenna equilateral triangle array, the physically feasible range of signal incidence angles is 0° to 360° for full coverage. This is because a three-antenna system eliminates the inherent front-back ambiguity of a two-antenna system. The system expresses this angle range constraint as two inequalities: θ - 360° ≤ 0 and 0° - θ ≤ 0, ensuring that the estimated angle remains within the valid range.

[0042] The physical constraint on arrival time difference stems from a fundamental fact: electromagnetic waves propagate at a finite speed. Therefore, the time difference between any two antennas receiving the same signal cannot exceed the time required for the signal to propagate the distance between the two antennas. In this embodiment, the lengths of the three baselines are all 10 meters, and the speed of light c = 3×10 8 m / s, so the maximum possible arrival time difference is 10 / (3×10 8 )≈33.3 nanoseconds. The system imposes a constraint on the arrival time difference of each observation: |τ 12 |≤33.3ns, |τ 13 |≤33.3ns, |τ 23 |≤33.3ns. These constraints apply not only to the actual observations, but also indirectly constrain the angle estimates through the theoretical model function h(θ).

[0043] The unique feature of the three-antenna system is its geometric closure property. When the signal is incident from a certain direction, the three arrival time differences must satisfy a vector closure relationship. Specifically, the time difference from antenna 1 to antenna 3 must be equal to the time difference from antenna 1 to antenna 2 plus the time difference from antenna 2 to antenna 3, that is, τ 13 =τ 12 +τ 23 This relationship holds true for measurements in the time, frequency, and phase domains. Ideally, this equality holds strictly; however, in actual measurements, slight deviations may occur due to various errors. By incorporating this relationship as an equality constraint into the optimization framework, the system enforces this geometric consistency in the estimated results.

[0044] To uniformly handle these constraints within the mathematical optimization framework, the system organizes all inequality constraints into a vector function g(θ). For the angular range constraint, g1(θ) = θ - 360° and g2(θ) = 0° - θ. For the arrival time difference constraint, since the theoretical model h(θ) calculates the theoretical time difference between each baseline, the constraint is expressed as g3(θ) = |h1(θ)| - 33.3 × 10 -9 、g4(θ)=|h2(θ)|-33.3×10 -9 、g5(θ)=|h3(θ)|-33.3×10 -9 The complete inequality constraint vector g(θ)=[g1(θ),g2(θ),g3(θ),g4(θ),g5(θ)] T , requiring g(θ)≤0.

[0045] The equality constraint is expressed by the function c(θ). In a three-antenna system, the geometric closure constraint requires h3(θ)-h1(θ)-h2(θ)=0. After expansion, it becomes: (d 23 / c)×sin(θ-120°)-(d 12 / c)×sin(θ)-(d 13 / c) × sin(θ - 60°) = 0. Since the three baselines in this embodiment are of equal length, this constraint simplifies to: sin(θ - 120°) - sin(θ) - sin(θ - 60°) = 0. This trigonometric equation ensures the geometric consistency of the estimated angle and the three time difference measurements.

[0046] In actual implementation, the system also takes into account the accuracy of numerical calculations. Due to the rounding error of floating-point operations, the strict equality constraint c(θ)=0 may never be accurately satisfied numerically. Therefore, the system sets a very small tolerance ε_c=10 -8 , when |c(θ)|<ε_c, the equality constraint is considered to be satisfied. Similarly, for the inequality constraint, when g_i(θ)<ε_g=10-8 The constraint is considered active when . It's important to note that the introduction of constraints improves the reliability of angle estimation. For example, when an observation is severely disturbed, resulting in an abnormal measurement value, without constraints, the least squares estimate may produce physically impossible results, such as angles exceeding 360° or time differences exceeding physical limits. With these constraints, the optimization algorithm automatically finds a solution that both closely fits the majority of normal observations and satisfies all physical and geometric constraints, providing a more robust estimate.

[0047] S104, constructing a constrained least squares optimization problem by combining the observation equation, the observation weight matrix, and the constraint conditions; It should be noted that after completing the establishment of the observation equation, the construction of the weight matrix and the setting of the constraints, all the aforementioned elements are integrated into a unified mathematical optimization framework to form a solvable constrained least squares problem. Specifically, the weighted least squares objective function f(θ) is first constructed. This function quantifies the degree of weighted deviation between the actual observation value and the theoretical prediction value at a given angle θ. In specific implementation, the system calculates the observation residual vector zh(θ), where z contains 9 actual measured arrival time differences, and h(θ) is the corresponding theoretical calculated value. For example, when θ=45°, the theoretical time difference of baseline 1-2 h1(45°)=(10 / 3×10 8 )×sin(45°)≈23.6 nanoseconds. If the actual measured value z1=24.1 nanoseconds, the residual of this observation is 0.5 nanoseconds.

[0048] Objective function f(θ)=(zh(θ)) T W(zh(θ)) performs a weighted summation of the squares of the residuals using the weight matrix W. In this embodiment, assuming that the phase-domain observation z7 has a weight of w7 = 10000 (corresponding to high-precision σ7 = 0.01), and the time-domain observation z1 has a weight of w1 = 156.25 (corresponding to lower-precision σ1 = 0.08), a 0.1 nanosecond phase-domain residual contributes as much to the objective function as a 0.8 nanosecond time-domain residual. This weighting mechanism ensures that the high-precision observation plays a dominant role in angle estimation.

[0049] However, simply minimizing the objective function f(θ) can lead to unreasonable results. For example, if all observations are affected by systematic biases in the same direction, unconstrained optimization may produce angle estimates that are outside the physical range. Therefore, the system incorporates constraints into the optimization framework through the Lagrangian method. The Lagrangian function L(θ,λ,μ)=f(θ)+λ T c(θ)+μ T g(θ) cleverly transforms the constrained optimization problem into finding the stable point of the function.

[0050] In this function, λ is the Lagrange multiplier associated with the equality constraint c(θ). For the geometric closure constraint c(θ)=h3(θ)-h1(θ)-h2(θ)=0 for a three-antenna system, the corresponding multiplier λ reflects the "strength" of the constraint. When the observed data naturally satisfies the closure relation well, the value of λ is small. Conversely, when the observations are inconsistent, the optimization process automatically adjusts the value of λ to strike a balance between satisfying the constraint and fitting the observations.

[0051] Similarly, the vector μ=[μ1,μ2,μ3,μ4,μ5] T Contains Lagrange multipliers corresponding to the five inequality constraints. These multipliers have a special property: when a constraint is not activated (e.g. g3(θ)=|h1(θ)|-33.3×10 -9 <0, that is, the time difference is much smaller than the physical limit), the corresponding μ3=0; only when the constraint approaches or reaches the boundary, the corresponding multiplier will take a positive value and affect the optimization process.

[0052] The KKT (Karush-Kuhn-Tucker) optimality condition established by the system provides the necessary and sufficient conditions for judging whether the solution is the optimal solution to the constrained optimization problem. Gradient condition ∇θL=∇f+A T λ+B T μ=0 requires that at the optimal point, the linear combination of the gradient of the objective function and the gradient of the constraint function is balanced. T W(zh(θ)), J is the Jacobian matrix of the observation function. For this embodiment, the first row of J is ∂h1 / ∂θ=(10 / 3×10 8 )×cos(θ), which reflects the effect of angle change on the baseline 1-2 time difference.

[0053] The Jacobian matrix of the equality constraint A=∂c / ∂θ describes how the closure constraint changes with angle. By using the chain rule, we get A=∂h3 / ∂θ-∂h1 / ∂θ-∂h2 / ∂θ=(10 / 3×10 8 )×[cos(θ-120°)-cos(θ)-cos(θ-60°)]. Similarly, the inequality constraint Jacobian matrix B contains the gradient information of each inequality constraint.

[0054] Complementary relaxation condition μ T g(θ)=0 embodies a subtlety of optimization theory: for each inequality constraint, either the constraint is inactive (g_i(θ)<0) and the corresponding multiplier is zero (μᵢ=0), or the constraint is on the boundary (g_i(θ)=0) and the multiplier may be positive. This ensures that only the constraints that actually matter affect the final solution.

[0055] S105 , using a constrained optimization algorithm to solve the constructed optimization problem, and iteratively calculating to obtain an angle estimation value that satisfies all constraints and minimizes the objective function.

[0056] In this example, assume the system is at iteration 5, with the current angle estimate θ5 = 32.46°, the equality constraint multiplier λ5 = -0.015, and the inequality constraint multiplier vector μ5 = [0, 0, 0.002, 0, 0]ᵀ (indicating that the third inequality constraint is nearing activation). Based on these current values, the system constructs a quadratic programming subproblem to determine the next search direction.

[0057] The construction process of the quadratic subproblem embodies the core idea of ​​sequential quadratic programming: approximate the original nonlinear objective function with a quadratic function at the current point, and approximate the nonlinear constraint with a linear constraint. Specifically, the system calculates the objective function gradient ∇f5=-2J5 at the current point T W(zh(θ5)), where J5 is the Jacobian matrix of the observation function at θ5. Assuming that the main component of the weighted residual vector W(zh(θ5)) comes from the phase domain observation, the gradient value is calculated as ∇f5=-3.2×10 -7 .

[0058] The Hessian matrix H5 approximates the second-order derivative information of the Lagrangian function, which determines the curvature characteristics of the quadratic subproblem. The system uses the BFGS method to maintain and update this matrix. From the 4th to the 5th iteration, the parameter changes s4=θ5-θ4=32.46°-32.38°=0.08°, and the gradient of the Lagrangian function changes from ∇L4=-2.8×10 -7 Change to ∇L5=-3.5×10 -7 , so y4=-0.7×10 -7 Applying the BFGS update formula, the system calculates the rank-two correction term: the first term (y4y4 T ) / (y4 T s T ) adds curvature information along the y4 direction, the second term - (H4s4s4 T H4) / (s4 T H4s4) removes outdated curvature information to obtain a more accurate Hessian matrix approximation H5.

[0059] Linearize the constraints at the current point. Equality constraint c(θ5)=h3(θ5)-h1(θ5)-h2(θ5)=2.1×10 -10 (minor violation), its Jacobian A5 = -0.85×10 -8 In the inequality constraint, since μ3>0, it means g3(θ5)=|h1(θ5)|-33.3×10 -9The linearized constraints require that the search direction d satisfies A5d+c(θ5)=0 (making the equality constraint tend to be satisfied) and B5d+g(θ5)≤0 (keeping the inequality constraint satisfied).

[0060] Solve the quadratic program min(1 / 2)dᵀH5d+∇f5 T d, the system uses the interior point method or the active set method to obtain the search direction d5 = 0.18°. This direction comprehensively considers the multiple requirements of objective function descent, equality constraint satisfaction, and inequality constraint feasibility.

[0061] However, taking a full step along the direction d5 may result in constraint violation or an increase in the objective function. The system determines the appropriate step size α5 through line search. The merit function merit(θ)=f(θ)+ρ||c(θ)||1+ρΣmax(0,g i (θ)) comprehensively evaluates the quality of the solution, where the penalty parameter ρ = 150 ensures that constraint violations are appropriately penalized. The system checks the value of the merit function at different step sizes: when α = 1, merit(θ5 + d5) = 0.0125, when α = 0.5, merit(θ5 + 0.5d5) = 0.0098, and when α = 0.25, merit(θ5 + 0.25d5) = 0.0089. Since merit(θ5) = 0.0092, α5 = 0.5 is selected to satisfy the Armijo criterion merit(θ5 + 0.5d5) ≤ merit(θ5 + 10 -4 ×0.5×(-0.012).

[0062] The updated value is θ6=θ5+0.5×0.18°=32.55°, and the system continues to iterate. After several iterations, at the 12th iteration, the system checks the convergence condition: the gradient norm of the Lagrangian function ||∇L 12 ||=7.3×10 -7 <10 -6 ; parameter change ||θ 13 -θ 12 ||=5.6×10 -9 <10 -8 ; The objective function remains basically unchanged |f(θ 13 )-f(θ 12 )|=8.9×10 -11 <10 -10 ; Equality constraint || c(θ 13 )||=6.2×10 -9 <10 -8 ; All inequality constraints max(g(θ 13 ))=-0.008<10 -8All convergence conditions are met at the same time, the system determines convergence and outputs the optimal angle estimate θ`=32.5847°.

[0063] In one possible implementation of the present invention: After obtaining the optimal angle estimate θ' and the corresponding Jacobian matrix J and weight matrix W, first pass F=J T WJ-∑λ i ∇ 2 c i -∑μ j ∇ 2 g j Calculate the Fisher information matrix F, where J is the partial derivative matrix of the observation equation with respect to θ, W is the observation weight matrix, and λ i 、μ j are the equality constraints c i and the inequality constraint g j The Lagrange multiplier of ∇ 2 c i 、∇ 2 g j is the second-order derivative matrix; then invert F to obtain the parameter covariance matrix C_θ=F -1 , and extract the standard error σ_θ=√diag(C_θ); under the normal distribution assumption, with the critical value z_{α / 2} at the significance level α, the (1–α) confidence interval is calculated as [θ`–z_{α / 2}·σ_θ,θ`+z_{α / 2}·σ_θ]; if the uncertainty of the angle function f(θ) needs to be evaluated, its standard error is calculated according to the error propagation formula σ_f=|∂f / ∂θ|·σ_θ.

[0064] In one possible implementation of the present invention: When performing statistical analysis on the optimal angle estimation result θ`, the residual vector z–h(θ`) is first standardized using the square root of the observation weight matrix W to obtain the standardized residual r=W 1 / 2 (z–h(θ`)), theoretically, it should approximately obey the standard normal distribution; then, by calculating the chi-square statistic χ 2 =(z–h(θ`))ᵀW(z–h(θ`)) and perform a chi-square test with degrees of freedom df=m–p–q (where m, p, and q represent the number of observations, the number of parameters to be estimated, and the number of equality constraints, respectively) to test the goodness of fit of the model to the overall observations. In order to identify the influence of a single observation point on the estimation result, this method further constructs a hat matrix H=J(J T WJ) -1 J T W, whose diagonal element h iiMeasures the leverage effect of the i-th observation on the entire fitting, and is based on the standardized residual r i and h ii Calculate Cook distance Di = (r i 2 / p)·\[h ii / (1–h ii ) 2 ], when D i When the threshold value of 4 / m is exceeded, the observation is determined to be an abnormal influence point, prompting the subsequent removal or recalibration of the point.

[0065] In one possible implementation of the present invention: In the specific implementation, this method further performs robust processing on the residual and matrix numerical stability: after completing the standardized residual r=W in each iteration 1 / 2 After calculating (z–h(θ`)), the residuals outside the range of 3σ are first marked as outliers and the Huber loss function is used to calculate the ρ(r)={r 2 / 2,|r|≤k;k|r|-k 2 / 2,|r|>k} (where the adjustment parameter k is 1.345) replaces the traditional quadratic loss and reweights the residuals by the influence function ψ, thereby updating the weights in the iterative reweighted least squares (IRLS) process w i (k+1) =ψ(r i (k) / σ) / (r i ( k) / σ) If at any iteration, the Jacobian matrix condition number cond(JᵀWJ) exceeds 10 12 , then introduce the Tikhonov regularization term λI (λ is the regularization parameter) into the normal equation and construct (J T WJ+λI)θ=J T Wz To ensure numerical stability; when the number of iterations exceeds 100 and the preset convergence condition is still not met, it will automatically switch to the trust region subproblem min f(θ k ) + ∇f k T d + (1 / 2)d T H k d By dynamically adjusting the trust region radius ||d|| ≤ Δ k Ensure the convergence of the algorithm and ultimately obtain accurate angle estimation that is both resistant to outlier interference and numerically stable.

[0066] A second embodiment of the present invention provides an AIS signal constrained least squares angle estimation device, comprising: A quality assessment unit is used to receive AIS signal arrival time difference observation data from multiple observation domains and perform quality assessment on each observation data to obtain the measurement accuracy of each observation data; an observation model construction unit, configured to establish a unified observation equation based on the arrival time difference observation data, and to construct an observation weight matrix based on the measurement accuracy of each observation data, wherein the observation equation is used to express the functional relationship between the observation value and the arrival angle of the signal to be estimated; A constraint condition establishing unit, configured to establish constraint conditions based on the physical characteristics and geometric configuration of the antenna array, wherein the constraint conditions include physical constraints that limit the physically feasible range of the estimated angle and geometric constraints for ensuring consistency of multi-antenna measurement; An optimization problem construction unit, configured to construct a constrained least squares optimization problem by combining the observation equation, the observation weight matrix, and the constraint conditions; The angle estimation value calculation unit is used to solve the constructed optimization problem using a constrained optimization algorithm, and iteratively calculate the angle estimation value that satisfies all constraints and minimizes the objective function.

[0067] A third embodiment of the present invention provides an AIS signal constrained least squares angle estimation device, including a memory and a processor, wherein the memory stores a computer program, and the computer program can be executed by the processor to implement an AIS signal constrained least squares angle estimation method as described in any one of the above.

[0068] A fourth embodiment of the present invention provides a computer-readable storage medium, characterized in that it stores a computer program, and the computer program can be executed by a processor of the device where the computer-readable storage medium is located to implement an AIS signal constrained least squares angle estimation method as described in any one of the above items.

[0069] The present invention provides an angle estimation method, apparatus, device, and storage medium for AIS signal constrained least squares. By collecting TDOA observation data in the time, frequency, and phase domains, the accuracy of each observation channel is evaluated based on signal-to-noise ratio, goodness-of-fit, and phase stability. The multi-source observations are then mapped into a unified nonlinear observation equation, and a weighting matrix is ​​constructed based on the accuracy of each observation, giving higher weight to the high-precision channel. A set of equality and inequality constraints is then formed by combining the physical incident angle range constraints of the dual / triple antenna array, the arrival time difference range constraints determined by the baseline length and propagation velocity, and the three-antenna closed TDOA geometric constraints. Finally, an optimization objective is constructed in the form of constrained least squares and iteratively solved using a sequential quadratic programming (SQP) algorithm. The optimal angle estimate is output while satisfying all physical and geometric constraints, achieving high-precision estimation of the AIS signal arrival angle in complex electromagnetic environments.

[0070] For example, the computer programs described in the third and fourth embodiments of the present invention may be divided into one or more modules, which are stored in the memory and executed by the processor to implement the present invention. The one or more modules may be a series of computer program instruction segments capable of performing specific functions, which describe the execution process of the computer program in the device for implementing constrained least squares angle estimation for AIS signals. For example, the apparatus described in the second embodiment of the present invention.

[0071] The processor may be a central processing unit (CPU), other general-purpose processors, digital signal processors (DSP), application-specific integrated circuits (ASIC), field-programmable gate arrays (FPGA), other programmable logic devices, discrete gate or transistor logic devices, discrete hardware components, etc. The general-purpose processor may be a microprocessor or any conventional processor. The processor serves as the control center of the AIS signal constrained least squares angle estimation method, and utilizes various interfaces and circuits to connect the various components of the AIS signal constrained least squares angle estimation method.

[0072] The memory can be used to store the computer programs and / or modules. The processor implements various functions of an AIS signal constrained least squares angle estimation method by running or executing the computer programs and / or modules stored in the memory and accessing data stored in the memory. The memory may primarily include a program storage area and a data storage area. The program storage area may store an operating system and at least one application required for a function (such as audio playback or text conversion); the data storage area may store data generated based on the use of the mobile phone (such as audio data and text message data). Furthermore, the memory may include high-speed random access memory (RAM) and non-volatile memory, such as a hard disk, internal memory, a plug-in hard disk, a Smart Media Card (SMC), a Secure Digital (SD) card, a Flash Card, at least one disk storage device, a flash memory device, or other volatile solid-state memory device.

[0073] If the implementation module is implemented as a software functional unit and sold or used as an independent product, it can be stored in a computer-readable storage medium. Based on this understanding, the present invention can also implement all or part of the process steps in the above-mentioned method embodiments by using a computer program to instruct the relevant hardware. The computer program can be stored in a computer-readable storage medium. When executed by a processor, the computer program can implement the steps of each of the above-mentioned method embodiments. The computer program includes computer program code, which can be in source code form, object code form, executable file, or some intermediate form. The computer-readable medium can include: any entity or device capable of carrying the computer program code, recording medium, USB flash drive, mobile hard drive, magnetic disk, optical disk, computer memory, read-only memory (ROM), random access memory (RAM), electric carrier signal, telecommunication signal, and software distribution medium. It should be noted that the content of the computer-readable medium can be appropriately increased or decreased based on the requirements of legislation and patent practice in a jurisdiction. For example, in some jurisdictions, according to legislation and patent practice, computer-readable media does not include electric carrier signals and telecommunication signals.

[0074] It should be noted that the device embodiments described above are merely illustrative, wherein the units described as separate components may or may not be physically separated, and the components displayed as units may or may not be physical units, that is, they may be located in one place, or they may be distributed across multiple network units. Some or all of the modules may be selected according to actual needs to achieve the purpose of the present embodiment. In addition, in the drawings of the device embodiments provided by the present invention, the connection relationship between the modules indicates that there is a communication connection between them, which may be specifically implemented as one or more communication buses or signal lines. A person of ordinary skill in the art can understand and implement the present invention without inventive effort.

[0075] The above description is merely a preferred embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any changes or substitutions that can be easily conceived by a person skilled in the art within the technical scope disclosed in the present invention should be included in the scope of protection of the present invention. Therefore, the scope of protection of the present invention should be based on the scope of protection of the claims.

Claims

1. An AIS signal constrained least squares angle estimation method, characterized in that: include: Receive AIS signal arrival time difference observation data from multiple observation domains, and perform quality assessment on each observation data to obtain the measurement accuracy of each observation data; Establishing a unified observation equation based on the arrival time difference observation data, and constructing an observation weight matrix based on the measurement accuracy of each observation data, wherein the observation equation is used to express the functional relationship between the observation value and the arrival angle of the signal to be estimated; Establishing constraints based on the physical characteristics and geometric configuration of the antenna array, wherein the constraints include physical constraints that limit the physically feasible range of the estimated angle and geometric constraints for ensuring consistency of multi-antenna measurements; Combining the observation equation, the observation weight matrix, and the constraint conditions, constructing a constrained least squares optimization problem; The constrained optimization algorithm is used to solve the constructed optimization problem, and the angle estimation value that satisfies all constraints and minimizes the objective function is obtained through iterative calculation.

2. The AIS signal constrained least squares angle estimation method according to claim 1, characterized in that: The unified observation equation is established based on the arrival time difference observation data, specifically: The m observation values ​​are combined into an observation vector z=[z1,z2,...,z m ] T , where each observation zᵢ is the arrival time difference τ obtained by time domain cross correlation i _time, or arrival time difference τ obtained by frequency domain fitting i _freq or carrier phase to obtain the arrival time difference τ i _phase; Establish the theoretical model function h(θ)=[h1(θ),h2(θ),...,h m (θ)] T , where θ is the angle of arrival of the signal to be estimated; For a dual-antenna system, the theoretical model function is expressed as h i (θ)=(d / c)×sin(θ), where d is the antenna baseline length and c is the signal propagation speed; For a three-antenna system, the theoretical model function is expressed as h1(θ)=(d 12 / c)×sin(θ),h2(θ)=(d 13 / c)×sin(θ-α 13 )、h3(θ)=(d 23 / c)×sin(θ-α 23 ), where d ij is the baseline length between the i-th antenna and the j-th antenna, α ij is the angle between the corresponding baseline and the reference direction; Establish the observation equation z=h(θ)+ε, where ε is the observation error vector.

3. The angle estimation method of an AIS signal constrained least squares method according to claim 1, characterized in that: The observation weight matrix is ​​constructed according to the measurement accuracy of each observation data, specifically: The quality of the observation data from the time domain, frequency domain and phase domain is evaluated separately: the measurement accuracy σ of the time domain observation data is i The measurement accuracy σ of the frequency domain observation data is determined according to the signal-to-noise ratio SNR and the interpolation quality. i The measurement accuracy σ of the phase domain observation data is determined based on the goodness of fit and coherence. i According to the phase stability and ambiguity confidence, the inverse square of the precision of each observation data is used to construct the initial weight matrix: W=diag([1 / σ1 2 ,1 / σ2 2 ,...,1 / σ m 2 ]), when there is correlation between observations, the weight matrix is ​​modified to W=R -1 , where R is the observation covariance matrix; During the iteration process, the weight w is adaptively adjusted through the robust weight function according to the size of the residual ᵢ i =w i0 ×robust_weight(residual i ).

4. The angle estimation method of an AIS signal constrained least squares method according to claim 1, characterized in that: The constraints include: Angle range constraints: -90°≤θ≤90° for a two-antenna system, 0°≤θ<360° for a three-antenna system; Time difference of arrival range constraint: The absolute value of the time difference of arrival between any two antennas does not exceed the ratio of the baseline length to the signal propagation speed, which is expressed as |τ ij |≤d ij / c; Three-antenna closure constraint: The time difference from the first antenna to the third antenna is equal to the time difference from the first antenna to the second antenna plus the time difference from the second antenna to the third antenna, which is expressed as τ 13 =τ 12 +τ 23 ; The inequality constraint is expressed as g(θ)≤0, where g(θ)=[θ-θ max ,θ min -θ,|h(θ)|-τ max ] T ; The equality constraint is expressed as c(θ)=0, and for a three-antenna system, c(θ)=h3(θ)-h1(θ)-h2(θ).

5. The AIS signal constrained least squares angle estimation method according to claim 1, characterized in that: The constrained least squares optimization problem is constructed as follows: Define the weighted least squares objective function f(θ)=(zh(θ)) T W(zh(θ)), represents the weighted sum of squares of the differences between the observed values ​​and the theoretical model values; Construct the Lagrangian function L(θ,λ,μ)=f(θ)+λ T c(θ)+μ T g(θ), where λ is the Lagrange multiplier vector for equality constraints and μ is the Lagrange multiplier vector for inequality constraints; Establish KKT optimality conditions: gradient condition ∇θL=∇f+A T λ+B T μ=0, equality constraint c(θ)=0, inequality constraint g(θ)≤0, non-negativity constraint μ≥0, complementary slack condition μ T g(θ)=0, where A=∂c / ∂θ is the Jacobian matrix for equality constraints, and B=∂g / ∂θ is the Jacobian matrix for inequality constraints.

6. The AIS signal constrained least squares angle estimation method according to claim 5, characterized in that: The constrained optimization algorithm is used to solve the constructed optimization problem, and the angle estimation value that satisfies all constraints and minimizes the objective function is obtained by iterative calculation, specifically: At the kth iteration, estimate θ based on the current angle k and its multiplier λ k 、μ k Construct a quadratic subproblem: min(1 / 2)d T H k d+∇f k T d, the constraint condition is A k d+c(θ k )=0 and B k d+g(θ k )≤0, where d=θ k+1 -θ k is the search direction, H k is the Lagrangian function L(θ,λ,μ) at θ k The Hessian matrix approximation matrix at ∇f k is the objective function gradient; Update the Hessian matrix using the BFGS formula: H k+1 =H k +(y k y k T ) / (y k T s k )-(H k s k s k T H k ) / (s k T H k s k ), among them k =θ k+1 -θ k ,y k =∇L(θ k+1 ,l k+1 ,m k+1 )-∇L(θ k ,l k ,m k ) Determine the step size α using the Armijo criterion for line search k , so that merit_function(θ k +α k d k )≤merit_function(θ k )+c1α k ∇merit_function k T d k , where c1 is the line search parameter, merit(·) is the comprehensive evaluation function combining the objective function and the penalty term; When the following convergence conditions are met at the same time, the iteration terminates and the optimal estimate is output: Gradient norm ||∇L||<ε0; Parameter increment ||θ k+1 -θ k ||<ε θ ; Objective function change |f(θ k+1 )-f(θ k )|<ε f ; Constraint satisfaction ||c(θ k )||<ε c And max(g(θ k ))<ε g . Among them, ε0, ε θ , ε f , ε c and ε g are the preset convergence thresholds respectively.

7. An AIS signal constrained least squares angle estimation device, characterized in that: include: A quality assessment unit is used to receive AIS signal arrival time difference observation data from multiple observation domains and perform quality assessment on each observation data to obtain the measurement accuracy of each observation data; an observation model construction unit, configured to establish a unified observation equation based on the arrival time difference observation data, and to construct an observation weight matrix based on the measurement accuracy of each observation data, wherein the observation equation is used to express the functional relationship between the observation value and the arrival angle of the signal to be estimated; A constraint condition establishing unit, configured to establish constraint conditions based on the physical characteristics and geometric configuration of the antenna array, wherein the constraint conditions include physical constraints that limit the physically feasible range of the estimated angle and geometric constraints for ensuring consistency of multi-antenna measurement; An optimization problem construction unit, configured to construct a constrained least squares optimization problem by combining the observation equation, the observation weight matrix, and the constraint conditions; The angle estimation value calculation unit is used to solve the constructed optimization problem using a constrained optimization algorithm, and iteratively calculate the angle estimation value that satisfies all constraints and minimizes the objective function.

8. An AIS signal constrained least squares angle estimation device, characterized in that: The invention comprises a memory and a processor, wherein the memory stores a computer program, and the computer program can be executed by the processor to implement the angle estimation method of the AIS signal constrained least squares according to any one of claims 1 to 6.

9. A computer-readable storage medium, characterized in that A computer program is stored, and the computer program can be executed by a processor of the device where the computer-readable storage medium is located to implement the angle estimation method of AIS signal constrained least squares according to any one of claims 1 to 6.

Images