Can I speed this code up, looking for similarity between two 3-d matrices.
1 visualización (últimos 30 días)
Mostrar comentarios más antiguos
Is it possible to vectorize the following code in order to speed it up?
This little routine gets called tens of thousands of times in some code I'm working on, and it is listed as the major bottle neck in profile viewer.
What it does is shift a small 3D matrix A through a large 3D matrix B looking for the best match in B to A. My current thoughts are that it can be vectorized by reshaping and replicating A and B, but I can't figure out how to do it.
temp_position = zeros(3,1);
new_position = zeros(3,1);
ad_current = Inf;
for ii = 1:size(B,1)-size(A,1)+1
for jj = 1:size(B,2)-size(A,2)+1
for kk = 1:size(B,3)-size(A,3)+1
ad_new = sum(reshape(abs(B(ii:ii+size(A,1)-1,jj:jj+size(A,2)-1,kk:kk+size(A,3)-1) - A),[],1));
if ad_new < ad_current
ad_current = ad_new;
temp_position = [ii,jj,kk];
end
end
end
end
new_position = ... something + temp_position
2 comentarios
Jan
el 19 de Feb. de 2011
Is "B(jj:jj+size(A,1)-1,ii:ii+size(A,2)-1, ..." a typo? Do you mean "ii"<->"jj" here?
Respuesta aceptada
Bruno Luong
el 20 de Feb. de 2011
The trick here is loop on A (small size) and vectorized on B
% Test data
B=rand(50,60,100);
A=rand(2,3,4);
%%Engine
szA = size(A);
szB = size(B);
szC = szB-szA+1;
C = zeros(szC);
szBs = ones(1,6);
szBs([1 3 5]) = szA;
AA = reshape(A, szBs);
for n=1:numel(A)
first = cell(1,3);
[first{:}] = ind2sub(szA,n);
first = cat(2,first{:});
nb = floor((szB-first+1)./szA);
lgt = nb.*szA;
last = first + lgt - 1;
iB = arrayfun(@(i,j) i:j, first, last, 'Unif', false);
Bs = B(iB{:});
szBs([2 4 6]) = nb;
Bs = reshape(Bs, szBs);
BmA = bsxfun(@minus, Bs, AA);
d = sum(sum(sum(abs(BmA),1),3),5);
iC = arrayfun(@(i,s,j) i:s:j, first, szA, last, 'Unif', false);
C(iC{:}) = reshape(d, nb);
end
[mindiff loc] = min(C(:));
[ii jj kk] = ind2sub(szC, loc);
loc = [ii jj kk];
% Check
disp(mindiff)
disp(loc)
0 comentarios
Más respuestas (6)
Jan
el 20 de Feb. de 2011
Just some marginal changes:
sizeA = size(A);
sA3m1 = sizeA(3) - 1;
sizeB = size(B);
ad_current = Inf;
for ii = 1:sizeB(1)-sizeA(1)+1
Q = B(ii:ii+sizeA(1)-1, :, :);
for jj = 1:sizeB(2)-sizeA(2)+1
P = Q(:, jj:jj+sizeA(2)-1, :);
for kk = 1:sizeB(3)-sA3m1
ad_new = sum(reshape(abs( ...
P(:, :, kk:kk+sA3m1) - A),[],1));
if ad_new < ad_current
ad_current = ad_new;
temp_position = [ii,jj,kk];
end
end
end
end
EDITED: Q(:,:, kk:...) -> P(:, :, kk:...) Thanks Bruno
0 comentarios
Doug Hull
el 18 de Feb. de 2011
Is there a way that convn can be used to do this? I don't know the details, but I think it might help.
0 comentarios
tlawren
el 21 de Feb. de 2011
1 comentario
Bruno Luong
el 21 de Feb. de 2011
Are the As same size?
I don't understand how you run multiple A in your code. If you run multiple As with a for loop, I can't see how you function can be suddenly faster.
When asking a question, attach a little code is better than 1000 words.
tlawren
el 21 de Feb. de 2011
2 comentarios
Bruno Luong
el 21 de Feb. de 2011
Something is fishy. If you loop over A, the compute time is just proportional with the number of A. It did not make sense to me why your code is more efficient when running with many As.
Ver también
Categorías
Más información sobre Graphics Performance en Help Center y File Exchange.
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!