elec_sphere_project

 elec_sphere_project - calculates spherical projections of electrode positions

 Usage: [r,x,y,z] = elec_sphere_project(X,Y,Z,xo,yo,zo,[estim],[plot])

 Notes:    The general formula for a sphere, with radius r is given by:

           (x - xo)^2  +  (y - yo)^2  +  (z - zo)^2  =  r^2

           This function takes arguments for cartesian co-ordinates
           of X,Y,Z (assume Z > 0) and the center of the sphere (xo,yo,zo).
           If (xo,yo,zo) are not provided, they are assumed (0,0,0).

           Returned values are the fitted radius 'r' (constant)
           and the (x,y,z) Cartesian co-ordinates on the projected sphere.

           Options: estim - echo spherical radius estimates (default),
                    plot  - plot the input/projected xyz on a sphere.

0001 function [r,x,y,z] = elec_sphere_project(X,Y,Z,xo,yo,zo,varargin)
0002 
0003 % elec_sphere_project - calculates spherical projections of electrode positions
0004 %
0005 % Usage: [r,x,y,z] = elec_sphere_project(X,Y,Z,xo,yo,zo,[estim],[plot])
0006 %
0007 % Notes:    The general formula for a sphere, with radius r is given by:
0008 %
0009 %           (x - xo)^2  +  (y - yo)^2  +  (z - zo)^2  =  r^2
0010 %
0011 %           This function takes arguments for cartesian co-ordinates
0012 %           of X,Y,Z (assume Z > 0) and the center of the sphere (xo,yo,zo).
0013 %           If (xo,yo,zo) are not provided, they are assumed (0,0,0).
0014 %
0015 %           Returned values are the fitted radius 'r' (constant)
0016 %           and the (x,y,z) Cartesian co-ordinates on the projected sphere.
0017 %
0018 %           Options: estim - echo spherical radius estimates (default),
0019 %                    plot  - plot the input/projected xyz on a sphere.
0020 %
0021 
0022 % $Revision: 1.3 $ $Date: 2005/07/12 22:16:48 $
0023 
0024 % Licence:  GNU GPL, no express or implied warranties
0025 % History:  02/2002, Darren.Weber_at_radiology.ucsf.edu
0026 %                    adapted from elec_fit_sph
0027 %
0028 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0029 
0030 % initialise centroid, unless input parameters defined
0031 if ~exist('xo','var') xo = 0; else, if isempty(xo) xo = 0; end, end
0032 if ~exist('yo','var') yo = 0; else, if isempty(yo) yo = 0; end, end
0033 if ~exist('zo','var') zo = 0; else, if isempty(zo) zo = 0; end, end
0034 
0035 if (nargin>6),        estim = varargin{1};
0036 else                  estim = 1;
0037 end
0038 if (nargin>7),        plot = varargin{2};
0039 else                  plot = 0;
0040 end
0041 
0042 eegversion = '$Revision: 1.3 $';
0043 fprintf('ELEC_SPHERE_PROJECT [v %s]\n',eegversion(11:15));
0044 fprintf('...projecting electrodes to spherical points\n');
0045 tic;
0046 
0047 % Initialise r0 as a rough guess at the sphere radius
0048 rX = (max(X) - min(X)) / 2;
0049 rY = (max(Y) - min(Y)) / 2;
0050 rZ =  max(Z) - zo;
0051 r0 = mean([ rX rY rZ ]);
0052 
0053 if isequal(estim,1),
0054     fprintf('...initial spherical radius   = %f\n', r0); end
0055 
0056 % perform least squares estimate of spherical radius (r)
0057 options = optimset('fminsearch');
0058 [r,fval,exitflag,output] = fminsearch('elec_sphere_fit_optim',...
0059                            r0, options, X, Y, Z, xo, yo, zo);
0060 %fprintf('\n%s%f\n', 'Iterations = ', output.iterations);
0061 %fprintf('\n%s%d\n', 'Exit = ', exitflag);
0062 
0063 if isequal(estim,1),
0064     fprintf('...estimated spherical radius = %f\n', r); end
0065 
0066 
0067 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0068 % Find the projection point of X,Y,Z to the fitted sphere radius r
0069 
0070 % Convert Cartesian X,Y,Z to spherical (radians)
0071 theta = atan2( (Y-yo), (X-xo) );
0072 phi = atan2( sqrt( (X-xo).^2 + (Y-yo).^2 ), (Z-zo) );
0073 % do not recalc: r = sqrt( (X-xo).^2 + (Y-yo).^2 + (Z-zo).^2);
0074 
0075 %   Recalculate X,Y,Z for constant r, given theta & phi.
0076 R = ones(size(phi)) * r;    
0077 x = R .* sin(phi) .* cos(theta);
0078 y = R .* sin(phi) .* sin(theta);
0079 z = R .* cos(phi);
0080 
0081 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0082 % Plot the input & projected electrode positions on a sphere
0083 if isequal(plot,1),
0084     figure('NumberTitle','off','Name','Electrode Placements');
0085     set(gca,'Projection','perspective');
0086     set(gca,'DataAspectRatio',[1 1 1]);
0087 
0088     hold on
0089 
0090     plot3(x,y,z,'b.');
0091     plot3(X,Y,Z,'ro');
0092     legend('input xyz','fitted sphere','Location','BestOutside');
0093 
0094     [Xs,Ys,Zs]=sphere;
0095     Xs = Xs * r;
0096     Ys = Ys * r;
0097     Zs = Zs * r;
0098 
0099     surf(Xs,Ys,Zs,...
0100         'EdgeColor','none',...
0101         'FaceColor',[0.7 0.7 0.7],...
0102         'FaceAlpha',0.4);
0103     view(2);
0104     rotate3d;
0105 
0106     axis tight; hold off;
0107 end
0108 
0109 t = toc; fprintf('...done (%6.2f sec).\n\n',t);
0110 
0111 return