【问题标题】:General method to find submatrix in matlab matrixmatlab矩阵中求子矩阵的一般方法
【发布时间】:2013-09-20 18:42:44
【问题描述】:

我正在寻找一种在更大矩阵(任意维数)中找到矩阵(模式)的“好”方法。

例子:

total = rand(3,4,5);
sub = total(2:3,1:3,3:4);

现在我希望这发生:

loc = matrixFind(total, sub)

在这种情况下,loc 应该变为 [2 1 3]

现在我只对找到一个点(如果存在)感兴趣,并不担心舍入问题。可以假设sub“适合”total


这是我如何在 3 个维度上做到这一点,但感觉好像有更好的方法:

total = rand(3,4,5);
sub = total(2:3,1:3,3:4);
loc = [];
for x = 1:size(total,1)-size(sub,1)+1
    for y = 1:size(total,2)-size(sub,2)+1
        for z = 1:size(total,3)-size(sub,3)+1
            block = total(x:x+size(sub,1)-1,y:y+size(sub,2)-1,z:z+size(sub,3)-1);
            if isequal(sub,block)
                loc = [x y z]
            end
        end
    end
end

我希望为任意数量的维度找到一个可行的解决方案。

【问题讨论】:

  • 不确定它是否有助于解决问题,但可以将ndims(sub) 假定为等于ndims(total) 吗?
  • 请注意:对于仅二维的情况,Matlab File Exchange 中的函数findsubmat 有很好的实现思路(和代码 cmets)。
  • 与@ojdo 的第一个问题有点相关:我认为你必须更精确地定义你想要的输出,以防ndim(sub) < ndims(total)。我想在这种情况下,您的 for 循环可能找不到所有可能的解决方案。要求 ndims(sub) == ndims(total) 可能会简化一些事情。
  • @BasSwinckels 如果我尝试sub = total(2,1,3) 似乎没问题。 ndims(sub) 是 2,ndims(total) 是 3。即使维度是单例,定义 size(sub,3) 也会有所帮助。无论哪种方式,我的代码中都可能存在错误,但我希望目标很明确:我只希望它在 subs 'fits' in total 时工作。
  • @DennisJaheruddin 一个问题是 Matlab 没有实现真正的 N 维数组,而是做了一些奇怪的事情:它省略了尾随单维数(size(zeros(3,2,1)) == [3,2],而size(zeros(1,2,3)) == [1,2,3])和标量始终是二维的 (size(1) == [1,1])。但我的主要问题是与哪个维度进行比较,以防sub 的维度小于total

标签: matlab matrix find


【解决方案1】:

这是低性能,但(据说)任意维度的函数。它使用find 创建total 中潜在匹配位置的(线性)索引列表,然后检查total 的适当大小的子块是否与sub 匹配。

function loc = matrixFind(total, sub)
%matrixFind find position of array in another array

    % initialize result
    loc = [];

    % pre-check: do all elements of sub exist in total?
    elements_in_both = intersect(sub(:), total(:));
    if numel(elements_in_both) < numel(unique(sub))
        % if not, return nothing
        return
    end

    % select a pivot element
    % Improvement: use least common element in total for less iterations
    pivot_element = sub(1);

    % determine linear index of all occurences of pivot_elemnent in total
    starting_positions = find(total == pivot_element);

    % prepare cell arrays for variable length subscript vectors
    [subscripts, subscript_ranges] = deal(cell([1, ndims(total)]));


    for k = 1:length(starting_positions)
        % fill subscript vector for starting position
        [subscripts{:}] = ind2sub(size(total), starting_positions(k));

        % add offsets according to size of sub per dimension
        for m = 1:length(subscripts)
            subscript_ranges{m} = subscripts{m}:subscripts{m} + size(sub, m) - 1;
        end

        % is subblock of total equal to sub
        if isequal(total(subscript_ranges{:}), sub)
            loc = [loc; cell2mat(subscripts)]; %#ok<AGROW>
        end
    end
end

【讨论】:

  • if numel(elements_in_both) &lt; numel(sub) 如果 sub 中有双元素则失败。也许替换为numel(unique(sub))
  • @Gelliant:是的,很好,我已将您的建议添加到答案中。谢谢!
【解决方案2】:

