LA NAVIGATION INTÉGRÉE OPTIMALE PAR
FILTRAGE DE KALMAN

Philippe Allard

► To cite this version:

Philippe Allard. LA NAVIGATION INTÉGRÉE OPTIMALE PAR FILTRAGE DE KALMAN. [Rapport de recherche] Ecole Nationale de la Marine Marchande. 1984. hal-02025723

HAL Id: hal-02025723

https://hal.science/hal-02025723

Submitted on 19 Feb 2019

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

Public Domain
LA NAVIGATION INTÉGRÉE OPTIMALE

PAR FILTRAGE DE KALMAN

<figure><img src="image_1.png" /></figure>

Ph.Allard
1984
MEMOIRE

LA NAVIGATION INTEGRÉE OPTIMALE

PAR FILTRAGE DE KALMAN

par

Philippe ALLARD
AVANT-PROPOS

Ce mémoire se fonde sur une étude que nous avons menée au
sein du Service Hydrographique et Océanographique de la Marine (1982)

Il a pour but essentiel de communiquer l'expérience ainsi
acquise sur une méthode puissante de traitement de l'information et
sur son application à la navigation.

Pour l'exposé des notions théoriques nécessaires à sa lecture
nous avons largement emprunté à l'ouvrage "Introduction au Filtrage
Numérique" de j.c Radix.
## TABLE DES MATIERES

### Chapitre I NAVIGATION INTEGREE OPTIMALE

1 Définition, principe et objectifs

2 Choix des capteurs

3 Filtrage de Kalman et radionavigation integree optimale

4 Utilisation et interet du filtrage de Kalman

### Chapitre II FILTRAGE

1 Influence des bruits sur les systèmes physiques et sur les mesures

2 Le filtrage

3 Les grandes familles de filtres

### Chapitre III METHODE DES MOINDRES CARRES

1 Principes généraux

2 Application "classique" au filtrage

3 Le filtre de Wiener

### Chapitre IV METHODE DE KALMAN

1 Eléments sur les systèmes linéaires

2 Conditions d'application du filtre de Kalman

3 Description et interprétation du filtre de Kalman

4 Extensions de la méthode

5 Conclusion

### Chapitre V MODELE DE RADIONAVIGATION INTEGREE OPTIMALE

1 Modèle continu

2 Modèle discretisé

### Chapitre VI DEVELOPPEMENT DU LOGICIEL

1 Principes généraux

2 Simulation du modèle de radionavigation intégrée
# Appendice I ELEMENTS DE GEODESIE

1. Rayons de courbure de l'ellipsoïde

2. Géodésique sur l'ellipsoïde

3. Phases théoriques d'un système de radiolocalisation hyperbolique

# Appendice II ELEMENTS DE CALCUL DES PROBABILITES

1. Définition des probabilités

2. Calcul des probabilités

3. Variables aléatoires

4. Suites de variables aléatoires

5. Fontions aléatoires

6. Analyse spectrale des signaux aléatoires
# NAVIGATION INTEGREE OPTIMALE

## 1. DEFINITION, PRINCIPE ET OBJECTIFS
DE LA NAVIGATION INTEGREE

Si l'on cherche à définir en termes généraux la navigation intégrée, on peut dire qu'elle consiste à utiliser au mieux un ensemble d'équipements ou de capteurs de navigation, notamment en combinant de façon adéquate les informations qu'ils fournissent. Cette combinaison est habituellement réalisée par une technique de filtrage qui a pour but de tirer le meilleur parti de toute l'information disponible. Le principal objectif est ainsi d'obtenir un résultat global plus performant, qu' précision par exemple, que celui fourni par l'un quelconque des équipements.

## 2. CHOIX DES CAPTEURS DE NAVIGATION

Le succès dans la réalisation de l'objectif précédant repose pour commencer sur le choix des équipements de navigation. Le premier critère de selection doit être la complémentarité des informations fournies par les différents équipements. Par complémentarité il faut d'abord entendre celle de leurs performances et plus précisément celle des caractéristiques de répartition spectrale de leurs erreurs respectives.

Un bon choix sera d'associer un équipement ayant une erreur de type "haute fréquence" (un équipement de radionavigation par exemple) à un équipement ayant une erreur de type "basse fréquence" (une erreur ou une centrale inertielle par exemple). De plus la complémentarité doit être souvent comprise dans un sens plus large et concerner des aspects tels que la disponibilité temporelle (information continue ou discrète avec une periodicité plus ou moins grande), la disponibilité géographique (domaine de couverture des équipements non autonomes) et la nature des informations (position, vitesse, distance, relèvement, cap, etc.).

## 3. FILTRAGE DE KALMAN ET
RADIONAVIGATION INTEGREE OPTIMALE

La figure 1 montre le schéma de principe du système de navigation intégrée développé au chapitre V et utilisant un équipement de radio-localisation hyperbolique associé aux capteurs de vitesse et de cap.
I.2

<figure><img src="image_2.png" /><figcaption>fig I.1 système de navigation intégrée "radiolocalisation-estime"</figcaption></figure>

La combinaison des informations est réalisée par filtrage des différences entre les informations de la radiolocalisation et celles déduites des autres capteurs. Ces différences entre les phases représentent les différences de leurs erreurs respectives. Dans les systèmes intégrés performants et modernes la technique de filtrage qui s'impose, comme nous le verrons plus loin, est celle du filtrage "optimal" de Kalman. Ce filtrage fournit la "meilleure" estimee possible des erreurs.

Dans notre modèle ces estimées sont utilisées pour corriger la position "estime", pour renseigner sur l'environnement météo-océanique à travers la dérive globale due au courant et au vent, et pour étalonner le loch et le gyrocompas.

## 4. INTERET DU FILTRAGE DE KALMAN

La présentation détaillée du filtrage de Kalman fait l'objet du chapitre IV. Le but de ce paragraphe est d'indiquer, de façon simple, le principe du filtrage optimal de Kalman tel qu'il est utilisé pour la navigation intégrée et son intérêt pour ce type d'application.

<figure><img src="image_3.png" /><figcaption>fig I.2 principe du filtrage optimal de Kalman</figcaption></figure>
I.3

La figure 2 représente le schéma de principe simplifié de ce type d'utilisation du filtrage de kalman. Elle montre le caractère bouché du filtre. Les échanges représentés sont généralement multi-dimensionnel. Comme nous l'avons vu précédemment, les entrées Z du filtre sont des différences entre les informations délivrées par les capteurs c'est à dire les différences de leurs erreurs. Plus généralement ce sont des combinaisons linéaires de leurs erreurs. On appelle ces entrées des observations. Ces observations sont comparées aux combinaisons linéaires Z correspondantes des erreurs X estimées par le filtre. Les différences obtenues Z-Z, multipliées par des gains "optimaux" K, servent à corriger les erreurs estimées X, dont l'évolution est entretenue par un modèle de leurs dynamiques.

La partie 1a plus complexe du filtre concerne le calcul en temps réel des gains optimaux. Ce calcul fait notamment intervenir les modèles dynamiques statistiques des erreurs des capteurs. Il est généralement volumineux et délicat (intégration en temps réel, précision numérique,...).

La technique du filtrage de Kalman s'applique aux systèmes linéaires ou linéarisés. Elle a l'avantage sur les technique classiques de filtrage linéaire d'être plus générale et notamment beaucoup mieux adaptée au traitement de dynamiques non stationnaires, de perturbations statistiques non stationnaires et de systèmes multi-dimensionnels. Ces caractéristiques sont indissociables des systèmes de navigation intégrée performants.

Outre cette plus grande généralité, le filtrage de Kalman est une
technique optimale, dans le sens qu'elle minimise l'intégrale de l'er-
reur quadratique et conduit en ce sens aux meilleurs performances
possibles. Ceci n'est vrai que si les modèles utilisés représentent
exactement la réalité. En pratique, la complexité des phénomènes phy-
siques et les simplifications nécessaires à la réalisation conduisent
à utiliser des modèles approchés. Ces approximations et ces simplifi-
cations tendent à réduire le caractère "optimal" des performances.
Leur étude délicate fait partie de celle des compromis du develop-
pement. De plus, lors de la réalisation du filtre, la validite des in-
formations des différents capteurs, leur conformité aux modèles rete-
nus doivent être vérifiées par des tests. L'utilisation d'informations
non valides ou non conformes doit être interdite, sinon de graves dé-
fauts de fontionnement sont à craindre.
# FILTRAGE

## 1. INFLUENCE DES BRUITS SUR LES SYSTEMES PHYSIQUES ET SUR LES MESURES

D'une manière générale, l'étude expérimentale des systèmes matériels est gênée par l'influence des bruits, phénomènes fluctuants dont l'évolution ne peut être prévue exactement. Ces bruits interviennent à double titre :

- d'une part, dans le système physique étudié,

- d'autre part, dans les mesures que le physicien effectue pour étudier le système.

Cette situation est schématisée sur la figure 1 :

<figure><img src="image_4.png" /><figcaption>fig. II.4</figcaption></figure>

Le physicien se trouve ainsi en présence de signaux de la forme :

$$ m = \begin{array}{c} \uparrow \\ \text{valeur mesurée} \end{array} + \begin{array}{c} \downarrow \\ \text{sigmal utile} \end{array} + b \quad (1) $$

Les bruits, par définition, présentent un caractère aléatoire. L'étude de leur influence-tant sur les systèmes que sur les mesures-peut néanmoins être faite grâce à l'utilisation des méthodes du calcul des probabilités. En effet, les mathématiciens ont défini et étudié les fonctions aléatoires dont l'évolution est soumise au hasard; ces fonctions constituent d'excellents modèles représentatifs des bruits.

Le caractère des signaux utiles est plus ambigu. En effet, ils comprennent en général :

- une partie certaine pouvant normalement être calculée à l'avance en fonction des conditions initiales et d'une loi d'évolution;

- une partie aléatoire due à l'influence des bruits agissant sur le système.

Ceci nous conduit à décomposer tout signal utile de la manière suivante :

$$ x = x_c + x_a \quad (2) $$
## 2. LE FILTRAGE

Notre problème peut maintenant être précisé: connaissant la valeur mesurée $m$, il s'agit de faire une estimation (si possible la meilleur nous verrons plus loin ce qu'il faut entendre par là) de la valeur réelle $x$.

Pour avoir une chance de résoudre un tel problème, il nous faut :

1) Disposer d'un moyen qui nous permette de distinguer le signal utile $x$ du bruit de mesure. Ce moyen d'identification sera fourni par une étude statistique des habitudes du signal utile et du bruit de mesure.

2) Préciser ce que l'on entend par "la meilleure approximation" d'un paramètre : ce sera le but de critères d'optimalisation.

On appelle filtrage cette opération consistant à extraire le signal utile $x$ de la spmme $m=x+b$, seule information disponible. Ce terme est emprunté au langage des electroniciens.

Mais il existe une autre méthode pour traiter un message :

- échantillonner $m, c$ est à dire isoler des valeurs prises par $m$ à des instants $t_0, t_1, \dots, t_n$ suffisamment raprochés;

- digitaliser ou disjonctiser ces valeurs, c'est à dire exprimer ces nombres à l'aide de chiffres;

- traiter la suite des nombres ainsi obtenue par une méthode mathématique appropriée à l'aide d'un ordinateur.

Bien entendu, les filtres ainsi obtenus n'ont pas d'existence physique comme c'est le cas des filtres analogiques. Mais la puissance des méthodes mathématiques modernes et l'efficacité des ordinateurs ont donné à ce procédé (le filtraje numérique) une souplesse remarquable.

Toutes les méthodes dérivent d'un principe commun : distinguer le signal utile du bruit en utilisant deux catégories d'informations :

- informations a priori : lois (déterministes ou statistiques) auxquelles obéissent les signaux utiles et les bruits;

- information a posteriori : résultats de mesures.

<table><thead><tr><th rowspan="2">Méthodes</th><th colspan="2">Informations a priori</th></tr><tr><th>sur le signal x</th><th>sur le bruit b</th></tr></thead><tbody><tr><td>Moyens Caractéristiques</td><td>Forme analytique</td><td></td></tr><tr><td>Wiener</td><td>Auto-correlation et<br>correlation avec b</td><td>Auto-correlation<br>et corrélation avec x</td></tr><tr><td>Kalman</td><td>Équations d'état<br>Équipes statistiques<br>des perturbations</td><td>Propriétés<br>statistiques</td></tr></tbody></table>
## 3. LES GRANDES FAMILLES DE FILTRES

### 3.1 Filtre linéaire

On appelle linéaire un filtre vérifiant les relations suivantes :

$$ \hat{y}_1 + y_2 = \hat{y}_1 + \hat{y}_2 \quad (3) \quad \text{et} \quad \hat{a}_y = \hat{a}_y \quad (4) $$

($\hat{y}$ représente y filtré et a est une constante ).

La majorité des filtres numériques utilisés actuellement sont linéaires.

### 3.2 Filtre en temps réel

Un filtre est dit en temps réel lorsqu'il détermine une valeur filtrée $\hat{y}(t)$ en fonction d'informations y(t') avec t'≤t.

### 3.3 Filtre à mémoire finie, filtre récursif

Un filtre est dit à mémoire finie lorsqu'il détermine $\hat{y}(t)$ en fonction d'informations d'anciènneté au plus égale à T, T étant une durée déterminée. Par exemple un filtre de la forme

$$ \hat{y}(t) = \sum_{t'-n}^{n} A_i y(t+t') \quad (5) $$

est à mémoire finie (on a $T=N/\Delta t$).

Par contre un filtre est dit à mémoire infinie lorsqu'il détermine y(t) en fonction d'informations d'anciènneté quelconque. Par exemple un filtre de la forme

$$ \hat{y}(t) = Ay(t) + By(t-T) \quad (6) $$

est à mémoire infinie.

Les filtres numériques faisant intervenir les valeurs filtrées précédentes sont dits récursifs. Par exemple (6).

### 3.4 Filtre déterministe

Un filtre est dit déterministe lorsque son étude théorique ne fait pas intervenir le calcul des probabilités. Dans le cas contraire il est dit statistique ou stochastique.
# METHODE DES MOINDRES CARRES

## 1. PRINCIPES GENERAUX

La méthode des moindres carrés consiste à minimaliser une forme quadratique représentant l'écart entre une solution obtenue et une solution idéale.

Le choix d'une forme quadratique comme indice de qualité peut sembler a priori bien arbitraire : on aurait pu aussi décider de minimaliser une somme de modules d'écarts, ou encore une somme de modules d'écarts élevés à une puissance quelconque (3,4,...). En fait le succès considérable de la méthode des moindres carrés depuis près de deux siècles semble être dû aux raisons suivantes :

- la minimalisation d'une forme quadratique conduit souvant à un système d'équations linéaires;

- les sommes de carrés sont fréquemment l'objet de certaines interprétations : par exemple, en physique les energies sont exprimées par des sommes de carrés; de même les variances en calcul des probabilités.

L'application de la méthode des moindres carrés aux problèmes de filtrage peut être menée de différentes manières :

1) En admettant une forme analytique bien déterminée pour le signal et éventuellement pour le bruit; ceci constitue le principe d'application classique de cette méthode.

