CN105550497A - High-precision ballistic correction method - Google Patents

Abstract

The invention discloses a high-precision ballistic correction method. The method is characterized by comprising the following steps: estimating and correcting the state noise and measurement noise by utilizing a noise estimator; completing the one-step prediction process by utilizing the noise estimator; and completing the ballistic correction process on the basis of the noise estimator and a fading factor. According to the high-precision ballistic correction method, in the quasi-linear optimum filtering process, the state value is predicted by utilizing the measurement data on one hand, and a filtering gain matrix is adjusted through the self-adaptive estimation of a process noise estimator and a measurement noise estimator and the characteristic parameters of the corrected noises on the other hand, so that the influences on the ballistic correction result due to the imprecise noise parameter estimation can be reduced; and at the same time, the fading factor is introduced to carry out reasonable fading on the old data, and the weight of the new measurement data is increased, so that the aims of reducing the filtering estimation error and restraining the filtering divergence caused by imprecise system model and sudden change of state are achieved.

Description

A kind of high-precision projectile correction method

Technical field

The present invention relates to a kind of high-precision projectile correction method, specifically one relate in projectile correction process suppress disperse and put forward high-precision method, especially non-linear, the factor such as time variation and noise of model trajectory is considered, be intended to suppress cause in ballistic data makeover process due to the modeling error of model trajectory and the evaluated error of noise statistics disperse, improve trajectory formed precision.The present invention effectively can ensure the stability of projectile correction process and the accuracy of ballistic trajectory formation, and the ballistic trajectory being applicable to Gun Position Radar forms system.

Background technology

Gun Position Radar utilizes the flight parameter of the bullet obtained to carry out projectile correction, forms ballistic trajectory to carry out target following.The bullet flown in large spatial domain is subject to the impact of various factors, and measured value and calculated value out of true etc. as the height of projectile flight, Mach number, aerodynamic force parameter can cause the model parameter of bullet to alter a great deal.In addition, also certain difference can be there is between the experience estimated value arranged in projectile correction process and the parameter of actual projectile flight.When using pseudo-linear optimal smoothing filtering algorithm to carry out projectile correction, the error that model trajectory discretize and linearization bring, model trajectory describe out of true, the process noise adopted and the statistical property of measurement noise and the deviation of actual conditions is large etc. that factor all may cause state estimation and bullet practical flight data to depart from, cause filtering accuracy to decline, the requirement of Practical Project cannot be met.The enchancement factors such as actual environment change also can cause systematic procedure noise and measurement noise ANOMALOUS VARIATIONS, the mistake that can produce gain matrix and predicated error variance in filtering is estimated, in this case gain matrix will lose suitable weighting effect, the adjustment correcting action of new metric data to state estimation will reduce gradually, thus there is filtering divergence phenomenon, cause projectile correction failure, cannot tracking target be continued.

Summary of the invention

Technical matters to be solved by this invention is, provides a kind of high-precision projectile correction method utilizing noise estimator to estimate state-noise and measurement noise and correct; Further, the invention provides one utilizes noise estimator to complete one-step prediction, and in reduction pseudo-linear optimal smoothing filtering, the factor such as uncertain noise statistics and other random disturbance is on the high-precision projectile correction method of the impact of projectile correction precision; Further, the invention provides one and complete projectile correction based on noise estimator and fading factor, fading factor carries out adaptive adjustment to the covariance of the one-step prediction error of state estimation and filter gain matrix, suppresses the high-precision projectile correction method of the Divergent Phenomenon caused due to state mutation in model trajectory out of true and filtering.

For solving the problems of the technologies described above, the technical solution used in the present invention is:

A kind of high-precision projectile correction method, it is characterized in that, comprising: utilize the method that noise estimator is estimated state-noise and measurement noise and corrected, this method comprises the following steps:

S01, model trajectory:

The state equation that model trajectory adopts is:

$\frac{d}{dt}X\left(t\right)=f\left(X\right(t\left)\right)+\mathrm{\Γ}\left(t\right)W\left(t\right)---\left(1\right)$

In formula (1), X (t) represents the state of bullet, and W (t) is process noise, and Γ (t) is that noise drives battle array;

The measurement equation that model trajectory adopts is:

Z(t)＝h(X(t))+V(t)(2)

In formula (2), Z (t) represents measurement variable, and V (t) is measurement noise;

Formula (1), formula (2) are obtained after carrying out discretize:

X _k+1＝G(X _k)+Γ(X _k)W _k(3)

Z _k+1＝h(X _k+1)+V _k+1(4)

In formula (3), formula (4), $G\left(X_k\right)=X_k+f\left(X_k\right)T+A\left(X_k\right)f\left(X_k\right)\frac{T^2}{2},$ X _k=X (t) | _t=kT, T=t _k+1-t _kfor sampling interval, t _kfor a kth sampling time; w _k=W (t) | _t=kT, V _k+1=V (t) | _{t=(k+1) T};

S02, the calculation of initial value of projectile correction:

When revising ballistic trajectory, definition status vector is:

X＝[xyzv _xv _yv _zca _c](5)

