坐标方位角算法在GIS中的应用(最后修改)

基金颁发部门:国家自然科学基金委；项目名称：基于剖面的三维拓扑地质建模研究；编号:40572165；基金申请人:侯恩科坐标方位角算法在GIS中的应用王晓芳1,侯恩科2，董汉军3(1陕西职业技术学院计算机系，西安710100；2西安科技大学地质与环境工程学院,西安710054)摘要：在利用坐标反算方位角的反正弦、反余弦和反正切公式的基础上，给出了反正弦公式的C语言程序，又利用方位角的特性，构造了方位角排序算法，即半圆周排序和圆周排序算法。然后，将坐标反算求取方位角公式应用到了空间方向描述之中，并且将半圆周排序算法应用在了多边形面积计算中，圆周排序法应用在Voronoi图的构建之中。实验证明该坐标反算公式和两种排序算法在GIS图形算法研究方面很有实用价值。关键词：方位角；坐标反算；Voronoi图；空间方向描述ApplicationofAzimuthAngleAlgorithmontheGISWangXiaofang1,HouEnke2,DongHanjun3（1.computerdepartmentofshan’xioccupationalCollege,xi’an710100.China2.Dept.ofGeologyandEnvironmentEngineering,Xi’anUniversityofScienceandTechnology,Xi’an710054，China）Abstract：Basingonthedetailedintroductionofthecoordinate-inversecalculationformulas,suchasarc-sineformula,arc-cosineformulaandarc-tangentformula,thispaperhasshowntheClanguageofarc-sineformula,Inthesametime,itconstructsthesortalgorithm,namelyhalf-circlesortalgorithmandcirclesortalgorithm,byusingthespecialityofazimuthangle.Thentheauthorsapplytheformulasofcoordinate-inversecalculationtothespatialdirectiondescription.thehalf–circlesortalgorithmtocalculationofpolygon’sareaaswell.,AndcirclesortalgorithmisusedtodrawVoronoidiagram.IthasprovedthatthemethodsofcalculationandtwosortalgorithmsarepricelessinthestudyofGISgraphs.Keywords：azimuthangle;coordinate-inversecalculation;Voronoidiagram;spatialdirectiondescription1引言方位角和方向角都是用来定量描述方向关系的角。从测量学的角度上讲，方位角是从正北方向起，顺时针旋转到某一位置时所经过的角度。从几何学的角度上讲，方向角就是在二维平面直角坐标系下，与x轴正向逆时针方向形成的角。从本质上讲，方位角和方向角是等同的。很多学者在方位角的计算和应用方面作了很多研究。谢振波等给出了用BASIC语言编写的坐标方位角反算算法的代码，[1]罗来恩给出反算坐标方位角的算法步骤，但他们的代码作者简介：王晓芳(1978-)，女，山东郓城人，硕士，陕西职业技术学院计算机系讲师，主要研究方向为计算机图像处理和地理信息系统。侯恩科（1960-），男(汉)，博士生导师，西安科技大学地质与环境学院，主要从事地理信息系统与三维可视化方面的研究。基金项目：国家自然科学基金项目(40572165).都不够完善，没有考虑特殊情况的处理。[2]陈德标]对坐标方位角的计算给出较详细的方法，对特殊情况的坐标方位角做出了判断，它主要是在文献2的基础上进行了改进。文献2和文献3都具体给出了计算坐标方位角的公式。[3]涂群生则将坐标反算方位角应用到了可编程计算器中。[4]冯大幅直接用判断象限的方法来反算坐标方位角，[5]相比之下，陈德标对坐标反算求取方位角从原理和方法上都较为系统、全面。2坐标反算方位角及其公式坐标正算就是根据已知点坐标、已知边长和坐标方位角来求解未知点的坐标，而坐标反算则是根据两个已知点的坐标来反算其边长与坐标方位角。[6]假设A点为后视点，则首先判断∆x=xB-xA，∆y=yB-yA。式1给出了反算后方位角的取值范围。0,0)...arctan(1800)...arctan(1800,0)...arctan(3600,0)...arctan(0,0...2700,0...90yxxyxxyyxxyyxxyyxyxAB（1）参考文献1~3有坐标反算的不同公式，用反正弦函数和反余弦函数来求解坐标方位角，如公式（2）（反正弦函数公式）和（3）（反余弦函数公式）所示。其中坐标增量取值范围∆x，∆y不同时为0，通式的值域为[0°,360°]，即为“真”坐标方位角。很显然，用（2）式或（3）式编程求解方位角要比用（1）式简便很多。式中sgn()为取符号函数，sgn(x)=│x│/x，当然在实际编程中以x的值来确定，x=0，sgn(x)=0；x0，sgn(x)=1；x0，sgn(x)=-1。笔者用c语言编写了上面式（2）和式（3），通过实验，所得结果都是正确的。)](sgn2)sgn()(sgn)sgn([)](sgn)([sgn90)]sgn()[sgn()]sgn()sgn(2sgn[2180)sgn()arcsin(22222))sgn(1(yyxxyxyxabsyxxSyxAB（2）)](sgn2)sgn()(sgn)sgn([)](sgn)([sgn90)]sgn()[sgn()]sgn(1sgn[360))(sgn()sgn()arccos(2222yyxxyxyxabsyxabsySxAB（3）下面是用公式（2）（即反正弦函数公式）反算方位角的c程序。值得注意的是，公式（2）、和（3）计算出来的单位是度，而下面的c程序得出的单位是弧度。intsgn(doublex){inttemp;if(x0){temp=1;}elseif(x==0){temp=0;}elseif(x0){temp=-1;}returntemp;}doubleGetFWJsin(doublex1,doubley1,doublex2,doubley2){doubletemp;doubles=sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));intsgnx,sgny;if(x2-x1==0&&y2-y1==0){exit(0);}sgnx=sgn(y2-y1);sgny=sgn(x2-x1);temp=asin((y2-y1)/s)*sgnx+pi*pow(2,(1+sgnx)/2.0)*sgn(2-sgnx-sgny)*abs(sgnx*sgny)+pi/2.0*(sgnx*sgnx-sgny*sgny)*((-1)*sgnx+sgnx*sgnx+sgny-2*sgny*sgny);returntemp;}3方位角排序算法平面上各个定点都有确定的（x，y）。为此，我们可以通过坐标反算来求取各个点的坐标方位角，并按方位角的大小进行平面上的顺时针或逆时针排序。3.1半圆周排序算法半圆周排序就是方位角的取值在[0，180°]内，如图1所示，A，B等是平面内的点，为直观期间，做一独立坐标系x′o′y′。算法设计：1）通过比较各个点的x坐标值，寻找到x坐标值最小的点，如图1中A点，把A点作为后视点；2）利用上面的公式（2）或公式（3）求出各个点与A点所构成的坐标方位角（取值范围[0，180°]）；3）通过比较方位角的大小，对各个点与A点所构成的坐标方位角进行排序，从小到大排序，得出逆时针方向的顺序，如果从大到小，就是顺时针方向的顺序。ABCDEFGOO'xx'y'ABCDEFGOO'xx'yy'图1半圆周排序图2圆周排序Fig1Half-circlesortFig2Circlesort3.2圆周排序算法类似于半圆周排序，不过圆周排序中方位角在[0，360°]之间取值。如图2所示，A，B等是平面内的点，算法设计：1）首先要确定一中间点，如图2中的o′点，将它作为后视点；2）利用上面的公式（2）或公式（3）求出各个点与A点所构成的坐标方位角（取值范围[0，360°]）；3）通过比较方位角的大小，对各个点与A点所构成的坐标方位角进行排序，从小到大排序，得出逆时针方向的顺序，如果从大到小，就是顺时针方向的顺序。4方位角及方位角排序算法的应用4.1方位角公式在空间方向描述中的应用空间方向关系与空间距离、空间拓扑一样是GIS的重要理论问题，在GIS的空间分析、空间查询和空间建模方面有着重要的应用。空间方向的定性描述是用若干主方向粗略地描述空间方向。而定量描述则是用方位角来量测空间目标之间的方向关系，因此方位角是空间方向描述的一个重要手段。通过方位角可以知道线段AB相对于北方向的角度，只要有精确的参照系统，空间方向就可以确定下来。经纬仪测图、放样等就是基于此理。因此，式（2）和（3）无疑在空间方向描述中有着重要的应用，为空间方向关系描述和地理空间推理的研究奠定了良好的基础。4.2方位角算法在Voronoi图构建中的应用利用D-TIN可以构建平面点集Voronoi图，[8]其算法如下：（1）先构建离散平面点集的D-TIN；（2）求取每个三角形的外接圆的圆心；（3）对每一个离散点，利用方位角排序算法来连接三角形的外接圆心。下面详细描述方位角算法来对离散点进行顺时针或逆时针排序。假设A’、B’、C’等为外界圆的圆心，如图3所示。此时Voronoi多边形所包含D-TIN的离散点已知。运用式（2）或（3）求出每个圆心点与D-TIN的离散点（图3中o′点）所构成的坐标方位角，此时所求的坐标方位角的取值范围为[0°，360°]，利用圆周排序算法求得方位角序列，顺序连接方位角相应的圆心点就可以得到一个Voronoi多边形（如图3所示）。将各离散点的Voronoi图形成集合，即得到平面点集的Voronoi图。在图3中，O，A，B，C等为构成D-TIN的离散点，A’，B’，C’等为各三角形的外接圆圆心，则构成Voronoi多边形A’B’C’D’E’的顶点坐标A’、B’、C’等与O点的坐标方位角，依次有αOA′αOB′αOC′αOD′αOE′。从而可以确定连接次序为A′、B′、C′、D′、E′。AEDBCE'A'B'C'D'OOxyABCDE图3D-TIN和Voronoi图图4多边形面积求解图Fig3D-TINandVoronoidiagramFig4Theareasolutionofpolygon4．3方位角排序算法在多边形面积计算中的应用公式（5）式可用来计算多边形的面积。此公式只适用于凸多边形，对于非凸多边形，我们可以把非多边形分割成若干个凸多边形来进行计算。当然，前提条件是对这些点必须按顺时针或逆时针进行排序，而后才能计算。[7]具体算法如下：首先判断多边形是否凸多边形，如果是凸多边形，则应用方位角半圆周法进行排序，然后运用公式（5）进行计算；如果不是凸多边形，将其分割成若干凸多边形，然后对每一个凸多边形的顶点应用方位角半圆周法进行排序，再对每一个顶点排过序的凸多边形分别运用公式（5）进行计算，再求和就可以得到凸多边形的面积