CN105843040A

CN105843040A - Method and device for identifying unmanned helicopter kinetic parameters

Info

Publication number: CN105843040A
Application number: CN201610214448.2A
Authority: CN
Inventors: 蒋本忠; 白勍; 杨峥; 李明; 沙俊汀
Original assignee: Shenyang Shangbo Zhituo Technology Co Ltd
Current assignee: Shenyang Shangbo Zhituo Technology Co Ltd
Priority date: 2016-04-06
Filing date: 2016-04-06
Publication date: 2016-08-10
Anticipated expiration: 2036-04-06
Also published as: CN105843040B

Abstract

The invention discloses a method and a device for identifying unmanned helicopter kinetic parameters. The method comprises steps that initial estimates of to-be-identified parameters are acquired, a first system matrix and a first control matrix of the initial estimates containing the to-be-identified parameters are acquired, a response function of the system measurement data is acquired according to the acquired data, a system matrix containing the optimal solution and a coefficient matrix are acquired by utilizing a state space mode, the response function and a cost function according to the flight data, the first system matrix and the first control matrix, in combination with the update strategy, a second system matrix and a second control matrix after update are acquired; whether the second system matrix and the second control matrix satisfy the preset condition is determined, if not, the second system matrix and the second control matrix are taken as the first system matrix and the first control matrix to re-actuate the previous steps till the preset condition is satisfied, and present values of to-be-identified parameters in the second system matrix and the second control matrix are taken as identification accomplishment target values.

Description

The discrimination method of depopulated helicopter kinetic parameter and device

Technical field

It relates to unmanned machine control field, particularly relate to the discrimination method of a kind of depopulated helicopter kinetic parameter And device.

Background technology

Depopulated helicopter is estimated and controls based on accurate kinetics reference model.Unmanned straight from the end of the eighties Rise machine concept and propose until in the research of nearly 30 years of today, although propose substantial amounts of control method, and solved outstanding Stop with low-speed operations under the conditions of calm and tracking problem, but the potential maneuverability of depopulated helicopter is not opened completely Sending out, far apart with the maneuverability that flight control personnel are obtained by artificial remote control, therefore, from principal mode, nobody goes straight up at present Machine flight performance cannot meet it in civilian and military aspect application demand.

Advanced Flight Control Algorithm will be based on rational kinetic model, and for depopulated helicopter, by Complicated in its kinetics, it is the non-linear strongly coupled system of quiet instability, accurate kinetic model cannot obtain, system dynamic Learn Model Distinguish complicated with reference model mismatch problems be cause existing advanced control method cannot the major reason of actual application, Also it is the of paramount importance restraining factors that cannot promote of depopulated helicopter maneuverability.

Unmanned helicopter system modeling needs to carry out the reverse matching of kinetic model based on state of flight and control data, It it is sufficiently complex process.Although prior art exists, the model structure of depopulated helicopter is divided into rotor aerodynamics and body Kinetics two parts, and pass through Analysis on Mechanism with model structure it is assumed that establish parameterized linear model structure, then by single Enter the frequency domain optimal estimation method singly gone out and based on system, the response data that frequency sweep inputs is carried out transmitting the parameter optimum of function Solve.But, frequency domain optimal estimation method needs the selection of initial parameter, if initial parameter select it cannot be guaranteed that model prediction without Partially, then can cause directly dissipating of frequency domain optimal estimation method, but the discrimination method of traditional frequency domain is appropriate only for singly entering list Go out system, need transfer matrix result of calculation is carried out the substep identification of single-input single-output, the most in actual applications, for unmanned Helicopter System needs the substantial amounts of empirical trial carried out repeatedly, could select initial parameter and iterations, and amount of calculation is huge Greatly, identification efficiency is low.

Summary of the invention

The disclosure provides discrimination method and the device of a kind of depopulated helicopter kinetic parameter, for solving traditional frequency domain Discrimination method to be appropriate only for the amount of calculation that single-input single-output system causes huge, the inefficient problem of identification.

To achieve these goals, the disclosure provides the discrimination method of a kind of depopulated helicopter kinetic parameter, described side Method includes:

A. obtain the initial estimate of parameter to be identified, obtain comprising the first of the initial estimate of described parameter to be identified Sytem matrix and first controls matrix；

B. control matrix according to described the first system matrix and described first, obtain the receptance function of system measurement data, And determine the characteristic frequency point of described receptance function；

C. according to the flying quality of the described depopulated helicopter gathered, described the first system matrix and described first controls square State-space model, described receptance function and default cost function that battle array, utilization are preset obtain described parameter to be identified Optimal solution, obtains comprising sytem matrix and the coefficient matrix of the optimal solution of described parameter to be identified；

D. comprise sytem matrix and the coefficient matrix of the optimal solution of described parameter to be identified described in basis, utilize and preset more New Policy obtains the second system matrix after updating and second and controls matrix；

E. judge that described second system matrix and described second controls whether matrix meets pre-conditioned, described second be System matrix and described second control matrix be unsatisfactory for described pre-conditioned time, by described second system matrix and described second control Matrix controls matrix as described the first system matrix and described first and again performs step b～e, until described second system square Battle array and described second control matrix meet described pre-conditioned till, by described second system matrix and described second system matrix In the currency of parameter described to be identified that comprises as the desired value completing identification.

Optionally, described step c includes:

C1. gather the flying quality of described depopulated helicopter, and according to the described flying quality gathered, utilize described response Function obtains actual spectrum response；

C2. described the first system matrix and described first are controlled matrix as in described state-space model being System matrix and control matrix calculus spectral response estimated value；

C3. respond according to described spectral response estimated value and described actual spectrum and obtain error of spectrum；

C4. according to described error of spectrum, described cost function is utilized to obtain the epicycle iteration of described parameter to be identified Excellent solution；

C5. described the first system matrix and described first is updated according to the optimal solution of the epicycle iteration of described parameter to be identified Control matrix, and will update after the first system matrix and update after first control matrix as described state space mould Sytem matrix in type and control matrix, obtain the state-space model that epicycle iteration obtains；

C6. judge that the spectral response that the state-space model utilizing epicycle iteration to obtain calculates is estimated by time domain checking Whether value is more than error threshold with the error of described actual spectrum response；

When the spectral response estimated value that calculates of state-space model judging to utilize epicycle iteration to obtain and described reality When the error of spectral response is not more than described error threshold, using the optimal solution of described epicycle iteration as described parameter to be identified Optimal solution, obtains comprising sytem matrix and the coefficient matrix of the optimal solution of described parameter to be identified；When judging to utilize epicycle iteration The spectral response estimated value that the state-space model obtained calculates is more than described error with the error of described actual spectrum response During threshold value, update described the first system matrix and described first and control matrix and again perform step c2～c5, until utilizing this The spectral response estimated value that the state-space model that wheel iteration obtains calculates is not more than with the error of described actual spectrum response Till described error threshold.

Optionally, described flying quality includes the longitudinal velocity of described depopulated helicopter, lateral velocity, vertical velocity, indulges To acceleration, transverse acceleration, vertical acceleration, roll angle, the angle of pitch, course angle, roll angle speed, pitch rate with And course angle speed.

Optionally, before described step c2, also include:

Obtain based on described longitudinal velocity, described lateral velocity, described vertical velocity, described roll angle, the described angle of pitch, Described course angle, described roll angle speed, described pitch rate and the first differential model of described course angle speed；

According to described longitudinal velocity, described lateral velocity, described vertical velocity, described roll angle speed, the described angle of pitch Speed and the first differential model of described course angle speed, obtain respectively described longitudinal velocity, described lateral velocity, described hang down To speed, described roll angle speed, described pitch rate and described course angle speed with system mode and controlled quentity controlled variable phase Close linearly thinks incremental model in absolute terms；

Obtain the main rotor of described depopulated helicopter and winglet waves equation model；

By by described wave equation model respectively with described longitudinal velocity linearly think in absolute terms incremental model, described laterally Speed linearly think in absolute terms incremental model, described roll angle speed linearly think incremental model, described pitch rate in absolute terms Linear absolutization incremental model carries out coupling and obtains longitudinal acceleration model, transverse acceleration model, and roll angle acceleration Model, angle of pitch Fast track surgery；

Obtain the vertical passage of described depopulated helicopter and the kinetic model of course passage；

By the kinetic model of described vertical passage and course passage respectively with described vertical velocity linearly think increasing in absolute terms The linear absolutization incremental model of amount model and described course angle speed carries out coupling and obtains vertical acceleration model, course angle Fast track surgery；

According to described longitudinal acceleration model, transverse acceleration model, and roll angle Fast track surgery, the angle of pitch accelerate Degree model, described vertical acceleration model, course angle Fast track surgery and the control feedback quantity model of course-stability gyroscope Determine state equation matrix；

