function c = cartprod(varargin)
%c = allcomb(v1, v2, ....) create a matrix of the cartesian product of all vector
%v1, v2,... must be numeric vectors
%c(r, :) is a unique combination of the elements of v1, v2, ...
%c is size numel(v1)*numel(v2)*... x N (N: number of input vectors).
%TODO: add input validation
c = cell(size(varargin)); %create cell array to store output of ndgrid
[c{:}] = ndgrid(varargin{:}); %expand varargin into c-s-l and store output in c. C{i} is N-d array (N = number of vectors)
c = cat(numel(varargin) + 1, c{:}); %concatenate the N-d arrays into a N+1-d array
c = reshape(c, [], numel(varargin)); %and reshape into ?xN matrix
end
