一、概述
Geopoly 模块是R-Tree 扩展的替代接口,它使用GeoJSON符号 ( RFC-7946 ) 来描述二维多边形。Geopoly 包括用于检测一个多边形何时包含在另一个多边形中或与另一个多边形重叠、计算多边形包围的面积、对多边形进行线性变换、将多边形渲染为 SVG以及其他类似操作的功能。
Geopoly 的源代码包含在合并中,但不包含在库中,除非使用-DSQLITE_ENABLE_GEOPOLY编译时选项。
Geopoly 在“简单”多边形上运行 - 即边界不与自身相交的多边形。Geopoly 因此扩展了R-Tree 扩展的功能,它只能处理矩形区域。另一方面,R-Tree 扩展能够处理 1 到 5 个坐标维度,而 Geopoly 仅限于二维形状。
Geopoly 模块中的每个多边形都可以与任意数量的辅助数据字段相关联。
1.1. GeoJSON
GeoJSON 标准是使用 JSON 交换地理空间信息的语法。GeoJSON 是一个丰富的标准,可以描述几乎任何类型的地理空间内容。
Geopoly 模块只理解 GeoJSON 的一小部分,但却是一个关键的子集。特别是,GeoJSON 理解描述简单多边形的 JSON 顶点数组。
多边形由其顶点定义。每个顶点都是一个包含两个数值的 JSON 数组,这两个数值是顶点的 X 和 Y 坐标。多边形是至少包含四个顶点的 JSON 数组,因此是数组的数组。数组中的第一个和最后一个顶点必须相同。多边形遵循右手法则:当从一个顶点到下一个顶点追踪一条线时,线右侧的区域在多边形之外,左侧区域在多边形内部。换句话说,顶点的净旋转是逆时针的。
例如,以下 JSON 描述了一个位于 X 轴上且面积为 0.5 的等腰三角形:
[[0,0],[1,0],[0.5,1],[0,0]]
一个三角形有三个顶点,但三角形的 GeoJSON 描述有 4 个顶点,因为第一个和最后一个顶点是重复的。
1.2. 二进制存储格式
在内部,Geopoly 以二进制格式存储多边形 - SQL BLOB。下面给出了二进制格式的详细信息。所有 Geopoly 接口都能够接受 GeoJSON 格式或二进制格式的多边形。
2.使用 Geopoly 扩展
geopoly 表创建如下:
CREATE VIRTUAL TABLE newtab USING geopoly(a,b,c);
上面的语句创建了一个名为“newtab”的新 geopoly 表。每个 geopoly 表都包含一个内置整数“rowid”列和一个“_shape”列,其中包含与该表的该行关联的多边形。上面的示例还定义了三个名为“a”、“b”和“c”的辅助数据列,它们可以存储应用程序需要与每个多边形关联的任何附加信息。如果不需要存储辅助信息,可以省略辅助列列表。
使用普通的 INSERT 语句将新多边形存储在表中:
INSERT INTO newtab(_shape) VALUES('[[0,0],[1,0],[0.5,1],[0,0]]');
UPDATE 和 DELETE 语句的工作方式类似。
2.1. 查询
要使用索引地理空间搜索查询 geopoly 表,请在 WHERE 子句中使用函数 geopoly_overlap() 或 geopoly_within() 作为布尔函数,并将“_shape”列作为函数的第一个参数。例如:
SELECT * FROM newtab WHERE geopoly_overlap(_shape, $query_polygon);
前面的示例将返回 _shape 与 $query_polygon 参数中的多边形重叠的每一行。geopoly_within() 函数的工作方式类似,但只返回 _shape 完全包含在 $query_polygon 中的行。
WHERE 子句包含裸 geopoly_overlap() 或 geopoly_within() 函数的查询(以及 DELETE 和 UPDATE 语句)利用底层 R*Tree 数据结构进行快速查找,只需检查中的行的子集桌子。当然,检查的行数取决于 $query_polygon 的大小。大的 $query_polygons 通常需要查看比小的更多的行。
对 geopoly 表的 rowid 的查询也非常快,即使对于具有大量行的表也是如此。但是,辅助数据列都不是索引,因此针对辅助数据列的查询将涉及全表扫描。
三、特殊功能
geopoly 模块定义了几个可用于处理多边形的新 SQL 函数。这些函数的所有多边形参数可以是 GeoJSON 格式或内部二进制格式。
3.1. geopoly_overlap(P1,P2) 函数
如果 P1 和 P2 都是多边形,则 geopoly_overlap(P1,P2) 函数在 P1 和 P2 之间存在任何重叠时返回一个非零整数,或者在 P1 和 P2 完全不相交时返回零。如果 P1 或 P2 不是多边形,则此例程返回 NULL。
geopoly_overlap(P1,P2) 函数的特殊之处在于 geopoly 虚拟表知道如何使用 R*Tree 索引来优化其中 WHERE 子句使用 geopoly_overlap() 作为布尔函数的查询。只有 geopoly_overlap(P1,P2) 和 geopoly_within(P1,P2) 函数具有此功能。
3.2. geopoly_within(P1,P2) 函数
如果 P1 和 P2 都是多边形,则如果 P1 完全包含在 P2 中,则 geopoly_within(P1,P2) 函数返回一个非零整数,或者如果 P1 的任何部分在 P2 之外,则返回零。如果 P1 和 P2 是同一个多边形,则此例程返回非零值。如果 P1 或 P2 不是多边形,则此例程返回 NULL。
geopoly_within(P1,P2) 函数的特殊之处在于 geopoly 虚拟表知道如何使用 R*Tree 索引来优化其中 WHERE 子句使用 geopoly_within() 作为布尔函数的查询。只有 geopoly_within(P1,P2) 和 geopoly_overlap(P1,P2) 函数具有此功能。
3.3. geopoly_area(P) 函数
如果 P 是多边形,则 geopoly_area(P) 返回该多边形所包围的面积。如果 P 不是多边形,则 geopoly_area(P) 返回 NULL。
3.4. geopoly_blob(P) 函数
如果 P 是多边形,则 geopoly_blob(P) 返回该多边形的二进制编码作为 BLOB。如果 P 不是多边形,则 geopoly_blob(P) 返回 NULL。
3.5. geopoly_json(P) 函数
如果 P 是多边形,则 geopoly_json(P) 返回该多边形的 GeoJSON 表示形式作为文本字符串。如果 P 不是多边形,则 geopoly_json(P) 返回 NULL。
3.6. geopoly_svg(P,...) 函数
如果 P 是一个多边形,则 geopoly_svg(P,...) 返回一个文本字符串,它是该多边形的 可缩放矢量图形 (SVG) 表示形式。如果有多个参数,则第二个和后续参数将作为属性添加到每个 SVG 字形。例如:
SELECT geopoly_svg($polygon,'class="poly"','style="fill:blue;"');
如果 P 不是多边形,则 geopoly_svg(P,...) 返回 NULL。
请注意,geopoly 使用传统的右手笛卡尔坐标系,原点在左下角,而 SVG 使用左手坐标系,原点在左上角。geopoly_svg() 例程不尝试变换坐标系,因此显示的图像以镜像显示并旋转。如果不需要,可以使用 geopoly_xform() 例程将输出从笛卡尔坐标转换为 SVG 坐标,然后再将多边形传递到 geopoly_svg()。
3.7. geopoly_bbox(P) 和 geopoly_group_bbox(P) 函数
如果 P 是多边形,则 geopoly_bbox(P) 返回一个新的多边形,它是完全包围 P 的最小(轴对齐)矩形。如果 P 不是多边形,则 geopoly_bbox(P) 返回 NULL。
geopoly_group_bbox(P) 函数是 geopoly_bbox(P) 的聚合版本。geopoly_group_bbox(P) 函数返回将包含聚合期间看到的所有 P 值的最小矩形。
3.8. geopoly_contains_point(P,X,Y) 函数
如果 P 是多边形,则 geopoly_contains_point(P,X,Y) 返回一个非零整数当且仅当坐标 X,Y 在多边形 P 的内部或边界上。如果 P 不是多边形,则 geopoly_contains_point( P,X,Y) 返回 NULL。
3.9. geopoly_xform(P,A,B,C,D,E,F) 函数
geopoly_xform(P,A,B,C,D,E,F) 函数返回一个新的多边形,它是多边形 P 的仿射变换,其中变换由值 A、B、C、D、E、F 定义. 如果 P 不是有效的多边形,则此例程返回 NULL。
变换根据以下公式转换多边形的每个顶点:
x1 = A*x0 + B*y0 + E y1 = C*x0 + D*y0 + F
因此,例如,要将多边形移动一定量 DX、DY 而不改变其形状,请使用:
geopoly_xform($polygon, 1, 0, 0, 1, $DX, $DY)
围绕点 0, 0 将多边形旋转 R 弧度:
geopoly_xform($polygon, cos($R), sin($R), -sin($R), cos($R), 0, 0)
请注意,翻转多边形的变换可能会导致顶点的顺序颠倒。换句话说,变换可能导致顶点以顺时针顺序而不是逆时针顺序循环。这可以通过在转换后通过geopoly_ccw()函数发送结果来纠正。
3.10. geopoly_regular(X,Y,R,N) 函数
geopoly_regular(X,Y,R,N) 函数返回一个具有 N 条边、以 X、Y 为中心且外接半径为 R 的凸、简单、规则、等边、等角多边形。或者,如果 R 为负或如果 N小于 3,函数返回 NULL。N 值的上限为 1000,因此即使 N 值大于 1000,例程也不会渲染边数超过 1000 的多边形。
例如,下图:
由这个脚本生成:
SELECT '<svg width="600" height="300">'; WITH t1(x,y,n,color) AS (VALUES (100,100,3,'red'), (200,100,4,'orange'), (300,100,5,'green'), (400,100,6,'blue'), (500,100,7,'purple'), (100,200,8,'red'), (200,200,10,'orange'), (300,200,12,'green'), (400,200,16,'blue'), (500,200,20,'purple') ) SELECT geopoly_svg(geopoly_regular(x,y,40,n), printf('style="fill:none;stroke:%s;stroke-width:2"',color)) || printf(' <text x="%d" y="%d" alignment-baseline="central" text-anchor="middle">%d</text>',x,y+6,n) FROM t1; SELECT '</svg>';
3.11. geopoly_ccw(J) 函数
geopoly_ccw(J) 函数返回逆时针 (CCW) 旋转的多边形 J。
RFC-7946要求多边形使用 CCW 旋转。但该规范还指出,许多遗留 GeoJSON 文件不遵循该规范,并且包含顺时针 (CW) 旋转的多边形。geopoly_ccw() 函数对于读取遗留 GeoJSON 脚本的应用程序很有用。如果 geopoly_ccw() 的输入是格式正确的多边形,则不会进行任何更改。但是,如果输入多边形的循环是向后的,则 geopoly_ccw() 会反转循环顺序,使其符合规范,以便与 Geopoly 模块一起正常工作。
4.实施细节
geopoly 模块是R-Tree 扩展的扩展。Geopoly 使用与R-Tree 扩展相同的底层逻辑和影子表。Geopoly 只是提供了一个不同的接口,并提供了一些额外的逻辑来计算多边形解码、重叠和包含。
4.1. 多边形的二进制编码
Geopoly 使用二进制格式在内部存储所有多边形。二进制多边形包含一个 4 字节的标头,后面是一个坐标对数组,其中每个坐标的每个维度都是一个 32 位浮点数。
标头的第一个字节是标志字节。标志字节的最低有效位决定了标头后面的坐标对是以大端还是小端存储的。最低有效位的值为 0 表示大端,值为 1 表示小端。头中第一个字节的其他位保留用于将来的扩展。
标头中接下来的三个字节将多边形中的顶点数记录为大端整数。因此每个多边形有大约 1600 万个顶点的上限。
标题后面是坐标对数组。每个坐标都是一个 32 位浮点数。使用 32 位浮点值作为坐标意味着地球表面上的任何点都可以以大约 2.5 米的分辨率绘制。如果地图仅限于单个大陆或国家,则当然可以使用更高的分辨率。请注意,geopoly 模块中的坐标分辨率在幅度上与潮汐力引起的地球表面点的日常移动相似。
二进制格式的坐标列表不包含冗余。最后一个坐标不是第一个坐标的重复,因为它与 GeoJSON 一样。因此,与 GeoJSON 表示相比,多边形的二进制表示总是少一对坐标。
4.2. 影子表
geopoly 模块建立在R-Tree 扩展之上,并使用相同的底层影子表和算法。出于索引目的,每个多边形在影子表中表示为矩形边界框。底层 R-Tree 实现使用边界框来限制搜索空间。然后 geoploy_overlap() 和/或 geopoly_within() 例程进一步细化搜索以得到准确的答案。