【问题标题】:Gridding scattered data without interpolation网格化分散数据而不进行插值
【发布时间】:2019-02-11 20:38:27
【问题描述】:

我有三个向量,XYZ,它们表示在矩阵 V(即V = f(X,Y,Z))中找到的一些数据的 3d 坐标。以下是该数据的一些特征(假设/断言,如果您愿意的话):

  • XY 的长度相同,但(通常)不同于 Z 的长度。

    numel(X) == numel(Y);
    
  • X 中唯一元素的数量通常与 Y 中的唯一元素数量不同。

  • Z 中的所有值都是唯一的。

    numel(unique(Z)) == numel(Z);
    
  • V 的大小为 [numel(X), numel(Z)]

过去,我没有区分XY,而是使用单个索引来引用它们,这样可以轻松绘制具有“XY组合编号”的图表(即1:numel(X))在一个轴上,Z 在另一个轴上,所以V 的大小效果很好。

目前,我想分别查看XY 的效果,这就是为什么我想使用scatter3contourslice 的混合在3d 中可视化它。

分散部分很简单,因为我可以 repmat XY 沿着它们的单一维度 numel(Z) 次,对于 Z 也是如此(使用 numel(X))。这导致:

至于轮廓,这些需要以 3d 数组的形式提供数据——这意味着我必须将点放置在结构化网格中。使用meshgrid 创建网格坐标XXYYZZ 很简单:

[ XX, YY, ZZ ] = meshgrid( unique(X), unique(Y), unique(Z) );

我正在努力创建 3d 数组 VV

从上图我们可以看出,大部分卷不包含点 - 我非常希望保持这种状态。换句话说,理想的VV应该只包含与原始数据集对应的点,其余空间应该包含NaNs。

griddatainterp3 之类的函数执行插值,这将“填充”点云内部的漏洞 - 这是非常不希望的。我认为这里可以使用索引来使用来自V 的值填充VV,但我想不出办法。

我的问题是:如何生成不包含任何插值数据的VV

这是一个最小的例子:

%% Generate some data:
X = randi(10,100,1);
Y = randi(15,100,1) - 5;
Z = 1:50;
V = X./Y.*Z;

%% Scatter plot:
nXY = numel(X); nZ = numel(Z);
figure();
scatter3( reshape( repmat(X,[1, nZ]),  [], 1), ...
          reshape( repmat(Y,[1, nZ]),  [], 1), ...
          reshape( repmat(Z,[nXY, 1]), [], 1), ...
          [], V(:), '.');

%% Contours:
% Create the 3d grid:
[XX, YY, ZZ] = meshgrid( unique(X), unique(Y), unique(Z) );

% Preallocate VV:
VV = NaN(size(XX));

% Populate VV: <--------------------------------------------- Help needed with this stage
ind = randperm( numel(XX), numel(V) ); % PLACEHOLDER 
VV(ind) = V;

% Plot:
hold on; contourslice(XX, YY, ZZ, VV, X(2), Y(3), Z(10) );

【问题讨论】:

    标签: matlab matrix multidimensional-array indexing vectorization


    【解决方案1】:

    我相信通过将所有网格化索引与所有线性数据点匹配可能是一种过度杀伤力的方法。为此,我们需要注入一些维度,以便将 3d 数组 XX 等与 2d 数组 X 等元素进行比较:

    Xbc = reshape(X, [1,1,1,size(X)]);
    Ybc = reshape(Y, [1,1,1,size(Y)]);
    Zbc = reshape(Z, [1,1,1,size(Z)]);
    

    这些数组被重新整形,以便它们使用大小为[N,M,K] 的数组XX 等进行广播(“bc”代表广播)。所以元素比较现在有效:

    match = reshape((XX == Xbc) & (YY == Ybc) & (ZZ == Zbc), [size(XX), numel(V)]);
    

    如果V 的大小为[P,Q],则此逻辑数组的大小为[N,M,K,P,Q]。它包含任意数量的trues:

    >> sum(match(:)) == numel(V)
    
    ans =
    
      logical
    
       1
    

    所以现在我们需要沿着前三个维度挑选出对应的索引,并将它们与V 的正确元素配对。我们需要一些线性到多指数的肘部润滑脂:

    [ii,jj,kk,ll] = ind2sub(size(match), find(match));
    

    现在左边的所有数组的大小都是[numel(V), 1];前三个为您提供XX 等的索引,最后一个为您提供V 的索引。

    V_inds = ll;
    VV_inds = sub2ind(size(VV), ii, jj, kk);
    
    VV(VV_inds) = V(V_inds);
    

    现在,由于某种原因,我在结果中只看到 5000 个中的 3750 个唯一索引:

    >> numel(VV_inds)           
    
    ans =
    
            5000
    
    >> numel(unique(VV_inds))
    
    ans =
    
            3750
    

    除了XY 的值重复导致您的一些原始数据点重叠之外,我似乎找不到任何其他原因,因此您实际上无法将它们表示为唯一点的 3d 网格(因为某些 3d 点包含多个数据点)。我相信以下证明了这一点:

    >> size(unique([X,Y], 'rows'))
    
    ans =
    
        75     2
    
    >> size([X,Y])
    
    ans =
    
       100     2
    

    有 100 个 (x,y) 对,但只有 75 个唯一的。无论您如何将这些与正交z 点集结合起来,最终都会出现重复点。因此,您要么必须剔除原始数据中的冗余,要么需要找到另一种表示形式(或为冲突值取平均值)。


    我想我也有一个更高效的版本,使用unique 在运行期间生成的索引。请注意,我假设您使用meshgrid 而不是ndgrid 来生成网格,因此生成的数组(以及VV)的尺寸对应于X、@分别为 987654349@ 和 Z

    % take the indices
    [uX, ~, iX] = unique(X);
    [uY, ~, iY] = unique(Y);
    [uZ, ~, iZ] = unique(Z);
    
    % generate mesh and allocate result
    [XX, YY, ZZ] = ndgrid(uX, uY, uZ);
    VV = NaN(size(XX));
    
    % switch from `iX`, `iY` and `iZ` to a 2d mesh of size `[P,Q]` where `iX` and `iY` are of size `[P,1]` and `iZ` is of size `[Q,1`]:
    % a.k.a. lazy repmat
    iXbig = iX + 0*iZ.';
    iYbig = iY + 0*iZ.';
    iZbig = iZ.' + 0*iX;
    
    % turn 3d indices into linear index into VV
    VV_inds = sub2ind(size(VV), iXbig, iYbig, iZbig);
    
    % profit
    VV(VV_inds) = V;
    

    【讨论】:

      猜你喜欢
      • 1970-01-01
      • 2014-03-26
      • 1970-01-01
      • 1970-01-01
      • 2012-02-19
      • 1970-01-01
      • 1970-01-01
      • 2015-07-21
      • 2019-10-27
      相关资源
      最近更新 更多