>>> from env_helper import info; info()
页面更新时间: 2023-04-15 20:01:51
运行环境:
    Linux发行版本: Debian GNU/Linux 12 (bookworm)
    操作系统内核: Linux-6.1.0-7-amd64-x86_64-with-glibc2.36
    Python版本: 3.11.2

8.4. 开始使用 SpatiaLite

8.4.1. 熟悉Geometry

在开始之前,先运行:

>>> import sqlite3 as sqlite
>>> conn = sqlite.connect("spalite.db")
>>> conn.enable_load_extension(True)
>>> conn.execute('SELECT load_extension("mod_spatialite.so.7")')
<sqlite3.Cursor at 0x7f62a40a0640>
>>> cursor = conn.cursor()

现在执行第一个查询。

>>> sql = 'SELECT name , AsText(Geom) from pcapital limit 5'
>>> cursor.execute(sql)
---------------------------------------------------------------------------

OperationalError                          Traceback (most recent call last)

Cell In [4], line 2
      1 sql = 'SELECT name , AsText(Geom) from pcapital limit 5'
----> 2 cursor.execute(sql)


OperationalError: no such table: pcapital
>>> for x in cursor:
>>>     print(x)

AsText() 是SpatiaLite 函数, 它返回 Geometry 字段的 WKT 值。 在前面, 我们使用 HEX() 函数返回无法查阅的二进制数据。 现在, AsText() 函数返回有用,且是易于理解的字符串。

点是最简单的 Geometry 类,且它只由一对 [X, Y] 坐标构成。

下面我们使用不同的方式来执行前面的查询。

>>> sql = 'SELECT name , X(Geom), Y(Geom) FROM pcapital limit 5;'
>>> cursor.execute(sql)
>>> for x in cursor:
>>>     print(x)

SpatiaLite的 X() 函数,返回点的 X 坐标。 Y() 函数返回点的 Y 坐标。

使用下面的 Geometry格式转换函数:

>>> sql = "SELECT HEX(GeomFromText('POINT(10 20) '));"
>>> cursor.execute(sql)
>>> for rec in cursor:
>>>     print(rec)
>>> sql = "SELECT HEX(AsBinary(GeomFromText('POINT(10 20) ')));"
>>> cursor.execute(sql)
>>> for rec in cursor:
>>>     print(rec)
>>> sql = "SELECT AsText(GeomFromWKB(X'010100000000000000000024400000000000003440'));"
>>> cursor.execute(sql)
>>> for rec in cursor:
>>>     print(rec)
  • SpatiaLite的函数 GeomFromText() 返回使用内部 BLOB 方式表示的几何图形;

  • 函数 AsBinary() 返回WKB(Well Known Binary)表示的几何图形;

  • 函数 GeomFromWKB() 将WKB的值,转换成对应的内部BLOB的值。

WKB是实现了OpenGIS规范的一种表示方法

GEOMETRY 类

开始之前,先运行:

>>> import sqlite3 as sqlite
>>> conn = sqlite.connect("spalite.db")
>>> conn.enable_load_extension(True)
>>> conn.execute('SELECT load_extension("mod_spatialite.so.7")')
>>> cursor = conn.cursor()

这一部分主要是探讨 OpenGIS 中的定义,SpatiaLite支持的不同 Geometry 类。 简单来讲,任何的 Geometry 类都是一种几何类型。

LINESTRING类型

了解了POINT 类,下面我们来熟悉一下 LINESTRING 类型的 简单的查询。

>>> sql = "SELECT ogc_fid, AsText(Geometry) FROM hyd2_4l limit 5"
>>> cursor.execute(sql)
>>> [print(rec) for rec in cursor]

LINESTRING 是另一个 GEOMETRY 类,并由许多点组成。 在上面你得到一个非常简单的 LINESTRING ,使用四个顶点来表示线。 在实际的GIS数据中,为数以千计的顶点组成的线计数并不简单。

为了更进一步了解,我们来看下面的例子。

>>> sql = '''SELECT ogc_fid, NumPoints(Geometry), GLength(Geometry) ,Dimension(Geometry),
>>>   GeometryType(Geometry) FROM hyd2_4l ORDER BY NumPoints(Geometry) DESC LIMIT 5'''
>>> cursor.execute(sql)
>>> [print(x) for x in cursor]

SpatiaLite 的 NumPoint() 函数返回线几何要素的顶点的数目。 GLength() 函数返回以地图单位计算的线几何类型的几何长度。 Dimension() 函数返回任何一种几何类的维度(对线来讲是1)。 GeometryType() 函数返回任何 Geometry 类型的值。

下面我们来看更多的函数。

>>> sql = '''SELECT ogc_fid, NumPoints(Geometry),
>>>         AsText(StartPoint(Geometry)), Y(PointN(Geometry, 2))
>>>         FROM hyd2_4l ORDER BY NumPoints(Geometry) DESC LIMIT 5'''
>>> cursor = cursor.execute(sql)
>>> [print(x) for x in cursor]

SpatiaLite 的 StartPoint() 函数返回线状几何要素的第一个点。 EndPoint() 函数返回线状几何要素的最后一个点。 PointN() 函数返回选中的顶点和返回点状要素,每个点通过相对索引来标志; 第一个点通过索引值1来标识, 第二个点通过索引值2来标识,以下类推。

你可以通过GEOMETRY类将内部函数的返回值作为参数传递给外部函数,进而查看大量的 SpatiaLite 函数。

多边形类型

>>> sql = 'SELECT name, AsText(Geometry) FROM region_popu WHERE ogc_fid=52'
>>> cursor.execute(sql)
>>> [print(x) for x in cursor]

POLYGON 是另外一个GEOMETRY类。 POLYGON 是一个非常简单的多边形,且只有此部的环(没有内部的洞)。 要注意,POLYGON可以包含何意数据的洞, 通过内部环分隔开来。

外部环是一个简单的 LINESTRING (内部洞也是 LINESTRING )。 注意 POLYGON 是一个闭合的几何类型, POLYGON的第一个点与最后一个点的位置是完全相同的。

>>> sql = '''SELECT Area(Geometry), AsText(Centroid(Geometry)),
>>>              Dimension(Geometry), GeometryType(Geometry) FROM region_popu
>>>              ORDER BY Area(Geometry) DESC LIMIT 5'''
>>> cursor.execute(sql)
>>> [print(x) for x in cursor]

在前面我们使用了 SpatiaLite 的 Dimension() 函数与 GeometryType() 函数。 对于 POLYGON 类型,这两个函数的意义与其他都是一样的。 SpatiaLite的 Area() 函数返回多边形的几何面积, Centroid() 函数返回多边形几何要素的质心。

>>> sql = '''SELECT ogc_fid, NumInteriorRings(Geometry),
>>>         NumPoints(ExteriorRing(Geometry)),
>>>         NumPoints(InteriorRingN(Geometry, 1))
>>>         FROM region_popu ORDER BY NumInteriorRings(Geometry) DESC LIMIT 5'''
>>> cursor.execute(sql)
>>> [print(x) for x in cursor]

SpatiaLite 的 ExteriorRing() 函数返回给定几何要素的外部线环。 任何有效的多边形要素必须要有一个外部线环,并且这个线环是闭合的。 SpatiaLite 的 NumInteriorRings() 函数返回多边形中内部洞的数目可以是一个有效的多边形可以有一些洞,也可以没有。

SpatiaLite 的 InteriorRingN() 函数以 LINESTRING 的格式返回第 N 个内部洞。 每个洞都以相对索引来标识:第一个的索引值是1, 第二个的索引值是2,其余依次类推。

ring 一个 LINESTRING, 所以我们可以使用 NumPoints()函数来获取相关的顶点的数目。 在无效的要素上,则返回结果为NULL

查看多边形的坐标

>>> sql = '''SELECT AsText(InteriorRingN(Geometry, 1)),
>>>     AsText(PointN(InteriorRingN(Geometry, 1), 4)),
>>>     X(PointN(InteriorRingN(Geometry, 1), 5))
>>>     FROM region_popu WHERE ogc_fid=2'''
>>> cursor.execute(sql)
>>> [print(x) for x in cursor]

我们已经在前面遇到过用法了。

对于POLYGON来说,它变得更加乏味,但仍然容易理解。 例如,为了获得最后一列,我们使用了InteriorRingN()来获取第一个内部环, 然后通过PointN()获得第五个顶点。

最后我们可以调用Y()来获取坐标值。

更多的类型

点、线、面是 GEOMETRY 中的基本类。 但是 GEOMETRY 也支持复杂类。

  1. MULTIPOINT 是属于同一个实体的两个或更多点的集合。

  2. MULTILINESTRING 是两个或更多线状要素。

  3. MULTIPOLYGON 是两个或更多多边形要素。

  4. GEOMETRYCOLLECTION 是包含多种要素类型的集合。

对于上面这些类型就不进行更多说明了。总体来讲, 这个类型比基本的类型会多一些性质。

  1. SpatiaLite 的 NumGeometries() 函数返回集合中元素的数目。

  2. GeometryN() 函数返回集合中的第 N 个元素。

  3. GLength() 函数由 MULTILINESTRING 集合中所有线状要素组成的各单独长度的和。

  4. Area() 函数返回 MULTIPOLYGON 集中中所有多边形要素的单独面积的和。

  5. Centroid() 函数返回 MULTIPOLYGON 的平均质心。

GEOMETRY 外包矩形

在开始之前,先运行:

>>> import sqlite3 as sqlite
>>> conn = sqlite.connect("spalite.db")
>>> conn.enable_load_extension(True)
>>> conn.execute('SELECT load_extension("mod_spatialite.so.7")')
>>> cursor = conn.cursor()

下面是所有 GEOMETRY 类的基本属性:

>>> sql = '''SELECT Name, AsText(Envelope(Geometry)) FROM region_popu LIMIT 5'''
>>> cursor.execute(sql)
>>> for x in cursor:
>>>     print(x)

SpatiaLite 的 Envelope() 函数总是返回给定多边形的最小边界矩形。 (Minimum Bounding Rectangle - MBR)。 最小边界矩形由5个点组成(第一个点与最后一个点是相同的)。 需要注意, 在不同的投影方式,同一多边形的最小边界矩形是不一样的(包括形状、大小与方向)。

MBR 也被称为矩形边界框(bounding boxes, BBOX)

有效的点如下所示。

  • POINT #1: minX,minY

  • POINT #2: maxX,minY

  • POINT #3: maxX,maxY

  • POINT #4: minX,maxY

  • POINT #5: minX,minY

MBR 在使用中是比较特别的。 通过 MBR, 可以对多边形的空间关系进行简单与大致的分析。 由于 MBR 计算起来非常快, 所以在提高数据处理速度中得到了广泛的应用。

到目前,通过前面的实例,对于空间数据处理的基本核心有了基本的了解。

8.4.2. 针对数据表的操作

开始运行

>>> import os,shutil,stat
>>> import sqlite3 as sqlite
>>> tmp_db = '/tmp/xx_new_db.sqlite'
>>> if os.path.exists(tmp_db):
>>>     os.remove(tmp_db)
>>> shutil.copy("/gdata/test-2.3.sqlite", tmp_db)
>>> os.chmod(tmp_db, stat.S_IRUSR + stat.S_IWUSR)
>>> conn = sqlite.connect(tmp_db)
>>> conn.enable_load_extension(True)
>>> conn.execute('SELECT load_extension("mod_spatialite.so.7")')
>>> cursor = conn.cursor()

创建表

创建一个新表,并在其中插入一些记录:

>>> cursor.execute('CREATE TABLE MyTable (name TEXT NOT NULL, geom BLOB NOT NULL)')
>>> cursor.execute("INSERT INTO MyTable (name, geom) VALUES ('one', GeomFromText('POINT(1 1)'))")
>>> cursor.execute("INSERT INTO MyTable (name, geom) VALUES ('two', GeomFromText('POINT(2 2)'))")
>>> cursor.execute("INSERT INTO MyTable (name, geom) VALUES ('three', GeomFromText('POINT(3 3)'))")
>>> cursor.execute("SELECT name, AsText(geom) FROM MyTable;")
>>> [print(rec) for rec in cursor]

CREATE TABLE 语句是用来创建一个新表的,创建的时候要指明它要包含的属性列。 属性列的定义可以在后期添加、删除,或修改类型都可以。

还应该了解,一个地理空间属性列的类型是 BLOBINSERT INTO 语句可以在数据表中存储新行。 使用 GeoMFromText() ,你可以创建新的几何对象数据。 最后,使用 SELECT 语句,执行了一个最终的检查。

更新表

更新行与插入或删除它们是一样简单的:

>>> cursor.execute("SELECT name, AsText(geom) FROM mytable")
>>> [print(rec) for rec in cursor]
  • UPDATE SQL语句允许您修改任何列值。

  • 你只需要SET列名称和新的值即可替换当前的。你可以像往常一样使用GeomFromText()函数获取几何值。 PK_UID列是一个特殊的函数,作为这个表的主要键。

  • 在每个表中,PRIMARY KEY列都保证您唯一 一行可能包含一个选定的值,从而确保一致性。

如果你的一张表格,由于任何原因,没有主要的钥匙,别担心; SQLite自动管理表中每一行都有用,每个行的名称都是ROWID; 您可以直接引用SQL表达式,就像普通的列表一样。

还要注意,在这个例子中,你显然没有使用任何交易。(您既不使用BEGIN也不使用COMMIT)

  • 一个新的事务隐式地为每个SQL语句被开始

  • 并且在处理后自动执行隐式COMMIT SQL语句。

  • 这可能很容易导致或多或少严重的性能下降; 到目前为止,SQLite更倾向于您明确的BEGIN和COMMIT 交易,特别是当你执行许多和许多 连续的INSERT和/或UPDATE。

选择数据创建新表格

以下示例需要原始的my_new_db.sqlite数据库。 因此,避免我们执行的DELETE和UPDATE导致任何错位,现在重建它更好。

SQLite的SQL语法为您提供了一种直观的方式来创建一个新的表,并同时使用从另一个表格中选出的数据进行填充。

8.4.3. 管理坐标参考与坐标转换

每个坐标值都是完全限定的域名,需要一个有明确标识。

此值标识了几何相关描述、几何对象定义的坐标空间。

作为一种普遍的规则,只有当它们的坐标是在同样的空间参考系统,不同的几何对象才能进行有意义的互操作。

在GIS中,一种常用的操作,是地理坐标重投影(coordinate reprojection), 是将不同的GIS数据转换成唯一的的空间参考系统,然后进行一些互操作和集成。

当您尝试使用属于不同坐标参考系统的两个GIS数据集时,实体下降得非常远,因为坐标的数值显然在两个不同的空间中。

但是,如果您为一个数据集应用一些时机,则将其放在另一个数据集的相同坐标参考系统中,实体将按预期正确重叠。

欧洲石油调查组[EPSG]维护和分配了一个大的数据集,它描述任何坐标系统和坐标变换都应用于世界各地的GIS数据。

SpatiaLite为任何类型的几何类实现SRID,并支持EPSG数据集来识别坐标参考系。

>>> import sqlite3 as sqlite
>>> db = sqlite.connect(':memory:')
>>> db.enable_load_extension(True)
>>> db.execute('SELECT load_extension("mod_spatialite.so.7")')
>>> cursor = db.cursor()
>>> cursor.execute('BEGIN')
>>> cursor.execute('SELECT InitSpatialMetaData();')

脚本为了初始化Spatial MetaData

开始新的spatialite.exe会话:

>>> cursor = cursor.execute('SELECT * FROM spatial_ref_sys LIMIT 5;')
>>> [print(rec) for rec in cursor]

spatial_ref_sys表是一个空间元数据, 通过执行init  _spatialite-2.3.sql脚本进行初始化 坐标参考系。auth _name和auth _srid分别标识权限 生成这些数据和原始Srid。如果你想了解更多,你可以 阅读以下网站内容:http://en.wikipedia.org/wiki/Geodeticftp://ftp.remotesensing.org/proj/OF90-284.pdf

>>> cursor.close()
>>> db.close()

查看数据库信息

首先我们使用新的数据库,来查看信息:

>>> import os,shutil,stat
>>> import sqlite3 as sqlite
>>> tmp_db = '/tmp/xx_new_db.sqlite'
>>> if os.path.exists(tmp_db):
>>>     os.remove(tmp_db)
>>> shutil.copy("spalite.db", tmp_db)
>>> os.chmod(tmp_db, stat.S_IRUSR + stat.S_IWUSR)
>>> conn = sqlite.connect(tmp_db)
>>> conn.enable_load_extension(True)
>>> conn.execute('SELECT load_extension("mod_spatialite.so.7")')
>>> cursor = conn.cursor()
>>> cursor.execute("select name from sqlite_master where type='table' order by name;")
>>> [print(rec) for rec in cursor]

查看数据库中的空间参考:

>>> cursor = cursor.execute('SELECT * FROM spatial_ref_sys LIMIT 5;')
>>> [print(rec) for rec in cursor]

进一步查看数据库的信息:

>>> cursor = cursor.execute('SELECT DISTINCT Srid(geom) FROM pcapital;')
>>> [print(rec) for rec in cursor]
>>> cursor = cursor.execute('''SELECT DISTINCT SRID(pcapital.geom), spatial_ref_sys.ref_sys_name FROM pcapital,
>>>         spatial_ref_sys WHERE SRID(pcapital.geom) = spatial_ref_sys.srid;''')
>>> [print(rec) for rec in cursor]
  1. SpatiaLite SRID()函数允许您识别标识任何类型的几何srid值。

  2. 你可以简单的运行,然后你会发现srid 32632真正的意思。

创建

您可以通过SQL查询,来完成这个非常复杂的任务。

>>> # cursor.execute('BEGIN')
>>> cursor.execute("SELECT AddGeometryColumn('Towns', 'wgs84', 4326, 'POINT', 2)")
>>> cursor.execute("UPDATE pcapital SET wgs84 = Transform(geom, 4326);")
>>> conn.commit()
>>> cursor.execute('SELECT AsText(geometry), Srid(geometry),AsText(wgs84), Srid(wgs84) FROM pcapital LIMIT 5;')
>>> [print(rec) for rec in cursor]
  1. WGS 84的srid是4326。

  2. 根据需要应用SpatiaLite Transform()函数,从原始函数获取一个新的几何体。

  3. NewTowns表中有两个替代几何体:

    • 原始的geom列包含32623中的几何。

    • UTM区32N坐标参考系。

    • 新的wgs84列包含4326中的几何。

    • WGS84坐标参考系。

    • 您现在可以根据需要和适当的方式随意使用一个或另一个。