PSD29007-Engtelecom(2019-2) - Prof. Marcos Moecke

De MediaWiki do Campus São José
Ir para: navegação, pesquisa

Registro on-line das aulas

Unidade 1

Unidade 1
Aula 1 (29 jul)
Aula 2 e 3 (1 e 5 ago)
Aula 4 (8 ago)
  • Revisão de Sinais e Sistemas no tempo discreto em Matlab:
  • Explorar a interface do Matlab.
  • Funções de visualização das variáveis no workspace.
  • Execução de instruções passo a passo.
  • Escrita de script .m
  • Uso da execução das seções de um script.
  • Incremento de valor e execução.
EXEMPLOS:
  • Leia com atenção e execute o exemplo (Moving-Avarage Filter) na página de help da função filter.
  • Revisão de Sinais e Sistemas no tempo discreto em Matlab:
  • Leia com atenção o help Using FFT, abra o script clicando no botão [Open this Example]. Execute o script seção após seção. Note o uso da fft para determinar a frequência das manchas solares.
  • Para melhorar o desempenho no Matlab recomendo que leiam a pagina do Help, . Também gostaria de lembra-los que a tecla F9 executa o código destacado no Help. Programação com scripts .m.
  • Leia sobre manchas solares para entender o que são os dados do segundo exemplo.
Sinais no dominio do tempo e dominio da frequencia. Uso da função fft
Exemplo de uso da FFT
%% Signal in Time Domain 
% Use Fourier transforms to find the frequency components of a signal buried in noise.
% Specify the parameters of a signal with a sampling frequency of 1 kHz and a signal duration of 1.5 seconds
Fs = 1000;            % Sampling frequency                    
T = 1/Fs;             % Sampling period       
L = 1500;             % Length of signal
t = (0:L-1)*T;        % Time vector

% Form a signal containing a 50 Hz sinusoid of amplitude 0.7 and a 120 Hz sinusoid of amplitude 1.
S = 0.7*sin(2*pi*50*t) + sin(2*pi*120*t);

% Corrupt the signal with zero-mean white noise with a variance of 4.
X = S + 2*randn(size(t));

% Plot the noisy signal in the time domain. It is difficult to identify the frequency components by looking at the signal X(t).
subplot(211);
plot(1000*t(1:200),X(1:200))
title('Signal Corrupted with Zero-Mean Random Noise')
xlabel('t (milliseconds)')
ylabel('X(t)')

%% Signal in Frequency Domain
% Compute the Fourier transform of the signal.
Y = fft(X);

% Compute the two-sided spectrum P2. Then compute the single-sided spectrum P1 based on P2 and the even-valued signal length L.
P2 = abs(Y/L);
P1 = P2(1:L/2+1);
P1(2:end-1) = 2*P1(2:end-1);

% Define the frequency domain f and plot the single-sided amplitude spectrum P1. 
% The amplitudes are not exactly at 0.7 and 1, as expected, because of the added noise. 
% On average, longer signals produce better frequency approximations.
f = Fs*(0:(L/2))/L;
subplot(212);
plot(f,P1)
ylim([0 1.05]) 
title('Single-Sided Amplitude Spectrum of X(t)')
xlabel('f (Hz)')
ylabel('|P1(f)|')

% Now, take the Fourier transform of the original, uncorrupted signal and retrieve the exact amplitudes, 0.7 and 1.0.
Y = fft(S);
P2 = abs(Y/L);
P1 = P2(1:L/2+1);
P1(2:end-1) = 2*P1(2:end-1);

plot(f,P1) 
title('Single-Sided Amplitude Spectrum of S(t)')
xlabel('f (Hz)')
ylabel('|P1(f)|')
  • Amostragem de Sinais (Experimento 1.2)
  • Relembrar teorema da amostragem. Efeito da amostragem abaixo da frequência de Nyquist. Aliasing.
  • Notar que as amostras de um sinal (3 Hz) e um sinal (7 Hz) são idênticas quando amostrado com um sinal de 10 Hz.
Experimento 1.2
%  Exemplos e Experimentos baseados no livro:
% DINIZ, P. S. R., DA SILVA, E. A. B., e LIMA NETTO, S. Processamento Digital de Sinais: Projeto e Análise de Sistemas. 2. ed. Porto Alegre: Bookman, 2014. 976 p. ISBN 978-8582601235.
%% Experimento 1.2
fs = 10; % frequencia (Hz) de amostragem dos sinais
Ts = 1/fs; fase = 0;
time = 0:Ts:(1-Ts);
f1 = 3; % frequencia (Hz) do sinal s_1
f2 = 7; % frequencia (Hz) do sinal s_2
s_1 = cos(2*pi*f1*time+fase);
s_2 = cos(2*pi*f2*time+fase);
fsa = 1000; % frequência auxiliar de amostragem usada apenas para representação dos sinais originais
Tsa = 1/fsa;
time_aux = 0:Tsa:(1-Tsa);
figure(1);
stem(time,s_1,'ob');
hold on;
plot(time_aux, cos(2*pi*f1*time_aux+fase),'--k');
stem(time,s_2,'+r');
plot(time_aux, cos(2*pi*f2*time_aux+fase),'--m');
hold off;
legend('s_1 discreto','s_1 contínuo','s_2 discreto','s_2 contínuo')
DICAS:
  • No help on-line da Mathworks, usando o botão [Try This Example > Try in your browser], permite executar o código no próprio browser sem ter nenhuma instalação do Matlab. Para verificar que o código realmente é executado mude a amplitude do ruído randômico para 0.1 ou 0.5, insira o comando close all antes da primeira linha, e execute todo o código [Run All]
  • No help do Matlab, usando o botão [Open this Example], é possível executar o código seção a seção.
Aula 5 (12 ago)
  • Revisão de Sinais e Sistemas no tempo discreto em Matlab:
