PostGIS地理类型入门：地理坐标与笛卡尔坐标的区别及应用

PostGIS地理类型入门：地理坐标与笛卡尔坐标的区别及应用

PostGIS中的地理类型（geography）用于处理地理坐标数据，与传统的笛卡尔坐标系不同，地理坐标系采用球面坐标系，能够更准确地计算地球表面的距离和方向。本文将详细介绍地理类型在PostGIS中的使用方法及其优势。

18. 地理

坐标为“地理”或“纬度/经度”的数据很常见。与 Mercator、UTM 或 Stateplane 中的坐标不同，地理坐标不是笛卡尔坐标。地理坐标并不表示在平面上绘制的距原点的直线距离。相反，这些球坐标描述了地球上的角坐标。在球坐标中 ，点由相对于参考子午线（经度）的旋转角度和相对于赤道（纬度）的角度来指定。

您可以将地理坐标视为近似笛卡尔坐标并继续进行空间计算。然而，距离、长度和面积的测量是没有意义的。由于球坐标测量角距离，因此单位为“度”。此外，索引和真/假测试（例如相交和包含）的近似结果可能会变得非常错误。随着接近两极或国际日期变更线等问题区域，点之间的距离会变得更大。

例如，这是洛杉矶和巴黎的坐标。

• 洛杉矶:
  POINT(-118.4079 33.9434)

• 巴黎:
  POINT(2.3490 48.8533)

以下使用标准 PostGIS 笛卡尔坐标计算洛杉矶和巴黎之间的距离。请注意，SRID 4326 声明了地理空间参考系统。

SELECT ST_Distance(
 'SRID=4326;POINT(-118.4079 33.9434)'::geometry, -- Los Angeles (LAX)
 'SRID=4326;POINT(2.5559 49.0083)'::geometry -- Paris (CDG)
);

121.898285970107

啊哈！ 122！但是，这是什么意思？

空间参考4326的单位是度。所以我们的答案是122度。但是（再次），这是什么意思？

在球体上，一个“平方度”的大小变化很大，随着远离赤道而变得更小。想象一下，当您向两极移动时，地球上的经线（垂直线）会变得更加接近。所以，122度的距离没有任何意义。这是一个无意义的数字。

为了计算有意义的距离，我们必须将地理坐标视为真正的球面坐标，而不是近似的笛卡尔坐标。我们必须将点之间的距离测量为球体（大圆的一部分）上的真实路径。

PostGIS通过geography类型提供此功能。

注解
不同的空间数据库有不同的“地理处理”方法

• 当SRID是地理位置时，Oracle试通过穿透式的地理坐标系计算来掩盖差异。
• SQL Server使用两种空间类型，“STGeometry”用于笛卡尔数据，“STGeography”用于地理。
• Informix Spatial是Informix的纯笛卡尔扩展，而Informix Geodetic是纯地理扩展。
• 与SQL Server类似，PostGIS使用两种类型，“几何体”和“地理”。

用geography而不是geometry类型，让我们再次尝试测量洛杉矶和巴黎之间的距离。

SELECT ST_Distance(
 'SRID=4326;POINT(-118.4079 33.9434)'::geography, -- Los Angeles (LAX)
 'SRID=4326;POINT(2.5559 49.0083)'::geography -- Paris (CDG)
);

9124665.27317673

一个很大的数字！计算的所有返回值geography 均以米为单位，因此我们的答案是 9125km。

旧版本的 PostGIS 支持使用ST_Distance_Spheroid(point, point, measurement)函数对球体进行非常基本的计算。然而，ST_Distance_Spheroid实质上是有限的。该函数仅适用于点，不支持跨极点或国际日期变更线的索引。

当提出“从洛杉矶到巴黎的航班到冰岛有多近？”这样的问题时，支持非点几何的需求就变得非常明显

在笛卡尔平面（紫色线）上使用地理坐标确实会产生非常错误的答案！使用大圆路线（红线）给出了正确的答案。如果我们将我们的LAX-CDG航班转换为字符串，并使用geography计算到冰岛某个点的距离，我们将得到以米为单位的正确答案（召回）。

