dSp first solution PDF

Preview

Summary

good...

Description

Lab P-3: Introduction to Complex Exponentials – Direction Finding 電機三 電機三乙 乙 405250175 詹舒閔 Pre-Lab (a) Enter the complex numbers z1 and z2 in MATLAB. Plot them with zvect(), and print them with zprint() . z1=2*exp(j*pi/3);z1 = ej𝛑/𝟑 z2=-sqrt(2)+5j; z2 = −√2+j5 p=zvect([z1,z2]);把複數向量印出來 set( p , 'linewidth' , 2 )加粗線 條 title('(a)-zvect');加上標題

zprint([z1,z2]);顯示複數數值(在command window) Z=

X 1 -1.414

+

jY Magnitude 1.732 2 5 5.196

Phase 1.047 1.846

Ph/pi Ph(deg) 0.333 60.00 0.588 105.79

(b) Compute the conjugate z_ and the inverse 1=z for both z1 and z2 and plot the results. In MATLAB, see help conj. Display the results numerically with zprint.

z3 = conj(z2);把Z2取共厄複數 z4 = conj(z1);把Z1取共厄複數 q=zvect([z4,z3]); 把複數向量印出來 set( q , 'linewidth' , 2 ) 加粗線 條 title('(b)-共厄複數'); 加上標題

z5 = 1./z1;把Z1取倒數 z6 = 1./z2; 把Z2取倒數 r=zvect ([z5,z6]); 把複數向量印出來 set( r , 'linewidth' , 2 ) 加粗線 條 title('(b)-倒數') 加上標題

zprint([z3,z4,z5,z6]); 顯示複數數 值(在command window)

Z=

X + -1.414 1 0.25 -0.05238

jY

Magnitude Phase Ph/pi Ph(deg) -5 5.196 -1.846 -0.588 -105.79 -1.732 2 -1.047 -0.333 -60.00 -0.433 0.5 -1.047 -0.333 -60.00 -0.1852 0.1925 -1.846 -0.588 -105.79

(c) The function zcat() can be used to plot vectors in a “head-to-tail” format. Execute the statement zcat([1+j,-2+j,1-2j]); to see how zcat() works when its input is a vector of complex numbers.

r=zcat([1+j,-2+j,1-2j]);把向量頭 街尾接起來 set( r , 'linewidth' , 2 ) 加粗 線條 title('(c)-zcat'); 加上標題

(d) Compute z1Cz2 and plot the sum using zvect(). Then use zcat() to plot z1 and z2 as 2 vectors head-to-tail, thus illustrating the vector sum. Use hold on to put all 3 vectors on the same plot. Ifyou want to see the numerical value of the sum, use zprint() to display it. z7=z1+z2;把Z1Z2相加 r=zvect(z7); 把複數向量印出來 set( r , 'linewidth' , 2 ) 加粗線 條 title('(d)-相加'); 加上標題 r=zcat ([z1,z2]) ;把Z1Z1頭尾連在一 起 set( r , 'linewidth' , 2 ) 加粗線 條 title('(d)-zcat'); 加上標題 s=zcat ([z1,z2],'r') ;把2張圖連在 一起 set( s , 'linewidth' , 2 ) 加粗線條 hold on; 先不要把圖印出來 r=zvect(z7); 把複數向量印出來 set( r , 'linewidth' , 2 ) 加粗線條 hold off ;把剛剛Hold on的圖和之後的圖顯示在同一張圖片 title('(d)-hold'); 加上標題 legend([s(1),r(1)],'[z1,z2]','z1+z2')標示向量顏色、代表意思

zprint([z7,z1,z2]); 顯示複數數值(在command window) Z=

X + -0.4142 1 -1.414

jY Magnitude 6.732 6.745 1.732 2 5 5.196

Phase 1.632 1.047 1.846

Ph/pi Ph(deg) 0.520 93.52 0.333 60.00 0.588 105.79

