convolução matlab

Prévia do material em texto

Análise de Sinais – 2007/2008 
António Espírito Santo 1
 
Análise de Sinais 
(5001) 
Ficha Prática Laboratorial VII 
Curso: Eng. Electrotécnica 
Duração prevista: 2 horas 
 
Convolução de Sinais 
I. Sinais analógicos 
Operações matemáticas realizadas em sinais analógicos exigem que os sinais sejam definidos no 
mesmo intervalo em que a operação é realizada. Algumas das características dos sinais são referidos de 
seguida. 
Função impulso de Dirac 
Esta função é caracterizada pela sua duração nula, amplitude infinita, mas de área unitária. 
Matematicamente podem ser definidas as seguintes propriedades. 
1. ( ) 0 0t t∂ = ⇐ ≠ 
2. ( ) 1t dt
+∞
−∞
∂ =∫ 
A função de Dirac pode ser usada sempre que se pretenda formular a aquisição de um sinal, usando 
para isso a expressão 
( ) ( ) ( )x t x t dτ τ τ
+∞
−∞
= ∂ −∫ 
É assumido que a função x(t) é contínua e finita. 
Função degrau unitário 
A função é nula até ao momento t(0), nesse ponto apresenta uma descontinuidade e toma o valor 
unitário. 
1 0
( )
0 0
se t
u t
se t
≥⎧= ⎨ <⎩ 
Análise de Sinais – 2007/2008 
António Espírito Santo 2
Funções trigonométricas 
Este tipo de funções é central para a análise de sinais. São definidas pela frequência angular 2 fω π= , 
onde f é a frequência do sinal. 
1. função coseno: ( ) cos( )x t tω= 
2. função seno: ( ) sin( )x t tω= 
Função exponencial complexa 
Esta função pode ser escrita como a soma de duas partes: uma real; e a outra imaginária. 
( ) cos( ) sin( )i tx t e t i tω ω ω= = + 
A relação de Euler pode aqui ser introduzida com vantagens 
1. cos( )
2
i t i te et
ω ω
ω
−+= 
2. cos( )
2
i t i te et
i
ω ω
ω
−−= 
II. Sinais Discretos 
Regra geral, um sinal discreto é uma sequência de números reais ou complexos. Sendo cada elemento 
da sequência identificado por um índice. Por exemplo, xn-2, xn-1, xn, xn+1, xn+2. Podemos facilmente definir 
operações matemáticas nesses sinais. Assim, a soma, a multiplicação por uma constante ou por outro sinal, 
e a operação de atraso são conseguidas usando a identificação dos índices. 
Dadas as sequências [xn ] = [0.5 2.4 3.2 4.5] e [yn] = [0.0 2.2 7.2 5.5]. 
Podemos definir a adição dos sinais discretos como sendo 
[zn] = [xn] + [yn] = [0.5 4.6 10.4 10.0] 
A multiplicação por uma constante c = 2 resultará em: 
[wn] = 2 *[xn] = [1.0 4.8 6.4 9.0] 
Impulso de dirac 
O impulso unitário é sempre nulo com a excepção de na origem dos tempos ser unitário. Matematicamente 
pode ser definido por 
1 0
0 0n
n
n
⇐ =⎧∂ = ⎨ ⇐ ≠⎩ com ] [,n∈ −∞ +∞ 
A função de impulso é muito útil e pode ser usada como o elemento constituinte de qualquer sinal 
discreto. Seja xn uma sequência, podemos definir o sinal xn como se segue. 
 
Análise de Sinais – 2007/2008 
António Espírito Santo 3
n k n k
k
x x
+∞
−
=−∞
= ∂∑ 
Degrau unitário 
Este sinal é nulo para valores de índice negativo, e unitário para valores de índice não nulo. 
1 0
0 0n
n
u
n
⇐ ≥⎧= ⎨ ⇐ <⎩ com ] [,n∈ −∞ +∞ 
Funções trigonométricas 
A função trigonométrica discreta é definida como se segue. Sendo n o índice da sequência e ω a 
frequência angular. 
1. função coseno: cos( )nx nω= 
2. função seno: sin( )nx nω= 
Funções exponencial complexa 
Esta função pode ser escrita como a soma de duas partes: uma real; e a outra imaginária. 
( ) i nx t e ω= 
A relação de Euler pode aqui ser introduzida com vantagens 
1. cos( )
2
i n i ne en
ω ω
ω
−+= 
2. cos( )
2
i n i ne en
i
ω ω
ω
−−= 
III. Convolução de sinais contínuos 
A operação de convolução é definida em sistemas LTI (linear and Time-Invariant), ou seja, lineares e 
invariantes no tempo. Tendo presente o facto anterior podemos estabelecer que: 
- A resposta de um sistema a uma entrada a impulso possui também uma resposta na forma de impulso; 
 
- A resposta de um impulso deslocado no tempo leva a uma resposta também deslocada no tempo, 
devido à característica do sistema ser invariante no tempo; 
Análise de Sinais – 2007/2008 
António Espírito Santo 4
 
- Se o impulso for afectado de um factor de escala, a resposta será também afectada desse mesmo facto 
de escala, justificada pela linearidade da multiplicação por um escalar; 
 
- A soma de um número de entradas de impulsos afectados de um factor de escala é igual à soma das 
respostas aos impulsos afectadas dos mesmos factores de escala, já que a operação de multiplicação é 
linear; 
 
Podemos reconhecer que a soma anterior corresponde a um integral, ou seja podemos converter os dois 
lados, entrada e saídas, em dois integrais. Reconhecendo que a entrada é a função f(t), também 
reconhecemos que a saída é o integral de convolução. 
Como já foi referido, o integral de convolução fornece um método matemático de expressar a saída de 
um sistema LTI baseada num sinal arbitrário de entrada x(t), e na sua resposta a impulso h(t). O integral de 
convolução é expresso por 
( ) ( ) ( )y t x h t dτ τ τ+∞−∞= −∫ 
A operação de convolução pode ser alternativamente representado pelo operador * 
( ) ( )* ( )y t x t y t= 
Esta operação matemática possui as propriedades: 
1. Associativa - [ ] [ ]1 2 3 1 2 3( )* ( ) * ( ) ( )* ( )* ( )f t f t f t f t f t f t= 
2. Comutativa - 1 2 2 1( )* ( ) ( )* ( )f t f t f t f t= 
3. Distributiva - [ ]1 2 3 1 2 1 3( )* ( ) ( ) ( )* ( ) ( )* ( )f t f t f t f t f t f t f t+ = + 
4. Deslocamento - 1 2 1 2( )* ( ) ( )* ( )f t f t f t f tτ τ− = − 
5. Convolução com impulso - 1 1( )* ( ) ( )f t t f t∂ = 
6. Duração do resultado de convolução – se duração de 1( )f t é T1 e a duração de 2 ( )f t é T2. Então a 
duração do resultado da convolução é T1+T2. 
Análise de Sinais – 2007/2008 
António Espírito Santo 5
7. Casualidade – Se 1( )f t e 2 ( )f t são ambos sinais casuais então o resultado da convolução é 
também casual. 
IV. Convolução de sinais discretos 
O conceito da convolução de sinais discretos é em tudo idêntico ao dos sinais contínuos. Nesta situação, o 
integral é substituído por um somatório expresso por 
[ ] [ ] [ ]
k
y n x k h n k
+∞
=−∞
= −∑ 
Tal como no caso da convolução contínua, é usado o operador *, de onde resulta a notação 
[ ] [ ]* [ ]y n x n y n= 
A operação de convolução contínua observa o mesmo conjunto de propriedades que a convolução 
contínua. Havendo a registar que agora a duração do resultado de convolução, no caso de 1( )f n possuir 
N1 amostras e 2 ( )f n possuir N2 amostras, ser igual é N1+N2 -1. 
V. Exemplo de convolução contínua 
‰ Exemplo I 
Um sistema LTI é descrito por uma resposta a impulso dada por h(t) = exp(-t)u(t) (em que u(t) é um sinal 
em degrau). Uma entrada x(t) é aplicada ao sistema, esta entrada possui um valor de 0.6 no intervalo 
-1 < t < 0.5, um valor de 0.3 no intervalo 0.5 < t < 3, e um valor nulo no tempo restante (Fig. 1). 
 
Fig. 1 – Resposta a impulso h(t) e entrada x(t) Fig. 2 – Resposta a impulso h(t) e entrada x(t) 
invertida no tempo 
Para realizar a operação de convolução um dos sinais deve ser invertido no tempo. Neste caso optou-se 
por inverter no tempo o sinal x(t) (Fig. 2). 
Análise de Sinais – 2007/2008 
António Espírito Santo 6
O próximo passo é deslocar temporalmente o sinal de entrada. Como os gráficos são feitos em ordem à 
variável τ, o sinal de entrada invertido e deslocado é dado por x(t-τ). O valor de t é a quantidade da qual o 
sinal invertido x(-τ) é deslocado no tempo. Correspondendo em simultâneo ao valor da saída y(t) fornecido 
por 
( ) ( ) ( )y t x t h dτ τ τ+∞−∞= −∫ 
Valor negativos de τ deslocam o sinal para a esquerda, enquanto que valores positivos deslocam o sinal 
para a direita. 
Como o valor da saída num instante particular de tempo t corresponde à área do produtode h(τ) e 
x(t-τ), os limites de integração correspondem às definições dos sinais. Como x(t) está definido em 
segmentos de tempo, o integral pode ser fraccionado nesses pequenos intervalos de tempo. 
O primeiro segmento corresponde a um deslocamento suficientemente para a esquerda para que não 
ocorra qualquer sobreposição dos sinais (Fig. 3), e logo o resultado da multiplicação é nula e também o 
valor de y(t). O intervalo de tempo no qual a saída do sistema apresenta um valor nulo pode ser obtido 
observando os sinais. A entrada x(τ) é nula para τ < -1, então x(-τ) é nulo para τ > 1. A transição do sinal 
 x(-τ) de 0.6 para 0 ocorre sempre para valores de τ = t+1. De igual forma, a transição de 0.3 para 0.6 
ocorre em τ = t-0.5, e a transição de 0 para 0.3 ocorre em τ = t-3. Podemos então concluir que o intervalo 
de tempo no decorrer do qual não existe sobreposição corresponde a τ <= 0, ou seja, para t+1 < 0. Logo 
y(t) = 0 para valores de t < -1. 
 
Fig. 3 – Entrada deslocada para a esquerda sem 
sobreposição 
Fig. 4 – Sobreposição do primeiro segmento 
Para valores de t > -1 passa a haver sobreposição (Fig. 4). Nesta primeira situação haverá sobreposição 
com x(t) = 0.6. Isto ocorre para t+1 > 0 e para t-0.5 < 0. Os limites de integração serão o intervalo de 
sobreposição. O que corresponde o intervalo de τ = 0 até τ = t +1. O intervalo de tempo para definir y(t) é 
dado pelo integral de -1< t < 0.5, o que corresponde a saída y(t) = 0.6*[1 - 0.3679*exp(-t)]. 
Análise de Sinais – 2007/2008 
António Espírito Santo 7
 
Fig. 5 – Sobreposição do segundo segmento Fig. 6 – Sobreposição de ambos os segmentos 
A próxima situação corresponde a x(t-τ) completamente sobreposto e x(t-τ) = 0.3 parcialmente 
sobreposto a h(τ) (Fig. 5). Isto ocorre quando t-0.5 > 0 e t-3 < 0. Os limites de integração são de 0 a t-0.5 
com x(t-τ) = 0.3 e de t-0.5 a t+1 com x(t-τ) = 0.6. A expressão obtida para y(t) neste intervalo de tempo é 
y(t) = 0.3*[1 + 0.9130*exp(-t)]. Sendo válida para o intervalo de tempo 0.5 < t < 3. 
A última situação corresponde a ter uma sobreposição completa de todos os segmentos (Fig. 6). Isto 
ocorre para t-3 > 0. Uma vez que h(t) está definida e é não nulo para todos os valores finitos de tempo, o 
tipo de sobreposição não varia. Os limites de integração são de t-3 até t-0.5 para x(t-τ) = 0.3 e para t-05 até 
t+1 para x(t-τ) = 0.6. O intervalo é definido para t > 3, e a expressão resultante para a saída 
y(t)= 6.2995*exp(-t). 
 
Fig. 7 – Saída y(t) 
O problema está resolvido usando o código em Matlab listado no anexo I. 
‰ Exemplo II 
A resposta a impulso h(t) e o sinal de entrada x(t) para um sistema LTI são fornecidas em seguida. 
Ambos os sinais são casuais já que são nulos para tempos negativos. 
 
O sinal de saída pode ser obtido efectuando a convolução analítica dos dois sinais. O que devolveria para a 
saída y(t) 
Análise de Sinais – 2007/2008 
António Espírito Santo 8
 
 
Este resultado pode também ser obtido usando a função conv do MatLab. Mas devemos ter em atenção 
dois factores. 
O primeiro diz respeito a um factor de escala que devemos usar. Esta função aproxima o integral de 
convolução usando um somatório. Assumindo que os dados dos vectores x(t) e y(t) são amostras dos 
sinais no domínio contínuo, com um período de amostragem de dt segundo, o resultado obtido ao usar a 
função conv deve, neste caso, ser multiplicado por dt = 0.01 segundos. 
O segundo factor a que devemos prestar atenção diz respeito ao facto de que, quer o sinal de entrada 
x(t), quer o sinal de saída h(t), são representados por vectores de comprimento N. A aplicação da função 
conv devolve um vector de resultados com 2N-1 elementos, espaçados no tempo em dt segundos. A 
interpretação dos resultados deve ser a seguinte. 
- Os primeiros elementos correspondem ao resultado correcto, já que nesse intervalo de tempo 
possuímos valores que caracterizam a entrada e a resposta a impulso. 
- Os restantes N-1 elementos são calculados na assumpção de que ambos os sinais, x(t) e h(t), são 
nulos nos instante de tempo fora do intervalo representado pelos vectores. 
No nosso exemplo, ambos os sinais existem para os tempos t > 0, mas no MatLab apenas são 
representados nos primeiros 10 segundos. Como x(t) e h(t) não são nulos para além dos 10 segundos, 
apenas podemos considerar os primeiros N elementos do resultado vectorial devolvido pela função conv. 
 
 
Fig. 8 – Entrada x(t), resposta a impulso h(t), e 
saída y(t) 
Fig. 9 – Resultado da multiplicação para t = 2.6 seg 
O código do MatLab que resolve este problema pode ser encontrado no anexo II. 
VI. Exemplo de convolução discreta 
Uma aproximação gráfica da operação discreta pode ajudar a compreender a operação. Para isso 
vamos assumir que o sinal de entrada x[n] é casual e de comprimento m, e que h[n] é causal e de 
comprimento k. Sendo o sistema LTI, podemos escrever a soma finita como 
Análise de Sinais – 2007/2008 
António Espírito Santo 9
1
0
[ ] [ ] [ ]
m
l
y n x l h n l
−
=
= −∑ 
De observar que para qualquer índice de amostra n, temos uma soma do produto da amostra x[l] com a 
amostra h[l-]. Ou seja, isto corresponde a multiplicar as amostras de x[n] pelas amostras de h[n] invertidas 
no tempo seguida da soma de todos os termos. O resultado final a obter pode ser representado 
graficamente por 
 
O processo de cálculo, passo a passo, pode ser representado graficamente através dos passos 
seguintes. 
 
 
Análise de Sinais – 2007/2008 
António Espírito Santo 10
 
 
Este processo pode ser observado nas applets situadas em http://www.jhu.edu/~signals/ 
 
 
Análise de Sinais – 2007/2008 
António Espírito Santo 11
Anexo I 
 
% ***** MATLAB Code Starts Here ***** 
% 
%convolution_01_mat 
% 
fig_size = [50 50 800 600]; 
th1 = linspace(0,10,1001); th = [0 th1]; % Impulse response h(t) 
h1 = exp(-th1); h = [0 h1]; 
% 
tx = [-1 -1 0.5 0.5 3 3]; % Input signal x(t) 
x = [0 0.6 0.6 0.3 0.3 0]; 
% 
figure(1),plot(th,h,tx,x),grid,xlabel('Time (s)'),ylabel('Amplitude'),... 
title('Impulse Response h(t) and Input Signal x(t)'),... 
axis([-5 10 0 1]),text(0.5,0.85,'h(t)'),text(3.2,0.25,'x(t)'),set(gcf,'Position',fig_size) 
% 
figure(2),plot(th,h,-tx,x),grid,xlabel('Time (s)'),ylabel('Amplitude'),... 
title('Impulse Response h(t) and Time-Reversed Input Signal x(-t)'),... 
axis([-5 10 0 1]),text(0.5,0.85,'h(t)'),text(-2.5,0.32,'x(-t)'),set(gcf,'Position',fig_size) 
% 
% Shifting x(t-tau) left 2 seconds (t = -2) 
% 
figure(3),plot(th,h,-tx-2,x),grid,xlabel('Tau (s)'),ylabel('Amplitude'),... 
title('Impulse Response h(tau) and Reversed & Shifted Input Signal x(t-tau)'),... 
axis([-5 10 0 1]),text(0.5,0.85,'h(tau)'),text(-2.5,0.62,'x(t-tau)'),text(6,0.65,'t = -2 s'),... 
set(gcf,'Position',fig_size) 
% 
% Shifing x(t-tau) right 1.5 seconds (t = -0.5) 
% 
figure(4),plot(th,h,-tx-0.5,x),grid,xlabel('Tau (s)'),ylabel('Amplitude'),... 
title('Impulse Response h(tau) and Reversed & Shifted Input Signal x(t-tau)'),... 
axis([-5 10 0 1]),text(0.5,0.85,'h(tau)'),text(-3,0.32,'x(t-tau)'),text(6,0.65,'t = -0.5 s'),... 
set(gcf,'Position',fig_size) 
% 
% Shifting x(t-tau) right 2.5 seconds (t = 2) 
% 
figure(5),plot(th,h,-tx+2,x),grid,xlabel('Tau (s)'),ylabel('Amplitude'),... 
title('Impulse Response h(tau) and Reversed & Shifted Input Signal x(t-tau)'),... 
axis([-5 10 0 1]),text(0.5,0.85,'h(tau)'),text(1.5,0.62,'x(t-tau)'),text(6,0.65,'t = +2 s'),... 
set(gcf,'Position',fig_size) 
% 
% Shifting x(t-tau) right 2 seconds (t = 4) 
% 
figure(6),plot(th,h,-tx+4,x),grid,xlabel('Tau (s)'),ylabel('Amplitude'),... 
title('Impulse Response h(tau) and Reversed & Shifted Input Signal x(t-tau)'),... 
axis([-5 10 0 1]),text(0.5,0.85,'h(tau)'),text(3.5,0.62,'x(t-tau)'),text(6,0.65,'t = +4 s'),... 
set(gcf,'Position',fig_size)% 
ty = [-1:0.01:0.5]; % 1st segment of y(t) 
y = 0.6*[1 - 0.3679*exp(-ty)]; 
% 
ty1 = [0.5:0.01:3]; % 2nd segment of y(t) 
y1 = 0.3*[1 + 0.9130*exp(-ty1)]; 
% 
ty2 = [3:0.01:10]; % 3rd segment of y(t) 
y2 = 6.2995*exp(-ty2); 
% 
ty = [-5 ty ty1 ty2]; % Combining the segments 
y = [0 y y1 y2]; 
Análise de Sinais – 2007/2008 
António Espírito Santo 12
% 
figure(7),plot(ty,y),grid,xlabel('Time (s)'),ylabel('Amplitude'),... 
title('Output Signal y(t)'),axis([-5 10 0 0.5]),set(gcf,'Position',fig_size) 
% 
% 
% ***** MATLAB Code Stops Here ***** 
 
Análise de Sinais – 2007/2008 
António Espírito Santo 13
 
Anexo II 
 
% ***** MATLAB Code Starts Here ***** 
% 
%CONVOLUTION_02_MAT 
% 
fig_size = [232 84 774 624]; 
t = 0:0.01:10; % time vector 
nt = length(t); 
dt = t(2); 
h = 2*t.*exp(-t) + exp(-2*t) - exp(-3*t); % impulse response 
x = 1 - exp(-1.5*t); % input signal 
% 
figure(1),clf,plot(t,x,t,h),grid,xlabel('Time (s)'),ylabel('Amplitude'),... 
title('Input x(t) and Impulse Response h(t)'),... 
text(2.5,0.95,'x(t)'),text(2.5,0.45,'h(t)'),set(gcf,'Position',fig_size) 
% 
figure(2),clf,plot(t,x,2.5-t,h),grid,xlabel('Time Tau (s)'),ylabel('Amplitude'),... 
title('Input x(tau) and Impulse Response h(t-tau), t = 2.5 s'),... 
text(2.5,0.95,'x(tau)'),text(2.25,0.65,'h(t-tau)'),set(gcf,'Position',fig_size) 
% 
% Calculation of the product of x(tau) & h(t-tau) for t = 2.5 s 
% 
tau = t(1:251); 
tmtau = 2.5 - tau; 
hx = (1 - exp(-1.5*tau)) .* ( 2*tmtau.*exp(-tmtau) + exp(-2*tmtau) - exp(-3*tmtau) ); 
% 
figure(3),clf,plot(t(1:251),hx),grid,xlabel('Time Tau (s)'),ylabel('Amplitude'),... 
title('Product of Input x(tau) and Impulse Response h(t-tau), t = 2.5 s'),... 
text(0.75,0.25,'Value of y(t) at t = 2.5 s equals the area under this curve'),... 
set(gcf,'Position',fig_size) 
% 
[numh,denh] = residue([0;2;1;-1;],[-1;-1;-2;-3],0); % calculation of transfer function 
[numx,denx] = residue([1; -1],[0; -1.5],0); % calculation of input transform 
[numy,deny] = series(numx,denx,numh,denh); % calculation of output transform 
% 
[resy,poly,ky] = residue(numy,deny); % partial fraction expansion for output 
% 
y1 = lsim(numh,denh,x,t); % calculation of output using "lsim" 
% 
figure(4),clf,plot(t,y1,t,x,t,h,t(251),y1(251),'r*'),grid,xlabel('Time 
(s)'),ylabel('Amplitude'),... 
title('Input x(t),Impulse Response h(t), and Output y(t)'),... 
text(2.5,0.9,'x(t)'),text(3.5,0.3,'h(t)'),text(3.5,1.6,'y(t)'),... 
text(2.7,y1(251),sprintf('y(t=2.5) = %g',y1(251))),set(gcf,'Position',fig_size) 
% 
y2 = conv(x,h) * dt; % calculation of output using "conv" 
y2 = y2(1:nt); 
% 
y3 = 0; 
figure(5),clf,plot(t,y1,'g-',t,y2,'r-.'),grid,xlabel('Time (s)'),ylabel('Amplitude'),... 
title('Output y(t) from LSIM and from CONV'),set(gcf,'Position',fig_size) 
% 
w = logspace(-2,3,501); % frequency vector 
% 
[magx,phx] = bode(numx,denx,w); % frequency domain data 
[magh,phh] = bode(numh,denh,w); 
[magy,phy] = bode(numy,deny,w); 
% 
Análise de Sinais – 2007/2008 
António Espírito Santo 14
figure(6),clf,subplot(211),semilogx(w,20*log10(magx),w,20*log10(magh),... 
w,20*log10(magy)),grid,xlabel('Frequency (r/s)'),ylabel('Magnitude (db)'),... 
title('Frequency-Domain Magnitudes of Input Signal, System, and Output Signal'),... 
subplot(212),semilogx(w,phx,w,phh,w,phy),grid,xlabel('Frequency (r/s)'),ylabel('Phase (deg)'),... 
title('Frequency-Domain Phases of Input Signal, System, and Output Signal'),... 
set(gcf,'Position',fig_size) 
% 
% ***** MATLAB Code Stops Here *****