function  [v,w]= housemap(a);
%% function  [v,w]= housemap(a);
%% To get w so (I-ww^t)a=v, v=[alpha,0,...,0]^t,  |alpha|= ||a||_2
%%  INPUT:
%%    a - a column vector;
%%  OUTPUT:
%%    v - a column vector of the form v=[alpha,0,...,0]^t,  |alpha|= ||a||_2
%%    w - a column vector 
%%
n=length(a);
if a(1)>=0
   alpha=norm(a,2);
else
   alpha = -norm(a,2);
end;
v= [alpha,zeros(1,n-1)]';
r=sqrt((alpha^2-a(1)*alpha)/2);
w(2:n,1)=a(2:n)/(2*r);w(1,1)=(a(1)-alpha)/(2*r);