2) En admettant certaines propriétés statistiques concernant le signal et le bruit; ceci constitue le point de départ de la méthode de Wiener.

3) En admettant à la fois certaines propriétés analytiques sur le signal, et certaines propriétés statistiques sur le signal et sur le bruit; ceci constitue le point de départ de la méthode de Kalman, forme moderne de la méthode des moindres carrés.

## 2. APPLICATION "CLASSIQUE" AU FILTRAGE

### 2.1 Définition des formes analytiques

Nous devons traiter un message

$$m = \text{signal } s + \text{bruit } b$$
Le critère d'identification qui va nous permettre de séparer s et b doit être leur forme analytique. Le filtre est donc déterministe. En toute rigueur il serait souhaitable de disposer :

- d'une forme analytique pour le signal, notée $s_{th}$;

- d'une forme analytique pour le bruit, notée $b_{th}$.

La méthode des moindres carrés consisterait alors à ajuster les paramètres entrant dans la définition de $s_{th}$ et $b_{th}$, pour minimaliser une somme (ou une intégrale) des carrés des écarts

$$m - \bar{m} = (s+b) - (s_{th} + b_{th}) \quad (1)$$

La résolution des systèmes d'équations traduisant cette minimalisation fournirait une estimation des paramètres en question, d'ou une expression analytique du signal et du bruit.

## 2.2 Optimalisation

Nous disposons :

d'un message $m(t)=\text{signal utile } s(t) + \text{bruit } b(t)$

d'une forme analytique pour $s(t)$ : $s_{th}=s_{th}(t;a_1,a_2,\dots,a_k)$

t est une variable indépendante représentant le temps, ou une longueur, un numéro d'ordre...

$a_1,a_2,\dots,a_k$ sont des paramètres dont il va falloir faire l'estimation.

Cette estimation, et par suite celle de $s$, va être effectuée en minimalisant la quantité $Q$ définie par

$$Q = \int_{a}^{b} (m(t) - s_{th}(t))^2 dt \quad (2)$$

si t est continu, ou par

$$Q = \sum_{n=n_1}^{n_2} (m(t_n) - s_{th}(t_n))^2 \quad (3)$$

si t est discret.

Ce résultat peut en général être obtenu en résolvant le système d'équations :

$$\frac{\partial Q}{\partial a_1} = 0, \quad \frac{\partial Q}{\partial a_2} = 0, \quad \dots, \quad \frac{\partial Q}{\partial a_k} = 0 \quad (4)$$

Le nombre des équations (4) est égale au nombre $k$ des inconnues $a_i$.

A priori, les relations (4) traduisent seulement la stationnarité de $Q$; cependant si la fonction $s_{th}$ est bien choisie, les relations (4) conduisent effectivement à un minimum; en particulier ceci se produit lorsque le système (4) est linéaire.
## 2.3 Résolution d'un système linéaire redondant

Considérons un système de n équations linéaires à deux inconnues

$$ 
\begin{aligned}
& a_1 x + b_1 y + c_1 = 0 \\
& \cdot \\
& a_n x + b_n y + c_n = 0
\end{aligned} 
\qquad (5) 
$$

Ce système est redondant dès que n est supérieur à deux et il n'admet alors aucune solution, car dans le cas général les droites représentant ces équations ne sont pas concourantes.

Considérons la quantité :

$$ Q = \sum_{i=1}^{n} (a_i x + b_i y + c_i)^2 = Q(x, y) \qquad (6) $$

6'est la somme des carrés des résidus de chaque équation pour un couple $(x,y)$. Elle serait nulle pour la solution du système (5) si elle existait. La solution convenant "au mieux" à l'ensemble des équations (5) est alors celle qui minimise la quantité $Q$. Elle donc solution du système :

$$ \frac{\partial Q}{\partial x} = 0 \qquad \frac{\partial Q}{\partial y} = 0 \qquad (7) $$

qui est un système de deux équations linéaires à deux inconnues :

$$ 
\begin{aligned}
& \left(\sum_{i=1}^{n} a_i^2\right) x + \left(\sum_{i=1}^{n} a_i b_i\right) y + \left(\sum_{i=1}^{n} a_i c_i\right) = 0 \\
& \left(\sum_{i=1}^{n} a_i b_i\right) x + \left(\sum_{i=1}^{n} b_i^2\right) y + \left(\sum_{i=1}^{n} b_i c_i\right) = 0
\end{aligned} 
\qquad (8) 
$$

admettant toujours une solution et une seule.

Ce problème apparaît lorsqu'on calcule une estimation de la position d'un navire à partir de plusieurs lieux assimilables localement à des droites.

## 3. LE FILTRE DE WIENER

Vers 1949 Ni Wiener a posé explicitement le problème du filtrage d'un signal aléatoire stationnaire. Cette manière de voir le traitement de l'information était alors tout à fait originale; elle conduisit Wiener à proposer des filtres linéaires minimalisant une intégrale (ou une somme) de carrés d'éarts entre le signal utile et le signal estimé. Cette méthode diffère cependant de la méthode classique des moidres carrés en ce sens qu'elle ne caractérise pas le signal et le bruit par leur forme analytique, mais par leurs propriétés statistiques résumées par :

- les valeurs moyennes, supposées nulles,

- les fonctions corrélations (auto-corrélation et corrélations mutuelles).

Le filtre de Wiener est donc un filtre statistique.

Dans le cas de signaux échantillonnés de la forme :

$$ m_k = x_k + b_k \qquad (9) $$
on recherche le filtre numérique dont les coefficients de pondération $A_i$ minimalisent la quantité :

$$Q = \lim_{K \to \infty} \frac{1}{2k+1} \sum_{k=-K}^{+K} (x_k - \hat{m}_k)^2 \quad (10)$$

avec

$$\hat{m}_k = \sum_{i=-N/2}^{N/2} A_i m_{k+i} \quad (11)$$

Ces coefficients $A_i$ doivent vérifier le système linéaire

$$\frac{\partial Q}{\partial A_i} = 0 \quad (i=-N/2 \text{ à } +N/2) \quad (12)$$

En remplaçant $\hat{m}_k$ par sa valeur (11) dans l'expression (10), et en effectuant le calcul des dérivées partielles de (12), nous obtenons :

$$\sum_{j=-N/2}^{N/2} A_j R_{mm} (i-j) = R_{mx} (i) \quad (i=-N/2 \text{ à } +N/2) \quad (13)$$

système de $(N'+N+1)$ équations à $(N'+N+1)$ inconnues permettant le calcul explicite des $A_i$ en fonction de

$$R_{mm}(1) = \lim_{K \to \infty} \frac{1}{2k+1} \sum_{-K}^{K} m_h^2 m_{h+1}^2 \quad (14)$$

(auto-corrélation du message), et de

$$R_{mx}(1) = \lim_{K \to \infty} \frac{1}{2k+1} \sum_{-K}^{K} m_h^x m_{h+1}^x \quad (15)$$

(corrélation mutuelle du signal et du message).

(On admet que $x$ et $m$ sont aléatoires, stationnaires, ergodiques et de valeur moyenne nulle.)

La connaissance des fonctions $R_{mm}$ et $R_{mx}$ permet la résolution de (13), et finalement la détermination du filtre optimal de Wiener. Ainsi, les coefficients $A_i$ ne dépendent pas spécialement de telle ou telle valeur de $x$ (ou de $m$), mais simplement de leurs propriétés statistiques.

Un cas limite est celui où le signal et le message ont une corrélation mutuelle nulle : les "seconds membres" $R_{mx}(i)$ et par suite tous les coefficients de pondération sont nuls. Dans ce cas spécialement défavorable, la connaissance du message $m$ n'apprend rien sur le signal utile!
# METHODE DE KALMAN

La méthode de kalman version moderne et probabiliste de la méthode des moindres carrés, concerne essentiellement les systèmes linéaires. Ceci nous conduit à faire quelques rappels à leur sujet, au paragraphe suivant. L'étude effective du filtre de Kalman commencera au §2.

## 1. ELEMENTS SUR LES SYSTEMES LINEAIRES

### 1.1 Variables d'état. Equations d'état

Nous considérons un système dont l'état est défini par un ensemble fini de variables d'état $x_1(t), x_2(t), \dots, x_n(t)$ obéissant à un système d'équations de la forme suivante :

$$ 
\begin{aligned}
\frac{dx_1}{dt} &= A_{11}x_1 + A_{12}x_2 + \dots + A_{1n}x_n + B_1 \\
\frac{dx_2}{dt} &= A_{21}x_1 + A_{22}x_2 + \dots + A_{2n}x_n + B_2 \\
\vdots & & \vdots \\
\frac{dx_n}{dt} &= A_{n1}x_1 + A_{n2}x_2 + \dots + A_{nn}x_n + B_n
\end{aligned}
\qquad (1)
$$

La variable indépendante $t$ représente le temps.

Les seconds membres du système d'équations (1) sont des formes linéaires des variables d'état $x_i(t)$ : on dit que ce système est linéaire.

### 1.2 Forme matricielle des équations d'état

Il est très commode d'écrire (1) sous forme matricielle; dans ce but on considère :

Une matrice colonne définissant le vecteur d'état

$$ X(t) = \begin{bmatrix} x_1(t) \\ x_2(t) \\ \vdots \\ x_n(t) \end{bmatrix} \qquad (2) $$
Une matrice carrée définissant les coefficients du système :

$$A(t) = [A_{ij}(t)] \quad (3)$$

Une matrice colonne définissant l'excitation ou la commande du système :

$$B(t) = \begin{bmatrix} B_1(t) \\ B_2(t) \\ \vdots \\ B_n(t) \end{bmatrix} \quad (4)$$

Dans ces conditions le système (1) s'écrit simplement :

$$\dot{X}(t) = \frac{dX(t)}{dt} = A(t)X(t) + B(t) \quad (5)$$

## 1.3 Résolution des équations d'état. Matrice de transition

Etant donné une équation homogène :

$$\dot{X}(t) = A(t)X(t) \quad (6)$$

assortie des conditions initiales

$$X(t_o) = [X(t)]_{t=0}$$

Il est en général possible de mettre sa solution sous la forme :

$$X(t) = \phi(t, t_o)X(t_o) \quad (7)$$

$\phi(t, t_o)$ appelée matrice de transition, peut être considéré comme un opérateur linéaire transformant $X(t_o)$ en $X(t)$. Cette matrice carrée est solution de l'équation différentielle matricielle suivante :

$$\dot{\phi}(t) = A(t)\phi(t) \quad (8)$$

par sa définition même, l'opérateur $\phi(t, t_o)$, applicable à tous les instants $t_o$ et $t$, doit obéir aux relations :

$$\phi(t, t) = I \quad (9)$$

$$\phi(t_3, t_2)\phi(t_2, t_1) = \phi(t_3, t_1) \quad (10)$$

$$\phi(t_1, t_2) = \phi^{-1}(t_2, t_1) \quad (11)$$

La matrice de transition $\phi$, définie à propos de la résolution du système homogène (6), permet de résoudre le système complet (5); en effet, on démontre que la solution de (5) peut être mise sous la forme suivante :

$$X(t) = \phi(t, t_o)X(t_o) + \int_{t_o}^{t} \phi(t, s)B(s)ds \quad (12)$$
Le premier terme de (12) décrit l'influence de l'état initial du système.

Le deuxième terme décrit l'influence de l'excitation B(t) sur le système.

Jusqu'à présent, nous avons considéré des variables d'état variant d'une manière continue, Interessons nous maintenant aux suites de nobres obtenus par échantillonnage de ces variables d'état, à des instants successifs $t_0, t_1, \dots, t_k, \dots, t_N$.

Dans le cas d'une équation homogène, le passage de $X(t_{k-1})$ à $X(t_k)$ se fait exactement par application de (7); on a :

$$X(t_k) = \phi_{k, k-1} X(t_{k-1}) \quad (13)$$

avec

$$\phi_{k, k-1} = \phi(t_k, t_{k-1})$$

Dans le cas d'une équation complète (5), nous avons :

$$X(t_k) = \phi_{k, k-1} X(t_{k-1}) + F_{k-1} \quad (14)$$

avec

$$F_{k-1} = \int_{t_{k-1}}^{t_k} \phi(t_k, s) B(s) ds \quad (15)$$

L'équation (14) décrit l'évolution des grandeurs d'état échantillonnées du système étudié; elle constitue la forme digitalisée ou discretisée de l'équation (12).

## 2. CONDITIONS D'APPLICATION DU FILTRE DE KALMAN

L'application du filtre de Kalman à l'étude d'un système physique suppose que l'on ait réussi à représenter ce dernier par un système d'équations linéaires (S). Par ailleurs, le filtre lui même est également constitué d'un système d'équations linéaires (F). Les systèmes (S) et (F) peuvent présenter deux formes : continues, il s'agit alors d'équations differentielles; digitales, il s'agit alors d'équations aux differences. Le choix entre les formes continues et digitales est guidé par les circonstances suivantes :

- l'exploitation du filtre implique l'utilisation d'un ordinateur en raison de la complexité des opérations à effectuer et de la précision nécessaire à chaque opération;

