ACTIVITÉ ARDUINO/PYTHON : Etude d’un mouvement d’oscillations avec un accéléromètre MPU6050 (tracé de graphe en temps réel et modélisation)

Merci à l’auteur de cet article qui a été une ressource très précieuse pour la réalisation de cette activité :

http://gilles.thebault.free.fr/spip.php?article32

Objectifs et enjeux

Dans le cadre des nouveaux programmes CPGE, nous avons recherché un moyen de réaliser des acquisitions d’oscillations forcées avec un accéléromètre et un microcontrôleur. Nous souhaitons réaliser un graphe en temps réel et modéliser la courbe pour mesurer l’amplitude et la fréquence des oscillations.

Accéléromètre et carte Arduino

L’accéléromètre-gyroscope utilisé ici est le MPU 6050.

Les branchements sont relativement simples car il suffit de brancher les bornes Vcc, GND, SCL et SDA pour récupérer les 3 composantes de l’accélération et les 3 composantes de l’angle de rotation.

Il est nécessaire d’installer préalablement les bibliothèques I2C et MPU6050 développées par Jeff Rowberg du MIT :

Pour installer ces biblothèques :

  • Ouvrir le logiciel Arduino
  • À partir du menu [Croquis][Inclure une librairie][Ajouter la librairie .ZIP], installer la librairie « I2dev.zip» (disponible sur ce lien : I2dev.zip)
  • À partir du menu [Croquis][Inclure une librairie][Ajouter la librairie .ZIP], installer la librairie « MPU6050.zip» (disponible sur ce lien : MPU6050.zip)

Étalonnage de l’accéléromètre

Pour étalonner l’accéléromètre, téléverser le programme IMU_Zero à partir du menu [Fichier]Exemples][MPU6050]. On laissera l’accéléromètre au repos et à plat.

Ouvrir le moniteur série, puis attendre la fin de l’étalonnage (l’opération prend quelques minutes)

https://www.youtube.com/watch?v=FUc_gd_pTFA

Dispositif expérimental pour l’acquisition

Nous utilisons un vibreur de Melde alimenté par un GBF (modèle avec lequel il est possible de régler précisément la fréquence).

Il faudra utiliser un suiveur de puissance type Geneboost entre le GBF et le vibreur pour éviter les chutes de tension.

Matériel utilisé :

  • GBF Centrad GF265
  • Geneboost Pierron
  • Un ressort de 20 N/m coupé en deux (donc deux ressorts de 40 N/m !)
  • Accéléromètre MPU6050 fixé sur pièce imprimée en 3D qui contient une masse de 20 g.
  • Support imprimé en 3D sur masse de 500 g

Les fichiers STL des pièces 3D sont disponibles ici.

Téléversement du programme Arduino

Le programme téléversé dans la carte Arduino est le suivant, il est librement inspiré du programme de l’article disponible à cette adresse : http://gilles.thebault.free.fr/spip.php?article32

Certaines lignes ont été désactivées (transformées en commentaires avec \\) pour ne garder que les valeurs des 3 composantes de l’accélération

    #include "Wire.h"  // Arduino Wire library
    #include "I2Cdev.h"  //bibliothèque I2Cdev à installer
    #include "MPU6050.h" //bibliothèque MPU6050 à installer
    // AD0 low = 0x68 (default for InvenSense evaluation board)
    // AD0 high = 0x69
    MPU6050 accelgyro;
    long temps;
    int16_t ax, ay, az;  //mesures brutes
    int16_t gx, gy, gz;
     
    void setup() {
      Wire.begin();  // bus I2C
      Serial.begin(9600); // liaison série
      while (!Serial) {
        ; // wait for serial port to connect. Needed for native USB (LEONARDO)
      }
      accelgyro.initialize();  // initialize device
      temps = millis();
      
      
    }
     
    void loop() {
//      accelgyro.getMotion6(&ax, &ay, &az, &gx, &gy, &gz);
      temps = millis();
      // On peut aussi utiliser ces méthodes
      accelgyro.getAcceleration(&ax, &ay, &az);
      //accelgyro.getRotation(&gx, &gy, &gz);
     
      // Affichage accel/gyro x/y/z
//      Serial.print("a/g:\t");
//      Serial.print(ax); 
//      Serial.print("\t");
      Serial.print("t : ");
      Serial.print("\t");
      Serial.print(temps); 
      Serial.print("\t");
      Serial.print("ax : "); 
      Serial.print("\t");
      Serial.print(ax);//    
      Serial.print("\t");
      Serial.print("ay : "); 
      Serial.print("\t");
      Serial.print(ay);//    
      Serial.print("\t");
      Serial.print("az : "); 
      Serial.print("\t");  
      Serial.print(az);      
      Serial.println("\t");

//      Serial.print(gx); 
//      Serial.print("\t");
//      Serial.print(gy); 
//      Serial.print("\t");
//      Serial.print(gz); 
//      Serial.println("\t");
      delay(10);  
    }

