Matlab plot -- seasonly spatial lake surface temperature (fvcom/observation/difference)

fvcom

clear all;close all;clc;
filename=('interpolated_obs_HighRes.mat');

%% read and interpolate observation data and save interpolated data
plot_only=input('Plot only? Yes = 1, No = 0\n')

if (plot_only==0)
yr=2019;
day=365;
map=load('/Users/chuyanzhao/Desktop/FVCOM_output/Lake_Superior/boundary_superior.dat');

% read observation temperature
for dd=1:day
    dd
    ncload(['/Users/chuyanzhao/Desktop/gls_observation_data/obs_glsea_SST/2019_glsea_sst/' num2str(yr) '_' num2str(dd,'%.3i') '_glsea_sst.nc'],'sst');
    sst_glsea(dd,:,:)=sst;

end

% read observation grids
ncload(['/Users/chuyanzhao/Desktop/gls_observation_data/obs_glsea_SST/2019_glsea_sst/' num2str(yr) '_' num2str(1,'%.3i') '_glsea_sst.nc'],'lon','lat');

[xlon0,xlat0]=meshgrid(lon,lat);

% read fvcom grids
ncload(['/Users/chuyanzhao/Desktop/FVCOM_output/Lake_Superior/2019/raw_data/HighRes/gls_lst_2019_HighRes.nc'],'lon','lat','temp');
temp_mod=squeeze(temp(:,1,:));
xlon=lon-360;
xlat=lat;

% interpolate observation data to fvcom grids
for dd=1:day
    dd
    clear ind
    sst_glsea_tmp=squeeze(sst_glsea(dd,:,:));
    ind=find(sst_glsea_tmp>=0);% export the location of lake
    clear tmp
    tmp=interp_scatter_quick(xlon0(ind),xlat0(ind),sst_glsea_tmp(ind),xlon,xlat,0.5,1e-5,2,2);
    temp_obs(dd,:)=tmp(:)';

end

save(filename,'temp_obs','temp_mod','yr','xlon','xlat','map','filename');

end


%% plot seasonly fvcom spatial lake surface temperature

% read interpolated observation data
load(filename);

% seasonly fvcom temperature
size(temp_mod)
lst_mod=mean_by_segments(temp_mod,[0 2160 4344 6552 8761],1);
size(lst_mod)

% plot
cmin=0;cinc=1;cmax=20;
ncc=round([cmax-cmin]/cinc);
ccmap=flipud(spectral(ncc));

fontsize=8;