- la théorie du filtre est moins difficile à établir sous la forme digitale.

Ceci nous incite évidemment à ne considérer que la forme digitale du filtre de Kalman; cependant, il se trouve que dans de nombreux cas le modèle mathématique représentant le système à étudier est consétitué d'un système différentiel auquel devrait donc être normalement
appliqué le filtre continu de Kalman; la solution habituellement adoptée consiste à digitaliser le modèle mathématique du système, c'est à dire le remplacer par un système d'équations aux differences équivalent, auquel peut être commoément appliqué le filtre de Kalman digital.

Ainsi, l'objet des paragraphes suivants va être tout d'abord de préciser la forme du modèle mathématique décrivant le système physique étudié, puis de donner le moyen d'en effectuer la digitalisation.

## 2.1 Modèle mathématique adopté

L'état du système est représenté par un vecteur d'état $X(t)$, obéissant à une équation de la forme suivante :

$$ \frac{dX(t)}{dt} = A(t)X(t) + F(t) + W(t) \quad (16) $$

avec

* A matrice représentant la dynamique du système

* F commande

* W bruit du système (bruit blanc)

D'après (12), la solution de cette équation est :

$$ X(t) = \phi(t, t_o)X(t_o) + \int_{t_o}^{t} \phi(t,s)F(s)ds + \int_{t_o}^{t} \phi(t,s)W(s)ds \quad (17) $$

Par ailleurs, l'application du filtre de Kalman suppose que le système est soumis à des mesures portant sur certaines combinaisons liéaires des variables d'état. Ces mesures sont les composantes d'un vecteur $Z(t)$ obéissant à la relation suivante :

$$ Z(t) = H(t)X(t) + V(t) \quad (17) $$

avec

* H matrice de mesure

* V bruit de mesure (bruit blanc)

Le nombre m des mesures est en général différent du nombre n de variables d'état; la matrice H est alors rectangulaire.

## 2.2 Digitalisation du modèle mathématique

### a) digitalisation de l'équation de transition

Nous posons

$$ X(t_k) = X_k \quad (19) $$

$$ \phi(t_k, t_{k-1}) = \phi_{k, k-1} \quad (20) $$

Il en résulte que l'équation de transition (17) peut être mise sous la forme récurrente suivante :

$$ X_k = \phi_{k, k-1} X_{k-1} + F_{k-1} + W_{k-1} \quad (21) $$
les termes $\phi_{k,k-1}^X$ et $F_k$ sont obtenus par application des formules établies $\phi_{k,k-1}^{\text{au Sf}}$ le terme aléatoire $W_{k-1}$, d'après (15), admet pour expression

$$W_{k-1} = \int_{t_{k-1}}^{t_k} \phi(t_k, s) W(s) ds \quad (22)$$

ce qui implique

$$EW_{k-1} = \int_{t_{k-1}}^{t_k} \phi(t_k, s) (EW(s)) ds = 0 \quad (23)$$

la valeur moyenne du terme $W_{k-1}$ est donc nulle.

En ce qui concerne le calcul de la variance des composantes de $W_k$ (nécessaires au calcul de filtrage) l'application de la formule (23) posent quelques problèmes dus en particulier à l'intégration du vecteur aléatoire non stationnaire $\phi(t_k, s)W(s)$ et à l'obligation d'obtenir des $W_k$ succéssifs non corrélatifs, condition nécessaire à l'application du filtre de Kalman digital.

Cette question peut être résolue apprимативноm en assimilant les composantes $w_i(t)$ de $W(t)$ à des bruits pseudo-blancs échantillonnés à la cadence $1/(t_k-t_{k-1})$ dont on connaît les densités spectrales de puissance moyenne $S_{ij}^*$ et les densités spectrales mutuelles $S_{ij}$. Il est possible d'en déduire une valeur approchée de $Q_k$ de la manière suivante :

Dans (22), $\phi(t_k, s)$ est voisinne de la matrice unité I; on obtient ainsi

d'ou

$$W_{k-1} \neq \int_{t_{k-1}}^{t_k} I W(s) ds \\ EW_{k-1} = 0 \quad (24)$$

et

$$Q_k = EW_{k-1} W_{k-1}^T \# E \left( \int_{t_{k-1}}^{t_k} W(s) ds \int_{t_{k-1}}^{t_k} (W(s))^T ds \right) \\ = \left[ E \left( \int_{t_{k-1}}^{t_k} w_i(s) ds \right) \left( \int_{t_{k-1}}^{t_k} w_j(s) ds \right) \right] \quad (25)$$

D'après AII(35,36), cette expression peut approximativement être remplacée par

$$[\varphi_{ij}(0)] (t_k - t_{k-1})^2 \quad (26)$$

$\varphi_{ij}(s)$ étant la corrélation mutuelle de $w_i(t)$ et de $w_j(t)$.

Mais d'après AII(51,53), nous avons

$$\varphi_{ij}(0) = E(w_i(t) w_j(t)) = \frac{1}{2\pi} \int_{-\pi/\Delta t}^{\pi/\Delta t} S_{ij}^* du = \frac{S_{ij}^*}{\Delta t} \quad (27)$$

en faisant $\Delta t = t_k - t_{k-1}$, nous obtenons finalement

$$Q_k = [S_{ij}^*] (t_k - t_{k-1}) \quad (27\text{bis})$$
b) digitalisation de l'équation de mesure

Nous posons

$$H(t_k) = H_k \quad (28)$$

$$Z(t_k) = Z_k \quad (29)$$

Dans ces conditions, l'équation (18) devient

$$Z_k = H_k X_k + V_k \quad (30)$$

les $V_k$ formant une suite de vecteurs aléatoires respectant

$$EV_k = 0 \quad (31)$$

$$EV_k V_j^T = R_k \delta_{kj} \quad (\text{les } V_k \text{ succéssifs sont indépendants}) \quad (32)$$

Ces relations traduisent sous forme digitale, les propriétés d'un bruit pseudo-blanc.

En pratique, il arrive fréquemment que les propriétés statistiques des erreurs $V_k$ soient simplement définies par une matrice de covariance, que l'on égale à $R_k$. Par exemple, si les erreurs $(v_{k1}, v_{k2}, \dots, v_{km})$, à l'instant $t_k$, sont indépendantes et d'écarts types $\sigma_{k1}^2, \sigma_{k2}^2, \dots, \sigma_{km}^2$, on a :

$$R_k = \begin{bmatrix} \sigma_{k1}^2 & 0 & \dots & 0 \\ 0 & \sigma_{k2}^2 & \dots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \dots & \sigma_{km}^2 \end{bmatrix} \quad (33)$$

Le processus de digitalisation est résumé sur la figure 1.

<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /></figure>
## 2.3 Entrées, Sorties

Le but du filtre de Kalman est d'obtenir une estimation du vecteur d'état $X = X(t)$, compte tenu des informations disponibles à un instant $t_q$ ($t_q \le t \le t_p$). Cette estimation est notée de la manière suivante : $\hat{X}_{p/q} = P_{p/q}$.

Le filtre fourni de plus les matrices de covariance des erreurs d'estimation

$$E(X_p \hat{X}_{p/q}) (X_p \hat{X}_{p/q})^T = P_{p/q} \quad (34)$$

Ces grandeurs sont calculées pap à pas en fonction des données suivantes:

- matrice de transition $\Phi_{k,k-1}$

- matrice de mesure $H_k$

- estimation initiale du vecteur d'état, notée $\hat{X}_{o/o}$

- covariance de l'erreur d'estimation initiale, notée $P_{o/o}$

- propriétés statistiques des bruits :

- valeurs moyennes nulles

$$EV_k = 0 \quad (31)$$

$$EW_k = 0 \quad (23)$$

- matrices de covariance

$$EW_j W_k^T = 0 \quad (\text{si } j \neq k) \quad (35)$$

$$EW_k W_k^T = Q_k \quad (36)$$

$$EV_j V_k^T = 0 \quad (\text{si } j \neq k) \quad (37)$$

$$EV_k V_k^T = R_k \quad (38)$$

- indépendance des bruits du système et de mesure

$$EW_k V_j^T = 0 \quad (39)$$

- indépendance du bruit de mesure et de l'erreur d'estimation initiale

$$EV_k (X_o - \hat{X}_{o/o})^T = 0 \quad (40)$$

# 3. DESCRIPTION ET INTERPRETATION DU FILTRE DE KALMAN

Nous allons successivement étudier l'estimation de divers types de signaux, de complexité croissante, ce qui nous fera finalement à boutir au filtre de Kalman.
## 3.1 Prévision d'un signal analytique

Considérons tout d'abord un système idéal qui ne serait soumis à aucun bruit ($W_k=0$) et dont le vecteur d'état serait l'objet de mesures parfaites, également sans bruit ($W_k=0$). L'évolution du système serait intégralement fixée, dans ce cas idéal, par son équation de transition (12). La prévision de cette évolution consisterait simplement à appliquer cette équation en utilisant un jeu de mesures commes conditions initiales. C'est la situation que rencontrent les astronomers dans le calcul de la position des planètes.

## 3.2 Prévision d'un signal aléatoire

Considérons maintenant le cas d'un signal Y aléatoire, de valeur moyenne $E(Y)$ nulle, de variance $E(Y^2)$ déterminée. Comment prédire un tel signal? Le problème peut sembler absurde. Il existe cependant une réponse simple: lameiller estimation est simplement

$$ \hat{Y} = 0 \quad (41) $$

Ca résulte, quelque peut décevant, est conforme au bon sens; en effet, d'après ce que l'on sait de Y, on a aucune raison de prédire une valeur positive ou une valeur négative. De plus, l'estimation (41) peut être justifiée en notant qu'elle conduit à une variance d'erreur minimale. On a en effet :

$$ (Y-\hat{Y})=E(Y-\hat{Y})^2=EY^2-2\hat{Y}EY+\hat{Y}^2=EY^2+\hat{Y}^2 \quad (42) $$

La variance $EY^2$ étant fixée, cette quantité est effectivement minimale quand $\hat{Y}=0$.

## 3.3 Estimation d'un signal mixte

Nous en arrivons au vecteur d'état d'un système linéaire régi par des équations de la forme (21,30), soumis à un bruit $W_k$ et objet de mesures $Z_k$.

Supposons que nous disposions d'une estimation $\hat{X}_{k-1/k-1}$ de X relative à l'instant $t_{k-1}$. Que pouvons nous en déduire relativement à l'instant $t_k$? Deux informations :

1) Une estimation de $X(t_{k-1})$, par application de l'équation de transition (14), à l'estimation $\hat{X}_{k-1/k-1}$, conformement à la méthode définie au §3.1. Nous obtenons ainsi :

$$ \hat{X}_{k/k-1} = \hat{\phi}_{k,k} \hat{X}_{k-1/k-1} + F_{k-1} \quad (43) $$

Cette formule, à l'intuition, peut être déduite de (21) et de (41), si l'on admet que l'opération d'estimation est linéaire.

2) Une estimation des résultats de mesures que l'on doit obtenir à l'instant $t_k$ :

$$ \hat{Z}_{k/k-1} = H_k \hat{X}_{k/k-1} \quad (44) $$

Cette formule, elle aussi conforme à l'intuition, traduit qu'une estimation correcte $\hat{X}_{k/k-1}$ doit conduire à une estimation correcte $\hat{Z}_{k/k-1}$.

Le problème se pose maintenant d'améliorer l'estimation $\hat{X}_{k/k-1}$, compte tenu du résultat de mesure $Z_k$. Par souci de simplicité, nous adoptons une correction linéaire de la forme suivante :
$$ \hat{\mathbf{X}}_{k/k} = \hat{\mathbf{X}}_{k/k-1} + K_k (\mathbf{Z}_k - \mathbf{H}_k \hat{\mathbf{X}}_{k/k-1}) \quad (45) $$

La matrice $K_k$ permettant le calcul du terme correctif en fonction de l'écart ($\mathbf{Z}_k - \mathbf{H}_k \hat{\mathbf{X}}_{k/k-1}$), esr appelé matrice de correction. Si l'on veut donner au filtré un caractère véritablement linéaire, cette matrice doit pouvoir être calculée indépendamment des résultats de mesure. Cela sera vérifié au §3.5.

Remplaçons $\mathbf{X}_{k/k-1}$ par son expression (43) dans la relation (45);
nous obtenons:

$$ \hat{\mathbf{X}}_{k/k} = (I - K_k \mathbf{H}_k) (\emptyset_{k,k-1} \hat{\mathbf{X}}_{k-1/k-1} + \mathbf{F}_{k-1}) + K_k \mathbf{Z}_k \quad (46) $$

Cette relation est applicable aux estimations successives $\hat{\mathbf{X}}_{1/1}, \hat{\mathbf{X}}_{2/2}, \dots, \hat{\mathbf{X}}_{k/k}$; l'hypothèse de linéarité traduite par (45) impose donc les conséquences suivantes :

- les estimations succéssives sont des formes linéaires des conditions initiales $\hat{\mathbf{X}}_{0/0}$ et des résultats de mesure $Z_1 Z_2 \dots Z_k$

- pour le calcul de $\hat{\mathbf{X}}_{k/k}$, l'estimation précédente $\hat{\mathbf{X}}_{k/k-1}$ et la mesure nouvelle $Z_k$ constituent un résumé suffisant de $\hat{\mathbf{X}}_{k/k-1}$ et de $Z_1 Z_2 \dots Z_k$. Il en résulte une économie de mémoire, dans le cas d'un calcul au-tomatique.

La forme (45) de la relation de correction étant admise, il reste à expliciter les matrices $K_k$ en fonction des propriétés du système (traduites par les matrices de transition $\emptyset_{k+1,k}$), des propriétés des bruits (traduites par les matrices de covariance $Q_k$ et $R_k$), et enfin de la définition des mesures (traduite par la matrice $\mathbf{H}_k$).

Ce travail sera effectué au §3.5.

### 3.4 Interprétation probabiliste de la méthode d'estimation

Approfondissons le caractère aléatoire des $\mathbf{X}$.