X, y and z represent the directive distance of bullet, height and lateral deviation respectively in formula (5); v _x, v _yand v _zrepresent three components of velocity of shot respectively; C and a _crepresent ballistic coefficient and projectile correction coefficient respectively;

If radar records the two initial point sampling values of shell for (x ₁, y ₁, z ₁) and (x ₂, y ₂, z ₂), obtain the initial valuation of state vector X thus for:

${\hat{X}}_0=\[\begin{array}{cccccccc}{x}_2& y_2& z_2& \frac{x_2-x_1}{T}& \frac{y_2-y_1}{T}& \frac{z_2-z_1}{T}& c& a_c\end{array}\]---\left(6\right)$

Initial prediction error variance matrix is:

In formula (7), for the variance of system stochastic error, it is the empirical estimating value of ballistic coefficient and projectile correction coefficient mean square deviation.

The present invention also comprises the process utilizing noise estimator to complete one-step prediction, and this process comprises the following steps:

S03, projectile correction:

S03-1, calculates one-step prediction smooth value:

S03-1-1, utilizes the calculating correlation parameter that predicts the outcome of initial value and formula (7),

Order ${\hat{X}}_k^{*}={\hat{X}}_k,(k=0,1,2,...)---\left(8\right)$

Computing mode function:

$f\left({\hat{X}}_k^{*}\right)=f\left(X\right(t\left)\right){|}_{X\left(t\right)={\hat{X}}_k^{*}}=\left[\begin{array}{c}{v}_x\\ v_y\\ v_z\\ -{\mathrm{cN}}_1{x}_{4r}+a_c{v}_x/v_r\\ -{\mathrm{cN}}_1{v}_y-g+a_c{v}_y/v_r\\ -{\mathrm{cN}}_1{x}_{6r}+B_Z+a_c{v}_z/v_r\\ 0\\ 0\end{array}\right]{|}_{X\left(t\right)={\hat{X}}_k^{*}}---\left(9\right)$

X in formula (9) _4r=v _x-w _xand x _6r=v _z-w _zrepresent bullet respectively in x and z direction by the actual speed after air speed influence, w _xand w _zfor wind speed is at the component in x and z direction, namely ρ is atmospheric density, and to be that bullet is maximum block area to S d is caliber, and m is Shell body quality; C _x(Ma) being coefficient of air resistance, is the function of Mach number Ma, and g is the acceleration of gravity at height above sea level y place, B _zfor side direction lift acceleration;

The partial derivative of computing mode function

$A\left({\hat{X}}_k^{*}\right)=\frac{\∂f\left(X\right(t\left)\right)}{\∂X\left(t\right)}{|}_{X\left(t\right)={\hat{X}}_k^{*}}---\left(10\right)$

Calculate one step state transition matrix

$\mathrm{\φ}\left({\hat{X}}_k^{*}\right)=\frac{\∂G\left(X\right)}{\∂X}{|}_{X={\hat{X}}_k^{*}}---\left(11\right)$

The variance matrix of computation process noise and the variance matrix of measurement noise:

${\hat{Q}}_k=E\left\{W_k{W}_k^H\right\}---\left(12\right)$

${\hat{R}}_k=E\left\{V_k{V}_k^H\right\}---\left(13\right)$

In formula (12), (13), H represents conjugate transpose;

S03-1-2, carries out predicted estimate:

Calculate one-step prediction estimated value:

$\begin{array}{c}{\hat{X}}_{k+1/k}=G\left({\hat{X}}_k^{*}\right)+\mathrm{\φ}\left({\hat{X}}_k^{*}\right)({\hat{X}}_k-{\hat{X}}_k^{*})\\ ={\hat{X}}_k^{*}+f\left({\hat{X}}_k^{*}\right)T+A\left({\hat{X}}_k^{*}\right)f\left({\hat{X}}_k^{*}\right)\frac{T^2}{2}+\mathrm{\φ}\left({\hat{X}}_k^{*}\right)({\hat{X}}_k-{\hat{X}}_k^{*})\end{array}---\left(14\right)$

Calculate one-step prediction error covariance:

${\hat{P}}_{k+1/k}=\mathrm{\φ}\left({\hat{X}}_k^{*}\right)P_k{\mathrm{\φ}}^H\left({\hat{X}}_k^{*}\right)+Q_k^{*}---\left(15\right)$

Calculate and measure transition matrix:

$B\left({\hat{X}}_{k+1/k}\right)=\frac{\∂h\left(X\right)}{\∂X}{|}_{X={\hat{X}}_{k+1/k}}---\left(16\right)$

Calculated gains matrix:

$K_k={\hat{P}}_{k+1/k}{B}^H\left({\hat{X}}_{k+1/k}\right){\[B\left({\hat{X}}_{k+1/k}\right){\hat{P}}_{k+1/k}{B}^H\left({\hat{X}}_{k+1/k}\right)+{\hat{R}}_k^{*}\]}^{-1}---\left(17\right)$

S03-1-3, utilizes noise estimator correction noise:

The variance matrix of process noise is revised:

