Vistas de página en total

jueves, 20 de agosto de 2015

FILTRO DE MEDIA MOVIL

  Introducción

   En esta entrada voy a resolver un típico problema que nos surge al leer sensores e intentar mostrar o tomar los datos.
   Normalmente los sensores, tienen un error suficientemente grande como para distorsionarnos la lectura. Por ejemplo, un termopar tipo k puede tener un error de +-2,2ºC. Suficiente, para volvernos locos la lectura.
   Pero eso no significa que mida mal, sino que va a estar oscilando sobre la medida real 2,2 grados.
   En el caso de otros sensores como son presostatos podremos encontrar igualmente un error sobre el fondo de escala de un 1%.
 Posíblemente el error de la sonda y la precisión de la medida no sean iguales. Si leemos un sensor con un 1% de error con un ADC de más de  7 bits estaremos viendo oscilaciones en la medida que no son realmente oscilaciones de presión sino errores en las medidas del sensor.
   Por lo tanto según esto ¿Para qué usar entradas analógicas  de 16 bits para leer un sensor, si todo lo  vamos a leer error?.
   Cuando leamos con 10 o 16 bits, vemos unas oscilaciones como ruido, que nos impedirán determinar exactamente donde se encuentra la sonda. Si visualizamos el dato en un display y lo hacemos rápido, veremos oscilar constantemente el dígito de menos peso, y si  lo hacemos cada segundo estaremos leyendo datos que aunque nos parezcan correctos, tienen un error.

   Planteamiento

  ¿ Entonces que podemos hacer?, ¿Usar solo resoluciones bajas?. Pues no.
    Aunque la sonda tiene un error, por ejemplo 2,2º en el termopar. Este error siempre oscila alrededor de la medida real. por lo tanto el error es como un ruido blanco aleatorio en radio.
   Por ello podremos promediarlo con un filtro de MEDIA MÓVIL, y obtener el valor real, con precisión.

   Si las medidas se realizan de forma que se procesen posteriormente, lo más correcto es tomar las medidas reales y almacenarlas tal cual, para poder aplicar el filtro posteriormente con un programa, de forma que pueda variarse las condiciones del filtro según se desee e incluso probar con varias opciones.
   Pero si lo que queremos es mostrarlo en un display, tenemos un pequeño problema.
  
    En un filtro de media móvil de largo M, la salida actual consiste en el promedio de las últimas M muestras de la entrada.
   La función de transferencia de este filtro es:

  

     M> significa más suavizado, y M< menos suavizado

   Para tomar el valor filtrado de un punto de la onda debemos conocer el valor real de la últimas M medidas.
 

   Solución

   Para que podamos utilizar los datos del display estos deben refrescarse con una cierta cadencia. Es decir debemos mostrar el dato con unos 200ms mínimo de cadencia. Más rápido es tontería, porque no veremos nada.
  Aprovechando este tiempo, normalmente se realiza un bucle de M lecturas para calcular la media y mostrarla en el display, Esto realmente genera un retraso de M lecturas y si el sensor es lento , como es el caso de un sensor de temperatura  DS18B20, o un sensor MAX6675,  la respuesta en el visor es muy lenta.
  Por ejemplo el MAX6675 debemos leerlo cada 200ms si queremos asegurarnos que de la lectura sea correcta. Si queremos un filtro de orden 5 (M=5), esto nos dará un bucle de 1000 ms. Por lo que no podremos refrescar la pantalla en intetervalos inferiores a un segundo. Si el orden es mayor el refresco está más comprometido aún.

  Para evitar este retraso voy a hacer una función que realice la meda sobre una cola. Así podré refrescar el dato en cada lectura del sensor, aun aplicando el filtro. Dado que los datos están en una pila, solo iré añadiendo la última lectura real a la pila en cada lectura al tiempo que realizaré la media con las M-1 lecturas anteriores que ya están en la pila. Así puedo conseguir un refresco de los datos cada 200ms independientemente del número de muestras que tome.
  El único retraso que veré al aumentar M  será el  retraso a  la respuesta a un escalón, que irá aumentando al aumentar M.
         

    Operatividad

   Para realizar el filtro, vamos a realizar una pila de M datos. Iremos empujando los datos que vayamos leyendo sobre la pila. En cada lectura realizaremos el cálculo del valor medio a mostrar. Trasladando los valores de la formula a la pila, y considerando que siempre almaceno el nuevo dato en x[0] empujando los demás hacia el fondo, el valor más antiguo está en x[M] y el más actual en x[0].

  Hay que tener en cuenta que un filtro genera un retraso. El filtro provoca una respuesta más lenta ante escalones bruscos. Esto es una desventaja y una ventaja.  El fin precisamente del filtro es eliminar los picos y escalones de frecuencia alta, que no dan datos, sino ruido y tomar una media con ellos para mostrar un valor más real. Pero si necesitamos mostrar  datos en una señal con rápidos cambios, puede que el filtro nos genere un retraso adicional. Este retraso se verá en el encendido inicial o cuando la sonda se conecte o desconecte, pues en estos casos el salto es muy elevado.
  Como el valor mostrado es la media, esta variará de forma progresiva desde el valor inicial al valor final al producirse un escalón. Si El escalón resultase muy alto, aparecerán medidas de valores intermedias durante algunos segundos. Pero esto ya es algo que estamos acostumbrados a ver en los sistemas comerciales.
   En cambio veremos la medida con mucha más precisión y estabilidad, eliminando reduciendo drásticamente el error del sensor.
   Por lo tanto, hay que valorar el parámetro M del filtro. Cuanto mayor sea más estable se verá la medida, y mejor filtrará el ruido, pero solo servirá para sistemas con alta inercia. Si el sistema puede variar sus parámetros de forma rápida habrá que reducir M a 3 o menos, para que el sistema tenga una respuesta más adecuada.
   En este gráfico puede verse la diferencia entre M altos y bajos.

  Gráfico tomado de "Introducción a los Filtros Digitales clase 10"

