Matrix operations in three dimensions

I'm creating synthetic surfaces by summing multiple frequency components of random phase and orientation. Is there a faster way to do this? It seems that I could build the 3D components matrix without a loop, but I'm not sure how.

% Create domain
[X,Y] = meshgrid(0:.1:50,0:.1:50); 
f = logspace(0,3,10); 
theta = 90*rand(size(f)); 
A = f/1000; 

% Build one surface for each frequency: 
components = NaN(size(X,1),size(X,2),length(f)); 
for k = 1:length(f)
    components(:,:,k) = A(k) * sin(f(k)*(X*sind(theta(k))+Y*cosd(theta(k))) + rand);
end

% Sum all components to create surface: 
Z = sum(components,3);

Answer 1

[James Tursa]

You can build components without a loop by reshaping the vectors A, f, theta as 1x1xN and then using bsxfun with @times to get the 3rd dimension slices. But the code would look quite ugly and I don't know if the speed improvement (if any) will be worth it to you. E.g.,

N = numel(f);
f = reshape(f,1,1,N);
theta = reshape(theta,1,1,N);
A = reshape(A,1,1,N);
Components = bsxfun(@times,A,sin(bsxfun(@plus,bsxfun(@times,X,f.*sind(theta)) + bsxfun(@times,Y,f.*cosd(theta)),rand(1,1,N))));

[Chad Greene]

Thanks James. For the work I'm doing your bsxfun method seems to be about 10% faster than looping.