对于任意数量的维度,您可以尝试convn

C = convn(total,reshape(sub(end:-1:1),size(sub)),'valid'); % flip dimensions of sub to be correlation
[~,indmax] = max(C(:));
% thanks to Eitan T for the next line
cc = cell(1,ndims(total)); [cc{:}] = ind2sub(size(C),indmax); subs = [cc{:}]

感谢 Eitan T 建议使用逗号分隔的列表作为广义的 ind2sub。

最后,您应该使用isequal 测试结果,因为这不是归一化的互相关,这意味着局部子区域中较大的数字会夸大相关值,从而可能导致误报。如果您的total 矩阵与大值区域非常不均匀,您可能需要在C 中搜索其他最大值。

【讨论】:

  • 这看起来很有希望,但有两点:1. 不能保证峰值与子矩阵的“左上角”相对应。 2、可以通过逗号分隔的列表获取ind2sub的输出参数(以here为例)。
  • 如果我们将sub2*total 进行比较,它仍然会返回匹配项。或者如果我们在total 中有一些元素比sub 中的值大一些数量级,那么这些区域将是可能的匹配项。如果您按照匹配滤波器的思路进行思考,则需要考虑正在匹配的信号的功率。
  • @Mohsen Nosratinia,你是对的——因此是我原始答案的最后一段。标准化的互相关是最好的。
  • 这个解决方案出奇的紧凑,但是我真的不明白它能做什么,不能做什么。例如,对于简单的情况它会失败:total = [1 2 3]; subs = 2。有没有办法使用这个解决方案而不用担心误报或否定?
  • 当模板/sub 矩阵较小时,误报的可能性变高。此外,正如 Mohsen 和我都指出的那样,您很可能在明显大于矩阵其余部分的 total 矩阵区域中得到误报。所以,这当然不是一个完整的解决方案,除非你实现 N 维归一化互相关 (NCC),它没有这个缺点。这对于典型的 FFT 解决方案更难做到,并且 非常 用加窗求和减慢了幼稚的方式 - 与在每个位置进行简单的等价测试相比,速度太慢了,不值得做!
【解决方案3】:

这是基于对原始矩阵 total 进行所有可能的移位,并将移位后的 total 的左上角等子矩阵与寻找的模式 subs 进行比较。使用字符串生成移位,并使用circshift 应用。

大部分工作都是矢量化完成的。只使用了一层循环。

该函数查找所有匹配项,而不仅仅是第一个。例如:

>> total = ones(3,4,5,6);
>> sub = ones(3,3,5,6);
>> matrixFind(total, sub)
ans =

     1     1     1     1
     1     2     1     1

函数如下:

function sol = matrixFind(total, sub)

nd = ndims(total);
sizt = size(total).';
max_sizt = max(sizt);
sizs = [ size(sub) ones(1,nd-ndims(sub)) ].'; % in case there are
% trailing singletons

if any(sizs>sizt)
    error('Incorrect dimensions')
end

allowed_shift = (sizt-sizs);
max_allowed_shift = max(allowed_shift);
if max_allowed_shift>0
    shifts = dec2base(0:(max_allowed_shift+1)^nd-1,max_allowed_shift+1).'-'0';
    filter = all(bsxfun(@le,shifts,allowed_shift));
    shifts = shifts(:,filter); % possible shifts of matrix "total", along 
    % all dimensions
else
    shifts = zeros(nd,1);
end

for dim = 1:nd
    d{dim} = 1:sizt(dim); % vectors with subindices per dimension
end
g = cell(1,nd);
[g{:}] = ndgrid(d{:}); % grid of subindices per dimension
gc = cat(nd+1,g{:}); % concatenated grid
accept = repmat(permute(sizs,[2:nd+1 1]), [sizt; 1]); % acceptable values
% of subindices in order to compare with matrix "sub"
ind_filter = find(all(gc<=accept,nd+1));

sol = [];
for shift = shifts
    total_shifted = circshift(total,-shift);
    if all(total_shifted(ind_filter)==sub(:))
        sol = [ sol; shift.'+1 ];
    end
end

【讨论】:

    猜你喜欢
    • 2013-06-18
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2012-03-16
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2015-05-20
    相关资源
    最近更新 更多