function [hit,s] = sphere(R,xs,x0,u,inside)

% Calculates the distance to a sphere of radius R centered at (xs,ys,zs)

% Inputs
%%%%%%%%
%  double
%     R :  Radius of the sphere (scalar)
%     xs:  coordinates of the center of the sphere (row vector)
%     x0:  coordinates of the particle (row vector)
%     u :  direction cosines of the particle (row vector)
%  logical
%     inside: == true   particle thinks it is inside
%           : == false  particle thinks it is outside
%      

% Outputs
%%%%%%%%%
%  logical
%     hit: == true   particle would hit  the surface
%        : == false  particle would miss the surface
%
%  double
%     s: distance to the surface (if hit)
%     
% Internal variables
%  double
%     xr: relative coordinates of the particle wrt the sphere (row vector)
%     A:  Quadratic coefficient A (communicated to quadric) (scalar)
%     B:  Quadratic coefficient B (communicated to quadric) (scalar)
%     C:  Quadratic coefficient C (communicated to quadric) (scalar)

% Compute relative coordinates

xr = x0 - xs;

% Compute the quadratic coefficients
A = 1; % i.e. assuming u**2 + v**2 + w**2 = 1
B = u*xr';
C = xr*xr' - R.^2;

% Call the quadratic solver 
[hit,s] = quadric(A,B,C,inside);

return