From 0d94818479457549d5e635bb082d5fc79950db04 Mon Sep 17 00:00:00 2001
From: lkugler <lukas.kugler@gmail.com>
Date: Tue, 18 Oct 2022 11:08:57 +0200
Subject: [PATCH] new

---
 dartwrf/evaluate_posterior.py | 27 +++++++++++++++++
 dartwrf/evaluate_prior.py     | 31 +++++++++++++++++++
 dartwrf/obsseq_2dim.py        | 56 +++++++++++++++++++++++++++++++++++
 3 files changed, 114 insertions(+)
 create mode 100755 dartwrf/evaluate_posterior.py
 create mode 100755 dartwrf/evaluate_prior.py
 create mode 100755 dartwrf/obsseq_2dim.py

diff --git a/dartwrf/evaluate_posterior.py b/dartwrf/evaluate_posterior.py
new file mode 100755
index 0000000..e1e24a4
--- /dev/null
+++ b/dartwrf/evaluate_posterior.py
@@ -0,0 +1,27 @@
+import os, sys, shutil, warnings
+import time as time_module
+import datetime as dt
+import numpy as np
+
+from config.cfg import exp, cluster
+from dartwrf import assim_synth_obs as aso
+
+
+
+if __name__ == "__main__":
+
+     assim_time = dt.datetime.strptime(sys.argv[1], "%Y-%m-%d_%H:%M")
+
+     print(" computing posterior observations ")
+     aso.write_list_of_inputfiles_posterior(assim_time)
+
+    # prepare an obsseq without rejected observations
+     if exp.use_existing_obsseq:  # from another exp
+          oso_input = assim_time.strftime(exp.use_existing_obsseq)
+     else:  # from same exp
+          oso_input = cluster.archivedir+'/obs_seq_out' + assim_time.strftime("/%Y-%m-%d_%H:%M_obs_seq.out-beforeQC")
+          if not os.path.isfile(oso_input):
+               oso_input = cluster.archivedir+'/obs_seq_out' + assim_time.strftime("/%Y-%m-%d_%H:%M_obs_seq.out")
+          shutil.copy(oso_input, cluster.dartrundir+'/obs_seq.out')
+
+     aso.evaluate(assim_time, output_format="%Y-%m-%d_%H:%M_obs_seq.final-eval_posterior_allobs")
\ No newline at end of file
diff --git a/dartwrf/evaluate_prior.py b/dartwrf/evaluate_prior.py
new file mode 100755
index 0000000..0e0ec66
--- /dev/null
+++ b/dartwrf/evaluate_prior.py
@@ -0,0 +1,31 @@
+import os, sys, shutil, warnings
+import time as time_module
+import datetime as dt
+import numpy as np
+
+from config.cfg import exp, cluster
+from dartwrf.utils import symlink, copy, sed_inplace, append_file, mkdir, try_remove, print, shell
+from dartwrf import assim_synth_obs as aso
+from osselyze.utils import get_prior_config
+
+if __name__ == "__main__":
+
+     assim_time = dt.datetime.strptime(sys.argv[1], "%Y-%m-%d_%H:%M")
+
+     print("> computing posterior observations ...")
+     print(">> prepare prior ensemble")
+     prior_exp, prior_init = get_prior_config(cluster.archive_base, exp.expname, init=assim_time)
+     prior_valid_time = assim_time
+     aso.prepare_prior_ensemble(assim_time, prior_init, prior_valid_time, cluster.archive_base+'/'+prior_exp)
+
+    # prepare an obsseq without rejected observations
+     if exp.use_existing_obsseq:  # from another exp
+          oso_input = assim_time.strftime(exp.use_existing_obsseq)
+     else:  # from same exp
+          oso_input = cluster.archivedir+'/obs_seq_out' + assim_time.strftime("/%Y-%m-%d_%H:%M_obs_seq.out-beforeQC")
+          if not os.path.isfile(oso_input):
+               oso_input = cluster.archivedir+'/obs_seq_out' + assim_time.strftime("/%Y-%m-%d_%H:%M_obs_seq.out")
+          shutil.copy(oso_input, cluster.dartrundir+'/obs_seq.out')
+
+
+     aso.evaluate(assim_time, output_format="%Y-%m-%d_%H:%M_obs_seq.final-eval_prior_allobs")
\ No newline at end of file
diff --git a/dartwrf/obsseq_2dim.py b/dartwrf/obsseq_2dim.py
new file mode 100755
index 0000000..c61d865
--- /dev/null
+++ b/dartwrf/obsseq_2dim.py
@@ -0,0 +1,56 @@
+"""Create obs_seq.out files with collapsed vertical dimension
+Specifically, one observation per column which is the maximum of the column
+"""
+
+from copy import copy
+import os, sys, shutil, warnings
+import time as time_module
+import datetime as dt
+import numpy as np
+
+from config.cfg import exp, cluster
+from dartwrf import assim_synth_obs as aso
+from dartwrf import obsseq
+
+
+if __name__ == "__main__":
+
+     assim_time = dt.datetime.strptime(sys.argv[1], "%Y-%m-%d_%H:%M")
+
+    # prepare an obsseq without rejected observations
+     if exp.use_existing_obsseq:  # from another exp
+          oso_input = assim_time.strftime(exp.use_existing_obsseq)
+
+     # only assured to work with single obstype
+     if len(exp.observations) > 1:
+          raise NotImplementedError()
+     n_obs = exp.observations[0]['n_obs']
+
+     # existing obsseq with multi levels
+     oso = obsseq.ObsSeq(oso_input)
+
+     nlev = len(oso.df)/n_obs
+     if nlev - int(nlev) != 0:
+          raise RuntimeError()
+     nlev = int(nlev)  # levels per obs
+     
+     # copy will be modified
+     output = copy(oso)
+     output.df = output.df.iloc[0::nlev]  #  every nth level = first level
+
+     #print(output.df, oso.df)
+
+     # iterate through, set value to max
+     for i_obs in range(0, n_obs):  # go through n_obs (all columns)
+
+          i_obs_subset = i_obs*nlev  # jumps by nlev (from one to next column)
+          column = oso.df.loc[0 + i_obs_subset:nlev + i_obs_subset, :]  # select column
+
+          output.df.loc[i_obs_subset, ('observations')] = float(column['observations'].max())
+          output.df.loc[i_obs_subset, ('truth')] = float(column['truth'].max())
+
+     print(output.df) #, 'observations'], output.df.loc[i_obs, 'observations'])
+
+     fout = cluster.archivedir + assim_time.strftime("/obs_seq_out/%Y-%m-%d_%H:%M_obs_seq.out")
+     os.makedirs(cluster.archivedir+'/obs_seq_out', exist_ok=True)
+     output.to_dart(fout)
-- 
GitLab