Variação do Experimento 2.2
%  Exemplos e Experimentos baseados no livro:
% DINIZ, P. S. R., DA SILVA, E. A. B., e LIMA NETTO, S. Processamento Digital de Sinais: Projeto e Análise de Sistemas. 2. ed. Porto Alegre: Bookman, 2014. 976 p. ISBN 978-8582601235.
%% Experimento 2.2
% Resposta em frequencia usando a função freqz
N = 1;
num = [1 0 0 0];
den = poly([0.8 0.2])
%den = [1 0.6 -0.16];
% modo 1
%[H,w]=freqz(num,den,[0:pi/100:N*pi-pi/100]);
%plot(w/pi, abs(H));
% modo 2
%[H,w]=freqz(num,den);
%plot(w/pi, abs(H));
% modo 3
%[H,w]=freqz(num, den, 'whole');
%plot(w/pi, abs(H));
% modo 4
freqz(num, den, 'whole');
figure(2);
zplane(num,den);

%% Resposta em frequencia substituindo z -> e^(jw)
syms z
Hf(z) = symfun(z^2/(z-0.2)/(z+0.8),z);
pretty(Hf)
latex(Hf)
N = 1;
w = [0:pi/100:N*pi-pi/100];
plot(w/pi,abs(Hf(exp(1i*w))))
%title(['$' latex(Hf) '$'],'interpreter','latex')
text(0.2,2,['H(z) = ' '$$' latex(Hf) '$$'],'interpreter','latex')
xlabel(['w/' '$$' '\pi' '$$'],'interpreter','latex')
  1. Verifique a diferença entre os tipos de plots comentados no código.
  2. substitua o denominador de H(z) por dois polos em [-0.8 -0.8].
  3. verifique o que ocorre se forem utilizados polos complexos conjugados [0.3-0.4i 0.3+0.4i 0.1]
  4. verifique o que ocorre se forem utilizados polos complexos não conjugados [0.3-0.4i 0.3+0.8i]
  5. verifique o que ocorre se os polos estiverem fora do circulo unitário [1.2 -0.2]. Interprete este resultado

Aula 1 EAD

Aula 6 (15 ago EAD)
  • Revisão de Sinais e Sistemas no tempo discreto em Matlab:
  • Filtros Digitais
ATIVIDADE 1

Execute no Matlab o código abaixo, e analise os 3 filtros implementados através dos seus zeros e polos. Busque tirar conclusões sobre a influência da posição dos polos e zeros (ver o gráfico do plano z) e correlacione com a resposta de frequência em magnitude (gráfico do freqz).

%% Experimento 2.3 - Filtros Digitais
% Exemplos e Experimentos baseados no livro:
% DINIZ, P. S. R., DA SILVA, E. A. B., e LIMA NETTO, S. Processamento Digital de Sinais: Projeto e Análise de Sistemas. 2. ed. Porto Alegre: Bookman, 2014. 976 p. ISBN 978-8582601235.
% FILE: Exp2_3.m
 
%% 1º filtro
p1 = 0.9*exp(1j*pi/4);
Z = [1 -1 ]'; P = [p1 p1']';
[num,den] = zp2tf(Z,P,1);
[h,w] = freqz(num,den);
figure(1); plot(w,abs(h)/max(abs(h)));
figure(2); zplane(num,den);
 
