具有约 2000 万个采样点和千兆字节数据的交互式大型绘图

问题描述：

我这里遇到了一个问题（我的 RAM）：它无法容纳我想要绘制的数据。我有足够的 HD 空间。有什么解决方案可以避免我的数据集被“遮蔽”吗？

具体来说，我处理数字信号处理，并且必须使用高采样率。我的框架（GNU Radio）以二进制形式保存值（以避免使用太多磁盘空间）。我将其解压缩。之后我需要绘图。我需要可缩放且可交互的绘图。这就是问题所在。

这是否有任何优化潜力，或者是否有其他可以处理更大数据集的软件/编程语言（如 R 等）？实际上，我希望在我的图中有更多的数据。但我没有使用其他软件的经验。GNUplot 失败了，方法与以下类似。我不了解 R（jet）。

import matplotlib.pyplot as plt
import matplotlib.cbook as cbook
import struct

"""
plots a cfile

cfile - IEEE single-precision (4-byte) floats, IQ pairs, binary
txt - index,in-phase,quadrature in plaintext

note: directly plotting with numpy results into shadowed functions
"""

# unpacking the cfile dataset
def unpack_set(input_filename, output_filename):
    index = 0   # index of the samples
    output_filename = open(output_filename, 'wb')

    with open(input_filename, "rb") as f:

        byte = f.read(4)    # read 1. column of the vector

        while byte != "":
        # stored Bit Values
            floati = struct.unpack('f', byte)   # write value of 1. column to a variable
            byte = f.read(4)            # read 2. column of the vector
            floatq = struct.unpack('f', byte)   # write value of 2. column to a variable
            byte = f.read(4)            # next row of the vector and read 1. column
            # delimeter format for matplotlib 
            lines = ["%d," % index, format(floati), ",",  format(floatq), "
"]
            output_filename.writelines(lines)
            index = index + 1
    output_filename.close
    return output_filename.name

# reformats output (precision configuration here)
def format(value):
    return "%.8f" % value            

# start
def main():

    # specify path
    unpacked_file = unpack_set("test01.cfile", "test01.txt")
    # pass file reference to matplotlib
    fname = str(unpacked_file)
    plt.plotfile(fname, cols=(0,1)) # index vs. in-phase

    # optional
    # plt.axes([0, 0.5, 0, 100000]) # for 100k samples
    plt.grid(True)
    plt.title("Signal-Diagram")
    plt.xlabel("Sample")
    plt.ylabel("In-Phase")

    plt.show();

if __name__ == "__main__":
    main()

类似 plt.swap_on_disk() 的东西可以将内容缓存在我的 SSD 上；)

解决方案 1：

对 Ubuntu 上具有 1000 万点散点图基准的开源交互式绘图软件的调查

受到以下描述的用例的启发：https://stats.stackexchange.com/questions/376361/how-to-find-the-sample-points-that-have-statistically-meaningful-large-outlier-r我使用完全相同的输入文件对一些绘图程序进行了基准测试。

基本上，我想要：

• 绘制多维数据的 XY 散点图，希望以 Z 为点颜色

• 使用鼠标从图中交互地选择一些有趣的点

• 查看所选点的所有维度（至少包括 X、Y 和 Z），尝试了解为什么它们是 XY 散点图中的异常值

该问题可以通过以下简化的测试数据来表示：

python -c 'for i in range(10000000): print(f"{i},{i*2},{i*4}")' > 10m1.csv
echo 5000000,20000000,-1 >> 10m1.csv

（~239 MB）的前几行10m1.csv如下所示：

10m1.csv

0,0,0
1,2,4
2,4,8
3,6,12
4,8,16

最后一个，即第 1000 万个，是异常值，如下所示：

5000000,20000000,-1

所以我们基本上有：

• 倾角为 2 的一条线，其上有 1000 万个点

• 加上线外的一个异常点，位于图的顶部中心

类似于：

Y

^
|
|
|       +       +
|
|             +
|
|           +
|
|         +
|
|       +
|
|     +
|
|   +
|
| +
|
+-------------------> X

