C# + NetTopologySuite读取Shp存入PostGIS数据库

之前文章中，我们下载了西安路网数据Shp。分享一个下载各省市边界、路网等数据的地址-CSDN博客

今天我们将Shp数据存入PostGIS空间数据库。

VS项目Nuget引入NetTopologySuite，NetTopologySuite.IO.Esri.Shapefile，NetTopologySuite.IO.PostGis库

<PackageReference Include="NetTopologySuite" Version="2.5.0" />
<PackageReference Include="NetTopologySuite.IO.Esri.Shapefile" Version="1.0.0" />
<PackageReference Include="NetTopologySuite.IO.PostGis" Version="2.1.0" />

创建一个model，用来存储路网数据中的属性数据和空间数据：

public class RoadModel
{
    public long OsmId { get; set; }

    public int RoadTypeCode { get; set; }

    public string RoadType { get; set; }

    public string RoadName { get; set; }

    public string RoadRef { get; set; }

    public string Oneway { get; set; }

    public int MaxSpeed { get; set; }

    public string Bridge { get; set; }

    public string Tunnel { get; set; }

    /// <summary>
    /// geometry: wkb
    /// </summary>
    public byte[] RoadShape { get; set; }
}

读取Shp文件的所有feature（涉及到Geometry，WKB的转换）：

private const string roadShpFilePath = @"D:\workdata\xian\610100_new.shp";

private List<RoadModel> GetAttributesFromShapeFile()
{
    var features = NetTopologySuite.IO.Esri.Shapefile.ReadAllFeatures(roadShpFilePath);
    var pgWriter = new PostGisWriter();
    var roadModels = new List<RoadModel>();

    foreach (var feature in features)
    {
        //此时通过NTS库读取的Geometry SRID为0，没有坐标系
        var geometry = feature.Geometry;
        
        //这一步为Geometry添加4326坐标系
        var geometryForWGS84 = GeometryFactory.Default.WithSRID(4326).CreateGeometry(geometry);
        
        //如果这一步直接使用NTS库转换Geometry为WKB，此时得到的WKB还是没有坐标系的
        //var wkbForWGS84 = geometryForWGS84.ToBinary();
        
        //需要通过NTS PG库将含有坐标系的Geometry转成含有坐标系的EWKB
        var wkbForWGS84 = pgWriter.Write(geometryForWGS84);

        roadModels.Add(new RoadModel
        {
            OsmId = Convert.ToInt64(feature.Attributes.GetOptionalValue("osm_id")),
            RoadTypeCode = Convert.ToInt32(feature.Attributes.GetOptionalValue("code")),
            RoadType = feature.Attributes.GetOptionalValue("fclass").ToString(),
            RoadName = feature.Attributes.GetOptionalValue("name").ToString(),
            RoadRef = feature.Attributes.GetOptionalValue("ref").ToString(),
            Oneway = feature.Attributes.GetOptionalValue("oneway").ToString(),
            MaxSpeed = Convert.ToInt32(feature.Attributes.GetOptionalValue("maxspeed")),
            Bridge = feature.Attributes.GetOptionalValue("bridge").ToString(),
            Tunnel = feature.Attributes.GetOptionalValue("tunnel").ToString(),
            RoadShape = wkbForWGS84,
        });
    }
    return roadModels;
}

在PostGIS数据库中创建表road_xian并赋予坐标系：

CREATE TABLE public.road_xian (
	osm_id int8 NOT NULL,
	road_type_code int4 NULL,
	road_type varchar NULL,
	road_name varchar NULL,
	road_ref varchar NULL,
	oneway varchar NULL,
	max_speed int4 NULL,
	bridge varchar NULL,
	tunnel varchar NULL,
	shape public.geometry NULL,
	CONSTRAINT road_xian_pk PRIMARY KEY (osm_id)
);
--此时表的shape字段坐标系为0，所以执行下面这句sql为shape字段添加上4326 WGS84坐标系

SELECT UpdateGeometrySRID('road_xian','shape',4326);

PostGIS数据库的安装可参考之前的文章Linux 安装PostGIS 3.1-CSDN博客

VS项目Nuget再引入Npgsql，Npgsql.LegacyPostgis库：

<PackageReference Include="Npgsql" Version="8.0.1" />
<PackageReference Include="Npgsql.LegacyPostgis" Version="5.0.17" />

通过NpgsqlCommand将数据插入db，shape(Geometry)字段需要用DbType.Object类型：

private int InsertIntoPg(List<RoadModel> roadModels)
{
    var allNum = 0;
    var connString = "Host=192.168.137.133;Username=postgres;Password=postgres;Database=gis_lp";
    var connection = new NpgsqlConnection(connString);
    connection.Open();

    var sql = $"INSERT INTO road_xian (osm_id, road_type_code, road_type, road_name, road_ref, oneway, max_speed, bridge, tunnel, shape) " +
                                    $"VALUES (@osm_id, @road_type_code, @road_type, @road_name, @road_ref, @oneway, @max_speed, @bridge, @tunnel, @shape)";

    foreach (var model in roadModels)
    {
        using (var command = new NpgsqlCommand(sql, connection))
        {
            var parameters = ConstructParameters(model);
            command.Parameters.AddRange(parameters.ToArray());
            var num = command.ExecuteNonQuery();
            allNum += num;
        }
    }

    return allNum;
}

private List<NpgsqlParameter> ConstructParameters(RoadModel roadModel)
{
    var osmIdPar = new NpgsqlParameter("osm_id", NpgsqlDbType.Bigint)
    {
        Value = roadModel.OsmId
    };
    var roadTypeCodePar = new NpgsqlParameter("road_type_code", NpgsqlDbType.Integer)
    {
        Value = roadModel.RoadTypeCode
    };
    var roadTypePar = new NpgsqlParameter("road_type", NpgsqlDbType.Varchar)
    {
        Value = roadModel.RoadType
    };
    var roadNamePar = new NpgsqlParameter("road_name", NpgsqlDbType.Varchar)
    {
        Value = roadModel.RoadName
    };
    var roadRefPar = new NpgsqlParameter("road_ref", NpgsqlDbType.Varchar)
    {
        Value = roadModel.RoadRef
    };
    var onewayPar = new NpgsqlParameter("oneway", NpgsqlDbType.Varchar)
    {
        Value = roadModel.Oneway
    };
    var maxSpeedPar = new NpgsqlParameter("max_speed", NpgsqlDbType.Integer)
    {
        Value = roadModel.MaxSpeed
    };
    var bridgePar = new NpgsqlParameter("bridge", NpgsqlDbType.Varchar)
    {
        Value = roadModel.Bridge
    };
    var tunnelPar = new NpgsqlParameter("tunnel", NpgsqlDbType.Varchar)
    {
        Value = roadModel.Tunnel
    };
    var geoPar = new NpgsqlParameter("shape", DbType.Object)
    {
        Value = roadModel.RoadShape
    };
    return new List<NpgsqlParameter>
    {
        osmIdPar,roadTypeCodePar,roadTypePar,roadNamePar,roadRefPar,onewayPar,maxSpeedPar,bridgePar,tunnelPar,geoPar
    };
}

启动项目，大约18s后，29902条数据成功插入db。

 

打开QGIS，在左侧PostgreSQL建立新的连接：

 

选择我们刚才创建的road_xian表，添加到当前地图项目中预览，与之前下载的Shp数据重叠一致，数据OK！