vb_current_vector

 Calculate current direction
   [V,xx] = vb_current_vector(V,xx,Lmode)
 --- INPUT
 V         : Cortex vertex coordinate [NV , 3]
 xx        : Normal direction of cortex [NV , 3]
 V(n,:)    : 3-D coordinate   at n-th vertex
 xx(n,:)   : Normal direction at n-th vertex
 Lmode     : number of independent curent direction
       = 1 :  current vextor = xx(n,:)
       = 2 :  current vextor = orthogonal vectors to V(n,:)
       = 3 :  current vextor = xx(n,:) & orthogonal vectors to xx(n,:)

 --- OUTPUT
  V( NV*Lmode , 3) = repmat(V, [Lmode 1])
 xx( NV*Lmode , 3)

 Lmode = 1 :
   xx = xx(input) : Normal direction

 Lmode = 2 : 
   xx(n,:)      = 1st orthogonal vectors to V(n,:)
   xx(n + NV,:) = 2nd orthogonal vectors to V(n,:)

 Lmode = 3 :
    xx(n       ,:) =  xx(input) : Normal direction
    xx(n + NV  ,:) =  1st orthogonal vectors to xx(n,:)
    xx(n + 2*NV,:) =  2nd orthogonal vectors to xx(n,:)
               
 if xx == [] &  Lmode = 3 
    xx(n       ,:) =  [1 0 0]
    xx(n + NV  ,:) =  [0 1 0]
    xx(n + 2*NV,:) =  [0 0 1]
   
 written by M. Sato  2005-8-8
 modified by M. Sato  2006-7-14

 Copyright (C) 2011, ATR All Rights Reserved.
 License : New BSD License(see VBMEG_LICENSE.txt)

0001 function    [V,xx] = vb_current_vector(V,xx,Lmode)
0002 % Calculate current direction
0003 %   [V,xx] = vb_current_vector(V,xx,Lmode)
0004 % --- INPUT
0005 % V         : Cortex vertex coordinate [NV , 3]
0006 % xx        : Normal direction of cortex [NV , 3]
0007 % V(n,:)    : 3-D coordinate   at n-th vertex
0008 % xx(n,:)   : Normal direction at n-th vertex
0009 % Lmode     : number of independent curent direction
0010 %       = 1 :  current vextor = xx(n,:)
0011 %       = 2 :  current vextor = orthogonal vectors to V(n,:)
0012 %       = 3 :  current vextor = xx(n,:) & orthogonal vectors to xx(n,:)
0013 %
0014 % --- OUTPUT
0015 %  V( NV*Lmode , 3) = repmat(V, [Lmode 1])
0016 % xx( NV*Lmode , 3)
0017 %
0018 % Lmode = 1 :
0019 %   xx = xx(input) : Normal direction
0020 %
0021 % Lmode = 2 :
0022 %   xx(n,:)      = 1st orthogonal vectors to V(n,:)
0023 %   xx(n + NV,:) = 2nd orthogonal vectors to V(n,:)
0024 %
0025 % Lmode = 3 :
0026 %    xx(n       ,:) =  xx(input) : Normal direction
0027 %    xx(n + NV  ,:) =  1st orthogonal vectors to xx(n,:)
0028 %    xx(n + 2*NV,:) =  2nd orthogonal vectors to xx(n,:)
0029 %
0030 % if xx == [] &  Lmode = 3
0031 %    xx(n       ,:) =  [1 0 0]
0032 %    xx(n + NV  ,:) =  [0 1 0]
0033 %    xx(n + 2*NV,:) =  [0 0 1]
0034 %
0035 % written by M. Sato  2005-8-8
0036 % modified by M. Sato  2006-7-14
0037 %
0038 % Copyright (C) 2011, ATR All Rights Reserved.
0039 % License : New BSD License(see VBMEG_LICENSE.txt)
0040 
0041 % Number of cortex points
0042 NV = size(V,1);
0043 
0044 switch    Lmode,
0045 case    1,
0046     return;
0047 case    2,
0048     % Transform into spherical coordinate
0049     [phi,theta,r] = cart2sph(V(:,1),V(:,2),V(:,3));
0050     
0051     % Calculate tangential vector to radial direction
0052     [x1,y1,z1] = sph2cart(phi+0.5*pi, zeros(NV,1),  1); 
0053     [x2,y2,z2] = sph2cart(  phi,      theta+0.5*pi, 1); 
0054     xx           = [ x1,y1,z1 ; x2,y2,z2];
0055 
0056 case    3,
0057     if isempty(xx),
0058         % Current direction vector
0059         xx = zeros( 3*NV, 3 ); 
0060         ix = 1:NV;
0061         xx(ix       ,1) = 1; % x : [1 0 0]
0062         xx(ix + NV  ,2) = 1; % y : [0 1 0]
0063         xx(ix + 2*NV,3) = 1; % z : [0 0 1]
0064     else
0065         % Transform into spherical coordinate
0066         [phi,theta,r] = cart2sph(xx(:,1),xx(:,2),xx(:,3));
0067         
0068         % Calculate tangential vector to normal direction xx
0069         [x1,y1,z1] = sph2cart(phi+0.5*pi, zeros(NV,1),  1); 
0070         [x2,y2,z2] = sph2cart(  phi,      theta+0.5*pi, 1); 
0071         xx           = [ xx; x1,y1,z1 ; x2,y2,z2];
0072     end
0073 end;
0074 
0075 % Replicate Vertex
0076 V  = repmat(V, [Lmode 1]);