此基准测试的目标是在图形图上找到点 (5000000,20000000)，然后从中确定第三列的值，这是-1我们的测试。

当我第一次写这个答案时，我使用了以下生成的 10.csv：

python -c 'for i in range(10000000): print(f"{i},{i*2},{i*4}")' > 10m.csv

没有异常值。虽然这测试了性能，但它没有测试选择能力，因此目标是当我有动力时将每个测试迁移到 10m1.csv。

我还做了一个 10 点 + 异常值示例，以防我想评估某些无法处理 10m 点数的工具的可用性：

i=0;
while [ "$i" -lt 10 ]; do
  echo "$i,$((2 * i)),$((4 * i))"; i=$((i + 1));
done > 11.csv
echo 5,20,-1 >> 11.csv

为了增加趣味性，我还准备了一个更大的 10 亿点数据集，以防任何程序可以处理 1000 万个点！CSV 文件有点不稳定，所以我转到了 HDF5：

#!/usr/bin/env python3

import h5py
import numpy

size = 1000000000

with h5py.File('1b.hdf5', 'w') as f:
    x = numpy.arange(size + 1)
    x[size] =  size / 2
    f.create_dataset('x', data=x, dtype='int64')
    y = numpy.arange(size + 1) * 2
    y[size] =  3 * size / 2
    f.create_dataset('y', data=y, dtype='int64')
    z = numpy.arange(size + 1) * 4
    z[size] = -1
    f.create_dataset('z', data=z, dtype='int64')

这将生成一个~23 GiB 文件，类似于10m1.csv包含以下内容：

• 10 亿个点在一条直线上，就像10m.csv

• 图表中心顶部有一个异常点

我还创建了 10m1.csv 的 SQLite 版本，因为这可能是实践中最合理的格式之一，因为它允许很好理解的 SQL 查询、明确的索引控制和二进制数字数据：

f=10m.sqlite
rm -f "$f"
n=10000000
time sqlite3 "$f" 'create table t(x integer, y integer, z integer)'
time sqlite3 "$f" "insert into t select value as x, value * 2 as y, value * 2 as z from generate_series(0, $((n - 1)))"
time sqlite3 "$f" "insert into t values ($((n/2)), $((3*n/2)), -1)"
time sqlite3 "$f" 'create index txy on t(x, y)'

我还用 n = 10 亿运行了该代码来生成一个1b.sqlite。generate_series这是我迄今为止找到的最快的插入方法：使用 Python 将大量数据批量插入 SQLite生成时间：插入耗时 100 秒，创建索引耗时 290 秒。磁盘大小：40 GB。

我通过 (x, y) 进行索引，因为这可能加快查看器工具尝试获取给定 xy 矩形中的所有点的查询速度。生成的 10m1.sqlite 大约为 367 MB，由于索引，它比 CSV 大。

除非本节另有说明，否则测试均在 Ubuntu 18.10 中进行，测试使用的笔记本电脑为 ThinkPad P51，配备 Intel Core i7-7820HQ CPU（4 核/8 线程）、2x Samsung M471A2K43BB1-CRC RAM（2x 16GiB）、NVIDIA Quadro M1200 4GB GDDR5 GPU。

结果摘要

考虑到我非常具体的测试用例，并且我是许多所评测软件的首次用户，这是我观察到的结果：

它能处理 1000 万点吗：

Vaex 2.0.2

https://github.com/vaexio/vaex

安装并让 hello world 运行，如下所示：如何在 Vaex 中进行交互式 2D 散点图缩放/点选择？

我用高达 10 亿个点测试了 vaex，并且它起作用了，太棒了！

它是“Python 脚本优先”，具有很好的可重复性，并允许我轻松地与其他 Python 东西交互。

Jupyter 设置有几个活动部件，但一旦我使用 virtualenv 运行它，它就变得非常棒了。

要在 Jupyter 中加载我们的 CSV 运行：

import vaex
df = vaex.from_csv('10m.csv', names=['x', 'y', 'z'],)
df.plot_widget(df.x, df.y, backend='bqplot')

