如何使用给定数据绘制 3D 数密度散点图?
how to do 3D number-density scatter plot using given data?
我有一个立方体盒子,每个方向的大小为 300。我有一些数据(包含 X、Y、Z 坐标)代表此框中近 100 万个数据点的分布。我想为它们的数字密度指定一种颜色(在这种情况下,它是一个密集的数量,用于描述可数对象的集中程度数据点)。换句话说,使用颜色来说明哪个部分在数据点方面比其他部分更浓缩。最终图像中颜色条的索引应表示使用该颜色指定的数据点的百分比。
我试图通过将立方体中的整个 space 分成 100 万个更小的立方体(每个立方体在所有方向上的长度均为 3)来做到这一点。通过计算这些立方体中的粒子数量,我将知道它们在盒子中的分布情况以及其中存在的数据点数量。然后我可以为他们指定一种我没有成功计算和指定的颜色。任何建议表示赞赏。
%reading the files
[FileName,PathName,FilterIndex] = uigetfile('H:\*.txt','MultiSelect','on');
numfiles = size(FileName,2);%('C:\final1.txt');
j=1;
X=linspace(0,300,100);
for ii = 1:numfiles
FileName{ii}
entirefile = fullfile(PathName,FileName{ii});
a = importdata(entirefile);
x = a(:,2);
y = a(:,3);
z = a(:,4);
%% I DON'T KNOW HOW TO CREAT THIS LOOP TO COUNT FOR THE NUMBER OF PARTICLES WITHIN EACH DEFINED CUBE %%
for jj = 2:size(X,2)
%for kk=1:m
if x(:)<X(jj) & y(:)<X(jj) & z(:)<X(jj)
x;
end
%end
end
h=figure(j);
scatter3(x, y, z, 'filled', 'MarkerSize', 20);
cb = colorbar();
cb.Label.String = 'Probability density estimate';
end
我需要得到类似下图的结果。但我需要每种颜色指定的数据点百分比。提前致谢。
这里是link一个采样数据。
这是一种计算点云 3D 密度的方法。尽管我担心您提供的示例数据不会产生与您的示例图像相同的 3D 分布。
为了计算密度,方法分为几个步骤:
- 计算
[X,Y]
平面中的二维密度:计算每个 (x,y)
bin 中的点数。但是,在这个阶段,这个点数包含给定 bin 的所有 Z
列。
- 对于每个 非空
(x,y)
bin,计算沿 Z
列的分布。我们现在有落在每个 (x,y,z)
bin 中的点数。计算 density/percentage 只需将每个计数除以总点数即可。
- 现在,对于每个 非空
(x,y,z)
bin,我们确定属于该 bin 的点的线性索引。然后,我们将 bin 值(颜色、百分比或与此 bin 关联的任何值)分配给所有已识别的点。
- 显示结果。
在代码中,它是这样的:
%% Import sample data
entirefile = '1565015520323.txt' ;
a = importdata(entirefile);
x = a(:,1);
y = a(:,2);
z = a(:,3);
npt = numel(x) ; % Total Number of Points
%% Define domain and grid parameters
nbins = 100 ;
maxDim = 300 ;
binEdges = linspace(0,maxDim,nbins+1) ;
%% Count density
% we start counting density along in the [X,Y] plane (Z axis aglomerated)
[Nz,binEdges,~,binX,binY] = histcounts2(y,x,binEdges,binEdges) ;
% preallocate 3D containers
N3d = zeros(nbins,nbins,nbins) ; % 3D matrix containing the counts
Npc = zeros(nbins,nbins,nbins) ; % 3D matrix containing the percentages
colorpc = zeros(npt,1) ; % 1D vector containing the percentages
% we do not want to loop on every block of the domain because:
% - depending on the grid size there can be many
% - a large number of them can be empty
% So we first find the [X,Y] blocks which are not empty, we'll only loop on
% these blocks.
validbins = find(Nz) ; % find the indices of non-empty blocks
[xbins,ybins] = ind2sub([nbins,nbins],validbins) ; % convert linear indices to 2d indices
nv = numel(xbins) ; % number of block to process
% Now for each [X,Y] block, we get the distribution over a [Z] column and
% assign the results to the full 3D matrices
for k=1:nv
% this block coordinates
xbin = xbins(k) ;
ybin = ybins(k) ;
% find linear indices of the `x` and `y` values which are located into this block
idx = find( binX==xbin & binY==ybin ) ;
% make a subset with the corresponding 'z' value
subZ = z(idx) ;
% find the distribution and assign to 3D matrices
[Nz,~,zbins] = histcounts( subZ , binEdges ) ;
N3d(xbin,ybin,:) = Nz ; % total counts for this block
Npc(xbin,ybin,:) = Nz ./ npt ; % density % for this block
% Now we have to assign this value (color or percentage) to all the points
% which were found in the blocks
vzbins = find(Nz) ;
for kz=1:numel(vzbins)
thisColorpc = Nz(vzbins(kz)) ./ npt * 100 ;
idz = find( zbins==vzbins(kz) ) ;
idx3d = idx(idz) ;
colorpc(idx3d) = thisColorpc ;
end
end
assert( sum(sum(sum(N3d))) == npt ) % double check we counted everything
%% Display final result
h=figure;
hs=scatter3(x, y, z, 3 , colorpc ,'filled' );
xlabel('X'),ylabel('Y'),zlabel('Z')
cb = colorbar ;
cb.Label.String = 'Probability density estimate';
正如我一开始所说,结果与您的示例图像略有不同。此样本集产生以下分布:
如果你想 "double check" 结果不是垃圾,你可以查看每个轴上的二维密度结果,并检查它是否与你的点的表观分布相匹配:
%% Verify on 3 axis:
Nz = histcounts2(y,x,binEdges,binEdges) ./ npt *100 ;
Nx = histcounts2(z,y,binEdges,binEdges) ./ npt *100 ;
Ny = histcounts2(x,z,binEdges,binEdges) ./ npt *100 ;
figure
ax1=subplot(1,3,1) ; bz = plotDensity(Nz,ax1) ; xlabel('X'),ylabel('Y') ;
ax2=subplot(1,3,2) ; bx = plotDensity(Nx,ax2) ; xlabel('Y'),ylabel('Z') ;
ax3=subplot(1,3,3) ; by = plotDensity(Ny,ax3) ; xlabel('Z'),ylabel('X') ;
点击图片放大:
plotDensity.m
的代码:
function hp = plotDensity(Ndist,hax)
if nargin<2 ; hax = axes ; end
hp = bar3(Ndist,'Parent',hax) ;
for k = 1:length(hp)
zdata = hp(k).ZData;
hp(k).CData = zdata;
hp(k).FaceColor = 'interp';
end
shading interp
我有一个立方体盒子,每个方向的大小为 300。我有一些数据(包含 X、Y、Z 坐标)代表此框中近 100 万个数据点的分布。我想为它们的数字密度指定一种颜色(在这种情况下,它是一个密集的数量,用于描述可数对象的集中程度数据点)。换句话说,使用颜色来说明哪个部分在数据点方面比其他部分更浓缩。最终图像中颜色条的索引应表示使用该颜色指定的数据点的百分比。
我试图通过将立方体中的整个 space 分成 100 万个更小的立方体(每个立方体在所有方向上的长度均为 3)来做到这一点。通过计算这些立方体中的粒子数量,我将知道它们在盒子中的分布情况以及其中存在的数据点数量。然后我可以为他们指定一种我没有成功计算和指定的颜色。任何建议表示赞赏。
%reading the files
[FileName,PathName,FilterIndex] = uigetfile('H:\*.txt','MultiSelect','on');
numfiles = size(FileName,2);%('C:\final1.txt');
j=1;
X=linspace(0,300,100);
for ii = 1:numfiles
FileName{ii}
entirefile = fullfile(PathName,FileName{ii});
a = importdata(entirefile);
x = a(:,2);
y = a(:,3);
z = a(:,4);
%% I DON'T KNOW HOW TO CREAT THIS LOOP TO COUNT FOR THE NUMBER OF PARTICLES WITHIN EACH DEFINED CUBE %%
for jj = 2:size(X,2)
%for kk=1:m
if x(:)<X(jj) & y(:)<X(jj) & z(:)<X(jj)
x;
end
%end
end
h=figure(j);
scatter3(x, y, z, 'filled', 'MarkerSize', 20);
cb = colorbar();
cb.Label.String = 'Probability density estimate';
end
我需要得到类似下图的结果。但我需要每种颜色指定的数据点百分比。提前致谢。
这里是link一个采样数据。
这是一种计算点云 3D 密度的方法。尽管我担心您提供的示例数据不会产生与您的示例图像相同的 3D 分布。
为了计算密度,方法分为几个步骤:
- 计算
[X,Y]
平面中的二维密度:计算每个(x,y)
bin 中的点数。但是,在这个阶段,这个点数包含给定 bin 的所有Z
列。 - 对于每个 非空
(x,y)
bin,计算沿Z
列的分布。我们现在有落在每个(x,y,z)
bin 中的点数。计算 density/percentage 只需将每个计数除以总点数即可。 - 现在,对于每个 非空
(x,y,z)
bin,我们确定属于该 bin 的点的线性索引。然后,我们将 bin 值(颜色、百分比或与此 bin 关联的任何值)分配给所有已识别的点。 - 显示结果。
在代码中,它是这样的:
%% Import sample data
entirefile = '1565015520323.txt' ;
a = importdata(entirefile);
x = a(:,1);
y = a(:,2);
z = a(:,3);
npt = numel(x) ; % Total Number of Points
%% Define domain and grid parameters
nbins = 100 ;
maxDim = 300 ;
binEdges = linspace(0,maxDim,nbins+1) ;
%% Count density
% we start counting density along in the [X,Y] plane (Z axis aglomerated)
[Nz,binEdges,~,binX,binY] = histcounts2(y,x,binEdges,binEdges) ;
% preallocate 3D containers
N3d = zeros(nbins,nbins,nbins) ; % 3D matrix containing the counts
Npc = zeros(nbins,nbins,nbins) ; % 3D matrix containing the percentages
colorpc = zeros(npt,1) ; % 1D vector containing the percentages
% we do not want to loop on every block of the domain because:
% - depending on the grid size there can be many
% - a large number of them can be empty
% So we first find the [X,Y] blocks which are not empty, we'll only loop on
% these blocks.
validbins = find(Nz) ; % find the indices of non-empty blocks
[xbins,ybins] = ind2sub([nbins,nbins],validbins) ; % convert linear indices to 2d indices
nv = numel(xbins) ; % number of block to process
% Now for each [X,Y] block, we get the distribution over a [Z] column and
% assign the results to the full 3D matrices
for k=1:nv
% this block coordinates
xbin = xbins(k) ;
ybin = ybins(k) ;
% find linear indices of the `x` and `y` values which are located into this block
idx = find( binX==xbin & binY==ybin ) ;
% make a subset with the corresponding 'z' value
subZ = z(idx) ;
% find the distribution and assign to 3D matrices
[Nz,~,zbins] = histcounts( subZ , binEdges ) ;
N3d(xbin,ybin,:) = Nz ; % total counts for this block
Npc(xbin,ybin,:) = Nz ./ npt ; % density % for this block
% Now we have to assign this value (color or percentage) to all the points
% which were found in the blocks
vzbins = find(Nz) ;
for kz=1:numel(vzbins)
thisColorpc = Nz(vzbins(kz)) ./ npt * 100 ;
idz = find( zbins==vzbins(kz) ) ;
idx3d = idx(idz) ;
colorpc(idx3d) = thisColorpc ;
end
end
assert( sum(sum(sum(N3d))) == npt ) % double check we counted everything
%% Display final result
h=figure;
hs=scatter3(x, y, z, 3 , colorpc ,'filled' );
xlabel('X'),ylabel('Y'),zlabel('Z')
cb = colorbar ;
cb.Label.String = 'Probability density estimate';
正如我一开始所说,结果与您的示例图像略有不同。此样本集产生以下分布:
如果你想 "double check" 结果不是垃圾,你可以查看每个轴上的二维密度结果,并检查它是否与你的点的表观分布相匹配:
%% Verify on 3 axis:
Nz = histcounts2(y,x,binEdges,binEdges) ./ npt *100 ;
Nx = histcounts2(z,y,binEdges,binEdges) ./ npt *100 ;
Ny = histcounts2(x,z,binEdges,binEdges) ./ npt *100 ;
figure
ax1=subplot(1,3,1) ; bz = plotDensity(Nz,ax1) ; xlabel('X'),ylabel('Y') ;
ax2=subplot(1,3,2) ; bx = plotDensity(Nx,ax2) ; xlabel('Y'),ylabel('Z') ;
ax3=subplot(1,3,3) ; by = plotDensity(Ny,ax3) ; xlabel('Z'),ylabel('X') ;
点击图片放大:
plotDensity.m
的代码:
function hp = plotDensity(Ndist,hax)
if nargin<2 ; hax = axes ; end
hp = bar3(Ndist,'Parent',hax) ;
for k = 1:length(hp)
zdata = hp(k).ZData;
hp(k).CData = zdata;
hp(k).FaceColor = 'interp';
end
shading interp