Analyse fréquentielle #7: Le fenêtrage – Implémentation avec Arduino

Objectifs

Savoir implémenter le fenêtrage avec Arduino
Savoir supprimer la composante DC d’un signal
Savoir la forme mathématique des fenêtres basiques
Comprendre la notion du fenêtrage
Savoir l’utilité du fenêtrage
Savoir les différents types des fenêtres
Introduire la carte Portenta H7
Savoir implémenter la DFT en C/Arduino
Test de la DFT sur Arduino Mega/Due
Analyse des performances temporelles du code sur Mega/Due
Analyse de la DFT d’un signal sinusoïdal
Analyse de la DFT d’une entrée réelle
Prendre consciente du problème du sur-échantillonnage
Etc.

Voir le tuto pour plus de détails

La fonction getData()

float getData(float *dataout,int nn, unsigned long ts_us, int analogpin, float gain, int noDC)

Nous avons optimisé la fonction getData() utilisée dans le tuto afin d’obtenir une mesure précise de la fréquence. Elle permet de lire un buffer de données de taille N dans un pin analogique. La fonction peut supprimer ou non la composante DC du signal, ainsi l’ajustement du gain du signal acquis. Elle retourne la valeur approchée de la fréquence d’échantillonnage ainsi le tableau des données. Ci-dessous la définition de la fonction.

float getData(float *dataout,int nn, unsigned long ts_us, int analogpin, float gain, int noDC)
{
  // Lecture des données brutes  
  float ts_f=0.0;
  unsigned long t11; 
  for(int i=0; i<nn; i++)
  {
    t11=micros();
    dataout[i]=(float)analogRead(54+analogpin); 
    delayMicroseconds(ts_us);
    ts_f+=(float)micros()-(float)t11; 
  }
  ts_f=ts_f/(float)nn; 
  
  // Conversion en V
  for(int i=0; i<nn; i++)
    dataout[i]=gain*dataout[i]*3.3/4095.0;
  
  // Suppression de la valeur DC 
  if(noDC==0) 
  {
    float Vm=0.0; 
    for(int i=0; i<nn;i++) Vm+=dataout[i]; 
    Vm/=(float)nn; 
    for(int i=0; i<nn;i++) dataout[i]=dataout[i]-Vm; 
  } 
  return 1E6/ts_f;
}

La fonction setWindow()

float setWindow(float *datain, float *dataout, int n, int type)

La fonction setWindow() permet d’appliquer une fenêtre à un buffer de données. Elle prend en entrée le buffer et sa taille, ainsi le type de la fenêtre :

Type=0: Fenêtre de Hann
Type=1: Fenêtre de Hamming
Type=2: Fenêtre de Blackman
Type=3: Fenêtre Gaussienne
Type=4: Fenêtre Sin²
Type=X: Sans Fenêtre

Les amplitudes sont corrigées par rapport à un coefficient constant. Ci-dessous la définition de la fonction.

float getData(float *dataout,int nn, unsigned long ts_us, int analogpin, float gain, int noDC)
{
  // Lecture des données brutes  
  float ts_f=0.0;
  unsigned long t11; 
  for(int i=0; i<nn; i++)
  {
    t11=micros();
    dataout[i]=(float)analogRead(54+analogpin); 
    delayMicroseconds(ts_us);
    ts_f+=(float)micros()-(float)t11; 
  }
  ts_f=ts_f/(float)nn; 
  
  // Conversion en V
  for(int i=0; i<nn; i++)
    dataout[i]=gain*dataout[i]*3.3/4095.0;
  
  // Suppression de la valeur DC 
  if(noDC==0) 
  {
    float Vm=0.0; 
    for(int i=0; i<nn;i++) Vm+=dataout[i]; 
    Vm/=(float)nn; 
    for(int i=0; i<nn;i++) dataout[i]=dataout[i]-Vm; 
  } 
  return 1E6/ts_f;
}