Avant que toute mesure soit faite les résultats de mesures ($Z_1 Z_2 \dots Z_k$) doivent être considérés comme aléatoires. Par suite, les estimations $\hat{\mathbf{X}}$, combinaisons linéaires des $Z_k$, doivent aussi être considérées comme aléatoires a priori. Nous allons d'ailleurs vérifier que si les estimations initiales sont sans biais (1) les estimations suivantes $\hat{\mathbf{X}}_{k/k-1}, \hat{\mathbf{X}}_{k/k-2}, \dots$ sont également sans biais. Ceci implique :

$$ \hat{\mathbf{EX}}_{k/k-1} = \hat{\mathbf{EX}}_{k/k-1} = \hat{\mathbf{EX}}_{k/k} \quad (47) $$

Procédons par récurrence et montrons que

$$ (48) \quad \hat{\mathbf{EX}}_{k-1} = \hat{\mathbf{EX}}_{k-1/k-1} \quad \Rightarrow \quad \hat{\mathbf{EX}}_{k} = \hat{\mathbf{EX}}_{k/k-1} = \hat{\mathbf{EX}}_{k/k} \quad (49) $$

Commençons par vérifier que la valeur moyenne de $(\mathbf{Z}_k - \mathbf{H}_k \hat{\mathbf{X}}_{k/k-1})$ est nulle; nous avons :

$$ \mathbf{E}(\mathbf{Z}_k - \mathbf{H}_k \hat{\mathbf{X}}_{k/k-1}) = \mathbf{E}(\mathbf{H}_k \mathbf{X}_k + \mathbf{V}_k - \mathbf{H}_k \hat{\mathbf{X}}_{k/k-1}) = \mathbf{H}_k \mathbf{E}(\mathbf{X}_k - \hat{\mathbf{X}}_{k/k-1}) $$

mais

$$ \mathbf{X}_k = \emptyset_{k,k-1} \mathbf{X}_{k-1} + \mathbf{F}_{k-1} + \mathbf{W}_{k-1} $$

<footnote>(1) Une estimation $\hat{\mathbf{a}}$ de $\mathbf{a}$ est dite sans biais, lorsqu'elle vérifie $\hat{\mathbf{a}} = \mathbf{a}$.</footnote>
$$ \hat{X}_{k/k-1} = \phi_{k, k-1} \hat{X}_{k-1/k-1} + F_{k-1} $$

par soustraction des deux égalités précédentes, nous obtenons :

$$ X_k - \hat{X}_{k/k-1} = \phi_{k, k-1} (X_{k-1} - \hat{X}_{k-1/k-1}) + W_{k-1} $$

$$ E(X_k - \hat{X}_{k/k-1}) = \phi_{k, k-1} (EX_{k-1} - \hat{EX}_{k-1/k-1}) $$

parsuite, d'après l'hypothèse (48) :

$$ EX_k = \hat{EX}_{k/k-1} \quad (50) $$

et

$$ E(Z_k - H_k \hat{X}_{k/k-1}) = 0 \quad (51) $$

En prenant la valeur moyenne des deux membres de (45), nous obtenons compt tenu de (50) :

$$ \hat{EX}_k/k = \hat{EX}_{k/k-1} = EX_k \quad (c.q.f.d) $$

Dans ces conditions, les écarts

$$ (X_k - \hat{X}_{k/k-1}) \text{ et } (X_k - \hat{X}_{k/k}) $$

sont aléatoires et de valeur moyenne nulle; ceci permet d'exprimer leur matrice de covariance par les expressions suivantes :

$$ P_{k/k-1} = E(X_k - \hat{X}_{k/k-1}) (X_k - \hat{X}_{k/k-1})^T \quad (52) $$

$$ P_{k/k} = E(X_k - \hat{X}_{k/k}) (X_k - \hat{X}_{k/k})^T \quad (53) $$

Ces matrices concernant des écarts a priori doivent pouvoir être calculées indépendamment des résultats de mesure, cela sera vérifié au §3.5.

De plus, la considérationde de $P_{k/k}$ peut nous aider à apprécier l'intérêt d'une méthode d'estimation.Il est en effet possible de tra- duire la valeur d'une telle méthode par la somme des variances ddes composantes de $(X_k - \hat{X}_{k/k})$; meilleur est la méthode, plus faible doit être cette somme. Cette dernière peut d'ailleurs être mise sous la forme :

$$ \text{trace de } P_{k/k} = E(X_k - \hat{X}_{k/k}) (X_k - \hat{X}_{k/k})^T \quad (54) $$

LA METHODE DE KALMAN CONSISTE JUSTEMENT A CHOISIR LA MATRICE DE CORRECTION $K_k$ DE MANIERE A MINIMALISER (54). LA METHODE DE KALMAN EST DONC OPTIMALE, EN CE SENS QU'EL- LE CONDUIT A UNE VARIANCE D'ERREUR MINIMALE;

Le choix des matrices de correction $K_k$ peut aussi être considéré comme un compromis entre :

- la confiance portée à l'estimation précédente $\hat{X}_{k/k-1}$, représentée par $tr P_{k/k-1}$; plus cette confiance est grande (c'est à dire plus $tr P_{k/k-1}$ est petite), moins on se laisse impressionner par un éventuel écart $(Z_k - \hat{Z}_{k/k-1})$, qui peut simplement être dû au bruit de mesure $V_k$. Ainsi, une confiance élevée en $\hat{X}_{k/k-1}$ tend à faire choisir de faibles valeurs pour $K_k$;
- l'intérêt porté aux nouvelles mesures, dû à ce que l'on veut tenir compte :

a) de l'imperfection des estimations précédentes, traduites par tr P<sub>k/k-1</sub>

b) de l'influence du bruit $W_k$ agissant sur le système, d'une manière imprévisible, entre les instants $t_{k-1}$ et $t_k$.

Ce point de vue tend à faire choisir de fortes valeurs pour $K_k$.

Ainsi, le choix explicite de $K_k$ doit être fonction de :

$P_{k/k-1}$, traduisant la confiance dans l'ancienne estimation

$Q_k$, traduisant l'influence du bruit sur le système

$R_k$, traduisant l'influence du bruit sur les mesures.

Ces grandeurs, et par suite $K_k$, sont indépendantes des résultats de mesure.

Calculons la matrice de covariance de l'écart $Z_k - \hat{Z}_{k/k-1}$ entre les valeurs mesurées et les valeurs prévues pour les mesures. Dans le cas d'estimations initiales non biaisées (ce que nous admettons), la valeur moyenne de cet écart est nulle (d'après (51)). D'autre part, en sous-trayant (30), (44) membre à membre, nous obtenons :

$$ (Z_k - \hat{Z}_{k/k-1})^T = (X_k - \hat{X}_{k/k-1})^T H_k^{-1} + V_k $$

en tenant compte de (52,38), nous en tirons :

$$ U_{k/k-1} = E(Z_k - \hat{Z}_{k/k-1}) (Z_k - \hat{Z}_{k/k-1})^T = H_k P_{k/k-1} H_k^{-1} + R_k \quad (55) $$

on obtiendrait de même

$$ U_{k/k} = E(Z_k - \hat{Z}_{k/k}) (Z_k - \hat{Z}_{k/k})^T = H_k P_{k/k} H_k^{-1} + R_k \quad (56) $$

Ces résultats permettent une attitude critique vis-à-vis des modalités d'utilisations du filtre.

En effet, à chaque cycle de calcul, il est possible de comparer chaque composante

$$ \delta_{ki} = Z_{ki} - \hat{Z}_{k/k-1,i} \quad \text{de} \quad Z_k - \hat{Z}_{k/k-1} $$

à l'écart type $\delta_{ki}$ déduit de la matrice $U_{k/k-1}$. Tout écart $\delta$ z supérieur à $\delta_{ki}$ (n étant égal, par exemple, à 3,4,5) traduit l'apparition d'un résultat de mesure a priori très imprable; ce dépassement constitue un signal d'alarme dont il convient de tenir compte; il peut être dû :

- à une mesure aberrante, elle même imputable à la non-transmission d'une donnée, à une panne dans un appareil de mesure, etc... Dans ce cas si le phénomène aberrant ne dure pas, ou n'intéresse qu'une faible partie des mesures, le filtrage peut très bien continuer, à condition, bien entendu que l'on élimine les mesures incriminées;

- à un manque d'exactitude du modèle mathématique utilisé pour représenter le fonctionnement du système étudié. Une telle circonstance peut se produire, en particulier, lorsque l'on a adopté un modèle simple plifié, afin de réduire la complexité des calculs de filtrage. L'écart entre les comportements du modèle et du système peut alors se creuser au fur et à mesure de l'écoulement du temps; il en résulte l'apparition de plus en plus insistante de mesures aberrantes, bien que la qualité des mesures effectuées soit conforme aux prévisions statistiques. On
dit alors que le filtre diverge. Dans ce cas un palliatif consiste à remplacer les matrices $P_k/k$ par des matrices dont les éléments diagonaux (positifs) sont plus élevés, ou encore, à reprendre périodiquement une nouvelle série d'estimations.

Le calcul des $P_k/k$ (pouvant être effectué a priori, c'est-à-dire indépendamment des résultats de mesure) présente à lui seul un intérêt considérable. En effet, il constitue un calcul statistique d'ereur applicable à la méthode d'estimation choisie : les termes diagonaux de $P_k/k$ représentent la confiance que l'on doit accorder aux estimations $\hat{X}_k/k$, compte tenu des perturbations prévues et de la connaissance initiale que l'on a du système. De plus l'évolution des termes diagonaux, en fonction du temps, nous avertit des conditions de fonctionnement du filtre. Par exemple, lorsque l'on recherche la valeur d'une variable d'état constante, la décroissance d'un terme diagonal de $P_k/k$ montre que l'ensemble des informations fournies au filtre est suffisant pour lui permettre d'améliorer les estimations en cours. Au contraire, la non décroissance d'un tel terme montre que le filtre fonctionne dans de mauvaises conditions (voir fig 2).

<figure><img src="image_7.png" /><figcaption>fig. 2</figcaption></figure>

<figure><img src="image_8.png" /></figure>

Cas a) : $\sigma_{\delta Z}$ décroit rapidement : les conditions de fonctionnement du filtre sont satisfaisantes.

Cas b) : $\sigma_{\delta Z}$ ne décroit pratiquement pas : les mesures $Z_k$ ne permettent pas au filtre d'améliorer les estimations initiales $\hat{X}_k/k$ et $P_{k/0}$ (cette circonstance peut d'ailleurs être simplement due à une sous estimation de $P_{k/0}$)

Il est possible de controler la signification et la valeur des matrices $P_k/k$ en effectuant un certain nombre d'applications du filtre à des données calculées par simulation (méthode de Monte-Carlo). A cet effet, des mesures sont calculées de la manière suivante :

$$Z_k \text{ calculé} = H_k \hat{X}_k + V_k$$

$\hat{X}_k$ est obtenu par un calcul théorique déduit du modèle mathématique auquel est supposé obéir le système étudié.

$V_k$ est obtenu à partir d'une table de nombre aléatoires respectant les propriétés prévues pour les bruits de mesure.

A chaque application du filtre, on obtient des erreurs d'estimation $\hat{x}_k$ (composantes de $\hat{X}_k - \hat{X}_{k-1}/k$) qui doivent rester compatibles avec les variances correspondantes $\hat{\sigma}_{\delta x_k}$ (déduites de la diagonale de $P_k/k$) (fig 3).

Dans le cas d'une application à un système réel, il est impossible de comparer les $\hat{\delta}_{x_k}$ aux $\hat{\sigma}_{\delta x_k}$, car on ne connaît pas $\hat{x}_k$. Par contre, on
connait $Z_k$, ce qui permet de comparer les $\delta z_k$ aux variances $\sigma_{z_k}^2$ déduites de la diagonale de $U_k^{-1}$ (fig 4). Cette comparaison constitue un test de validité de l'application du filtre.

<figure><img src="image_9.png" /><figcaption>fig. 3</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>fig. 4
Cas d'une divergence</figcaption></figure>

Nous constatons ainsi que, dans l'exploitation des résultats fournis par le filtre de Kalman, les résultats a priori sont aussi importants que les résultats a posteriori, ces derniers étant constitués par les estimations de $X_k$.

## 3.5 Formules de filtrage

Les formules utilisées sont de type récurrent; à chaque cycle de calcul, on part de $\hat{X}_{k-1/k-1}$ et de $P_{k-1/k-1}$, et on aboutit à $\hat{X}_k/k$ et $P_k/k$. L'organisation du calcul est représentée schématiquement sur le diagramme 5.

<figure><img src="image_11.png" /><figcaption>fig. 5</figcaption></figure>
Examinons les différentes phases :

a) Passage de $P_{k-1/k-1}$ à $P_{k/k-1}$

En soustrayant membre à membre (21) et (43), nous obtenons :

$$X_k - \hat{X}_{k/k-1} = \phi_{k,k-1} (X_{k-1} - \hat{X}_{k-1/k-1}) + W_{k-1}$$

En remplaçant $(X_k - \hat{X}_{k/k-1})$ par cette dernière expression dans (52) nous en tirons :

$$P_{k/k-1} = \phi_{k,k-1} E(X_{k-1} - \hat{X}_{k-1/k-1}) (X_{k-1} - \hat{X}_{k-1/k-1})^T \phi_{k,k-1}^T + EW_{k-1} W_{k-1}^T$$

L'indépendance de $W_{k-1}$ et de $(X_{k-1} - \hat{X}_{k-1/k-1})$ a fait disparaître les valeurs moyennes de leurs produits.

Par définition des matrices P et Q, cette dernière égalité peut être mise sous la forme :

$$P_{k/k-1} = \phi_{k,k-1} P_{k-1/k-1} \phi_{k,k-1}^T + Q_{k-1} \quad (57)$$

L'évolution de la matrice $P_{k/k-1}$ depuis l'instant $t_{k-1}$ jusqu'à l'instant $t_k$, est donc fonction :

- de la dynamique du sytème représentée par $\phi_{k,k-1}$

- du bruit agissant sur le système représenté par $Q_{k-1}$

b) Passage de $P_{k/k-1}$ à $P_{k/k}$

D'après (30) et (45), nous avons :