%% 2º filtro
z1 = exp(1j*pi/8);
z2 = exp(1j*3*pi/8);
p1 = 0.9*exp(1j*pi/4);
Z = [1 -1 z1 z1' z2 z2']';
P = [p1 p1' p1 p1' p1 p1']';
[num,den] = zp2tf(Z,P,1);
[h,w] = freqz(num,den);
figure(1); plot(w,abs(h)/max(abs(h)));
figure(2); zplane(num,den);
 
%% 3º filtro
z1 = exp(1j*pi/8);
z2 = exp(1j*3*pi/8);
p1 = 0.99*exp(1j*pi/4);
p2 = 0.9*exp(1j*pi/4 - 1j*pi/30);
p3 = 0.9*exp(1j*pi/4 + 1j*pi/30);
Z = [1 -1 z1 z1' z2 z2']';
P = [p1 p1' p2 p2' p3 p3']';
[num,den] = zp2tf(Z,P,1);
[h,w] = freqz(num,den);
figure(1); plot(w,abs(h)/max(abs(h)));
figure(2); zplane(num,den);


ATIVIDADE 2
  • A filtragem de sinais digitais pode ser realizada de diferentes formas:
  • convolução (y = conv(x,h)), onde x(n) é o sinal de entrada e h(n) é a resposta ao impulso do filtro (sistema linear invariante no tempo),
  • filtragem no domínio do tempo (y = a1.x(n)+ a2.x(n-1)+ .. ak.x(n-k));
  • no domínio da frequência (y = ifft(fft(x)fft(h))
Variação do Experimento 3.1
%% Variação do Experimento 3.1 do livro:
% DINIZ, P. S. R., DA SILVA, E. A. B., e LIMA NETTO, S. Processamento Digital de Sinais: Projeto e Análise de Sistemas. 2. ed. Porto Alegre: Bookman, 2014. 976 p. ISBN 978-8582601235.
% FILE: Ex3_1.m
% Exemplificando as possiveis formas de realizar a filtragem de um sinal x(n)

clc; clear all; close all;
%% Definindo valores iniciais
Nh = 10; Nx = 20;
%Nh = 400; Nx = 10000;
x = ones(1,Nx);
% A resposta ao impulso de um sistema h(n) 
% no filtro FIR aos coeficientes b(n) = h(n) 
h = [1:Nh]; b = h;
%% Filtrando o sinal e medindo tempos

% OPÇÃO 1 - Filtragem utilizando a convolução
% NOTE: length(y) = length(x) + length(h) -1

tic;  % iniciar a contagem do tempo
y1 = conv(x,h); 
t(1) = toc; % terminar a contagem e mostrar tempo no console

% OPÇÃO 2 - filtragem utilizando a equação recursiva
% NOTE: length(y) = length(x)

tic;
y2 = filter(b,1,x);
t(2) = toc;

% OPÇÃO 3 - filtragem utilizando a equação recursiva 
% aumentando o tamanho de x para que length(y3) = length(y1)
x3 = [x zeros(1,length(h)-1)];

tic;
y3 = filter(h,1,x3); 
t(3) = toc;

length_y = length(x) + length(h) - 1;

% OPÇÃO 4 - filtragem utilizando a FFT 
% a y = IFFT(FFT(x)*FFT(h))

tic;
X = fft(x,length_y);
H = fft(h,length_y);
Y4 = X.*H;
y4 = ifft(Y4);
t(4) = toc;

% OPÇÃO 5 - filtragem utilizando a função fftfilt
% a y = IFFT(FFT(x)*FFT(h))

tic
y5 = fftfilt(h,x3);
t(5) = toc;

disp('Comprimento do vetor de saída length(y)')
disp(['    ' num2str([length(y1) length(y2) length(y3) length(y4) length(y5)])])
disp('Tempo usado na filtragem em micro segundos')
disp(['    ' num2str(t*1e6) ' us'])

%%  Plotando o gráfico
subplot(411);stem(y1);
hold on;
stem(y2,'xr');
stem(y3,'+m');
legend('y1', 'y2', 'y3')
hold off
subplot(412);stem(y1, 'ob');legend('y1')
subplot(413);stem(y2, 'xr'); hold on; stem(zeros(size(y1)),'.w');hold off; legend('y2')
subplot(414);stem(y3, '+m');legend('y3')
  • Verificar as funções tic e toc
  • Notar a diferença de tempo de processamento entre os processos de filtragem.
  • A situação pode ser muito diferente conforme muda o tamanho do sinal e ordem do filtro (h(n)). Modifique a resposta ao impulso e o sinal de entrada modificando os valores das variáveis de tamanho: Nh = 10, 100, 1000; Nx = 20, 1000, 10000;
  • Em função do sistema operacional e reserva de memória para as variáveis é importante desprezar a primeira medida de tempo. Realize 3 medidas de tempo para cada uma das 5 opções de filtragem, com pelo menos duas combinações de comprimento Nh e Nx. Relate os resultados obtidos no comentário da atividade no Moodle.
Aula 7 (19 ago)
  • Exercício - Sinal DTMF com ruído
  • Verifique se o Matlab está reproduzindo corretamente o som.
%% Carregando o som
clear, close, clc
load handel;

%% Reproduzindo o som 
sound(y,Fs)
 
% Reproduzindo o som 
%soundsc(y,Fs)
 
% Reproduzindo o som 
%player = audioplayer(y, Fs);
%play(player);
  • Usando o Matlab (ou Audacity) para gerar um sinal DTMF correspondente a um número N e adicionar um ruido ao sinal. Opcionalmente utilize um sinal DTMF gravado
  • Utilizar uma frequência de amostragem de 8000Hz de fazer a duração do sinal igual a 2 segundos.
Sinal 1234567890*#
  • Para adicionar o ruído utilize a função y = awgn(x,snr), ou y = x + nivel*randn(n).
  • Observe este sinal no domínio do tempo (DT) e domínio da frequência (DF).
%% Carregando o som
clear, close, clc
[y,Fs] = audioread('DTMF_8kHz.ogg');

%% Reproduzindo o som 
sound(y,Fs)

%% Visualizando o som no DT
time = [0:length(y)-1]'/Fs;
plot(time',y'); xlabel('segundos');
xlim([0 time(end)]), ylim([-1 1]);

%% Visualizando o som no DF
Nfreq = length(y);
freq = linspace(0,2*pi,Nfreq)'*Fs/pi/2;
Y = fft(y,Nfreq)/Nfreq;
plot(freq,abs(Y)); xlabel('Hertz');
xlim([0 Fs/2]);
  • Utilizar no Audacity um sinal DTMF ("1234567890") com fa= 8kHz
  • Visualizar no domínio do tempo e frequência.
  • Realizar a filtragem passa-baixas com fc = 1066 Hz, (selecionar a maior atenuação permitida)
  • Realizar a filtragem passa-faixa com f0 = 770 Hz e B = 70Hz (selecionar a maior ordem permitida)
  • Repetir o procedimento anterior para um sinal de ruído branco.
  • Consulte a documentação do Matlab sobre
     fft, ifft, fftshift, randn
    
  • Consulte a documentação do Matlab sobre
     plot, grid, subplot, hold, xlabel, ylabel, title, legend, xlim, ylim, log10, log
    
  • Consulte a documentação do Matlab sobre text, zp2tf, tf2zp, fftfilt, awgn
  • Ver pag. 141 a 145 e 230 a 235 de [1]

Unidade 2

Unidade 2
Aula 8 (22 ago)
  • Filtros Analógicos:
  • Função de transferência
  • Resposta em frequência: para obter a resposta em frequência é necessário avaliar
  • O projeto de filtros analógicos é realizado em 2 etapas:
  1. projeto de um filtro passa baixas (LP) protótipo normalizado com frequência de passagem
  2. transformação em frequência para o tipo de filtro (LP, HP, BP ou BS)
  • Análise básica de filtros analógicos com Matlab.
Dado um sistema linear invariante no tempo, representado pela função de transferência , obter a resposta de frequência do sistema (Magnitude e Fase).
b = [1 1];
a = [1 1 5];
[z1,p1,k]=tf2zp(b,a)
z2 = roots(b);
p2 = roots(a);
zplane(b,a);
%%
freqs(b,a);
%%
syms s  w
H(s) = (s+1)/(s^2 + s + 5);
pretty(H(1j*w))
latex(H(1j*w))
%%
ws = logspace(-2, 1, 1000);
h = H(1j*ws);
subplot(211)
semilogx(ws,abs(h)); grid on;
subplot(212)
semilogx(ws,angle(h)/pi*180); grid on;
  • Projeto de filtros analógicos do tipo Butterworth
  • A aproximação de magnitude de filtros analógicos pode ser realizado usando as aproximações de Butterworth, Chebyshev (tipo 1 ou 2) e Cauer.

TiposFiltrosHs.png

Proposta de exercício
  • Use os polinômios de Butterworth com ordens de 1 a 10 mostrados na tabela abaixo para obter os filtros .
n Fatores Polinomiais de
1
2
3
4
5
6
7
8
9
10
  • Escolha uma ordem n (entre 5 e 10)
  • Plote a resposta em frequência em escala log da amplitude (em dB) e da fase (em rad/pi).
  • Qual é o ganho do filtro na banda passante?
  • Qual é a frequência de corte (-3dB) do filtro.
  • Qual é o salto de de fase que ocorre em algumas frequências?
  • Qual é o fator de atenuação em dB/decada após a frequência de corte?
  • Faça o diagrama de polos e zeros desse filtro.
  • Procure observar o que ocorre com a posição dos polos do filtro.
  • Calcule o valor do módulo dos pólos.


Aula 9, 10 (16, 26 ago)
  • Projeto de filtros analógicos passa baixas (low pass - LP) do tipo Butterworth, considerando: é a frequência de passagem, é a atenuação em dB na frequência de passagem, é a frequência de stopband, é a atenuação em dB na frequência de stopband.

MascaraFiltroLP.png

  • Escalando as frequências em relação a , teremos que , e são as frequências de passagem e stopband do filtro protótipo , que tem ganho unitário e frequência de passagem 1.
  • Se considerarmos o caso particular em que na frequência de passagem o ganho (em escala linear) deve ser , que corresponde a um ganho (em escala log) , ou atenuação .
  • Considere que , teremos
  • Para projetar o filtro é necessário:
1) determinar a ordem do filtro:
2) obter os polos do filtro:
3) obter a função de transferência:
, onde
  • No caso de um filtro LP é necessário ainda obter a função de transferência do filtro especificado fazendo a transformação de frequência


