CN103809439B - A kind of be applicable to remote control system controling parameters pass rank discrimination method - Google Patents

Abstract

The present invention relates to a kind of be applicable to remote control system controling parameters pass rank discrimination method, the remote operating process simplification of sampling different by input/output signal is a pair of rate sampling system, utilize the cycle submodel of method establishment one based on minimum sampling of resampling, and the intermediateness value of resampling system is estimated with this, the resampling cycle is greater than owing to exporting the sampling period, the output estimation difference of system can not directly revise intermediateness value, so the method that there is employed herein state change difference and difference is revised the state estimation in Kalman Prediction process, the present invention can carry out identification to the two rate remote control system parameters in the known and unknown two kinds of situations of state value, the feasibility of simulation results show the method and validity.

Description

A kind of be applicable to remote control system controling parameters pass rank discrimination method

Technical field

The invention belongs to system control parameters identification field, be specifically related to a kind of be applicable to remote control system controling parameters pass rank discrimination method.

Background technology

Two rate system is the general designation to the input sample frequency system different with output sampling frequency rate.Two rate system is not only widely used in Chemical Engineering Process Control field, also have great importance in spationautics application, such as in Space teleoperation process, by the input control signal of main side, the space manipulator from end is controlled, input control signal is position and the torque signals of main side hand controller, output feedback signal comprises the position of remote operating robotic arm joint, speed and torque signals etc., generally, main side control signal frequency is higher than the frequency receiving feedback signal, this is relevant with the ability of the processing signals of spaceborne computer on the one hand, also relevant with information capacity with the sample frequency of sensor on the other hand, can by a pair of rate that the is reduced to system of this remote operating process according to input control signal and output feedback signal, by carrying out identification to the controling parameters of this system thus drawing the control output relation of remote operating process.

The discrimination method being directed to many rates system and two rate system at present mainly contains: strong fraction transfer function model discrimination method, frequency-domain model discrimination method, wavelet model discrimination method, Stute space model identification method etc.Wherein, state-space method is an important method in two rate identification, it mainly can not be directly used in identification for rate data two in separate manufacturing firms model, single rate system model constant when two rate sampling system being converted into by " lift technique ", and with reference to MOESP, the subspace state space system identifications such as N4SID and CAV carry out identification to systematic parameter, although LPTV is converted into LTI system by lift technique, the direct application of conveniently recognized method, but also can cause the increase of number of parameters to be identified and calculated amount and the undesirable of identification precision, in addition, framework identified parameters cycle length of " lifting " once, the convergent cycle of identified parameters can be caused so long, speed of convergence is slower.Therefore, the method that needs one are new solves and utilizes lift technique to carry out parameter identification institute problems faced to two rate system.

Summary of the invention

The technical matters solved

In order to avoid the deficiencies in the prior art part, the present invention propose a kind of be applicable to remote control system controling parameters pass rank discrimination method.

Technical scheme

What be applicable to remote control system controling parameters passs a rank discrimination method, it is characterized in that step is as follows:

Step 1: gather remote operating input control signal and output feedback signal, described input control signal is position and the torque signals of main side hand controller; Described output feedback signal comprises the position of remote operating robotic arm joint, speed and torque signals;

Be a pair of rate system model by remote operating process simplification, its state-space model is:

