PostGIS

创建日期:2024-06-21
更新日期:2025-01-19

PostGIS官方网站:http://postgis.net/

PostGIS在线文档:http://postgis.net/documentation/

PostGIS在线手册:http://postgis.net/docs/

pgRouting官方网站:https://pgrouting.org/

pgRouting在线手册:https://docs.pgrouting.org/3.3/en/index.html

PostGIS离线安装:http://postgis.net/docs/postgis_installation.html#install_short_version

安装PostGIS

1、打开PostgreSQL 14,Application Stack Builder。勾选Spatial Extensions下的PostGIS插件,点击下一步安装。

PostGIS 3.2.1 bundle includes PostGIS 3.2.1 w GDAL 3.4.2 (SQLite 3.30.1, OpenJPEG 2.4.0, Expat 2.4.8, FreeXL 1.0.6)  , GEOS 3.10.2, Proj 7.2.1, SFCGAL 1.4.1, pgRouting 3.3.1,  osm2pgrouting 2.3.8, ogr_fdw 1.1.1 spatial foreign data wrapper extension, and pgPointcloud 1.2.0.

创建GIS数据库

1、使用navicat连接到PostgreSQL,创建数据库nyc。

2、新建查询,输入以下代码并执行。

CREATE EXTENSION IF NOT EXISTS plpgsql;
CREATE EXTENSION postgis;
CREATE EXTENSION postgis_raster; -- OPTIONAL
CREATE EXTENSION postgis_topology; -- OPTIONAL

3、验证PostGIS是否启用。

WITH city AS (
  SELECT 'Gotham' AS name,
       ST_Buffer(ST_Point(0,0), 10) AS geom
)
 , superhero(name,geom) AS (
    VALUES
      ('Bat Boy', ST_Point(0.1,0))
    , ('Bat Girl', ST_Point(1,1) )
)

SELECT superhero.name
FROM city INNER JOIN superhero
       ON ST_Contains(city.geom, superhero.geom);

几何类型

Point:POINT (1 2)

LineString:LINESTRING (1 2, 3 4, 5 6)

LinearRing:LINEARRING (0 0 0, 4 0 0, 4 4 0, 0 4 0, 0 0 0)

Polygon:POLYGON ((0 0 0,4 0 0,4 4 0,0 4 0,0 0 0),(1 1 0,2 1 0,2 2 0,1 2 0,1 1 0))

MultiPoint:MULTIPOINT ( (0 0), (1 2) )

MultiLineString:MULTILINESTRING ( (0 0,1 1,1 2), (2 3,3 2,5 4) )

MultiPolygon:MULTIPOLYGON ~(~(~(1 5, 5 5, 5 1, 1 1, 1 5)), ((6 5, 9 1, 6 1, 6 5~)~)~)

GeometryCollection:GEOMETRYCOLLECTION ( POINT(2 3), LINESTRING(2 3, 3 4))

Triangle:TRIANGLE ((0 0, 0 9, 9 0, 0 0))

TIN:TIN Z ( ((0 0 0, 0 0 1, 0 1 0, 0 0 0)), ((0 0 0, 0 1 0, 1 1 0, 0 0 0)) )

WKT

POINT(0 0)

POINT Z (0 0 0)

POINT ZM (0 0 0 0)

POINT EMPTY

LINESTRING(0 0,1 1,1 2)

LINESTRING EMPTY

POLYGON((0 0,4 0,4 4,0 4,0 0),(1 1, 2 1, 2 2, 1 2,1 1))

MULTIPOINT((0 0),(1 2))

MULTIPOINT Z ((0 0 0),(1 2 3))

MULTIPOINT EMPTY

MULTILINESTRING((0 0,1 1,1 2),(2 3,3 2,5 4))

MULTIPOLYGON~(~(~(0 0,4 0,4 4,0 4,0 0),(1 1,2 1,2 2,1 2,1 1)), ((-1 -1,-1 -2,-2 -2,-2 -1,-1 -1~)~)~)

GEOMETRYCOLLECTION(POINT(2 3),LINESTRING(2 3,3 4))

GEOMETRYCOLLECTION EMPTY

EWKT

POINT(0 0 0) ~-~- XYZ

SRID=32632;POINT(0 0) ~-~- XY with SRID

POINTM(0 0 0) ~-~- XYM

POINT(0 0 0 0) ~-~- XYZM

SRID=4326;MULTIPOINTM(0 0 0,1 2 1) ~-~- XYM with SRID

MULTILINESTRING((0 0 0,1 1 0,1 2 1),(2 3 1,3 2 1,5 4 1))

POLYGON((0 0 0,4 0 0,4 4 0,0 4 0,0 0 0),(1 1 0,2 1 0,2 2 0,1 2 0,1 1 0))