https://www.youtube.com/watch?v=Wb6xFCufiT4

Récupération des données avec Python

Une fois ce code Arduino téléversé, nous pouvons récupérer les valeurs avec Python (par biais d’un IDE comme Pyzo, Spyder, IDLE, Sublime Text,..). Il est alors possible de tracer un graphe en temps réel, traiter les données, modéliser les courbes etc…

Dans cet article , je ne rentrerai pas dans les détails pour les étapes de récupération de données. pour cela, je vous invite à consulter l’article qui explique les différentes étapes sur ce lien :

Récupération des données d’une carte Arduino avec Python

Voici le script Python à exécuter :


 #importation des modules
 import serial
 import serial.tools.list_ports # pour la communication avec le port série
 import matplotlib.pyplot as plt  # pour le tracé de graphe
 from matplotlib import animation # pour la figure animée
 import time # gestion du temps
 import numpy as np # numpy pour l'importation des donnees en format txt
 from scipy.optimize import curve_fit
 

 liste_a = [] # liste pour stocker les valeurs de distance
 liste_t = []
 t_acquisition = 10.0 # en s
 amax =2 # en g
 amin= 0 # en g
 

 dt=0.1
 

 

 #pour le graphe en teamax= 3 # en temps réel
 def animate(i):
     line1 = Data.readline()
     print (line1)
     # on retire les caractères d'espacement en début et fin de chaîne
     listeDonnees = line1.strip()
     # on sépare les informations reçues séparées par les espaces et on stocke ces informations dans une liste pour chacune de lignes
     listeDonnees = line1.split()
     print (listeDonnees)
 

 

     if len(listeDonnees) == 12 : # parfois des lignes de données vides peuvent être envoyées, il faut les "écarter"
         accelx = (float(listeDonnees[5].decode()))/16834 # après consulation des données, nous choisissons le 6 ème élément de listeDonnees, on convertit l'accélération en g
         accely = (float(listeDonnees[8].decode()))/16834 # après consulation des données, nous choisissons le 6 ème élément de listeDonnees, on convertit l'accélération en g
         accelz = (float(listeDonnees[11].decode()))/16834 # après consulation des données, nous choisissons le 6 ème élément de listeDonnees, on convertit l'accélération en g
         accel =np.sqrt(accelx**2+accely**2 +accelz**2)
         temps = (float(listeDonnees[2].decode()))/1000.0 # après consulation des données, nous choisissons le 1er élément de listeDonnees
 

         while temps <= t_acquisition:
             liste_a.append(accel)
             print("a = %f"%(accel), " g") # affichage de la valeur de la distance
             liste_t.append(temps)
             print("temps = %f"%(temps), " s") # affichage de la valeur du temps en partant de 0
             line.set_data(liste_t,liste_a)
             return line,
 

 

 

 

 

 

 # Fonction pour la récupération des données série venant de la carte Arduino
 def recup_port_Arduino() :
     ports = list(serial.tools.list_ports.comports())
     for p in ports:
         if 'Arduino' in p.description :
             mData = serial.Serial(p.device,9600)
     print(mData.is_open) # Affiche et vérifie que le port est ouvert
     print(mData.name) # Affiche le nom du port 
     return mData
 

 

 

 

 

 

 Data =recup_port_Arduino() #récupération des données
 

 # Création figure
 fig=plt.figure()
 line, = plt.plot([],[])
 plt.xlim(0, t_acquisition)
 plt.ylim(amin,amax)
 plt.xlabel('temps en s')
 plt.ylabel('a en g')
 plt.grid()
 

 

 #Animation
 ani = animation.FuncAnimation(fig, animate, frames=2000,  interval=20,repeat=False)
 

 plt.show()
 

 plt.close(fig)
 Data.close()
 

 

 

 

 #Ecriture dans un fichier txt
 lines=['t\ta\n'] #première ligne du fichier txt
 for i in range (len (liste_a)):
     line = str(liste_t[i]) +'\t'+ str(liste_a[i])+'\n'
     lines.append(line)
 

 fichier = open('data_accelerometre.txt', 'w').writelines(lines) #création d'un nouveau fichier texte
 

 

 

 t = np.array(liste_t) 
 acc = np.array(liste_a) 
 

 

 # Fonction d'estimation de la fréquence
 def estim_freq(y) : 
     compt = 0
     moy = np.mean(y)
     etat_old = False
     etat_new = False
     for i in range (len(y)) :
         if  y[i] < moy :
             etat_new = True
         else :
             etat_new = False
         if etat_old != etat_new :
             etat_old = etat_new
             compt += 1
 

     return (compt/(2*t_acquisition))
 

 # Fonction d'estimation des valeurs des paramètres de la modélisation
 def get_p0(x, y): 
     
     A0 = (np.max(y)-np.min(y))/2
     f0 =estim_freq(y)
     phase0 =0
     offset0 = np.mean(y)
     
     
     return [A0, f0, phase0,offset0]
 

 def f(x,a,b,c,d):
     return (a*np.sin(2.*np.pi*b*x+c)+d)
 

 popt,pcov = curve_fit (f,t,acc,p0=get_p0(t,acc))
 

 # popt,pcov = curve_fit (f,t,acc)
 

 texte = 'Accélération = '+str(round(float(popt[0]),2))+' sin (2pi*'+str(round(float(popt[1]),2))+'*t+'+str(round(float(popt[2]),2))+') + '+str(round(float(popt[3]),2))+'\n' +'A = '+str(round(float(popt[0]),2))+'; f = '+str(round(float(popt[1]),2))+' ; phase ='+str(round(float(popt[2]),2))+' ; offset = '+str(round(float(popt[3]),2))
 

 # afficher points avec croix rouges. Inserer texte (titre, nom des axes,…)
 plt.figure()
 plt.scatter(t, acc, c = 'red', marker = '+')
 plt.plot(t,f(t,*popt),'g--',label = texte)
 plt.xlabel("t en s")
 plt.ylabel("a en g")
 plt.legend()   # pour afficher les légendes (label)
 plt.show()
 

 

 print (texte)
 
