function varargout=rnd2plm(lmcosi,meth)
% lmcosip=RND2PLM(lmcosi,meth)
%
% Randomizes a spherical-harmonic coefficient matrix
%
% INPUT:
%
% lmcosi    Input matrix listing l, m, Ccos and Csin
% meth      'kevin' a la Kevin (faster for low degrees)
%           'frederik' a la Frederik  (faster for high degrees)
%
% OUTPUT:
%
% lmcosip   Output matrix listing l, m, Ccos and Csin
%
% SEE ALSO: 
%
% PLM2RND
% 
% REQUIRES:
%
% RandOrthMat, from
% http://www.mathworks.com/matlabcentral/fileexchange/11783-randorthmat
%
% Written by Kevin Lewis, 02/21/2010
% Last modified by fjsimons-at-alum.mit.edu, 02/21/20100

% Supply the default, i.e. a low-resolution terrestrial gravity field
defval('lmcosi',...
       kindeks(rindeks(fralmanac('EGM96','SHM'),1:addmup(255)-3),1:4))
defval('meth','kevin')

switch meth 
  case 'kevin'
   tic
   l=lmcosi(:,1);
   m=lmcosi(:,2);
   co=lmcosi(:,3);
   si=lmcosi(:,4);
   for i=min(l):max(l)
     ind=find(l==i);
     lcoeffs=RandOrthMat(2*i+1)*[flipud(si(ind(2:end))) ; co(ind)];
     si2(ind(2:end))=flipud(lcoeffs(1:i));
     co2(ind)=lcoeffs(i+1:end);
   end
   lmcosip=[l m co2(:) si2(:)];
   toc
case 'frederik'
 tic
 lmcosip=[lmcosi(:,1:2) zeros(length(lmcosi),2)];
 % Loop over the degrees
 for l=lmcosi(1,1):lmcosi(end,1)
   % Extract the COSINE coefficients and their indices at increasing degree
   [Ccos,b,e]=shcos(lmcosi,l);
   % Extract the SINE coefficients and return them at decreasing degree
   Csin=lmcosi(e:-1:b+1,4);
   % Randomize by orthogonal multiplication
   lcoeffs=RandOrthMat(2*l+1)*[Csin ; Ccos];
   % Assign the SINE coefficients (m>0) to the right spots in the larger matrix
   % Note the trick with the sum to convert the empty to a zero if needed
   lmcosip(e:-1:b+1,4)=sum(lcoeffs(1:l),2);
   % Assign the COSINE (m<=0) coefficients to the right spots
   lmcosip(b:e,3)=lcoeffs(l+1:end);
 end
 toc
end

% Prepare output
varns={lmcosip};
varargout=varns(1:nargout);