[go: up one dir, main page]

CN111193528A - Gaussian filtering method based on nonlinear network system under non-ideal conditions - Google Patents

Gaussian filtering method based on nonlinear network system under non-ideal conditions Download PDF

Info

Publication number
CN111193528A
CN111193528A CN201911400580.2A CN201911400580A CN111193528A CN 111193528 A CN111193528 A CN 111193528A CN 201911400580 A CN201911400580 A CN 201911400580A CN 111193528 A CN111193528 A CN 111193528A
Authority
CN
China
Prior art keywords
time
formula
state
matrix
model
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201911400580.2A
Other languages
Chinese (zh)
Other versions
CN111193528B (en
Inventor
宋申民
赵凯
张秀杰
谭立国
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Harbin Institute of Technology Shenzhen
Original Assignee
Harbin Institute of Technology Shenzhen
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Harbin Institute of Technology Shenzhen filed Critical Harbin Institute of Technology Shenzhen
Priority to CN201911400580.2A priority Critical patent/CN111193528B/en
Publication of CN111193528A publication Critical patent/CN111193528A/en
Application granted granted Critical
Publication of CN111193528B publication Critical patent/CN111193528B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • HELECTRICITY
    • H04ELECTRIC COMMUNICATION TECHNIQUE
    • H04BTRANSMISSION
    • H04B1/00Details of transmission systems, not covered by a single one of groups H04B3/00 - H04B13/00; Details of transmission systems not characterised by the medium used for transmission
    • H04B1/69Spread spectrum techniques
    • H04B1/707Spread spectrum techniques using direct sequence modulation
    • H04B1/7097Interference-related aspects
    • H04B1/7103Interference-related aspects the interference being multiple access interference
    • H04B1/7105Joint detection techniques, e.g. linear detectors
    • HELECTRICITY
    • H04ELECTRIC COMMUNICATION TECHNIQUE
    • H04BTRANSMISSION
    • H04B17/00Monitoring; Testing
    • H04B17/30Monitoring; Testing of propagation channels
    • H04B17/309Measuring or estimating channel quality parameters
    • HELECTRICITY
    • H04ELECTRIC COMMUNICATION TECHNIQUE
    • H04LTRANSMISSION OF DIGITAL INFORMATION, e.g. TELEGRAPHIC COMMUNICATION
    • H04L43/00Arrangements for monitoring or testing data switching networks
    • H04L43/08Monitoring or testing based on specific metrics, e.g. QoS, energy consumption or environmental parameters
    • H04L43/0823Errors, e.g. transmission errors
    • H04L43/0829Packet loss
    • HELECTRICITY
    • H04ELECTRIC COMMUNICATION TECHNIQUE
    • H04LTRANSMISSION OF DIGITAL INFORMATION, e.g. TELEGRAPHIC COMMUNICATION
    • H04L43/00Arrangements for monitoring or testing data switching networks
    • H04L43/08Monitoring or testing based on specific metrics, e.g. QoS, energy consumption or environmental parameters
    • H04L43/0852Delays
    • HELECTRICITY
    • H04ELECTRIC COMMUNICATION TECHNIQUE
    • H04LTRANSMISSION OF DIGITAL INFORMATION, e.g. TELEGRAPHIC COMMUNICATION
    • H04L43/00Arrangements for monitoring or testing data switching networks
    • H04L43/50Testing arrangements

Landscapes

  • Engineering & Computer Science (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Signal Processing (AREA)
  • Environmental & Geological Engineering (AREA)
  • Quality & Reliability (AREA)
  • Physics & Mathematics (AREA)
  • Electromagnetism (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

The invention relates to a Gaussian filtering method of a nonlinear network system based on a Gaussian filtering method of the nonlinear network system under a non-ideal condition. The invention aims to solve the problems that the existing method does not consider the related noise, one-step random delay measurement and data packet loss which may occur in a nonlinear network system, and the problem that the estimation accuracy of a filter is reduced or even diverged due to model-based linear approximation or neglect of delay measurement. The Gaussian filtering method based on the nonlinear network system under the nonideal condition comprises the following steps: firstly, establishing a system model and a sensor measurement model; step two, providing hypothesis and lemma; thirdly, designing a Gaussian filter based on the second step; and step four, based on a third-order sphere diameter volume rule, approximating the Gaussian weighted integral in the step three to obtain a numerical form of the designed filter. The invention can be applied to the technical field of spacecraft and aircraft navigation.

Description

Gaussian filtering method based on non-linear network system under non-ideal condition
Technical Field
The invention relates to a Gaussian filtering method of a nonlinear network system, in particular to a state estimation method for a nonlinear system with correlated noise, one-step random delay measurement and data packet loss.
Background
In recent years, the estimation problem of network systems has attracted extensive attention[1-3]([1]L.Schenato,“Optimalestimation in networked control systems subject to random delay and packetdrop,”IEEE transactions on automatic control,vol.53,no.5,pp.1311,2008.
[2]W.A.Zhang,L.Yu,G.Feng,“Optimal linear estimation for networkedsystems with communication constraints,”Automatica,vol.47,no.9,pp.1992-2000,2011.
[3]R.Caballero-
Figure RE-GDA0002450881680000011
Hermoso-Carazo, J.Linares-P rez, "Optimal station estimation for network Systems with random parameter matrices, calibrated responses and delayed measurements," International Journal of General Systems, vol.44, No.2, pp.142-154, 2015.). Kalman Filter (KF)[4](R.E.Kalman, "A new approach filtering and prediction schemes," Journal of basic Engineering, vol.82, No.1, pp.35-45, 1960.) is an effective method for solving the state estimation problem. However, it is only suitable for ideal linear systems. In fact, non-linearity, correlated noise, random delay and packet loss are ubiquitous. In this context, the estimation problem of a non-linear network system with synchronous correlated noise and one-step random delay measurement and multi-packet loss is considered.
For non-linear systems, there are many estimation methods. For systems with weak non-linearity, the Extended KF (EKF) can achieve acceptable accuracy based on a first order Taylor series expansion[5](Y.Bar-Shalom,X.R.Li,T.Kirubarajan,“Estimation with applications to tracking and navigation:theory algorithmsand software,” John Wiley&Sons, 2004.). Differential Filters (DDF) unlike EKF algorithms which require computation of the Jacobian matrix[6](M.
Figure RE-GDA0002450881680000012
N.k. poulsen, o.ravn, "New definitions in statistics for nonlinear systems," Automatica, vol.36, No.11, pp.1627-1638, 2000.) the problem of the EKF algorithm possibly falling into local linearization is overcome by approximating the nonlinear function using stirling interpolation. Based on the fact that it is easier to approximate a probability distribution than to approximate any non-linear function or transformation, a series of sigma point filters, such as Particle Filters (PF) have been proposed[7](A.Doucet, S.Godsill, C.Andrieu, "On sequential Monte Carlo sampling methods for Bayesian filtering," statics and computing, vol.10, No.3, pp.197-208, 2000.), Unscented KF (UKF)[8](S.J.Julie, J.K.Uhlmann, "unknown filtering and nonlinear estimation," Proceedings of the IEEE, vol.92, No.3, pp.401-422, 2004.) and the volume KF (CKF)[9](I.Arasaratnam, S.Haykin, "Cubauture kalman filters," IEEEtransformations on automatic control, vol.54, No.6, pp.1254-1269, 2009.). Meanwhile, in order to adapt to different application scenarios, many improved versions of the above algorithm are proposed. For example, EKF based on second order taylor expansion, unscented PF, DDF with second order stirling interpolation, adaptive UKF, higher order CKF, square root CKF, etc. Because the UKF and the CKF have simpler forms and the numerical stability of the CKF algorithm is much higher than that of the UKF, the CKF algorithm is widely applied to the field of practical engineering, such as target tracking and mobile terminal positioning.
For systems with correlated noise, in[10](X.X.Wang, Y.Liang, Q.Pan, et al, "A Gaussian adaptive regenerative filter for nonlinear systems with coated chemicals," Automatica, vol.48, No.9, pp.2290-2297, September, 2012.), a new pseudo process equation was reconstructed in which the process noise is uncorrelated with the observation noise. In that[11](X.X.Wang, Y.Liang, Q.Pan, et., "Design and implementation of Gaussian filter for nonlinear system with linear deleted measures and correlated lipids," Applied Mathematics and combinatorial, vol.232, pp. 1011-1024, 2014.), a baseIn the projection theorem, a GASF of a nonlinear system is proposed. Then, [11 ]]The method of (1) for designing a Gaussian filter for a nonlinear system with random delay measurements and associated noise[12](G.B. Chang, "Comments on" A Gaussian adaptive reconstruction filters for nonlinear systems with corrected colors, "Automatica, vol.50, No.2, pp.655-656, February, 2014.). In that[13](G.B.Chang,“Alternative formulation of the Kalmanfilter for correlated process and observation noise,”IET Science,Measurement&Technology, vol.8, No.5, pp.310-318, September, 2014.) [13]And [14]The two filtering algorithms in (a) have proven to be theoretically equivalent in a linear system. In that[14](H.Yu,X.J.Zhang,S.Wang,et al.,“Alternative framework of the Gaussian filter for non-linear systems withsynchronously correlated noises,”IET Science, Measurement&Technology, vol.10, No.4, pp.306-315, July, 2016.) two alternative frameworks were developed based on state enhancement and conditional Gaussian distributions and demonstrated [14 []And [10-11]The equivalence of the algorithm (c). Then, at [15 ]](K.Zhao,S.M.Song,“Alternativeframework of the Gaussian filter for non-linear systems with randomly delayedmeasurements and correlated noises,”IET Science, Measurement&Technology, vol.12, No.2, pp.161-168, 2017.) and [16](S.L.Sun, L.Xie, W.Xiao, et al, "optimal estimation for systems with multiple packet routes," Automatica, vol.44, No.5, pp.1333-1342, May,2008.) the alternative framework of conditional Gaussian distributions is extended to nonlinear systems with correlated noise and one-step delay measurements.
For systems with one-step random delay measurement, the delay is typically described using a bernoulli-distributed random variable. For systems with multiple packet losses, zero input compensation, hold input compensation, and prediction compensation are separately proposed. For the problem of observation packets sent from the sensor to the filter, one or more packets arriving at the filter within one sampling period are discussed in [17] (S.L. Sun, "Optimal line filters for discrete-time Systems with distributed delay and lost measures with/without time stages," IEEETransactions on Automatic Control, vol.58, No.6, pp.1551-1556, June,2013.) and [18] (S.L. Sun, G.H. Wang, "model and estimation for networked Systems with distributed delay and packets," Systems & Systems, vol.73, pp.6-16, 20146, 201416), respectively. In [19] (S.L. Sun, J.Ma, "Linear timing for network control systems with random transmission delay and drop routes," Information Sciences, vol.269, pp.349-365, June,2014.), a more general case is considered in which data packets are transmitted from the sensor and controller to the filter and actuator, respectively. In the above case, when packet loss occurs, the measurement is invalid. To avoid this, a hold-in compensation mechanism was developed. In [20] (y.liang, t.w.chen, q.pan, "Optimal linear state estimator with multiple packet tdroutes," IEEE Transactions on Automatic Control, vol.55, No.6, pp.1428-1433, June,2010.), in the case where there is a plurality of packets lost, input compensation is kept for the best linear estimation problem. In [21] (j.ma, s.l.sun, "Information Fusion estimators for systems with multiplex sensors of differential packet drop rates," Information Fusion, vol.12, No.3, pp.213-222, July,2011.), the same approach is considered in the centralized and distributed Fusion estimation problem in the Linear Minimum Variance (LMV) sense. In [22] (S.L.Sun, L.Xie, W.Xiao, et al, "Optimal linearity for systems with multiple packet routes," Automatica, vol.44, No.5, pp.1333-1342, May,2008.), an Optimal linear estimator was developed based on a novel packet loss model. However, keeping the input compensation does not take into account the latest information of the system, and therefore a new compensation mechanism is proposed in which the latest observed predicted value is used as compensation. In [23] (J.Ding, S.L.Sun, J.Ma, N.Li, "Fusion optimization for Multi-Sensor network Systems with Packet Loss Compensation," Information Fusion, vol.45, pp.138-149, January,2019.), a centralized and distributed Fusion estimator was designed based on a new Compensation mechanism. In [24] (e.i. silver, m.a. solis, "An alternative look at both constant-gain Kalman filter for state estimation over channels," IEEE Transactions on Automatic Control, vol.58, No.12, pp.3259-3265, Decumber, 2013.), a type of handover estimator is considered in which missing data is replaced by the best estimate.
It is worth mentioning that in [25-27] ([25] J.Ma, S.L.Sun, "Linear estimators for network systems with one-step random delay and multiple packet drop based compression compensation," IET Signal Processing, vol.11, No.2, pp.197-204, April,2017.
[26]C.Zhu,Y.Xia,L.Xie,et al.,“Optimal linear estimation for systemswith transmission delays and packet dropouts,”IET signal Processing,vol.7,no.9,pp.814-823,December, 2013.
[27] L.l.sun, "Optimal linear estimators for discrete-time systems with one-step random delays and multiple packet routes," Acta automatic Sinica, vol.38, No.3, pp.349-354, March,2012 "), a one-step random delay measurement and multiple packet losses were considered simultaneously in three different models. In [27], possible packet loss, delay measurement and compensation values are described by introducing two uncorrelated bernoulli random variables, which ensure that the system can receive one of the three quantities as a measurement value at any time in order to maintain the accuracy of the filtering algorithm. However, in [27], it does not consider the case where the delay measurement and the real-time measurement arrive in synchronization. In [26], the measurement model is changed because adjacent measurements may arrive at the same time, where the measurement received at the most recent time in the data processing center is used as a compensation value for the current epoch, and the superiority of the algorithm in [26] compared to [27] is shown in [26 ]. Further, in [25], a new measurement model is proposed by changing the compensation mechanism in which the measurements of the current epoch are used as a compensator instead of the one-step prediction of the latest measurements received by the data processing center. It has been demonstrated that the algorithm proposed in [25] has higher estimation accuracy and less computational burden than [26], since the model in [25] always uses the latest measurement information.
Because the existing method does not consider that the data processing center can receive two adjacent measurements simultaneously in the nonlinear system, information loss may occur when the actual system processes such problems, thereby affecting the estimation accuracy of the system, and even possibly causing filter divergence.
Disclosure of Invention
The invention aims to solve the problems that the existing method does not consider the related noise, one-step random delay measurement and data packet loss which may occur in a nonlinear network system, and the problem that the estimation accuracy of a filter is reduced or even diverged due to model linear approximation or neglect of delay measurement, and provides a Gaussian filtering method based on the nonlinear network system under the nonideal condition.
The Gaussian filtering method based on the non-linear network system under the non-ideal condition comprises the following specific processes:
firstly, establishing a system model and a sensor measurement model;
step two, providing hypothesis and lemma;
thirdly, designing a Gaussian filter based on the second step;
and step four, based on a third-order sphere diameter volume rule, approximating the Gaussian weighted integral in the step three to obtain a numerical form of the designed filter.
The invention has the beneficial effects that:
the invention provides a state estimation method of a nonlinear network system with synchronous correlated noise, one-step random delay and a plurality of data packet losses, which considers that the system is a general nonlinear system, designs a Gaussian recursive filtering algorithm aiming at the synchronous correlated noise, one-step random delay measurement and data packet losses which may occur in the system, can ensure that the system obtains a high-precision estimation value, avoids the divergence of the system and ensures the stability of the system.
Drawings
FIG. 1 is a diagram of the estimation of the state in the UNGM model by the algorithm of the present invention and the algorithm of document [23 ];
FIG. 2 is a diagram of the estimated root mean square error of the state in the UNGM model for the algorithm of the present invention and the algorithm in document [23 ];
FIG. 3 is a graph of the estimated error of the state in the UNGM model for the algorithm of the present invention and the algorithm of document [23 ];
FIG. 4 is a diagram of the RMS error estimation of the algorithm of the present invention and the algorithm of document [23] with respect to states in a strongly non-linear model;
FIG. 5 is a diagram of the estimation error of the algorithm of the present invention and the algorithm of document [23] with respect to a state in a strongly non-linear model.
Detailed Description
The first embodiment is as follows: the Gaussian filtering method based on the non-linear network system under the non-ideal condition in the embodiment comprises the following specific processes:
firstly, establishing a system model and a sensor measurement model;
step two, providing hypothesis and lemma;
thirdly, designing a Gaussian filter based on the second step;
and step four, based on a third-order sphere diameter volume rule, approximating the Gaussian weighted integral in the step three to obtain a numerical form of the designed filter.
The second embodiment is as follows: the first embodiment is different from the first embodiment in that a system model and a sensor measurement model are established in the first step; the specific process is as follows:
establishing a nonlinear discrete time system model with correlated noise:
xk+1=f(xk)+ωk(7)
establishing a general nonlinear measurement model:
zk=h(xk)+υk(8)
in the formula, xk+1System state at time k +1, xkSystem state at time k, xk,xk+1∈Rn,RnIs an n-dimensional real number space; z is a radical ofkIs the sensor model at time k, zk∈Rm,RmIs m-dimensional real number space; f (-) and h (-) are known non-linear functions; omegak∈RnAnd upsilonk∈RmIs correlated zero mean white Gaussian noise and has a covariance of
Figure RE-GDA0002450881680000061
In the formula, deltaklIs the Kronecker delta function, QkAnd RkRespectively process noise and measurement noise covariance, SkIs the cross-covariance, l is the time l, ωl∈RnAnd upsilonl∈RmIs correlated zero mean gaussian white noise;
in the present invention, the following is considered: one step of random delay and packet loss may occur during data transmission; the data packet is sent only once; the estimator may receive a maximum of two measurement data simultaneously. That is, the estimator may receive zero, one or two packets. In practice, the measurement equation can be described by the following model. Considering communication bandwidth, delay measurement and data packet loss, a general nonlinear measurement model is further established as follows:
Figure RE-GDA0002450881680000062
in the formula, zkIs the sensor model at time k; z is a radical ofk-1Is the sensor model at time k; z is a radical ofk|k-1Is when z iskThe compensation amount when the compensation is lost is predicted for the measurement value at the k moment in one step; gamma raykand ηkIs an uncorrelated Bernoulli distribution variable and satisfies
Figure RE-GDA0002450881680000063
Figure RE-GDA0002450881680000064
P is the probability;
Figure RE-GDA0002450881680000065
is an intermediate variable; y iskIs a k-time sensor model with measurement skew and packet loss.
Other steps and parameters are the same as those in the first embodiment.
The third concrete implementation mode: the second embodiment is different from the first or second embodiment in that assumptions and reasoning are given in the second step; the specific process is as follows:
corresponding hypothesis and lemma are given, wherein the hypothesis is the premise of filter design, and the lemma is used for facilitating filter derivation;
hypothesis 1. hypothesis ωkkkand ηkAnd x0Is not related, and x0Satisfy the requirement of
Figure RE-GDA0002450881680000066
In the formula, x0Is the initial value of the number of the first,
Figure RE-GDA0002450881680000067
an estimated value of the initial value, E [ ]]To meet expectations (·)TIs composed of
Figure RE-GDA0002450881680000068
T is the transpose of the first image,
Figure RE-GDA0002450881680000069
the initial value corresponds to the covariance;
theory 1.A ═ aij]n×nIs a real valued matrix, B ═ diag { B }1,…,bnC and C ═ diag { C }1,…,cnIs a diagonal random matrix, defining
Figure RE-GDA0002450881680000071
In the formula, aijIs the ith row and jth column element of the A matrix, [ a ]ij]n×nIs the ith row and the jth column element is aijN × n matrix of (d), diag denotes arranging the following elements thereof as a diagonal matrix, b1Is the 1 st element of the diagonal of the matrix, bnIs the nth element of the diagonal of the matrix, c1Is the 1 st element of the diagonal of the matrix, cnFor the nth element of the diagonal of the matrix,E{BACTIs a first pair matrix BACTmaking a multiplication operation and then asking for an expectation, an Hadamard product defining [ M ⊙ N for arbitrary matrices M and N of the same dimension]ij=Mij·NijIt is clear that in the corresponding equation (12), even though A contains uncertainty, (12) still satisfies that A needs to be E [ A ] as long as A is uncorrelated with B, C]Replacing;
for the convenience of derivation hereinafter, the augmentation system is
Figure RE-GDA0002450881680000072
Then, there are
Figure RE-GDA0002450881680000073
Figure RE-GDA0002450881680000074
In the formula, Xk+1In order to augment the system with a wide range of systems,
Figure RE-GDA0002450881680000075
phi and
Figure RE-GDA0002450881680000076
are all intermediate variables;
wherein,
Figure RE-GDA0002450881680000077
and is
Figure RE-GDA0002450881680000078
Figure RE-GDA0002450881680000079
Wherein I is an identity matrix, 0 is a zero matrix,
Figure RE-GDA00024508816800000710
is the expectation of phi, E phi]For a period of phiThe physician can watch the disease,
Figure RE-GDA00024508816800000711
is composed of
Figure RE-GDA00024508816800000712
In the expectation that the position of the target is not changed,
Figure RE-GDA00024508816800000713
is a pair of
Figure RE-GDA00024508816800000714
Calculating expectation;
Figure RE-GDA00024508816800000715
a ⊥ b means a does not relate to b. (a) (.)TWhere represents the same amount as a. T represents the rank when it is the upper-label. Y isk=L{yk,yk-1,…,y1Where L { · } represents a linear space formed by a · leaf. 0 and I represent the zero matrix and identity matrix of the appropriate dimensions. E [. C]Where E represents the mathematical expectation. E [ a | b ]]Representing the condition expectation of a under the condition of b. P represents a covariance matrix. N (-) represents the distribution function.
Integral operation is represented by ^ integral.
Other steps and parameters are the same as those in the first or second embodiment.
The fourth concrete implementation mode: the difference between this embodiment and one of the first to third embodiments is that a gaussian filter is designed in the third step; the specific process is as follows:
the one-step prediction mean and covariance matrix are given as follows
Xk+1|k=Xk+1|k-1+Kkεk(1)
Figure RE-GDA0002450881680000081
In the formula, Xk+1|k-1In order to perform the two-step prediction,
Figure RE-GDA0002450881680000082
predicting the covariance matrix, ε, for two stepskIn order to be a new message,
Figure RE-GDA0002450881680000083
is an innovation covariance matrix, T is transposed, KkFor the gain matrix, the following is defined
Figure RE-GDA0002450881680000084
In the formula,
Figure RE-GDA0002450881680000085
is X under the measurement condition of time k-1k+1And epsilonkCross covariance matrix of (a);
the mean and covariance matrices after measurement modification are as follows:
Xk+1|k+1=Xk+1|k+Mk+1εk+1(4)
Figure RE-GDA0002450881680000086
in the formula,
Figure RE-GDA0002450881680000087
as an innovation covariance matrix, εk+1To be new, Mk+1Is a gain matrix;
Figure RE-GDA0002450881680000088
in the formula,
Figure RE-GDA0002450881680000089
is X under the measurement condition of time kk+1And epsilonk+1Cross covariance matrix of (a);
A. prediction
Theorem 1 for systems (13) - (14), the one-step prediction mean and covariance are as follows
Xk+1|k=Xk+1|k-1+Kkεk(17)
Figure RE-GDA00024508816800000810
Xk+1|k-1And
Figure RE-GDA00024508816800000811
the compounds are given in (20) and (24);
new message epsilonkSum covariance
Figure RE-GDA00024508816800000812
Has been given at the last moment;
and (3) proving that:
the certification process is divided into two parts:
1) the two-step prediction mean and covariance are calculated.
According to GASF, there are
Xk+1|k=Xk+1|k-1+Kkεk(19)
Because of omegak⊥L{yk-1,…,y1Thus, E [ omega ]k|Yk-1]Is equal to 0, so
Figure RE-GDA0002450881680000091
Bringing (20) into (19) to obtain (17);
then calculate
Figure RE-GDA0002450881680000092
From (19) to
Figure RE-GDA0002450881680000093
A covariance-based definition sum (21) of
Figure RE-GDA0002450881680000094
Rewriting Xk+1|k-1And with Xk+1Are subtracted to obtain
Figure RE-GDA0002450881680000095
According to the definition of covariance, there are
Figure RE-GDA0002450881680000096
Here, the
Figure RE-GDA0002450881680000097
Substituting (24) into (22) to obtain (18).
2) A gain matrix is calculated.
The gain matrix is defined as
Figure RE-GDA0002450881680000101
Here, ,
Figure RE-GDA0002450881680000102
then calculate
Figure RE-GDA0002450881680000103
And
Figure RE-GDA0002450881680000104
first, calculate
Figure RE-GDA0002450881680000105
Using ε in (46)k+1For εkIs provided with
Figure RE-GDA0002450881680000106
Here, the
Figure RE-GDA0002450881680000107
Figure RE-GDA0002450881680000111
Figure RE-GDA0002450881680000112
Figure RE-GDA0002450881680000113
In (31), N1As follows
Figure RE-GDA0002450881680000114
In N1In, remove
Figure RE-GDA0002450881680000115
All but a few of the amounts are known,
Figure RE-GDA0002450881680000116
the calculation is as follows,
Figure RE-GDA0002450881680000117
here, ωk-1|k-1And upsilonk-1|k-1The calculation is as follows.
yk-1And ωk-1The joint distribution of (A) can be as follows
Figure RE-GDA0002450881680000118
According to the conditional Gaussian distribution theorem, obtaining
Figure RE-GDA0002450881680000119
Figure RE-GDA00024508816800001110
For upsilonk-1Available as above
Figure RE-GDA0002450881680000121
Figure RE-GDA0002450881680000122
Of note is εk-1=yk-1-yk-1|k-2And
Figure RE-GDA0002450881680000123
has been calculated at time k-1;
then calculate
Figure RE-GDA0002450881680000124
Figure RE-GDA0002450881680000125
Here, ,
Figure RE-GDA0002450881680000126
Figure RE-GDA0002450881680000127
Figure RE-GDA0002450881680000128
Figure RE-GDA0002450881680000129
bringing (27) and (38) into (26) and (26) into (25), K can be obtainedk
The certification is over.
B. Correction
Theorem 2 for systems (13) - (14), the state corrections for the mean and covariance of the Gaussian filter are as follows:
Xk+1|k+1=Xk+1|k+Mk+1εk+1(43)
Figure RE-GDA00024508816800001210
wherein,
Figure RE-GDA00024508816800001211
Figure RE-GDA0002450881680000131
in the formula, yk+1For the k +1 moment sensor model with measurement lag and packet loss, N (-) is a distribution function, xk+1|kState x at time k +1k+1The one-step prediction value of (1),
Figure RE-GDA0002450881680000132
is a state xk+1One step of predicting covariance, xk|kIs state x at time kkIs determined by the estimated value of (c),
Figure RE-GDA0002450881680000133
is state x at time kkCovariance matrix of (v)k|kMeasuring noise upsilon for time kkAn estimated value of (d);
the certification process is divided into two parts:
1) calculating innovation epsilonk+1Sum covariance
Figure RE-GDA0002450881680000134
According to new message epsilonk+1Is defined as
εk+1=yk+1-yk+1|k(47)
Here, the
Figure RE-GDA0002450881680000135
In (48), xk+1|k,
Figure RE-GDA0002450881680000136
xk|kAnd
Figure RE-GDA0002450881680000137
bring (48) into (47), get (46).
Overwrite yk+1|kIs provided with
Figure RE-GDA0002450881680000138
Here upsilonk|kThe calculation process is as follows (36)k-1|k-1. Then, the user can use the device to perform the operation,
Figure RE-GDA0002450881680000141
is obviously provided with
Figure RE-GDA0002450881680000142
Figure RE-GDA0002450881680000143
Thus, it is possible to provide
Figure RE-GDA0002450881680000144
To calculate
Figure RE-GDA0002450881680000145
According to introduction 1, define
E[φAφT]=Φ⊙E[A](52)
Figure RE-GDA0002450881680000146
Figure RE-GDA0002450881680000147
Figure RE-GDA0002450881680000148
Order to
Figure RE-GDA0002450881680000149
Then the
Figure RE-GDA0002450881680000151
Figure RE-GDA0002450881680000152
Figure RE-GDA0002450881680000153
Figure RE-GDA0002450881680000154
Through a series of algebraic operations, obtain
Figure RE-GDA0002450881680000155
Hereinafter, the term on the right side is calculated (51).
Item1.
Figure RE-GDA0002450881680000156
Here, ,
Figure RE-GDA0002450881680000157
Figure RE-GDA0002450881680000158
Figure RE-GDA0002450881680000159
Figure RE-GDA00024508816800001510
Item2.
Figure RE-GDA00024508816800001511
here, ,
Figure RE-GDA00024508816800001512
and is
Figure RE-GDA00024508816800001513
Can obtain the product
Figure RE-GDA00024508816800001514
and
Figure RE-GDA00024508816800001515
Figure RE-GDA0002450881680000161
Where N is2And N1In the same form, but N1K in (1) needs to be replaced by k +1,
Item3.
Figure RE-GDA0002450881680000162
Item4.
Figure RE-GDA0002450881680000163
Item5.
Figure RE-GDA0002450881680000164
Item6.
Figure RE-GDA0002450881680000165
Item7.
Figure RE-GDA0002450881680000166
Item8.
Figure RE-GDA0002450881680000167
Item9.
Figure RE-GDA0002450881680000168
Item10.
because of upsilonk+1And upsilonkIs not related, available
Figure RE-GDA0002450881680000171
Item11.
Figure RE-GDA0002450881680000172
Item12.
Figure RE-GDA0002450881680000173
Item13.
Figure RE-GDA0002450881680000174
Bring (61), (66), (68) - (78) into (5)1) Is obtained by
Figure RE-GDA0002450881680000175
2) Computing a gain matrix Mk+1.
From (45), only calculation is needed
Figure RE-GDA0002450881680000176
Figure RE-GDA00024508816800001711
Here, ,
Figure RE-GDA0002450881680000178
can be similar to
Figure RE-GDA0002450881680000179
And (6) obtaining.
Figure RE-GDA00024508816800001710
Here, ,
Figure RE-GDA0002450881680000181
Figure RE-GDA0002450881680000182
Figure RE-GDA0002450881680000183
Figure RE-GDA0002450881680000184
bringing (82) - (84) into (81) can obtain
Figure RE-GDA0002450881680000185
Will be provided with
Figure RE-GDA0002450881680000186
And
Figure RE-GDA0002450881680000187
bring in (79), can obtain
Figure RE-GDA0002450881680000188
The certification is complete.
Other steps and parameters are the same as those in one of the first to third embodiments.
The fifth concrete implementation mode: the difference between this embodiment and the first to the fourth embodiment is that the gaussian weighted integral in step four is approximated based on the third-order sphere diameter volume rule to obtain the numerical form of the designed filter; the specific process is as follows:
based on the sphere diameter volume rule algorithm, the numerical implementation of the algorithm provided by the invention is as follows.
And (3) prediction:
1. decomposition of
Figure RE-GDA0002450881680000189
Figure RE-GDA00024508816800001810
In the formula,
Figure RE-GDA00024508816800001811
one-step predictive covariance matrix, M, for the state at time kk|k-1Is the intermediate variable(s) of the variable,
Figure RE-GDA00024508816800001812
augmenting system for time k
Figure RE-GDA00024508816800001813
The one-step prediction of the covariance matrix,
Figure RE-GDA00024508816800001814
is an intermediate variable;
in N1In the middle, let
Figure RE-GDA00024508816800001815
Figure RE-GDA00024508816800001816
Figure RE-GDA00024508816800001817
In the formula, N1Is a Gaussian distribution, Γk,k-1|k-1Being intermediate variables, is the state x of the constructkAnd noise vk-1Estimation based on measurements at time k-1, xk|k-1One-step prediction of state at time k, vk-1|k-1For measuring the estimated value of the noise at the time k-1, Πk,k-1|k-1Is the intermediate variable(s) of the variable,
Figure RE-GDA0002450881680000191
the covariance matrix is predicted for the state at time k in one step,
Figure RE-GDA0002450881680000192
is a state xkAnd noise vk-1Based on the cross covariance matrix measured at time k-1,
Figure RE-GDA0002450881680000193
covariance matrix, M, for the measured noise at time k-1k,k-1|k-1Is an intermediate variable;
2. calculating volume points
xi,k|k-1=Mk|k-1ξi+xk|k-1,i=1,…2n (88)
Figure RE-GDA0002450881680000194
Figure RE-GDA0002450881680000195
In the formula, xi,k|k-1volumetric point, ξ, for time k with respect to a one-step predictori、ζi
Figure RE-GDA00024508816800001913
Sigma points and dimensions of 2n, 4n and 2(n + m), respectively;
Figure RE-GDA0002450881680000197
for component x in joint estimationk|k-1The volume point of (a) is,
Figure RE-GDA0002450881680000198
for component x in joint estimationk-1|k-1Volume point of (1), Xi,k|k-1Augmenting system federation states for time k
Figure RE-GDA0002450881680000199
One step of predicting the volume point, δi,k,k-1|k-1To correspond to gammak,k-1|k-1In xk|k-1The ith volume point of the part, gammai,k,k-1|k-1To correspond to gammak,k-1|k-1Is on vk-1|k-1Ith volume point of the part, Mk,k-1|k-1Being an intermediate variable, Γk,k-1|k-1Is an intermediate variable, n is a system state dimension, and m is a measurement noise dimension;
3. calculating propagated volume points
Figure RE-GDA00024508816800001910
Figure RE-GDA00024508816800001911
Figure RE-GDA00024508816800001912
In the formula,
Figure RE-GDA0002450881680000201
is xi,k|k-1The volume points after propagation through the system model,
Figure RE-GDA0002450881680000202
is xi,k|k-1The volume point, h (-) after propagation through the metrology model is a known non-linear function,
Figure RE-GDA0002450881680000203
is composed of
Figure RE-GDA0002450881680000204
The volume points, f (-) after propagation through the system model are known non-linear functions,
Figure RE-GDA0002450881680000205
is composed of
Figure RE-GDA0002450881680000206
The volume points after propagation through the metrology model,
Figure RE-GDA0002450881680000207
is deltai,k,k-1|k-1Volume points after propagation through the system model;
definition of
Figure RE-GDA0002450881680000208
Figure RE-GDA0002450881680000209
Figure RE-GDA00024508816800002010
Figure RE-GDA00024508816800002011
Figure RE-GDA00024508816800002012
Figure RE-GDA00024508816800002013
Figure RE-GDA00024508816800002014
In the formula I1、l2、l3、l4、l5、l6、l7、l8、l9、l10、l11、l12、l13、l14、l15Is an intermediate variable;
4. calculating the corresponding quantities in equations 1 and 2
Figure RE-GDA00024508816800002015
Figure RE-GDA00024508816800002016
Wherein,
Figure RE-GDA00024508816800002017
Figure RE-GDA00024508816800002018
Figure RE-GDA0002450881680000211
where ω isk-1|k-1And upsilonk-1|k-1Given in (34) and (36);
bring in (103) - (104)
Figure RE-GDA0002450881680000212
To obtain
Figure RE-GDA0002450881680000213
For epsilon obtained at the last momentkAnd
Figure RE-GDA0002450881680000214
carry in (1) - (2) to obtain X directlyk+1|kAnd
Figure RE-GDA0002450881680000215
in the formula, xk|k-1The mean is predicted for one step of the state at time k,
Figure RE-GDA0002450881680000216
the covariance is predicted for one step of the state at time k,
Figure RE-GDA0002450881680000217
is a state x under the measurement condition at the time k-1k+1The cross-covariance matrix with the innovation at time k,
Figure RE-GDA0002450881680000218
is a state x under the measurement condition at the time k-1kCross covariance matrix with innovation at time k, QkProcess noise at time k;
and (3) correction:
1. decomposition of
Figure RE-GDA0002450881680000219
Figure RE-GDA00024508816800002110
Figure RE-GDA00024508816800002111
In the formula,
Figure RE-GDA00024508816800002112
one-step prediction for state at time k +1Covariance matrix, Mk+1|kIs the intermediate variable(s) of the variable,
Figure RE-GDA00024508816800002113
estimating a covariance matrix, M, for a state at time kk|kIs the intermediate variable(s) of the variable,
Figure RE-GDA00024508816800002114
augmenting system joint states for k +1 moments
Figure RE-GDA00024508816800002115
The one-step prediction of the covariance matrix,
Figure RE-GDA00024508816800002116
is an intermediate variable;
in N2In the middle, let
Figure RE-GDA00024508816800002117
Figure RE-GDA00024508816800002118
In the formula, N2Is a Gaussian distribution, Γk+1,k|kIs an intermediate variable, Πk+1,k|kIs gammak+1,k|kThe corresponding covariance matrix is then used as a basis,
Figure RE-GDA00024508816800002119
the cross covariance matrix is predicted for one step at time k +1,
Figure RE-GDA00024508816800002120
is a state x under the measurement condition of time kk+1And measure noise vkThe cross-covariance matrix of (a) is,
Figure RE-GDA00024508816800002121
measuring a posterior covariance matrix of the noise at the time k;
2. calculating volume points
xi,k+1|k=Mk+1|kξi+xk+1|k,i=1,…2n (109)
xi,k|k=Mk|kξi+xk|k,i=1,…2n (110)
Figure RE-GDA0002450881680000221
Figure RE-GDA0002450881680000222
In the formula, xi,k+1|kOne-step prediction of the volume point, x, for the state quantity at the time k +1k+1|kFor one-step prediction of state quantities at the time k +1, xi,k|kCorresponding product point, x, for the estimated value of state quantity at time kk|kIs an estimate of the state quantity at time k, Xi,k+1|kAugmenting system joint state quantities for k +1 moments
Figure RE-GDA0002450881680000223
The volume point is predicted in one step,
Figure RE-GDA0002450881680000224
augmenting system joint state quantities for k +1 moments
Figure RE-GDA0002450881680000225
Corresponding to xk+1|kThe volume point of the part of the volume,
Figure RE-GDA0002450881680000226
augmenting system joint state quantities for k +1 moments
Figure RE-GDA0002450881680000227
Corresponding to xk|kThe volume point of the part of the volume,
Figure RE-GDA0002450881680000228
is an intermediate variable, Xk+1|kIs a joint state quantity at the time of k +1
Figure RE-GDA0002450881680000229
One step prediction, δi,k+1,k|kIs a combined quantity
Figure RE-GDA00024508816800002210
Corresponding to xk+1|kPartial volume point, gammai,k+1,k|kIs a combined quantity
Figure RE-GDA00024508816800002211
Corresponding upsilonk|kPartial volume points, Mk+1,k|kBeing an intermediate variable, Γk+1,k|kIs an intermediate variable;
3. calculating propagated volume points
Figure RE-GDA00024508816800002212
Figure RE-GDA00024508816800002213
Figure RE-GDA00024508816800002214
Figure RE-GDA00024508816800002215
In the formula,
Figure RE-GDA00024508816800002216
is xi,k+1|kThe volume points propagated through the system model,
Figure RE-GDA00024508816800002217
is xi,k|kThe volume points propagated through the system model,
Figure RE-GDA00024508816800002218
is xi,k+1|kThe volume points propagated through the metrology model,
Figure RE-GDA00024508816800002219
is xi,k|kThe volume points propagated through the metrology model,
Figure RE-GDA00024508816800002220
is composed of
Figure RE-GDA00024508816800002221
The volume points propagated through the metrology model,
Figure RE-GDA00024508816800002222
is composed of
Figure RE-GDA00024508816800002223
The volume points propagated through the metrology model,
Figure RE-GDA00024508816800002224
is deltai,k+1,k|kVolume points propagated through the metrology model;
definition of
Figure RE-GDA0002450881680000231
Figure RE-GDA0002450881680000232
Figure RE-GDA0002450881680000233
Figure RE-GDA0002450881680000234
Figure RE-GDA0002450881680000235
Figure RE-GDA0002450881680000236
Figure RE-GDA0002450881680000237
Figure RE-GDA0002450881680000238
Of formula (II) to'1、l′2、l′3、l′4、l′5、l′6、l′7、l′8、l′9、l′10、l′11、l′12、l′13、l′14、l′15、l′16、l′17Is an intermediate variable;
4. calculating the corresponding quantities in equations 4, 5 and 6
Figure RE-GDA0002450881680000239
In the formula, yk+1A k +1 moment sensor model with measurement time lag and data packet loss exists;
Figure RE-GDA0002450881680000241
Figure RE-GDA0002450881680000242
Figure RE-GDA0002450881680000243
order to
Figure RE-GDA0002450881680000244
Then there is
Figure RE-GDA0002450881680000245
Figure RE-GDA0002450881680000246
Figure RE-GDA0002450881680000247
Figure RE-GDA0002450881680000248
In the formula,
Figure RE-GDA0002450881680000249
an innovation covariance matrix at the time k +1, phi is an intermediate variable, Hadamard product, psi is an intermediate variable, and xi is an intermediate variable,
Figure RE-GDA00024508816800002410
is a state x measured at time kkThe cross-covariance matrix with the innovation at time k +1,
Figure RE-GDA00024508816800002411
a cross-covariance matrix of the state at time k +1 and the innovation at time k +1 measured at time k, SkFor the cross-correlation noise at time k, omegak|kMeasuring process noise omega for time kkAn estimated value of (d);
bringing (127) - (128) into
Figure RE-GDA0002450881680000251
To obtain
Figure RE-GDA0002450881680000252
From (125) - (126), epsilon is obtainedk+1And
Figure RE-GDA0002450881680000253
will epsilonk+1And
Figure RE-GDA0002450881680000254
carry in (43) - (44) to obtain Xk+1|k+1And
Figure RE-GDA0002450881680000255
Xk+1|k+1=Xk+1|k+Mk+1εk+1(43)
Figure RE-GDA0002450881680000256
other steps and parameters are the same as in one of the first to fourth embodiments.
The first embodiment is as follows:
numerical simulation analysis
To verify the validity of the control algorithm designed by the present invention, two non-linear models were used to test the validity of the proposed algorithm. Wherein, 'original', 'thisapaper' and '23' represent reference signals, the estimation results of the algorithm proposed herein and the estimation results of the algorithm in [23] are extended to a nonlinear system based on EKF.
1) Single variable non-static growth model (UNGM)
UNGM is represented as follows:
Figure RE-GDA0002450881680000257
Figure RE-GDA0002450881680000258
where ω iskAnd upsilonkIs zero mean white Gaussian noise and the covariance is satisfied
Figure RE-GDA0002450881680000259
Initial value of x0Filter initial value is x ═ 0.30|00 and P0|0Probability satisfied by random variables of Bernoulli distribution
Figure RE-GDA00024508816800002510
And
Figure RE-GDA00024508816800002511
300 independent shadesA tecalol simulation is performed. And the Error at time k (Error) and the Root Mean Square Error (RMSE) are defined as follows
Figure RE-GDA00024508816800002512
Figure RE-GDA00024508816800002513
Here, the
Figure RE-GDA00024508816800002514
And
Figure RE-GDA00024508816800002515
the true and estimated values of nth Monte Carlo guidelines at time k are shown, N is the total number of simulations, FIGS. 1-3 are the corresponding simulation results, where FIG. 1 illustrates that the estimation results herein are closer to the true values and have better results than the literature [23]]Fig. 2 and 3 show the algorithm and document [23]]The advantages of the algorithm are further illustrated by the estimation error and the estimation root mean square error of the algorithm, the algorithm fluctuation is small, and the algorithm robustness is good.
2) Strongly non-linear model
The strong nonlinear model is given as follows:
Figure RE-GDA0002450881680000261
zk=cos(x1,k)+x2,kx3,kk(135)
here, xkAnd zkIs the system state and quantity measurement, ωkAnd upsilonkIs zero mean white Gaussian noise and the covariance is satisfied
Figure RE-GDA0002450881680000262
Initial state is x0=[-0.7 1 1]T. The initial value of the filter is x0|0=[0 0 0]TAnd P0|0=I3
Figure RE-GDA0002450881680000263
And
Figure RE-GDA0002450881680000264
as described above.
Fig. 4-5 are corresponding simulation results, where fig. 4 and 5 are respectively an estimation error and an estimated root mean square error of the algorithm in this document and the algorithm in document [23], and the simulation results show that the proposed algorithm still has a strong processing capability even for a strong nonlinear system, while the algorithm in document [23] hardly obtains a satisfactory true value, which illustrates that the algorithm has a stronger capability of processing problems when applied to an actual system.
The present invention is capable of other embodiments and its several details are capable of modifications in various obvious respects, all without departing from the spirit and scope of the present invention.

