Medidor de energía doméstico simple Arduino


Michaele Klements nos apunta que para obtener mediciones perfectamente precisas, necesitamos monitorear tanto el voltaje de suministro como la corriente, pero para un monitoreo doméstico simple que brinde estimaciones de costos al mínimo de centavos, ¿por qué no simplificar las cosas?

Entonces, este simple medidor propuesto por Michaele mide la corriente de suministro a su hogar a través de un CT (transformador de corriente) y luego hace un par de cálculos para darle su corriente, potencia, potencia máxima y kilovatios hora consumidos. También es muy fácil agregar su tarifa local y mostrar el costo de la electricidad utilizada hasta la fecha.

Es evidente además que este circuito puede ser la base para volcar los datos a una BBDD o simplemente usando una plataforma de IoT como por ejemplo Cayenne , poder consultas las medidas online , obtener el histórico, realizar acciones correctoras , etc

Los mínimos componentes que necesitará para un medidor de energía doméstico de este tipo son los siguientes:

  • Arduino Uno ( o similar ) 
  • Pantalla LCD (o pantalla LCD-  
  • CT – Talema AC1030
  • Resistencia de carga de 56Ω
  • Condensador de 10 µF
  • 2 resistencias divisoras de 100K

En el siguiente video podemos ver el montaje del circuito así como una pequeña demostración del proyecto en funcionamiento

Cómo hacer el medidor de energía

Primero, debe comenzar ensamblando los componentes en el CT o en su placa de pruebas para crear su sensor actual que produce una señal que su Arduino puede entender. Un Arduino solo tiene entradas de voltaje analógicas que miden 0-5 V CC, por lo que debe convertir la salida de corriente del CT en una referencia de voltaje y luego escalar la referencia de voltaje en un rango de 0-5 V.

Ensamblaje

Si va a instalar su medidor de potencia en algún lugar de forma permanente, es posible que desee soldar las resistencias y el condensador directamente en el TC para que no se suelten. Si simplemente está probando este proyecto por diversión, entonces una placa de pruebas es perfecta.

El circuito básico para la conexión del CT al Arduino se muestra en el siguiente diagrama:

energy-meter-wiring-diagram

A continuación se muestra un diseño de circuito de tablero. Tenga en cuenta que TinkerCAD no es compatible con un transformador de corriente. Por tanto, se ha utilizado un generador de señales para generar una señal de ejemplo:

Arduino home energy meter

El blindaje de la pantalla LCD ya se activa en las entradas analógicas, pero el blindaje solo utiliza A0. Simplemente suelde los tres cables de su sensor de corriente en los encabezados de los pines en el escudo y use A1 como su entrada de sensor como se muestra a continuación.

current-sensor-connections

Una vez que haya conectado todos sus componentes, debe conectar su sensor a lo que desea monitorear. Si desea monitorear un par de aparatos, entonces debe conectar el CT al cable de entrada de un enchufe múltiple, todo lo que conecte al enchufe múltiple luego se contará.

Alternativamente, puede conectar el CT directamente a la red eléctrica de su hogar y monitorear el uso de toda la casa como se ha hecho aquí. De cualquier manera, debe colocar el CT alrededor de uno de los cables de suministro, preferiblemente el cable rojo «vivo». Asegúrese de ponerlo solo alrededor de 1, ya que no funcionará si está alrededor de ambos y no se puede conectar alrededor del cable de tierra (cable amarillo, verde pelado) ya que la energía no se extrae a través de este cable. Si lo está conectando a la red, conéctelo a uno de los cables de salida después del interruptor principal como se muestra a continuación.

NB: tenga cuidado al conectar el medidor de potencia a la red eléctrica de su hogar y asegúrese de que la alimentación de su placa esté apagada antes de hacer nada en la caja de red. No retire ningún cable ni quite ningún tornillo antes de consultar con su autoridad local, puede requerir que un electricista certificado le instale el CT.

Elegir diferentes componentes

Básicamente, hay cuatro componentes que deben elegirse o dimensionarse correctamente para su medidor de energía.

Elegir un transformador de corriente

El primero es el CT o transformador de corriente. El que se utiliza aquí es el Talema AC1030 que puede detectar una corriente nominal de 30 A y una corriente máxima de 75 A. A 220 VCA, teóricamente puede detectar hasta 16,5 kW durante cortos períodos de tiempo, pero está dimensionado para detectar continuamente 6,6 kW, lo que es adecuado para un hogar pequeño. Para calcular cuántos amperios necesita detectar el suyo, tome la potencia continua máxima que espera detectar y divídala por su voltaje (generalmente 110 V o 220 V según su país).

Dimensionamiento de la resistencia de carga

A continuación, debe dimensionar su resistencia de carga R3, esto convierte su corriente CT en una referencia de voltaje. Comience dividiendo su corriente primaria (el máximo como se usó anteriormente) por la relación de vueltas de su CT (disponible en la hoja de datos). Esto debería estar alrededor de 500-5000 a 1. Este artículo trabajó en 42A con una relación de vueltas de 1000: 1 dando una corriente secundaria de 0.042A o 42mA. Su voltaje de referencia analógica para el Arduino es de 2.5V, por lo que para determinar la resistencia usa R = V / I – R = 2.5 / 0.042 = 59.5Ω. El valor de resistencia estándar más cercano es 56 Ω, por lo que se utilizó este.

Aquí hay algunas opciones en diferentes TC y sus resistencias de carga ideales (en tamaños estándar):

  • Murata 56050C – 10A – 50: 1 – 13Ω
  • Talema AS-103 – 15A – 300: 1 – 51Ω
  • Talema AC-1020 – 20A – 1000: 1 – 130Ω
  • Alttec L01-6215 – 30A – 1000: 1 – 82Ω
  • Alttec L01-6216 – 40A – 1000: 1 – 62Ω
  • Talema ACX-1050 – 50A – 2500: 1 – 130Ω
  • Alttec L01-6218 – 60A – 1000: 1 – 43Ω
  • Talema AC-1060 – 60A – 1000: 1 – 43Ω
  • Alttec L01-6219 – 75A – 1000: 1 – 33Ω
  • Alttec L01-6221 – 150A – 1000: 1 – 18Ω
  • CTYRZCH SCT-013-000 – 100A – Resistencia de carga incorporada – 
  • TOOGOO SCT-013-000 – 100A – 

El condensador utilizado es de 10 µF, que debería ser suficiente para la mayoría de los rangos de TC para aplicaciones domésticas.

Finalmente, necesita dos resistencias divisorias para obtener el voltaje de referencia de 2.5V del Arduino. Deben tener el mismo valor, por lo que R1 = R2 y no necesitamos mucha corriente, por lo que este artículo usa dos resistencias de 100K.

firmware


#include <LiquidCrystal.h>

int currentPin = 1; // Asignar entrada CT al pin 1
doubles kilos= 0;
int PeakPower = 0;
LiquidCrystal lcd (8, 9, 4, 5, 6, 7); // Asignar pines de pantalla LCD, según los requisitos de pantalla LCD

void setup () 
{ 
  lcd.begin (16,2); // columnas filas. utilice 16,2 para una pantalla LCD de 16x2, etc.
  lcd.clear ();
  lcd.setCursor (0,0); // coloca el cursor en la columna 0, fila 0 (la primera fila)
  lcd.print ("En ejecución");
}

void loop () 
{ 
  int actual = 0;
  int maxCurrent = 0;
  int minCurrent = 1000;
  for (int i = 0; i <= 200; i ++) // Supervisa y registra la entrada de corriente durante 200 ciclos para determinar                                                                     la corriente máxima y mínima
  {
    current = analogRead (currentPin); // Lee la entrada actual y registra la corriente máxima y mínima
    if(actual> = maxCurrent)
      maxCurrent = actual;
    else if (actual <= minCurrent)
      minCurrent = actual;
  }
  if (maxCurrent <= 517)
  {
    maxCurrent = 516;
  }
  doble RMSCurrent = ((maxCurrent - 516) * 0,707) /11,8337; // Calcula la corriente RMS según el valor máximo
  int RMSPower = 220 * RMSCurrent; // Calcula la potencia RMS asumiendo un voltaje de 220 VCA, cambie a 110 VCA en consecuencia
  if (RMSPower> peakPower)
  {
    PeakPower = RMSPower;
  }
  kilos = kilos + (RMSPower * (2.05 / 60/60/1000)); // Calcule los kilovatios hora utilizados
  delay (2000);
  lcd.clear ();
  lcd.setCursor (0,0); // Muestra todos los datos actuales
  lcd.print (RMSCurrent);
  lcd.print ("A");
  lcd.setCursor (10,0);
  lcd.print (RMSPower);
  lcd.print ("W");
  lcd.setCursor (0,1);
  lcd.print (kilos);
  lcd.print ("kWh");
  lcd.setCursor (10,1);
  lcd.print (PeakPower);
  lcd.print ("W");
}

Aquí está el enlace para descargar el  código del medidor .

Debido a que su configuración, CT, resistencias y voltaje de entrada pueden ser diferentes, hay un factor de escala en el esquema que deberá cambiar antes de obtener resultados precisos; consulte la calibración a continuación. Si su pantalla LCD está conectada a los mismos pines que se utilizan aquí y su CT está conectado al mismo pin de entrada, al menos debería llenar la pantalla con algunas cifras, aunque lo más probable es que sean incorrectas y algunas pueden ser negativas.

Si no desea usar o no tiene una pantalla LCD, también puede modificar el boceto para enviarlo a la ventana serial del Arduino IDE como se muestra a continuación.

// Michael Klements
// La vida del bricolaje
// 27 de octubre de 2014

int currentPin = 1; // Asignar entrada CT al pin 1
double kilos  = 0;
int PeakPower = 0;

void setup () 
{ 
  Serial.begin (9600); // Iniciar la comunicación en serie
  Serial.println ("En ejecución");
}

void loop () 
{ 
  int actual = 0;
  int maxCurrent = 0;
  int minCurrent = 1000;
  for (int i = 0; i <= 200; i ++) // Supervisa y registra la entrada de corriente durante 200 ciclos para determinar la corriente máxima y mínima
  {
    current = analogRead (currentPin); // Lee la entrada actual y registra la corriente máxima y mínima
    if (actual> = maxCurrent)
      maxCurrent = actual;
    else if (actual <= minCurrent)
      minCurrent = actual;
  }
  if (maxCurrent <= 517)
  {
    maxCurrent = 516;
  }
  double RMSCurrent = ((maxCurrent - 516) * 0,707) /11,8337; // Calcula la corriente RMS según el valor máximo
  int RMSPower = 220 * RMSCurrent; // Calcula la potencia RMS asumiendo un voltaje de 220 VCA, cambie a 110 VCA en consecuencia
  if (RMSPower> peakPower)
  {
    PeakPower = RMSPower;
  }
  kilos = kilos + (RMSPower * (2.05 / 60/60/1000)); // Calcule los kilovatios hora utilizados
  delay (2000);
  Serial.print (RMSCurrent);
  Serial.println ("A");
  Serial.print (RMSPower);
  Serial.println ("W");
  Serial.print (kilos);
  Serial.println ("kWh");
  Serial.print (peakPower);
  Serial.println ("W");
}

Aquí está el enlace para descargar el código de salida serial del medidor .

Actualización de código

El código original del Energy Meter hacía uso de un período de tiempo fijo para calcular los kilovatios hora consumidos, esto se basaba en un tiempo de ciclo de 2050 ms y era bastante preciso.

Desde entonces, el código ha sido modificado para hacer uso de la función millis () incorporada que calcula el tiempo de ciclo exacto para cada ciclo con el fin de mejorar la precisión. Solo mejora alrededor del medio por ciento en la precisión del cálculo, pero es la mejor manera de hacerlo.

Aquí está el código mejorado:


#include <LiquidCrystal.h>

int currentPin = 1; // Asignar entrada CT al pin 1
duobles kilos = 0;
int PeakPower = 0;
unsigned long startMillis;
unsigned long endMillis;
LiquidCrystal lcd (8, 9, 4, 5, 6, 7); // Asignar pines de pantalla LCD, según los requisitos de pantalla LCD

void setup () 
{ 
  lcd.begin (16,2); // columnas filas. utilice 16,2 para una pantalla LCD de 16x2, etc.
  lcd.clear ();
  lcd.setCursor (0,0); // coloca el cursor en la columna 0, fila 0 (la primera fila)
  lcd.print ("Arduino");
  lcd.setCursor (0,1); // coloca el cursor en la columna 0, fila 1 (la segunda fila)
  lcd.print ("Medidor de energía");
  startMillis = millis ();
}

void loop () 
{ 
  int actual = 0;
  int maxCurrent = 0;
  int minCurrent = 1000;
  for (int i = 0; i <= 200; i ++) // Supervisa y registra la entrada de corriente durante 200 ciclos para determinar la corriente máxima y mínima
  {
    current = analogRead (currentPin); // Lee la entrada actual y registra la corriente máxima y mínima
    if (actual> = maxCurrent)
      maxCurrent = actual;
    else if (actual <= minCurrent)
      minCurrent = actual;
  }
  if (maxCurrent <= 517)
  {
    maxCurrent = 516;
  }
  double RMSCurrent = ((maxCurrent - 516) * 0,707) /11,8337; // Calcula la corriente RMS según el valor máximo
  int RMSPower = 220 * RMSCurrent; // Calcula la potencia RMS asumiendo un voltaje de 220 VCA, cambie a 110 VCA en consecuencia
  if (RMSPower> peakPower)
  {
    PeakPower = RMSPower;
  }
  endMillis = millis ();
  tiempo_largo_sin_firmar = endMillis - startMillis;
  kilos = kilos + ((doble) RMSPower * ((doble) tiempo / 60/60/1000000)); // Calcule los kilovatios hora utilizados
  startMillis = millis ();
  delay (2000);
  lcd.clear ();
  lcd.setCursor (0,0); // Muestra todos los datos actuales
  lcd.print (RMSCurrent);
  lcd.print ("A");
  lcd.setCursor (10,0);
  lcd.print (RMSPower);
  lcd.print ("W");
  lcd.setCursor (0,1);
  lcd.print (kilos);
  lcd.print ("kWh");
  lcd.setCursor (10,1);
  lcd.print (PeakPower);
  lcd.print ("W");
}

Aquí está el enlace para descargar el  código actualizado de  Millis Meter .

Para aquellos de ustedes que han leído que la función millis () se desborda después de aproximadamente 49 días, el código trata el rollover automáticamente haciendo uso de la variable larga sin firmar. Por ejemplo, si el desbordamiento ocurre en 10000, el milis inicial fue 9987 y el milis final fue 2031, la diferencia sería 2031-9987 = -7956 pero el valor no puede ser negativo ya que no está firmado, por lo que se convierte en -7956+ 10000 = 2044 que es la duración correcta.

Calibrar la lectura actual

Como se mencionó anteriormente, debido a que su configuración, CT, resistencias y voltaje de entrada pueden ser diferentes, hay un factor de escala en el esquema que deberá cambiar antes de obtener resultados precisos.

Para calibrar su medidor de energía, debe estar seguro de que la corriente que su medidor dice que se está extrayendo es la que espera que realmente se extraiga. Para hacer esto con precisión, necesita encontrar una carga calibrada. Estos no son fáciles de conseguir en un hogar normal, por lo que deberá encontrar algo que use una cantidad de energía establecida y constante. Usé un par de bombillas incandescentes y focos, estos vienen en una variedad de tamaños y su consumo es bastante cercano a lo que se indica en la etiqueta, es decir, una bombilla de 100W usa muy cerca de 100W de potencia real, ya que es casi completamente una carga puramente resistiva.

Enchufe una bombilla pequeña (100 W más o menos) y vea qué carga se muestra. Ahora deberá ajustar los usos del factor de escala en la línea de cálculo:

doble RMSCurrent = ((maxCurrent – 516) * 0,707) / 11,8337

En este caso fue 11.8337, puede ser mayor o menor dependiendo de su aplicación. Utilice una escala lineal para calcular esta cifra o, si no es bueno con las matemáticas, juegue con diferentes valores hasta que la carga que ha conectado se muestre en la pantalla del medidor de energía.

Una vez que haya calibrado su medidor de energía, lo reinicia y deja que haga su trabajo. A continuación se muestran dos imágenes en uso, ambas con una entrada de baja potencia y una entrada de alta potencia.

arduino-energy-meter-low-consumption

El primer número que se muestra es la corriente instantánea seguida de la potencia instantánea. En la línea inferior, los kilovatios hora utilizados desde el reinicio y luego la potencia máxima registrada desde el reinicio.

arduino-energy-meter-high-consumption

Introducción a python científico


Muchas áreas de carácter científico-técnico la adecuada elección del software y/o lenguaje de programación empleado es determinante, de cara a la potencia, versatilidad, facilidad de uso y acceso por parte de todos los usuarios en sus propios dispositivos, de manera generalizada y gratuita.

Dentro del software libre, uno de los que últimamente ha tenido una mejora sustancial, con la inclusión de potentes y versátiles nuevos módulos de cálculo simbólico (SymPy), numérico (NumPy, SciPy) y gráfico (PyPlot y Matplotlib) ha sido sin duda Python, y de ahí su vertiginosa evolución y expansión a nivel mundial, no sólo en el ámbito académico, sino también en el científico e industrial. De hecho, basta con echar un vistazo a las numerosas propuestas, tanto de comunidades de desarrolladores como de empresas privadas, surgidas a raíz de la versión de base inicial de Python, como por ejemplo IPython (interface interactivo de fácil uso, que gracias a Jupyter Notebook permite una versión HTML similar a los notebooks de Mathematica o Mapple) o Spyder (entorno integrado para cálculo científico parecido al de Matlab u Octave).

Por otro lado existen versiones completas de desarrollo, integrando Python como soporte de cálculo, pero con editores avanzados de texto para la programación y la depuración de código, ventanas de gráficos y datos, etc. La mayoría de estas plataformas integradas están disponibles para los distintos sistemas operativos Linux, MacOS X y Windows. Entre ellas cabría destacar Enthought Python Distribution (EPD), PyCharm y Anaconda CE (de Continuum Analytics).

Aunque no podamos abarcar todos los aspectos «básicos» del python científico, intentaremos en este resumen dar una idea de las principales librerías un funciones que podemos usar para NILM (Non-Intrusive Load Monitoring) sin olvidar los fundamentos de :Matplotlib y Numpy

Matplotlib: visualización con Python 

Matplotlib es una biblioteca completa para crear visualizaciones estáticas, animadas e interactivas en Python haciendo que las cosas fáciles sean fáciles y las difíciles posibles.

Nos permite crear :

  • Desarrollando gráficos de calidad de publicación con solo unas pocas líneas de código
  • Utilizando figuras interactivas que puedan hacer zoom, desplazarse, actualizar …

Personalizar

  • Tomando el control total de los estilos de línea, las propiedades de la fuente, las propiedades de los ejes …
  • Exportando e incrustando en varios formatos de archivo y entornos interactivos

Ampliar

  • Explorando la funcionalidad personalizada proporcionada por paquetes de terceros
  • Obteniendo más información sobre Matplotlib a través de los numerosos recursos de aprendizaje externos

Matplotlib es en resumen la librería de python para dibujar (equivalente al plot en matlab).

Puede encontrar mas información en el sitio oficial https://matplotlib.org/

Numpy

NumPy es una biblioteca para el lenguaje de programación Python que da soporte para crear vectores y matrices grandes multidimensionales, junto con una gran colección de funciones matemáticas de alto nivel para operar con ellas.

Numpy es pues una librería especializada para operaciones con matrices y vectores

Puede encontrar mas información en l sitio oficial https://numpy.org/

La imagen tiene un atributo ALT vacío; su nombre de archivo es image-11.png

Primeros pasos

Primero, es necesario importarlas al workspace

import numpy as np
import matplotlib.pyplot as plt

Opciones de visualizacion de matplotlib para un notebook

%matplotlib inline
plt.rcParams['figure.figsize'] = (13, 6)
plt.style.use('ggplot')

Otras importaciones:

import warnings
warnings.filterwarnings('ignore')

Crear arrays en python es muy sencillo y se puede hacer de forma nativa usando un tipo list. Sin embargo, aquí consideramos arrays del tipo numpy pues esto arrays incluyen funciones que facilitan las operaciones matemáticas y su manipulación

v=[1,2,3] # tipo list
v=np.array([1,2,3]) # array numpy
print (v)
print ("Dimensiones: " + str(v.ndim)) # numero de dimensiones
print ("Elementos: " + str(v.size)) # numero de elementos
print ("Longitud de las dimensiones: " + str(v.shape)) # longitud de cada dimensión
[1 2 3]
Dimensiones: 1
Elementos: 3
Longitud de las dimensiones: (3,)

Crear una matriz de 2 x 3:

v=np.array([[1,2,3], [4,5,6]])
print (v)
print ('Dimensiones: ' + str(v.ndim)) # numero de dimensiones
print ('Elementos: '+str(v.size)) # numero de elementos
print ('Longitud de las dimensiones: '+str(v.shape)) # longitud de cada dimensión
[[1 2 3]
 [4 5 6]]
Dimensiones: 2
Elementos: 6
Longitud de las dimensiones: (2, 3)

Crear una Matriz triple de 2 x 3 x 2 :

v=np.array([[[1,2], [3,4]],[[5,6],  [7,8]]])
print (v)
print ("Dimensiones: " + str(v.ndim)) # numero de dimensiones
print ("Elementos: "+str(v.size)) # numero de elementos
print ("Longitud de las dimensiones: "+str(v.shape) )# longitud de cada dimensión

[[[1 2], 
      [3 4]],

 [[5 6], 
      [7 8]]]
Dimensiones: 3
Elementos: 8
Longitud de las dimensiones: (2, 2, 2)

Utilizamos la función reshape para redimensionar los arrays

1 dimension

print (v.reshape(8,))
[1 2 3 4 5 6 7 8]

2 dimensiones

print (v.reshape(2,4))
[[1 2 3 4]
 [5 6 7 8]]

Matriz Identidad de 5×5

print (np.identity(5))
[[1. 0. 0. 0. 0.]
 [0. 1. 0. 0. 0.]
 [0. 0. 1. 0. 0.]
 [0. 0. 0. 1. 0.]
 [0. 0. 0. 0. 1.]]

Matriz de unos de 5×5

print ( np.ones([5,5]))
[[1. 1. 1. 1. 1.]
 [1. 1. 1. 1. 1.]
 [1. 1. 1. 1. 1.]
 [1. 1. 1. 1. 1.]
 [1. 1. 1. 1. 1.]]

Matriz de ceros de 5×5:

print (np.zeros([5,5]))
[[0. 0. 0. 0. 0.]
 [0. 0. 0. 0. 0.]
 [0. 0. 0. 0. 0.]
 [0. 0. 0. 0. 0.]
 [0. 0. 0. 0. 0.]]

Las operaciones por definición son elementwise

a=np.arange(5)
b=10*np.ones(5)
print ("vector a: "+str(a))
print ("vector b: "+str(b))
print ("suma b+a: "+str(b-a))
print ("resta b-a: "+str(b+a))
print ("producto b*a: "+str(b*a))
vector a: [0 1 2 3 4]
vector b: [10. 10. 10. 10. 10.]
suma b+a: [10.  9.  8.  7.  6.]
resta b-a: [10. 11. 12. 13. 14.]
producto b*a: [ 0. 10. 20. 30. 40.]

El producto de los vectores es:

a.dot(b)
100.0

Para las matrices tenemos que:

a=np.identity(3)
b=np.array([[1,2,3],[4,5,6],[7,8,9]])
print ("matriz a:\n"+str(a))
print ("matriz b:\n"+str(b))
print ("producto a*b:\n"+str(a.dot(b)))
print ("producto elementwise a.*b:\n"+str(a*b))
matriz a:
[[1. 0. 0.]
 [0. 1. 0.]
 [0. 0. 1.]]
matriz b:
[[1 2 3]
 [4 5 6]
 [7 8 9]]
producto a*b:
[[1. 2. 3.]
 [4. 5. 6.]
 [7. 8. 9.]]
producto elementwise a.*b:
[[1. 0. 0.]
 [0. 5. 0.]
 [0. 0. 9.]]

Vector formado por un rango de valores:

print ("De 0 a 10: " + str(np.arange(10)))
print ("De 10 a 20 de paso 0.5: "+str(np.arange(10,20,0.5)))
De 0 a 10: [0 1 2 3 4 5 6 7 8 9]
De 10 a 20 de paso 0.5: [10.  10.5 11.  11.5 12.  12.5 13.  13.5 14.  14.5 15.  15.5 16.  16.5
 17.  17.5 18.  18.5 19.  19.5]

Función linspace:

np.linspace(0,2*np.pi,10) # de 0 a 2*pi en 10 puntos equidistantes
array([0.        , 0.6981317 , 1.3962634 , 2.0943951 , 2.7925268 ,
       3.4906585 , 4.1887902 , 4.88692191, 5.58505361, 6.28318531])

función random:

np.random.rand(10)
array([0.63623588, 0.83924558, 0.35833155, 0.33835148, 0.53247758,
       0.0950348 , 0.2805706 , 0.47285484, 0.8696919 , 0.78361161])

Dibujar una función seno

t = np.arange(0.0, 2.0, 0.01)
s = np.sin(2*np.pi*t)
plt.plot(t, s)
plt.xlabel('time (s)')
plt.ylabel('voltage (mV)')
plt.title('Sinusoidal')
plt.grid(True)

Dibujar una función chirp

x=np.linspace(0,3*np.pi,500)
plt.plot(x,np.sin(x**2))
plt.title("A simple chirp")
Text(0.5, 1.0, 'A simple chirp')

Definir una función en Python

En python, las funciones pueden estar definidas en cualquier parte pero siempre antes de su llamada. En python, las anidaciones (bucles, if conditions, functions, etc.) se realizan mediante indentación, no existe el statement end. Las funciones se definen así:

def funcion_suma(x): 
    suma=0
    for i in x: 
        suma=suma+i
    return suma 
v=np.arange(10)
print (funcion_suma(v))
45

Aunque, como hemos dicho antes, numpy facilita las operaciones matemáticas y ya incluye una serie de operaciones:

print (v.sum())
print (v.cumsum())
print (v.mean())
45
[ 0  1  3  6 10 15 21 28 36 45]
4.5

Para saber más sobre numpy:

https://docs.scipy.org/doc/numpy-dev/user/quickstart.html

http://www.sam.math.ethz.ch/~raoulb/teaching/PythonTutorial/intro_numpy.html

(O simplemente «googleando»: numpy tutorial)

Pandas

Pandas (Python Data Analysis Library) es una librería de python para el análisis y manipulación de una gran cantidad de datos. También facilita el uso de «timeseries»

La llamada a la librería es:

import pandas as pd

Dado un archivo csv, la función read_csv carga los datos en un dataframe

# El parámetro parse_dates indica a pandas que al cargar este csv la primera columna [0] es de tipo datetime
df=pd.read_csv('data/events.csv',parse_dates=[0])

Las primeras N filas del dataFrame se puede visualizar de la siguiente forma

N=4
df.head(N)
timestamplabelphase
02011-10-20 12:22:01.473111A
12011-10-20 12:37:40.507111A
22011-10-20 13:23:55.390111A
32011-10-20 13:39:08.157111A

Y las N últimas columnas

df.tail(N)
timestamplabelphase
24812011-10-27 12:57:17.079111A
24822011-10-27 13:10:45.112111A
24832011-10-27 13:54:08.862111A
24842011-10-27 14:07:21.612111A

Podemos filtar por un cierto valor

df[df.phase=='B'].head()
timestamplabelphase
122011-10-20 15:45:54.590204B
132011-10-20 15:47:31.223204B
142011-10-20 16:09:00.424204B
182011-10-20 17:42:00.657155B
192011-10-20 17:42:04.407157B

Y hacer agrupaciones

df2=df.sort_values(['label','phase']).groupby(['label','phase']).count()
df2
timestamp
labelphase
101B26
102B25
103B24
108A16
111A619
1125A1
1126A1
1127A1
1200A1
1201A1

161 rows × 1 columns

Las nuevas columnas se crean fácilmente. Compatible con numpy.

df['x']=25*np.random.rand(len(df))
df['y']=100*np.sin(2*np.pi*np.linspace(0,2*np.pi,len(df)))
df['z']=df.x+df.y
df.head(5)
timestamplabelphasexyz
02011-10-20 12:22:01.473111A0.0167760.0000000.016776
12011-10-20 12:37:40.507111A5.6400571.5892417.229298
22011-10-20 13:23:55.390111A5.6322523.1780818.810333
32011-10-20 13:39:08.157111A9.2871414.76611914.053259
42011-10-20 14:25:51.473111A9.3135696.35295215.666521

Para dibujar sólo necesitamos usar la función plot

df.z.plot()
<matplotlib.axes._subplots.AxesSubplot at 0x7f6e3ba88080>

Existen ciertas funciones predefinidas que facilitan los cálculos

df.z.cumsum().plot()
<matplotlib.axes._subplots.AxesSubplot at 0x7f6e3b801e80>

Y se pueden concatenar

# Si integramos y derivamos obtenemos la misma señal
df.z.cumsum().diff().plot()
<matplotlib.axes._subplots.AxesSubplot at 0x7f6e32625320>

Unas de las herramientas más potentes de pandas es la manipulación de timeseries

df.index=df.timestamp
df.z.plot()
<matplotlib.axes._subplots.AxesSubplot at 0x7f6e325e6518>

Podemos filtar por fechas

d1='2011-10-21'
d2='2011-10-23'
df[(df.index>d1)&(df.index<d2)].z.plot()
<matplotlib.axes._subplots.AxesSubplot at 0x7f6e325d5588>

Existe una gran flexibilidad a la hora de resamplear un dataframe

# Cada día 
df.resample('1D',how='sum')
labelxyz
timestamp
2011-10-20233182044.60316112047.56801914092.171180
2011-10-21507894364.168081-4531.116500-166.948419
2011-10-221356456971.7315141811.9787858783.710300
2011-10-231021845445.726210-4029.5292421416.196968
2011-10-24462593554.5631237284.89716310839.460286
2011-10-25525523933.715959-4712.890577-779.174618
2011-10-26590843503.958137-7682.254604-4178.296467
2011-10-27181921325.9193697454.6146868780.534055
# Cada 6 horas
df.resample('6H',how='count')
timestamplabelphasexyz
timestamp
2011-10-20 12:00:00373737373737
2011-10-20 18:00:00135135135135135135
2011-10-21 00:00:00160160160160160160
2011-10-21 06:00:00676767676767
2011-10-21 12:00:00262626262626
2011-10-21 18:00:00828282828282
2011-10-22 00:00:00242424242424
2011-10-22 06:00:00107107107107107107
2011-10-22 12:00:00215215215215215215
2011-10-22 18:00:00203203203203203203
2011-10-23 00:00:00636363636363
2011-10-23 06:00:00646464646464
2011-10-23 12:00:00737373737373
2011-10-23 18:00:00237237237237237237
2011-10-24 00:00:00363636363636
2011-10-24 06:00:00393939393939
2011-10-24 12:00:00494949494949
2011-10-24 18:00:00162162162162162162
2011-10-25 00:00:00333333333333
2011-10-25 06:00:00919191919191
2011-10-25 12:00:00373737373737
2011-10-25 18:00:00152152152152152152
2011-10-26 00:00:00232323232323
2011-10-26 06:00:00616161616161
2011-10-26 12:00:00161616161616
2011-10-26 18:00:00196196196196196196
2011-10-27 00:00:00363636363636
2011-10-27 06:00:00555555555555
2011-10-27 12:00:00666666

Para aprender más sobre pandas:

http://pandas.pydata.org/pandas-docs/stable/tutorials.html

http://pandas.pydata.org/pandas-docs/stable/10min.html

Detector de eventos

Vamos a crear un detector de eventos.

Dado el consumo eléctrico de una vivienda (voltage y corriente) queremos detectar en que momento se produce una conexión de un dispositivo. Para ello, filtraremos la señal sinusoidal obteniendo el valor eficaz de la corriente cada cierto intervalo. Los cambios en el valor eficaz van a determinar las conexiones y desconexiones de los distintos dispositivos. Derivando este valor eficaz, obtenemos picos en los que existe un cambio en el valor eficaz y, por lo tanto, posibles candidatos a eventos de conexión/desconexión. Finalmente, usando un detector de picos filtraremos los eventos reales del resto.

Mediremos nuestros resultados usando métricas estándar de NILM.

Paso por paso

Importar pandas, numpy y matplotlib tal y como se ha visto anteriormente

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt

Definir una funcion llamada rms_function que devuelva un valor rms y que tenga como parámetro de entrada un vector de valores

# función rms 
def rms_function(x): 
    return np.sqrt(np.mean(np.square(x)))

Usar el siguiente path para cargar los datos en un dataframe df de pandas. Como parámetros: el índice es la columna 0 (index_col) y la fecha está en la columna 1 (parse_dates)

path='data/smart_meter_data.csv'
df= ...
path='data/smart_meter_data.csv'
df=pd.read_csv(path, parse_dates=[1],index_col=[0])

Mostrar las 5 primeras columnas del dataframe

df.head(5)
datetimeivlabelappl_namephase
02011-10-20 12:21:58.9730000.444955159.194375111RefrigeratorA
12011-10-20 12:21:58.9730830.402501160.677554111RefrigeratorA
22011-10-20 12:21:58.9731660.444955161.845163111RefrigeratorA
32011-10-20 12:21:58.9732491.102993163.107443111RefrigeratorA
42011-10-20 12:21:58.9733321.952074164.243495111RefrigeratorA

Imprimir mínimo y máximo de datetime y la diferencia de ambos

print (df.datetime.min())
print (df.datetime.max())
print (df.datetime.max()-df.datetime.min())
2011-10-20 12:21:58.973000
2011-10-20 12:23:03.713996
0 days 00:01:04.740996

Seleccionar datetime como índice del dataframe df

df.index=df.datetime

Periodo y frequencia de muestreo

# frecuencia
ts=df.datetime.diff().mean().total_seconds()
print (str(ts)+' seconds')
fs=1/ts
print ( str(fs)+' Hz')
8.3e-05 seconds
12048.192771084337 Hz

Dibujar Voltage (v) haciendo zoom en el intervalo de 100ms (6 periodos aproximadamente)

d1='2011-10-20 12:22:29.9'
d2='2011-10-20 12:22:30'
df[(df.index>d1)&(df.index<d2)].v.plot()
<matplotlib.axes._subplots.AxesSubplot at 0x7f6e2b9f86d8>
df.i.plot()
<matplotlib.axes._subplots.AxesSubplot at 0x7f6e2e274d30>

Resamplear mediante la función resample de pandas a 50ms (’50L’). La función rms_function se pasará como parámetro para obtener cada valor del resampleado. El resultado debe de guardarse en un dataframe nuevo llamado rms . Dibujar el resultado.

rms=pd.DataFrame(df.i.resample(....))
rms=pd.DataFrame(df.i.resample('50L',how=rms_function))
rms.plot()
<matplotlib.axes._subplots.AxesSubplot at 0x7f6e2f8848d0>

Hacer la derivada del dataframe rms y guardar el resultado en rms_diff.

rms_diff=rms.diff()

Dibujar el resultado (rms_diff)

rms_diff.plot()
<matplotlib.axes._subplots.AxesSubplot at 0x7f6e2f7f9fd0>

Guardar los valores de la columna «i» en una variable «y» en forma de array

y=rms_diff.i.values
Modifica los parámetros th_noise y dist de la función detect_peaks para obener los índices de los eventos y evaluar las métricas. Realizar el proceso 3 veces. ¿ Con qué valores de th_noise y dist se obtienen mejores resultados en las métricas?
th_noise=5
dist=5
from detect_peaks import detect_peaks
indexes=detect_peaks(y,mph=th_noise,mpd=dist)
dates=rms_diff.ix[indexes].index

Cuantos eventos hemos detectado

print (str(len(indexes))+' eventos detectados')
8 eventos detectados

Dibujamos los eventos y la corriente en una misma gráfica

plt.vlines(dates,-80,80)
df.i.plot()
<matplotlib.axes._subplots.AxesSubplot at 0x7f6e2f79ca58>

Métricas

from IPython.display import Image
Image(filename='metricas1.png')
Image(filename='metricas2.png')

Obtener las métricas: recall, precision y F1

FP=0.
P=9.
N=len(df)-P
TN=N-FP
P=9.
N=len(df)-P
TP=8.
FP=0.
FN=1.
TN=N-FP
recall=TP/(TP+FN)
precision=TP/(TP+FP)
F1=2*precision*recall/(precision+recall)
print (recall)
print (precision)
print (F1)
0.8888888888888888
1.0
0.9411764705882353

*Parámetros optimizados: * th_noise=0.1 y dist=5

*Con esto obtenemos : * recall=1, precision=1 y F1=1;