使用 SpatiaLite 构建位置到时区 API
SQLite 的 SpatiaLite 扩展增加了大量用于地理空间分析的函数,可以与 Datasette 结合使用来构建 GIS(地理信息系统)应用程序。
本教程将展示如何结合使用 SpatiaLite 和 Datasette 来创建一个 JSON API,该 API 可以返回地球上特定经度和纬度点的时区。
我们将构建什么
您可以在此处试用该 API。提供经度和纬度,它将返回相应的时区 ID:https://timezones.datasette.io/timezones/by_point
一些示例
- 英国布莱顿 位于 Europe/London 时区(JSON 格式)
- 美国旧金山 位于 America/Los_Angeles 时区(JSON 格式)
- 日本东京 位于 Asia/Tokyo 时区(JSON 格式)
设置开发环境
您需要为本教程准备好两件事
- 系统上已安装 SpatiaLite
- 允许
sqlite3
模块加载额外扩展的 Python 安装
推荐:使用 GitHub Codespaces
GitHub Codespaces 可以为您提供一个免费的开发环境,通过您的网络浏览器即可访问,并且预装了您所需的所有工具。
这是完成本教程最简单的方法。
访问此链接 创建一个新的 Codespace,其中包含了本教程其余部分所需的一切。
使用 Mac
在 Mac 上,您可以使用 Homebrew 安装 SpatiaLite 和 Datasette,如下所示
brew install spatialite-tools datasette
如果您使用的是 Mac,可能会发现您的 Python 安装无法加载外部 SQLite 模块。您可以通过运行以下命令进行检查
datasette --load-extension spatialite
如果您收到关于 enable_load_extension
的错误消息,请查阅 此页面 获取解决问题的提示。
构建数据库
要构建此项目,我们首先需要全球所有时区的地理形状数据。
timezone-boundary-builder 是 Evan Siroky 的一个项目,它使用来自 OpenStreetMap 的数据创建详细的时区多边形,然后以 GeoJSON 和 Shapefile 格式发布。结果数据根据 开放数据共享开放数据库许可证 (ODbL) 提供。
首先从最新版本下载 timezones-with-oceans.shapefile.zip
文件
wget https://github.com/evansiroky/timezone-boundary-builder/releases/download/2022g/timezones-with-oceans.shapefile.zip
(如果您正在使用 Codespaces,应在 Codespaces 界面的“终端”选项卡中运行所有这些命令。)
一个 Shapefile 是一组描述地理要素集合的文件。无需解压此 zip 文件——我们的工具可以直接使用它。
我们可以使用 shapefile-to-sqlite
命令将文件加载到新的 SpatiaLite 数据库中。
如果您在 Codespaces 中运行本教程,该工具已安装。否则,您可以按如下方式安装它
pip install shapefile-to-sqlite
要将 Shapefile 加载到数据库中,请使用以下命令
shapefile-to-sqlite timezones.db \
timezones-with-oceans.shapefile.zip \
--table timezones \
--spatial-index
这将创建一个名为 timezones.db
的新数据库文件,并将 Shapefile 中的形状加载到一个名为 timezones
的新表中。它还将在 timezones
表上设置一个空间索引,稍后会对此进行描述。
该命令将显示一个进度条,如下所示
zip://timezones-with-oceans.shapefile.zip
[########----------------------------] 22% 00:00:38
这可能需要几分钟才能完成。
在 Datasette 中浏览数据库
我们现在可以针对该文件启动 Datasette 来浏览数据
datasette timezones.db --load-extension spatialite
您应该会看到以下输出
INFO: Started server process [5385]
INFO: Waiting for application startup.
INFO: Application startup complete.
INFO: Uvicorn running on http://127.0.0.1:8001 (Press CTRL+C to quit)
如果您使用的是自己的计算机,现在可以访问 http://127.0.0.1:8001 查看 Datasette 界面。
如果您在 Codespaces 中运行,该工具应该会提供一个“在浏览器中打开”按钮。
您现在可以浏览数据库了。应该有一个 timezones
表,包含三列:id
、tzid
和 geometry
不过目前没什么可看的!geometry 列只是一些大的二进制数据块。
在地图上查看时区形状
要在地图上查看时区几何图形,我们可以安装一个 Datasette 插件。
Chris Amico 开发的 datasette-geojson-map 插件增加了直接在地图上渲染 GeoJSON 和 SpatiaLite 几何图形的功能。
在终端中按 Ctrl+C 停止 Datasette 服务器,然后运行以下命令
datasette install datasette-geojson-map
现在再次启动 Datasette
datasette timezones.db \
--load-extension spatialite \
--setting default_page_size 10
我们在这里添加了一个额外的设置,将默认分页大小设置为 10。这是因为有些时区多边形非常大,默认分页大小 100 可能需要很长时间才能渲染。
再次访问 timezones
表,您应该会看到类似这样的内容
查找点的时区
现在我们有了包含时区几何图形的完整表,我们可以构建一个 SQL 查询来显示地球上任何特定点的时区。
我们可以使用 SpatiaLite 的 within()
函数来实现这一点,该函数接受两个几何图形并检查一个是否包含在另一个之中。
首先,我们需要一个表示特定经纬度点的几何图形。
41.798, -87.696
是芝加哥市内的一个点 - 在 Google 地图上查看。
我们可以使用 SpatiaLite 的 MakePoint(longitude, latitude)
函数(请注意这里的经度和纬度是反过来的)来创建一个几何图形,然后可以将其与 within()
函数一起使用
select
tzid
from
timezones
where
within(
MakePoint(-87.696, 41.798),
timezones.Geometry
) = 1
如果第一个几何图形包含在第二个几何图形中,within(geom1, geom2)
返回 1
,否则返回 0
。
果然,这个查询返回了 America/Chicago
。
在 SQL 查询中使用参数
Datasette 有一个功能,SQL 查询可以包含命名参数 :like_this
,这些参数将被转换为表单字段并用于向查询提供新值。
尝试使用以下查询
select
tzid
from
timezones
where
within(
MakePoint(:longitude, :latitude),
timezones.Geometry
) = 1
这看起来应该可行...但如果您使用先前的坐标进行尝试,会发现它没有返回任何结果。
这是因为 MakePoint()
需要接收浮点数值,但 Datasette 将所有参数作为字符串传递。
如果您传递了无效类型,MakePoint()
返回 null
- 那么 within(null, geometry)
返回 -1
表示结果无效。
我们可以通过将字符串转换为浮点数来解决这个问题,如下所示
select
tzid
from
timezones
where
within(
MakePoint(cast(:longitude as float), cast(:latitude as float)),
timezones.Geometry
) = 1
这个查询达到了预期的效果:给定地球上任何一个点的经纬度,它将返回正确的时区。
使用索引加速
这个查询还有一个遗留问题:它相对较慢。在我的测试中,查询运行时间在 150ms 到 750ms 之间,这是因为需要将该点与数据库中的所有 455 个多边形进行比较。
我们第一次导入数据时使用了 --spatial-index
选项。以下是如何利用这个空间索引来加速的方法
select tzid
from
timezones
where
within(
MakePoint(cast(:longitude as float), cast(:latitude as float)),
timezones.Geometry
) = 1
and rowid in (
select
rowid
from
SpatialIndex
where
f_table_name = 'timezones'
and search_frame = MakePoint(cast(:longitude as float), cast(:latitude as float))
)
这里的技巧是额外的 and rowid in (...)
子查询。
SpatialIndex
是一个特殊的表 - 它是 SpatiaLite 创建的一个虚拟表。
您可以通过传递与其关联空间索引的另一个表的名称(在本例中是 timezones
)来查询该表。然后,您将一个几何图形作为 search_frame
传入,该表将返回一个 rowid
值列表,代表与您传入的几何图形的粗略边界框重叠的任何多边形。
请注意,这并非精确比较:您获得的一些行 ID 可能不会与几何图形精确相交。
但这已经足够近似了。您可以将这些值与更精确的 within()
检查结合起来,这样就只需对多边形总集的一个子集运行完整计算。
在我的本地测试中,这将查询所需的时间从 150ms 降低到不到 8ms - 这是显著的加速!
添加国家多边形
哪些时区与特定国家相关?
目前我们的数据无法告诉我们这一点。我们可以只过滤以 America/
开头的时区,但这将涵盖北美和南美的所有内容。如果我们能按国家浏览时区,那岂不是很棒?
我们可以尝试通过加载世界各国的多边形来回答这个问题。
datahub.io/core/geo-countries 包含我们所需的数据,该数据源自 Natural Earth 并根据 PDDL 许可证发布。
我们可以像这样下载一个包含国家数据的 GeoJSON 文件
wget https://datahub.io/core/geo-countries/r/countries.geojson
现在我们可以使用 geojson-to-sqlite 将其加载到 countries
数据库表中,并同样创建一个空间索引
pip install geojson-to-sqlite
geojson-to-sqlite timezones.db countries \
countries.geojson --spatial-index
Datasette 应该能识别新的 countries
表,并且 datasette-geojson-map
现在会显示这些国家的渲染轮廓。
要查找与特定国家相交的时区,我们可以使用以下查询
select
tzid,
geometry
from
timezones
where intersects(
timezones.geometry, (
select geometry from countries where admin = 'United Kingdom'
) = 1
)
我们在这里使用了一些额外的技巧。
intersects()
函数是 SpatiaLite 的一个函数,用于检查两个几何图形是否以任何方式相互相交。
我们还使用了子查询来访问我们感兴趣的国家的几何图形
select geometry from countries where admin = 'United Kingdom'
countries
表中的 admin
列包含国家名称 - 这里我们获取的是英国的几何图形。
以下是此查询的结果。
简化多边形
尝试将国家名称更改为“United States of America”,您可能会遇到一个问题:与美国相交的时区几何图形非常大,可能无法在地图上渲染!
我们可以使用另一个 SpatiaLite 函数来帮助解决这个问题:simplify()
。此函数应用 Douglas–Peucker 算法 将多边形简化为较少数量的点,这使得渲染变得容易得多。
simplify()
接受一个几何图形和一个精度值。经过反复试验,我发现精度值设置为 0.05
对于这些时区多边形效果很好
select
tzid,
simplify(geometry, 0.05) as geometry
from
timezones
where intersects(
timezones.geometry, (
select geometry from countries where admin = 'United States of America'
) = 1
)
这里必须使用 simplify(...) as geometry
,因为 datasette-geojson-map
插件会寻找一个名为 geometry
的输出列以便在地图上渲染。
以下是此查询针对美国的结果。
使用索引加速
我们可以使用与上面列出的时区查询类似的方式来使用空间索引
with country as (
select
geometry
from
countries
where
admin = :country
)
select
timezones.tzid,
simplify(timezones.geometry, 0.05) as geometry
from
timezones,
country
where
timezones.id in (
select
SpatialIndex.rowid
from
SpatialIndex,
country
where
f_table_name = 'timezones'
and search_frame = country.geometry
)
and intersects(timezones.geometry, country.geometry) = 1
这使用了与之前相同的技巧 - 一个 where timezones.id in (...)
子查询,它从 SpatialIndex
虚拟表中返回一个可能的 rowid
值列表。
如果一个表有一个整数主键 - 例如我们这里的 timezones
表 - SQLite 会将 rowid
设置为相同的值,这就是为什么我们可以在这个查询中比较 timezones.id
和 rowid
。
我们在这里使用了另一个技巧。为了避免两次运行查询来选择国家几何图形,我们转而将其捆绑到一个 CTE(公用表表达式)中,放在查询的开头。
with country as (...)
部分使得该查询的结果在 SQL 查询期间可以作为名为 country
的临时表使用。
使用预设查询定义元数据
现在我们已经弄清楚了驱动 API 所需的所有查询,是时候将它们捆绑到一个可以部署到互联网的配置中了。
Datasette 的元数据系统可用于提供有关数据库的额外信息,也可用于配置预设查询 - 带有名称的 SQL 查询,例如本教程开头展示的 by_point
示例。
这是一个 metadata.yml
文件,其中定义了 by_point
查询以及一个按国家代码查找时区的查询。
title: Time zones API
description: |
An API for looking up time zones by latitude/longitude
about: simonw/timezones-api
about_url: https://github.com/simonw/timezones-api
license: ODbL
license_url: http://opendatacommons.org/licenses/odbl/
source: timezone-boundary-builder
source_url: https://github.com/evansiroky/timezone-boundary-builder
allow_sql: false
databases:
timezones:
tables:
countries:
source: Natural Earth
source_url: https://www.naturalearthdata.com/
license: Open Data Commons Public Domain Dedication and License (PDDL) v1.0
license_url: https://opendatacommons.org/licenses/pddl/1-0/
about: geo-countries
about_url: https://datahub.io/core/geo-countries
queries:
by_point:
title: Find time zone by lat/lon
sql: |
select tzid
from
timezones
where
within(
MakePoint(cast(:longitude as float), cast(:latitude as float)),
timezones.Geometry
) = 1
and rowid in (
select
rowid
from
SpatialIndex
where
f_table_name = 'timezones'
and search_frame = MakePoint(cast(:longitude as float), cast(:latitude as float))
)
by_country:
title: Find time zones that intersect a country
sql: |
with country as (
select
geometry
from
countries
where
admin = :country
)
select
timezones.tzid,
simplify(timezones.geometry, 0.05) as geometry
from
timezones,
country
where
timezones.id in (
select
SpatialIndex.rowid
from
SpatialIndex,
country
where
f_table_name = 'timezones'
and search_frame = country.geometry
)
and intersects(timezones.geometry, country.geometry) = 1
除了提供名为 by_name
和 by_country
的预设查询外,此文件还包含显示我们用于数据库的数据来源的元数据。
它还设置了另一个重要选项
allow_sql: false
此选项阻止用户对我们发布的数据库执行自己的自定义 SQL 查询。只有我们定义的预设查询可用。
我们在这里使用该选项是因为 SpatiaLite 有大量函数,其中一些可能会导致底层的 Datasette 实例崩溃。
我们可以通过像这样启动 Datasette 来测试新的 metadata.yml
文件
datasette timezones.db --load-extension spatialite -m metadata.yml
将应用程序部署到 Fly
datasette publish 命令可用于将 Datasette 实例部署到各种不同的托管提供商。
Fly 是托管此应用程序的绝佳选择,因为该 API 可能会吸引大量流量。
Google Cloud Run 根据实例持续接收的流量量收费,这对该应用程序来说可能变得昂贵。
Fly 对实例收取固定的月费,外加带宽的额外费用。
此处可查看当前的 Fly 定价。截至撰写本文时,一个拥有 256MB RAM 的实例 - 足以舒适地托管此 API - 每月费用为 1.94 美元。
您需要安装 Fly CLI 工具
curl -L https://fly.io/install.sh | sh
然后运行以下命令向 Fly 进行身份验证
flyctl auth login
您还应该安装 datasette-publish-fly 插件
pip install datasette-publish-fly
准备好所有这些后,您可以像这样部署应用程序
datasette publish fly timezones.db \
--app timezones-api \
--setting default_page_size 10 \
--install datasette-geojson-map \
--metadata metadata.yml \
--spatialite
您需要一个唯一的 --app
名称(我已为此演示占用了 timezones-api
)。
--spatialite
标志确保为部署的应用程序配置了 SpatiaLite。
部署命令可能需要几分钟才能完成。完成后,您可以访问 https://your-app-name.fly.dev/
查看完成的应用程序,它已在线上。
以下是 timezones.datasette.io,就是使用此命令部署的。
源代码
本教程中所有内容的源代码都可以在 GitHub 上的 simonw/timezones-api 仓库中找到。