控制测量学课程设计讲解

控制测量学课程设计
学校：河南城建学院
班级：测绘工程（1）班
姓名：刘英坤
学号： 061412130
组别：第一组
指导教师：张杰
实习时间： 2015.05.11--2015.05.23
目录
一、实习目的及意义 (1)
1.1 实习目的 (1)
1.2 实际意义 (1)
三、实习过程 (2)
3.1 手算结果 (2)
3.2 Excel表格程序设计 (2)
3.2.1 总体规划设计 (2)
3.2.2 详细过程 (3)
3.2.3 程序调试与运行结果 (5)
3.3 MATLAB程序设计 (5)
3.3.1 总体规划设计 (5)
3.3.2 详细过程 (6)
3.3.3 程序调试与运行结果 (6)
四、实习总结与结论 (7)
五、参考文献 (8)
附录 (9)
一、实习目的及意义
1.1 实习目的
《控制测量学》是测绘学中一门重要的基础学科，控制测量的主要服务对象包括各种工程建设，城镇建设和土地规划与管理等工作。

在控制测量中内业包括数据处理、控制网布设技术设计等工作，外业包括控制点布设工作。

控制测量学课程设计旨在通过Excel表格、C/C++或MATLAB等汇编语言编写测量学或测量平差中一些基本的公式程序，巩固和掌握测量学或测量平差中基本概念和基本技能，提高实际动手能力，并通过实际编程加深对控制测量学的理解和掌握。

1.2 实际意义
测量数据处理是测量工作中一项核心内容，在实际工作中，面对上百点的水准网或导线网，通过手工计算不仅需要很大的人力和时间，而且在平差过程中牵扯到几何、三角函数的运算，大大增加了平差的难度和精度，将同等规模的网型平差数据放到电脑上，只需几分钟便可以完成全部计算，极大的提高了工作效率，更是测量数据处理的必然趋势。

