METHOD AND SYSTEM FOR ESTIMATING AERODYNAMIC ANGLES OF A FLYING BODY

Abstract

A computer implemented method for estimating angles of attack and/or an angle of sideslip of a flying body includes the steps of: providing data or measurements representative of a time derivative of a module of a speed of the flying body with respect to an air; providing data or measurements representative of corresponding projections of coordinate accelerations on a Cartesian reference frame fixed to the flying body; and calculating the angles of attack and/or the sideslip on a basis of an mathematical relationship of the angles of attack and the sideslip with the coordinate accelerations and the time derivative of the module of the speed of the flying body with respect to the air.

Claims

1. A method for estimating an angle α of attack and/or an angle β of sideslip RPM of a flying body comprising the following steps of: providing data or measurements representative of a time derivative of a module of a speed {dot over (V)} of the flying body with respect to an air; providing data or measurements representative of corresponding projections of coordinate accelerations (a.sub.x, a.sub.y, a.sub.z) on a Cartesian reference frame fixed to the flying body; and calculating the angle α of attack and/or the angle β of sideslip on a basis of a mathematical relationship of the angle α of attack and the angle β of sideslip with the coordinate accelerations and the derivative of the module of the speed of the flying body with respect to the air.

2. The method according to claim 1, wherein the step of calculating is not based on a model of six-degree-of-freedom dynamics, of a geometry and/or a mass and/or aerodynamics and/or propulsion forces and moments comprising buoyancy forces and/or flight control commands, comprising models in a linear or non-linear state space form, of the flying body.

3. The method according to claim 2, further comprising the step of: providing measurements representative of the corresponding projections on the Cartesian reference frame fixed to the flying body of a wind acceleration {dot over (w)}.sub.B; and the mathematical relationship also comprises the wind acceleration.

4. The method according to claim 3, wherein, after linearization, a mathematical relationship solved for the angle α.sub.lin of attack is as follows: ${\alpha }_{\mathrm{lin}}=\frac{{\dot{V}}_{\infty }-a_X\cos \beta -a_Y\sin \beta }{a_Z\cos \beta }.$ and a mathematical relationship solved for the angle β.sub.lin of sideslip is as follows: $.Math.{\beta }_{\mathrm{lin}}=\frac{\stackrel{.}{V}-\left(a_X\cos \alpha \right)-\left(a_Z\sin \alpha \right)}{a_Y}.$

5. The method according to claim 4, comprising the step of providing the data or the measurements representative of the corresponding projections on the Cartesian reference frame fixed to the flying body of angular velocities of the flying body as follows: (Ω=f(w)=f p, q, r)); and solving a system of equations comprising the mathematical relationship and a further mathematical formula expanded in a time domain expressed via the coordinate accelerations, the angle of attack and the angle of sideslip, the speed V.sub.∞ and the time derivative {dot over (V)}.sub.∞ of the speed with respect to the air, and the angular velocities; wherein unknown parameters of the system having mathematical formulas are the angles of attack and sideslip at a time t, and remaining parameters are either at the time t or at a time t-Δt and the remaining parameters are either measurable, when at the time t, or already calculated and therefore known, when referred to the time t-Δt.

6. The method according to claim 5, further comprising the following steps: measuring directly or indirectly the angle of attack and/or the angle of sideslip by a measuring device; calculating a solution of the mathematical relationship; comparing a measured angle with the solution; and generating a warning signal when a difference between the measured angle and the solution is compatible with a predetermined rule.

7. The method according to claim 1, wherein the steps of providing comprises a step of detecting a relative direct or indirect measurement of the time derivative of the speed with respect to the air and of the coordinate accelerations by a respective on-board measuring device of the flying body and the step of calculating is performed on a basis of the step of measuring either in real time or based on a step of storing measurement data following the step of measuring.

8. The method according to claim 7, wherein a device for detecting directly or indirectly the coordinate accelerations comprises a sensor unit and the sensor unit comprises exclusively an inertial measurement device having a plurality of accelerometers and gyroscopes.