float setWindow(float *datain, float *dataout, int n, int type)
{
  float Hi, t=0.0,t0=0.0, sig=0.0; 
  const float gain[6]={0.5581,0.5934,0.4765,0.6456,0.5581,1.0000};
    
  switch(type)
  {
    //Fenêtre de Hann
    case 0:
      for(int i=0;i<n; i++)
      {
        t=(float)i/float(n);  
        Hi=0.5-0.5*cos(2.0*PI*t); 
        dataout[i]=Hi*datain[i];
      }
      return gain[type];

    //Fenêtre de Hamming
    case 1:
      for(int i=0;i<n; i++)
      {
        t=(float)i/float(n);  
        Hi=0.54-0.46*cos(2.0*PI*t); 
        dataout[i]=Hi*datain[i];
      }
      return gain[type];

    // Fenêtre de Blackman
    case 2:
      for(int i=0;i<n; i++)
      {
        t=(float)i/float(n);  
        Hi=0.42-0.5*cos(2.0*PI*t)+0.08*cos(4.0*PI*t); 
        dataout[i]=Hi*datain[i];
      }
      return gain[type];
      
    //Fenêtre de Gauss
    case 3:
      for(int i=0;i<n; i++)
      {
        t=((float)i/float(n))-0.5;  
        t0=0.5;
        sig=t0/2.0;  
        Hi=exp(-0.5*t*t/(sig*sig)); 
        dataout[i]=Hi*datain[i];
      }
      return gain[type];

    //Fenêtre sin²
    case 4:
      for(int i=0;i<n; i++)
      {
        t=(float)i/float(n);  
        Hi=sin(2.0*0.5*PI*t)*sin(2.0*0.5*PI*t); // sin(x)²=(1-cos(2x))/2
        dataout[i]=Hi*datain[i];
      }
      return gain[type];

    default:
      for(int i=0;i<n; i++)
      {
        dataout[i]=datain[i];
      } 
      return gain[5];  
  }
}

Programme Complet

#define squareOut   8       // Sortie Carrée

#define N           256     // Données N=2^n, TFD Nf=2^(n-1)
                            // Agit sur la Résolution fréquentielle!!!!!!!! 
                            // df= Fs/N
                             
#define Ts_us       100      // Ts=1/Fs: Approchée (voir le programme) 
                            // Agit sur la fréquence MAX du spéctre!!!!!!!! 
                            // FMXA = Fs/2
/*
  Type=0: Fenêtre de Hann
  Type=1: Fenêtre de Hamming
  Type=2: Fenêtre de Blackman
  Type=3: Fenêtre Gaussienne
  Type=4: Fenêtre Sin²
  Type=X: Sans Fenêtre 
*/


bool State=false; 
unsigned long T1=0, T2=0, T3=0; 
float Ts=0.0, Fs=0.0; 
float Data[N],Data1[N], Nfft=N>>1; 
float dF=0.0; 
float MaxIdMax[2]={0.0,0.0}; 
float Am, ReIm[2]={0.0,0.0};


const float q=1.0/4096.0; 


void setup()
{
  // ADC 
  analogReadResolution(12);       // Ajustable <= 12 (4096 Niveaux) 

  // Sortie 
  analogWrite(squareOut,128);     // 1 KHz
  
  // Affichage 
  Serial.begin(250000); 
}


void loop()
{ 
  /*
   * Analyse Temporelle des Fenêtres
  */
   
  // 1.1 Lecture des données
  Fs=getData(Data, N, Ts_us, 0, 1000.0,0); 
  //Serial.println(Fs); return; 

  // 1.2 Fenêtrage 
  int N0=6;
  float Datai[N0][N], Coefs[N0];
  for(int i=0;i<N0; i++)
    Coefs[i]=setWindow(Data, Datai[i], N, i);
  
  // 1.3 Visualisation
  static bool state=false;
  if(state==false)
  {
    Serial.println("");
    for (int i=0;i<N;i++)
    {
      Serial.print(Datai[0][i]); Serial.print("\t");
      Serial.print(Datai[1][i]); Serial.print("\t");
      Serial.print(Datai[2][i]); Serial.print("\t");
      Serial.print(Datai[3][i]); Serial.print("\t"); 
      Serial.print(Datai[4][i]); Serial.print("\t");
      Serial.println(Datai[5][i]); 
    }
    Serial.println("");
    state=!state;
  }
  return; 



  /*
   * Analyse Fréquentielle des Fenêtres
  */

  // 2.1 DFT 
  float RealFFT[N0][N], ImagFFT[N0][N];
  for(int i=0;i<N0; i++)
    dF=DFT(Datai[i], N, RealFFT[i], ImagFFT[i], Fs);
    
  // 2.2 Affichage
  /*
  static bool state=false;
  float A; 
  if(state==false)
  {
    Serial.println("");
    for (int i=1;i<N>>1;i++)
    {
      Serial.print((i-1)*dF); Serial.print("\t");
      
      A=sqrt(RealFFT[0][i]*RealFFT[0][i]+ImagFFT[0][i]*ImagFFT[0][i])/(float)N; 
      Serial.print(A/Coefs[0]); Serial.print("\t");
      
      A=sqrt(RealFFT[1][i]*RealFFT[1][i]+ImagFFT[1][i]*ImagFFT[1][i])/(float)N; 
      Serial.print(A/Coefs[1]); Serial.print("\t");

      A=sqrt(RealFFT[2][i]*RealFFT[2][i]+ImagFFT[2][i]*ImagFFT[2][i])/(float)N; 
      Serial.print(A/Coefs[2]); Serial.print("\t");

      A=sqrt(RealFFT[3][i]*RealFFT[3][i]+ImagFFT[3][i]*ImagFFT[3][i])/(float)N; 
      Serial.print(A/Coefs[3]); Serial.print("\t");

      A=sqrt(RealFFT[4][i]*RealFFT[4][i]+ImagFFT[4][i]*ImagFFT[4][i])/(float)N; 
      Serial.print(A/Coefs[4]); Serial.print("\t");
      
      A=sqrt(RealFFT[5][i]*RealFFT[5][i]+ImagFFT[5][i]*ImagFFT[5][i])/(float)N; 
      Serial.println(A/Coefs[5]);
    }
    Serial.println("");
    state=!state;
  }
  */
}


