CN101846519B - Method for predicting flight technical error of lateral flight path control system - Google Patents

Abstract

The invention discloses a method used for predicting a flight technical error of a lateral flight path control system. The method comprises the following steps of: acquiring a closed loop transfer function of the lateral flight path control system of a pointed aircraft type or minimal state-space realization thereof; obtaining a side wind section and a wind speed at 20ft according to different heights (less than 1,000ft, between 1,000ft and 2,000ft, and greater than 2,000ft) of a considered area or by a relevant meteorological department, and further calculating to obtain a turbulence intensity; or obtaining the turbulence intensity by using MIL-8785C; calculating power spectral density of an input disturbing signal; calculating the FTE (flight technical error) of the lateral flight path control system and an FTE evaluation boundary; and if turbulence intensity data cannot be obtained, calculating the FTE values at three different atmospheric turbulence intensities, namely a light atmospheric turbulence intensity, a medium atmospheric turbulence intensity and a high atmospheric turbulence intensity, and then calculating an expected value of FTE. Due to the adoption of the method, accurate prediction before flight and short-term prediction in the fight of TSE can be performed on flight operation under PBN (performance-based navigation).

Description

A kind of Flight technical error Forecasting Methodology that is used for the side direction flight path control system

Technical field

The present invention relates to a kind of Flight technical error Forecasting Methodology based on the flight path control system master gain.

Background technology