Determine described state-space model according to described state equation matrix, the sytem matrix of described state-space model and Control matrix and include that described parameter to be identified, described parameter to be identified are according to described longitudinal acceleration model, transverse acceleration Model, and roll angle Fast track surgery, angle of pitch Fast track surgery, described vertical acceleration model, course angle acceleration mould The parameter determination to be identified controlled in feedback quantity model of type and course-stability gyroscope.

Optionally, described based on described longitudinal velocity, described lateral velocity, described vertical velocity, described roll angle, described The angle of pitch, described course angle, described roll angle speed, described pitch rate and the first differential mould of described course angle speed Type includes:

Wherein, u, v, w represent described longitudinal velocity, described lateral velocity and described vertical velocity respectively, and p, q, r are respectively Represent described roll angle speed, described pitch rate and described course angle speed,θ, ψ represent described roll angle, institute respectively State the angle of pitch and described course angle, I_xx, I_yy, I_zzRepresenting the rotary inertia of X-axis, Y-axis and Z axis respectively, m represents described unmanned straight The quality of the machine of liter, L, M, N are three axle moments of torsionIn parameter, X, Y, Z are for acting on described depopulated helicopter The force vector of center of gravityIn parameter；

Described longitudinal velocity, described lateral velocity, described vertical velocity, described roll angle speed, described pitch rate And the linear absolutization incremental model relevant to system mode and controlled quentity controlled variable of described course angle speed includes:

Wherein, u₀, v₀, w₀, p₀, q₀, r₀Representing u respectively, the original state of v, w, p, q, r, δ represents relative to original state Absolute increment；

The main rotor of described depopulated helicopter and the equation model of waving of winglet include:

Wherein, τ_f, A_b, B_b, A_lon, B_lat, K_cFor parameter to be identified, a, b represent that the longitudinal direction of described main rotor waves angle respectively Laterally waving the first harmonic component at angle, c, d represent that the longitudinal direction of described winglet is waved angle and laterally waves angle once respectively Harmonic component, δ_latRepresent crosswise joint input quantity, δ_lonRepresent longitudinally controlled input quantity；

Described longitudinal acceleration model and described transverse acceleration model include:

Wherein, X_u, X θ, X_a, Y_v,Y_bFor parameter to be identified；

Described roll angle Fast track surgery and described angle of pitch Fast track surgery include:

Wherein, L_u, L_v, L_b, L_w, M_u, M_v, M_a, M_wFor parameter to be identified；

The vertical passage of described depopulated helicopter includes with the kinetic model of course passage:

Wherein, N_v, N_r, N_ped, K_rfb, K_rFor parameter to be identified, δ_pedRepresent Heading control input quantity, δ_colRepresent vertical control Input quantity processed；

Described vertical acceleration model, course angle Fast track surgery include:

Wherein, Z_a, Z_b, Z_w, Z_r, Z_colFor parameter to be identified.

Optionally, described state-space model includes:

Wherein, Θ represents unknown parameter to be identified, and X is the system mode vector of n × 1 dimension, and U represents the control that r × 1 is tieed up Moment matrix, Y represents that system output is maintained in p × 1, and A is that n × n ties up sytem matrix, and B represents that n × r dimension controls matrix, and C represents that p × n ties up Calculation matrix；

Wherein, U includes: U=(δ_lat δ_lon δ_ped δ_col)^T, δ_latRepresent crosswise joint input quantity, δ_lonRepresent longitudinally control Input quantity processed, δ_pedRepresent Heading control input quantity, δ_colRepresent that vertical control inputs；Parameter in A and B is described ginseng to be identified Number.

The disclosure also provides for the device for identifying of a kind of depopulated helicopter kinetic parameter, and described device includes:

Initial estimation module, obtains the initial estimate of parameter to be identified for performing a., obtains comprising described to be identified The first system matrix of the initial estimate of parameter and first controls matrix；

RESPONSE CALCULATION module, is used for performing b. and controls matrix according to described the first system matrix and described first, obtain system The receptance function of system measurement data, and determine the characteristic frequency point of described receptance function；

Iterative processing module, for performing the c. flying quality according to the described depopulated helicopter gathered, described first is System matrix and described first controls matrix, utilizes state-space model, described receptance function and the default cost letter preset Number obtains the optimal solution of described parameter to be identified, obtains comprising sytem matrix and the coefficient square of the optimal solution of described parameter to be identified Battle array；

More new module, for perform d. according to described in comprise described parameter to be identified the sytem matrix of optimal solution and coefficient Matrix, utilizes the more New Policy preset to obtain the second system matrix after updating and second and controls matrix；

For e., judge module, judges that described second system matrix and described second controls whether matrix meets default bar Part, described second system matrix and described second control matrix be unsatisfactory for described pre-conditioned time, by described second system square Battle array and described second controls matrix and again performs step b～e, directly as described the first system matrix and described first control matrix To described second system matrix and described second control matrix meet described pre-conditioned, by described second system matrix and The currency of the parameter described to be identified comprised in described second system matrix is as the desired value completing identification.

Optionally, described iterative processing module, including:

Data acquisition and RESPONSE CALCULATION submodule, gather the flying quality of described depopulated helicopter, and root for performing c1. According to the described flying quality gathered, described receptance function is utilized to obtain actual spectrum response；

Submodule is estimated in response, is used for performing c2. and described the first system matrix and described first control matrix is made respectively For the sytem matrix in described state-space model and control matrix calculus spectral response estimated value；

Error obtains submodule, is used for performing c3. and obtains according to described spectral response estimated value and the response of described actual spectrum Take error of spectrum；

Calculating sub module, is used for performing c4. according to described error of spectrum, utilizes described cost function to obtain described to be identified The optimal solution of the epicycle iteration of parameter；

Model modification submodule, updates institute for performing c5. according to the optimal solution of the epicycle iteration of described parameter to be identified State the first system matrix and described first control matrix, and will update after the first system matrix and update after first control square Battle array is respectively as the sytem matrix in described state-space model and controls matrix, obtains the state space mould that epicycle iteration obtains Type；

Checking submodule, the state-space model utilizing epicycle iteration to obtain for performing c6. to be judged by time domain checking Whether the spectral response estimated value calculated is more than error threshold with the error of described actual spectrum response；

When the spectral response estimated value that calculates of state-space model judging to utilize epicycle iteration to obtain and described reality When the error of spectral response is not more than described error threshold, using the optimal solution of described epicycle iteration as described parameter to be identified Optimal solution, obtain comprising sytem matrix and the coefficient matrix of the optimal solution of described parameter to be identified；When judging to utilize epicycle repeatedly The spectral response estimated value that the state-space model that generation obtains calculates is more than described mistake with the error of described actual spectrum response During difference limen value, update described the first system matrix and described first and control matrix and again perform step c2～c5, until utilizing The spectral response estimated value that the state-space model that epicycle iteration obtains calculates is little with the error that described actual spectrum responds Till described error threshold.

Optionally, described flying quality includes the longitudinal velocity of described depopulated helicopter, lateral velocity, vertical velocity, indulges To acceleration, transverse acceleration, vertical acceleration, roll angle, the angle of pitch, course angle, roll angle speed, pitch rate and Course angle speed.

Optionally, also include: model obtains submodule, for before described step c2:

By the kinetic model of described vertical passage and course passage respectively with described vertical velocity linearly think increasing in absolute terms The linear absolutization incremental model of amount model and described course angle speed carries out coupling and obtains vertical acceleration model, course Angular acceleration model；

Wherein, X_u, X θ, X_a, Y_v,Y_bFor parameter to be identified；

Described vertical acceleration model, course angle Fast track surgery include:

Wherein, Z_a, Z_b, Z_w, Z_r, Z_colFor parameter to be identified.

Optionally, described state-space model includes:

In sum, the discrimination method of the depopulated helicopter kinetic parameter that the disclosure is provided and device, by bag Sytem matrix containing multiple parameter state spatial models to be identified and control matrix and carry out identification, can the simultaneously multiple ginseng of identification Number, it is achieved that multiple-input, multiple-output, it is possible to reduce amount of calculation, improves the identification efficiency of parameter.

Other feature and advantage of the disclosure will be described in detail in detailed description of the invention part subsequently.

Accompanying drawing explanation

Accompanying drawing is used to provide further understanding of the disclosure, and constitutes the part of description, with following tool Body embodiment is used for explaining the disclosure together, but is not intended that restriction of this disclosure.In the accompanying drawings:

Fig. 1 is the flow chart of the discrimination method of a kind of depopulated helicopter kinetic parameter that the disclosure one embodiment provides；

Fig. 2 is the flow chart of the discrimination method of a kind of depopulated helicopter kinetic parameter that the disclosure one embodiment provides；

Fig. 3 is the schematic diagram of the body coordinate system of the depopulated helicopter that the disclosure one embodiment provides；

Fig. 4 is the structural representation of the main rotor system of the depopulated helicopter that the disclosure one embodiment provides；

Fig. 5 is the main rotor relative velocity schematic diagram that the disclosure one embodiment provides；