9. The method according to claim 8, wherein a device for detecting directly or indirectly the time derivative of the module of the speed with respect to the air comprises the sensor unit and the sensor unit exclusively comprises a dynamic pressure sensor.

10. The method according to claim 9, wherein a device for detecting directly or indirectly the time derivative of the module of the speed with respect to the air comprises the sensor unit and the sensor unit exclusively comprises the dynamic pressure sensor and an air temperature sensor.

11. The method according to claim 1, wherein the angle of attack and/or the angle of sideslip are not calculated based on a vertical speed and/or inertial acceleration as such of the flying body.

12. The method according to claim 1, comprising the step of showing the angle of attack and/or the angle of sideslip on a display.

13. A computer implemented method for a flight simulation comprising the method according to claim 12.

14. A flying body comprising an avionic system, wherein the avionic system is programmed to implement the computer implemented method according to claim 13.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

[0037] The invention is described below on the basis of non-limiting examples illustrated by way of example in the accompanying figures, which refer respectively to:

[0038] FIG. 1: Cartesian reference frames useful for the mathematical definitions of parameters used in the method according to the invention;

[0039] FIG. 2: a wind variation diagram;

[0040] FIG. 3A: graphs showing the overlap between the theoretical results (solid line) and the numerical results according to the present invention (circles) for a first solution form;

[0041] FIG. 3B: graphs showing the overlap between the theoretical results (solid line) and the numerical results according to the present invention (circles) for a second solution form;

[0042] FIG. 4A: graphs showing the overlap between the theoretical results (solid line) and the numerical results according to the present invention (circles) for a first solution form; and

[0043] FIG. 4B: graphs showing the overlap between the theoretical results (solid line) and the numerical results according to the present invention (circles) for a second solution form.

DETAILED DESCRIPTION OF THE EMBODIMENTS

[0044] The present invention considers the following definitions, in which a bold lowercase letter refers to a vector quantity, a bold capital letter refers to a matrix quantity and a lowercase or capital letter to a scalar:

[0045] inertial velocity vector in the inertial reference frame: v.sub.I=[u.sub.i, v.sub.i, w.sub.i].sup.T

[0046] velocity vector with respect to the wind, i.e. surrounding air, in the ‘body’ reference frame: v.sub.B=[u, v, w].sup.T

[0047] wind velocity vector in the inertial reference frame: w.sub.I=[u.sub.wi, v.sub.wi, w.sub.wi].sup.T

[0048] With reference to the Euler angles, ψ, θ, φ, the following rotation matrices are defined to switch from an inertial reference system to a Body reference frame on board the aircraft:

[00003]$C_{\psi }=\left[\begin{array}{ccc}\cos \left(\psi \right)& \sin \left(\psi \right)& 0\\ -\sin \left(\psi \right)& \cos \left(\psi \right)& 0\\ 0& 0& 1\end{array}\right],C_{\theta }=\left[\begin{array}{ccc}\cos \left(\theta \right)& 0& -\sin \left(\theta \right)\\ 0& 1& 0\\ \sin \left(\theta \right)& 0& \cos \left(\theta \right)\end{array}\right],C_{\varphi }=\left[\begin{array}{ccc}1& 0& 0\\ 0& \cos \left(\varphi \right)& \sin \left(\varphi \right)\\ 0& -\sin \left(\varphi \right)& \cos \left(\varphi \right)\end{array}\right]$

[0049] and the following known relation:

[00004]$C_{I2B}=\text{⁠}C_{\varphi }{C}_{\theta }{C}_{\psi }==\left[\text{⁠}\begin{array}{ccc}\cos \theta \cos \psi & \cos \theta \sin \psi & -\sin \theta \\ \sin \mathrm{\varphi sin}\mathrm{\theta cos}\psi -\cos \mathrm{\varphi sin}\psi & \sin \mathrm{\varphi sin}\mathrm{\theta sin\psi }+\cos \mathrm{\varphi cos}\psi & \sin \varphi \cos \theta \\ \cos \mathrm{\varphi sin}\theta \cos \psi +\sin \mathrm{\varphi sin\psi }& \cos \mathrm{\varphi sin}\mathrm{\theta sin}\psi -\sin \mathrm{\varphi cos}\psi & \cos \varphi \cos \theta \end{array}\right]$

