Bài tập lớn xử lý số tín hiệu

Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 1 Project 1 COMPUTING THE DTFT: FINITE-LENGTH SIGNALS Trong project này chúng ta sẽ tập trung xử lý tín hiệu có độ dài hữu hạn. Điều này sẽ sử dụng hàm DTFT cho ở dưới đây. Project làm việc với những tín hiệu dạng xung và phổ tần số DTFT của chúng vì đây là những ví dụ dễ nhất cho tính toán phổ tần số của tín hiệu. Chú ý: Chúng ta cần 2 hàm để tính DTFT. Hàm MATLAB freqz sẽ dùng để xác định trong trường hợp tín hiệu infinite-length, nhưng một hàm mới sẽ cần để tính DTFT của tín hiệu finite-length. Hàm được gọi là DTFT(h,N) function [H,W]=DTFT(h,N) % DTFT calculate DTFT at N equally spaced frequencies % usage: H=DTFT(h,N) % h: finite-length input vector, whose length is L % N: number of frequencies for evaluation over [-pi,pi] % ==> constraint N>=L % H: DTFT values (complex) % W: (2nd output) vector of freqs where DTFT is computed % N=fix(N); L=length(h); h=h(:); if(N<L) error('DTFT: # data samles cannot exceed # freq samples') end W=(2*pi/N)*[0:(N-1)]'; mid=ceil(N/2)+1; W(mid:N)=W(mid:N)-2*pi; W=fftshift(W); H=fftshift(fft(h,N)); Chúng ta không cần phải đưa vào biến độ dài L, vì nó sẽ được xác định bằng tính độ dài vector h. Ngoài ra, khi tín DTFT là hàm tuần hoàn vùng từ là vùng tần số âm, nên phải sử dụng hàm fftshift để sắp xếp lại miền tần số theo thứ tự. Sử dụng DTFT tính vector H, vẽ đồ thị trên miền [- ] với chú ý là H(1) là tần số lấy mẫu ứng với 1.1 DTFT một xung Cho tín hiệu độ rời rạc độ dài hữu hạn sau 1 with 0 n L [ ] 0 elsewhere r n £ £ì = í î a. Tìm phổ tần số của tín hiệu trên: Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 2 2 2 2 -1 2 2 2 2 -1 2 --1 - - - - 0 - - 1 -2 - - 1 2 1 -2 1 2 2 2 1- ( ) [ ]* 1- sin( )( - ) sin( )( - ) sin( ) ( ) sin( ) sin( ) sin ( , ) sin( ) L L L L L iLL i in in i n i i i i i i i i ii L e R e r n e e e Le e e e e e e L R e e a c L q q q q q q q q q q q q qq w w q q q q w +¥ = ¥ = = = = = = = = å å Vậy phổ tần số của tín hiệu đã cho là -1 2 1 -2 1 2 sin( ) ( ) sin( ) Lii LR e e ww w w = Trong công thức trên ta có hàm số hay được dùng trong biến đổi DTFT gọi là hàm aliased sinc : 2 2 sin( ) sin ( , ) sin( ) L a c L w w w = b. Dùng hàm DTFT tính DTFT của tín hiệu xung 12 điểm. Sau đó vẽ đồ thị theo tần số góc ω trong khoảng từ -π tới π. x=ones(1,12); [H,W]=DTFT(x,1024); plot(W,abs(H)); grid on,xlabel('NORMALIZED FREQUENCY'),ylabel('abs(H(w))') đồ thị biên độ c. Tính hàm DTFT với số điểm lẻ L=15 x=ones(1,15); Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 3 [H,W]=DTFT(x,1024); plot(W,abs(H)); grid on,xlabel('NORMALIZED FREQUENCY'),ylabel('abs(H(w))') Qua hai đồ thị biên độ ta thấy các điểm không của hàm asinc(ω,L) là cách đều nhau và biên độ tại đó bằng 0. Đỉnh của biên độ đạt được là 15 d. Xác định các điểm không của hàm asinc(ω,L) Các điểm không thỏa mãn hệ sau : 1 sin( ) 0 2 1 sin( ) 0 2 Lw w ì =ïï í ï ¹ ïî Suy ra 2k L pw = do đó các điểm 0 là cách đều nhau và giá trị của nó được cho bởi công thức trên. 2.2 M-file for asinc , Xây dựng M-file cho hàm asinc Xây dựng hàm asinc(w ,L) để tính giá trị hàm trực tiếp theo công thức giải tích. Hàm gồm có 2 biến đầu vào là : độ dài L và vector tần số . Tại điểm ta sẽ cho giá trị hàm số bằng L. Hàm tính asinc(w ,L) function W = asinc(w,L) L = fix(L); len = length(w); for i = 1:1:len if(w(i)==0) Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 4 W(i) = L; %#ok else W(i) = sin(1/2*L*w(i))./ sin(1/2*w(i)); %#ok end end Sử dụng hàm asinc(w ,L) ta tính trực tiếp phổ tần số của tín hiệu r[n] ban đầu. Sau đó vẽ đồ thị biên độ và so sánh đồ thị thu được với đồ thị thu được từ sử dụng hàm DTFT trên. Tính giá trị DTFT theo công thức trực tiếp: x=ones(1,15); [H,W]=DTFT(x,1024);%lấy vector tần số lấy mẫu. H2=asinc(W,15).*exp(-j*W'*(15-1)/2); plot(W,abs(H)); grid on,xlabel('NORMALIZED FREQUENCY'),ylabel('abs(H(w))') Qua hai đồ thị biên độ tính theo hai cách trực tiếp và qua hàm DTFT ta thấy kết quả thu được từ 2 phương pháp là đồng nhất. Project 2 DTFT SYMMETRIE, DTFT ĐỐI XỨNG Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 5 Các tín hiệu finite-length thường được sử dụng có tính chất đối xứng để dễ tính DTFT hơn. Ví dụ các xung hình chữ nhật là chẵn xung quanh điểm n=0 thì tính toán DTFT của nó sẽ là hàm thực ,chẵn theo đối số . Trong project 2, sáu loại đối xứng chính sẽ được xem xét và minh họa bởi các ví dụ khác nhau. Chú ý : Cần thiết phải thay đổi hàm DTFT dể có thể chấp nhân một biến mới xác định thời điểm bắt đầu của tín hiệu x[n]. Thông thường các tín hiệu bắt đầu được tính tại thời điểm n=0 bởi hàm DTFT trong Project 1. Một kĩ năng cần thiết nữa trong Project 2 là vẽ đồ thị DTFT để biểu diễn tính đối xứng. Thông thương, các tần số góc được xác định trong khoảng từ 0 tới 2p . Tuy nhiên dể nghiên cứu tính đối xứng của đồ thị DTFT chúng ta cần phải đưa chúng các tần số góc về khoảng từ p p- ® . Khi kiểm tra các đặc tính thực và ảo của biến đổi có thể phần tín hiệu là bằng không lại không thực sự bằng 0. Điều đó có thể xảy ra do một số các giá trị vô cùng nhỏ trong quá trình làm tròn. Nếu các giá trị đó nằm trong khoảng từ thì có thể bỏ qua và kết luận là do làm tròn số gây ra. Tuy nhiên nếu cỡ thì chưa chắc đã là do làm tròn, trừ khi đã có một số lượng lớn các phép tính đã được thực hiện. Khi kiểm tra biến đổi là chẵn hay lẻ, nên sử dụng một chương trình để dịch chuyên biến đổi. Nghĩa là ta cần phải có một hàm MATLAB mới để thực hiện tương tự với biểu diễn toán hạng của dịch chuyển tần sô ( ) ( )i iY e X ew w-= . Hàm số sẽ thực hiện biến đổi với cặp vector (H,W). function [G, Wflipped ]= flipDTFT(H,W) % FlipDTFT flip the DTFT G(w)=H(-w) % usage: % [ G,Wflipped] =flipDTFT(H,W) % H=DTFT values complex % W=frequency samples % G= DTFT values % Wflipped= flipped frequency domain lies within [-pi,pi) % N=length(H); Wflipped=-W(N:-1:1); G=H(N:-1:1); %----- % now get everything back into the [-pi,pi) interval % assume that W was monotonically increasing % so Wflipped is also increasing % jkl=find(Wflipped(:)'<-pi); if(~isempty(jkl)) kk=[(length(jkl)+1):N jkl]; Wflipped(jkl)=Wflipped(jkl)+2*pi; Wflipped=Wflipped(kk); G=G(kk); end jkl=find(Wflipped(:)'>=(pi-100*eps)); if(~isempty(jkl)) kk=[jkl 1:(jkl(1)-1)]; Wflipped(jkl)=Wflipped(jkl)-2*pi; Wflipped=Wflipped(kk); G=G(kk); end 2.1 Zero—Phase Signals, Tín hiệu pha 0 Xử lý các tín hiệu zero-phase là một điều khó khăn trong MATLAB vì hàm DTFT trong Project 1 cho rằng tín hiệu bắt đầu từ n=0 trong khi các tín hiệu zero-phase phải có tính đối xứng xung quanh điểm n=0. Một cách để xử lý vấn đề này là ta thay đổi hàm DTFT bằng cách thêm vào 1 Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 6 biến mới xác định thời điểm bắt đầu của tín hiệu. Thời điểm ban đầu này sẽ làm thay đổi giá trị hàm DTFT theo hệ thức 0 0[ ] ( ) i nDTFT ix n n e X ew w-- ¬¾¾® a. Xây dựng hàm DTFT_n0(x,n0,N) để thực hiện sự dịch chuyển thời gian của tín hiệu bằng biến n0. function [H,W]=DTFT_n0(h,n0,N) %#ok %-- % usage [H,W]=DTFT_n0(h,n0,N) % h finite-length input vector % n0 xac dinh thoi diem bat dau tin hieu % N number of freuencies for evaluation over [-pi,pi] % H DTFT_n0 values % W (2nd output) vector of frequencies for evaluation over % [-pi,pi] % N=fix(N); [X,W]=DTFT(h,N); H=exp(i*n0*W).*X; %#ok b. Kiểm tra DTFT_n0 bằng cách tính DTFT của một xung 21 điểm bắt đầu tại n=-10. Kết quả sẽ là một hàm thực chẵn Ta sử dụng tín hiệu 1 with -1 n 1[ ] 0 elsewhere r n £ £ì = í î x=ones(1,21); [H,W]=DTFT_n0(x,10,25); H Kết quả tính với r[n] ,n0=10,N=25 như sau : n=1 0.2492 + 0.0000i n=13 21.0000 n=2 -0.6969 - 0.0000i n=14 3.8438 - 0.0000i n=3 1.0000 + 0.0000i n=15 -3.3951 + 0.0000i n=4 -1.0856 - 0.0000i n=16 2.7111 + 0.0000i n=5 0.9126 + 0.0000i n=17 -1.8782 + 0.0000i n=6 -0.4778 - 0.0000i n=18 1.0000 - 0.0000i n=7 -0.1831 - 0.0000i n=19 -0.1831 - 0.0000i n=8 1.0000 + 0.0000i n=20 -0.4778 - 0.0000i n=9 -1.8782 - 0.0000i n=21 0.9126 - 0.0000i n=10 2.7111 + 0.0000i n=22 -1.0856 + 0.0000i n=11 -3.3951 - 0.0000i n=23 1.0000 - 0.0000i n=12 3.8438 + 0.0000i n=24 -0.6969 - 0.0000i n=25 0.2492 + 0.0000i Theo bảng số liệu ta thấy phổ tần số là hoàn toàn thực và giá trị phần ảo là bằng 0. Đồng thời các giá trị thu được cũng đối xứng qua điểm n=0 nen là hàm chẵn.Điều đó phù hợp với lý thuyết. c. Vẽ đồ thị của DTFT trên Đồ thị phần thực Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 7 Đồ thị biên độ Vậy từ đồ thị ta lại thấy rằng DTFT của tín hiệu r[n] là chẵn và đối xứng, phù hợp với lý thuyết. Giá trị phần ảo cỡ không hoàn toàn bằng không nhưng đó là do sai số làm tròn trong quá trình tính toán. Chúng ta hoàn toàn có thể bỏ qua nó mà không ảnh hưởng tới kết quả thu được d. Chứng minh phần ảo là bằng 0 và pha biến thien từ 0 tới p 1 1 ( ) [ ]* 1 ( ) 1 2 cos( ) L i in in L L L in in R e r n e e e e n w w w w w w +¥ - - -¥ - - = = = + + = + å å å å Từ chứng minh trên ta thấy phần ảo là bằng 0. Pha của ( )iR e w biến thiên từ 0 của thành phần cố định 1 tới L . Mặt khác Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 8 2* k N pw = 1... 1k N= - Do dó pha của ( )iR e w biến thiên từ 0 tới Cũng cần phải chú ý là sự đối xứng này chỉ có thể quan sát khi số xung là lẻ, nếu số xung là chẵn thì còn phải kể tới một sự trễ pha của nửa tín hiệu. 2 .2 Triangular Pulse Cho một tín hiệu đối xứng khác dạng tam giác như sau L - n with 0 n< [ ] L+n w -L<n<0 0 es L n ith lewhere £ì ïD = í ï î Độ dài của tín hiệu là 2L-1 và có thể xác định bằng tích chập của 2 xung hình chữ nhật có L điểm. Kết quả của DTFT là một hàm asinc được làm vuông, pha bằng 0 vì đây là một tín hiệu đối xứng. a. Vẽ đồ thị 21 điểm xung trong khoảng 20 20n- £ £ . Sau đó tính DTFT của nó bằng hàm DTFT_n0 và vẽ đồ thị kết quả trong khoảng p w p- £ < Đồ thị: Đồ thị phần ảo của DTFT Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 9 Đồ thị phần thực của DTFT Đồ thị biên độ của DTFT Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 10 2.3 Symmetries in th DTFT Có nhiều các đặc tính đối xứng trong miền thời gian và tần số. Một nhóm các đặc tính dùng cho các hàm thuần thực hoặc thuần ảo, những thuộc tính khác cho hàm chẵn hoặc hàm lẻ. Ví dụ, DTFT của hàm chẵn là hàm chẵn. Mỗi ví dụ dưới đây sẽ tập trung vào một loại đối xứng. Đối xứng trong miền tần số có thể được thể hiện bằng vẽ đồ thị phần thực và phần ảo (hoặc biên độ và pha) va sử dụng thêm hàm flipDTFT để kiểm tra . a. DTFT của một tín hiệu thực là một đối xứng liên hợp, , nghĩa là liên hợp của DTFT giống như là dịch chuyển nó. Ví dụ với tín hiệu x[n]= với . Vẽ đồ thị biên độ và pha của DTFT và chứng tỏ rằng nó giống liên hợp phức đối xứng. So sánh biên độ a=0.9; nn=[0:1:21];x=a.^nn.*cos(2*pi*nn/sqrt(31)); %#ok [H,W]=DTFT(x,1024); [Y,Wflip]=flipDTFT(H,W); subplot(211), plot(W,abs(H)) xlabel('NORMALIZED FREQUENCY'),ylabel('Abs(H(w))'),grid on subplot(212),plot(W,abs(Y)) xlabel('NORMALIZED FREQUENCY'),ylabel('Abs of flipDTFT (H(w))'),grid on So sánh pha a=0.9; nn=[0:1:21];x=a.^nn.*cos(2*pi*nn/sqrt(31)); %#ok Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 11 [X,W]=DTFT(x,1024); [Y,Wflip]=flipDTFT(X,W); subplot(211), plot(W,angle(X)) xlabel('NORMALIZED FREQUENCY'),ylabel('Arg(H(w))'),grid on subplot(212),plot(W,angle(Y)) xlabel('NORMALIZED FREQUENCY'),ylabel('Arg of flipDTFT (H(w))'),grid on Dựa vào đồ thị ta có thể thấy đồ thị biên độ của DTFT và đồ thị liên hợp đối xứng là bằng nhau, đồ thị pha của DTFT và của liên hợp phức đối xứng là đối nhau qua trục tung . Đó chính là tính chất cần tìm. b. Nếu tín hiệu là thuần ảo, DTFT của nó sẽ là liên hợp phản đối xứng. Sử dụng hàm x[n] từ phần trên, định nghĩa hàm y[n] =jx[n]. So sánh biên độ a=0.9; nn=[0:1:21];x=j*a.^nn.*cos(2*pi*nn/sqrt(31)); %#ok [X,W]=DTFT(x,1024); [Y,Wflip]=flipDTFT(X,W); subplot(211), plot(W,abs(X)) xlabel('NORMALIZED FREQUENCY'),ylabel('Abs(H(w))'),grid on subplot(212),plot(W,abs(Y)) xlabel('NORMALIZED FREQUENCY'),ylabel('Abs of flipDTFT (H(w))'),grid on Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 12 So sánh phase a=0.9; nn=[0:1:21];x=j*a.^nn.*cos(2*pi*nn/sqrt(31)); %#ok [X,W]=DTFT(x,1024); [Y,Wflip]=flipDTFT(X,W); subplot(211), plot(W,angle(X)) xlabel('NORMALIZED FREQUENCY'),ylabel('Arg(H(w))'),grid on subplot(212),plot(W,angle(Y)) xlabel('NORMALIZED FREQUENCY'),ylabel('Arg of flipDTFT (H(w))'),grid on Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 13 Ta thấy biên độ của DTFT và của liên hợp phức phản đối xứng là bằng nhau. pha của chúng đối nhau. Đó là tính chất cần tìm. c. Một hàm chẵn trong miền thời gian cũng cho một hàm chẵn trong miền tần số. sừ dụng tín hiệu trong miền là chẵn. Tính DTFT của hàm trên và chứng minh nó chẵn với Đồ thị phần thực nn=[-29:1:29]; x=exp(i*2*pi*nn.^2/25); [H,W]=DTFT_n0(x,30,1024); plot(W,real(H)) xlabel('NORMALIZED FREQUENCY'),ylabel('Real(H(w))'),grid on Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 14 Đồ thị phần ảo nn=[-29:1:29]; x=exp(i*2*pi*nn.^2/25); [H,W]=DTFT_n0(x,30,1024); plot(W,imag(H)) xlabel('NORMALIZED FREQUENCY'),ylabel('Im(H(w))'),grid on Từ đồ thị ta thấy phần thực và ảo của DTFT đều đối xứng qua trục tung do đó đây là hàm chẵn. phù hợp với các tính toán lý thuyết d. Với tín hiệu lẻ trong khoảng , thì DTFT của nó sẽ là hàm lẻ theo Đồ thị phần thực x=[-19:1:19]; Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 15 [H,W]=DTFT_n0(x,19,1024); plot(W,real(H)) xlabel('NORMALIZED FREQUENCY'),ylabel('Real(H(w))'),grid on Đồ thị phần ảo x=[-19:1:19]; [H,W]=DTFT_n0(x,19,1024); plot(W,imag(H)) xlabel('NORMALIZED FREQUENCY'),ylabel('Im(H(w))'),grid on Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 16 Từ đồ thị ta thấy phần thực của DTFT có thể coi là bằng 0 còn phần ảo thì đối xứng qua gốc tọa độ. Chứng tỏ DTFT là hàm lẻ theo . e. Các tính chất chẵn/lẻ và thực/ảo có thể kết hợp lại. xét tín hiệu là hàm ảo ,lẻ. DTFT của nó sẽ là một hàm phản liên hợp và lẻ tương ứng. Kiểm tra bằng đồ thị. nn=[-20:1:20]; x=j*nn.*cos(nn); [H,W]=DTFT_n0(x,20,1024); plot(W,real(H)),grid on xlabel('NORMALIZED FREQUENCY'),ylabel('Real(H(w))'),grid on plot(W,imag(H)),grid on xlabel('NORMALIZED FREQUENCY'),ylabel('Im(H(w))'),grid on Đồ thị biên độ Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 17 Đồ thị pha Sử dụng hàm flipDTFT kiểm tra nn=[-20:1:20]; x=j*nn.*cos(nn); [H,W]=DTFT_n0(x,20,1024); [H,G]=flipDTFT(X,W); plot(G,real(H)),grid on xlabel('NORMALIZED FREQUENCY'),ylabel('Real of (H(w))'),grid on plot(G,imag(H)),grid on xlabel('NORMALIZED FREQUENCY'),ylabel('Im of (H(w))'),grid on Đồ thị biên độ Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 18 Đồ thị pha Từ đồ thị ta thấy 1. Biên độ DTFT và biên độ flipDTFT là bằng nhau, pha DTFT đối xứng với phần ảo flipDTFT qua trục tung. Đây là tính chất của hàm thuần ảo (câu C) 2. Biên độ của DTFT và của flipDTFT là bằng nhau,pha của chúng là đối nhau. Đó là tính chất của hàm lẻ -thời gian (câu d) Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 19 Project 3 : DTFT TDTFT OF INFINITE-LENGTH SIGNAL 3.1 Tín hiệu mũ. Sử dụng hàm freqz để tính toán DTFT ( )jX e w của tín hiệu [n] = (0.9) [n]nx u . 3.1.a Vẽ đồ thị biên độ, pha của theo tần số w trong khoảng [ ],p p- - Trước hết ta lấy vector giá trị biến đổi HH, tần số lấy mẫu WW (trên đường tròn đơn vị [ ]0,2p ) dung hàm freqz(b,a,n,’whole’); [n] = [n]nx a u DTFT¬¾¾® 1 ( ) 1 j j X e ae w w-= - (3.1) a = [1,-0.9]; b = [1]; N = 512; [HH,WW] = freqz(b,a,N,'whole'); subplot (221);plot(WW,abs(HH)); grid, title('MAGNITUDE RESPONSE BY FREQZ'),xlabel('NORMALIZED FREQUENCY'), ylabel('|X(e^-jw)|') subplot(222);plot(WW,angle(HH)*180/pi); grid, title('PHASE RESPONSE BY FREQZ' ),xlabel('NORMALIZED FREQUENCY'), ylabel('PHASE (DEGREE)') - Đồ thị biên độ, pha theo tần số w trong khoảng [ ]0,2p - Biến đổi HH, WW : chuyển vector tần số lấy mẫu từ khoảng [ ]0,2p về [ ],p p- mid = ceil (N/2) + 1; WW(mid:N)= WW(mid:N) - 2*pi; WW = fftshift(WW); Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 20 HH2 = [HH(mid : N)' HH(1 : mid -1)']; HH = HH2'; - Vẽ lại đồ thị biên độ và đò thị pha : subplot (211);plot(WW,abs(HH)); grid, title('MAGNITUDE RESPONSE BY FREQZ'),xlabel('NORMALIZED FREQUENCY'), ylabel('|X(e^-jw)|') subplot(212);plot(WW,angle(HH)*180/pi); grid, title('PHASE RESPONSE BY FREQZ' ),xlabel('NORMALIZED FREQUENCY'), ylabel('PHASE (DEGREE)') 3.1.b Biên độ là hàm chẵn và pha là hàm lẻ của w : 1 ( ) 1 j jX e ae w v-= - ÞBiên độ 2 1 ( ) 1 2 cos( ) jX e a a w w é ù =ë û - + ; pha : a sin ( ( )) arctan( ) 1 cos jX e a w wf w = - - 3.1.c Ta vẽ lại đồ thị biên độ, pha của tín hiệu ra theo tần số sử dụng công thức 3.1.b AA = 1./sqrt(1 - 2*0.9*cos(WW) + 0.9^2); PP = atan((-0.9*sin(WW))./(1-0.9*cos(WW))); subplot (223);plot(WW,AA); grid, title('MAGNITUDE RESPONSE BY FORMULAS'),xlabel('NORMALIZED FREQUENCY'), ylabel('|H(e^-jw)|') subplot (224);plot(WW,PP*180/pi); grid, title('PHASE RESPONSE BY FORMULAS'),xlabel('NORMALIZED FREQUENCY'), ylabel('PHASE (DEGREE)') Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 21 - Kết quả thu được hoàn toàn giống với kết quả sử dụng freqz 3.2 Số mũ phức. Nếu cho hệ số a trong công thức (3.1) bởi số phức 0 ja z re q= = thì sự biến đổi vẫn không đổi. Điều này rất có ý nghĩa trong việc xác định tác động của biên độ r và pha q của số phức lên DTFT. 3.2.a Cho 3 11 0 0.95 j z e p = , 0[n] = z [n] nx u với 0 n 30£ £ . Ta vẽ song song 2 đồ thị phần thực & phần ảo của [n]x theo n trên 2 panel: N = [0 : 30]; z = 0.95 * exp(sqrt(-1) * 3 * pi /11); Z = z.^N.*(N>0); subplot(211); plot(N,real(Z)); grid, title('REAL PART OF x[n] VESUS n'), xlabel ('n'), ylabel('REAL(x[n])') subplot(212); plot(N,imag(Z)); grid, title('IMAGINARY PART OF x[n] VESUS n'), xlabel ('n'), ylabel('IMAG(x[n])') Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 22 3.2.b Tiếp tục với tín hiệu 0[n] = z [n]nx u , 3 11 0 0.95 j z e p = , ta tính DTFT và vẽ đồ thị biên độ tín hiệu ra theo tần số w z = 0.95 * exp(sqrt(-1) * 3 * pi /11); a = [1,-z]; b = [1]; N = 512; [HH,WW] = freqz(b,a,N,'whole'); plot(WW,abs(HH)); grid, title('MAGNITUDE RESPONSE BY FREQZ'),xlabel('NORMALIZED FREQUENCY'), ylabel('|H(e^-jw)|') Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 23 - Trên Matlab ta có thể xác định được tọa độ đỉnh của đồ thị : %peak location of magnitude MaxMagnitue = max(abs(HH)) PeakLocation = WW(find(abs(HH) == MaxMagnitue)) - Kết quả thu được : MaxMagnitue = 19.9811 and PeakLocation = 0.8590. - Ta có : ( ) 1 1 [n] = ( ) [n] ( ) 1 1 DTFTj n j j j jx re u X e re e re q w q v q w- -¬¾¾® = =- - Biên độ : 2 1 ( ) 1 2 cos( ) jX e r r w w q = - - + (3-2-b)lớn nhất bằng : 1 1 r- khi có w q= Với 3 110 0.95 j z e p = thì ( ( ) )jmax X e w = 20 đạt được khi 0( ) 0.8590angle zw = = . Điêu này phù hợp với kết quả tính DTFT bằng freqz. 3.2.c Với 3 50 0.95 j z e p = , từ công thức (3-2-b) ta có thể mô tả đồ thị biên độ của tín hiệu ra sau khi biến đổi DTFT như sau: - Đồ thị có đỉnh ( ( ) )jmax X e w = 20 tại 0( ) 1.8850angle zw = = . - Tại các giá trị tấn số ngoài lân cận q biên độ giảm hẳn, độ lớn xâp xỉ 0.5 - Điều này được kiểm chứng băng việc vẽ đồ thị biên độ tín hiệu ra từ hàm freqz z = 0.95 * exp(sqrt(-1) * 3 * pi/5); a = [1,-z]; b = [1]; N = 512; [HH,WW] = freqz(b,a,N,'whole'); plot(WW,abs(HH)); grid, title('MAGNITUDE RESPONSE BY FREQZ'),xlabel('NORMALIZED FREQUENCY'), ylabel('|H(e^-jw)|') Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 24 %peak location of magnitude MaxMagnitue = max(abs(HH)) PeakLocation = WW(find(abs(HH) == MaxMagnitue)) MaxMagnitue = 19.9091 and PeakLocation = 1.8890 3.2.d Thay đổi biên độ của số phức r = 0.975, 0.95, 0.9, 0.8 , ta thực hiện vẽ đồ thị biên độ tín hiệu ra từ biến đổi DTFT z = 0.975 * exp(sqrt(-1) * 3 * pi /5);a = [1,-z];b = [1];N = 512; [HH,WW] = freqz(b,a,N,'whole'); subplot(221);plot(WW,20*log10(abs(HH))); grid, title('MAGNITUDE RESPONSE WITH r = 0.975'),xlabel('FREQUENCY (Hertz)'), ylabel('MAGNITUDE (dB)') z = 0.95 * exp(sqrt(-1) * 3 * pi /5);a = [1,-z];b = [1];N = 512; [HH,WW] = freqz(b,a,N,'whole'); subplot(222);plot(WW,20*log10(abs(HH))); grid, title('MAGNITUDE RESPONSE WITH r = 0.95'),xlabel('FREQUENCY (Hertz)'), ylabel('MAGNITUDE (dB)') z = 0.9 * exp(sqrt(-1) * 3 * pi /5);a = [1,-z];b = [1];N = 512; [HH,WW] = freqz(b,a,N,'whole'); subplot(223);plot(WW,20*log10(abs(HH))); grid, title('MAGNITUDE RESPONSE WITH r = 0.9'),xlabel('FREQUENCY (Hertz)'), ylabel('MAGNITUDE (dB)') z = 0.8 * exp(sqrt(-1) * 3 * pi /5);a = [1,-z];b = [1];N = 512; [HH,WW] = freqz(b,a,N,'whole'); subplot(224);plot(WW,20*log10(abs(HH))); grid, title('MAGNITUDE RESPONSE WITH r = 0.8'),xlabel('FREQUENCY (Hertz)'), ylabel('MAGNITUDE (dB)') Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 25 - Độ rộng tại -3dB: + r = 0.975 , wD = 3.485 – 0.2823 = 3.2027; + r = 0.95 , wD = 3.51 - 0.2577 = 3.2523; + r = 0.9 , wD = 3.559 - 0.2086 = 3.3504; + r = 0.8 , wD = 3.682 - 0.0859 = 3.5961; - Xây dựng công thức lien hệ độ rộng tại -3dB với r : 2 2 2 0.3 2 0.3 1 3 20lg( ) 1 2 cos( ) lg(1 2 cos( ) ) 0.3 1 10 os( ) 2 1 10 2arccos( ) 2 dB r r r r r c r r r w q w q w q w - = - - + Û - - + = + - Û - = + - Þ D = 3.3 Decaying Signusoid. Tính DTFT của tín hiệu thực : 2[ ] = 3(0.95) cos( ) [ ]7 3 n ny n u np p+ bằng hàm freqz và tính liên hợp phức ; 3.3.a Ta đã có biến đổi DTFT : 2 2 cos cos( ) [ ]= Ar cos( ) [ ] ( ) 1 2 cos j DTFTn j j j r e y n n u n Y e A r e r e w w w w f q f q f q - - - æ ö- - + ¬¾¾® = ç ÷- +è ø (3-6) Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 26 Sử dụng hàm [Y,W] = freqz(b,a,n) để lấy ra vector tần số mẫu W. Với A = 3,r = 0.95, 2 / 7q p= , / 3f p= ta tính được b = [3*cos(pi/3) , -3*0.95*cos(2*pi/7 - pi/3)] = [1.5 1.0412] a = [1, -2*0.95*cos(2*pi/7), 0.95^2] = [1.0000 -1.1846 0.9025] b = [3*cos(pi/3) , -3*0.95*cos(2*pi/7 + pi/3)]; a = [1, -2*0.95*cos(2*pi/7), 0.95^2] ; n = 512; [Y,W] = freqz(b,a,n); - Đồ thị biên độ, pha của tín hiệu ra theo tần số w subplot (211);plot(W,abs(Y)); grid, title('MAGNITUDE RESPONSE BY FREQZ'),xlabel('NORMALIZED FREQUENCY'), ylabel('|Y(e^-jw)|') subplot(212);plot(W,angle(Y)*180/pi); grid, title('PHASE RESPONSE BY FREQZ' ),xlabel('NORMALIZED FREQUENCY'), ylabel('PHASE (DEGREE)') 3.3.b Ta có - Tách tín hiệu /3 2 /7 /3 2 /73 32[ ] 3(0.95) cos( ) [ ] (0.95 ) [ ] (0.95 ) [ ]7 3 2 2 n j j n j j nny n u n e e u n e e u np p p pp p - -= + = + = ½( [n]x + *[n]x ) * *[n] ( ) [n] ( )DTFT DTFTj jx X e x X ew w-¬¾¾® Þ ¬¾¾® (3-3); Þ * 1 ( ) ( ) ( ) 2 j j jY e X e X ew w w-é ù= +ë û - Sử dụng freqz tính DTFT của tín hiệu [n]x : [XX,W1] = freqz([3*exp(sqrt(-1)*pi/3)], [1, -0.95 * exp(2*pi*sqrt(-1)/7)],n*2,'whole'); Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 27 X1 = XX(1:n); %lấy ra tín hiệu ra ( )jX e w trong khoảng tần số [0, ]p X2 = XX(n+1:2*n); %lấy ra tín hiệu ra *( )jX e w- trong khoảng tần số [0, ]p Temp = real(X2) - imag(X2)*sqrt(-1); X2 = Temp; Y2 = 1/2*(X1 + X2(n:-1:1)); subplot (211);plot(W,abs(Y2)); grid, title('MAGNITUDE RESPONSE BY VIA THE CONJUGATE’),xlabel(‘FREQUENCY'), ylabel('|Y(e^-jw)|') subplot(212);plot(W,angle(Y2)*180/pi); grid, title('PHASE RESPONSE BY VIA THE CONJUGATE’),xlabel(‘FREQUENCY'), ylabel('PHASE (DEGREE)') Đồ thị biên độ và đồ thị pha thu được từ 2 cách tính DTFT là giống nhau. Về số liệu sử dụng matlab kiểm tra : abs(Y2)-abs(Y), angle(Y2)-angle(Y) Thì sai số thu được từ 2 phương pháp tính là không đáng kể. Project 4 WINDOWING FOR DTFT Trong project này, có hai đặc tính của DTFT được minh họa : tính điều chế và đặc tính cửa sổ. Tính chất điều biến, có hiệu lực trong 1 trường hợp đặc biệt của hàm cửa sổ ở chỗ miền tần số chập lại tới 1 tần số đơn giản. 4.1 Modulation Property Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 28 Nhiều đặc tính DTFT có ý nghĩa và ứng dụng hữu ích, một trong số đó là biến điệu biên độ, có thể tìm thấy trong truyền thông và radar. Nếu tín hiệu x[n] được nhân với một hàm sin phức , thì kết quả biến đổi là một sự dịch chuyển ; trở thành . a. Chứng minh khẳng định trên với tín hiệu xung hình chữ nhật với độ dài L=21 và tần số lọc . Vẽ đồ thị kết quả nhờ DTFT. Kiểm tra đỉnh của biên độ DTFT đã được dịch chuyển một đoạn . Thử các giá trị để biểu diễn tính tuần hoàn của DTFT Xét tín hiệu xung chữ nhật định nghĩa là r[n]= Có phổ tần số là : -12 1 -2 1 2 sin( ) ( ) sin( ) lii lR e e ww w w = Xét tín hiệu 0[ ] [ ]* inx n r n e w= Có phổ tần số là : 0 0 0 0 ( )1 1 ( ) ( ) 0 0 1 ( ) [ ]* * 1 ill l in ini in in i e X e x n e e e e e w w w w ww w w w w - -+¥ - - - -- - - - -¥ - = = = = -å å å 0 0 0 0 0 0 0 ( ) ( ) ( ) 0 ( 1)( )2 2 2 2 ( ) ( ) ( ) 02 2 2 ( ) sin( ) 2* ( ) sin( ) 2 il il il i l i i i l e e e e e e e w w w w w w w w w w w w w w w w w w - - - - - - - - - - - - - - = = - - Vậy ta đã có 0( )( ) ( )iiX e R e w ww -= Với tín hiệu đang xét ta có ,L=21 (l=20) : Đồ thị biên độ ( )iR e w r=[0:20/1024:20]; l=21;omega0=2*pi/sqrt(31);[H,W]=DTFT(r,1025); R=sin(l*(W-omega0)/2)./sin((W-omega0)/2).* exp(j*(l-1)*(W-omega0)/2); plot(G,abs(R) Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 29 từ đồ thị ta thấy đỉnh của DTFT tín hiệu [ ]x n đã dịch sang phải một khoảng 0 1.129w = Với ta có đồ thị mới r=[0:20/1024:20]; l=21;omega0=2*pi/sqrt(31);[H,W]=DTFT(r,1025); W=W+2*pi; R=sin(l*(W-omega0)/2)./sin((W-omega0)/2).* exp(j*(l-1)*(W-omega0)/2); plot(G,abs(R)),grid on Từ đồ thị ta thấy đồ thị của DTFT tín hiệu cũng giống đồ thị của DTFT tín hiệu chứng tỏ DTFT tuần hoàn với chu kì Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 30 b. Lặp lại nhưng nhân với tín hiệu cosin cùng tần số x[n]= r[n]*cos(n) 1 1 0 0 1 ( ) [ ]* cos( )* ( )* 2 l l i in in in in inX e x n e n e e e ew w w w +¥ - - - - - - -¥ = = = +å å å ( 1) ( 1)1 ( 1) ( 1) ( 1) ( 1) 0 1 1 1 1 ( ) ( ) 2 2 1 1 il ill in in i i e e e e e e w w w w w w - - - +- - - - + - - - + - - = + = + - -å 1 1 ( 1) ( 1) 2 2 ( 1) ( 1) sin( sin(1 2 2( ) 1 12 sin( ) sin( ) 2 2 l l i i l l e e l l w w w w w w - - - - - + - + = + - - Đây chính là điều chế biên độ hai dải biên. Theo kết quả thu được ta dự đoán đồ thị sẽ là đối xứng nhau qua trục tung và có biên độ bằng 0.5 so với đồ thị DTFT của r[n]. Ngoài ra còn có hiện tượng điều dịch chuyển tần số về 2 phía trục hoành một khoảng là 1. Đồ thị biểu diễn biên độ. 4.2 Windowing Gives Frequency-Domain Convolution Đặc tính cửa sổ của DTFT của tích hai tín hiệu trong miền thời gian là tích chập trong miền tần số tuần hoàn của tích phân Fourier của chúng. a. Không thể tính được tích chập trong miền tần số theo công thức trên vì các giá trị của chúng ta là rời rạc. Tuy nhiên, một trường hợp mà kết quả trong miền tần số có thể xác định như một công thức – khi các tín hiệu không bị tác động của hàm cửa sổ là một hàm sin phức (ví dụ x[n]= . Khi đó là một xung trong tần số.,sự nhân chập suy ra . Đó là tính chất điều chế. Điều này dẫn tới biểu diễn DTFT của hàm cửa sổ dịch chuyển theo tần số. Tạo hàm cửa sổ hình sin với : Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 31 ở đây r[n] là xung hình chữ nhật có độ dài là L=32. Xung hình chữ nhật này có thể tạo nhờ hàm ones hoặc boxcar. Vẽ DTFT và chú ý rằng đỉnh bị dịch chuyển tới . Đồ thị biên độ Đồ thị biên độ teta=2*pi/sqrt(31); x=ones(1,32); y=[0:1:31]; z=exp(i*teta*y).*x; [H,W]=DTFT(z,1024); plot(W,abs(H)),grid on; xlabel('NORMALIZED FREQUENCY'),ylabel('Abs(H(w))'),grid on Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 32 Từ đồ thị ta thấy đỉnh biên độ đã dịch chuyển một khoảng . Tín hiệu xung hình sin tác động vào tín hiệu w[n]=r[n] đã cho ta tín hiệu y[n] mà DTFT của nó đã dịch chuyển một khoảng phù hợp với công thức . Đồng thời biên độ vẫn như nhau vì b. Cho hàm cửa sổ sau gọi là cửa sổ Hann. Áp dụng hàm cửa sổ ở trên vào 1 hàm sin có tần số là . Vẽ đồ thị thời gian và tính DTFT của nó sau đó vẽ đồ thị biên độ. Đồ thị trong miền thời gian Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 33 Đồ thị biên độ DTFT của tín hiệu w[n] nn=[0:1:31]; x=1/2-1/2*cos(2*pi/32*nn); [H,W]=DTFT(x,1024); plot(W,abs(H)); xlabel('NORMALIZED FREQUENCY'),ylabel('Abs(H(w))'),grid on Đồ thị DTFT của tín hiệu x[n]=w[n]. nn=[0:1:31]; x=(1/2-1/2*cos(2*pi/32*nn)).*sin(2*pi/sqrt(31)*nn); [H,W]=DTFT(x,1024); Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 34 plot(W,abs(H)); xlabel('NORMALIZED FREQUENCY'),ylabel('Abs(H(w))'),grid on Do tác động vào tín hiệu w[n] một tín hiệu dạng sin nên tín hiêu thu được sẽ bị dịch chuyển về 2 phía của trục hoành một khoảng bằng tần số của tín hiệu hình sin. Đồng thời biên độ của tín hiệu cũng bị giảm đi một nửa vì c. DTFT của hàm cửa số Hann có thể được viết bởi 3 hàm asinc. Điều này được thực hiện bằng cách coi như một cửa sổ chữ nhật tác động vào tín hiệu và sau đó sử dụng đặc tính điều biến. Ta thấy đồ thị biên độ được tạo ra từ đồ thị của tín hiệu xung chữ nhật như sau: Ta chia tín hiệu thành 3 tín hiệu khi đó w[n]=w1[n]+w2[n]+w3[n]. Do đó . Biên độ của sẽ bằng ½ biên độ của tín hiệu xung chữ nhật còn pha thì giống nhau, trong khi đó biên độ của và bằng ¼ của và pha của chúng thì lần lượt bị dịch chuyển một khoảng . Ta biểu diễn chúng dưới dạng số phức như sau : , , , Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 35 . Khi đó 2 2 2 ( ) ( ) ( ) (1 os( )) 2 4 4 2 i R R R RR e c L L L w p p py y y y= Ð + Ð - + Ð + = + Ð Cho ta thấy ( )iR e w cùng pha với DTFT của tín hiệu w[n]. xem hình trang 6 d. Khi chúng ta cho tín hiệu Sự thay đổi dấu trừ thành dấu cộng ở đây là quan trọng. Vẽ đồ thị trong miền thời gian theo n với ở đây L=32. Đồ thị biên độ DTFT Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 36 4.3 Convergence to True DTFT Còn một trường hợp mà chưa được đề cập đến : trường hợp dùng DTFT cho1 phần tử hữu hạn của một tín hiệu vô hạn. Trường hợp này thường được gặp trong thực tế vì chúng ta thường ghi lại (hoặc lưu trữ) một phần nhỏ của tín hiệu để phân tích. Dù sao, chúng ta cũng muốn suy ra tín hiệu đúng từ 1 đoạn tín hiệu giới hạn. Hàm DTFT là đủ để biểu diễn DTFT của mọi tín hiệu hữu hạn. Một ví dụ đơn giản của cửa sổ nhận được bằng việc lấy ra L điểm đầu tiên của lũy thừa dài vô han . Khi độ dài tăng lên thì kết quả sẽ hội tụ đến dạng giải tích được đưa ra trước đó cho . Với a=0.977, vẽ log của độ lớn DTFT cho vài độ dài khác nhau L=32,64,128 và 256. a=0.977;L=32; nn=[0:1:L-1]; x=a.^nn; [X,W]=DTFT(x,2048); subplot(221),plot(W,abs(X)),grid on xlabel('L=32'); L=64; nn=[0:1:L-1]; x=a.^nn; [X,W]=DTFT(x,2048); subplot(222),plot(W,abs(X)),grid on xlabel('L=64'); L=128; nn=[0:1:L-1]; x=a.^nn; [X,W]=DTFT(x,2048); subplot(223),plot(W,abs(X)),grid on, xlabel('L=128'); L=256; nn=[0:1:L-1]; x=a.^nn; [X,W]=DTFT(x,2048*4); subplot(224),plot(W,abs(X)),grid on, xlabel('L=256'); Xử Lý Tín Hiệu Số & Hệ Thống Nguyễn Duy Bình-Kỹ sư chất lượng cao K53 37 Ta thấy khi số lượng điểm được tăng lên thì đồ thị càng tiến gần đến đồ thi DTFT thu được từ công thức tính trực tiếp.