Un experimento aleatorio es un experimento en que el resultado varia de forma no predecible cuando se repite bajo las mismas condiciones.
Un experimento aleatorio se especifica dado el procedimiento experimental y el grupo de uno o mas medidas y observaciones.
Un solo procedimiento, pero con diferentes en las observaciones
puede tener mas de una medida u observación
las medidas pueden ser valores contínuos
pueden ser secuenciales, o una secuencia de experimentos más simples o subexperimentos
Espacio Muestral
Se define como resultado o muestra ζ (zeta griega) de un experimento aleatorio como el resultado que no puede ser descompuesto a partir de otros resultados, o una combinación de ellos.
El espacio muestral S de un experimento aleatorio se define como el grupo de todos los resultados posibles, que se especifican como:
Lista de todos los resultados a={0,1,2,3}
Mediante una propiedad de los resultados A={x: 0≤x≤3, x es entero}
Un espacio muestral discreto es contable, ejemplo las caras de un dado.
un espacio muestral contínuo es contable infinito o no contable, ejemplo los valores de voltaje de descarga de un capacitor, o el tiempo hasta la primera falla de un equipo.
Algunas veces es conveniente incluir en el espacio muestral el resultado imposible. Ejemplo, incluir que el espacio muestral es un número real positivo para el tiempo antes de la falla de un dispositivo, aunque conocemos que no tiene tiempo de operación infinito.
Eventos
Normalmente no es de interés un resultado, muestra u observación específica, el interés se enfoca más en que ocurra un evento, lo que requiere usar subgrupos de S. Por ejemplo, que el voltaje de una señal sea negativo {señal: -∞ < señal < 0}
Dos eventos de interés especial son: evento de certeza S, y el evento imposible o nulo ∅, que no resultados o nunca ocurre.
Eventos Clase
La teoría de probabilidades usa también lo que se conoce como una clase F con los eventos de interés, que son aquellos a los que se les asigna valores de probabilidad.
La teoría de probabilidad inicia con un grupo de axiomas que indican las propiedades que los valores de probabilidad deben satisfacer.
Supone que se definieron:
el experimento aleatorio con todos los posibles valores en el espacio muestral S
las clases y subgrupos llamados eventos
cada evento A tiene asignado un valor P[A]
con lo que se satisfacen los AXIOMAS:
Axioma I : 0 ≤ P[A] ≤ 1
Axioma II : P[S] = 1
Axioma III : Si A∩B=∅, P[A∪B] = P[A] + P[B]
Lo que se expresa como, que la probabilidad siempre es no negativa, que la probabilidad del espacio muestral es 1, que incluye la probabilidad P[∅] que a veces se le denomina «imposible» e igual a 0.
Para que un modelo sea útil, se debe poder realizar predicciones del comportamiento del sistema en el futuro, es decir debe tener un comportamiento regular.
Muchos modelos de probabilidad en ingeniería se fundamentan en los promedios obtenidos en las secuencias de resultados, así las muchas repeticiones o intentos del experimento aleatorio son consistentes y tienden aproximadamente al mismo valor.
Esta es la propiedad de regularidad estadítica.
En el experimento de la urna y pelotas con S={0, 1, 2}, se pueden procesar las muestras del archivo: aleatoriosurna100.txt
Para seleccionar los valores únicos de las muestras y contar el número de veces de cada uno se puede usar la instrucción de numpy np.unique
los valores unicos: [ 0. 1. 2.]
la cuenta de unicos: [33 31 36]
Nk:
[[ 0. 33.]
[ 1. 31.]
[ 2. 36.]]
Instrucciones en Python
# Muestras de una urna con pelotas S=[0,1,2]# propuesta: edelros@espol.edu.ecimport numpy as np
# INGRESO# archivotxt=input('nombre del archivo.txt: ')
archivotxt = 'aleatoriosurna100.txt'# PROCEDIMIENTO
muestras = np.loadtxt(archivotxt)
n = len(muestras)
# cuenta los únicos
[unicos,cuenta] = np.unique(muestras,return_counts=True)
# une en tabla N, y aplica T transpuesta
N = np.asarray((unicos, cuenta)).T
# SALIDA print('los eventos unicos: ', unicos)
print('la cuenta de unicos: ', cuenta)
print('Nk:')
print(N)
El número de veces que aparece la pelota 0 en el experimento es N0(n)=33,la pelota N1(n)=31 y N2(n)=36.
Frecuencia relativa
La frecuencia relativa de un resultado k, se define como:
fk(n)=nNk(n)
y lo que la regularidad estadística quiere decir es que fk(n) varía siempre acercandose a una constante pk cuando n tiende a ser muy grande.
Repita el experimento de la urna con n=1000 intentos, y determine las frecuencias relativas de cada pelota. Use el algoritmo de la hoja «modelo de probabilidad» para generar el archivo de aleatorios
Ejercicios
Referencia: Gubner 1.1 p.18
Se lanza una moneda 100 veces, y se anota la frecuencia relativa de las caras (0) y sellos (1). La experiencia nos indica que la frecuencia relativa debería estar alrededor de 0.5. Realice el experimento simulado y compare los resultados.
Muchos modelos de señales muestran variaciones impredecibles y un comportamiento aleatorio o probabilístico.
En un experimento aleatorio, repetido sobre las mismas condiciones, entrega resultados que varian de forma poco predecible. Si consideramos realizar un modelos determinístico, encontramos que no es el más apropiado al predecir siempre el mismo resultado.
Experimento aleatorio clásico
En un ánfora se colocan tres pelotas de igual tamaño numeradas con 0, 1 y 2, se agita el ánfora para que no exista un «orden» predeterminado.
Al sacar una pelota cualquiera, no se podría predecir cual será, luego se anota el número de la pelota y se la devuelve al ánfora.
El resultado del experimento, cada vez que se realiza un intento, es un número que pertenece al grupo S={0,1,2}, que formalmente se lo denomina espacio muestral
El resultado del experimento para 20 intentos, se obtiene un vector con los valores siguientes:
Simulación en Python
Para simular el experimento con python, se puede usar el generador de enteros aleatorios np.random.randint().
# Simulador de Urna con pelotas 0,1,2 # propuesta: edelros@espol.edu.ecimport numpy as np
import matplotlib.pyplot as plt
# INGRESO
intentos = 20
pelotas = 3
# PROCEDIMIENTO
aleatorios = np.random.randint(pelotas,size=intentos)
# SALIDAprint(aleatorios)
# grafica
plt.plot(aleatorios,'o-')
plt.xlabel('intentos')
plt.ylabel('aleatorios')
plt.margins(0.1)
plt.show()
Si diseña un equipo para transmitir el tono «la» de 440Hz, debe tomar en cuenta que no siempre es posible obtener una expresión matemática «sencilla» de una señal.
Escuche los siguientes sonidos, vea sus gráficas y realice sus observaciones:
Nota: si el archivo.wav tiene etiquetas, se descartan para el ejercicio. por lo que se presenta una advertencia «WavFileWarning».
c:\python34\lib\site-packages\scipy\io\wavfile.py:179: WavFileWarning: Chunk (non-data) not understood, skipping it.
WavFileWarning)
Se observa que:
todos los «sonidos» son monofónicos,
las frecuencias de muestreo de cada archivo son diferentes
las dimensiones de cada «sonido» o cantidad de muestras indican una duración diferente
Se recomienda escuchar nuevamente cada sonido y decidir el intervalo de observación para las gráficas. Se propone usar una ventana de observación o fragmento desde de los 1.2 segundos por 11 milisegundos.
Se extraen los fragmentos de cada sonido y se determinan los tiempos a los que corresponden basados en cada dt, con lo que se grafican los fragmentos.
Es una señal que puede ser escrita de forma determinística como:
x(t)=cos(2πft)
Siendo f=440 Hz
La señal de 440 Hz es usada como la frecuencia referencia para la afinación de todos los instrumentos musicales desde 1936, adoptada por ISO en 1955 y reafirmado en 1975 como ISO16.
Algoritmo en Python – Graficar en dominio del tiempo
El audio del archivo.wav se puede escuchar con windows media player. Para visualizar la forma de señal del archivo.wav, se puede usar python para leer el archivo y graficar la señal en el dominio del tiempo.
El archivo.wav debe encontrarse en el mismo directorio que el archivo .py de python.
# Señales analógicasimport numpy as np
import matplotlib.pyplot as plt
import scipy.io.wavfile as waves
# INGRESO# archivo = input('archivo.wav a leer: ')
archivo = '440Hz_44100Hz_16bit_05sec.wav'
[muestreo, sonido] = waves.read(archivo)
# SALIDA - Observación intermediaprint('frecuencia de muestreo: ', muestreo)
print('dimensiones de matriz: ', np.shape(sonido))
frecuencia de muestreo: 44100
dimensiones de matriz: (220500,)
Nota: si el archivo.wav tiene etiquetas, se descartan para el ejercicio. por lo que se presenta una advertencia «WavFileWarning».
c:\python34\lib\site-packages\scipy\io\wavfile.py:179: WavFileWarning: Chunk (non-data) not understood, skipping it. WavFileWarning)
Existen 220500 muestras en el archivo del ejemplo, que en una sola gráfica puede resultar muy denso y poco observable. La gráfica se limitará a una ventana de tiempo del archivo definida en el rango de tiempo por inicio y termina, con valores en el rango de milisegundos.
El rango de la ventana puede cambiarse si es de interés usar otro rango de tiempo.
De las dimensiones de la matriz (220500,), se determina que el audio es de tipo monofónico al tener una sola dimensión. El audio estéreo tiene dimensiones de (n,2).
Aunque la señal es simple de ver, y oir, la riqueza del sonido de un instrumento corresponde a otros componentes de frecuencia que los hace únicos al oido estando en la misma nota musical.
Los archivos para ser leidos deben copiarse al directorio donde se encuentra el algoritmo.py.
Los parametros obtenidos del archivo.txt son:
ΔY = deltaY = datos[1]
Δt = datos[0]
muestras en el rango de observación: k como tamaño del arreglo yentrada.
% matplotlib inline
# Modulacion Sigma-Delta Decodificador
# entrada y[n], salida: x[n]
# propuesta:edelros@espol.edu.ec
import numpy as np
import matplotlib.pyplot as plt
import scipy.io.wavfile as waves
# INGRESO - Lectura de datos desde archivos
# archivoparam = intput('archivo.txt de parametros: ')
# archivodatos = input('archivo.txt de datos: ')
# arhivoaudio = input('archivo.wav a grabar: ')
archivoparam = 'deltasigma_parametros.txt'
archivodatos = 'deltasigma_datos.txt'
archivoaudio = 'sigmadeltaaudio.wav'
param = np.loadtxt(archivoparam,dtype=float)
yentrada = np.loadtxt(archivodatos,dtype=int)
Decodificar la señal +1 y -1 en el vector yentrada consiste en acumular los valores de la secuencia en xdigital usando escalones de tamaño deltaY.
Se procede de forma semejante para los tiempos td usados para el eje de las abscisas usando pasos deltaT.
# PROCEDIMIENTO
deltaT = param[0] # Tamaño delta en eje tiempo
deltaY = param[1] # Tamaño delta en eje Y
k = len(yentrada) # número de muestras
xdigital = np.zeros(k, dtype='int16')
punto = np.zeros(k, dtype=int) # número de muestras
td = np.zeros(k, dtype=float) # tiempo muestreo
# DECOdifica Sigma-Delta
xdigital[0] = yentrada[0]
punto[0] = 0
td[0] = 0
for i in range(1,k):
punto[i] = i
td[i] = deltaT*i
xdigital[i] = xdigital[i-1]+yentrada[i]*deltaY
El resultado puede ser observado en los vectores, en una gráfica de los vectores y en un archivoaudio.wav.
que es tipo monofónica. La modulación delta-sigma se aplica usando los parámetros:
ΔY = deltaY = 0.3 * max(sonido)
rango de observación [0, 0.002] segundos para la gráfica
muestras en el rango de observación: k
Nota: En el caso de sonido estéreo, con dimensión kx2, se usa solo un canal (sonido[:,0])
Hay que tomar en cuenta que para graficar, la misma se satura con un número grande de muestras . Para escuchar el sonido el tiempo de ‘termina’ será de al menos 1 segundo.
Algoritmo en Python
A continuación se describe el algoritmo para codificar en Python, que el el bloque de ingreso especifica el rango de observación en segundos.
% matplotlib inline
# Modulacion Delta Sigma - Codificador
# entrada x(t), salida: y[n]
# propuesta:edelros@espol.edu.ec
import numpy as np
import matplotlib.pyplot as plt
import scipy.io.wavfile as waves
# INGRESO
# archivo = input('archivo de sonido:' )
archivo = '440Hz_44100Hz_16bit_05sec.wav'
muestreo, sonido = waves.read(archivo)
# rango de observación en segundos
inicia = 0
termina = 0.002
canal = 0 #Usar un canal en caso de estereo
c:\python34\lib\site-packages\scipy\io\wavfile.py:179: WavFileWarning: Chunk (non-data) not understood, skipping it. WavFileWarning)
Si el archivo.wav contiene etiquetas, el proceso de lectura emite una advertencia (WavFileWarning) que no afecta al proceso de codificación.
Se extrae solo una porción del sonido que contiene las k muestras en el intervalo de observación y se envía graficar.
La señal xdigital y ysalida se calcula con las diferencias entre el valor de cada muestra analógica y el valor xdigital anterior; se toma en cuenta solo el signo para acumular y generar el vector ysalida.
# PROCEDIMIENTO
# Codificar Sigma-Delta
deltaY = 0.1*np.max(sonido)
deltaT = 1/muestreo
# Usar un canal en caso de estereo
canales=sonido.shape
cuantos=len(canales)
canal = 0
if (cuantos==1): # Monofónico
uncanal=sonido[:]
if (cuantos==2): # Estéreo
uncanal=sonido[:,canal]
# Extrae solo una porcion del sonido
donde = int(inicia/deltaT)
# tiempo muestreo de la señal analógica
t = np.arange(inicia,termina,deltaT)
k = len(t)
muestra = np.copy(uncanal[donde:donde+k])
# Señal Digital
xdigital = np.zeros(k, dtype=float)
ysalida = np.zeros(k, dtype=int)
for i in range(1,k):
diferencia = muestra[i]-xdigital[i-1]
if (diferencia>0):
signo = 1
else:
signo = -1
xdigital[i] = xdigital[i-1]+signo*deltaY
ysalida[i] = signo
parametros=np.array([deltaT,deltaY,k])
Los resultados de pueden observar de dos maneras, en los archivos.txt de parámetros y datos, o en una gráfica.
para graficar y no saturar el gráfico, se recomienda observar solo una porción de los datos determinada por ‘verdesde’ y ‘verhasta’ muestras desde [0,k-1]
verdesde=0
verhasta=90
# Graficar
plt.figure(1) # define la grafica
plt.suptitle('Codificador Delta-Sigma')
plt.subplot(211) # grafica de 2x1 y subgrafica 1
plt.ylabel('x(t), x[n]')
plt.plot(t[verdesde:verhasta],muestra[verdesde:verhasta], 'g')
plt.step(t[verdesde:verhasta],xdigital[verdesde:verhasta], where='post',color='m') # Puntos x[n]
plt.subplot(212) # grafica de 2x1 y subgrafica 2
plt.ylabel('y[n]')
#plt.plot(ysalida, 'bo') # Puntos y[n]
plt.axis((verdesde,verhasta,-1.1,1.1))
puntos=np.arange(verdesde,verhasta,1) #posicion eje x para escalon
plt.step(puntos[verdesde:verhasta],ysalida[verdesde:verhasta], where='post')
plt.show()
Morse – Decodificador de un mensaje Morse desde el sonido en archivo.wav
Para analizar el sonido de un mensaje Morse, se procede de forma semejante al decodificador de un tono morse, revisando la duración de cada símbolo: punto o espacio tienen duración de una marca, la raya tiene una duración de 3 marcas, los simbolos se separan por una pausa con un intervalo de una marca.
Se divide el sonido del archivo.wav en partes o ventanas, usando la funcion separaventanas().
Cada ventana se analiza con la transformada de Fourier para encontrar la frecuencia del sonido mas fuerte o marcas de puntos '.', '-', ó ' '. El resultado de análisis de las ventanas es un vector con la marca de frecuencia más fuerte para cada una de las ventanas. (función marcasdeventanas())
La frecuencia o ‘tono frecuente’ se usa como una referencia para diferenciar si existe una señal para un punto o raya Morse. El tono frecuente se obtiene mediante la «moda» en el vector marcas usando scipy.stats.itemfreq() y usando la fila donde donde está la mayor cuenta usando np.argmax().
Para facilitar el análisis de un tono en una marca, se realiza una conversión a una secuencia binaria, usando '1' ó ,'0' si se encuentra el ‘tono frecuente’ en cada ventana.
Para el ejemplo, las ventanas se convierten al vector:
[1 1 1 1 1 1 0 0]
Al contar la duracion de cada intevalo, se determina si el sonido corresponde a un punto, una raya o una espacio:
Ejemplo:
. ... .--. --- .-..
Algoritmo en python
El sonido se obtiene de un archivo.wav del generador de tonos para un punto '.' o raya '-', del cual se leen los datos de sonido, y frecuencia de muestreo.
Las funciones para cada uno de los procesos descritos se encuentran a continuación
# Código Morse - Determina un mensaje desde sonido.wav# propuesta: edelros@espol.edu.ecimport numpy as np
import matplotlib.pyplot as plt
import scipy.io.wavfile as waves
import scipy.fftpack as fourier
import scipy.stats as stats
defseparaventanas(sonido, muestreo, duracion = 0.04, partes = 8 ):
# Divide el sonido en partes o ventanas para estudio# duracion = 0.04 # segundos de un punto# partes = 8 # divisiones de un punto# Extraer ventanas para espectro por partes
tventana = duracion/partes
mventana = int(muestreo * tventana) # muestras de una ventana# Ajuste de muestras de ventanas para matriz
anchosonido = int(len(sonido)/mventana)*mventana
sonidoajuste = np.resize(sonido,anchosonido)
ventanas = np.reshape(sonidoajuste,(-1,mventana))
# -1 indica que calcule las filasreturn(ventanas)
defmarcasdeventanas(ventanas, muestreo):
# Analiza con FFT todas las ventanas del sonido
filas, columnas = np.shape(ventanas)
marcas = np.zeros(filas,dtype=int)
# Espectro de Fourier de cada ventana# frecuencias para eje frq = fourier.fftfreq(mventana, 1/muestreo)
frq = fourier.fftfreq(columnas, 1/muestreo)
for f inrange(0,filas,1):
xf = fourier.fft(ventanas[f])
xf = np.abs(xf) # magnitud de xf
tono = np.argmax(xf) # tono, frecuencia mayor
marcas[f] = frq[tono]
return(marcas)
defsecuenciabinaria(marcas):
# Busca el tono frecuente mayor que cero
tonosventana = stats.itemfreq(marcas)
donde = np.argmax(tonosventana[1:,1])
tonopunto = tonosventana[donde+1,0]
# Convierte los tonos a secuencia binaria
secuencia = ''for valor in marcas:
if (valor==tonopunto):
secuencia = secuencia + '1'else:
secuencia = secuencia + '0'return(secuencia)
defduracionmarcas(secuencia):
# Duración de cada marca
conteo = []
caracter = '1'if (len(secuencia)>0):
caracter = secuencia[0]
i = 0
k = 0
while (i<len(secuencia)):
if (secuencia[i]==caracter):
k = k + 1
else:
conteo.append([int(caracter),k])
if (caracter=='1'):
caracter = '0'else:
caracter = '1'
k = 1
i = i+1
conteo = np.array(conteo)
return(conteo)
defmarcasmorse(conteo):
# Determina la base de un tono
veces = stats.itemfreq(conteo[:,1])
donde = np.argmax(veces[:,1])
base = veces[donde,0]
# Genera el código Morse
tolera = 0.2 # Tolerancia en relacion
bajo = 1-tolera
alto = 1+tolera
morse = ''for j inrange(0,len(conteo)):
relacion = conteo[j,1]/base
simbolo = conteo[j,0]
if (simbolo==1):
if (relacion>(1*bajo) and relacion<(1*alto)):
morse = morse + '.'if (relacion>(3*bajo) and relacion<(3*alto)):
morse = morse + '-'if simbolo==0 :
if (relacion>(3*bajo) and relacion<(3*alto)):
morse = morse + ' 'if (relacion>(7*bajo) and relacion<(7*alto)):
morse=morse+' 'return(morse)
# INGRESO # archivo = input('nombre del archivo: ')
archivo = 'morsetonoESPOL.wav'
muestreo, sonido = waves.read(archivo)
# PROCEDIMIENTO
ventanas = separaventanas(sonido, muestreo)
marcas = marcasdeventanas(ventanas, muestreo)
secuencia = secuenciabinaria(marcas)
conteo = duracionmarcas(secuencia)
morse = marcasmorse(conteo)
# SALIDAprint('codigo en morse: ')
print(morse)
Para facilitar el proceso de detección morse se convierten las marcas a una secuencia binaria, usando como referencia la frecuencia del tono de un punto determinada en la sección anterior.
Cambiar la secuencia a simbolos morse, consiste en determinar la relación entre las veces que aparecen los ‘1’s y ‘0’s, el resultado se puede guardar en un arreglo.
La base de cuántas marcas corresponden a un punto de estima como la «moda» de las veces en que se repite un uno o un cero.
una pausa se intercala al final del cada símbolo convertido a un tono.
Para la detección de un tono, se revisa la frecuencia activa en varios intervalos de tiempo o ventanas.
Al contar el número de intevalos (marcas) que aparece un tono vs. el número de intervalos de pausa, se podrá determinar si el sonido corresponde a un punto, una raya o un espacio. Es el proceso contrario a lo realizado para crear el tono, contanto las marcas, la relacion entre marcas y pausas indicará el símbolo recibido.
Algoritmo en python
El sonido para el análisis se obtiene de un archivo.wav, creado con el generador de tonos para un punto '.' o raya '-'. La lectura del archivo proporciona los datos de sonido, y frecuencia de muestreo.
Para el análisis, se supondrá que la ventana tiene la mitad de la duración de un punto, para el ejercicio es 0.04 segundos. Así, una ventana tiene la mitad de muestras que una marca'.' o pausa, obteniendo al menos dos ventanas por cada punto o pausa y tres ventanas por cada raya. Se cuentan las ventanas que tienen sonido o pausa y se comparan para determinar el símbolo.
Que existan al menos dos ventanas por cada punto tienen relación con el muestreo de Nyquist.
Lo expuesto se prueba primero analizando solo una ventana:
# Código Morse - Determina un simbolo a partir de UN sonido# propuesta: edelros@espol.edu.ecimport numpy as np
import matplotlib.pyplot as plt
import scipy.io.wavfile as waves
import scipy.fftpack as fourier
# INGRESO # archivo = input('nombre del archivo: ')# archivo = 'morsetonopunto.wav'
archivo = 'morsetonoraya.wav'
muestreo, sonido = waves.read(archivo)
# PROCEDIMIENTO
duracion = 0.04 # segundos de un punto
partes = 2 # divisiones de un punto# Extraer ventanas para espectro por partes
tventana = duracion/partes
mventana = int(muestreo * tventana) # muestras de una ventana
# Observacion intermedia
unaventana = sonido[0:mventana]
# SALIDA GRAFICAprint('muestreo:', muestreo)
print('muestras en sonido', len(sonido))
print('intervalo de ventana: ', tventana)
print('muestras por ventana: ', mventana)
plt.figure(1)
plt.subplot(211)
plt.plot(sonido, label='sonido')
plt.xlabel('muestras')
plt.ylabel('sonido y ventanas')
k = 1
for i inrange(0,int(len(sonido)/mventana)):
plt.axvline(x=k*mventana, color='r', linestyle='--')
k = k + 1
plt.subplot(212)
plt.plot(unaventana)
plt.xlabel('muestras')
plt.ylabel('unaventana')
plt.margins(0)
plt.show()
muestreo: 11025
muestras en sonido 1764
intervalo de ventana: 0.02
muestras por ventana: 220
Análisis de una ventana
Para determinar si existe un tono en una ventana, se analiza la existencia de una señal observada en el dominio de la frecuencia.
En el dominio de la frecuencia se requiere usar la «Transformada de Fourier» de la señal en la ventana. La Transformada rápida de Fourirer fft() se encuentra disponible en las librerias «scipy» y el rango de frecuencias para la gráfica se obtiene con fftfreq().
En el resultado se busca la frecuencia (np.argmax()) con mayor magnitud (np.abs()) y se determina si corresponde al «tono» morse esperado.
# Espectro de Fourier de una ventana [0:nventana]
xf = fourier.fft(unaventana)
xf = np.abs(xf) # magnitud de xf
tono = np.argmax(xf) # tono en donde# frecuencias para eje
frq = fourier.fftfreq(mventana, 1/muestreo )
# SALIDA GRAFICA /Observacion intermediaprint('frecuencia tono (Hz): ',frq[tono])
plt.figure(2)
plt.xlabel('Frecuencia Hz')
plt.ylabel('Magnitud')
plt.stem(frq,xf)
plt.show()
frecuencia tono (Hz): 451.022727273
Analisis del archivo.wav de un símbolo Morse
Para analizar todo el archivo de sonido, los datos se segmentan por «ventanas» en una matriz, cada fila corresponde a las muestras de una ventana de tiempo.
La matriz debe ser rectangular, por lo que es necesario confirmar que el número de muestras por ventana sean iguales. En caso de no ser así, se realiza un ajuste del ancho del sonido para que el numero de muestras en el archivo sea múltiplo del número de muestras por ventana.
Observar el resultado de crear varias ventanas de tiempo, es semejante a tener un cuaderno en el que cada hoja tiene una imagen con los valores de la señal en cada ventana.
# Ajuste de muestras de ventanas para matriz
antes = len(sonido)
anchosonido = int(len(sonido)/mventana)*mventana
sonido = np.resize(sonido,anchosonido)
ventanas = np.reshape(sonido,(-1,mventana))
# -1 indica que calcule las filas# SALIDAprint('muestras en sonido: ', antes)
print('muestras recortadas a: ', anchosonido)
print('ventanas: ')
print(ventanas)
Para observar el resultado se usa una gráfica 3d de superficie.
# Observacion intermediafrom mpl_toolkits.mplot3d import Axes3D
from matplotlib import cm
nfilas,ncolumnas = np.shape(ventanas)
x = np.arange(0,ncolumnas) # valores para cada eje
y = np.arange(0,nfilas)
X, Y = np.meshgrid(x, y) # valores para punto.
Z = ventanas # Matriz# Salida, gráfico 3d
fig = plt.figure(3,figsize=plt.figaspect(0.5))
ax = fig.add_subplot(1, 2, 1, projection='3d')
surf = ax.plot_surface(X, Y, Z, rstride=1, cstride=1,
cmap=cm.coolwarm, linewidth=0,
antialiased=False)
plt.xlabel('muestras')
plt.ylabel('ventana')
plt.show()
El proceso de analisis se repite para cada ventana (una fila) en la matriz, guardando cada resultado en un vector «marcas«.
# Analiza con FFT todas las ventanas del sonido
filas, columnas = np.shape(ventanas)
marcas = np.zeros(filas,dtype=float)
# Espectro de Fourier de cada ventana# frecuencias para eje
frq = fourier.fftfreq(mventana, 1/muestreo)
for fila inrange(0,filas,1):
xf = fourier.fft(ventanas[fila])
xf = np.abs(xf) # magnitud de xf
tono = np.argmax(xf) # tono, frecuencia mayor
marcas[fila] = frq[tono]
Observamos el resultado en una gráfica 2D.
# SALIDA GRAFICA /Observacion intermediaprint(marcas)
plt.figure(4)
plt.stem(marcas)
plt.xlabel('número de marca')
plt.ylabel('frecuencia del tono mas fuerte.')
plt.show()
Para determinar si es punto '.', raya'-' o espacio ' ', se cuentan los noceros y los ceros del arreglo de marcas, la relación entre ellos determina el símbolo que representa el sonido del archivo.wav.
Nota: Observe que cualquier frecuencia en una ventana marca un punto '.'. El detalle es tratado en «Morse – Decodificador del sonido de un mensaje» en : # Busca frecuencia del tono punto.
Programa resumido
El resultado del programa resumido se muestra a continuación:
# Código Morse - Determina un simbolo a partir de UN sonido# propuesta: edelros@espol.edu.ecimport numpy as np
import matplotlib.pyplot as plt
import scipy.io.wavfile as waves
import scipy.fftpack as fourier
# INGRESO # archivo = input('nombre del archivo: ')# archivo = 'morsetonopunto.wav'
archivo = 'morsetonoraya.wav'
muestreo, sonido = waves.read(archivo)
# PROCEDIMIENTO
duracion = 0.04 # segundos de un punto
partes = 2 # divisiones de un punto# Extraer ventanas para espectro por partes
tventana = duracion/partes
mventana = int(muestreo * tventana) # muestras de una ventana# Ajuste de muestras de ventanas para matriz
antes = len(sonido)
anchosonido = int(len(sonido)/mventana)*mventana
sonido = np.resize(sonido,anchosonido)
ventanas = np.reshape(sonido,(-1,mventana))
# -1 indica que calcule las filas# Analiza con FFT todas las ventanas del sonido
filas, columnas = np.shape(ventanas)
marcas = np.zeros(filas,dtype=float)
# Espectro de Fourier de cada ventana# frecuencias para eje
frq = fourier.fftfreq(mventana, 1/muestreo)
for f inrange(0,filas,1):
xf = fourier.fft(ventanas[f])
xf = np.abs(xf) # magnitud de xf
tono = np.argmax(xf) # tono, frecuencia mayor
marcas[f] = frq[tono]
# determina el simbolo
noceros = np.count_nonzero(marcas)
ceros = len(marcas)-noceros
relacion = int(noceros/ceros)
if (relacion==3):
simbolo = '-'if (relacion==1):
simbolo = '.'if (relacion==0):
simbolo = 'espacio'# SALIDAprint('relacion de marcas:', relacion)
print('simbolo morse detectado: ', simbolo)