我们立刻就能看到：

现在，我们可以用鼠标缩放、平移和选择点，更新速度非常快，全部在 10 秒内完成。在这里，我放大了图像以查看一些单独的点，并选择了其中几个（图像上较浅的矩形）：

使用鼠标选择后，这与使用df.select()方法具有完全相同的效果。因此，我们可以通过在 Jupyter 中运行来提取选定的点：

df.to_pandas_df(selection=True)

输出格式如下：

        x       y        z   index
0 4525460 9050920 18101840 4525460
1 4525461 9050922 18101844 4525461
2 4525462 9050924 18101848 4525462
3 4525463 9050926 18101852 4525463
4 4525464 9050928 18101856 4525464
5 4525465 9050930 18101860 4525465
6 4525466 9050932 18101864 4525466

由于 10M 点数运行良好，我决定尝试 1B 点数……而且它也运行良好！

import vaex
df = vaex.open('1b.hdf5')
df.plot_widget(df.x, df.y, backend='bqplot')

为了观察原始图上不可见的异常值，我们可以按照如何在 vaex 交互式 Jupyter bqplot plot_widget 中更改点样式以使单个点更大且可见？并使用：

df.plot_widget(df.x, df.y, f='log', shape=128, backend='bqplot')

生成结果：

选择点后：

我们获得异常值的完整数据：

   x          y           z
0  500000000  1500000000  -1

以下是创作者提供的演示，其中包含更有趣的数据集和更多功能：https://www.youtube.com/watch?v =2Tt0i823-ec&t=770

但不幸的是，没有内置的 sqlite 支持：https://github.com/vaexio/vaex/issues/864

在 Ubuntu 19.04 中测试。

访问 2.13.3

网站：https：//wci.llnl.gov/simulation/computer-codes/visit

许可证：BSD

由劳伦斯利弗莫尔国家实验室开发，该实验室是国家核安全局的一个实验室，所以你可以想象，如果我能让它工作起来，10m 个点对它来说不算什么。（查尔斯·J·默里 (Charles J. Murray) 的《超人：西摩·克雷的故事》(1997)很好地展示了在制造第一颗氢弹时，这些实验室对计算能力的消耗有多么大，因为你不能随意用核武器进行实验，即使你这样做了，你也无法真正测量你想要的东西，因为它爆炸得太快太热：计算机模型是必须的。他们认为，一群带着计算器的物理学家的妻子不会像早期的洛斯阿拉莫斯裂变弹那样成功。当以色列购买了他们的一台电脑时，每个人都立即认为这是用来制造核武器的。）

安装：没有 Debian 软件包，只需从网站下载 Linux 二进制文件即可。无需安装即可运行。另请参阅： https: //askubuntu.com/questions/966901/installing-visit

基于VTK，它是许多高性能图形软件使用的后端库。用 C 编写。

经过 3 个小时的使用 UI 后，我确实让它工作了，它确实解决了我的用例，详情如下：https://stats.stackexchange.com/questions/376361/how-to-find-the-sample-points-that-have-statistically-meaningful-large-outlier-r

以下是本文测试数据的样子：

以及一些精选的缩放：

这是选择窗口：

性能方面，VisIt 非常好：每个图形操作要么只花费很少的时间，要么立即完成。当我不得不等待时，它会显示一条“正在处理”消息，其中包含剩余工作的百分比，并且 GUI 不会冻结。

由于 10m 点工作得很好，我也尝试了 100m 点（一个 2.7G 的 CSV 文件），但不幸的是它崩溃了/进入了一种奇怪的状态，我看到htop4 个 VisIt 线程占用了我所有的 16GiB RAM 并且很可能由于 malloc 失败而死亡。

最初的开始有点痛苦：

• 如果你不是核弹工程师，许多默认设置会让你感觉很糟糕？例如：

  • 默认点大小为 1px（与显示器上的灰尘混淆）

  • 轴的比例从 0.0 到 1.0：如何在访问绘图程序上显示实际轴的数值，而不是从 0.0 到 1.0 的分数？

  • 多窗口设置，选择数据点时会出现令人讨厌的多个弹出窗口

  • 显示您的用户名和绘图日期（通过“控制”>“注释”>“用户信息”删除）

  • 自动定位默认值很糟糕：图例与轴冲突，找不到标题自动化，因此必须手动添加标签并重新定位所有内容