Fig. 6 is the mechanical construction drawing of the stable winglet of the depopulated helicopter that the disclosure one embodiment provides；

Fig. 7 is the structural frames of the device for identifying of a kind of depopulated helicopter kinetic parameter that the disclosure one embodiment provides Figure；

Fig. 8 is the structural frames of the device for identifying of the another kind of depopulated helicopter kinetic parameter that the disclosure one embodiment provides Figure.

Detailed description of the invention

It is described in detail below in conjunction with accompanying drawing detailed description of the invention of this disclosure.It should be appreciated that this place is retouched The detailed description of the invention stated is merely to illustrate and explains the disclosure, is not limited to the disclosure.

Fig. 1 is the flow chart of the discrimination method of a kind of depopulated helicopter kinetic parameter that the disclosure one embodiment provides, Seeing Fig. 1, the method includes:

Step 101, obtains the initial estimate of parameter to be identified, obtains comprising the initial estimate of described parameter to be identified The first system matrix and first control matrix；

Step 102, controls matrix according to described the first system matrix and described first, obtains the response of system measurement data Function, and determine the characteristic frequency point of described receptance function；

Step 103, according to the flying quality of the described depopulated helicopter gathered, described the first system matrix and described first Control matrix, utilize the state-space model preset, described receptance function and the acquisition of default cost function described to be identified The optimal solution of parameter, obtains comprising sytem matrix and the coefficient matrix of the optimal solution of described parameter to be identified；

Step 104, according to sytem matrix and the coefficient matrix of the described optimal solution comprising described parameter to be identified, utilizes pre- If more New Policy obtain update after second system matrix and second control matrix；

Step 105, it is judged that described second system matrix and described second controls whether matrix meets pre-conditioned, described Second system matrix and described second control matrix be unsatisfactory for described pre-conditioned time, by described second system matrix and described the Two control matrix performs step 102～105 again as described the first system matrix and described first control matrix, until described Second system matrix and described second control matrix meet described pre-conditioned till, by described second system matrix and described the The currency of the parameter described to be identified comprised in two system matrix is as the desired value completing identification.It is thus possible to determine described The final structure of state-space model, and then the state-space model that can be determined by is to control depopulated helicopter, further The Control System Design that can also be applied to depopulated helicopter.

Fig. 2 is the flow chart of the discrimination method of a kind of depopulated helicopter kinetic parameter that the disclosure one embodiment provides, Seeing Fig. 2, this step step 103 may include that

Step 1031, gather the flying quality of described depopulated helicopter, and according to the described flying quality gathered, utilize institute State receptance function and obtain actual spectrum response；

Step 1032, using described the first system matrix and described first control matrix as described state-space model In sytem matrix and control matrix calculus spectral response estimated value；

Step 1033, respond according to described spectral response estimated value and described actual spectrum and obtain error of spectrum；

Step 1034, according to described error of spectrum, utilize described cost function to obtain the epicycle of described parameter to be identified repeatedly The optimal solution in generation；

Step 1035, optimal solution according to the epicycle iteration of described parameter to be identified update described the first system matrix and institute State the first control matrix, and the first system matrix after updating and the first control matrix after updating are as described state Sytem matrix in spatial model and control matrix, obtain the state-space model that epicycle iteration obtains；

Step 1036, judge that the frequency spectrum that calculates of state-space model utilizing epicycle iteration to obtain is loud by time domain checking Answer the estimated value error with the response of described actual spectrum whether more than error threshold；

When the spectral response estimated value that calculates of state-space model judging to utilize epicycle iteration to obtain and described reality When the error of spectral response is not more than described error threshold, using the optimal solution of described epicycle iteration as described parameter to be identified Optimal solution, obtains comprising sytem matrix and the coefficient matrix of the optimal solution of described parameter to be identified；When judging to utilize epicycle iteration The spectral response estimated value that the state-space model obtained calculates is more than described error with the error of described actual spectrum response During threshold value, update described the first system matrix and described first and control matrix and again perform step 1032～1035, until sharp The error that the spectral response estimated value that the state-space model obtained by epicycle iteration calculates and described actual spectrum respond is not Till described error threshold.

Wherein, before performing above-mentioned steps 1032～1035, described default state-space model first should be obtained, Its process can be as follows:

Firstly, for can the abstract depopulated helicopter body for rigid body, can be according to the depopulated helicopter shown in Fig. 3 Body coordinate system, sets up newton--and Eulerian equation is as follows:

Wherein, m represents the quality of depopulated helicopter,Represent unmanned respectively Three-dimensional line speed under helicopter body coordinate system and angle of revolution speed, wherein u, v, w represent longitudinally (also referred to as forward direction) speed respectively Degree, laterally (the most lateral) speed and vertical velocity, p, q, r represent roll angle speed, pitch rate and course angle respectively Speed, and Then it is respectively the force vector acting on center of gravity under body coordinate system With three axle torque vector, this driving force and moment are both from rotor and stablize winglet and jointly act on.Nobody shown in Fig. 3 goes straight up to In the body coordinate system of machine, X-axis is heading, and Y-axis points to body dextrad, and it is downward that Z axis is perpendicular to XOY plane, and initial point O is nothing People's helicopter center of gravity.

The differential equation of first order of Eulerian angles and body angular speed is can get by coordinate system rotation relationship:

Wherein,θ, ψ represent roll angle, the angle of pitch and course angle respectively.So can set up nine measurable variables (described Longitudinal velocity, described lateral velocity, described vertical velocity, described roll angle, the described angle of pitch, described course angle, described rolling Angular speed, described pitch rate and described course angle speed) first differential model, as follows:

Wherein, I_xx, I_yy, I_zzRepresenting the rotary inertia of X-axis, Y-axis and Z axis respectively, L, M, N are three axle moments of torsionIn parameter, X, Y, Z are the force vector acting on described depopulated helicopter center of gravity In parameter.

Then, need above-mentioned dynamic mechanical equation is carried out linearization process, concrete: so that the complexity set up Nonlinear model is easy to follow-up model parameter estimation and control, needs to be operated model the linearisation near a little.Here Choose hovering mode, carry out linearisation.For hovering linearisation original state (X₀,U₀)=(0,0), i.e. linearisation point is zero point (coordinate origin).So, to modelDo linearisation and obtain following form:

Wherein,Γ represents location parameter.Wherein, A represents that sytem matrix, B represent Controlling matrix, U represents controlled quentity controlled variable.

So, its linearisation can be obtained for the kinetic part (the first six differential equation) in foregoing model (3) absolute Delta state spatial model, as follows:

Wherein, u₀, v₀, w₀, p₀, q₀, r₀Representing u respectively, the original state of v, w, p, q, r, δ represents relative to original state Absolute increment, such as: δ u=u-u₀, δ v=v-v₀Etc., by that analogy.It is assumed herein that external force suffered by helicopter body and Moment is the continuous function of controlled quentity controlled variable U and system mode X, makes the controlled quentity controlled variable be:

U=(δ_lon δ_lat δ_col δ_ped)^T

Wherein, δ_latRepresent crosswise joint input quantity, δ_lonRepresent longitudinally controlled input quantity, δ_pedRepresent Heading control input Amount, δ_colRepresent that vertical control inputs, then can obtain the expansion of power as follows and moment:

It should be noted that in above-mentioned obtained formula (6), not all derivative term all can work, Er Qiegong In formula (6), all of rotor power is all that instantaneous entrance body is dynamic (dynamical) with moment, and its size is directly proportional to controlled quentity controlled variable, so The most sufficiently consider the kinetic effect of main rotor.Launch so power should be given up in follow-up handling process with moment The derivative term that in formula, amplitude is less, the characteristic of rotor aerodynamics system to be taken into full account itself, show controlled quentity controlled variable and pass through After rotor system response, the most just act on body movement.

It is illustrated in figure 4 the structural representation (Three Degree Of Freedom sketch) of the main rotor system of depopulated helicopter, wherein Ω, β, ξ, Θ represent the speed of mainshaft respectively, wave angle, angle of oscillation and total square angle.

As shown above, operator can drive steering wheel to change the meansigma methods of total square angle Θ by total square, the most permissible By changing longitudinally and horizontal controlled quentity controlled variable, make total square angle periodic variation, i.e. cyclic pitch.As such, it is possible to use expression formula Total square angle variation relation with position of rotation Ψ is described:

Θ (Ψ)=Θ₀-B_latδ_latcosΨ-A_lonδ_lonsinΨ (7)

Stress in view of rotor is directly related in the speed of air corresponding thereto, sets rotor here relative to air Speed is U, U_TAnd U_PBe respectively the speed of air be U to rotor tangential and vertically-oriented component, its relation is concrete such as Fig. 5 Shown in shown main rotor relative velocity schematic diagram.And the relation shown in Fig. 5 can be passed through formula (8) and represent:

