在 ClickHouse，我们经常发现自己需要下载和探索越来越大的数据集，这不可避免地需要快速的互联网连接。当我想搬家和换国，并且需要确保自己拥有快速的互联网时，我发现了庞大的Ookla 数据集，它提供了全球所有位置的互联网性能数据。该数据集代表了全球最大的众包网络测试来源，涵盖了固定和移动互联网性能，覆盖了全球大部分陆地面积。
虽然 Ookla 已经提供了优秀的工具来使用该数据集按位置探索互联网性能，但这为我们提供了探索 ClickHouse 的 Apache Iceberg 支持以及不断扩展的地理函数的机会。
进一步限制一下，我很好奇仅使用 SQL 将地理多边形转换为 SVG 来可视化地理数据的可行性。这将带我们了解 h3 索引、多边形字典、墨卡托投影、颜色插值以及使用 UDF 计算质心。如果您好奇我们如何仅使用 SQL 得到以下结果，请继续阅读…
注意：ClickHouse 支持许多适合此数据集的可视化工具，包括但不限于 Superset。我们仅使用 SQL 的限制仅仅是为了好玩！
在深入研究数据集之前，我们将快速回顾一下 Apache Iceberg 的价值及其对 ClickHouse 的支持和相关性。对于熟悉此主题的人，可以随意跳到我们开始查询的部分！
Apache Iceberg
Apache Iceberg 在近年来获得了广泛关注，它是“数据湖”概念向“湖仓”演进的关键。通过为数据提供高性能的表格式，Iceberg 为存储在数据湖中的结构化数据带来了类似表的语义。更具体地说，这种开放的表格式与供应商无关，并允许数据文件（通常为结构化的 Parquet）作为表公开，并支持模式演变、删除、更新、事务，甚至 ACID 合规性。这代表了从基于 Hadoop 和 Hive 等技术的更传统的数据湖的逻辑演进。
也许最重要的是，这种数据格式允许像 Spark 和 Trino 这样的查询引擎安全地同时处理相同的表，同时还提供诸如分区之类的功能，这些功能可能会被查询计划器利用来加速查询。
ClickHouse 支持
虽然 ClickHouse 远不止一个简单的查询引擎，但我们认识到这些开放标准对于允许用户避免供应商锁定以及作为其他工具的合理数据交换格式的重要性。正如最近一篇博文中所讨论的，“云数据仓库的解耦”，我们认为湖仓是现代数据架构中的一个基本组件，与实时分析数据库完美互补。在这种架构中，湖仓提供经济高效的冷存储和一个共同的真相来源，并支持 ad hoc 查询。
但是，对于更苛刻的使用案例，例如需要低延迟查询和支持高并发性的面向用户的密集型数据应用程序，可以将子集移动到 ClickHouse。这允许用户从 Iceberg 的开放交换格式中受益，用于其冷存储和长期保留，同时使用 ClickHouse 作为实时数据仓库，以满足查询性能至关重要的需求。
使用 ClickHouse 查询 Iceberg
实现这一愿景需要能够从 ClickHouse 查询 Iceberg 文件。作为一种湖仓格式，这目前通过一个专用的表iceberg 函数来支持。这假设数据托管在兼容 S3 的服务上，例如 AWS S3 或 Minio。
在我们的示例中，我们已将 Ookla 数据集提供在以下公共 S3 存储桶中，供用户进行实验。
s3://datasets-documentation/ookla/iceberg/
我们可以使用带有 iceberg 函数（替换表名）的 DESCRIBE 查询来描述 Ookla 数据的模式。与依赖于实际数据采样的 Parquet 不同，模式是从 Iceberg 元数据文件中读取的。
DESCRIBE TABLE iceberg('https://datasets-documentation.s3.eu-west-3.amazonaws.com/ookla/iceberg/')
SETTINGS describe_compact_output = 1
┌─name────────────┬─type─────────────┐
│ quadkey │ Nullable(String) │
│ tile │ Nullable(String) │
│ avg_d_kbps │ Nullable(Int32) │
│ avg_u_kbps │ Nullable(Int32) │
│ avg_lat_ms │ Nullable(Int32) │
│ avg_lat_down_ms │ Nullable(Int32) │
│ avg_lat_up_ms │ Nullable(Int32) │
│ tests │ Nullable(Int32) │
│ devices │ Nullable(Int32) │
│ year_month │ Nullable(Date) │
└─────────────────┴──────────────────┘
10 rows in set. Elapsed: 0.216 sec.
类似地，可以使用标准 SQL（使用表函数替换表名）来查询行。在下面，我们对一些行进行采样并计算数据集的总大小。
SELECT *
FROM iceberg('https://datasets-documentation.s3.eu-west-3.amazonaws.com/ookla/iceberg/')
LIMIT 1
FORMAT Vertical
Row 1:
──────
quadkey: 1202021303331311
tile: POLYGON((4.9163818359375 51.2206474303833, 4.921875 51.2206474303833, 4.921875 51.2172068072334, 4.9163818359375 51.2172068072334, 4.9163818359375 51.2206474303833))
avg_d_kbps: 109291
avg_u_kbps: 15426
avg_lat_ms: 12
avg_lat_down_ms: ᴺᵁᴸᴸ
avg_lat_up_ms: ᴺᵁᴸᴸ
tests: 6
devices: 4
year_month: 2021-06-01
1 row in set. Elapsed: 2.100 sec. Processed 8.19 thousand rows, 232.12 KB (3.90 thousand rows/s., 110.52 KB/s.)
Peak memory usage: 4.09 MiB.
SELECT count()
FROM iceberg('https://datasets-documentation.s3.eu-west-3.amazonaws.com/ookla/iceberg/')
┌───count()─┐
│ 128006990 │
└───────────┘
1 row in set. Elapsed: 0.701 sec. Processed 128.01 million rows, 21.55 KB (182.68 million rows/s., 30.75 KB/s.)
Peak memory usage: 4.09 MiB
我们可以看到，Ookla 使用WKT 格式将地球表面划分为多边形。我们将在博客的后半部分深入探讨这一点；现在，我们将使用聚合计算每个多边形的平均设备数量以及下载速度的中位数、第 90 百分位数、第 99 百分位数和第 99.9 百分位数来完成我们的简单示例。
SELECT
round(avg(devices), 2) AS avg_devices,
arrayMap(m -> round(m / 1000), quantiles(0.5, 0.9, 0.99, 0.999)(avg_d_kbps)) AS download_mbps
FROM iceberg('https://datasets-documentation.s3.eu-west-3.amazonaws.com/ookla/iceberg/')
┌─avg_devices─┬─download_mbps────┐
│ 5.64 │ [51,249,502,762] │
└─────────────┴──────────────────┘
1 row in set. Elapsed: 4.704 sec. Processed 128.01 million rows, 3.75 GB (27.21 million rows/s., 797.63 MB/s.)
Peak memory usage: 22.89 MiB.
以上提供了一些关于其实时数据仓库如何运作的简单示例：湖仓中的冷 Iceberg 存储以及 ClickHouse 中的热查询以进行实时分析。
一些限制
以上只是一些非常简单的查询。在我们深入探讨如何可视化这些数据之前，值得强调 ClickHouse 的 Iceberg 支持需要进一步改进（并已计划）的一些领域。
- 目前，ClickHouse 不会读取具有基于行的删除 (v2) 或模式已发生演变的 Iceberg 文件。除非用户明确请求忽略它，否则他们会在这些情况下收到读取异常。
- 在查询评估期间，ClickHouse 目前不会利用分区来修剪扫描的数据。这代表了改进查询执行的重大可能性，并与最近对 Parquet 做出的改进相一致。
- 当前的实现绑定到 S3。将来，我们预计这种耦合将被移除，并将添加对 ClickHouse 支持的任何对象存储的支持。
- 与 ClickHouse 的排序键类似，Iceberg 支持对数据进行排序。通过确保 Iceberg 表中的数据文件被拆分和排序，查询可以在执行期间跳过文件。这有可能进一步减少扫描的数据量，从而提高性能。
探索数据
在完成了基础知识之后，让我们稍微调整一下这篇博文，更详细地探讨一下 Ookla 数据，目标是仅使用 SQL 来可视化互联网速度。
读取和丰富
如前所述，Ookla 数据将地球划分为大小不一的多个多边形，从中收集了互联网性能统计数据。这些多边形在WKT 格式中定义，这是一种用于表示空间对象的标记语言。由于底层的 Parquet 文件没有原生的多边形表示，因此这种格式对 Iceberg 和数据分发很有意义。
POLYGON((-160.02685546875 70.6435894914449, -160.021362304688 70.6435894914449, -160.021362304688 70.6417687358462, -160.02685546875 70.6417687358462, -160.02685546875 70.6435894914449))
但是，查询此数据需要转换为数值表示。对于 ClickHouse，可以使用readWKTPolygon 函数执行此操作。如下所示，这会将文本表示转换为Polygon 类型，我们稍后可以在该类型上使用ClickHouse 的大量几何函数。请注意，这些多边形包含经度和纬度值。
SELECT readWKTPolygon(assumeNotNull(tile)) AS poly, toTypeName(poly) as type
FROM iceberg('https://datasets-documentation.s3.eu-west-3.amazonaws.com/ookla/iceberg/')
LIMIT 1
FORMAT Vertical
Row 1:
──────
poly: [[(4.9163818359375,51.2206474303833),(4.921875,51.2206474303833),(4.921875,51.2172068072334),(4.9163818359375,51.2172068072334),(4.9163818359375,51.2206474303833)]]
type: Polygon
1 row in set. Elapsed: 1.920 sec. Processed 8.19 thousand rows, 232.12 KB (4.27 thousand rows/s., 120.91 KB/s.)
Peak memory usage: 4.10 MiB.
对于多边形和解析 WKT 数据的入门介绍，我们推荐观看我们自己的“Data with Mark” 中的以下视频。
请注意，以上视频使用了 Ookla 数据集中最新的文件，其中提供了质心。
虽然较新版本的 Ookla 数据集包含每个多边形的预计算质心（从 2023 年第二季度开始），但该数据集可以追溯到 2019 年。正如我们稍后将展示的那样，此质心（经纬度对）在渲染数据时非常有用。幸运的是，由于我们的多边形是非自交且封闭的，因此计算起来相对简单。对于好奇的读者，下面显示了一些 Python 代码示例（来源）。
def centroid(vertices):
x, y = 0, 0
n = len(vertices)
signed_area = 0
# computes the signed area and y and x accumulators - signed_area UDF function
for i in range(len(vertices)):
x0, y0 = vertices[i]
x1, y1 = vertices[(i + 1) % n]
# shoelace formula - maps to area function below
area = (x0 * y1) - (x1 * y0)
signed_area += area
x += (x0 + x1) * area
y += (y0 + y1) * area
# final step of centroid function
signed_area *= 0.5
x /= 6 * signed_area
y /= 6 * signed_area
return x, y
x, y = centroid([(-160.037841796875, 70.6363054807905), (-160.032348632812, 70.6363054807905),(-160.032348632812, 70.6344840663086), (-160.037841796875, 70.6344840663086), (-160.037841796875, 70.6363054807905)])
print(x,y)
python centroid.py
-160.03509089358298 70.63539286610323
如果我们将问题分解成几个函数，那么将其转换为我们在查询时可以使用的 SQL UDF 将是最简单的。
CREATE FUNCTION area AS (v1,v2) -> ((v1.1 * v2.2) - (v2.1 * v1.2))
CREATE FUNCTION signed_area AS polygon ->
(
arrayFold((acc, x) -> (
(acc.1) + ((((x.1).1) + ((polygon[((x.2) + 1) % (length(polygon) + 1)]).1)) * area(x.1, polygon[((x.2) + 1) % (length(polygon) + 1)])),
(acc.2) + ((((x.1).2) + ((polygon[((x.2) + 1) % (length(polygon) + 1)]).2)) * area(x.1, polygon[((x.2) + 1) % (length(polygon) + 1)])),
(acc.3) + area(x.1, polygon[((x.2) + 1) % (length(polygon) + 1)])), arrayZip(polygon, range(1, length(polygon) + 1)), (0.0, 0.0, 0.0)
)
)
CREATE FUNCTION centroid as polygon -> ((signed_area(polygon).1) / ((6 * (signed_area(polygon).3)) * 0.5), (signed_area(polygon).2) / ((6 * (signed_area(polygon).3)) * 0.5))
快速测试证实了我们的 UDF 给出的结果与 Python 实现相同。
SELECT centroid([(-160.037841796875, 70.6363054807905), (-160.032348632812, 70.6363054807905), (-160.032348632812, 70.6344840663086), (-160.037841796875, 70.6344840663086), (-160.037841796875, 70.6363054807905)]) AS centroid
┌─centroid────────────────────────────────┐
│ (-160.03509089358298,70.63539286610323) │
└─────────────────────────────────────────┘
1 row in set. Elapsed: 0.017 sec.
使用上述函数，我们现在可以读取我们的多边形数据并计算由经度和纬度组成的质心。正如我们将展示的那样，在可视化数据时，我们可以利用这些质心。
加载数据
虽然以下所有查询都可以与 iceberg 表函数一起使用，但这需要为每个查询下载数据，通常会将性能限制在我们的网络带宽上。由于我们将要运行许多查询，因此将数据加载到 ClickHouse 中以享受前面突出显示的性能优势是有意义的。我们修改表模式以包含转换后的多边形（Ring 类型 - 见下文）和质心（Point 类型），并使用上面的函数提取这些数据。
CREATE TABLE default.ookla
(
`quadkey` Nullable(String),
`polygon` Ring,
`centroid` Point,
`avg_d_kbps` Nullable(Int32),
`avg_u_kbps` Nullable(Int32),
`avg_lat_ms` Nullable(Int32),
`avg_lat_down_ms` Nullable(Int32),
`avg_lat_up_ms` Nullable(Int32),
`tests` Nullable(Int32),
`devices` Nullable(Int32),
`year_month` Nullable(Date)
)
ENGINE = MergeTree ORDER BY tuple()
INSERT INTO ookla SELECT
quadkey,
readWKTPolygon(assumeNotNull(tile))[1] AS polygon,
centroid(polygon) AS centroid,
avg_d_kbps,
avg_u_kbps,
avg_lat_ms,
avg_lat_down_ms,
avg_lat_up_ms,
tests,
devices,
year_month
FROM iceberg('https://datasets-documentation.s3.eu-west-3.amazonaws.com/ookla/iceberg/')
0 rows in set. Elapsed: 549.829 sec. Processed 128.01 million rows, 3.75 GB (232.81 thousand rows/s., 6.82 MB/s.)
Peak memory usage: 2.16 GiB.
读者可能会注意到，我们使用了 readWKTPolygon 函数返回的数组中的第一个元素。多边形可以包含多个环，每个环由一个封闭空间组成。第一个条目描述外部空间，后续数组条目指定任何孔。由于我们的多边形是封闭的，没有内部孔，因此我们只有一个条目 - 允许我们将多边形映射到更简单的 Ring 类型。
我们的示例仅包含固定设备的数据。出于好奇，移动设备的数据也存在。
多边形到 SVG
如前所述，我们现在可以使用多种可视化工具来渲染这些数据。但是，这些工具确实会带来开销 - 安装以及更重要的是获得熟悉度和专业知识所需的时间。鉴于 SVG 格式本身就是一种能够渲染多边形的基于矢量的格式，也许我们可以将这些多边形转换为 SVG 以提供数据的二维可视化表示作为地图？
这需要一种将每个多边形转换为SVG 多边形条目 的方法。理想情况下，我们还可以能够设置每个多边形的样式以传达度量标准，例如下载速度。一如既往，社区已经想到了这个问题，ClickHouse 提供了SVG 函数！
这允许我们将 ClickHouse 中的多边形转换为 SVG 多边形，并传递样式。
SELECT SVG([(-160.037841796875, 70.6363054807905), (-160.032348632812, 70.6363054807905), (-160.032348632812, 70.6344840663086), (-160.037841796875, 70.6344840663086), (-160.037841796875, 70.6363054807905)], 'fill: rgb(255, 55,38);') AS svg FORMAT Vertical
Row 1:
──────
svg: <polygon points="-160.038,70.6363 -160.032,70.6363 -160.032,70.6345 -160.038,70.6345 -160.038,70.6363" style="fill: rgb(255, 55,38);"/>
1 row in set. Elapsed: 0.002 sec.
点到像素
虽然以上内容很有希望，但如果要将其可视化为二维地图，我们不能简单地使用经纬度作为坐标。坐标表示三维表面上的位置，即地球表面，而不是二维表面。许多称为投影的技术存在于将经纬度转换为二维表面上的点，例如地图或计算机屏幕。由于不可能在不拉伸或压缩其表面某种方式的情况下使球体变平，因此此过程不可避免地会涉及对地球表面表示的某种妥协或扭曲。
常用于地图的一种流行投影是墨卡托投影。它有许多优点；主要是，它在局部尺度上保持角度和形状，使其成为导航的绝佳选择。它还有一个优点是恒定方位（行进方向）是直线，使导航变得简单。最后，虽然墨卡托投影准确地表示赤道附近的区域，但它会随着人们向两极移动而显着扭曲大小和形状。与赤道地区相比，格陵兰和南极洲等陆地面积看起来比实际大得多。
虽然 ClickHouse 没有内置的墨卡托函数，但我们可以使用另一个 SQL UDF 来解决此问题。在这里，我们将我们的点投影到一个 1024x1024 尺寸的像素空间中。
-- coord (lon, lat) format
CREATE OR REPLACE FUNCTION mercator AS (coord, width, height) -> (
((coord.1) + 180) * (width / 360),
(height / 2) - ((width * ln(tan((pi() / 4) + ((((coord.2) * pi()) / 180) / 2)))) / (2 * pi()))
)
SELECT mercator((-160.037841796875, 70.6363054807905), 1024, 1024) AS pixel_coords
┌─pixel_coords──────────────────┐
│ (56.78125,223.79687511938704) │
└───────────────────────────────┘
1 row in set. Elapsed: 0.001 sec.
要将投影应用于整个多边形，我们可以利用arrayMap 函数。
SELECT polygon, arrayMap(p -> mercator(p, 1024, 1024), polygon) AS pixels
FROM ookla
LIMIT 1
FORMAT Vertical
Row 1:
──────
polygon: [(-51.30615234375,-30.0358110426678),(-51.3006591796875,-30.0358110426678),(-51.3006591796875,-30.0405664305846),(-51.30615234375,-30.0405664305846),(-51.30615234375,-30.0358110426678)]
pixels: [(366.0625,601.6406251005949),(366.078125,601.6406251005949),(366.078125,601.6562501003996),(366.0625,601.6562501003996),(366.0625,601.6406251005949)]
1 row in set. Elapsed: 0.012 sec.
大小限制
有了我们的 SVG 和墨卡托函数，我们就可以生成第一个 SVG 了。以下查询迭代所有 1.28 亿个多边形，并将它们转换为一个 svg 元素。我们还输出包含背景颜色和viewBox 坐标的开始和结束 <svg> 标记。
对于此初始查询，我们希望确认此渲染技术有效，因此，使用固定样式，目前不考虑任何度量标准。我们将此查询的输出重定向到一个文件。请注意，我们使用输出格式CustomSeparated 来生成我们的 svg 开始和结束标记。
clickhouse client --query "SELECT SVG(arrayMap(p -> mercator(p, 1024, 1024), polygon), 'fill: rgb(255, 55,38);')
FROM ookla
FORMAT CustomSeparated
SETTINGS format_custom_result_before_delimiter = '<svg style=\"background-color: rgb(17,17,17);\" viewBox=\"0 0 1024 1024\" xmlns:svg=\"http://www.w3.org/2000/svg\" xmlns=\"http://www.w3.org/2000/svg\">', format_custom_result_after_delimiter = '</svg>', format_custom_escaping_rule = 'Raw'" > world.svg
对于许多人来说，可能可以预测的是，1.28 亿个多边形会产生一个非常大的 svg 文件！
(base) clickhouse@ClickHouse-MacBook-Pro ookla % ls -lh world.svg
-rw-r--r-- 1 dalemcdiarmid wheel 31G 2 Feb 17:04 world.svg
在 31GB 的情况下，大多数机器上的浏览器或工具都无法渲染它，甚至使用ImageMagick 等工具将其转换为压缩格式也是不现实的。如果您有关于如何实现此目标的想法，欢迎提出建议！
使用字典进行汇总
我们的第一个尝试可能是将数据汇总成更大的多边形，表示地理区域。为此，我们需要一种方法来确定我们的较小多边形位于哪个较大的多边形中。这可以通过polygonsWithinCartesian 等函数实现，但对于 1.28 亿个多边形来说，计算成本很高。
幸运的是，ClickHouse 有Polygon 字典 来解决此问题。这允许用户有效地搜索包含特定点的多边形。这需要我们的每个较小多边形都由一个点表示，例如质心 :)
例如，假设我们要将多边形汇总到国家/地区。我们将目标多边形存储在 countries 表中。
CREATE TABLE countries
(
`name` String,
`coordinates` MultiPolygon
)
ENGINE = MergeTree
ORDER BY name
INSERT INTO countries
SELECT name, coordinates
FROM (
SELECT JSONExtractString(JSONExtractString(json, 'properties'), 'ADMIN') as name,
JSONExtractString(JSONExtractRaw(json, 'geometry'), 'type') as type,
if(type == 'Polygon', [JSONExtract(JSONExtractRaw(JSONExtractRaw(json, 'geometry'), 'coordinates'), 'Polygon')], JSONExtract(JSONExtractRaw(JSONExtractRaw(json, 'geometry'), 'coordinates'), 'MultiPolygon')) as coordinates
FROM (SELECT arrayJoin(JSONExtractArrayRaw(json, 'features')) as json
FROM url('https://datasets-documentation.s3.eu-west-3.amazonaws.com/countries/countries.geojson', JSONAsString)))
我们国家/地区的多边形以 geoJSON 表示，可在此处获得。上面的 INSERT SELECT 使用 JSONExtract* 函数从 JSON 中提取这些多边形，并将它们存储在 countries 表中。
我们还有一个关于此主题的视频，涵盖了查询和后续字典。
CREATE DICTIONARY country_polygons
(
`coordinates` MultiPolygon,
`name` String
)
PRIMARY KEY coordinates
SOURCE(CLICKHOUSE(TABLE 'countries'))
LIFETIME(MIN 0 MAX 0)
LAYOUT(POLYGON(STORE_POLYGON_KEY_COLUMN 1))
此字典按多边形的坐标索引，并允许我们使用dictGet 函数检索与特定经纬度关联的多边形。例如，
SELECT dictGet(country_polygons, 'name', (-2.72647, 55.599621)) AS location
┌─location───────┐
│ United Kingdom │
└────────────────┘
1 row in set. Elapsed: 0.012 sec.
请注意这里的检索速度。除了在我们以后的可视化中很有用之外，此字典还允许我们快速了解不同国家/地区的互联网速度。例如，我们可以快速确定南乔治亚岛可能是居住在拥有快速互联网的地方的最佳地点！
SELECT dictGet(country_polygons, 'name', centroid) AS country,
avg(avg_d_kbps) AS download_speed,
sum(devices) AS total_devices
FROM ookla
GROUP BY country
ORDER BY download_speed DESC
LIMIT 10
┌─country──────────────────────────────────┬─────download_speed─┬─total_devices─┐
│ South Georgia and South Sandwich Islands │ 198885 │ 4 │
│ Singapore │ 197264.31481017102 │ 3591446 │
│ Hong Kong S.A.R. │ 166375.523948868 │ 3891226 │
│ China │ 158978.10242450528 │ 53841880 │
│ Thailand │ 155309.49656953348 │ 27506456 │
│ Iceland │ 152290.62935414084 │ 281514 │
│ United States of America │ 149532.95964270137 │ 246966470 │
│ South Korea │ 149115.5647068738 │ 2129690 │
│ United Arab Emirates │ 148936.4285346194 │ 8380900 │
│ Jersey │ 146030.67963711882 │ 100902 │
└──────────────────────────────────────────┴────────────────────┴───────────────┘
10 rows in set. Elapsed: 5.467 sec. Processed 256.01 million rows, 6.66 GB (46.83 million rows/s., 1.22 GB/s.)
Peak memory usage: 129.74 MiB.
最后，正如我们将在下面看到的，能够根据国家/地区的名称有效地检索其多边形也可能会有所帮助。这可以通过使用相同 countries 表的简单字典来实现。
CREATE DICTIONARY default.country_names
(
`coordinates` MultiPolygon,
`name` String
)
PRIMARY KEY name
SOURCE(CLICKHOUSE(TABLE 'countries'))
LIFETIME(MIN 0 MAX 0)
LAYOUT(COMPLEX_KEY_HASHED())
SELECT dictGet(country_names, 'coordinates', 'Yemen') AS coords
[[[(53.308239457000155,12.11839453700017),(53.310267307000146,12.11144367200015)..
1 row in set. Elapsed: 0.007 sec.
有关字典的更多详细信息，我们推荐以下博文。
全球 SVG
有了我们的字典、质心、墨卡托 UDF 和 SVG 函数，我们就可以生成第一个（可用的）SVG 了。
上面的 SVG 为 9.2MB，由以下查询生成，该查询在 11 秒内完成。如密钥所示，红色表示互联网速度慢（下载速率），而绿色表示连接速度更快。
WITH
-- settings determining height, width and color range
1024 AS height, 1024 AS width, [234, 36, 19] AS lc, [0, 156, 31] AS uc,
country_rates AS
(
-- compute the average download speed per country using country_polygons dictionary
SELECT
dictGet(country_polygons, 'name', centroid) AS country,
avg(avg_d_kbps) AS download_speed
FROM ookla
-- exclude 'Antarctica' as polygon isnt fully closed
WHERE country != 'Antarctica'
GROUP BY country
),
(
-- compute min speed over all countries
SELECT min(download_speed)
FROM country_rates
) AS min_speed,
(
-- compute max speed over all countries
SELECT max(download_speed)
FROM country_rates
) AS max_speed,
country_colors AS
(
SELECT
country,
download_speed,
(download_speed - min_speed) / (max_speed - min_speed) AS percent,
-- compute a rgb value based on linear interolation between yellow and red
format('rgb({},{},{})', toUInt8((lc[1]) + (percent * ((uc[1]) - (lc[1])))), toUInt8((lc[2]) + (percent * ((uc[2]) - (lc[2])))), toUInt8((lc[3]) + (percent * ((uc[3]) - (lc[3]))))) AS rgb,
-- get polygon for country name and mercator project every point
arrayMap(p -> arrayMap(r -> arrayMap(x -> mercator(x, width, height), r), p), dictGet(country_names, 'coordinates', country)) AS pixel_poly
FROM country_rates
)
-- convert country polygon to svg
SELECT SVG(pixel_poly, format('fill: {};', rgb))
FROM country_colors
-- send output to file
INTO OUTFILE 'countries.svg'
FORMAT CustomSeparated
SETTINGS format_custom_result_before_delimiter = '<svg style="background-color: rgb(17,17,17);" viewBox="0 0 1024 1024" xmlns:svg="http://www.w3.org/2000/svg" xmlns="http://www.w3.org/2000/svg">', format_custom_result_after_delimiter = '</svg>', format_custom_escaping_rule = 'Raw'
246 rows in set. Elapsed: 1.396 sec. Processed 384.02 million rows, 8.06 GB (275.15 million rows/s., 5.78 GB/s.)
Peak memory usage: 39.54 MiB.
尽管只是组合了上述概念，但这里发生了很多事情。我们对 SQL 进行了注释，希望对您有所帮助，但总而言之
- 我们在第一个 CTE 中声明常量，例如颜色渐变（lc = 黄色，uc = 红色）和图像高度。
- 对于每个国家/地区，我们使用
country_polygons字典计算平均速度，类似于我们之前的查询。 - 对于以上值，我们确定最大和最小互联网速度。这些值用于使用颜色范围之间简单的线性插值 生成每个国家/地区的颜色。
- 对于每个国家/地区，使用其名称，我们使用反向字典
country_names获取多边形。这些多边形很复杂，许多内部孔以环的形式表示。因此，我们需要使用嵌套的arrayMap函数使用墨卡托投影函数对其进行投影。 - 最后，我们使用 SVG 函数和 RGB 值将项目国家/地区多边形转换为 SVG。这使用 INTO OUTFILE 子句输出到本地文件。
使用 h3 索引进行汇总
虽然上述查询证明了从我们的地理数据渲染 SVG 的概念，但在其汇总方面过于粗糙。按国家/地区可视化互联网速度很有趣，但这仅仅是在之前的查询的基础上略有改进，该查询按国家/地区列出了互联网速度。
我们可以尝试使用公共来源 中更大的多边形填充我们的 country_polygons 字典。但是，这似乎有点限制，理想情况下，我们将能够按指定的解析度生成矢量，从而控制我们渲染的多边形的大小。
H3 是 Uber 开发的一种地理索引系统，其中地球表面被划分为一个均匀的六边形单元格网格。此系统是分层的，即顶层上的每个六边形（“父”）可以分成七个更小的六边形（“子”），依此类推。使用索引系统时，我们在 0 到 15 之间定义一个解析度。解析度越高，每个六边形越小，因此覆盖地球表面所需的六边形数量越多。
解析度对多边形数量的影响可以通过使用h3NumHexagons 函数的简单 ClickHouse 查询来说明，该函数返回每个解析度的六边形数量。
SELECT h3NumHexagons(CAST(number, 'UInt8')) AS num_hexagons
FROM numbers(0, 16)
┌─res─┬────num_hexagons─┐
│ 0 │ 122 │
│ 1 │ 842 │
│ 2 │ 5882 │
│ 3 │ 41162 │
│ 4 │ 88122 │
│ 5 │ 2016842 │
│ 6 │ 14117882 │
│ 7 │ 98825162 │
│ 8 │ 691776122 │
│ 9 │ 4842432842 │
│ 10 │ 33897029882 │
│ 11 │ 237279209162 │
│ 12 │ 1660954464122 │
│ 13 │ 11626681248842 │
│ 14 │ 81386768741882 │
│ 15 │ 569707381193162 │
└─────┴─────────────────┘
16 rows in set. Elapsed: 0.001 sec.
ClickHouse 通过一系列函数支持 h3，包括通过 geoToH3 函数将经度和纬度转换为 h3 索引的功能，并指定分辨率。使用此函数，我们可以将我们的中心点转换为 h3 索引，然后在此索引上聚合并计算统计信息，例如平均下载速度。聚合后，我们可以使用函数 h3ToGeoBoundary 反转此过程，将 h3 索引及其六边形转换为多边形。最后，我们使用此墨卡托投影和 SVG 函数生成我们的 SVG。我们将在下面可视化此过程
上面，我们展示了分辨率为 6 的结果，它生成了 14,117,882 个多边形，文件大小为 165MB。这需要 Chrome 在 Macbook M2 上渲染大约 1 分钟。通过修改分辨率，我们可以控制多边形的大小和可视化中的细节级别。如果我们将分辨率降低到 4，我们的文件大小也会降低到 9.2MB，因为我们的可视化变得更粗糙，Chrome 几乎可以立即渲染它。
我们在这里的查询
Our query here:
WITH
1024 AS height, 1024 AS width, 6 AS resolution, 30 AS lc, 144 AS uc,
h3_rates AS
(
SELECT
geoToH3(centroid.1, centroid.2, resolution) AS h3,
avg(avg_d_kbps) AS download_speed
FROM ookla
– filter out as these are sporadic and don’t cover the continent
WHERE dictGet(country_polygons, 'name', centroid) != 'Antarctica'
GROUP BY h3
),
(
SELECT min(download_speed) FROM h3_rates
) AS min_speed,
(
SELECT quantile(0.95)(download_speed) FROM h3_rates
) AS max_speed,
h3_colors AS
(
SELECT
-- sqrt gradient
arrayMap(x -> (x.2), arrayConcat(h3ToGeoBoundary(h3), [h3ToGeoBoundary(h3)[1]])) AS longs,
sqrt(download_speed - min_speed) / sqrt(max_speed - min_speed) AS percent,
-- oklch color with gradient on hue
format('oklch(60% 0.23 {})', toUInt8(lc + (percent * (uc - lc)))) AS oklch,
arrayMap(p -> mercator((p.2, p.1), width, height), h3ToGeoBoundary(h3)) AS pixel_poly
FROM h3_rates
-- filter out points crossing 180 meridian
WHERE arrayExists(i -> (abs((longs[i]) - (longs[i + 1])) > 180), range(1, length(longs) + 1)) = 0
)
SELECT SVG(pixel_poly, format('fill: {};', oklch))
FROM h3_colors
INTO OUTFILE 'h3_countries.svg' TRUNCATE
FORMAT CustomSeparated
SETTINGS format_custom_result_before_delimiter = '<svg style="background-color: rgb(17,17,17);" viewBox="0 0 1024 1024" xmlns:svg="http://www.w3.org/2000/svg" xmlns="http://www.w3.org/2000/svg">', format_custom_result_after_delimiter = '</svg>', format_custom_escaping_rule = 'Raw'
55889 rows in set. Elapsed: 10.101 sec. Processed 384.02 million rows, 8.06 GB (38.02 million rows/s., 798.38 MB/s.)
Peak memory usage: 23.73 MiB.
除了将中心点转换为 h3 索引及其聚合之外，此查询中还有其他一些有趣的部分
- 我们将经度提取到一个单独的数组中，并将第一个条目也附加到末尾。即
arrayMap(x -> (x.2), arrayConcat(h3ToGeoBoundary(h3), [h3ToGeoBoundary(h3)[1]]))。这使我们可以检查连续的经度点对，并过滤掉跨越 180 子午线的点（这些点会渲染成线，破坏我们的视觉效果），即WHERE arrayExists(i -> (abs((longs[i]) - (longs[i + 1])) > 180), range(1, length(longs) + 1)) = 0 - 我们选择使用 Oklch 而不是使用 RGB 颜色渐变。这提供了几个比 RGB（和 HSL）更好的优势，提供了更符合预期的结果和更平滑的渐变，同时更容易生成，因为我们的调色板只需要更改一个维度（色相）。在上面，我们的调色板从红色（低带宽）到绿色（高带宽），在 30 到 144 色相之间使用平方根渐变。这里的
sqrt为我们提供了更均匀的值分布，如果为线性，这些值自然会聚集在一起，从而使我们能够区分区域。我们对 144 使用第 95 个百分位数在这里也有帮助。
栅格化作为替代方案
上面提到的 h3 索引方法对于渲染基于高达约 7 的分辨率的图像非常灵活且实用。高于此的值会生成大多数浏览器无法渲染的文件，例如，7 和 8 会生成 600MB 的文件大小（几乎可以渲染！）和 1.85GB。某些地理区域（例如非洲）的数据缺失在这些更高分辨率下也变得更加明显，因为每个中心点对应的多边形变得更小，即更多的空间是空的。
此方法的替代方法是将图像栅格化，为每个中心点生成一个像素，而不是一个多边形。这种方法可能也更灵活，因为当我们指定所需的分辨率（通过墨卡托）函数时，这将直接映射到用于表示中心点的像素数量。对于相同数量的多边形，数据表示在磁盘上的大小也应该更小。
这种方法的原理很简单。我们使用墨卡托函数投影每个中心点，指定分辨率并对 x 和 y 值进行四舍五入。每个 x 和 y 对代表一个像素，我们可以使用简单的 GROUP BY 为其计算统计信息，例如平均下载速度。以下使用此方法计算每个 x,y 像素的平均下载速度（对于 2048 x 2048 的分辨率），并按行和列顺序返回像素。
SELECT
(y * 2048) + x AS pos,
avg(avg_d_kbps) AS download_speed
FROM ookla
GROUP BY
round(mercator(centroid, 2048, 2048).1) AS x,
round(mercator(centroid, 2048, 2048).2) AS y
ORDER BY pos ASC WITH FILL FROM 0 TO 2048 * 2048
使用上面描述的相同方法，以下查询使用下载速度的最小值和最大值范围计算颜色渐变，为每个像素输出 RGBA 值。下面，我们将 alpha 通道固定为 255，并对我们的颜色渐变使用线性函数，其中第 95 个百分位数作为我们的上限（以帮助压缩我们的渐变）。在这种情况下，我们使用查询创建一个视图 ookla_as_pixels 以简化可视化此数据的下一步。请注意，我们使用 (17, 17, 17, 255) 的背景色来表示下载速度为 null 的像素，即没有数据点，并确保我们的 RGBA 值为 UInt8。
CREATE OR REPLACE VIEW ookla_as_pixels AS
WITH [234, 36, 19] AS lc, [0, 156, 31] AS uc, [17,17,17] As bg,
pixels AS
(
SELECT (y * 2048) + x AS pos,
avg(avg_d_kbps) AS download_speed
FROM ookla
GROUP BY
round(mercator(centroid, 2048, 2048).1) AS x,
round(mercator(centroid, 2048, 2048).2) AS y
ORDER BY
pos WITH FILL FROM 0 TO 2048 * 2048
),
(
SELECT min(download_speed)
FROM pixels
) AS min_speed,
(
SELECT quantile(0.95)(download_speed)
FROM pixels
) AS max_speed,
pixel_colors AS (
SELECT
least(if(isNull(download_speed),0, (download_speed - min_speed) / (max_speed - min_speed)), 1.0) AS percent,
if(isNull(download_speed),bg[1], toUInt8(lc[1] + (percent * (uc[1] - lc[1])))) AS red,
if(isNull(download_speed),bg[2], toUInt8(lc[2] + (percent * (uc[2] - lc[2])))) AS green,
if(isNull(download_speed),bg[3], toUInt8(lc[3] + (percent * (uc[3] - lc[3])))) AS blue
FROM pixels
) SELECT red::UInt8, green::UInt8, blue::UInt8, 255::UInt8 as alpha FROM pixel_colors
将此数据转换为图像最容易使用Canvas 元素和一些简单的 js 来完成。我们实际上每个像素有 4 个字节，总共有 2048 x 2048 x 4 个字节，大约 16MB。使用 putImageData，我们将此数据直接渲染到画布中。
<!doctype html>
<html>
<head>
<meta charset="utf-8">
<link rel="icon" href="favicon.png">
<title>Simple Ookla Visual</title>
</head>
<body>
<div id="error"></div>
<canvas id="canvas" width="2048" height="2048"></canvas>
<script>
async function render(tile) {
const url = `https://127.0.0.1:8123/?user=default&default_format=RowBinary`;
const response = await fetch(url, { method: 'POST', body: 'SELECT * FROM ookla_as_pixels' });
if (!response.ok) {
const text = await response.text();
let err = document.getElementById('error');
err.textContent = text;
err.style.display = 'block';
return;
}
buf = await response.arrayBuffer();
let ctx = tile.getContext('2d');
let image = ctx.createImageData(2048, 2048);
let arr = new Uint8ClampedArray(buf);
for (let i = 0; i < 2048 * 2048 * 4; ++i) {
image.data[i] = arr[i];
}
ctx.putImageData(image, 0, 0, 0, 0, 2048, 2048);
let err = document.getElementById('error');
err.style.display = 'none';
}
const canvas = document.getElementById("canvas");
render(canvas).then((err) => {
if (err) {
let err = document.getElementById('error');
err.textContent = text;
} else {
err.style.display = 'none'
}
})
</script>
</body>
</html>
此代码由Alexey Milovidov 提供，他通过他的ADS-B 数据的交互式可视化和分析将此提升到了一个全新的水平。
上述方法利用了 ClickHouse 的HTTP 接口和以RowBinary 格式返回结果的能力。创建视图大大简化了我们的 SQL 查询，只需 SELECT * FROM ookla_as_pixels 即可。我们的结果图像。
这代表了我们 SVG 方法的一个不错的替代方案，即使我们确实不得不使用比仅仅 SQL 多一点的东西 :)
结论
这篇文章介绍了 Iceberg 数据格式，并讨论了 ClickHouse 中当前的支持状态，并举例说明了如何查询和导入数据。我们探索了以 Iceberg 格式分发的 Ookla 互联网速度数据集，并借此机会探索了 h3 索引、多边形字典、墨卡托投影、颜色插值和使用 UDF 计算中心点——所有这些都旨在仅使用 SQL 来可视化数据。
虽然我们主要关注下载速度和固定设备，但我们希望看到其他人使用类似的方法来探索其他指标。或者，我们欢迎您提出改进上述方法的建议！