Claims (5)

1. The Gaussian filtering method based on the nonlinear network system under the nonideal condition is characterized by comprising the following steps: the method comprises the following specific processes:
firstly, establishing a system model and a sensor measurement model;
step two, providing hypothesis and lemma;
thirdly, designing a Gaussian filter based on the second step;
and step four, based on a third-order sphere diameter volume rule, approximating the Gaussian weighted integral in the step three to obtain a numerical form of the designed filter.
2. The Gaussian filtering method based on the non-linear network system under the non-ideal condition as claimed in claim 1 is characterized in that: establishing a system model and a sensor measurement model in the first step; the specific process is as follows:
establishing a nonlinear discrete time system model with correlated noise:
xk+1=f(xk)+ωk(7)
establishing a general nonlinear measurement model:
zk=h(xk)+υk(8)
in the formula, xk+1System state at time k +1, xkSystem state at time k, xk,xk+1∈Rn,RnIs an n-dimensional real number space; z is a radical ofkIs the sensor model at time k, zk∈Rm,RmIs m-dimensional real number space; f (-) and h (-) are known non-linear functions; omegak∈RnAnd upsilonk∈RmIs correlated zero mean white Gaussian noise and has a covariance of
Figure RE-FDA0002450881670000011
In the formula, deltaklIs the Kronecker delta function, QkAnd RkRespectively process noise and measurement noise covariance, SkIs the cross-covariance, l is the time l, ωl∈RnAnd upsilonl∈RmIs correlated zero mean gaussian white noise;
considering communication bandwidth, delay measurement and data packet loss, a general nonlinear measurement model is further established as follows:
Figure RE-FDA0002450881670000012
in the formula, zkIs the sensor model at time k; z is a radical ofk-1Is the sensor model at time k; z is a radical ofk|k-1Is when z iskThe compensation amount when the compensation is lost is predicted for the measurement value at the k moment in one step; gamma raykand ηkIs an uncorrelated Bernoulli distribution variable and satisfies
Figure RE-FDA0002450881670000013
Figure RE-FDA0002450881670000014
P is the probability;
Figure RE-FDA0002450881670000016
Figure RE-FDA0002450881670000015
is an intermediate variable; y iskIs a k-time sensor model with measurement skew and packet loss.
3. The Gaussian filtering method based on the nonlinear network system under the non-ideal condition as claimed in claim 1 or 2, characterized in that: giving assumptions and lemmas in the second step; the specific process is as follows:
hypothesis 1. hypothesis ωkkkand ηkAnd x0Is not related, and x0Satisfy the requirement of
Figure RE-FDA0002450881670000021
In the formula, x0Is the initial value of the number of the first,
Figure RE-FDA0002450881670000022
an estimated value of the initial value, E [ ]]To meet expectations (·)TIs composed of
Figure RE-FDA0002450881670000023
T is the transpose of the first image,
Figure RE-FDA0002450881670000024
the initial value corresponds to the covariance;
theory 1.A ═ aij]n×nIs a real valued matrix, B ═ diag { B }1,…,bnC and C ═ diag { C }1,…,cnIs a diagonal random matrix, defining
Figure RE-FDA0002450881670000025
In the formula, aijIs the ith row and jth column element of the A matrix, [ a ]ij]n×nIs the ith row and the jth column element is aijN × n matrix of (d), diag denotes arranging the following elements thereof as a diagonal matrix, b1Is the 1 st element of the diagonal of the matrix, bnIs the nth element of the diagonal of the matrix, c1Is the 1 st element of the diagonal of the matrix, cnFor the nth element of the diagonal of the matrix, E { BACTIs a first pair matrix BACTa multiplication operation is performed and then an expectation is asserted, an Hadamard product;
the augmentation system is
Figure RE-FDA0002450881670000026
Then, there are
Figure RE-FDA0002450881670000027
Figure RE-FDA0002450881670000028
In the formula, Xk+1In order to augment the system with a wide range of systems,
Figure RE-FDA0002450881670000029
phi and
Figure RE-FDA00024508816700000210
are all intermediate variables;
wherein,
Figure RE-FDA00024508816700000211
and is
Figure RE-FDA0002450881670000031
Figure RE-FDA0002450881670000032
Wherein I is an identity matrix, 0 is a zero matrix,
Figure RE-FDA0002450881670000033
is the expectation of phi, E phi]In order to make the expectation for phi,
Figure RE-FDA0002450881670000034
is composed of
Figure RE-FDA0002450881670000035
In the expectation that the position of the target is not changed,
Figure RE-FDA0002450881670000036
is a pair of
Figure RE-FDA0002450881670000037
And (4) making expectations.
4. The Gaussian filtering method based on the nonlinear network system under the non-ideal condition as claimed in claim 3 is characterized in that: designing a Gaussian filter in the third step; the specific process is as follows:
the one-step prediction mean and covariance matrix are given as follows
Xk+1|k=Xk+1|k-1+Kkεk(1)
Figure RE-FDA0002450881670000038
In the formula, Xk+1|k-1In order to perform the two-step prediction,
Figure RE-FDA0002450881670000039
predicting the covariance matrix, ε, for two stepskIn order to be a new message,
Figure RE-FDA00024508816700000310
is an innovation covariance matrix, T is transposed, KkFor the gain matrix, the following is defined
Figure RE-FDA00024508816700000311
In the formula,
Figure RE-FDA00024508816700000312
is X under the measurement condition of time k-1k+1And epsilonkCross covariance matrix of (a);
the mean and covariance matrices after measurement modification are as follows:
Xk+1|k+1=Xk+1|k+Mk+1εk+1(4)
Figure RE-FDA00024508816700000313
in the formula,
Figure RE-FDA00024508816700000314
as an innovation covariance matrix, εk+1To be new, Mk+1Is a gain matrix;
Figure RE-FDA00024508816700000315
in the formula,
Figure RE-FDA00024508816700000316
is X under the measurement condition of time kk+1And epsilonk+1Cross covariance matrix of (2).
5. The Gaussian filtering method based on the nonlinear network system under the non-ideal condition as recited in claim 4, characterized in that: based on the third-order sphere diameter volume rule in the fourth step, the Gaussian weighted integral in the third step is approximated to obtain the numerical form of the designed filter; the specific process is as follows:
and (3) prediction:
1. decomposition of
Figure RE-FDA0002450881670000041
Figure RE-FDA0002450881670000042
In the formula,
Figure RE-FDA0002450881670000043
one-step predictive covariance matrix, M, for the state at time kk|k-1Is the intermediate variable(s) of the variable,
Figure RE-FDA0002450881670000044
augmenting system for time k
Figure RE-FDA0002450881670000045
The one-step prediction of the covariance matrix,
Figure RE-FDA0002450881670000046
is an intermediate variable;
in N1In the middle, let
Figure RE-FDA0002450881670000047
Figure RE-FDA0002450881670000048
Figure RE-FDA0002450881670000049
In the formula, N1Is a Gaussian distribution, Γk,k-1|k-1Being intermediate variables, is the state x of the constructkAnd noise vk-1Estimation based on measurements at time k-1, xk|k-1One-step prediction of state at time k, vk-1|k-1For measuring the estimated value of the noise at the time k-1, Πk,k-1|k-1Is the intermediate variable(s) of the variable,
Figure RE-FDA00024508816700000410
the covariance matrix is predicted for the state at time k in one step,
Figure RE-FDA00024508816700000411
is a state xkAnd noise vk-1Based on the cross covariance matrix measured at time k-1,
Figure RE-FDA00024508816700000412
covariance matrix, M, for the measured noise at time k-1k,k-1|k-1Is an intermediate variable;
2. calculating volume points
xi,k|k-1=Mk|k-1ξi+xk|k-1,i=1,…2n (88)
Figure RE-FDA00024508816700000413
Figure RE-FDA00024508816700000414
In the formula, xi,k|k-1volumetric point, ξ, for time k with respect to a one-step predictori、ζi
Figure RE-FDA00024508816700000415
Sigma points and dimensions of 2n, 4n and 2(n + m), respectively;
Figure RE-FDA0002450881670000051
for component x in joint estimationk|k-1The volume point of (a) is,
Figure RE-FDA0002450881670000052
for component x in joint estimationk-1|k-1Volume point of (1), Xi,k|k-1Augmenting system federation states for time k
Figure RE-FDA0002450881670000053
One step of predicting the volume point, δi,k,k-1|k-1To correspond to gammak,k-1|k-1In xkk-1The ith volume point of the part, gammai,k,k-1|k-1To correspond to gammak,k-1|k-1Is on vk-1|k-1Ith volume point of the part, Mk,k-1|k-1Being an intermediate variable, Γk,k-1|k-1Is an intermediate variable, n is a system state dimension, and m is a measurement noise dimension;
3. calculating propagated volume points
Figure RE-FDA0002450881670000054
Figure RE-FDA0002450881670000055
Figure RE-FDA0002450881670000056
In the formula,
Figure RE-FDA0002450881670000057
is xi,k|k-1The volume points after propagation through the system model,
Figure RE-FDA0002450881670000058
is xi,k|k-1The volume point, h (-) after propagation through the metrology model is a known non-linear function,
Figure RE-FDA0002450881670000059
is composed of
Figure RE-FDA00024508816700000510
The volume point, f (-) after propagation through the system model is a known non-linear functionThe number of the first and second groups is,
Figure RE-FDA00024508816700000511
is composed of
Figure RE-FDA00024508816700000512
The volume points after propagation through the metrology model,
Figure RE-FDA00024508816700000513
is deltai,k,k-1|k-1Volume points after propagation through the system model;
definition of
Figure RE-FDA00024508816700000514
Figure RE-FDA00024508816700000515
Figure RE-FDA00024508816700000516
Figure RE-FDA00024508816700000517
Figure RE-FDA00024508816700000518
Figure RE-FDA0002450881670000061
Figure RE-FDA0002450881670000062
In the formula I1、l2、l3、l4、l5、l6、l7、l8、l9、l10、l11、l12、l13、l14、l15Is an intermediate variable;
4. calculating the corresponding quantities in equations 1 and 2
Figure RE-FDA0002450881670000063
Figure RE-FDA0002450881670000064
Wherein,
Figure RE-FDA0002450881670000065
Figure RE-FDA0002450881670000066
Figure RE-FDA0002450881670000067
bring in (103) - (104)
Figure RE-FDA0002450881670000068
To obtain
Figure RE-FDA0002450881670000069
For epsilon obtained at the last momentkAnd
Figure RE-FDA00024508816700000610
carry in (1) - (2) to obtain X directlyk+1|kAnd
Figure RE-FDA00024508816700000611
in the formula, xk|k-1The mean is predicted for one step of the state at time k,
Figure RE-FDA00024508816700000612
the covariance is predicted for one step of the state at time k,
Figure RE-FDA00024508816700000613
is a state x under the measurement condition at the time k-1k+1The cross-covariance matrix with the innovation at time k,
Figure RE-FDA00024508816700000614
is a state x under the measurement condition at the time k-1kCross covariance matrix with innovation at time k, QkProcess noise at time k;
and (3) correction:
1. decomposition of
Figure RE-FDA00024508816700000615
Figure RE-FDA00024508816700000616
Figure RE-FDA00024508816700000617
In the formula,
Figure RE-FDA0002450881670000071
one-step prediction of covariance matrix, M, for state at time k +1k+1|kIs the intermediate variable(s) of the variable,
Figure RE-FDA0002450881670000072
estimating a covariance matrix, M, for a state at time kk|kIs the intermediate variable(s) of the variable,
Figure RE-FDA0002450881670000073
augmenting system joint states for k +1 moments
Figure RE-FDA0002450881670000074
The one-step prediction of the covariance matrix,
Figure RE-FDA0002450881670000075
is an intermediate variable;
in N2In the middle, let
Figure RE-FDA0002450881670000076
Figure RE-FDA0002450881670000077
In the formula, N2Is a Gaussian distribution, Γk+1,k|kIs an intermediate variable, Πk+1,k|kIs gammak+1,k|kThe corresponding covariance matrix is then used as a basis,
Figure RE-FDA0002450881670000078
the cross covariance matrix is predicted for one step at time k +1,
Figure RE-FDA0002450881670000079
is a state x under the measurement condition of time kk+1And measure noise vkThe cross-covariance matrix of (a) is,
Figure RE-FDA00024508816700000710
measuring a posterior covariance matrix of the noise at the time k;
2. calculating volume points
xi,k+1|k=Mk+1|kξi+xk+1|k,i=1,…2n (109)
xi,k|k=Mk|kξi+xk|k,i=1,…2n (110)
Figure RE-FDA00024508816700000711
Figure RE-FDA00024508816700000712
In the formula, xi,k+1|kOne-step prediction of the volume point, x, for the state quantity at the time k +1k+1|kFor one-step prediction of state quantities at the time k +1, xi,k|kCorresponding product point, x, for the estimated value of state quantity at time kk|kIs an estimate of the state quantity at time k, Xi,k+1|kAugmenting system joint state quantities for k +1 moments
Figure RE-FDA00024508816700000713
The volume point is predicted in one step,
Figure RE-FDA00024508816700000714
augmenting system joint state quantities for k +1 moments
Figure RE-FDA00024508816700000715
Corresponding to xk+1|kThe volume point of the part of the volume,
Figure RE-FDA00024508816700000716
augmenting system joint state quantities for k +1 moments
Figure RE-FDA00024508816700000717
Corresponding to xk|kThe volume point of the part of the volume,
Figure RE-FDA00024508816700000718
is an intermediate variable, Xk+1|kIs a joint state quantity at the time of k +1
Figure RE-FDA00024508816700000719
One step prediction, δi,k+1,k|kIs a combined quantity
Figure RE-FDA00024508816700000720
Corresponding to xk+1|kPartial volume point, gammai,k+1,k|kIs a combined quantity
Figure RE-FDA0002450881670000081
Corresponding upsilonk|kPartial volume points, Mk+1,k|kBeing an intermediate variable, Γk+1,k|kIs an intermediate variable;
3. calculating propagated volume points
Figure RE-FDA0002450881670000082
Figure RE-FDA0002450881670000083
Figure RE-FDA0002450881670000084
Figure RE-FDA0002450881670000085
In the formula,
Figure RE-FDA0002450881670000086
is xi,k+1|kThe volume points propagated through the system model,
Figure RE-FDA0002450881670000087
is xi,k|kThe volume points propagated through the system model,
Figure RE-FDA0002450881670000088
is xi,k+1|kThe volume points propagated through the metrology model,
Figure RE-FDA0002450881670000089
is xi,k|kThe volume points propagated through the metrology model,
Figure RE-FDA00024508816700000810
is composed of
Figure RE-FDA00024508816700000811
The volume points propagated through the metrology model,
Figure RE-FDA00024508816700000812
is composed of
Figure RE-FDA00024508816700000813
The volume points propagated through the metrology model,
Figure RE-FDA00024508816700000814
is deltai,k+1,k|kVolume points propagated through the metrology model;
definition of
Figure RE-FDA00024508816700000815
Figure RE-FDA00024508816700000816
Figure RE-FDA00024508816700000817
Figure RE-FDA00024508816700000818
Figure RE-FDA00024508816700000819
Figure RE-FDA00024508816700000820
Figure RE-FDA00024508816700000821
Figure RE-FDA0002450881670000091
Of formula (II) to'1、l′2、l′3、l′4、l′5、l′6、l′7、l′8、l′9、l′10、l′11、l′12、l′13、l′14、l′15、l′16、l′17Is an intermediate variable;
4. calculating the corresponding quantities in equations 4, 5 and 6
Figure RE-FDA0002450881670000092
In the formula, yk+1A k +1 moment sensor model with measurement time lag and data packet loss exists;
Figure RE-FDA0002450881670000093
Figure RE-FDA0002450881670000094
Figure RE-FDA0002450881670000095
in the formula,
Figure RE-FDA0002450881670000096
an innovation covariance matrix at the time k +1, phi is an intermediate variable, Hadamard product, psi is an intermediate variable, and xi is an intermediate variable,
Figure RE-FDA0002450881670000097
is a state x measured at time kkThe cross-covariance matrix with the innovation at time k +1,
Figure RE-FDA0002450881670000098
a cross-covariance matrix of the state at time k +1 and the innovation at time k +1 measured at time k, SkFor the cross-correlation noise at time k, omegak|kMeasuring process noise omega for time kkAn estimated value of (d);
bringing (127) - (128) into
Figure RE-FDA0002450881670000099
To obtain
Figure RE-FDA00024508816700000910
From (125) - (126), epsilon is obtainedk+1And
Figure RE-FDA00024508816700000911
will epsilonk+1And
Figure RE-FDA00024508816700000912
carry in (43) - (44) to obtain Xk+1k+1And
Figure RE-FDA00024508816700000913
Xk+1|k+1=Xk+1|k+Mk+1εk+1(43)
Figure RE-FDA0002450881670000101
CN201911400580.2A 2019-12-30 2019-12-30 Gaussian filtering method based on non-linear network system under non-ideal condition Active CN111193528B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201911400580.2A CN111193528B (en) 2019-12-30 2019-12-30 Gaussian filtering method based on non-linear network system under non-ideal condition

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201911400580.2A CN111193528B (en) 2019-12-30 2019-12-30 Gaussian filtering method based on non-linear network system under non-ideal condition