Wherein, what β represented wing currently waves angle, U_∞Represent the free-stream velocity of air, α_DRepresent that air freely flows Dynamic speed and the angle of propeller hub plane.Thus the expression formula that can obtain rotor and speed U that is involutory of air is:

The expression formula being apparent from air force influent stream angle Φ according to relation shown in Fig. 5 is:

Φ=arctan (U_P/U_T) (10)

According to relation shown in Fig. 5, and infinitesimal lift knowable to linear aerodynamic theory is that expression formula is:

Wherein, ρ represents atmospheric density, C_lαRepresenting that wing rises curvature, c represents that rotor chord length, α=Θ-Φ represent wing The angle of attack.The least in view of Φ, then it is believed that lift infinitesimal has an approximate expression:

dF_Z=dL (12)

Wherein, F_ZRepresent lift.

Furthermore, it is contemplated that the feature of small-sized heligyro own, can make the following assumptions:

1) paddle is rigidity, and it is the least to wave angle beta, and ratio μ=U/ Ω R≤0.3；

2) ignore the inlet loss of blade tip end, and have U_P/U_T＜＜ 1；

3) only consider that rotor is waved impact by angular acceleration, angular velocity and linear acceleration；

4) the lead-lag corner ξ impact on rotor stress is ignored.

The most then there is air force dF_aero=dF_z, centrifugal force dF_cent=mdy Ω²Y β, inertia forceWave torque reaction M_k=k_ββ.So for waving axle, momental equation can be obtained as follows:

If definitionThen above formula (13) can turn to Following expression:

It should be noted that the main rotor for unmanned helicopter system, the most do not wave restriction, i.e. k_β=0 and w_β =Ω, now to wave angular frequency identical with anglec of rotation frequency, and according to the form of above-mentioned equation (14) it can be seen that wave angle Kinetics equation be a second-order system meeting pouring weight spring-damper pattern (Occur in dF_zIn).So may be used Should be a continuous function with 2 π as cycle with solution β (Ψ) finding out the above-mentioned differential equation, it is possible to use fourier level The form of number represents, shown in formula specific as follows:

β (Ψ)=β₀-β_1ccosΨ-β_1ssinΨ-β_2ccos2Ψ-β_2ssin2Ψ-....... (15)

In view of the amplitude contribution of harmonic component more than secondary and secondary not over 10%, therefore neglect high order Harmonic wave, only retains first harmonic, and replaces β to express simplicity a_1c, replace β with b_1s, then waving angle expression formula can be with such as Lower expression formula describes:

β (Ψ)=β₀-acosΨ-bsinΨ (16)

If here, definition vectorSimultaneously take account of and obtained above wave the second order that angle is obeyed Expression formula (formula 14), brings the form (formula 15) solved into the differential equation, and available expression formula is as follows:

\overset{\cdot\cdot}{\overset{&RightArrow;}{a}} + D \overset{\cdot}{\overset{&RightArrow;}{a}} + K \overset{&RightArrow;}{a} = F - - - (17)

Wherein, D, K, F are the coefficient of second order differential equation.Theoretical research for wing tip plane finds, β₀Naturally frequency Rate is λ_βWith rotor to wave frequency consistent, much larger than the motion frequency of fuselage, it can thus be assumed that be that high-frequency signal is by hardware " eat up ", it is impossible to be coupled to body movement, can not consider.And first harmonic amount a, its natural frequency of b is λ_β-1 and body Motion frequency is close, therefore can be coupled to body movement.Therefore make in following presentationD, K, F all do corresponding fall Dimension.So in the case of meeting above-mentioned condition, dynamics research finds D, and K, F meet following formula:

Wherein,For rotor coefficient, A₁=B_latδ_lat, B₁=A_lonδ_lon, they are brought into above formula, and neglect The Derivative Terms slightly waved and body angular acceleration item (they impacts are the least after deliberation), the most just can obtain waving equation Differential equation of first order expression formula:

In like manner, consider further that the frame for movement of stable winglet, be illustrated in figure 6 the machinery of the stable winglet of depopulated helicopter Structure chart.By according to shown in Fig. 6, if that copies main rotor waves equation, c and d is made to be respectively the longitudinal direction of winglet with horizontal Wave angle first harmonic component, neglect horizontal stroke, longitudinal coupling simultaneously, winglet wave the expression that differential equation of first order should be following Form:

According to the mechanical driving structure in Fig. 6, stablize functioning as in amount horizontal, longitudinally controlled of winglet direct Superposition winglet waves angle, and the expression formula of such controlled quentity controlled variable can be rewritten as:

It addition, consider further that main rotor waves the controlled quentity controlled variable coupling terms not considered in modeling, ignore small quantity simultaneouslyWithThen can get complete main rotor and winglet to wave equation model as follows:

Wherein, above-mentioned various in, τ_f, A_b, B_b, A_lon, B_lat, C_lon, D_lat, K_cBeing parameter to be identified, a, b are respectively Representing that the longitudinal direction of described main rotor is waved angle and laterally waves the first harmonic component at angle, c, d represent the vertical of described winglet respectively To waving angle and laterally waving the first harmonic component at angle.Above-mentioned main rotor and winglet to wave equation model (formula 24) abundant Reflect control steering wheel driving be applied to rotor after, the response that rotor is made.Then, just above-mentioned rotor can be responded direct coupling Close the laterally and longitudinally equation of body movement, obtain the expression formula of linear acceleration and angular acceleration, wherein by by described Wave equation model respectively with described longitudinal velocity linearly think in absolute terms incremental model, described lateral velocity linearly think increasing in absolute terms Amount model, described roll angle speed linearly think in absolute terms incremental model, described pitch rate linearly think incremental model in absolute terms Carry out coupling and obtain longitudinal acceleration model, transverse acceleration model, and roll angle Fast track surgery, angle of pitch acceleration mould Type, wherein:

Wherein, X_u, X_θ, X_a, Y_v,Y_bIt is parameter to be identified；

Wherein, L_u, L_v, L_b, L_w, M_u, M_v, M_a, M_wIt is parameter to be identified.Consider further thatWithSo Parameterized linear differential equation has fully described the laterally and longitudinally kinetics of low speed and depopulated helicopter during hovering.

On the other hand, for vertical and the model of course passage, it is first-order system to have tested proof, and can recognize Tail-rotor for depopulated helicopter has the response characteristic being exceedingly fast, and is i.e. directly controlled quentity controlled variable to be applied to body movement, so Consider further that the x wire speed pair that impact and the vertical fin of course passage are caused by airborne stabilizing gyroscope angular velocity single order feedback The effect in course, the vertical passage that can obtain described depopulated helicopter is as follows with the kinetic model of course passage:

So, identical with the method for transverse and longitudinal modeling, consider further that the unmodel parts of coupling effect is compensated for, therefore boat To the kinetic model last with vertical channel, i.e. vertical acceleration model, course angle Fast track surgery can be expressed as:

Wherein, N_v, N_r, N_ped, K_rfb, K_r, Z_a, Z_b, Z_w, Z_r, Z_colFor parameter to be identified, r_fbFor course-stability gyroscope (often For steady course passage control electronic equipment, such as GV-1 etc.) control feedback quantity.

In sum, a total of 13 first increment differential equations of kinetic model obtained above (i.e. δ u, δ v, δ p, δ Q,δ θ, a, b, δ w, δ r, δ r_fb, the first increment differential equation of c, d), it is fully described depopulated helicopter operating point [u₀ v₀ w₀] ≈ [0 0 0], [p₀ q₀ r₀] dynamics near ≈ [0 0 0], i.e. dynamics under hovering flight mode, The form of upstate equation matrix is stated:

Parameter matrix therein and state vector are defined as:

U=(δ_lat δ_lon δ_ped δ_col)^T；

Wherein, X is 13 × 1 system mode vectors, and A is 13 × 13 sytem matrixes, and B is 13 × 4 control matrixes, and U is 4 × 1 Controlled quentity controlled variable matrix enters, the unknown parameters representated by symbol in above-mentioned matrix A and B, is parameter to be identified.

It addition, it is worth mentioning that the measurement equation of system can be carried out according to the sensor configuration of unmanned helicopter platform Selected, for general unmanned helicopter platform, three axial speed, acceleration, body angular speed, attitude angle all may be used To carry out on-line measurement by airborne sensor, example, measure equation and may include that

Wherein, I_m×mIt is m × m unit matrix, 0_m×nIt it is m × n null matrix.Above-mentioned equation (31) comprises 13 system shapes State variable and four control variable, and sytem matrix A and 44 unknown parameters controlling in matrix B.So, for controlling institute For the reference model used, matrix calculus process is sufficiently complex, simultaneously because some parameter cannot directly be measured, therefore can not Input under particular flight state (hovering, cruise etc.) obtains with output data identification.This phenomenon is referred to as parametrization and surpasses More, it is necessary to unknown parameter is increased corresponding constraint and just can carry out system identification.Meanwhile, strong between unknown parameter phase Guan Xing, and respectively control passage, as longitudinally with laterally, seriously coupled, cause initial parameter value to be difficult to select, be so directed to Equation uses discrimination method based on frequency domain sufficiently complex and calculation cost big, and the collection of test data, except hovering mode Outward, flying as before low speed, fly before high speed, high for handling the remote control skill requirement of hands, handling hands often cannot be at completion system It is maintained at while excitation in model of flight to be identified.