${\hat{Q}}_{k+1}=(1-d_{k+1}){\hat{Q}}_k+d_{k+1}(K_k{\tilde{Z}}_{k+1}{\tilde{Z}}_{k+1}^T{K}_k^T+P_{k+1/k}-{\mathrm{\Φ}}_k{P}_k{\mathrm{\Φ}}_k^T)---\left(18\right)$

The variance matrix of measurement noise is revised:

${\hat{R}}_{k+1}=(1-d_{k+1}){\hat{R}}_k+d_{k+1}({\tilde{Z}}_{k+1}{\tilde{Z}}_{k+1}^T-B_{k+1/k}{\hat{P}}_{k+1/k}{B}_{k+1/k}^T)---\left(19\right)$

In formula (18), (19), residual error newly ceases be calculated as:

${\tilde{Z}}_{k+1}=Z_{k+1}-H\left(X_{k+1/k}\right)---\left(20\right)$

In formula (18), (19), d _k=(1-α)/(1-α ^k), 0 < α < 1 is forgetting factor;

S03-1-4, carries out self-adaptive smooth filtering:

Calculation of filtered smooth value:

$\begin{array}{c}{\hat{X}}_{k/k+1}={\hat{X}}_k+P_k{\mathrm{\&phi;}}_k\left({\hat{X}}_k^{*}\right)B^T\left({\hat{X}}_{k+1/k}\right)\\ \&CenterDot;{\&lsqb;B\left({\hat{X}}_{k+1/k}\right){\hat{P}}_{k+1/k}{B}^T\left({\hat{X}}_{k+1/k}\right)+R_{k+1}\&rsqb;}^{-1}\\ \&CenterDot;\&lsqb;Z_{k+1}-H\left({\hat{X}}_{k+1/k}\right)\&rsqb;\end{array}---\left(21\right)$

S03-1-5, iterative computation

Order repeat S03-1-1 to S03-1-5 computation process at least three times, can try to achieve after so carrying out at least three iteration approximate value, complete the calculating of one-step prediction smooth value.

The present invention also comprises the process completing projectile correction based on noise estimator and fading factor, and this process comprises:

S03-2, utilizes fading factor correction filter result:

S03-2-1, calculates fading factor, makes fading factor