Prediction and navigation prediction a middle or short term before the enforcement of PBN (based on the navigation of performance) need be navigated to TSE (system's total error); And FTE (Flight technical error) and NSE (Navigation system error) are two chief components of TSE, therefore will directly influence the execution of PBN to the accurately predicting algorithm of FTE.The PBN navigation is based on the navigation concepts of multisensor; By RNP (Required Navigation Performance; Required navigation performance) and RNAV (Regional NAVigation; Area navigation) forms, mainly rely on the characteristic such as autonomous, complete of the high precision, high coverage rate of GNSS (Global Navigation Satellite System, GPS), round-the-clock and inertial navigation; Comprehensive other navigational system (like radio navigation system etc.) optimal combination, reach possible optimum navigation performance.

In September, 2007, International Civil Aviation Organization (ICAO) requires each contracting party before the end of the year 2009, to formulate to accomplish the PBN implementary plan, accomplishes whole implementations in 2016, carries out the transition to PBN with the whole world one coordinated mode of making peace fully from traditional offline mode.PBN uses airport construction, navigational facility layout and the spatial domain of China civil aviaton will produce significant impact; To the lasting safety of effective promotion civil aviaton, improve flight quality, increase airspace capacity, minimizing surface facility input and energy-saving and emission-reduction etc. and all have significant positive role.

FTE is one of chief component of TSE; And more become main TSE source (because multi-sensor combined navigation locate the earth to the utmost reduced NSE) advancing nearly stage F TE; Therefore very important to the accurate prediction of FTE, it calculates with measuring in real time and is all accomplished by Flight Management Computer (FMS).Abroad except that Boeing provides the statistical value of FTE of a small amount of type, there is not any related data; Domestic also still do not have any related data.

The present invention predicts for the boat of the side direction FTE in each stage of navigation is preceding and navigation prediction a middle or short term provides exact algorithm and boundary estimation solution.

Summary of the invention

Technology of the present invention is dealt with problems: the deficiency that overcomes prior art; Provide a kind of Flight technical error that is used for the side direction flight path control system to confirm method, this method makes the flight operation under the PBN navigation can carry out preceding prediction of accurate boat and the underway short-term forecasting of TSE.

Technical solution of the present invention, a kind of Flight technical error that is used for the side direction flight path control system is confirmed method, performing step is following:

(1) obtains the side direction automatic flight control system (Automatic Flight ControlSystem) of the aircraft type that is directed against

Closed loop transfer function, matrix G (s) is shown below: G (s)=C (sI-A) ^-1B+D, wherein s is the Laplace domain variable, I is and the unit matrix of matrix A with line number same column number; Or a minimum state space that obtains G (s) realizes, shape be shown below (general piecemeal notation)

Wherein A, B, C, D are four constant matricess,

(2) obtain crosswind section and the wind speed W of 20ft place according to the flying height difference of being concerned about the zone or through relevant meteorological department ₂₀, and then calculate turbulence intensity; Or by U.S. army mark MIL-8785C acquisition turbulence intensity; Said flying height is divided into height smaller or equal between 1000ft, the 1000ft-2000ft, more than or equal to 2000ft,

As flying height h≤1000ft:

$\left\{\begin{array}{c}{\mathrm{\&sigma;}}_w=0.1{\mathrm{<figure-callout\; id="20"\; label="W"\; filenames="HSA00000110114400012.png,HSA00000110114400023.png"\; state="\{\{state\}\}">W</figure-callout>}}_{20}\\ \frac{{\mathrm{\&sigma;}}_u}{{\mathrm{\&sigma;}}_w}=\frac{{\mathrm{\&sigma;}}_v}{{\mathrm{\&sigma;}}_w}=\frac{1}{{(0.177+0.000823h)}^{0.4}}\end{array}\right.$

σ wherein _u, σ _v, σ _wBe respectively the turbulence intensity of the aircraft fuselage longitudinal axis, side direction, vertical direction, u wherein, v, w be aircraft along x, y, the speed component of z axle, the ft/s of unit;

As flying height h >=2000ft:

The turbulence intensity of aircraft fuselage side direction, vertical direction is found from the icon that MIL-8785C provides, and this icon provides the function of the probability that turbulence intensity surmounted as height and this turbulence intensity, and the relation of turbulence intensity is: σ _v=σ _w

As flying height 1000ft＜h＜2000ft:

Turbulence intensity σ _v, σ _wFunction σ as height _u(h) with the σ at 1000ft place and 2000ft place _v, σ _wLinear interpolation obtains;

With the turbulence intensity on the u direction is example, is shown below:

$\frac{{\mathrm{\&sigma;}}_u\left(h\right)-{\mathrm{\&sigma;}}_u\left(1000\right)}{{\mathrm{\&sigma;}}_u\left(2000\right)-{\mathrm{\&sigma;}}_u\left(1000\right)}=\frac{h-1000}{2000-1000}$

σ wherein _u(1000), σ _u(2000) represent turbulence intensity on the u direction at 1000ft and 2000ft place respectively, h representes flying height.

(3) calculate the power spectrum density of input disturbance signal according to the side direction turbulence intensity of step (2);

${\mathrm{\&Phi;}}_v\left(\mathrm{\&Omega;}\right)={\mathrm{\&sigma;}}_v^2\frac{L_v}{\mathrm{\&pi;}}\frac{1+3{\left(L_v\mathrm{\&Omega;}\right)}^2}{{[1+{\left(L_v\mathrm{\&Omega;}\right)}^2]}^2}$

ω＝ΩV

Φ wherein _u, Φ _v, Φ _wBe u, v, the power spectrum density of w, the ft of unit ³/ s ²L _u, L _v, L _wBe along u, v, the space scale of the power spectrum of w, the ft of unit, V are aircraft airspeed, Ω is a spatial frequency;

(4) G that obtains (s) and the step (3) according to step (1) obtains power spectrum density, and calculation side is to flight path control system Flight technical error FTE and FTE valuation border;

$E\left\{y^{T}y\right\}\&le;\frac{1}{2\mathrm{\&pi;}}{\left[\underset{\mathrm{\&omega;}\&Element;B_d}{\sup }\right\{\stackrel{\&OverBar;}{\mathrm{\&sigma;}}\left[G\left(s\right)\right]\left\}\right]}^2\underset{i}{\mathrm{\&Sigma;}}\underset{\mathrm{\&omega;}\&Element;B_d}{\&Integral;}{\mathrm{\&Phi;}}_{u_i{u}_i}\left(\mathrm{\&omega;}\right)---\left(1\right)$

$E\left\{y^{T}y\right\}\&GreaterEqual;\frac{1}{2\mathrm{\&pi;}}{\left[\underset{\mathrm{\&omega;}\&Element;B_d}{\sup }\right\{\underset{\&OverBar;}{\mathrm{\&sigma;}}\left[G\left(s\right)\right]\left\}\right]}^2\underset{i}{\mathrm{\&Sigma;}}\underset{\mathrm{\&omega;}\&Element;B_d}{\&Integral;}{\mathrm{\&Phi;}}_{u_i{u}_i}\left(\mathrm{\&omega;}\right)---\left(2\right)$

Wherein σ[G (s)], Be side direction flight path control closed-loop system master gain, promptly

Be the maximum singular value of G (s), σ[G (s)] is the minimum singular value of G (s), B _dBe main disturbing signal (being that amplitude is 10 times of other signals or above disturbing signal in the input vector) spectral bandwidth 2-3 doubly,

It is the power spectrum density of i component among the side direction flight path control closed-loop system input vector u;

Above-mentioned flight path control system Flight technical error FTE comes down to a stochastic process of obeying the zero-mean Gaussian distribution, thereby has only the variance on the statistical significance meaningful, and formula (1), (2) have provided the border of valuation up and down of FTE variance; Wherein the variance of Flight technical error is represented on the inequality left side, though because the Flight technical error one-component among the output vector y just, therefore 10 times of its magnitudes, can use E{y to other component (referring to Fig. 4) ^TY} is similar to σ ²(FTE).In addition, because variance and system's singular value represent with σ all that traditionally in order not cause confusion, the formula left side is without σ ²(FTE);

(5) if the turbulence intensity data in the step (2) can't obtain, then need calculate light, in, weigh three kinds of FTE values under the different atmospheric turbulence intensities, calculate the FTE expectation value again; Said light strength of vortex is smaller or equal to 15knots, and middle strength of vortex is greater than 15, and less than 45knots, heavy strength of vortex is that formula is shown below more than or equal to 45knots:

$E\left({\mathrm{\&sigma;}}_{\mathrm{FTE}}^2\right)=P_l\&times;{\mathrm{\&sigma;}}_l^2\left(d\right)+P_m\&times;{\mathrm{\&sigma;}}_m^2\left(d\right)+P_s\&times;{\mathrm{\&sigma;}}_s^2\left(d\right)$

Wherein,

The expectation value of expression FTE variance, P _l, P _m, P _sBe respectively light, in, weigh the probability that the turbulent flow of intensity in three occurs, it calculates with MIL-8785C is basis, concrete numerical value is seen table-2.

Principle of the present invention and derivation:

The present invention is based on lineary system theory, adopts the analytical approach of theory of random processes, is more specifically: covariance analysis and vector power spectral density analytical approach.Facing central principle and derivation down concisely sets forth.

The input vector of side direction flight path control system and the power spectrum density of output vector are shown below:

Wherein

The expression Fourier transform, E{u (t) u ^T(t+ τ) } and E{y (t) y ^T(t+ τ) } be respectively the covariance matrix of input vector and output vector, G (s) is a side direction flight path control system closed loop transfer function, matrix, then knows have following formula to set up according to theory of random processes:

Φ _yy(ω)＝G(jω)Φ _uu(ω)G ^T(-jω)

Notice that again the autocovariance of output vector is obtained by following formula:

$E\left\{y^{T}y\right\}=\frac{1}{2\mathrm{\&pi;}}{\&Integral;}_{-\&infin;}^{\&infin;}\mathrm{tr}\left[{\mathrm{\&Phi;}}_{\mathrm{yy}}\left(\mathrm{\&omega;}\right)\right]\mathrm{d\&omega;}$

Tr [Φ wherein _Yy(ω)] be Φ _YyMark (ω).If G (s) is stable, then have in addition:

$E\left\{y^{T}y\right\}=\frac{1}{2\mathrm{\&pi;}}{\&Integral;}_{-\&infin;}^{\&infin;}\mathrm{tr}\left[{\mathrm{\&Phi;}}_{\mathrm{yy}}\left(\mathrm{\&omega;}\right)\right]\mathrm{d\&omega;}$

$=\frac{1}{2\mathrm{\&pi;}}{\&Integral;}_{-\&infin;}^{\&infin;}\mathrm{tr}\left[G\left(\mathrm{j\&omega;}\right){\mathrm{\&Phi;}}_{\mathrm{uu}}\left(\mathrm{\&omega;}\right)G^T(-\mathrm{j\&omega;})\right]\mathrm{d\&omega;}$

$=\frac{1}{2\mathrm{\&pi;}}{\&Integral;}_{-\&infin;}^{\&infin;}\underset{i}{\mathrm{\&Sigma;}}{\mathrm{\&sigma;}}_i^2\left({\mathrm{\&Phi;}}_{\mathrm{uu}}^{1/2}\left(\mathrm{\&omega;}\right)G\left(\mathrm{j\&omega;}\right)\right)\mathrm{d\&omega;}$

σ in the following formula _iRepresent i singular value,

Notice:

${\underset{\&OverBar;}{\mathrm{\&sigma;}}}^2\left(G\left(\mathrm{j\&omega;}\right)\right)\&le;\frac{\mathrm{tr}\left({\mathrm{\&Phi;}}_{\mathrm{yy}}\left(\mathrm{\&omega;}\right)\right)}{\mathrm{tr}\left({\mathrm{\&Phi;}}_{\mathrm{uu}}\left(\mathrm{\&omega;}\right)\right)}\&le;{\stackrel{\&OverBar;}{\mathrm{\&sigma;}}}^2\left(G\left(\mathrm{j\&omega;}\right)\right)$

Wherein σ(G (j ω)) and

Be the master gain of G (j ω) (that is: be respectively G (j ω) minimum, maximum singular value), then further can obtain:

$E\left\{y^{T}y\right\}\&le;\frac{1}{2\mathrm{\&pi;}}{\&Integral;}_{-\&infin;}^{\&infin;}{\stackrel{\&OverBar;}{\mathrm{\&sigma;}}}^2\left(G\left(\mathrm{j\&omega;}\right)\right)\mathrm{tr}\left[{\mathrm{\&Phi;}}_{\mathrm{uu}}\left(\mathrm{\&omega;}\right)\right]\mathrm{d\&omega;}$

$E\left\{y^{T}y\right\}\&GreaterEqual;\frac{1}{2\mathrm{\&pi;}}{\&Integral;}_{-\&infin;}^{\&infin;}{\underset{\&OverBar;}{\mathrm{\&sigma;}}}^2\left(G\left(\mathrm{j\&omega;}\right)\right)\mathrm{tr}\left[{\mathrm{\&Phi;}}_{\mathrm{uu}}\left(\mathrm{\&omega;}\right)\right]\mathrm{d\&omega;}$

Consider that again the concentration of energy of the power spectrum density overwhelming majority of turbulent perturbation, is derived as follows so can top border estimate equation further be done referring to Fig. 3 in a more limited frequency domain scope:

$E\left\{y^{T}y\right\}\&le;\frac{1}{2\mathrm{\&pi;}}{\&Integral;}_{-\&infin;}^{\&infin;}{\stackrel{\&OverBar;}{\mathrm{\&sigma;}}}^2\left(G\left(\mathrm{j\&omega;}\right)\right)\mathrm{tr}\left[{\mathrm{\&Phi;}}_{\mathrm{uu}}\left(\mathrm{\&omega;}\right)\right]\mathrm{d\&omega;}$

$=\frac{1}{2\mathrm{\&pi;}}{\&Integral;}_{-\&infin;}^{\&infin;}{\stackrel{\&OverBar;}{\mathrm{\&sigma;}}}^2\left(G\left(\mathrm{j\&omega;}\right)\right)\underset{i}{\mathrm{\&Sigma;}}{\mathrm{\&Phi;}}_{u_i{u}_i}\left(\mathrm{\&omega;}\right)\mathrm{d\&omega;}$

$\&ap;\frac{1}{2\mathrm{\&pi;}}{\left[\underset{\mathrm{\&omega;}\&Element;B_d}{\sup }\right\{\stackrel{\&OverBar;}{\mathrm{\&sigma;}}\left[(G\left(s\right)\right]\left\}\right]}^2\underset{i}{\mathrm{\&Sigma;}}\underset{\mathrm{\&omega;}\&Element;B_d}{\&Integral;}{\mathrm{\&Phi;}}_{u_i{u}_i}\left(\mathrm{\&omega;}\right)$

$E\left\{y^{T}y\right\}\&GreaterEqual;\frac{1}{2\mathrm{\&pi;}}{\&Integral;}_{-\&infin;}^{\&infin;}{\underset{\&OverBar;}{\mathrm{\&sigma;}}}^2\left(G\left(\mathrm{j\&omega;}\right)\right)\mathrm{tr}\left[{\mathrm{\&Phi;}}_{\mathrm{uu}}\left(\mathrm{\&omega;}\right)\right]\mathrm{d\&omega;}$

$=\frac{1}{2\mathrm{\&pi;}}{\&Integral;}_{-\&infin;}^{\&infin;}{\underset{\&OverBar;}{\mathrm{\&sigma;}}}^2\left(G\left(\mathrm{j\&omega;}\right)\right)\underset{i}{\mathrm{\&Sigma;}}{\mathrm{\&Phi;}}_{u_i{u}_i}\left(\mathrm{\&omega;}\right)\mathrm{d\&omega;}$

$\&ap;\frac{1}{2\mathrm{\&pi;}}{\left[\underset{\mathrm{\&omega;}\&Element;B_d}{\sup }\right\{\underset{\&OverBar;}{\mathrm{\&sigma;}}\left[G\left(s\right)\right]\left\}\right]}^2\underset{i}{\mathrm{\&Sigma;}}\underset{\mathrm{\&omega;}\&Element;B_d}{\&Integral;}{\mathrm{\&Phi;}}_{u_i{u}_i}\left(\mathrm{\&omega;}\right)$

Wherein

Expression G (s) is belonging to B _dFrequency range on the supremum of minimum singular value, subscript u _iI component of signal in the expression input vector.

The present invention's advantage compared with prior art is following:

(1) utilizes master gain and turbulence power spectrum concentration of energy zone to estimate the FTE border, under the prerequisite of assurance, simplify effectively and calculate based on the conservative property of security consideration.

(2) utilize the closed loop transfer function, matrix of crabbing flight path control system or its state space to realize, embody that FTE receives from rudders pneumatic power parameter, flying quality and automatic flight control system combined influence; Thereby realized numerous FTE source parameters is effectively embodied.

(3) the present invention is the first systems approach that can carry out the FTE prediction; Because the proposition of method is based on thorough, the deep analysis that FTE is generated physical mechanism; And the establishment of its border estimation formulas pushes away card based on the theory of strictness, makes the correctness of method and validity that abundant assurance arranged.

Description of drawings

Fig. 1 advances closely for side direction flight path control system of the present invention/air route side direction FTE algorithm flow chart;

Fig. 2 is the master gain curve of aircraft side direction flight path control system closed loop transfer function;

Fig. 3 is an atmospheric turbulence disturbance spectral density;

Fig. 4 is the aircraft side direction flight path control system of having considered under the true atmospheric disturbance, to the response curve of the turbulent perturbation of three kinds of (light, in, heavy) intensity;

The icon that Fig. 5 provides for American army mark MIL-8785C.

Embodiment

Present embodiment of the present invention has adopted the lateral linear gasification movable model of a large transport airplane, and based on analyzing with the side direction flight path control system of GLQG/LTR robust control system designing method.Since the approach flight stage be in all legs to the highest part of requirement such as security and error precision, present embodiment is in the state of flight of laggard nearly leg.Flying height is 900ft, and air speed is 229.67ft/s.Correspond respectively to light, in, weigh three types of atmospheric turbulence intensity W ₂₀Be respectively 15knots, 30knots or 45knots.

Step 1:

(1) obtains the side direction automatic flight control system closed loop transfer function, matrix G (s) of the type that is directed against, be shown below: G (s)=C (sI-A) ^-1B+D

Or realize in a minimum state space that obtains G (s); Shape is shown below (general piecemeal notation): wherein A, B, C, D are four constant matricess, when the closed loop configuration of the side direction flight path control system that obtains certain type, have just obtained the information of these four matrixes.

$\left(\begin{array}{cc}A& B\\ C& D\end{array}\right)$

The value of A in the present embodiment, B, C, D is shown in form.

The value of table 1 matrix A (1-10 row)

The value of table 1 (continuing) matrix A (11-20 row)

The value of table 2 matrix B

The value of table 3 Matrix C

The value of table 4 matrix D

Step 2:

Meteorological department of authorities obtains the crosswind section from the airport, and is extracted in the mean wind speed W of 20ft eminence ₂₀

A) if flying height is lower than and equals 1000ft, then obtain σ by (1) formula _v:

$\left\{\begin{array}{c}{\mathrm{\&sigma;}}_w=0.1{\mathrm{<figure-callout\; id="20"\; label="W"\; filenames="HSA00000110114400012.png,HSA00000110114400023.png"\; state="\{\{state\}\}">W</figure-callout>}}_{20}\\ \frac{{\mathrm{\&sigma;}}_v}{{\mathrm{\&sigma;}}_w}=\frac{1}{{(0.177+0.000823h)}^{0.4}}\end{array}\right.---\left(1\right)$

If flying height is higher than and equals 2000ft, then find σ by the probability tables that surmounts among the MIL-8785C (as shown in Figure 5) _vValue; When flying height is in (1000,2000), then through σ to 1000ft and 2000ft place _vCarry out linear interpolation and obtain the turbulence intensity standard deviation of desired height, wherein σ _u, σ _v, σ _wBe respectively the turbulence intensity of the aircraft fuselage longitudinal axis, side direction, vertical direction.

B) calculate turbulent perturbation yardstick information according to aircraft place flying height,, then obtain L by (2) formula if flying height is lower than 1000ft _vIf flying height is higher than then L of 2000ft _v=1750ft;

$\left\{\begin{array}{c}{L}_w=h\\ L_v=\frac{h}{{(0.177+0.000823h)}^{1.2}}\end{array}\right.---\left(2\right)$

When flying height is in (1000,2000), then through L to 1000ft and 2000ft place _vCarry out linear interpolation and obtain the scale of turbulence information of desired height.

Step 3:

A) obtain aircraft airspeed value V, can obtain the analytical function of Dryden turbulent flow pulsation frequency spectrum by (3), (4) formula.

ω＝ΩV (3)

${\mathrm{\&Phi;}}_v\left(\mathrm{\&Omega;}\right)={\mathrm{\&sigma;}}_v^2\frac{L_v}{\mathrm{\&pi;}}\frac{1+3{\left(L_v\mathrm{\&Omega;}\right)}^2}{{[1+{\left(L_v\mathrm{\&Omega;}\right)}^2]}^2}---\left(4\right)$

Φ wherein _v, Φ _wBe v, the power spectrum density of w, the ft of unit ³/ s ²L _v, L _wBe along v, the space scale of the power spectrum of w, the ft of unit; V is an aircraft airspeed, and Ω is a spatial frequency, and the Dryden turbulent flow pulsation frequency spectrum under the moderate turbulent perturbation is shown in figure-3; This frequency spectrum is an even function, because of its curve about longitudinal axis symmetry, so only need to represent monolateral spectrum curve.

B) formed filter (Forming Filter) is σ by standard deviation shown in (5) formula _WnWhite noise drive to generate turbulent perturbation.(6) in, V is an aircraft airspeed, Dt be sample time of presetting at interval, reduce with the increase of aircraft speed.

$F_v\left(s\right)=\frac{1}{1+L_{v}s}---\left(5\right)$

Dx＝V·Dt (6)

σ _wn＝σ _v(2L _v/Dx) ^1/2 (7)

Step 4:

A) calculate input signal vector u=[δ according to (8) formula _aδ _rn _w] spectral density.

