嵌套if语句的循环向量化

3

问题

我正在尝试优化我的代码运行时间,并且之前已经提出过一个类似的问题,其中包含了几个嵌套的 if 语句。向量化嵌套的 if 语句

由于我在那里发布的代码有些冗长,而且我仍然在努力实现嵌套循环的向量化,所以我想再次发问,这次我将提供一些更简单的代码:

代码

NB_list_all=zeros(length(BML),4);
for NB=1:length(BML);
    NB_list=zeros(4,1);
    %in +x direction
    if isempty(find(BML(:,2)==BML(NB,2)+x_blockdimension & BML(:,3)==BML(NB,3), 1));
        NB_list(1,1)=0;
    else
        NB_list(1,1)=find(BML(:,2)==(BML(NB,2)+x_blockdimension)& BML(:,3)==BML(NB,3));

    end
NB_list_z(NB,:)=NB_list;
end

% BML(:,2) stores the x-coordinate
% BML(:,3) stores the y-coordinate

一些示例数据

BML=
    1 1005  115
    2 1100  115
    3 1419  120
    4 1424  120
    5 660   115
    6 655   115

注意,BML的大小为170,000 x 7。
代码描述:
我尝试使用这段代码来查找我的点云中距离“x_blockdimension”一定距离的下一个点。如果没有找到任何点,则将该条目设置为零。由于这需要大量时间处理1800万个点(而且我不仅仅是在一个方向上查找),因此我正在寻找一种通过使用向量化或逻辑索引来优化此过程的方法。如果有其他方法可以改善运行时间,我会很高兴得到任何提示。
我尝试过:
if isempty(find(BML(:,2)==BML(:,2)+x_blockdimension & BML(:,3)==BML(:,3), 1));
    NB_list(1,1)=0;
else
    NB_list(1,1)=find(BML(:,2)==(BML(:,2)+x_blockdimension)& BML(:,3)==BML(:,3));

end

但它并没有真正做我想要的事情。 我希望能得到一些帮助!
2个回答

1
如果我正确理解输入的格式,您可以使用广播和 bsxfun 实现矢量化解决方案,示例如下 -
% Perform broadcasted comparison corresponding to the iterative comparison
% in original code to get a boolean array/mask. 
% Get the row, col indices for the mask, which will correspond to the index
% values and positions where those values are to be stored in the output. 
[R,C] = find(bsxfun(@eq,BML(:,2),BML(:,2).'+x_blockdimension) & ...
                                     bsxfun(@eq,BML(:,3),BML(:,3).'));

% Setup output array and store the indices at respcetive positions. 
NB_list_z_out = zeros(size(BML,1),numel(NB_list));
NB_list_z_out(C,1) = R;

请注意,似乎输出仅编辑输出数组中的第一列元素,因此在最后一步使用NB_list_z_out(C,1)进行索引。
另一种方法可以建议关注内存效率,并且还可以提高性能,获取 RC,以后可以像之前列出的方法一样使用。实现应该是这样的 -
% Filter out with "bsxfun(@eq,BML(:,3),BML(:,3).'))".
[~,~,idx] = unique(BML(:,3),'stable');
vidx = find(ismember(idx,find(accumarray(idx(:),1)>1)));

% Filter out on remaining ones with valid indices (vidx)
[R1,C1] = find(bsxfun(@eq,BML(vidx,2),BML(vidx,2).'+x_blockdimension));
R = vidx(R1);
C = vidx(C1);

@KiW,根据您的实际情况,“BML”的大小是多少? - Divakar
@ Divakar http://s000.tinyupload.com/?file_id=05161766586627414815 这是"BML"文件,x_blockdimension只是等于5。 - KiW
如果我想添加第三个维度(即x,y,z),我是否可以像您刚才编辑的那样以同样的方式添加它,还是必须有两个稳定的索引?...说实话,您写的那些行对我来说很难理解:D - KiW
@KiW 取决于您在 z 上使用条件的方式。但是,bsxfun 部分将扩展以包括 BML(vidx,4) - Divakar
在处理大型数据集的内存问题上,您可以计算vidx,它基本上意味着有效的索引来遍历行。因此,您可以坚持使用循环代码并使用vidx。因此,for NB=1:length(BML)可以替换为for NB=1:length(vidx)并继续使用您的循环代码,这将为您提供R1C1的等效项。在所有操作结束时,我们需要使用vidx(R1)等进行索引。 - Divakar
显示剩余18条评论

1
如果您知道BML中每行只有0或1个匹配项,那么您可以对它们进行排序并使用diff,而不是使用循环:
%%  Find matches for x dimension

% sort on x dimension using sortrows, then split the matrix again
BMLx= sortrows(BML(:,[2:end,1]));
sorted_xx = BMLx(:,1:end-1);
idx = BMLx(:,end);

diff_ = diff(sorted_xx);
x_diff_match = find(diff_(:,1)==x_blockdimension & (diff_(:,2:end)==0));
% or you might want to use abs(a-b)<told

% assign all zeros as default
NB_list_x=zeros(length(BML),1);
% alocate matches
NB_list_x(idx(x_diff_match)) = idx(x_diff_match+1)

嘿,谢谢你的回答。我知道对于相同的y值,x + x-维度只有一个值,但总共有更多满足条件的值。如果我看你的代码,似乎没有包含这个条件,或者是我太慢了看不到? - KiW
抱歉,误解了问题。现在检查其他差异(other_dimesions==0)。 - Dave Rayner
我尝试实现您的解决方案,看起来是个好主意,特别是对于大数据集,但我只得到了零值作为邻居。老实说,我不知道为什么...如果您有时间,我在这里上传了数据集:s000.tinyupload.com/?file_id=05161766586627414815...可能更容易看到结果。我猜你的意思是´(diff_(:,2)==0)´,否则矩阵维度不一致。但如果您有任何想法为什么结果只是零,我会很高兴的:/ 感谢您的帮助! - KiW
如果您有解决方案,我会非常高兴,因为您的方法似乎更快,但我只是想不出问题所在... - KiW
嗨KiW,我刚看了一下你的tinyupload文件。问题可能是tinyupload文件的第一列包含大整数,而你帖子中的示例使用的是1:end。如果这些是其他索引号码,那么把它们替换为1:length(BML)是否可行?否则,如果第一列实际上是x值,则在第一列之前连接一个1:n向量。 - Dave Rayner

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接