在处理地理信息系统(GIS)数据时,经常需要将多个地图段合并为一个单一的Shapefile,并将其坐标系统从OSGB1936转换为WGS84。本文将介绍如何使用MapWinGIS实现这一过程,并提供相应的代码示例。
MapWinGIS是一个用于GIS数据操作的库,它支持多种地图数据格式,包括Shapefile。Shapefile是一种常见的GIS数据格式,它包含地理数据的几何形状和属性信息。每个Shapefile由多个Shape实体和一组字段组成,字段与每个形状相关联,并且可以自定义数量和名称。
Shapefile中的每个形状由多个点组成,这些点定义了形状的坐标。常见的形状类型包括多边形。此外,点可以被分组到多个部分中,部分通过添加信息到形状上来定义,这些信息指定了每个部分的起始点索引。每个点都有与之相关的参数,包括坐标(x, y, z)和测量值(M)。每个形状也有为Shapefile定义的字段的值。
本文提供的代码是用VB.NET编写的,它涵盖了构建Shapefile的基础,包括合并点和字段以创建单独的形状,然后将形状合并到Shapefile中。需要注意的是,坐标转换函数Convert_OSGB36_to_WGS384_HF没有包含在代码中,但可以在Hannah Fry的网站找到Python版本的转换方法。
主函数Main接受一个包含文件名(包括完整路径)的集合和一个输出Shapefile的名称作为字符串。它处理集合中的每个Shapefile,并尝试将所有Shapefile的内容合并到一个输出Shapefile中。
为了创建新字段,需要指定字段的名称、类型和宽度。以下是VB.NET代码示例:
Dim Newfield As New MapWinGIS.Field
Newfield.Name = sfIn.Field(FieldIndex).Name
Newfield.Type = sfIn.Field(FieldIndex).Type
Newfield.Width = sfIn.Field(FieldIndex).Width
If Not sfOut.EditInsertField(Newfield, FieldIndex) Then
Debug.Print("Failed to add new Field. Error was " & sfOut.ErrorMsg(sfOut.LastErrorCode))
Stop
End If
创建新形状时,需要指定类型,例如多边形。在合并Shapefile时,从源Shapefile中提取类型。以下是VB.NET代码示例:
Dim ShapeType As Integer = sfIn.Shape(ShapeIndex).ShapeType
Dim NewShape As New MapWinGIS.Shape
If Not NewShape.Create(ShapeType) Then
Debug.Print("New shape creation failed, aborting import of " & _
sfIn.Filename & " at ShapeIndex = " & ShapeIndex)
Exit Sub
End If
在转换坐标时,需要逐个复制每个点,并在过程中进行转换。以下是VB.NET代码示例:
For PointIndex As Integer = 0 To sfIn.Shape(ShapeIndex).numPoints - 1
Dim X As Double = sfIn.Shape(ShapeIndex).Point(PointIndex).x
Dim Y As Double = sfIn.Shape(ShapeIndex).Point(PointIndex).y
Dim Z As Double = sfIn.Shape(ShapeIndex).Point(PointIndex).Z
Dim M As Double = sfIn.Shape(ShapeIndex).Point(PointIndex).M
Convert_OSGB36_to_WGS384_HF(Lat, Lon, Y, X)
Dim NewPoint As New MapWinGIS.Point
NewPoint.x = Lon
NewPoint.y = Lat
NewPoint.Z = Z
NewPoint.M = M
If Not NewShape.InsertPoint(NewPoint, PointIndex) Then
Debug.Print("Point add Failed, aborting import of " & _
sfIn.Filename & " at ShapeIndex = " & ShapeIndex)
Exit Sub
End If
Next
每个形状可以有一个或多个部分。部分只是点的子集。点之间的关系决定了它们的显示方式,顺时针部分被填充,逆时针部分为空。以下是VB.NET代码示例:
Dim FirstPointInPartIndex As Integer
For PartIndex As Integer = 0 To sfIn.Shape(ShapeIndex).NumParts - 1
FirstPointInPartIndex = sfIn.Shape(ShapeIndex).Part(PartIndex)
NewShape.InsertPart(FirstPointInPartIndex, PartIndex)
Next
复制字段值相对复杂,因为并非所有源Shapefile都有完全相同的字段列,因此需要按名称匹配字段,即根据字段名称查找正确的列。以下是VB.NET代码示例:
Dim OutFieldIndex As Integer
Dim AddedField As Boolean
For FieldIndex As Integer = 0 To sfIn.NumFields - 1
Dim InFieldName As String = sfIn.Field(FieldIndex).Name
AddedField = False
For OutFieldIndex As Integer = 0 To sfOut.NumFields - 1
If sfOut.Field(OutFieldIndex).Name = InFieldName Then
If sfOut.EditCellValue(OutFieldIndex, sfOut.NumShapes - 1, _
sfIn.CellValue(FieldIndex, ShapeIndex)) Then
AddedField = True
Else
Debug.Print("Failed to add cell value from import of " & _
sfIn.Filename & " at ShapeIndex = " & ShapeIndex)
End If
End If
Next
If Not AddedField Then
Debug.Print("Failed to find Filed Heading " & InFieldName & _
" in output file, skipping entry for ShapeIndex = " & ShapeIndex)
End If
Next
在使用MapWinGIS时,最大的挑战是要知道事情的顺序或顺序,例如,不能在启用编辑之前编辑形状:
ShapeFile.StartEditingShapes()
ShapeFile.StartEditingTable()
Debug.print(ShapeFile.ErrorMsg(ShapeFile.LastErrorCode))