Realizaré el código en C para Pic.


Código 1

 

////////////////////////////////////////////////////////////////////
// Librer¡a Filtro Media Movil                                    //
// compilador PIC CCS                                             //
// Realizado por Jose Angel Moneo                                 //
// Funciones:                                                     //
//  float MediaMovil(int ntabla,float valor)  
//    devuelve el valor medio calculado al almacenar un nuevo dato //
//    El valor se calcula y maneja con la tabla ntabla            //
////////////////////////////////////////////////////////////////////


 //*********************** Definición de pila*********************

//Estas constantes deben ser definidas en el programa principal, pues podrá manejarse varias pilas por el programa.


//#define NSensor_MediaMovil 1 //defineimos el número de sensores

//#define Muestras_MediaMovil 10 //definimos el número de muestras de la media movil


//En caso de no haber definidio las tablas se define 1 con 10 muestras

#ifndef NSensor_MediaMovil
   #define NSensor_MediaMovil 1
#endif


#ifndef Muestras_MediaMovil
  #define Muestras_MediaMovil 5
#endif

//Pilas de datos para el cálculo
float  PilaSensor[NSensor_MediaMovil][Muestras_MediaMovil];


//Devuelve un nuevo valor de la media al añadir un nuevo dato a la pila
float MediaMovil(int ntabla,float valor)
{ int i;
  float Ynuevo;
  Ynuevo=0;
   //Hace hueco moviendo los datos hacia el fondo y suma los datos para la media
  for(i=Muestras_MediaMovil-1;i>0;i--)
       {PilaSensor[ntabla][i] = PilaSensor[ntabla][i-1];
        Ynuevo+=PilaSensor[ntabla][i];}
   //Empuja el nuevo dato en la pila y suma el nuevo dato
   PilaSensor[ntabla][0]=valor;
   Ynuevo+=PilaSensor[ntabla][0];
   //Calcula el nuevo valor medio
   Ynuevo/=Muestras_MediaMovil;  //calcula la media
   return Ynuevo;

}


Código 2. El mismo algoritmo pero con cola circular mediante punteros (más rapido)