MULTIPOLYGON~(~(~(0 0 0,4 0 0,4 4 0,0 4 0,0 0 0),(1 1 0,2 1 0,2 2 0,1 2 0,1 1 0)),((-1 -1 0,-1 -2 0,-2 -2 0,-2 -1 0,-1 -1 0~)~)~)

GEOMETRYCOLLECTIONM( POINTM(2 3 9), LINESTRINGM(2 3 4, 3 4 5) )

MULTICURVE( (0 0, 5 5), CIRCULARSTRING(4 0, 4 4, 8 4) )

POLYHEDRALSURFACE( ((0 0 0, 0 0 1, 0 1 1, 0 1 0, 0 0 0)), ((0 0 0, 0 1 0, 1 1 0, 1 0 0, 0 0 0)), ((0 0 0, 1 0 0, 1 0 1, 0 0 1, 0 0 0)), ((1 1 0, 1 1 1, 1 0 1, 1 0 0, 1 1 0)), ((0 1 0, 0 1 1, 1 1 1, 1 1 0, 0 1 0)), ((0 0 1, 1 0 1, 1 1 1, 0 1 1, 0 0 1)) )

TRIANGLE ((0 0, 0 9, 9 0, 0 0))

TIN( ((0 0 0, 0 0 1, 0 1 0, 0 0 0)), ((0 0 0, 0 1 0, 1 1 0, 0 0 0)) )

常用查询

1、查询所有GIS扩展。

SELECT name, default_version,installed_version
FROM pg_available_extensions WHERE name LIKE 'postgis%' or name LIKE 'address%';

2、启用GIS扩展。

CREATE EXTENSION postgis;
CREATE EXTENSION postgis_raster;
CREATE EXTENSION postgis_sfcgal;
CREATE EXTENSION fuzzystrmatch; --needed for postgis_tiger_geocoder;
--optional used by postgis_tiger_geocoder, or can be used standalone
CREATE EXTENSION address_standardizer;
CREATE EXTENSION address_standardizer_data_us;
CREATE EXTENSION postgis_tiger_geocoder;
CREATE EXTENSION postgis_topology;

3、查询版本。

select postgis_full_version();

4、创建表,插入和查询GIS数据。(WKT)

CREATE TABLE geotable (id int4, name varchar(255), geom geometry);

insert into geotable (id, name, geom)
values (1, 'A Place', ST_GeomFromText('POINT(-126.4 45.32)', 4326));

select id, name, ST_AsText(geom) geom from geotable;

5、插入和查询GIS数据。(EWKT)

insert into geotable (id, name, geom)

values (2, 'A City', ST_GeomFromEWKT('SRID=4326;POINT(-127.4 32.32)'));

select id, name, ST_AsEWKT(geom) geom from geotable;

6、创建地理表并插入数据。

CREATE TABLE global_points (

id SERIAL PRIMARY KEY,

name VARCHAR(64),

location geography(POINT,4326)

);

INSERT INTO global_points (name, location) VALUES ('Town', 'SRID=4326;POINT(-110 30)');

INSERT INTO global_points (name, location) VALUES ('Forest', 'SRID=4326;POINT(-109 29)');

INSERT INTO global_points (name, location) VALUES ('London', 'SRID=4326;POINT(0 49)');

  1. 查询所有地理列。

select * from geography_columns;

  1. 创建索引。

CREATE INDEX global_points_gix ON global_points USING GIST ( location );

  1. 获取所有支持的SRID。

select * from spatial_ref_sys;

  1. 查询距离点1000km内的要素。

注意:输入输出单位都是米。

SELECT name FROM global_points WHERE ST_DWithin(location, 'SRID=4326;POINT(-110 29)'::geography, 1000000);

  1. 计算点到线的距离。(球面,正确)

SELECT ST_Distance('LINESTRING(-122.33 47.606, 0.0 51.5)'::geography, 'POINT(-21.96 64.15)'::geography);

  1. 计算点到线的距离。(平面,错误)

SELECT ST_Distance('LINESTRING(-122.33 47.606, 0.0 51.5)'::geometry, 'POINT(-21.96 64.15)'::geometry);

  1. 创建空间表。

CREATE TABLE roads (

~ id SERIAL PRIMARY KEY,

~ name VARCHAR(64),

~ geom geometry(LINESTRING,3005)

~ );

CREATE TABLE geoms(gid serial PRIMARY KEY, geom geometry );

CREATE TABLE pts(gid serial PRIMARY KEY, geom geometry(POINT) );

CREATE TABLE pts(gid serial PRIMARY KEY, geom geometry(POINTZ,3005) );

CREATE TABLE lines(gid serial PRIMARY KEY, geom geometry(LINESTRINGZM) );

