domingo, 15 de octubre de 2017

Convertir un modelo de Simulink a código VHDL con HDL Coder [Ejemplo 1]

Para este primer ejemplo, generarémos el código VHDL para una función booleana sencilla de la forma X = AB xor C desde un diagrama en Simulink.

Bloques requeridos
HDL Coder > Commonly Used Blocks > In1
HDL Coder > Commonly Used Blocks > Out1
HDL Coder > Logic and Bit Operations > Logical Operators

Diagrama
Una vez terminado el bloque debemos irnos a la pestaña Code > HDL Code > Generate HDL. Si no hay ningún error deberá aparecer lo siguiente en la ventana de comandos:

### Generating HDL for 'Ejemplo_HDL'.
### Starting HDL check.
### Begin VHDL Code Generation for 'Ejemplo_HDL'.
### Working on Ejemplo_HDL as hdlsrc\Ejemplo_HDL\Ejemplo_HDL.vhd.
### Creating HDL Code Generation Check Report Ejemplo_HDL_report.html
### HDL check for 'Ejemplo_HDL' complete with 0 errors, 3 warnings, and 0 messages.
### HDL code generation complete. 


Por defecto el código generado aparecerá en el directorio  hdlsrc\. El código generado para este ejemplo es el siguiente:

martes, 26 de septiembre de 2017

Jojutla

De niño crecí en un lugar que, junto con la formación académica me mis papás, moldeó en gran parte mi amor por la ciencia. Vivía en una unidad habitacional que tiene una vista muy amplia hacia el valle del suroeste de Morelos con dos volcanes en el fondo. El paso de las estaciones es muy notorio en esa zona de México gracias a los cambios en el color de la vegetación en los cerros, el ir y venir de los glaciares de los volcanes y las "lluvias" de tisne del final de la zafra azucarera. La cerros también me daban una referencia para darme cuenta que el Sol salía en un lugar diferente a lo largo del año, efecto que noté muy bien durante la formación de las mañanas en la secundaria. Las fumarolas del Popocatepetl eran muy comunes y siempre visibles por el clima. Las tormentas eléctricas también son particularmente intensas en esa región, algo que sólo noté hasta que me mudé a otros estados. Hasta la fecha no he estado en un lugar de México en dónde haya tantas especies de plantas y animales en un mismo patio como en ese lugar. Jojutla siempre me recordó que la Tierra es un lugar vivo y de constante cambio. Este mes lo ha hecho más que nunca y nos ha dado a muchos un recordatorio de lo frágiles que somos los humanos ante el poder de la naturaleza

* * *

El 26 de septiembre de 2007 publiqué la primera estrada en este, mi primer y único blog continuamente activo. Me hubiera gustado que mi entrada de X aniversario tuviera un tema menos sombrío. El de 19 de septiembre de este mes ocurrió uno de los sismos más devastadores en la historia de mi país. Mis primeras entradas hablaban de la región, ahora en buena parte en ruinas, en la que pasé toda mi infancia. Nunca nos pasó por la mente el vivir en Jojutla un desastre tan grande. Los sismos siempre fueron frecuentes pero nunca pasaron un breve susto. De niños crecimos haciendo colectas de vivieres hacia lugares que nos parecían muy lejanos. Este fue un tema recurrente entre mis amigos de la infancia con quienes estuve apoyado como brigadista de remoción de escombros desde el martes pasado hasta el domingo. Caminar por aquella ciudad ciudad tan querida reducida a una zona de guerra me provocó una sensación extraña de irrealidad. Al mismo tiempo el ver el compromiso de tanta gente me animó de forma inesperada. Por ahora las emergencias ya están atendidas y solo queda trabajar en las demoliciones y reconstrucción. Aún hay un largo camino por trabajar. Me hubiera quedado más tiempo pero desafortunadamente tenía acceso a internet sólo por breves momentos y no me dí cuenta que la BUAP reanudaba labores hasta el 2 de octubre. Sólo tenía que regresar por mis tramites de titulación. Posiblemente regrese a apoyar este jueves.  

jueves, 7 de septiembre de 2017

Clasificador k-NN supervisado con Scikit-Learn

