function varargout=jxm2plm(fname)
% lmcosi=jxm2plm(fname)
% [hlm,dels,dems]=jxm2plm(fname)
%
% Converts a spherical harmonics file from Jerry X. Mitrovica 
%
% INPUT:
%
% fname     A complete file name string
%
% OUTPUT:
%
% lmcosi    A matrix of expansion coefficients in a format suitable for
%           plotting and analysis by PLOTPLM, PLM2XYZ, etc., OR
%
% hlm       The coefficients in a linear vector [indeks(lmcosi(:,3:4)',':')]
% dels      The spherical harmonic degrees [this would be lmcosi(:,1)]
% dems      The spherical harmonic orders [this would be lmcosi(:,2)]
%
% Last modified by fjsimons-at-alum.mit.edu, 07/23/2010

% Load the file
fid=fopen(fname,'r');
% Read the header
head=fgetl(fid);
LT=sscanf(head,'%f',2);
% Spherical harmonic bandwidth
L=LT(1);
% Time
T=LT(2);
% Read the coefficients
COF=fscanf(fid,'%f');
fclose(fid);
% Check that the length is as I expect it
difer(length(COF)-2*addmup(L)-2*(L/2+1),[],[],NaN)

% Now rearrange
mods=(L-[0:L]+1);
% The indices of the blank pairs
modu=cumsum(mods+mod(mods,2));
modb=modu(~~mod(mods,2));
blanx=sort([2*modb-1 2*modb]);
% Check these are indeed blanks
difer(COF(blanx),[],[],NaN)
% Remaining are the non-blanks
COF=skip(COF,blanx);
% Which should be the same number as expected
difer(length(COF)-2*addmup(L),[],[],NaN)

% Now split into the real and imaginary parts
C=COF(1:2:end);
S=COF(2:2:end);
% The first L+1 should be zero
difer(S(1:L+1),[],[],NaN)

% Create a blank array with FJS format
[dems,dels,mz,lmcosi]=addmon(L);
% Fake arrival to align with readshds
hlm=lmcosi(:,3:4);
modm=[cumsum([0 mods(1:end-1)])+1 addmup(L)+1];
% Populate this with the JXM coefficients
for m=0:L
  % Stick in the "cosine" coefficients
  hlm(addmup(m-1:L-1)+m+1,1)=C(modm(m+1):modm(m+2)-1);
  % Stick in the "sine" coefficients
  hlm(addmup(m-1:L-1)+m+1,2)=S(modm(m+1):modm(m+2)-1);
end

% Some conventions need to be realigned
CSC=(-1).^dems;
dom=sqrt(2-(dems==0));
hlm(:,1)=hlm(:,1).*CSC.*dom;
hlm(:,2)=hlm(:,2).*CSC.*dom;

% Preserve the natural order of things
lmcosi(:,3:4)=hlm;

% Fake arrival to align with readshds
hlm=hlm';
hlm=hlm(:);

if nargout==1
  varargout={lmcosi};
end
if nargout==3
  varargout{1}=hlm;
  varargout{2}=dels;
  varargout{3}=dems;
end