eeg_interp_sph_spline_g

 [Gx] = eeg_interp_sph_spline_g(COS)

 COS is the cosine matrix from elec_cosines

 Gx is the solution to Eq. 3 of Perrin et al. (1989)

 g(x) = 1/(4*pi) * (for n=1:inf, sum = sum + ( ( (2*n+1)/(n^m * (n+1)^m) ) * Pn(x) ) );

 where x is the cosine between two points on a sphere,
 m is a constant > 1 and Pn(x) is the nth degree Legendre
 polynomial.  Perrin et al. (1989) evaluated m=1:6 and 
 recommend m=4, for which the first 7 terms of Pn(x) 
 are sufficient to obtain a precision of 10^-6 for g(x)

 Notes:    This function solves g(x), Eq. 3 in
           Perrin et al (1989).  Electroenceph. & Clin. 
             Neurophysiology, 72: 184-187. Corrigenda (1990),
             Electroenceph. & Clin. Neurophysiology, 76: 565.
             (see comments in the .m file for details).

 see elec_cosines, eeg_interp_sph_spline, LegendreP

0001 function [Gx] = eeg_interp_sph_spline_g(COS)
0002 
0003 % [Gx] = eeg_interp_sph_spline_g(COS)
0004 %
0005 % COS is the cosine matrix from elec_cosines
0006 %
0007 % Gx is the solution to Eq. 3 of Perrin et al. (1989)
0008 %
0009 % g(x) = 1/(4*pi) * (for n=1:inf, sum = sum + ( ( (2*n+1)/(n^m * (n+1)^m) ) * Pn(x) ) );
0010 %
0011 % where x is the cosine between two points on a sphere,
0012 % m is a constant > 1 and Pn(x) is the nth degree Legendre
0013 % polynomial.  Perrin et al. (1989) evaluated m=1:6 and
0014 % recommend m=4, for which the first 7 terms of Pn(x)
0015 % are sufficient to obtain a precision of 10^-6 for g(x)
0016 %
0017 % Notes:    This function solves g(x), Eq. 3 in
0018 %           Perrin et al (1989).  Electroenceph. & Clin.
0019 %             Neurophysiology, 72: 184-187. Corrigenda (1990),
0020 %             Electroenceph. & Clin. Neurophysiology, 76: 565.
0021 %             (see comments in the .m file for details).
0022 %
0023 % see elec_cosines, eeg_interp_sph_spline, LegendreP
0024 
0025 % $Revision: 1.2 $ $Date: 2005/07/13 06:03:50 $
0026 
0027 % Licence:  GNU GPL, no implied or express warranties
0028 % History:  08/2001 Darren.Weber_at_radiology.ucsf.edu, with
0029 %                   mathematical advice from
0030 %                   Dr. Murk Bottema (Flinders University of SA)
0031 %           10/2003 Darren.Weber_at_radiology.ucsf.edu, with
0032 %                   mathematical advice and LegendreP function from
0033 %                   Dr. Tom.Ferree_at_radiology.ucsf.edu
0034 %
0035 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0036 
0037 eegversion = '$Revision: 1.2 $';
0038 fprintf('EEG_INTERP_SPH_SPLINE_G [v %s]\n',eegversion(11:15)); tic
0039 fprintf('...calculating Legendre function of the cosine matrix\n');
0040 
0041 m = 4;
0042 N = 7;
0043 
0044 n = [1:N]';
0045 Series = (2*n + 1) ./ (n.^m .* (n+1).^m);
0046 %fprintf('%12.10f  ',Series) % note how Series diminishes quickly
0047 %0.1875000000  0.0038580247  0.0003375772  0.0000562500
0048 %0.0000135802  0.0000041778  0.0000015252
0049 
0050 
0051 % Perrin et al. (1989) recommend tabulation of g(x) for
0052 % x = linspace(-1,1,2000) to be used as a lookup given actual
0053 % values for cos(Ei,Ej).
0054 if isempty(COS),
0055   error('...cosine matrix is empty.\n');
0056   %msg = sprintf('...cosine matrix empty, generating COS.\n');
0057   %warning(msg);
0058   %COS = linspace(-1,1,2000)';
0059 end
0060 
0061 Gx = zeros(size(COS));
0062 
0063 CONST = 1/(4*pi);
0064 
0065 for i = 1:size(COS,1),
0066   for j = 1:size(COS,2),
0067     
0068     % P7x is an 8x1 column array, starting at P0x (we only need P1x - P7x)
0069     P7x = LegendreP(N,COS(i,j));
0070     
0071     Gx(i,j) = CONST * dot( Series, P7x(2:N+1) );
0072     
0073   end
0074 end
0075 
0076 t=toc; fprintf('...done (%6.2f sec)\n\n',t);
0077 
0078 return