NSDT工具推荐： Three.js AI纹理开发包 - YOLO合成数据生成器 - GLTF/GLB在线编辑 - 3D模型格式在线转换 - 可编程3D场景编辑器 - REVIT导出3D模型插件 - 3D模型语义搜索引擎 - AI模型在线查看 - Three.js虚拟轴心开发包 - 3D模型在线减面 - STL模型在线切割 - 3D道路快速建模
postgis-users 邮件列表中的一位社区成员最近提出了一个问题：
我有一张高程点表，我想计算出穿过这些点的航线的高程剖面图。 如何？
这个问题很好地展示了我最喜欢的一些空间工具，包括索引、球体上的点对点距离查询和最近邻查询。 原问题的作者很好地分享了他的高程数据，所以我可以给你看一些图片和一种方法。
1、每2KM高程点
数据集为欧洲1.5M高程点，每2km一个点。 当绘制在地图上时，你可以看到它非常密集。
表结构很简单：
Column | Type | Collation | Nullable | Default
-----------+------------------+-----------+----------+---------
latlng | geography | | |
lat | double precision | | |
lng | double precision | | |
elevation | real | | |我们所需要的只是 latlng 列上的空间索引，以使我们的最近邻搜索更快。
CREATE INDEX elevation_latlng_x ON elevation USING GIST (latlng)2、航线
航线不同于陆上航线，因为它们不仅“像乌鸦飞过一样”，而且“像乌鸦飞过一个球体”。 在更长的距离和更高的纬度上，效果更明显，但大圆路线与点对点路线并不完全相同。
例如，伦敦和华沙之间的直线可以这样构造：
SELECT ST_MakeLine(ST_MakePoint(-0.061,51.482), -- London
ST_MakePoint(21.252,52.150)) -- Warsaw但是直线（平面上的最短路径）不是我们想要的，我们想要一个大圆（球体上的最短路径）。 因此，我们将线转换为地理类型，它使用大圆在顶点之间进行插值，然后对线进行 ST_Segmentize，以在端点之间添加新顶点。
SELECT ST_Segmentize(geography(
ST_MakeLine(
ST_MakePoint(-0.061,51.482),
ST_MakePoint(21.252,52.150))
), 4000)红线是大圆，蓝线是原来的笛卡尔线。
3、最近点
对于最终查询，我们将构建马德里和瓜达拉哈拉之间的飞行路径，这两个城市距离足够近，我们可以看到实际的高程点。 我们将沿着该路径每 4KM 放置一个顶点，并找到离每个顶点最近的高程点。
要找到最近点，我们使用 <-> “索引有序距离”运算符。 当放在查询的 ORDER BY 子句中时 <-> 不仅计算左右操作数之间的距离，它还返回从可能的最小距离开始的距离。
因此，查询离马德里最近的高程点 (-3.6996,40.4009) 是：
SELECT
e.elevation,
ST_AsText(e.latlng, 3) AS coords,
geography(ST_MakePoint(-3.6996,40.4009)) <-> e.latlng AS dist
FROM elevation e
ORDER BY dist
LIMIT 1返回结果如下（在 2ms 内，因为有索引辅助）：
elevation | coords | dist
-----------+----------------------+--------------------
650.8142 | POINT(-3.707 40.409) | 1045.41856429102934、最近的飞行路线
索引辅助最近邻运算符 <-> 是一个奇怪的运算符，因为它需要运算符一侧的文字值。 它不能接受两个变量参数，那么我们如何将具有多个点的飞行路径连接到高程网格？
答案是使用 LATERAL JOIN 来控制迭代并将飞行路径点（一次一个）提供给最近的邻居计算。
以块为单位，查询首先构建飞行路径。
WITH path AS (
SELECT ST_Segmentize(geography(
ST_MakeLine(
ST_MakePoint(-3.6996,40.4009), -- Madrid
ST_MakePoint(-3.1675,40.6582)) -- Guadalajara
), 4000) AS geog
)然后它使用 ST_DumpPoints 将路径分解为每个顶点的一条记录：
points AS (
SELECT (ST_DumpPoints(geog::geometry)).geom FROM path
)最后，它对高程点进行横向连接，为每个路径点找到最近的高程。
SELECT
d.elevation,
round(d.dist) AS distance,
ST_AsText(d.latlng, 3) AS elevation_pt,
ST_AsText(points.geom, 3) AS path_pt
FROM points
CROSS JOIN LATERAL (
SELECT
e.elevation, e.latlng,
points.geom::geography <-> e.latlng AS dist
FROM elevation e
ORDER BY dist
LIMIT 1
) d最终结果大约需要 14ms 运行，并找到飞行路线上每个点的高度。
elevation | distance | elevation_pt | path_pt
-----------+----------+----------------------+----------------------
650.8142 | 1045 | POINT(-3.707 40.409) | POINT(-3.7 40.401)
630.8268 | 563 | POINT(-3.66 40.415) | POINT(-3.666 40.417)
701.6845 | 728 | POINT(-3.641 40.436) | POINT(-3.633 40.433)
632.26556 | 882 | POINT(-3.595 40.443) | POINT(-3.6 40.449)
611.845 | 733 | POINT(-3.575 40.464) | POINT(-3.567 40.465)
607.7102 | 701 | POINT(-3.533 40.488) | POINT(-3.534 40.482)
644.5825 | 1057 | POINT(-3.509 40.491) | POINT(-3.501 40.498)
637.69604 | 166 | POINT(-3.467 40.515) | POINT(-3.467 40.514)
666.1685 | 1335 | POINT(-3.447 40.536) | POINT(-3.434 40.53)
673.436 | 373 | POINT(-3.401 40.543) | POINT(-3.401 40.546)
717.8294 | 938 | POINT(-3.358 40.567) | POINT(-3.368 40.562)
696.11383 | 910 | POINT(-3.335 40.57) | POINT(-3.334 40.578)
709.38855 | 723 | POINT(-3.292 40.594) | POINT(-3.301 40.594)
716.1025 | 697 | POINT(-3.273 40.615) | POINT(-3.268 40.61)
655.39514 | 862 | POINT(-3.226 40.621) | POINT(-3.234 40.626)
659.81757 | 501 | POINT(-3.207 40.642) | POINT(-3.201 40.642)
654.524 | 943 | POINT(-3.164 40.666) | POINT(-3.168 40.658)5、结束语
总结一下——我推荐的计算高程剖面和飞行路线的方法是：
- 创建空间索引以加快最近邻搜索
- 使用地理类型在线之间查询，因此最短路径在球体上而不是直线上
- 使用索引有序距离运算符查找两点之间的最短距离
- 使用最近邻运算符与高程点的横向连接。
原文链接：Elevation Profiles and Flightlines with PostGIS
BimAnt翻译整理，转载请标明出处
