Traitement des signaux : filtre à moyenne mobile exponentielle (EMA)

Précédemment, dans l'article sur l'initiation au traitement des signaux, nous avons présenté deux classes de filtres : à réponse impulsionnelle finie (RIF) ou à réponse impulsionnelle infinie (RII). Nous avons vu comment le filtre à moyenne mobile peut être exprimé sous les formes RIF et RII, mais quels sont les avantages de l'un par rapport à l'autre ?

Si l'on reprend l'exemple utilisé dans mon blog précédent, le filtre RIF développé a la forme suivante :

y[5] = (x[5]+x[4]+x[3]+x[2]+x[1]+x[0])/5,

Ici, nous constatons que nous avons besoin de :

1. 5 opérations de multiplication et
2. 4 opérations de sommation.

Les opérations de multiplication sont particulièrement coûteuses en termes de calcul. Par conséquent, si l'on examine à nouveau la forme RII, on constate qu'elle ne nécessite que :

1. 3 opérations de multiplication et
2. 2 opérations de sommation

y[6]=(x[6]+y[5]-x[1])/5

Cela réduit considérablement le coût de calcul ! C'est un avantage pour les dispositifs embarqués tels que les microcontrôleurs, car ils dépensent moins de ressources à chaque pas de temps discret pour effectuer des calculs.

Par exemple, lorsque j'utilise la fonction Python 'time.time' pour le filtre à moyenne mobile à 11 points sous formes RIF et RII, avec tous les paramètres (taille de fenêtre, fréquence d'échantillonnage, taille d'échantillon, etc.) identiques, j'obtiens les temps d'exécution suivants : 51 ms, 27 ms, respectivement.

Exemple de filtre RII à temps discret

Maintenant que nous savons pourquoi les filtres RII sont plus performants sur les microcontrôleurs, examinons un exemple de projet utilisant une carte Arduino UNO et une unité de mesure inertielle (IMU) MPU6050 de DFRobot (Figure 1). Nous appliquerons le filtre à moyenne mobile exponentielle (EMA) aux données IMU afin d'observer les différences entre les données brutes et les données lissées.

Figure 1 : Schéma fonctionnel de la connexion entre le MPU6050 et l'Arduino Uno. (Source de l'image : Mustahsin Zarif)

Figure 2 : Connexion entre le MPU6050 et l'Arduino Uno. (Source de l'image : Mustahsin Zarif)

Le filtre à moyenne mobile exponentielle est de forme récursive :

y[n] = α*x[n] + (1- α)*y[n-1]

Il est récursif parce que toute sortie actuelle que nous mesurons dépend également des sorties précédentes, c'est-à-dire que le système a une mémoire.

La constante alpha (𝞪) détermine le poids que nous voulons donner à l'entrée actuelle par rapport aux sorties précédentes. Pour plus de clarté, développons l'équation pour obtenir :

y[n] = α*x[n] + (1- α )*(α*x[n−1]+(1−α)*y[n−2])

y[n] = α*x[n] + (1- α )*x[n−1]+α*(1−α)2*x[n−2])+ ...

y[n] = k=0nα*(1−α)k*x[n−k]

Nous constatons que plus la valeur alpha est élevée, plus l'entrée actuelle affecte la sortie actuelle. C'est une bonne chose, car si le système évolue, les valeurs passées ne sont pas aussi représentatives du système actuel. En revanche, cela serait négatif si, par exemple, le système subissait un changement soudain et momentané par rapport à la normale. Dans ce cas, nous voudrions que notre sortie suive la tendance des sorties précédentes.

Sans plus attendre, voyons comment le code d'un filtre EMA fonctionnerait pour le MPU6050.

Code de filtre EMA :

Copier#include <wire.h>
#include <mpu6050.h>

MPU6050 mpu;

#define BUFFER_SIZE 11  // Window size

float accelXBuffer[BUFFER_SIZE];
float accelYBuffer[BUFFER_SIZE];
float accelZBuffer[BUFFER_SIZE];
int bufferCount = 0;  