• 功能太多了，所以很难找到你想要的

• 该手册非常有用，但它是一份 386 页的 PDF 巨著，日期不祥地标明为“2005 年 10 月版本 1.5”。我想知道他们是否用它开发了Trinity！这是一个很棒的 Sphinx HTML，在我最初回答这个问题后创建

• 没有 Ubuntu 软件包。但是预构建的二进制文件确实可以工作。

我将这些问题归因于：

• 它已经存在很长时间了，并且使用了一些过时的 GUI 理念

• 你不能直接点击情节元素来改变它们（例如轴、标题等），而且有很多功能，所以很难找到你想要的功能

我也很喜欢 LLNL 基础设施泄漏到该存储库中的方式。例如，请参阅docs/OfficeHours.txt和该目录中的其他文件！我很同情 Brad，他是“星期一早上的人”！哦，答录机的密码是“Kill Ed”，别忘了这一点。

Paraview 5.9.0

网站：https://www.paraview.org/

许可证：BSD

测试平台：Ubuntu 20.10。

安装：

sudo apt install paraview

或者从网站下载预构建版本来获取最新版本。这是我为本次评测所做的，因为 apt 版本只有 5.7.0。我下载了ParaView-5.9.0-MPI-Linux-Python3.8-64bit.tar.gz。

由Kitware和洛斯阿拉莫斯国家实验室开发，后来由桑迪亚国家实验室（其他两个 NNSA 实验室）开发，因此我们再次期待它能够轻松处理数据。它还基于 VTK 并用 C++ 编写，这更加令人期待。

然而，我感到很失望：出于某种原因，10m 点数让 GUI 变得非常缓慢且反应迟钝，无法使用。每当我点击某个东西（例如隐藏线条）时，都要花几十秒的时间。我认为在某个时候它只是出了故障，完全停止了响应。

我对广告中宣传的“我正在工作，请稍等”这一控制功能很满意，但当这种情况发生时 GUI 会冻结吗？不可接受。

htop 显示 Paraview 正在使用 8 个线程和 3GB RAM，因此 CPU 和内存都没有达到最大容量。

从 GUI 角度来看，Paraview 非常漂亮且现代，在没有卡顿的情况下比 VisIt 好得多。

自从10m1.csv杀死它之后，我进行了测试，11.csv看看除了性能之外，我是否能够解决我的问题，答案是肯定的：

• paraview 11.csv

• 从弹出窗口中选择 CSV 阅读器

• 属性 属性 在左侧应用

• 右键单击管道浏览器上的 CSV

• 添加过滤器 > 按字母顺序 > 绘制数据。为什么要绘制过滤器？对于初次使用的用户来说，这不是很直观，相关：paraview：从 csv 文件绘制数据我相信这是一旦你进一步了解过滤器可以做什么的概括后就会明白的事情之一，但仍然如此。

• 属性>应用

• 取消选择“使用索引作为 x 轴”

• X 数组名称：字段 0

• 系列参数删除字段 0 和字段 2

• 选择字段 1 并：

  • 线条样式：无

  • 标记样式: 十字

  • 标记尺寸：根据需要增加或减少

• 图上方的“矩形选择”图标

• 选择异常值（点突出显示）

• 为绘图过滤器添加另一个过滤器：“提取选择”

• 申请

最后！我得到一个仅包含所选异常值的表格，并显示“字段 2”的值作为 -1：

是的，虽然这不完全是轻松的事，但我最终还是做到了。

另一个缺点是，与 VisIt 相比，Paraview 缺少一些功能，例如：

• 我找不到如何根据第三列设置散点图的颜色：如何像 gnuplot 调色板一样在 Paraview 中根据第三列的值为散点图点着色？

Mayavi 4.6.2

网址：https://github.com/enthought/mayavi