$$\hat{X}_{k/k} = \hat{X}_{k/k-1} + K_k (H_k X_k + V_k - H_k \hat{X}_{k/k-1})$$

en soustrayant $X_k$ de chacun des membres de cette équation, et après transposition nous obtenons :

$$(\hat{X}_{k/k} - X_k)^T = (\hat{X}_{k/k-1} - X_k)^T (I - K_k H_k)^T + V_k^T K_k^T$$

en prenant la valeur moyenne du produit des deux expressions précédentes nous obtenons :

$$P_{k/k} = (I - K_k H_k) P_{k/k-1} (I - K_k H_k)^T + K_k R_k K_k^T \quad (58)$$

c) Détermination de $K_k$

Kalman a montré que l'adoption du critère somme de variances minimale (soit $tr P_{k/k}$ minimale) implique l'utilisation d'une matrice de correction donnée par

$$K_k = P_{k/k-1} H_k^T (H_k P_{k/k-1} H_k^T + R_k)^{-1} \quad (59)$$

il en résulte que l'expression (58) se réduit à :

$$P_{k/k} = P_{k/k-1} - K_k H_k P_{k/k-1} \quad (60)$$
Les formules (59) et (60) confirment la possibilité de calculer à l'avance les matrices $P_{k/k-1}$, $P_{k/k}$ et $K_k$ indépendamment des résultats de mesure.

De plus la formule (60) implique l'inégalité :

$$ \text{tr } P_{k/k} \le \text{tr } P_{k/k-1} $$

On peur en effet vérifier que le produit $K_k H_k P_{k/k-1}$ est toujours de la forme $BB^T$; par suite, dans l'expression de $\text{tr } P_{k/k}$, ses éléments diagonaux sont tous positifs ou nuls et viennent se déduire de ceux de $P_{k/k-1}$. La relation (60) met ainsi en évidence le gain en précision qu'on attend de l'exploitation de la mesure $Z_k$, pour l'estimation de $X_k$.

La formule (59) peut être écrite sous la forme :

$$ K_k = P_{k/k} H_k^{-1} \qquad (59\text{bis}) $$

mettant en évidence la proportionnalité des termes de $K_k$ à ceux de $P_{k/k}$; interprétation : plus la covariance de $X_k - \hat{X}_k$ (calculable à priori) esr élevée, plus le filtre est "nerveux", pour un $R_k$ donné.

d) Démarrage du calcul

L'application des formules récurrentes précédentes nécessite le choix de valeurs initiales $X_0$ et $P_0/o_0$. Habituellement $X_0$ est considéré comme aléatoire, et on choisit :

$$ \hat{X}_0 = EX_0 \qquad P_{0/o_0} = E(X_0 - \hat{X}_0)(X_0 - \hat{X}_0)^T \qquad (61,62) $$

Les estimations ultérieures $\hat{X}_{k/k-1}$ et $\hat{X}_{k/k}$ sont alors sans biais.

Le fait de choisir des valeurs inexactes pour $P_{0/o_0}$ fausse les évaluations des matrices $P_{k/k-1}$ et $P_{k/k}$; en effet ces dernières ne peuvent être corrigées par les mesures, puisqu'elles en sont indépendantes.

Par contre l'influence d'un mauvais choix pour $X_0$ peut se trouver peu à peu neutralisé par l'exploitation des mesures nouvelles $Z_k$.

# 4. EXTENSIONS DE LA METHODE

## 4.1 Application aux systèmes non linéaires

Le filtre optimal de Kalman est applicable à tout système dont l'évolution est fixée par des équations differentielles linéaires. Ceci lui confère une assez grande généralité d'emploi car il n'existe pas de limitation théorique au nombre de variables d'état utilisées de plus les coefficients des équations du système peuvent varier en fonction du temps.

Cependant il arriver fréquemment que les systèmes à étudier n'entrent pas dans ce cadre. Il est néanmoins souvent possible d'appliquer le filtre de Kalman à un système non linéaire, à condition d'en connaître une solution approchée, dite "nominale" ou de "référence".
En effet, admettons que l'équation d'état puisse être mise sous la forme :

$$ \dot{X} = f(X, B, t) \quad (63) $$

X est le vecteur d'état

B est le vecteur de commande.

Supposons d'autre part que nous nous intéressions à une solution de (63) correspondant aux conditions initiales X dont nous ne connaissons qu'une valeur approchée notée $X_{\text{ref}}$. Nous pouvons obtenir, par résolution de (63), une solution de référence obéissant à

$$ X_{\text{ref}} = f(X_{\text{ref}}, B, t) \quad (64) $$

$$ (X_{\text{ref}})_{t=0} = X_{\text{ref}} \quad (65) $$

Introduisons l'écart entre le vecteur d'état réel et le vecteur d'état de référence :

$$ \Delta X = X - X_{\text{ref}} \quad (66) $$

par soustraction membre à membre de (63) et (64) nous obtenons :

$$ \dot{X} - \dot{X}_{\text{ref}} = f(X, B, t) - f(X_{\text{ref}}, B, t) \quad (67) $$

Si la valeur approchée $X$ est suffisamment voisinne de $X$, nous pouvons espérer que $X_{\text{ref}}(t)$ va rester suffisamment voisin de $X(t)$ pour que l'on puisse faire un développement limité du second membre de (67) et en déduire, compte tenu de (66) :

$$ \frac{d}{dt} \Delta X = \left[ \frac{\partial f}{\partial X} \right]_{\Delta X} \Delta X \quad (68) $$

Nous obtenons ainsi un système d'équation relatif à $\Delta X$, nouveau vecteur d'état.

Ces équations sont linéaires, car les coefficients de la matrice $[\partial f/\partial X]$ sont indépendants de $\Delta X$ (ils ne dépendent que de la solution de référence).

Le filtre de Kalman peut donc être appliqué à ce système, ce qui permet de faire des estimations de $\Delta X$ et finalement de $X$ :

$$ \hat{X} = X_{\text{ref}} + \Delta \hat{X} \quad (69) $$

nb: La solution de référence doit être choisie de manière que la relation (68) reste valable; cette containte peut limiter sa durée d'utilisation.

## 4.2 Détermination de paramètres constants

Il est possible de déterminer la valeur d'un paramètre constant intervenant dans le fonctionnement du système en introduisant :

- une variable d'état supplémentaire : $x_{n+1} = a$

- une équation d'état supplémentaire

$$ \frac{d}{dt} x_{n+1} = 0 \quad (70) $$

La méthode de Kalman peut alors être appliquée comme précédemment.
## 4.3 Cas des bruits colorés

La théorie du filtre de Kalman suppose des bruits blancs, dans le cas continu, ou pseudo-blanc dans le cas digital. Il est cependant possible de traiter des problèmes d'estimation dans lesquels interviennent des bruits colorés, en adjoignant au système étudié des filtres de mise en forme effectuant la transformation

bruit blanc ou pseudo-blanc → bruit coloré

Il en résulte un nouveau système auquel est appliqué le filtre de Kalman (fig. 6).

<figure><img src="image_12.png" /><figcaption>fig. 6</figcaption></figure>

L'application de cet artifice peut augmenter sensiblement la complexité des calculs de filtrage. En effet, le nombre de variables d'équation et le nombre d'équations d'état sont augmentés d'une unité au moins par composante de bruit coloré.

Précisons ces idées par un exemple simple, correspondant aux conditions suivantes :

- les composantes des vecteurs bruit sont indépendantes les unes des autres, elles peuvent donc être traitées séparément;

- ces composantes sont aléatoires et stationnaires; leur auto-corrélation est de la forme

$$ \sigma^2 e^{-|s|/T} $$

$\sigma^2$ étant leur variance,

T étant leur constante de temps de corrélation,

- les intervalles d'échantillonnage ($t_k - t_{k-1}$) sont égaux à $\Delta t$.

Dans ces conditions, d'après les relations AII (44,46), chacune de ces composantes peut être considérée comme obtenue à la sortie d'un filtre ($f$ de t. $1/1+jwF_1$) dont l'entrée est excitée par un bruit blanc de densité spectrale $2T\sigma^2$. L'équation temporelle d'un tel filtre est

$$ c(t) + T \frac{d}{dt} c(t) = b(t) \quad (71) $$

$b(t)$ bruit blanc de densité spectrale $2T\sigma^2$

$c(t)$ bruit coloré.

L'introduction de n bruits colorés implique donc :

- la création de n variables d'état $c_i$

- l'établissement de n nouvelles équations d'état de la forme (71)

- l'application du filtre de Kalman au nouveau système ainsi défini. Lequel est soumis à un bruit blanc de composantes $b(t)$.
## 5. CONCLUSION

La méthode de Kalman et le formalisme qui lui est associé présentent de nombreux avantages :

- la forme récurrente des équations à utiliser facilite leur exploitation par un ordinateur

- la forme des équations de base (21,30) se prête souvent à une exploitation satisfaisante des informations dont on dispose; d'une part, la dynamique du système et la nature des mesures y sont représentées par les matrices $\emptyset$ et $H$ respectivement; d'autre part, les perturbations du système et les bruits de mesure y sont représentés par leurs propriétés statistiques (valeurs moyennes, matrices de covariance)

- les optimisations effectuées ont une signification statistique tout à fait claire : les variances des erreurs sont minimales

- la méthode fournit non seulement une estimation des variables d'état, mais encore, par l'intermédiaire de la matrice de covariance d'erreur $P_k$, une estimation de la confiance que l'on peut accorder à ses résultats.

En contre partie, il faut néanmoins reconnaitre certainnes difficultés pouvant être rencontrées lors de l'application pratique de cette méthode :

- le volume des calculs à effectuer, dû à la croissance du nombre des variables d'état

- l'écart pouvant exister entre le modèle mathématique utilisé et les propriétés réelles du système à étudier, l'influence de cet écart, à la longue, peut conduire à la divergence du filtre, c'est-à-dire des écarts ($Z_k - 2_{k/k}$) invraisemblables, compte tenu des écarts-types déduits de la matrice $U_k$.

Ces deux aspects du problème sont liés; en effet, une représentation plus exacte de la réalité se traduit en général par une augmentation du nombre $n$ de variables d'état utilisées. Mais l'application du filtre de Kalman à l'aide d'un ordinateur conduit à l'utilisation d'un nombre de mémoires croissant comme $n^2$ et à une durée de chaque cycle de calcul croissant comme $n^3$. Il est donc, en pratique du plus haut intérêt de réduire $n$. Cette réduction peut être envisagée de diverses manières :

- en simplifiant la représentation des bruits (passage de bruits colorés à des bruits blancs, suppression du bruit du système)

- en laissant de côté certains phénomènes secondaires, ce qui conduit à l'élimination des variables d'état correspondantes.

On se trouve ainsi en présence d'un compromis à respecter entre la complexité admissible pour les calculs et l'exactitude nécessaire au modèle mathématique adopté. Les ainsi obtenus sont qualifiés de sous-optimaux, bien qu'ils conduisent souvent à de bons résultats.

De plus, en pratique, l'incertitude règne souvent sur les propriétés statistiques des bruits et des conditions initiales à traduire par les matrices $(Q_k, R_k, P_{k/0})$. Il en résulte la mise en service de filtres qualifiés également de sub-optimaux, pouvant être, par exemple, trop timides $(R_k)$ trop élevé, $(P_{k/0})$ et $(Q_k)$ trop faibles) ou trop entre-
prenants($R_k$ trop faible, $P_{o/o}$ et $Q_k$ trop élevés) (fig.7).

<figure><img src="image_13.png" /><figcaption>fig. 7 - Influence de $R_k$ et $Q_k$ sur la recherche d'une constante</figcaption></figure>

Bien entendu, un affichage fantaisiste de $Q_k, R_k, P_{o/o}$ conduit à des estimations $P_{k/k}$ tout aussi fantaistes. Cependant :

1) L'affichage des matrices $P_{o/o}, Q_k, R_k$ peut tolérer une marge d'erreur appréciable (par exemple : 50 %) sans perturber gravement les performances du filtre, en ce qui concerne l'estimation des variables d'état.

2) L'affichage de termes diagonaux trop importants dans $P_{o/o}$ est progressivement corrigé par le filtre lui-même.

3) L'affichage d'estimations ($\hat{X} = EX$) inexactes (par exemple : zéro) est souvent sans importance réelle : le filtre les corrige rapidement en fonction des $Z_k$.
# MODELE DE RADIONAVIGATION INTEGREE OPTIMALE

## 1. MODELE CONTINU

Nous supposons disposer des informations suivantes :

- la vitesse $V_1$ mesurée par un loch électromagnétique

- le cap $C_g$ mesuré avec un compas gyroscopique

- les phases reçues dans un système de radiolocalisation à n lieux.

### 1.1 Correction de l'estime

Une position "estime" $L_e, G_e$ est calculée simplement à partir des informations cap et vitesses :

$$ \frac{d}{dt} G_e = \frac{1}{R_e} V_1 \cdot \sin C_g \quad (1) $$

$$ \frac{d}{dt} L_e = \frac{1}{R_n} V_1 \cdot \cos C_g \quad (2) $$

ou $R_e$ est le rayon de courbure du parallèle en $L_e$, sur l'ellipsoïde, et $R_n$ celui du méridien.

En fait, les valeurs réelles de la vitesse et du cap sont $V_r$ et $C_r$, de plus le navire subit une dérive de composantes Est et Nord notées respectivement $D_e$ et $D_n$. La position réelle évolue donc suivant :

$$ \frac{d}{dt} Gr = \frac{1}{R_e} (V_r \cdot \sin Cr + D_e) \quad (3) $$

$$ \frac{d}{dt} Lr = \frac{1}{R_n} (V_r \cdot \cos Cr + D_n) \quad (4) $$

La correction de position pour passer de l'estime à la réalité est :

$$ \mathcal{S} G = Gr - Ge \quad (5) $$

$$ \mathcal{S} L = Lr - Le \quad (6) $$

d'ou l'évolution de cette correction :

$$ \frac{d}{dt} \mathcal{S} G = \frac{1}{R_e} (V_r \cdot \sin Cr - V_l \cdot \sin C_g + D_e) \quad (7) $$

