[ reconstrucción ] [ sobre-muestreo ] [ sub-muestreo ] [ algoritmo ] [ gráfica ]
..
1. Muestreo y Reconstrucción
La salida de un convertidor contínuo a discreto C-D es una señal discreta en tiempo con un número infinito de alias. Al observar un gráfico de espectro de frecuencias ω es necesario recordar que existen más señales fuera del intervalo mostrado.
Considere una señal x(t) simple, de una sola frecuencia, cuyo espectro contiene dos líneas a ±ω0 con amplitudes de (1/2 A e±ω0)
x(t) = A cos\Big(\omega_0 t + \varphi \Big) x[n] = x(n/f_s) A cos((\omega_0 /f_s)*n + \varphi = \frac{1}{2}A e^{j\varphi}e^{j(\omega_0/f_s)n} + \frac{1}{2}A e^{-j\varphi}e^{j(-\omega_0/f_s)n}también tiene dos frecuencias ω en ±ω0/fs, con todos los alias en
±ω0/fs + 2πk, siendo k= ±1, ±2, ±3, …
[ reconstrucción ] [ sobre-muestreo ] [ sub-muestreo ] [ algoritmo ] [ gráfica ]
..
2. Sobre-muestreo
Referencia: McClellan ejercicio 4.5 p133
Realice el espectro contínuo de la señal x(t), la señal muestreada a fs=500 y el espectro discreto resultante. Observe que la frecuencia de muestreo es superior a la frecuencia de la señal x(t) para evitar el problema de aliasing. El proceso se conoce como sobre-muestreo. Para el ejercicio se usa una frecuencia 2.5 veces la recomendada por Nyquist.
x(t) = cos\Big(2\pi (100) t + \pi /3 \Big)
Los convertidores D-C trasforman el espectro discreto en tiempo a una salida de espectro contínuo en el tiempo, seleccionando solo un par de las líneas de todas las posibilidades mostradas.
Para ser consistentes con la operación Digital hacia Analógico (D-A) se tomará siempre la frecuencia mas baja de cada grupo de alias o alias principal con |ω|<π.
El resultado siempre se encontrará entre [-fs/2, +fs/2].
En resumen, en sobre-muestreo la frecuencia original f0 es menor que fs/2 permite la reconstrucción mas exacta.
Para el análisis con el algoritmo se obtienen los siguientes parámetros para realizar la gráfica.
señal(t): cos(100*t*(2*pi) + pi/3)
espectro:
freq : [-100. 100.]
ampl : [1/2 1/2]
fase : [-pi/3 pi/3]
etiq : ['$\\frac{1}{2}$$ e^j(- \\frac{\\pi}{3})$'
'$\\frac{1}{2}$$ e^j(\\frac{\\pi}{3})$']
freq_max : 100.0
freq_min : 100.0
x_expand : cos(2*pi*n/5 + pi/3)
muestreo_tipo: sobre-muestreo
muestreo: 11
tamaño ki: 11 ; tamaño ti: 121
fs: 500 ; dt: 0.002 ; fs_veces: 12
x(t): cos(100*t*(2*pi) + pi/3)
x[n]: cos(n*(2*pi)/5 + pi/3)
señal[n]: cos(n*(2*pi)/5 + pi/3)
espectro_n:
freq : [-2.4 -1.6 -0.4 0.4 1.6 2.4]
ampl : [1/2 1/2 1/2 1/2 1/2 1/2]
fase : [-pi/3 pi/3 -pi/3 pi/3 -pi/3 pi/3]
etiq : ['$\\frac{1}{2}$$ e^j(\\frac{\\pi}{3})$'
'$\\frac{1}{2}$$ e^j(- \\frac{\\pi}{3})$'
'$\\frac{1}{2}$$ e^j(- \\frac{\\pi}{3})$'
'$\\frac{1}{2}$$ e^j(\\frac{\\pi}{3})$'
'$\\frac{1}{2}$$ e^j(- \\frac{\\pi}{3})$'
'$\\frac{1}{2}$$ e^j(\\frac{\\pi}{3})$']
freq_factor : pi
freq_max : 2.4
freq_min : 0.4
BW : 2.0
freq_alias : [ True True False False True True]
señal reconstruida con alias principal:
cos(100*t*(2*pi) + pi/3)
Las instrucciones del bloque de ingreso para el algoritmo son:
# INGRESO # usar np.pi para evitar frecuencia simplificada por sympy. x = sym.cos(DosPi*100*t + pi/3) x_etiqueta = 'x(t)' ; x_etiqueta_n = 'x[n]' # señal muestreada fs = 500 # muestreo discreto, freq_sampling fs_veces = 12 # sobremuestreo para x(t) muestrasN = 10 + 1 # para x[n] aliasn = 1 # alias añadidos a x[n] tolera = 1e-9 # tolerancia a error, números iguales
[ reconstrucción ] [ sobre-muestreo ] [ sub-muestreo ] [ algoritmo ] [ gráfica ]
..
3. Sub-muestreo y aliasing
Cuando fs < 2f0 la señal se encuentra sub-muestreada y se produce el aliasing. Para el ejercicio anterior, si el ejercicio lo realiza con fs = 80Hz, siendo f0=100 Hz, la tasa de muestreo es menor que la frecuencia de Nyquist y se presenta la señal de alias.
las instrucciones para el bloque de ingreso del algoritmo son las mismas que en el ejercicio anterior, excepto que fs = 80 Hz
# señal muestreada fs = 80 # muestreo discreto, freq_sampling
el resultado con el algoritmo es una señal reconstruida de tan solo 20Hz
x_0(t) = cos\Big(2\pi (20) t + \pi /3 \Big)señal(t): cos(100*t*(2*pi) + pi/3)
espectro:
freq : [-100. 100.]
ampl : [1/2 1/2]
fase : [-pi/3 pi/3]
etiq : ['$\\frac{1}{2}$$ e^j(- \\frac{\\pi}{3})$'
'$\\frac{1}{2}$$ e^j(\\frac{\\pi}{3})$']
freq_max : 100.0
freq_min : 100.0
x_expand : cos(5*pi*n/2 + pi/3)
muestreo_tipo: sub-muestreo
muestreo: 11
tamaño ki: 11 ; tamaño ti: 121
fs: 80 ; dt: 0.0125 ; fs_veces: 12
x(t): cos(100*t*(2*pi) + pi/3)
x[n]: cos(5*n*(2*pi)/4 + pi/3)
señal[n]: cos(5*n*(2*pi)/4 + pi/3)
espectro_n:
freq : [-4.5 -2.5 -0.5 0.5 2.5 4.5]
ampl : [1/2 1/2 1/2 1/2 1/2 1/2]
fase : [-pi/3 -pi/3 -pi/3 pi/3 pi/3 pi/3]
etiq : ['$\\frac{1}{2}$$ e^j(\\frac{\\pi}{3})$'
'$\\frac{1}{2}$$ e^j(- \\frac{\\pi}{3})$'
'$\\frac{1}{2}$$ e^j(- \\frac{\\pi}{3})$'
'$\\frac{1}{2}$$ e^j(- \\frac{\\pi}{3})$'
'$\\frac{1}{2}$$ e^j(\\frac{\\pi}{3})$'
'$\\frac{1}{2}$$ e^j(\\frac{\\pi}{3})$']
freq_factor : pi
freq_max : 4.5
freq_min : 0.5
BW : 4.0
freq_alias : [ True True False False True True]
señal reconstruida con alias principal:
cos(20*t*(2*pi) + pi/3)
[ reconstrucción ] [ sobre-muestreo ] [ sub-muestreo ] [ algoritmo ] [ gráfica ]
..
4. Folding o plegado de señal por sub-muestreo
Para una señal x(t) de 100 Hz se aplica muestreo insuficiente en el intervalo [f0,2f0] se presenta aliasing denominado «folding» o plegado de señal. Para el ejercicio de prueba y con fs = 125, se tiene:
x(t) = cos\Big(2\pi (100) t + \pi /3 \Big)
Los componentes de frecuencia entre ±π también son ±0.4π, sin embargo el componente en +0.4π es un alias del componente de la frecuencia negativa -1.6π, que genera el «plegado» o «folding». La reconstrucción de la señal se realiza a f = 0.4(fs/2π) = fs/5 = 25 Hz. También se observa una inversión de fase. El resultado es semejante a una señal muestreada de 25 Hz con inversión de fase.
señal(t): cos(100*t*(2*pi) + pi/3)
espectro:
freq : [-100. 100.]
ampl : [1/2 1/2]
fase : [-pi/3 pi/3]
etiq : ['$\\frac{1}{2}$$ e^j(- \\frac{\\pi}{3})$'
'$\\frac{1}{2}$$ e^j(\\frac{\\pi}{3})$']
freq_max : 100.0
freq_min : 100.0
x_expand : cos(8*pi*n/5 + pi/3)
muestreo_tipo: sub-muestreo y folding
muestreo: 11
tamaño ki: 11 ; tamaño ti: 121
fs: 125 ; dt: 0.008 ; fs_veces: 12
x(t): cos(100*t*(2*pi) + pi/3)
x[n]: cos(4*n*(2*pi)/5 + pi/3)
señal[n]: cos(4*n*(2*pi)/5 + pi/3)
espectro_n:
freq : [-3.6 -1.6 -0.4 0.4 1.6 3.6]
ampl : [1/2 1/2 1/2 1/2 1/2 1/2]
fase : [-pi/3 -pi/3 pi/3 -pi/3 pi/3 pi/3]
etiq : ['$\\frac{1}{2}$$ e^j(\\frac{\\pi}{3})$'
'$\\frac{1}{2}$$ e^j(- \\frac{\\pi}{3})$'
'$\\frac{1}{2}$$ e^j(- \\frac{\\pi}{3})$'
'$\\frac{1}{2}$$ e^j(- \\frac{\\pi}{3})$'
'$\\frac{1}{2}$$ e^j(\\frac{\\pi}{3})$'
'$\\frac{1}{2}$$ e^j(\\frac{\\pi}{3})$']
freq_factor : pi
freq_max : 3.6
freq_min : 0.3999999999999999
BW : 3.2
freq_alias : [ True True False False True True]
señal reconstruida con alias principal:
cos(25.0*t*(2*pi) - (pi/3))
Otro caso a observar es cuando fs=100Hz, el resultado es una constante, pues siempre se obtiene la muestra de la señal con la misma magnitud.
5. Algoritmo en Python
El algoritmo contiene componentes desarrollados en secciones anteriores. Para el espectro de frecuencias en dominio del tiempo se basa en lo descrito en Espectro – Operaciones en dominio de tiempo y frecuencia.
El desarrollo del muestreo usando la señal x(t) y fs toma como referencia lo descrito en Muestreo y aliasing.
La construcción del espectro discreto a partir de x[n] usa lo desarrollado en Espectro x[n] – Señal discreta en tiempo
# ejercicio 4.5 p133 Muestreo y Reconstrucción de señales sinusoidales # telg1034 DSP fiec-espol edelros@espol.edu.ec import numpy as np import matplotlib.pyplot as plt import sympy as sym import telg1034 as dsp # variables continuas from telg1034 import t,A,w,f,p,pi,DosPi,I,equivalentes # variables discretas from telg1034 import n # INGRESO # usar np.pi para evitar frecuencia simplificada por sympy. x = sym.cos(DosPi*100*t + pi/3) x_etiqueta = 'x(t)' ; x_etiqueta_n = 'x[n]' # señal muestreada fs = 500 # muestreo discreto, freq_sampling fs_veces = 12 # sobremuestreo para x(t) muestrasN = 10 + 1 # para x[n] aliasn = 1 # alias añadidos a x[n] tolera = 1e-9 # tolerancia a error, números iguales # PROCEDIMIENTO ---------------- # operaciones para espectro x(t) # Ref: Blog|Espectro-Operaciones en dominio de tiempo y frecuencia x_term = dsp.x_list_term_Add(x) Xe_term = dsp.cos_to_euler_one_term(x_term) x_term_expand = dsp.euler_to_cos_list(Xe_term) un_espectro = dsp.cos_spectrum_list(x_term) # espectro de cada señal freq = un_espectro['freq'] ampl = un_espectro['ampl'] fase = un_espectro['fase'] etiqueta = un_espectro['etiq'] # freq_max y freq_min freq_max = float(max(freq)) if len(freq[freq>0])>0: freq_min = float(min(freq[freq>0])) else: freq_min = 0 un_espectro['freq_max'] = freq_max un_espectro['freq_min'] = freq_min # revisa muestreo Nyquist if freq_max<=(fs/2): muestreo_tipo = 'sobre-muestreo' if freq_max>(fs): muestreo_tipo = 'sub-muestreo' if freq_max>(fs/2) and freq_max<=fs: muestreo_tipo = 'sub-muestreo y folding' # Muestreo de señales sinusoidales # Ref: blog|2.1 Muestreo y aliasing xn = x.subs(t,n/fs) # muestreo x[n] dt = 1/fs ki = np.arange(0,muestrasN,1,dtype=int) # muestreo x(t) dtc = dt/fs_veces # para x(t) ti = np.arange(0,ki[muestrasN-1]*dt + dtc, dtc) tki = ki*dt xt = sym.lambdify(t,x, modules=equivalentes) xki = xt(tki) xti = xt(ti) # operaciones para espectro x[n] # Ref: blog|2.2 Espectro x[n] – Señal discreta en tiempo x_term = dsp.x_list_term_Add(xn) Xe_term = dsp.cos_to_euler_one_term(x_term) x_term_expand = dsp.euler_to_cos_list(Xe_term) un_espectro_n = dsp.cos_spectrum_list(x_term) # espectro de cada señal freq_n = np.array(un_espectro_n['freq']) ampl_n = np.array(un_espectro_n['ampl']) fase_n = np.array(un_espectro_n['fase']) etiq_n = np.array(un_espectro_n['etiq']) un_espectro['x_expand'] = x_term_expand # Añadir alias a señal freqn_conteo = len(freq_n) freqn_alias = np.zeros(freqn_conteo,dtype=bool) for i in range(0,freqn_conteo,1): freq_i = freq_n[i] ampl_i = ampl_n[i] fase_i = fase_n[i] etiq_i = etiq_n[i] for k in range(1,aliasn+1,1): # añade freqm[i] + 2*pi*k freq_n = np.concatenate([freq_n,[freq_i+2*pi*k]]) ampl_n = np.concatenate([ampl_n,[ampl_i]]) fase_n = np.concatenate([fase_n,[fase_i]]) etiq_n = np.concatenate([etiq_n,[etiq_i]]) freqn_alias = np.concatenate([freqn_alias,[True]]) # añade freqm[i] - 2*pi*k freq_n = np.concatenate([freq_n,[-(freq_i.evalf()+2*np.pi*k)]]) ampl_n = np.concatenate([ampl_n,[ampl_i]]) fase_n = np.concatenate([fase_n,[-fase_i]]) texto = '$' + sym.latex(ampl_i)+'$' if fase_i!=sym.S.Zero: texto = texto+f'$ e^j('+sym.latex(fase_i)+')$' etiq_n = np.concatenate([etiq_n,[texto]]) freqn_alias = np.concatenate([freqn_alias,[True]]) # ordena frecuencias orden = np.argsort(freq_n) freq_n = freq_n[orden] ampl_n = ampl_n[orden] fase_n = fase_n[orden] etiq_n = etiq_n[orden] freqn_alias = freqn_alias[orden] # revisa factor pi en freq_n unfactor = sym.S.One ; factor_pi = False freqn_conteo = len(freq_n) for i in range(0,freqn_conteo,1): if freq_n[i].has(pi): factor_pi = True if factor_pi: freq_n = np.array(freq_n/pi,dtype=float) unfactor = pi # actualiza espectro de señal con factor pi,alias,orden un_espectro_n['freq_factor'] = unfactor un_espectro_n['freq'] = freq_n un_espectro_n['ampl'] = ampl_n un_espectro_n['fase'] = fase_n un_espectro_n['etiq'] = etiq_n # ancho de banda, freq_max y freq_min freq_posi = freq_n[freq_n>0] freq_max_n = float(max(freq_posi)) freq_min_n = 0 if len(freq_posi)>0: freq_min_n = float(min(freq_posi)) freq_centro_n = (freq_max_n + freq_min_n)/2 un_espectro_n['freq_max'] = freq_max_n un_espectro_n['freq_min'] = freq_min_n un_espectro_n['BW'] = freq_max_n - freq_min_n # aliasing, Revisa espectro de frecuencias freq_conteo = len(freq_n) freq_alias = np.zeros(freq_conteo,dtype=bool) for i in range (0,freq_conteo,1): unafreq = freq_n[i] k = 1 unalias = unafreq + k*2 while unalias<=freq_max_n: for j in range(i+1,freq_conteo,1): if abs(unalias - freq_n[j])<tolera: if abs(freq_n[i])>abs(freq_n[j]): freq_alias[i] = True else: freq_alias[j] = True k = k+1 unalias = unafreq + k*2 freq_alias_conteo = np.sum(freq_alias) un_espectro_n['freq_alias'] = freq_alias # reconstruye señal sin alias freq_alias0 = np.invert(freq_alias) freq_x0 = freq_n[freq_alias0] ampl_x0 = ampl_n[freq_alias0] fase_x0 = fase_n[freq_alias0] x0_cuenta = len(freq_x0) x0t = sym.S.Zero for i in range(0,x0_cuenta,1): if freq_x0[i]>=0: una_freq = freq_x0[i]*unfactor*(fs/(2*pi)) una_freq = dsp._float_is_int(una_freq) # si es entero x0t = x0t + 2*ampl_x0[i]*sym.cos(DosPi*una_freq*t+sym.UnevaluatedExpr(fase_x0[i])) x0tr = x0t.subs(pi,np.pi) x0tr = sym.lambdify(t,x0tr,modules=equivalentes) x0_ti = x0tr(ti) # SALIDA print('señal(t): ',x) print(' espectro:') for parametro in un_espectro: print(' ',parametro,':',un_espectro[parametro]) print('\nmuestreo_tipo: ',muestreo_tipo) print('muestreo:',muestrasN) print(' tamaño ki:',len(ki),'; tamaño ti:',len(ti)) print('fs:',fs,' ; dt:',dt, ' ; fs_veces:',fs_veces) print('x(t):',x) print('x[n]:',xn) print('\nseñal[n]: ',xn) print(' espectro_n:') for parametro in un_espectro_n: print(' ',parametro,':',un_espectro_n[parametro]) print('\nseñal reconstruida con alias principal:') print(x0t)
[ reconstrucción ] [ sobre-muestreo ] [ sub-muestreo ] [ algoritmo ] [ gráfica ]
..
6. Gráfica
Con los componentes descritos en la parte de algoritmos, se ajustan las gráficas de cada sección para cada sub-gráfica en la figura.
# GRAFICAS ---------------------------------------- # Grafica de espectro de frecuencias # Ref: blog|2.2 Espectro x[n] – Señal discreta en tiempo graf_dx = 0.12 fig_espectro = plt.figure() ampl_max = float(max(ampl)) freq_max = float(max(max(freq),fs))*(1+graf_dx/2) graf_fasor = fig_espectro.add_subplot(311) if freq_max!=0: graf_fasor.set_xlim([-freq_max*(1+graf_dx), freq_max*(1+graf_dx)]) else: graf_fasor.set_xlim([-1,1]) graf_fasor.set_ylim([0,ampl_max*(1+graf_dx*2)]) graf_fasor.grid() graf_fasor.axhline(0,color='black') graf_fasor.axvline(0,linestyle='dotted',color='grey') graf_fasor.stem(freq,ampl) # grafica magnitud for k in range(0,len(freq),1): # etiquetas de fasor graf_fasor.annotate(etiqueta[k],xy=(freq[k],ampl[k]), xytext=(0,5),textcoords='offset points',ha='center') # fs frecuencia de muestreo graf_fasor.stem(fs,ampl_max/2,linefmt = 'C3:', ) # fs graf_fasor.annotate('fs',xy=(fs,ampl_max/2), xytext=(0,5),textcoords='offset points',ha='center') graf_fasor.set_ylabel('Magnitud '+x_etiqueta) graf_fasor.set_xlabel('freq Hz', labelpad=0,) texto = '$' + sym.latex(x)+'$' +' ; fs='+str(fs) texto = texto + ' ; alias='+str(aliasn) graf_fasor.set_title(texto) # GRAFICA x(t), x[n], x(t)_alias0 graf_xt = fig_espectro.add_subplot(312) graf_xt.plot(ti,xti, color='gray',label='x(t)') graf_xt.plot(ti,x0_ti, color='blue',label='x(t)_alias0', linestyle='dotted') graf_xt.stem(tki,xki,linefmt = 'C0:',label='x[n]') graf_xt.set_xlabel('t') graf_xt.set_ylabel('amplitud') graf_xt.legend() graf_xt.grid() # GRAFICAS de espectro de frecuencias discretas--------- # Ref: blog|2.2 Espectro x[n] – Señal discreta en tiempo freq_n = un_espectro_n['freq'] ampl_n = un_espectro_n['ampl'] etiq_n = un_espectro_n['etiq'] freq_alias = un_espectro_n['freq_alias'] ampl_max = float(max(ampl_n))*(1+graf_dx*2) freq_max = float(max(un_espectro_n['freq_max'],1))*(1+graf_dx/2) # grafica graf_fasorn = fig_espectro.add_subplot(313) graf_fasorn.set_xlim([-freq_max,freq_max]) graf_fasorn.set_ylim([min([min(ampl_n),0]),ampl_max]) graf_fasorn.grid() graf_fasorn.axhline(0,color='black') graf_fasorn.axvline(0,linestyle='dotted',color='grey') # grafica magnitud freq_noalias = np.invert(freq_alias) graf_fasorn.stem(freq_n[freq_noalias],ampl_n[freq_noalias], linefmt = 'C0--') graf_fasorn.stem(freq_n[freq_alias],ampl_n[freq_alias], linefmt = 'C1--',label = 'alias') if un_espectro_n['freq_factor'] != sym.S.One: unfactor = r'$'+sym.latex(un_espectro_n['freq_factor'])+'$' freq_etiq = [] for un_num in freq_n: freq_etiq.append(f''+str(round(un_num,4))+unfactor) graf_fasorn.set_xticks(ticks=freq_n,labels=freq_etiq) for k in range(0,len(freq_n),1): # etiquetas de fasor graf_fasorn.annotate(etiq_n[k],xy=(freq_n[k],ampl_n[k]), xytext=(0,5),textcoords='offset points', ha='center') graf_fasorn.set_ylabel('Magnitud '+x_etiqueta_n) graf_fasorn.set_xlabel('freq rad') graf_fasorn.legend() plt.tight_layout() plt.show()
[ reconstrucción ] [ sobre-muestreo ] [ sub-muestreo ] [ algoritmo ] [ gráfica ]