void setup() {
  Serial.begin(115200);
  Wire.begin();

  mpu.initialize();

  if (!mpu.testConnection()) {
    Serial.println("MPU6050 connection failed!");
    while (1);
  }

  int16_t ax, ay, az;
  for (int i = 0; i < BUFFER_SIZE; i++) {
    mpu.getMotion6(&ax, &ay, &az, NULL, NULL, NULL);
    accelXBuffer[i] = ax / 16384.0;
    accelYBuffer[i] = ay / 16384.0;
    accelZBuffer[i] = az / 16384.0;
  }
  bufferCount = BUFFER_SIZE;
}

void loop() {
  int16_t accelX, accelY, accelZ;

  mpu.getMotion6(&accelX, &accelY, &accelZ, NULL, NULL, NULL);

  float accelX_float = accelX / 16384.0;
  float accelY_float = accelY / 16384.0;
  float accelZ_float = accelZ / 16384.0;

  if (bufferCount < BUFFER_SIZE) {
    accelXBuffer[bufferCount] = accelX_float;
    accelYBuffer[bufferCount] = accelY_float;
    accelZBuffer[bufferCount] = accelZ_float;
    bufferCount++;
  } else {
    for (int i = 1; i < BUFFER_SIZE; i++) {
      accelXBuffer[i - 1] = accelXBuffer[i];
      accelYBuffer[i - 1] = accelYBuffer[i];
      accelZBuffer[i - 1] = accelZBuffer[i];
    }
    accelXBuffer[BUFFER_SIZE - 1] = accelX_float;
    accelYBuffer[BUFFER_SIZE - 1] = accelY_float;
    accelZBuffer[BUFFER_SIZE - 1] = accelZ_float;
  }

//calculate EMA using acceleration values stored in buffer
  float emaAccelX = accelXBuffer[0];
  float emaAccelY = accelYBuffer[0];
  float emaAccelZ = accelZBuffer[0];
  float alpha = 0.2;

  for (int i = 1; i < bufferCount; i++) {
    emaAccelX = alpha * accelXBuffer[i] + (1 - alpha) * emaAccelX;
    emaAccelY = alpha * accelYBuffer[i] + (1 - alpha) * emaAccelY;
    emaAccelZ = alpha * accelZBuffer[i] + (1 - alpha) * emaAccelZ;
  }

  Serial.print(accelX_float); Serial.print(",");
  Serial.print(accelY_float); Serial.print(",");
  Serial.print(accelZ_float); Serial.print(",");
  Serial.print(emaAccelX); Serial.print(",");
  Serial.print(emaAccelY); Serial.print(",");
  Serial.println(emaAccelZ);

  delay(100);
}
</mpu6050.h></wire.h>

Lorsque nous exécutons ce code et vérifions le traceur série, nous pouvons observer des lignes irrégulières et lisses par paires pour les accélérations sur les axes x, y et z, avec une taille de fenêtre de 11 et une valeur alpha de 0,2 (Figures 3 à 5).

Figure 3 : Valeurs d'accélération brutes et filtrées dans la direction x. (Source de l'image : Mustahsin Zarif)

Figure 4 : Valeurs d'accélération brutes et filtrées dans la direction y. (Source de l'image : Mustahsin Zarif)

Figure 5 : Valeurs d'accélération brutes et filtrées dans la direction z. (Source de l'image : Mustahsin Zarif)

Rendre le code plus intelligent

Nous comprenons maintenant pourquoi les filtres RII sont mieux adaptés aux contrôleurs que les filtres RIF, en raison du nombre nettement inférieur de calculs de sommation et de multiplication requis. Cependant, lors de l'implémentation de ce code, la sommation et la multiplication ne sont pas les seuls calculs effectués : nous devons décaler les échantillons chaque fois qu'un nouvel échantillon de temps arrive, et ce processus, en arrière-plan, requiert de la puissance de calcul. Par conséquent, au lieu de décaler tous les échantillons à chaque intervalle de temps d'échantillonnage, nous pouvons utiliser des tampons circulaires.

Voici comment nous procédons : nous avons un pointeur qui mémorise l'index de l'échantillon de données reçu. Ensuite, chaque fois que le pointeur pointe sur le dernier élément dans le tampon, il pointe sur le premier élément du tampon suivant, et les nouvelles données remplacent les données qui étaient stockées ici auparavant, puisqu'il s'agit maintenant des données les plus anciennes dont nous n'avons plus besoin (Figure 6). Par conséquent, cette méthode nous permet de garder la trace de l'échantillon le plus ancien dans le tampon et de le remplacer sans avoir à décaler les échantillons à chaque fois pour placer les nouvelles données dans le dernier élément du tableau.