Publications (2)

Publication Number Publication Date
CN111193528A true CN111193528A (en) 2020-05-22
CN111193528B CN111193528B (en) 2021-09-17

Family

ID=70709614

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201911400580.2A Active CN111193528B (en) 2019-12-30 2019-12-30 Gaussian filtering method based on non-linear network system under non-ideal condition

Country Status (1)

Country Link
CN (1) CN111193528B (en)

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112034713A (en) * 2020-09-07 2020-12-04 山东大学 Optimal state estimation method and system for moving target in non-ideal network environment
CN112653410A (en) * 2020-12-15 2021-04-13 中国地质大学(武汉) Design method of Gaussian filter of noise correlation and measurement one-step random delay system
CN113011475A (en) * 2021-01-29 2021-06-22 深圳信息职业技术学院 Distributed fusion algorithm considering correlated noise and random parameter matrix
CN113985737A (en) * 2021-10-27 2022-01-28 湘潭大学 Research on networked control system with time delay and packet loss

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5995539A (en) * 1993-03-17 1999-11-30 Miller; William J. Method and apparatus for signal transmission and reception
CN106372646A (en) * 2016-08-30 2017-02-01 上海交通大学 Multi-target tracking method based on SRCK-GMCPHD filtering
CN106871905A (en) * 2017-03-02 2017-06-20 哈尔滨工业大学 Gaussian filtering substitutes framework composition air navigation aid under a kind of non-ideal condition
CN107045490A (en) * 2017-05-09 2017-08-15 衢州学院 A kind of method for estimating state of nonlinear system
CN109345875A (en) * 2018-09-28 2019-02-15 哈尔滨工程大学 An Estimation Method to Improve the Measurement Accuracy of Ship Automatic Identification System

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5995539A (en) * 1993-03-17 1999-11-30 Miller; William J. Method and apparatus for signal transmission and reception
CN106372646A (en) * 2016-08-30 2017-02-01 上海交通大学 Multi-target tracking method based on SRCK-GMCPHD filtering
CN106871905A (en) * 2017-03-02 2017-06-20 哈尔滨工业大学 Gaussian filtering substitutes framework composition air navigation aid under a kind of non-ideal condition
CN107045490A (en) * 2017-05-09 2017-08-15 衢州学院 A kind of method for estimating state of nonlinear system
CN109345875A (en) * 2018-09-28 2019-02-15 哈尔滨工程大学 An Estimation Method to Improve the Measurement Accuracy of Ship Automatic Identification System

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
于晗: "一类非理想条件下非线性系统的高斯滤波算法及其应用研究", 《中国博士学位论文全文数据库·工程科技Ⅱ辑》 *

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112034713A (en) * 2020-09-07 2020-12-04 山东大学 Optimal state estimation method and system for moving target in non-ideal network environment
CN112034713B (en) * 2020-09-07 2021-10-19 山东大学 Optimal state estimation method and system for moving target in non-ideal network environment
CN112653410A (en) * 2020-12-15 2021-04-13 中国地质大学(武汉) Design method of Gaussian filter of noise correlation and measurement one-step random delay system
CN113011475A (en) * 2021-01-29 2021-06-22 深圳信息职业技术学院 Distributed fusion algorithm considering correlated noise and random parameter matrix
CN113985737A (en) * 2021-10-27 2022-01-28 湘潭大学 Research on networked control system with time delay and packet loss