[0050] Referring to FIG. 1, the following definition of the aircraft velocity vector applies, in which the subscript ‘I’ refers to an inertial Cartesian reference frame:

[0051] The wind speed and direction w.sub.I is ideally estimated for example from meteorological data or calculated through an on-board system, e.g. known air data system with vanes and with the speed measured or obtained in an inertial frame, e.g. using data from a satellite global positioning system (GPS).

[0052] In the present description the wind is measured in the inertial Cartesian reference frame. Nothing prohibits to consider the wind in Cartesian body reference frame and therefore using the w.sub.B notation. The transition from one to the other is obtained by means of a rotation matrix that uses the Euler angles. In order not to burden the following discussion, we prefer to consider the wind measured or estimated in inertial axes, w.sub.I.

[0053] Furthermore, the following definitions of acceleration vectors in the Body and Inertial reference frames are valid:

[00005]$a_B={\left(\frac{{\mathrm{dv}}_I}{\mathrm{dt}}\right)}_B=C_{I2B}{\stackrel{.}{v}}_I={\stackrel{.}{v}}_B+{\Omega }_B{v}_B+C_{I2B}{\stackrel{.}{w}}_I\mathrm{and}{a}_I={\left(\frac{{\mathrm{dv}}_I}{\mathrm{dt}}\right)}_I={\left(\frac{{\mathrm{dC}}_{B2I}{v}_B}{\mathrm{dt}}\right)}_I+{\left(\frac{{\mathrm{dw}}_I}{\mathrm{dt}}\right)}_I.Math.{\stackrel{.}{v}}_I=C_{B2I}{\stackrel{.}{v}}_B+{\stackrel{.}{C}}_{B2I}{v}_B+{\stackrel{.}{w}}_I.Math.{\stackrel{.}{v}}_I=C_{B2I}\left({\stackrel{.}{v}}_B+C_{I2B}{\stackrel{.}{C}}_{B2I}{v}_B\right)+{\stackrel{.}{w}}_I.Math.{\stackrel{.}{v}}_I=C_{B2I}\left({\stackrel{.}{v}}_B+{\Omega }_B{v}_B\right)+{\stackrel{.}{w}}_I\mathrm{where}{\Omega }_B=\left[\begin{array}{ccc}0& -r& q\\ r& 0& -p\\ -q& p& 0\end{array}\right]$

[0054] where p, q and r are respectively the components in the body reference frame of the aircraft angular velocity measured for example by the onboard gyroscopes.

[0055] Simplified Solution for Time t

[0056] The following relationship is also valid, introducing the angles of attack, α, and sideslip, β:

v.sub.B=V.sub.∞[cos βcos α,sin β,cos β sin α].sup.T=V.sub.∞î.sub.B

[0057] Considering the equations of speed and acceleration in the body reference frame and

[00006]${\dot{V}}_{\infty }=\frac{{\nu }_B^T{\dot{\nu }}_B}{V_{\infty }},$

i.e. a time derivative of an airspeed, the result is:

[00007]${\dot{V}}_{\infty }{V}_{\infty }=v_B^T{\stackrel{.}{v}}_B=v_B^T\left(a_B-\Omega v_B-C_{I2B}{\stackrel{.}{w}}_I\right)==v_B^T{a}_B--v_B^T{C}_{I2B}{\stackrel{.}{w}}_I==V_{\infty }{\hat{i}}_B^T\left(a_B-C_{I2B}{\stackrel{.}{w}}_I\right).Math..Math.{\stackrel{.}{V}}_{\infty }={\hat{i}}_B^T\left(a_B-C_{I2B}{\stackrel{.}{w}}_I\right).$