////////////////////////////////////////////////////////////////////
// Librer¡a Filtro Media Movil  V2                                //
// compilador PIC CCS                                             //
// Realizado por Jose Angel Moneo                                 //
// Funciones:                                                     //
//  float MediaMovil(int ntabla,float valor)  
//    devuelve el valor medio calculado al almacenar un nuevo dato //
//    El valor se calcula y maneja con la tabla ntabla            //
////////////////////////////////////////////////////////////////////



 //*********************** Definición de pila*********************

//Estas constantes deben ser definidas en el programa principal, pues podrá manejarse varias pilas por el programa.


//#define NSensor_MediaMovil 1 //defineimos el número de sensores

//#define Muestras_MediaMovil 3 //definimos el número de muestras de la media movil


//En caso de no haber definidio las tablas se define 1 con 10 muestras

#ifndef NSensor_MediaMovil
   #define NSensor_MediaMovil 1
#endif


#ifndef Muestras_MediaMovil
  #define Muestras_MediaMovil 3
#endif



//Pilas de datos para el cálculo
int16  ColaSensor[NSensor_MediaMovil][Muestras_MediaMovil];  //Tabla historico últimos valores
int16 VADC[NSensor_MediaMovil];                  //Suma Valores de las entradas analógicas pasadas por filtro media movil
int antiguo=0;
int nuevo=0;


//Devuelve un nuevo valor de la media al añadir un nuevo dato a la pila
int16 MediaMovil(int ntabla,int16 valor)
{    // Calcula la nueva suma
  VADC[NSensor_MediaMovil]-=VADC[ntabla]-ColaSensor[ntabla][antiguo]; //Resta el valor más antiguo
  VADC[NSensor_MediaMovil]+=valor;  //Calcula la nueva suma

    ColaSensor[ntabla][nuevo]=valor;  //Almacena el nuevo valor en la cola

   //Actualiza los punteros de cola    
    if(antiguo++>Muestras_MediaMovil)  antiguo=0;//Hace la cola circular  
    if(nuevo++>Muestras_MediaMovil)  nuevo=0;//Hace la cola circular    
   //Devuelve el nuevo valor medio
   return (VADC[NSensor_MediaMovil]/Muestras_MediaMovil);


}







 Implementación


   Para implementarlo en el programa principal donde se toman los datos, basta con llamar a la función con la lectura actual.

     Ejemplo de implementación.

     Voy a usar el filtro para estabilizar y mejorar la lectura de tres sensor termopar tipo k leído con un  MAX6675 con PIC.
      Para leer el sensor a través del MAX6675 lo haremos por SPI por lo que la instrucción de lectura usaremos la función SPI_XFER
    Configuraremos la lectura de SPI de los tres sensores mediante

  #use spi(MASTER, CLK=PIN_C7, DI=PIN_C6, ENABLE=PIN_C1, MODE=0, BITS=16, STREAM=Sensor1)

    Cambiando el ENABLE de cada MAX6675 y el STREAM por el nombre que deseemos.

   Para leer el valor del sensor usaríamos SPI_XFER(Sensor1,0)>>3, despreciando los bits de control, y para obtener su valor en grados lo multiplicaríamos por 0.25.
  Por lo tanto, leer su valor en grados sería (SPI_XFER(Sensor1,0)>>3)*0.25.
   Esto lo podríamos mostrar ya en el display, pero veremos como el valor puede oscilar más de un grado en cada lectura.

   Para corregirlo vamos a implementar la función de media móvil. 
  Primero incluimos la función en el programa y definimos los parámetros que deseamos.