The upper and lower border of covariance of b) calculating the output signal vector by (9) formula, (10) formula respectively.

$E\left\{y^{T}y\right\}\&le;\frac{1}{2\mathrm{\&pi;}}{\left[\underset{\mathrm{\&omega;}\&Element;B_d}{\sup }\right\{\stackrel{\&OverBar;}{\mathrm{\&sigma;}}\left[G\left(s\right)\right]\left\}\right]}^2\underset{i}{\mathrm{\&Sigma;}}\underset{\mathrm{\&omega;}\&Element;B_d}{\&Integral;}{\mathrm{\&Phi;}}_{u_i{u}_i}\left(\mathrm{\&omega;}\right)---\left(9\right)$

$E\left\{y^{T}y\right\}\&GreaterEqual;\frac{1}{2\mathrm{\&pi;}}{\left[\underset{\mathrm{\&omega;}\&Element;B_d}{\sup }\right\{\underset{\&OverBar;}{\mathrm{\&sigma;}}\left[G\left(s\right)\right]\left\}\right]}^2\underset{i}{\mathrm{\&Sigma;}}\underset{\mathrm{\&omega;}\&Element;B_d}{\&Integral;}{\mathrm{\&Phi;}}_{u_i{u}_i}\left(\mathrm{\&omega;}\right)---\left(10\right)$