float getData(float *dataout,int nn, unsigned long ts_us, int analogpin, float gain, int noDC)
{
  // Lecture des données brutes  
  float ts_f=0.0;
  unsigned long t11; 
  for(int i=0; i<nn; i++)
  {
    t11=micros();
    dataout[i]=(float)analogRead(54+analogpin); 
    delayMicroseconds(ts_us);
    ts_f+=(float)micros()-(float)t11; 
  }
  ts_f=ts_f/(float)nn; 
  
  // Conversion en V
  for(int i=0; i<nn; i++)
    dataout[i]=gain*dataout[i]*3.3/4095.0;
  
  // Suppression de la valeur DC 
  if(noDC==0) 
  {
    float Vm=0.0; 
    for(int i=0; i<nn;i++) Vm+=dataout[i]; 
    Vm/=(float)nn; 
    for(int i=0; i<nn;i++) dataout[i]=dataout[i]-Vm; 
  } 
  return 1E6/ts_f;
}

float setWindow(float *datain, float *dataout, int n, int type)
{
  float Hi, t=0.0,t0=0.0, sig=0.0; 
  const float gain[6]={0.5581,0.5934,0.4765,0.6456,0.5581,1.0000};
    
  switch(type)
  {
    //Fenêtre de Hann
    case 0:
      for(int i=0;i<n; i++)
      {
        t=(float)i/float(n);  
        Hi=0.5-0.5*cos(2.0*PI*t); 
        dataout[i]=Hi*datain[i];
      }
      return gain[type];

    //Fenêtre de Hamming
    case 1:
      for(int i=0;i<n; i++)
      {
        t=(float)i/float(n);  
        Hi=0.54-0.46*cos(2.0*PI*t); 
        dataout[i]=Hi*datain[i];
      }
      return gain[type];

    // Fenêtre de Blackman
    case 2:
      for(int i=0;i<n; i++)
      {
        t=(float)i/float(n);  
        Hi=0.42-0.5*cos(2.0*PI*t)+0.08*cos(4.0*PI*t); 
        dataout[i]=Hi*datain[i];
      }
      return gain[type];
      
    //Fenêtre de Gauss
    case 3:
      for(int i=0;i<n; i++)
      {
        t=((float)i/float(n))-0.5;  
        t0=0.5;
        sig=t0/2.0;  
        Hi=exp(-0.5*t*t/(sig*sig)); 
        dataout[i]=Hi*datain[i];
      }
      return gain[type];

    //Fenêtre sin²
    case 4:
      for(int i=0;i<n; i++)
      {
        t=(float)i/float(n);  
        Hi=sin(2.0*0.5*PI*t)*sin(2.0*0.5*PI*t); // sin(x)²=(1-cos(2x))/2
        dataout[i]=Hi*datain[i];
      }
      return gain[type];

    default:
      for(int i=0;i<n; i++)
      {
        dataout[i]=datain[i];
      } 
      return gain[5];  
  }
}

float DFT(float *data, int Ndata, float *real, float *imag, float fs)
{
  int nfft=floor((float)Ndata/2.0); 
  float df=(fs/2.0)/(float)nfft; 
  float somme_r=0.0, somme_i=0.0;  
  float wt=0.0; 
  
  for (int i=0; i<nfft; i++)
  {
      for (int j=0; j<Ndata; j++)
      {
        wt=-2.0*PI*(float)i*(float)j/Ndata;
        somme_r+=data[j]*cos(wt); 
        somme_i+=data[j]*sin(wt); 
      }      
      real[i]= somme_r; somme_r=0.0;
      imag[i]= somme_i; somme_i=0.0;
  }
  return df; 
}

Objectifs

La fonction getData()

La fonction setWindow()

Programme Complet

Laisser un commentaire Annuler la réponse