function f = matfactorial(a)

N = a(:);
a(N>170) = 171;
m = max([1;a(:)]);
N = [1 1 cumprod(2:m)];
f = N(a+1);

% b = a(:);
% B = length(b);
% for i = 1:B
%     f(i) = factorial(b(i));
% end
% f = reshape(f,size(a));