Therefore for the model structure of reduced parameter so that it is while reducing amount of calculation, degree of accuracy is sufficient for control Device design requirement processed, and suppress MODAL TRANSFORMATION OF A for the impact of control performance in the design of controller, reach to use single The model of flight of easy identification, such as hovering mode, corresponding model carries out the purpose of full envelope flight control, here, first come Carry out the modal model structure abbreviation that hovers.Hovering mode under, depopulated helicopter laterally and longitudinally kinetics can with Approximate Decoupling, Then can be following three part by above-mentioned kinetic model STRUCTURE DECOMPOSITION: horizontal dynamic, longitudinal dynamics, and vertical and The Coupled Dynamics in course, only controls matrix and has coupling parameter, and ignore the coupling ginseng in sytem matrix between every part Number, then can increase simultaneously wait to distinguish by coupling parameter zero setting in the longitudinal direction sytem matrix with interconnection in controlling matrix The parameter known compensates Coupled Dynamics.So, every part relates only to 2 and controls passage and at most 5 state variables, distinguishes Knowing complexity to be substantially reduced, this three department patterns structure is represented by following three groups of state equations, wherein horizontal dynamic models Including:

y_lon=(I_3×3 0_3×2)δX_lon=C_lonδX_lon

Wherein, M_lon, M_lat, X_lon, X_lat, C_lon, C_lat, X_lon, y_lonIt is parameter to be identified.

Longitudinal Dynamic Model includes:

y_lat=(I_3×3 0_3×2)δX_lat=C_latδX_lat

Wherein, D_lon, D_lat, y_latIt is parameter to be identified.

Vertical and course Coupling Dynamic Model includes:

y_yaw-heave=(I_2×2 0_2×1)δX_yaw-heave=C_yaw-heaveδX_yaw-heave

Wherein, A_yaw-heave, X_yaw-heave, B_yaw-heave, u_yaw-heave, y_yaw-heave, X_yaw-heaveIt is parameter to be identified, its In_yaw-heaveRefer to course-vertical.

It addition, for the aerial mission that depopulated helicopter is conventional, the transformation of model of flight generally includes hovering (speed Less than 5 meter per seconds), low cruise (speed more than 5 less than 15 meter per seconds), high-performance cruise (speed is more than 15 meter per seconds) and landing Mode (near-earth distance, less than 2 times of main rotor diameters, imitates impact with existing), also has the transient state in these four kinds of mutual transition of mode, And wherein only hovering mode defines the clearest and the most definite, other mode and transient state are the definition of empirical form mostly, it is difficult to according to number According to accurately judging.Therefore under floating state, identification system kinetic model simplicity the most is the most feasible, simultaneously under this mode Depopulated helicopter, operating point range width, carry out less action just can excitation system, data acquisition danger is minimum, and this is also The application use it as the major reason of Reference model for control system.

About parameter frequency domain identification method, the input of system can be been described by as follows with output by frequency response function:

Y (j ω)=H (j ω) U (j ω) (36)

Wherein, U (j ω) and Y (j ω) is system input u (t) and the Fourier Tranform of system output y (t), H (j respectively ω) it is the Fourier transform function of system impulse response.It is the discrete of finite length for discrete sampling system, u (t) and y (t) Sample, then can pass through discrete Fourier transform (DFT) and describe the relation of input and output at frequency domain, be shown below:

Wherein, w_k=k Ω_sIt is discrete point in frequency, Ω_s=2 π/(NT_s) it is frequency sampling interval, T_sIt is sampling time interval, N=T_d/T_sIt is sampling number, T_dIt it is sampled data segment length.The auto-correlation density function G of input_uu(j ω) and input and output are assisted Correlation density function G_uy(j ω) can be estimated with following formula:

Then discrete point in frequency ω_kFrequency response be estimated as follows；

Association's correlation function of amplitude is defined at thisCorrelation degree between u (t) and y (t) described:

Wherein, G_yy(j ω) is output auto-correlation spectrum density function.

According to above-mentioned formula (39～41) may determine that frequency domain estimation difference majorized function (minimizing cost function) as Under:

Wherein ω_iIt is Frequency point, ε (ω_i, Θ) and it is that current estimated value based on model parameter vector Θ to be identified records in advance To system output amplitude deviation and the vector that formed of phase deviation (determining according to above-mentioned formula 39～41), W (ω_i) It is to optimize weight matrix.Parameter Θ to be identified can obtain by minimizing cost function (formula 43).Additionally frequency domain optimizes Cost function is complex nonlinear form, the numerical method iterative Θ that can be searched for by cosecant.

About the parameter frequency domain identification method of multiple-input, multiple-output (Multiple Input Multiple Output, MIMO), The state-space model that acquisition MIMO goes out:

Wherein, Θ represents unknown parameter to be identified, and X is the system mode vector of n × 1 dimension, and U represents the control that r × 1 is tieed up Moment matrix, Y represents that system output is maintained in p × 1, and A is that n × n ties up sytem matrix, and B represents that n × r dimension controls matrix, and C represents that p × n ties up Calculation matrix；Wherein, U includes: U=(δ_lat δ_lon δ_ped δ_col)^T, δ_latRepresent crosswise joint input quantity, δ_lonRepresent longitudinally control Input quantity processed, δ_pedRepresent Heading control input quantity, δ_colRepresent that vertical control inputs；Parameter in A and B is described ginseng to be identified Number.The most just complete the process of the described default state-space model of above-mentioned acquisition.

For mimo system, aforesaid cost function (formula 43) is the most applicable, it is therefore desirable to design new can describing All MIMO input and optimizing index functions exported.For MIMO linear parameterization state-space model, its impulse response square Battle array can be expressed as:

T (jw, Θ)=C (Θ) [jwI-A (Θ)]^-1B(Θ)+D(Θ) (45)

So, if having chosen discrete characteristic frequency point (ω₁,ω₂......ω_n), then, based on former optimizing index The expression-form of function (i.e. company 43), can define new optimizing index function as follows:

Wherein Δ | |_lIt is that the l impulse response function is at ω_kThe output predicted magnitude of point and actual output response measurement The deviation of amplitude, Δ ∠_lIt is that the l impulse response function is at ω_kThe output predicted phase of point and actual output response measurement phase The deviation of position；W_γIt is the weight matrix being decided by amplitude association correlation function, W_gAnd W_pIt is amplitude and the phase estimation difference of two squares respectively Weight matrix, n_TIt it is T (jw, Θ) entry of a matrix prime number amount.

But, can there are following two problems in actual applications in the identification flow process of said method: owing to using cosecant Searching method, the initial estimate of parameter to be identified needs artificial selected, and this job demand enters based on substantial amounts of flying quality The craft that pedestrian is is attempted, and does not has fixing selection rule；Moreover, iterative process does not has convergence criteria, it is also desirable to people For selected iterations, then making repeated attempts and just can search optimized parameter estimated value, empirical factor is big, and often obtain is Suboptimal solution, it is impossible to obtain the unbiased esti-mator of parameter, even causes nonlinear iteration search numerical value instability problem to cause iteration to send out Dissipate.For two above application problem, choose from initial value separately below and provide solution in terms of iterative process convergence criteria The convergence theoretical proof of the identification flow process after the scheme of determining and improvement.

In order to obtain the initial estimate of parameter to be identified, it is defined as follows one-step prediction deviation criterion at this:

Wherein,It is the system mode one-step prediction of t, when dt is sampling Between.If setting column vectorRepresent the parameter to be identified of the i-th row in matrix [A (θ) B (θ)], then can obtain:

Wherein φ_i,tIt it is tCorresponding coefficient vector, F_i,tIt it is one-step prediction vectorKnown portions.

IfThen have:

Wherein,F_i=(F_i,1F_i,2...F_i,N)

In order to minimize target function (formula 49), i.e.ThenIt must is fulfilled forThen have Then pass through linear regression optimization and can solve the initial estimate of parameter to be identifiedAs follows:

Therefore, the initial estimate obtaining parameter to be identified in step 101 can pass through above-mentioned formula (47,50) acquisition.

For the convergence criteria of successive iterations process, the judgement inequality that can be defined as follows:

||A⁽ⁱ⁾ _(s+1)-A⁽ⁱ⁾ _(s)||+||B⁽ⁱ⁾ _(s+1)-B⁽ⁱ⁾ _(s)| | ＜ ε (51)