Figure 6 : Exemple d'illustration d'un tampon circulaire. (Source de l'image : Mustahsin Zarif)

Voici à quoi ressemble le code pour la mise en œuvre d'un filtre EMA avec des tampons circulaires. Pouvez-vous essayer de l'exécuter pour un gyroscope plutôt qu'un accéléromètre ? Jouez également avec les coefficients !

Code pour filtre EMA avec tampons circulaires :

Copier#include <wire.h>

#include <mpu6050.h>

MPU6050 mpu;

#define BUFFER_SIZE 11  // Window size

float accelXBuffer[BUFFER_SIZE];

float accelYBuffer[BUFFER_SIZE];

float accelZBuffer[BUFFER_SIZE];

int bufferIndex = 0;  

void setup() {

  Serial.begin(115200);

  Wire.begin();
 

  mpu.initialize();


  if (!mpu.testConnection()) {

    Serial.println("MPU6050 connection failed!");

    while (1);

  }

  int16_t ax, ay, az;

  for (int i = 0; i < BUFFER_SIZE; i++) {

    mpu.getMotion6(&ax, &ay, &az, NULL, NULL, NULL);

    accelXBuffer[i] = ax / 16384.0;

    accelYBuffer[i] = ay / 16384.0;

    accelZBuffer[i] = az / 16384.0;

  }

}

void loop() {

  int16_t accelX, accelY, accelZ;

  mpu.getMotion6(&accelX, &accelY, &accelZ, NULL, NULL, NULL);

  float accelX_float = accelX / 16384.0;

  float accelY_float = accelY / 16384.0;

  float accelZ_float = accelZ / 16384.0;

  accelXBuffer[bufferIndex] = accelX_float;

  accelYBuffer[bufferIndex] = accelY_float;

  accelZBuffer[bufferIndex] = accelZ_float;

  bufferIndex = (bufferIndex + 1) % BUFFER_SIZE; //circular buffer implementation 

  float emaAccelX = accelXBuffer[bufferIndex];

  float emaAccelY = accelYBuffer[bufferIndex];

  float emaAccelZ = accelZBuffer[bufferIndex];

  float alpha = 0.2;

  for (int i = 1; i < BUFFER_SIZE; i++) {

    int index = (bufferIndex + i) % BUFFER_SIZE;

    emaAccelX = alpha  accelXBuffer[index] + (1 - alpha)  emaAccelX;

    emaAccelY = alpha  accelYBuffer[index] + (1 - alpha)  emaAccelY;

    emaAccelZ = alpha  accelZBuffer[index] + (1 - alpha)  emaAccelZ;

  }

  Serial.print(accelX_float); Serial.print(",");

  Serial.print(emaAccelX); Serial.print(",");

  Serial.print(accelY_float); Serial.print(",");

  Serial.print(emaAccelY); Serial.print(",");

  Serial.print(accelZ_float); Serial.print(",");

  Serial.println(emaAccelZ);

  delay(100);

}
</mpu6050.h></wire.h>

Résumé

Dans ce blog, nous avons abordé la différence entre les filtres RII et RIF en mettant l'accent sur leur efficacité de calcul. En prenant un petit exemple de réduction du nombre d'opérations requises de RIF à RII, nous pouvons imaginer l'efficacité des filtres RII lorsque les applications sont mises à l'échelle, ce qui est important pour les applications en temps réel avec une puissance matérielle limitée.

Nous avons également examiné un projet d'exemple utilisant une carte Arduino Uno et une unité IMU MPU6050, où nous avons déployé un filtre à moyenne mobile exponentielle pour réduire le bruit dans les données des capteurs tout en capturant le comportement du signal sous-jacent. Enfin, dans un souci d'efficacité, nous avons également vu un exemple de code plus intelligent en utilisant des tampons circulaires au lieu de décaler les données à chaque intervalle de temps.

Dans le prochain blog de cette série, nous utiliserons la fonctionnalité FPGA de Red Pitaya pour implémenter un circuit numérique de filtre RIF à 4 prises !