(e) Compute z1z2 and z2=z1 and plot the answers using zvect() to show how the angles of z1 and z2 determine the angles of the product and quotient. Use zprint() to display the results numerically. z8 = z1 .* z2;Z1Z2 r=zvect ([z8]); 把複數向量印出來 set( r , 'linewidth' , 2 ) 加粗 線條 title('(e)-z1*z2');加上標題 z9 = z2 ./ z1;Z2/Z1 r=zvect ([z9]); 把複數向量印出來 set( r , 'linewidth' , 2 ) 加粗 線條 title('(e)-z2/z1'); 加上標題

a=dot(z1,z2)/norm(z1,2)/norm(z2,2); 算夾角 rad2deg(acos(a)); 把弧度變成角度 r=zvect ([a]); 把複數向量印出來 set( r , 'linewidth' , 2 ) 加粗線條 title('(e)-z1和z2的夾角'); 加上標題 c=dot(z8,z9)/norm(z8,2)/norm(z9,2); 算夾角 rad2deg(acos(c));把弧度變成角度 r=zvect ([c]); 把複數向量印出來 set( r , 'linewidth' , 2 ) 加粗線條 title('(e)-z1*z2和z2/z1的夾角'); 加上 標題

zprint([a,c]); 顯示複數數值(在command window) Z= X + jY Magnitude Phase Ph/pi Ph(deg) 0.6973 0.7168 1 0.799 0.254 45.79 -0.5 -0.866 1 -2.094 -0.667 -120.00 (f) Make a 2 _ 2 subplot that displays four plots in one window: similar to the four operations done previously: (i) z1, z2, and the sum z1Cz2 on a single plot; (ii) z2 and z_2 on the same plot; (iii) z1 and 1=z1 on the same plot; and iv) z1z2. Add a unit circle and x-y axis to each plot for reference.zprint([a,c]); subplot(2,2,1);在2x2的四個圖中，畫第一個小圖 alpha=0:pi/20:2*pi; R=1; x=R*cos(alpha); y=R*sin(alpha); plot(x,y,'k'); axis equal hold on;畫半徑為1的圓 s=zcat (z1,'r') ;用紅色畫Z1 set( s , 'linewidth' , 2 ) 加粗線條 hold on;保留圖形，先不印 d=zcat (z2,'g') ; set( d , 'linewidth' , 2 ) 加粗線條 hold on; 保留圖形，先不印 r=zvect(z7);畫出Z7向量 set( r , 'linewidth' , 2 ) 加粗線條 hold off ;將圖形顯示出來 title('(f)-(i)'); 加上標題 legend([s(1),d(1),r(1)],'z1','z2','z1+z2') 標示向量顏色、代表意 思

subplot(2,2,2); 在2x2的四個圖中，畫第二個小圖 alpha=0:pi/20:2*pi; R=1; x=R*cos(alpha); y=R*sin(alpha); plot(x,y,'k'); axis equal畫半徑為1的圓 hold on; 保留圖形，先不印 d=zcat (z2,'g') ;用綠色畫Z2 set( d , 'linewidth' , 2 ) 加粗線條 hold on; 保留圖形，先不印 z3 = conj(z2);Z3=Z2的共厄負數 r=zvect(z3); 畫出Z3向量 set( r , 'linewidth' , 2 ) 加粗線條 hold off ; 將圖形顯示出來 title('(f)-(ii)'); 加上標題 legend([d(1),r(1)],'z2','z2*') 標示向量顏色、代表意思 subplot(2,2,3); 在2x2的四個圖中，畫第三個小圖 alpha=0:pi/20:2*pi; R=1; x=R*cos(alpha); y=R*sin(alpha); plot(x,y,'k'); axis equal hold on; 畫半徑為1的圓 s=zcat (z1,'r') ;用紅色畫出Z1 set( s , 'linewidth' , 2 ) 加粗線條 hold on; 保留圖形，先不印 d=zcat (1./z1,'g') ; 用綠色畫出1/Z1 set( d , 'linewidth' , 2 ) 加粗線條 hold off ; 將圖形顯示出來 title('(f)-(iii)'); 加上標題 legend([s(1),d(1)],'z1','1/z1') 標示向量顏色、代表意思

