function [xi] = symtri(N)

u  = rand(N,1);
xi = zeros(N,1);

for ii = 1:N,
  if (u(ii) < 0.5),
    xi(ii) = sqrt(u(ii)/2);
  else,
    xi(ii) = 1 - sqrt((1-u(ii))/2);
  end
end