Wherein, σ[G (s)],

Be side direction flight path control closed-loop system master gain, i.e. minimum, the maximum singular value of G (s), master gain is the function of frequency, and is as shown in Figure 2, B _dBe main disturbing signal spectral bandwidth 2-3 doubly.Be similar to the analysis in the step 5, formula (9), (10) are the upper and lower boundary algorithm for estimating of side direction FTE.

Control system Flight technical error FTE comes down to a stochastic process of obeying the zero-mean Gaussian distribution, thereby has only the variance on the statistical significance meaningful, below (1), (2) formula provided the border of valuation up and down of FTE variance; Wherein the variance of Flight technical error is represented on the inequality left side, though because the Flight technical error one-component among the output vector y just, therefore 10 times of its magnitudes, can use E{y to other component (referring to Fig. 4) ^TY} is similar to σ ²(FTE).In addition, because variance and system's singular value represent with σ all that traditionally in order not cause confusion, the formula left side is without σ ²(FTE).

The side direction flight path system dynamic system that has merged turbulent perturbation spectral shaping wave filter is driven by the white noise with corresponding turbulence intensity; Can obtain the FTE variance under the varying strength according to formula (9), (10) prediction; Once as shown in Figure 4 based on the emulation of True Data; Shown among the figure light, in, it is heavy that (light strength of vortex is for smaller or equal to 15knots, and middle strength of vortex is greater than 15, less than 45knots; Heavy strength of vortex is more than or equal to 45knots) under three kinds of turbulent perturbations, in the end advance the Flight technical error curve of nearly leg aircraft.Because the essence of FTE is stochastic process, therefore each curve is respectively the once realization under the corresponding turbulent perturbation intensity among the figure.

