Processamento Digital de Sinais – Laboratório

Prof. Eduardo Parente Ribeiro

Roteiro para o  Octave e Matlab
  1. Introdução: O matlab é um ambiente computacional para manipulação de matrizes e visualização dos dados. A versão atual é a 5.1, porem a versão que estamos utilizando é a 4.2. O Octave, atualmente chamado de GNU Octave, visto que é distribuído gratuitamente sob a licença GNU (acompanha os fontes em C++). O Octave funciona em ambiente unix enquanto o Matlab atualmente roda tanto em Windows quanto Unix. O Octave foi desenvolvido por universidades americanas no curso de Química. Ele utiliza uma notação idêntica ao matlab, diferindo apenas em pequenos aspectos. De maneira breve, podemos dizer que o Octave incorpora alguns conceitos que o Matlab não utiliza, como estrutura de dados (na ultima versão do Matlab isto já foi incorporado), e o Matlab possui a parte gráfica mais desenvolvida (O Octave utiliza o Gnuplot para a geração dos gráficos).
  2. Obtendo Ajuda: O primeiro comando a se aprender, e o mais importante é o comando ‘help’. Ele mostra uma lista dos tópicos é que existe ajuda disponível. Para obter ajuda sobre o tópico pode se usar ‘help tópico’. Exemplos:

  3. > help
    > help inv
    No Octave, o manual online pode ser acessado com ‘help –i’, ou ‘help –i tópico’. O programa ‘info’ é invocado e para navegar nesta programa pode-se usar os seguintes comandos:
    Espaço: Próxima Pag.
    DEL: Pag. Anterior

    b: inicio da seção (beggining)

    e: final da seção (ending)

    TAB: pula p/ proxima referencia

    ENTER: Entra na referencia

    u: volta pro nivel anterior (up)

    n: Entra no proximo tópico (next)

    p: Entra no tópico anterior (previous)
     

  1.  Manipulações Básicas

  2. Formando um vetor:

    a = [ 1 2 3 7 6 5]

    Resulta em:

    a =

    1 2 3

    Formando uma Matriz

    b = [1 2 3; 4 5 6; 7 8 9]

    Resulta em :

    b =

    1 2 3

    4 5 6

    7 8 9

    Cada elemento do vetor ou matriz pode ser uma expressão:

    c = [3 –4 2.3 sqrt(2) (4+5)/2]

    Resulta em:

    c =

    3 -4 2.3 1.41 4.5

Transpondo um vetor, ou uma matriz:

Usa-se o operador ‘ (apóstrofe que é a aspa simples. Não confundir com o acento agudo nem o grave).

c = c’

resulta em:

3

-4

2.3

1.41

4.5

b = b’

1 4 7

2 5 8

3 6 9

de novo.

b = b’

resulta em:

1 2 3

4 5 6

7 8 9

O vetor nada mais é do que um caso particular da matriz. O vetor linha é uma matriz de uma linha e várias colunas enquanto o vetor coluna é uma matriz de uma coluna e várias linhas.

Formando um vetor de contagem.

Trata-se de um vetor como outro qualquer porem pode ser gerado com uma notação simplificada. Por ex.:

i=1:10

resulta em:

i =

1 2 3 4 5 6 7 8 9 10

ou então:

x = 0:.1:1

resulta em

x =

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

Acessando os elementos de um vetor ou matriz:

aproveitando o vetor x acima temos:

x(2)

resulta em

0.1

que é o segundo elemento do vetor x. A numeração sempre começa de 1 e não de zero como na linguagem "C’.

De forma análoga, para a matriz b definida anteriormente:

b(2,3)

resulta em:

6

que é o elemento da segunda linha terceira coluna. Sempre (linha, coluna).

Tambem pode-se se acessar uma sub-matriz:

b(1:2,1:2)

resulta em:

1 2

4 5

Isto é equivalente á b([1 2],[1 2]), ou seja estamos usando vetores para indexar uma matriz (ou outro vetor). Outro exemplo:

b([1 3],[1 3])

resulta em:

1 3

7 9

Selecionando uma coluna inteira (o mesmo pode ser usado para linha)

b(:,3)

resulta em:

2

5

8

Modificando um elemento de uma matriz ou vetor.

Da mesma forma que acessamos para leitura um elemento, o mesmo pode ser modificado, colocando-o a esquerda do operador de atribuição "=".

b(2,2)=0

resulta em:

1 2 3

4 0 6

7 8 9

Pode-se tambem operar em sub-matriz, mas deve-se sempre observar a seguinte regra: O tamanho da matriz do lado esquerdo do operador "=" deve ser igual ao tamanho da matriz do lado direito (mesmo número de linhas e colunas).

b([1 3],[2 3])=[0.2 0.3; 0.8 0.9]

resulta em:

1 0.2 0.3

4 0 6

7 0.8 0.9
 
 

Tratando uma matriz como vetor.

Às vezes é interessante operar com todos os elementos (ou alguns elementos) de uma matriz como se fosse um vetor. Por exemplo para achar o maior elemento de um vetor usa-se a função max().

