An ellipsoid obeys the inequality (x-c)'*Q*(x-c)<=1 for some appropriate 3x3 matrix Q and 3x1 vector c (which is also the center of the ellipsoid). If you have it expressed this way, generation of the mask is very simple:
[dx,dy,dz]=ndgrid((1:nx)-c(1), (1:ny)-c(2)), (1:nz)-c(3) );
X=[dx(:),dy(:),dz(:)]';
mask = sum(X.*(Q*X))<=1;
mask=reshape(mask,nx, ny,nz);
