From e5850dfbfeca361cb6b19d10f36db81b349477b8 Mon Sep 17 00:00:00 2001
From: lkugler <lukas.kugler@gmail.com>
Date: Tue, 23 Nov 2021 21:08:37 +0100
Subject: [PATCH] new initialconditions routines

---
 scripts/prep_IC_prior.py        |  60 ++++++++++++++++
 scripts/prepare_namelist.py     |   9 +--
 scripts/prepare_wrf_initials.py | 118 --------------------------------
 scripts/update_IC.py            |  58 ++++++++++++++++
 4 files changed, 120 insertions(+), 125 deletions(-)
 create mode 100755 scripts/prep_IC_prior.py
 delete mode 100755 scripts/prepare_wrf_initials.py
 create mode 100755 scripts/update_IC.py

diff --git a/scripts/prep_IC_prior.py b/scripts/prep_IC_prior.py
new file mode 100755
index 0000000..a4d2a73
--- /dev/null
+++ b/scripts/prep_IC_prior.py
@@ -0,0 +1,60 @@
+import os, sys, warnings, glob
+import datetime as dt
+
+from config.cfg import exp, cluster
+from utils import copy, clean_wrfdir, try_remove
+
+"""
+Sets initial condition data (wrfinput/wrfrst file) in the run_WRF directory for each ensemble member 
+
+You have 2 options to restart:
+1) using wrfout files  (function create_wrfinput_from_wrfout)
+2) using wrfrst files  (function create_wrfrst_in_WRF_rundir)
+
+Ad 1: copy wrfout from prior to archivedir
+
+Ad 2: copies wrfrst to run_WRF directory
+
+"""
+
+def create_wrfrst_in_WRF_rundir(time, prior_init_time, prior_path_exp):
+    """copies wrfrst to run_WRF directory (for next WRF run)
+    """
+    for iens in range(1, exp.n_ens+1):
+        clean_wrfdir(cluster.wrf_rundir(iens))
+    
+        prior_wrfrst = prior_path_exp + prior_init_time.strftime('/%Y-%m-%d_%H:%M/') \
+                       +str(iens)+time.strftime('/wrfrst_d01_%Y-%m-%d_%H:%M:%S')
+        wrfrst = cluster.wrf_rundir(iens) + time.strftime('/wrfrst_d01_%Y-%m-%d_%H:%M:%S')
+        print('copy prior (wrfrst)', prior_wrfrst, 'to', wrfrst)
+        copy(prior_wrfrst, wrfrst)
+        
+        # remove all wrfrst (but not the one used)
+        files_rst = glob.glob(prior_path_exp + prior_init_time.strftime('/%Y-%m-%d_%H:%M/'+str(iens)+'/wrfrst_*'))
+        files_rst.remove(prior_wrfrst)
+        for f in files_rst:
+            print('removing', f)
+            try_remove(f)
+
+def create_updated_wrfinput_from_wrfout(time, prior_init_time, prior_path_exp):
+    """Same as create_wrfout_in_archivedir, but output is `wrfinput` in WRF run directory"""
+
+    print('writing updated wrfout to WRF run directory as wrfinput')
+    for iens in range(1, exp.n_ens+1):
+        prior_wrfout = prior_path_exp + prior_init_time.strftime('/%Y-%m-%d_%H:%M/') \
+                       +str(iens)+time.strftime('/wrfout_d01_%Y-%m-%d_%H:%M:%S')
+        post_wrfout = cluster.wrf_rundir(iens) + '/wrfinput_d01' 
+        copy(prior_wrfout, post_wrfout)
+        print(post_wrfout, 'created.')
+
+
+if __name__ == '__main__':
+    prior_path_exp = sys.argv[1]
+    prior_init_time = dt.datetime.strptime(sys.argv[2], '%Y-%m-%d_%H:%M')
+    prior_valid_time = dt.datetime.strptime(sys.argv[3], '%Y-%m-%d_%H:%M')
+
+    use_wrfout_as_wrfinput = False
+    if use_wrfout_as_wrfinput:
+        create_updated_wrfinput_from_wrfout(prior_valid_time, prior_init_time, prior_path_exp)
+    else:
+        create_wrfrst_in_WRF_rundir(prior_valid_time, prior_init_time, prior_path_exp)
diff --git a/scripts/prepare_namelist.py b/scripts/prepare_namelist.py
index 716c079..c71cd9b 100755
--- a/scripts/prepare_namelist.py
+++ b/scripts/prepare_namelist.py
@@ -42,15 +42,8 @@ def run(iens, begin, end, hist_interval=5, radt=5, archive=True,
         os.makedirs(archdir, exist_ok=True)
     else:
         archdir = './'
-    print('namelist for run from', begin, end, 'output to', archdir)
     sed_inplace(rundir+'/namelist.input', '<archivedir>', archdir)
 
-    if rst_inname:  # rst_inname did not work -> link file to directory
-        # sed_inplace(rundir+'/namelist.input', '<initdir>', rst_inname) 
-        fname = begin.strftime('/wrfrst_d01_%Y-%m-%d_%H:%M:%S')
-        symlink(rst_inname+fname, rundir+fname)
-        print('linked', rst_inname+fname, 'to rundir')
-
     # set times
     for k, v in {'<y1>': '%Y', '<m1>': '%m', '<d1>': '%d',
                  '<HH1>': '%H', '<MM1>': '%M'}.items():
@@ -59,6 +52,8 @@ def run(iens, begin, end, hist_interval=5, radt=5, archive=True,
                  '<HH2>': '%H', '<MM2>': '%M'}.items():
         sed_inplace(rundir+'/namelist.input', k, end.strftime(v))
 
+    print(rundir+'/namelist.input created')
+    print('runtime:', begin, end, 'output to', archdir)
     #########################
     if archive:
         
diff --git a/scripts/prepare_wrf_initials.py b/scripts/prepare_wrf_initials.py
deleted file mode 100755
index 3c00abc..0000000
--- a/scripts/prepare_wrf_initials.py
+++ /dev/null
@@ -1,118 +0,0 @@
-print('loading modules')
-import os, sys, warnings, glob
-import datetime as dt
-import netCDF4 as nc
-
-from config.cfg import exp, cluster
-from utils import symlink, copy, mkdir, clean_wrfdir, try_remove
-
-"""
-Sets initial condition data (wrfinput/wrrst file) in the run_WRF directory for each ensemble member 
-
-You have 2 options to restart:
-1) using wrfout files  (function create_updated_wrfinput_from_wrfout)
-2) using wrfrst files  (function create_wrfrst_in_WRF_rundir)
-
-Ad 1: copy wrfout from prior to archivedir
-      overwrite DA-updated variables with DART output
-
-Ad 2: copies wrfrst to run_WRF directory
-      overwrites DA-updated variables with DART output fields
-
-# assumes T = THM (dry potential temperature as prognostic variable)
-"""
-update_vars = ['Times',]
-update_vars.extend(exp.update_vars) # 'U', 'V', 'T', 'PH', 'MU', 'QVAPOR', 'QCLOUD', 'QICE', 'TSK', 'CLDFRA']
-updates = ','.join(update_vars)
-
-def create_wrfrst_in_WRF_rundir(time, prior_init_time, exppath_firstguess):
-    """copies wrfrst to run_WRF directory (for next WRF run)
-    """
-    for iens in range(1, exp.n_ens+1):
-        clean_wrfdir(cluster.wrf_rundir(iens))
-    
-        prior_wrfrst = exppath_firstguess + prior_init_time.strftime('/%Y-%m-%d_%H:%M/') \
-                       +str(iens)+time.strftime('/wrfrst_d01_%Y-%m-%d_%H:%M:%S')
-        wrfrst = cluster.wrf_rundir(iens) + time.strftime('/wrfrst_d01_%Y-%m-%d_%H:%M:%S')
-        print('copy prior (wrfrst)', prior_wrfrst, 'to', wrfrst)
-        copy(prior_wrfrst, wrfrst)
-    
-        filter_out = cluster.dartrundir+'/filter_restart_d01.'+str(iens).zfill(4)
-        print('update assimilated variables => overwrite', updates, 'in', wrfrst, 'from', filter_out) 
-        os.system(cluster.ncks+' -A -v '+updates+' '+filter_out+' '+wrfrst)
-    
-        print('writing T into THM of wrfrst')  # assumes T = THM (dry potential temperature as prognostic variable)
-        with nc.Dataset(filter_out, 'r') as ds_filter:
-            with nc.Dataset(wrfrst, 'r+') as ds_wrfrst:
-                ds_wrfrst.variables['THM_1'][:] = ds_filter.variables['T'][:]
-                ds_wrfrst.variables['THM_2'][:] = ds_filter.variables['T'][:]
-        print(wrfrst, 'created.')
-        
-        # remove all wrfrst (but not the one used)
-        files_rst = glob.glob(exppath_firstguess + prior_init_time.strftime('/%Y-%m-%d_%H:%M/'+str(iens)+'/wrfrst_*'))
-        files_rst.remove(prior_wrfrst)
-        for f in files_rst:
-            try_remove(f)
-
-def create_updated_wrfinput_from_wrfout(time, prior_init_time, exppath_firstguess):
-    """Same as create_wrfout_in_archivedir, but output is `wrfinput` in WRF run directory"""
-
-    print('writing updated wrfout to WRF run directory as wrfinput')
-    for iens in range(1, exp.n_ens+1):
-        prior_wrfout = exppath_firstguess + prior_init_time.strftime('/%Y-%m-%d_%H:%M/') \
-                       +str(iens)+time.strftime('/wrfout_d01_%Y-%m-%d_%H:%M:%S')
-        post_wrfout = cluster.wrf_rundir(iens) + '/wrfinput_d01' 
-        copy(prior_wrfout, post_wrfout)
-
-        if True: # overwrite DA updated variables
-            filter_out = cluster.dartrundir+'/filter_restart_d01.'+str(iens).zfill(4)
-            os.system(cluster.ncks+' -A -v '+updates+' '+filter_out+' '+post_wrfout)
-
-            # need to overwrite THM manually
-            with nc.Dataset(filter_out, 'r') as ds_filter:
-                with nc.Dataset(post_wrfout, 'r+') as ds_wrfout:
-                    ds_wrfout.variables['THM'][:] = ds_filter.variables['T'][:]
-        print(post_wrfout, 'created.')
-
-def create_wrfout_in_archivedir(time, prior_init_time, exppath_firstguess):
-    """Put updated wrfout in archive dir (because wrf restart writes no 0 minute wrfout)
-
-    probably not necessary, if using the wrf namelist option `write_hist_at_0h_rst`
-    """
-    print('writing updated wrfout to archive (for verification)')
-    for iens in range(1, exp.n_ens+1):
-        prior_wrfout = exppath_firstguess + prior_init_time.strftime('/%Y-%m-%d_%H:%M/') \
-                       +str(iens)+time.strftime('/wrfout_d01_%Y-%m-%d_%H:%M:%S')
-        post_wrfout_archive = cluster.archivedir +time.strftime('/%Y-%m-%d_%H:%M/') \
-                       +str(iens)+time.strftime('/wrfout_d01_%Y-%m-%d_%H:%M:%S')
-        # copy template wrfout (including cycled variables)
-        os.makedirs(os.path.dirname(post_wrfout_archive), exist_ok=True)
-        copy(prior_wrfout, post_wrfout_archive) 
-
-        # overwrite DA updated variables
-        filter_out = cluster.dartrundir+'/filter_restart_d01.'+str(iens).zfill(4)
-        #filter_out = cluster.archivedir +time.strftime('/%Y-%m-%d_%H:%M/assim_stage0/filter_restart_d01.'+str(iens).zfill(4))
-        os.system(cluster.ncks+' -A -v '+updates+' '+filter_out+' '+post_wrfout_archive)
-        
-        # need to overwrite THM manually
-        with nc.Dataset(filter_out, 'r') as ds_filter:
-            with nc.Dataset(post_wrfout_archive, 'r+') as ds_wrfout:
-                ds_wrfout.variables['THM'][:] = ds_filter.variables['T'][:]
-        print(post_wrfout_archive, 'created.')
-
-if __name__ == '__main__':
-    time = dt.datetime.strptime(sys.argv[1], '%Y-%m-%d_%H:%M')
-    prior_init_time = dt.datetime.strptime(sys.argv[2], '%Y-%m-%d_%H:%M')
-    exppath_firstguess = str(sys.argv[3])
-
-
-    use_wrfout_as_wrfinput = False
-    if use_wrfout_as_wrfinput:
-        create_updated_wrfinput_from_wrfout(time, prior_init_time, exppath_firstguess)
-    else:
-        create_wrfrst_in_WRF_rundir(time, prior_init_time, exppath_firstguess)
-
-        # this is done with the right WRF namelist entry (write first wrfout)
-        #create_wrfout_in_archivedir(time, prior_init_time, exppath_firstguess)
-
-
diff --git a/scripts/update_IC.py b/scripts/update_IC.py
new file mode 100755
index 0000000..8158cc9
--- /dev/null
+++ b/scripts/update_IC.py
@@ -0,0 +1,58 @@
+import os, sys, warnings
+import datetime as dt
+import netCDF4 as nc
+
+from config.cfg import exp, cluster
+"""
+Updates initial condition (wrfinput/wrfrst files) in the run_WRF directories with ./filter output
+
+# assumes T = THM (dry potential temperature as prognostic variable)
+"""
+
+use_wrfrst = True  # recommended to be True
+if use_wrfrst:
+    initials_fmt = '/wrfrst_d01_%Y-%m-%d_%H:%M:%S'
+else:
+    initials_fmt = '/wrfinput_d01' 
+
+
+def update_initials_in_WRF_rundir(time):
+    """updates wrfrst in run_WRF directory with posterior state from ./filter
+
+    Args:
+        time (dt.datetime):     time of assimilation (directory preceeding ./assim_stage0/...)
+    """
+
+    # which WRF variables will be updated?
+    update_vars = ['Times',]
+    update_vars.extend(exp.update_vars)
+    updates = ','.join(update_vars)
+
+    for iens in range(1, exp.n_ens+1):
+        ic_file = cluster.wrf_rundir(iens) + time.strftime(initials_fmt)
+        if not os.path.isfile(ic_file):
+            raise IOError(ic_file+' does not exist, updating impossible!')
+        else:
+            # overwrite DA updated variables
+            filter_out = cluster.archivedir+time.strftime('/%Y-%m-%d_%H:%M/assim_stage0/filter_restart_d01.'+str(iens).zfill(4))
+            print('update assimilated variables => overwrite', updates, 'from', filter_out) 
+            os.system(cluster.ncks+' -A -v '+updates+' '+filter_out+' '+ic_file)
+        
+            # assumes T = THM (dry potential temperature as prognostic variable)
+            print('writing T into THM')  
+            with nc.Dataset(filter_out, 'r') as ds_filter:
+
+                if use_wrfrst:
+                    with nc.Dataset(ic_file, 'r+') as ds_wrfrst:
+                        ds_wrfrst.variables['THM_1'][:] = ds_filter.variables['T'][:]
+                        ds_wrfrst.variables['THM_2'][:] = ds_filter.variables['T'][:]
+                else:
+                    with nc.Dataset(ic_file, 'r+') as ds_wrfout:
+                        ds_wrfout.variables['THM'][:] = ds_filter.variables['T'][:]
+
+            print(ic_file, 'updated.')
+
+
+if __name__ == '__main__':
+    time = dt.datetime.strptime(sys.argv[1], '%Y-%m-%d_%H:%M')
+    update_initials_in_WRF_rundir(time)
-- 
GitLab