开发者：Enthought

安装：

sudo apt-get install libvtk6-dev
python3 -m pip install -u mayavi PyQt5

VTK Python 版本。

Mayavi 似乎非常专注于 3D，我找不到如何在其中绘制 2D 绘图，因此遗憾的是它不适合我的用例。

但是，为了检查性能，我改编了https://docs.enthought.com/mayavi/mayavi/auto/example_scatter_plot.html中的示例，其中包含 1000 万个点，并且它运行良好且没有滞后：

import numpy as np
from tvtk.api import tvtk
from mayavi.scripts import mayavi2

n = 10000000
pd = tvtk.PolyData()
pd.points = np.linspace((1,1,1),(n,n,n),n)
pd.verts = np.arange(n).reshape((-1, 1))
pd.point_data.scalars = np.arange(n)

@mayavi2.standalone
def main():
   from mayavi.sources.vtk_data_source import VTKDataSource
   from mayavi.modules.outline import Outline
   from mayavi.modules.surface import Surface
   mayavi.new_scene()
   d = VTKDataSource()
   d.data = pd
   mayavi.add_source(d)
   mayavi.add_module(Outline())
   s = Surface()
   mayavi.add_module(s)
   s.actor.property.trait_set(representation='p', point_size=1)
main()

输出：

但是我无法放大到足以看清单个点，近处的 3D 平面太远了。也许有办法？

Mayavi 的一个很酷的方面是，开发人员付出了很多努力，让您能够从 Python 脚本中很好地启动和设置 GUI，就像 Matplotlib 和 gnuplot 一样。似乎在 Paraview 中也可以做到这一点，但文档至少没有那么好。

总体来说，它的功能不如 VisIt / Paraview 丰富。例如，我无法直接从 GUI 加载 CSV：如何从 Mayavi GUI 加载 CSV 文件？

Gnuplot 5.2.2

网站：http://www.gnuplot.info/

当我需要快速而粗略地操作时，gnuplot 确实很方便，它总是我首先尝试的东西。

安装：

sudo apt-get install gnuplot

对于非交互式使用，它可以相当好地处理10m点：

#!/usr/bin/env gnuplot
set terminal png size 1024,1024
set output "gnuplot.png"
set key off
set datafile separator ","
plot "10m1.csv" using 1:2:3:3 with labels point

耗时 7 秒：

但如果我尝试与

#!/usr/bin/env gnuplot
set terminal wxt size 1024,1024
set key off
set datafile separator ","
plot "10m.csv" using 1:2:3 palette

和：

gnuplot -persist main.gnuplot

那么初始渲染和缩放感觉太慢了。我甚至看不到矩形选择线！

另请注意，对于我的用例，我需要使用超文本标签，如下所示：

plot "10m.csv" using 1:2:3 with labels hypertext

但标签功能存在性能错误，包括非交互式渲染。但我报告了这个问题，Ethan 在一天内就解决了：https://groups.google.com/forum/#!topic/ comp.graphics.apps.gnuplot/qpL8aJIi9ZE

不过，我必须说，对于异常值选择，有一个合理的解决方法：只需为所有点添加带有行 ID 的标签！如果附近有很多点，您将无法读取标签。但对于您关心的异常值，您可能会！例如，如果我将一个异常值添加到我们的原始数据中：

cp 10m.csv 10m1.csv
printf '2500000,10000000,40000000
' >> 10m1.csv

并将绘图命令修改为：

#!/usr/bin/env gnuplot
set terminal png size 1024,1024
set output "gnuplot.png"
set key off
set datafile separator ","
plot "10.csv" using 1:2:3:3 palette with labels

这显著减慢了绘图速度（上述修复后 40 分钟！！！），但产生了合理的输出：

所以通过一些数据过滤，我们最终可以到达那里。

Matplotlib 1.5.1、numpy 1.11.1、Python 3.6.7

网站：https://matplotlib.org/

当我的 gnuplot 脚本开始变得太疯狂时，我通常会尝试 Matplotlib。

numpy.loadtxt单独这件事花了大约 10 秒钟，所以我知道这不会顺利进行：

