function [fx,fy,fz]=hopfmap(a,b,c)

proj=2;             % 1 - stereographic, 2- another projection
phi=0:pi/100:2*pi;
theta=atan2(b,-a)-phi;
alpha=sqrt((1+c)/2);
beta=sqrt((1-c)/2);

% Inverse Hopf Map
x=beta*cos(phi);
y=beta*sin(phi);
w=alpha*cos(theta);
z=alpha*sin(theta);

if proj==1
    % Stereographic Projection
    fx=x./(1-w);
    fy=y./(1-w);
    fz=z./(1-w);
else
    % Different Projection
    r=acos(w)/pi;
    %r=1;
    fx=r./sqrt(1-w.^2).*x;
    fy=r./sqrt(1-w.^2).*y;
    fz=r./sqrt(1-w.^2).*z;
end