Scikit-Learn es un módulo para Python que incluye varias rutinas de clasificación, regresión y clustering entre otras herramientas matemáticas utilizadas en machine learning y minería de datos. En esta entrada únicamente trataré un ejemplo sencillo que permitirá introducir a la idea central detrás del clustering o análisis de grupos.

Imaginemos que existe una especie de conejos en la que las hembras tienden a ser grandes y tener pelaje de color gris claro mientras que los machos tienden a ser más pequeños y tener un pelaje más oscuro. Crearemos un programa que implemente un clasificador k-NN (k-nearest neighbors) utilizando el modulo scikit-learn. Este algoritmo de clasificación es muy sencillo. Su trabajo es dividir un espacio de n dimensiones (dónde cada una de ellas en la práctica representa un atributo o característica de un objeto) en N regiones (dónde cada una representa una clase de objetos). El modelo de división del espacio se basa en la razón de clases de los 'k' datos vecinos más cercanos a un punto de dicho espacio. Por ejemplo, para el caso en dónde solo existen n = 2 características y N = 2 clases (triángulos y cuadrados):
Para k = 3 (circulo solido), el elemento desconocido en verde será clasificado dentro de la categoría de los triángulos (2 triángulos vs 1 cuadrado). Para un k = 5 (circulo punteado), el elemento sera clasificado dentro de la categoría de los cuadrados (3 cuadrados vs 2 triángulos). El modelo de clasificación debe poder hacer predicciones para todos los posibles datos que caigan en cualquier punto del espacio de características. El proceso de generación de este modelo se conoce como entrenamiento. Existen dos caminos: supervisado, cuando el conjunto de datos de entrenamiento esta clasificado desde un principio para el algoritmo y no supervisado, cuando no se le dice algoritmo a que clase pertenece cada elemento del conjunto de entrenamiento. En nuestro ejemplo de los conejos utilizaremos un entrenamiento supervisado a partir de un conjunto de datos de entrenamiento que generé con la siguiente distribución (conejos.csv):
El color de los conejos está especificado por valores en escala de grises de 8-bits (0-255). El código se puede dividir en 4 secciones: carga y acondicionado de los datos de entrenmiento, instanciación del clasificador desde el módulo, entrenamiento del clasificador y finalmente la predicción de la clase de un elemento arbitrario. Una vez entrenado el clasificador, debe poder retornar una predicción de la clase a la cual puede pertenecer el vector [color,tamaño]. El código es el siguiente:

Para más detalles recomiendo revisar la documentación de sklearn.neighbors.KNeighborsClassifier.

jueves, 27 de julio de 2017

De finales y comienzos

Nunca pensé que terminar un tratamiento psicológico se sintiera como terminar un libro, juego o una serie que te gustó mucho. Te libraste del cruel villano después tantas batallas, finalmente, está hecho, pero ahora sabes que no volverás a ver a quienes te acompañaron en esta larga campaña. Quizá haya un nuevo enemigo y los antiguos combatientes deban ser reunidos de nuevo en el futuro. Pero por ahora camino sobre las ruinas del régimen vencido con una mezcla de libertad y soledad en el pecho... Pero me si me detengo a pensarlo, ¿no es este el inicio de una historia completamente nueva?

lunes, 17 de julio de 2017

He Has Left Us Alone But Shafts of Light Sometimes Grace the Corners of Our Rooms


Tenía bastante que no disfrutaba tanto un disco de Post-Rock. El disco debut de Thee Silver Mt. Zion Memorial Orchestra. Poniéndome un poco en contexto, acabo de emborracharme solo con una Calavera Witbier mientras lo escuchaba. La lluvia de verano en Puebla aún suena allá afuera. Puedo ahorrarme una larga reseña de por qué este disco en una obra de arte si únicamente les recomiendo salir a buscar unas buenas cervezas para acompañarlo. Acompañarlo en soledad y con el mejor equipo de sonido que tengan. Descargar

jueves, 13 de julio de 2017

Ejemplo de programación multithread en Python

Como anexo a la entrada anterior, pongo un ejemplo del mismo programa pero en Python utilizando el módulo threadings. En este programa hago un cambio importante, estoy utilizando un producto punto del módulo numpy, una operación que está optimizada y escrita en C lo que me da tiempos de ejecución menores que los del programa anterior [Ref]:

Ejecutamos y medimos el tiempo de ejecución con:

time (python MatMulThreads.py)

