clear

fig_ps=[752 56 1537 1038];
cm=brewermap(64,'*Spectral'); 

% NSIDC SIarea
nnx=304; nny=448;
nlat=readbin('NSIDC_IceCon/psn25lats.dat',[nnx nny],1,'int32')/100000;
nlon=readbin('NSIDC_IceCon/psn25lons.dat',[nnx nny],1,'int32')/100000;
ix=find(nlon>180);
nlon(ix)=nlon(ix)-360;
% read southern data grid coordinates
snx=316; sny=332;
slat=readbin('NSIDC_IceCon/pss25lats.dat',[snx sny],1,'int32')/100000;
slon=readbin('NSIDC_IceCon/pss25lons.dat',[snx sny],1,'int32')/100000;
ix=find(slon>180);
slon(ix)=slon(ix)-360;

%fn='NSIDC_IceCon/NSIDC_LLC_SIarea.mat';
fn='NSIDC_IceCon/NSIDC_LLC_SIarea_2020.mat';
load(fn)


yr2012=2012;
for mn=[3 9 12]
if mn==12
day1=1; day2=365;
else
day1=datenum(yr2012,mn,1)-datenum(yr2012,1,1)+1;
day2=datenum(yr2012,mn+1,1)-datenum(yr2012,1,1);
end
dd=day1:day2;

%sst_mur=nanmean(nic_nsidc(:,:,dd),3);
%sst_llc=nanmean(nic_llc  (:,:,dd),3);
%sst_DIF=nic_llc(:,:,dd)-nic_nsidc(:,:,dd);
sst_mur=nanmean(sic_nsidc(:,:,dd),3);
sst_llc=nanmean(sic_llc  (:,:,dd),3);
sst_DIF=sic_llc(:,:,dd)-sic_nsidc(:,:,dd);
sst_dif=nanmean(sst_DIF,3);
sst_rms=sqrt(nanmean(sst_DIF.^2,3));


figure(mn)
set(gcf,'Position',fig_ps)

subplot(221)
%m_proj('stereo','lat',90,'lon',0,'rad',40)
%m_pcolor(nlon,nlat,sst_mur)
m_proj('stereo','lat',-90,'lon',0,'rad',40)
m_pcolor(slon,slat,sst_mur)

caxis([0 1])
shading flat,thincb(1);
%m_grid('xtick',-180:30:180,'ytick',40:10:80)
m_grid('xtick',-180:30:180,'ytick',-80:10:-40,'Xaxislocation','top','Yaxislocation','middle')
title(['NSIDC SIarea'])

subplot(222)
%m_proj('stereo','lat',90,'lon',0,'rad',40)
%m_pcolor(nlon,nlat,sst_llc)
m_proj('stereo','lat',-90,'lon',0,'rad',40)
m_pcolor(slon,slat,sst_llc)
caxis([0 1])
shading flat,thincb(1);
%m_grid('xtick',-180:30:180,'ytick',40:10:80)
m_grid('xtick',-180:30:180,'ytick',-80:10:-40,'Xaxislocation','top','Yaxislocation','middle')
title(['LLC2160 SIarea'])

cx=subplot(223);
%m_proj('stereo','lat',90,'lon',0,'rad',40)
%m_pcolor(nlon,nlat,sst_dif)
m_proj('stereo','lat',-90,'lon',0,'rad',40)
m_pcolor(slon,slat,sst_dif)
caxis([-1 1]*1)
shading flat,thincb(1);
%m_grid('xtick',-180:30:180,'ytick',40:10:80)
m_grid('xtick',-180:30:180,'ytick',-80:10:-40,'Xaxislocation','top','Yaxislocation','middle')
colormap(cx,bluewhitered)
title(['mean LLC2160 - NSIDC SIarea'])

cx=subplot(224);
%m_proj('stereo','lat',90,'lon',0,'rad',40)
%m_pcolor(nlon,nlat,sst_rms)
m_proj('stereo','lat',-90,'lon',0,'rad',40)
m_pcolor(slon,slat,sst_rms)
caxis([0 1]);thincb(1);
colormap(cx,jet)
%m_grid('xtick',-180:30:180,'ytick',40:10:80)
m_grid('xtick',-180:30:180,'ytick',-80:10:-40,'Xaxislocation','top','Yaxislocation','middle')
title(['rms LLC2160 - NSIDC SIarea'])

%print('-dpng',['fig_NSIDC_SIarea'])

end