Re: Vectorizing loops

[c.]

On 5 Nov 2011, at 07:57, andrewcd wrote:

> Hi All,
> 

> My code (which works) with all of the loops is:
> G=[]
> for j = 1:columns(N)
>               disp(j)
>               disp ('of 259')
> for i = 1:rows(N)
> for k = 1:size(N,3)
> G(size(G,1)+1,4) = N(i,j,k);
> G(size(G,1),1) = Nxy(1,j,k);
> G(size(G,1),2) = Nxy(i,1,k);
> G(size(G,1),3) = k;
> end end end

the following seems to work for me:

## construct fake data
Nxy = reshape (1:(30*29*8), [30, 29, 8]);
N   = Nxy;

## original code
G = [];

for j = 1:columns(N)

disp (j)
disp ('of 259')

for i = 1:rows(N)
for k = 1:size(N,3)
G(size(G,1)+1,4) = N(i,j,k);
G(size(G,1),1) = Nxy(1,j,k);
G(size(G,1),2) = Nxy(i,1,k);
G(size(G,1),3) = k;
end end end

## vectorized code
X = Nxy(1, :, :);
X = repmat (X, [rows(N), 1, 1]);
Y = Nxy(:, 1, :);
Y = repmat (Y, [1, columns(N), 1]); 

T = repmat ([1:size(N, 3)](:), [rows(N)*columns(N), 1]);
X = permute (X, [3, 1, 2]);
Y = permute (Y, [3, 1, 2]);
N = permute (N, [3, 1, 2]);

G1 = [X(:), Y(:), T, N(:)];

## comparison
assert (G1, G)

HTH,
c.