Estos fueron los tiempos de ejecución en una Raspberry Pi 3 para matrices de 2000x2000:

Debe tomarse en cuenta que a pesar de que ambos programas son semanticamente iguales, su implementación no lo es. Es de notarse que en este caso no hay tanta mejoría entre 2 y 4 threads.

miércoles, 12 de julio de 2017

Programación multithread en C (pthread.h)

Un thread o hilo de ejecución es el subconjunto mínimo de instrucciones dentro de un programa que pueden ser ejecutadas secuencialmente [existe el caso particular en el que un solo hilo de ejecución representa al programa entero]. Este subconjunto de instrucciones viene asociado con los registros y locaciones de memoria que son usados durante su ejecución. Si un programa está divido en varios hilos de ejecución, cada uno de ellos puede ejecutarse de forma concurrente (o de forma paralela si se cuenta con múltiples procesadores). Es importante tener claros los conceptos de proceso e hilo de ejecución cuando se va hacer programación paralela. Recomiendo ver los videos sobre hilos de ejecución y gestión de procesos de la UCAM dónde resumen muy bien los temas en menos de 10 minutos (español). De manera conceptual, podemos representar la relación entre los hilos de ejecución (threads) y un proceso de la siguiente manera:
Se conoce como Pthreads o POSIX threads a un modelo de ejecución estandar (independiente del lenguaje de programación) para sistemas Unix. Esta definido como un conjunto de rutinas y tipos en C implementados en el header: pthread.h.

Si un programa en C/C++ no crea nuevos hilos explícitamente, será ejecutado como un proceso de hilo único. Para la crear un nuevo thread con la librería pthreads utilizamos la función pthread_create:

int pthread_create(thread_ID,att,rutina_concurrente,argumento) 

thread_ID: apuntador al identificador de hilo (tipo pthread_t)
att: apuntador a la estructura que contiene los atributos del hilo. NULL para atributos por defecto.
rutina_concurrente: puntero a la función que contiene la rutina que va a ejecutar el hilo.
argumento: argumento opcional a pasar hacia rutina_concurrente (tipo void*). NULL si ninguno es requerido.

Pueden encontrar una descripción más detalla del resto de las funciones de la librería aquí.

Ejemplo: Multiplicación matricial 

En esta entrada quiero hacer énfasis en la creación de hilos y ejecución. Como ejemplo de esto, probaremos un programa que realice una multiplicación de matrices cuadradas de forma paralela. Consideremos el algoritmo general para la multiplicación de matrices:
 Una forma de paralelizar el producto matricial es dividir las operaciones entre filas y columnas entre los hilos (no es la forma más eficiente pero será muy útil para aprender a dividir un programa en hilos independientes). El número total de productos punto entre n filas y  m columnas es igual a n*m. Si usamos matrices de 16x16, el total de operaciones entre filas y columnas sería de 256. Si dividimos estas operaciones en partes iguales entre 4 hilos de ejecución, cada hilo llevará a cabo 64 operaciones. Podemos hacer esto dividiendo las filas de la matriz A a lo largo del índice 'i'. Un ejemplo para matrices 8x8 repartido en dos hilos:
Cada hilo en este caso realiza 8 operaciones de dot(fila,columna) del total de 16. Con 4 hilos, cada uno accedería a una única fila de la matriz A realizado 4 operaciones. La rutina que pasaremos como argumento a la función pthread_create() es la siguiente:

void *matmul(void* id_arg){
  int i,j,k;
  long  id = (long) id_arg;
  int rows_per_thr = col/num_of_threads;
  int start_index = id*rows_per_thr;
  int final_index = (id+1)*rows_per_thr;

  for(i=start_index;i
< final_index;i++){
   for(j=0;j
< col;j++){
    for(k=0;k
< row;k++){
      C[i][j] += A[i][k]*B[k][j]; 

    }
   }
  }
}


Todos los hilos ejecutarán exactamente la misma rutina pero los datos a los que accederán (elementos de 'A','B' y 'C') dependieran del identificador de cada hilo. Ya que los hilos tienen acceso a la memoria compartida del proceso principal podrán leer directamente las locaciones de memoria de las variables globales del programa. Sólo un identificador de hilo será pasado como argumento para la rutina paralela. La creación de threads se hace dentro de main() de la siguiente manera:

  pthread_t tid[num_of_threads];
  long rank;

  //Creación de threads
  for (rank = 0; rank < num_of_threads; rank++)
     pthread_create(&tid[rank], NULL,matmul , (void*) rank);