subplot(2,2,4); 在2x2的四個圖中，畫第四個小圖 alpha=0:pi/20:2*pi; R=1; x=R*cos(alpha); y=R*sin(alpha); plot(x,y,'k'); axis equal

hold on; 畫半徑為1的圓 s=zcat (z1.*z2,'r') ; 用紅色畫出Z1Z2 set( s , 'linewidth' , 2 ) 加粗線條 hold off ; 將圖形顯示出來 title('(f)-z1*z2'); 加上標題

Lab 05 (a) Generate the signal −2ej50πt -ej50π (t−0.02) +(2 − j3)ej50πt make a plot versus t.Use the syn sin function and take a range for t that will cover 3 periods. Include the MATLAB code and theplot with your report.

syn sin function [xx,tt] = syn_sin(fk, Xk, fs, dur, tstart) %SYN_SIN Function to synthesize a sum of cosine waves % usage: % [xx,tt] = syn_sin(fk, Xk, fs, dur, tstart) % fk = vector of frequencies % (these could be negative or positive) % Xk = vector of complex amplitudes: Amp*e?(j*phase) % fs = the number of samples per second for the time axis % dur = total time duration of the signal % tstart = starting time (default is zero, if you make this input optional)

% xx = vector of sinusoidal values % tt = vector of times, for the time axis % % Note: fk and Xk must be the same length. % Xk(1) corresponds to frequency fk(1), % Xk(2) corresponds to frequency fk(2), etc. tt= tstart : 1/fs : dur; xx= zeros(1,length(tt)); M=length(fk); for k=1:M xx = xx + real(Xk(k)*exp(j*2*pi*fk(k)*tt)); end end [x,t] = syn_sin([25,25,25],[-2,-1*exp(-pi),(23*j)],10000,0.12,0); −2ej50πt -ej50π(t−0.02) +(2 − j3)ej50πt (3個週期，0.12s) q= plot(t,x); xlabel('time'); 標示Ｘ軸 ylabel('amplitude'); 標示Ｙ軸 title('x(t)'); 加上標題 set( q , 'linewidth' , 2 )

(b) From the plot of x.t/ versus t, measure the frequency, phase and amplitude of the sinusoidal signal by hand.Show annotations on the plots to indicate how these measurements were made and what the values are. Compareto the calculation in part (c).t =0.04s f = 1/0.04 = 25Hz Amplitude = 3; Φ = -ωt t = 0.01s Φ = -2*180*25*0.01 = -90o

(c) Use the phasor addition theorem and MATLAB to determine the magnitude and phase of x.t/ z1 = -2*exp(j*50*pi*t); −2ej50πt z2 = -1*exp(j*50*pi*(t-0.02));ej50π (t−0.02) z3 = (2-3*j)*exp(j*50*pi*t); (2 − j3)ej50πt z = z1+z2+z3; zprint(z); Z= X + 1

jY

Magnitude -3 3.162

Phase -1.249

Ph/pi -0.398

Ph(deg) -71.57

Lab-exercise (c) The received signals at the receivers, called x1.t/ and x2.t/, are delayed copies of the transmitter signalwhere is the transmitted (sinusoidal) signal.3 Assume that the source signal s.t/ is a zero-phase sinusoid at f0 D 400 Hz; and also assume that the amplitudeof the transmitted signal is 1000. Make a plot of x1.t/ and x2.t/ when xv D 100 meters. Use subplot to put both signals on the figure. Plot only 3 periods and then measure the relative time-shift between the two received signals by comparing peak locations. clear all; f=400; yr =100; d =0.4; c=1000/3; xv=100; (a)給的參數 path1= sqrt(xv*xv+yr*yr) ;計算path1的距離 path2= sqrt((d-xv)*(d-xv)+yr*yr) ; 計算path2的距離 t1=path1./c;距離/速率=時間 t2=path2./c; 距離/速率=時間 subplot(3,1,1);把圖片分成3份，第一部分 t = 0 : 0.0001 : 0.0075 ; 3個週期，一次增加0.0001 y = 1000 * exp(j*800*pi*t) ;帶入1000ej800πt yy=plot (t,y,'r') ;畫圖 hold on ;