二、实习任务
通过Excel和其他编程语言，以小组为单位，以抽签的方式对下列选题进行程序设计：
任务1：测量学教材P193附合导线（第9题）的平差计算；
任务2：附合导线的平差计算；（如图2.1）；
任务3：测量学教材P193闭合导线（第8题）的平差计算；
任务4：闭合导线的平差计算（如图2.2）；
任务5：坐标的正反算P13 ；测边交会
P171测角交会P172;
任务6：坐标的正反算P13 ；边角交会P174;
任务7：测量平差教材P99水准网（条件平差）的平差计算；(要求MATLAB调用EXCEL或文本）
图2.1 附和导线计算算例图2.2 闭合导线计算算例
任务8：测量平差教材P99水准网（间接平差）的平差计算； (要求MATLAB调用EXCEL或文本）
三、实习过程
3.1 手算结果
通过手动推算该闭合导线，一方面是为了重新理清闭合导线近似坐标平差方法的基本方法和公式，另一方面通过计算结果，可以检验程序运行结果是否正确。

计算结果如下：
图3.1 手动计算结果
3.2 Excel表格程序设计
3.2.1 总体规划设计
（1）程序设计流程图
（2）数据结构设计
图3.2.1 Excel程序设计流程图
图3.2.2 Excel表格结构设计
3.2.2 详细过程
（1）数据输入
将已知数据和观测数据以此输入到表格的对应位置。

（2）公式编制
1）观测角度之和
B22=INT(SUM(B8:B20)+SUM(C8:C20)/60+INT(SUM(D8:D20)/3600))
C22=IF((SUM(C8:C16)+INT(SUM(D8:D16)/60))<60,(SUM(C8:C16)+INT(SUM(D 8:D16)/60)),SUM(C8:C16)+INT(SUM(D8:D16)/60)，INT(SUM(C8:C16)+INT(SUM(D8:
D16)/60)))
D22=IF(SUM(D8:D16)<60,SUM(D8:D16),SUM(D8:D16)-(INT(SUM(D8:D16)/60) )*60)
2）角度闭合差计算
E25=SUM(B22*3600,C22*60,D22)-(C2-2)*180*3600;
3）角度改正数计算
E8,E10,E12,E14==IF(MOD(540*3600-B22*3600-C22*60D22,C2)=0,(540*3600 *B22*3600-C22*60-D22)/C2);
4）坐标方位角计算
坐标方位角计算比较复杂，首先判断S2输入的是左角还是右角，以左角为例，
方位角计算公式为：α
后=α
前
+β
测
-180°，如果计算结果大于或等于360°，应减
去360*n（n表示计算结果整除360°的值），如果计算结果小于0，应加360°。

这是需要用到IF（）函数，具体公式如下：
I9=IF(INT((SUM(I7*3600+J7*60+K7,F8*3600+G8*60+H8)180*3600)/3600)>3 60,INT((SUM(I7*3600+J7*60+K7,F8*3600+G8*60+H8)-180*3600)/3600)-INT(INT ((SUM(I7*3600+J7*60+K7,F8*3600+G8*60+H8)-180*3600)/3600)/360),
IF(INT((SUM(I7*3600+J7*60+K7,F8*3600+G8*60+H8)-180*3600)/3600)<0,
INT((SUM(I7*3600+J7*60+K7,F8*3600+G8*60+H8)-180*3600)/3600)+360,
INT((SUM(I7*3600+J7*60+K7,F8*3600+G8*60+H8)-180*3600)/3600)));
J9=INT((INT(SUM(I7*3600+J7*60+K7,F8*3600+G8*60+H8)-180*3600)-I9*36 00)/60);
K9==(INT(SUM(I7*3600+J7*60+K7,F8*3600+G8*60+H8)-180*3600)-I9*3600-J9*60);
I11,I13,I15计算公式与I9雷同；J11,J13,J15计算公式与J9雷同；K11,K13, K15与K9计算公式雷同；只需改变计算公式中对应表格的行号即可。

5）纵横坐标增量计算
根据平差后的坐标方位角可计算出坐标增量。

纵坐标增量（ΔX）计算公式为：M7=L7*COS((I7+J7/60+K7/3600)*PI()/180); M9,M11，M13，M15与M7的计算公式类似。

横坐标增量（ΔY）计算公式为：N7=L7*SIN((I7+J7/60+K7/3600)*PI()/180); N9,N11,N13,N15与N7的计算公式类似。

6）辅助计算
边长之和计算公式：L22=SUM(L7:L19);
计算公式：M22=SUM(M7:M19);
纵坐标闭合差f
x
横坐标闭合差f
计算公式：N22=SUM(N7:N19);
y
导线全长闭合差f计算公式：J29=SQRT(SUM(POWER(ABS(M22),2),
POWER(ABS(N22),2)));
导线全长相对闭合差T计算公式：导线全长相对闭合差T为分子为1的分式所以只需计算分母即可，计算公式为：O27==INT(L22/J29);
7）改正后纵横坐标增量计算
纵坐标增量改正数计算：O7=IF(M22<=0,ABS(M22)*L7/L22,-M22*L7/L22);
改正后纵坐标增量计算：O8 =M7+O7；O9,O11,O13,O15的计算公式与O7相似；O10,O12,O14,O16的计算公式与O8相似。

横坐标增量改正数计算：P7=IF(N22<=0,ABS(N22)*L7/L22,-N22*L7/L22);
改正后纵坐标增量计算：P8 =N7+P7;P9,P11,P13,P15的计算公式与P7相似；P10,P12,P14,P16的计算公式与O8相似。

8）导线纵横坐标计算
导线纵坐标计算：Q8=SUM(Q6,O8);
导线横坐标计算：R8==SUM(P8,R6);
3.2.3 程序调试与运行结果
图3.2.3 Excel程序调试结果
3.3 MATLAB程序设计
3.3.1 总体规划设计
（1）程序设计流程图
图3.3.1 MATLAB程序设计流程图
（2）数据结构设计通过Excel表格制作的起算数据（如图所示）：
图3.3.2 MATLAB程序设计起算数据设计
3.3.2 详细过程
（1）读取数据
MATLAB内置xlsread函数，可以读取Excel表格中的数据，但对于闭合导线网型不同，数据格式略有变化，大致分为三种情况：已知一点坐标和一个起算方位角、已知两点坐标且都位于网型内、已知两点坐标且有连接角，对应三种情况，闭合环模式以此为0、1、2,对于任务四属于第一种模式：Model=0。

（2）数据计算
在闭合导线网近似坐标推算过程中，可以分为七项工作包括角度闭合差、坐标方位角、坐标增量计算及闭合差调整、坐标计算等，主要的计算公式可以参照[5]，在此不过多赘述。

在角度闭合差调整过程中，按照“反其符号，平均分配”的原则对个导线转折角的观测值进行改正。

在MATLAB中，角度闭合差计算较为简单，首先需要将角度观测值转化为以度为单位的值，由于转化后的观测值为一行（或列）的矩阵，可直接利用sum()函数将其求和，然后判断角度闭合差是否大于限差，如果小于限差执行改正观测值，否则退出退出程序。

上述详细代码见附录Angelsure子函数。

由于坐标方位角的起算数据为度分秒格式，所以需要将其转化为十进制度形式，此时需要用到fix（）取整函数，通过Inputdlg输入对话框，提示输入左角（“L/l”）或右角（“R/r”），通过判断左角还是右角调用相应的计算公式。

详细代码见附录renum子函数和主程序代码。

由于MATLAB可以直接进行矩阵运算，按照上述公式直接编写代码，省去了过多的循环语句，然后判断导线全长相对闭合差是否小于限差，当导线相对闭合差在允许范围以内时，可将坐标增量依照“反其符号，按比例分配”的原则对纵、横坐标增量进行改正。

上述详细代码见附录CoordValue子函数。

（3）数据输出
在MATLAB中，plot函数可以将计算数据图形形式输出，fprintf函数可以将计算数据文本文件格式输出。

在图形输出编辑中，首先建立一个figure界面，然后通过plot函数通过设置其中参数，改变点、线的显示情况。

比如在测绘标准中，已知点一般用正上三角形表示，可以通过”^”参数改变该点的视图。

详细代码见附录drawskhMap子函数。

3.3.3 程序调试与运行结果
（1）程序调试及改正
坐标方位角计算错误：在坐标方位角推算过程中，对应的内角进行计算，而在起算数据中观测角是从第一个或第二个一直进行输入的，因此需要将观测角数
据的第一行移至最后一行。

（2）程序运行结果
四、实习总结与结论
导线平差计算是测绘专业技术人员必须具备的能力，Excel是一款优秀的电子表格软件，几乎在每一台电脑都有安装该办公软件，同时在科学计算中，作为专业的电子表格Excel，也常被用作数据处理，它可以实现制表、计算、统计、逻辑分析等强大功能，无论闭合导线，还是附和导线网平差计算，需要的正是这种信息化电子产品。

在利用Excel进行导线网平差过程中，它提供给了数据和公式的自动填充、表格格式的自动套用、自动求和、自动计算、记忆式输入、选择列表、自动更正、拼写检查、审核、排序和筛选等众多功能。

通过Excel进行闭合导线处理，可以制作出与手工计算一样的表格，可以清晰、直观的显示计算的全部过程，而且便于检核和修改。

在Excel中进行控制网平差处理，不依赖于计算机语言，也不需要任何汇编语言，在运算时间上不逊于编程上的计算，通过该课程设计，意识到曾经被我们忽视的Excel，竟有如此强大的功能，所以今后无论对那款软件，其功能都不不可小觑。

MATLAB是一种解释执行的语言，灵活方便，且比Visual Basic和C/C++语言更加接近我们平时书写计算公式的思维方式，用MATLAB语言编写程序犹如在演算纸上排列出公式并求解问题，以此，MATLAB易学易懂，编程效率高。

在控制网平差中，以模块的形式编写代码，不仅可以使用与任何网型，而且可以重复调用，省去了大量时间。

通过该课程设计，全面系统的理解了编译原理程序构造的一般原理和基本实现方法。

把学过的计算机平差原理以及相关课程的知识强化，能够把课堂上学的知识通过自己设计的程序表示出来，加深了对理论知识的理解和对自己动手能力的提高。

以前对与计算机平差系统概念上的认识是模糊的，现在通
图3.3.3 MATLAB程序设计结果
过自己动手做实验，从实践上认识了操作系统是如何处理命令的，如何协调计算机内部各个部件运行记忆是怎么运行的，对计算机平差原理的认识更加深刻。

在这次课程设计中，通过自己动手，加深理解编程系统的内部功能及内部实现，培养实践动手能力和程序开发能力的目的。

把死板的课本知识变得生动有趣，激发了学习的积极性。

通过两周的时间，我们分别运用Excel表格和MATLAB编程的方法，对任务四中的闭合导线网平差算例进行平差处理，其中我们发现无论那种方法，首先需要掌握牢靠的理论知识外，其次需要一定的工作经验，更需要有耐心的调试代码，并仔细排查错误。

比较两种平差方法，各有千秋，在将该参加工作之后，都可以根据工作情况稍作修改，便可运用到实际工作中去。

五、参考文献
[1] 董云.在Excel表中进行单一导线简易平差的方法[J].地矿测绘.2004年2月，第20卷第2期:21-25
[2] 王芳.Excel在控制测量学中的应用[J]. 内江师范学院学报，2013年8月，第28卷第8期：96-100
[3] 杨冲.如何用EXCEL进行附和导线近似平差的计算[J].山西建筑.2012年5月，第38卷第14期：238-239
[4] 董春来.《MATLAB语言及测绘数据处理应用》[M].成都：西南交通大学出版社，2011年8月:161-175
[5] 顾孝烈等.《测量学（第四版）》[M].上海：同济大学出版社，2011年2月：151-169
[6] **.《基于MATLAB的导线网坐标计算》[D].江西抚州：东华理工大学长江学院
附录
%主函数--------------------------------------
function Closenetwork
[Angle,Disnce,Coord,Azimu,Model,SideNum] = Getdata;
prompt = { '请输入观测角方向：'};
Title = '参数输入';
Line = 1;
type = char(inputdlg(prompt,Title,Line));
[rAzimu rAngel disAngel] = Angelsure(Angle,Model,Azimu,Coord,type);
[rCoordX rCoordY] = CoordValue(Disnce,rAzimu,Coord,Model);
rCoord = [rCoordX rCoordY];
FigTitle = '闭合导线网概略图';
hWnd = drawskhMap(rCoordX,rCoordY,FigTitle,Model);
pause(0.6);
ButtonName = questdlg('请问是否需要导出所算数据！', '系统提示', ...
'导出', '取消','导出');
switch ButtonName
case '导出',
[filename pathname] = uiputfile({'*.txt','文本文档(*.txt)';...
'*.xls','Excel(*.xls)'},' 导出平差后数据');
if isequal(filename,0) || isequal(pathname,0)
return
end
pathfilename = strcat(pathname,filename);
Preser(rAzimu,rAngel, disAngel,rCoord,pathfilename);
close(hWnd);
case '取消',
close(hWnd);
end
------------------------------------------------------------------
%读取数据
function [Angle,Disnce,Coord,Azimu,Model,SideNum] = Getdata
[FileName PathName] = uigetfile({'*.pnw','闭合网数据(*.pnw）';'*.xls','Excel(*.xls)'},' 读取数据');
format long
fid = fopen (strcat(PathName,FileName),'r+');
if fid ==0
return;
end
SideNum = xlsread(FileName,1,'A2:A2');
Model = xlsread(FileName,1,'B2:B2');
Angle = xlsread(FileName,1,strcat('C2:',['E',num2str(SideNum+1)]));
if Model ==0
Coord = xlsread(FileName,1,'G2:H2');
Azimu = xlsread(FileName,1,'I2:I2');
Disnce = xlsread(FileName,1,strcat('F2:',['F',num2str(SideNum+1)]));
elseif Model ==1
Coord = xlsread(FileName,1,'G2:H3');
Disnce = xlsread(FileName,1,strcat('F2:',['F',num2str(SideNum)]));
Azimu = 0;
else
Coord = xlsread(FileName,1,'G2:H3');
Disnce = xlsread(FileName,1,strcat('F2:',['F',num2str(SideNum+1)]));
Azimu = 0;
end
---------------------------------------------------------------------
%角度闭合差计算、角度改正和坐标方位角计算
function [rAzimu rAngel disAngel] = Angelsure(Angle,Model,Azimu,Coord,type) %format long
[row col] = size(Angle);
cAngle = Angle(1,:);
Angle(1,:) = [];
Angle(row,:) = cAngle;
cAngel = Angle(1:row,1) + Angle(1:row,2)./60 + Angle(1:row,3)./3600; disAngel = sum(cAngel)- 180*(row-2);
if abs(disAngel) >= (40*sqrt(row))
msgbox('角度闭合差超限,请重新测量！');
return;
end
rAngel = cAngel - disAngel/row;
if Model ==0
Azimulen = length(Azimu);
for i=1:Azimulen
dAzimu(i) = fix(Azimu(i));
mAzimu(i) = fix((Azimu(i)-dAzimu(i))*100);
sAzimu(i) = ((Azimu(i)-dAzimu(i))*100-mAzimu(i))*100;
cAzimu(i) = dAzimu(i)+mAzimu(i)/60+sAzimu(i)/3600;
end
%rAzimu=Azimusure(cAzimu);
else
cAzimu =zeros(1,row+1);
cdCoord = Coord(1,1)- Coord(2,1);
clCoord = Coord(1,2)-Coord(2,2);
if cdCoord ==0 && clCoord >0
cAzimu(1) = 90;
elseif cdCoord ==0 && clCoord <0
cAzimu(1) = 270;
end
if cdCoord >0 && clCoord >=0
cAzimu(1) = atan(abs(clCoord/cdCoord));
elseif cdCoord <0 && clCoord >0
cdCoord = abs(cdCoord);
cAzimu(1) =(pi-atan(abs(clCoord/cdCoord)))*180/pi;
elseif cdCoord <0 && clCoord <=0
cdCoord = abs(cdCoord);
clCoord = abs(clCoord);
cAzimu(1) =(pi+atan(abs(clCoord/cdCoord)))*180/pi;
elseif cdCoord >0 && clCoord <0
clCoord = abs(clCoord);
cAzimu(1) =(2*pi-atan(abs(clCoord/cdCoord)))*180/pi;
end
cAzimu(1) = cAzimu(1)+180;
end
if strcmp(type,'L') || strcmp(type,'l')
for i=1:row
cAzimu(i+1) = cAzimu(i) +rAngel(i)-180;
if cAzimu(i+1)>=360
cAzimu(i+1) = cAzimu(i+1)-fix(cAzimu(i+1)/360)*360;
elseif cAzimu(i+1)< 0
cAzimu(i+1) =360 + cAzimu(i+1);
end
end
end
if strcmp(type,'R') || strcmp(type,'r')
for i=1: row
cAzimu(i+1) = cAzimu(i) -rAngel(i) +180;
if cAzimu(i+1)>=360
cAzimu(i+1) = cAzimu(i+1)-fix(cAzimu(i+1)/360)*360;
elseif cAzimu(i+1)< 0
cAzimu(i+1) =360 + cAzimu(i+1);
end
end
end
rAzimu =cAzimu;
SECode = abs(rAzimu(1)-rAzimu(row+1));
if abs(SECode)/3600 > 0.00001
errordlg('计算坐标方位角未闭合，请检查数据是否存在粗差！','系统提示','on'); end
rAzimu(row+1) =[];
----------------------------------------------------------------------
%坐标增量计算及调整
function [rCoordX rCoordY] = CoordValue(Disnce,rAzimu,Coord,Model)
%format long
len = length(Disnce);
if Model==0 || Model==2
for i=1:len
cIncrenX(i) = Disnce(i)*cos(rAzimu(i)*pi/180);
cIncrenY(i) = Disnce(i)*sin(rAzimu(i)*pi/180);
end
Disncesure = sqrt((sum(cIncrenX))^2 + (sum(cIncrenY))^2);
if Disncesure/sum(Disnce) > 1/2000
msgbox('距离闭合差超限，请重新测量距离！');
return
end
[row col] =size(Disnce);
rCoordX = Coord(:,1);
rCoordY = Coord(:,2);
for i=1:row
rIncrenX(i)=vpa(cIncrenX(i)-sum(cIncrenX)*Disnce(i)/sum(Disnce),6);
rCoordX(i+1)=rIncrenX(i) + rCoordX(i);
rIncrenY(i)= vpa(cIncrenY(i)-sum(cIncrenY)*Disnce(i)/sum(Disnce),6);
rCoordY(i+1) = rCoordY(i) + rIncrenY(i);
end
end
if Model==1
for i=1:len
cIncrenX(i) = Disnce(i)*cos(rAzimu(i+1)*pi/180);
cIncrenY(i) = Disnce(i)*sin(rAzimu(i+1)*pi/180);
end
sum(cIncrenX)
sum(cIncrenY)
Disncesure = sqrt((sum(cIncrenX))^2 + (sum(cIncrenY))^2);
if Disncesure/sum(Disnce) > 1/2000
msgbox('距离闭合差超限，请重新测量距离！');
return
end
[row col] =size(Disnce);
rCoordX(1,1:col) = Coord(:,1);
rCoordY(1,1:col) = Coord(:,2);
for i=1:col
rIncrenX(i)=vpa(cIncrenX(i)-sum(cIncrenX)*Disnce(i)/sum(Disnce),6);
rCoordX(i+1)=rIncrenX(i) + rCoordX(i);
rIncrenY(i)= vpa(cIncrenY(i)-sum(cIncrenY)*Disnce(i)/sum(Disnce),6);
rCoordY(i+1) = rCoordY(i) + rIncrenY(i);
end
end
----------------------------------------------------------------------
%导线概略图显示
function hWnd=drawskhMap (x,y,Title,Model)
ScreenSize = get(0,'ScreenSize');
len = length(x);
hWnd = figure('Color',[1,1,0.5]',...
'Position',[(ScreenSize(3)-650)/2,(ScreenSize(4)-500)/2,650,500],...
'Menu','none',...
'NumberTitle','off',...
'Name',Title);
if Model ==0
h=plot(y(1),x(1),'^',y(2:len-1),x(2:len-1),'o',y(len),x(len),y,x,'Color','g','LineWidth',2,'MarkerSize',1 2);
else
h=plot(y(1:2),x(1:2),'^',y(3:len-1),x(3:len-1),'o',y(len),x(len),y,x,'Color','g','LineWidth',2,'MarkerSiz e',12);
end
pause(0.2);
----------------------------------------------------------------------
%计算成果保存
function Preser(rAzimu,rAngel,disAngel,rCoord,pathfilename)
if ~ischar(pathfilename)
error('该参数为字符串，请输入字符串！');
end
filelen = length(pathfilename);
str = pathfilename(filelen-3:filelen);
%dstr = length(data);
fid = fopen(pathfilename,'w+');
len = length(rCoord);
reCoord = reshape(rCoord,len/2,2);
r=length(rAngel);
if strcmp(str,'.txt')
fprintf(fid,'观测改正后值为：\r\n');
fprintf(fid,' %5f \r\n',rAngel);
fprintf(fid,'角度闭合差为：\n');
fprintf(fid, [num2str(3600*disAngel),'″']);
fprintf(fid,' 角度闭合差限差为：\n');
fprintf(fid,[num2str(40*sqrt(r)),' ″\r\n']);
fprintf(fid,'坐标方位角为：\r\n');
fprintf(fid,' %5f \r\n',rAzimu);
fprintf(fid,'待定点坐标为：\r\n');
fprintf(fid,' 横坐标（X/m）纵坐标（Y/m）\r\n');
fprintf(fid,' %3f %3f\r\n',reCoord');
end
if strcmp(str,'.xls')
xlswrite(pathfilename,'观测改正后值','A1:A1');
xlswrite(pathfilename,rAngel,['A2:A',num2str(r+2)]);
PA = ['A',num2str(r+4),':A',num2str(r+4)];
xlswrite(pathfilename,'角度闭合差',PA);
PB = ['B',num2str(r+4),':B',num2str(r+4)];
xlswrite(pathfilename,num2str(3600*disAngel),PB);
PB = ['C',num2str(r+4),':C',num2str(r+4)];
xlswrite(pathfilename,'角度闭合差限差',PC);
PD =['D',num2str(r+4),':D',num2str(r+4)];
xlswrite(pathfilename,num2str(40*sqrt(r)),PD);
xlswrite(pathfilename,'坐标方位角','B1:B1');
xlswrite(pathfilename,rAzimu,['B2:B',num2str((dstr-1)/3)]);
xlswrite(pathfilename,'待定点坐标','D1:D1');
xlswrite(pathfilename,reCoord',['D2:F',num2str((dstr-1)/3)]);
end
--------------------------------------------------------------------
%转化十进制角度
function result = renum(c,n)
if nargin<=1
error('参数输入错误，缺少必要元素！')
end
if c>=0
cresult = (c-fix(c))*10^n;
result = fix(c)+round(cresult); elseif c<0
c=abs(c);
cresult = (c-fix(c))*10^n;
result = fix(c)+round(cresult);
result = -1*result;
end。