
%% ================================================================
% John Kerl
% kerl.john.r@gmail.com
% 2007-02-08
%
% A Matlab script to display the two-sphere along with a tangent plane.
%% ================================================================

%% ----------------------------------------------------------------
ntheta=20;
nphi=20;
ns = 20;
nt = 20;

thetamin=0;
thetamax=2*pi;

phimin=0;
phimax=pi;

smin = -1;
smax =  1;

tmin = -1;
tmax =  1;

dtheta=(thetamax-thetamin)/ntheta;
dphi=(phimax-phimin)/nphi;
ds=(smax-smin)/ns;
dt=(tmax-tmin)/nt;

%% ----------------------------------------------------------------
clf;   % Clear the figure in case this is a redisplay
rotate3d on; % Permit me to rotate the graph with click and drag
hold on; % All plot3 output goes on the same graph

%% ----------------------------------------------------------------
% Draw latitude lines
for phi = phimin : dphi : phimax
	theta = [ thetamin : dtheta : thetamax ];
	x = cos(theta) * sin(phi);
	y = sin(theta) * sin(phi);
	z = cos(phi) * ones(1, ntheta+1);
	plot3(x,y,z,'g')
end

% Draw longitude lines
for theta = thetamin : dtheta : thetamax
	phi = [ phimin : dphi : phimax ];
	x = cos(theta) * sin(phi);
	y = sin(theta) * sin(phi);
	z = cos(phi);
	plot3(x,y,z,'g')
end

% Draw the equator
phi = pi/2;
theta = [ thetamin : dtheta : thetamax ];
x = cos(theta) * sin(phi);
y = sin(theta) * sin(phi);
z = cos(phi) * ones(1, ntheta+1);
plot3(x,y,z,'r')

% Draw the Greenwich meridian
theta = 0
phi = [ phimin : dphi : phimax ];
x = cos(theta) * sin(phi);
y = sin(theta) * sin(phi);
z = cos(phi);
plot3(x,y,z,'k')

% Matlab line-style codes:
%           b     blue          .     point              -     solid
%           g     green         o     circle             :     dotted
%           r     red           x     x-mark             -.    dashdot 
%           c     cyan          +     plus               --    dashed   
%           m     magenta       *     star             (none)  no line
%           y     yellow        s     square
%           k     black

%% ----------------------------------------------------------------
% Point q is (1/2, 1/2, 1/sqrt 2).  This has spherical coordinates
% theta = 45 deg and phi = 45 deg.
qx = 1/2;
qy = 1/2;
qz = 1/sqrt(2);

q = [qx; qy; qz];

qtheta = pi/4;
qphi   = pi/4;

%% - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
% Tangent vector 1 is d/dx using graph coordinates
% Tangent vector 2 is d/dy using graph coordinates

%  q =
% [ x ]
% [ y ]
% [ (1 - x^2 - y^2)^1/2 ]

% d/dx = 
% [ 1 ]
% [ 0 ]
% [ -x/z ]

% d/dy = 
% [ 0 ]
% [ 1 ]
% [ -y/z ]

ux = 1;
uy = 0;
uz = -qx/qz;

vx = 0;
vy = 1;
vz = -qy/qz;

%% - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
% Tangent vector 1 is d/dtheta using spherical coordinates
% Tangent vector 2 is d/dphi   using spherical coordinates

%  q =
% [ cos(theta) sin(phi)]
% [ sin(theta) sin(phi)]
% [ cos(phi)]

% d/dtheta =
% [-sin(theta) sin(phi)]
% [ cos(theta) sin(phi)]
% [ 0]

% d/dphi =
% [ cos(theta) cos(phi)]
% [ sin(theta) cos(phi)]
% [-sin(phi)]

% d/dtheta at q =
% [-1/2]
% [ 1/2]
% [ 0]

% d/dphi at q =
% [ 1/2 ]
% [ 1/2 ]
% [-1/sqrt(2)]

%ux = -sin(qtheta) * sin(qphi);
%uy =  cos(qtheta) * sin(qphi);
%uz =  0;
%
%vx =  cos(qtheta) * cos(qphi);
%vy =  sin(qtheta) * cos(qphi);
%vz = -sin(qphi);

u = [ux; uy; uz];
v = [vx; vy; vz];

%q
%u
%v
%u' * v
%q' * u
%q' * v


% Draw t lines
for s = smin : ds : smax
	t = [ tmin : dt : tmax ];
	x = qx + s * ux + t * vx;
	y = qy + s * uy + t * vy;
	z = qz + s * uz + t * vz;
	plot3(x,y,z,'b')
end

% Draw s lines
for t = tmin : dt : tmax
	s = [ smin : ds : smax ];
	x = qx + s * ux + t * vx;
	y = qy + s * uy + t * vy;
	z = qz + s * uz + t * vz;
	plot3(x,y,z,'b')
end