CREATE TABLE polys(gid serial PRIMARY KEY, geom geometry(POLYGON,4267) );

  1. 获取所有几何列。

select * from geometry_columns;

  1. 判断是否有效。

select  ST_IsValid('LINESTRING(0 0, 1 1)');

  1. 限制表的几何体有效。

ALTER TABLE mytable

~ ADD CONSTRAINT geometry_valid_check

~ CHECK (ST_IsValid(geom));

常用函数

地理数据类型支持的函数:

ST_Area – 返回面积

ST_AsBinary – 返回WKB

ST_AsEWKT – 返回EWKT

ST_AsGML – 返回GML

ST_AsGeoJSON – 返回GeoJSON

ST_AsKML – 返回KML

ST_AsSVG – 返回SVG的Path数据

ST_AsText – 返回WKT

ST_Azimuth – 返回方位角

ST_Buffer – 计算缓冲区

ST_Centroid – 返回几何中心

ST_CoveredBy – 是否被包含

ST_Covers – 是否包含

ST_DWithin – 返回几何体内部要素

ST_Distance – 计算两要素距离

ST_GeogFromText – 根据WKT生成几何数据

ST_GeogFromWKB – 根据WKB生成几何数据

ST_GeographyFromText – 从WKT生成地理几何数据

~= - 判断两几何体是否相同

ST_Intersection – 计算交集

ST_Intersects – 判断是否相交

ST_Length – 计算长度

ST_Perimeter – 计算周长

ST_Project - 返回从起点按距离和方位角(方位角)投影的点。

ST_Segmentize - 返回不超过给定距离的段的修改几何/地理。

ST_Summary – 返回几何体综述。

<-> - 返回A和B的距离。

&& - 返回true如果A和B的包围盒相交。

其他函数:

ST_IsSimple:判断是否是简单几何体

ST_IsValid:判断几何体是否有效。

空间索引

B-tree:为属性数据创建索引,对空间索引无效

GiST(广义搜索树):将数据按不能方面划分,PostGIS使用R-Tree进行索引。

BRIN(区块范围指数):通过汇总表记录空间范围。

SP-GiST(空间划分的广义搜索树):四叉树、kd树、radix(基数)树。

GiST索引

创建2D索引:CREATE INDEX [indexname] ON [tablename] USING GIST ( [geometryfield] );

创建n维索引:CREATE INDEX [indexname] ON [tablename] USING GIST ([geometryfield] gist_geometry_ops_nd);

生产环境并行创建索引:CREATE INDEX CONCURRENTLY [indexname] ON [tablename] USING GIST ( [geometryfield] );

让PostgreSQL收集表格统计信息,优化查询计划:VACUUM ANALYZE [table_name] [(column_name)];

BRIN索引

优点:创建快,占用空间小,性能高

缺点:需要按空间列排序,不会自动更新索引,查询速度比GiST慢

排序方法:CREATE TABLE table_sorted AS

~ SELECT * FROM table  ORDER BY geom;

排序方法2:CREATE INDEX idx_temp_geohash ON table

~ USING btree (ST_GeoHash( ST_Transform( geom, 4326 ), 20));

CLUSTER table USING idx_temp_geohash;

地理数据类型支持BRIN索引:CREATE INDEX [indexname] ON [tablename] USING BRIN ( [geog_col] );

SP-GiST索引

创建索引:CREATE INDEX [indexname] ON [tablename] USING SPGIST ( [geometryfield] );

创建3D索引:CREATE INDEX [indexname] ON [tablename] USING SPGIST ([geometryfield] spgist_geometry_ops_3d);

空间查询

空间关系:ST_Contains, ST_Crosses, ST_Disjoint, ST_Equals, ST_Intersects, ST_Overlaps, ST_Touches, ST_Within

示例:

SELECT city.name, state.name, city.geom

FROM city JOIN state ON ST_Intersects(city.geom, state.geom);

SELECT ST_Relate( 'LINESTRING (1 1, 5 5)',

~ 'POLYGON ((3 3, 3 7, 7 7, 7 3, 3 3))' );

SELECT a.id

FROM roads a, roads b

WHERE a.id != b.id

~ AND a.geom && b.geom

~ AND ST_Relate(a.geom, b.geom, '1*1~*~1~*~*');**

SELECT geom

FROM geom_table

WHERE ST_Distance( geom, 'SRID=312;POINT(100000 200000)' ) < 100

SELECT geom

FROM geom_table

WHERE ST_DWithin( geom, 'SRID=312;POINT(100000 200000)', 100 )

SELECT sum(ST_Length(the_geom))/1000 AS km_roads FROM bc_roads;

SELECT

~ ST_Area(the_geom)/10000 AS hectares

FROM bc_municipality

WHERE name = 'PRINCE GEORGE';