Para qualquer
  • Teremos
  • Para projetar o filtro é necessário:
1) determinar a ordem do filtro:
2) obter os polos do filtro:
3) obter a função de transferência:
, onde e .
NOTA: o valor também pode ser obtido a partir de , pois corresponde ao último termo do polinômio .
  • No caso de um filtro LP é necessário ainda obter a função de transferência do filtro especificado fazendo a transformação de frequência


Aula 11 e 12 (29 ago e 2 set)
  • Projeto de filtros analógicos do tipo Chebyshev I.
  • Polinômios de Chebyshev:

Os polinômios de Chebyshev de primeira ordem são definidos pela relação recursiva:

Os primeiros cinco polinômios de Chebyshev de primeira ordem são:

  • Determine a ordem mínima necessária considerando: é a frequência de passagem do filtro LP, é a atenuação em dB na frequência de passagem, é a frequência de stopband do filtro, é a atenuação em dB na frequência de stopband, , , são as frequências de passagem e stopband do filtro protótipo.
  • Em seguida obter os polos do filtro:
, onde
  • Para obter a função de transferência:
, onde
onde
se n é par
se n é impar
é o último termo do denominador D(p)
  • Projeto de filtros analógicos do tipo Butterworth, Chebyshev I e II e Cauer (eliptico) usando funções do Matlab.
%% Projeto de filtro passa-baixas usando funções do Matlab  
%% Especificações do filtro 
Wp =16000; Ws = 20000; Ap = 0.3; As = 20; G0= 3;
% Para analisar o filtro projetado, use fvtool(b,a) para observar plano s, resposta em magnitude, fase e atraso de grupo

%% Butterworth
[n,Wn] = buttord(Wp, Ws, Ap, As,'s')
[b,a] = butter(n,Wn, 's');

%% Chebyshev I
n = cheb1ord(Wp, Ws, Ap, As,'s')
[b,a] = cheby1(n,Ap, Wp, 's');

%% Chebyshev II
n = cheb2ord(Wp, Ws, Ap, As,'s')
[b,a] = cheby2(n,As, Ws, 's');

%% Elliptic - Cauer
[n, Wn] = ellipord(Wp, Ws, Ap, As,'s')
[b,a] = ellip(n,Ap,As, Wn, 's');
  • Ver pag. 204 a 208 de [2]


Aula 13 e 14 (5 e 9 set)
  • Filtros Analógicos:
  • Transformações de frequência de filtros analógicos
  • passa-baixas () -> passa-baixas ()
  • Cálculo do protótipo com
  • Substituição de variáveis
  • passa-baixas () -> passa-altas ()
  • Cálculo do protótipo com
  • Substituição de variáveis
  • passa-baixas () -> passa-faixa ( e )
  • Cálculo do protótipo com
  • Substituição de variáveis
onde e
  • passa-baixas () -> rejeita-faixa ( e )
  • Cálculo do protótipo com
  • Substituição de variáveis
onde e
  • Ver pag. 208 a 218 de [2]
  • Exemplos de Filtros Analógicos:
  • Exemplo 1: Filtro passa-baixas ( = 941Hz, = 1209 Hz, = 1 dB, = 20 dB)
  • Exemplo 2: Filtro passa-altas ( = 1209 Hz, = 941Hz, = 1 dB, = 20 dB)
  • Exemplo 3: Filtro passa-faixa ( = 811 Hz, = 895,5 Hz = 770 Hz, = 1209 Hz, = 1 dB, = 30 dB)
  • Exemplo 4: Filtro rejeita-faixa ( = 53 Hz, = 58 Hz, = 62 Hz = 67 Hz, = 2 dB, = 25 dB)
NOTA:
  • No calculo do filtro lembre-se de usar as frequências angulares para , , , .
  • onde () é a frequência de passagem em Hz (rad/s), () é a frequência de rejeição em Hz (rad/s), () é a frequência central em Hz (rad/s), () é a largura de banda em Hz (rad/s).
  • Confira os projetos dos filtros plotando as respostas em frequência dos filtros protótipo H(p) e filtro final H(s) de cada um dos exemplos.
Aula 15 (16 set)
  • Filtros Digitais: Filtros IIR: transformações do tempo contínuo no tempo discreto
  • Obter a especificação do filtro em angulo entre 0 e 1, onde 1 corresponde a metade da frequência de amostragem
  • Obter o valor desse angulo predistorcido para compensar a distorção na frequência causada pela transformação bilinear , onde
  • passa-baixas () -> passa-baixas ()
  • Substituição de variáveis
  • Cálculo do protótipo com
  • passa-baixas () -> passa-altas ()
  • Substituição de variáveis
  • Cálculo do protótipo com
  • passa-baixas () -> passa-faixa ( e )
  • Substituição de variáveis
  • Cálculo do protótipo com