SELECT ST_Distance(
 ST_GeographyFromText('LINESTRING(-118.4079 33.9434, 2.5559 49.0083)'), -- LAX-CDG
 ST_GeographyFromText('POINT(-22.6056 63.9850)') -- Iceland (KEF)
);

502454.906643729

因此，在LAX-CDG航线上，离冰岛最近的距离（从其国际机场测量）是相对较小的502公里。

对于跨越国际日期变更线的特征，处理地理坐标的笛卡尔方法完全失效。从洛杉矶到东京最短的大圆航线要穿越太平洋。最短的笛卡尔路线穿过大西洋和印度洋。

SELECT ST_Distance(
 ST_GeometryFromText('Point(-118.4079 33.9434)'), -- LAX
 ST_GeometryFromText('Point(139.733 35.567)')) -- NRT (Tokyo/Narita)
 AS geometry_distance,
ST_Distance(
 ST_GeographyFromText('Point(-118.4079 33.9434)'), -- LAX
 ST_GeographyFromText('Point(139.733 35.567)')) -- NRT (Tokyo/Narita)
 AS geography_distance;

 geometry_distance | geography_distance
-------------------+--------------------
 258.146005837336 | 8833954.76996256

18.1. 使用地理

为了将几何数据加载到地图表中，首先需要将几何投影到EPSG:4326(经度/纬度)中，然后需要将其更改为地理。ST_Transform(geometry,srid)函数将坐标转换为地理位置，而Geography(geometry)函数或 ::geography 后缀将坐标“强制转换”为地理位置。

CREATE TABLE nyc_subway_stations_geog AS
SELECT
 ST_Transform(geom,4326)::geography AS geog,
 name,
 routes
FROM nyc_subway_stations;

在地理表中创建空间索引与在几何表中创建空间索引完全相同:

CREATE INDEX nyc_subway_stations_geog_gix
ON nyc_subway_stations_geog USING GIST (geog);

不同之处在于:地理索引可以正确处理覆盖两极或国际日期变更线的查询，而几何索引则不能。

这里有一个查询，可以找到帝国大厦500米范围内的所有地铁站。

WITH empire_state_building AS (
 SELECT 'POINT(-73.98501 40.74812)'::geography AS geog
)
SELECT name,
 ST_Distance(esb.geog, ss.geog) AS distance,
 degrees(ST_Azimuth(esb.geog, ss.geog)) AS direction
FROM nyc_subway_stations_geog ss,
 empire_state_building esb
WHERE ST_DWithin(ss.geog, esb.geog, 500);

地理类型只有少量本地函数:

• ST_AsText(geography)返回text类型
• ST_GeographyFromText(text)返回geography类型
• ST_AsBinary(geography)返回bytea类型
• ST_GeogFromWKB(bytea)返回geography类型
• ST_AsSVG(geography)返回text类型
• ST_AsGML(geography)返回text类型
• ST_AsKML(geography)返回text类型
• ST_AsGeoJson(geography)返回text类型
• ST_Distance(geography, geography)返回double类型
• ST_DWithin(geography, geography, float8)返回boolean类型
• ST_Area(geography)返回double类型
• ST_Length(geography)返回double类型
• ST_Covers(geography, geography)返回boolean类型
• ST_CoveredBy(geography, geography)返回boolean类型
• ST_Intersects(geography, geography)返回boolean类型
• ST_Buffer(geography, float8)返回 geography类型
• ST_Intersection(geography, geography)返回geography类型

18.2. 创建地理表

用于创建包含地理列的新表的SQL与用于创建几何表的SQL非常相似。但是，地理包括在创建表时直接指定对象类型的能力。例如:

CREATE TABLE airports (
 code VARCHAR(3),
 geog GEOGRAPHY(Point)
 );
INSERT INTO airports
 VALUES ('LAX', 'POINT(-118.4079 33.9434)');
INSERT INTO airports
 VALUES ('CDG', 'POINT(2.5559 49.0083)');