#define NSensor_MediaMovil 3
#define Muestras_MediaMovil 3
#include "Filtro_MediaMovil.c"

   En este caso definimos 3 tablas, una por sensor y um M de 3, para que la respuesta a cambios sea rápida.

   Para tomar el dato vamos a pasar todas las lecturas por la media movil. Para ello usaremos la instrucicón:

 T_Sensor1= MediaMovil(0,(float)(SPI_XFER(Sensor1,0)>>3)*0.25); 

  Esto nos leerá el sensor, lo convertirá en grados y lo almacenará en la tabla de la media móvil ya en coma flotante. Así mismo la función MediaMovil, nos devolverá el valor filtrado en coma flotante, el cual almacenamos en una variable en coma flotante llamada T_Sensor1. que posteriormente usaremos es para mostrar en el display.
  Esto mismo lo haremos con todos los sensores.

 T_Sensor1= MediaMovil(0,(float)(SPI_XFER(Sensor1,0)>>3)*0.25); 
 T_Sensor2= MediaMovil(0,(float)(SPI_XFER(Sensor2,0)>>3)*0.25); 
 T_Sensor3= MediaMovil(0,(float)(SPI_XFER(Sensor3,0)>>3)*0.25); 


Ya está... ahora lo mostraremos donde queramos.



3 comentarios :

  1. Existe una mejora substancial al código en dos aspectos:
    1. Si el vector es importante en micros de poca RAM eso se complica.
    2. En el mismo orden con una cantidad de muestras importantes, se sigue consumiendo tiempo.

    Cuando yo estudié el promedio movil lo vi desde el punto de vista matemático, tal como tu lo has planteado.
    Entonces la suma de elementos puede considerarse como la suma de cada elemento divido la cantidad de muestras. Llamaré PS a PilaSensor y M a Muestras_MediaMovil.
    Bastaría en todo momento con guardar PS[0]/M para en la siguiente iteracción, tomar una nueva muestra, sumarla ponderada y restar la primera de la cola que como dije llamé PS[0]/M pero solo fue a titulo indicativo porque ni siquiera requiero el vector por lo que libero la memoria usada.
    Asi que en todo momento llevo la sumatoria, guardo el ultimo primer item de la cola ponderado y agrego el nuevo valor muestreado y lo pondero.
    Reducción de tiempo y memoria fantásticos a costo 0 con toda la ventaja del promedio móvil.

    ResponderEliminar
    Respuestas
    1. Muchas gracias Ricardo. No conocía la media ponderada.... Un problema de no haber dado estadística en mis tiempos en la carrera, y de haberme dedicado a otros campos... pero siempre hay algo que aprender.
      Muchas gracias. La usaré a partir de ahora. Gracias a tu comentario he conocido la Media movil ponderada y la media ponderada exponencial.... Aunque en los dos casos seguir´na precisando de una cola de valores almacenados. De todas formas he añadido una versión más rápida de código mediante el uso de punteros, que permite usar la cola sin mover los datos, para conocer el más nuevo y el más antiguo.
      No entiendo de todas formas tu comentario, pues por lo que he visto, en la ponderación sigues necesitando la cola, para conocer y descontar el valor más antiguo. Tu simplemente usas M como valor de ponderación igual en todos los casos. Al almacenar el valor ponderado, en vez del raw, puedes perder precisión en la suma, si lo haces sin coma flotante. En muchos casos yo quiero hacer una media de los valores raw, de las señales analógicas, almacenadas en int. Así ahorro espacio de cola. Ese es el segundo algoritmo.
      Por otro lado lo interesante de la media móvil ponderada puede ser precisamente dar más peso a los valores más cercanos.
      Lo plantearé ahora que conozco estos algoritmos. Gracias. ;-)

      Eliminar
    2. Por lo tanto tu algoritmo y el segundo que he puesto deberían de ser idénticos. y para poder realizar una media movil ponderada o una media movil ponderada exponencial, habría que recurrir a un bucle de sumas ponderando en cada caso, pues los valores de ponderación deberá ir adaptándose al a posición que cada valor ocupa en la tabla. Por lo que es más simple partir del primer código y en la suma incluida en el código añadir el factor de ponderación correspondiente. Dado que en este momento, lo que menos importa es la potencia del procesador... por 5€ puedes comprar un PIC32, realmente el algoritmo no es crítico. Es más crítico la capacidad del algoritmo, y en este caso reconozco que una algoritmo de Media Movil ponderado con la ponderación puesta a i/Muestras_MediaMovil, donde i es el indice de la tabla es mucho más potente que una media movil simple, pues permite reaccionar más deprisa a los cambios.

      Eliminar