#!/usr/bin/env python3

import numpy
import matplotlib.pyplot as plt

x, y, z = numpy.loadtxt('10m.csv', delimiter=',', unpack=True)
plt.figure(figsize=(8, 8), dpi=128)
plt.scatter(x, y, c=z)
# Non-interactive.
#plt.savefig('matplotlib.png')
# Interactive.
plt.show()

首先，非交互式尝试给出了良好的输出，但耗时 3 分 55 秒......

然后，交互式版本在初始渲染和缩放时花费了很长时间。不可用：

请注意此屏幕截图中的缩放选择，它应该立即缩放并消失，但在等待缩放计算时却在屏幕上停留了很长时间！

由于某种原因，我不得不注释掉plt.figure(figsize=(8, 8), dpi=128)交互式版本才能工作，否则它会崩溃：

RuntimeError: In set_size: Could not set the fontsize

散景 1.3.1

https://github.com/bokeh/bokeh

Ubuntu 19.04 安装：

python3 -m pip install bokeh

然后启动 Jupyter：

jupyter notebook

现在，如果我绘制 1m 个点，一切都会完美运行，界面非常棒而且速度很快，包括缩放和悬停信息：

from bokeh.io import output_notebook, show
from bokeh.models import HoverTool
from bokeh.transform import linear_cmap
from bokeh.plotting import figure
from bokeh.models import ColumnDataSource
import numpy as np

output_notebook()
N = 1000000
source = ColumnDataSource(data=dict(
    x=np.random.random(size=N) * N,
    y=np.random.random(size=N) * N,
    z=np.random.random(size=N)
))
hover = HoverTool(tooltips=[("z", "@z")])
p = figure()
p.add_tools(hover)
p.circle(
    'x',
    'y',
    source=source,
    color=linear_cmap('z', 'Viridis256', 0, 1.0),
    size=5
)
show(p)

初始视图：

缩放后：

如果我上升到 10m，它就会阻塞，htop表明 Chromium 有 8 个线程在不可中断 IO 状态下占用了我所有的内存。

这询问了如何引用点：如何引用选定的散景数据点

芘維茲

https://pyviz.org/

TODO 评估。

集成 Bokeh + datashader + 其他工具。

演示 1B 数据点的视频：https://www.youtube.com/watch?v =k27MJJLJNT4 “PyViz：用 30 行 Python 代码可视化 10 亿个数据点的仪表板”，由“Anaconda, Inc.”于 2018-04-17 发布。

西伯恩

https://seaborn.pydata.org/

TODO 评估。

已经有一个关于如何使用 seaborn 可视化至少 5000 万行数据的QA 。

sqlitebrowser 3.12.2

https://github.com/sqlitebrowser/sqlitebrowser

我尝试了这个，看看它是否可以处理 10m1.sqlite，但不幸的是它不能。真可惜！

但它能够直接绘制查询结果，这非常酷。

它看起来是这样的：

在此图中，我将 10m1.sqlite 加载到工具中，然后开始浏览数据。

但它只绘制了为浏览而加载的数据。

您可以单击图表右下角的按钮“加载所有数据并重新绘制图表”，但这会打开一个进度条，每 3 秒增加 1%，所以看起来不太有希望，我放弃了。

在 Ubuntu 23.04 上测试。

SQL 直方图查询

我想知道为什么我无法轻松找到一个使用此工具作为后端的交互式 UI 工具。索引数据库上的 SQL 直方图似乎是处理事情的最合理方法。例如，使用 10 个步骤并忽略空箱：

div=10
x=0
y=0
x2=10000000
y2=20000000

dx=$(((x2 - x) / div))
dy=$(((y2 - y) / div))
time sqlite3 10m1.sqlite --cmd '.mode csv' <<EOF
select
  floor(x/$dx)*$dx as x,
  floor(y/$dy)*$dy as y,
  count(*) as cnt
from t
where
  x >= $x and x < $x2 and
  y >= $y and y < $y2
group by 1, 2
order by 1, 2
EOF

我们达到：