En este punto del código los hilos comienzan a ejecutar independiente la rutina que se les ha asignado. Sin embargo aquí hay un problema: el programa principal continuara con su ejecución sin esperar a que los hilos terminen sus tareas. Para resolver esto es necesario indicarle explícitamente al programa principal esperar a cada uno de los hilos con la función pthread_join():

  //Unión de threads
  for (rank = 0; rank < num_of_threads; rank++)
      pthread_join(tid[rank], NULL);


Nota: en este ejemplo se considera que la mayoría de las implementaciones de Pthreads los hilos tienen el atributo de joinable por defecto. Para mayor portabilidad se recomienda activar este atributo explícitamente. Pueden encontrar un ejemplo de esto aquí.

Para simplificar el código, las matrices serán leídas en nuestro programa desde un archivo de texto. Creamos este archivo con dos matrices de enteros aleatorios (y una matriz de resultado para verificar la salida del programa en C) con el siguiente código en Python:

import numpy as np
n= 16
A = np.random.randint(0,9,size = (n,n))
B = np.random.randint(0,9,size = (n,n))
C = np.matmul(A,B)
namein = "matext"+str(n)+"x"+str(n)+".txt"
nameres = "resultado"+str(n)+"x"+str(n)+".txt"
np.savetxt(namein,np.concatenate((A,B),0),fmt = "%d",delimiter = " ")
np.savetxt(nameres,C,fmt="%d",delimiter= " ")

El código completo en C es el siguiente:

Compilamos y ejecutamos desde terminal con:

 gcc MatMulThreads.c -lpthread ; ./a.out | less 

Para verificar la reducción en el tiempo de ejecución, corrí el programa en una Raspbery Pi 3 que tiene un quad-core ARM Cortex-A53. Generé dos matrices de 1000x1000 y estos fueron los resultados variando el número de threads:
Este programa tiene varias limitaciones, una de ellas es que sólo admite números pares de hilos. Se deja como ejercicio modificarlo para cualquier número de hilos.  

Bien, hasta aquí se han cubierto los aspectos de creación y unión de threads. Existen aún varios temas importantes más como la sincronización de hilos, especialmente en los accesos de memoria (en el ejemplo de la multiplicaciones de matrices no existen conflictos de memoria ya que cada hilo trabaja con una parte diferente de los arreglos). Estos temas son tratados de forma consista en este artículo y en estas diapositivas de la Universidad de Buenos Aires. 

[Escribí el mismo ejemplo de esta entrada en Python utilizando el módulo threading. Pueden revisarlo aquí.]

Referencias
POSIX Threads Programming, Blaise Barney, Lawrence Livermore National Laboratory

martes, 4 de julio de 2017

Instalación de GNU Scientific Library (Debian)

GNU Scientific Library (GSL) es una librería numérica para C/C++, al parecer, poco conocida fuera del ámbito académico especializado. Tiene un repertorio de rutinas comparable con el popular módulo SciPy para Python. Encuentro muy útil GSL para sistemas con Linux embebido en robótica o IoT, especialmente en aplicaciones que requieran la operación de estos sistemas de tiempo real. Recomendaría echarle un vistazo si trabajan en proyectos serios con Raspberry Pi's o UDO's que requieran cálculos numéricos avanzados. Su documentación está bien estructurada y es muy amigable. El único inconveniente es que no hay muchos ejemplos y ayuda ahí afuera (ni siquiera en Stack Overflow) así que habrá que leer cuidadosamente el manual de referencia.

Para instalar GSL en Debian y derivados:

sudo apt-get install libgsl0ldbl libgsl0-dev 

Probaremos la instalación con el siguiente ejemplo con operaciones de números complejos (revisar documentación):

#include <stdio.h>
#include 
<gsl/gsl_complex.h>
#include 
<gsl/gsl_complex_math.h>