Wherein A⁽ⁱ⁾ _(s)It is to estimate for the s time to A (θ) the i-th row unknown element, B⁽ⁱ⁾ _(s)It is the s to B (θ) the i-th row element Secondary estimation, ε ∈ R is can the threshold value of manual setting.In order to ensure stablizing of non-linear optimal searching algorithm, it is necessary to avoid parameter Sudden change.Following more New Policy of use at this:

Wherein α ∈ (0,1),WithIt is by minimizing the current estimated value of parameter obtained by cost function (formula 47), A_(s)It is parameter matrix A the estimated value of the s time, B_(s)Thing parameter matrix B is the estimated value of the s time.Thus just obtained step More New Policy described in 101.

Therefore the flow process shown in Fig. 1 can be understood as a kind of enhancement mode frequency domain identification flow process, step 101 in this flow process～ 105 can being summarized as with example:

1) flow process initial set s=1, then utilizes formula (50) and minimizes cost function (formula 47) and calculate and wait to distinguish Know the initial estimate of parameter, obtain A₍₀₎And B₍₀₎；

2) according to A₍₀₎And B₍₀₎, calculate the receptance function of system measurement data based on formula (37～41), and choose object Characteristic frequency point (ω₁,ω₂......ω_n)；

3) according to A₍₀₎And B₍₀₎, described receptance function, characteristic frequency point (ω₁,ω₂......ω_n), formula 44 and Littleization cost function (46) performs workflow management parametric optimal solution shown in Fig. 2, and wherein matrix A and the B of formula 44 replaces with respectively A₍₀₎And B₍₀₎, obtain

4) according to the more New Policy shown in formula 52～53 and calculating A_(s+1)And B_(s+1)；

If being unsatisfactory for pre-conditioned, then take s=s+1, by A_(s+1)And B_(s+1)Replace original A₍₀₎And B₍₀₎, and according to Need to adjust the value of some of which parameter, then by current A_(s+1)And B_(s+1)The value of middle parameter as next round iteration wait distinguish Know the estimated value of parameter, then return step 2) again perform above-mentioned flow process；If meeting pre-conditioned, iteration terminates, and obtainsAndWhereinAndIn the value of each parameter to be identified complete exactly The desired value of identification.Wherein, above-mentioned pre-conditioned can be the condition shown in formula 51.Further it will be understood that this Shen Parameter to be identified involved in please be sytem matrix A above and controls the unknown parameter that matrix B includes, in the application Involved other unknown parameters being not included in this sytem matrix A and control matrix B can be understood as determining system square Intermediate parameters during battle array A and control matrix B, it is not necessary to carry out identification.

Proof procedure be presented herein below:

By minimizing parameter initial estimate obtained by cost functionIt is the unbiased of actual value θ and consistent estimation, front The enhancement mode frequency domain identification flow process that literary composition is proposed meets:

Therefore, it was demonstrated that:

Again

Because θ=(θ₁...θ_n), thereforeIt it is the unbiased consistent Estimation of θ.For the cost being made up of amplitude and phase angle quadratic form Function, it has proved that and if only if iteration initial valueIt is to have during the unbiased consistent Estimation of θ Wherein N is the number of selected characteristic frequency point.

So, consider further that optimum results updates step (formula 52 and formula 53), and for ith iteration, have:

Wherein α is to minimize the estimates of parameters obtained by cost function (formula 47), then the estimation of the s time iteration is nothing Inclined and consistent, consider further that and select different characteristic frequency points to estimate then have for ith iteration:

Card is finished.

It addition, about being described as follows of collection of flying quality:

For flight every time, flight control hands will control passage by remote control equipment to the four of rotor flying robot In (vertical passage, interconnection, vertical passage and course passage) one applying frequency sweep sequence inputting, use simultaneously other three Rotor flying robot is maintained at hovering flight mode (speed is less than 5m/s) by the controlled quentity controlled variable of individual passage.Inputted by frequency sweep, The sine wave applying approximation 1～20Hz carrys out the response of excitation system different frequency, provides wider frequency band for Model Distinguish.Root According to the characteristic frequency region that the rotor of rotor flying robot is different with fuselage, the low-frequency data less than 3Hz is used for fuselage power The department of the Chinese Academy of Sciences divides the identification of (formula 25～30), and the part higher than 10Hz is for the kinetics of rotor system part (formula 21～24) Identification, the most just can carry out parameter identification by the frequency domain estimation method proposed above.Additionally in the middle of flight experiment, all Controlled quentity controlled variable can the state of flight variable of examining system quantity of state be recorded with the sample frequency of 50Hz.The aforementioned flight collected Data can be passed through band elimination filter (-3dB, at 10Hz) and filter mechanism's vibrations noise.It addition, to each control passage record The frequency sweep input of 30s is used for identification with output response data.It addition, the hovering flight data (non-frequency sweep) acquiring 5s are used for distinguishing Know the simulating, verifying of result.

Finally, in order to verify the parameter identification degree of accuracy of model structure and discrimination method further, it is defined as follows at this Deviation root-mean-square criterion:

Wherein y_i(t) andBeing measured value and the model predication value of i-th output variable, N is the dimension of output vector. So, criterion V directly reflects the mean error of model, is shown the ginseng of frequency domain estimation method that the disclosure used by emulation The model structure of half decoupling of number estimated result and simplification can describe rotor flying robot under hovering mode accurately Dynamics, therefore, the model of gained can apply to control depopulated helicopter and Control System Design.

Fig. 7 is the structural frames of the device for identifying of a kind of depopulated helicopter kinetic parameter that the disclosure one embodiment provides Figure, sees Fig. 7, and described device for identifying 700 includes:

Initial estimation module 710, for perform a. obtain parameter to be identified initial estimate, comprised described in wait to distinguish The first system matrix and first knowing the initial estimate of parameter controls matrix；

RESPONSE CALCULATION module 720, is used for performing b. and controls matrix according to described the first system matrix and described first, obtain The receptance function of system measurement data, and determine the characteristic frequency point of described receptance function；

Iterative processing module 730, for performing the c. flying quality according to the described depopulated helicopter gathered, described first Sytem matrix and described first controls matrix, utilizes state-space model, described receptance function and the default cost preset Function obtains the optimal solution of described parameter to be identified, obtains comprising sytem matrix and the coefficient of the optimal solution of described parameter to be identified Matrix；

More new module 740, for perform d. according to described in comprise described parameter to be identified optimal solution sytem matrix with Coefficient matrix, utilizes the more New Policy preset to obtain the second system matrix after updating and second and controls matrix；

For e., judge module 750, judges that described second system matrix and described second controls whether matrix meets default Condition, described second system matrix and described second control matrix be unsatisfactory for described pre-conditioned time, by described second system Matrix and described second controls matrix and again performs step b～e as described the first system matrix and described first control matrix, Until described second system matrix and described second control matrix meet described pre-conditioned till, by described second system matrix With the currency of the parameter described to be identified comprised in described second system matrix as the desired value completing identification.

Optionally, described iterative processing module 730, including:

Optionally, as shown in Figure 8, this device 700 also includes: model obtains submodule 760, in described step c2 Before:

Wherein, u₀, v₀, w₀, p₀, q₀, r₀Representing u respectively, the original state of v, w, p, q, r, δ represents relative to initial shape The absolute increment of state；

Wherein, X_u, X θ, X_a, Y_v,Y_bFor parameter to be identified；

Wherein, N_v, N_r, N_ped, K_rfb, K_rFor parameter to be identified, r_fbRepresent the control feedback quantity of course-stability gyroscope, δ_ped Represent Heading control input quantity, δ_colRepresent vertical control input quantity；

Described vertical acceleration model, course angle Fast track surgery include:

Wherein, Z_a, Z_b, Z_w, Z_r, Z_colFor parameter to be identified.

Optionally, described state-space model includes:

Wherein Θ represents unknown parameter to be identified, and X is the system mode vector of n × 1 dimension, and U represents the control that r × 1 is tieed up Moment matrix, Y represents that system output is maintained in p × 1, and A is that n × n ties up sytem matrix, and B represents that n × r dimension controls matrix, and C represents that p × n ties up Calculation matrix；

The preferred implementation of the disclosure is described in detail above in association with accompanying drawing, but, the disclosure is not limited to above-mentioned reality Execute the detail in mode, in the technology concept of the disclosure, multiple letter can be carried out with technical scheme of this disclosure Monotropic type, these simple variant belong to the protection domain of the disclosure.

It is further to note that each the concrete technical characteristic described in above-mentioned detailed description of the invention, at not lance In the case of shield, can be combined by any suitable means, in order to avoid unnecessary repetition, the disclosure to various can The compound mode of energy illustrates the most separately.

Additionally, combination in any can also be carried out between the various different embodiment of the disclosure, as long as it is without prejudice to this Disclosed thought, it should be considered as disclosure disclosure of that equally.

Claims

1. the discrimination method of a depopulated helicopter kinetic parameter, it is characterised in that described method includes:

A. obtain the initial estimate of parameter to be identified, obtain comprising the first system of the initial estimate of described parameter to be identified Matrix and first controls matrix；

