s2Eva2016TII_T4 resolver en dominio de frecuencia

Ejercicio: 2Eva2016TII_T4 resolver en dominio de frecuencia

Partiendo de,

y(t) = x(t) \circledast h(t) Y(\omega) = X(\omega) H(\omega) Y\Big(\frac{\omega}{3}\Big) = X\Big(\frac{\omega}{3}\Big) H\Big(\frac{\omega}{3}\Big)

se tiene también que:

g(t) = x(3t) \circledast h(3t) G(\omega) = \frac{1}{3}X\Big(\frac{\omega}{3}\Big)\text{ } \frac{1}{3}H\Big(\frac{\omega}{3}\Big) = \frac{1}{9}X\Big(\frac{\omega}{3}\Big)H\Big(\frac{\omega}{3}\Big) = \frac{1}{9}Y\Big(\frac{\omega}{3}\Big) G(\omega) = \frac{1}{3} \Bigg[\frac{1}{3}Y\Big(\frac{\omega}{3}\Big)\Bigg] g(t) = \frac{1}{3} y(3t)

A = 1/3
B = 3

 

s2Eva2016TII_T3 LTI DT sistemas en serie

Ejercicio: 2Eva2016TII_T3 LTI DT sistemas en serie

literal a

S1: w[n] = \frac{1}{2}w[n-1] + x[n] S2: y[n] = \alpha y[n-1] + \beta w[n] SG: y[n] = -\frac{1}{8}y[n-2] + \frac{3}{4}y[n-1]+x[n]

usando transformada z:

S1:

W(z) = \frac{1}{2} z^{-1} W(z) + X(z) W(z) - \frac{1}{2} z^{-1} W(z) = X(z) W(z)\Big[1 - \frac{1}{2} z^{-1} \Big] = X(z) W(z) = \frac{ X(z)}{\Big[1 - \frac{1}{2} z^{-1} \Big]}

S2:

Y(z) = \alpha z^{-1} Y(z) + \beta W(z) Y(z) - \alpha z^{-1} Y(z) = \beta W(z) Y(z) \Big[1 - \alpha z^{-1} \Big] = \beta W(z)

sustituyendo la ecuacion de S1 para W(z)

Y(z) \Big[1 - \alpha z^{-1} \Big] = \beta \frac{ X(z)}{\Big[1 - \frac{1}{2} z^{-1} \Big]} Y(z) \frac{1}{\beta} \Big[1 - \alpha z^{-1} \Big]\Big[1 - \frac{1}{2} z^{-1} \Big] = X(z) Y(z) \Big[ \frac{1}{\beta} - \frac{1/2+\alpha}{\beta} z^{-1}+ \frac{1}{2} \frac{\alpha}{\beta} z^{-2}\Big]= X(z)

SG:

Y(z) = -\frac{1}{8} z^{-2} Y(z) + \frac{3}{4} z^{-1}Y(z)+X(z) Y(z)+\frac{1}{8} z^{-2} Y(z) - \frac{3}{4} z^{-1}Y(z) = X(z) Y(z)\big[ 1 - \frac{3}{4} z^{-1} +\frac{1}{8} z^{-2} \Big] = X(z)

comparando con la ecuación de S2

\frac{1}{\beta} = 1 \beta = 1 \frac{1}{2} \frac{\alpha}{\beta} = \frac{1}{8} \alpha = \frac{2\beta}{8} = \frac{1}{4}

comprobar con

- \frac{1/2+\alpha}{\beta} = \frac{3}{4}

se confirma que  α = 1/4 y β=1

la función de transferencia es:

Y(z)\Big[ 1 - \frac{3}{4} z^{-1} +\frac{1}{8} z^{-2} \Big] = X(z) \frac{Y(z)}{X(z)} = \frac{1}{\Big[ 1 - \frac{3}{4} z^{-1} +\frac{1}{8} z^{-2} \Big]} H(z) = \frac{z^2}{z^2 - \frac{3}{4} z +\frac{1}{8}} \frac{H(z)}{z} = \frac{z}{z^2 - \frac{3}{4} z +\frac{1}{8}}

usando las raices para:

z^2 - \frac{3}{4} z +\frac{1}{8} = \Big[ z-\frac{1}{4}\Big] \Big[z-\frac{1}{2}\Big]

y la parte derecha de la ecuación:

\frac{z}{\Big[ z-\frac{1}{4}\Big] \Big[z-\frac{1}{2}\Big]} = \frac{C_1}{z-\frac{1}{4}}+\frac{C_2}{z-\frac{1}{2}}

despejando para C1 y haciendo z=1/4,

C_1 = \frac{z}{\Big[z-\frac{1}{2}\Big]} = \frac{\frac{1}{4}}{\Big[\frac{1}{4}-\frac{1}{2}\Big]} = -1

despejando para C2 y haciendo z=1/2,

C_2 = \frac{z}{\Big[z-\frac{1}{4}\Big]} = \frac{\frac{1}{2}}{\Big[\frac{1}{2}-\frac{1}{4}\Big]} = 2

se H(z) se resume en,

H(z) = - \frac{z}{z - \frac{1}{4}} +2\frac{z}{z-\frac{1}{2}}

Para obtener h[n] usando la antitransformada,

h[n]=z^{-1} \Big[ H(z) \Big] h[n]=z^{-1} \Big[- \frac{z}{z - \frac{1}{4}} +2\frac{z}{z-\frac{1}{2}} \Big] h[n]= 2 \Big[ \frac{1}{2}\Big]^n \mu [n] + \Big[ \frac{1}{4}\Big]^n \mu [n] h[n]= \Bigg[ 2 \Big[ \frac{1}{2}\Big]^n - \Big[ \frac{1}{4}\Big]^n \Bigg] \mu [n]

siendo la forma de la respuesta un impulso, es un sistema IIR.


Algoritmo en Python

Usando la expresión H(z) se obtiene:

 Hz:
      2     
     z      
------------
 2   3*z   1
z  - --- + -
      4    8

 Hz en fracciones parciales
     z        2*z  
- ------- + -------
  z - 1/4   z - 1/2

 Hz en factores
          2         
         z          
--------------------
(z - 0.5)*(z - 0.25)

 {Q_polos:veces}: {1/2: 1, 1/4: 1}
 {P_ceros:veces}: {0: 2}

estabilidad asintótica en z:
circ1_dentro : 2
circ1_repetidos : 0
circ1_sobre : 0
circ1_fuera : 0
unicos : 2
repetidos : 0
asintota : estable

 h[n]:
/      n        n\             
\- 0.25  + 2*0.5 /*Heaviside(n)
>>> 

2Eva2016TII_T3 graf Hz polos

2Eva2016TII_T3 graf Hz polos02

añadiendo instrucciones para graficar h[n] se obtiene

señal discreta h[n]
n   : [0. 1. 2. 3. 4. 5. 6. 7. 8. 9.]
h[n]: [1.         0.75       0.4375     0.234375   0.12109375
       0.06152344 0.03100586 0.01556396 0.00779724 0.00390244]

2Eva2016TII_T3 graf Hzpolos03

Instrucciones en Python

strong>Nota: cuando se produzca el siguiente error con Numpy para evaluar una expresión con exponente negativo,

Traceback (most recent call last):
  File "D:\MATG1052Ejemplos\Transformadaz\ejercicio....py", line 93, in 
    fi  = f_n(ki)
  File "", line 2, in _lambdifygenerated
    return (9/59)*4**(-n)*Heaviside(n, 1/2)
ValueError: Integers to negative integer powers are not allowed.

proceda actualizando los valores a evaluar como tipo real (dtype float), tan solo usando en la línea de ki con lo siguiente:

ki  = np.arange(0,muestras_fn,1.0)

quedando las instrucciones de la siguiente forma, que si evalua valores para realizar gráficas.

# Transformada z- Fracciones parciales
# http://blog.espol.edu.ec/telg1001/lti-dt-transformada-z-xz-fracciones-parciales-con-python/
import numpy as np
import sympy as sym
import matplotlib.pyplot as plt
import telg1001 as fcnm
#sym.SYMPY_DEBUG=True

# INGRESO
z = sym.Symbol('z')
n = sym.Symbol('n', real=True)

# coeficientes como racional en dominio 'ZZ' enteros
a0 = sym.Rational(3,4)
a1 = sym.Rational(1,8)

Pz = z**2
Qz = z**2-a0*z+a1

#Pz = z*z**2
#Qz = (z-1)*(z**2-(a0)*z+a1)

F = Pz/Qz

# para graficar
f_nombre = 'H'    # nombre de función[z]: H,X,Y, etc
muestras_fn = 10  # muestras para f[n]

# PROCEDIMIENTO
Fz  = fcnm.apart_z(F)
Fz_factor = sym.factor(F.evalf())
Fz_factor = fcnm._round_float_is_int(Fz_factor)

# polos y ceros de Hz
[P,Q] = Fz.as_numer_denom()
P = sym.poly(P,z)
Q = sym.poly(Q,z)
P_ceros = sym.roots(P)
Q_polos = sym.roots(Q)

estable_z = fcnm.estabilidad_asintotica_z(Q_polos)

# Inversa de transformada z
fn = 0*n ; Fz_revisar = [] ; Qz2_term =[]
term_sum = sym.Add.make_args(Fz)
for term_k in term_sum:
    term_kn = fcnm.inverse_z_transform(term_k,z,n)
    if type(term_kn)==tuple:
        fn = fn + term_kn[0]
    elif term_kn is not None:
        fn = fn + term_kn
    elif term_kn is None:
        f_noeval = f_noeval + term_k
    Qz2 = fcnm.Q_cuad_z_parametros(term_k)
    if Qz2:
        Qz2_term.append(Qz2)
fn = fn.collect(sym.Heaviside(n))
fn = fn.collect(sym.DiracDelta(n))
fn = fcnm._round_float_is_int(fn)

# SALIDA
print('\n '+f_nombre+'z:')
sym.pprint(F)
print('\n '+f_nombre+'z en fracciones parciales')
sym.pprint(Fz)
print('\n '+f_nombre+'z en factores')
sym.pprint(Fz_factor)
print('\n {Q_polos:veces}:',Q_polos)
print(' {P_ceros:veces}:',P_ceros)
if len(Qz2_term)>0:
    print('\nparametros cuadraticos: ')
    for i in range(0,len(Qz2_term),1):
        for unterm in Qz2_term[i]:
            print(' termino:',unterm)
            fcnm.print_resultado_dict(Qz2_term[i][unterm])
print('\nestabilidad asintótica en z:')
fcnm.print_resultado_dict(estable_z)
print('\n '+f_nombre.lower()+'[n]:')
sym.pprint(fn)
if len(Fz_revisar)>0:
    print('revisar terminos sin transformada de tabla:')
    for un_term in Fz_revisar:
        print(un_term)

# # GRAFICA  -----------
fig_ROC = fcnm.graficar_Fz_polos(Fz_factor,Q_polos,P_ceros,
                      muestras=101,f_nombre=f_nombre)

fig_Fz = fcnm.graficar_Fs(Fz_factor,Q_polos,P_ceros,
                     muestras=101,
                     f_nombre=f_nombre)

# graficar f[n] -------
f_n = sym.lambdify(n,fn.expand(),modules=fcnm.equivalentes)
ki  = np.arange(0,muestras_fn,1.0)
fi  = f_n(ki)

print('\nseñal discreta '+f_nombre.lower()+'[n]')
print('n   :',ki)
print(f_nombre.lower()+'[n]:',fi)

# graficar f[n]
fig_fn, grafxn = plt.subplots()
plt.axvline(0,color='grey')
plt.stem(ki,fi)
plt.grid()
plt.xlabel('n')
plt.ylabel(f_nombre.lower()+'[n]')
etiqueta = r''+f_nombre.lower()+'[n]= $'+str(sym.latex(fn))+'$'
plt.title(etiqueta)

plt.show()

 


Literal c

revisando los polos y ceros:

ceros : z = 0  y z = 0
polos: z=1/4 y z=1/2

Dado que todos los polos se encuentran dentro del círculo de radio unitario , el sistema es asintóticamente estable, por lo que es BIBO o EASA estable

literal d

S(z) = \Big[ \frac{z}{z-1} \Big] \frac{z^2}{z^2-\frac{3}{4}z + \frac{1}{8}} S(z) = \frac{z^3}{(z-1)(z-1/4)(z-1/2)}

aplicando el mismo método anterior, se tiene:

# coeficientes como racional en dominio 'ZZ' enteros
a0 = sym.Rational(3.4)
a1 = sym.Rational(1,8)

Pz = z*z**2
Qz = (z-1)*(z**2-(a0)*z+a1)
F = Pz/Qz

# para graficar
f_nombre = 'S'    # nombre de función[z]: H,X,Y, etc
muestras_fn = 10  # muestras para f[n]

con resultado:

 Sz:
           3          
          z           
----------------------
        / 2   3*z   1\
(z - 1)*|z  - --- + -|
        \      4    8/

 Sz en fracciones parciales
     z          2*z        8*z   
----------- - ------- + ---------
3*(z - 1/4)   z - 1/2   3*(z - 1)

 Sz en factores
              3             
             z              
----------------------------
(z - 1)*(z - 0.5)*(z - 0.25)

 {Q_polos:veces}: {1: 1, 1/2: 1, 1/4: 1}
 {P_ceros:veces}: {0: 3}

 s[n]:
/    n             \             
|0.25         n   8|             
|----- - 2*0.5  + -|*Heaviside(n)
\  3              3/             
S(z) = \frac{1}{3}\frac{z}{z-1/4} -2\frac{z}{z-1/2} + \frac{8}{3}\frac{z}{z-1}

aplicando la transformada inversa

s(n)= \Bigg[ \frac{1}{3}\Big[ \frac{1}{4} \Big]^n -2\Big[\frac{1}{2} \Big]^n + \frac{8}{3} \Bigg] \mu [n]
señal discreta s[n]
n   : [0 1 2 3 4 5 6 7 8 9]
h[n]: [1.         1.75       2.1875     2.421875   
 2.54296875 2.60449219 2.63549805 2.65106201 2.65885925
 2.66276169]

2Eva2016TII_T3 graf Hz polos03

2Eva2016TII_T3 graf Hz polos05

2Eva2016TII_T3 graf Hzpolos04

 

s2Eva2016TII_T2 LTI CT Circuito RC respuesta de frecuencia H(ω), impulso h(t)

Ejercicio: 2Eva2016TII_T2 LTI CT Circuito RC respuesta de frecuencia H(ω), impulso h(t)

literal a

v_1 (t) = v_R (t) +v_C (t) v_1 (t) = R i(t) +v_2 (t) i(t) = C \frac{\delta v_2 (t)}{\delta t} v_2(t) = RC \frac{\delta v_2 (t)}{\delta t} +v_2 (t) V_1 (\omega) = j \omega RC V_2 (\omega) + V_2(\omega) = V_2(\omega) [1+j\omega RC] H(\omega) = \frac{V_2 (\omega)}{V_1(\omega)} = \frac{1}{1+j\omega RC} = \frac{1}{1+j\frac{\omega}{\omega_c}} \omega _c = \frac{1}{RC} H(\omega) =\frac{V_2(\omega)}{V_1(\omega)} = \frac{1}{1+j\frac{\omega}{2}} \begin{cases} |H(\omega) = \frac{1}{\sqrt{1+\big( \frac{\omega}{2}\big)^2}} \\ \theta_{H(\omega)} = -tg^{-1} \big( \frac{\omega}{2}\big) \end{cases} \omega _c = \frac{1}{RC}

literal c

La respuesta impulso del filtro LPF, se obtiene mediante:

h(t) = \mathcal{F}^{-1} \Big[ H(\omega) \Big] = \mathcal{F}^{-1} \Big[ \frac{1}{1+j(\omega /2)} \Big] =2\mathcal{F}^{-1} \Big[ \frac{1}{j\omega +2} \Big] h(t)=2e^{2t}\mu (t)

literal d

Método 1: usando la respuesta de frecuencia

\begin{cases} |H(\omega) = \frac{1}{\sqrt{1+\big( \frac{\omega}{2}\big)^2}} \\ \theta_{H(\omega)} = -tg^{-1} \big( \frac{\omega}{2}\big) \end{cases} \begin{cases} |H(50) = \frac{1}{\sqrt{1+\big( \frac{50}{2}\big)^2}} = \frac{1}{\sqrt(626)} \\ \theta_{H(\omega)} = -tg^{-1} \big( \frac{50}{2}\big) = -87.7\end{cases} v_2(t) = |H(50)| \sin \big(50t+\theta_{H(50)} \big) v_2(t) = \frac{1}{\sqrt(626)} \sin \big(50t-87.7) \big)

método 2:

V_2(\omega ) = V_1(\omega) H(\omega) V_1 (\omega) = \mathcal{F}[v_1(t) ] = \mathcal{F} [ \sin (50t) ] = j \pi \delta (\omega +50) - j \pi \delta (\omega-50) V_2 (\omega) = V_1(\omega) H(\omega) V_2 (\omega) = \Big[ j \pi \delta (\omega +50) - j \pi \delta (\omega-50) \Big] \Big[ \frac{1}{1+j(\omega/2)} \Big] = j \pi \Big[ \delta (\omega + 50)\frac{1}{1+j(\omega /2)} - \delta (\omega - 50)\frac{1}{1+j(\omega /2)}\Big] = j \pi \Big[ \delta (\omega + 50)\frac{1}{1+j(50 /2)} - \delta (\omega - 50)\frac{1}{1+j(50 /2)}\Big] = j \pi \Big[ \delta (\omega + 50)\frac{1+j25}{626} - \delta (\omega - 50)\frac{1-j25}{626}\Big] = j \frac{\pi}{626} \Big[ \delta (\omega + 50) - \delta (\omega - 50) + j25 \delta (\omega + 50) + j25 \delta (\omega - 50)\Big] = \frac{1}{626} \Big[ j\pi \delta (\omega + 50) - j \pi \delta(\omega - 50)\Big] - \frac{25}{626} \Big[ \pi \delta (\omega + 50) +\pi \delta(\omega - 50)\Big] v_2(t) = \mathcal{F}^{-1} [V_2 (\omega)] = \frac{1}{626} \mathcal{F}^{-1}\Bigg[ \Big[ j\pi \delta (\omega + 50) - j \pi \delta(\omega - 50)\Big] - 25 \Big[ \pi \delta (\omega + 50) +\pi \delta(\omega - 50)\Big] \Bigg] v_2(t) = \frac{1}{626} \Big[ \sin (50t) - 25 \cos (50t) \Big]

usando fasores:

v_2(t) = \frac{1}{626} \cos (50t-177.709) = \frac{1}{626} \sin (50t-87.70)

En la salida, existe un factor de atenuación de 0.04 y un retardo de 87.70°.
Como la señal de entrada se reproduce de manera exacta en su salida a pesar tener amplitud diferente y un retardo en el tiempo.

s2Eva2012TI_T3 LTI DT causal, coeficientes de respuesta impulso h[n]

Ejercicio: 2Eva2012TI_T3 LTI DT causal, coeficientes de respuesta impulso h[n]

La ecuación de diferencias descrita se usa para crear el diagrama de bloques.

y[n] - \frac{5}{4} y[n-1] + \frac{1}{36} y[n-2] + \frac{1}{18} y[n-3] = x[n] - \frac{1}{2} x[n-1]

literal a. Respuesta impulso h[n]

La forma para usar transformadas z requiere que se desplazar tres unidades.

y[n+3] - \frac{5}{4} y[n+2] + \frac{1}{36} y[n+1] + \frac{1}{18} y[n] = x[n+3] - \frac{1}{2} x[n+2] z^3 Y[z] - \frac{5}{4} z^2 Y[z] + \frac{1}{36} z Y[z] + \frac{1}{18} Y[z] = z^3X[z] - \frac{1}{2} z^2 X[z] (z^3 - \frac{5}{4} z^2 + \frac{1}{36} z + \frac{1}{18}) Y[z] = (z^3 - \frac{1}{2} z^2) X[z] H[z] = \frac{Y[z]}{X[z]} = \frac{z^3 - \frac{1}{2} z^2}{z^3 - \frac{5}{4} z^2 + \frac{1}{36} z + \frac{1}{18}}

las raíces del denominador obtenidas con Sympy-Python:

>>> sym.roots(Pz)
{1/2: 1, 0: 2}
>>> sym.factor(Pz)
z**2*(2*z - 1)/2
>>> sym.factor(Qz)
(4*z - 1)*(9*z**2 - 9*z - 2)/36
>>> sym.roots(Qz)
{1/4: 1, 1/2 - sqrt(17)/6: 1, 1/2 + sqrt(17)/6: 1}
>>> 

Observación: existe una raíz con distancia al origen superior al radio=1, por lo que habría un termino creciente no acotado.

2Eva2012TI_T3 graf Hz polos01

que se observa también como:

2Eva2012TI_T3 graf Hz polos02

El modelo se puede plantear como:

H[z] = \frac{z^2\Big(z-\frac{1}{2}\Big)}{\Big(z-\frac{1}{4}\Big)\Big(z^2-z-\frac{2}{9}\Big)}

para usar fracciones parciales modificadas, se multiplica cada lado por 1/z:

\frac{H[z]}{z} = \frac{z\Big(z-\frac{1}{2}\Big)}{\Big(z-\frac{1}{4}\Big)\Big(z^2-z-\frac{2}{9}\Big)}

Se puede plantear un modelo de respuesta por cada raíz como:

\frac{H[z]}{z} = \frac{k1}{z-\frac{1}{4}}+\frac{k_2 z +k_3}{z^2-z-\frac{2}{9}}

usando el método de Heaviside,

k_1 = \frac{z\Big(z-\frac{1}{2}\Big)}{\cancel{\Big(z-\frac{1}{4}\Big)}\Big(z^2-z-\frac{2}{9}\Big)}\Big|_{z=\frac{1}{4}} k_1 = \frac{\frac{1}{4}\Big(\frac{1}{4}-\frac{1}{2}\Big)}{\Big(\frac{1}{4}^2-\frac{1}{4}-\frac{2}{9}\Big)} k_1 = -\frac{1}{16}\frac{1}{\frac{(4)(9)-(16)(9)-(2)(16)(4)}{(16)(4)(9)}} =\frac{36}{236}

con lo que ahora en la expresión se convierte en:

\frac{H[z]}{z} = \frac{z\Big(z-\frac{1}{2}\Big)}{\Big(z-\frac{1}{4}\Big)\Big(z^2-z-\frac{2}{9}\Big)} =\frac{\frac{36}{236}}{z-\frac{1}{4}}+\frac{k_2 z +K_3}{z^2-z-\frac{2}{9}} =

Usando el método de los factores cuadráticos, se multiplica ambos lados por z y z→∞

\frac{H[z]}{z} z = \frac{z^2\Big(z-\frac{1}{2}\Big)}{\Big(z-\frac{1}{4}\Big)\Big(z^2-z-\frac{2}{9}\Big)} =\frac{\frac{36}{236}z}{z-\frac{1}{4}}+z\frac{k_2 z +K_3}{z^2-z-\frac{2}{9}} =

se divide numerador y denominador del lado izquierdo para 1/z3 y el lado derecho el primer termino 1/z y el segundo termino 1/z2

\frac{\Big(1-\frac{1}{2}\frac{1}{z}\Big)}{\Big(1-\frac{1}{4}\frac{1}{z}\Big)\Big(1-\frac{1}{z}-\frac{2}{9}\frac{1}{z^2}\Big)} =\frac{\frac{36}{236}}{1-\frac{1}{4}\frac{1}{z}}+\frac{k_2 +K_3\frac{1}{z}}{1-\frac{1}{z}-\frac{2}{9}\frac{1}{z}}

y cuando z→∞

\frac{\Big(1-\frac{1}{2}(0)\Big)}{\Big(1-\frac{1}{4}(0)\Big)\Big(1-(0)-\frac{2}{9}(0)\Big)} =\frac{\frac{36}{236}}{1-\frac{1}{4}(0)}+\frac{k_2 +K_3(0)}{1-(0)-\frac{2}{9}(0)} 1 =\frac{36}{236}+k_2 k_2 = 1-\frac{36}{236}=\frac{200}{236} = \frac{50}{59}

con lo que K2=50/59 , para encontrar K3 se usa un valor conveniente de z=0

H[0] = \frac{(0)\Big((0)-\frac{1}{2}\Big)}{\Big((0)-\frac{1}{4}\Big)\Big((0)^2-(0)-\frac{2}{9}\Big)} =\frac{\frac{36}{236}}{(0)-\frac{1}{4}}+\frac{\frac{50}{59}(0) +K_3}{(0)^2-(0)-\frac{2}{9}} = 0 =\frac{\frac{36}{236}}{-\frac{1}{4}}+\frac{K_3}{-\frac{2}{9}} = K_3 = \frac{2}{9}\frac{\frac{36}{236}}{-\frac{1}{4}} = -4\frac{2}{9}\frac{36}{236} = -\frac{8}{59} \frac{H[z]}{z} = \frac{36}{236}\frac{1}{z-\frac{1}{4}}+\frac{50/59 z -8/59}{z^2-z-\frac{2}{9}} = \frac{H[z]}{z} = \frac{36}{236}\frac{1}{z-\frac{1}{4}}+\frac{2}{59}\frac{25 z -4}{z^2-z-\frac{2}{9}}

y restaurando en fracciones parciales al multiplicar por z cada lado

H[z] = \frac{9}{59}\frac{z}{z-\frac{1}{4}}+\frac{2}{59}\frac{z(25 z -4)}{z^2-z-\frac{2}{9}}

los valores no se ajustan al modelo planteado en el enunciado, y existe un polo con radio>1, por lo que el sistema es creciente,y otro polo en posición r<0.

h[n] = a \alpha^n \mu [n] + b \beta^n \mu[n] + c \gamma^n \mu [n]

obtenga entonces los valores pertinentes.

a= 9/59 = 0.1525 b= … creciente sin cota, con polo mayor a 1 c=… con polo negativo, lo que es un término con signo alternado.
α= 1/4 = 0.25 β = 1/2 + np.sqrt(17)/6 = 1.1871842709362768 γ = 1/2 – np.sqrt(17)/6 = -0.18718427093627676

Resultados iniciales con el algoritmo

Se encuentra que hay un término que no se puede usar para encontrar la transformada z.

 Hz:
            2      
       3   z       
      z  - --      
           2       
-------------------
        2          
 3   5*z    z    1 
z  - ---- + -- + --
      4     36   18

 Hz en fracciones parciales
50*z*(z - 4/25)       9*z     
--------------- + ------------
   / 2       2\   59*(z - 1/4)
59*|z  - z - -|               
   \         9/               

 Hz en factores
               2                       
              z *(z - 0.5)             
---------------------------------------
           / 2                        \
(z - 0.25)*\z  - z - 0.222222222222222/

 {Q_polos:veces}: {1/4: 1, 1/2 - sqrt(17)/6: 1, 1/2 + sqrt(17)/6: 1}
 {P_ceros:veces}: {1/2: 1, 0: 2}

parametros cuadraticos: 
 termino: 50*z*(z - 4/25)/(59*(z**2 - z - 2/9))
r : 120.72788283210394
gamma : None
beta : None
theta : None

estabilidad asintótica en z:
circ1_dentro : 2
circ1_repetidos : 0
circ1_sobre : 0
circ1_fuera : 1
unicos : 3
repetidos : 0
asintota : inestable

 h[n]:
   -n             
9*4  *Heaviside(n)
------------------
        59        
revisar terminos sin transformada de tabla:
50*z*(z - 4/25)/(59*(z**2 - z - 2/9))

señal discreta h[n]
n   : [0. 1. 2. 3. 4. 5. 6. 7. 8. 9.]
h[n]: [1.52542373e-01 3.81355932e-02 9.53389831e-03 2.38347458e-03
 5.95868644e-04 1.48967161e-04 3.72417903e-05 9.31044756e-06
 2.32761189e-06 5.81902973e-07]

Instrucciones en Python

Desarrollo con las expresiones iniciales, sin corrección sobre los parámetros que se indican en el enunciado.

Nota: cuando se produzca el siguiente error con Numpy para evaluar una expresión con exponente negativo,

Traceback (most recent call last):
  File "D:\MATG1052Ejemplos\Transformadaz\ejercicio03.py", line 93, in 
    fi  = f_n(ki)
  File "", line 2, in _lambdifygenerated
    return (9/59)*4**(-n)*Heaviside(n, 1/2)
ValueError: Integers to negative integer powers are not allowed.

proceda actualizando los valores a evaluar como tipo real (dtype float), tan solo usando en la línea de ki con lo siguiente:

ki  = np.arange(0,muestras_fn,1.0)

quedando las instrucciones de la siguiente forma, que si evalúa valores para realizar gráficas.

# Transformada z- Fracciones parciales
# http://blog.espol.edu.ec/telg1001/lti-dt-transformada-z-xz-fracciones-parciales-con-python/
import numpy as np
import sympy as sym
import matplotlib.pyplot as plt
import telg1001 as fcnm
#sym.SYMPY_DEBUG=True

# INGRESO
z = sym.Symbol('z')
n = sym.Symbol('n', real=True)

# coeficientes como racional en dominio 'ZZ' enteros
a0 = sym.Rational(1,2)
a1 = sym.Rational(5,4)
a2 = sym.Rational(1,36)
a3 = sym.Rational(1,18)

Pz = z**3-a0*z**2
Qz = z**3-a1*z**2+a2*z+a3

F = Pz/Qz

# para graficar
f_nombre = 'H'    # nombre de función[z]: H,X,Y, etc
muestras_fn = 10  # muestras para f[n]

# PROCEDIMIENTO
Fz  = fcnm.apart_z(F)
Fz_factor = sym.factor(F.evalf())
Fz_factor = fcnm._round_float_is_int(Fz_factor)

# polos y ceros de Hz
[P,Q] = F.as_numer_denom()
P = sym.poly(P,z)
Q = sym.poly(Q,z)
P_ceros = sym.roots(P)
Q_polos = sym.roots(Q)

estable_z = fcnm.estabilidad_asintotica_z(Q_polos)

# Inversa de transformada z
fn = 0*n ; f_noeval = 0*n ; Qz2_term =[]
term_sum = sym.Add.make_args(Fz)
for term_k in term_sum:
    term_kn = fcnm.inverse_z_transform(term_k,z,n)
    if type(term_kn)==tuple:
        fn = fn + term_kn[0]
    elif term_kn is not None:
        fn = fn + term_kn
    elif term_kn is None:
        f_noeval = f_noeval + term_k
    Qz2 = fcnm.Q_cuad_z_parametros(term_k)
    if Qz2:
        Qz2_term.append(Qz2)
fn = fn.collect(sym.Heaviside(n))
fn = fn.collect(sym.DiracDelta(n))

# SALIDA
print('\n '+f_nombre+'z:')
sym.pprint(F)
print('\n '+f_nombre+'z en fracciones parciales')
sym.pprint(Fz)
print('\n '+f_nombre+'z en factores')
sym.pprint(Fz_factor)
print('\n {Q_polos:veces}:',Q_polos)
print(' {P_ceros:veces}:',P_ceros)
if len(Qz2_term)>0:
    print('\nparametros cuadraticos: ')
    for i in range(0,len(Qz2_term),1):
        for unterm in Qz2_term[i]:
            print(' termino:',unterm)
            fcnm.print_resultado_dict(Qz2_term[i][unterm])
print('\nestabilidad asintótica en z:')
fcnm.print_resultado_dict(estable_z)
print('\n '+f_nombre.lower()+'[n]:')
sym.pprint(fn)
if not f_noeval==sym.S.Zero:
    print('revisar terminos sin transformada de tabla:')
    print(f_noeval)

# # GRAFICA  -----------
fig_ROC = fcnm.graficar_Fz_polos(Fz_factor,Q_polos,P_ceros,
                      muestras=101,f_nombre=f_nombre)

fig_Fz = fcnm.graficar_Fs(Fz_factor,Q_polos,P_ceros,
                     muestras=101,
                     f_nombre=f_nombre)

# graficar f[n] -------
f_n = sym.lambdify(n,fn.expand(),modules=fcnm.equivalentes)
ki  = np.arange(0,muestras_fn,1.0)
fi  = f_n(ki)

print('\nseñal discreta '+f_nombre.lower()+'[n]')
print('n   :',ki)
print(f_nombre.lower()+'[n]:',fi)

# graficar f[n]
fig_fn, grafxn = plt.subplots()
plt.axvline(0,color='grey')
plt.stem(ki,fi)
plt.grid()
plt.xlabel('n')
plt.ylabel(f_nombre.lower()+'[n]')
etiqueta = r''+f_nombre.lower()+'[n]= $'+str(sym.latex(fn))+'$'
plt.title(etiqueta)

plt.show()

 

s2Eva2011TI_T2 LTI DT Determinar H[z] desde bloques

Ejercicio: 2Eva2011TI_T2 LTI DT Determinar H[z] desde bloques

literal a. expresión de la función de transferencia

El diagrama de bloques del enunciado se reordena de la siguiente forma:

El nuevo diagrama muestra que el sistema tiene dos sub-componentes en paralelo.

H[z] = -\frac{(11/2)z +7}{z^2-1z-2} +\frac{z}{2} -\frac{9}{2}

como existen varios componentes, se pueden tratar de forma separada.

H[z] = H_1[z] +\frac{z}{2} -\frac{9}{2} H_1[z] = -\frac{(11/2)z +7}{z^2-z-2}

los polos y ceros de H1

 {Q_polos:veces}: {2: 1, -1: 1}
 {P_ceros:veces}: {-14/11: 1}
H_1[z] = -\frac{(11/2)z +7}{(z-2)(z+1)}

usando fracciones parciales modificadas

\frac{H_1[z]}{z} = -\frac{(11/2)z +7}{z(z-2)(z+1)} = \frac{k_1}{z}+\frac{k_2}{z-2}+\frac{k_3}{z+1}

usando el algoritmo se tiene:

\frac{H_1[z]}{z} = \frac{\frac{7}{2}}{z}-\frac{3}{z-2}-\frac{\frac{1}{2}}{z+1}

restaurando a fracciones parciales, al multiplicar cada lado por z, se obtiene la función de transferencia H[z] completa como componentes en paralelo junto con H1

H_1[z] = \frac{7}{2}-\frac{3z}{z-2}-\frac{1}{2}\frac{z}{z+1}

usando la tabla de transformadas z

h_1[n] = \frac{7}{2}\delta[n]-3(2)^n \mu[n] -\frac{1}{2}(-1)^n \mu[n]

Usando el algoritmo de la sección X[z] Fracciones parciales modificadas con Python para H1[z]

con entrada:

# coeficientes como racional en dominio 'ZZ' enteros
a0 = sym.Rational(11,2)

Pz = -(a0*z+7)
Qz = z**2-z-2

se obtiene:

 Hz:
  11*z    
- ---- - 7
   2      
----------
 2        
z  - z - 2

 Hz en fracciones parciales
      z        3*z    7
- --------- - ----- + -
  2*(z + 1)   z - 2   2

 Hz en factores
-3.5*(0.785714285714286*z + 1) 
-------------------------------
      (0.5*z - 1)*(z + 1)      

 {Q_polos:veces}: {2: 1, -1: 1}
 {P_ceros:veces}: {-14/11: 1}

estabilidad asintótica en z:
circ1_dentro : 0
circ1_repetidos : 0
circ1_sobre : 1
circ1_fuera : 1
unicos : 2
repetidos : 0
asintota : inestable

 h[n]:
/      n       \                               
|  (-1)       n|                7*DiracDelta(n)
|- ----- - 3*2 |*Heaviside(n) + ---------------
\    2         /                       2       

señal discreta h[n]
n   : [0. 1. 2. 3. 4. 5. 6. 7. 8. 9.]
h[n]: [    0.     -5.5   -12.5   -23.5   -48.5
         -95.5  -192.5  -383.5  -768.5 -1535.5]

resultado que se completan con los términos de los otros componentes, al incorporar la expresión con los elementos en paralelo como se desarrolla en el siguiente literal.

2Eva2011TI_T2 graf Hz polos01

2Eva2011TI_T2 graf Hz polos02

Instrucciones en Python

# Transformada z- Fracciones parciales
# http://blog.espol.edu.ec/telg1001/lti-dt-transformada-z-xz-fracciones-parciales-con-python/
import numpy as np
import sympy as sym
import matplotlib.pyplot as plt
import telg1001 as fcnm
#sym.SYMPY_DEBUG=True

# INGRESO
z = sym.Symbol('z')
n = sym.Symbol('n', real=True)

# coeficientes como racional en dominio 'ZZ' enteros
a0 = sym.Rational(11,2).limit_denominator(100)

Pz = -(a0*z+7)
Qz = z**2-z-2
F = Pz/Qz

# para graficar
f_nombre = 'H'    # nombre de función[z]: H,X,Y, etc
muestras_fn = 10  # muestras para f[n]

# PROCEDIMIENTO
Fz  = fcnm.apart_z(F)
Fz_factor = sym.factor(F.evalf())

# polos y ceros de Hz
[P,Q] = F.as_numer_denom()
P = sym.poly(P,z)
Q = sym.poly(Q,z)
P_ceros = sym.roots(P)
Q_polos = sym.roots(Q)

estable_z = fcnm.estabilidad_asintotica_z(Q_polos)

# Inversa de transformada z
fn = 0*n ; Fz_revisar = [] ; Qz2_term =[]
term_sum = sym.Add.make_args(Fz)
for term_k in term_sum:
    term_kn = fcnm.inverse_z_transform(term_k,z,n)
    if type(term_kn)==tuple:
        fn = fn + term_kn[0]
    elif term_kn is not None:
        fn = fn + term_kn
    elif term_kn is None:
        f_noeval = f_noeval + term_k
    Qz2 = fcnm.Q_cuad_z_parametros(term_k)
    if Qz2:
        Qz2_term.append(Qz2)
fn = fn.collect(sym.Heaviside(n))
fn = fn.collect(sym.DiracDelta(n))
#fn = fcnm._round_float_is_int(fn)

# SALIDA
print('\n '+f_nombre+'z:')
sym.pprint(F)
print('\n '+f_nombre+'z en fracciones parciales')
sym.pprint(Fz)
print('\n '+f_nombre+'z en factores')
sym.pprint(Fz_factor)
print('\n {Q_polos:veces}:',Q_polos)
print(' {P_ceros:veces}:',P_ceros)
if len(Qz2_term)>0:
    print('\nparametros cuadraticos: ')
    for i in range(0,len(Qz2_term),1):
        for unterm in Qz2_term[i]:
            print(' termino:',unterm)
            fcnm.print_resultado_dict(Qz2_term[i][unterm])
print('\nestabilidad asintótica en z:')
fcnm.print_resultado_dict(estable_z)
print('\n '+f_nombre.lower()+'[n]:')
sym.pprint(fn)
if len(Fz_revisar)>0:
    print('revisar terminos sin transformada de tabla:')
    for un_term in Fz_revisar:
        print(un_term)

# # GRAFICA  -----------
fig_ROC = fcnm.graficar_Fz_polos(Fz_factor,Q_polos,P_ceros,
                      muestras=101,f_nombre=f_nombre)

fig_Fz = fcnm.graficar_Fs(Fz_factor,Q_polos,P_ceros,
                     muestras=101,
                     f_nombre=f_nombre)

# graficar f[n] -------
f_n = sym.lambdify(n,fn.expand(),modules=fcnm.equivalentes)
ki  = np.arange(0,muestras_fn,1.0)
fi  = f_n(ki)

print('\nseñal discreta '+f_nombre.lower()+'[n]')
print('n   :',ki)
print(f_nombre.lower()+'[n]:',fi)

# graficar f[n]
fig_fn, grafxn = plt.subplots()
plt.axvline(0,color='grey')
plt.stem(ki,fi)
plt.grid()
plt.xlabel('n')
plt.ylabel(f_nombre.lower()+'[n]')
etiqueta = r''+f_nombre.lower()+'[n]= $'+str(sym.latex(fn))+'$'
plt.title(etiqueta)

plt.show()

literal b. Ecuación de diferencias

Se realiza la conversión por la suma de cada componente (en paralelo):

Para y1[n]:

Y_1[z] = H_1[z]X[z] = \Big[-\frac{(11/2)z +7}{z^2-z-2}\Big] X[z] Y_1[z] [z^2-z-2]= -(11/2)zX[z] -7X[z] y_![n+2] - y_1[n+1] -2y_1[n] = -\frac{11}{2} x[n+1] -7x[n]

Para y2[n]:

Y_2[z] = H_2[z]X[z] = \frac{z}{2}X[z] y_2[n] = \frac{1}{2}x[n+1]

Para y3[n]:

Y_3[z] = H_3[z]X[z] = -\frac{9}{2} X[z] y_3[n] = -\frac{9}{2} x[n]

se suman las expresiones obtenidas de Y1[z] + Y2[z]+Y3[z]

y[n+2] - y[n+1] -\cancel{2y[n]} +\cancel{y[n]} + \cancel{y[n]}= -\frac{11}{2} x[n+1] -7x[n] +\frac{1}{2}x[n+1]-\frac{9}{2} x[n]

la ecuación de diferencias simplificada es;

y[n+2] - y[n+1] = -5 x[n+1] -\frac{23}{2} x[n]

El sistema global se puede reescribir nuevamente en z como

z^2Y[z] - zY[z] = -5 zX[z] -\frac{23}{2} X[z] z(z-1)Y[z] = (-5 z -\frac{23}{2}) X[z]

se tiene que el sistema tienen polos en 0 y 1, que se encuentran dentro del radio 1 del plano z.

literal d. Respuesta al impulso H[z]

H[z] se obtiene a partir de la última ecuación

H[z] = \frac{Y[z]}{H[z]} = -\frac{5 z +\frac{23}{2}}{z(z-1)}

aplicando fracciones parciales modificadas:

\frac{H[z]}{z} = -\frac{5 z +\frac{23}{2}}{z^2(z-1)} \frac{H[z]}{z} = \frac{33}{2} \frac{1}{z} +\frac{23}{2}\frac{1}{z^2} -\frac{33}{2} \frac{1}{z-1}

restaurando fracciones parciales al multiplicar por z

H[z] = \frac{33}{2} +\frac{23}{2} \frac{1}{z} -\frac{33}{2} \frac{z}{z-1}

usando la tabla de transformada z se tiene:

h[n] = \frac{33}{2}\delta [n] +\frac{23}{2} \delta [n-1] - \frac{33}{2} \mu[n]

usando el algoritmo con entrada:

# coeficientes como racional en dominio 'ZZ' enteros
a0 = sym.Rational(11,2)
a1 = sym.Rational(9,2)
a2 = sym.Rational(23,2)
#F = -(a0*z+7)/(z**2-z-2) + z/2 -a1
F = -(5*z+a2)/(z*(z-1))

se tiene como resultado:

 Hz:
-5*z - 23/2
-----------
 z*(z - 1) 

 Hz en fracciones parciales
  16.5*z          11.5
- ------ + 16.5 + ----
  z - 1            z  

 Hz en factores
-11.5*(0.434782608695652*z + 1) 
--------------------------------
           z*(z - 1)            

 {Q_polos:veces}: {1: 1, 0: 1}
 {P_ceros:veces}: {-2.30000000000000: 1}

estabilidad asintótica en z:
circ1_dentro : 1   circ1_repetidos : 0 
circ1_sobre : 1  circ1_fuera : 0
unicos : 2  repetidos : 0
asintota : marginalmente estable

 h[n]:
16.5*DiracDelta(n) + 11.5*DiracDelta(n - 1) - 16.5*Heaviside(n)

señal discreta h[n]
n   : [0. 1. 2. 3. 4. 5. 6.]
h[n]: [  0.   -5.  -16.5 -16.5 -16.5 -16.5 -16.5]

2Eva2011TI_T2 graf Hz polos03

2Eva2011TI_T2 graf Hz polos04

2Eva2011TI_T2 graf Hz polos05

s2Eva2010TI_T3 LTI CT simplificar H(s) por Bloques

Ejercicio: 2Eva2010TI_T3 LTI CT simplificar H(s) por Bloques

Se ubican algunos puntos de referencia sobre el diagrama de bloques para plantear las ecuaciones, dejando para el último el bloque en serie del exponencial e(-3s) o retraso en tiempo

W_1(s) = X(s) W_2(s) = W_1(s)-\frac{6}{s+5} W_4(s) W_2(s) = X(s)- -\frac{6}{s+5} W_4(s) W_3(s) = \frac{1}{s} W_2(s) = \frac{1}{s} \Big[ X(s) -\frac{6}{s+5} W_4(s) \Big] W_3(s) = \frac{1}{s} X(s) -\frac{6}{s(s+5)} W_4(s) W_4(s) = W_3(s) +\frac{1}{s+4}X(s)

Se plantea la simplificación de H(s) =W4(s)/X(s), dejando el término exponencial o retraso de tiempo para el final. Por lo que se despeja W4(s)

W_4(s) = \frac{1}{s} X(s) -\frac{6}{s(s+5)} W_4(s) +\frac{1}{s+4}X(s) W_4(s) +\frac{6}{s(s+5)} W_4(s) = \frac{1}{s} X(s) +\frac{1}{s+4}X(s) \Big[1 +\frac{6}{s(s+5)}\Big] W_4(s) = \Big[ \frac{1}{s} +\frac{1}{s+4} \Big] X(s) \frac{s(s+5)+6}{s(s+5)} W_4(s) = \frac{(s+4)+s}{s(s+4)} X(s) \frac{s^2+5s+6}{s+5} W_4(s) = \frac{2(s+2)}{s+4} X(s) \frac{(s+2)(s+3)}{s+5} W_4(s) = \frac{2(s+2)}{s+4} X(s) \frac{s+3}{s+5} W_4(s) = \frac{2}{s+4} X(s) \frac{W_4(s)}{X(s)} = \frac{2(s+5)}{(s+4)(s+3)}

La respuesta al impulso H1(s) tiene polos en s=-3 y s=-4, que se encuentran en el lado izquierdo del plano imaginario s. Por lo que sus componentes en el dominio del tiempo son exponenciales decrecientes, el sistema es asintóticamente estable.

El grado del polinomio P del numerador es menor al grado del polinomio Q del numerador.

Separando en fracciones parciales para H1(s) :

H_1(s) =\frac{2(s+5)}{(s+4)(s+3)} H_1(s) = \frac{k_1}{s+4} +\frac{k_2}{s+3} = k_1 =\frac{2(s+5)}{\cancel{(s+4)}(s+3)}\Bigg|_{s=-4} = \frac{2(-4+5)}{(-4+3)} = \frac{2(1)}{-1} = -2 k_2 =\frac{2(s+5)}{(s+4)\cancel{(s+3)}} \Bigg|_{s=-3} = \frac{2(-3+5)}{(-3+4)} = \frac{2(2)}{1} = 4 H_1(s) = -\frac{2}{s+4} +\frac{4}{s+3} =

Finalmente, añadiendo el término exponencial, que es el retraso en tiempo del sistema:

H(s) = \Big[-\frac{2}{s+4} +\frac{4}{s+3} \Big] e^{-3s}

s2Eva2010TI_T3 polos

se observa el comportamiento de H(s) junto a los polos en la parte real e imaginaria:

s2Eva2010TI_T3 polos H(s)

el resultado con el algoritmo para el literal a es:

 H(s) = P(s)/Q(s):
/    2       4  \  -3*s
|- ----- + -----|*e    
\  s + 4   s + 3/      
 H(s) en factores:
           -3*s
2*(s + 5)*e    
---------------
(s + 3)*(s + 4)

 h(t) :
/   9  -3*t      12  -4*t\                 
\4*e *e     - 2*e  *e    /*Heaviside(t - 3)

polosceros:
exp(-3*s) : {'Q_polos': {-3: 1, -4: 1}, 'P_ceros': {-5: 1}, 'Hs_k': 2*(s + 5)/((s + 3)*(s + 4))}
Q_polos : {-3: 1, -4: 1}
P_ceros : {-5: 1}

Estabilidad de H(s):
 n_polos_real : 2
 n_polos_imag : 0
 enRHP : 0
 unicos : 0
 repetidos : 0
 asintota : estable

literal b. mostrar h(t)

H(s) = \Big[-\frac{2}{s+4} +\frac{4}{s+3} \Big] e^{-3s}

Usando la tabla de transformadas de Laplace, para la línea 5 y la propiedad de desplazamiento en t, se obtiene:

h(t) = \Big[-2 e^{-4(t-3)} +4e^{-3(t-3)} \Big] \mu(t-3)

literal c. Y(s) con entrada exponencial decreciente

La señal de salida Y(s)=H(s)X(s) ante una entrada X(s)=1/(s+5) que es la transformada de Laplace de x(t)=e-5t μ(t)

Y(s) = H(s)x(s) = \Big[\frac{2(s+5)}{(s+4)(s+3)}\Big]\Big[\frac{1}{s+5}\Big] e^{-3s} Y(s) = \frac{2}{(s+4)(s+3)}e^{-3s}

Tarea: Realizar el desarrollo analítico y revisar los resultados del algoritmo

al realizar las fracciones parciales, se obtiene:

Y(s) = \Big[-\frac{2}{s+4} +\frac{2}{s+3} \Big] e^{-3s}

Con transformada inversa de Laplace, usando nuevamete la tabla de transformadas de Laplace, para la línea 5 y la propiedad de desplazamiento en t, se obtiene:

y(t) = \Big(-2 e^{-4(t-3)} +2e^{-3(t-3)}\Big) \mu (t-3)

s2Eva2010TI_T3 y(t)

El resultado del algoritmo para el literal c del ejercicio es:

 X(s): 
  1  
-----
s + 5

Respuesta entrada cero ZIR H(s) y condiciones iniciales
term_cero : 0
ZIR :
0
yt_ZIR :
0

 ZSR respuesta estado cero:
ZSR :
/    2       2  \  -3*s
|- ----- + -----|*e    
\  s + 4   s + 3/      
yt_ZSR :
/   9  -3*t      12  -4*t\                 
\2*e *e     - 2*e  *e    /*Heaviside(t - 3)

 Y(s)_total = ZIR + ZSR:
/    2       2  \  -3*s
|- ----- + -----|*e    
\  s + 4   s + 3/      

 y(t)_total = ZIR + ZSR:
/   9  -3*t      12  -4*t\                 
\2*e *e     - 2*e  *e    /*Heaviside(t - 3)
>>>

Instrucciones en Python

Usando los bloques desarrollados en la Unidad 4 Sistemas LTI – Laplace  y las funciones resumidas como telg1001.py que pueden ser usados en cada pregunta.

# Y(s) Respuesta total con entada cero y estado cero
# Qs Y(s) = Ps X(s) ; H(s)=Ps/Qs
# http://blog.espol.edu.ec/telg1001/
import sympy as sym
import matplotlib.pyplot as plt
import telg1001 as fcnm

# INGRESO
s = sym.Symbol('s')
t = sym.Symbol('t', real=True)
d = sym.DiracDelta(t)
u = sym.Heaviside(t)

# Literal a. H(s) y estabilidad
Hs = 2*(s+5)/((s+4)*(s+3)) * sym.exp(-3*s)
#Hs = 1+0*s cuando es constante

# literal c. X(s) Señal de entrada
xt = sym.exp(-5*t)*u

# condiciones iniciales, [y'(0),y(0)] orden descendente
t0 = 0
cond_inicio = [0, 0] # estado cero no se usan

# Grafica, intervalo tiempo [t_a,t_b]
t_a = 0 ; t_b = 5
muestras = 101  # 51 resolucion grafica

# PROCEDIMIENTO
Hs = fcnm.apart_s(Hs) # fracciones parciales
Hs_fc = fcnm.factor_exp(Hs) # en factores
Hs_Qs2 = fcnm.Q_cuad_s_parametros(Hs_fc)

polosceros = fcnm.busca_polosceros(Hs)
Q_polos = polosceros['Q_polos']
P_ceros = polosceros['P_ceros']

estable = fcnm.estabilidad_asintotica_s(Q_polos)

# H(t) respuesta al impulso
ht = 0*s
term_suma = sym.Add.make_args(Hs)
for term_k in term_suma:
    ht_k = sym.inverse_laplace_transform(term_k,s,t)
    # simplifica log(exp()) ej: e**(-2s)/(s**2)
    if ht_k.has(sym.log):
        ht_k = sym.simplify(ht_k,inverse=True)
    ht  = ht + ht_k
lista_escalon = ht.atoms(sym.Heaviside)
ht = sym.expand(ht,t) # terminos suma
ht = sym.collect(ht,lista_escalon)

# PROCEDIMIENTO Respuesta ZIR, ZSR
Xs = fcnm.laplace_transform_suma(xt)

# ZIR_s respuesta entrada cero de s
sol_ZIR = fcnm.respuesta_ZIR_s(Hs,cond_inicio)
ZIR = sol_ZIR['ZIR']
yt_ZIR = sol_ZIR['yt_ZIR']

# ZSR respuesta estado cero, Y(s) a entrada X(s)
sol_ZSR = fcnm.respuesta_ZSR_s(Hs,Xs)
ZSR = sol_ZSR['ZSR']
yt_ZSR = sol_ZSR['yt_ZSR']

# Respuesta total Y(s) y y(t)
Ys = ZIR + ZSR
Ys = fcnm.apart_s(Ys)
yt = yt_ZIR + yt_ZSR
lista_escalon = yt.atoms(sym.Heaviside)
yt = sym.collect(yt,lista_escalon)

# SALIDA
print(' H(s) = P(s)/Q(s):')
sym.pprint(Hs)
print(' H(s) en factores:')
sym.pprint(Hs_fc)
if len(Hs_Qs2)>0:
    print('\nH(s) parámetros cuadraticos:')
    fcnm.print_resultado_dict(Hs_Qs2)

print('\n h(t) :')
sym.pprint(ht)

print('\npolosceros:')
fcnm.print_resultado_dict(polosceros)

print('\nEstabilidad de H(s):')
for k in estable:
    print('',k,':',estable[k])

print('\n X(s): ')
sym.pprint(Xs)
print('\nRespuesta entrada cero ZIR H(s) y condiciones iniciales')

if not(sol_ZIR == sym.nan): # existe resultado
    fcnm.print_resultado_dict(sol_ZIR)
else:
    print(' insuficientes condiciones iniciales')
    print(' revisar los valores de cond_inicio[]')

print('\n ZSR respuesta estado cero:')
fcnm.print_resultado_dict(sol_ZSR)

print('\n Y(s)_total = ZIR + ZSR:')
sym.pprint(Ys)
print('\n y(t)_total = ZIR + ZSR:')
sym.pprint(yt)

# Graficas polos, H(s), con polos h(t) --------
muestras_H = 201
figura_s  = fcnm.graficar_Fs(Hs_fc,Q_polos,P_ceros,f_nombre='H',solopolos=True)
figura_Hs = fcnm.graficar_Fs(Hs_fc,Q_polos,P_ceros,muestras=muestras_H,f_nombre='H')
figura_ht = fcnm.graficar_ft(ht,t_a,t_b,muestras,f_nombre='h')
# GRAFICAS y(t),x(t),h(t) ---------------------
figura_ft = fcnm.graficar_xh_y(xt,ht,yt,t_a,t_b,muestras)
plt.show()

s2Eva2010TI_T1 LTI DT Bloques de H[z] para ecuación de diferencias

Ejercicio: 2Eva2010TI_T1 LTI DT Bloques de H[z] para ecuación de diferencias

2E2010TI_T1 LTID diagrama1

litera a. La función de transferencia H(z)

El diagrama mostrado se reordena y añade anotaciones para facilitar la escritura de la ecuación de diferencias:

Dado que hay dos bloques (1/z) o de retraso hacia abajo, se identifica que los polinomios son de segundo orden (N=2).

Los coeficientes de retraso se leen de arriba hacia abajo como a0,-a1,-a2 ordenados en la siguiente tabla:

coeficientes
retraso adelanto
a0 = 1 b0 = 0
a1 = -(1.6) = -1.6 b1 = 4
a2 = -(-0.63) = 0.63 b2 = -4

Se aplica lo mismo para los coeficientes de adelanto y se completa el modelo para los polinomios del numerador y denominador:

H[z] = \frac{b_0 z^N +b_1 z^{N-1} + \text{ ... } + b_{N-1} z + b_N}{z^N + a_1 z^{N-1} +\text{ ... } + a_{N-1}z +a_N} H[z] = \frac{4z-4}{z^2 -1.6 z + 0.63}

Con el polinomio del denominador Qz  se encuentran las raíces

polos: {0.70: 1, 0.90: 1}

se aplican fracciones parciales modificadas al dividir cada lado para z,

\frac{H[z]}{z} = \frac{4}{z}\frac{z-1}{(z-0.7)(z-0.9)} = \frac{k_1}{z-0.7} +\frac{k_2}{z-0.9}+\frac{k_3}{z}

aplicando el método de Heaviside:

k_1 = \frac{4}{z}\frac{z-1}{\cancel{(z-0.7)}(z-0.9)}\Big|_{z=0.7} = \frac{4}{0.7}\frac{(0.7)-1}{((0.7)-0.9)} = 8.5714 k_2 = \frac{4}{z}\frac{z-1}{(z-0.7)\cancel{(z-0.9)}}\Big|_{z=0.9} = \frac{4}{0.9}\frac{(0.9)-1}{((0.9)-0.7)} = -2.2222 k_3 = \frac{4}{\cancel{z}}\frac{z-1}{(z-0.7)(z-0.9)}\Big|_{z=0} = 4\frac{(0)-1}{((0)-0.7)((0)-0.9)} = -6.3492 \frac{H[z]}{z} = \frac{8.5714}{z-0.7} -\frac{2.2222}{z-0.9}-\frac{6.3492}{z}

restaurando a fracciones parciales al multiplicar por z

H[z] = \frac{8.5714z}{z-0.7} -\frac{2.2222z}{z-0.9}-6.3492

usando la tabla de transformadas z

h[n] = 8.5714 (0.7)^n \mu[n] -2.2222(0.9)^n \mu[n] -6.3492 \delta[n]

Algoritmo en Python

Usando el algoritmo de la sección X[z] Fracciones parciales modificadas con Python

con entrada:

Pz = 4*z-4
Qz = z**2-1.6*z+0.63

y resultado:

 Hz:
     4*z - 4     
-----------------
 2               
z  - 1.6*z + 0.63

 Hz en fracciones parciales
8.57142857142857*z   2.22222222222222*z                   
------------------ - ------------------ - 6.34920634920635
     z - 0.7              z - 0.9                         

 Hz en factores
     4*(z - 1)     
-------------------
(z - 0.9)*(z - 0.7)

 {Q_polos:veces}: {0.700000000000000: 1, 0.900000000000000: 1}
 {P_ceros:veces}: {1: 1}

estabilidad asintótica en z:
circ1_dentro : 2
circ1_repetidos : 0
circ1_sobre : 0
circ1_fuera : 0
unicos : 2
repetidos : 0
asintota : estable

 h[n]:
/                    n                       n\                               
\8.57142857142857*0.7  - 2.22222222222222*0.9 /*Heaviside(n) 

  - 6.34920634920635*DiracDelta(n)

señal discreta h[n]
n   : [0 1 2 3 4 5 6 7 8 9]
h[n]: [ 0.   4.   2.4   1.32   0.6   0.1284
       -0.17256    -0.356988   -0.462468   -0.51504636]

2Eva2010TI_T1 Hz polos01

2Eva2010TI_T1 h[n]

Observaciones:

Las raíces características o frecuencias naturales del sistema se encuentran dentro del círculo de radio unitario. El sistema es asintóticamente estable, que implica que es BIBO estable.

literal c. Ecuación de diferencias

H[z] = \frac{4z-4}{z^2 -1.6 z + 0.63} Y[z](z^2 -1.6 z+ 0.63) = (4z-4)X[z] y[n+2] - 1.6 y[n+1] + 0.63y[n] = 4x[n+1]-4x[n]

literal d. entrada x[n]

x(t) = cos(1500t) = cos (2πf) con intervalo de muestreo Ts = 0.0015 = 1/fs.

De la señal de entrada se tiene que 2πf=1500, por lo que f = 1500/(2π)

Para un buen muestreo de la señal de entrada, se debe cumplir la frecuencia de Nyquist, es decir fs>2f

fs = \frac{1}{0.0015} = 666.66 \text{ ; } f= \frac{1500}{2 \pi}=238.73 fs > 2*f = 2(238.73) = 477.46

que cumple con la frecuencia mínima de muestreo y no se presentan inconvenientes para el funcionamiento de la señal por sub-muestreo

Instrucciones en Python

# Transformada z- Fracciones parciales
# http://blog.espol.edu.ec/telg1001/lti-dt-transformada-z-xz-fracciones-parciales-con-python/
import numpy as np
import sympy as sym
import matplotlib.pyplot as plt
import telg1001 as fcnm
#sym.SYMPY_DEBUG=True

# INGRESO
z = sym.Symbol('z')
n = sym.Symbol('n', real=True)

# coeficientes como racional en dominio 'ZZ' enteros
Pz = 4*z-4
Qz = z**2-1.6*z+0.63

F = Pz/Qz

# para graficar
f_nombre = 'H'    # nombre de función[z]: H,X,Y, etc
muestras_fn = 10  # muestras para f[n]

# PROCEDIMIENTO
Fz  = fcnm.apart_z(F)
Fz_factor = sym.factor(F.evalf())
Fz_factor = fcnm._round_float_is_int(Fz_factor)

# polos y ceros de Hz
[P,Q] = F.as_numer_denom()
P = sym.poly(P,z)
Q = sym.poly(Q,z)
P_ceros = sym.roots(P)
Q_polos = sym.roots(Q)
estable_z = fcnm.estabilidad_asintotica_z(Q_polos)

# Inversa de transformada z
fn = 0*n ; Fz_revisar = [] ; Qz2_term =[]
term_sum = sym.Add.make_args(Fz)
for term_k in term_sum:
    term_kn = fcnm.inverse_z_transform(term_k,z,n)
    if type(term_kn)==tuple:
        fn = fn + term_kn[0]
    elif term_kn is not None:
        fn = fn + term_kn
    elif term_kn is None:
        f_noeval = f_noeval + term_k
    Qz2 = fcnm.Q_cuad_z_parametros(term_k)
    if Qz2:
        Qz2_term.append(Qz2)
fn = fn.collect(sym.Heaviside(n))
fn = fn.collect(sym.DiracDelta(n))
fn = fcnm._round_float_is_int(fn)

# SALIDA
print('\n '+f_nombre+'z:')
sym.pprint(F)
print('\n '+f_nombre+'z en fracciones parciales')
sym.pprint(Fz)
print('\n '+f_nombre+'z en factores')
sym.pprint(Fz_factor)
print('\n {Q_polos:veces}:',Q_polos)
print(' {P_ceros:veces}:',P_ceros)
if len(Qz2_term)>0:
    print('\nparametros cuadraticos: ')
    for i in range(0,len(Qz2_term),1):
        for unterm in Qz2_term[i]:
            print(' termino:',unterm)
            fcnm.print_resultado_dict(Qz2_term[i][unterm])
print('\nestabilidad asintótica en z:')
fcnm.print_resultado_dict(estable_z)
print('\n '+f_nombre.lower()+'[n]:')
sym.pprint(fn)
if len(Fz_revisar)>0:
    print('revisar terminos sin transformada de tabla:')
    for un_term in Fz_revisar:
        print(un_term)

# # GRAFICA  -----------
fig_ROC = fcnm.graficar_Fz_polos(Fz_factor,Q_polos,P_ceros,
                      muestras=101,f_nombre=f_nombre)

fig_Fz = fcnm.graficar_Fs(Fz_factor,Q_polos,P_ceros,
                     muestras=101,
                     f_nombre=f_nombre)

# graficar f[n] -------
f_n = sym.lambdify(n,fn.expand(),modules=fcnm.equivalentes)
ki  = np.arange(0,muestras_fn,1.0)
fi  = f_n(ki)

print('\nseñal discreta '+f_nombre.lower()+'[n]')
print('n   :',ki)
print(f_nombre.lower()+'[n]:',fi)

# graficar f[n]
fig_fn, grafxn = plt.subplots()
plt.axvline(0,color='grey')
plt.stem(ki,fi)
plt.grid()
plt.xlabel('n')
plt.ylabel(f_nombre.lower()+'[n]')
etiqueta = r''+f_nombre.lower()+'[n]= $'+str(sym.latex(fn))+'$'
plt.title(etiqueta)

plt.show()

 

s2Eva2009TII_T2 LTI DT Dado h[n], y[n] determine x[n]

Ejercicio: 2Eva2009TII_T2 LTI DT Dado h[n], y[n] determine X[n]

literal a. Expresión de x[n]

Con los datos de h[n] y y[n], se obtienen las transformadas H[z] y Y[z],

h[n] = 2 \Big( \frac{1}{3} \Big)^n \mu [n-1] = 2 \Big( \frac{1}{3} \Big)^1\Big( \frac{1}{3} \Big)^{-1}\Big( \frac{1}{3} \Big)^n \mu [n-1] h[n] = \frac{2}{3} \Big( \frac{1}{3} \Big)^{n-1} \mu [n-1]

la transformada z es:

H[z] = \frac{2}{3} \frac{1}{z-\frac{1}{3}}

Para la señal de salida Yc[z] conocida:

y_c[n] = (-2)^{n} \mu [n-1] = (-2)^1(-2)^{n-1} \mu [n-1] y_c[n] = -2 (-2)^{n-1} \mu [n-1]

la transformada z es:

Y_c[z] = -2\frac{1}{z+2}

La señal de entrada x[n]

x[n] = a\delta [n] + b (c)^{n-1} \mu [n-1] X[z] = a + b \frac{1}{z-c}

La señal de salida  Y[z] esperada en dominio z se obtiene como Y[z] = H[z]X[z], la expresión se escribe como:

y_e [z] = \Bigg[\frac{2}{3} \frac{1}{z-\frac{1}{3}}\Bigg] \Bigg[a + b \frac{1}{z-c}\Bigg] = \frac{2}{3} \Bigg[\frac{a}{z-\frac{1}{3}} + \frac{b}{(z-c)(z-\frac{1}{3})}\Bigg] = \frac{2}{3} \Bigg[\frac{a(z-c) + b}{(z-c)(z-\frac{1}{3})}\Bigg] y_e [z] = \frac{2}{3} \Bigg[\frac{a(z-c) + b}{(z-c)(z-\frac{1}{3})}\Bigg] = \frac{2}{3} \Bigg[\frac{az - ac + b}{(z-c)(z-\frac{1}{3})}\Bigg]

igualando las expresiones con Y[z] conocida con Y[z] esperada:

Y_c[z] = Y_e[z] -2\frac{1}{z+2} = \frac{2}{3} \Bigg[ \frac{az - ac + b}{(z-c)(z-\frac{1}{3})} \Bigg] -3\frac{1}{z+2} = \frac{az - ac + b}{(z-c)(z-\frac{1}{3})}

para que el denominador quede (z+2), se iguala el término con (z-c), con lo que c=-2

-3\frac{1}{z+2} = \frac{az-a(-2) + b}{(z+2)(z-\frac{1}{3})} -3 = \frac{az+2a + b}{(z-\frac{1}{3})} -3 (z-\frac{1}{3}) = az +2a + b -3 z + 1 = az + 2a + b

comparando el término z, se tiene que a=-3, quedando solo la parte constante para determinar el valor de b

2a + b = 1 2(-3) + b = 1

se tiene que b= 7

teniendo la expresión de la entrada como:

x[n] = -3\delta [n] + 7 (-2)^{n-1} \mu [n-1]

literal b. Ecuacion de diferencias H[z]

Se inicia con los datos de H[z]

H[z] = \frac{X[z]}{Y[z]} = \frac{2}{3} \frac{1}{z-\frac{1}{3}} \Big(z-\frac{1}{3} \Big) Y[z] = \frac{2}{3}X[z] zY[z]-\frac{1}{3} Y[z] = \frac{2}{3}X[z] y[n+1]-\frac{1}{3} y[n] = \frac{2}{3}x[n]

literal c. Diagrama de bloques

Para el diagrama de bloques se desplaza para despejar y[n]

y[n]-\frac{1}{3} y[n-1] = \frac{2}{3}x[n-1] y[n]= \frac{1}{3} y[n-1] + \frac{2}{3}x[n-1]

H[z] = \frac{X[z]}{Y[z]} = \frac{2}{3} \frac{1}{z-\frac{1}{3}}

Observaciones:

Las raíces características o frecuencias naturales del sistema H[z]  se encuentran dentro del círculo de radio unitario. El sistema es asintóticamente estable, que implica que es BIBO estable.

 

 

 

s1Eva2016TII_T4 LTI DT Sistema 2do orden

Ejercicio: 1Eva2016TII_T4 LTI DT Sistema 2do orden

literal a. Ecuación de diferencias de coeficientes constantes

Al diagrama del ejercicio se añaden las referencias de y[n], y[n-1] y y[n-2] de acuerdo a los bloques de retraso (delay) para escribir la expresión:

y[n] = x[n] + \frac{3}{4} y[n-1] - \frac{1}{8} y[n-2] y[n] - \frac{3}{4} y[n-1] + \frac{1}{8} y[n-2] = x[n]

Para usar los operadores ‘E’ y encontrar el polinomio característico, dado que el sistema es LTI – DT, se  puede desplazar 2 unidades:

y[n+2] - \frac{3}{4} y[n+1] + \frac{1}{8} y[n] = x[n+2]

La expresión usando notación de operadores ‘E’ es:

\Big[ E^2 - \frac{3}{4} E + \frac{1}{8} \Big] y[n] = E^2 x[n]

donde el numerador P[E] = E2
y el denominador Q[E] = E2 – (3/4)E +(1/8)

H(E) = \frac{P[E]}{Q[E]} = \frac{E^2}{E^2-\frac{3}{4}E+\frac{1}{8}}

literal b. Respuesta a impulso

A partir de la notación de operadores’E’,  se busca los modos característicos en el polinomio del denominador Q(E). También se expresa como y[n] como resultado de entrada cero x[n]=0

Q[E] = E^2 - \frac{3}{4} E + \frac{1}{8} \gamma^2 - \frac{3}{4} \gamma + \frac{1}{8} = 0

usando la fórmula o el algoritmo se busca las raices del polinomio Q(E):

(\gamma - \frac{1}{4})(\gamma - \frac{1}{2}) = 0 \gamma_1 = \frac{1}{4} ; \gamma_2 = \frac{1}{2}

Siendo raices reales y no repetidas, la respuesta del sistema tiene la forma:

y_c[n] = c_1 \Big( \frac{1}{4} \Big) ^n +c_2 \Big( \frac{1}{2} \Big) ^n

se convierte a la forma:

h[n] =\frac{b_n}{a_n} \delta [n] + y_c[n] \mu [n] h[n] = \Bigg[ c_1 \Big( \frac{1}{4} \Big) ^n +c_2 \Big( \frac{1}{2} \Big) ^n \Bigg] \mu [n]

Para encontrar los coeficientes c1 y c2, se requieren dos valores iniciales. En el literal a se indica que el sistema es causal, en consecuencia «No es posible obtener una salida antes que se aplique la entrada» y tenemos que: y[-1] = 0; y[-2] = 0.

La ecuación original de y[n] para respuesta a impulso tiene entrada x[n]=δ[t]

y[n] = \frac{3}{4} y[n-1] - \frac{1}{8} y[n-2] + x[n] y[n] = \frac{3}{4} y[n-1] - \frac{1}{8} y[n-2] + \delta [n]

donde el impulso tiene valor de 1 solo en n=0, haciendo equivalente h[n] equivale a y[n] , entonces:

h[n] = \delta [n] + \frac{3}{4} h[n-1] - \frac{1}{8} h[n-2]

n = 0

h[0] = \delta [0] + \frac{3}{4} h[-1] - \frac{1}{8} h[-2] h[0] = \delta [0] + \frac{3}{4} (0) - \frac{1}{8} (0) = 1

n=1

h[1] = \delta [1] + \frac{3}{4} h[0] - \frac{1}{8} h[-1] h[1] = (0) + \frac{3}{4} (1) - \frac{1}{8} (0) = \frac{3}{4}

con lo que es posible encontrar c1 y c2,

h[0] = 1 = c_1 \Big( \frac{1}{4} \Big) ^0 \mu [0] + c_2 \Big( \frac{1}{2} \Big) ^0 \mu (0) c_1 + c_2 = 1 h[1] = \frac{3}{4} = c_1 \Big( \frac{1}{4} \Big) ^1 \mu [1] + c_2 \Big( \frac{1}{2} \Big) ^1 \mu (1) \frac{1}{4} c_1 + \frac{1}{2} c_2 = \frac{3}{4}

resolviendo:

c_1 = -1 c_2 = 2 h[n] = - \Big( \frac{1}{4} \Big) ^n \mu [n] + 2 \Big( \frac{1}{2} \Big) ^n \mu [n] h[n] = \Bigg[ 2 \Big( \frac{1}{2} \Big) ^n - \Big( \frac{1}{4} \Big) ^n \Bigg] \mu [n]

Observaciones:

Las raíces características o frecuencias naturales del sistema se encuentran dentro del círculo de radio unitario. El sistema es asintóticamente estable, que implica que es BIBO estable.

h[n] no es de la forma k δ[n], por lo que el sistema global es con memoria.

La forma de respuesta al impulso se vuelve evidente que el sistema es IIR.

Algoritmo en Python

Usando el algoritmo de LTID – Respuesta impulso. Ejercicio con Python:

 respuesta impulso: 
hi: [1.   0.75]
Matriz: 
[[1.   1.  ]
 [0.5  0.25]]
Ch:  [ 2. -1.]
n>=0, hn:  
      n          n
- 0.25  + 2.0*0.5 
>>> 

Instrucciones Python

# Sistema LTID. Respuesta impulso
# QE con raices Reales NO repetidas Lathi ejemplo 3.13 pdf271
# Lathi ejemplo 3.18 y 3.19 pdf278,
# blog.espol.edu.ec/telg1001
import numpy as np
import sympy as sym
import matplotlib.pyplot as plt

# INGRESO
# coeficientes E con grado descendente
QE = [1.0, -3/4, 1/8]
PE = [1., 0., 0.]
# condiciones iniciales ascendente ...,y[-2],y[-1]
inicial = [0, 0.]

tolera = 1e-6  # casi_cero
muestras = 10  # para grafica

# PROCEDIMIENTO
# Respuesta a ENTRADA CERO
# raices, revisa numeros complejos
gamma = np.roots(QE)
revisaImag = np.iscomplex(gamma)
escomplejo = np.sum(revisaImag)

# coeficientes de ecuacion
m_q = len(QE)-1
Ac = np.zeros(shape=(m_q,m_q),dtype=float)

# revisa si parte compleja <tolera o casi_cero 
if escomplejo>0:
    for i in range(0,m_q,1):
        valorimag = np.imag(gamma[i])
        if np.abs(valorimag)<tolera:
            gamma[i] = float(np.real(gamma[i]))
    sumaj = np.sum(np.abs(np.imag(gamma)))
    if sumaj <tolera:
        print(sumaj)
        gamma = np.real(gamma)
        escomplejo = 0

# revisa repetidos
unicoscuenta = np.unique(gamma,return_counts=True)
repetidas = np.sum(unicoscuenta[1]-1)

# Respuesta impulso h[n]
ki = np.arange(-m_q,m_q,1)
hi = np.zeros(m_q, dtype=float)
xi = np.zeros(2*m_q, dtype=float)
xi[m_q] = 1 # impulso en n=0
Ah = np.zeros(shape=(m_q,m_q),dtype=float)

# h[n] iterativo
p_n = len(PE)
for i in range(0,m_q,1):
    for k in range(m_q,2*m_q,1):
        hi[i] = hi[i] - QE[k-m_q]*hi[m_q-k+1]
    for k in range(0,p_n,1):
        hi[i] = hi[i] + PE[k]*xi[m_q+i]
Bh = np.copy(hi[:m_q])

# coeficientes de ecuacion
for i in range(0,m_q,1):
    for j in range(0,m_q,1):
        Ah[i,j] = gamma[j]**(i)
ch = np.linalg.solve(Ah,Bh)

# ecuacion hn
n = sym.Symbol('n')
hn = 0*n
for i in range(0,m_q,1):
    hn = hn + ch[i]*(gamma[i]**n)

# SALIDA
if escomplejo == 0: 
    print('\n respuesta impulso: ')
    print('Bh:',Bh)
    print('Matriz: ')
    print(Ah)
    print('Ch: ',ch)
    print('n>=0, hn:  ')
    sym.pprint(hn)
else:
    print(' existen raices con números complejos.')
    print(' usar algoritmo de la sección correspondiente.')
    

# grafica datos
ki = np.arange(-m_q,muestras,1)
hi = np.zeros(muestras+m_q)

if escomplejo == 0: 
    # evaluación de h[n]
    h_n = sym.lambdify(n,hn)
    hi[m_q:] = h_n(ki[m_q:])

    # grafica h[n]
    plt.stem(ki,hi,label='h[n]',
             markerfmt ='C1o',
             linefmt='C2--')
    plt.legend()
    plt.grid()
    plt.ylabel('h[n]')
    plt.xlabel('ki')
    plt.title('h[n] = '+str(hn))
    plt.show()


literal c.  respuesta de paso s[n]

s[n] = \sum_{k=-\infty}^{n} h[k] s[n \rightarrow \infty] = \sum_{k=-\infty}^{\infty} h[k] = \sum_{k=-\infty}^{\infty} \Bigg[ 2 \Big( \frac{1}{2} \Big) ^k - \Big( \frac{1}{4} \Big) ^k \Bigg] \mu [k] = \sum_{k=0}^{\infty} \Bigg[ 2 \Big( \frac{1}{2} \Big) ^k - \Big( \frac{1}{4} \Big) ^k \Bigg]

considerando que:

\sum_{k=0}^{\infty} \alpha ^k = \frac{1}{1-\alpha} ; |\alpha|<1

se reemplaza con:

s[n \rightarrow \infty] = \sum_{k=0}^{\infty} \Bigg[ 2 \Big( \frac{1}{2} \Big) ^k - \Big( \frac{1}{4} \Big) ^k \Bigg] = 2\frac{1}{1-\frac{1}{2}} - \frac{1}{1-\frac{1}{4}} = 2(2) - \frac{4}{3} s[n \rightarrow \infty] =\frac{8}{3}

un método mas detallado es:

s[n] = \sum_{k=-\infty}^{n} h[k] s[n] = \sum_{k=-\infty}^{n} \Bigg[ 2 \Big( \frac{1}{2} \Big) ^k - \Big( \frac{1}{4} \Big) ^k \Bigg] \mu [k] s[n] = \sum_{k=0}^{n} \Bigg[ 2 \Big( \frac{1}{2} \Big) ^k - \Big( \frac{1}{4} \Big) ^k \Bigg] s[n] = 2 \sum_{k=0}^{n} \Big( \frac{1}{2} \Big) ^k - \sum_{k=0}^{n} \Big( \frac{1}{4} \Big) ^k

considerando que:

\sum_{k=0}^{n} \alpha ^k = \frac{1-\alpha ^{n+1}}{1-\alpha}

se tiene lo siguiente,

s[n] = 2 \frac{1-\frac{1}{2}^{n+1}}{1-\frac{1}{2}} - \frac{1-\frac{1}{4}^{n+1}}{1-\frac{1}{4}} s[n] = 2 \frac{1-\frac{1}{2} \Big(\frac{1}{2} \Big)^{n}}{\frac{1}{2}} - \frac{1-\frac{1}{4} \Big( \frac{1}{4} \Big)^{n}}{\frac{3}{4}} s[n] = 4 \Bigg[ 1-\frac{1}{2} \Big(\frac{1}{2}\Big)^{n} \Bigg]- \frac{4}{3} \Bigg[ 1- \frac{1}{4} \Big(\frac{1}{4} \Big)^{n} \Bigg] s[n] = 4-2 \Big(\frac{1}{2}\Big)^{n} -\frac{4}{3}+ \frac{1}{3} \Big(\frac{1}{4} \Big)^{n} lim _ {n \to \infty}s[n] = lim _ {n \to \infty} \Bigg[ \frac{8}{3}-2 \Big(\frac{1}{2}\Big)^{n} + \frac{1}{3} \Big(\frac{1}{4} \Big)^{n} \Bigg] lim _ {n \to \infty}s[n] = \frac{8}{3}

que es el mismo resultado que con el método presentado al inicio del literal.



Solución alterna: Usando transformada z

Revisar: LTID Transformada z – X[z] Fracciones parciales modificadas con Python

A partir de la ecuación de diferencias:

y[n+2] - \frac{3}{4} y[n+1] + \frac{1}{8} y[n] = x[n+2]

en transformada z (unidad 7), que es semejante a operador ‘E’,

z^2 Y[z]- \frac{3}{4} z Y[z] + \frac{1}{8}Y[z] = z^2 X[z] \Big[ z^2 - \frac{3}{4} z + \frac{1}{8}\Big] Y[z] = z^2 X[z]

la función de transferencia o respuesta al impulso se expresa como:

H[z] = \frac{X[z]}{Y[z]} = \frac{z^2}{z^2 - \frac{3}{4} z + \frac{1}{8}}

usandolas raices del polinomio Q(E)

H[z] = \frac{z^2}{(z - \frac{1}{4})(z - \frac{1}{2})}

para facilitar las operaciones se usan fracciones parciales modificadas,

\frac{H[z]}{z} = \frac{z}{(z - \frac{1}{4})(z - \frac{1}{2})} = \frac{k_1}{z - \frac{1}{4}} + \frac{k_2}{z - \frac{1}{2}}

Usando el método de «cubrir» de Heaviside:

k_1 = \frac{z}{\cancel{(z - \frac{1}{4})}(z - \frac{1}{2})} \Big|_{z=1/4} = \frac{1/4}{1/4 - \frac{1}{2}} = -1 k_2 = \frac{z}{(z - \frac{1}{4})\cancel{(z - \frac{1}{2})}} \Big|_{z=1/2} = \frac{1/2}{1/2 - \frac{1}{4}} = 2

reemplazando k1, k2 y restaurando a fracciones parciales al multiplicar ambos lados por z,

H[z] = \frac{-z}{z - \frac{1}{4}} + \frac{2z}{z - \frac{1}{2}}

se puede usar la tabla de transformadas z para obtener el mismo resultado que el desarrolo anterior con operador ‘E’

h[n] = - \Big( \frac{1}{4} \Big) ^n \mu [n] + 2 \Big( \frac{1}{2} \Big) ^n \mu [n]

s1Eva2016TII_T3 LTI DT Sistema 1er orden con constante a

Ejercicio: 1Eva2016TII_T3 LTI DT Sistema 1er orden con constante a

A partir del diagrama se tiene que:

literal a.1 ecuación que relaciona entrada-salida

El sistema es discreto de primer orden,

y[n] = x[n] + a y[n-1] y[n] - a y[n-1] = x[n]

literal a.2 Respuesta al impulso

Para la respuesta de impulso se escribe la expresión en notación de operador E, con x[n]=δ[n]

y[n+1] - a y[n] = x[n+1] Ey[n] - a y[n] = E x[n] (E -a) y[n] = E x[n]

para determinan las raíces del Q[E]

Q[E] = (\gamma-a) = 0 \gamma = a

que establece la ecuación

y_c[n] = c_1 a^n

Para encontrar el valor de la constante se evalua la expresión del sistema con  n=0, teniendo en cuenta que la ecuación general de respuesta a impulso:

h[n] = \frac{b_N}{a_N} \delta (n) + y_c[n] \mu [n] h[n] = y_c[n] \mu [n] = c_1 a^n \mu [n] y[n] = x[n] +a y[n-1] h[n] = \delta [n] +a h[n-1] h[0] = \delta [0] +a h[0-1]

El sistema es causal si para n<0 los valores son cero, h[-1] =0

h[0] = 1 +a (0) = 1

con lo que  par encontrar la constante c1,

h[0] = 1 = c_1 a^0 \mu [0] =c_1 (1)(1) = c_1

la constante c1 =1, quedando com respuesta a impulso:

h[n] = a^n \mu [n]

literal a.3 Respuesta de paso

s[n] = \sum_{k=-\infty}^{n} h[k] = \sum_{k=-\infty}^{n} a^k \mu [k]

Siendo:

\sum_{k=0}^{n} \alpha^k = \frac{1-\alpha^{n+1}}{1-\alpha} \mu [n]

se tiene que:

s[n] = \frac{1-a^{n+1}}{1-a} \mu [n]

Solución alterna, El mismo resultado se obtiene mediante:

s[n] = h[n] \circledast \mu [n]

usando la tabla de convolucion discreta:

= \sum_{k=- \infty}^{\infty} a^k \mu[k] \mu[n-k] = \sum_{k=0}^{n} a^k s[n] = \frac{1-a^{n+1}}{1-a} \mu [n]

literal b.

dado que h[n] no es de la forma k δ[n], el sistema es con memoria.

Si las raíces características, valores característicos o frecuencias naturales de un sistema, se encuentran dentro del círculo de radio unitario, el sistema es asintóticamente estable e implica que es BIBO estable.

h[n] = an μ[n] el sistema es de tipo IIR.


literal c. Sistema inverso

h[n] \circledast h_{inv} = \delta[n] y[n] = x[n] + a y[n-1] h[n] = \delta[n] + a h[n-1] h[n] - a h[n-1] = \delta[n] h[n] \circledast \delta [n] - a h[n] \circledast \delta [n-1] = \delta[n] h[n] \circledast (\delta [n] - a \delta [n-1]) = \delta[n] h_{inv}[n] = (\delta [n] - a \delta [n-1]) = \delta[n]

h[n] no es de la forma k δ[n], por lo que el sistema inverso es con memoria.

h[n]=0 para n<0, se entiende que el sistema inverso es Causal.

La respuesta impulso del sistema inverso es absolutamente sumable o convergente, por lo que el sistema es BIBO estable.

la forma de hinv[n] el sistema es de tipo FIR.


literal d. Diagrama de bloques del sistema inverso

El diagrama del sistema inverso sería: