一种数字地形图空间坐标系与投影的全信息无损转换方法与流程

1.本发明涉及大地测量学和地理信息制图领域，尤其涉及一种数字地形图空间坐标系与投影的全信息无损转换方法。


背景技术：

2.在测绘工作和地理信息制图时，根据工作的需要，采用的制图数学基础各不相同。单就坐标系统来说，我们可能采用的坐标系统就有1954年北京坐标系、1980年西安坐标系、wgs-84坐标系、国家2000坐标系；就投影方式来说，可能采用的有高斯-克吕格(gauss-kruger)投影、墨卡托(mercator)投影、utm（universal transverse mercator）投影、兰勃特等角投影(lambert conformal conic)等。
3.坐标系统和投影的多样性，工程上的应用多数采用城市坐标系或者工程独立坐标系下的高斯-克吕格投影（gauss
ꢀ‑
kruger projection）、国家基本比例尺地形图多采用cgcs2000坐标系下的高斯投影高斯-克吕格投影（gauss
ꢀ‑
kruger projection）或者兰勃特投影（lambert projection），航海中多采用wgs-84坐标系下的墨卡托投影（mercator projection），采用何种投影是由地图所在地理区域、比例尺大小和用途等因素所决定的，不同的投影方式有各自的优点和缺点。当我们将具有不同的坐标体统和不同的投影方式的资料进行综合利用时，需要将这些数学基础不同的数字地图统一到同一个坐标系统和投影方式下，这时候就需要进行数字地形图在不同坐标系和投影方式下的变换。
4.若要完成数字地形图在不同坐标系和投影方式下的变换这个工作，一般有三个途径：其一是将原始测量数据进行坐标系转换和投影变换，然后将变换后的数据重新进行制图；其二是将数字地形图转换输入到某些专业的制图软件中，利用这些软件的投影变换功能实现数字地形图在不同坐标系和投影方式下的变换；三是对原图形进行平移、旋转和缩放，使其转换为新的坐标系和投影。第一种方法是一种严密的方法，但将耗费大量的人力，费时费力，成本高；第二种方法，成本高昂，而且会引起很多原始地物属性信息的丢失；第三种方法将产生不均匀变形，严重影响地形图的转换精度。


技术实现要素：

5.本发明的目的在于针对上述现有技术的不足，提供了一种数字地形图空间坐标系与投影的全信息无损转换方法，足够高的精度完成数字地形图在不同的空间坐标系和投影方式下的全信息和无损的转换。
6.为实现上述目的，本发明采用了如下技术方案：本发明提供了一种数字地形图空间坐标系与投影的全信息无损转换方法，包括以下步骤：s1、数据准备，将原数据文件，通过初步平移、旋转处理后，使数字地形图中显示的地物坐标与实际坐标一致，保存为dxf格式的数字地形图；s2、采用计算机高级语言编程，读取所述dxf格式的数字地形图；
s3、对读取的所述dxf格式的数字地形图的数据进行分析判断，根据读取的信息采取不同的处理方式，具体为：s31、将实体段以外的数据按顺序输出到目标文件；s32、判断所述dxf格式的数字地形图中的实体段的位置；s33、按顺序读取实体段每个实体对象的所有坐标信息；s34、对每个所述实体对象的所有坐标进行坐标转换；s35、将处理后的信息重新组合输出到目标转换文件；s4、数据处理完毕，关闭所述原数据文件和所述目标转换文件。
7.进一步，所述s2中，读取所述dxf格式的数字地形图的具体步骤为：s21、识别实体对象所在数字地形图数据文件中的段，一个完整的dxf文件由四个段和一个文件结尾组成，分别为标题段、标段、块段、实体段和文件结束段；其中，实体段包含了所有的实体对象，并记录了每个实体的名称、图层名、线型名、颜色、基面高度、厚度、坐标的属性数据；s22、数据处理软件需要识别并提取每个所述实体对象包含的坐标信息，并进行坐标配对，然后按照选定的转换方法和输入的转换参数进行坐标转换；s23、非坐标信息保持不变。
8.进一步，所述数字地形图中对实体对象坐标的判别为：寻找标识符为10～19的行数据，并记为第一行，第二行必然是坐标，第三行必然是标识符为20～29的行数据，第四行必然是坐标，将所述第二行数据和第四行数据组成坐标对。
9.进一步，设为所述坐标对在新坐标系中的坐标，转换方法为：s51、根据原数字地形图的投影方式，将所述坐标对的平面坐标通过墨卡托反算、高斯反算或兰勃特反算，转换为原坐标系中的大地经纬度坐标；s52、将所述原坐标系中的大地经纬度坐标转换为原坐标系中的空间直角坐标；s53、将所述原坐标系中的空间直角坐标转换为新坐标系中的空间直角坐标系；（1）
其中，其中，其中，其中，为三个平移参数，单位为米；为三个旋转角度参数，单位为弧度；为尺度参数，无量纲；s54、将所述新坐标系中的空间直角坐标，转换为新坐标系中的大地经纬度坐标；s55、通过墨卡托正算、高斯正算或者兰勃特正算公式将所述新坐标系中大地经纬度坐标转换为新坐标系和投影方式下的平面直角坐标。
10.进一步，设为所述坐标对在新坐标系中的坐标，转换方法为：（2）其中，表示原坐标系下的平面直角坐标，单位为米；表示新坐标系下的平面直角坐标，单位为米；
表示平移参数，单位为米；表示旋转参数，单位为弧度；表示尺度参数，无量纲。
11.进一步，每个所述坐标对转换完毕后，顺序读取并组成下一个坐标对并进行转换、数据输出，直至所述实体段中所有的实体对象坐标转换结束。
12.进一步，所述实体段中所有的实体对象转换结束后，需要重新组合，并按照原来各自在所述数字地形图中出现的顺序输出到新的转换结果文件中。
13.本发明的有益效果为：该方法针对不同的坐标系和投影的数字地形图，可采用空间七参数模型或平面四参数模型、多项式拟合模型等进行转换，易于采用计算机编程实现，转换速度快，转换精度高，转换完成后的数据文件不丢失任何信息，同时保持了各地形图要素之间的拓扑关系。
14.采用本方法进行数据处理，能显著提升数字地形图制图的效率，保证地形图的转换精度，节省大量的成本。
附图说明
15.图1 为本发明一种数字地形图空间坐标系与投影的全信息无损转换方法流程图；图2为实施例一转换前图形属性显示图；图3为实施例一转换后图形属性显示图。
具体实施方式
16.为了使本发明的目的、技术方案及优点更加清楚明白，下面结合附图，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅用以解释本发明，并不用于限定本发明。
17.请参阅图1，一种数字地形图空间坐标系与投影的全信息无损转换方法,包括以下步骤： s1、数据准备，将原数据文件，通过初步平移、旋转处理后，使数字地形图中显示的地物坐标与实际坐标一致，保存为dxf格式的数字地形图；s2、采用计算机高级语言编程，读取所述dxf格式的数字地形图；s3、对读取的所述dxf格式的数字地形图的数据进行分析判断，根据读取的信息采取不同的处理方式，具体为：s31、将实体段以外的数据按顺序输出到目标文件；s32、判断所述dxf格式的数字地形图中的实体段的位置；s33、按顺序读取实体段每个实体对象的所有坐标信息；s34、对每个所述实体对象的所有坐标进行坐标转换；s35、将处理后的信息重新组合输出到目标转换文件；s4、数据处理完毕，关闭所述原数据文件和所述目标转换文件。
18.所述s2中，读取所述dxf格式的数字地形图的具体步骤为：s21、识别实体对象所在数字地形图数据文件中的段，一个完整的dxf文件由四个段和一个文件结尾组成，分别为标题段、标段、块段、实体段和文件结束段；
其中，实体段包含了所有的实体对象，并记录了每个实体的名称、图层名、线型名、颜色、基面高度、厚度、坐标的属性数据；s22、数据处理软件需要识别并提取每个所述实体对象包含的坐标信息，并进行坐标配对，然后按照选定的转换方法和输入的转换参数进行坐标转换；s23、非坐标信息保持不变。
19.所述数字地形图中对实体对象坐标的判别为：寻找标识符为10～19的行数据，并记为第一行，第二行必然是 坐标，第三行必然是标识符为20～29的行数据，第四行必然是坐标，将所述第二行数据和第四行数据组成坐标对。
20.设为所述坐标对在新坐标系中的坐标，转换方法为：s51、根据原数字地形图的投影方式，将所述坐标对的平面坐标通过墨卡托反算、高斯反算或兰勃特反算，转换为原坐标系中的大地经纬度坐标；s52、将所述原坐标系中的大地经纬度坐标转换为原坐标系中的空间直角坐标；s53、将所述原坐标系中的空间直角坐标转换为新坐标系中的空间直角坐标系；（1）其中，
为三个平移参数，单位为米；为三个旋转角度参数，单位为弧度；为尺度参数，无量纲；s54、将所述新坐标系中的空间直角坐标，转换为新坐标系中的大地经纬度坐标；s55、通过墨卡托正算、高斯正算或者兰勃特正算公式将所述新坐标系中大地经纬度坐标转换为新坐标系和投影方式下的平面直角坐标。
21.设为所述坐标对在新坐标系中的坐标，转换方法为：（2）其中，表示原坐标系下的平面直角坐标，单位为米；表示新坐标系下的平面直角坐标，单位为米；表示平移参数，单位为米；表示旋转参数，单位为弧度；表示尺度参数，无量纲。
22.每个所述坐标对转换完毕后，顺序读取并组成下一个坐标对并进行转换、数据输出，直至所述实体段中所有的实体对象坐标转换结束。
23.所述实体段中所有的实体对象转换结束后，需要重新组合，并按照原来各自在所
述数字地形图中出现的顺序输出到新的转换结果文件中。
24.实施例一某地形图数据文件中有包含3个节点的多段线（polyline）对象（主要属性为：线宽为0.5m，颜色为red红色），和一个“河流”text对象（主要属性为：颜色为green绿色），原坐标系为1954年北京坐标系，高斯正形投影，投影中央子午线为120度，需要转换到wgs-84坐标系下的墨卡托投影（mercator projection），投影中央子午线为121度、参考纬度为31度。
25.采用空间七参数模型转换，假定转换参数为：,,,，，，“河流”text对象在原坐标系中的坐标为（550000，3520000），多段线（polyline）对象在原坐标系中的起点坐标为（550020，3520000）该数字地形图保存为dxf格式，共4744行，其中，entities实体段位于1675行至1800行，各行数据如下所示：0section2entities0text51233301f100acdbentity80623100acdbtext1055000020352000030040
2.51界河5115100acdbtext0circle512c3301f100acdbentity80100acdbcircle10550026.7562203520017.187300400.20circle51443301f100acdbentity80100acdbcircle
10550020203520000300400.20circle514d3301f100acdbentity80100acdbcircle10550038.7203203520027.875300400.20lwpolyline514e3301f100acdbentity8062
1100acdbpolyline903700430.51055002020352000010550026.7562203520017.18710550038.7203203520027.8750endsec采用计算机高级语言编制程序，完成以下步骤：（1）读取该数据文件，同时建立输出文件；（2）识别“entities”标识符所在的位置和该位置后面的第一个endsec的位置，这两个位置之间的数据即为需要处理的数据；（3）从逐行读取entities段的数据，寻找标识符为10、20、30的位置，这三个标识符后分别为实体对象的东坐标、北坐标和高程，该段中有1个acdbtext标识符、3个acdbcircle和1个acdbpolyline标识符，分别表示1个text文字注记对象、3个circle圆对象和一个polyline多段线对象，其后紧跟的10、20标识符位置后面的数字就是各自对象的坐标；（4）程序读取原数字地形图的dxf格式的数据文件时，按顺序先会读取到acdbtext标识符，并获取到（550000，3520000）坐标对并进行以下步骤的转换计算：第一步：将（550000，3520000）按照程序输入的投影参数转换为1954年北京坐标系下的大地经纬度坐标，通过高斯反算，得到的结果为=（31
°
48
′
04.5669
″
，120
°
31
′
40.7677
″
）；第二步：将1954年北京坐标系中的大地经纬度坐标转换为1954年北京坐
标系中的空间直角坐标系=（-2756086.7215，4673684.1258，3341787.9707）；第三步：将1954年北京坐标系中的空间直角坐标系按照程序输入的投影参数转换为wgs-84坐标系中的空间直角坐标系=（-2756055.7949，4673866.2764，3341729.9651）；第四步，将wgs-84坐标系中的空间直角坐标系，转换为wgs-84坐标系中的大地经纬度坐标=（31
°
48
′
0.6293
″
， 120
°
31
′
36.2377
″
）；第五步，通过墨卡托正算公式，将wgs-84坐标系中大地经纬度坐标按照程序输入的投影参数转换为平面直角坐标=（-45199.0442，3186860.292）。
26.（5）重复前面的操作，按顺序转换所有的坐标；（6）实体段中所有的实体对象坐标转换结束后，需要与其它的信息重新组合，并按照原来各自在数字地形图中出现的顺序输出到新的转换结果文件中。
27.（7）数据处理完毕后，关闭原数据文件和目标转换文件，转换结果宣告完成。
28.转换后的数据文件依旧为4744行，其中entities(实体)段位于1675行至1800行，各行数据如下所示：0section2entities0text51233301f100acdbentity80623100acdbtext10-45199
203186860300402.51界河5115100acdbtext0circle512c3301f100acdbentity80100acdbcircle10-45172203186877300400.20circle51443301f100
acdbentity80100acdbcircle10-45178.9203186860300400.20circle514d3301f100acdbentity80100acdbcircle10-45159.9203186888300400.20lwpolyline514e3301f
100acdbentity80621100acdbpolyline903700430.510-45178.920318686010-4517220318687710-45159.92031868880endsec数字地形图转换前后的图形各属性显示如图2、图3所示。
29.从上图2、图3所示的转换前后的地形图属性来看，坐标系及投影转换的结果正确无误，转换对象的主要属性保持不变， 真正做到了全信息的无损转换，同时保持了各地形图要素之间的拓扑关系。
30.需要强调的是，采用计算机高级语言编制的处理程序，输入模型对应的转换参数，然后按照本发明的思路对数据文件进行处理，为了便于理解，前面的例子很简单，只列举了空间七参数转换模型，只有少量对象的坐标需要转换，实际上，我们需要转换的数字地形图可能有上千上万个对象，也就有成千上万个坐标需要转换，坐标转换模型也可根据实际情况，采用平面四参数转换模型、多项式拟合模型等，为了便于计算机编程计算，可以先把原地形图数据进行预分析，将所有需要转换的坐标赋值给一个数组，然后进行转换
计算。但无论坐标数量是多少，坐标转换选用哪种计算模型，以上计算的原理和步骤都是类似的。
31.以上所述实施例仅表达了本发明的实施方式，其描述较为具体和详细，但并不能因此而理解为对本发明专利范围的限制。应当指出的是，对于本领域的普通技术人员来说，在不脱离本发明构思的前提下，可以采用不同的计算机高级语言编程。
32.以上所述实施例仅表达了本发明的实施方式，其描述较为具体和详细，但并不能因此而理解为对本发明专利范围的限制。应当指出的是，对于本领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干变形和改进，这些都属于本发明的保护范围。因此，本发明专利的保护范围应以所附权利要求。