%-------------------------------------------------------------------------
% ReadTest_SpectralTopographyModel
%-------------------------------------------------------------------------
% Christian Hirt, TU Munich, & Michael Kuhn, Curtin University
%
% This archive contains the (2pi) fully-normalized spherical harmonic 
% coefficients (SHCs) of the topography (shape) model of the Moon
%
% MoonTopo2600pa_shape
%
% taken from Mark Wieczorek's model archive. The file is part of our 
% archive because it was used as input topography model both in the
% * spatial-domain gravity forward modelling and 
% * spectral-domain gravity forward modelling 

% The MoonTopo2600pa_shape is based on LOLA topography, given in the
% principal axis coordinate system of the Moon. The file provides radii 
% of the Moon's surface relative to the selenocentre.

% Further details on the use of the  MoonTopo2600pa_shape model to derive 
% gravity from topography are found in the publication
%
% Hirt, C. and M. Kuhn (2017), Convergence and divergence in 
% spherical harmonic series of the gravitational field generated by 
% high-resolution planetary topography – a case study for the Moon. 
% Journal of Geophysical Research – Planets, 122, doi:10.1002/2017JE005298.
%
% and the development of the MoonTopo2600pa_shape model is described in 
%
% Wieczorek,  M.A.  (2015),  Gravity  and  Topography  of  the  Terrestrial  
% Planets.  In:  Schubert  G  (ed)  Treatise  on Geophysics, 2nd edition, 
% Oxford, 153?193, doi:10.1016/B978?0?444?53802?4.00169?X.
%
% ------------------------------------------------------------------------
% Christian Hirt, 02 August 2017
% ------------------------------------------------------------------------

clear all
close all

% six spectral gravity forward models ------------------------------------
fileName = 'MoonTopo2600pa_shape.bshc'
     
% read the SHCs from the file --------------------------------------------
[n m C S]      = read_SHCs_binary_bshc_format2014(fileName);
% now the SHCs are loaded into Matlab and can be further processed.

% compute dimensionless potential degree variances -----------------------
[n_vec sigma2] = degree_variances2012(n,m,C,S,2600);

% plot the degree variances ----------------------------------------------
figure(1);
plot(n_vec,sigma2);
grid on
set(gca,'XMinorTick','off','YMinorTick','off');
set(gca,'XMinorGrid','off','YMinorGrid','off');
set(gca,'XScale','linear','YScale','log');
xlabel('Harmonic degree');
ylabel('Degree variance');
title(['Degree variances of  SHCs from file ' strrep(fileName,'_','-')] )
print(strcat('Degree_variances_',fileName,'.png'),'-r300','-dpng')
 
    
  