Step 5:

If can't obtain crosswind section reliably, perhaps need predict side direction FTE generally, then need to calculate the expectation value of FTE according to (11) formula according to the probability that surmounts among the MIL-8785C.

$E\left({\mathrm{\&sigma;}}_{\mathrm{FTE}}^2\right)=P_l\&times;{\mathrm{\&sigma;}}_l^2\left(d\right)+P_m\&times;{\mathrm{\&sigma;}}_m^2\left(d\right)+P_s\&times;{\mathrm{\&sigma;}}_s^2\left(d\right)---\left(11\right)$

Present embodiment is to listed target type in the table-1; Suppose to obtain to be concerned about regional 20ft eminence mean wind speed; Then must calculate three kinds of turbulent perturbation frequency spectrums under the varying strength; And because of it plays the FTE value that disturbance produces, and according to the expectation value that Probability p robability calculates the FTE variance that surmounts among the MIL-8785C, shown in 2 tables.

Table is-2 light, in, severe (Light; Moderate, severe) the FTE variance upper bound under the turbulent perturbation

$E\left({\mathrm{\&sigma;}}_{\mathrm{FTE}}^2\right)=P_l\&times;{\mathrm{\&sigma;}}_l^2\left(d\right)+P_m\&times;{\mathrm{\&sigma;}}_m^2\left(d\right)+P_s\&times;{\mathrm{\&sigma;}}_s^2\left(d\right)---\left(12\right)$