subplot(3,1,2); 把圖片分成3份，第二部分

y1 = 1000 * exp(j*800*pi*(t-t1)); 帶入1000ej800π(t−t1) y11=plot (t,y1,'g'); 畫圖 hold on ; subplot(3,1,3); 把圖片分成3份，第三部分 y2 = 1000 * exp(j*800*pi*(t-t2)); 帶入1000ej800π(t−t2) y22=plot (t,y2,'b'); 畫圖 hold off ;將圖片顯示 xlabel('time'); x軸 ylabel('amplitude'); y軸 title('c'); 加上標題 legend([yy(1),y11(1),y22(1)],'t','t1','t2'); 標示向量顏色、代表 意思

在第一次波峰時，t2的波鋒在最左，左邊表示超前，所以t2收到訊號的時間最早。 (d) How do we convert relative time-shift into the direction ? The answer is the following equation which can besolved for For the relative time-shift obtained in the previous part, calculate . In addition, use geometry and the valuesof xv and yr to figure out what the “true value” of should be. Verify that your calculated value of is very close to the true value. rad = ((asin(t1-t2)*c)/d);算出θ的弧度(套公式算) sin=xv/path1; rad1=asin(sin); 算出θ的弧度(自己算) deg= (rad*360)/(2*pi);轉成角度 deg1=(rad1*360)/(2*pi); rad = 0.7064 rad1 =

(公式) (數學運算)

0.7854 deg = 40.4736 deg1 = 45.0000 (e) The objective in the rest of this lab is to write a MATLAB function that will process the received signals to finddirection. To do this, the received signals will be given as complex amplitudes, and a MATLAB function calledDF genis supplied for generating the receiver signals.4DF genis supplied for generating the receiver signals. function [X1,X2,theta] = DF_gen(xx) %DF_GEN generate complex amplitudes at the two receivers % for the Direction Finding Lab % usage [X1,X2,theta] = DF_gen; % % X1 = complex amplitude at Receiver #1 % X2 = complex amplitude at Receiver #2 % theta = the TRUE value of the "direction" in DEGREES % % alternate usage: % [X1,X2,theta] = DF_gen(xx); % % xx = vector of "x positions" of the vehicle % then X1, X2 and theta are vectors f=400; yr =100; d =0.4; c=1000/3; omg=2*pi*f; xx=2*yr*randn; p1 = j*yr - xx; p2 = j*yr + d - xx; R1 = abs(p1); R2 = abs(p2); X0 = 1000*exp(j*2*pi*rand); X1 = X0*exp(j*omg*R1/c)./R1; X2 = X0*exp(j*omg*R2/c)./R2; theta = (180/pi)*angle(j*yr - xx) - 90; end

(f) Debug the MATLAB program from the previous part by using the MATLAB function DF gen to test individual cases. Then run your function for the vehicle moving from -400 meters to +500 meters, in steps of one meter. Compare the computed value of to the true value of (given by DF gen) by plotting both on the same graph. xx=-400:1:500;車子位移 [x1,x2,theta]=DF_gen(xx);套用公式 p = angle(x1./x2);x1/x2的角度 theta2 = 180/pi*angle(100+j*xx);實際角度 s=plot(xx,theta , 'r');模擬角度畫圖 set( s , 'linewidth' , 2 ) hold on ; q=plot(xx, theta2, 'g'); 模擬角度畫圖 set( q , 'linewidth' , 2 ) title('value and true value'); 加上標題 legend([s(1),q(1)],'value','true value'); xlabel('time'); ylabel('amplitude'); hold off;

