clc
clear all
close all

s2s_forcing_path = '/discover/nobackup/projects/gmao/geos_carb/fzeng/FORECASTS/GEOS5/CLIM/GEOS5v2/BCSD_Final/6-hourly/2007/jan01/ens1/';
s2s_forcing_file = [s2s_forcing_path 'GEOS5.200701.nc4'];

swgdn_data = ncread(s2s_forcing_file,'SWGDN');
lat_data = ncread(s2s_forcing_file,'latitude');
lon_data = ncread(s2s_forcing_file,'longitude');


exp_tile_coord_file = '/discover/nobackup/jkolassa/Eunjee_test/IDS_fcst_CTRL_M09_2007_ens1_SD/output/SMAP_EASEv2_M09_GLOBAL/rc_out/IDS_fcst_CTRL_M09_2007_ens1_SD.ldas_tilecoord.bin';
exp_tile_coord = read_tilecoord_le(exp_tile_coord_file);

exp1_swdown_mat = NaN*ones(exp_tile_coord.N_tile,48*4);
exp2_swdown_mat = NaN*ones(exp_tile_coord.N_tile,48*4);

for d = 1%:4

d_str = num2str(d,'%02i');

exp1_forcing_path = '/discover/nobackup/jkolassa/Eunjee_test/IDS_fcst_CTRL_M09_2007_ens1_SD/output/SMAP_EASEv2_M09_GLOBAL/cat/ens0000/Y2007/M01/';
exp1_forcing_file = [exp1_forcing_path 'IDS_fcst_CTRL_M09_2007_ens1_SD.tavg24_1d_lfs_Nt.200701' d_str '.nc4'];
exp1_swdown_data = ncread(exp1_forcing_file,'SWdown');
exp1_time_data = ncread(exp1_forcing_file,'time_stamp');
exp1_time_data = exp1_time_data';

exp_date_vec = NaN*ones(length(exp1_time_data),1);

for i = 1:size(exp1_time_data,1)
year_vec = str2double(exp1_time_data(i,1:4));
month_vec = str2double(exp1_time_data(i,5:6));
day_vec = str2double(exp1_time_data(i,7:8));
hour_vec = str2double(exp1_time_data(i,10:11));
min_vec = str2double(exp1_time_data(i,12:13));

exp_date_vec(i) = datenum(year_vec,month_vec,day_vec,hour_vec,min_vec,0);
end

exp2_forcing_path = '/discover/nobackup/jkolassa/Eunjee_test/IDS_fcst_CTRL_M09_2007_ens1_SD_2/output/SMAP_EASEv2_M09_GLOBAL/cat/ens0000/Y2007/M01/';
exp2_forcing_file = [exp2_forcing_path 'IDS_fcst_CTRL_M09_2007_ens1_SD_2.tavg24_1d_lfs_Nt.200701' d_str '.nc4'];
exp2_swdown_data = ncread(exp2_forcing_file,'SWdown'); 

force_date_hh = [3 9 15 21];
for i = 1:length(force_date_hh)
    force_date_vec(i) = datenum(year_vec,month_vec,day_vec,force_date_hh(i),0,0);
end

force_plot_vec = swgdn_data(127,251,1:4);
force_plot_vec = NaN*ones(size(exp1_swdown_data(1,:)));
one_index = exp_date_vec < (datenum(year_vec,month_vec,day_vec,6,0,0));
force_plot_vec(one_index) = swgdn_data(127,251,1);
two_index = (exp_date_vec < (datenum(year_vec,month_vec,day_vec,12,0,0))) & ...
(exp_date_vec >= (datenum(year_vec,month_vec,day_vec,6,0,0)));
force_plot_vec(two_index) = swgdn_data(127,251,2);
three_index = (exp_date_vec < (datenum(year_vec,month_vec,day_vec,18,0,0))) & ...
(exp_date_vec >= (datenum(year_vec,month_vec,day_vec,12,0,0)));
force_plot_vec(three_index) = swgdn_data(127,251,3);
four_index = (exp_date_vec < (datenum(year_vec,month_vec,day_vec,24,0,0))) & ...
(exp_date_vec >= (datenum(year_vec,month_vec,day_vec,18,0,0)));
force_plot_vec(four_index) = swgdn_data(127,251,4);

figure
plot(exp_date_vec,exp1_swdown_data(1,:),'-r','LineWidth',2)
hold on 
plot(exp_date_vec,exp2_swdown_data(1,:),'-b','LineWidth',2)
hold on
plot(exp_date_vec,force_plot_vec,'-k','LineWidth',2)
grid on
datetick('x',15)
xlabel('time [HH:MM]')
ylabel('SWdown [W/m^2]')
legend('experiment 1','experiment 2','S2S forcing')
saveas(gcf,'/home/jkolassa/Projects/Eunjee_test/SWdown_point','png')



end