[0058] The above equation, for the purpose of a solution for the calculation of the aerodynamic angles included in the term î.sub.B.sup.T, can be further simplified by assuming that in the time interval the acceleration of the wind is zero or by assuming a pre-defined wind speed variation on the basis of estimates of the same aerodynamic angles relative to previous time instants for each considered time interval. Such assumptions are quite common and allow for relatively high precision calculations in a generic flight condition. The hypothesis of negligible wind variation, {dot over (w)}.sub.I≈0, provides exact results in constant (or null) wind conditions and acceptable results in moderately variable wind conditions.

[0059] Following this assumption, {dot over (w)}.sub.I≈0, the following formula valid in a generic time t is obtained:

a.sub.x cos βcos α+a.sub.y sin β+a.sub.z cos β sin α={dot over (V)}.sub.∞  EQ. 1

[0060] which provides an explicit correlation between the angles of attack and sideslip to directly or indirectly measurable quantities, i.e. aircraft time derivative of an airspeed with respect to the air (calculable through a conventional air data system) and the coordinate accelerations expressed in the body reference frame (measured with respect to an inertial reference system, e.g. GPS or calculated by means of an inertial measurement device and the angles of Euler obtained starting from the measurements of the gyroscopes of the inertial measurement device). For the sake of simplicity, the time derivative of the true airspeed is taken into consideration below, since according to their respective assumptions, the time derivative of indicated, measured and equivalent airspeeds respectively obtained from the IAS, CAS and EAS speeds can also be considered. Examples of calculation of the time derivative are:

[00008]${\stackrel{.}{V}}_{\infty }=\frac{V_{\infty ,t}-V_{\infty ,t-\Delta t}}{\Delta t}{\stackrel{.}{V}}_{\infty }=\frac{V_{\infty ,t+\Delta t}-V_{\infty ,t-\Delta t}}{2\Delta t}$

[0061] According to a first embodiment of the present invention, it is possible to calculate one of the two angles α and β when the other is known, for example when the angle of attack is measured with known systems.

[0062] In fact, expressing the sine and cosine through the known parametric formulas.

[00009]$\sin \beta =\frac{2t}{1+t^2}$

[00010]$\cos \beta =\frac{1-t^2}{1+t^2}$

[0063] Where

[00011]$t=\tan \frac{\beta }{2}$

with β≠π+2kπ, k∈

[0064] It is therefore possible to write:

[00012]$\left(a_x\cos \alpha +a_z\sin \alpha \right)\cos \beta +a_y\sin \beta =A\cos \beta +B\sin \beta =\frac{1-t^2}{1+t^2}A+\frac{2t}{1+t^2}B={\dot{V}}_{\infty }=C.Math.\left(1+t^2\right)C=\left(1-t^2\right)A+2Bt.Math.\left(C+A\right)t^2-2Bt+\left(C-A\right)=0$

[0065] Out of the two possible solutions of such second order equations, the most realistic one is selected; if there are two close solutions, the selected one is that minimizing the time derivative with respect to the preceding step.

[0066] The possible analytic solution of angle α from a, {dot over (V)}.sub.∞ e β is similar to what is shown previously.

[0067] For instance, the formula to calculate the angle of attack from the angle of sideslip is:

[00013]$\alpha =2{\tan }^{-1}\frac{a_Z\cos \beta \pm \sqrt{{\left(a_Z\cos \beta \right)}^2-\left[{\left({\stackrel{.}{V}}_{\infty }-a_Y\sin \beta \right)}^2-{\left(a_X\cos \beta \right)}^2\right]}}{{\stackrel{.}{V}}_{\infty }-a_Y\sin \beta +a_X\cos \beta }$

[0068] Furthermore, it is possible to approximate the sine with the angle and the cosine to 1 and

[00014]${\alpha }_{\mathrm{lin}}=\frac{{\stackrel{.}{V}}_{\infty }-a_X\cos \beta -{\alpha }_Y\sin \beta }{a_Z\cos \beta }$