$$ \frac{d}{dt} \mathcal{S} L = \frac{1}{R_n} (V_r \cdot \cos Cr - V_l \cdot \sin C_g + D_n) \quad (8) $$

Les écarts entre les valeurs réelles et mesurées du cap et de la vitesse sont supposés suffisamment petits pour pouvoir linéariser les équations précédentes. Nous obtenons donc au premier ordre :

$$ \frac{d}{dt} \mathcal{S} G = \frac{1}{R_e} \left[ (V_r - V_l) \cdot \sin C_g + (C_r - C_g) \cdot V_l \cdot \cos C_g + D_e \right] \quad (9) $$
$$ \frac{d}{dt} \delta L = \frac{1}{Rn} \left[ (V_r - V_l) \cos C_g - (C_r - C_g) V_l \sin C_g + D_n \right] \quad (10) $$

## 1.2 Corrections du cap et de la vitesse

Les corrections instantanées aux mesures du cap et de la vitesse sont constituées d'une partie constante, la correction proprement di- te de chaque capteur, et d'une partie aléatoire de moyenne nulle que nous supposons stationnaire. C'est-à-dire :

$$ \begin{aligned} \delta C_g &= \delta C + B_c && (11) \\ \delta V_l &= \delta V + B_v && (12) \end{aligned} $$

avec

$$ \begin{array}{l@{\quad}l@{\quad}l} \delta C = \text{cte} & (13) & \text{corrections du gyro} \\ \delta V = \text{cte} & (14) & \text{et du loch} \\ \delta B_c = 0 & (15) & \\ \delta B_v = 0 & (16) & \text{bruits du gyro} \\ \sigma_{B_c}^2 = \text{cte} = \sigma_C^2 & (17) & \text{et du loch} \\ \sigma_{B_v}^2 = \text{cte} = \sigma_V^2 & (18) & \end{array} $$

Ces bruits sont indépendants entre eux, mais étant dus à l'élec- tronique des capteurs et surtout à la mécanique de leurs principes, pour chacun d'eux les valeurs succéssives sont certainement corrélées. L'hypothèse la plus simple mais retenant cet aspect physique consis- te à admettre une constante de temps de corrélation. Leur évolution est alors décrite par les équations suivantes :

$$ \frac{d}{dt} B_c = - \frac{1}{T_c} B_c + \frac{1}{T_c} W_c \quad (19) $$

$$ \frac{d}{dt} B_v = - \frac{1}{T_v} B_v + \frac{1}{T_v} W_v \quad (20) $$

avec

$T_c$ constante de temps de corrélation du gyro
$T_v$ " " " loch
$W_c$ bruit blanc de densité spectrale $2T_c\sigma_C^2$
$W_v$ " " " $2T_v\sigma_V^2$

## 1.3 Evaluation de la dérive

La dérive est la résultante d'un courant marin et de l'effet du vent :

- le courant marin est supposé constant, ce qui revient à supposer que le temps de réponse du filtre est petit par rapport à celui des variations du courant;

- la dérive due au vent fluctue avec celui-ci autour d'une valeur moyenne qui dépend de la force moyenne du vent et de son gisement par rapport au navire.

Si le navire maintient une même route et si le vent est station-naire, en force et en direction, à l'échelle du temps de réponse du
filtre, alors la dérive peut se décomposer en

$$ \mathrm{De} = \mathrm{Ce} + \mathrm{Fe} \quad (21) $$

$$ \mathrm{Dn} = \mathrm{Cn} + \mathrm{Fn} \quad (22) $$

avec

$$ \begin{aligned} \mathrm{Ce}=\mathrm{cte} & \quad (23) \\ \mathrm{Cn}=\mathrm{cte} & \quad (24) \end{aligned} $$ 
 résultante du courant marin et de la valeur moyenne de la dérive

$$ \begin{aligned} \mathrm{EFe}=0 & \quad (25) \\ \mathrm{EFn}=0 & \quad (26) \end{aligned} $$ 
 fluctuations, supposées stationnaires, de la dérive due au vent

$$ \begin{aligned} \sigma^2_{\mathrm{Fe}} &= \mathrm{cte} = \sigma^2_{\mathrm{e}} & (27) \\ \sigma^2_{\mathrm{Fn}} &= \mathrm{cte} = \sigma^2_{\mathrm{n}} & (28) \end{aligned} $$

Ces fluctuations de la dérive sont les réactions hydrodynamiques du navire aux fluctuations du vent. Si ces dernières sont stationnaires et corrélées, il en est de même des fluctuations de la dérive. Et en admettant encore l'hypothèse la plus simple, leur évolution est décrite par les équations :

$$ \frac{d}{dt} \mathrm{Fe} = - \frac{1}{\mathrm{Te}} \mathrm{Fe} + \frac{1}{\mathrm{Te}} \mathrm{We} \quad (29) $$

$$ \frac{d}{dt} \mathrm{Fn} = - \frac{1}{\mathrm{Tn}} \mathrm{Fn} + \frac{1}{\mathrm{Tn}} \mathrm{Wn} \quad (30) $$

avec

Te et Tn constantes de temps de corrélation, fonctions du vent et du navire

We et Wn bruits blancs de densités spectrales respectives $2\mathrm{Te}\sigma^2_{\mathrm{e}}$ et $2\mathrm{Tn}\sigma^2_{\mathrm{n}}$.

## 1.4 Équations d'évolution

Le système est finalement décrit par les dix variables d'état :

$\delta G, \delta L, \mathrm{Ce}, \mathrm{Cn}, \delta C, \delta V, \mathrm{Fe}, \mathrm{Fn}, \mathrm{Bc}, \mathrm{Bv}$

qui vérifient le système d'équations :

$$ \frac{d}{dt} \delta G = \frac{1}{\mathrm{Re}} \left( \mathrm{Ce} + \delta C \cdot \mathrm{Vl} \cdot \cos \gamma_g + \delta V \cdot \sin \gamma_g + \mathrm{Fe} + \mathrm{Bc} \cdot \mathrm{Vl} \cdot \cos \gamma_g + \mathrm{Bv} \cdot \sin \gamma_g \right) \quad (31) $$

$$ \frac{d}{dt} \delta L = \frac{1}{\mathrm{Rn}} \left( \mathrm{Cn} + \delta C \cdot \mathrm{Vl} \cdot \sin \gamma_g + \delta V \cdot \cos \gamma_g + \mathrm{Fn} - \mathrm{Bc} \cdot \mathrm{Vl} \cdot \sin \gamma_g + \mathrm{Bv} \cdot \cos \gamma_g \right) \quad (32) $$

$$ \frac{d}{dt} \mathrm{Ce} = 0 \quad (33) $$

$$ \frac{d}{dt} \mathrm{Cn} = 0 \quad (34) $$

$$ \frac{d}{dt} \delta C = 0 \quad (35) $$

$$ \frac{d}{dt} \delta V = 0 \quad (36) $$
$$ \frac{d}{dt} Fe = -Ke.Fe + Ke.We \qquad Ke=1/Te \qquad (37) $$

$$ \frac{d}{dt} Fn = -Kn.Fn + Kn.Wn \qquad Kn=1/Tn \qquad (38) $$

$$ \frac{d}{dt} Bc = -Kc.Bc + Kc.Wc \qquad Kc=1/Tc \qquad (39) $$

$$ \frac{d}{dt} Bv = -Kv.Bv + Kv.Wv \qquad Kv=1/Tv \qquad (40) $$

## 1.5 Conditions initiales

Nous supposons partir d'une position estime initiale exacte (port, recalage Transit) et en l'absence d'information sur les courants et sur les biais des capteurs nous admettons que ces quantités sont nulles. Les bruits et les fluctuations, de moyennes nulles, sont initialisés à zéro. Le vecteur estimation initiale est alors le vecteur nul.

Si un étalonnage antérieur du loch et du gyro a déjà fourni des corrections, nous adoptons ces valeurs comme initiales.

Si nous disposons d'une valeur très probable du courant, grace à des documents nautiques par exemple, nous l'affichons aussi.

## 1.6 Équations de mesures

La correction de position est déduite de l'observation des phases: la position estime Le, Ge fournit pour chaque phase une valeur théorique $\varphi_e$, la position réelle $Lr, Gr$ donne une valeur réelle $\varphi_r$, en fait les mesures donnent une valeur $\varphi_m$ affectée d'un bruit $\sqrt{\varphi}$.

L'écart entre la phase mesurée et la phase théorique est :

$$ \delta\varphi = \varphi_m - \varphi_e = (\varphi_r + \sqrt{\varphi}) - \varphi_e \qquad (41) $$

soit en linéarisant

$$ -\delta\varphi = \frac{\partial\varphi_e}{\partial G} \delta G + \frac{\partial\varphi_e}{\partial L} \delta L + \sqrt{\varphi} \qquad (42) $$

Nous supposons que l'étalonnage du système de radiolocalisation, dans la zone de navigation, a fourni les variances $\sigma_{\varphi}^2$ des bruits ainsi que les corrections à apporter aux phases mesurées pour que les bruits soient de moyennes nulles.

Le nombre de phases observées est quelconque mais au moins égal à deux pour pouvoir estimer les corrections en longitude et en latitude. Le système des équations de mesure est alors constitué de n≥2 équations de la forme (42).
## 2. MODELE DISCRETISE

### 2.1 Équations d'évolution

Les variables $C_k$ et $V_k$ sont échantillonnées à la période $\Delta t$. Les valeurs correspondant à l'instant $t_k$ sont notées $C_k$ et $V_k$. Les équations differentielles du §1.4 deviennent les équations aux differences suivantes :

$$ \delta G_{k+1} = \delta G_k + \frac{\Delta t}{Re} (C_{ek} + \delta C_k \cdot V_k \cdot \cos C_k + \delta V_k \cdot \sin C_k + \\ \qquad Fe_k + Bc_k \cdot V_k \cdot \cos C_k + Bv_k \cdot \sin C_k) \quad (31bis) $$

$$ \delta L_{k+1} = \delta L_k + \frac{\Delta t}{Rn} (C_{nk} - \delta C_k \cdot V_k \cdot \sin C_k + \delta V_k \cdot \cos C_k + \\ \qquad Fn_k - Bc_k \cdot V_k \cdot \sin C_k + Bv_k \cdot \cos C_k) \quad (32bis) $$

$$ C_{ek+1} = C_{ek} \quad (33bis) $$

$$ C_{nk+1} = C_{nk} \quad (34bis) $$

$$ \delta C_{k+1} = \delta C_k \quad (35bis) $$

$$ \delta V_{k+1} = \delta V_k \quad (36bis) $$

$$ Fe_{k+1} = (1-Ke \cdot \Delta t) \cdot Fe_k + Ke \cdot \Delta t \cdot We_k \quad (37bis) $$

$$ Fn_{k+1} = (1-Kn \cdot \Delta t) \cdot Fn_k + Kn \cdot \Delta t \cdot Wn_k \quad (38bis) $$

$$ Bc_{k+1} = (1-Kc \cdot \Delta t) \cdot Bc_k + Kc \cdot \Delta t \cdot Wc_k \quad (39bis) $$

$$ Bv_{k+1} = (1-Kv \cdot \Delta t) \cdot Bv_k + Kv \cdot \Delta t \cdot Wv_k \quad (40bis) $$

soit matriciellement :

$$ X_{k+1} = \phi_{k+1/k} X_k + W_k $$

avec

$$ X_k^T = \begin{bmatrix} \delta G_k, \delta L_k, C_{ek}, C_{nk}, \delta C_k, \delta V_k, Fe_k, Fn_k, Bc_k, Bv_k \end{bmatrix} $$
$$ W_k^T = \begin{bmatrix} 0, 0, 0, 0, 0, 0, Ke \cdot \Delta t \cdot We_k, \dots \end{bmatrix} $$

Le bruit blanc $W_k$ qui engendre la fluctuation $Fe$ doit avoir pour densité spectrale $2\pi e^{-\frac{(\delta G_k)^2 + (\delta L_k)^2 + (\delta C_k)^2 + (\delta V_k)^2 + (\delta Fe_k)^2 + (\delta Fn_k)^2 + (\delta Bc_k)^2 + (\delta Bv_k)^2}{2}}$ (§1.3). Le bruit $Ke \cdot W_k$ est donc lui aussi blanc et sa densité spectrale vaut $Ke$ fois la précédente (AII §6.1) c'est-à-dire $2Ke$. Il en est de même pour les autres bruits $Kn \cdot W_n, Kc \cdot Wc$ et $Kv \cdot Wv$.

D'après IV(27bis) nous avons une valeur approchée de la matrice de covariance du bruit $W_k$, qui donne dans notre cas particulier :
$$Q_k = \Delta t \begin{bmatrix} 0 & \cdots & 0 \\ 2K_e\sigma_e^2 & 2K_n\sigma_n^2 & 0 \\ 2K_c\sigma_c^2 & 2K_v\sigma_v^2 & 0 \end{bmatrix}$$

L'expression de la matrice $\phi_{k+1/k}$ est donnée à la fin du chapitre.

## 2.2 Équations de mesure

A l'instant $t_k$ on dispose pour chaque phase $\varphi$ de la valeur mesurée $\varphi_m$ et de la valeur théorique $\varphi_e$, on en déduit leur différence $\delta\varphi$ et le vecteur de mesure :

$$Z_k^T = [\delta \varphi_{1k}, \dots, \delta \varphi_{nk}]$$

qui est lié au vecteur d'état par :

$$Z_k = H_k \cdot X_k + V_k$$

où $H_k$ est la matrice de mesure déduite des équations (42) :

$$H_k = \begin{bmatrix} \frac{\partial \varphi_1}{\partial G} & \frac{\partial \varphi_1}{\partial L} & 0 & \dots & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ \frac{\partial \varphi_n}{\partial G} & \frac{\partial \varphi_n}{\partial L} & 0 & \dots & 0 \end{bmatrix}$$

et $V_k$ est le vecteur bruit :

$$V_k^T = [V \varphi_1 k, \dots, V \varphi_n k]$$

dont la matrice de covariance est d'après IV(33) :