Lab S-2: Direction Finding: Time-Difference or Phase Difference 2.1 Direction from Phase Difference Use the LabDFgen2015(UserID,1) function with its second argument set to 1 to generate two phase measurements, along with an intersensor istance d, and the propagation velocity c. These values are displayed in the MATLAB command window. Use the theory of direction finding to determine the direction to the source . Write your method and results (and explanation) on the verification page [r1,r2,fs] = LabDFgen2015('405250175',1); f=780; Sinusoid Frequency of the Source omg=2*pi*f;ω phase1=-2.6+2*pi;把相位變正的 phase2=-1.4+2*pi; 把相位變正的 t1=phase1/omg;相位/角度=時間 t2=phase2/omg; 相位/角度=時間 d1=t1*350; 時間*速率=距離 d2=t2*350; 時間*速率=距離 sin=d2/(d1*(1-(2*0.449)/d1));套公式算sin(θ) theta=asin(sin);算出θ theta = -0.5814 =============Version-2 Tues 22-Jan 8am============================ For GT ID = 405250175 the parameters for Direction from Phase Difference are: Propagation Velocity = 350 meters/sec Sinusoid Frequency of the Source = 780 Hz Intersensor Spacing, d = 44.9 cm (CENTIMETERS) Phase at Receiver #1 = -2.6 rads Phase at Receiver #2 = -1.4 rads 2.2 Ambiguous Direction Start from the phases generated previously, i.e., LabDFgen2015(UserID,1). Since the phase of a sinusoid might not be unique, the same is true of the phase difference. By changing ฀', show that there is at least one more direction in the first or fourth quadrant that satisfies the equations used in the previous part. lambda=350/780;算波長 dp1=phase1-phase2;∆φ

dp2=phase1-phase2+3*lambda;改變∆φ deg1=dp1*180/pi;換成角度 deg2=dp2*180/pi; phase3=phase2+3*lambda; t1=phase1/omg;算時間 t3=phase3/omg; d1=t1*350; d3=t3*350;算距離 sin2=d3/(d1*(1-(2*0.449)/d1));代公式算角度 theta2=asin(sin2); theta2 = 0.0035 (第四象限) deg1 = -68.7549 (第一象限) deg2 = 8.3740 2.4 Direction from a Speech Signal SourceIn this part, use the [r1,r2,fs]=LabDFgen2015(UserID,2) function with its second argument set to 2 to generate two time-shifted speech signals, along with parameter values for the intersensor distance d, the propagation velocity c, and the sampling rate fs. The values of distance d and propagation velocity c are displayed in the MATLAB command window. [r1,r2,fs]=LabDFgen2015('405250175',2); =============Version-2 Tues 22-Jan 8am============================ For GT ID = 405250175 the parameters for Direction from Time-Shifted Speech are: Propagation Velocity = 350 meters/sec Sampling Frequency of the Signal = 11000 samples/sec Intersensor Spacing, d = 150 cm (CENTIMETERS)

Use the [r1,r2,fs]=LabDFgen2015(UserID,3)function with its second argument set to 3 to generate twotime-shifted noise signals, along with an intersensor distance d, the propagation velocity c, and samplingrate fs. The values of distance d and propagation velocity c are displayed in the MATLAB commandwindow. Use your template matching function (based on the AMDF) to obtain an estimate for the TDOAfrom which you can determine the direction ./. Write your method and results on the Lab-HW page. In order to use your template matching function it is necessary to take a small part of one signal andtry to find it within the other signal. One suggested approach would be to take the “center-half” of the first signal and match it to the second signal. Center-half means drop the first 25% and last 25% of the points in the vector. Once you get the best match, the shift index must be compensated for the center-half operation. Finally, the shift index must be converted to a time (in secs) using the sampling rate. [r1,r2,fs] = LabDFgen2015('405250175',3); =============Version-2 Tues 22-Jan 8am============================ For GT ID = 405250175 the parameters for Direction from Time-Shifted Noise are: Propagation Velocity = 350 meters/sec Sampling Frequency of the Signal = 9500 samples/sec Intersensor Spacing, d = 260 cm (CENTIMETERS)...