From e1495b9611749006242962dae3d021d5f1056c1e Mon Sep 17 00:00:00 2001
From: Lukas Kugler <lukas.kugler@univie.ac.at>
Date: Thu, 27 Feb 2025 12:36:07 +0100
Subject: [PATCH] unfinished
---
config/jet.py | 4 +-
dartwrf/evaluate_obs_space.py | 4 +-
dartwrf/evaluate_prior.py | 4 +-
docs/source/tutorial1.rst | 4 +-
free_forecast2.py | 49 ++++++++++
multiple_exps.py | 173 ++++++++++++++++++++--------------
tests/test_assim.py | 2 +-
7 files changed, 162 insertions(+), 78 deletions(-)
create mode 100644 free_forecast2.py
diff --git a/config/jet.py b/config/jet.py
index ed06145..5b30a22 100755
--- a/config/jet.py
+++ b/config/jet.py
@@ -2,7 +2,7 @@
cluster_defaults = dict(
max_nproc = 20,
- max_nproc_for_each_ensemble_member = 16,
+ max_nproc_for_each_ensemble_member = 9,
use_slurm = True,
# binaries
@@ -20,7 +20,7 @@ cluster_defaults = dict(
# paths used as input
dir_wrf_src = '/jetfs/home/lkugler/data/compile/WRF-4.3/run',
- dir_dart_src = '/jetfs/home/lkugler/data/compile/DART/DART-10.8.3_10pct/models/wrf/work/',
+ dir_dart_src = '/jetfs/home/lkugler/data/compile/DART/DART-10.8.3/models/wrf/work/',
dir_rttov_src = '/jetfs/home/lkugler/data/compile/RTTOV13/rtcoef_rttov13/',
dir_dartwrf_dev = '/jetfs/home/lkugler/DART-WRF/',
diff --git a/dartwrf/evaluate_obs_space.py b/dartwrf/evaluate_obs_space.py
index 05d8b67..b772111 100755
--- a/dartwrf/evaluate_obs_space.py
+++ b/dartwrf/evaluate_obs_space.py
@@ -21,9 +21,9 @@ def evaluate_one_time(assim_time, valid_time, use_other_obsseq=False):
aso.prepare_prior_ensemble(valid_time, prior_init_time=assim_time, prior_valid_time=valid_time, prior_path_exp=cluster.archivedir)
dart_nml.write_namelist()
- if os.path.isfile(exp.use_existing_obsseq):
+ if os.path.isfile(exp.assimilate_existing_obsseq):
# use the existing obs_seq.out file
- shutil.copy(exp.use_existing_obsseq, cluster.dart_rundir+'/obs_seq.out')
+ shutil.copy(exp.assimilate_existing_obsseq, cluster.dart_rundir+'/obs_seq.out')
else:
# is there a pre-existing obs file from assimilation before?
f_oso = valid_time.strftime(cluster.pattern_obs_seq_out)
diff --git a/dartwrf/evaluate_prior.py b/dartwrf/evaluate_prior.py
index d218751..8f7ed49 100755
--- a/dartwrf/evaluate_prior.py
+++ b/dartwrf/evaluate_prior.py
@@ -21,8 +21,8 @@ if __name__ == "__main__":
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)
+ if exp.assimilate_existing_obsseq: # from another exp
+ oso_input = assim_time.strftime(exp.assimilate_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):
diff --git a/docs/source/tutorial1.rst b/docs/source/tutorial1.rst
index ed73fa4..25f5332 100644
--- a/docs/source/tutorial1.rst
+++ b/docs/source/tutorial1.rst
@@ -61,7 +61,7 @@ In case you want to generate new observations, like for an observing system simu
.. code-block:: python
- exp.use_existing_obsseq = False
+ exp.assimilate_existing_obsseq = False
in this case, you need to set the path to WRF nature run files from where DART can generate observations:
@@ -77,7 +77,7 @@ You can use pre-existing observation files with
.. code-block:: python
- exp.use_existing_obsseq = '/usr/data/sim_archive/exp_ABC/obs_seq_out/%Y-%m-%d_%H:%M_obs_seq.out'
+ exp.assimilate_existing_obsseq = '/usr/data/sim_archive/exp_ABC/obs_seq_out/%Y-%m-%d_%H:%M_obs_seq.out'
where time-placeholders (``%Y-%m-%d_%H:%M``) are filled in later, depending on the assimilation time.
diff --git a/free_forecast2.py b/free_forecast2.py
new file mode 100644
index 0000000..70182c1
--- /dev/null
+++ b/free_forecast2.py
@@ -0,0 +1,49 @@
+import sys
+import datetime as dt
+import pandas as pd
+from dartwrf.workflows import WorkFlows
+from dartwrf.utils import Config
+
+
+# import default config for jet
+from config.jet import cluster_defaults
+from config.defaults import dart_nml, CF_config
+
+
+ensemble_size = 1
+t0 = dt.datetime(2008, 7, 30, 8)
+
+for c_s in [0.15, 0.2]:
+
+ id = None
+ cfg = Config(name='exp_nat250m_Cs'+str(c_s),
+ model_dx = 2000,
+ ensemble_size = ensemble_size,
+ geo_em_forecast = '/jetfs/home/lkugler/data/sim_archive/geo_em.d01.nc.2km_200x200',
+ time = t0,
+ c_s = str(c_s),
+
+ input_profile = '/jetfs/home/lkugler/data/sim_archive/nat_250m_1600x1600x100/2008-07-30_08:00/1/input_sounding',
+ verify_against = 'nat_250m_blockavg2km',
+ **cluster_defaults, # type: ignore
+ )
+
+
+ w = WorkFlows(cfg)
+ w.prepare_WRFrundir(cfg)
+ id = w.run_ideal(cfg, depends_on=id)
+
+ timedelta_integrate = dt.timedelta(hours=10)
+ restart_interval = 9999
+
+ cfg.update( WRF_start=t0,
+ WRF_end=t0+timedelta_integrate,
+ restart=False,
+ restart_interval=restart_interval,
+ hist_interval_s=300,
+ )
+ # only verification depends on the output
+ id_long = w.run_WRF(cfg, depends_on=id)
+ id_long = w.run_RTTOV(cfg, depends_on=id_long)
+
+ w.verify(cfg, depends_on=id_long)
\ No newline at end of file
diff --git a/multiple_exps.py b/multiple_exps.py
index 44019cf..6aa2a94 100644
--- a/multiple_exps.py
+++ b/multiple_exps.py
@@ -1,3 +1,4 @@
+import sys
import datetime as dt
import pandas as pd
from dartwrf.workflows import WorkFlows
@@ -6,89 +7,103 @@ from dartwrf.utils import Config
# import default config for jet
from config.jet import cluster_defaults
-from config.defaults import dart_nml
+from config.defaults import dart_nml, CF_config, vis
# test multiple assimilation windows (11-12, 12-13, 13-14, )
timedelta_btw_assim = dt.timedelta(minutes=15)
-assim_times_start = pd.date_range('2008-07-30 11:00', '2008-07-30 13:00', freq='1h')
ensemble_size = 40
dart_nml['&filter_nml'].update(num_output_state_members=ensemble_size,
ens_size=ensemble_size)
-vis = dict(
- kind='MSG_4_SEVIRI_BDRF', sat_channel=1,
- km_between_obs=12, skip_border_km=8.0,
- error_generate=0.03, error_assimilate=0.06,
- loc_horiz_km=20,
- height=6000, loc_vert_km=6,
- )
+# which DART version to use?
+assimilate_cloudfractions = True
+cf1 = dict(kind='CF192km', loc_horiz_km=9999)
+cf2 = dict(kind='CF96km', loc_horiz_km=96)
+cf3 = dict(kind='CF48km', loc_horiz_km=48)
+cf4 = dict(kind='CF24km', loc_horiz_km=24)
+cf5 = dict(kind='CF12km', loc_horiz_km=12)
-for t0 in assim_times_start:
+assimilate_these_observations = [cf5]
- id = None
- cfg = Config(name=t0.strftime('exp_nat250_VIS_obs12_loc20_oe2_inf4-0.5_%H%M'),
- model_dx = 2000,
- ensemble_size = ensemble_size,
- dart_nml = dart_nml,
- geo_em_forecast = '/jetfs/home/lkugler/data/sim_archive/geo_em.d01.nc.2km_200x200',
- time = t0,
-
- use_existing_obsseq = False,
- nature_wrfout_pattern = '/jetfs/home/lkugler/data/sim_archive/nat_250m_1600x1600x100/*/1/wrfout_d01_%Y-%m-%d_%H_%M_%S',
- #geo_em_nature = '/jetfs/home/lkugler/data/sim_archive/geo_em.d01.nc.2km_200x200',
- geo_em_nature = '/jetfs/home/lkugler/data/sim_archive/geo_em.d01.nc.250m_1600x1600',
-
- update_vars = ['U', 'V', 'W', 'THM', 'PH', 'MU', 'QVAPOR', 'QCLOUD', 'QICE', 'QSNOW', 'PSFC'],
- input_profile = '/jetfs/home/lkugler/data/sim_archive/nat_250m_1600x1600x100/2008-07-30_08:00/1/input_sounding',
- nature_exp_verif = 'nat_250m_blockavg2km',
- assimilate_these_observations = [vis,],
- **cluster_defaults, # type: ignore
+if assimilate_cloudfractions:
+ cluster_defaults.update(
+ dir_dart_src = '/jetfs/home/lkugler/data/compile/DART/DART-v10.9.2-andrea/models/wrf/work/',
)
+else:
+ for obs in assimilate_these_observations:
+ if obs['kind'] == 'MSG_4_SEVIRI_BDRF':
+ cluster_defaults.update(
+ dir_dart_src = '/jetfs/home/lkugler/data/compile/DART/DART-10.8.3_10pct/models/wrf/work/',
+ )
+
- w = WorkFlows(cfg)
- w.prepare_WRFrundir(cfg)
- #id = w.run_ideal(cfg, depends_on=id)
+
+time0 = dt.datetime(2008, 7, 30, 11)
+
+id = None
+cfg = Config(name='exp_nat250_VIS_SO12_loc12_oe2_inf4-0.5',
+ model_dx = 2000,
+ ensemble_size = ensemble_size,
+ dart_nml = dart_nml,
+ geo_em_forecast = '/jetfs/home/lkugler/data/sim_archive/geo_em.d01.nc.2km_200x200',
+ time = time0,
+
+ assimilate_these_observations = assimilate_these_observations,
+ assimilate_cloudfractions = assimilate_cloudfractions,
+ CF_config=CF_config,
+
+ assimilate_existing_obsseq = False,
+ #nature_wrfout_pattern = '/jetfs/home/lkugler/data/sim_archive/nat_250m_1600x1600x100/*/1/wrfout_d01_%Y-%m-%d_%H_%M_%S',
+ #geo_em_nature = '/jetfs/home/lkugler/data/sim_archive/geo_em.d01.nc.2km_200x200',
+ # geo_em_nature = '/jetfs/home/lkugler/data/sim_archive/geo_em.d01.nc.250m_1600x1600',
- # assimilate at these times
- time0 = cfg.time
- assim_times = pd.date_range(time0, time0 + dt.timedelta(hours=1), freq=timedelta_btw_assim)
- last_assim_time = assim_times[-1]
+ update_vars = ['U', 'V', 'W', 'THM', 'PH', 'MU', 'QVAPOR', 'QCLOUD', 'QICE', 'QSNOW', 'PSFC'],
+ #input_profile = '/jetfs/home/lkugler/data/sim_archive/nat_250m_1600x1600x100/2008-07-30_08:00/1/input_sounding',
+ verify_against = 'nat_250m_blockavg2km',
+ **cluster_defaults)
- # loop over assimilations
- for i, t in enumerate(assim_times):
-
- if i == 0:
- if t == dt.datetime(2008, 7, 30, 11):
- prior_init_time = dt.datetime(2008, 7, 30, 8)
- else:
- prior_init_time = t - dt.timedelta(minutes=15)
-
- cfg.update(time = t,
- prior_init_time = prior_init_time,
- prior_valid_time = t,
- prior_path_exp = '/jetfs/home/lkugler/data/sim_archive/exp_nat250m_noDA/',)
- else:
- cfg.update(time = t,
- prior_init_time = assim_times[i-1],
- prior_valid_time = t,
- prior_path_exp = cfg.dir_archive,)
-
-
- id = w.assimilate(cfg, depends_on=id)
-
- # 1) Set posterior = prior
- id = w.prepare_IC_from_prior(cfg, depends_on=id)
-
- # 2) Update posterior += updates from assimilation
- id = w.update_IC_from_DA(cfg, depends_on=id)
-
- if t == last_assim_time:
- restart_interval = 9999
- timedelta_integrate = dt.timedelta(hours=4)
+
+w = WorkFlows(cfg)
+w.prepare_WRFrundir(cfg)
+#id = w.run_ideal(cfg, depends_on=id)
+
+# assimilate at these times
+assim_times = pd.date_range(time0, time0 + dt.timedelta(hours=4), freq=timedelta_btw_assim)
+last_assim_time = assim_times[-1]
+
+# loop over assimilations
+for i, t in enumerate(assim_times):
+
+ if i == 0:
+ if t == dt.datetime(2008, 7, 30, 11):
+ prior_init_time = dt.datetime(2008, 7, 30, 8)
else:
- restart_interval = timedelta_btw_assim.total_seconds()/60 # in minutes
- timedelta_integrate = dt.timedelta(minutes=15)
+ prior_init_time = t - dt.timedelta(minutes=15)
+
+ cfg.update(time = t,
+ prior_init_time = prior_init_time,
+ prior_valid_time = t,
+ prior_path_exp = '/jetfs/home/lkugler/data/sim_archive/exp_nat250m_noDA/',)
+ else:
+ cfg.update(time = t,
+ prior_init_time = assim_times[i-1],
+ prior_valid_time = t,
+ prior_path_exp = cfg.dir_archive,)
+
+
+ id = w.assimilate(cfg, depends_on=id)
+
+ # 1) Set posterior = prior
+ id = w.prepare_IC_from_prior(cfg, depends_on=id)
+
+ # 2) Update posterior += updates from assimilation
+ id = w.update_IC_from_DA(cfg, depends_on=id)
+
+ if t != last_assim_time:
+ # integrate until next assimilation
+ timedelta_integrate = dt.timedelta(minutes=15)
+ restart_interval = timedelta_btw_assim.total_seconds()/60 # in minutes
cfg.update( WRF_start=t,
WRF_end=t+timedelta_integrate,
@@ -98,4 +113,24 @@ for t0 in assim_times_start:
)
# 3) Run WRF ensemble
- id = w.run_WRF(cfg, depends_on=id)
\ No newline at end of file
+ id = w.run_WRF(cfg, depends_on=id)
+ id = w.run_RTTOV(cfg, depends_on=id)
+
+ if t.minute == 0 and i != 0:
+ # full hour but not first one
+ # make long forecasts without restart files
+ timedelta_integrate = dt.timedelta(hours=4)
+ restart_interval = 9999
+
+ cfg.update( WRF_start=t,
+ WRF_end=t+timedelta_integrate,
+ restart=True,
+ restart_interval=restart_interval,
+ hist_interval_s=300,
+ )
+ # only verification depends on the output
+ # but we only have one run directory, so we need to wait until the other forecast is done
+ id = w.run_WRF(cfg, depends_on=id)
+ id = w.run_RTTOV(cfg, depends_on=id)
+
+w.verify(cfg, depends_on=id)
\ No newline at end of file
diff --git a/tests/test_assim.py b/tests/test_assim.py
index 875d670..c68ed35 100644
--- a/tests/test_assim.py
+++ b/tests/test_assim.py
@@ -32,7 +32,7 @@ def test_dartwrf_cycled_da():
geo_em_forecast = '/jetfs/home/lkugler/data/sim_archive/geo_em.d01.nc.2km_200x200',
time = dt.datetime(2008, 7, 30, 11),
- use_existing_obsseq = False,
+ assimilate_existing_obsseq = False,
nature_wrfout_pattern = '/jetfs/home/lkugler/data/sim_archive/nat_250m_1600x1600x100/*/1/wrfout_d01_%Y-%m-%d_%H_%M_%S',
#geo_em_nature = '/jetfs/home/lkugler/data/sim_archive/geo_em.d01.nc.2km_200x200',
geo_em_nature = '/jetfs/home/lkugler/data/sim_archive/geo_em.d01.nc.250m_1600x1600',
--
GitLab