int main (void){
 
  //Declaración de variables complejas
  gsl_complex z,w,u;
 
  //Asignación de valores
  GSL_SET_COMPLEX(&z,3,4);  //z = 3+4i
  GSL_SET_COMPLEX(&w,4,5);  //w = 4+5i

  //ejemplos de operaciones
  u = gsl_complex_add(z,w);
  printf("z + w = %2.f + i%.2f\n",u.dat[0],u.dat[1]);
  u = gsl_complex_mul(z,w);
  printf("z*w = %.2f + i%.2f\n",u.dat[0],u.dat[1]);
  u = gsl_complex_log(z);
  printf("log(z) = %.2f + i%.2f\n",u.dat[0],u.dat[1]);
  u = gsl_complex_exp(gsl_complex_add(z,w));
  printf("exp(z+w) = %.2f + i%.2f\n",u.dat[0],u.dat[1]);
  return 0;
}


Compilamos y ejecutamos desde terminal con:

gcc GSL_Complex_demo.c -lgsl -lgslcblas ; ./a.out

martes, 27 de junio de 2017

Segmentación de imágenes por umbral con scikit-image

Segmentar una imagen a través de umbrales es probablemente la forma más sencilla de identificar regiones de interés en una imagen. Esto hace que este método sea vulnerable a ruido, sombras y otros factores más a la hora de utilizarlo en la práctica. Sin embargo es una forma sencilla de introducir a las funciones de segmentación del módulo scikit-image. Primero, para este ejemplo tomaremos la siguiente pintura de Kazimir Malevich:
 Separaremos la imagen en secciones para las cuales obtendremos un array con los índices de cada pixel para cada una de ellas. Para ello cargamos primero la imagen:

#Cargar imagen
im = skimage.io.imread("malevich-map-style.png")
imGS = skimage.color.rgb2gray(im) #Escala de grises


 y creamos una imagen binarizada a partir de un umbral en la imagen en escala de grises:

#Umbral
thresh = 0.7
bw = binary_erosion(closing(imGS < thresh, np.ones((3,3))))


La variable bw es la imagen binarizada corregida y es obtenida aplicando la función closing seguida de binary_erosion (ambas contenidas en skimage.morphology). Closing() realiza una operación de dilatación seguida de una erosión, lo que permite cerrar "agujeros" oscuros en una imagen. El primer argumento que toma en este ejemplo es la imagen binarizada resultante de aplicar el umbral imGS < 0.7 y el segundo es una matriz que determina el tamaño de la dilatación/erosión. Binary_erosion es solo una operación de erosión que elimina las falsas regiones debidas a las zonas brillantes del fondo. Teniendo la imagen binarizada limpia, obtenermos los segmentos:

#Segmentación
labeled_im =
skimage.measure.label(bw)
lab_img = label2rgb(labeled_im, image=imGS,bg_label=0)


dónde labeled_im es una nueva imagen obtenida con la función label()en dónde cada pixel de la imagen en escala de grises original esta clasificado con su número de región . Si esto no es muy claro, pueden visualizar este arreglo de datos con el explorador de variables de Spyder dando doble click en el nombre de la variable:
La imagen lab_img es una visualización gráfica de las etiquetas de región sobre la imagen original que al mostrarla con Matplotlib y utilizando la función plt.text() junto con los datos contenidos en labeled_im, el resultado final es el siguiente:
 Código completo:

# -*- coding: utf-8 -*-
"""
Created on Tue Jun 27 13:27:53 2017

@author: rodolfo
"""

import numpy as np
from matplotlib import pyplot as plt
import skimage
from skimage.morphology import closing, binary_erosion
from skimage.color import label2rgb

#Cargar imagen
im =
skimage.io.imread("malevich-map-style.png")
imR = im[:,:,0]
imG = im[:,:,1]
imB = im[:,:,2]
imGS = skimage.color.rgb2gray(im)

#Umbral
thresh = 0.7
bw = binary_erosion(closing(imGS < thresh, np.ones((3,3))))

#Segmentación
labeled_im =
skimage.measure.label(bw)
lab_img = label2rgb(labeled_im, image=imGS,bg_label=0)

#Array de arrays de indices por segmento
reg_index = []
for i in range(0, np.max(labeled_im)):
    u = np.where(labeled_im == i)
    reg_index.append(zip(u[0],u[1]))


#Plots
R = []
cr = []
for i in range(0,np.max(labeled_im)):
    R.append(np.where(labeled_im==i+1))
    cr.append((np.mean(R[i][0]),np.mean(R[i][1])))