B. according to described the first system matrix and described first control matrix, the receptance function of acquisition system measurement data, and really The characteristic frequency point of fixed described receptance function；

C. according to the flying quality of described depopulated helicopter gathered, described the first system matrix and described first control matrix, The state-space model preset, described receptance function and default cost function is utilized to obtain the optimum of described parameter to be identified Solve, obtain comprising sytem matrix and the coefficient matrix of the optimal solution of described parameter to be identified；

D. comprise sytem matrix and the coefficient matrix of the optimal solution of described parameter to be identified described in basis, utilize the renewal plan preset Slightly obtain the second system matrix after updating and second and control matrix；

E. judge that described second system matrix and described second controls whether matrix meets pre-conditioned, at described second system square Battle array and described second control matrix be unsatisfactory for described pre-conditioned time, by described second system matrix and described second control matrix As described the first system matrix and described first control matrix again perform step b～e, until described second system matrix and Described second control matrix meet described pre-conditioned till, will in described second system matrix and described second system matrix wrap The currency of the parameter described to be identified contained is as the desired value completing identification.

The discrimination method of depopulated helicopter kinetic parameter the most according to claim 1, it is characterised in that described step c Including:

C1. gather the flying quality of described depopulated helicopter, and according to the described flying quality gathered, utilize described receptance function Acquisition actual spectrum responds；

C2. described the first system matrix and described first are controlled matrix as the system square in described state-space model Battle array and control matrix calculus spectral response estimated value；

C4. according to described error of spectrum, the optimal solution of the epicycle iteration of the described cost function described parameter to be identified of acquisition is utilized；

C5. update described the first system matrix and described first according to the optimal solution of the epicycle iteration of described parameter to be identified to control Matrix, and the first system matrix after updating and the first control matrix after updating are as in described state-space model Sytem matrix and control matrix, obtain epicycle iteration obtain state-space model；

C6. the spectral response estimated value that calculates of state-space model that judges to utilize epicycle iteration to obtain by time domain checking and Whether the error of described actual spectrum response is more than error threshold；

When the spectral response estimated value that calculates of state-space model judging to utilize epicycle iteration to obtain and described actual spectrum The error of response is when being not more than described error threshold, using the optimal solution of described epicycle iteration as the optimum of described parameter to be identified Solve, obtain comprising sytem matrix and the coefficient matrix of the optimal solution of described parameter to be identified；When judging to utilize epicycle iteration to obtain The error of the spectral response estimated value that calculates of state-space model and the response of described actual spectrum more than described error threshold Time, update described the first system matrix and described first and control matrix and again perform step c2～c5, until utilizing epicycle repeatedly The spectral response estimated value that the state-space model that generation obtains calculates is not more than described with the error of described actual spectrum response Till error threshold.

Method the most according to claim 1, it is characterised in that described flying quality includes the longitudinal direction of described depopulated helicopter Speed, lateral velocity, vertical velocity, longitudinal acceleration, transverse acceleration, vertical acceleration, roll angle, the angle of pitch, course angle, Roll angle speed, pitch rate and course angle speed.

The discrimination method of depopulated helicopter kinetic parameter the most according to claim 3, it is characterised in that in described step Before c2, also include:

According to described longitudinal velocity, described lateral velocity, described vertical velocity, described roll angle speed, described pitch rate And the first differential model of described course angle speed, obtain described longitudinal velocity, described lateral velocity, described vertical speed respectively Degree, described roll angle speed, described pitch rate and described course angle speed relevant to system mode and controlled quentity controlled variable Linearly think incremental model in absolute terms；

By by described wave equation model respectively with described longitudinal velocity linearly think incremental model, described lateral velocity in absolute terms Linearly think in absolute terms incremental model, described roll angle speed linearly think in absolute terms incremental model, described pitch rate linear Absolutization incremental model carries out coupling and obtains longitudinal acceleration model, transverse acceleration model, and roll angle Fast track surgery, Angle of pitch Fast track surgery；

By the kinetic model of described vertical passage and course passage respectively with described vertical velocity linearly think difference module in absolute terms The linear absolutization incremental model of type and described course angle speed carries out coupling and obtains vertical acceleration model, course angle acceleration Degree model；

According to described longitudinal acceleration model, transverse acceleration model, and roll angle Fast track surgery, angle of pitch acceleration mould The control feedback quantity model of type, described vertical acceleration model, course angle Fast track surgery and course-stability gyroscope determines State equation matrix；

Described state-space model, the sytem matrix of described state-space model and control is determined according to described state equation matrix Matrix includes that described parameter to be identified, described parameter to be identified are according to described longitudinal acceleration model, transverse acceleration model, And roll angle Fast track surgery, angle of pitch Fast track surgery, described vertical acceleration model, course angle Fast track surgery and The parameter determination to be identified controlled in feedback quantity model of course-stability gyroscope.

The discrimination method of depopulated helicopter kinetic parameter the most according to claim 4, it is characterised in that described based on institute State longitudinal velocity, described lateral velocity, described vertical velocity, described roll angle, the described angle of pitch, described course angle, described rolling The first differential model of corner speed, described pitch rate and described course angle speed includes:

Wherein, u, v, w represent described longitudinal velocity, described lateral velocity and described vertical velocity respectively, and p, q, r represent respectively Described roll angle speed, described pitch rate and described course angle speed,θ, ψ represent described roll angle respectively, described in bow The elevation angle and described course angle, I_xx, I_yy, I_zzRepresenting the rotary inertia of X-axis, Y-axis and Z axis respectively, m represents described depopulated helicopter Quality, L, M, N are three axle moments of torsionIn parameter, X, Y, Z are for acting on described depopulated helicopter weight The force vector of the heartIn parameter；

δ \overset{\cdot}{u} = (- w_{0} δ q - {δwq}_{0} + v_{0} δ r + {δvr}_{0}) + δ X / m

δ \overset{\cdot}{v} = (- u_{0} δ r - {δur}_{0} + w_{0} δ p + {δwp}_{0}) + δ Y / m

δ \overset{\cdot}{w} = (- v_{0} δ p - {δvp}_{0} + u_{0} δ q + {δuq}_{0}) + δ Z / m

δ \overset{\cdot}{p} = (- q_{0} δ r - {δqr}_{0}) (I_{y y} - I_{z z}) / I_{x x} + δ L / I_{x x}

δ \overset{\cdot}{q} = (- p_{0} δ r - {δpr}_{0}) (I_{z z} - I_{x x}) / I_{y y} + δ M / I_{y y}

δ \overset{\cdot}{r} = (- p_{0} δ q - {δpq}_{0}) (I_{x x} - I_{y y}) / I_{z z} + δ N / I_{z z}

Wherein, u₀, v₀, w₀, p₀, q₀, r₀Representing u respectively, the original state of v, w, p, q, r, it is exhausted that δ represents relative to original state To increment；

