考虑一个矩阵X
。我必须用矩阵 Z(大小为row2-row1+1
,col2-col1+1
)更新X
、X(row1:row2, col1:col2)
的子矩阵,但仅限于逻辑矩阵L
(大小为row2-row1+1
,col2-col1+1
)为真的位置。
例如,如果
X=[ 1 2 3 4 5 6
11 12 13 14 15 16
21 22 23 24 25 26
31 32 33 34 34 36]
Z=[31 41
32 42]
L=[ 1 0
0 1]
row1 = 2; row2 = 3; col1 = 3; col2 = 4
然后在更新后我应该得到:
X=[ 1 2 3 4 5 6
11 12 31 14 15 16
21 22 23 42 25 26
31 32 33 34 34 36]
目前我做以下工作:
Y = X(row1:row2, col1:col2);
Y(L) = Z(L);
X(row1:row2, col1:col2) = Y;
这段代码处于一个紧密循环中,根据 Matlab (v2019a) 的分析器是我程序的主要瓶颈。在实际代码中X
是一个 2000x1500x3 的立方体; 循环中的row1
、row2
、col1
、col2
、Z
和L
变化。
问题是它是否可以重写为单个/更快的作业。
谢谢。
老实说,在没有看到您的实际代码的情况下,我感觉您的解决方案可能会尽可能快。我之所以这么说,是因为我通过创建一些更接近您的实际问题的随机样本数据来测试几种不同的解决方案。我假设X
是大小为 2000 x 1500-x3 的uint8
型图像,Z
大小为N
xN
(即我们只会修改X
的第一页),L
是一个N
×N
的逻辑数组,行和列索引是随机选择的:
X = randi([0 255], 2000, 1500, 3, 'uint8');
N = 20; % Submatrix size
Z = randi([0 255], N, N, 'uint8');
L = (rand(N, N) > 0.5);
row1 = randi([1 2000-N]);
row2 = row1+N-1
col1 = randi([1 1500-N]);
col2 = col1+N-1;
然后,我测试了 3 种不同的解决方案:您的原始解决方案、使用find
和sub2ind
为X
创建线性索引的解决方案,以及为X
创建逻辑索引的解决方案:
% Original solution:
Y = X(row1:row2, col1:col2, 1);
Y(L) = Z(L);
X(row1:row2, col1:col2, 1) = Y;
% Linear index solution:
[rIndex, cIndex] = find(L);
X(sub2ind(size(X), rIndex+row1-1, cIndex+col1-1)) = Z(L);
% Logical index solution
[R, C, ~] = size(X);
fullL = false(R, C);
fullL(row1:row2, col1:col2) = L;
X(fullL) = Z(L);
我使用timeit
随机生成的样本数据反复测试了这些,发现您的原始解决方案始终是最快的。线性指数解决方案非常接近,但速度稍慢。逻辑索引解决方案需要两倍以上的时间。
让我们定义一些示例数据:
X = randi(9,5,6);
Y = 10+X;
row1 = 2;
row2 = 4;
col1 = 3;
col2 = 4;
L = logical([0 1; 0 0; 1 1]);
然后:
ind_subm = bsxfun(@plus, (row1:row2).',size(X,1)*((col1:col2)-1));
% linear index for submatrix
ind_subm_masked = ind_subm(L);
% linear index for masked submatrix
X(ind_subm_masked) = Y(ind_subm_masked);
示例结果:
之前
X
:X = 6 2 1 7 9 6 3 3 3 5 5 7 6 3 8 6 5 4 7 4 1 3 3 4 2 5 9 5 5 9
L
:L = 3×2 logical array 0 1 0 0 1 1
X
之后:X = 6 2 1 7 9 6 3 3 3 15 5 7 6 3 8 6 5 4 7 4 11 13 3 4 2 5 9 5 5 9