$$R_k = \begin{bmatrix} \sigma_e^2 \varphi_1^2 & & & \\ & \sigma_n^2 \varphi_n^2 & & \\ & & \sigma_e^2 \varphi_1^2 & \\ & & & \sigma_n^2 \varphi_n^2 \end{bmatrix}$$
<table><tbody><tr><td>1</td><td>0</td><td>Δt / Re</td><td>0</td><td>Δt / Re cosCk</td><td>Δt / Re sinCk</td><td>Δt / Re</td><td>0</td><td>Δt / Re cosCk</td><td>Δt / Re sinCk</td></tr><tr><td>0</td><td>1</td><td>0</td><td>Δt / Rn</td><td>Δt / Rn sinCk</td><td>Δt / Rn cosCk</td><td>0</td><td>Δt / Rn</td><td>Δt / Rn sinCk</td><td>Δt / Rn cosCk</td></tr><tr><td>0</td><td>0</td><td>1</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td></tr><tr><td>0</td><td>0</td><td>0</td><td>1</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td></tr><tr><td>0</td><td>0</td><td>0</td><td>0</td><td>1</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td></tr><tr><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>1</td><td>0</td><td>0</td><td>0</td><td>0</td></tr><tr><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>1-KcΔt</td><td>0</td><td>0</td><td>0</td></tr><tr><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>1-KnΔt</td><td>0</td><td>0</td></tr><tr><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>1-KcΔt</td><td>0</td></tr><tr><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>1-KvΔt</td></tr></tbody></table>
# DEVELOPPEMENT DU LOGICIEL

## 1. PRINCIPES GENERAUX

Le developpement d'un système performant de navigation intégrée esr une entreprise délicate par la complexité des études et de la réalisation du logiciel. Il exige naturellement une expérience de la programmation et de la validation des logiciels temps réel. Pour être réussi le développement d'un tel logiciel necessite l'application d'une méthodologie et l'utilisation des moyens de simulation et de validation appropriés.

La figure suivante représente les différentes étapes d'une telle méthodologie :

<figure><img src="image_14.png" /></figure>

Dans le chapitre précédent nous avons réalisé les deux premieres étapes sur l'exemple d'un système de radionavigation.

Au chapitre IV paragraphe 3.4, nous avons déjà signalé l'information disponible avant toute application du filtre à des données. Elle consiste dans le calcul a priori des matrices de covariance $P_{k/k}$ dont les éléments diagonaux représentent la confiance que l'on doit accorder aux estimations $\hat{x}_{k/k}$. En particulier leur variation indique si les mesures permettront au filtre d'améliorer l'estimation initiale. Cette démarche constitue la simulation préliminaire des performances.
L'étude empirique complète des performances et de la sensibilité du filtre aux erreurs d'affichage ($P_o/\sigma^2, R_k, Q_k$), est réalisée en simulant des données à l'entrée du filtre $\delta$ (méthode de Monte-Carlo) et en comparant les écarts entre estimation et valeur vraie de chaque variable d'état à sa variance.

Elle consiste dans le cas particulier de notre exemple à réaliser une simulation dynamique complète, en temps réel, de l'environnement du calculateur embarqué. C'est-à-dire de simuler les sorties des différents capteurs (loch, gyro, récepteur Oméga,...) dans leur forme réelle et pour des trajectoires types.

Enfin l'expérimentation sur des données réelles mais dans des situations connues (par exemple courant mesuré par un courantomètre, loch et gyro étalonnés par une autre méthode) établit ou non l'adéquation du modèle et des propriétés statistiques des bruits à la réalité.

## 2. SIMULATION DU SYSTEME DE RADIONAVIGATION INTEGREE

### 2.1 Route surface réelle

Le choix d'une trajectoire surface (tronçon rectiligne, gyration, ...) et d'une période d'acquisition $\Delta t$ détermine pour chaque instant $t_k=k \cdot \Delta t$ les valeurs réelles du cap et de la vitesse, par exemple :

$$ \text{tronçon rectiligne} \quad Cr_k=cte \quad Vr_k=cte \quad (1) $$

$$ \text{gyration} \quad Cr_k=cte+k \cdot \Delta C \quad Vr_k=cte \quad (2) $$

### 2.2 Route surface mesurée

Pour que les capteurs nécessitent les corrections $\delta C$ et $\delta V$ il faut qu'ils fournissent les valeurs $Cr-\delta C$ et $Vr-\delta V$. À cette route surface biaisée il convient encore d'ajouter le bruit des capteurs.

Nous avons déjà dit que ce bruit est certainement corrélé mais pour simplifier notre modèle nous avons supposé une auto-corrélation à constante de temps ($V.\S 1.2$). Nous devons donc élaborer un bruit de cette forme. Nous allons traiter en détail, à titre d'exemple, le cas du bruit du gyro.

Pour obtenir un bruit $B_c$ de variance $\sigma^2$ et corrélé avec la constante de temps $T_c$ nous devons appliquer la formule $V(39)$ à un bruit blanc $W_c$ de densité spectrale $2T_c\sigma^2$, puis échantillonner $B_c$, ou appliquer la formule $V(39\text{bis})$ à un bruit pseudo-blanc $W_c$ de même densité spectrale.

D'après AII(55), le bruit pseudo-blanc qui a cette densité spectrale doit avoir pour variance :

$$ \sigma^2 = 2 \cdot T_c \cdot \sigma^2 / \Delta t \quad (3) $$
En général un ordinateur fournit seulement une sequence aléatoire uniforme entre 0 et 1 :

$$0 \le u \le 1 \qquad Eu=1/2 \qquad \sigma_u^2=1/12 \qquad (4)$$

le bruit pseudo-blanc est alors obtenu par :

$$Wc_k = \sigma \sqrt{12}(u_k - 1/2) \qquad (5)$$

Finalement le cap mesuré est donné par :

$$C_{g_k} = C_{r_k} - \delta C_k + B_{c_k} \qquad (6)$$

avec

$$B_{c_k+1} = (1-K_c \Delta t)B_{c_k} + 2\sqrt{6} t.K_c \sigma_c(u_k - 1/2) \qquad (7)$$

de même la vitesse mesurée est générée suivant :

$$V_{l_k} = V_{r_k} - \delta V_k + B_{v_k} \qquad (8)$$

avec

$$B_{v_k+1} = (1-K_v \Delta t)B_{v_k} + 2\sqrt{6} t.K_v \sigma_v(u_k - 1/2) \qquad (9)$$

le même symbole $u_k$ désignant un autre tirage de la variable aléatoire uniforme $u$.

## 2.3 Dérive réelle

Le choix d'un courant constant est exprimé par les composantes $C_e$ et $C_n$.

La fluctuation de la dérive est elle aussi supposée corrélée, avec des constantes de temps ($V$ §1.3). Le raisonnement développé au paragraphe précédent est applicable. Le dérive réelle à l'instant $t_k$ est générée suivant :

$$D_{e_k} = C_e + F_{e_k} \qquad (10)$$

$$D_{n_k} = C_n + F_{n_k} \qquad (11)$$

avec

$$F_{e_k+1} = (1-K_e \Delta t)F_{e_k} + 2\sqrt{6} t.K_e \sigma_e(u_k - 1/2) \qquad (12)$$

$$F_{n_k+1} = (1-K_n \Delta t)F_{n_k} + 2\sqrt{6} t.K_n \sigma_n(u_k - 1/2) \qquad (13)$$

## 2.4 Position réelle

La connaissance de la route surface réelle et de la dérive réelle permet d'intégrer les équations $V(3,4)$, d'ou :

$$Gr_{k+1} = Gr_k + \frac{\Delta t}{Re}(V_{r_k} \cdot \sin Cr_k + D_{e_k}) \qquad (14)$$

$$Lr_{k+1} = Lr_k + \frac{\Delta t}{Rn}(V_{r_k} \cdot \cos Cr_k + D_{n_k}) \qquad (15)$$
## 2.5 Phases réelles

La connaissance de la position réelle et des caractéristiques du réseau de radiolocalisation permet de calculer sur l'ellipsoïde de référence les phases théoriques (AI §3).

## 2.6 Phases mesurées

Pour chaque phase $\varphi$ dont nous avons la valeur réelle $\varphi_{rk}$ à l'instant $t_k$, l'étalonnage a fourni une correction $C\varphi$ valable pour la zone de navigation. La phase mesurée subit une bruit $V_0$ supposé blanc de variance $\sigma_\varphi^2$ (V §1.6). La valeur mesurée à l'instant $t_k$ est donc générée suivant :

$$ \varphi_{mk} = \varphi_{rk} - C\varphi + \sigma_\varphi \sqrt{12(u_k - 1/2)} $$
# ELEMENTS DE GEODESIE

## 1. Rayons de courbure de l'ellipsoïde

Un ellipsoïde de révolution est caractérisé par le demi-grand axe a et l'aplatissement

$$ \alpha = (a-b)/a \quad (1) $$

ou b est le demi-petit axe.

Par exemple pour l'ellipsoïde Fischer (1960) :

$$ a = 6378.166\,\text{m} \qquad \alpha = 1/298.3 $$

On en déduit le carré de l'excentricité :

$$ e^2 = (a^2 - b^2)/a^2 = 2\alpha - \alpha^2 \quad (2) $$

Pour un point de l'ellipsoïde, de latitude L, on appelle grande normale la quantité

$$ N = \frac{a}{\sqrt{1 - e^2 \sin^2 L}} \quad (3) $$

c'est la longueur du segment de normale allant du point à l'axe de révolution.

Les rayons de courbure du méridien et du parallèle, du point, sont alors donnés par :

$$ R_n = N \qquad R_e = N \cos L \quad (4) $$

## 2. Géodésique sur l'ellipsoïde

Soient A(La, Ga) et B(Lb, Gb) deux points sur un ellipsoïde de demi-grand axe a et d'aplatissement. La longueur de la géodésique reliant les deux points est donnée approximativement par la formule d'Andoyer.

On calcule d'abord la distance sur la sphère de rayon a :

$$ d = \text{Arc cos}(\sin La \cdot \sin Lb + \cos La \cdot \cos Lb \cdot \cos (Ga - Gb)) \quad (5) $$

puis l'on corrige de l'effet de l'aplatissement :

$$ AB = a \left( d - \frac{\pi}{4} \left( \frac{d - 3 \sin d}{1 + \cos d} \right)^2 + \frac{d + 3 \sin d}{1 - \cos d} \sin^2 d \right) \quad (6) $$

Pour des besoins plus précis il convient d'utiliser la formule de Sodano ($10^{-8}$ près).
### 3. Phases théoriques d'un système de radiolocalisation hyperbolique

La phase théorique en un point M de l'ellipsoïde relative à un couple de stations A-B s'exprime en fonction des distances géodésiques à ces stations.

Par exemple, pour le système Oméga, en repérant le point et les stations dans le système géodésique Mercury utilisant l'ellipsoïde Fischer (1960), la valeur théorique est donnée par la formule :

$$ \varphi_{AB} = \frac{MA - MB}{\lambda} + 900 \quad (7) $$

avec $\lambda=29 468m$, soit une vitesse de propagation de 300 574km/s pour une fréquence de 10,2kHz.
# ELEMENTS DE CALCUL DES PROBABILITES

## 1. DEFINITION DES PROBABILITES

Un événement est dit aléatoire lorsqu'il peut se produire de diverses manières, sans que l'on soit capable de prédire suivant laquelle il se réalisera.

Le fait qu'un événement soit aléatoire n'implique pas automatiquement que l'on ne puisse rien dire à son sujet. Les probabilités ont justement pour but d'exprimer d'une manière quantitative la connaissance que l'on a d'un phénomène aléatoire.

Bien entendu, les probabilités ainsi définies sont fonctions des renseignements dont on dispose; c'est pourquoi, en toute rigueur, il faudrait remplacer l'expression

"la probabilité de l'événement E", notée P(E), par

"la probabilité affectée à l'événement E, compte tenu des renseignements R dont on dispose", notée P(E/R).

## 2. CALCUL DES PROBABILITES

A priori, il peut sembler bien audacieux de vouloir représenter par des nombres notre savoir sur des événements aléatoires. Les créateurs du calcul des probabilités y sont cependant parvenus en admettant les bases suivantes :

- la probabilité d'un événement certain est égale à l'unité;

- la probabilité d'un événement impossible est égale à zéro;

- axiome des probabilités totales: si deux événements a et b sont aléatoires, on a

$$P(a \text{ ou } b) = P(a) + P(b) - P(a \text{ et } b) \quad (1)$$

P(a ou b) est la probabilité affectée à l'événement "a se produit, ou bien b se produit mais a et b ne se produisent pas tous les deux"

P(a et b) est la probabilité affectée à l'événement "a et b se produisent"

- axiome des probabilités composées:

$$P(a \text{ et } b) = P(a) \cdot P(b/a) \quad (2)$$

P(b/a) est la probabilité de b lorsque l'on sait que a s'est produit.
### 3. VARIABLES ALEATOIRES

On appelle variable aléatoire une variable X dont la valeur x ne peut être prédite exactement.

La correspondance entre les valeurs x et les probabilités associées $P(X=x)$ est la loi de probabilité de X.

Le calcul des probabilités a pour but de déterminer les lois de probabilité des v.a. En particulier il arrive fréquemment que l'on puisse déterminer pour une v.a réelle continue X, la probabilité pour qu'elle soit comprise entre x et x+dx. On définit ainsi la densité de probabilité de X, $f(x)$, par la relation :

$$f(x)dx = P(x<X<x+dx) \quad (3)$$

cette densité vérifie la relation

$$\int_a^b f(x)dx = 1 \quad (4)$$

(a,b) définissant la plage des valeurs possibles de X.

Lorsque l'on connaît la loi de probabilité d'une variable aléatoire X, on peut en général calculer :

- sa valeur moyenne ou espérance mathématique

$$EX = \int_{-\infty}^{+\infty} xf(x)dx \quad (5)$$

- sa variance

$$E(X-EX)^2 = \int_{-\infty}^{+\infty} (x-EX)^2 f(x)dx = \sigma_X^2 \quad (6)$$

### 4. SUITES DE VARIABLES ALEATOIRES

Considérons une suite d'événements aléatoires

