第八講測量程序設計在測繪中應用初步性

1、測量程序設計應用測量程序設計應用1.1.測量數據處理：進行平差計算（近似、測量數據處理：進行平差計算（近似、精確）；廣義測量平差計算（濾波、配精確）；廣義測量平差計算（濾波、配置、抗差估計等）；其它測量中有關計置、抗差估計等）；其它測量中有關計算算2.2.從事科學研究所必要的計算（讀研）從事科學研究所必要的計算（讀研）等等測量程序設計流程測量程序設計流程一、一、Matlab用于經典平差計算用于經典平差計算 1.條件平差條件平差函數模型函數模型11010rrnrnALA0WAV隨機模型隨機模型12020PQD minPVVT平差準則平差準則 條件平差就是在滿足條件平差就是在滿足r個條件方程式的條

2、件下，個條件方程式的條件下，求使函求使函 數數 最小的最小的 值。值。PVVTV1）、列平差值條件方程：）、列平差值條件方程：11010rrnrnALA000022110221102211rLrLrLrbLbLbLbaLaLaLannnnnn), 2 , 1(,nirbaiii條件方程系數條件方程系數000,rba常數項常數項2）、條件方程：）、條件方程：11010rrnrnALA 代入平差值條件方程中，得到代入平差值條件方程中，得到將將VLL11110rrnrWVA為條件方程閉合差為條件方程閉合差1rW0AALW3）、改正數方程：）、改正數方程：按求函數條件極值的方法引入聯系系數向量按求函數

3、條件極值的方法引入聯系系數向量TrbarkkkK),(1組成拉氏函數，并計算其最小值：組成拉氏函數，并計算其最小值：)(2WAVKPVVTT 將將對對V求一階導數并令其為零求一階導數并令其為零 022AKPVVTT 則：則：KAPVTKQAKAPVTT14）、法方程：）、法方程： 由條件方程由條件方程 AV+W=0和改正數方程和改正數方程 ，則，則得：得： 0WKAQAT記作：記作：1110rrrrraaWKNNaa為滿秩方陣，為滿秩方陣， rARAQARNRTaa)()()(由于由于)()(011AALAQAWNKTaaKQAVTVLLKQAVT5）、檢核和精度計算）、檢核和精度計算檢查檢查

4、AV+W=0 ?或或?00ALAAQNQAQQrPVVaaTLT120,例例1 1 教材誤差理論與測量平差基礎教材誤差理論與測量平差基礎P74P74 設對下圖中的三個內角作同精度觀測，得觀測值：設對下圖中的三個內角作同精度觀測，得觀測值：L L1 1=42=42o o12122020，L L2 2=78=78o o09090909，L L3 3=59=59o o38384040，試按條件平差求三個內角得平差值。試按條件平差求三個內角得平差值。clcclcdisp(disp(條件平差示例條件平差示例1)1)disp(disp(三角形內角觀測值三角形內角觀測值) ) L1 = 42 12 20L1

5、 = 42 12 20; ; % %注意格式注意格式L2 = 78 9 9L2 = 78 9 9; ;L3 = 59 38 40L3 = 59 38 40; ;L = L1; L2; L3L = L1; L2; L3; ;disp(disp(度分秒轉換為弧度度分秒轉換為弧度)LL = dms2rad(mat2dms(L)LL = dms2rad(mat2dms(L)A = 1 1 1A = 1 1 1; ;w = sum(LL(:) w = sum(LL(:) pi pi ; ; % %閉合差閉合差w = dms2mat(rad2dms(ww = dms2mat(rad2dms(w); ; %

6、 %弧度化為度分秒弧度化為度分秒P = eye(3);P = eye(3);Naa = ANaa = A* *inv(P)inv(P)* *AA; ;K = -inv(Naa)K = -inv(Naa)* *w w; ; V =inv(PV =inv(P) )* * A A* *K K; ; LP = L + VLP = L + V LP = LP = dms2rad(mat2dms(LP)dms2rad(mat2dms(LP) sumLP = sum(LP) sumLP = sum(LP) if(sum(LP) = pi) if(sum(LP) = pi) disp( disp(檢核正確檢核

7、正確) else else disp( disp(檢核錯誤檢核錯誤) end end例例2 2 教材誤差理論與測量平差基礎教材誤差理論與測量平差基礎P75P75在下圖中，、為已知水準點，其高程為在下圖中，、為已知水準點，其高程為HA=12.013m, HA=12.013m, HB = 10.013mHB = 10.013m。為了確定點及點的高程，共觀測了。為了確定點及點的高程，共觀測了四個高差，高差觀測值及相應的水準路線的距離為：四個高差，高差觀測值及相應的水準路線的距離為：h1 = -1.004m, S1 = 2km;h1 = -1.004m, S1 = 2km;h2 = 1.516m, S

8、2 = 1km;h2 = 1.516m, S2 = 1km;h3 = 2.512m, S3 = 2km;h3 = 2.512m, S3 = 2km;h4 = 1.520m, S4 = 1.5kmh4 = 1.520m, S4 = 1.5km試求和點高程的平差值。試求和點高程的平差值。clcclcclearclearh1 = -1.004;h1 = -1.004;h2 = 1.516;h2 = 1.516;h3 = 2.512;h3 = 2.512;h4 = 1.520;h4 = 1.520;HA = 12.013HA = 12.013HB = 10.013HB = 10.013L = h1 h

9、2 h3 h4L = h1 h2 h3 h4s1 = 2; s1 = 2; s2 = 1;s2 = 1;s3 = 2;s3 = 2;s4 = 1.5;s4 = 1.5;s=s1 s2 s3 s4s=s1 s2 s3 s4A = 1 1 -1 0; 0 1 0 -1A = 1 1 -1 0; 0 1 0 -1w1 = h1 + h2 - h3 + HA - HB;w1 = h1 + h2 - h3 + HA - HB;w2 = h2 - h4;w2 = h2 - h4;w = w1; w2w = w1; w2; ;P = diag(1./s)P = diag(1./s); ;Naa = ANaa

10、 = A* *inv(P)inv(P)* *AA; ;K = -inv(Naa)K = -inv(Naa)* *w w; ;V = inv(P)V = inv(P)* *AA* *K K; ;LP= L + V;LP= L + V;if if LP(1,1)+LP(2,1)-LP(3,1)+HA-HB=0 & LP(1,1)+LP(2,1)-LP(3,1)+HA-HB=0 & LP(2,1)-LP(4,1)=0LP(2,1)-LP(4,1)=0 disp( disp(檢核正確檢核正確)elseelse disp( disp(檢核錯誤檢核錯誤)endend disp( disp(

11、平差后的高程值平差后的高程值) HC = HA + LP(1,1) HC = HA + LP(1,1) HD = HA + LP(1,1) + LP(4,1) HD = HA + LP(1,1) + LP(4,1)例例3 如圖如圖3.1所示水準網，所示水準網，A、B兩點為高程已知，各兩點為高程已知，各觀測高差及路線長度如表觀測高差及路線長度如表3.1所列。用所列。用條件平差法條件平差法計計算求知點的高程平差值及算求知點的高程平差值及p2和和p3之間平差后高差值之間平差后高差值的中誤差。的中誤差。高差觀測值高差觀測值(m)對應線路長對應線路長度度(km)已知點高程已知點高程（m）h1 = 1.3

12、59h2 = 2.009h3 = 0.363h4 =-0.640h5 = 0.657h6 = 1.000h7 = 1.6501122112H1= 35.000H2= 36.000表表3.1解解 n=7,t=3,r=n-t=4.平差值條件方程為：平差值條件方程為：000031643765521BAHhhHhhhhhhhhh改正數條件方程為：改正數條件方程為：0403070731643765521vvvvvvvvvvv令令C=1，觀測值的權倒數為：，觀測值的權倒數為：21122111PKAPVWNKAAPNTT111,-)(mmVT2.0 1.2- 3.8- 0.3- 4.4- 2.8 0.43-

13、得：得：7h76543216570110000hhhhhhhhhh0.7416AQfNQAfQffQfTTTffT1,01100005.1mm7ffThQrPVV下面求平差后下面求平差后的中誤差：的中誤差：中誤差為中誤差為利用利用Matlab編寫上題的代碼為：編寫上題的代碼為：clear A=-1 1 0 0 -1 0 0;0 0 0 0 1 1 -1;0 0 1 -1 0 -1 0;1 0 -1 0 0 0 0;q=1 1 2 2 1 1 2;Q=diag(q); %生成對角矩陣生成對角矩陣QP=inv(Q); %inv()為矩陣求逆運算為矩陣求逆運算W=-7;7;3;-4;N=A*Q*A;

14、 K=-inv(N)*W;V=Q*A*KL=1.359; 2.009; 0.363; -0.640; 0.657; 1.000; 1.650;LP=L+(V/1000)f=0;0;0;0;1;1;0;sigma=sqrt(V*P*V)/4) Qff=f*Q*f-f*Q*A*inv(N)*A*Q*f;Dff=sigma*sqrt(Qff)vfunction Lv,V,d0=szwtjpc(L,A,W,P) vnl=length(L); vma, na=size(A);vnw=length(W);vif nl=nav disp( 觀測值個數與改正數個數不相等觀測值個數與改正數個數不相等! )v r

15、eturnvendvif ma=nwv warning(所列條件方程有誤所列條件方程有誤) v returnvendv mp,np=size(P);v if mp=npv warning(權陣應為方陣權陣應為方陣) v returnvendvQ=inv(P);vV=-Q*A*inv(A*Q*A)*WvLv=L+V/1000vd0=sqrt(V*P*V/nw)v%運行時必須把函數放在運行時必須把函數放在G盤中，更改路徑盤中，更改路徑2、間接平差、間接平差 在一個控制網中，設有在一個控制網中，設有t t個獨立參數，將每一個個獨立參數，將每一個觀測值都表達成所選參數的函數，以此為基礎進行平觀測值都表

16、達成所選參數的函數，以此為基礎進行平差，求得參數的估計值。差，求得參數的估計值。選擇參數應做到選擇參數應做到足數足數（參數的個數等于必要觀測（參數的個數等于必要觀測數）和數）和獨立獨立（參數間不存在函數關系）。利用參（參數間不存在函數關系）。利用參數將觀測值表示為數將觀測值表示為dXAL其中其中L L為觀測值，為觀測值，為誤差，或者表示為為誤差，或者表示為lXA其中其中l(wèi) lL Ld.d.設誤差設誤差和參數和參數X X的估計值分別為的估計值分別為V V 和和X則有則有l(wèi)XAVxXX0若則誤差方程表示為則誤差方程表示為lxAV其中常數項為其中常數項為)(0dAXLl由最小二乘準則，所求參數的改正

17、數應該滿足由最小二乘準則，所求參數的改正數應該滿足minPVVT目標函數對目標函數對x x求一階導數，并令其為零求一階導數，并令其為零022)(PAVxddVPVxdPVVdTTT轉置后得到轉置后得到0PVAT把誤差方程代入上式后得到把誤差方程代入上式后得到0PlAxPAATT設設PAANTAPlAWT則法方程為則法方程為0WxNA由此求得參數改正數的唯一解為由此求得參數改正數的唯一解為WNxA1將其代入誤差方程，可求得改正數將其代入誤差方程，可求得改正數V V，最后得到觀測，最后得到觀測值得平差值及參數的平差值分別為值得平差值及參數的平差值分別為VLLxXX0精度評定精度評定單位權中誤差的估

18、計公式為單位權中誤差的估計公式為tnPVVT0由協因數傳播定律得參數的權逆陣為由協因數傳播定律得參數的權逆陣為1 AxxXXNQQ參數的中誤差為參數的中誤差為jjjXXXQ0若設參數的函數為若設參數的函數為),(21tXXXfF以以xXX0代入，并按泰勒級數展開，取一次項代入，并按泰勒級數展開，取一次項ttxfxffF110記為記為xfxfxfTttF11由此得到參數函數的權逆陣為由此得到參數函數的權逆陣為fNffQfQATXXTFF1中誤差為中誤差為FFFQ0改正數為改正數為lIPAANlxAVTA)(1TAVVAANQQ1平差值及其協因數陣分別為平差值及其協因數陣分別為VLLTAVVVLL

19、VLLLLAANQQQQQ1例例1 1 教材教材P124P124 A A和和B B是已知高程的水是已知高程的水準點，準點，C C、D D、E E是待定點，是待定點，A A、B B高程、觀測高差和高程、觀測高差和相應的水準路線長度見相應的水準路線長度見下表，試按間接平差求下表，試按間接平差求各待定點的高程平差值。各待定點的高程平差值。線路編號觀測高差（m）線路長度（km）已知高程（m）1234567+1.359+2.009+0.363+1.012+0.657+0.238-0.5951.11.72.32.72.41.42.6HA=5.016HB=6.016程序代碼如下：程序代碼如下：disp(水準

20、網間接平差示例水準網間接平差示例)disp(已知高程已知高程)Ha = 5.015 % 已知點高程已知點高程Hb = 6.016 % 已知點高程已知點高程disp(觀測高差，單位觀測高差，單位m)L = 1.359;2.009; 0.363; 1.012; 0.657; -0.357 disp(系數矩陣系數矩陣B)B = 1 0; 0 1; 1 0; 0 1; -1 1; -1 0l = 0; 0; 4; 3; 7; 2;S = 1.1, 1.7, 2.3, 2.7, 2.4, 4.0;P = C./S ; P = diag(P); % 定義權陣定義權陣disp(參數的解參數的解)x = in

21、v(B*P*B)*B*P*ldisp(誤差誤差V(mm), 各待定點的高程平差值各待定點的高程平差值L1（m）)V = B*x - l ; % 誤差方程誤差方程(mm)L1 = L + V/1000 % 觀測值的平差值，觀測值的平差值，disp(精度評定精度評定)n = 6; % 觀測值的個數觀測值的個數t = 2; % 必要觀測數必要觀測數delta = sqrt(V*P*V/(n t)disp(C是單位權觀測高差的線路公里數，是單位權觀測高差的線路公里數，S是線路長度是線路長度)C = ones(1,6);disp(水準網間接平差示例水準網間接平差示例)disp(已知高程已知高程)Ha =

22、 5.015 % 已知點高程已知點高程Hb = 6.016 % 已知點高程已知點高程disp(觀測高差觀測高差)L = 1.359;2.009; 0.363; 1.012; 0.657; -0.357 function jiesuan2( )fid=fopen(E:szjj.txt,r) %注意路徑注意路徑temp=fscanf(fid,%f%f%f%f,4,6)fclose(fid); temp=temp; B=temp(:,1,2); l=temp(:,3); S=temp(:,4); S=S; disp(C是單位線路公里數，是單位線路公里數，S是線路長度是線路長度) C = ones(1

23、,6)P = C./S ; P = diag(P) % 定義權陣定義權陣disp(參數的解參數的解)x=inv(B*P*B)*B*P*l disp(誤差誤差V(mm), 各待定點的高程平差值各待定點的高程平差值L1（m）) V=B*x-l % 誤差方程誤差方程(mm)L1 = L + V/1000 % 觀測值的平差值，觀測值的平差值，disp(精度評定精度評定)n = 6; % 觀測值的個數觀測值的個數t = 2; % 必要觀測數必要觀測數delta = sqrt(V*P*V/(n t)例例2. 圖圖4.3中中,A,B,C是已知點是已知點, 為待定點為待定點,網中網中觀測了觀測了12個角度和個

24、角度和6條邊長。已知測角中誤差為條邊長。已知測角中誤差為 ，邊長測量中誤差為邊長測量中誤差為 cm，起算數據及觀測值分，起算數據及觀測值分別列表于表別列表于表4.1和表和表4.2。21PP、5 . 1 0 . 2ABC 點號點號 坐坐 標標 (m) 坐標方位角坐標方位角邊邊 長長(m)xy4899.8468781.9454548.795130.8121099.4437572.622 14 00 35.77123 10 57.974001.1177734.443表表4.1編號編號觀測角觀測角編號編號觀測角觀測角編號編號觀測邊觀測邊（m）12345684 07 38.237 46 34.958 0

25、5 44.133 03 03.2126 01 55.720 55 02.378910111274 18 16.877 27 59.128 13 43.255 21 09.972 22 25.852 16 20.51314151617182463.943414.715216.236042.945085.085014.99表表4.2.2211TYXYXX )./(5625. 025 . 1222222cmpSS秒.11814418118lxBV.)(1PlBPBBxTT解：本題解：本題n =18,即有即有18個誤差方程，其中有個誤差方程，其中有12個角度誤差個角度誤差方程和方程和4個邊長誤差方程。

26、必要觀測數個邊長誤差方程。必要觀測數t=4.采用間接平差采用間接平差，現取待定點坐標的平差值為參數，即現取待定點坐標的平差值為參數，即以測角中誤差為該邊角網的單位權觀測值中誤差，則測邊的權為：以測角中誤差為該邊角網的單位權觀測值中誤差，則測邊的權為：解題的關鍵是列誤差方程解題的關鍵是列誤差方程在此基礎上解出在此基礎上解出完成精度評定。完成精度評定。解題具體過程如下。解題具體過程如下。21APAP和.,73312LLLABAPABAP.sin,cos11110101APAPAAPAPASYYSXX（一）。首先計算待定點的（一）。首先計算待定點的近似近似坐標坐標。利用觀測值和已知坐標求出。利用觀測

27、值和已知坐標求出的近似坐標方位角，具體計算公式為的近似坐標方位角，具體計算公式為進一步可由坐標增量公式進一步可由坐標增量公式求得待定的近似坐標為：求得待定的近似坐標為：)556.2475,928.5656(),099.2944,875.663(0201PP點近似坐標類似計算2PkjkjkkjkjkjjkjkjjkjkjkySXxSYySXxSY)()()()(200200200200（二）（二）.計算坐標方位角改正數方程的系數。計算坐標方位角改正數方程的系數。按按計算。計算計算。計算時時以米為單位，以米為單位，而而因數值較小，一般采用厘米為單位，此時因數值較小，一般采用厘米為單位，此時有關計算

28、見表有關計算見表-1-1。1 x1 y2 x2 y方向方向的系數（秒的系數（秒/cm）AP1AP2BP2CP1CP2P1 P2-0.2244000.26140-0.0384-0.337900-0.21940-0.40950-0.7966-0.243400.38640.038400.2572-0.552800.0840.4095表表-1000jkjkjkYXS、kkjjyxyx,65.2062 （三）（三）.計算測邊的誤差方程的系數。即計算測邊的誤差方程的系數。即ikjkjkkjkjkjjkjkjjkjkilySYxSXySYxSXv00000000的系數。的系數。1 x1 y2 x2 y方向方

29、向邊長誤差方程系數邊長誤差方程系數AP1AP2BP2CP1CP2P1 P2 -0.833 0-0.6430-0.996-0.9960.55300-0.76600.09300.307-0.91500.2120.99600.9520.4030-0.977-0.093表表-2)(00jhjkiiLl0jkiiSLll（四）（四）.確定誤差方程。在計算確定誤差方程。在計算或或時，利用近似坐標求出每個邊的近似坐標方位角和近似邊長。時，利用近似坐標求出每個邊的近似坐標方位角和近似邊長。從而可確定從而可確定的系數。進一步可得出誤差方程的系數項，見表的系數。進一步可得出誤差方程的系數項，見表-3。1 x1 y

30、2 x2 ylpiiS表表-3編號編號1 12 23 34 45 56 67 78 89 91010111112120 00 00 00 00 00 0-0.2244-0.22440.03840.03840.1860.1860.29980.2998-0.0384-0.0384-0.2614-0.26140 00 00 00 00 00 0-0.3379-0.33790.40950.40950.07160.07160.19010.1901-0.4095-0.40950.21940.21940.55320.55320.24340.2434-0.7966-0.7966-0.2434-0.24340.

31、62980.6298-0.3864-0.38640.79660.7966-0.835-0.8350.03840.0384-0.0384-0.0384-0.348-0.348 0.3864 0.3864-0.8055-0.80550.55280.55280.25720.2572-0.5528-0.55280.63680.6368-0.084-0.084-0.2572-0.2572-0.1523-0.15230.40950.4095-0.4095-0.40950.32550.32550.0840.084-2.6-2.6-0.11-0.11-0.02-0.020.480.480.300.300.42

32、0.420.020.020.050.05-0.97-0.970.240.241.051.05-5.09-5.091 11 11 11 11 11 11 11 11 11 11 11 11313141415151616171718180 00 00 0-0.643-0.643-0.833-0.833-0.996-0.9960 00 00 0-0.766-0.7660.5530.5530.0930.0930.3070.307-0.915-0.9150.2120.2120 00 00.9960.9960.9520.9520.4030.403-0.977-0.9770 00 0-0.093-0.093

33、0.10.112.112.19.79.710.710.70 00.10.10.56250.56250.56250.56250.56250.56250.56250.56250.56250.56250.56250.5625可得法方程為可得法方程為（五）（五）.平差值計算與精度評定。由法方程可解得參數改正數為：平差值計算與精度評定。由法方程可解得參數改正數為：2965. 00249. 00481. 00074. 02709. 00775. 01583. 09888. 00137. 07970. 0,8997. 03647. 38065. 65297. 312211對稱bbNyxyx02148. 19

34、415. 81620. 67469. 24194. 32866. 01437. 00229. 02866. 03174. 43363. 08608. 01437. 03363. 00459. 10861. 00229. 08608. 00861. 04270. 12211yxyxx lxBV21,PP由上式可解出由上式可解出，再根據，再根據得各個改正數為得各個改正數為V= 1.9 -1.6 1.4 -3.2 0.7 -4.8 1.3 4.3 -2.5 0 6.1 -2.8 -10.7 -8.2 -3.4 -1.8 2.3T進而可求出觀測值的平差值。進而可求出觀測值的平差值。單位權中誤差單位權中

35、誤差的樅橫坐標的協因數和的樅橫坐標的協因數和點位中誤差點位中誤差分別為分別為2967. 30tnPVVT2965. 0,2709. 0,9888. 0,7970. 022221111YYXXYYXXQQQQ2.4833(cm)4.1622(cm)222202111101YYXXPYYXXPQQQQ利用利用Matlab編寫上題的代碼為：編寫上題的代碼為：clearB=0 0 0.5532 -0.8055;0 0 0.2434 0.5528;0 0 -0.7966 0.2572;0 0 -0.2434. -0.5528;0 0 0.6298 0.6368;0 0 -0.3864 -0.084;-0

36、.2244 -0.3379 0.7966. -0.2572;0.0384 0.4095 -0.835 -0.1523;0.186 0.0716 0.0384 0.4095;.0.2998 0.1901 -0.0384 -0.4095;-0.0384 -0.4095 -0.348 0.3255;.-0.2614 0.2194 0.3864 0.084;0 0 0.307 0.952;0 0 -0.915 0.403;0 0 0.212. -0.977;-0.643 -0.766 0 0;-0.833 0.553 0 0;-0.996 0.093 0.996 -0.093;l=-2.6;-0.11

37、;-0.02;0.48;0.30;0.42;0.02;0.05;-0.97;0.24;1.05;-5.09;0.1;12.1.;9.7;10.7;0;0.1;X0=663.875,2944.099,5656.928,2475.556;p=1 1 1 1 1 1 1 1 1 1 1 1 0.5625 0.5625 0.5625 0.5625 0.5625 0.5625;P=diag(p);N=B*P*B;W=B*P*l;x=inv(N)*W;XP=X0+x/100invN=inv(N) v=B*x-l;sigma=sqrt(l*P*l-W*x)/14)sigmaP1=sigma*sqrt(inv

38、N(1,1)+invN(2,2)sigmaP2=sigma*sqrt(invN(3,3)+invN(4,4)v二、二、MATLABMATLAB在數字測量技術中簡單應用在數字測量技術中簡單應用v1.1.坐標正、反算坐標正、反算v1 1）坐標正算）坐標正算v 坐標正算即已知坐標正算即已知A A點坐標點坐標 ,A ,A點至點至B B點水平點水平距離距離 以及坐標方位角以及坐標方位角 ，求，求B B點坐標點坐標 的過程。計算公式如下：的過程。計算公式如下：v ABSABABAABABABABAABABAByySyyxxSxxsincosvax=input(請輸入請輸入A的縱坐標的縱坐標:)vay=in

39、put(請輸入請輸入A點的橫坐標點的橫坐標:)vs=input(請輸入請輸入AB長度長度:)vad=input(請輸入請輸入AB的坐標方位角度的坐標方位角度:)vam=input(請輸入請輸入AB的坐標方位角分的坐標方位角分:)vas=input(請輸入請輸入AB的坐標方位角秒的坐標方位角秒:)va=ad+am/60+as/3600;va=a*pi/180;vbx = ax + s * cos(a) vby = ay + s * sin(a) v2）坐標反算坐標反算v設設A A、B B兩已知點的坐標分別為兩已知點的坐標分別為 和和 v，則直線，則直線ABAB的坐標方位角的坐標方位角 和水平距離

40、和水平距離 為：為：v v ),(AAyx),(BByx2222)()(arctanarctanABABABABABABABABABAByxyyxxSxyxxyyvclearvclcvx=2507.69 2299.82 2192.43 2179.72 2166.72;vy=1215.63 1303.80 1383.98 1556.41 1757.29;vplot(x,y)vhold onvfor i=1:5v plot(x(i),y(i),o)vendvfor i=1:4 a(i)=x(i+1)-x(i); b(i)=y(i+1)-y(i); s(i)=sqrt(a(i)2+b(i)2);en

41、dvs %給出各邊長給出各邊長vfor i=1:4v if a(i)=0 m(i)=atand(abs(b(i)/a(i); v %結果是度，得到是與坐標縱軸結果是度，得到是與坐標縱軸x軸構成的銳角軸構成的銳角 k(i)=m(i); if a(i)0 & b(i)0 c(i)=k(i); elseif a(i)0 c(i)=180-k(i); elseif a(i)0 & b(i)0v c(i)=90;v elsev c(i)=270;v endv endv d=fix(c)v e=fix(c-d)*60)v f=(c-d)*60-e)*60v g=d e fv v 2.可以進

42、行地形圖繪制（若有可以進行地形圖繪制（若有11個碎部點坐標）個碎部點坐標）v 115.14v129.56v138.66v2112.55v2210.51v232.70v3122.03v3222.97v3332.94v4255.58v4324.55vfid=fopen(E:dengzhixian.txt,r);vdata = fscanf(fid,%f %f %f,3,inf);v%按行取，按列存按行取，按列存vdata=data;vfclose(fid);vx =data(:,1);y =data(:,2);z =data(:,3);vxb=min(x); xe=max(x); %得到最小值最大

43、值得到最小值最大值vyb=min(y); ye=max(y); %得到最小值最大值得到最小值最大值vXl=linspace(xb,xe,200); vYl=linspace(yb,ye,200); vX,Y=meshgrid(Xl,Yl); %坐標矩陣坐標矩陣vmethod=cubic; %插值方法插值方法vZ =griddata(x,y,z,X,Y,method);vfigure; %準備圖形窗口準備圖形窗口vhold on; vcontour(X,Y,Z);%畫等高線畫等高線vcolormap cool; vbox on; %加邊框加邊框vxlabel(x);ylabel(y);vaxis

44、 equal;vtitle(等值線圖等值線圖:(單位單位:scale);vhold offv3.導線網近似平差計算導線網近似平差計算v 規(guī)范規(guī)定，對于二級及以下等級的圖根導線規(guī)范規(guī)定，對于二級及以下等級的圖根導線允許對單一導線、單結點導線網采用近似平差方允許對單一導線、單結點導線網采用近似平差方法進行計算。導線近似平差的基本思路是將角度法進行計算。導線近似平差的基本思路是將角度觀測值和邊長觀測值觀測值和邊長觀測值分別分別進行平差計算。進行平差計算。v附合導線附合導線坐標近似平差計算步驟如下：坐標近似平差計算步驟如下：v（1 1）方位角計算；）方位角計算；v（2 2）由轉角公式計算其角度閉合差；

45、）由轉角公式計算其角度閉合差；v（3 3）將角度閉合差除以測邊數分配到各觀測角中；將角度閉合差除以測邊數分配到各觀測角中；v（4 4）計算調整后的導線邊的方位角；）計算調整后的導線邊的方位角；v（5）計算坐標增量；）計算坐標增量；v（6）計算坐標增量閉合差并分配；）計算坐標增量閉合差并分配；v（7）用調整后的方位角計算各點坐標。）用調整后的方位角計算各點坐標。v對于閉合導線和支導線的情況可類似計算對于閉合導線和支導線的情況可類似計算v注意：在進行控制網精平時，一般采用間接注意：在進行控制網精平時，一般采用間接平差，誤差方程形式統(tǒng)一、規(guī)律性強、便于平差，誤差方程形式統(tǒng)一、規(guī)律性強、便于程序開發(fā)；

46、而使用條件平差等，其方程形式程序開發(fā)；而使用條件平差等，其方程形式多樣、規(guī)律性差，不利于程序開發(fā)。所以，多樣、規(guī)律性差，不利于程序開發(fā)。所以，一般建議使用間接平差方法進行程序編寫。一般建議使用間接平差方法進行程序編寫。 如下圖所示，已知四控制點平面坐標分如下圖所示，已知四控制點平面坐標分別為：別為：A（599.951,224.856）B（704.816,141.165）、）、C（747.166,572.726）、）、D(889.339,622.134)，各觀測角數據如，各觀測角數據如下：下：1 =741030，2 =2790512，3=675529，4 =2761011，5 =802346，其導，其導線各邊長分別為：線各邊長分別為：S1=143.825m ，S2=124.777m,S3=188.950m,S4=117.338m,vA Dv 2 4v s1v 1 s2 s3 s4 5v B Cv 3vn=input(請輸入導線邊數請輸入導線邊數); %觀測數據輸入觀測數據輸入vfor i=1:nv a=num2str(i);v string1=(char(請輸入第請輸入第),