%---------plot figure3---------------
clc,clear,close all

load('./fig3.mat')

%% plot
month_name = {'Jan.','Feb.','Mar.','Apr.','May','Jun.',...
    'Jul.','Aug.','Sep.','Oct.','Nov.','Dec.'};

figure('color','w')
set(gcf,'position',[100,70,1400,700])

for m=1:12
    subplot(3,4,m)
    
    Ns_monthly_mean_2018 = Ns_monthly_mean(1,m,:);
    Ns_monthly_mean_2019 = Ns_monthly_mean(2,m,:);
    Ns_monthly_mean_2020 = Ns_monthly_mean(3,m,:);
    Ns_monthly_mean_2021 = Ns_monthly_mean(4,m,:);
    
    Ns_std_2018 = Ns_std(1,m,:);
    Ns_std_2019 = Ns_std(2,m,:);
    Ns_std_2020 = Ns_std(3,m,:);
    Ns_std_2021 = Ns_std(4,m,:);
    
    hold on
    
    plot(Ns_monthly_mean_2018(:),'ro-','MarkerSize',5);
    plot(Ns_monthly_mean_2019(:),'bo-','MarkerSize',5);
    plot(Ns_monthly_mean_2020(:),'o-','color',[0.4660 0.6740 0.1880],'MarkerSize',5);
    plot(Ns_monthly_mean_2021(:),'o-','color',[0.9290 0.6940 0.1250],'MarkerSize',5);
    
    errorbar(1:12,Ns_monthly_mean_2018(:),Ns_std_2018(:),'ro-','MarkerSize',5);
    errorbar(1:12,Ns_monthly_mean_2019(:),Ns_std_2019(:),'bo-','MarkerSize',5);
    errorbar(1:12,Ns_monthly_mean_2020(:),Ns_std_2020(:),'o-','color',[0.4660 0.6740 0.1880],'MarkerSize',5);
    errorbar(1:12,Ns_monthly_mean_2021(:),Ns_std_2021(:),'o-','color',[0.9290 0.6940 0.1250],'MarkerSize',5);
    
    hold off
    grid on
    set(gca,'Layer','bottom')
    set(gca,'xminortick','on');
    set(gca,'yminortick','on');
    
    xticks(1:2:13);
    xticklabels({'18','20','22',...
                 '00','02','04','06'})
    xlim([0,13])
    ylim([0,520])
    
    box on
    ax_tmp = gca;
    text(ax_tmp.XLim(1)+0.82*(ax_tmp.XLim(2)-ax_tmp.XLim(1)),...
        ax_tmp.YLim(1)+0.93*(ax_tmp.YLim(2)-ax_tmp.YLim(1)),...
        month_name{m},'FontSize',12)
    text(ax_tmp.XLim(1)+0.02*(ax_tmp.XLim(2)-ax_tmp.XLim(1)),...
        ax_tmp.YLim(1)+0.92*(ax_tmp.YLim(2)-ax_tmp.YLim(1)),...
        ['($',char('a'+m-1),'$)'],'interpreter','latex','fontweight','bold','FontSize',12)
    clear ax_tmp
    
    if m>=9
        xlabel('Local Time (hour)','fontweight','bold','FontSize',12)
    end
    if m==1 || m==5 || m==9
        ylabel('Occurence Number','interpreter','none','fontweight','bold','FontSize',12)
    end
end
clear m

legend('2018','2019','2020','2021',...
    'interpreter','latex',...
    'Orientation','horizontal','Box','off',...
    'Location',[0.40,0.00,0.20,0.03],...
    'fontweight','normal','FontSize',12)