$E_1$, auquel est attachée la v.a $X_1$

E_n " " " X_n

Deux cas peuvent se présenter :

1) Les v.a $X_1...X_n$ sont indépendantes en probabilité. Cela signifie que la loi de probabilité d'une v.a $X_i$ est indépendante de la valeur prise par une v.a $X_j$ (avec $j \neq i$). On a dans ce cas :

$$P(X_i=x_i) = P(X_i=x_i/X_j=x_j) \quad (7)$$

avec $i \neq j$
2) Les v.a $X_1...X_n$ ne sont pas indépendantes en probabilité. Dans ce cas on peut avoir :

$$P(X_i=x_i) \neq P(X_i=x_i/X_j=x_j) \quad (8)$$

Les notions de lois de probabilités, de densités de probabilité, de moyennes, de variance, s'étendent aux ensembles de N v.a. Par exemple, dans le cas de deux v.a X et Y, on définit :

- la densité du couple X,Y,f(x,y), par la relation :

$$f(x,y)dxdy = P(x<X<x+dx;y<Y<y+dy) \quad (9)$$

- la densité marginale de X :

$$g(x) = \int_{-\infty}^{+\infty} f(x,y)dy \quad (10)$$

- la densité marginale de Y :

$$h(y) = \int_{-\infty}^{+\infty} f(x,y)dx \quad (11)$$

On démontre les propriétés suivantes :

$$1) \text{ Si k est un nombre certain, on a : } E(kX)=kE(X) \quad (12)$$

$$2) \text{ On a toujours : } E(X+Y)=E(X)+E(Y) \quad (13)$$

L'opérateur E faisant passer d'une v.a à sa valeur moyenne est donc linéaire.

$$3) \text{ On a toujours : } \left|E(XY)\right| \leq \sqrt{EX^2 EY^2} \quad (14)$$

$$EXY=EYX \quad (15)$$

$$4) \text{ Si X et Y sont indépendantes, on a : } EXY=EX.EY \quad (16)$$

5) Si X et Y ont des valeurs moyennes nulles, leur coefficient de corrélation

$$r_{XY} = \frac{E(XY)}{\sqrt{EX^2.EY^2}} \quad (17)$$

admet les propriétés suivantes :

$$- r \leq 1 \quad (\text{d'après } (14)) \quad (18)$$

$$- X \text{ et } Y \text{ indépendantes } \Rightarrow r=0 \quad (19)$$

$$- Y=kX \Rightarrow r=1 \quad (20)$$

Ce coefficient de corrélation représente donc, dans une certaine mesure, la liaison en probabilité entre X et Y.

Etant donné une suite de n v.a $X_1...X_n$ de valeur moyenne nulle, on peut lui associer un vecteur aléatoire de dimension n :

$$X^T = [X_1, ..., X_n] \quad (21)$$
On définit alors la matrice de covariance de X :

$$ E(XX^T) = \begin{bmatrix} EX_1^2 & EX_1X_2 & \dots & EX_1X_n \\ EX_2^2 & EX_2^2 & \dots & EX_2X_n \\ \vdots & \vdots & \ddots & \vdots \\ EX_n^2 & EX_nX_2 & \dots & EX_nX_n \end{bmatrix} \quad (22) $$

Cette matrice est symétrique.

Signalons la relation :

$$ \text{Trace de } EXX^T = EX^T X = EX_1^2 + EX_2^2 + \dots + EX_n^2 \quad (23) $$

# 5. FONCTIONS ALEATOIRES

## 5.1 Définition

Lorsqu'une v.a dépend d'un paramètre réel continu t, elle prend le nom de fonction aléatoire (f.a).

D'un point de vue théorique, il serait souhaitable de pouvoir considérer la f.a comme un tout (t variant de - à +); il faudrait connaitre sa loi de probabilité, c'est-à-dire la probabilité d'appré- rition de chaque fonction f(t). Ceci représente une tache considéra- ble; nous limiterons donc notre ambition à l'étude des paramètres suivants, dont nous admettons l'existence, et qui décrivent partiel- lement le f.a :

- valeur moyenne pour chaque valeur de t :

$$ m(t) = E f(t) \quad (24) $$

- fonction auto-corrélation

$$ \varphi(t,s) = E f(t) f(t+s) \quad (25) $$

Cette fonction traduit la dépendance en probabilité entre f(t) et f(t+s). On a en particulier

$$ \varphi(t,0) = E(f^2(t)) \quad (26) $$

qui est la variance de la v.a f(t) lorsqu'elle de moyenne nulle.

## 5.2 F.a stationnaire

Il arrive fréquemment que le phénomène conduisant à l'étude de f(t) présente un caractère de permanence suffisant pour que la valeur moyenne et l'auto-corrélation de f(t) ne dépendent pas de t; on a dans ce cas :

$$ E f(t) = m = cte \quad (27) $$
$$Ef(t)f(t+s)=\varphi(s) \quad (28)$$

la f.a $f(t)$ est dite alors stationnaire du second ordre. Son auto-corrélation vérifie les relations suivantes :

$$\varphi(-s)=\varphi(s) \quad (29)$$

$$\varphi(s) \le \varphi(0) \quad (30)$$

Etant donné deux fonctions aléatoires stationnaires du 2ᵉ ordre, on définit la corrélation mutuelle :

$$\varphi_{12}(s) = \mathrm{E}f_1(t)f_2(t+s) \quad (31)$$

### 5.3 Intégration d'une f.a

Etant donné une f.a $f(t)$, nous admettons qu'il est possible de considérer son intégrale $f(u)$ du comme une f.a $F(t)$.

La permutabilité des opérateurs $E$ (valeur moyenne) et (somme), d'après (13), nous conduit à admettre celle des opérateurs $E$ et $du$. Nous obtenons ainsi la valeur moyenne de $F(t)$ :

$$E(F(t)) = E \int_{0}^{t} f(u) du = \int_{0}^{t} (E f(u)) du \quad (32)$$

Si nous partons d'une f.a stationnaire du second ordre de valeur moyenne nulle alors son intégrale a pour variance :

$$E(F(t))^2 = 2 \int_{0}^{t} (t-u) \varphi(u) du \quad (33)$$

Pour deux f.a stationnaires, $f_1(t)$ et $f_2(t)$, de valeur moyenne nulle et de corrélation mutuelle $\varphi_{12}(s)$, nous avons de même :

$$E(F_1(t)F_2(t)) = 2 \int_0^t (t-u) \varphi_{12}(u) du \quad (34)$$

Lorsque les corrélations sont continues à l'origine et que la durée $t=\Delta t$ est très petite, les expressions (33,34) se réduisent respectivement à :

$$E(F(\Delta t))^2 \approx \varphi(0)(\Delta t)^2 \quad (35)$$

et

$$E(F_1(\Delta t)F_2(\Delta t)) \approx \varphi_{12}(0)(\Delta t)^2 \quad (36)$$

## 6. ANALYSE SPECTRALE DES SIGNAUX ALEATOIRES

### 6.1 Densité spectrale d'un signal continu

Etant donné un signal aléatoire stationnaire $x(t)$ dont la valeur
moyenne est nulle, introduisons la transformée de Fourier $S(\omega)$ de l'auto-corrélation $\varphi_x(s)$ de $x(t)$; nous avons :

$$S(\omega) = \int_{-\infty}^{+\infty} \varphi_x(s) e^{-j\omega s} ds \qquad (37)$$

$$\varphi_x(s) = \frac{1}{2\pi} \int_{-\infty}^{+\infty} S(\omega) e^{j\omega s} d\omega \qquad (38)$$

Nous nous limitons au cas où ces intégrales ont un sens; la fonction $S(\omega)$ ainsi obtenue, jamais négative, est la densité spectrale de puissance moyenne de $x(t)$ : $S(\omega)dw$ représente la puissance de $x(t)$ dans la bande $(\omega, \omega+dw)$. Elle définit une analyse spectrale du signal aléatoire $x(t)$.

Cette propriété est confirmée par les résultats suivants :

A) En faisant $s=0$ dans (38), nous obtenons :

$$\varphi_x(0) = \int_{-\infty}^{+\infty} S(f) df \quad \text{avec } \omega = 2\pi f \qquad (39)$$

mais, par définition de $\varphi_x(s)$, et d'après (26) :

$$\varphi_x(0) = \sigma_x^2 \qquad (40)$$

Ceci est à l'origine de cet abus de langage habituel consistant à appeler $\sigma_x^2$ puissance moyenne de $x(t)$ et de considérer $S(\omega)$ comme une densité de puissance moyenne.

B) Si l'on applique un signal aléatoire stationnaire $x_e(t)$ à un filtre linéaire de f de t $F(\omega)$, il en résulte un signal $x_s(t)$ dont la densité spectrale est

$$S_s(\omega) = S_e(\omega) \left| F(j\omega) \right|^2 \qquad (41)$$

$S_e(\omega)$ étant la densité spectrale de $x_e(t)$.

Signalons enfin que l'on appelle densité spectrale mutuelle de deux signaux aléatoires stationnaires $x_1(t)$ et $x_2(t)$, la transformée de Fourier de leur corrélation mutuelle $\varphi_{12}(s)$; on a :

$$S_{12}(\omega) = \int_{-\infty}^{+\infty} \varphi_{12}(s) e^{-j\omega s} ds \qquad (42)$$

d'ou

$$E x_1(t)x_2(t) = \varphi_{12}(0) = \frac{1}{2\pi} \int_{-\infty}^{+\infty} S_{12}(\omega) d\omega \qquad (43)$$
## 6.2 Notion de bruit blanc

On appelle bruit blanc un signal, idéal, de valeur moyenne nulle dont la densité spectrale de puissance moyenne est une constante positive. Un tel bruit ne peut évidemment pas exister dans un système physique car sa puissance, d'après (39), serait infinie. Par opposition, on dit qu'un bruit est coloré si sa densité spectrale est variable en fonction de .

Le concept de bruit blanc présente néanmoins un intérêt, car il permet de représenter d'une manière simple l'application à un système physique (P) d'un bruit dont la bande de fréquence (pour laquelle $S(\omega)$ conserve une valeur constante) est très grande vis-à-vis de la bande passante de (P) :

<figure><img src="image_15.png" /></figure>

Précisons ces idées sur un exemple :

On applique un bruit blanc (densité spectrale=1) à l'entrée d'un filtre passe-bas du 1er ordre, de f de t

$$F(\omega) = \frac{1}{1 + j\omega\Gamma} \qquad (44)$$

il en résulte à la sortie un signal aléatoire stationnaire $x_s(t)$ caractérisé par :

$$- \text{une densité spectrale} \quad \frac{1}{1+w^2 T^2} \qquad (45)$$

- une auto-corrélation

$$\varphi_s(s) = \frac{1}{2\pi} \int_{-\infty}^{+\infty} \frac{e^{jws}}{1+w^2 T^2} dw = \frac{|s|}{2T} \qquad (46)$$

- une puissance

$$\sigma_s^2 = \varphi_s(0) = \frac{1}{2T} \qquad (47)$$
valeur finie, trés voisine de celle qu'aurait provoquée l'application d'un signal d'entrée $x_e(t)$ à très large bande de fréquence ($\gg 1/T$).

Enfin, on montre que l'auto-corrélation du bruit blanc de densité spectrale $S$ est une impulsion de Dirac d'amplitude $S$ :

$$ \varphi(s) = S\delta(s) \qquad (48) $$

Cette relation signifie que, pour un bruit blanc, deux valeurs $x(t)$ et $x(t+s)$ ne sont pas corrélées, aussi petit que soit $s$.

## 6.3 Densité spectrale d'un signal échantillonné

Soit $x(t)$ un signal continu de densité spectrale $S(\omega)$, le signal échantillonné $x(t_k)$, à la période $\Delta t$, a pour densité spectrale :

$$ S^*(\omega) = \sum_{k=-\infty}^{+\infty} S(\omega + k\omega_e) \qquad (49) $$

avec

$$ \omega_e = 2\pi/\Delta t \qquad (50) $$

Il est possible d'étendre aux signaux aléatoires échantillonnés $x(t_k)$ la relation entre variance et densité spectrale; nous avons en effet :

$$ E(x(t_k))^2 = \sigma_{\text{ech}}^2 = \frac{1}{2\pi} \int_{-\pi/\Delta t}^{\pi/\Delta t} S^*(\omega) d\omega \qquad (51) $$

Pour deux signaux aléatoires stationnaires échantillonnés $x_1(t)$ et $x_2(t)$, on définit la densité spectrale mutuelle par

$$ S_{12}^*(\omega) = \sum_{k=-\infty}^{+\infty} S_{12}^*(\omega + k\omega_e) \qquad (52) $$

Nous avons alors :

$$ \mathrm{E}_{1}(t_{k})x_{2}(t_{k}) = \varphi_{12}(0) = \frac{1}{2\pi} \int_{-\pi/\Delta t}^{\pi/\Delta t} S_{12}^{*}(\omega) d\omega \qquad (53) $$

Considérons maintenant un filtre numérique de f de t $F_d(\omega)$, auquel on applique un signal aléatoire échantillonné $x^*(t)$ de densité spectrale $S^*(\omega)$. Il en résulte à la sortie un signal aléatoire échantillonné $x_s^*(t)$ de densité spectrale

$$ S_s^*(\omega) = S_e^*(\omega) |F_d(\omega)|^2 \qquad (54) $$

## 6.4 Bruit pseudo-blanc

Supposons qu'à chaque instant : $t_k = t_o + k\Delta t$ soit associée une v.a $x_k$ vérifiant : $\mathrm{E}_{x_k}=0$ $\mathrm{E}_{x_k x_j} = \sigma_{\delta_{jk}}^2$ (les $x_k$ sont indépendants les uns des autres).
La suite ainsi définie constitue, dans le domaine des signaux échantillonnés, l'équivalence des bruits blancs étudiés au §6.2 d'ou son appellation : bruit pseudo-blanc.

Connaissant la variance $\sigma^2$ des $x_k$, et admettant pour eux une densité spectrale constante $S^*$, l'application de (51) nous conduit à :

$$ \sigma^2 = \frac{S^*}{\Delta t} \qquad (55) $$