thus to obtain a linearized form:

[0069] FIGS. 3A and 3B compare the numerical results obtained by applying the two formulas, i.e. nonlinear (left) and linearized (right) of the angle of attack shown above, where the angle of sideslip is known from a conventional system of measurement of the latter angle. The comparison data (in continuous line) are obtained from the solution of linearized longitudinal dynamic differential equations of a known aircraft.

[0070] Referring to the angle of sideslip, formulas obtained similarly to those of the angle of attack are:

[00015]$\beta =2{\tan }^{-1}\frac{a_Y\sqrt{a_Y^2-\left({\stackrel{.}{V}}_{\infty }^2-a_X^2{\cos }^2\alpha -a_Z^2{\sin }^2\alpha \right)}}{{\stackrel{.}{V}}_{\infty }+a_X\cos \alpha +a_Z\sin \alpha }\mathrm{and}.Math.{\beta }_{\mathrm{lin}}=\frac{{\stackrel{.}{V}}_{\infty }-\left(a_X\cos \alpha \right)-\left(a_Z\sin \alpha \right)}{a_Y}$

[0071] As indicated by the above formulas, the aerodynamic angles can only be calculated in the case of non-zero coordinate acceleration. In case of hovering, i.e. with zero coordinated acceleration, a respective value of the angles calculated in the last instant of time in which the coordinate acceleration is different from zero is stored and kept until the coordinated acceleration returns to be different from zero.

[0072] It is also possible to solve the equation for both angles.

[0073] General Time Domain Solution

[0074] Preferably, considering a term as an integration over time of an acceleration to obtain the speed:

[00016]$v_B\left(t\right)=v_{B,t-\Delta t}+{\int }_{t-\Delta t}^t{\stackrel{.}{v}}_B\mathrm{dt}=v_{B,t-\Delta t}+{\int }_{t-\Delta t}^t\left(a_B-{\Omega }_B{v}_B-C_{I2B}{\stackrel{.}{w}}_I\right)\mathrm{dt}==v_{B,t-\Delta t}+{\int }_{t-\Delta t}^t{a}_B\mathrm{dt}-{\int }_{t-\Delta t}^t{\Omega }_B{v}_B\mathrm{dt}-{\int }_{t-\Delta t}^t{C_{I2B}\left(\frac{d{w}_I}{\mathrm{dt}}\right)}_I\mathrm{dt}$

[0075] and remembering the complete formula of the scalar derivative expressed in time t−Δt:

[00017]$V_{\infty ,t-\Delta t}{\stackrel{.}{V}}_{\infty ,t-\Delta t}=={\left[v_B\left(t\right)-{\int }_{t-\Delta t}^t{a}_B\mathrm{dt}+{\int }_{t-\Delta t}^t{\Omega }_B{v}_B\mathrm{dt}+{\int }_{t-\Delta t}^t{C_{I2B}\left(\frac{d{w}_I}{\mathrm{dt}}\right)}_I\mathrm{dt}\right]}^T{\left(a_B-C_{I2B}{\stackrel{.}{w}}_I\right)}_{t-\Delta t}.Math..Math.V_{\infty ,t-\Delta t}{\stackrel{.}{V}}_{\infty ,t-\Delta t}+{\left[{\int }_{t-\Delta t}^t{a}_B\mathrm{dt}-{\int }_{t-\Delta t}^t{C_{I2B}\left(\frac{d{w}_I}{\mathrm{dt}}\right)}_I\mathrm{dt}\right]}^T{\left(a_B-C_{I2B}{\stackrel{.}{w}}_I\right)}_{t-\Delta t}=={\left[v_B\left(t\right)+{\int }_{t-\Delta t}^t{\Omega }_B{v}_B\mathrm{dt}\right]}^T{\left(a_B-C_{I2B}{\stackrel{.}{w}}_I\right)}_{t-\Delta t}$