plt.figure(1)
plt.imshow(im)

plt.figure(2)
plt.imshow(lab_img,cmap="hot")
for i in range(0,np.max(labeled_im)):
    plt.text(cr[i][1],cr[i][0],str(i+1),color='w')


sábado, 20 de mayo de 2017

Experimentando con autómatas celulares

Encontré por accidente el libro The New Kind of Science de Stephen Wolfram mientras buscaba libros de estadística en la biblioteca de la FCFM de la BUAP. Después de leer los primeros capítulos no pude evitar probar en Python algunos autómatas celulares. Es increíble el como con reglas tan simples se puede generar tanta complejidad. Sólo observen, tomemos una una linea de n casillas de largo. Cada casilla puede contener un '1' o un '0'. Supongamos que inicialmente la linea contiene solo '0' exceptuando la casilla n/2 que inicializamos con '1'. Nuestro pequeño autómata se encargará de generar una nueva linea de n casillas debajo de la inicial rellenando las casillas de la nueva linea de acuerdo a una función booleana cuyas entradas son casillas de la linea anterior. El autómata tiene un tamaño de 3 casillas y va rellenado las casillas de la nueva linea moviéndose una casilla a la vez hasta terminar todo el array. Un ejemplo para un autómata sería:

 dónde las casillas blancas representan '0' y las negras '1'. Consideremos un autómata un poco más complicado definido por el siguiente conjunto de reglas:

 Podemos expresar la regla anterior como una función booleana de la siguiente forma:
la cual podemos reducir mediante factorización a la forma:
El código en Python para este autómata celular queda:

import numpy as np
import matplotlib.pyplot as plt

n = 500
iline = np.zeros(n)
iline[n/2] = 1;
lines= [iline]

def rule(A,B,C):
    X = A and not(C) or not(A) and C
    return X

for i in range(0,n):
    line = lines[i]
    nline = np.zeros(n)
    for j in range(0,n-2):
        if rule(line[j],line[j+1],line[j+2]) :
            nline[j+1] = 1
    lines.append(nline)
   
lines = np.array(lines)

plt.imshow(lines, cmap = 'gray')


Y el espectacular resultado graficado es:
La figura resultante es un fractal conocido como Triángulo de Sierpinski y el autómata que lo genera se conoce como Regla 90.

viernes, 14 de abril de 2017

Listas de interes

Hoy en el foro de soporte de Facebook confirmé mis sospechas: eliminaron las listas de interés desde noviembre del año pasado. Por suerte había guardado los enlaces a mis listas por lo que aún puedo verlas, pero ya no se tiene la opción de crear más. Ok, sé que no era de las opciones más populares y de no ser por las protestas en el foro seguiría creyendo que era el único que las utilizaba, pero eran una herramienta muy útil para organizar los likes y también para evitar cerrarse en nuestras burbujas de gustos (gracias al machine learning que le da prioridad de aparición a lo que más se ajusta a nuestros gustos). Una de mis listas más útiles recopila más de 80 fanpages de centros culturales en la ciudad de Puebla (si no la borran en el futuro, pueden verla aquí). Mi lista favorita, y arma secreta para cazar convocatorias de escuelas de verano, estancias de investigación, etc, tiene más de 100 fanpages de centros de investigación, universidades y programas de vinculación en México y unos cuantos en el extranjero. Si no la eliminan, lo que espero jamás pase, pueden seguirla aquí. Es una lastima que hayan quitado una de la herramientas que más me han servido para vivir el arte y la ciencia en el mundo real.

viernes, 24 de marzo de 2017

Spoon - Kill The Moonlight (2002)


Este es sin dudas un disco que podría recomendarle a cualquier persona con un sello de garantia y, considerando mis gustos, esto no es algo que podría decir de cualquier otro. Llevo ya algunas semanas escuchándolo. Me parece que tiene un sonido muy accesible y al mismo tiempo una calidad que bien podría ponerse a un lado de un Yankee Hotel Foxtrot o un The Moon & Antartica. Me costó bastante encontrarlo pero la ardua búsqueda entre enlaces muertos valió mucho la pena. Era justo ponerlo al alcance de todos por aquí.

sábado, 18 de marzo de 2017

Monterrey Vol. 2