onde e
  • passa-baixas () -> rejeita-faixa ( e )
  • Substituição de variáveis
  • Cálculo do protótipo com
onde e
  • Realizar os projetos dos exemplos anteriores, considerando uma frequência de amostragem de 8 kHz.
  • Ver pag. 219 a 229 de [2]
  • Ver pag. 403 a 415 e 434 a 435 de [1]
Aula 16 (19 set)
O projeto dos filtros digitais IIR baseados na transformada bilinear no Matlab é realizada em dois passos: (1) Determinação da ordem do filtro; (2) Determinação dos coeficientes do numerador e denominador de .
fa = 200;
fN = fa/2;
wo = 60/fN;  bw = 10/fN;
[b,a] = iirnotch(wo,bw);
fvtool(b,a);
syms z;
N(z) = poly2sym(b,z);
D(z) = poly2sym(a,z);
H(z) = N(z)/D(z);
pretty(vpa(H(z),3))
fa = 8000;
fN = fa/2;
wo = 941/fN;  bw = 100/fN;
[b,a] = iirpeak(wo,bw);
fvtool(b,a);
syms z;
N(z) = poly2sym(b,z);
D(z) = poly2sym(a,z);
H(z) = N(z)/D(z);
pretty(vpa(H(z),3))
fa = 8000; fN = fa/2;
fo = 1000;  bw = 20/fN;
[b,a] = iircomb(fa/fo,bw,'peak');  % ou use a flag 'notch'
fvtool(b,a);
syms z;
N(z) = poly2sym(b,z);
D(z) = poly2sym(a,z);
H(z) = N(z)/D(z);
pretty(vpa(H(z),3))

Unidade 3

Unidade 3 - Filtros FIR
Aula 16 (23 set)
  • Filtros Digitais: Filtros FIR
  • Filtros de fase linear: simétricos e antisimétricos (Tipo 1, 2, 3 e 4)
  • Filtros de fase linear: propriedades (respostas em frequência possíveis, distribuição dos zeros em simetria quadrantal)
  • Observar as propriedades dos FIR tipo 1, 2, 3 e 4. Observe a resposta de frequência, fase, atraso de grupo, coeficientes e diagrama de polos e zeros dos filtros abaixo.
N = 10;
bi = 2*(rand(1,N)-0.5)
%% Tipo I - LP, HP, BS, BP
b = [bi (2*rand(1,1)-0.5) flip(bi)];
...
%% Tipo II - LP, BP
% tem um zero em -1 
b = [bi  flip(bi)];
...

%% Tipo III - BP
% tem um zero em 1  e -1 
b = [bi  0 -flip(bi)];
...

%% Tipo IV - BP, HP
% tem um zero em 1   
b = [bi  -flip(bi)];
...

FIR tipo1.png
Figura 1 - Propriedades do filtro FIR de fase linear (Tipo 1)

Aula 17 e 18 (26 e 30 set)
  • Coeficientes da série de Fourier de filtros ideias: LP, HP, BP, BS
  • Passa-baixas (Low-pass)
  • Passa-altas (High-pass)
  • Passa-faixa (Band-pass)
  • Rejeita-banda (Band-stop)
  • Janela retangular, fenômeno de Gibbs
  • Uso de funções de janelamento temporal no projeto de filtros digitais.
  • Tipos de janelas temporais usadas no projeto de filtros digitais.
  • Retangular
  • Bartlett
  • Hanning
  • Hamming
  • Blackman
  • em todas as janelas quando
onde é para par e para impar
ver também apodization function
L = 64; 
wvtool(rectwin(L), triang(L), bartlett(L), hann(L), hamming(L), blackman(L), blackmanharris(L), nuttallwin(L));
Tabela 5.1
Janela
Retangular 13.3 20.33 0.92/M
Triangular 26.6 27.41
Bartlett 26.5 27.48
Hann 31.5 44.03 3.11/M
Bartlett-Hanning 35.9 40.77
Hamming 42.5 54.08 3.32/M
Bohman 46.0 51.84 7.01/M
Parzen 53.1 56.89
Blackman 58.1 75.25 5.56/M
Flat Top 88.0 106.3
Blackman-Harris 92.1 108.8
Nutfall 93.8 109.7
  • Dados acima obtidos para um filtro passa baixas de ordem N = 64 com
  • Projeto de filtro FIR utilizando janelas temporais fixas.
  • Exemplo de projeto
Projetar um filtro passa baixas usando uma janela temporal fixa (verificar a janela que atende a especificação)
wp = 0.2*pi; Ap = 0.2 dB; Gp = 0 dB
ws = 0.3*pi; As = 60 dB;
  • Informar qual o tipo de janela, a ordem obtida, e o valor de wc do projeto final
  • Exemplo de projeto
Projetar um filtro LP usando uma janela temporal fixa (hamming, bartlett-hanning, hanning).
wp = 0.4*pi; Ap = 1 dB; Gp = 0 dB
ws = 0.6*pi; As = 40 dB;
  • Comparar os 3 tipos de janela, a ordem obtida, e o valor de wc em cada projeto.
Use como uma estimativa inicial os valores da Tabela 5.1 pag. 268.
  • PASSO 1 - Escolher o tipo de janela de acordo com a atenuação do lóbulo lateral Asl e As.
  • PASSO 2 - Estimar a ordem N1 do filtro considerando os parâmetros Dw
  • PASSO 3 - Calcule os coeficientes clp do filtro LP , calcule os valores da janela w e obtenha a resposta ao impulso do filtro h = clp * w.
  • PASSO 4 - Ajuste o ganho do filtro para que a resposta na banda de passagem fique abaixo da especificação maxima.
  • PASSO 5 - Verifique o valor real de Dwr = wAs-wAp, e faça a correção da ordem do filtro em função do desvio constatado. N2 = N*Dwr/Dw.
  • PASSO 6 - Corrija o valor de projeto dos coeficientes Clp do filtro ideal, a janela e a resposta ao impulso.
  • Repita o PASSO 3 até 6, até obter um filtro que atenda as especificações de Dw.
  • PASSO 7 - Desloque a frequência de corte wc de modo a obter o valor correto de wp. wc2 = wp + (wp-wAp).
  • Projeto de filtro FIR.
  • Projete os dois filtros projetados anteriormente como IIR, utilizando 3 janelas diferentes. Compare os filtros obtidos com os filtros IIR.