0,0,1000000
1000000,2000000,1000000
2000000,4000000,1000000
3000000,6000000,1000000
4000000,8000000,1000000
5000000,10000000,1000000
5000000,20000000,1
6000000,12000000,1000000
7000000,14000000,1000000
8000000,16000000,1000000
9000000,18000000,1000000

查询需要 6 秒，因此它可以处理 10m 个点，但不能扩展到 1B。

由于我们已经只有​​ 1 个点，因此我们可以在该范围内进行完整列出：

x=5000000
y=20000000
x2=6000000
y2=40000000
time sqlite3 10m1.sqlite --cmd '.mode csv' <<EOF
select *
from t
where
  x >= $x and x < $x2 and
  y >= $y and y < $y2
order by x, y
EOF

这立即给出了最终的期望：

5000000,20000000,-1

因此，此 GUI 可能具有点的最大限制，其中：

• 如果超出限制，则使用热图

• 否则，查询该箱中的完整单个点，并在图中绘制单个点

如何将 SQL 扩展到 1B 行：R 树索引

要扩展到 1B，我们需要 r 树/空间索引，这使我们能够有效地对多个列进行不等式运算。SQLite 有它们，但使用起来有点烦人：

• 所有内容都转换为 32 位浮点数，没有 64 位也没有整数：https://sqlite-users.sqlite.narkive.com/mom9X2r3/sqlite-is-it-possible-to-support-64-bit-value-in-rtree-module

• 无法重复使用现有表中的数据

• 无法存储点，只能存储矩形：https://sqlite-users.sqlite.narkive.com/PMLhd3md/r-tree-storage-optimization-for-points

尽管存在这些限制，我最终还是进行了 100 米点测试，使用重复的 x/y 列，创建时间为：30 分钟，文件大小为：5.9 GB

然后是封顶计数横扫：

max=100
div=10
x=0
y=0
x2=100000000
y2=200000000
dx=$(((x2 - x) / div))
dy=$(((y2 - y) / div))

cx=0
while [ $cx -lt $x2 ]; do
  cy=0
  while [ $cy -lt $y2 ]; do
    printf "$cx,$cy,"
    sqlite3 100mr.sqlite --cmd '.mode csv' <<EOF
select count(x) from (
  select x from t
  where
    x >= $cx and x < $((cx + dx)) and
    y >= $cy and y < $((cy + dy))
  limit $max
)
EOF
    cy=$((cy+dy))
  done
  cx=$((cx+dx))
done

仅用 0.2 秒就完成了，这真是太神奇了。如果不是因为生成时间太长，它可能会扩展到 1B。

不幸的是，PostgreSQL 索引创建速度也没有更快：如何使用 SQLite R 树将简单空间索引移植到 Postgres？尽管它至少支持点而不仅仅是矩形。

在 Ubuntu 23.04 上测试。

解决方案 2：

因此，您的数据并不大，而您在绘制数据时遇到困难这一事实表明工具存在问题。Matplotlib 有很多选项，输出也很好，但它占用大量内存，并且它从根本上假设您的数据很小。但还有其他选择。

举例来说，我使用以下命令生成了一个 20M 数据点文件“bigdata.bin”：

#!/usr/bin/env python
import numpy
import scipy.io.numpyio

npts=20000000
filename='bigdata.bin'

def main():
    data = (numpy.random.uniform(0,1,(npts,3))).astype(numpy.float32)
    data[:,2] = 0.1*data[:,2]+numpy.exp(-((data[:,1]-0.5)**2.)/(0.25**2))
    fd = open(filename,'wb')
    scipy.io.numpyio.fwrite(fd,data.size,data)
    fd.close()

if __name__ == "__main__":
    main()

这会生成一个大小约为 229MB 的文件，这个文件并不是很大；但您表示想要生成更大的文件，因此最终会达到内存限制。

首先让我们集中讨论非交互式绘图。首先要意识到的是，每个点都有字形的矢量图将是一场灾难——对于 2000 万个点中的每一个，其中大多数点无论如何都会重叠，尝试渲染小十字或圆圈或其他东西将是一场灾难，会生成巨大的文件并花费大量时间。我认为这就是默认情况下 matplotlib 沉没的原因。