Monterrey me pareció mucho más bonito en esta segunda visita. Sospecho que en aquella primera ocasión había un bias emocional que me hacía ver todo triste y polvoriento. Esta vez también fue más cómodo al tener apoyo de hospedaje en el hotel que está a dos estaciones de metro de la UANL. Vine, como en la primera vez, a presentar un trabajo al TASP. Mi poster de esto año trató sobre los resultados de las pruebas de la montura para mi radiotelescopio. Aunque casi no estuve en mi poster a la hora de las presentaciones, discutí mi trabajo con bastantes personas en los pasillos y en las comidas. Tuve muy buenos comentarios y salieron varias ofertas de trabajo. Esto fue algo muy reconfortante después de las interminables asoleadas y vueltas al taller para el maquinado de las piezas. Aprendí bastante sobre el trabajo de investigación que se realiza en México en materia de clima espacial. Me enteré de proyectos como TLALOCNet y RICE. Me doy cuenta que hace falta mucha más vinculación de ingenieros en el área de instrumentación geofísica pues parece estar creciendo rápido en el país. Pienso que sería buena idea ir comenzando a hacer divulgación aquí en la facultad de electrónica de la BUAP (antes de que me vaya). Como siempre, conocí a gente muy interesante, geofísicos y biólogos muy apasionados con su trabajo. Además de distraerme un rato (una semana entera), fue un viaje bien aprovechado.

domingo, 12 de marzo de 2017

Convertir coordenadas ecuatoriales a altazimutales (horizontales) con astropy

Algunos meses después de entender la geometría y hacer un programa para convertir entre estas coordenadas para mi tesis de licenciatura, descubrí que las conversiones podían hacerse con el módulo de astropy. No fue tiempo perdido pues entendí bastante bien el proceso. Para una conversión de muy buena precisión recomendaría usar astropy aunque sus funciones tienen el inconveniente de ser relativamente lentas (promedian una ejecución de 3 ms).

# -*- coding: utf-8 -*-
"""
Created on Sun Mar 12 20:28:24 2017

@author: rodolfo
"""
from datetime import datetime
import astropy.units as u
from astropy.time import Time
from astropy.coordinates import SkyCoord, EarthLocation, AltAz, Angle

#Locación
Puebla = EarthLocation(lat=Angle('19d02m36s'), 
                       lon=Angle('-98d21m6.94s'), 
                       height=2152*u.m)
#Fecha y hora
utc_time = Time(datetimeutcnow(), scale='utc')
time = utc_time 

#Ecuatoriales
Sirio = SkyCoord(Angle('101d17m18.7469s'), 
                 Angle('-16d42m47.3141s'), 
                 frame='icrs')
#Horizontales
Sirio = Sirio.transform_to(AltAz(obstime=time,location=Puebla))

print "Altura: ", Sirio.alt.degree
print "Azimut; ", Sirio.az.degree

jueves, 26 de enero de 2017

Días comunes

A principios de diciembre compré en la librería universitaria, en su semana remate, un libro de la editorial Tusquets (la de cuadritos) que tomé al azar. "Lo que queda por vivir" de John Updike (en $30, por cierto). Me está gustando bastante. Es un libro de cuentos con un estilo que me recuerda a Chejov, una narrativa donde no ocurre nada extraordinario y sin embargo es difícil detener su lectura. Estas últimas noches he leído a ambos escritores en paralelo y me ha parecido curioso que de todos estos relatos de lo "ordinario" me haya llegado una revelación tan obvia y al mismo tiempo tan extraordinaria: cada uno de nosotros está viviendo su propia historia, y es un historia sin pausas. Creemos saberlo, creemos que es obvio pero si lo piensan con cuidado, realmente no lo habíamos entendido. Tenemos la idea de que lo que llamamos "nuestra historia" es únicamente un conjunto de "milestones" (logros académicos, un nacimiento, un accidente...) y un conjunto de "snapshots" (viajes memorables, borracheras con buenos amigos, un noviazgo que nos marcó, una despedida en una estación/aeropuerto...). Estos elementos son verdaderamente claves en nuestra vida, pero pienso que el hecho de reducir el continuo de nuestra existencia en ellos nos da la falsa sensación de que todo lo ocurrido entre esas memorias es mero espacio vacío. No solemos darnos cuenta de que el más ordinario de nuestros días es también digno de un relato. Deberíamos empezar a ver esa normalidad también como un género literario de la vida, porque así, con todo y su brevedad, nos daríamos cuenta de que tiene más páginas de las que pensábamos.