N = <ordem>
h_fir = fir1(N,Wn,hamming(N+1));
[Hw,w] =freqz(h_fir);
plot(w/pi,20*log10(abs(Hw)))
title(['hamming N = ' num2str(N)])
%fvtool(h_fir,1)
  • Ver pag. 256 a 265 de [2]
  • Ver artigos:
Aula 18 (03 out)
  • Projetar os filtros LP, HP e BP de acordo com as especificações dadas.
Aula 19 (07 out)
  • Filtros Digitais: Filtros FIR
  • Projeto de filtro FIR utilizando janelas temporais ajustáveis
L = 64; 
r = 60;    % Chebyshev e Tukey
alpha = 3; % Gauss
betha = 8; % Kaiser
nbar = 10; % Taylor
wvtool(kaiser(L,betha), chebwin(L,r), gausswin(L,alpha),tukeywin(L,r), taylorwin(L,nbar,-r));

Para a janela de Kaiser, a estimação do fator e da ordem do filtro são obtidos por:

onde é a atenuação do lóbulo lateral e é a largura da banda de transição em rad/amostra.

A janela de Kaiser é definida por:

onde : é a função de Bessel de ordem zero [1]

Utilizando o Matlab é possível estimar esses valores utilizando a função kaiserord. Exemplo da obtenção de um filtro passa baixa com , , atenuação de 40 dB na "stopband"

fsamp = 8000;
fcuts = [1000 1500];
mags = [1 0];
devs = [0.01 0.01];
[n,Wn,beta,ftype] = kaiserord(fcuts,mags,devs,fsamp);

Com os parâmetros é possível projetar o filtro usando a função fir1, que utiliza o método da janela para o projeto do filtro.

h_fir = fir1(n,Wn,ftype,kaiser(n+1,beta),'noscale');
[Hw,w] =freqz(h_fir);
plot(w*fsamp/2/pi,20*log10(abs(Hw)))
title(['Kaiser filter N = ' num2str(n)])
%fvtool(h_fir,1)
  • Ver as funções fir1, kaiserord do Matlab.
  • Ver pag. 266 a 273 de [2]
  • Uso das funções window e fir1 do Matlab para projeto de filtro FIR


Aula 22 (10 out)
  • Filtros Digitais: Filtros FIR


Digital Filters with Linear Phase].

  • Exemplo do projeto de um filtro passa-baixas, com minima ordem (Filtro de Parks-McClellan) com frequência de passagem de 1000 Hz e frequência de rejeição de 1500 Hz, dada uma frequência de amostragem de 8000 Hz. Considere que a atenuação na banda de rejeição é de no mínimo 40 dB e o ripple máximo na banda passante é de 0.4 dB.
fa = 8000;

Ap = 0.4;          
Ar = 40;          
       
fp = 1000;
fr = 1500;

f = [fp fr];    
a = [1 0];        
dev = [(10^(Ap/20)-1)/(10^(Ap/20)+1)  10^(-Ar/20)]; 
[n,fo,ao,w] = firpmord(f,a,dev,fa);
b = firpm(n,fo,ao,w);
[h,w] = freqz(b,1,1024,fa);
plot(w, 20*log10(abs(h))); hold on;
plot([0 fr fr fa/2], [Ap/2 Ap/2 -Ar -Ar],':m')
plot([0 fp fp], [-Ap/2 -Ap/2 -(Ar+30)],':m');
ylim([-(Ar+30) Ap/2+10])


Unidade 4 - Realização de Filtros

Unidade 4 - Realização de Filtros
Aula 24 (17 out)
  • Realização de Filtros
  • Realização de filtros FIR: Forma Direta.
FIR FD MathWorks.png
Figura 1 - Realização de filtros FIR na Forma Direta
  • Realização de filtros FIR: Forma Transposta. A transposição consiste na inversão do fluxo de todos os sinais, substituição de nós de soma por derivações e as derivações por soma. A entrada e saída também devem ser invertidas. A realização da transposição não altera o sistema implementado.
FIR FDT MathWorks.png
Figura 2 - Realização de filtros FIR na Forma Transposta
FIR FDT2 MathWorks.png
Figura 3 - Realização de filtros FIR na Forma Transposta
  • Realização de filtros FIR de fase linear: simétrico tipo I e II e antissimétrico tipo III e IV.
FIR Sym2 MathWorks.png
Figura 4 - Realização de filtros FIR de fase linear Simétrico I
FIR Sym1 MathWorks.png
Figura 5 - Realização de filtros FIR de fase linear Simétrico II
FIR AntiSym3 MathWorks.png
Figura 6 - Realização de filtros FIR de fase linear Antisimétrico III
FIR AntiSym4 MathWorks.png
Figura 7 - Realização de filtros FIR de fase linear Antisimétrico IV
  • Realização de Filtros usando o comando realizemdl do MatLab
Fs = 30000;              % Sampling Frequency
Fpass = 12000;           % Passband Frequency
Fstop = 13000;           % Stopband Frequency
Dpass = 0.01;            % Passband Ripple
Dstop = 0.01;            % Stopband Attenuation
flag  = 'scale';         % Sampling Flag

% Calculate the order from the parameters using KAISERORD.
[N,Wn,BETA,TYPE] = kaiserord([Fpass Fstop]/(Fs/2), [1 0], [Dstop Dpass]);

% Calculate the coefficients using the FIR1 function.
b  = fir1(N, Wn, TYPE, kaiser(N+1, BETA), flag);

hFIR = dsp.FIRFilter;
hFIR.Numerator = b;

% Para definir diretamente os coeficientes
realizemdl(hFIR)

% Para definir os coeficientes através de uma matriz de entrada
realizemdl(hFIR,'MapCoeffsToPorts','on');
Aula 25 (21 out)
  • Uso do Simulink
Aula 27 (31 out)
  • Avaliação A1.