max(c)

resulta em:

4.5

Quando se usa a função max() em uma matriz, o que é retornado é um vetor contendo os maiores elementos de cada coluna (vide help max).

max(b)

resulta em:

7 0.8 6

Se quiser obter o maior elemento de toda a matriz, basta tratar a matriz como um vetor:

max(b(:))

resulta em:

7
 
 

  1. Processamento de Sinais
A função conv() calcula a convolução:

a=conv([1 1 1 1],[1 1 1 1 ])

Vamos criar um sinal x=[1 1 1 1 1] e uma resposta impulsional h=[1 0.5 0.25]

A saída do sistema linear e invariante no tempo (LTI) é dada pela convulção da entrada com a resposta imulsional:

y = conv(x,h);

Visualizando a exponencial complexa

n=0:99; % define um vetor n indo de 0 a 99

w=0.3; %define a frequência digital

d=exp(j*w*n);

Para visualizar utilize o comando subplot:

subplot(4,1,1) % divide a tela em 4 linhas x 1 coluna e seleciona o primeiro grafico

plot(real(d))

subplot(4,1,2) % seleciona o segundo grafico

plot(imag(d))

subplot(4,1,3)

plot(abs(d))

subplot(4,1,4)

plot(angle(d)) % ou então plot(unwrap(angle(d)))

verifique agora outras frequências: w=0.6; w=0.3+2*pi;

Visualizando a FFT

defina um sinal x:

t=0:31;

w=2*pi/32;

x=sin(w*t);

X=fft(x);

subplot(3,1,1); plot(t,x,’^’)

subplot(3,1,2); plot(t,abs(X),’^’)

subplot(3,1,3); plot(t,angle(X))

Agora em crie outra tela para plotar a frequência w=4*pi/100;

figure(2)

w=4*pi/32;

x=sin(w*t);

X=fft(x);

subplot(2,1,1); plot(t,x,‘^’)

subplot(2,1,2); plot(t,abs(X),’^’)

O mesmo para w=15*2*pi/32; na terceira tela. E w=8*2*pi/32 na quarta tela.

Efeito do Janelamento

Agora complete o sinal x com 64 zeros e plote a FFT :

x=[x zeros(1,64)];

tt=0:95;

X=fft(x);

subplot(2,1,1); plot(tt,x)

subplot(2,1,2); plot(t,abs(X),’-@’)
 
 

Verifique a propriedade de linearidade da fft, somando duas senoides (32 pontos) e verificando o espectro. w1=2*pi/32 e w2=10*pi/32; Ex. x3=x2+x1; X3=fft(x3); Depois faça o mesmo, mas usando w2=10.5*pi/32.

Interprete o gráfico.

Verifique a transformada de fourier de algumas janelas e plote-as no mesmo gráfico:

j1=[zeros(1,32) hamming(64)’ zeros(1,32)];

j2=[zeros(1,32) hanning(64)’ zeros(1,32)];

j3=[zeros(1,32) ones(1,64) zeros(1,32);

J1=fft(j1); J2=fft(j2); J3=fft(j3);

n=0:127;

subplot(2,1,1)

plot(n,j1,n,j2,n,j3,n,j4);

subplot(2,1,2)

plot(n,abs(J1),n,abs(J2),n,abs(J3));

Depois, dê um zoom, visualizando os pontos de 1 a 10:

n=0:9;

plot(n,abs(J1(1:10)),n,abs(J2(1:10)),n,abs(J3(1:10)));
 
 

Outras funções: sinc, diric (M), sawtooth (M)

Filtragem Digital

veja a função filter(): help filter

Calcule a Função de transferencia do filtro definido por b=[0.0186 0.0743 0.1114 0.0743 0.0186];

a=[1.0000 -1.5704 1.2756 -0.4844 0.0762]; por dois metodos diferentes.
 
 
 
 

Funções

Crie um função com o seu nome, que retorne a raiz quarta de um numero.

ex. crie um arquivo ‘nome.m’ contendo:

Function a=nome(b)

a = b ^ (1/4);

depois execute:

d = nome(16)
 
 

Verifique a diferença entre a convolução Linear e a convolução circular para as seguintes sequências: a=[1 2 3 4 5]; b=[1 1 1 1]. Sugestão: a função conv() calcula a convolução. Pelo teorema da convolução circular o produto das tranformadas equivale a convolução circular no tempo.

Agora prencha as sequencias com zero e calcule o produto das tranformadas e tire a inversa para obter a convolução linear.
 
 
 
 
 
 

Filtragem

Utilizando o filtro dado pelos coeficientes b=[0.0186 0.0743 0.1114 0.0743 0.0186];

a=[1.0000 -1.5704 1.2756 -0.4844 0.0762]; realize a filtragem do sinal:

t=0:31;

w1=2*pi/32; w2=8*pi/32; w3=12.5*pi/32;

x=sin(w1*t)+sin(w2*t)+sin(w3*t)+0.5*rand(1,32);

Usando três métodos: 1) função filter(), 2) função conv(); 3) multiplicação por H(w)