[0076] According to a first approximation (which ultimately leads to a numerically solvable definition of the problem):

∫.sub.t−Δt.sup.tΩ.sub.Bv.sub.Bdt=(Ω.sub.Bv.sub.B).sub.tΔt

[0077] the expression becomes:

[00018]$V_{\infty ,t-\Delta t}{\stackrel{.}{V}}_{\infty ,t-\Delta t}+{\left[\text{⁠}{\int }_{t-\Delta t}^t{a}_B\mathrm{dt}-{\int }_{t-\Delta t}^t{a}_B\mathrm{dt}-{\int }_{t-\Delta t}^t{C_{I2B}\left(\frac{d{w}_I}{\mathrm{dt}}\right)}_I\mathrm{dt}\right]}^T{\left(a_B-C_{I2B}{\stackrel{.}{w}}_I\right)}_{t-\Delta t}=={\left[v_{B,t}+{\Omega }_{B,t}{v}_{B,t}\Delta t\right]}^T{\left(a_B-C_{I2B}{\stackrel{.}{w}}_I\right)}_{t-\Delta t}=\left[v_{B,t}^T+{\left({\Omega }_{B,t}{v}_{B,t}\right)}^T\Delta t\right]{\left(a_B-C_{I2B}{\stackrel{.}{w}}_I\right)}_{t-\Delta t}=V_{\infty ,i}{\hat{i}}_{B,t}^T\left(I-{\Omega }_{B,t}\Delta t\right){\left(a_B-C_{I2B}{\stackrel{.}{w}}_I\right)}_{t-\Delta t}$

[0078] It is possible to apply even higher order schemes to approximate v.sub.B,t−τ.

[0079] For example, if the v.sub.B,t−τ=v.sub.B,t−{dot over (v)}.sub.B,t.sup.T integral is:

[00019]$\begin{array}{c}{\int }_{t-\Delta t}^t{\Omega }_B{v}_B\mathrm{dt}={\Omega }_{B,t}{\int }_{t-\Delta t}^t\left(v_{B,t}-{\stackrel{.}{v}}_{B,t}\tau \right)\mathrm{dt}\\ ={\Omega }_{B,t}{v}_{B,t}\Delta t-{\Omega }_{B,t}{\stackrel{.}{v}}_{B,t}\frac{\Delta t^2}{2}\end{array}$

[0080] and, substituting in the general expression:

[00020]$V_{\infty ,t-\Delta t}{\stackrel{.}{V}}_{\infty ,t-\Delta t}+{\left[\text{⁠}{\int }_{t-\Delta t}^t{a}_B\mathrm{dt}-{\int }_{t-\Delta t}^t{C_{I2B}\left(\frac{d{w}_I}{\mathrm{dt}}\right)}_I\mathrm{dt}\right]}^T{\left(a_B-C_{I2B}{\stackrel{.}{w}}_I\right)}_{t-\Delta t}=={\left[v_{B,t}+{\Omega }_{B,t}{v}_{B,t}\Delta t-{{\Omega }_{B,t}\left(a_B-{\Omega }_B{v}_B-C_{I2B}{\stackrel{.}{w}}_I\right)}_t\frac{\Delta t^2}{2}\right]}^T{\left(a_B-C_{I2B}{\stackrel{.}{w}}_I\right)}_{t-\Delta t}.Math..Math.V_{\infty ,t-\Delta t}{{\stackrel{.}{V}}_{\infty ,t-\Delta t}++\left[{\int }_{t-\Delta t}^t{a}_B\mathrm{dt}-{\int }_{t-\Delta t}^t{C_{I2B}\left(\frac{d{w}_I}{\mathrm{dt}}\right)}_I\mathrm{dt}+{{\Omega }_{B,t}\left(a_B-C_{I2B}{\stackrel{.}{w}}_I\right)}_t\frac{\Delta t^2}{2}\right]}^T{\left(a_B-C_{I2B}{\stackrel{.}{w}}_I\right)}_{t-\Delta t}=={\left[v_{B,t}+{\Omega }_{B,t}{v}_{B,t}\Delta t+{\Omega }_{B,t}{\Omega }_{B,t}{v}_{B,t}\frac{\Delta t^2}{2}\right]}^T{\left(a_B-C_{I2B}{\stackrel{.}{w}}_I\right)}_{t-\Delta t}==V_{\infty ,i}{\hat{i}}_{B,t}^T\left(I-{\Omega }_{B,t}\Delta t+{\Omega }_{B,t}^2\frac{\Delta t^2}{2}\right){\left(a_B-C_{I2B}{\stackrel{.}{w}}_I\right)}_{t-\Delta t}$

