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. Le but est d’uitliser Python pour 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

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

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);  
    }

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 d'accélération
 liste_t = []
 t_acquisition = 10.0 # en s
 amax =2 # en g
 amin= -2 # en g
 

 # dt=0.1
 

 

 #pour le graphe 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"
         accel = (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
 

         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 l'accélération
             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')
 fichier.writelines(lines) #création d'un nouveau fichier texte
 

 

 

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

 

 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))
 

 

 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)
 

 Xcalc = np.linspace(0,max(t) , 1024) # création de points pour le tracé du modèle : on crée 1024 points régulièrement espacés entre 0 et la valeur max de I
 

 

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

 # pop,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(Xcalc,f(Xcalc,*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)
 

Petite remarque sur les données provenant de l’accéléromètre

Il s’agit d’un accéléromètre 16 bits. Pour le calibre -/+2 g, il va pouvoir donner 2^16 valeurs comprises entre -32 768 et 32768.

Une valeur de 32768/2 =16384 correspond donc à 1 g.

Pour faire la conversion de l’accélération en g dans le script, il suffit donc de diviser la valeur mesurée par 16384.

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 l’équation du modèle :

 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)

 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 des données de la dernière acquisition (sans recommencer l’acquisition) à partir des valeurs sauvegardées dans un fichier txt (il peut être nécessaire d’indiquer le chemin d’accès au 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/(20)
 


 

 

 

 

 

 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)
 

 Xcalc = np.linspace(0,max(t) , 1024) # création de points pour le tracé du modèle : on crée 1024 points régulièrement espacés entre 0 et la valeur max de I
 

 
 

 Xcalc = np.linspace(0,max(t) , 1024) # création de points pour le tracé du modèle : on crée 1024 points régulièrement espacés entre 0 et la valeur max de I
 

 #
 # popt,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.scatter(t, acc, c = 'red', marker = '+')
 plt.plot(Xcalc,f(Xcalc,*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)
 

 

 

Une vidéo qui récapitule les différentes étapes :

Lien Github des scripts Python et Arduino (et aussi des fichiers stl des pièces 3D) :

https://github.com/jonasforlot/python-arduino/tree/main/Donn%C3%A9es%20s%C3%A9rie%20acc%C3%A9l%C3%A9rom%C3%A8tre%20Arduino

Complément : acquisition avec l’accéléromètre MMA8451

Il est nécessaire d’installer préalablement la bibliothèque Adafruit MMA8451 :

Pour installer cette biblothèque :

Téléversement du programme Arduino

#include <Wire.h>
#include <Adafruit_MMA8451.h>
#include <Adafruit_Sensor.h>

Adafruit_MMA8451 mma = Adafruit_MMA8451();
 long temps;

void setup(void) {
   Serial.begin(9600);
   temps = millis();
 //  Serial.println("Adafruit MMA8451 test!");
 mma.begin();
 mma.setRange(MMA8451_RANGE_2_G);
 }

 void loop() {
   // Read the 'raw' data in 14-bit counts
   temps = millis();
   mma.read();
   Serial.print("t : ");
   Serial.print("\t");
   Serial.print(temps); 
   Serial.print("\t");
   Serial.print("ax : "); 
   Serial.print("\t");
   Serial.print(mma.x); 
   Serial.print("\t");
   Serial.print("ay : "); 
   Serial.print("\t");
   Serial.print(mma.y); 
   Serial.print("\t");
   Serial.print("ay : "); 
   Serial.print("\t"); 
   Serial.print(mma.z); 
   Serial.println("\t");
 Serial.println();
   delay(10);
 }

Récupération des données avec Python

Voici le script Python à exécuter :

Le code Python est (quasiment !) identique

Il s’agit ici d’un accéléromètre 14 bits. Pour le calibre -/+2 g, il va pouvoir donner 2^14 valeurs comprises entre -8192 et 8192.

Une valeur de 8192/2 =4096 correspond donc à 1 g.

Pour faire la conversion de l’accélération en g dans le script, il suffit donc de diviser la valeur mesurée par 4096 (au lieu de 16384).

 #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 d'accélération
 liste_t = []
 t_acquisition = 10.0 # en s
 amax =2 # en g
 amin= -2 # en g
 

 # dt=0.1
 

 

 #pour le graphe 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"
         accel = (float(listeDonnees[5].decode()))/4096 # après consulation des données, nous choisissons le 6 ème élément de listeDonnees, on convertit l'accélération en g
 

         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 l'accélération
             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')
 fichier.writelines(lines) #création d'un nouveau fichier texte
 

 

 

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

 

 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))
 

 

 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)
 

 Xcalc = np.linspace(0,max(t) , 1024) # création de points pour le tracé du modèle : on crée 1024 points régulièrement espacés entre 0 et la valeur max de I
 

 

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

 # pop,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(Xcalc,f(Xcalc,*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)