Also Published As

Publication number Publication date
CN111193528B (en) 2021-09-17

Similar Documents

Publication Publication Date Title
CN111193528A (en) Gaussian filtering method based on nonlinear network system under non-ideal conditions
CN103141039B (en) Method for time synchronization in a communication network
CN108364014A (en) A kind of multi-sources Information Fusion Method based on factor graph
CN110376620B (en) Real-time clock difference prediction method, device and computer equipment
CN109754013B (en) Electric power system hybrid measurement fusion method based on unscented Kalman filtering
CN108562290B (en) Navigation data filtering method and device, computer equipment and storage medium
Caballero-Águila et al. Centralized, distributed and sequential fusion estimation from uncertain outputs with correlation between sensor noises and signal
CN102914308A (en) Anti-outlier federated filtering method based on innovation orthogonality
CN116088303B (en) A Time-varying Recursive Estimation Method for Uncertain Complex Dynamic Network State
CN107292410A (en) tunnel deformation prediction method and device
CN113011475A (en) Distributed fusion algorithm considering correlated noise and random parameter matrix
Hajikhani et al. A recursive method for clock synchronization in asymmetric packet-based networks
CN116224746A (en) High-stability time reference establishment method for satellite-ground atomic clock fusion
Zhao et al. Gaussian filter for nonlinear networked systems with synchronously correlated noises and one-step randomly delayed measurements and multiple packet dropouts
CN116595897B (en) A method and device for state estimation of nonlinear dynamic systems based on message passing
Caballero-Águila et al. Optimal fusion filtering in multisensor stochastic systems with missing measurements and correlated noises
Zhao et al. Gaussian filter for nonlinear stochastic uncertain systems with correlated noises
CN114445459B (en) Continuous-discrete maximum correlation entropy target tracking method based on variable decibel leaf theory
CN114745041A (en) Intelligent clock synchronization method for distributed satellite network
Shao Distributed consensus Kalman filtering for asynchronous multi-rate sensor networks
CN119535946A (en) A Beidou and long wave fusion timing method and related device
CN103607181A (en) Spatially distributed variable order adaptive system identification method
CN106871905B (en) Gaussian filtering substitution framework combined navigation method under non-ideal condition
Zhao et al. Fusion estimation for nonlinear multi-sensor networked systems with packet loss compensation and correlated noises
CN107886058B (en) Noise-dependent two-stage volumetric Kalman filter estimation method and system

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
CB03 Change of inventor or designer information

Inventor after: Song Shenmin

Inventor after: Tan Liguo

Inventor after: Zhao Kai

Inventor after: Zhang Xiujie

Inventor before: Song Shenmin

Inventor before: Zhao Kai

Inventor before: Zhang Xiujie

Inventor before: Tan Liguo

CB03 Change of inventor or designer information
GR01 Patent grant
GR01 Patent grant