Aula 28 e 29 (04 e 07 nov)
  • Uso do Simulink
  • Uso do [4] para projeto de filtro IIR, FIR equiripple e FIR com janela.
  • Realizar a montagem do modelo indicado na figura abaixo e fazer a simulação, usar fa = 8000Hz.

SimulationFilterPSD1.png

DICAS:

  • Ative a visualização das dimensões das portas [Display > Signal & Ports > Signal Dimensions] ou [ALT]+D+S+D+[ENTER]
  • Ative a visualização do tipo de dados das portas [Display > Signal & Ports > Port Data Types] ou [ALT]+D+S+D+D+[ENTER]
  • Ative a visualização da cor para indicar os tempos de amostragem sinais [Display > Sample Time > Colors] ou [ALT]+D+T+C+[ENTER]
  • Fazer a análise no DT e DF para um ruido uniforme (-1 a 1)

DTnoiseFilterPSD1.png DFnoiseFilterPSD1.png

  • Fazer a análise no DT e DF para um sinal de chirp (0 a 4000Hz)

DTchirpFilterPSD1.png DFchirpFilterPSD1.png

  • Diferença entre processamento por amostra e processamento por quadro (Sample- and Frame-Based Concepts).
  • Uso dos módulos buffer e unbuffer para converter entre fluxos de dados de diferentes taxas de frames, ou entre sample-based e frame-based.
  • Exemplos:

Para configurar o Simulink para sistemas discretos execute o comando dspstartup.m antes de abrir um novo modelo.

Aula 30 (11 nov)
  • Realização de filtros IIR de 2ª ordem: Forma Direta I e II, e Forma Transposta I e II.
  • Separando H(z) em dois blocos , e obtendo o sinal intermediário W(z) ou Y(z) dependendo da ordem dos blocos.
H1 H2 MathWorks.png
Figura 7 - Separação do filtro IIR H(z) em H1(z) e H2(z)
Com o ordenamento dos blocos e em ordem direta teremos a Forma Direta I:
Podemos obter a realização de na forma direta.
Para obter a realização de , é necessário reescrever a saída em função de e das saídas anteriores e :
IIR FD1 MathWorks.png
Figura 8 - Realização de filtros IIR na Forma Direta I
Com o ordenamento dos blocos e em ordem reversa teremos a Forma Direta II:
IIR FD2a MathWorks.png
Figura 9 - Realização de filtros IIR na Forma Direta II
Considerando que os sinais no centro são idênticos podemos simplificar e obter a Forma Direta II (Canônica):
IIR FD2b MathWorks.png
Figura 10 - Realização de filtros IIR na Forma Direta II Canônica
Considerando as regras de transposição podemos obter a forma transposta I e II. A transposição consiste na inversão do fluxo de todos os sinais, substituição de nós de soma por derivações e as derivações por soma. A entrada e saída também devem ser invertidas. A realização da transposição não altera o sistema implementado.
IIR FT1 MathWorks.png
Figura 11 - Realização de filtros IIR na Forma Transposta I
IIR FT2 MathWorks.png
Figura 12 - Realização de filtros IIR na Forma Transposta II
  • Realização de filtros IIR de ordem maior que 2: Forma Direta I e II, Transposta I e II, Cascata, Paralela
  • Os filtros IIR de ordem superior a 2 podem ser implementados nas FD I ou II e na FT I ou II. No entanto nessa configuração tendem a ficar instáveis ao terem os coeficientes quantizados, e também terem uma significativa alteração da resposta em frequência. Para reduzir esses problemas uma possível solução é a decomposição em filtros de 2ª ordem para serem associados na forma em Cascata ou Paralela.

ATUAL

Aula 31 (14 nov)
  • Filtros Digitais: Ferramentas do Matlab para projeto
  • Realização de filtros FIR: Cascata, Polifase
  • Vantagens do uso de filtro Polifase:
1) Quando o sinal será subamostrado (downsampling) de "D" amostras após a filtragem, a complexidade da implementação é reduzida de "D" vezes, pois apenas uma das "fases" precisa ser implementada.
2) Para reduzir o harware a ser implementado, é possível implementar apenas uma das "fases" do filtro e trocar "D" vezes os coeficientes.
  • Filtros Digitais: Quantização
Ver Use Filter Designer with DSP System Toolbox Software
x=-0.2;
% Word length = 8, fraction length = 7
q=quantizer([8,7]);
xq=quantize(q,x);
binxq=num2bin(q,xq)
% Word length = 16, fraction length = 15
q1=quantizer([16 15]); 
xq1 = quantize(q1,x);
binxq1=num2bin(q1,xq1)
  • Sobre ponto flutuante [7]


Unidade 5 - Projeto Final

Unidade 5 - Projeto Final
Aula XX (18 nov)
PF - Projeto de uma Unidade de Resposta Audível (URA) (Entrega e prazos ver Moodle)
  • Projeto de uma Unidade de Resposta Audível (URA). O projeto é constituído de um receptor DTMF, que discrimina os números DTMF "123456789*0#". A saída desse discriminador deverá ser usada para comutar o sinal de audio de entrada com (fa = 44,1 kHz) para a saída de audio (0 até 9). Uma URA teria ainda um sistema de controle das mensagens a ser enviadas depois de cada interação com o usuário, mas para simplificar é considerado que após a recepção do sinal DTMF a comutação já é realizada.
  • A frequência de amostragem fa2 do sinal de entrada no sistema mostrado abaixo é de 44,1/N kHz, no entanto o sinal gerado no AUDACITY é amostrado em fa1 = 44,1 kHz, portanto antes do circuito abaixo é necessário incluir um filtro antialiasing (low pass) com fc = (44,1/2)/N kHz e um circuito para subamostrar (downsampling) o sinal com fa1 = 44,1 kHz. Esse filtro deverá ser do tipo IIR, sendo de aproximação Butterworth ou Chebyshev tipo 1, para que tenha uma resposta de frequência monotonicamente decrescente após a frequência de corte Assim o sistema terá duas frequências de amostragem, fa1 = 44,1 kHz até o subamostrador e fa2 = 44,1/N kHz depois deste circuito. O fator de subamostragem do sinal de entrada a ser usado por cada equipe é mostrado na tabela abaixo.