\{\begin{matrix} τ_{f} \overset{\cdot}{a} = - a - τ_{f} q + A_{b} b + A_{l o n} (δ_{l o n} + K_{c} c) + A_{l a t} δ_{l a t} \\ τ_{f} \overset{\cdot}{b} = - b - τ_{f} p + B_{b} b + B_{l a t} (δ_{l a t} + K_{d} d) + B_{l o n} δ_{l o n} \end{matrix}

Wherein, τ_f, A_b, B_b, A_lon, B_lat, K_cFor parameter to be identified, a, b represent that the longitudinal direction of described main rotor waves angle and horizontal stroke respectively To waving the first harmonic component at angle, c, d represent that the longitudinal direction of described winglet is waved angle and laterally waves the first harmonic at angle respectively Component, δ_latRepresent crosswise joint input quantity, δ_lonRepresent longitudinally controlled input quantity；

Wherein, X_u, X_θ, X_a, Y_v,Y_bFor parameter to be identified；

\{\begin{matrix} δ \overset{\cdot}{p} = L_{u} δ u + L_{v} δ v + L_{b} b + L_{w} δ w \\ δ \overset{\cdot}{q} = M_{u} δ u + M_{v} δ v + M_{a} a + M_{w} δ w \end{matrix}

\{\begin{matrix} \overset{\cdot}{r} = N_{v} v + N_{r} r + N_{p e d} (δ_{p e d} - r_{f b}) \\ {\overset{\cdot}{r}}_{f b} = - K_{r f b} r_{f b} + K_{r} r \end{matrix}

\overset{\cdot}{w} = (- v_{0} p + u_{0} q) + Z_{w} w + Z_{c o l} δ_{c o l}

Wherein, N_v, N_r, N_ped, K_rfb, K_rFor parameter to be identified, r_fbRepresent the control feedback quantity of course-stability gyroscope, δ_pedRepresent Heading control input quantity, δ_colRepresent vertical control input quantity；

Described vertical acceleration model, course angle Fast track surgery include:

\{\begin{matrix} δ \overset{\cdot}{r} = N_{v} δ v + N_{p} δ p + N_{r} δ r + N_{p e d} (δ_{p e d} - δ r_{f b}) \\ δ {\overset{\cdot}{r}}_{f b} = - K_{r f b} δ r_{f b} + K_{r} δ r \end{matrix}

δ \overset{\cdot}{w} = (- v_{0} δ p + u_{0} δ q) + Z_{a} a + Z_{b} b + Z_{w} δ w + + Z_{r} δ r + Z_{c o l} δ_{c o l}

Wherein, Z_a, Z_b, Z_w, Z_r, Z_colFor parameter to be identified.

6. according to the discrimination method of the arbitrary described depopulated helicopter kinetic parameter of claim 1-4, it is characterised in that described State-space model includes:

\{\begin{matrix} \overset{\cdot}{X} = A (Θ) X + B (Θ) U \\ Y = C (Θ) X \end{matrix}

Wherein, Θ represents unknown parameter to be identified, and X is the system mode vector of n × 1 dimension, and U represents the controlled quentity controlled variable square tieed up r × 1 Battle array, Y represents that system output is maintained in p × 1, and A is that n × n ties up sytem matrix, and B represents that n × r dimension controls matrix, and C represents that p × n dimension is measured Matrix；

Wherein, U includes: U=(δ_lat δ_lon δ_ped δ_col)^T, δ_latRepresent crosswise joint input quantity, δ_lonRepresent longitudinally controlled input Amount, δ_pedRepresent Heading control input quantity, δ_colRepresent that vertical control inputs；Parameter in A and B is described parameter to be identified.

7. the device for identifying of a depopulated helicopter kinetic parameter, it is characterised in that described device includes:

Initial estimation module, obtains the initial estimate of parameter to be identified, obtains comprising described parameter to be identified for performing a. The first system matrix of initial estimate and first control matrix；

RESPONSE CALCULATION module, is used for performing b. and controls matrix according to described the first system matrix and described first, and the system that obtains is surveyed The receptance function of amount data, and determine the characteristic frequency point of described receptance function；

Iterative processing module, for performing the c. flying quality according to the described depopulated helicopter gathered, described the first system square Battle array and described first controls matrix, utilizes state-space model, described receptance function and the default cost function preset to obtain Take the optimal solution of described parameter to be identified, obtain comprising sytem matrix and the coefficient matrix of the optimal solution of described parameter to be identified；

More new module, for perform d. according to described in comprise described parameter to be identified the sytem matrix of optimal solution and coefficient square Battle array, utilizes the more New Policy preset to obtain the second system matrix after updating and second and controls matrix；

For e., judge module, judges that described second system matrix and described second controls whether matrix meets pre-conditioned, Described second system matrix and described second control matrix be unsatisfactory for described pre-conditioned time, by described second system matrix and institute State the second control matrix and again perform step b～e as described the first system matrix and described first control matrix, until described Second system matrix and described second control matrix meet described pre-conditioned till, by described second system matrix and described the The currency of the parameter described to be identified comprised in two system matrix is as the desired value completing identification.

The device for identifying of depopulated helicopter kinetic parameter the most according to claim 7, it is characterised in that at described iteration Reason module, including:

Data acquisition and RESPONSE CALCULATION submodule, gather the flying quality of described depopulated helicopter for performing c1., and according to adopting The described flying quality of collection, utilizes described receptance function to obtain actual spectrum response；

Submodule is estimated in response, is used for performing c2. and described the first system matrix and described first are controlled matrix as institute State the sytem matrix in state-space model and control matrix calculus spectral response estimated value；

Error obtains submodule, is used for performing c3. and responds acquisition frequency according to described spectral response estimated value and described actual spectrum Spectrum error；

Calculating sub module, is used for performing c4. according to described error of spectrum, utilizes described cost function to obtain described parameter to be identified The optimal solution of epicycle iteration；

Model modification submodule, updates described for performing c5. according to the optimal solution of the epicycle iteration of described parameter to be identified One sytem matrix and described first control matrix, and will update after the first system matrix and update after first control matrix divide Not as the sytem matrix in described state-space model and control matrix, obtain the state-space model that epicycle iteration obtains；

By time domain checking, checking submodule, judges that the state-space model utilizing epicycle iteration to obtain calculates for performing c6. Whether the spectral response estimated value gone out is more than error threshold with the error of described actual spectrum response；

The device for identifying of depopulated helicopter kinetic parameter the most according to claim 7, it is characterised in that described flight number According to including the longitudinal velocity of described depopulated helicopter, lateral velocity, vertical velocity, longitudinal acceleration, transverse acceleration, vertical adding Speed, roll angle, the angle of pitch, course angle, roll angle speed, pitch rate and course angle speed.

The device for identifying of depopulated helicopter kinetic parameter the most according to claim 9, it is characterised in that also include: mould Type obtains submodule, for before described step c2:

The device for identifying of 11. depopulated helicopter kinetic parameters according to claim 10, it is characterised in that described based on Described longitudinal velocity, described lateral velocity, described vertical velocity, described roll angle, the described angle of pitch, described course angle, described The first differential model of roll angle speed, described pitch rate and described course angle speed includes:

Wherein, u, v, w represent described longitudinal velocity, described lateral velocity and described vertical velocity respectively, and p, q, r represent respectively Described roll angle speed, described pitch rate and described course angle speed,θ, ψ represent described roll angle respectively, described in bow The elevation angle and described course angle, I_xx, I_yy, I_zzRepresenting the rotary inertia of X-axis, Y-axis and Z axis respectively, m represents described depopulated helicopter Quality, L, M, N are three axle moments of torsionIn parameter, X, Y, Z are for acting on described depopulated helicopter center of gravity Force vectorIn parameter；

δ \overset{\cdot}{u} = (- w_{0} δ q - {δwq}_{0} + v_{0} δ r + {δvr}_{0}) + δ X / m

δ \overset{\cdot}{v} = (- u_{0} δ r - {δur}_{0} + w_{0} δ p + {δwp}_{0}) + δ Y / m

δ \overset{\cdot}{w} = (- v_{0} δ p - {δvp}_{0} + u_{0} δ q + {δuq}_{0}) + δ Z / m

δ \overset{\cdot}{p} = (- q_{0} δ r - {δqr}_{0}) (I_{y y} - I_{z z}) / I_{x x} + δ L / I_{x x}

δ \overset{\cdot}{q} = (- p_{0} δ r - {δpr}_{0}) (I_{z z} - I_{x x}) / I_{y y} + δ M / I_{y y}

δ \overset{\cdot}{r} = (- p_{0} δ q - {δpq}_{0}) (I_{x x} - I_{y y}) / I_{z z} + δ N / I_{z z}

\{\begin{matrix} τ_{f} \overset{\cdot}{a} = - a - τ_{f} q + A_{b} b + A_{l o n} (δ_{l o n} + K_{c} c) + A_{l a t} δ_{l a t} \\ τ_{f} \overset{\cdot}{b} = - b - τ_{f} p + B_{b} b + B_{l a t} (δ_{l a t} + K_{d} d) + B_{l o n} δ_{l o n} \end{matrix}

Wherein, X_u, X_θ, X_a, Y_v,Y_bFor parameter to be identified；

\{\begin{matrix} δ \overset{\cdot}{p} = L_{u} δ u + L_{v} δ v + L_{b} b + L_{w} δ w \\ δ \overset{\cdot}{q} = M_{u} δ u + M_{v} δ v + M_{a} a + M_{w} δ w \end{matrix}

\{\begin{matrix} \overset{\cdot}{r} = N_{v} v + N_{r} r + N_{p e d} (δ_{p e d} - r_{f b}) \\ {\overset{\cdot}{r}}_{f b} = - K_{r f b} r_{f b} + K_{r} r \end{matrix}

\overset{\cdot}{w} = (- v_{0} p + u_{0} q) + Z_{w} w + Z_{c o l} δ_{c o l}

Described vertical acceleration model, course angle Fast track surgery include:

\{\begin{matrix} δ \overset{\cdot}{r} = N_{v} δ v + N_{p} δ p + N_{r} δ r + N_{p e d} (δ_{p e d} - δ r_{f b}) \\ δ {\overset{\cdot}{r}}_{f b} = - K_{r f b} δ r_{f b} + K_{r} δ r \end{matrix}

δ \overset{\cdot}{w} = (- v_{0} δ p + u_{0} δ q) + Z_{a} a + Z_{b} b + Z_{w} δ w + + Z_{r} δ r + Z_{c o l} δ_{c o l}

Wherein, Z_a, Z_b, Z_w, Z_r, Z_colFor parameter to be identified.

12. according to the device for identifying of the arbitrary described depopulated helicopter kinetic parameter of claim 7-11, it is characterised in that institute State state-space model to include:

\{\begin{matrix} \overset{\cdot}{X} = A (Θ) X + B (Θ) U \\ Y = C (Θ) X \end{matrix}