${\mathrm{\&lambda;}}_k=\frac{Tr\left(N_k\right)}{Tr\left(M_k\right)}=\left\{\begin{array}{cc}{\mathrm{\&lambda;}}_k,& {\mathrm{\&lambda;}}_k\&GreaterEqual;1\\ 1,& {\mathrm{\&lambda;}}_k<1\end{array}\right.---\left(22\right)$

In formula (22), Tr represents and asks matrix trace; $N_k=b_k-H\left({\hat{X}}_{k/k+1}\right)Q_k^{*}{H}^T\left({\hat{X}}_{k/k+1}\right)-l_k{R}_k;$ $M_k=H\left({\hat{X}}_{k/k+1}\right)\mathrm{\&phi;}\left({\hat{X}}_{k/k+1}\right){\hat{P}}_{k+1/k}{\mathrm{\&phi;}}^T\left({\hat{X}}_{k/k+1}\right)H^T\left({\hat{X}}_{k/k+1}\right);$ l _k＝1-d _k；

S03-2-2, carries out one-step prediction modified result,

Calculate the error co-variance matrix of the one-step prediction of band fading factor

$P_{k+1/k}={\mathrm{\&lambda;}}_k\mathrm{\&phi;}\left({\hat{X}}_k^{*}\right)P_k{\mathrm{\&phi;}}^T\left({\hat{X}}_k^{*}\right)+Q_k^{*}---\left(23\right)$

Calculated gains matrix

$K_{k+1}=P_{k+1/k}{B}^T\left({\hat{X}}_{k+1/k}\right){\&lsqb;B\left({\hat{X}}_{k+1/k}\right)P_{k+1/k}{B}^T\left({\hat{X}}_{k+1/k}\right)+R_{k+1}\&rsqb;}^{-1}---\left(24\right)$

Calculation of filtered estimated value

${\hat{X}}_{k+1}={\hat{X}}_{k+1/k}+K_{k+1}\&lsqb;Z_{k+1}-H\left({\hat{X}}_{k+1/k}\right)\&rsqb;---\left(25\right)$

Calculation of filtered error covariance matrix

$P_{k+1}=\&lsqb;I-K_{k+1}B\left({\hat{X}}_{k+1/k}\right)\&rsqb;P_{k+1/k}{\&lsqb;I-K_{k+1}B\left({\hat{X}}_{k+1/k}\right)\&rsqb;}^H+K_{k+1}{R}_{k+1}{K}_{k+1}^H---\left(26\right)$

S03-3, makes k+1 → k, repeats S03-1 to S03-3 step, completes projectile correction process.

The number of times of described iterative computation is three times.

The high-precision projectile correction method of one provided by the invention, the projectile correction method of the present invention's design is in pseudo-linear optimal filtering process, metric data is utilized to predict state value on the one hand, the characterisitic parameter passing through process noise estimator and the adaptive estimation of measurement noise estimator and correction noise on the other hand, to regulate filter gain matrix, reduces the impact that noise parameter estimates lax pair projectile correction result; Introduce fading factor reasonably to fade to outmoded data simultaneously, increase the weight of new metric data, by carrying out real-time online adjustment to the covariance matrix of predicated error and filter gain matrix, reach the object reducing filtering evaluated error, suppress the filtering divergence caused due to system model out of true and state mutation.

Figure of description

Fig. 1 is the projectile correction process flow diagram that the present invention is based on self-adaptation fading factor;

Fig. 2 is the Divergent Phenomenon figure of prior art based on pseudo-linear optimal smoothing filtering algorithm projectile correction;

The simulation result figure of high precision projectile correction method in Fig. 3 the present invention.

Embodiment

Below in conjunction with accompanying drawing, the present invention is further described.

As shown in FIG. 1 to 3, in order to verify the projectile correction method validity based on noise estimator and fading factor of the present invention, 1100m/s is elected as at shell initial velocity, directive angle is elected 25o condition Imitating as and is generated radar measurement data, adopt ballistic trajectory correction algorithm of the present invention and pseudo-linear optimal smoothing filtering algorithm respectively, Filtering Analysis is carried out to the shell flying quality that radargrammetry obtains, the three-dimensional artificial result of the ballistic trajectory formed is as Fig. 2, shown in Fig. 3, wherein " emulated data " is the radar measurement data of simulation, the ballistic trajectory that " filtering data " is estimated for filtering.As shown in Figure 2, in Fig. 2, there is filtering divergence phenomenon, and used the projectile correction method based on noise estimator and fading factor effectively to inhibit in Fig. 3 to disperse.

A kind of high-precision projectile correction method, it is characterized in that, comprising: utilize the method that noise estimator is estimated state-noise and measurement noise and corrected, this method comprises the following steps:

S01, model trajectory:

The state equation that model trajectory adopts is:

$\frac{d}{dt}X\left(t\right)=f\left(X\right(t\left)\right)+\mathrm{\&Gamma;}\left(t\right)W\left(t\right)---\left(1\right)$

In formula (1), X (t) represents the state of bullet, and W (t) is process noise, and Γ (t) is that noise drives battle array;

The measurement equation that model trajectory adopts is:

Z(t)＝h(X(t))+V(t)(2)

In formula (2), Z (t) represents measurement variable, and V (t) is measurement noise;

Formula (1), formula (2) are obtained after carrying out discretize:

X _k+1＝G(X _k)+Γ(X _k)W _k(3)

Z _k+1＝h(X _k+1)+V _k+1(4)

In formula (3), formula (4), $G\left(X_k\right)=X_k+f\left(X_k\right)T+A\left(X_k\right)f\left(X_k\right)\frac{T^2}{2},$ X _k=X (t) | _t=kT, T=t _k+1-t _kfor sampling interval, t _kfor a kth sampling time; w _k=W (t) | _t=kT, V _k+1=V (t) | _{t=(k+1) T};

S02, the calculation of initial value of projectile correction:

When revising ballistic trajectory, definition status vector is:

X＝[xyzv _xv _yv _zca _c](5)

X, y and z represent the directive distance of bullet, height and lateral deviation respectively in formula (5); v _x, v _yand v _zrepresent three components of velocity of shot respectively; C and a _crepresent ballistic coefficient and projectile correction coefficient respectively;

If radar records the two initial point sampling values of shell for (x ₁, y ₁, z ₁) and (x ₂, y ₂, z ₂), obtain the initial valuation of state vector X thus for:

${\hat{X}}_0=\&lsqb;\begin{array}{cccccccc}{x}_2& y_2& z_2& \frac{x_2-x_1}{T}& \frac{y_2-y_1}{T}& \frac{z_2-z_1}{T}& c& a_c\end{array}\&rsqb;---\left(6\right)$

Initial prediction error variance matrix is:

In formula (7), for the variance of system stochastic error, it is the empirical estimating value of ballistic coefficient and projectile correction coefficient mean square deviation.

The present invention also comprises the process utilizing noise estimator to complete one-step prediction, and this process comprises the following steps:

S03, projectile correction:

S03-1, calculates one-step prediction smooth value:

S03-1-1, utilizes the calculating correlation parameter that predicts the outcome of initial value and formula (7),

Order ${\hat{X}}_k^{*}={\hat{X}}_k,(k=0,1,2,...)---\left(8\right)$

Computing mode function:

$f\left({\hat{X}}_k^{*}\right)=f\left(X\right(t\left)\right){|}_{X\left(t\right)={\hat{X}}_k^{*}}=\left[\begin{array}{c}{v}_x\\ v_y\\ v_z\\ -{\mathrm{cN}}_1{x}_{4r}+a_c{v}_x/v_r\\ -{\mathrm{cN}}_1{v}_y-g+a_c{v}_y/v_r\\ -{\mathrm{cN}}_1{x}_{6r}+B_Z+a_c{v}_z/v_r\\ 0\\ 0\end{array}\right]{|}_{X\left(t\right)={\hat{X}}_k^{*}}---\left(9\right)$

X in formula (9) _4r=v _x-w _xand x _6r=v _z-w _zrepresent bullet respectively in x and z direction by the actual speed after air speed influence, w _xand w _zfor wind speed is at the component in x and z direction, namely ρ is atmospheric density, and to be that bullet is maximum block area to S d is caliber, and m is Shell body quality; C _x(Ma) being coefficient of air resistance, is the function of Mach number Ma, and g is the acceleration of gravity at height above sea level y place, B _zfor side direction lift acceleration;

The partial derivative of computing mode function

$A\left({\hat{X}}_k^{*}\right)=\frac{\&part;f\left(X\right(t\left)\right)}{\&part;X\left(t\right)}{|}_{X\left(t\right)={\hat{X}}_k^{*}}---\left(10\right)$

Calculate one step state transition matrix

$\mathrm{\&phi;}\left({\hat{X}}_k^{*}\right)=\frac{\&part;G\left(X\right)}{\&part;X}{|}_{X={\hat{X}}_k^{*}}---\left(11\right)$

The variance matrix of computation process noise and the variance matrix of measurement noise:

${\hat{Q}}_k=E\left\{W_k{W}_k^H\right\}---\left(12\right)$

${\hat{R}}_k=E\left\{V_k{V}_k^H\right\}---\left(13\right)$

In formula (12), (13), H represents conjugate transpose;

S03-1-2, carries out predicted estimate:

Calculate one-step prediction estimated value:

$\begin{array}{c}{\hat{X}}_{k+1/k}=G\left({\hat{X}}_k^{*}\right)+\mathrm{\&phi;}\left({\hat{X}}_k^{*}\right)({\hat{X}}_k-{\hat{X}}_k^{*})\\ ={\hat{X}}_k^{*}+f\left({\hat{X}}_k^{*}\right)T+A\left({\hat{X}}_k^{*}\right)f\left({\hat{X}}_k^{*}\right)\frac{T^2}{2}+\mathrm{\&phi;}\left({\hat{X}}_k^{*}\right)({\hat{X}}_k-{\hat{X}}_k^{*})\end{array}---\left(14\right)$

Calculate one-step prediction error covariance:

${\hat{P}}_{k+1/k}=\mathrm{\&phi;}\left({\hat{X}}_k^{*}\right)P_k{\mathrm{\&phi;}}^H\left({\hat{X}}_k^{*}\right)+Q_k^{*}---\left(15\right)$

Calculate and measure transition matrix:

$B\left({\hat{X}}_{k+1/k}\right)=\frac{\&part;h\left(X\right)}{\&part;X}{|}_{X={\hat{X}}_{k+1/k}}---\left(16\right)$

Calculated gains matrix:

$K_k={\hat{P}}_{k+1/k}{B}^H\left({\hat{X}}_{k+1/k}\right){\&lsqb;B\left({\hat{X}}_{k+1/k}\right){\hat{P}}_{k+1/k}{B}^H\left({\hat{X}}_{k+1/k}\right)+{\hat{R}}_k^{*}\&rsqb;}^{-1}---\left(17\right)$

S03-1-3, utilizes noise estimator correction noise:

The variance matrix of process noise is revised:

${\hat{Q}}_{k+1}=(1-d_{k+1}){\hat{Q}}_k+d_{k+1}(K_k{\tilde{Z}}_{k+1}{\tilde{Z}}_{k+1}^T{K}_k^T+P_{k+1/k}-{\mathrm{\&Phi;}}_k{P}_k{\mathrm{\&Phi;}}_k^T)---\left(18\right)$

The variance matrix of measurement noise is revised:

${\hat{R}}_{k+1}=(1-d_{k+1}){\hat{R}}_k+d_{k+1}({\tilde{Z}}_{k+1}{\tilde{Z}}_{k+1}^T-B_{k+1/k}{\hat{P}}_{k+1/k}{B}_{k+1/k}^T)---\left(19\right)$

In formula (18), (19), residual error newly ceases be calculated as:

${\tilde{Z}}_{k+1}=Z_{k+1}-H\left(X_{k+1/k}\right)---\left(20\right)$

In formula (18), (19), d _k=(1-α)/(1-α ^k), 0 < α < 1 is forgetting factor;

S03-1-4, carries out self-adaptive smooth filtering:

Calculation of filtered smooth value:

$\begin{array}{c}{\hat{X}}_{k/k+1}={\hat{X}}_k+P_k{\mathrm{\&phi;}}_k\left({\hat{X}}_k^{*}\right)B^T\left({\hat{X}}_{k+1/k}\right)\\ \&CenterDot;{\&lsqb;B\left({\hat{X}}_{k+1/k}\right){\hat{P}}_{k+1/k}{B}^T\left({\hat{X}}_{k+1/k}\right)+R_{k+1}\&rsqb;}^{-1}\\ \&CenterDot;\&lsqb;Z_{k+1}-H\left({\hat{X}}_{k+1/k}\right)\&rsqb;\end{array}---\left(21\right)$

S03-1-5, iterative computation

Order repeat S03-1-1 to S03-1-5 computation process at least three times, can try to achieve after so carrying out at least three iteration approximate value, complete the calculating of one-step prediction smooth value.

The present invention also comprises the process completing projectile correction based on noise estimator and fading factor, and this process comprises:

S03-2, utilizes fading factor correction filter result:

S03-2-1, calculates fading factor, makes fading factor

${\mathrm{\&lambda;}}_k=\frac{Tr\left(N_k\right)}{Tr\left(M_k\right)}=\left\{\begin{array}{cc}{\mathrm{\&lambda;}}_k,& {\mathrm{\&lambda;}}_k\&GreaterEqual;1\\ 1,& {\mathrm{\&lambda;}}_k<1\end{array}\right.---\left(22\right)$

In formula (22), Tr represents and asks matrix trace; $N_k=b_k-H\left({\hat{X}}_{k/k+1}\right)Q_k^{*}{H}^T\left({\hat{X}}_{k/k+1}\right)-l_k{R}_k;$ $M_k=H\left({\hat{X}}_{k/k+1}\right)\mathrm{\&phi;}\left({\hat{X}}_{k/k+1}\right){\hat{P}}_{k+1/k}{\mathrm{\&phi;}}^T\left({\hat{X}}_{k/k+1}\right)H^T\left({\hat{X}}_{k/k+1}\right);$ l _k＝1-d _k；

S03-2-2, carries out one-step prediction modified result,

Calculate the error co-variance matrix of the one-step prediction of band fading factor

$P_{k+1/k}={\mathrm{\&lambda;}}_k\mathrm{\&phi;}\left({\hat{X}}_k^{*}\right)P_k{\mathrm{\&phi;}}^T\left({\hat{X}}_k^{*}\right)+Q_k^{*}---\left(23\right)$

Calculated gains matrix

$K_{k+1}=P_{k+1/k}{B}^T\left({\hat{X}}_{k+1/k}\right){\&lsqb;B\left({\hat{X}}_{k+1/k}\right)P_{k+1/k}{B}^T\left({\hat{X}}_{k+1/k}\right)+R_{k+1}\&rsqb;}^{-1}---\left(24\right)$

Calculation of filtered estimated value

${\hat{X}}_{k+1}={\hat{X}}_{k+1/k}+K_{k+1}\&lsqb;Z_{k+1}-H\left({\hat{X}}_{k+1/k}\right)\&rsqb;---\left(25\right)$

Calculation of filtered error covariance matrix

$P_{k+1}=\&lsqb;I-K_{k+1}B\left({\hat{X}}_{k+1/k}\right)\&rsqb;P_{k+1/k}{\&lsqb;I-K_{k+1}B\left({\hat{X}}_{k+1/k}\right)\&rsqb;}^H+K_{k+1}{R}_{k+1}{K}_{k+1}^H---\left(26\right)$

S03-3, makes k+1 → k, repeats S03-1 to S03-3 step, completes projectile correction process.

The number of times of described iterative computation is three times.

The above is only the preferred embodiment of the present invention; be noted that for those skilled in the art; under the premise without departing from the principles of the invention, can also make some improvements and modifications, these improvements and modifications also should be considered as protection scope of the present invention.

Claims (4)

1. a high-precision projectile correction method, it is characterized in that, comprising: utilize the method that noise estimator is estimated state-noise and measurement noise and corrected, this method comprises the following steps:

S01, model trajectory:

The state equation that model trajectory adopts is:

$\frac{d}{dt}X\left(t\right)=f\left(X\right(t\left)\right)+\mathrm{\&Gamma;}\left(t\right)W\left(t\right)---\left(1\right)$

In formula (1), X (t) represents the state of bullet, and W (t) is process noise, and Γ (t) is that noise drives battle array;

The measurement equation that model trajectory adopts is:

Z(t)＝h(X(t))+V(t)(2)

In formula (2), Z (t) represents measurement variable, and V (t) is measurement noise;

Formula (1), formula (2) are obtained after carrying out discretize:

X _k+1＝G(X _k)+Γ(X _k)W _k(3)

Z _k+1＝h(X _k+1)+V _k+1(4)

In formula (3), formula (4), $G\left(X_k\right)=X_k+f\left(X_k\right)T+A\left(X_k\right)f\left(X_k\right)\frac{T^2}{2},$ X _k=X (t) | _t=kT, T=t _k+1-t _kfor sampling interval, t _kfor a kth sampling time; v _k+1=V (t) | _{t=(k+1) T};

S02, the calculation of initial value of projectile correction:

When revising ballistic trajectory, definition status vector is:

X＝[xyzv _xv _yv _zca _c](5)

X, y and z represent the directive distance of bullet, height and lateral deviation respectively in formula (5); v _x, v _yand v _zrepresent three components of velocity of shot respectively; C and a _crepresent ballistic coefficient and projectile correction coefficient respectively;

If radar records the two initial point sampling values of shell for (x ₁, y ₁, z ₁) and (x ₂, y ₂, z ₂), obtain the initial valuation of state vector X thus for:

${\hat{X}}_0=\&lsqb;\begin{array}{cccccccc}{x}_2& y_2& z_2& \frac{x_2-x_1}{T}& \frac{y_2-y_1}{T}& \frac{z_2-z_1}{T}& c& a_c\end{array}\&rsqb;---\left(6\right)$

Initial prediction error variance matrix is:

$P_0=\left[\begin{array}{cccccccc}{\mathrm{\&sigma;}}_x^2& & & & & & & \\ & {\mathrm{\&sigma;}}_y^2& & & & 0& & \\ & & {\mathrm{\&sigma;}}_z^2& & & & & \\ & & & 2{\mathrm{\&sigma;}}_x^2/T^2& & & & \\ & & & & 2{\mathrm{\&sigma;}}_y^2/T^2& & & \\ & & 0& & & 2{\mathrm{\&sigma;}}_z^2/T^2& & \\ & & & & & & {\mathrm{\&sigma;}}_c^2& \\ & & & & & & & {\mathrm{\&sigma;}}_{ac}^2\end{array}\right]---\left(7\right)$

In formula (7), for the variance of system stochastic error, it is the empirical estimating value of ballistic coefficient and projectile correction coefficient mean square deviation.

2. the high-precision projectile correction method of one according to claim 1, it is characterized in that: also comprise the process utilizing noise estimator to complete one-step prediction, this process comprises the following steps:

S03, projectile correction:

S03-1, calculates one-step prediction smooth value:

S03-1-1, utilizes the calculating correlation parameter that predicts the outcome of initial value and formula (7),

Order ${\hat{X}}_k^{*}={\hat{X}}_k(k=0,1,2,...)---\left(8\right)$

Computing mode function:

$f\left({\hat{X}}_k^{*}\right)=f\left(X\right(t\left)\right){|}_{X\left(t\right)={\hat{X}}_k^{*}}=\left[\begin{array}{c}{v}_x\\ v_y\\ v_z\\ -{\mathrm{cN}}_1{x}_{4r}+a_c{v}_x/v_r\\ -{\mathrm{cN}}_1{v}_y-g+a_c{v}_y/v_r\\ -{\mathrm{cN}}_1{x}_{6r}+B_Z+a_c{v}_z/v_r\\ 0\\ 0\end{array}\right]{|}_{X\left(t\right)={\hat{X}}_k^{*}}---\left(9\right)$

X in formula (9) _4r=v _x-w _xand x _6r=v _z-w _zrepresent bullet respectively in x and z direction by the actual speed after air speed influence, w _xand w _zfor wind speed is at the component in x and z direction, namely ρ is atmospheric density, and to be that bullet is maximum block area to S d is caliber, and m is Shell body quality; C _x(Ma) being coefficient of air resistance, is the function of Mach number Ma, and g is the acceleration of gravity at height above sea level y place, B _zfor side direction lift acceleration;

The partial derivative of computing mode function

$A\left({\hat{X}}_k^{*}\right)=\frac{\&part;f\left(X\right(t\left)\right)}{\&part;X\left(t\right)}{|}_{X\left(t\right)={\hat{X}}_k^{*}}---\left(10\right)$

Calculate one step state transition matrix

$\mathrm{\&phi;}\left({\hat{X}}_k^{*}\right)=\frac{\&part;G\left(X\right)}{\&part;X}{|}_{X={\hat{X}}_k^{*}}---\left(11\right)$

The variance matrix of computation process noise and the variance matrix of measurement noise:

${\hat{Q}}_k=E\left\{W_k{W}_k^H\right\}---\left(12\right)$

${\hat{R}}_k=E\left\{V_k{V}_k^H\right\}---\left(13\right)$

In formula (12), (13), H represents conjugate transpose;

S03-1-2, carries out predicted estimate:

Calculate one-step prediction estimated value:

$\begin{array}{c}{\hat{X}}_{k+1/k}=G\left({\hat{X}}_k^{*}\right)+\mathrm{\&phi;}\left({\hat{X}}_k^{*}\right)({\hat{X}}_k-{\hat{X}}_k^{*})\\ ={\hat{X}}_k^{*}+f\left({\hat{X}}_k^{*}\right)T+A\left({\hat{X}}_k^{*}\right)f\left({\hat{X}}_k^{*}\right)\frac{T^2}{2}+\mathrm{\&phi;}\left({\hat{X}}_k^{*}\right)({\hat{X}}_k-{\hat{X}}_k^{*})\end{array}---\left(14\right)$

Calculate one-step prediction error covariance:

${\hat{P}}_{k+1/k}=\mathrm{\&phi;}\left({\hat{X}}_k^{*}\right)P_k{\mathrm{\&phi;}}^H\left({\hat{X}}_k^{*}\right)+Q_k^{*}---\left(15\right)$

Calculate and measure transition matrix:

$B\left({\hat{X}}_{k+1/k}\right)=\frac{\&part;h\left(X\right)}{\&part;X}{|}_{X={\hat{X}}_{k+1/k}}---\left(16\right)$

Calculated gains matrix:

$K_k={\hat{P}}_{k+1/k}{B}^H\left({\hat{X}}_{k+1/k}\right){\&lsqb;B\left({\hat{X}}_{k+1/k}\right){\hat{P}}_{k+1/k}{B}^H\left({\hat{X}}_{k+1/k}\right)+{\hat{R}}_k^{*}\&rsqb;}^{-1}---\left(17\right)$

S03-1-3, utilizes noise estimator correction noise:

The variance matrix of process noise is revised:

${\hat{Q}}_{k+1}=(1-d_{k+1}){\hat{Q}}_k+d_{k+1}(K_k{\tilde{Z}}_{k+1}{\tilde{Z}}_{k+1}^T{K}_k^T+P_{k+1/k}-{\mathrm{\&Phi;}}_k{P}_k{\mathrm{\&Phi;}}_k^T)---\left(18\right)$

The variance matrix of measurement noise is revised:

${\hat{R}}_{k+1}=(1-d_{k+1}){\hat{R}}_k+d_{k+1}({\tilde{Z}}_{k+1}{\tilde{Z}}_{k+1}^T-B_{k+1/k}{\hat{P}}_{k+1/k}{B}_{k+1/k}^T)---\left(19\right)$

In formula (18), (19), residual error newly ceases be calculated as:

${\tilde{Z}}_{k+1}=Z_{k+1}-H\left(X_{k+1/k}\right)---\left(20\right)$

In formula (18), (19), d _k=(1-α)/(1-α ^k), 0 < α < 1 is forgetting factor;

S03-1-4, carries out self-adaptive smooth filtering:

Calculation of filtered smooth value:

$\begin{array}{c}{\hat{X}}_{k/k+1}={\hat{X}}_k+P_k{\mathrm{\&phi;}}_k\left({\hat{X}}_k^{*}\right)B^T\left({\hat{X}}_{k+1/k}\right)\\ \&CenterDot;{\&lsqb;B\left({\hat{X}}_{k+1/k}\right){\hat{P}}_{k+1/k}{B}^T\left({\hat{X}}_{k+1/k}\right)+R_{k+1}\&rsqb;}^{-1}\\ \&CenterDot;\&lsqb;Z_{k+1}-H\left({\hat{X}}_{k+1/k}\right)\&rsqb;\end{array}---\left(21\right)$

S03-1-5, iterative computation

Order repeat S03-1-1 to S03-1-5 computation process at least three times, can try to achieve after so carrying out at least three iteration approximate value, complete the calculating of one-step prediction smooth value.

3. the high-precision projectile correction method of one according to claim 2, it is characterized in that: also comprise the process completing projectile correction based on noise estimator and fading factor, this process comprises:

S03-2, utilizes fading factor correction filter result:

S03-2-1, calculates fading factor, makes fading factor

${\mathrm{\&lambda;}}_k=\frac{Tr\left(N_k\right)}{Tr\left(M_k\right)}=\left\{\begin{array}{cc}{\mathrm{\&lambda;}}_k,& {\mathrm{\&lambda;}}_k\&GreaterEqual;1\\ 1,& {\mathrm{\&lambda;}}_k<1\end{array}\right.---\left(22\right)$

In formula (22), Tr represents and asks matrix trace; $N_k=b_k-H\left({\hat{X}}_{k/k+1}\right)Q_k^{*}{H}^T\left({\hat{X}}_{k/k+1}\right)-l_k{R}_k;$

$M_k=H\left({\hat{X}}_{k/k+1}\right)\mathrm{\&phi;}\left({\hat{X}}_{k/k+1}\right){\hat{P}}_{k+1/k}{\mathrm{\&phi;}}^T\left({\hat{X}}_{k/k+1}\right)H^T\left({\hat{X}}_{k/k+1}\right);l_k=1-d_k;$

S03-2-2, carries out one-step prediction modified result,

Calculate the error co-variance matrix of the one-step prediction of band fading factor

$P_{k+1/k}={\mathrm{\&lambda;}}_k\mathrm{\&phi;}\left({\hat{X}}_k^{*}\right)P_k{\mathrm{\&phi;}}^T\left({\hat{X}}_k^{*}\right)+Q_k^{*}---\left(23\right)$

Calculated gains matrix

$K_{k+1}=P_{k+1/k}{B}^T\left({\hat{X}}_{k+1/k}\right){\&lsqb;B\left({\hat{X}}_{k+1/k}\right)P_{k+1/k}{B}^T\left({\hat{X}}_{k+1/k}\right)+R_{k+1}\&rsqb;}^{-1}---\left(24\right)$

Calculation of filtered estimated value

${\hat{X}}_{k+1}={\hat{X}}_{k+1/k}+K_{k+1}\&lsqb;Z_{k+1}-H\left({\hat{X}}_{k+1/k}\right)\&rsqb;---\left(25\right)$

Calculation of filtered error covariance matrix

$P_{k+1}=\&lsqb;I-K_{k+1}B\left({\hat{X}}_{k+1/k}\right)\&rsqb;P_{k+1/k}{\&lsqb;I-K_{k+1}B\left({\hat{X}}_{k+1/k}\right)\&rsqb;}^H+K_{k+1}{R}_{k+1}{K}_{k+1}^H---\left(26\right)$

S03-3, makes k+1 → k, repeats S03-1 to S03-3 step, completes projectile correction process.

4. the high-precision projectile correction method of one according to claim 2, is characterized in that: the number of times of described iterative computation is three times.