Procesadores digitales de
señales
- Práctica 1. Sistemas LTI y
análisis de Fourier. - Práctica 2.
Muestreo. - Práctica 3 La
Transformada Z. - Práctica 4. La DFT y la FFT.
y Práctca 5.Diseño de filtros
digitales.
Nota: Para la realización de
esta memoria se ha
utilizado el procesador de
textos "Microsoft Word
2002 XP"
PRACTICA 1. Sistemas LTI y
análisis de Fourier.
>> help conv
La función
conv(x,y) sirve para calcular la convolucion de x e y. El vector
que resulta es de longitud Length(a)+length(b)-1. Si A y B son
vectores de
coeficientes polinómicos, la convolución es
equivalente a multiplicar los dos polinómios.
CONV Convolution and polynomial
multiplication.
C = CONV(A, B) convolves vectors A and
B. The resulting
vector is length
LENGTH(A)+LENGTH(B)-1.
If A and B are vectors of polynomial
coefficients, convolving
them is equivalent to multiplying the
two polynomials.
See also DECONV, CONV2, CONVN, FILTER
and, in the Signal
Processing Toolbox, XCORR,
CONVMTX.
>> help filter
La función filter() permite
calcular la salida de un sistema LTI
causal cuando éste es definido mediante una
ecuación en diferencias de coeficientes
constantes.
Si x es un vector con la secuencia de
entrada y a y b son vectores con los coeficientes
de la ecuación , la expression y=filter(b,a,x) devuelve un
vector y con la respuesta del sistema LTI
causal.
FILTER One-dimensional digital filter.
Y = FILTER(B,A,X) filters the data in vector X with
the
filter described by vectors A and B to create the
filtered
data Y. The filter is a "Direct Form II
Transposed"
implementation of the standard difference
equation:
a(1)*y(n) = b(1)*x(n) + b(2)*x(n-1) + … +
b(nb+1)*x(n-nb)
– a(2)*y(n-1) – … – a(na+1)*y(n-na)
If a(1) is not equal to 1, FILTER normalizes the
filter
coefficients by a(1).
FILTER always operates along the first non-singleton
dimension,
namely dimension 1 for column vectors and non-trivial
matrices,
and dimension 2 for row vectors.
[Y,Zf] = FILTER(B,A,X,Zi) gives access to
initial and final
conditions, Zi and Zf, of the delays. Zi is a vector of
length
MAX(LENGTH(A),LENGTH(B))-1 or an array of such vectors,
one for
each column of X.
FILTER(B,A,X,[],DIM) or FILTER(B,A,X,Zi,DIM) operates
along the
dimension DIM.
See also FILTER2 and, in the Signal Processing Toolbox,
FILTFILT.
Ejercicio 1.2 La función
"conv"
c=1; % Control
if c==1
x=ones(1,6); % Secuencia
y=conv(x,x); % Convolución x*x
Nx=0; % Indices
Mx=5;
Ny=Nx+Nx;
My=Mx+Mx;
ny=Ny:My;
figure % Nueva figura
stem(ny,y); % Resultados
title('Figura 1.1')
end
Ejercicio 1.2.1 Promediador
Móvil
c=1; % Control
if c==1
n=-30:30; % Tiempo
discreto
h=1/7*ones(1,7); % Respuesta impulsional
x=sin(2*pi*n/13); % Secuencia
y=conv(h,x); % Convolución
V=zeros(1,81); % Vector de ceros
Nh=-3; % Indices h[n]
Mh=3;
Nx=-33; % Indices x[n]
Mx=33;
Ny=Nx+Nx; % Indices y[n]
My=Mx+Mx;
ny=Ny:My;
V(8:8+66)=y; % Salida ajustada
nV=-40:40; % Tiempo discreto
figure % Nueva figura
stem(nV,V); % Resultados
title('Figura 1.2')
end
Ejercicio 1.3 La Función Filter
c=1; % Control
if c==1
a=[1 2]; % Coeficientes y[n]
b=[1 -3]; % Coeficientes x[n]
n=0:5; % Tiempo discreto
x=n; % Secuencia x[n]
y=filter(b,a,x); % Filtrado de x[n]
figure % Nueva figura
stem(x,y); % Resultados
title('Figura 1.3')
end
Para ver el gráfico seleccione la
opción "Descargar" del menú superior
Ejercicio 1.3.1
c=1; % Control
if c==1
%— Punto 1
% Coeficientes y[n]
a=[1 (-1.8*cos(pi/16)) 0.81];
b=[1 1/2]; % Coeficientes x[n]
n=0:50; % Tiempo discreto
x=zeros(1,51); % Vector de ceros
x(1)=1;
y=filter(b,a,x); % Filtrado de x[n]
figure % Nueva figura
stem(n,y); % Resultados
title('Figura 1.4')
Para ver el gráfico seleccione la
opción "Descargar" del menú superior
%—
n=-10:50; % Tiempo discreto
x=zeros(1,61); % Vector de ceros
x(11)=1;
y=filter(b,a,x); % Filtrado de x[n]
figure % Nueva figura
stem(n,y); % Resultados
title('Figura 1.4 Modificada la x[n]')
Para ver el gráfico seleccione la
opción "Descargar" del menú superior
%— Punto 2
n=-10:50; % Tiempo discreto
x=zeros(1,61); % Vector de ceros
x(11:61)=1;
y=filter(b,a,x); % Filtrado de x[n]
figure % Nueva figura
stem(n,y); % Resultados
title('Figura 1.5')
Para ver el gráfico seleccione la
opción "Descargar" del menú superior
%—
n=-10:50; % Tiempo discreto
x=zeros(1,61); % Vector de ceros
x(11:61)=1;
y=filter(b,a,x); % Filtrado de x[n]
yc=cumsum(y,1); % "Filtrado" de y[n]
figure % Nueva figura
stem(n,yc); % Resultados
title('Figura 1.5 con Cumsum')
Para ver el gráfico seleccione la
opción "Descargar" del menú superior
%—
n=-10:50; % Tiempo discreto
x=zeros(1,61); % Vector de ceros
x(11:61)=1;
y=filter(b,a,x); % Filtrado de x[n]
y1=diff(y); % Recuperación
n1=-9:50; % Nuevo intervalo
figure % Nueva figura
stem(n1,y1); % Resultados
title('Respuesta al impulso con diff')
Para ver el gráfico seleccione la
opción "Descargar" del menú superior
%— Punto 3
a=[1]; % Coeficientes y[n]
b=1/7*ones(1,7); % Coeficientes x[n]
n=-30:30; % Tiempo discreto
x=sin(2*pi*n/13); % Secuencia
y=filter(b,a,x); % Filtrado de x[n]
figure % Nueva figura
stem(n,y); % Resultados
title('Figura 1.2 con función filter')
Para ver el gráfico seleccione la
opción "Descargar" del menú superior
end
Ejercicio 1.4 Cancelación de
ecos
c=1; % Control
if c==1
% Lectura
WAV
[voz,fs]=wavread('ho1.wav');
sound(voz,fs) % Reproducción
% [voz fs]=loadw16('ho1.wav');
pause % Control de
reproducción
%– Reproducción con ecos 1.4.1
N=1000; % Retardo
alfa=0.5; % Factor de atenuación
% Sistema
b=[1 zeros(1,999) alfa];
% Filtrado
ve=filter(b,1,voz);
sound(ve,fs) % Reproducción con eco
%—
figure % Nueva figura
stem(b); % Resultados
title('Figura 1.6')
pause % Control de reproducción
%– Sistema inverso 1.4.2
% Coeficientes y[n]
a=[1 zeros(1,999) alfa];
b=[1]; % Coeficientes x[n]
y=filter(b,a,ve); % Señal sin ecos
sound(y,fs) % Reproducción sin eco
%— Cálculo del sistema
inverso
n=0:4000; % Intervalo discreto
a=[1 zeros(1,999) alfa]; % Coeficientes de
y[n]
b=1; % Coeficientes de x[n]
x=zeros(1,4000); % Secuencia impulso
x(1)=1;
hi=filter(b,a,x); % Filtrado
figure % Nueva figura
subplot(211)
plot(n,hi) % Respuesta al impulso
title('Respuesta al impulso nuevo sistema')
b=[1 zeros(1,999) alfa]; % Sistema con eco
he=filter(b,1,x); % Filtrado
ht=conv(hi,he); % Convolución
subplot(212)
plot(n,ht) % Respuesta al impulso
title('Convolución sistema inverso y con
eco')
end
Ejercicio 1.5.1 Señales
de duración finita
c=1; % Control
if c==1
p=ones(1,15); % Pulso de duración 15
ep=fft(p,256); % Espectro del pulso
figure % Nueva figura
subplot(221)
plot(abs(ep)) % Módulo del espectro
title('Espectro del pulso ( 0 y 1 )')
ep=fftshift(ep); % Espectro del pulso
% Frecuencia
f=linspace(-0.5,0.5-(1/256),256);
subplot(222)
plot(f,abs(ep)) % Módulo del espectro
title('Espectro del pulso ( -0.5 y 0.5 )')
subplot(223)
plot(f,ep) % Módulo y fase del
espectro
title('Módulo y Fase del espectro')
subplot(224)
plot(f,imag(ep)) % Módulo y fase del
espectro
title('Fase del espectro')
Para ver el gráfico seleccione la
opción "Descargar" del menú superior
ep=ep.*exp(j*2*pi*f*7);
figure % Nueva figura
subplot(211)
plot(f,real(ep)) % Módulo del
espectro
title('Módulo del espectro')
subplot(212)
plot(f,imag(ep)) % Módulo del
espectro
title('Fase del espectro')
Para ver el gráfico seleccione la
opción "Descargar" del menú superior
Ejercicio 1.5.2 Ecuaciones en
diferencia
c=1; % Control
if c==1
% Espectro señal exponencial
[H,W]=freqz(1,[1 -0.9],256);
figure % Nueva figura
subplot(221)
plot(W/(2*pi),abs(H)) % Módulo del espectro
title('Módulo del espectro ( 0 y pi )')
% Espectro señal exponencial
[H,W]=freqz(1,[1 -0.9],256,'whole');
subplot(222)
plot(W/(2*pi),abs(H)) % Módulo del espectro
title('Módulo del espectro ( 0 y 2*pi )')
W=W-pi; % Espectro entre -0.5 y 0.5
H=fftshift(H);
subplot(223)
plot(W/(2*pi),abs(H)) % Módulo del espectro
title('Módulo del espectro ( -0.5 y 0.5 )')
subplot(224)
plot(W/(2*pi),angle(H)) % Fase del espectro
title('Fase del espectro ( -0.5 y 0.5 )')
Para ver el gráfico seleccione la
opción "Descargar" del menú superior
end
>> help fftshift
FFTSHIFT es útil para visualizar
la Transformada de Fourier con la componente de la frecuencia
cero en el centro del espectro. Por vectores, FFTSHIFT(X)
intercambia las mitades izquierdas y derechas de X. Por matrices,
FFTSHIFT(X) intercambia los primeros y terceros cuadrantes y los
segundos y cuartos cuadrantes.
FFTSHIFT Shift zero-frequency component
to center of spectrum.
For vectors, FFTSHIFT(X) swaps the left
and right halves of
X. For matrices, FFTSHIFT(X) swaps the
first and third
quadrants and the second and fourth
quadrants. For N-D
arrays, FFTSHIFT(X) swaps "half-spaces"
of X along each
dimension.
FFTSHIFT(X,DIM) applies the FFTSHIFT
operation along the
dimension DIM.
FFTSHIFT is useful for visualizing the
Fourier transform with
the zero-frequency component in the
middle of the spectrum.
See also IFFTSHIFT, FFT, FFT2,
FFTN.
Ejercicio 2.2 SEÑALES
SINUSOIDALES
c=1; % Control
if c==1
fs=8000; % Frecuencia de muestreo
T=1/fs; % Período de muestreo
f0=1000; % Frecuencia señal sinusoidal
n=0:7999; % Tiempo discreto
x=sin(2*pi*f0*n*T); % Secuencia x[n]
sound(x,fs) % Reproducción de x[n]
figure % Nueva figura
t=n*T; % Intervalo temporal
% Visualiza resultados
subplot(211), stem(n(1:20),x(1:20))
title('Señal sinusoidal x[n]'), grid
xlabel('n'), ylabel('Unidades Lineales')
subplot(212), plot(t(1:20),x(1:20))
title('Señal sinusoidal x(t)'), grid
xlabel('Tiempo (seg)'), ylabel('Unidades Lineales')
figure % Nueva figura
N=length(x); % Longitud de x[n]
% Espectro de x[n]
X=fftshift(fft(x,N))/N;
% Frecuencia discreta
f=linspace(-0.5,0.5-(1/N),N);
% Visualiza resultados
subplot(211), plot(f,abs(X))
title('Espectro discreto de x[n]'), grid
xlabel('fd=f/fs'), ylabel('Unidades Lineales')
subplot(212), plot(f*fs,abs(X))
title('Espectro continuo de x[n]'), grid
xlabel('Frecuencia (Hz)'), ylabel('Unidades
Lineales')
Para ver el gráfico seleccione la
opción "Descargar" del menú superior
end
Ejercicio 2.3 SEÑAL CHIRP
c=1; % Control
if c==1
fs=8000; % Frecuencia de muestreo
T=1/fs; % Período de muestreo
f0=3000; % Frecuencia señal
sinusoidal
beta=2000; % Beta
n=0:1*(7999); % Tiempo discreto
% Secuencia x[n]
x=sin(2*pi*f0*n*T+0.5*beta*(n*T).^2);
sound(x,fs) % Reproducción de x[n]
figure % Nueva figura
t=n*T; % Intervalo temporal
% Visualiza resultados
subplot(211), stem(n(1:100),x(1:100))
title('Señal chirp x[n]'), grid
xlabel('n'), ylabel('Unidades Lineales')
subplot(212), plot(t(1:100),x(1:100))
title('Señal chirp x(t)'), grid
xlabel('Tiempo (seg)'), ylabel('Unidades
Lineales')
Para ver el gráfico seleccione la
opción "Descargar" del menú superior
figure % Nueva figura
N=length(x); % Longitud de x[n]
% Espectro de x[n]
X=fftshift(fft(x,N))/N;
% Frecuencia discreta
f=linspace(-0.5,0.5-(1/N),N);
% Visualiza resultados
subplot(211), plot(f,abs(X))
title('Espectro discreto de x[n]'), grid
xlabel('fd=f/fs'), ylabel('Unidades
Lineales')
subplot(212), plot(f*fs,abs(X))
title('Espectro continuo de x[n]'), grid
xlabel('Frecuencia (Hz)'), ylabel('Unidades
Lineales')
Para ver el gráfico seleccione la
opción "Descargar" del menú superior
end
Ejercicio 2.4 CAMBIO DE LA
FRECUENCIA DE MUESTREO
Ejercicio 2.4.1 Diezmado
c=1; % Control
if c==1
n=0:124; % Tiempo discreto
x1=0.4*(n-62); % Secuencia x1[n]
x1=(sinc(x1)).^2;
x2=0.2*(n-62); % Secuencia x2[n]
x2=(sinc(x2)).^2;
figure % Nueva figura
% Visualiza resultados
subplot(211), plot(n,x1)
title('Secuencia x1[n]'), grid
xlabel('n'), ylabel('Unidades Lineales')
subplot(212), plot(n,x2)
title('Secuencia x2[n]'), grid
xlabel('n'), ylabel('Unidades Lineales')
figure % Nueva figura
% Espectro de x1[n]
X1=fftshift(fft(x1,2048));
% Espectro de x2[n]
X2=fftshift(fft(x2,2048));
% Frecuencia discreta
f=linspace(-0.5,0.5-(1/2048),2048);
% Visualiza resultados
subplot(211), plot(f,abs(X1))
title('Espectro discreto de x1[n]'), grid
xlabel('fd=f/fs'), ylabel('Unidades
Lineales')
subplot(212), plot(f,abs(X2))
title('Espectro discreto de x2[n]'), grid
xlabel('fd=f/fs'), ylabel('Unidades
Lineales')
%— Diezmado por un factor M
% Secuencia diezmada xd1[n]
xd1=x1(1:2:length(x1));
% Secuencia diezmada xd2[n]
xd2=x2(1:2:length(x2));
figure % Nueva figura
% Visualiza resultados
subplot(211), plot(xd1)
title('Secuencia diezmada xd1[n]'), grid
xlabel('n'), ylabel('Unidades Lineales')
subplot(212), plot(xd2)
title('Secuencia diezmada xd2[n]'), grid
xlabel('n'), ylabel('Unidades Lineales')
figure % Nueva figura
% Espectro de x1[n]
Xd1=fftshift(fft(xd1,2048));
% Espectro de x2[n]
Xd2=fftshift(fft(xd2,2048));
% Frecuencia discreta
f=linspace(-0.5,0.5-(1/2048),2048);
% Visualiza resultados
subplot(211), plot(f,abs(Xd1))
title('Espectro discreto de xd1[n]'), grid
xlabel('fd=f/fs'), ylabel('Unidades
Lineales')
subplot(212), plot(f,abs(Xd2))
title('Espectro discreto de xd2[n]'), grid
xlabel('fd=f/fs'), ylabel('Unidades
Lineales')
%— Filtro antialiasing
fc=0.25; % Frecuencia de corte 0.5/M
% Respuesta al impulso del filtro
h=2*fc*sinc(2*fc*(-32:32)).*(hamming(65)');
figure % Nueva figura
freqz(h,1) % Visualiza resultados
%— Filtrando x1[n]
xf1=conv(h,x1); % Secuencia filtrada xf2[n]
% Secuencia filtrada diezmada xfd2[n]
xfd1=xf1(1:2:length(xf1));
figure % Nueva figura
% Espectro de xfd2[n]
Xfd1=fftshift(fft(xfd1,2048));
% Frecuencia discreta
f=linspace(-0.5,0.5-(1/2048),2048);
% Visualiza resultados
plot(f,abs(Xfd1))
title('Espectro discreto de xfd1[n]'), grid
xlabel('fd=f/fs'), ylabel('Unidades
Lineales')
end
Ejercicio 2.4.2 Interpolación
c=1; % Control
if c==1
n=0:124; % Tiempo discreto
x1=0.4*(n-62); % Secuencia x1[n]
x1=(sinc(x1)).^2;
x2=0.2*(n-62); % Secuencia x2[n]
x2=(sinc(x2)).^2;
%— Interpolando por un factor L
figure % Nueva figura
% Visualiza resultados
stem(x1(52:72))
title('Secuencia x1[n]'), grid
xlabel('n'), ylabel('Unidades Lineales')
figure % Nueva figura
stem(x2(52:72))
title('Secuencia x2[n]'), grid
xlabel('n'), ylabel('Unidades Lineales')
% Secuencia interpolada xe1[n]
xe1=zeros(1,3*length(x1));
xe1(1:3:length(xe1))=x1;
% Secuencia interpolada xe2[n]
xe2=zeros(1,3*length(x2));
xe2(1:3:length(xe2))=x2;
figure % Nueva figura
% Visualiza resultados
stem(xe1(177:197))
title('Secuencia interpolada xe1[n]'), grid
xlabel('n'), ylabel('Unidades Lineales')
figure % Nueva figura
stem(xe2(177:197))
title('Secuencia interpolada xe2[n]'), grid
xlabel('n'), ylabel('Unidades Lineales')
figure % Nueva figura
% Espectro de xe1[n]
Xe1=fftshift(fft(xe1,2048));
% Espectro de xe2[n]
Xe2=fftshift(fft(xe2,2048));
% Frecuencia discreta
f=linspace(-0.5,0.5-(1/2048),2048);
% Visualiza resultados
subplot(211), plot(f,abs(Xe1))
title('Espectro discreto de xe1[n]'), grid
xlabel('fd=f/fs'), ylabel('Unidades
Lineales')
subplot(212), plot(f,abs(Xe2))
title('Espectro discreto de xe2[n]'), grid
xlabel('fd=f/fs'), ylabel('Unidades
Lineales')
%— Filtro interpolador
fc=0.16; % Frecuencia de corte 0.5/L
% Respuesta al impulso del filtro
h=(2*fc*sinc(2*fc*(-32:32)).*(hamming(65)'))*3;
figure % Nueva figura
freqz(h,1) % Visualiza resultados
%— Filtrando xe1[n] y xe2[n]
% Xfe1=conv(h,xe1); % Secuencia filtrada
xfe1[n]
% xfe2=conv(h,xe2); % Secuencia filtrada
xfe2[n]
xfe1=filter(h,1,xe1); % Secuencia filtrada
xfe1[n]
xfe2=filter(h,1,xe2); % Secuencia filtrada
xfe2[n]
figure % Nueva figura
% Visualiza resultados
stem(xfe1(209:229))
title('Secuencia interpolada filtrada xfe1[n]'),
grid
xlabel('n'), ylabel('Unidades Lineales')
figure % Nueva figura
stem(xfe2(209:229))
title('Secuencia interpolada filtrada xfe2[n]'),
grid
xlabel('n'), ylabel('Unidades Lineales')
figure % Nueva figura
% Espectro de xfe1[n]
Xfe1=fftshift(fft(xfe1,2048));
% Espectro de xfe2[n]
Xfe2=fftshift(fft(xfe2,2048));
% Frecuencia discreta
f=linspace(-0.5,0.5-(1/2048),2048);
% Visualiza resultados
subplot(211), plot(f,abs(Xfe1))
title('Espectro discreto de xfe1[n]'), grid
xlabel('fd=f/fs'), ylabel('Unidades
Lineales')
subplot(212), plot(f,abs(Xfe2))
title('Espectro discreto de xfe2[n]'), grid
xlabel('fd=f/fs'), ylabel('Unidades
Lineales')
end
Ejercicio 2.4.3 Cambio fs por un factor
racional
c=1; % Control
if c==0
% Carga archivo de
voz
[voz,fs]=wavread('ho1.wav');
% Versiones anteriores a MATLAB 5.3
% [voz fs]=loadw16('ho1.wav');
clc
disp('Pausa de voz… cualquier tecla para
continuar')
sound(voz,fs) % Reproducción archivo de
voz
pause
%— Diezmado por un factor M
% Voz diezmada
vozd=voz(1:3:length(voz));
disp('Pausa de voz diezmada… cualquier tecla para
continuar')
sound(vozd,fs/3) % Reproducción voz
diezmada
pause
%— Filtro antialiasing
fc=0.16; % Frecuencia de corte 0.5/M
% Respuesta al impulso del filtro
h=2*fc*sinc(2*fc*(-32:32)).*(hamming(65)');
%— Filtrando voz
% Voz filtrada
vozf=filter(h,1,voz);
% Voz diezmada
vozfd=vozf(1:3:length(vozf));
disp('Pausa de voz filtrada diezmada… cualquier tecla
para continuar')
sound(vozfd,fs/3) % Reproducción voz filtrada
diezmada
pause
%— Interpolando por un factor L
clc
disp('Pausa de voz… cualquier tecla para
continuar')
sound(voz,fs) % Reproducción archivo de
voz
pause
% Voz interpolada
voze=zeros(1,3*length(voz));
voze(1:3:length(voze))=voz;
disp('Pausa de voz interpolada… cualquier tecla para
continuar')
sound(voze,3*fs) % Reproducción voz
interpolada
pause
%— Filtro interpolador
fc=0.16; % Frecuencia de corte 0.5/L
% Respuesta al impulso del filtro
h=(2*fc*sinc(2*fc*(-32:32)).*(hamming(65)'))*3;
%— Filtrando voz
% Voz filtrada
vozfe=filter(h,1,voze);
disp('Pausa de voz filtrada expandida… cualquier tecla
para continuar')
sound(vozfe,3*fs) % Reproducción voz filtrada
expandida
pause
%— Cambio de la frecuencia de muestreo
fs'=fs*(L/M)
%— Solución 1: L=10, M=8
%— Solución 2: L=5, M=4
clc % Borra la pantalla
% Control del cambio de fs
L=5; % Factor de interpolación
M=4; % Factor de diezmado
%— Interpolando por un factor L
disp('Pausa de voz… cualquier tecla para
continuar')
sound(voz,fs) % Reproducción archivo de
voz
pause
% Voz interpolada
voze=zeros(1,L*length(voz));
voze(1:L:length(voze))=voz;
%— Filtro interpolador
fc=0.5/L; % Frecuencia de corte 0.5/L
% Respuesta al impulso del filtro
h=(2*fc*sinc(2*fc*(-32:32)).*(hamming(65)'))*L;
%— Filtrando voz
% Voz filtrada
vozfe=filter(h,1,voze);
%— Diezmado por un factor M
% Voz diezmada
vozfed=vozfe(1:M:length(vozfe));
% Reproducción voz filtrada, expandida y
diezmada
disp('Pausa de voz filtrada, expandida y diezmada…
cualquier tecla para continuar')
sound(vozfed,fs*(L/M))
pause
end
>> help
polar
La función polar() permite
hacer diagramas en
forma polar, es decir, dado un numero complejo, polar()
representa su módulo y fase.
POLAR(THETA, RHO) hace un diagrama
usando los coordenadas polares del angulo de la THETA , en
radianes, frente el radio
RHO.
POLAR Polar coordinate
plot.
POLAR(THETA, RHO) makes a plot
using polar coordinates of
the angle THETA, in radians,
versus the radius RHO.
POLAR(THETA,RHO,S) uses the
linestyle specified in string S.
See PLOT for a description of
legal linestyles.
See also PLOT, LOGLOG, SEMILOGX,
SEMILOGY.
>> help
roots
La función roots(C)
encuentra las raíces polinómicas del vector C,
donde los ceros y los polos son las raíces de los
polinomios del numerador y del denominador de la función
de transferencia, respectivemente.
ROOTS Find polynomial
roots.
ROOTS(C) computes the roots of the
polynomial whose coefficients
are the elements of the vector C.
If C has N+1 components,
the polynomial is C(1)*X^N + … +
C(N)*X + C(N+1).
See also POLY, RESIDUE,
FZERO.
>> help
poly
La función poly() construye
un polinomio a partir de sus raíces.
POLY Convert roots to
polynomial.
POLY(A), when A is an N by N
matrix, is a
row vector with
N+1 elements which are the
coefficients of the
characteristic polynomial,
DET(lambda*EYE(SIZE(A)) – A) .
POLY(V), when V is a vector, is a
vector whose elements are
the coefficients of the polynomial
whose roots are the
elements of V . For vectors, ROOTS
and POLY are inverse
functions of each other, up to
ordering, scaling, and
roundoff error.
ROOTS(POLY(1:20)) generates
Wilkinson's famous example.
See also ROOTS, CONV, RESIDUE,
POLYVAL.
Ejercicio 3.2 MANEJO DE FUNCIONES DE
TRANSFERENCIA
Ejercicio 3.2.1 REPRESENTACION GRAFICA DE POLOS
Y CEROS
c=1; % Control
if c==1
%— Punto 1
z=0.7*exp(j*pi/3); % Variable compleja z
figure % Nueva figura
% Visualiza resultados
polar(angle(z),abs(z),'o')
title('Variable compleja en forma polar'),
grid
z=0.5+j*0.5; % Variable compleja z
Para ver el gráfico seleccione la
opción "Descargar" del menú superior
figure % Nueva figura
% Visualiza resultados
polar(angle(z),abs(z),'o')
title('Variable compleja en forma cartesiana'),
grid
Para ver el
gráfico seleccione la opción "Descargar" del
menú superior
figure % Nueva figura
% Variables
complejas
z=[0.5+j*0.5 0.7*exp(j*pi/3)];
% Visualiza resultados
polar(angle(z),abs(z),'x')
title('Variables complejas'), grid
Para ver el gráfico seleccione la
opción "Descargar" del menú superior
%— Punto 2
clc % Borra la pantalla
disp('Raices del polinomio z^2-z-2')
A=[1 -1 -2]; % Polinomio z^2-z-2
r=roots(A) % Raices del polinomio z^2-z-2
%— Punto 3
disp('Polinomio asociado con las raices…')
A=poly([-1 2]) % Polinomio asociado con las
raices…
%— Punto 4
A=[1 -1 -2]; % Polinomio z^2-z-2
figure % Nueva figura
subplot(211)
zplane(A); % Diagrama de ceros y polos
title('Diagrama de ceros y polos para A,[]'),
grid
subplot(212)
zplane(A'); % Diagrama de ceros y polos
title('Diagrama de ceros y polos para A`,[]'),
grid
figure % Nueva figura
zplane(A); % Diagrama de ceros y polos
title('Diagrama de ceros y polos para A,[0]'),
grid
%—
% Sistema con H1(z)
bH1=roots([1]); % Ceros de H1(z)
% Polos de H1(z)
aH1=roots([1 -0.85*sqrt(2) 0.7225]);
% Sistema con H2(z)
% Ceros de H2(z)
bH2=roots([1 -0.85*sqrt(2) 0.7225]);
aH2=roots([1]); % Polos de H2(z)
% Sistema con H3(z)
bH3=roots([1]) % Ceros de H3(z)
% Polos de H3(z)
aH3=roots([1 -0.95*sqrt(2) 0.9025])
% Sistema con H4(z)
% Ceros de H4(z)
bH4=roots([1 -0.95*sqrt(2) 0.9025]);
aH4=roots([1]) % Polos de H4(z)
figure % Nueva figura
subplot(223)
zplane(bH1,aH1); % Diagrama de ceros y polos
title('Sistema con H1(z)'), grid
subplot(221)
zplane(bH2,aH2); % Diagrama de ceros y polos
title('Sistema con H2(z)'), grid
subplot(224)
zplane(bH3,aH3); % Diagrama de ceros y polos
title('Sistema con H3(z)'), grid
subplot(222)
zplane(bH4,aH4); % Diagrama de ceros y polos
title('Sistema con H4(z)'), grid
end
Ejercicio 3.2.2 RESPUESTA EN FRECUENCIA A PARTIR DE
LA TZ
c=1; % Control
if c==1
%— Punto 1
clc % Borra la pantalla
disp('Evaluación
del polinomio z^2-z-2 en z=-1')
A=[1 -1 -2]; % Polinomio z^2-z-2
z=-1; % Evaluación en z=-1
polyval(A,z) % Evaluación delpolinomio z^2-z-2
en z=-1
%—
disp('Evaluación [H(0), H(0.5)]…')
% Sistema con H1(z)
cH1=[1]; % Ceros de H1(z)
% Polos de H1(z)
pH1=[1 -0.85*sqrt(2) 0.7225];
z=1; % H(0)
H1=(polyval(cH1,z))/(polyval(pH1,z));
z=-1; % H(0.5)
H1=[H1 (polyval(cH1,z))/(polyval(pH1,z))];
% Sistema con H2(z)
% Ceros de H2(z)
cH2=[1 -0.85*sqrt(2) 0.7225];
pH2=[1]; % Polos de H2(z)
z=1; % H(0)
H2=(polyval(cH2,z))/(polyval(pH2,z));
z=-1; % H(0.5)
H2=[H2 (polyval(cH2,z))/(polyval(pH2,z))];
% Sistema con H3(z)
cH3=[1]; % Ceros de H3(z)
% Polos de H3(z)
pH3=[1 -0.95*sqrt(2) 0.9025];
z=1; % H(0)
H3=(polyval(cH3,z))/(polyval(pH3,z));
z=-1; % H(0.5)
H3=[H3 (polyval(cH3,z))/(polyval(pH3,z))];
% Sistema con H4(z)
% Ceros de H4(z)
cH4=[1 -0.95*sqrt(2) 0.9025];
pH4=[1]; % Polos de H4(z)
z=1; % H(0)
H4=(polyval(cH4,z))/(polyval(pH4,z));
z=-1; % H(0.5)
H4=[H4 (polyval(cH4,z))/(polyval(pH4,z))];
H1, H2, H3, H4 % Evaluación de H1, H2, H3 y
H4
%— Punto 2
% Respuesta en frecuencia
[H,W]=freqz([1 -1 -2],1,512,'whole');
figure % Nueva figura
subplot(211)
% Módulo respuesta en frecuencia
plot(W/(2*pi),abs(H))
title('Módulo de la respuesta en frecuencia'),
grid
ylabel('Unidades Lineales')
subplot(212)
% Fase respuesta en frecuencia
plot(W/(2*pi),angle(H))
title('Fase de la respuesta en frecuencia'),
grid
ylabel('Unidades Lineales')
%—
% Sistema con H1(z)
bH1=[1]; % Coeficientes b
% Coeficientes a
aH1=[1 -0.85*sqrt(2) 0.7225];
% Respuesta en frecuencia
[H1,W1]=freqz(bH1,aH1,512,'whole');
% Sistema con H2(z)
% Coeficientes b
bH2=[1 -0.85*sqrt(2) 0.7225];
aH2=[1]; % Coeficientes a
% Respuesta en frecuencia
[H2,W2]=freqz(bH2,aH2,512,'whole');
% Sistema con H3(z)
bH3=[1]; % Coeficientes b
% Coeficientes a
aH3=[1 -0.95*sqrt(2) 0.9025];
% Respuesta en frecuencia
[H3,W3]=freqz(bH3,aH3,512,'whole');
% Sistema con H4(z)
% Coeficientes b
bH4=[1 -0.95*sqrt(2) 0.9025];
aH4=[1]; % Coeficientes a
% Respuesta en frecuencia
[H4,W4]=freqz(bH4,aH4,512,'whole');
figure % Nueva figura
subplot(221)
% Módulo respuesta en frecuencia
plot(W1/(2*pi),abs(H1))
title('Módulo de H1'), grid
ylabel('Unidades Lineales')
subplot(223)
% Fase respuesta en frecuencia
plot(W1/(2*pi),angle(H1))
title('Fase de H1'), grid
ylabel('Unidades Lineales')
subplot(222)
% Módulo respuesta en frecuencia
plot(W2/(2*pi),abs(H2))
title('Módulo de H2'), grid
ylabel('Unidades Lineales')
subplot(224)
% Fase respuesta en frecuencia
plot(W2/(2*pi),angle(H2))
title('Fase de H2'), grid
ylabel('Unidades Lineales')
figure % Nueva figura
subplot(221)
% Módulo respuesta en frecuencia
plot(W3/(2*pi),abs(H3))
title('Módulo de H3'), grid
ylabel('Unidades Lineales')
subplot(223)
% Fase respuesta en frecuencia
plot(W3/(2*pi),angle(H3))
title('Fase de H3'), grid
ylabel('Unidades Lineales')
subplot(222)
% Módulo respuesta en frecuencia
plot(W4/(2*pi),abs(H4))
title('Módulo de H4'), grid
ylabel('Unidades Lineales')
subplot(224)
% Fase respuesta en frecuencia
plot(W4/(2*pi),angle(H4))
title('Fase de H4'), grid
ylabel('Unidades Lineales')
end
Ejercicio 3.2.3 RESPUESTA AL IMPULSO
Para ver el ejercicio y el
gráfico seleccione la opción "Descargar" del
menú superior
Ejercicio 3.2.4 DESCOMPOSICION EN FRACCIONES
SIMPLES
c=1; % Control
if c==1
%— Punto 1
a=[6 -5 1]; % Coeficientes a de y[n]
b=[1 -1]; % Coeficientes b de x[n]
[h,t]=impz(b,a,11); % h[n]
figure % Nueva figura
stem(t,h) % Respuesta al impulso del sistema
title('Respuesta al impulso del sistema'), grid
ylabel('Unidades Lineales')
%— Punto 2
clc % Borra la pantalla
disp('Expansión en fracciones simples…')
% Expansión en fracciones simples
[r,p,k]=residue(b,a)
%— Punto 3
n=0:10; % Tiempo discreto
% h[n]
h=-0.5*(0.5).^n+0.6667*(0.3333).^n;
figure % Nueva figura
stem(n,h) % Respuesta al impulso del sistema
title('Respuesta al impulso del sistema por fracciones
simples'), grid
ylabel('Unidades Lineales')
end
Ejercicio 3.3 PROPIEDADES EN EL DOMINIO
TRANSFORMADO
Ejercicio 3.3.1 FASE LINEAL
Para ver el
ejercicio y el gráfico seleccione la opción
"Descargar" del menú superior
PRACTICA 4. La DFT y la FFT.
y PRACTCA 5.Diseño
de filtros digitales.
>> help
fft
FFT(X) calcula la
Transformada de Fourier discreta (DFT) del vector
X.
FFT Discrete
Fourier transform.
FFT(X) is the
discrete Fourier transform (DFT) of vector X.
For
matrices, the FFT
operation is applied to each column. For N-D
arrays, the FFT
operation operates on the first non-singleton
dimension.
FFT(X,N) is the
N-point FFT, padded with zeros if X has less
than N points and
truncated if it has more.
FFT(X,[],DIM) or
FFT(X,N,DIM) applies the FFT operation across
the
dimension
DIM.
For length N input
vector x, the DFT is a length N vector X,
with
elements
N
X(k) = sum
x(n)*exp(-j*2*pi*(k-1)*(n-1)/N), 1 <= k <=
N.
n=1
The inverse DFT
(computed by IFFT) is given by
N
x(n) = (1/N) sum
X(k)*exp( j*2*pi*(k-1)*(n-1)/N), 1 <= n <=
N.
k=1
The relationship
between the DFT and the Fourier coefficients a and b
in
N/2
x(n) = a0 + sum
a(k)*cos(2*pi*k*t(n)/(N*dt))+b(k)*sin(2*pi*k*t(n)/(N*dt))
k=1
is
a0 = X(1)/N, a(k)
= 2*real(X(k+1))/N, b(k) = -2*imag(X(k+1))/N,
where x is a
length N discrete signal sampled at times t with spacing
dt.
See also IFFT,
FFT2, IFFT2, FFTSHIFT.
>> help
fftshift
Calcula la
componente de la frecuencia cero en el centro del
espectro.
FFTSHIFT Shift
zero-frequency component to center of spectrum.
For vectors,
FFTSHIFT(X) swaps the left and right halves of
X. For matrices,
FFTSHIFT(X) swaps the first and third
quadrants and the
second and fourth quadrants. For N-D
arrays,
FFTSHIFT(X) swaps "half-spaces" of X along each
dimension.
FFTSHIFT(X,DIM)
applies the FFTSHIFT operation along the
dimension
DIM.
FFTSHIFT is useful
for visualizing the Fourier transform with
the zero-frequency
component in the middle of the spectrum.
See also
IFFTSHIFT, FFT, FFT2, FFTN.
>> help
angle
La función
angle calcula la fase , en radianes, de una matriz
compleja.
ANGLE Phase
angle.
ANGLE(H) returns
the phase angles, in radians, of a matrix with
complex
elements.
See also ABS,
UNWRAP.
>> help
freqz
FREQZ Z calcula la
transformada de la respuesta en frecuencia digital del
filtro.
FREQZ Z-transform
digital filter frequency response. When N is an
integer, [H,W] =
FREQZ(B,A,N) returns the N-point frequency
vector W and the
N-point complex frequency response vector G
of the filter
B/A:
-1
-nb
jw B(z) b(1) +
b(2)z + …. + b(nb+1)z
H(e) = —- =
—————————-
-1
-na
A(z) 1 + a(2)z +
…. + a(na+1)z
given numerator
and denominator coefficients in vectors B and
A. The frequency
response is evaluated at N points equally
spaced around the
upper half of the unit circle. To plot
magnitude and
phase of a filter:
[h,w] =
freqz(b,a,n);
mag = abs(h);
phase = angle(h);
semilogy(w,mag),
plot(w,phase)
FREQZ(B,A,N,'whole') uses N points around the whole
unit circle.
FREQZ(B,A,W)
returns the frequency response at frequencies
designated in
vector W, normally between 0 and pi. (See
LOGSPACE to
generate W). See also YULEWALK, FILTER, FFT,
INVFREQZ, and
FREQS.
>> help
fir1
B = FIR1(N, Wn)
diseña un filtro digital de los "paso-bajo" de orden N y
devuelve los coeficientes del filtro en el vector B de longitud
N+1. La frecuencia de atajo Wn debe estar entre 0 < Wn <
1,0, donde 1,0 corresponde a la mitad del rango de la muestra. El
filtro B es real y tiene fase lineal. El aumento normalizado del
filtro en Wn es de -6 dB.
FIR1 FIR filter
design using the window method.
B = FIR1(N,Wn)
designs an N'th order lowpass FIR digital
filter
and returns the
filter coefficients in length N+1 vector B.
The cut-off
frequency Wn must be between 0 < Wn < 1.0, with
1.0
corresponding to
half the sample rate. The filter B is real and
has linear phase.
The normalized gain of the filter at Wn is
-6
dB.
>> help
Kaiser
La función
Kaiser() devuelve el número de elementos de la ventana
Kaiser.
KAISER
KAISER(N,beta) returns the BETA-valued N-point Kaiser
window.
Ejercicio 3.3 PROPIEDADES EN EL
DOMINIO TRANSFORMADO
Ejercicio 3.3.2
ESTABILDAD
Para ver el
ejercicio y el gráfico seleccione la opción
"Descargar" del menú superior
Ejercicio 3.3.3 SISTEMA PASO
TODO
Para ver el
ejercicio y el gráfico seleccione la opción
"Descargar" del menú superior
Ejercicio 3.4.2 RESPUESTA EN
FRECUENCIA DE POLOS
Para ver el
ejercicio y el gráfico seleccione la opción
"Descargar" del menú superior
Ejercicio 3.4.3 APLICACION:
DISEÑO DE FILTROS EN HENDIDURA
Para ver el
ejercicio y el gráfico seleccione la opción
"Descargar" del menú superior
David Savall
Climent
20 de mayo de
2003