Gnuplot 可以轻松处理这个问题：

gnuplot> set term png
gnuplot> set output 'foo.png'
gnuplot> plot 'bigdata.bin' binary format="%3float32" using 2:3 with dots

甚至 Matplotlib 也可以表现得更加谨慎（选择栅格后端，并使用像素来标记点）：

#!/usr/bin/env python
import numpy
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt

datatype=[('index',numpy.float32), ('floati',numpy.float32), 
        ('floatq',numpy.float32)]
filename='bigdata.bin'

def main():
    data = numpy.memmap(filename, datatype, 'r') 
    plt.plot(data['floati'],data['floatq'],'r,')
    plt.grid(True)
    plt.title("Signal-Diagram")
    plt.xlabel("Sample")
    plt.ylabel("In-Phase")
    plt.savefig('foo2.png')

if __name__ == "__main__":
    main()  

现在，如果您想要交互，您将不得不将数据放入绘图区，然后动态放大。我不知道有什么 Python 工具可以帮助您立即完成此操作。

另一方面，绘制大数据是一项相当常见的任务，而且有适合这项工作的工具。Paraview 是我个人最喜欢的，VisIt是另一个。它们都主要用于 3D 数据，但 Paraview 尤其适用于 2D 数据，并且非常具有交互性（甚至有一个 Python 脚本界面）。唯一的技巧是将数据写入Paraview可以轻松读取的文件格式。

解决方案 3：

一个较新项目在处理大型数据集方面具有强大潜力：Bokeh ，它正是出于这种考虑而创建的。

实际上，只有与绘图比例相关的数据才会被发送到显示后端。这种方法比 Matplotlib 方法快得多。

解决方案 4：

当然，您可以优化文件的读取：您可以直接将其读入 NumPy 数组，以充分利用 NumPy 的原始速度。您有几个选择。如果 RAM 是一个问题，您可以使用memmap，它将大部分文件保存在磁盘上（而不是 RAM 中）：

# Each data point is a sequence of three 32-bit floats:
data = np.memmap(filename, mode='r', dtype=[('index', 'float32'), ('floati','float32'), ('floatq', 'float32')])

如果 RAM 不是问题，您可以使用fromfile将整个数组放入 RAM 中：

data = np.fromfile(filename, dtype=[('index', 'float32'), ('floati','float32'), ('floatq', 'float32')])

然后可以使用 Matplotlib 的常用plot(*data)功能进行绘图，可能通过另一个解决方案中提出的“放大”方法。

解决方案 5：

我建议做一些复杂一点但应该可行的事情：以不同的分辨率、不同的范围构建图表。

以 Google Earth 为例。如果您以最大级别取消缩放以覆盖整个地球，则分辨率最低。当您放大时，图片会变得更详细，但仅限于您正在放大的区域。

所以基本上对于你的图（它是二维的吗？三维的？我假设它是二维的），我建议你构建一个覆盖整个 [0，n] 范围的大图，分辨率较低，2个较小的图，覆盖 [0，n/2] 和 [n/2 + 1，n]，分辨率是大图的两倍，4个较小的图，覆盖 [0，n/4] ... [3 * n / 4 + 1，n]，分辨率是上面2的两倍，依此类推。

不确定我的解释是否清楚。另外，我不知道是否有任何现有的绘图程序可以处理这种多分辨率图形。

解决方案 6：

我想知道通过加快查找点数是否会有好处？（一段时间以来，我一直对 R*（r 星）树很感兴趣。）

我想知道在这种情况下使用像 r* 树这样的东西是否可以解决问题。（缩小时，树中较高的节点可能包含有关更粗糙、缩小的渲染的信息，更靠近叶子的节点包含单个样本）

甚至可能将树（或最终使用的任何结构）映射到内存中，以保持较高的性能和较低的 RAM 使用率。（你将内存管理任务转移给内核）

希望这有意义......有点闲聊。已经很晚了！