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

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 monthly fvcom spatial lake surface temperature

% read interpolated observation data
load(filename);

% monthly fvcom temperature
clear lst_mod
size(temp_mod)
lst_mod=mean_by_segments(temp_mod,[0 744 1416 2160 2880 3624 4344 5088 5832 6552 7296 8016 8761],1);
size(lst_mod)

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

figure(101);clf;hold on;set_portrait;
    hh101=subplot(3,1,1);cla;hold on;box on;
    hh101_pp=get(hh101,'position');
    set(hh101,'position',hh101_pp+[0.1 0 -0.2 0]);
    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,'Monthly mean LST -- January');
    text(-92,46.4,'FVCOM');

    %
    hh102=subplot(3,1,2);cla;hold on;box on;
    hh102_pp=get(hh102,'position');
    set(hh102,'position',hh102_pp+[0.1 0.05 -0.2 0]);%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,'Monthly mean LST -- Feburary');
    text(-92,46.4,'FVCOM');

    %
    hh103=subplot(3,1,3);cla;hold on;box on;
    hh103_pp=get(hh103,'position');
    set(hh103,'position',hh103_pp+[0.1 0.1 -0.2 0]);
    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,'Monthly mean LST -- March');
    text(-92,46.4,'FVCOM');

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

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

figure(201);clf;hold on;set_portrait;
    hh201=subplot(3,1,1);cla;hold on;box on;
    hh201_pp=get(hh201,'position');
    set(hh201,'position',hh201_pp+[0.1 0 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_mod(4,:)',[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,'Monthly mean LST -- April');
    text(-92,46.4,'FVCOM');

    %
    hh202=subplot(3,1,2);cla;hold on;box on;
    hh202_pp=get(hh202,'position');
    set(hh202,'position',hh202_pp+[0.1 0.05 -0.2 0]);%set figure 
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_mod(5,:)',[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,'Monthly mean LST -- May');
    text(-92,46.4,'FVCOM');

    %
    hh203=subplot(3,1,3);cla;hold on;box on;
    hh203_pp=get(hh203,'position');
    set(hh203,'position',hh203_pp+[0.1 0.1 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_mod(6,:)',[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,'Monthly mean LST -- June');
    text(-92,46.4,'FVCOM');

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

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

figure(301);clf;hold on;set_portrait;
    hh301=subplot(3,1,1);cla;hold on;box on;
    hh301_pp=get(hh301,'position');
    set(hh301,'position',hh301_pp+[0.1 0 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_mod(7,:)',[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,'Monthly mean LST -- July');
    text(-92,46.4,'FVCOM');

    %
    hh302=subplot(3,1,2);cla;hold on;box on;
    hh302_pp=get(hh302,'position');
    set(hh302,'position',hh302_pp+[0.1 0.05 -0.2 0]);%set figure 
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_mod(8,:)',[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,'Monthly mean LST -- August');
    text(-92,46.4,'FVCOM');

    %
    hh303=subplot(3,1,3);cla;hold on;box on;
    hh303_pp=get(hh303,'position');
    set(hh303,'position',hh303_pp+[0.1 0.1 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_mod(9,:)',[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,'Monthly mean LST -- September');
    text(-92,46.4,'FVCOM');

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

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

figure(401);clf;hold on;set_portrait;
    hh401=subplot(3,1,1);cla;hold on;box on;
    hh401_pp=get(hh401,'position');
    set(hh401,'position',hh401_pp+[0.1 0 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_mod(10,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh401,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,'Monthly mean LST -- October');
    text(-92,46.4,'FVCOM');
    
    %
    hh402=subplot(3,1,2);cla;hold on;box on;
    hh402_pp=get(hh402,'position');
    set(hh402,'position',hh402_pp+[0.1 0.05 -0.2 0]);%set figure 
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_mod(11,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh402,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,'Monthly mean LST -- November');
    text(-92,46.4,'FVCOM');

    %
    hh403=subplot(3,1,3);cla;hold on;box on;
    hh403_pp=get(hh403,'position');
    set(hh403,'position',hh403_pp+[0.1 0.1 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_mod(12,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh403,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,'Monthly mean LST -- December');
    text(-92,46.4,'FVCOM');

    hh401.Colormap=ccmap;
    hh402.Colormap=ccmap;
    hh403.Colormap=ccmap;

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


%% plot monthly observation spatial lake surface temperature

clear lst_obs
size(temp_obs)
lst_obs=mean_by_segments(temp_obs,[0 31 59 90 120 151 181 212 243 273 304 334 365],1);
size(lst_obs)

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

figure(501);clf;hold on;set_portrait;
    hh501=subplot(3,1,1);cla;hold on;box on;
    hh501_pp=get(hh501,'position');
    set(hh501,'position',hh501_pp+[0.1 0 -0.2 0]);
    hh901_pp=get(hh501,'position');
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_obs(1,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh501,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,'Monthly mean LST -- January');
    text(-92,46.4,'OBS');

    hh502=subplot(3,1,2);cla;hold on;box on;
    hh502_pp=get(hh502,'position');
    set(hh502,'position',hh502_pp+[0.1 0.05 -0.2 0]);%set figure 
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_obs(2,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh502,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,'Monthly mean LST -- Feburary');
    text(-92,46.4,'OBS');

    %
    hh503=subplot(3,1,3);cla;hold on;box on;
    hh503_pp=get(hh503,'position');
    set(hh503,'position',hh503_pp+[0.1 0.1 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_obs(3,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh503,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,'Monthly mean LST -- March');
    text(-92,46.4,'OBS');

    hh501.Colormap=ccmap;
    hh502.Colormap=ccmap;
    hh503.Colormap=ccmap;

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

figure(601);clf;hold on;set_portrait;
    hh601=subplot(3,1,1);cla;hold on;box on;
    hh601_pp=get(hh601,'position');
    set(hh601,'position',hh601_pp+[0.1 0 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_obs(4,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh601,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,'Monthly mean LST -- April');
    text(-92,46.4,'OBS');
    % colormap(ax(1), ccmap)

    %
    hh602=subplot(3,1,2);cla;hold on;box on;
    hh602_pp=get(hh602,'position');
    set(hh602,'position',hh602_pp+[0.1 0.05 -0.2 0]);%set figure 
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_obs(5,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh602,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,'Monthly mean LST -- May');
    text(-92,46.4,'OBS');

    %
    hh603=subplot(3,1,3);cla;hold on;box on;
    hh603_pp=get(hh603,'position');
    set(hh603,'position',hh603_pp+[0.1 0.1 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_obs(6,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh603,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,'Monthly mean LST -- June');
    text(-92,46.4,'OBS');

    hh601.Colormap=ccmap;
    hh602.Colormap=ccmap;
    hh603.Colormap=ccmap;

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

figure(701);clf;hold on;set_portrait;
    hh701=subplot(3,1,1);cla;hold on;box on;
    hh701_pp=get(hh701,'position');
    set(hh701,'position',hh701_pp+[0.1 0 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_obs(7,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh701,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,'Monthly mean LST -- July');
    text(-92,46.4,'OBS');
    % colormap(ax(1), ccmap)

    %
    hh702=subplot(3,1,2);cla;hold on;box on;
    hh702_pp=get(hh702,'position');
    set(hh702,'position',hh702_pp+[0.1 0.05 -0.2 0]);%set figure 
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_obs(8,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh702,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,'Monthly mean LST -- August');
    text(-92,46.4,'OBS');

    %
    hh703=subplot(3,1,3);cla;hold on;box on;
    hh703_pp=get(hh703,'position');
    set(hh703,'position',hh703_pp+[0.1 0.1 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_obs(9,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh703,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,'Monthly mean LST -- September');
    text(-92,46.4,'OBS');

    hh701.Colormap=ccmap;
    hh702.Colormap=ccmap;
    hh703.Colormap=ccmap;

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

figure(801);clf;hold on;set_portrait;
    hh801=subplot(3,1,1);cla;hold on;box on;
    hh801_pp=get(hh801,'position');
    set(hh801,'position',hh801_pp+[0.1 0 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_obs(10,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh801,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,'Monthly mean LST -- October');
    text(-92,46.4,'OBS');
    % colormap(ax(1), ccmap)

    %
    hh802=subplot(3,1,2);cla;hold on;box on;
    hh802_pp=get(hh802,'position');
    set(hh802,'position',hh802_pp+[0.1 0.05 -0.2 0]);%set figure 
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_obs(11,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh802,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,'Monthly mean LST -- November');
    text(-92,46.4,'OBS');

    %
    hh803=subplot(3,1,3);cla;hold on;box on;
    hh803_pp=get(hh803,'position');
    set(hh803,'position',hh803_pp+[0.1 0.1 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_obs(12,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh803,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,'Monthly mean LST -- December');
    text(-92,46.4,'OBS');

    hh801.Colormap=ccmap;
    hh802.Colormap=ccmap;
    hh803.Colormap=ccmap;

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

%% plot monthly difference between observation and fvcom

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

%plot
color_scale1=[
73,0,0
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
7,0,73]/255;

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

figure(901);clf;hold on;set_portrait;
    hh901=subplot(3,1,1);cla;hold on;box on;
    hh901_pp=get(hh901,'position');
    set(hh901,'position',hh901_pp+[0.1 0 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_diff(1,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh901,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,'Monthly mean LST -- January');
    text(-92,46.4,'Difference');

    hh902=subplot(3,1,2);cla;hold on;box on;
    hh102_pp=get(hh902,'position');
    set(hh902,'position',hh102_pp+[0.1 0.05 -0.2 0]);%set figure 
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_diff(2,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh902,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,'Monthly mean LST -- Feburary');
    text(-92,46.4,'Difference');

    %
    hh903=subplot(3,1,3);cla;hold on;box on;
    hh903_pp=get(hh903,'position');
    set(hh903,'position',hh903_pp+[0.1 0.1 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_diff(3,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh903,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,'Monthly mean LST -- March');
    text(-92,46.4,'Difference');

    hh901.Colormap=ccmap;
    hh902.Colormap=ccmap;
    hh903.Colormap=ccmap;

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

figure(1001);clf;hold on;set_portrait;
    hh1001=subplot(3,1,1);cla;hold on;box on;
    hh1001_pp=get(hh1001,'position');
    set(hh1001,'position',hh1001_pp+[0.1 0 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_diff(4,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh1001,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,'Monthly mean LST -- April');
    text(-92,46.4,'Difference');
    % colormap(ax(1), ccmap)

    %
    hh1002=subplot(3,1,2);cla;hold on;box on;
    hh1002_pp=get(hh1002,'position');
    set(hh1002,'position',hh1002_pp+[0.1 0.05 -0.2 0]);%set figure 
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_diff(5,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh1002,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,'Monthly mean LST -- May');
    text(-92,46.4,'Difference');

    %
    hh1003=subplot(3,1,3);cla;hold on;box on;
    hh1003_pp=get(hh1003,'position');
    set(hh1003,'position',hh1003_pp+[0.1 0.1 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_diff(6,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh1003,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,'Monthly mean LST -- June');
    text(-92,46.4,'Difference');

    hh1001.Colormap=ccmap;
    hh1002.Colormap=ccmap;
    hh1003.Colormap=ccmap;

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

figure(1101);clf;hold on;set_portrait;
    hh1101=subplot(3,1,1);cla;hold on;box on;
    hh1101_pp=get(hh1101,'position');
    set(hh1101,'position',hh1101_pp+[0.1 0 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_diff(7,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh1101,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,'Monthly mean LST -- July');
    text(-92,46.4,'Difference');

    %
    hh1102=subplot(3,1,2);cla;hold on;box on;
    hh1102_pp=get(hh1102,'position');
    set(hh1102,'position',hh1102_pp+[0.1 0.05 -0.2 0]);%set figure 
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_diff(8,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh1102,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,'Monthly mean LST -- August');
    text(-92,46.4,'Difference');

    %
    hh1103=subplot(3,1,3);cla;hold on;box on;
    hh1103_pp=get(hh1103,'position');
    set(hh1103,'position',hh1103_pp+[0.1 0.1 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_diff(9,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh1103,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,'Monthly mean LST -- September');
    text(-92,46.4,'Difference');

    hh1101.Colormap=ccmap;
    hh1102.Colormap=ccmap;
    hh1103.Colormap=ccmap;

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

figure(1201);clf;hold on;set_portrait;
    hh1201=subplot(3,1,1);cla;hold on;box on;
    hh1201_pp=get(hh1201,'position');
    set(hh1201,'position',hh1201_pp+[0.1 0 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_diff(10,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh1201,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,'Monthly mean LST -- October');
    text(-92,46.4,'Difference');
    % colormap(ax(1), ccmap)

    %
    hh1202=subplot(3,1,2);cla;hold on;box on;
    hh1202_pp=get(hh1202,'position');
    set(hh1202,'position',hh1202_pp+[0.1 0.05 -0.2 0]);%set figure 
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_diff(11,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh1202,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,'Monthly mean LST -- November');
    text(-92,46.4,'Difference');

    %
    hh1203=subplot(3,1,3);cla;hold on;box on;
    hh1203_pp=get(hh1203,'position');
    set(hh1203,'position',hh1203_pp+[0.1 0.1 -0.2 0]);
    set(gca,'tickdir','out');
    scattercontourf(xlon,xlat,lst_diff(12,:)',[cmin:cinc:cmax]);
    caxis([cmin cmax]);
    colormap(hh1203,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,'Monthly mean LST -- December');
    text(-92,46.4,'Difference');

    hh1201.Colormap=ccmap;
    hh1202.Colormap=ccmap;
    hh1203.Colormap=ccmap;

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