function Z=join_tr(X,Y)
% join_tr - Cartesian product of two arrays, transpose of join
%
% join_tr :: [[M,N]], [[P,L]] -> [[M+P,N*L]].
%
% The result contains every pairing of rows from the two arguments

%Z=[repmat(X,size(Y,1),1) kron(Y,ones(size(X,1),1))];

% this is probably faster (no multiplications!)
m=size(X,2); n=size(Y,2);
Z=[repmat(X,1,n); Y(:,repmat((1:n),m,1))];
