Postgis教程
一、PostgreSQL与PostGIS的关系
PostgreSQL 是世界上技术最先进的开源数据库,其前身是1977年一个源于Berkeley名为Ingres的非关系型数据库,其项目领导人为Michael Stonebraker教授。1982年该教授商业化了Ingres;1985年,Michael Sto nebraker教授回到Berkeley,开始对新的数据库设计进行研究,并于次年在美国防务高级研究项目局(DARPA)、陆军研究办公室(ARO)、国家科学基金(NSF)以及ESL, I nc等机构的赞助下启动了Postgres (Post-Ingres)项目。
Postgres 在1987年形成第一个Demo,1989年发布第一个版本,直到1993年的4.2版本,由于外部用户过多,做技术支持和维护源代码的时间影响到了对数据库的研究,因此Berkeley中止了该项目。在此期间,Postgres项目就已经被使用在了一些GIS系统中。
Postgres项目并未就此消亡,在1994年两个Berkeley的研究生向Postgres中加入了SQL语言解释器,将之改名为Postgre95并发布到了互联网上。经过一些黑客的修改,1996年Postgres95再次更名为PostgreSQL,并采用BSD许可证发布了第一个开源版本。经过多年发展,PostgreSQL已经发展成为一个技术非常先进的开源数据库,其支持特性之多性能之强可与诸多高级商业数据库比肩。
这里提供了一个世界上主流数据库的特性比较,有兴趣的可以从中管窥PostgreSQL在数据库领域中的成就。
http://en.wikipedia.o rg/wiki/Co mparison_of_relatio nal_database_management_systems
PostGIS 则是PostgreSQL的一个扩展,目的是使PostgreSQL支持空间数据的存储和使用,其本质类似于ArcSDE和Oracle Spatial Extensio n。PostGIS是采用GPL许可发布的,完整地实现了OGC的《Simple Features Specificatio n for SQL》规范,并于2006年获得OGC认证。在此基础上,PostGIS还对规范进行了一些扩展,在后面的特性中我们可以慢慢了解到。
二、PostGIS中的几何类型
PostGIS支持所有OGC规范的“Simple Features”类型,同时在此基础上扩展了对3DZ、3DM、4D坐标的
支持。
1. OGC的WKB和WKT格式
OGC定义了两种描述几何对象的格式,分别是WKB(Well-Known Binary)和WKT(Well-Known Text)。
在SQL语句中,用以下的方式可以使用WKT格式定义几何对象:
POINT(0 0) ——点
LINESTRING(0 0,1 1,1 2) ——线
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) ——多点
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))) ——几何集合
以下语句可以使用WKT格式插入一个点要素到一个表中,其中用到的Geo mFromText等函数在后面会有详细介绍:
INSERT INTO table ( SHAPE, NAME )
VALUES ( Geo mFro mText('POINT(116.39 39.9)', 4326), '北京');
2. EWKT、EWKB和Canonical格式
EWKT和EWKB相比OGC WKT和WKB格式主要的扩展有3DZ、3DM、4D坐标和内嵌空间参考支持。
以下以EWKT语句定义了一些几何对象:
POINT(0 0 0) ——3D点
SRID=32632;POINT(0 0) ——内嵌空间参考的点
POINTM(0 0 0) ——带M值的点
POINT(0 0 0 0) ——带M值的3D点
SRID=4326;MULTIPOINTM(0 0 0,1 2 1) ——内嵌空间参考的带M值的多点
以下语句可以使用EWKT格式插入一个点要素到一个表中:
INSERT INTO table ( SHAPE, NAME )
VALUES ( Geo mFro mEWKT('SRID=4326;POINTM(116.39 39.9 10)'), '北京' )
Cano nical格式是16进制编码的几何对象,直接用SQL语句查询出来的就是这种格式。
3. SQL-MM格式
SQL-MM格式定义了一些插值曲线,这些插值曲线和EWKT有点类似,也支持3DZ、3DM、4D坐标,但是不支持嵌入空间参考。
以下以SQL-MM语句定义了一些插值几何对象:
CIRCULARSTRING(0 0, 1 1, 1 0) ——插值圆弧
COMPOUNDCURVE(CIRCULARSTRING(0 0, 1 1, 1 0),(1 0, 0 1)) ——插值复合曲线CURVEPOLYGON(CIRCULARSTRING(0 0, 4 0, 4 4, 0 4, 0 0),(1 1, 3 3, 3 1, 1 1)) ——曲线多边形MULTICURVE((0 0, 5 5),CIRCULARSTRING(4 0, 4 4, 8 4)) ——多曲线
MULTISURFACE(CURVEPOLYGON(CIRCULARSTRING(0 0, 4 0, 4 4, 0 4, 0 0),(1 1, 3 3, 3 1, 1 1)),((10 10, 14 12, 11 10, 10 10),(11 11, 11.5 11, 11 11.5, 11 11))) ——多曲面
三、PostGIS中空间信息处理的实现
1. spatial_ref_sys表
在基于PostGIS模板创建的数据库的public模式下,有一个spatial_ref_sys表,它存放的是OGC规范的空间参考。我们取我们最熟悉的4326参考看一下:
它的srid存放的就是空间参考的Well-Kno wn ID,对这个空间参考的定义主要包括两个字段,srtext存放的是以字符串描述的空间参考,proj4text存放的则是以字符串描述的PROJ.4 投影定义(PostGIS使用PROJ.4实现投影)。
4326空间参考的srtext内容:
GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS
84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],TOWGS84[0,0,0,0,0,0,0],AUTHORITY["E PSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.01745329251 994328,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]]
4326空间参考的proj4text内容:
+pro j=lo nglat +ellps=WGS84 +datum=WGS84 +no_defs
2. geometry_columns表
geo metry_columns表存放了当前数据库中所有几何字段的信息,比如我当前的库里面有两个空间表,在geo metry_columns表中就可以找到这两个空间表中几何字段的定义:
其中f_table_schema字段表示的是空间表所在的模式,f_table_name字段表示的是空间表的表名,
f_geo metry_co lumn字段表示的是该空间表中几何字段的名称,srid字段表示的是该空间表的空间参考。
3. 在PostGIS中创建一个空间表
在PostGIS中创建一个包含几何字段的空间表分为2步:第一步创建一个一般表,第二步给这个表添加几何字段。
以下先在test模式下创建一个名为cities的一般表:
create table test.cities (id int4, name varchar(20))
再给cities添加一个名为shape的几何字段(二维点):
select AddGeo metryColumn('t est', 'cities', 'shape', 4326, 'POINT', 2)
4. PostGIS对几何信息的检查
PostGIS可以检查几何信息的正确性,这主要是通过IsValid函数实现的。
以下语句分辨检查了2个几何对象的正确性,显然,(0, 0)点和(1,1)点可以构成一条线,但是(0, 0)点和(0, 0)点则不能构成,这个语句执行以后的得出的结果是TRUE,FALSE。
select IsValid('LINESTRING(0 0, 1 1)'), IsValid('LINESTRING(0 0,0 0)')
默认PostGIS并不会使用IsValid函数检查用户插入的新数据,因为这会消耗较多的CPU资源(特别是复杂的几何对象)。当你需要使用这个功能的时候,你可以使用以下语句为表新建一个约束:
ALTER TABLE cities
ADD CONSTRAINT geo metry_valid
CHECK (IsValid(shape))
这时当我们往这个表试图插入一个错误的空间对象的时候,会得到一个错误:
INSERT INTO test.cities ( shape, name )
VALUES ( Geo mFro mText('LINESTRING(0 0,0 0)', 4326), '北京');
ERROR: new row for relatio n "cities" v iolates check constraint "geo metry_valid"
SQL 状态: 23514
5. PostGIS中的空间索引
数据库对多维数据的存取有两种索引方案,R-Tree和GiST(Generalized Search Tree),在PostgreSQL 中的GiST比R-Tree的健壮性更好,因此PostGIS对空间数据的索引一般采用GiST实现。
以下的语句给sde模式中的cities表添加了一个空间索引shape_index_cities,在pgAdmin中也可以通过图形界面完成相同的功能。
CREATE INDEX shape_index_cities
ON sde.cities
USING gist
(shape);
另外要注意的是,空间索引只有在进行基于边界范围的查询时才起作用,比如“&&”操作。
四、PostGIS中的常用函数
以下内容包括比较多的尖括号,发布到blogger的时候会显示不正常,内容太多我也无暇一个个手动改代码,因此如有问题就去参考PostGIS官方文档。
首先需要说明一下,这里许多函数是以ST_[X]yyy形式命名的,事实上很多函数也可以通过xyyy的形式访问,在PostGIS的函数库中我们可以看到这两种函数定义完全一样。
1. OGC标准函数
管理函数:
添加几何字段AddGeo metryColumn(, , , , , )
删除几何字段DropGeo metryCo lumn(, , )
检查数据库几何字段并在geo metry_columns中归档Probe_Geo metry_Columns()
给几何对象设置空间参考(在通过一个范围做空间查询时常用)ST_SetSRID(geo metry, integer)
几何对象关系函数:
获取两个几何对象间的距离ST_Distance(geo metry, geometry)
如果两个几何对象间距离在给定值范围内,则返回TRUE ST_DWithin(geo metry, geo metry, float)
判断两个几何对象是否相等
(比如LINESTRING(0 0, 2 2)和LINESTRING(0 0, 1 1, 2 2)是相同的几何对象)ST_Equals(geo metry, geo metry)
判断两个几何对象是否分离ST_Disjoint(geo metry, geo metry)
判断两个几何对象是否相交ST_Intersects(geo metry, geo metry)
判断两个几何对象的边缘是否接触ST_Touches(geo metry, geo metry)
判断两个几何对象是否互相穿过ST_Crosses(geo metry, geo metry)
判断A是否被B包含ST_Within(geo metry A, geo metry B)
判断两个几何对象是否是重叠ST_Overlaps(geo metry, geo metry)
判断A是否包含B ST_Contains(geo metry A, geometry B)
判断A是否覆盖 B ST_Co v ers(geo metry A, geo metry B)
判断A是否被B所覆盖ST_Co veredBy(geo metry A, geo metry B)
通过DE-9IM 矩阵判断两个几何对象的关系是否成立ST_Relate(geo metry, geo metry,
intersectio nPatternMatrix)
获得两个几何对象的关系(DE-9IM矩阵)ST_Relate(geo metry, geo metry)
几何对象处理函数:
获取几何对象的中心ST_Centroid(geo metry)
面积量测ST_Area(geo metry)
长度量测ST_Length(geo metry)
返回曲面上的一个点ST_PointOnSurface(geo metry)
获取边界ST_Bo undary(geo metry)
获取缓冲后的几何对象ST_Buffer(geo metry, do uble, [integer])
获取多几何对象的外接对象ST_Co nvexHull(geo metry)
获取两个几何对象相交的部分ST_Intersectio n(geo metry, geo metry)
将经度小于0的值加360使所有经度值在0-360间ST_Shift_Longitude(geo metry)
获取两个几何对象不相交的部分(A、B可互换)ST_SymDifference(geo metry A, geo metry B) 从A去除和B相交的部分后返回ST_Difference(geometry A, geo metry B)
返回两个几何对象的合并结果ST_Unio n(geo metry, geo metry)
返回一系列几何对象的合并结果ST_Unio n(geo metry set)
用较少的内存和较长的时间完成合并操作,结果和ST_Unio n相同ST_MemUnion(geo metry set)
几何对象存取函数:
获取几何对象的WKT描述ST_AsText(geo metry)
获取几何对象的WKB描述ST_AsBinary(geo metry)
获取几何对象的空间参考ID ST_SRID(geo metry)
获取几何对象的维数ST_Dimensio n(geo metry)
获取几何对象的边界范围ST_Envelope(geo metry)
判断几何对象是否为空ST_IsEmpty(geo metry)
判断几何对象是否不包含特殊点(比如自相交)ST_IsSimple(geo metry)
判断几何对象是否闭合ST_IsClosed(geo metry)
判断曲线是否闭合并且不包含特殊点ST_IsRing(geo metry)
获取多几何对象中的对象个数ST_NumGeo metries(geo metry)
获取多几何对象中第N个对象ST_Geo metryN(geo metry,int)
获取几何对象中的点个数ST_NumPoints(geo metry)
获取几何对象的第N个点ST_PointN(geo metry,integer)
获取多边形的外边缘ST_ExteriorRing(geo metry)
获取多边形内边界个数ST_NumI nterio rRings(geo metry)
同上ST_NumI nterio rRing(geo metry)
获取多边形的第N个内边界ST_InteriorRingN(geometry,integer)
获取线的终点ST_EndPoint(geo metry)
获取线的起始点ST_StartPoint(geo metry)
获取几何对象的类型GeometryType(geo metry)
类似上,但是不检查M值,即POINTM对象会被判断为point ST_Geo metryType(geo metry)
获取点的X坐标ST_X(geo metry)
获取点的Y坐标ST_Y(geo metry)
获取点的Z坐标ST_Z(geo metry)
获取点的M值ST_M(geo metry)
几何对象构造函数:
参考语义:
Text:WKT
WKB:WKB
Geom:Geo metry
M:Multi
Bd:BuildArea
Coll:Co llectio n ST_Geo mFro mText(text,[])
ST_PointFro mText(text,[])
ST_LineFro mText(text,[])
ST_LinestringFro mText(text,[])
ST_PolyFro mText(text,[])
ST_Polygo nFro mText(text,[])
ST_MPointFro mText(text,[])
ST_MLineFro mText(text,[])
ST_MPolyFro mText(text,[])
ST_Geo mCollFro mText(text,[])
ST_Geo mFro mWKB(bytea,[])
ST_Geo metryFro mWKB(bytea,[])
ST_PointFro mWKB(bytea,[])
ST_LineFro mWKB(bytea,[])
ST_LinestringFro mWKB(bytea,[])
ST_PolyFro mWKB(bytea,[])
ST_Polygo nFro mWKB(bytea,[])
ST_MPointFro mWKB(bytea,[])
ST_MLineFromWKB(bytea,[])
ST_MPolyFro mWKB(bytea,[])
ST_Geo mCollFro mWKB(bytea,[])
ST_BdPolyFro mText(text WKT, integer SRID)
ST_BdMPolyFro mText(text WKT, integer SRID)
2. PostGIS扩展函数
管理函数:
删除一个空间表(包括geo metry_columns中的记录)DropGeo metryTable([], ) 更新空间表的空间参考UpdateGeo m etrySRID([], , , )
更新空间表的统计信息update_geo metry_stats([, ])
参考语义:
Geos:GEOS库
Jts:JTS库
Proj:PROJ4库postgis_versio n()
postgis_lib_version()
postgis_lib_build_date()
postgis_script_build_date()
postgis_scripts_installed()
postgis_scripts_released()
postgis_geos_versio n()
postgis_jts_version()
postgis_proj_version()
postgis_uses_stats()
postgis_full_version()
几何操作符:
A范围=B范围A = B
A范围覆盖B范围或A范围在B范围左侧A &<> B
A范围在B范围左侧A <<>> B
A范围覆盖B范围或A范围在B范围下方A &<| B A范围覆盖B范围或A范围在B范围上方A |&> B A范围在B范围下方A <<| B A范围在B范围上方A |>> B
A=B A ~= B
A范围被B范围包含A @ B
A范围包含B范围A ~ B
A范围覆盖B范围A && B
几何量测函数:
量测面积ST_Area(geo metry)
根据经纬度点计算在地球曲面上的距离,单位米,地球半径取值6370986米ST_distance_sphere(point, point)
类似上,使用指定的地球椭球参数ST_distance_spheroid(point, point, spheroid)
量测2D对象长度ST_length2d(geo metry)
量测3D对象长度ST_length3d(geo metry)
根据经纬度对象计算在地球曲面上的长度ST_length_spheroid(geo metry,spheroid)
ST_length3d_sphero id(geo metry,spheroid)
量测两个对象间距离ST_distance(geo metry, geo metry)
量测两条线之间的最大距离ST_max_distance(linestring,linestring) 量测2D对象的周长ST_perimeter(geo metry)
ST_perimeter2d(geo metry)
量测3D对象的周长ST_perimeter3d(geo metry)
量测两点构成的方位角,单位弧度ST_azimuth(geo metry, geo metry)
几何对象输出:
参考语义:
NDR:Little Endian
XDR:big-endian
HEXEWKB:Cano nical
SVG:SVG 格式
GML:GML 格式
KML:KML 格式
GeoJson:GeoJson 格式
ST_AsBinary(geo metry,{'N DR'|'XDR'})
ST_AsEWKT(geo metry)
ST_AsEWKB(geo metry, {'NDR'|'XDR'})
ST_AsHEXEWKB(geo metry, {'N DR'|'XDR'})
ST_AsSVG(geo metry, [rel], [precision])
ST_AsGML([versio n], geo metry, [precision])
ST_AsKML([version], geo metry, [precision])
ST_AsGeo Json([version], geo metry, [precision], [options])
几何对象创建:
参考语义:
Dump:转储ST_Geo mFro mEWKT(text)
ST_Geo mFro mEWKB(b ytea)
ST_MakePoint(, , [], [])
ST_MakePointM(, , )
ST_MakeBox2D(, )
ST_MakeBox3D(, )
ST_MakeLine(geo metry set)
ST_MakeLine(geo metry, geo metry)
ST_LineFro mMultiPoint(multipoint)
ST_MakePolygo n(linestring, [linestring[]])
ST_BuildArea(geo metry)
ST_Polygo nize(geo metry set)
ST_Co llect(geo metry set)
ST_Co llect(geo metry, geo metry)
ST_Dump(geo metry)
ST_DumpRings(geo metry)
几何对象编辑:
给几何对象添加一个边界,会使查询速度加快ST_AddBBOX(geo metry) 删除几何对象的边界ST_DropBBOX(geo metry)
添加、删除、设置点ST_AddPoint(linestring, point, [])
ST_Remo vePoint(linestring, offset)
ST_SetPoint(linestring, N, point)
几何对象类型转换ST_Force_collection(geo metry)
ST_Fo rce_2d(geometry)
ST_Fo rce_3dz(geo metry), ST_Force_3d(geo metry),
ST_Fo rce_3dm(geo metry)
ST_Fo rce_4d(geo metry)
ST_Multi(geo metry)
将几何对象转化到指定空间参考ST_Transform(geo metry,integer)
对3D几何对象作仿射变化ST_Affine(geo metry, float8, float8, flo at8, float8, float8, flo at8, flo at8, float8, float8, flo at8, flo at8, float8)
对2D几何对象作仿射变化ST_Affine(geo metry, flo at8, flo at8, float8, flo at8, flo at8, float8)
对几何对象作偏移ST_Translate(geometry, flo at8, float8, float8)
对几何对象作缩放ST_Scale(geo metry, float8, float8, flo at8)
对3D几何对象作旋转ST_RotateZ(geo metry, flo at8)
ST_RotateX(geo metry, float8)
ST_RotateY(geo metry, float8)
对2D对象作偏移和缩放ST_TransScale(geo metry, flo at8, flo at8, float8, float8)
反转ST_Reverse(geo metry)
转化到右手定则ST_ForceRHR(geo metry)
参考IsSimple函数
使用Douglas-Peuker算法ST_Simplify(geo metry, tolerance)
ST_Simplify PreserveTopology(geo metry, tolerance)
讲几何对象顶点捕捉到网格ST_SnapToGrid(geo metry, o riginX, originY, sizeX, sizeY)
ST_SnapToGrid(geo metry, sizeX, sizeY), ST_SnapToGrid(geo metry, size)
第二个参数为点,指定原点坐标ST_SnapToGrid(geo metry, geo metry, sizeX, sizeY, sizeZ, sizeM)
分段ST_Segmentize(geo metry, maxlength)
合并为线ST_LineMerge(geo metry)
线性参考:
根据location(0-1)获得该位置的点ST_line_interpo late_point(linestring, locatio n)
获取一段线ST_line_substring(linestring, start, end)
根据点获取location(0-1)ST_line_locate_point(LineString, Point)
根据量测值获得几何对象ST_locate_alo ng_measure(geo metry, flo at8)
根据量测值区间获得几何对象集合ST_locate_between_measures(geo metry, float8, float8)
杂项功能函数:
几何对象的摘要ST_Summary(geo metry)
几何对象的边界ST_bo x2d(geo metry)
ST_bo x3d(geo metry)
多个几何对象的边界ST_extent(geo metry set)
0=2d, 1=3dm, 2=3dz, 3=4d ST_zmflag(geo metry)
是否包含Bounding Bo x ST_HasBBOX(geo metry)
几何对象的维数:2、3、4 ST_ndims(geo metry)
子对象的个数ST_nrings(geo metry)
ST_npoints(geo metry)
对象是否验证成功ST_isvalid(geo metry)
扩大几何对象ST_expand(geo metry, float)
计算一个空间表的边界范围ST_estimated_extent([schema], table, geocolumn)
获得空间参考ST_find_srid(, , )
几何对象使用的内存大小,单位byte ST_mem_size(geo metry)
点是否在圆上ST_point_inside_circle(,,,)
获取边界的X、Y、Z ST_XMin(box3d)
ST_YMin(box3d)
ST_ZMin(box3d)
ST_XMax(bo x3d)
ST_YMax(bo x3d)
ST_ZMax(bo x3d)
构造一个几何对象的数组ST_Accum(geo metry set)
长事务支持:
启用/关闭长事务支持,重复调用无副作用EnableLo ngTransactions() DisableLongTransactions()
检查对行的update和delete操作是否已授权CheckAuth([],
下面我们希望PostGIS帮助我们实现一些空间分析。我们以污染发生的城市为起点,当地风向为主方向,构造一个30度开角的范围;这个范围将是污染扩散的主要方向,扩散的范围主要和风的强度有关;在构造这个区域以后,为了保险起见,我们在对其进行一定范围的缓冲,最后得到每个污染源可能扩散的范围。我们构造的SQL语句如下:
select *,ST_AsGeoJso n(ST_Buffer(ST_Po lygo nFro mText('POLY GO N(('||ST_X(shape)||' '||ST_Y(shape)||','
||ST_X(shape)+velocity*co s((directio n+15)*PI()/180)/20||' '||ST_Y(shape)+velo city*sin((directio n+15)*PI()/180)/20||','
||ST_X(shape)+velocity*co s((directio n-15)*PI()/180)/20||' '||ST_Y(shape)+velocity*sin((directio n-15)*PI()/180)/20||','
||ST_X(shape)||' '||ST_Y(shape)||'))'), velocity/50)) fro m sde.wind