Equipe Membros Fator subamostragem (N)
E1 Alisson & Guilherme 11
E2 Felipe & Yan 10
E3 Adonis & Victor 8
E4 Rafael & Suyan 6
  • As especificações do discriminador de frequência, mostrado na figura, são:
DiscriminadorDTMF.png

De acordo com ETSI ES 201 235-3 - Specification of Dual Tone Multi-Frequency (DTMF) Transmitters and Receivers; Part 3: Receivers as características dos componentes do receptor são: — low-pass filter F <= 960 Hz, no entanto adotaremos a frequência sqrt(941*1209); — high-pass filter F >= 1190 Hz, no entanto adotaremos a frequência sqrt(941*1209); — two threshold comparators; — eight channel filters (two pole filters: -1,5 dB band pass limits at the nominal frequency ± (1,5%+ 2 Hz)); — eight rectifiers; — eight single pole (RC) filters, time constant \tau = 4 a 5 ms. The thresholds and decision logic are regulated in order not to recognize character signals with a level below -30 dBm and to recognize valid character signals with a level in the range of -4 dBm to -30 dBm. </syntaxhighlight> De acordo com ETSI ES 201 235-2 - Specification of Dual Tone Multi-Frequency (DTMF) Transmitters and Receivers; Part 3: Transmitters as características dos componentes do transmissor são:

High group: signalling frequencies , which have nominal values of 1 209 Hz, 1 336 Hz, 1 477 Hz and 1 633 Hz Low group: signalling frequencies, which have nominal values of 697 Hz, 770 Hz, 852 Hz and 941 Hz The tolerances of the output frequencies shall be within ±1,5 % of their nominal values. The sending levels when the DTMF transmitter is terminated with the reference impedance ZR,shall be: - for the high frequency group: -9,0 dBV, +2,0 dB / -2,5 dB; - for the low frequency group: -11,0 dBV, +2,5 dB / -2,0 dB where dBV - Absolute voltage level expressed in decibels relative to 1 volt. The duration of any individual DTMF tone combination sent shall not be less than 65 ms. The duration of the pause between any individual DTMF tone combination shall not be less than 65 ms. </syntaxhighlight>

  • Os filtros das linhas e das colunas devem ser do tipo FIR. As diferenças de atraso entre os filtros devem ser compensadas usados elementos de atraso ().
  • Sugere-se utilizar para o filtro do detector de envoltório um filtro IIR do tipo Butterworth.
  • O circuito retificador deve se implementado pela função abs.
  • O circuito TC - Threshold Comparator não deve ser implementado, pois tem apenas a função de rejeitar sinais de entrada que estão acima de -3dBm ou abaixo de -40 dBm,
  • Como sinais de entrada serão utilizados:
  • 1) Sinais DTMF gerado através da soma das duas frequência correspondente a linha e coluna.
  • 2) Os sinais DTMF necessários para testar todas as linhas e colunas (7 frequencias DTMF) com duração de 0,5 segundos.
  • 3) Um sinal DTMF correspondente a sequência de teclas "123456789*0#" com duração de tom de 65ms e pausa de 65ms para cada tecla. (ver ETSI ES 201 235-2 - Specification of Dual Tone Multi-Frequency (DTMF) Transmitters and Receivers; Part 2: Transmitters)
  • 4) Os sinais do item (2) e (3) somados a um ruído branco, cuja relação sinal/ruído deve alterável entre 0 dB até 80 dB.
  • A seleção do sinal de entrada pode ser feita através de manual switch ou multiport switch.
  • Os discriminadores das linhas e colunas deverão ser feitos usando filtros centrados nas frequências centrais das linhas e colunas correspondentes.
  • As entradas DTMF deverão deverão gerar valores binário correspondente 1 => 0001 a 9 => 1001, 0 -> 1010, * => 1011, # => 1100.
  • Caso nenhuma frequência seja ativada o discriminador deverá indicar um código de "0000" (0); caso seja ativada apenas uma frequência o discriminador deverá indicar um código de ERRO 1 "1101" (13); caso sejam ativados 2 linhas e 1 coluna ou 1 linha e duas colunas o discriminador deverá indicar um código de ERRO 2 "1110" (14); caso contrário o discriminador deverá indicar um código de ERRO 3 "1111" (15);
  • Neste projeto os sinais de entrada a serem utilizados deverão ser gerados com frequência de amostragem especificada. Antes de realizar o processamento indicado no diagrama do receptor DTMF, deve ser feita uma filtragem passa-baixa com um filtro de no mínimo 2 polos do tipo Butterworth ou Chebychev 1, com fp = (f_s/2)/N kHz, seguido de uma subamostragem Downsampling de um fator N para a nova frequência f_s2.
  • A contante de tempo é equivalente a uma frequência de corte de
NORMAS DE DTMF
  • Ver as especificações DTMF em:
  • Outros links auxiliares:

Avaliações

  • Prova escrita A1
  • Entrega do Projeto Final. O projeto é avaliado nos quesitos:
PFe - Documento de Especificação (apresentado no relatório);
PFp - Implementação do Projeto;
PFr - Relatório do Projeto (excluído a especificação);
PFi - Avaliação individual do aluno no projeto (conceito subjetivo atribuído pelo professor a partir da observação e da apresentação do projeto).

Referências Bibliográficas

  1. 1,0 1,1 1,2 1,3 DINIZ, P. S. R., DA SILVA, E. A. B., e LIMA NETTO, S. Processamento Digital de Sinais: Projeto e Análise de Sistemas. 2. ed. Porto Alegre: Bookman, 2014. 976 p. ISBN 978-8582601235
  2. 2,0 2,1 2,2 2,3 2,4 2,5 2,6 2,7 2,8 SHENOI, B. A. Introduction to Digital Signal Processing and Filter Design. 1.ed. New Jersey: John Wiley-Interscience, 2006. 440 p. ISBN 978-0471464822
  3. LATHI, Bhagwandas P. Sinais e Sistemas Lineares. 2. ed. Porto Alegre: Artmed-Bookman, 2007. 856 p. ISBN 978-8560031139


Curso de Engenharia de Telecomunicações