3.5 Respuesta a impulsos – discreto

Referencia: Lathi 3.8 pdf193, Oppenheim 2.1.2 p77/pdf105, Schaum 2.6.C p62,

Una señal discreta de entrada x[n] se representa como una superposición de versiones escaladas de un conjunto de impulsos unitarios desplazados δ[n-k], cada uno con valor diferente de cero en un solo punto en el tiempo, especificado por el valor de k.

La respuesta de un sistema lineal y[n] a x[n] es la superposición de las respuestas escaladas del sistema a cada uno de estos impulsos desplazados.

Si usamos hk[n] como la respuesta del sistema lineal al impulso unitario desplazado por δ[n-k]. La respuesta y[n] del sistema lineal a la entrada x[n] en la ecuacion será la combinación lineal ponderada de las respuesta básicas.

y[n] = \sum_{k=-\infty}^{+\infty} x[k]h_{k}[n]

Si se conoce la respuesta de un sistema lineal al conjunto de impulsos unitarios desplazados, podemos construir una respuesta a una entrada arbitraria.

Si el sistema lineal también es invariante en el tiempo, entonces estas respuestas a impulsos unitarios desplazados en el tiempo son todas las versiones desplazadas en el tiempo unas de otras.

h[n] es la salida del sistem LTI cuando δ[n] es la entrada. Entonces para un sistema LTI la ecuación se vuelve.

y[n] = \sum_{k=-\infty}^{+\infty} x[k]h[n-k]

El resultado se conoce como la «suma de convolución» o «suma de superposición» y la operación miembro derecho de la ecuación se llama convolución de las secuencias x[n] y h[n] que se representa de manera simbólica como:

y[n] = x[n]*h[n]

Ejemplo

Referencia: Openheim Ejemplo 2.4 p85/pdf113

Considere una entrada x[n] y una respuesta al impulso unitario h[n] dada por:

x[n] = (u[n]-u[n-5])
h[n] = αn (u[n]-u[n-7])

con α>1.  Para el ejercicio, α=1.5.

Nota: si α es un entero, por ejemplo 2, usar α=2.0, para que la operación potencia se realice con números reales, como entero se puede saturar y dar error.


Al aplicar la suma de convolución obtiene:

Para aplicar el algoritmo se requiere definir u[n], por ser parte de las funciones x[n] y h[n]. Dado que la operación requiere valores fuera del rango muestreado para n, la sección suma convolución utiliza las funciones en lugar de los vectores xi, hi.

La función está definida en un intervalo simétrico, por lo que el rango de trabajo [a,b] se mantiene de la forma [-b,b] en las intrucciones.

# Ejemplo Oppenheim Ejemplo 2.3 p83/pdf111
import numpy as np

# INGRESO
# Rango [a,b], simétrico a 0
b = 15 ; a = -b

alfa = 1.5
u = lambda n: np.piecewise(n,n>=0,[1,0])

x = lambda n: u(n)-u(n-5)
h = lambda n: (alfa**n)*(u(n)-u(n-7))


# PROCEDIMIENTO
ni = np.arange(a,b+1,1)
xi = x(ni)
hi = h(ni)

# Suma de Convolucion x[n]*h[n]
muestras = len(xi)
yi = np.zeros(muestras, dtype=float)
for i in range(0,muestras):
    suma=0
    for k in range(0,muestras):
        suma = suma + x(ni[k])*h(ni[i]-ni[k])
    yi[i]= suma

# yi = np.convolve(xi,hi,'same')

# SALIDA - GRAFICA
import matplotlib.pyplot as plt

plt.figure(1)
plt.suptitle('Suma de Convolución x[n]*h[n]')

plt.subplot(311)
plt.stem(ni,xi,
         linefmt='b--',
         markerfmt='bo',
         basefmt='k-')
plt.ylabel('x[n]')

plt.subplot(312)
plt.stem(ni,hi,
         linefmt='b--',
         markerfmt='ro',
         basefmt='k-')
plt.ylabel('h[n]')

plt.subplot(313)
plt.stem(ni,yi,
         linefmt='g-.',
         markerfmt='mo',
         basefmt='k-')
plt.ylabel('x[n]*h[n]')
plt.xlabel('n')

plt.show()

La suma convolución se encuentra también disponible con Numpy en np.convolve(), la  sección de suma convoluacion se puede reemplazar y obtener los mismos resultados. Considere que para éste caso se usan los vectores xi y hi.

yi = np.convolve(xi,hi,'same')

el algoritmo se puede aplicar a otros ejercicio para comprobar los resultados.