$=6.7749\mathrm{<figure-callout\; id="1"\; label="e"\; filenames="HSA00000110114400021.png"\; state="\{\{state\}\}">e</figure-callout>}1$

D representes FTE value, P in the formula (12) _l, P _m, P _s, respectively expression light, in, the probability that occurs of severe turbulent perturbation.

Step 6:

Write code with airborne equipment (being generally Flight Management Computer FMC) or uphole equipment (being generally RNAV/RNP navigation prediction platform) and carry out FTE Boundary Prediction method; Can obtain one of major part among the TSE (another major part is NSE); With predict NSE value (the present invention only is used for the prediction of FTE) addition can be to ANP (Actual NavigationPerformance; Actual navigation performance) makes the preceding prediction of short-term or boat; And only learn whether ANP meets the RNP index of appointment, could implement the PBN navigation.Therefore the present invention is one of key issue in the prediction of PBN navigation performance.

Consider the sky high cost of practical flight; Adopt in the embodiment of the invention and FTE border method of estimation is verified based on the emulation of Live Flying condition and flight control system parameter; The figure that characterizes two key factors (flight control system, turbulent flow frequency spectrum) characteristic in the simulation process is Fig. 2 (the master gain curve of aircraft side direction flight path control system closed loop transfer function) and Fig. 3 (atmospheric turbulence disturbance spectral density).Simulation result is shown in Figure 4, has promptly considered the aircraft side direction flight path control system under the true atmospheric disturbance, to the response curve of the turbulent perturbation of three kinds of (light, in, heavy) intensity.D among Fig. 4 representes side direction FTE value; ψ representes course angle; These two variablees have constituted two component of signals in the output vector; Simulation result by Fig. 4 is visible, and the magnitude of side direction FTE is more than 10 times of another course angle magnitude, so simulation result has verified that vectorial covariance analysis method is approximately vectorial variance the correctness of the variance of side direction FTE.

The present invention does not set forth the known technology that part belongs to those skilled in the art in detail.

Claims (1)

1. Flight technical error Forecasting Methodology that is used for the side direction flight path control system is characterized in that performing step is following:

(1) obtains the side direction automatic flight control system closed loop transfer function, matrix G (s) of the aircraft type that is directed against, be shown below: G (s)=C (sI-A) ^-1B+D, wherein A, B, C, D are four constant matricess, and s is the Laplace domain variable, and I is and the unit matrix of matrix A with line number same column number; Or the minimum state space realization of acquisition G (s), promptly $\left(\begin{array}{cc}A& B\\ C& D\end{array}\right),$ Wherein A, B, C, D are four constant matricess;

(2) obtain crosswind section and the wind speed W20 of 20ft place according to the flying height difference of being concerned about the zone or through relevant meteorological department, and then calculate turbulence intensity; Or by U.S. army mark MIL-8785C acquisition turbulence intensity; Said flying height is divided into height smaller or equal between 1000ft, the 1000ft-2000ft, more than or equal to 2000ft,

When flying height h≤1000ft:

$\left\{\begin{array}{c}{\mathrm{\&sigma;}}_w=0.1W_{20}\\ \frac{{\mathrm{\&sigma;}}_v}{{\mathrm{\&sigma;}}_w}=\frac{1}{{(0.177+0.000823h)}^{0.4}}\end{array}\right.$

σ wherein _v, σ _wBe respectively the turbulence intensity of aircraft fuselage side direction, vertical direction, v wherein, w be aircraft along y, the speed component of z axle, the ft/s of unit;

When flying height h >=2000ft:

Turbulence intensity along aircraft fuselage side direction, vertical direction finds from the icon that MIL-8785C provides, and this icon provides the function of the probability that turbulence intensity surmounted as height and this turbulence intensity, and the turbulence intensity that relates to relation is: σ _v=σ _w

When flying height 1000ft＜h＜2000ft:

Turbulence intensity σ _v, σ _wFunction σ as height _v(h), σ _w(h) with the σ at 1000ft place and 2000ft place _v, σ _wLinear interpolation obtains;

(3) calculate the input disturbance signal respectively according to the side direction turbulence intensity of step (2), i.e. the power spectrum density of turbulent flow;

${\mathrm{\&Phi;}}_v\left(\mathrm{\&Omega;}\right)={\mathrm{\&sigma;}}_v^2\frac{L_v}{\mathrm{\&pi;}}\frac{1+3{\left(L_v\mathrm{\&Omega;}\right)}^2}{{[1+{\left(L_v\mathrm{\&Omega;}\right)}^2]}^2}$

ω＝ΩV

Φ wherein _v(ω) be the power spectrum density of v, with ω=Ω V substitution Φ _v(Ω) can obtain the ft of unit ³/ s ²Lv is the space scale along the power spectrum of v, and the ft of unit, V are aircraft airspeed, and Ω is a spatial frequency;

(4) G (s) and the step (3) that obtain according to step (1) obtain power spectrum density, and calculation side is to flight path control system FTE and FTE valuation border,

$E\left\{y^{T}y\right\}\&le;\frac{1}{2\mathrm{\&pi;}}{\left[\underset{\mathrm{\&omega;}\&Element;B_d}{\sup }\right\{\stackrel{\&OverBar;}{\mathrm{\&sigma;}}\left[G\left(s\right)\right]\left\}\right]}^2\underset{i}{\mathrm{\&Sigma;}}\underset{\mathrm{\&omega;}\&Element;B_d}{\&Integral;}{\mathrm{\&Phi;}}_{u_i{u}_i}\left(\mathrm{\&omega;}\right)---\left(1\right)$

$E\left\{y^{T}y\right\}\&GreaterEqual;\frac{1}{2\mathrm{\&pi;}}{\left[\underset{\mathrm{\&omega;}\&Element;B_d}{\sup }\right\{\underset{\&OverBar;}{\mathrm{\&sigma;}}\left[G\left(s\right)\right]\left\}\right]}^2\underset{i}{\mathrm{\&Sigma;}}\underset{\mathrm{\&omega;}\&Element;B_d}{\&Integral;}{\mathrm{\&Phi;}}_{u_i{u}_i}\left(\mathrm{\&omega;}\right)---\left(2\right)$

Wherein σ[G (s)],

Be side direction flight path control closed-loop system master gain, promptly Be the maximum singular value of G (s), σ[G (s)] is the minimum singular value of G (s), and Bd is main disturbing signal, and the 2-3 of spectral bandwidth times, said main disturbing signal is a turbulent perturbation, is that amplitude is 10 times of other signals or above disturbing signal in the input vector;

Be the power spectrum density of i component among the side direction flight path control closed-loop system input vector u, owing to only consider main disturbance, among the present invention

Be turbulent flow power spectrum density Φ _v(ω);

(5) if the turbulence intensity data in the step (2) can't obtain, then need calculate light, in, weigh three kinds of FTE values under the different atmospheric turbulence intensities, calculate the FTE expectation value again; Said light strength of vortex is smaller or equal to 15knots, and middle strength of vortex is greater than 15, and less than 45knots, heavy strength of vortex is that formula is shown below more than or equal to 45knots:

$E\left({\mathrm{\&sigma;}}_{\mathrm{FTE}}^2\right)=P_l\&times;{\mathrm{\&sigma;}}_l^2\left(d\right)+P_m\&times;{\mathrm{\&sigma;}}_m^2\left(d\right)+P_s\&times;{\mathrm{\&sigma;}}_s^2\left(d\right)$

Wherein, The expectation value of expression FTE variance, P _l, P _m, P _sBe respectively light, in, weigh the probability that the turbulent flow of intensity in three occurs.

Images