
%-------------------------------------------------------------------------
% read SHCs in binary format (bshc)
%-------------------------------------------------------------------------
% format is: nmin, nmax, whole C-block, whole S-block, with the SHCs
% within the blocks ascending in degree, then in order.
% further: storage as double (required for scientific numbers)
%          and in big endian

function [n m C S] = read_SHCs_binary_bshc_format2014(filename_in)

 fid = fopen(filename_in,'rb')
 [data, cnt] = fread(fid, inf,'double');
 fclose(fid); 
 
 disp('Min and max deg of bshc-file')
 % extract min and max SHC
 input_nmin = data(1)
 input_nmax = data(2)
 
 %work out number of SHC records for C and S
 
 n_down = ((input_nmin-1)+1)*((input_nmin-1)+2)/2;
 n_up = ( input_nmax+1)*( input_nmax+2)/2;
 n_rows = n_up-n_down;
 offset = 2; % for correct computation of indices to access C, S blocks
 
 % access the C and S data blocks
 C = data(offset+1:offset+n_rows);
 S = data(offset+n_rows+1:offset+2*n_rows);
 
 % build new n, m vectors
 
 % build matrices with degrees and orders
 a = 0:input_nmax;          % vector from min to to max_degree
 b = a'*ones(1,input_nmax+1); % extend this into matrix
 degree_matrix = tril(b);            % set upper triangle to zero
 order_matrix = tril(b');
 clear b n m

% disp('Build n m vectors')
 ix1 = 1; ix2 = 1+input_nmin;
 
 for i = input_nmin+1: input_nmax+1
   % i
    n(ix1:ix2) = degree_matrix(i,1:i);
    m(ix1:ix2) = order_matrix(i,1:i);
    ix1 = ix2 +1; ix2 = ix1 + i;
 end
end