https://youtu.be/W3Gki7TVplk

Quelques précisions concernant la modélisation : la fonction CURVE FIT de scipy

Les lignes de code concernant la modélisation mérite des éclaircissements !

Nous utilisons la fonction curve_fit de scipy (en l’important avec

 from scipy.optimize import curve_fit

Nous définissons une fonction pour le modèle à trouver :

 def f(x,a,b,c,d):
return (a*np.sin(2.*np.pi*b*x+c)+d)

Pour ajuster la courbe par rapport au modèle, il faut écrire :

popt,pcov = curve_fit (f,t,acc)

Seul popt nous intéresse, il s’agit des paramètres d’optimisation de la courbe par rapport au modèle (c’est à dire ici ces 4 valeurs: amplitude (a), fréquence(b), phase(c), offset(d) )

Même si pcov ne nous intéresse pas,nous sommes obligés d’écrire popt,pcov !

… malheureusement cela ne suffit pas 🙁 Il faut « aider » le programme en donnant une estimation de ces paramètres pour que l’optimisation se fasse correctement.

Pour cela il faut définir au préalable une fonction get_p0(x,y) pour estimer ces valeurs :

  • Pour l’amplitude, on divise par deux l’écart entre la valeur minimale et la valeur maximale des accélérations mesurées
  • Par défaut, on définit la phase égale à 0
  • Pour la tension de décalage (ou offset), on calcule la moyenne des valeurs d’accélération mesurées

 def get_p0(x, y): 
     
     A0 = (np.max(y)-np.min(y))/2
     f0 =estim_freq(y)
     phase0 =0
     offset0 = np.mean(y)

Et pour l’estimation de la fréquence ?

On définit (encore !) une petite fonction qui va compter les passages par la valeur moyenne (donc toutes les demi-périodes). Cette fonction va retourner la valeur nombre de comptages/temps d’acquisition … qu’il faut diviser par deux pour une période entière !

Voici cette fonction (qu’il faut définir avant get_p0) :

 def estim_freq(y) : 
     compt = 0
     moy = np.mean(y)
     etat_old = False
     etat_new = False
     for i in range (len(y)) :
         if  y[i] < moy :
             etat_new = True
         else :
             etat_new = False
         if etat_old != etat_new :
             etat_old = etat_new
             compt +=
     return (compt/(2*t_acquisition))

En petit bonus, un script Python permettant de faire la modélisation directement à partir de la dernière acquisition (sans recommencer l’acquisition) à partir des valeurs sauvegardées dans un fichier txt :

import matplotlib.pyplot as plt # pour les graphiques
 import numpy as np # numpy pour l'importation des donnees en format txt
 from scipy.optimize import curve_fit
 importation des donnees txt obtenues apres pointage en supprimant la premiere ligne dans le fichier texte
 lines = open('data_accelerometre.txt').readlines() #on lit les lignes du fichier texte
 open('data_new.txt', 'w').writelines(lines[1:]) #création d'un nouveau fichier texte sans la première ligne
 data = np.loadtxt('data_new.txt')# importation du nouveau fichier texte pour récupérer les valeurs det, x et y dans un tableau
 t = data[:,0] # selection de la premiere colonne
 acc = data[:,1] # selection de la deuxieme colonne
 def estim_freq(y) : 
     compt = 0
     moy = np.mean(y)
     etat_old = False
     etat_new = False
     for i in range (len(y)) :
         if  y[i] < moy :
             etat_new = True
         else :
             etat_new = False
         if etat_old != etat_new :
             etat_old = etat_new
             compt += 1
 return (compt/(2*10.0)) #temps d'acquisition de 10s
 Fonction d'estimation des valeurs des paramètres de la modélisation
 def get_p0(x, y): 
 A0 = (np.max(y)-np.min(y))/2 f0 =estim_freq(y) phase0 =0 offset0 = np.mean(y) return [A0, f0, phase0,offset0]
 def f(x,a,b,c,d):
     return (anp.sin(2.np.pibx+c)+d)
  
 pop,pcov = curve_fit (f,t,acc)
 popt,pcov = curve_fit (f,t,acc,p0=get_p0(t,acc))
 texte = 'Accélération = '+str(round(float(popt[0]),2))+' sin (2pi'+str(round(float(popt[1]),2))+'t+'+str(round(float(popt[2]),2))+') + '+str(round(float(popt[3]),2))+'\n' +'A = '+str(round(float(popt[0]),2))+'; f = '+str(round(float(popt[1]),2))+' ; phase ='+str(round(float(popt[2]),2))+' ; offset = '+str(round(float(popt[3]),2))
 afficher points avec croix rouges. Inserer texte (titre, nom des axes,…)
 plt.figure(1)
 plt.plot(t, acc, c = 'red', marker = '+')
 plt.plot(t,f(t,*popt),'g--',label = texte)
 plt.xlabel("t en s")
 plt.ylabel("a en g")
 plt.legend()   # pour afficher les légendes (label)
 plt.show()
 print (texte)
 print(get_p0(t,acc))
 estim_freq(acc)

Laisser un commentaire

Votre adresse e-mail ne sera pas publiée. Les champs obligatoires sont indiqués avec *