Continúe con el ejercicio 2.5 del libro.

LTI en dominio frecuencia

Oppenheim Ejemplo 4.15 (p.345 pdf)

Ejemplo: Considere un sistema continuo LTI con respuesta al impulso:

h(t) = δ(t-t0)

La respuesta en frecuencia es la transformada de Fourier :

H(jω) = e-jωt0

A cualquier entrada x(t) con transformada X(jω) , la salida en dominio de la frecuencia es:

Y(jω) = H(jω)X(jω) = e-jωt0 X(jω)

que usando la propiedad de desplazamiento:

y(t)=x(t-t0)

Oppenheim Ejemplo 4.18 (p.346 pdf)

El Filtrado selectivo en frecuencia se puede representar en un filtro PASA-BAJO ideal:

H(jω) = { 1    |ω|< ωc
        { 0    |ω|> ωc

cuya respuesta en el dominio de la frecuencia es:

h(t) = sen(ωc t)/πt

Como ejemplo: para una frecuencia fs=20 y usando la transformada para una función par se obtiene:

Observaciones: h(t) no es cero para t<0, en consecuencia el filtro pasabajos ideal no es causal. NO es fácil aproximarse mucho a un filtro ideal con componentes reales, por lo que se contruye ms facilmente un filtro no ideal.

Transformada Laplace

Para una señal X(t), Su transformada de Laplace esta definida como:

X(s) = \int_{-\infty}^{\infty} x(t) e^{-st} \delta t

la señal x(t) es la inversa de la transformada X(s), que se obtiene como:

x(t) = \frac{1}{2πj} \int_{c-j\infty}^{c+j\infty} x(t) e^{st} \delta t

donde c es una constante seleccionada para asegurar la convergencia de la integral.

Los pares de ecuaciones conocidos como «Pares de la transformada de Laplace» se escriben simbólicamente como:

X(s) = L[x(t)] x(t) = L^{-1}[X(s)]

y propiedades como

L[a_1 x_1(t) + a_2 x_2 (t)] = a_1 X_1(s) + a_2 X_2 (s)

que muestra la linealidad de la transformada

Región de convergencia (ROC)

También conocida como región de existencia, es el grupo de valores de s, en la region del plano complejo, para los cuales la integral de la transformada converge.

Ejemplo

referencia: Lathi ejemplo 4.1. pdf231

Para una señal x(t) = e-atu(t), encuentre la transformada X(s) y su región de convergencia (ROC)

X(s) = \int_{-\infty}^{\infty} e^{-at} \mu (t) e^{-st} \delta t

pero tomando en cuenta que u(t)=0 para t<0 y u(t)=1 para t≥0:

X(s) = \int_{0}^{\infty} e^{-at} \mu(t) e^{-st} \delta t X(s) = \int_{0}^{\infty} e^{-at} e^{-st} \delta t = \int_{0}^{\infty} e^{-(s+a)t} \delta t = \frac{1}{s+a} \Big|_{0}^{\infty}

Si s es un número complejo y t →∞, el término exponencial no necesariamente desaparece.

El número complejo se puede expresar como:

z = α + jβ
e-zt = e-(α + jβ)t = e-αte-jβt.

El valor |e-jβt| = 1 sin importar el valor de βt.
Entonces si t→∞, e-zt →0 solo si α>0
y e-zt →∞ en el caso que α<0.

La región de convergencia (ROC) de X(s) es que la parte real s>-a.

Tabla de transformadas de Laplace (wikipedia)

Tabla de transformadas y propiedades de Laplace (Universidad de Córdova)

Propiedad de diferenciación

\frac{\delta x}{\delta t} \Leftrightarrow sX(s) - x(0^{-}) \frac{d^2x}{dt^2} \Leftrightarrow s^2X(s) - sx(0^{-}) - x'(0^{-}) \frac{\delta^3x}{\delta t^3} \Leftrightarrow s^3 X(s) - s^2 x(0^{-}) - sx'(0^{-}) - x''(0^{-})

Transformada Laplace Ejercicio01

Ejemplo 1:

Referencia: Oppenheim ejemplo 9.37 p718/pdf746, semejante al ejercicio sistema-modelo LTI Lathi 1.8-1 pdf/p.80. Oppenheim problema 2.61c pdf/p.191

Suponga un sistema LTI causal que se describe mediante la ecuación diferencial:

\frac{d^2y(t)}{dt} + 3 \frac{dy(t)}{dt} + 2y(t) = x(t)

junto con la condición de reposo inicial:

H(s) = \frac{1}{s^2 + 3s +2}

Suponga que x(t) = α μ(t) , entonces la entrada del sistema es x(s) = α

Siempre que las señales de entrada son cero para t<0, la convolución es la multiplicación de los términos. (no aplica para valores t<0)

y(s) = H(s) X(s) = \Big[\frac{1}{s^2 + 3s +2} \Big]\Big[\frac{\alpha}{s}\Big] y(s) = \frac{\alpha}{2}\frac{1}{s} - \alpha\frac{1}{s+1} +\frac{\alpha}{2}\frac{1}{s+2}

por lo que aplicando la anti transformada:

y(t) = \alpha \Big [ \frac{1}{2} - e^{-t} + \frac{1}{2} e^{- 2t} \Big] u(t)

Usando Python y Sympy

Para crear la expresión polinómica simple, no usamos el factor α, que se añade luego de obtener el resultado, se obtiene:

    a         a      a 
--------- - ----- + ---
2*(s + 2)   s + 1   2*s

resultado al que hay que multiplicar por la constante α todos sus términos.

las instrucciones son:

# Oppenheim ejemplo 9.37 p718/pdf746
import sympy as sym

s = sym.Symbol('s')
a = sym.Symbol('a')
Hs = (1/(s**2+3*s+2))*(a/s)
separa = Hs.apart(s)

sym.pprint(separa)

Ejemplo 2

Referencia: Oppenheim ejemplo 9.38 p719/747

Considere el sistema caracterizado por la ecuación diferencial con condiciones iniciales.

\frac{d^2y(t)}{dt^2} + 3 \frac{dy(t)}{dt} + 2y(t) = x(t) y(0^{-}) = \beta \text{, } y'(0^{-}) = \gamma

sea x(t) = αμ(t)

Aplicando la transformada de laplace en ambos lados de la ecuación con condiciones iniciales y aplicando la propiedad de diferenciación:

\frac{d^2y(t)}{dt^2} = s^2Y(s) - y(0^{-})s - y'(0^{-}) = s^2Y(s) - \beta s - \gamma \frac{dy(t)}{dt} = sY(s) - \beta

resultados que se insertan en la ecuacion original:

[s^{2}Y(s) - \beta s - \gamma] + 3[sY(s) - \beta] + 2Y(s) = \frac{\alpha}{s}

para simplificar como:

[s^{2}+3s +2] Y(s) - \beta s - (\gamma + 3 \beta) = \frac{\alpha}{s} [s^{2}+3s +2] Y(s) = \frac{\alpha}{s}+ \beta s + (\gamma + 3 \beta) (s+1)(s+2)Y(s) = \frac{\alpha}{s}+ \beta s + (\gamma + 3 \beta) Y(s) = \frac{\alpha}{s(s+1)(s+2)}+ \frac{\beta s}{(s+1)(s+2)} + \frac{(\gamma + 3 \beta)}{(s+1)(s+2)}

se encuentra que:

Y(s) = \frac{\beta (s+3)}{(s+1)(s+2)} + \frac{\gamma}{(s+2)(s+2)} + \frac{\alpha}{s(s+2)(s+2)}

donde Y(s) es la transformada unilateral de Laplace de y(t)

el último término, representa la respuesta del sistema LTI causal y la condición de reposo inicial, conocido también como la respuesta en estado cero.

Una interpretación análoga se aplica a los primeros dos términos del miembro derecho de la ecuación, que representa la respuesta del sistema cuando la entrada es cero (α=0) , conocida también como respuesta a entrada cero.

Observe que la respuesta a entrada cero es una función lineal de los valores de las condiciones iniciales. Demostrando que la respuesta es la superposición del estado cero y respuesta a entrada cero.

si α=2, β = 3 y γ=-5, al realizar la expansión en fracciones parciales encontramos que:

Y(s) = \frac{1}{s} - \frac{1}{s+1} + \frac{3}{s+2}

con la aplicación de las tablas de transformadas se obtiene:

y(t) = \Big[ 1 - e^{-t} + 3 e^{-2t} \Big] \mu (t) \text{, para } t\gt 0

usando Python y Sympy

se puede trabajar la ecuación para obtener:

ecuacion de entrada
    2                                   a
Ys*s  + 3*Ys*s + 2*Ys - b*s - 3*b - g = -
                                        s
ecuacion Ys:  
        2               
 a + b*s  + 3*b*s + g*s 
[----------------------]
      / 2          \    
    s*\s  + 3*s + 2/    
fracciones parciales de Ys: 
 a    a - 2*b - 2*g   a - 2*b - g
--- + ------------- - -----------
2*s     2*(s + 2)        s + 1   
sustituye a=2, b= 3 y g=-5
  3       1     1
----- - ----- + -
s + 2   s + 1   s
>>>

con las siguientes instrucciones, presentadas por cada sección de la parte de desarrollo analítico, para ir comparando los resultados.

# Oppenheim ejemplo 9.38 p719/747
import sympy as sym

s,a,b,g  = sym.symbols('s a b g')
Ys = sym.Symbol('Ys')

d2y = (s**2)*Ys - b*s - g
dy  = s*Ys - b

ecuacion = sym.Eq(d2y + 3*dy + 2*Ys, a/s)
ecuacionYs = sym.solve(ecuacion,Ys)
parciales = ecuacionYs[0].apart(s)

solucion = parciales.subs({a:2,b:3,g:-5})

print('ecuacion de entrada')
sym.pprint(ecuacion)
print('ecuacion Ys:  ')
sym.pprint(ecuacionYs)
print('fracciones parciales de Ys: ')
sym.pprint(parciales)
print('sustituye a=2, b= 3 y g=-5')
sym.pprint(solucion)

La definición de la ecuacion se realiza usando las variables como símbolos, incluyendo Y(s) como Ys. En éste ejercicio no existe X(s), por lo que no se la incluyó.

Luego se busca la expresión solución para Ys, por simplicidad se la cambia a fracciones parciales, asi es más sencillo usar la tabla de transformadas de Laplace.

Para una solución con los valores de a,b y g, se usa un diccionario {} con las parejas de variables y sus valores. Compare los resultados con lo expuesto en la parte teórica.

FFT de una señal coseno

Señal de entrada: x(t)= cos(ω0 t)
ω0 = (2π fs) en radianes, para el ejemplo fs=20 Hz

x(ω) =π[δ(ω-ω0)+δ(ω+ω0)] con eje en radianes.
La gráfica se muestra en Hz. fs0/2π
Se muestran las partes:
– real e imaginaria de la transformada y
– magnitud y fase

Script de python para usar la transformada rapida de Fourier.
En la seccion de función de entrada, se muestra para cuatro funciones: coseno, constante, impulso y sinc. Para usar una de ellas se elimina el simbolo de comentario (##)

# Transformadas Rapida de Fourier
# entrada funcion par
# Propuesta: edelros@espol.edu.ec
import numpy as np
import scipy.fftpack as fourier
import matplotlib.pyplot as plt

# Definir la funcion de Entrada par
def entradax(t,fs):
    # función matemática CAMBIAR AQUI
    x=np.cos(2*np.pi*fs*t)

##    x=1  #constante
    
##    x=0   # impulso
##    if abs(t)<(1/fs):
##        x=1

##    if t==0:  #para sinc aumentar k=8
##        t=(1/fs)/100000
##    x=(np.sin(2*np.pi*fs*t)/(np.pi*t))

    # función matemática CAMBIAR AQUI
    return(x)

# PROGRAMA
# INGRESO
fs=20       #float(input('frecuencia Señal Hz:'))
ancho=1     #ancho del intervalo
t0=0        #intervalo de t
tn=ancho*1/fs
resolucion=1    #resolucion de transformada

# PROCEDIMIENTO
n= 1024*resolucion 
dt=(tn-t0)/n   # intervalo de muestreo
# Analógica Referencia para mostrar que es par
t=np.arange(-tn,tn,dt)  # eje tiempo analógica
m=len(t)
xanalog=np.zeros(m, dtype=float)
for i in range(0,m):
    xanalog[i]=entradax(t[i],fs)

# FFT: Transformada Rapida de Fourier
# Analiza la parte t>=0 de xnalog muestras[n:2*n]
# Transformada de Fourier
xf=fourier.fft(xanalog[n:2*n])
xf=fourier.fftshift(xf)
# Rango de frecuencia para eje
frq=fourier.fftfreq(n, dt)
frq=fourier.fftshift(frq)

# x[w] real
xfreal=(1/n)*np.real(xf)
# x[w] angulo
xfimag=(1/n)*np.imag(xf)
# x[w] magnitud
xfabs=(1/n)*np.abs(xf)
# x[w] angulo en radianes
xfangle=np.angle(xf)# np.unwrap(np.angle(xf))


#SALIDA
plt.figure(1)   # define la grafica
plt.suptitle('Señal de Entrada')
plt.ylabel('xanalog[t]')
plt.xlabel('tiempo')
plt.plot(t,xanalog)
plt.margins(0,0.05)
plt.grid()

plt.figure(2)   # define la grafica
plt.suptitle('Transformada Rápida Fourier FFT')
ventana=ancho*0.01 # ventana de frecuencia alrededor f=0
ra=int((n/2)*(1-ventana)) #rango de indices
rb=int((n/2)*(1+ventana))+1
a=frq[ra] #rango de ejes uniformes entre graficos
b=frq[rb]
c=np.min([-0.05,1.1*np.min(xfreal)])
d=1.1*np.max(xfabs)
ejes=[a,b,c,d]

plt.subplot(221)    # grafica de 2x2, subgrafica 1
plt.ylabel('x[f] real')
plt.xlabel(' frecuencia (Hz)')
plt.axis(ejes)
plt.stem(frq[ra:rb],xfreal[ra:rb])
plt.grid()

plt.subplot(223)    # grafica de 2x2, subgrafica 2
plt.ylabel('x[f] imag')
plt.xlabel(' frecuencia (Hz)')
c=np.min([-0.05,1.1*np.min(xfimag)])
ejes=[a,b,c,d]
plt.axis(ejes)
plt.stem(frq[ra:rb],xfimag[ra:rb])
plt.grid()

plt.subplot(222)    # grafica de 2x2, subgrafica 3
plt.ylabel('x[f] magnitud')
plt.xlabel(' frecuencia (Hz)')
c=np.min([-0.05,1.1*np.min(xfabs)])
ejes=[a,b,c,d]
plt.axis(ejes)
plt.stem(frq[ra:rb],xfabs[ra:rb])

plt.grid()

plt.subplot(224)    # grafica de 2x2, subgrafica 4
plt.ylabel('x[f] fase radianes')
plt.xlabel(' frecuencia (Hz)')
c=1.1*np.min(xfangle)
d=1.1*np.max(xfangle)
ejes=[a,b,c,d]
plt.axis(ejes)
plt.stem(frq[ra:rb],xfangle[ra:rb])
plt.grid()

plt.show()

En el caso del seno(ω0 t), las graficas resultantes son:


con transformada:

continuar el ejercicio con las otras funciones.

Periodicidad señales analógicas y digitales

Considere la exponencial compleja discreta con frecuencia (ω0 + 2π)

ej(ω0 + 2π)n = ej2πn e0n = e0n

Se puede ver que la exponencial con frecuencia ω0 + 2π es la misma que aquella con la frecuencia ω0.
Existe una situación diferente en el caso contínuo, en el cual las señales e0t son todas distintas para distintos valores de ω0.
En el caso discreto, éstas señales no son distintas, ya que la señal con frecuencia ω0 es idéntica a las señales con frecuencias ω0 ± 2π , ω0 ± 4π, … y las que le siguen. Por lo tanto, al considerar las exponenciales complejas, necesitamos tomar en cuenta solamente un intervalo de frecuencia de longitud 2π dentro del cual se escoge ω0.
ω0N = 2π m
o de otra forma
ω0/2π = m/N

Para el caso de m=1 y N=2

otra prueba con m=7 y N=4

con m=2 y N=1

El código en python para realizar las gráficas es:

# Señales discretas
# y(t)  a  y[t]
# propuesta:edelros@espol.edu.ec
import numpy as np
import math
import matplotlib.pyplot as plt

# Definir la funcion para el ejemplo
def ejemplo_funcion(f,t):
    # función matemática CAMBIAR AQUI
    y=np.cos(f*t)
    # función matemática CAMBIAR AQUI
    return(y)

# Programa para graficar la función
# INGRESO
rango=float(input('rangos en periodos de analogica:'))
fa=float(input('frecuencia de analógica:'))
fd=float(input('frecuencia digital:'))

# PROCEDIMIENTO
n=500  # puntos en eje t para el gráfico analogica
t0=-rango*2*np.pi*fa/2
tn=rango*2*np.pi*fa/2
delta=(tn-t0)/n
deltaD=(2*np.pi*fa)/(2*fd)
nd=int((tn-t0)//deltaD +1)   # número de muestras para gráfica
# Dimensiona los arreglos
t=np.zeros(n, dtype=float)  # para la analogica
yanalog=np.zeros(n, dtype=float)
td=np.zeros(nd, dtype=float)    # Para la digital
ydigital=np.zeros(nd, dtype=float)
# valores para analogica
for i in range(0,n):
    t[i]=t0+delta*i     # para el eje x
    yanalog[i]=ejemplo_funcion(fa,t[i])
#valores para digital
for i in range(0,nd):
    td[i]=t0+deltaD*i   #para el eje x
    ydigital[i]=ejemplo_funcion(fa,td[i])

# SALIDA
#Escala y para el grafico
ymax=np.max(yanalog)+0.1*np.max(yanalog)# rango en el eje y
ymin=np.min(yanalog)-0.1*np.max(yanalog)
plt.figure(1)       # define la grafica
plt.suptitle('Señal Analogica vs digital')
plt.subplot(211)    # grafica de 2x1 y subgrafica 1
plt.ylabel('y(t)')
plt.axis((t0,tn,ymin,ymax))
plt.plot(t,yanalog)
plt.subplot(212)    # grafica de 2x1 y subgrafica 2
plt.ylabel('Digital: y[t]')
plt.axis((t0,tn,ymin,ymax))
plt.plot(td,ydigital,'ro')
plt.show()

 

Fracciones Parciales – Separar con Python

Referencia: Oppenheim 9.4 p671/pdf699. Lathi B5 pdf21, Schaum/HSu 3.8.C p120

Se pueden separar en fraciones parciales del siguiente ejemplo:

H(s) = \frac{2s^2 +5}{s^2 +3s+2}

se requiere obtener los residuos, polos y términos directos al separar en fracciones parciales entre dos polinomios P(s)/Q(s).

\frac{P(s)}{Q(s)} = \frac{ceros(1)}{s - polos(1)}+\frac{ceros(2)}{s - polos(2)}+\text{ ... } \text{ ... }+\frac{ceros(n)}{s - polos(n)}+ganancia(s)

Se puede realizar de dos formas usando Python, en ambas se requiere definir los vectores P y Q.

Sympy – fórmulas simbólicas

La primera usando Sympy para crear las fórmulas simbólicas y luego formar una  fracción, separar las partes Sympy.apart()

Para los datos proporcionados en el ejemplo se tiene como resultado:

P = [2, 0, 5]
Q = [1, 3, 2]

Fracciones parciales:
2 - 13/(s + 2) + 7/(s + 1)

con las instrucciones sympy:

# Expansión en fracciones parciales
# P es numerador, Q es denominador

P = [2, 0, 5]
Q = [1, 3, 2]

import sympy as sym
s = sym.Symbol('s')
n = len(P)
m = len(Q)
Ps = 0
for i in range(0,n,1):
    Ps = Ps + P[i]*s**(n-i)
Qs = 0
for i in range(0,m,1):
    Qs = Qs + Q[i]*s**(m-i)
fraccion = Ps/Qs
parciales = sym.apart(fraccion)

# Salida - polinomio
print('Fracciones parciales:')
print(parciales)

Las fórmulas simbólicas permiten agrupar y realizar operaciones entre varios sistemas además de expandir, separar, y evaluar las expresiones. Util principalmente para la sección de transformadas de Laplace.

Un resumen de sympy se muestra en:

Fórmulas simbólicas – Sympy


Scipy – Procesar señales

La segunda forma es usando libreria scipy.signal.residue con los valores de P y Q:

ceros(n)   :  [-13.   7.]
polos(n)   :  [-2. -1.]
ganancia(s):  [ 2.]
>>>

Con los resultados el equivalente en Laplace es:

transformado al dominio del tiempo,
H(t) = (−13e−2t + 7e−t)u(t) + 2δ(t)

las instrucciones adicionales a la sección anterior son:

# Usando scipy.signal.residue

import scipy.signal as senal
ceros,polos,ganancia = senal.residue(P,Q)

# SALIDA - coeficientes
print()
print('ceros(n)   : ', ceros)
print('polos(n)   : ', polos)
print('ganancia(s): ', ganancia)

Diagramas de bloques con Z

Considere el sistema LTI causal descrito mediante:

y[n] - (1/4) y[n-1] = x[n]

Ejemplo Openheim 10.28. con función del sistema:

H(z)  = 1/ [1- (1/4) z-1]

Aqui z-1 esla función del sistema de un retardo unitario. el diagrama de bloques en la figura contiene un lazo de retroalimentación.

y[z] - (1/4) z-1 y[z] = x[z]

y[z]  = x[z] + (1/4) z-1 y[z]

Ejemplo Openheim 10.29. Considere un sistema LTI causal con función del sistema:

H(z) = [ 1 - 2z-1 ] / [ 1 - (1/4)z-1 ] 
     = {1/ [ 1 - (1/4)z-1 ]} [ 1 - 2z-1 ]

Diagramas de Bloques con s

Representaciones en diagramas de bloques para los sistemas LTI causales descritos por ecuaciones diferenciales y funciones racionales del sistema. (Oppenheim 9.8.2)

Ejemplo 1.

referencia: Oppenheim ejemplo 9.28 p708/pdf736

Considere un sistema LTI causal

\frac{\delta y(t)}{\delta t} + 3y(t) = x(t) Dy(t) + 3y(t) = x(t) Dy(t) = x(t) -3y(t) y(t) = \frac{1}{D}[x(t) - 3y(t)]

usando s en lugar de D

y(s) = \frac{1}{s} [x(s) - 3 y(s)]

usando la transformada de Laplace, cambiando s por d/dt desde la primera ecuación, la función del sistema es:

( s+3) y(s) = x(s)

por lo que y(s)/x(s) es:

H(s) = \frac{1}{s+3}

al cambiar el signo en el coeficiente de la multiplicación de la trayectoria de retroalimentación, obtenemos una forma de diagrama de bloques similar a la anterior, pero con la que se verifica que:

H(s) = \frac{\frac{1}{s}}{1+\frac{3}{s}} = \frac{1}{s+3}

que confirma que:

\frac{Y(s)}{X(s)} = H(s) = \frac{H_1(s)}{1+H_1(s)H_2(s)}

siendo H1 el bloque superior, y H2 el bloque inferior en el diagrama.


Ejemplo 2

referencia: Oppenheim ejemplo 9.30 p711/pdf739

Considere el sistema causal de segundo orden con la funcion del sistema:

H(s) = \frac{1}{(s+1)(s+2)} = \frac{1}{s^2 + 3s +2}

la entrada x(t) y la salida y(t) para este sistema satisfacen la ecuación diferencial:

\frac{\delta^2(t)}{\delta t^2} + 3 \frac{\delta y(t)}{dt} + 2 y(t) = x(t) (s^2 + 3s + 2)y(t) = x(t)

pero también se puede armar como:

H(s) = \frac{1}{(s+1)(s+2)} = \frac{k_1}{s+1} + \frac{k_2}{s+2}

siendo el numerador :

k_1s + 2k_1 + k_2s + k_2 = (k_1 + k_2)s + (2k_1 +k_2)

que se observa que debe dar:
(k_1 + k_2)s + (2k_1 +k_2) = 0s + 1

entonces

k_1 + k_2 = 0 2k_1 + k_2 = 1

con lo que

k1 = 1 y
K2 = -1
H(s) = \frac{1}{s+1} -\frac{1}{s+2}

Si

H(s) = \frac{1}{(s+1)(s+2)}

se escribe como

H(s) = \Big[\frac{1}{s+1}\Big]\Big[\frac{1}{s+2}\Big]

que se representa como la serie o cascada de dos sistemas de primer orden.

que por cierto, también es la convolución de dos sistemas en serie.

También es posible realizar las fracciones parciales con Sympy y se obtiene:

    1       1  
- ----- + -----
  s + 2   s + 1
>>>

usando las instrucciones:

# Oppenheim 9.30 p711/pdf739
import sympy as sym

s = sym.Symbol('s')
Hs = 1/((s+1)*(s+2))
separa = Hs.apart()

sym.pprint(separa)

Ejemplo 3.

referencia: Oppenheim ejemplo 9.31.p712/pdf739

Considere la funciín del sistema

H(s) = \frac{2s^2 + 4s - 6}{s^2 + 3s + 2}

se puede reescribir como:

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

y si las derivadas se toman del sistema anterior, se tiene una forma mas simplificada:

también, con un poco de trabajo, reordenando las ecuaciones, se podría convertir en:

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

o también en :

H(s) = 2 + \frac{6}{s+2} - \frac{8}{s+1}

que requiere un sistema en forma paralelo.
QUEDA COMO TAREA HACER EL DIAGRAMA.

Se presenta el ejercicio usando Sympy:

      6       8  
2 + ----- - -----
    s + 2   s + 1
>>>

que se obtiene usando las instrucciones:

# Oppenheim 9.31.p712/pdf739
import sympy as sym

s = sym.Symbol('s')
Hs = (2*s**2+4*s-6)/(s**2+3*s+2)
separa = Hs.apart()

sym.pprint(separa)

Sistemas – Invarianza en el tiempo

Referencia: Oppenheim 1.6.5 p50/pdf78, Lathi 1.7-2 pdf76, Schaum/Hsu p18

Si el comportamiento de un sistema y las características del mismo están fijos en el tiempo, se los llama invariantes en el tiempo.

Por ejemplo, para un circuito RC, los valores de la resistencia y capacitor fijos no varian si se realiza un experimento o medición hoy o mañana.

Expresando lo mismo como:
Un sistema es invariante en el tiempo si aun corrimiento de tiempo en la señal de entrada ocasiona un corrimiento en el tiempo en la señal de salida.

y(t) = x(t)
y(t-t0) = x(t-t0)

Ejemplo 1

Referencia: Oppenheim 1.14 p51/pdf78

Para muestra de lo indicado, considere un sistema definido como:

y(t) = \sin(x(t)) \text{, siendo x(t)=t} y(t-t_0) = \sin(x(t - t_0))

Por otro lado, si la respuesta de la señal se modifica en el tiempo, el sistema será VARIANTE en el tiempo:

Ejemplo 2

Referencia: Oppenheim 1.14 p51/pdf78. Schaum ejercicio 1.39 p48

Considere el sistema:

y(t) = x(2t)

que tiene escalamiento en el tiempo. y(t) es una versión comprimida de x(t). Los corrimientos/desplazamientos en el tiempo también serán comprimidos por un factor de 2, por lo que por ésta razón el sistema no es invariante en el tiempo.

Siguiendo las indicaciones en el problema presentado en Schaum, donde T se conoce como el «operador lineal«:

El sistema tiene una relación de entrada-salida dada por:

y(t) = T{x(t)} = x(2t)

y t0 es un desplazamiento en el tiempo.

Sea y2(t) la respuesta a x2(t) = x(t-t0) entonces:

y2(t) = T{x2(t)} = x2(2t) = x(2t-t0)

y(t-t0) = x(2(t-t0))

y(t-t0) ≠ y2(t), por lo que el sistema no es invariante en el tiempo. 
El sistema x(2t) se conoce como compresor, pues crea una secuencia de salida cada 2 valores de la secuencia de entrada.

Ejemplo 2

Referencia: Schaum/HSU ejercicio 1.38 p47

Desarrollar semejante al ejercicio anterior, en forma analítica y usando Python.

x(t) = t^2 y(t) = t x(t)

sist_invariatiempo02


Desarrollo en Python

Ejemplo 01

# Sistemas – Invariantes en el tiempo
# Oppenheim 1.14 p51/pdf78
import numpy as np

# INGRESO
a = 0; b = 2*np.pi ; dt=0.1
t0 = 1

xt = lambda t: t
yt = lambda t: np.sin(xt(t))

# PROCEDIMIENTO
t = np.arange(a,b,dt)

x1 = xt(t)
y1 = yt(t)


# Corrimiento en la entrada
t2x = t + t0
x2 = x1
# efecto en salida
t2y = t + t0
y2 = y1

# evaluación de tiempos desplazados
y3 = yt(t-t0)

# SALIDA - GRAFICA
import matplotlib.pyplot as plt
plt.figure(1)

plt.plot(t,y1, label='y1(t)')
plt.plot(t2y,y2, 'm.', label='y2(t)')
plt.plot(t,y3, label='y1(t-t0)')

plt.xlabel('t')
plt.legend()
plt.show()

Ejemplo 02

# Sistemas – Invariantes en el tiempo
# Oppenheim 1.16 p52/pdf80
import numpy as np

# INGRESO
a = -5; b = 5 ; dt=0.1
t0 = 2

u = lambda t: np.piecewise(t,t>=0,[1,0])
xt = lambda t: u(t+2)-u(t-2)
yt = lambda t: xt(2*t)

# PROCEDIMIENTO
t = np.arange(a,b,dt)

x1 = xt(t)
y1 = yt(t)

# Corrimiento en la entrada
t2x = t + t0
x2 = x1
# efecto compresor en salida
t2y = t + t0/2
y2 = y1

# evaluación de tiempos desplazados
x3 = xt(t-t0)
y3 = yt(t-t0)

# SALIDA - GRAFICA
import matplotlib.pyplot as plt
plt.figure(1)
plt.subplot(311)
plt.xlim([min(t),max(t)])
plt.plot(t,x1, label='x1(t)')
plt.plot(t,y1, label='y1(t)')
plt.legend()

plt.subplot(312)
plt.xlim([min(t),max(t)])
plt.plot(t2x,x2, label='x2(t)=x1(t-t0)')
plt.plot(t2y,y2, label='y2(t)')
plt.legend()

plt.subplot(313)
plt.xlim([min(t),max(t)])
plt.plot(t,y3,color ='orange', label='y1(t-t0)')
plt.xlabel('t')
plt.legend()
plt.show()