INSERT INTO airports
 VALUES ('KEF', 'POINT(-22.6056 63.9850)');

在表定义中，

GEOGRAPHY(Point)

将机场数据类型指定为点。新的地理字段不会在geometry_columns视图中注册。相反，它们被注册在一个名为geography_columns的视图中。

SELECT * FROM geography_columns;

 f_table_name | f_geography_column | srid | type
--------------------------+--------------------+------+----------
 nyc_subway_stations_geog | geog | 0 | Geometry
 airports | geog | 4326 | Point

注解
上面的输出中省略了一些列。

18.3. 转换为几何图形

虽然地理类型的基本函数可以处理许多用例，但有时您可能需要访问仅由几何类型支持的其他函数。幸运的是，您可以在地理和几何之间来回转换对象。

PostgreSQL的语法惯例是将::typename追加到您希望强制转换的值的末尾。因此，

2::text

将把数字2转换为文本字符串'2'。

'POINT(0 0)'::geometry

将点的文本表示转换为几何点。

ST_X(point)函数仅支持几何类型。我们如何从地理位置读取 X 坐标？

SELECT code, ST_X(geog::geometry) AS longitude FROM airports;

 code | longitude
------+-----------
 LAX | -118.4079
 CDG | 2.5559
 KEF | -21.8628

By appending

::geometry

to our geography value, we convert the object to a geometry with an SRID of 4326. From there we can use as many geometry functions as strike our fancy. But, remember -- now that our object is a geometry, the coordinates will be interpreted as Cartesian coordinates, not spherical ones.

18.4. 为什么(不)使用地理

地理坐标是被普遍接受的坐标——每个人都知道经纬度是什么意思，但很少有人知道UTM坐标是什么意思。为什么不一直用地理呢?

• 首先，如前所述，直接支持地理类型的可用功能(目前)要少得多。您可能需要花费大量时间来解决地理类型的限制。
• 其次，在球体上的计算比笛卡尔计算要昂贵得多。例如，距离的笛卡尔公式(毕达哥拉斯)涉及一次调用sqrt()。距离的球面公式(Haversine)涉及两次sqrt()调用，一次arctan()调用，四次sin()调用和两次cos()调用。三角函数是非常昂贵的，球面计算涉及到很多三角函数。

结论?

如果您的数据在地理上紧凑的(包含在一个州，县或市)，使用具有笛卡尔投影的几何类型，这对您的数据有意义。请参阅http://epsg.io网站，并键入您所在地区的名称，以选择可能的参考系统。

如果你需要用地理上分散的数据集(覆盖世界的大部分)来测量距离，使用地理类型。在geography中工作节省的应用程序复杂性将抵消任何性能问题。而转换为geometry 可以抵消大多数功能限制。

18.5. 函数列表

• ST_Distance(geometry, geometry): 对于几何类型，返回两个几何图形之间的二维笛卡尔最小距离（基于空间参考）（以投影单位表示）。对于地理类型，默认返回两个地理区域之间的球体最小距离（以米为单位）。
• ST_GeographyFromText(text):从众所周知的文本表示或扩展 (WKT) 返回指定的地理值。
• ST_Transform(geometry, srid):返回一个新的几何图形，其坐标转换为整数参数引用的 SRID。
• ST_X: 返回点的 X 坐标，如果不可用则返回 NULL。输入必须是一个点。
• ST_Azimuth(geography_A, geography_B): 以弧度为单位返回从A到B的方向。
• ST_DWithin(geography_A, geography_B, R):如果A距离B在R米以内返回true。

脚注
缓冲区和交集函数实际上是在转换为几何图形之上的包装器，并不是在球坐标中进行的。因此，对于无法清晰地转换为平面表示的范围非常大的对象，它们可能无法返回正确的结果。

例如，ST_Buffer(geography,distance)函数将地理对象转换为“最佳”投影，对其进行缓冲，然后将其转换回地理。如果没有“最佳”投影(对象太大)，则操作可能失败或返回格式错误的缓冲区。