$\left\{\begin{array}{c}x(K+1)=\mathrm{Ax}\left(K\right)+\mathrm{Bu}\left(K\right)+w\left(K\right)\\ y\left(K\right)=\mathrm{Cx}\left(K\right)+\mathrm{Du}\left(K\right)+v\left(K\right)\end{array}\right.;$

Wherein, A, B, C, D are systematic parameter, and K represents system the K time sampling, and x (K) is state value, and u (K) is input signal, and y (K) is output signal;

The control inputs sampling period of system is made to be T ₁with export sampling period be T ₂, h is T ₁and T ₂highest common factor, T ₁and T ₂meet: T ₁=ph, T ₂=qh and p and q two number is relatively prime, w (K) and v (K) represents process noise and measurement noises respectively, and the system model of system is: $\left\{\begin{array}{c}x(\mathrm{kh}+h)=\mathrm{ax}\left(\mathrm{kh}\right)+\mathrm{bu}\left(\mathrm{kh}\right)+r\left(\mathrm{kh}\right)\\ y\left(\mathrm{kh}\right)=\mathrm{cx}\left(\mathrm{kh}\right)+\mathrm{du}\left(\mathrm{kh}\right)+v\left(\mathrm{kh}\right)\end{array}\right.,$ Wherein, a, b, c, d are submodel parameter, and k is the kth time sampling of submodel, and r (kh) and v (kh) is respectively submodel process noise and measurement noises;

Parameters relationship between two models is: $a=A^{\frac{1}{p}},b=B/\underset{i=0}{\overset{p-1}{\mathrm{\Σ}}}{A}^{\frac{i}{p}},c=C,d=D,r\left(\mathrm{kh}\right)=w\left(\mathrm{kh}\right)/\underset{i=0}{\overset{p}{\mathrm{\Σ}}}{A}^{\frac{i}{p}};$

Step 2: make residual values o (kh) its expression formula be:

$\left(\mathrm{kh}\right)=(y\left(\mathrm{kh}\right)-\hat{c}\left(\mathrm{kh}\right)x\left(\mathrm{kh}\right)-\hat{d}\left(\mathrm{kh}\right)u\left(\mathrm{kh}\right))\underset{i=0}{\overset{q}{\mathrm{\Σ}}}{\hat{a}}^i$

Its update cycle is identical with the cycle qh of output variable y (kh), varivance matrix and error co-variance matrix expression formula are M (kh)=x (kh) x (kh-h), P (kh)=x (kh) x (kh), wherein x (kh)=x (kh)-x (kh), the estimated value that x (kh) is x (kh), order $\widehat{\mathrm{\θ}}\left(\mathrm{kh}\right)={[\hat{c}\left(\mathrm{kh}\right),d\left(\mathrm{kh}\right)]}^T,$ utilization pushes through recurrence equation:

Wherein, a (kh), d (kh) is respectively the estimated value of a (kh), b (kh), a (kh), d (kh), then the more new formula of Kalman filtering gain matrix L (kh), L1 (kh), P (kh) and M (kh) is:

$\begin{array}{c}L\left(\mathrm{kh}\right)=[a\left(\mathrm{kh}\right)(M\left(\mathrm{kh}\right){\hat{c}}^T(\mathrm{kh}-h)-P\left(\mathrm{kh}\right){\hat{c}}^T\left(\mathrm{kh}\right))+L1\left(\mathrm{kh}\right)\mathrm{Ro}][2\mathrm{Ro}+\hat{c}\left(\mathrm{kh}\right)P\left(\mathrm{kh}\right)\\ {\hat{c}}^T\left(\mathrm{kh}\right)-2\hat{c}\left(\mathrm{kh}\right)M\left(\mathrm{kh}\right){\hat{c}}^T(\mathrm{kh}-h)+\hat{c}(\mathrm{kh}-h)P(\mathrm{kh}-h){\hat{c}}^T(\mathrm{kh}-h){]}^{-1}\\ L1\left(\mathrm{kh}\right)=K1\hat{c}\left(\mathrm{kh}\right)\\ P(\mathrm{kh}-h)=a\left(\mathrm{kh}\right)P\left(\mathrm{kh}\right)a^T\left(\mathrm{kh}\right)+\mathrm{Rr}+L\left(\mathrm{kh}\right)[a\left(\mathrm{kh}\right)(P\left(\mathrm{kh}\right){\hat{c}}^T\left(\mathrm{kh}\right)\\ -M\left(\mathrm{kh}\right){\hat{c}}^T(\mathrm{kh}-h))+L1\left(\mathrm{kh}\right)\mathrm{Ro}]+L1\left(\mathrm{kh}\right){(L1\left(\mathrm{kh}\right)+2I)}^T\mathrm{Ro}\\ M(\mathrm{kh}+h)=a\left(\mathrm{kh}\right)P\left(\mathrm{kh}\right)+\mathrm{Rr}+L\left(\mathrm{kh}\right)(P\left(\mathrm{kh}\right){\hat{c}}^T\left(\mathrm{kh}\right)-M\left(\mathrm{kh}\right){\hat{c}}^T(\mathrm{kh}-h))\end{array}$

Wherein, cov (o (kh), o (kh) ^t)=Ro δ _ij, cov (r (kh), r (kh) ^t)=Rr δ _ij, wherein, cov () represents covariance computing, δ _ijfor kronecker-δ function, K1 is constant, and the state value obtaining submodel is further estimated:

$\begin{array}{c}x(\mathrm{kh}+h)=\hat{a}x\left(\mathrm{kh}\right)+\hat{b}u\left(\mathrm{kh}\right)+L\left(\mathrm{kh}\right)[y\left(\mathrm{kh}\right)-y(\mathrm{kh}-h)-(\hat{c}\left(\mathrm{kh}\right)-\hat{c}(\mathrm{kh}-h))x\left(\mathrm{kn}\right)\\ -(\hat{d}\left(\mathrm{kh}\right)-\hat{d}(\mathrm{kh}-h))u\left(\mathrm{kh}\right)]+L1\left(\mathrm{kh}\right)o\left(\mathrm{kh}\right)\end{array}$

The input and output state moment state estimation of original system is obtained: x (KT according to the parameter corresponding relation of two rate remote control system and submodel ₁)=x (Kph), x (KT ₂)=x (Kqh);

Step 3: pass rank parameter identification according to the state value that submodel dopes to two rate remote control system, expression is as follows:

Wherein, i, j represent the sampling instant of input and output signal respectively,

θ _x(iph)＝[A(iph),B(iph)] ^T，

θ _i(jqh)＝[C(jqh),D(jqh)] ^T，

X(iph+ph)＝x(iph+ph)，

A (iph), B (iph), C (jqh), D (jqh) represent the estimated value of A (iph), B (iph), C (jqh), D (jqh) respectively;

The parameter of next step submodel system is obtained again according to the result of passing rank identification:

$A^{\frac{1}{p}}\hat{b}\underset{i=0}{\overset{p-1}{\mathrm{\Σ}}}{A}^{\frac{i}{p}},\widehat{c=C,d=D;}$

Step 4: circulation step 2 and step 3 pass rank parameter identification process until complete all two rate remote control systems.

In described step 2, Ro and Rr is estimated by following formula:

$\mathrm{Rr}=\frac{1}{L}\underset{i=1}{\overset{L}{\mathrm{\Σ}}}r(\mathrm{kh}+\mathrm{ih})r^T(\mathrm{kh}+\mathrm{ih})$

$\mathrm{Ro}=\frac{1}{L}\underset{i=1}{\overset{L}{\mathrm{\Σ}}}o(\mathrm{kh}+\mathrm{ih})o^T(\mathrm{kh}+\mathrm{ih})$

Wherein, $r(\mathrm{kh}+\mathrm{ih})=x(\mathrm{kh}+(i+1)h)-a(\mathrm{kh}+\mathrm{ih})x(\mathrm{kh}+\mathrm{ih})-\hat{b}(\mathrm{kh}+\mathrm{ih})u(\mathrm{kh}+\mathrm{ih}),$ L represents the data length of data r (kh) or o (kh).

Beneficial effect

The present invention propose a kind of be applicable to remote control system controling parameters pass rank discrimination method, the remote operating process simplification of sampling different by input/output signal is a pair of rate sampling system, utilize the cycle submodel of method establishment one based on minimum sampling of resampling, and the intermediateness value of resampling system is estimated with this, the resampling cycle is greater than owing to exporting the sampling period, the output estimation difference of system can not directly revise intermediateness value, so the method that there is employed herein state change difference and difference is revised the state estimation in Kalman Prediction process, the present invention can carry out identification to the two rate remote control system parameters in the known and unknown two kinds of situations of state value, the feasibility of simulation results show the method and validity.

The invention has the beneficial effects as follows: be two rate sampling system model by remote operating process simplification, identified parameters is there is more in method for the identification again of conventional lift skill upgrading variable, the problems such as the large and identification precision of calculated amount is undesirable, the present invention can effective estimation procedure state value, and reduce calculated amount and computing time, compared to former method, calculation procedure is simple, is convenient to apply in engineering.Through simulating, verifying this method, there is good accuracy, robustness and validity.

Accompanying drawing explanation

What Fig. 1 was applicable to remote control system controling parameters passs rank discrimination method

The process flow diagram of Fig. 2 submodel state-updating

Fig. 3 Parameter identification process

Fig. 4 Parameter identification result

Embodiment

Now in conjunction with the embodiments, the invention will be further described for accompanying drawing:

As shown in Figure 1, it comprises following four steps to the process flow diagram of the embodiment of the present invention:

Step 1: remote operating input control signal and output feedback signal are gathered, input control signal is position and the torque signals of main side hand controller, output feedback signal comprises position, the speed and torque signals etc. of remote operating robotic arm joint, ignore the impact of system time delay and consider the frequency of the sample frequency of output feedback signal lower than input control signal, to being a pair of rate system model by remote operating process simplification after input/output signal process, its state-space model is:

$\left\{\begin{array}{c}x(K+1)=\mathrm{Ax}\left(K\right)+\mathrm{Bu}\left(K\right)+w\left(K\right)\\ y\left(K\right)=\mathrm{Cx}\left(K\right)+\mathrm{Du}\left(K\right)+v\left(K\right)\end{array}\right.---\left(1\right)$

Wherein, A, B, C, D are systematic parameter, and K is expression the K time sampling, and x (K) is state value, and u (K) is input signal, and y (K) is output signal, and the control inputs sampling period making system is T ₁with export sampling period be T ₂, h is T ₁and T ₂highest common factor, then T ₁and T ₂expression formula be T ₁=ph, T ₂=qh and p and q two number is relatively prime, w (K) and v (K) represents process noise and measurement noises respectively, and it is all stable zero-mean white noise sequence, meet: E (w (i))=0, E (v (i))=0cov (w (i), w (j) ^t)=Rw δ _ij, cov (v (i), v (j) ^t)=Rv δ _ij, cov () represents covariance computing, δ _ijfor kronecker-δ function.Set up and have mutually isostructural submodel with former pair of rate discrete system:

$\left\{\begin{array}{c}x(\mathrm{kh}+h)=\mathrm{ax}\left(\mathrm{kh}\right)+\mathrm{bu}\left(\mathrm{kh}\right)+r\left(\mathrm{kh}\right)\\ y\left(\mathrm{kh}\right)=\mathrm{cx}\left(\mathrm{kh}\right)+\mathrm{du}\left(\mathrm{kh}\right)+v\left(\mathrm{kh}\right)\end{array}\right.---\left(2\right)$

Wherein, a, b, c, d are submodel parameter, k is the kth time sampling of submodel, r (kh) and v (kh) is respectively submodel process noise and measurement noises, export sampling period T1 and be greater than input refresh cycle T2, meet p < q, equal with system (2) coefficient of correspondence under the framework cycle by system (1), submodel coefficient can be obtained:

$a=A^{\frac{1}{p}},b=B/\underset{i=0}{\overset{p-1}{\mathrm{\&Sigma;}}}{A}^{\frac{i}{p}},c=C,d=D,r\left(\mathrm{kh}\right)=w\left(\mathrm{kh}\right)/\underset{i=0}{\overset{p}{\mathrm{\&Sigma;}}}{A}^{\frac{i}{p}}.$

Step 2, submodel state estimation.Its process flow diagram is as shown in Figure 2:

First, calculate the residual values of submodel, the update cycle of residual values o (kh) is identical with the cycle qh of output variable y (kh), and its expression formula is:

$\left(\mathrm{kh}\right)=(y\left(\mathrm{kh}\right)-\hat{c}\left(\mathrm{kh}\right)x\left(\mathrm{kh}\right)-\hat{d}\left(\mathrm{kh}\right)u\left(\mathrm{kh}\right))\underset{i=0}{\overset{q}{\mathrm{\&Sigma;}}}{\hat{a}}^i---\left(3\right)$

Its update cycle is identical with the cycle qh of output variable y (kh), varivance matrix and error co-variance matrix expression formula are M (kh)=x (kh) x (kh-h), P (kh)=x (kh) x (kh), wherein x (kh)=x (kh)-x (kh), the estimated value that x (kh) is x (kh), order $\widehat{\mathrm{\&theta;}}\left(\mathrm{kh}\right)={[\hat{c}\left(\mathrm{kh}\right),d\left(\mathrm{kh}\right)]}^T,$ utilization pushes through recurrence equation:

Wherein, a (kh), d (kh) is respectively the estimated value of a (kh), b (kh), a (kh), d (kh), then the more new formula of Kalman filtering gain matrix L (kh), L1 (kh), P (kh) and M (kh) is:

$\begin{array}{c}L\left(\mathrm{kh}\right)=[a\left(\mathrm{kh}\right)(M\left(\mathrm{kh}\right){\hat{c}}^T(\mathrm{kh}-h)-P\left(\mathrm{kh}\right){\hat{c}}^T\left(\mathrm{kh}\right))+L1\left(\mathrm{kh}\right)\mathrm{Ro}][2\mathrm{Ro}+\hat{c}\left(\mathrm{kh}\right)P\left(\mathrm{kh}\right)\\ {\hat{c}}^T\left(\mathrm{kh}\right)-2\hat{c}\left(\mathrm{kh}\right)M\left(\mathrm{kh}\right){\hat{c}}^T(\mathrm{kh}-h)+\hat{c}(\mathrm{kh}-h)P(\mathrm{kh}-h){\hat{c}}^T(\mathrm{kh}-h){]}^{-1}\\ \end{array}---\left(5\right)$

$\begin{array}{c}L1\left(\mathrm{kh}\right)=K1\hat{c}\left(\mathrm{kh}\right)\end{array}----\left(6\right)$

$\begin{array}{c}P(\mathrm{kh}-h)=a\left(\mathrm{kh}\right)P\left(\mathrm{kh}\right)a^T\left(\mathrm{kh}\right)+\mathrm{Rr}+L\left(\mathrm{kh}\right)[a\left(\mathrm{kh}\right)(P\left(\mathrm{kh}\right){\hat{c}}^T\left(\mathrm{kh}\right)-M\left(\mathrm{kh}\right){\hat{c}}^T(\mathrm{kh}-h))\\ +L1\left(\mathrm{kh}\right)\mathrm{Ro}]+L1\left(\mathrm{kh}\right){(L1\left(\mathrm{kh}\right)+2I)}^T\mathrm{Ro}\\ \end{array}----\left(7\right)$

$\begin{array}{c}M(\mathrm{kh}+h)=a\left(\mathrm{kh}\right)P\left(\mathrm{kh}\right)+\mathrm{Rr}+L\left(\mathrm{kh}\right)(P\left(\mathrm{kh}\right){\hat{c}}^T\left(\mathrm{kh}\right)-M\left(\mathrm{kh}\right){\hat{c}}^T(\mathrm{kh}-h))\end{array}---\left(8\right)$

Wherein, cov (o (kh), o (kh) ^t)=Ro δ _ij, cov (r (kh), r (kh) ^t)=Rr δ _ij, wherein δ i _jfor kronecker-δ function, K1 is constant, and the state value obtaining submodel is further estimated:

$\begin{array}{c}x(\mathrm{kh}+h)=\hat{a}x\left(\mathrm{kh}\right)+\hat{b}u\left(\mathrm{kh}\right)+L\left(\mathrm{kh}\right)[y\left(\mathrm{kh}\right)-y(\mathrm{kh}-h)-(\hat{c}\left(\mathrm{kh}\right)-\hat{c}(\mathrm{kh}-h))x\left(\mathrm{kn}\right)\\ -(\hat{d}\left(\mathrm{kh}\right)-\hat{d}(\mathrm{kh}-h))u\left(\mathrm{kh}\right)]+L1\left(\mathrm{kh}\right)o\left(\mathrm{kh}\right)\end{array}---\left(9\right)$

The parameter corresponding relation of the two rate system of last basis and submodel obtains the input and output state moment state estimation of former pair of rate system: x (KT ₁)=x (Kph), x (KT ₂)=x (Kqh).

Step 3: two rate remote control system parameter pass rank identification.Adopt the method for passing rank identification to carry out identification for two rate remote control system parameter, detailed process is as follows:

Wherein, i, j represent the sampling instant of input and output signal respectively,

θ _x(iph)＝[A(iph),B(iph)] ^T，

θ _i(jqh)＝[C(jqh),D(jqh)] ^T，

X(iph+ph)＝x(iph+ph)，

A (iph), B (iph), C (jqh), D (jqh) represent the estimated value of A (iph), B (iph), C (jqh), D (jqh) respectively.And then the parameter of next step submodel system is obtained according to the result of passing rank identification:

$A^{\frac{1}{p}}\hat{b}\underset{i=0}{\overset{p-1}{\mathrm{\&Sigma;}}}{A}^{\frac{i}{p}},\widehat{c=C,d=D;}$

Step 4: the formula (3) in circulation step 2 and step 3 to formula (17) until complete two rate remote control system parameter pass rank identification process.

Ro and Rr is estimated by following formula in step 2:

$\mathrm{Rr}=\frac{1}{L}\underset{i=1}{\overset{L}{\mathrm{\&Sigma;}}}r(\mathrm{kh}+\mathrm{ih})r^T(\mathrm{kh}+\mathrm{ih})---\left(16\right)$

$\mathrm{Ro}=\frac{1}{L}\underset{i=1}{\overset{L}{\mathrm{\&Sigma;}}}o(\mathrm{kh}+\mathrm{ih})o^T(\mathrm{kh}+\mathrm{ih})---\left(17\right)$

Wherein, $r(\mathrm{kh}+\mathrm{ih})=x(\mathrm{kh}+(i+1)h)-a(\mathrm{kh}+\mathrm{ih})x(\mathrm{kh}+\mathrm{ih})-\hat{b}(\mathrm{kh}+\mathrm{ih})u(\mathrm{kh}+\mathrm{ih}),$ L represents the data length of data r (kh) or o (kh).

For method of the present invention, carry out simulating, verifying, the state space equation of the two rate remote control system of order is as follows:

$\left\{\begin{array}{c}x(k+1)=0.9x\left(k\right)+1.2u\left(k\right)+w\left(k\right)\\ y\left(k\right)=1.1x\left(k\right)-0.5u\left(k\right)+v\left(k\right)\end{array}\right.---\left(18\right)$

Wherein, input sample cycle T ₁=3h, the output sampling period is T ₂=4h, h=1s, framework period p qh=12h, input signal is quasi stationary sequence E [u (k)]=1.5, and system noise and measurement noises are respectively Rw=0.05 ², Rv=0.05 ², system (18) submodel data length to be identified is 10000.Carry out identification to the parameter of system, as shown in Figure 3, Figure 4, simulation result shows for identification process and result, and the present invention effectively can carry out identification to systematic parameter under the known and unknown situation of state value, and speed of convergence is very fast, and calculated amount is little, and identification result is accurate.

Claims (2)

1. what be applicable to remote control system controling parameters passs a rank discrimination method, it is characterized in that step is as follows:

Step 1: gather remote operating input control signal and output feedback signal, described input control signal is position and the torque signals of main side hand controller; Described output feedback signal comprises the position of remote operating robotic arm joint, speed and torque signals;

Be a pair of rate system model by remote operating process simplification, its state-space model is:

$\left\{\begin{array}{c}x(K+1)=Ax\left(K\right)+Bu\left(K\right)+w\left(K\right)\\ y\left(K\right)=Cx\left(K\right)+Du\left(K\right)+v\left(K\right)\end{array}\right.;$

Wherein, A, B, C, D are systematic parameter, and K represents system the K time sampling, and x (K) is state value, and u (K) is input signal, and y (K) is output signal;

The control inputs sampling period of system is made to be T ₁with export sampling period be T ₂, h is T ₁and T ₂highest common factor, T ₁and T ₂meet: T ₁=ph, T ₂=qh and p and q two number is relatively prime, w (K) and v (K) represents process noise and measurement noises respectively, and the system model of system is: $\left\{\begin{array}{c}x(kh+h)=ax\left(kh\right)+bu\left(kh\right)+r\left(kh\right)\\ y\left(kh\right)=cx\left(kh\right)+du\left(kh\right)+v\left(kh\right)\end{array}\right.,$ Wherein, a, b, c, d are submodel parameter, and k is the kth time sampling of submodel, and r (kh) and v (kh) is respectively submodel process noise and measurement noises;

Parameters relationship between two models is: $a=A^{\frac{1}{p}},$ $b=B/\underset{i=0}{\overset{p-1}{\&Sigma;}}{A}^{\frac{i}{p}},$ c＝C，d＝D， $r\left(kh\right)=w\left(kh\right)/\underset{i=0}{\overset{p}{\&Sigma;}}{A}^{\frac{i}{p}};$

Step 2: make residual values o (kh) its expression formula be:

$o\left(kh\right)=\left(y\right(kh)-\hat{c}(kh\left)x\right(kh)-\hat{d}(kh\left)u\right(kh\left)\right)/\underset{i=0}{\overset{q}{\&Sigma;}}{\hat{a}}^i$

Its update cycle is identical with the cycle qh of output variable y (kh), varivance matrix and error co-variance matrix expression formula are M (kh)=x (kh) x (kh-h), P (kh)=x (kh) x (kh), wherein estimated value, order $\widehat{\mathrm{\&theta;}}\left(kh\right)={\&lsqb;\hat{c}\left(kh\right),d\left(kh\right)\&rsqb;}^T,$ recurrence equation is utilized to obtain:

Wherein, be respectively the estimated value of a (kh), b (kh), a (kh), d (kh), then the more new formula of Kalman filtering gain matrix L (kh), L1 (kh), P (kh) and M (kh) is:

$L\left(kh\right)=\&lsqb;a\left(kh\right)\left(M\right(kh\left){\hat{c}}^T\right(kh-h)-P(kh\left){\hat{c}}^T\right(kh\left)\right)+L1\left(kh\right)Ro\&rsqb;\&lsqb;2Ro+\hat{c}\left(kh\right)P\left(kh\right)$

${\hat{c}}^T\left(kh\right)-2\hat{c}\left(kh\right)M\left(kh\right){\hat{c}}^T(kh-h)+\hat{c}(kh-h)P(kh-h){\hat{c}}^T(kh-h){\&rsqb;}^{-1}$

$L1\left(kh\right)=K1\hat{c}\left(kh\right)$

$\begin{array}{c}P(kh+h)=a\left(kh\right)P\left(kh\right)a^T\left(kh\right)+Rr+L\left(kh\right)\&lsqb;a\left(kh\right)\left(P\right(kh\left){\hat{c}}^T\right(kh)\\ -M\left(kh\right){\hat{c}}^T(kh-h))+L1\left(kh\right)Ro\&rsqb;+L1\left(kh\right){\left(L1\right(kh)+2I)}^{T}Ro\end{array}$

$M(kh+h)=a\left(kh\right)P\left(kh\right)+Rr+L\left(kh\right)\left(P\right(kh\left){\hat{c}}^T\right(kh)-M(kh\left){\hat{c}}^T\right(kh-h\left)\right)$

Wherein, cov (o (kh), o (kh) ^t)=Ro δ _ij, cov (r (kh), r (kh) ^t)=Rr δ _ij, wherein, cov () represents covariance computing, δ _ijfor kronecker-δ function, K1 is constant, and the state value obtaining submodel is further estimated:

$\begin{array}{c}x(kh+h)=\hat{a}x\left(kh\right)+\hat{b}u\left(kh\right)+L\left(kh\right)\&lsqb;y\left(kh\right)-y(kh-h)-\left(\hat{c}\right(kh)-\hat{c}(kh-h\left)\right)x\left(kh\right)\\ -\left(\hat{d}\right(kh)-\hat{d}(kh-h\left)\right)u\left(kh\right)\&rsqb;+L1\left(kh\right)o\left(kh\right)\end{array}$

The input and output state moment state estimation of original system is obtained: x (KT according to the parameter corresponding relation of two rate remote control system and submodel ₁)=x (Kph), x (KT ₂)=x (Kqh);

Step 3: pass rank parameter identification according to the state value that submodel dopes to two rate remote control system, expression is as follows:

Wherein, i, j represent the sampling instant of input and output signal respectively,

θ _x(iph)＝[A(iph),B(iph)] ^T，

θ _i(jqh)＝[C(jqh),D(jqh)] ^T，

X(iph+ph)＝x(iph+ph)，

represent the estimated value of A (iph), B (iph), C (jqh), D (jqh) respectively;

The parameter of next step submodel system is obtained again according to the result of passing rank identification:

$a=A^{\frac{1}{p}},$ $\hat{b}=B/\underset{i=0}{\overset{p-1}{\&Sigma;}}{A}^{\frac{i}{p}},$ $\hat{c}=C,$ d＝D；

Step 4: circulation step 2 and step 3 pass rank parameter identification process until complete all two rate remote control systems.

2. what be applicable to remote control system controling parameters according to claim 1 passs rank discrimination method, it is characterized in that:

In described step 2, Ro and Rr is estimated by following formula:

$Rr=\frac{1}{L}\underset{i=1}{\overset{L}{\&Sigma;}}r(kh+ih)r^T(kh+ih)$

$Ro=\frac{1}{L}\underset{i=1}{\overset{L}{\&Sigma;}}o(kh+ih)o^T(kh+ih)$

Wherein, $r(kh+ih)=x(kh+(i+1\left)h\right)-a(kh+ih)x(kh+ih)-\hat{b}(kh+ih)u(kh+ih),$

L represents the data length of data r (kh) or o (kh).