[0081] Considering a constant wind within the time interval, and thus a zero wind acceleration, the formulas respectively become:

[00021]$V_{\infty ,t-\Delta t}{\stackrel{.}{V}}_{\infty ,t-\Delta t}+{\left[\text{⁠}{\int }_{t-\Delta t}^t{a}_B\mathrm{dt}+{\Omega }_{B,t}{a}_{B,t}\frac{\Delta t^2}{2}\right]}^T{a}_{B,t-\Delta t}==V_{\infty ,t}{\hat{i}}_{B,t}^T\left(I-{\Omega }_{B,t}\Delta t+{\Omega }_{B,t}^2\frac{\Delta t^2}{2}\right)a_{B,t-\Delta t}\mathrm{and}{V}_{\infty ,t-\Delta t}{\stackrel{.}{V}}_{\infty ,t-\Delta t}+{\left[\text{⁠}{\int }_{t-\Delta t}^t{a}_B\mathrm{dt}\right]}^T{a}_{B,t-\Delta t}==V_{\infty ,t}{\hat{i}}_{B,t}^T\left(I-{\Omega }_{B,t}\Delta t\right)a_{B,t-\Delta t}$

[0082] It is therefore possible to define a system of two time domain equations and two unkowns i.e. angles of attack and sideslip included explicitly in the versor î.sub.B,t at time t, that considering the above mentioned first approximation, is:

[00022]$\mathrm{SYS}1\{\begin{array}{c}{\stackrel{.}{V}}_{\infty ,t}={\hat{i}}_{B,t}^T{a}_{B,t}\\ V_{\infty ,t-\Delta t}{\stackrel{.}{V}}_{\infty ,t-\Delta t}+{\left[\text{⁠}{\int }_{t-\Delta t}^t{a}_B\mathrm{dt}\right]}^T{a}_{B,t-\Delta t}=V_{\infty ,t}{\hat{i}}_{B,t}^T\left(I-{\Omega }_{B,t}\Delta t\right)a_{B,t-\Delta t}\end{array}$

[0083] Wherein all other parameters are either measurable at time t or are already calculated at time t−Δt. It is important to note that first formula in SYS 1 is EQ 1 written in a different form.

[0084] In the non-linear form, the system may for example be solved by a Levenberg Marquardt method or another equivalent method.

[0085] Alternatively, it is possible to linearize the trigonometric functions:

î.sub.B.sup.T≈[1,β,α].sup.T

[0086] In this manner, the above mentioned system can be approximated and expressed as:

[00023]$\{\begin{array}{c}{A}_0={\hat{i}}_{B,t}^T{b}_0\\ A_0={\hat{i}}_{B,t}^T{b}_1\end{array}\mathrm{wherein}:A_0={\stackrel{.}{V}}_{\infty ,t},A_1=\frac{V_{\infty ,t-\Delta t}{\stackrel{.}{V}}_{\infty ,t-\Delta t}+{\left[\text{⁠}{\int }_{t-\Delta t}^t{a}_B\mathrm{dt}\right]}^T{a}_{B,t-\Delta t}}{V_{\infty ,t}},b_0=a_{B,t}\mathrm{and}{b}_1=\left(I-{\Omega }_{B,t}\Delta t\right)a_{B,t-\Delta t}$

[0087] Expliciting the aerodinamic angles, the system becomes:

[00024]$\{\begin{array}{c}{A}_0-b_{0x}=b_{0z}\alpha +b_{0y}\beta \\ A_1-b_{1x}=b_{1z}\alpha +b_{1y}\beta \end{array}$

[0088] and solutions are:

[00025]$\{\begin{array}{c}\alpha =\frac{.Math.\begin{array}{cc}{A}_0-b_{0x}& b_{0y}\\ A_1-b_{1x}& b_{1y}\end{array}.Math.}{.Math.\begin{array}{cc}{b}_{0z}& b_{0y}\\ b_{1x}& b_{1y}\end{array}.Math.}\\ \beta =\frac{.Math.\begin{array}{cc}{b}_{0z}& A_0-b_{0x}\\ b_{1z}& A_1-b_{1x}\end{array}.Math.}{.Math.\begin{array}{cc}{b}_{0z}& b_{0y}\\ b_{1x}& b_{1y}\end{array}.Math.}\end{array}$

[0089] FIGS. 4A and 4B show on the left the overlap between the nonlinear solution for the angle of attack (solid line) with the solution of the non-linear system. On the right there is instead the overlap with the results of the linearized system: the error of a few tenths of a degree is widely acceptable because it is lower than the errors produced by the actual systems for measuring these angles by at least one order of magnitude. In both cases, the wind variation law illustrated in FIG. 2 is applied.

[0090] According to another embodiment of the present invention, the SYS can be expanded forward and/or backward in time adding as many equations as necessary to achieve the required estimation accuracy. It is therefore possible to define a system of n+1 equations and two unknowns, i.e. angles of attack and sideslip included explicitly in the versor î.sub.B,t at time t, that considering the above mentioned first approximation, is:

[00026]$\mathrm{SYS}2\{\begin{array}{c}{\stackrel{.}{V}}_{\infty ,t}={\widehat{\mathrm{.Math.}}}_{B,t}^T{a}_{B,t}\\ {\stackrel{.}{V}}_{\infty ,t-\Delta t}{V}_{\infty ,t-\Delta t}+{\left[\underset{t-\Delta t}{\stackrel{t}{\int }}{a}_B\mathrm{dt}\right]}^T{a}_{B,t-\Delta t}=V_{\infty ,t}{\widehat{\mathrm{.Math.}}}_{B,t}^T\left(I-{\Omega }_{B,t}\Delta t\right)a_{B,t-\Delta t}\\ .Math.\\ {\stackrel{.}{V}}_{\infty ,t-n\Delta t}{V}_{\infty ,t-n\Delta t}+{\left[\underset{t-n\Delta t}{\stackrel{t}{\int }}{a}_B\mathrm{dt}\right]}^T{a}_{B,t-n\Delta t}=\\ V_{\infty ,t}{\widehat{\mathrm{.Math.}}}_{B,t}^T\left(I-{\Omega }_{B,t}\Delta t\right)a_{B,t-n\Delta t}\end{array}$

[0091] Wherein all other parameters are either measurable at time t or are already calculated at time t−Δt, t−2Δt, . . . ,t−nΔt.

[0092] According to an embodiment of EQ 1 where coordinate acceleration and raw aerodynamic angles are given, it is possible to estimate {dot over (V)}.sub.∞ via EQ1:

{dot over (V)}.sub.∞=a.sub.X cos {tilde over (β)}cos {tilde over (α)}+a.sub.Y sin {tilde over (β)}+a.sub.Z cos {tilde over (β)}sin {tilde over (α)}

[0093] where {tilde over (α)} and {tilde over (β)} are raw measures of α and β respectively or {tilde over (α)} is the measure of α and {tilde over (β)} is an approximation of β or vice versa. Then the estimated {dot over (V)}.sub.∞ is used in SYS1 or SYS2 to refine the aerodynamic angle estimations.

[0094] Finally, it is clear that it is possible to make changes or variations to the method described and illustrated here without departing from the scope of protection as defined by the attached claims.