viernes, 6 de enero de 2017

Muestras I/Q: Un ejemplo con GNU Radio

En la entrada anterior resumía todo lo relevante sobre las señales analíticas (I/Q). Esta vez haremos un ejemplo que nos ayudará a visualizarlas mejor y hacernos de una intuición sobre su naturaleza. El ejemplo consiste en obtener las componentes I/Q de la siguiente señal:
Al ser una señal real, esperamos que tenga un espectro simétrico tal como podemos observar en el sink de frecuencia:
Para obtener las componentes I/Q de esta señal requerimos implementar un modulador en cuadratura con un filtro pasabajas para cada componente y posteriormente usar un bloque de conversión Float to Complex. No hará falta agregar simular un ADC ya que la señal generada ya es discreta. Si aún no están familiarizados con GNU Radio tengo un entrada introductoria [Asumo también cierto conocimiento en filtros digitales, aunque los bloques de estos son fáciles de usar y pueden consultar su documentación aquí] y tampoco estará de más tener muy claro el concepto de frecuencia en tiempo discreto. El flowgraph nos queda de la siguiente manera:
Click para agrandar
He agregado un slider gráfico que nos permite cambiar de forma interactiva la frecuencia de los generadores de señal que representan al oscilador local en cuadratura. Esto me parece un recurso didáctico muy bueno porque ayuda visualizar la relación que existe entre la frecuencia del oscilador local y la frecuencia que vemos en el centro (0 Hz) del espectro de la señal analítica. Basta asignar una frecuencia de 14.5 KHz para sintonizar a la componente de esa misma frecuencia que existe en nuestra señal de prueba:

jueves, 5 de enero de 2017

Entendiendo a las señales analíticas (muestras I/Q)

Las muestras I/Q son una representación matemática discreta de una señal de radiofrecuencia cuya característica más útil es tener un espectro de frecuencia negativa nulo (un espectro compacto y asimétrico). Dado que el procedimiento de obtención de estas muestras involucra una reducción de frecuencia, su uso esta muy extendido en en las aplicaciones de radio definido por software y comunicaciones digitales en general. Las dos características anteriores nos dan dos ventajas: la primera elimina la redundancia en el espectro y la segunda permite reducir en el ancho de banda.

¿Cómo se obtiene una representación compleja de una señal real? Para el caso continuo el proceso es el siguiente:
 Simplemente se multiplica a la señal de radiofrecuencia por un sinusoide complejo y se aplica un filtro pasa-bajas para eliminar la redundancia en el espectro. Esto es lo que se conoce como señal analítica. Este procedimiento se vuelve más claro en el dominio de la frecuencia:
Recordemos que una multiplicación en el dominio del tiempo se manifiesta como una convolución en el dominio de la frecuencia. La convolución entre el espectro asimétrico de un sinuosoide complejo produce un corrimiento en frecuencia, el cual está dado por la frecuencia del sinusoide complejo (que lleva la etiqueta de LO por Local Oscillator) .

Ahora bien, el proceso anterior es puramente matemático. ¿Cómo se hace en el mundo real? El problema se resuelve implementando un modulador en cuadratura como se muestra en el siguiente diagrama:
 Tomando en cuenta a la identidad de Euler, una multiplicación por un sinusoide complejo se emula como la multiplicación y filtrado separado de dos sinusoides reales en cuadratura (heterodino en cuadratura). Posteriormente, al digitalizar independientemente a cada una de las componentes, se consideran a las muestras del ADC en cuadratura como la parte imaginaria multiplica por -1. La I viene de In-phase y la Q de Quadrature. Prácticamente todos los dongles SDR que pueden conseguirse fácilmente en el mercado como los RTL's y el FunCube funcionan de esta manera. En la siguiente entrada haremos un ejemplo con GNU Radio para clarificar aún más las cosas.

Referencias:
"RTL-SDR: Inexpensive Software Defined Radio", EE123 University of California Berkeley
"I/Q Data for Dummies", Mikael Q Kuisma
"Complex Sinusoids", PDS Related