intTypePromotion=1
zunia.vn Tuyển sinh 2024 dành cho Gen-Z zunia.vn zunia.vn
ADSENSE

Bài tập Xử lý tín hiệu số 2: Thiết kế bộ lọc số IIR thông cao

Chia sẻ: Hoàng Văn Thuận | Ngày: | Loại File: DOCX | Số trang:9

266
lượt xem
20
download
 
  Download Vui lòng tải xuống để xem tài liệu đầy đủ

Nhằm giúp các bạn có thêm tài liệu phục vụ nhu cầu học tập và nghiên cứu, mời các bạn cùng tham khảo nội dung bài tập Xử lý tín hiệu số 2 "Thiết kế bộ lọc số IIR thông cao" dưới đây. Hy vọng nội dung tài liệu phục vụ hữu ích nhu cầu học tập và nghiên cứu.

Chủ đề:
Lưu

Nội dung Text: Bài tập Xử lý tín hiệu số 2: Thiết kế bộ lọc số IIR thông cao

  1. BÀI TẬP XỬ LÝ TÍN HIỆU SỐ 2 Họ và tên: Nguyễn Đình Quý                    Nguyễn Phước Thắng                   Trần Văn Duy Phước Lớp:          11DT3 Nhóm:       40 NỘI DUNG: Thiết kế bộ lọc số IIR thông cao thỏa mãn các yêu cầu đề bài.  YÊU CẦU THIẾT KẾ:           Stopband edge: 0.3π                                              Stopband attenuation: As = 15 dB                                              Passband edge: 0.5π                                              Passband ripple: Rp = 1 dB BÀI 1. Thiết kế bằng cách sử dụng bộ lọc Chebyshev­II, phép biến đổi song tuyến  tính( bilinear transform), hàm zmapping(để chuyển từ bộ lọc số thông thấp sang thông cao).  Phương Pháp Thiết Kế Bộ Lọc Số IIR: Thiết kế bộ lọc tương tự thông thấp. Áp dụng biến đổi song tuyến tính để thu được bộ lọc số thông thấp. Áp dụng biến đổi Frequency­band(dùng hàm zmapping) để chuyển từ bộ lọc số  thông thấp sang bộ lọc số thông cao.  Các bước thực hiện: Bước 1: Từ yêu cầu đề bài, ta có ws và wp của bộ lọc thông cao. Đầu tiên, ta  cần thiết kế bộ lọc thông thấp có các tần số wslp và wplp để sau đó chuyển  đổi thành bộ lọc thông cao. Ta chọn giá trị wplp bằng 1 giá trị hợp lý, chọn  wplp = 0.2*pi. Từ công thức chuyển đổi Frequency Band trong miền Z, ta có: 
  2. Thay wp và wplp vào công thức trên, ta thu được giá trị của alpha. Ta xác định  giá trị wslp từ công thức:  Với Z = exp(j*wslp) và z = exp(j*ws). Suy ra:  Wslp = angle(­(exp(­j*ws)+alpha)/(1+alpha*exp(­j*ws))) Bước 2: Thiết kế bộ lọc analog thông thấp Chebysev­II a.i.1. Chọn giá trị T  = 1. a.i.2. Prewarp các giá trị tần số cắt wplp và wslp. OmegaP = (2/T)*tan(wplp/2) OmegaS = (2/T)*tan(wslp/2) a.i.3. Thiết kế bộ lọc tương tự Chebysev­II với các thông số:  OmegaP, OmegaS, As, Rp. Ta thu được các đa thức ở tử và mẫu  của hàm truyền ở miền s. Bước 3: Áp dụng bilinear transform để chuyển bộ lọc thông thấp tương tự  thành bộ lọc thông thấp số. Ta thu được các đa thức ở tử và mẫu của hàm  truyền của bộ lọc thông thấp ở miền Z. Bước 4: Sử dụng hàm zmapping(Frequency­band transform) để chuyển đổi từ  bộ lọc thông thấp thành bộ lọc thông cao. Ta thu được các đa thức ở tử và mẫu  của hàm truyền của bộ lọc thông cao ở miền z.  Bước 5: Chuyển đổi các đa thức ở tử và mẫu của hàm truyền từ Direct form  sang Cascade form. Bước 6: Tìm hàm đáp ứng xung của bộ lọc số thông cao. Bước 7: Tính các hàm đáp ứng biên độ, đáp ứng biên độ tương đối ở dB, đáp  ứng pha, độ trễ nhóm của bộ lọc. Bước 8: Tính As và Rp của bộ lọc thực tế. Bước 9: Vẽ đồ thị và so sánh, đối chiếu bộ lọc thực tế với yêu cầu đề bài.  Thiết kế: Sử dụng các Function: Sử dụng các Function:
  3. Hàm u_chb1ap: function [b,a] = u_chb2ap(N,Rp,Omegac); % Unnormalized Chebyshev­2 Analog Lowpass Filter Prototype % ­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­ % [b,a] = u_chb1ap(N,Rp,Omegac); %      b = numerator polynomial coefficients %      a = denominator polynomial coefficients %      N = Order of the Elliptic Filter %     Rp = Passband Ripple in dB; Rp > 0 % Omegac = Cutoff frequency in radians/sec [z,p,k] = cheb2ap(N,Rp); a = real(poly(p)); aNn = a(N+1); p = p*Omegac; a = real(poly(p)); aNu = a(N+1);   k = k*aNu/aNn; b0 = k; B = real(poly(z)); b = k*B; Hàm afd_chb2: function [b,a] = afd_chb1(Wp,Ws,Rp,As); % Analog Lowpass Filter Design: Chebyshev­2 % ­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­ % [b,a] = afd_chb2(Wp,Ws,Rp,As); %  b = Numerator coefficients of Ha(s) %  a = Denominator coefficients of Ha(s) % Wp = Passband edge frequency in rad/sec; Wp > 0 % Ws = Stopband edge frequency in rad/sec; Ws > Wp > 0 % Rp = Passband ripple in +dB; (Rp > 0) % As = Stopband attenuation in +dB; (As > 0) if Wp 
  4. ep = sqrt(10^(Rp/10)­1); A = 10^(As/20); OmegaC = Wp; OmegaR = Ws/Wp; g = sqrt(A*A­1)/ep; N = ceil(log10(g+sqrt(g*g­1))/log10(OmegaR+sqrt(OmegaR*OmegaR­1))); fprintf('\n*** Chebyshev­2 Filter Order = %2.0f \n',N) [b,a]=u_chb2ap(N,Rp,OmegaC); Hàm zmapping: function [bz,az] = zmapping(bZ,aZ,Nz,Dz) % Frequency band Transformation from Z­domain to z­domain % ­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­ % [bz,az] = zmapping(bZ,aZ,Nz,Dz) % performs: % b(z) b(Z)| % ---- = ----| N(z) % a(z) a(Z)|@Z = ---- % D(z) % bzord = (length(bZ)­1)*(length(Nz)­1); azord = (length(aZ)­1)*(length(Dz)­1); bz = zeros(1,bzord+1); for k = 0:bzord     pln = [1]; for l = 0:k­1         pln = conv(pln,Nz); end     pld = [1]; for l = 0:bzord­k­1         pld = conv(pld,Dz); end     bz = bz+bZ(k+1)*conv(pln,pld); end az = zeros(1,azord+1); for k = 0:azord     pln = [1]; for l = 0:k­1         pln = conv(pln,Nz); end     pld = [1]; for l = 0:azord­k­1         pld = conv(pld,Dz);
  5. end     az = az+aZ(k+1)*conv(pln,pld); end az1 = az(1); az = az/az1; bz = bz/az1; Hàm freqz_m: function [db,mag,pha,grd,w] = freqz_m(b,a); % Modified version of freqz subroutine % ­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­ % [db,mag,pha,grd,w] = freqz_m(b,a); %  db = Relative magnitude in dB computed over 0 to pi        radians % mag = absolute magnitude computed over 0 to pi radians  % pha = Phase response in radians over 0 to pi radians % grd = Group delay over 0 to pi radians %   w = 501 frequency samples between 0 to pi radians %   b = numerator polynomial of H(z)   (for FIR: b=h) %   a = denominator polynomial of H(z) (for FIR: a=[1]) % [H,w] = freqz(b,a,1000,'whole'); H = (H(1:1:501))'; w = (w(1:1:501))';  mag = abs(H);  db = 20*log10((mag+eps)/max(mag)); pha = angle(H); grd = grpdelay(b,a,w); Hàm dir2cas: function [b0,B,A] = dir2cas(b,a); % DIRECT­form to CASCADE­form conversion (cplxpair version) % ­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­­ % [b0,B,A] = dir2cas(b,a) % b0 = gain coefficient %  B = K by 3 matrix of real coefficients containing bk's %  A = K by 3 matrix of real coefficients containing ak's %  b = numerator polynomial coefficients of DIRECT form %  a = denominator polynomial coefficients of DIRECT form % compute gain coefficient b0 b0 = b(1); b = b/b0; a0 = a(1); a = a/a0; b0 = b0/a0; % M = length(b); N = length(a);
  6. if N > M     b = [b zeros(1,N­M)]; elseif M > N     a = [a zeros(1,M­N)]; N = M; else     NM = 0; end % K = floor(N/2); B = zeros(K,3); A = zeros(K,3); if K*2 == N;     b = [b 0];     a = [a 0]; end broots = cplxpair(roots(b)); aroots = cplxpair(roots(a)); for i=1:2:2*K     Brow = broots(i:1:i+1,:);     Brow = real(poly(Brow));     B(fix((i+1)/2),:) = Brow;     Arow = aroots(i:1:i+1,:);     Arow = real(poly(Arow));     A(fix((i+1)/2),:) = Arow; end Code chương trình: clc; clear; close all; wp = 0.5*pi; ws = 0.3*pi; As = 15; Rp = 1; % Determine digital lowpass cutoff frequencies(wp & ws) wplp = 0.2*pi; alpha = -(cos((wplp+wp)/2))/(cos((wplp-wp)/2)); wslp = angle(-(exp(-j*ws)+alpha)/(1+alpha*exp(-j*ws))); T = 1; Fs = 1/T; OmegaP = (2/T)*tan(wplp/2); OmegaS = (2/T)*tan(wslp/2); %Chebysev Prototype Lowpass Filter % ve do thi thong thap analog [b1,a1]=afd_chb2(wplp,wslp,Rp,As);
  7. [db1,mag1,pha1,w1]=freqs_m(b1,a1,0.5*pi) [ha1,x1,t1]=impulse(b1,a1) figure(1) subplot(2,2,1); plot(w1/pi,mag1);title('Magnitude Response');grid; xlabel('Frequency in pi units');ylabel('|H|'); subplot(2,2,2); plot(w1/pi,db1);title('Magnitude Response in db');grid; xlabel('Frequency in pi units');ylabel('decibels'); subplot(2,2,3); plot(w1/pi,pha1);title('Phase Response');grid; xlabel('Frequency in pi units');ylabel('Radians'); subplot(2,2,4); plot(t1,ha1);title('Impulse Response');grid; xlabel('n');ylabel('h(n)'); %Bilinear transformation [cs,ds] = afd_chb2(OmegaP,OmegaS,Rp,As); [blp,alp] = bilinear(cs,ds,Fs); %Transform digital lowpass into highpass filter Nz = -[alpha,1]; Dz = [1,alpha]; [b,a] = zmapping(blp,alp,Nz,Dz); %From direct to cascade form [C,B,A] = dir2cas(b,a); %Impulse response [h,t] = impz(b,a); %Get Magnitude/Phase Response in db, Group Delay [db,mag,pha,grd,w] = freqz_m(b,a); %Calculate actual Rp,As delta_w = 2*pi/1000; Rp = -(min(db(wp/delta_w+1:1:501))) As = -round(max(db(1:1:ws/delta_w+1))) %Plot figure(3) subplot(3,2,1); plot(w/pi,mag);axis tight;title('Magnitude Response');grid; xlabel('Frequency in pi units');ylabel('|H|'); subplot(3,2,2); plot(w/pi,db);axis([0,1,-30,0]);title('Magnitude Response in db');grid; xlabel('Frequency in pi units');ylabel('decibels'); subplot(3,2,3); plot(w/pi,pha);title('Phase Response');grid; xlabel('Frequency in pi units');ylabel('Radians'); subplot(3,2,4); plot(t,h);title('Impulse Response');grid; xlabel('n');ylabel('h(n)'); subplot(3,2,5); plot(w/pi,grd);axis([0 1 0 10]);title('Group Delay');grid; xlabel('Frequency in pi units');ylabel('Samples');
  8.  Kết quả: Đồ thị:
ADSENSE

CÓ THỂ BẠN MUỐN DOWNLOAD

 

Đồng bộ tài khoản
2=>2