clc
clear all
close all

%---------------------------------------------------
% EXPERIMENT DEFINITION
%---------------------------------------------------
exp_path = '/discover/nobackup/jkolassa/test_dir/';
exp_name = 'GEOSldas_CN40_local';
exp_resolution = 'SMAP_EASEv2_M09_GLOBAL';
exp_start_year = 1981;
exp_start_month = 1;
model_name = 'Catchment-CN'; % 'Catchment' or 'Catchment-CN'
%---------------------------------------------------
% END EXPERIMENT DEFINITION
%---------------------------------------------------


if strcmp(model_name,'Catchment-CN')
   model_str = 'catchcnclm40';
elseif strcmp(model_name,'Catchment')
   model_str = 'catch';
end

y_str = num2str(exp_start_year,'%04i');
m_str = num2str(exp_start_month,'%02i');

% define file names

rst_in_file = [exp_path '/' exp_name '/output/' exp_resolution '/rs/ens0000/Y' y_str '/M' m_str '/' exp_name '.' model_str '_internal_rst.' y_str m_str '01_0000']; 
cat_param_file = [exp_path '/' exp_name '/output/' exp_resolution '/rc_out/Y' y_str '/M' m_str '/' exp_name '.ldas_catparam.' y_str m_str '01_0000z.bin'];
tile_coord_file = [exp_path '/' exp_name '/output/' exp_resolution '/rc_out/' exp_name '.ldas_tilecoord.bin'];

% read coordinate data

tile_coord = read_tilecoord( tile_coord_file);
N_tile = tile_coord.N_tile;

% read parameters

[cat_param, cat_param_units] = read_catparam( cat_param_file, N_tile, 0 );

% create new rst file and remove parameters that need to be overwritten

save_file = [rst_in_file '_new'];
cmd_1 = ['ncks -x -v CDCR1,CDCR2,VGWMAX' ' ' rst_in_file ' ' save_file];
[status,cmdout] = system(cmd_1)

% write updated parameters to new file

nccreate(save_file,"CDCR1","Dimensions",{"tile",size(cat_param.cdcr1,1)}, ...
                   "FillValue","disable");
nccreate(save_file,"CDCR2","Dimensions",{"tile",size(cat_param.cdcr2,1)}, ...
                   "FillValue","disable");
nccreate(save_file,"VGWMAX","Dimensions",{"tile",size(cat_param.vgwmax,1)}, ...
                   "FillValue","disable");

ncwrite(save_file,'CDCR1',single(cat_param.cdcr1));
ncwriteatt(save_file,'CDCR1',"long_name","moisture_threshold");
ncwriteatt(save_file,'CDCR1',"units","kg m-2");

ncwrite(save_file,'CDCR2',single(cat_param.cdcr2));
ncwriteatt(save_file,'CDCR2',"long_name","max_water_content");
ncwriteatt(save_file,'CDCR2',"units","kg m-2");

ncwrite(save_file,'VGWMAX',single(cat_param.vgwmax));
ncwriteatt(save_file,'VGWMAX',"long_name","max_rootzone_water_content");
ncwriteatt(save_file,'VGWMAX',"units","kg m-2");

% move original rst file

cmd_2 = ['mv ' rst_in_file ' ' rst_in_file '.orig'];
[status,cmdout] = system(cmd_2)

% move new file to file name that model run looks for

cmd_3 = ['mv ' save_file ' ' rst_in_file];
[status,cmdout] = system(cmd_3)



