How do I convert an array into a set of indices using colon notation?

I have two 4D matrices. Each is approximately 150 x 150 x 150 x 3 (the "3" is exact in both matrices). The 4th dimension in Matrix1 specifies a location from which to extract data in Matrix2.
So, for example, if Matrix1(20,30,25,:) = [22 32 15], this is the location in Matrix2 which we want to use to extract data for that 3D location. So, the result ShiftPart1(20,30,25,:) = Matrix2(22,32,15,:).
I can calculate this in a nested loop as follows, but it is extremely inefficient, take ~ 2 minutes. I need to run these calculations thousands of times, so it must take < 1 sec. Is there a vectorized way of running this to avoid the nested loops? Thanks!
function [ShiftPart1] = getIndexShift(Matrix1, Matrix2)
size1 = size(Matrix1);
for j = 1:size1(1),
for k = 1:size1(2),
for sl = 1:size1(3),
ShiftPart1(j,k,sl,1:3) = squeeze(double(Matrix2(Matrix1(j,k,sl,1:3))));
end,
end,
end

Answers (1)

[m, n, k, ~] = size(Matrix1);
[d,c,b,a] = ndgrid(1:3,1:k,1:n,1:m);
ii = sub2ind([m,n,k,3],a(:),b(:),c(:),d(:));
out = reshape(Matrix2(Matrix1(ii)),[m,n,k,3]);

1 Comment

Thanks, but this doesn't work. This gives the same result as if I simply used out = Matrix2(Matrix1);
I just realized that the looping structure which I provided will not give the right answer either! That's because, as in my example, Matrix2(Matrix1(20,30,25,1:3)) is NOT the same as Matrix2(22,32,15,1:3). The former is Matrix2([22 32 15]) and is NOT the correct answer. Here is the correct looping structure:
for j = 1:size1(1), for k = 1:size1(2), for sl = 1:size1(3), ShiftPart1(j,k,sl,1:3) = squeeze(double(Matrix2(Matrix1(j,k,sl,1),Matrix1(j,k,sl,2),Matrix1(j,k,sl,3),:))); end, end, end
The issue with doing it the simple way may be due to the fact that there are non-unique elements, i.e. multiple Matrix1 elements can point to the same point in Matrix2. Thus, there may be no simplified form which solves this, but thanks for any further insights you can provide! Best, Mark

Sign in to comment.

Asked:

on 22 Jun 2017

Commented:

on 22 Jun 2017

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!