如何绘制(x,y,z)点并显示它们的密度？

Question

我的数据文件是一组围绕轴原点定位的（x，y，z）点。它们代表了一种测量失败的点。这些点在此链接中。

Gnuplot可以将它们绘制出来。

set encoding iso_8859_1
set term postscript eps enhanced color size 4.7in,4in
set xlabel "X"
set ylabel "Y"
set zlabel "Z"
set output "test_gp.eps"
set style line 1 lc rgb '#0060ad' pt 7 ps 0.5 lt 1 lw 0.5 # --- blue
set style fill  transparent solid 0.15 noborder
splot "data.dat" u 1:2:3 w p ls 1 title "P_{error}"

通过生成的图

问题在于图示无法显示错误最可能发生的位置。 如果可能的话，我想展示点密度的变化情况。

我不知道是否可以改变点的颜色或使它们透明，并尝试用最密集的点来表示位置。

如果无法表示3D点密度，则另一种方法是在三个平面（x=0，y，z）、（x，y=0，z）和（x，y，z=0）中投影点密度。

问候

编辑： 我可以用不同的颜色绘制实验成功（0）和错误（1）的点位置。此link文件有第四列包含所有数据样本（0和1）。

splot图

splot "data_all.dat" u 1:2:3:4 w points ls 1 palette title "P_{error}"

是

但是这个图形并没有显示点的密度。例如，在Mathematica中，这些数据样本的密度图为：

我该如何使用 Gnuplot 绘制密度图？很可能Mathematica正在对中间的点进行插值，并赋予它们0到1之间的值，但我不知道如何用Gnuplot实现。

- user1993416

---

使用您的代码，数据点不会透明。相反，例如将线条颜色设置为rgb '#cc0060ad'，根据方案#aarrggbb。 aa将是alpha通道。关于计算密度，我不知道Mathematica如何做或如何在gnuplot中完成，但基本上它是在某个小体积中计算点数并分配颜色。这让我想起了“平滑频率”和分箱，但是在3D中。 - theozh

@theozh 我认为在Gnuplot中这是不可能的。 我可以通过制作矩阵来制作2D密度图。矩阵的每个位置都是一个xy bin（2D bin），其值是落入该2D bin中的xy数字的计数。在链接https://stackoverflow.com/questions/46159983/count-number-of-points-in-2d-bins?noredirect=1#comment79374603_46159983中，该解决方案展示了如何从xy数据样本生成那些数字的矩阵。有了矩阵，就可以像http://gnuplot.sourceforge.net/demo/heatmaps.html的第一个示例一样制作热力图。 - user1993416

@theozh 这里的问题类似，但是生成3D箱更加困难，并且没有使用Gnuplot制作3D热图的示例。因此，我认为Gnuplot没有实现这种方法。一种解决方法是为三个坐标平面（x，y，z=0），（x，y=0，z），（x=0，y，z）制作热图。然后在XYZ坐标系中将这三个热图一起表示（如投影）。我可以生成这三个热图（虽然我还没有尝试过，但我相信是可能的），但我不知道如何将它们三个连接在坐标系中。 - user1993416

Answer 1

@user1993416，我猜你可以用gnuplot做些事情。您可能需要调整参数Radius来确定在该半径内围绕某个点的点数，并计算密度。在我的8年老电脑上，大约需要2-3分钟处理1000个点。

更新：重新访问这篇旧文章后，我注意到测试数据创建中存在不必要的循环，通过仅使用数组而不是写入表格（似乎比数组慢得多）可以使脚本运行更快。现在，处理500、1000和2000个点大约需要2、7和26秒。因此，运行时间为O(n^2)。

脚本：（适用于gnuplot>=5.2.0）

### 3D density plot
reset session
set term wxt
TimeStart = time(0.0)

# create some random test data
set table $Data
    set samples 1000
    plot '+' u (invnorm(rand(0))):(invnorm(rand(0))):(invnorm(rand(0))) with table
unset table

# put the datafile/dataset into arrays
stats $Data u 0 nooutput
N = STATS_records
array X[N]
array Y[N]
array Z[N]
array C[N]
stats $Data u (X[$0+1]=$1, Y[$0+1]=$2, Z[$0+1]=$3) nooutput

# look at each datapoint and its sourrounding
Radius = 2
V      = pi*4/3*Radius**3   # sphere volume
p0     = 0;  dp = 10        # progress steps
do for [i=1:N] {
    if (p=real(i)/N*100, p>=p0 ? p0=(floor(p/dp)+1)*dp : 0) {
        print sprintf("Progress: %.0f%%", p)
    } 
    C[i] = 0
    stats $Data u ( C[i] = C[i] + \
          (sqrt((X[i]-$1)**2 + (Y[i]-$2)**2 + (Z[i]-$3)**2) <= Radius )) nooutput
}
print sprintf("Time elapsed: %.3f sec",time(0.0)-TimeStart)

set key noautotitle
set xyplane relative 0
set view equal xyz
set view 65,45,1.2
set palette rgb 33,13,10

splot X u (X[$0+1]):(Y[$0+1]):(Z[$0+1]):(C[$0+1]/V) w p ps 1 pt 7 lc palette z

set terminal gif animate delay 20
set output "SO53659762.gif"
    a0     = 40
    a1     = 60
    Frames = 24
    do for [i=0:Frames] {
        set view 65,(sin(2*pi*i/(Frames+1))*(a1-a0)+a0)
        replot
    }
set output
### end of script

结果：（输出到gif终端）

- theozh

---

谢谢。我需要用我的样本来测试你的解决方案。我有一些XYZ点，其中一些是错误的（'1'），其他数据是成功的（'0'）。我想到我应该只测试一类数据（例如错误），看看输出是否提供了错误可能发生的位置的信息。然而，我认为只考虑一种类型的数据不会给我所有的信息。我会尝试一下。 - user1993416

我的Gnuplot 5.2没有“wxt”终端。我应该尝试哪个终端？我想测试gif，但最终，我只想在报告中使用其中一个视图。 - user1993416

此外，是否可以在三维表示中添加XY平面上点的热力图？ - user1993416

你应该能够使用标准终端而不是wxt。我插入这个只是因为我第一次使用splot...时使用了我的标准终端并将动画输出到了gif终端。我相信两个数据类可以被处理，还有在XY平面上的热力图...需要更多的编码。 - theozh

我想要添加的新编码与我在这里询问的不同。因此，我将开一个新问题来询问编码扩展。如果允许的话，我想邀请您来看一下。谢谢。 - user1993416