【正文】 os(pi/6)。y3 = d*sin(pi/6)。plot(0,0,39。go39。,x2,y2,39。go39。,x3,y3,39。go39。)。 %畫圖hold onn=randn(1,3)。 %噪聲產生程序n=n/std(n)。n=nmean(n)。r_mean=1。r_var=2。n=r_mean+sqrt(r_var)*n。i = 0。 %移動臺的位置坐標隨機生成M = 1。while i M x = [2*D*rand(1,2)1*D]。 if(abs(x(1))+abs(x(2))/sqrt(3))=D amp。amp。 abs(x(2))=D*sqrt(3)/2 i = i+1。 hold on plot(x(1)+X,x(2)+Y,39。r+39。)。 endendx = x + [X Y]。K2=x2^2+y2^2。K3=x3^2+y3^2。H=[x2 y2。x3 y3]。r1=sqrt(sum(x.^2))+n(1)。r2=sqrt(sum((x[x2 y2]).^2))+n(2)。r3=sqrt(sum((x[x3 y3]).^2))+n(3)。b=1/2.*[K2r2^2+r1^2。K3r3^2+r1^2]。x_est=inv(H39。*H)*H39。*b。 %定位出的橫坐標Eucdis=x_estx39。 %求定位點與真實的的差XXXXXXX 學士學位論文23 Eucdis=sqrt(sum(Eucdis.^2))。x_est=x_est39。hold on。plot(x_est(1),x_est(2),39。ro39。)。axis square。Eucdisoriginal=[x(1),x(2)]estimate=[x_est(1),x_est(2)]附錄二注釋：此程序主要基于TOA技術運用最小二乘法（LS）來減少NLOS傳播產生的誤差。%============= TOA calculation and MS Position Location ===================% 基站 (BS) 的數目 = 3% 移動臺（MS）的數目 = 1clear allclcc_speed = 3e8。 % 電波的傳播速度為光速n = 3。 % Number of BSBS = [0 0。 70 20。 10 70]。 % 假設基站的坐標，其中必有一基站坐標位于坐標原點MS_ini = [30 30]。 % 移動臺的位置的真實坐標% 為接收信號提供隨機的延時，以仿真現實的環(huán)境sigma = 1e8。 % 延時的時間是x*10 ns。 picture_TOA = 1。 % 打開圖片顯示% 計算兩坐標點的的距離的公式是：% Distance = sqrt((X x)^2 + (Y y^2))。% 計算基站到移動臺的真實距離 (simple Euclidian distance)% 當我們計算出真正的歐式距離后我們可以給與隨機的延時for i = 1:nactual_distance_MS_BS(i) = sqrt((MS_ini(1)BS(i,1))^2 + (MS_ini(2)BS(i,2))^2)。endif picture_TOA == 1%///////////////// Plot IDEAL position of MS ///////////////////////% 說明: circle(c, r, n) figure (1) TITLE(39。Ideal position of MS39。) for i =1:n c = BS(i,:)。 r = actual_distance_MS_BS(i)。 nsides = 128。 col = [0 0 1]。 nsides = max(round(nsides),3)。 % 這里我們要確定基站的數目要= 3 a = [0:2*pi/nsides:2*pi]。XXXXXXX 學士學位論文24 line(r*cos(a)+c(1), r*sin(a)+c(2), 39。color39。, col)。 hold on plot(c(1), c(2))。 hold on end plot(MS_ini(1), MS_ini(2), 39。k39。, 39。MarkerFaceColor39。,39。g39。,39。MarkerSize39。,10)。 grid。 hold offend % 對傳輸路徑產生隨機的延時 % delay = randn * sigma。這是此算法的基本原理 %求出每個移動臺信號到達基站的時間 for i = 1:n TOA(i) = actual_distance_MS_BS(i)/c_speed + abs(randn * sigma)。%在這里我們增加了一個隨機的延遲噪聲，目的是為了仿真出真正的NLOS的室外環(huán)境。 % as follows: L(i) = TOA(i) * c_speed。 end if picture_TOA == 1 figure (2) TITLE(39。Found position of MS39。) for i = 1:n c = BS(i,:)。 r = L(i)。 nsides = 128。 col = [1 0 0]。 nsides = max(round(nsides),3)。 % 這里我們要確定基站的數目要= 3 a = [0:2*pi/nsides:2*pi]。 line(r*cos(a)+c(1), r*sin(a)+c(2), 39。color39。, col)。 hold on plot(c(1), c(2))。 hold on end grid on。 end % 這里我們使用最小二乘法（ LEAST SQUARES） 來估計移動臺的位置 for i = 2:n P(i1,1) = BS(i,1)。 P(i1,2) = BS(i,2)。 end B(1,:) = (BS(2,1)^2 + BS(2,2)^2) L(2)^2 + L(1)^2。 B(2,:) = (BS(3,1)^2 + BS(3,2)^2) L(3)^2 + L(1)^2。XXXXXXX 學士學位論文25 k = (P39。*P)^(1)*P39。*B*。 if picture_TOA == 1 plot(k(1), k(2), 39。o39。, 39。MarkerFaceColor39。,39。g39。,39。MarkerSize39。,5)。 hold off end Error = sqrt((MS_ini(1) k(1))^2 + (MS_ini(2) k(2))^2) [k(1), k(2)] MS_ini附錄三注釋：基于到達時間（TOA）的三種算法仿真程序%初始化程序clc。 % 清屏clear all。 % 清除工作空間close all。 % 關閉所有圖形warning off。% ====仿真節(jié)點布置====so = [0 0 10 10。0 10 0 10]。 %基站坐標（BS） uo = [2,3]39。 % 移動臺坐標（MS）M = size(so,2)。 % 基站的數目N = size(so,1)。 % 定位的維度（二維）figure(1)。 %畫出各節(jié)點的布置圖plot(uo(1),uo(2),39。r^39。)。hold on。grid on。for i = 1 : M plot(so(1,i),so(2,i),39。bo39。)。end。axis([5,15,5,15])。xlabel(39。xlabel39。)。ylabel(39。ylabel39。)。 d=(sqrt(sum((uo*ones(1,M)so).^2,1)))39。 % BS與MS的真實距離sigma2 = 。 %高斯噪聲的方差 r = d + randn(M,1)*sigma2。 %測量距離的方程 iter = 30。 %迭代的次數% ==== NewtonRaphson 算法====u = [4,5]39。 %隨機生成移動臺的位置坐標 for i = 1:iter H = hessian_nls(so,u,r)。 g = grad_nls(so,u,r)。 u = u inv(H)*g。 uo_nr(:,i)=u。end。% ====NLS GaussNewton 算法====XXXXXXX 學士學位論文26u = [3,2]39。 for i = 1:iter G = jacob(so,u)。 f_TOA = sqrt(sum((ones(M,1)*u39。so39。).^2,2))。 u= u + inv(G39。*G)*G39。*(rf_TOA)。 uo_gn(:,i) = u。end。% ====最速下降法====u = [3,2]39。mu = 。for i=1:iter g = grad_nls(so,u,r)。 u = umu*g。 uo_sd(:,i) = u。end。 % ====畫出仿真結果圖====figure(2)。plot(1:iter,uo_nr(1,:),39。k.39。)。hold on。plot(1:iter,uo_gn(1,:),39。ko39。)。hold on。plot(1:iter,uo_sd(1,:),39。k+39。)。hold on。legend(39。NewtonRaphson39。,39。GaussNetwon39。,39。steepest descent39。)。ylabel(39。estimate of x39。)xlabel(39。number of iterations39。)figure(3)。plot(1:iter,uo_nr(2,:),39。k.39。)。hold on。plot(1:iter,uo_gn(2,:),39。ko39。)。hold on。plot(1:iter,uo_sd(2,:),39。k+39。)。hold on。legend(39。NewtonRaphson39。,39。GaussNetwon39。,39。steepest descent39。)。ylabel(39。estimate of y39。)xlabel(39。number of iterations39。)souo 子函數一 function H = hessian_nls(so,u,r)% NLS Hessian 矩陣計算% so: 基站位置矩陣% u ：移動臺位置矩陣% r ：TOA測量的距離向量 M = size(so,2)。 % 基站的數目N = size(so,1)。 % 空間維度A = 0。B = 0。XXXXXXX 學士學位論文27C = 0。for i = 1:M A = A + (u(1)so(1,i))^2/(norm(uso(:,i)))^2 ((r(i)norm(uso(:,i)))*(u(2)so(2,i))^2)/(norm(uso(:,i)))^3。 B = B + (r(i)*(u(1)so(1,i))*(u(2)so(2,i)))/(norm(uso(:,i)))^3。 C = C + (u(2)so(2,i))^2/(norm(uso(:,i)))^2 ((r(i)norm(uso(:,i)))*(u(1)so(1,i))^2)/(norm(uso(:,i)))^3。end。H = 2.*[A B。 B C]。子函數二function G = jacob(so,u)% 雅克比矩陣計算M = size(so,2)。N = size(so,1)。dx = zeros(M,1)。dy = zeros(M,1)。for i = 1:M dx(i) = (u(1)so(1,i))/norm(uso(:,i))。 dy(i) = (u(2)so(N,i))/norm(uso(:,i))。 end。G = [dx,dy]。子函數三function g = grad_nls(so,u,r)% NLS梯度矩陣計算M = size(so,2)。 N = size(so,1)。 A = 0。B = 0。for i = 1:M A = A + ((r(i)norm(uso(:,i)))*(u(1)so(1,i)))/norm(uso(:,i))。 B = B + ((r(i)norm(uso(:,i)))*(u(2)so(2,i)))/norm(uso(:,i))。end。g = 2.*[A,B]39。