figure(101);clf;hold on;set_portrait;
    hh101=subplot(4,1,1);cla;hold on;box on;
    hh101_pp=get(hh101,'position');
    set(hh101,'position',hh101_pp+[0.1 0 -0.375 -0.0077]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_mod(1,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh101,ccmap);
    plot_fvcom_obc(map,[210 180 140]/255);
    axis([-92.2 -84.4 46.2 49.2]);
    lcb=colorbar('location','eastoutside');
    ylabel(lcb,'LST (\circC)');
    text(-92,49,'Seasonly mean LST -- Winter','FontSize',fontsize);
    text(-92,46.4,'FVCOM','FontSize',fontsize);
    set(gca,'FontSize',fontsize);

    %
    hh102=subplot(4,1,2);cla;hold on;box on;
    hh102_pp=get(hh102,'position');
    set(hh102,'position',hh102_pp+[0.1 0.04 -0.375 -0.0077]);%set figure 
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_mod(2,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh102,ccmap);
    plot_fvcom_obc(map,[210 180 140]/255);
    axis([-92.2 -84.4 46.2 49.2]);
    lcb=colorbar('location','eastoutside');
    ylabel(lcb,'LST (\circC)');
    text(-92,49,'Seasonly mean LST -- Spring','FontSize',fontsize);
    text(-92,46.4,'FVCOM','FontSize',fontsize);
    set(gca,'FontSize',fontsize);

    %
    hh103=subplot(4,1,3);cla;hold on;box on;
    hh103_pp=get(hh103,'position');
    set(hh103,'position',hh103_pp+[0.1 0.08 -0.375 -0.0077]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_mod(3,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh103,ccmap);
    plot_fvcom_obc(map,[210 180 140]/255);
    axis([-92.2 -84.4 46.2 49.2]);
    lcb=colorbar('location','eastoutside');
    ylabel(lcb,'LST (\circC)');
    text(-92,49,'Seasonly mean LST -- Summer','FontSize',fontsize);
    text(-92,46.4,'FVCOM','FontSize',fontsize);
    set(gca,'FontSize',fontsize);

    %
    hh104=subplot(4,1,4);cla;hold on;box on;
    hh104_pp=get(hh104,'position');
    set(hh104,'position',hh104_pp+[0.1 0.12 -0.375 -0.0077]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_mod(4,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh104,ccmap);
    plot_fvcom_obc(map,[210 180 140]/255);
    axis([-92.2 -84.4 46.2 49.2]);
    lcb=colorbar('location','eastoutside');
    ylabel(lcb,'LST (\circC)');
    text(-92,49,'Seasonly mean LST -- Fall','FontSize',fontsize);
    text(-92,46.4,'FVCOM','FontSize',fontsize);
    set(gca,'FontSize',fontsize);

    hh101.Colormap=ccmap;
    hh102.Colormap=ccmap;
    hh103.Colormap=ccmap;
    hh104.Colormap=ccmap;

    saveas(gcf,['/Users/chuyanzhao/Desktop/FVCOM_output/Lake_Superior/2019/figure/HighRes/fig_lake_lst_seanonly_mod_' num2str(yr) '.png'],'png');  

%% plot seasonly observation spatial lake surface temperature

% seansonly observation temperature
clear lst_obs
size(temp_obs)
lst_obs=mean_by_segments(temp_obs,[0 90 181 273 365],1);
size(lst_obs)

cmin=0;cinc=1;cmax=20;
ncc=round([cmax-cmin]/cinc);
ccmap=flipud(spectral(ncc));

% plot
figure(201);clf;hold on;set_portrait;
    hh201=subplot(4,1,1);cla;hold on;box on;
    hh201_pp=get(hh201,'position');
    set(hh201,'position',hh201_pp+[0.1 0 -0.375 -0.0077]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_obs(1,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh201,ccmap);
    plot_fvcom_obc(map,[210 180 140]/255);
    axis([-92.2 -84.4 46.2 49.2]);
    lcb=colorbar('location','eastoutside');
    ylabel(lcb,'LST (\circC)');
    text(-92,49,'Seasonly mean LST -- Winter');
    text(-92,46.4,'OBS');
    set(gca,'FontSize',fontsize);

    hh202=subplot(4,1,2);cla;hold on;box on;
    hh202_pp=get(hh202,'position');
    set(hh202,'position',hh202_pp+[0.1 0.04 -0.375 -0.0077]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_obs(2,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh202,ccmap);
    plot_fvcom_obc(map,[210 180 140]/255);
    axis([-92.2 -84.4 46.2 49.2]);
    lcb=colorbar('location','eastoutside');
    ylabel(lcb,'LST (\circC)');
    text(-92,49,'Seasonly mean LST -- Spring');
    text(-92,46.4,'OBS');
    set(gca,'FontSize',fontsize);

    %
    hh203=subplot(4,1,3);cla;hold on;box on;
    hh203_pp=get(hh203,'position');
    set(hh203,'position',hh203_pp+[0.1 0.08 -0.375 -0.0077]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_obs(3,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh203,ccmap);
    plot_fvcom_obc(map,[210 180 140]/255);
    axis([-92.2 -84.4 46.2 49.2]);
    lcb=colorbar('location','eastoutside');
    ylabel(lcb,'LST (\circC)');
    text(-92,49,'Seasonly mean LST -- Summer');
    text(-92,46.4,'OBS');
    set(gca,'FontSize',fontsize);

    %
    hh204=subplot(4,1,4);cla;hold on;box on;
    hh204_pp=get(hh204,'position');
    set(hh204,'position',hh104_pp+[0.1 0.12 -0.375 -0.0077]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_obs(4,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh204,ccmap);
    plot_fvcom_obc(map,[210 180 140]/255);
    axis([-92.2 -84.4 46.2 49.2]);
    lcb=colorbar('location','eastoutside');
    ylabel(lcb,'LST (\circC)');
    text(-92,49,'Seasonly mean LST -- Fall');
    text(-92,46.4,'OBS');
    set(gca,'FontSize',fontsize);

    hh201.Colormap=ccmap;
    hh202.Colormap=ccmap;
    hh203.Colormap=ccmap;
    hh204.Colormap=ccmap;


    saveas(gcf,['/Users/chuyanzhao/Desktop/FVCOM_output/Lake_Superior/2019/figure/HighRes/fig_lake_lst_seasonly_obs_' num2str(yr) '.png'],'png');  

%% plot seasonly difference between observation and fvcom

% calculate the difference between observation and fvcom
lst_diff=lst_obs-lst_mod;

% plot
color_scale1=[
103,0,31
178,24,43
214,96,77
244,165,130
253,219,199
247,247,247
247,247,247
209,229,240
146,197,222
67,147,195
33,102,172
5,48,97]/255;

cmin=-6;cinc=1;cmax=6;
ccmap=flipud(color_scale1);

figure(301);clf;hold on;set_portrait;
    hh301=subplot(4,1,1);cla;hold on;box on;
    hh301_pp=get(hh301,'position');
    set(hh301,'position',hh301_pp+[0.1 0 -0.375 -0.0077]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_diff(1,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh301,ccmap);
    plot_fvcom_obc(map,[210 180 140]/255);
    axis([-92.2 -84.4 46.2 49.2]);
    lcb=colorbar('location','eastoutside');
    ylabel(lcb,'LST (\circC)');
    text(-92,49,'Seasonly mean LST difference -- Winter');
    text(-92,46.4,'Difference');
    set(gca,'FontSize',fontsize);

    hh302=subplot(4,1,2);cla;hold on;box on;
    hh302_pp=get(hh302,'position');
    set(hh302,'position',hh302_pp+[0.1 0.04 -0.375 -0.0077]); 
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_diff(2,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh302,ccmap);
    plot_fvcom_obc(map,[210 180 140]/255);
    axis([-92.2 -84.4 46.2 49.2]);
    lcb=colorbar('location','eastoutside');
    ylabel(lcb,'LST (\circC)');
    text(-92,49,'Seasonly mean LST difference -- Spring');
    text(-92,46.4,'Difference');
    set(gca,'FontSize',fontsize);

    %
    hh303=subplot(4,1,3);cla;hold on;box on;
    hh303_pp=get(hh303,'position');
    set(hh303,'position',hh303_pp+[0.1 0.08 -0.375 -0.0077]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_diff(3,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh303,ccmap);
    plot_fvcom_obc(map,[210 180 140]/255);
    axis([-92.2 -84.4 46.2 49.2]);
    lcb=colorbar('location','eastoutside');
    ylabel(lcb,'LST (\circC)');
    text(-92,49,'Seanonly mean LST difference -- Summer');
    text(-92,46.4,'Difference');
    set(gca,'FontSize',fontsize);

    %
    hh304=subplot(4,1,4);cla;hold on;box on;
    hh304_pp=get(hh304,'position');
    set(hh304,'position',hh304_pp+[0.1 0.12 -0.375 -0.0077]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_diff(4,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh304,ccmap);
    plot_fvcom_obc(map,[210 180 140]/255);
    axis([-92.2 -84.4 46.2 49.2]);
    lcb=colorbar('location','eastoutside');
    ylabel(lcb,'LST (\circC)');
    text(-92,49,'Seanonly mean LST difference -- Fall');
    text(-92,46.4,'Difference');
    set(gca,'FontSize',fontsize);

    hh301.Colormap=ccmap;
    hh302.Colormap=ccmap;
    hh303.Colormap=ccmap;
    hh304.Colormap=ccmap;

    saveas(gcf,['/Users/chuyanzhao/Desktop/FVCOM_output/Lake_Superior/2019/figure/HighRes/fig_lake_lst_seanonly_obs-mod-diff_' num2str(yr) '.png'],'png');