首页 文章

体素邻域索引 - 使用线性索引检测26邻居访问中的“越界”

提问于
浏览
4

好吧,我不知道怎么用 Headers 来描述我的问题,我希望我得到的那个是正确的 .

我有一个矩阵(在下面的例子中是 M ),这是一个3D图像,在这种情况下,由11x11x11体素组成(我认为它只是为了方便,大小只是一个例子) .

在我的代码中,我需要到达一些体素的26个邻居,为此我使用了一些奇特的线性索引:http://www.mathworks.com/matlabcentral/answers/86900-how-to-find-all-neighbours-of-an-element-in-n-dimensional-matrix

问题是,如果 point 位于"boundary"的"boundary"中,则尝试访问某些超出范围的值,这将产生错误 .

要解决这个问题,一个好的方法是在 M 周围创建一个边界,使其在每个维度上都有2个大小,并用零填充,但是 I really would like to avoid changing M ,因为我的代码比示例中的代码复杂得多 .

我找不到任何办法,我有点卡在这里 . 有什么建议吗?

EDIT: @Dan答案有效,但是我想看看是否有使用这种线性索引方法的可能解决方案 .

% Example data
M=round(randn(11,11,11))~=0;

% Fancy way of storing 26 neigh indices for good accesing 
s=size(M);
N=length(s);
[c1{1:N}]=ndgrid(1:3);
c2(1:N)={2};
neigh26=sub2ind(s,c1{:}) - sub2ind(s,c2{:});

point=[5 1 6];

% This will work unless the point is in the boundary (like in this example)
neighbours=M(sub2ind(s,point(1),point(2),point(3))+neigh26)

2 回答

  • 3

    Basic theory: 使用 NaNs 将输入数组扩展为 left-rightup-downone more on each sides of the third dimension . 这将允许我们使用统一的 3x3x3 网格,然后使用那些 NaNs 来检测超出输入数组边界的元素,因此将被丢弃 .

    Code

    %// Initializations
    sz_ext = size(M)+2; %// Get size of padded/extended input 3D array
    M_ext = NaN(sz_ext); %// Initialize extended array
    M_ext(2:end-1,2:end-1,2:end-1) = M; %// Insert values from M into it
    
    %// Important stuff here : Calculate linear offset indices within one 3D slice
    %// then for neighboring 3D slices too
    offset2D = bsxfun(@plus,[-1:1]',[-1:1]*sz_ext(1)); %//'
    offset3D = bsxfun(@plus,offset2D,permute([-1:1]*sz_ext(1)*sz_ext(2),[1 3 2]));
    
    %// Get linear indices for all points
    points_linear_idx = sub2ind(size(M_ext),point(:,1)+1,point(:,2)+1,point(:,3)+1);
    %// Linear indices for all neighboring elements for all points; index into M_ext
    neigh26 = M_ext(bsxfun(@plus,offset3D,permute(points_linear_idx,[4 3 2 1])))
    

    How to use: 因此, 4th 维度中的每个切片表示 27 元素(邻近加上元素本身)为 3x3x3 数组 . 因此, neigh26 将是 3x3x3xN 数组,其中 N 是点数组中的点数 .

    Example: 例如,让我们假设 MPoint 中的一些随机值 -

    M=rand(11,11,11);
    point = [
        1 1 4;
        1 7 1]
    

    在使用这些输入运行早期代码时,我得到这样的东西 -

    neigh26(:,:,1,1) =
           NaN       NaN       NaN
           NaN    0.5859    0.4917
           NaN    0.6733    0.6688
    neigh26(:,:,2,1) =
           NaN       NaN       NaN
           NaN    0.0663    0.5544
           NaN    0.3440    0.3664
    neigh26(:,:,3,1) =
           NaN       NaN       NaN
           NaN    0.3555    0.1257
           NaN    0.4424    0.9577
    neigh26(:,:,1,2) =
       NaN   NaN   NaN
       NaN   NaN   NaN
       NaN   NaN   NaN
    neigh26(:,:,2,2) =
           NaN       NaN       NaN
        0.7708    0.3712    0.2866
        0.7088    0.3743    0.2326
    neigh26(:,:,3,2) =
           NaN       NaN       NaN
        0.4938    0.5051    0.9416
        0.1966    0.0213    0.8036
    
  • 1

    那线性索引的东西是必不可少的吗?因为处理边界条件非常容易,所以你使用下标索引和 minmax 这样:

    p = [5, 1, 6];
    
    neighbourhood = M(max(1,p(1)-1)):min(p(1)+1,end),
                      max(1,p(2)-1)):min(p(2)+1,end),
                      max(1,p(3)-1)):min(p(3)+1,end))
    
    %// Get rid of the point it self (i.e. the center)
    neighbours = neighbourhood([1:13, 15:end])
    

    这样,如果您想要更广泛的社区,您也可以轻松地概括这一点:

    p = [5, 1, 6];
    n = 2;
    neighbourhood = M(max(1,p(1)-n)):min(p(1)+n,end),
                      max(1,p(2)-n)):min(p(2)+n,end),
                      max(1,p(3)-n)):min(p(3)+n,end))
    
    %// Get rid of the point it self (i.e. the center)
    mid = ceil(numel(neigbourhood)/2);
    neighbours = neighbourhood([1:mid-2, mid+1:end])
    

    或者如果您喜欢保持立方体形状,那么可能:

    neighbours = neighbourhood;
    neighbours(mid) = NaN;
    

    如果你想在代码中多次使用它,最好将它重构为只返回索引的m文件函数:

    function ind = getNeighbours(M,p,n)
        M = zeros(size(M));
        M(max(1,p(1)-n)):min(p(1)+n,end), max(1,p(2)-n)):min(p(2)+n,end), max(1,p(3)-n)):min(p(3)+n,end)) = 1;
        M(p(1), p(2), p(3)) = 0;
        ind = find(M);
    end
    

相关问题