Skip to content
Snippets Groups Projects
Commit 3f42ea3b authored by lkugler's avatar lkugler
Browse files

corrections

parent 34d1bf9c
No related branches found
No related tags found
No related merge requests found
Show whitespace changes
Inline Side-by-side
docs/source/notebooks/tutorial1.ipynb
+ 6
− 1
View file @ 3f42ea3b
%% Cell type:markdown id:fd5c3005-f237-4495-9185-2d4d474cafd5 tags: %% Cell type:markdown id:fd5c3005-f237-4495-9185-2d4d474cafd5 tags:
# Tutorial 1: The assimilation step # Tutorial 1: The assimilation step
DART-WRF is a python package which automates many things like configuration, saving configuration and output, handling computing resources, etc. DART-WRF is a python package which automates many things like configuration, saving configuration and output, handling computing resources, etc.
The data for this experiment is accessible for students on the server srvx1. The data for this experiment is accessible for students on the server srvx1.
%% Cell type:markdown id:93d59d4d-c514-414e-81fa-4ff390290811 tags: %% Cell type:markdown id:93d59d4d-c514-414e-81fa-4ff390290811 tags:
### Configuring the experiment ### Configuring the experiment
Firstly, you need to configure the experiment in `config/cfg.py`. Firstly, you need to configure the experiment in `config/cfg.py`.
Let's go through the most important settings: Let's go through the most important settings:
- expname should be a unique identifier and will be used as folder name - expname should be a unique identifier and will be used as folder name
- model_dx is the model resolution in meters - model_dx is the model resolution in meters
- n_ens is the ensemble size - n_ens is the ensemble size
- update_vars are the WRF variables which shall be updated by the assimilation - update_vars are the WRF variables which shall be updated by the assimilation
- filter_kind is 1 for the EAKF (see the DART documentation for more) - filter_kind is 1 for the EAKF (see the DART documentation for more)
- prior and post_inflation defines what inflation we want (see the DART docs) - prior and post_inflation defines what inflation we want (see the DART docs)
- sec is the statistical sampling error correction from Anderson (2012) - sec is the statistical sampling error correction from Anderson (2012)
```python ```python
exp = utils.Experiment() exp = utils.Experiment()
exp.expname = "test_newcode" exp.expname = "test_newcode"
exp.model_dx = 2000 exp.model_dx = 2000
exp.n_ens = 10 exp.n_ens = 10
exp.update_vars = ['U', 'V', 'W', 'THM', 'PH', 'MU', 'QVAPOR', 'QCLOUD', 'QICE', 'PSFC'] exp.update_vars = ['U', 'V', 'W', 'THM', 'PH', 'MU', 'QVAPOR', 'QCLOUD', 'QICE', 'PSFC']
exp.filter_kind = 1 exp.filter_kind = 1
exp.prior_inflation = 0 exp.prior_inflation = 0
exp.post_inflation = 4 exp.post_inflation = 4
exp.sec = True exp.sec = True
``` ```
In case you want to generate new observations like for an observing system simulations experiment, OSSE), set In case you want to generate new observations like for an observing system simulations experiment, OSSE), set
```python ```python
exp.use_existing_obsseq = False`. exp.use_existing_obsseq = False`.
``` ```
`exp.nature` defines which WRF files will be used to draw observations from, e.g.: `exp.nature` defines which WRF files will be used to draw observations from, e.g.:
```python ```python
exp.nature = '/users/students/lehre/advDA_s2023/data/sample_nature/' exp.nature = '/users/students/lehre/advDA_s2023/data/sample_nature/'
``` ```
`exp.input_profile` is used, if you create initial conditions from a so called wrf_profile (see WRF guide). `exp.input_profile` is used, if you create initial conditions from a so called wrf_profile (see WRF guide).
```python ```python
exp.input_profile = '/doesnt_exist/initial_profiles/wrf/ens/raso.fc.<iens>.wrfprof' exp.input_profile = '/doesnt_exist/initial_profiles/wrf/ens/raso.fc.<iens>.wrfprof'
``` ```
For horizontal localization half-width of 20 km and 3 km vertically, set For horizontal localization half-width of 20 km and 3 km vertically, set
```python ```python
exp.cov_loc_vert_km_horiz_km = (3, 20) exp.cov_loc_vert_km_horiz_km = (3, 20)
``` ```
You can also set it to False for no vertical localization. You can also set it to False for no vertical localization.
#### Single observation #### Single observation
Set your desired observations like this. Set your desired observations like this.
```python ```python
t = dict(plotname='Temperature', plotunits='[K]', t = dict(plotname='Temperature', plotunits='[K]',
kind='RADIOSONDE_TEMPERATURE', kind='RADIOSONDE_TEMPERATURE',
n_obs=1, # number of observations n_obs=1, # number of observations
obs_locations=[(45., 0.)], # location of observations obs_locations=[(45., 0.)], # location of observations
error_generate=0.2, # observation error used to generate observations error_generate=0.2, # observation error used to generate observations
error_assimilate=0.2, # observation error used for assimilation error_assimilate=0.2, # observation error used for assimilation
heights=[1000,], # for radiosondes, use range(1000, 17001, 2000) heights=[1000,], # for radiosondes, use range(1000, 17001, 2000)
cov_loc_radius_km=50) # horizontal localization half-width cov_loc_radius_km=50) # horizontal localization half-width
exp.observations = [t,] # select observations for assimilation exp.observations = [t,] # select observations for assimilation
``` ```
#### Multiple observations #### Multiple observations
To generate a grid of observations, use To generate a grid of observations, use
```python ```python
vis = dict(plotname='VIS 0.6µm', plotunits='[1]', vis = dict(plotname='VIS 0.6µm', plotunits='[1]',
kind='MSG_4_SEVIRI_BDRF', sat_channel=1, kind='MSG_4_SEVIRI_BDRF', sat_channel=1,
n_obs=961, obs_locations='square_array_evenly_on_grid', n_obs=961, obs_locations='square_array_evenly_on_grid',
error_generate=0.03, error_assimilate=0.03, error_generate=0.03, error_assimilate=0.03,
cov_loc_radius_km=20) cov_loc_radius_km=20)
``` ```
But caution, n_obs should only be one of the following: But caution, n_obs should only be one of the following:
- 22500 for 2km observation density/resolution - 22500 for 2km observation density/resolution
- 5776 for 4km; - 5776 for 4km;
- 961 for 10km; - 961 for 10km;
- 256 for 20km; - 256 for 20km;
- 121 for 30km - 121 for 30km
For vertically resolved data, like radar, n_obs is the number of observations at each observation height level. For vertically resolved data, like radar, n_obs is the number of observations at each observation height level.
%% Cell type:markdown id:16bd3521-f98f-4c4f-8019-31029fd678ae tags: %% Cell type:markdown id:16bd3521-f98f-4c4f-8019-31029fd678ae tags:
### Configuring the hardware ### Configuring the hardware
In case you use a cluster which is not supported, configure paths inside `config/clusters.py`. In case you use a cluster which is not supported, configure paths inside `config/clusters.py`.
### Assimilate observations ### Assimilate observations
We start by importing some modules: We start by importing some modules:
```python ```python
import datetime as dt import datetime as dt
from dartwrf.workflows import WorkFlows from dartwrf.workflows import WorkFlows
``` ```
To assimilate observations at dt.datetime `time` we set the directory paths and times of the prior ensemble forecasts: To assimilate observations at dt.datetime `time` we set the directory paths and times of the prior ensemble forecasts:
```python ```python
prior_path_exp = '/users/students/lehre/advDA_s2023/data/sample_ensemble/' prior_path_exp = '/users/students/lehre/advDA_s2023/data/sample_ensemble/'
prior_init_time = dt.datetime(2008,7,30,12) prior_init_time = dt.datetime(2008,7,30,12)
prior_valid_time = dt.datetime(2008,7,30,12,30) prior_valid_time = dt.datetime(2008,7,30,12,30)
assim_time = prior_valid_time assim_time = prior_valid_time
``` ```
Finally, we run the data assimilation by calling To set up the experiment, call
```python ```python
w = WorkFlows(exp_config='cfg.py', server_config='srvx1.py') w = WorkFlows(exp_config='cfg.py', server_config='srvx1.py')
```
It will also create the output folders and backup the configuration files and scripts.
Finally, we run the data assimilation by calling
```python
w.assimilate(assim_time, prior_init_time, prior_valid_time, prior_path_exp) w.assimilate(assim_time, prior_init_time, prior_valid_time, prior_path_exp)
``` ```
Congratulations! You're done! Congratulations! You're done!
%% Cell type:code id:82e809a8-5972-47f3-ad78-6290afe4ae17 tags: %% Cell type:code id:82e809a8-5972-47f3-ad78-6290afe4ae17 tags:
``` python ``` python
``` ```
......
docs/source/notebooks/tutorial2.ipynb
+ 28
− 26
View file @ 3f42ea3b
%% Cell type:markdown id:fd5c3005-f237-4495-9185-2d4d474cafd5 tags: %% Cell type:markdown id:fd5c3005-f237-4495-9185-2d4d474cafd5 tags:
# Tutorial 2: Forecast after DA # Tutorial 2: Forecast after DA
**Goal**: To run a cycled data assimilation experiment. **Goal**: To run an ensemble of forecasts.
[`cycled_exp.py`](https://github.com/lkugler/DART-WRF/blob/master/generate_free.py) contains an example which will be explained here: [free_forecast.py ](https://github.com/lkugler/DART-WRF/blob/master/free_forecast.py) contains examples.
Now there are two options: Initialize the forecast with either (1) or (2). Run the forecast with (3)
1) To start a forecast from an existing forecast, i.e. from WRF restart files 1) Initialize a forecast from defined profiles of temperature, humidity and wind, i.e. from a `wrf_profile` (see WRF guide)
2) To start a forecast from defined thermodynamic profiles, i.e. from a `wrf_profile` 2) Initialize a forecast from an existing forecast, i.e. from WRF restart files - optionally with updates from data assimilation.
3) Run the forecast
### Restart a forecast ### 1) Initialize from sounding profiles
tutorial missing
### 2) Initialize a forecast from a previous forecast
To run a forecast from initial conditions of a previous forecasts, we import these modules To run a forecast from initial conditions of a previous forecasts, we import these modules
```python ```python
import datetime as dt import datetime as dt
from dartwrf.workflows import WorkFlows from dartwrf.workflows import WorkFlows
``` ```
Let's say you want to run a forecast starting at 9 UTC until 12 UTC. Let's say you want to run a forecast starting at 9 UTC until 12 UTC.
Initial conditions shall be taken from a previous experiment in `/user/test/data/sim_archive/exp_abc` which was initialized at 6 UTC and there are WRF restart files for 9 UTC. Initial conditions shall be taken from a previous experiment in `/user/test/data/sim_archive/exp_abc` which was initialized at 6 UTC and there are WRF restart files for 9 UTC.
Then the code would be Then the code would be
```python ```python
prior_path_exp = '/user/test/data/sim_archive/exp_abc' prior_path_exp = '/user/test/data/sim_archive/exp_abc'
prior_init_time = dt.datetime(2008,7,30,6) prior_init_time = dt.datetime(2008,7,30,6)
prior_valid_time = dt.datetime(2008,7,30,9) prior_valid_time = dt.datetime(2008,7,30,9)
w = WorkFlows(exp_config='cfg.py', server_config='srvx1.py') w = WorkFlows(exp_config='cfg.py', server_config='srvx1.py')
begin = dt.datetime(2008, 7, 30, 9) begin = dt.datetime(2008, 7, 30, 9)
end = dt.datetime(2008, 7, 30, 12) end = dt.datetime(2008, 7, 30, 12)
w.prepare_WRFrundir(begin) w.prepare_WRFrundir(begin)
w.prepare_IC_from_prior(prior_path_exp, prior_init_time, prior_valid_time) w.prepare_IC_from_prior(prior_path_exp, prior_init_time, prior_valid_time)
w.run_ENS(begin=begin, # start integration from here
end=end, # integrate until here
output_restart_interval=9999, # do not write WRF restart files
)
``` ```
Note that we use predefined workflow functions like `run_ENS`.
#### 2b) Optional: Update posterior with increments from assimilation
### Forecast run after Data Assimilation In order to continue a forecast after assimilation you need the posterior = prior (1) + increments (2)
In order to continue after assimilation you need the posterior = prior (1) + increments (2)
1. Set posterior = prior 1. Prepare initial conditions from a prior forecast:
```python ```python
id = w.prepare_IC_from_prior(prior_path_exp, prior_init_time, prior_valid_time, depends_on=id) id = w.prepare_IC_from_prior(prior_path_exp, prior_init_time, prior_valid_time, depends_on=id)
``` ```
2) Update posterior with increments from assimilation
After this, the wrfrst files are updated with assimilation increments from DART output and copied to the WRF's run directories so you can continue to run the forecast ensemble. 2. Update the initial conditions from data assimilation:
```python ```python
id = w.update_IC_from_DA(time, depends_on=id) id = w.update_IC_from_DA(time, depends_on=id)
``` ```
3) Define how long you want to run the forecast and when you want WRF restart files. Since they take a lot of space, we want as few as possible. After this, the wrfrst files are updated with assimilation increments from DART output and copied to the WRF's run directories so you can continue to run the forecast ensemble.
### 3) Run the forecast
Define how long you want to run the forecast and when you want WRF-restart files. Since they take a lot of space, we want as few as possible.
```python ```python
timedelta_integrate = dt.timedelta(hours=5) timedelta_integrate = dt.timedelta(hours=5)
output_restart_interval = 9999 # any value larger than the forecast duration
w.run_ENS(begin=begin, # start integration from here
end=time + timedelta_integrate, # integrate until here
output_restart_interval=9999, # do not write WRF restart files
)
``` ```
If you run DA in cycles of 15 minutes, it will be If you want to assimilate in 15 minutes again, use
```python ```python
timedelta_integrate = dt.timedelta(hours=5) timedelta_integrate = dt.timedelta(hours=5)
timedelta_btw_assim = dt.timedelta(minutes=15) timedelta_btw_assim = dt.timedelta(minutes=15)
output_restart_interval = timedelta_btw_assim.total_seconds()/60 output_restart_interval = timedelta_btw_assim.total_seconds()/60
```
3) Run WRF ensemble
```python
id = w.run_ENS(begin=time, # start integration from here id = w.run_ENS(begin=time, # start integration from here
end=time + timedelta_integrate, # integrate until here end=time + timedelta_integrate, # integrate until here
output_restart_interval=output_restart_interval, output_restart_interval=output_restart_interval
depends_on=id) )
``` ```
%% Cell type:code id:400244f1-098b-46ea-b29d-2226c7cbc827 tags: %% Cell type:code id:400244f1-098b-46ea-b29d-2226c7cbc827 tags:
``` python ``` python
``` ```
......
docs/source/notebooks/tutorial3.ipynb
+ 15
− 3
View file @ 3f42ea3b
%% Cell type:markdown id:fd5c3005-f237-4495-9185-2d4d474cafd5 tags: %% Cell type:markdown id:fd5c3005-f237-4495-9185-2d4d474cafd5 tags:
# Tutorial 3: Cycle forecast and assimilation # Tutorial 3: Cycle forecast and assimilation
**Goal**: To run a cycled data assimilation experiment. **Goal**: To run a cycled data assimilation experiment.
[`cycled_exp.py`](https://github.com/lkugler/DART-WRF/blob/master/generate_free.py) contains an example which will be explained here:
For example, your experiment can look like this
The script [cycled_exp.py](https://github.com/lkugler/DART-WRF/blob/master/cycled_exp.py) contains an example which will be explained here.
In this example, we assume that our computing jobs run on a cluster with a "job scheduler", so that the jobs don't run immediately, but when there are free resources. This means that we tell each job to wait on another job's completion by using the keyword `depends_on=id`.
```python ```python
w = WorkFlows(exp_config='cfg.py', server_config='jet.py')
timedelta_integrate = dt.timedelta(minutes=15)
timedelta_btw_assim = dt.timedelta(minutes=15)
prior_path_exp = '/jetfs/home/lkugler/data/sim_archive/exp_v1.19_P2_noDA' prior_path_exp = '/jetfs/home/lkugler/data/sim_archive/exp_v1.19_P2_noDA'
init_time = dt.datetime(2008, 7, 30, 13) init_time = dt.datetime(2008, 7, 30, 13)
time = dt.datetime(2008, 7, 30, 14) time = dt.datetime(2008, 7, 30, 14)
last_assim_time = dt.datetime(2008, 7, 30, 14) last_assim_time = dt.datetime(2008, 7, 30, 14)
forecast_until = dt.datetime(2008, 7, 30, 14, 15) forecast_until = dt.datetime(2008, 7, 30, 14, 15)
w.prepare_WRFrundir(init_time) w.prepare_WRFrundir(init_time)
id = w.run_ideal(depends_on=id) id = w.run_ideal()
prior_init_time = init_time prior_init_time = init_time
prior_valid_time = time prior_valid_time = time
while time <= last_assim_time: while time <= last_assim_time:
# usually we take the prior from the current time # usually we take the prior from the current time
# but one could use a prior from a different time from another run # but one could use a prior from a different time from another run
# i.e. 13z as a prior to assimilate 12z observations # i.e. 13z as a prior to assimilate 12z observations
prior_valid_time = time prior_valid_time = time
id = w.assimilate(time, prior_init_time, prior_valid_time, prior_path_exp, depends_on=id) id = w.assimilate(time, prior_init_time, prior_valid_time, prior_path_exp, depends_on=id)
# 1) Set posterior = prior # 1) Set posterior = prior
id = w.prepare_IC_from_prior(prior_path_exp, prior_init_time, prior_valid_time, depends_on=id) id = w.prepare_IC_from_prior(prior_path_exp, prior_init_time, prior_valid_time, depends_on=id)
# 2) Update posterior += updates from assimilation # 2) Update posterior += updates from assimilation
id = w.update_IC_from_DA(time, depends_on=id) id = w.update_IC_from_DA(time, depends_on=id)
# How long shall we integrate? # How long shall we integrate?
timedelta_integrate = timedelta_btw_assim timedelta_integrate = timedelta_btw_assim
output_restart_interval = timedelta_btw_assim.total_seconds()/60 output_restart_interval = timedelta_btw_assim.total_seconds()/60
if time == last_assim_time: #this_forecast_init.minute in [0,]: # longer forecast every full hour if time == last_assim_time: #this_forecast_init.minute in [0,]: # longer forecast every full hour
timedelta_integrate = forecast_until - last_assim_time # dt.timedelta(hours=4) timedelta_integrate = forecast_until - last_assim_time # dt.timedelta(hours=4)
output_restart_interval = 9999 # no restart file after last assim output_restart_interval = 9999 # no restart file after last assim
# 3) Run WRF ensemble # 3) Run WRF ensemble
id = w.run_ENS(begin=time, # start integration from here id = w.run_ENS(begin=time, # start integration from here
end=time + timedelta_integrate, # integrate until here end=time + timedelta_integrate, # integrate until here
output_restart_interval=output_restart_interval, output_restart_interval=output_restart_interval,
depends_on=id) depends_on=id)
# as we have WRF output, we can use own exp path as prior # as we have WRF output, we can use own exp path as prior
prior_path_exp = cluster.archivedir prior_path_exp = cluster.archivedir
id_sat = w.create_satimages(time, depends_on=id) id_sat = w.create_satimages(time, depends_on=id)
# increment time # increment time
time += timedelta_btw_assim time += timedelta_btw_assim
# update time variables # update time variables
prior_init_time = time - timedelta_btw_assim prior_init_time = time - timedelta_btw_assim
w.verify_sat(id_sat) w.verify_sat(id_sat)
w.verify_wrf(id) w.verify_wrf(id)
w.verify_fast(id) w.verify_fast(id)
``` ```
#### Job scheduling status #### Job scheduling status
If you work on a server with a queueing system, the script submits jobs into the SLURM queue with dependencies so that SLURM starts the jobs itself as soon as resources are available. Most jobs need only a few cores, but model integration is done across many nodes. You can look at the status with If you work on a server with a queueing system, the script submits jobs into the SLURM queue with dependencies so that SLURM starts the jobs itself as soon as resources are available. Most jobs need only a few cores, but model integration is done across many nodes. You can look at the status with
```bash ```bash
$ squeue -u `whoami` --sort=i $ squeue -u `whoami` --sort=i
JOBID PARTITION NAME USER ST TIME NODES NODELIST(REASON) JOBID PARTITION NAME USER ST TIME NODES NODELIST(REASON)
1710274 mem_0384 prepwrfr lkugler PD 0:00 1 (Priority) 1710274 mem_0384 prepwrfr lkugler PD 0:00 1 (Priority)
1710275 mem_0384 IC-prior lkugler PD 0:00 1 (Dependency) 1710275 mem_0384 IC-prior lkugler PD 0:00 1 (Dependency)
1710276 mem_0384 Assim-42 lkugler PD 0:00 1 (Dependency) 1710276 mem_0384 Assim-42 lkugler PD 0:00 1 (Dependency)
1710277 mem_0384 IC-prior lkugler PD 0:00 1 (Dependency) 1710277 mem_0384 IC-prior lkugler PD 0:00 1 (Dependency)
1710278 mem_0384 IC-updat lkugler PD 0:00 1 (Dependency) 1710278 mem_0384 IC-updat lkugler PD 0:00 1 (Dependency)
1710279 mem_0384 preWRF2- lkugler PD 0:00 1 (Dependency) 1710279 mem_0384 preWRF2- lkugler PD 0:00 1 (Dependency)
1710280_[1-10] mem_0384 runWRF2- lkugler PD 0:00 1 (Dependency) 1710280_[1-10] mem_0384 runWRF2- lkugler PD 0:00 1 (Dependency)
1710281 mem_0384 pRTTOV-6 lkugler PD 0:00 1 (Dependency) 1710281 mem_0384 pRTTOV-6 lkugler PD 0:00 1 (Dependency)
1710282 mem_0384 Assim-3a lkugler PD 0:00 1 (Dependency) 1710282 mem_0384 Assim-3a lkugler PD 0:00 1 (Dependency)
1710283 mem_0384 IC-prior lkugler PD 0:00 1 (Dependency) 1710283 mem_0384 IC-prior lkugler PD 0:00 1 (Dependency)
1710284 mem_0384 IC-updat lkugler PD 0:00 1 (Dependency) 1710284 mem_0384 IC-updat lkugler PD 0:00 1 (Dependency)
1710285 mem_0384 preWRF2- lkugler PD 0:00 1 (Dependency) 1710285 mem_0384 preWRF2- lkugler PD 0:00 1 (Dependency)
1710286_[1-10] mem_0384 runWRF2- lkugler PD 0:00 1 (Dependency) 1710286_[1-10] mem_0384 runWRF2- lkugler PD 0:00 1 (Dependency)
1710287 mem_0384 pRTTOV-7 lkugler PD 0:00 1 (Dependency) 1710287 mem_0384 pRTTOV-7 lkugler PD 0:00 1 (Dependency)
``` ```
%% Cell type:code id:400244f1-098b-46ea-b29d-2226c7cbc827 tags: %% Cell type:code id:400244f1-098b-46ea-b29d-2226c7cbc827 tags:
``` python ``` python
``` ```
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment