clear all addpath(genpath('/data/home/yxzhang/DarwinFuture/library')); addpath(genpath('/data/home/yxzhang/NASA')); addpath(genpath('/data/home/ppwu/draw')); load("grid.mat"); load("mycmap.mat"); lon_144x90=lon_144x90-180; landareas=shaperead('ne_10m_land.shp'); geodir = '/data/home/yxzhang/IGSM/regrid/'; ar =ncread(strcat(geodir,'grid_igsm.nc'),'rA'); % grid area wet=ncread(strcat(geodir,'grid_igsm.nc'),'HFacC'); % fraction of water in grid (vertically) delz=[10, 15, 20, 30, 40, ... 55, 70, 90, 110, 135, ... 160, 190, 225, 265, 310, ... 360, 415, 475, 540, 610, ... 685, 765]; tvolume=zeros(144,90,22); for level=1:22 tvolume(:,:,level)=ar*delz(level); end wvolume=tvolume.*wet; depth={'10m','25m','45m','75m','115m', ... '170m','240m','330m','440m','575m', ... '735m','925m','1150m','1415m','1725m', ... '2085m','2500m','2975m','3515m','4125m', ... '4810m','5575m'}; density=[950,1048,1146,950,1048,1146,950,1048,1146, ... 950,1048,1146,950,1048,1146,950,1048,1146, ... 900,1048,1196,900,1048,1196,900,1048,1196, ... 900,1048,1196,900,1048,1196,900,1048,1196, ... 1410,1410,1410,1410,1410,1410, ... 1050,1050,1050,1050,1050,1050, ... 1050,1050,1050,1050,1050,1050, ... 550,550,550,550,550,550]; Dir_Prefix = '/data/home/ympeng/MITgcm-Plastic/verification/global_hg_igsm/run_Covid19_optimum_microp_macrop/'; yy=[1,5,80]; ww=0.22; ll=0.185; x1=0.11; x2=0.35; x3=0.59; y1=0.19; y2=0.365; y3=0.54; year=linspace(349920,520560,80); year1=linspace(2021,2100,80); % % -------------- data reading and processing -------------- % % (kg m-3) for yr=yy %% Surface (kg m-3) data_s=zeros(144,90,22,40); for n=[1:12 19:30 37:40 43:46 49:52 55:58] dataf_s=rdmds(strcat(Dir_Prefix,'PTRtave',num2str(n,'%02d')),year(yr)); dataf_s(dataf_s<0.0)=0.0; data_s(:,:,:,n) = dataf_s.*wvolume; end data2=zeros(144,90); for n=[1:12 19:30 37:40 43:46 49:52 55:58] if(density(n)<1048) data2(:,:)=data2(:,:)+data_s(:,:,1,n); else data2(:,:)=data2(:,:)+data_s(:,:,1,n)*0.05; end end dataf2=data2./ar*1000000*1000; data3=log10(dataf2); data3(data3==-Inf)=NaN; u=zeros(144,90); u(:,:)=data3(:,:); pdata=zeros(90,144); pdata=u(:,:)'; pdata2(:,1:72)=pdata(:,73:144); pdata2(:,73:144)=pdata(:,1:72); % p1=figure(1); set(gcf,'Color',[1 1 1]); if yr==1 subplot(3,3,1) pcolor(lon_144x90,lat_144x90,pdata2); set(gca,'position',[x1 y3 ww ll]) elseif yr==5 subplot(3,3,2) pcolor(lon_144x90,lat_144x90,pdata2); set(gca,'position',[x2 y3 ww ll]) elseif yr==80 subplot(3,3,3) pcolor(lon_144x90,lat_144x90,pdata2); set(gca,'position',[x3 y3 ww ll]) end shading flat colormap(widl) axesm('mercator','Frame','off','Grid','off',... 'origin',[0 180 0]); mapshow(landareas,'FaceColor',[0.45 0.45 0.45],... 'EdgeColor',[0.45 0.45 0.45]); grid off caxis([-4,2]); hcb1=colorbar('position',[0.82 0.558 0.01 0.15],... 'FontWeight','Bold','FontName',... 'Liberation Serif','FontSize',8, ... 'Ticks',-4:2:2,'YTickLabel',... {'10^{-4}','10^{-2}','10^{0}','10^{2}'}); set(get(hcb1,'Title'),'string','g km^{-2}',... 'FontWeight','Bold') %% Beach (kg m-3) data_b=zeros(144,90,40); for n=[1:12 19:30 37:40 43:46 49:52 55:58] dataf_bc=rdmds(strcat(Dir_Prefix,'diags/Beach',num2str(n,'%02d')),year(yr)); % kg m-3 dataf_bc(dataf_bc<0.0)=0.0; data_b(:,:,n) = dataf_bc.*delz(1); end data2=squeeze(nansum(data_b*1000000*1000,3)); data3=log10(data2); data3(data3==-Inf)=NaN; u=zeros(144,90); u(:,:)=data3(:,:); pdata=zeros(90,144); pdata=u(:,:)'; pdata2(:,1:72)=pdata(:,73:144); pdata2(:,73:144)=pdata(:,1:72); % figure set(gcf,'Color',[1 1 1]); if yr==1 subplot(3,2,4) pcolor(lon_144x90,lat_144x90,pdata2); set(gca,'position',[x1 y2 ww ll]) elseif yr==5 subplot(3,3,5) pcolor(lon_144x90,lat_144x90,pdata2); set(gca,'position',[x2 y2 ww ll]) elseif yr==80 subplot(3,3,6) pcolor(lon_144x90,lat_144x90,pdata2); set(gca,'position',[x3 y2 ww ll]) end shading flat colormap(widl) axesm('mercator','Frame','off','Grid','off',... 'origin',[0 180 0]); mapshow(landareas,'FaceColor',[0.45 0.45 0.45],... 'EdgeColor',[0.45 0.45 0.45]); grid off caxis([-2,4]); hcb2=colorbar('position',[0.82 0.381 0.01 0.15],... 'FontWeight','Bold','FontName',... 'Liberation Serif','FontSize',8, ... 'Ticks',-2:2:4,'YTickLabel',... {'10^{-2}','10^{0}','10^{2}','10^{4}'}); %% Bottom (kg m-2) data_bt=zeros(144,90,40); for n=[1:12 19:30 37:40 43:46 49:52 55:58] dataf_bt=rdmds(strcat(Dir_Prefix,'diags/Bottom',num2str(n,'%02d')),year(yr)); % kg m-2 dataf_bt(dataf_bt<0.0)=0.0; data_bt(:,:,n) = dataf_bt; end data2=squeeze(nansum(data_bt*1000000*1000,3)); data3 = log10(data2); data3(data3==-Inf)=NaN; u=zeros(144,90); u(:,:)=data3(:,:); pdata=zeros(90,144); pdata=u(:,:)'; pdata2(:,1:72)=pdata(:,73:144); pdata2(:,73:144)=pdata(:,1:72); % figure set(gcf,'Color',[1 1 1]); if yr==1 subplot(3,3,7) pcolor(lon_144x90,lat_144x90,pdata2); set(gca,'position',[x1 y1 ww ll]) elseif yr==5 subplot(3,3,8) pcolor(lon_144x90,lat_144x90,pdata2); set(gca,'position',[x2 y1 ww ll]) elseif yr==80 subplot(3,3,9) pcolor(lon_144x90,lat_144x90,pdata2); set(gca,'position',[x3 y1 ww ll]) end shading flat colormap(widl) axesm('mercator','Frame','off','Grid','off',... 'origin',[0 180 0]); mapshow(landareas,'FaceColor',[0.45 0.45 0.45],... 'EdgeColor',[0.45 0.45 0.45]); grid off caxis([-14,6]); hcb3=colorbar('position',[0.82 0.208 0.01 0.15],... 'FontWeight','Bold','FontName',... 'Liberation Serif','FontSize',8, ... 'Ticks',-14:5:6,'YTickLabel',... {'10^{-14}','10^{-9}','10^{-4}','10^{1}','10^{6}'}); end hAxis=axes('visible','off'); h1=text(-0.045,0.64,'Surface'); set(h1,'fontsize',10,'FontName','Liberation Serif',... 'FontWeight','Bold','rotation',90,... 'HorizontalAlignment','center') h2=text(-0.045,0.43,'Beach'); set(h2,'fontsize',10,'FontName','Liberation Serif',... 'FontWeight','Bold','rotation',90,... 'HorizontalAlignment','center') h3=text(-0.045,0.21,'Sediment'); set(h3,'fontsize',10,'FontName','Liberation Serif',... 'FontWeight','Bold','rotation',90,... 'HorizontalAlignment','center') t1=text(-0.012,0.58,'a'); set(t1,'fontsize',9.5,'fontname','Liberation Serif',... 'FontWeight','Bold'); t2=text(0.295,0.58,'b'); set(t2,'fontsize',9.5,'fontname','Liberation Serif',... 'FontWeight','Bold'); t3=text(0.605,0.58,'c'); set(t3,'fontsize',9.5,'fontname','Liberation Serif',... 'FontWeight','Bold'); t4=text(-0.012,0.36,'d'); set(t4,'fontsize',9.5,'fontname','Liberation Serif',... 'FontWeight','Bold'); t5=text(0.295,0.36,'e'); set(t5,'fontsize',9.5,'fontname','Liberation Serif',... 'FontWeight','Bold'); t6=text(0.605,0.36,'f'); set(t6,'fontsize',9.5,'fontname','Liberation Serif',... 'FontWeight','Bold'); t7=text(-0.012,0.155,'g'); set(t7,'fontsize',9.5,'fontname','Liberation Serif',... 'FontWeight','Bold'); t8=text(0.295,0.155,'h'); set(t8,'fontsize',9.5,'fontname','Liberation Serif',... 'FontWeight','Bold'); t9=text(0.605,0.155,'i'); set(t9,'fontsize',9.5,'fontname','Liberation Serif',... 'FontWeight','Bold'); t10=text(0.19,0.15,'2021'); set(t10,'fontsize',9.5,'fontname','Liberation Serif',... 'FontWeight','Bold'); t11=text(0.495,0.15,'2025'); set(t11,'fontsize',9.5,'fontname','Liberation Serif',... 'FontWeight','Bold'); t12=text(0.807,0.15,'2100'); set(t12,'fontsize',9.5,'fontname','Liberation Serif',... 'FontWeight','Bold'); tite('microplastic','fontsize',12,'fontname','Liberation Serif',... 'FontWeight','Bold'); print(gcf,'-dpng','-r300',... strcat('/data/home/ympeng/MITgcm-Plastic/draw/Fig5.png'));