clc; clear all;
% % % =====================================================================
% % % load the data
save('data_trig_int_example','lon4interp','lat4interp','den_inp','R1inp','R2inp','dSp','G','fntl','ngl','ilays','nres','Rh','tlon','tlat','tg','ief','iefs','layss','csys','units');    

% % % =====================================================================
% % % plot the triangle elements
S3m=max(dSp(:));

fg=figure('visible','off');
scatter(lon4interp(1:end),lat4interp(1:end),6,dSp/S3m*100,'filled'); colorbar;
% title(['Normalized area of each spherical triangle, in %']);
set(findall(gcf,'Type','text'),'FontSize',10);
set(findall(gcf,'Type','axes'),'FontSize',10);
grid on;
set(gcf, 'PaperUnits', 'inches', 'PaperPosition', [0 0 1600 900]/100); % 'renderer','painters'
print (fg, '-dpng','-r100', 'Trig_grid_area');
close all;  

% % % =====================================================================
% % % interate 2 layers in L8
[v,vx,vy,vz,vxx,vxy,vxz,vyy,vyz,vzz,slay,sss,lonv,latv] = newton_trig_int_ver1(lon4interp,lat4interp,den_inp,R1inp,R2inp,dSp,G,fntl,...
                                   ngl,ilays,nres,Rh,tlon,tlat,tg,ief,iefs,layss,csys,units);
[lont,latt]= meshgrid(lonv,latv);

% % % =====================================================================
% % % plot the result
fg=figure('visible','on'); h=pcolor(lont,latt,vz(:,:,1));colorbar;set(h,'linestyle','none');hcl=colorbar;% set(gca,'clim',[-1 1]);%pp=1e8; colormap(cmap); shading flat;