This web page is intended to serve as a guide through the practice exercises of this
tutorial. Exercises are split into seven main sections, each of which focuses on
a particular aspect of using the MPAS-JEDI data assimilation system.
In case you would like to refer to any of the lecture slides from previous days,
you can open the Tutorial Agenda in another window.
The test dataset can be downloaded from Here.
Note that Slurm-based job scripts used in this tutorial are included in mpasjedi_tutorial_slurmjob.tar.
You can proceed through the sections of this practical guide at your own pace. It is highly recommended
to go through the exercises in order, since later exercises may require the output of earlier ones.
Clicking the grey headers will expand each section or subsection.
The practical exercises in this tutorial have been tailored to work on
the Taiwania-3 system.
Taiwania-3 is an HPC cluster that provides most of the libraries needed by MPAS-JEDI and its pre- and post-processing
tools through modules. In general, before compiling and running MPAS-JEDI on your own system,
you will need to install spack-stack.
However, this tutorial does not cover the installation of spack-stack, which was pre-installed on Taiwania-3.
Logging onto Taiwania-3 from your laptop or desktop by
$ ssh -Y [account]@twnia3.nchc.org.tw ('-Y' is on Mac, should be '-X' for non-Mac laptop)
then using password and DUO two-factor authentication to get onto Taiwania-3. Note that
'-Y' is necessary in order to use X11-forwarding for direct graphics display from Taiwania-3.
It is recommened to at least login onto Taiwania-3 with two terminals for different tasks.
It is assumed that each user is under bash when doing these exercises.
First, copying the pre-prepared 'mpas_jedi_tutorial' directory to your own disk space:
$ cd /work/$USER
$ cp -r /work/gpsarc207/mpas_jedi_tutorial .
$ ls -l mpas_jedi_tutorial
total 18
drwxr-xr-x 3 gpsarc207 TRI1122359 4096 Sep 26 09:18 background
drwxr-xr-x 3 gpsarc207 TRI1122359 4096 Sep 26 09:18 background_120km
drwxr-xr-x 5 gpsarc207 TRI1122359 4096 Sep 26 09:18 B_Matrix
drwxr-xr-x 6 gpsarc207 TRI1122359 4096 Sep 25 15:09 bufr_lib
drwxr-xr-x 2 gpsarc207 TRI1122359 4096 Sep 26 09:18 crtm_coeffs_v2.3
drwxr-xr-x 3 gpsarc207 TRI1122359 4096 Sep 26 09:18 ensemble
-rw-r-xr-- 1 gpsarc207 TRI1122359 235 Oct 2 14:12 gnu-openmpi-taiwania3.sh
drwxr-xr-x 6 gpsarc207 TRI1122359 4096 Sep 27 14:34 graphics
drwxr-xr-x 2 gpsarc207 TRI1122359 4096 Sep 26 09:18 localization_pregenerated
drwxr-xr-x 4 gpsarc207 TRI1122359 4096 Oct 2 15:24 mpas_bundle_v2_DP
drwxr-xr-x 4 gpsarc207 TRI1122359 4096 Oct 2 14:23 mpas_bundle_v2_SP
drwxr-xr-x 3 gpsarc207 TRI1122359 4096 Oct 2 14:22 MPAS_JEDI_yamls_scripts
drwxr-xr-x 2 gpsarc207 TRI1122359 4096 Oct 2 08:21 MPAS_namelist_stream_physics_files
drwxr-xr-x 6 gpsarc207 TRI1122359 4096 Sep 25 15:09 obs2ioda_prebuild
drwxr-xr-x 3 gpsarc207 TRI1122359 4096 Sep 26 09:18 obs_bufr
drwxr-xr-x 3 gpsarc207 TRI1122359 4096 Sep 26 09:18 obs_ioda_pregenerated
drwxr-xr-x 4 gpsarc207 TRI1122359 4096 Sep 26 13:34 omboma_from2experiments
This copy will take some time as the size of the whole directory is ~33GB!
You may pay special attention to the following directories:
- mpas_bundle_v2_SP - Pre-build mpas-jedi with single-precision mode,
which will be used by most of DA exercises in this tutorial.
- obs2ioda_prebuild - Pre-build bufr2ioda converter as a backup plan in case you fail to compile it yourself.
- obs_ioda_pregenerated - Pre-converted ioda obs as a backup plan in case you fail the obs conversion exercise.
- localization_pregenerated - Pre-generated BUMP localization files
as a backup plan in case you fail the localization generation exercise.
Next, you may check what shell you are using by 'echo $SHELL'. If your are not under
bash, you should consider to change to bash by typing 'bash'. Note that
the exercises should work in principal with other shells
(if you know the difference of different shells), but we did not test them under other shell
environment.
Taiwania-3 uses the LMOD package to manage the software development.
Running module list to see what modules are loaded by default
right after you log in.
It should print something similar to below:
$ module list
No modules loaded
We will run all of the practical exercises in this tutorial in your work directories,
i.e., /work/$USER/mpas_jedi_tutorial.
Running jobs on Taiwania-3 requires the submission of a job script to a batch queueing system,
which will allocate requested computing resources to your job when they become available.
In general, it's best to avoid running any compute-intensive jobs on the login nodes, and
the practical instructions to follow will guide you in the process of submitting jobs when
necessary.
As a first introduction to running jobs on Taiwania-3, there are several key commands worth
noting:
- sbatch - This command submits a job script, which describes a job
to be run on one or more batch nodes.
- squeue -u $USER - This command tells you the status of your pending and running jobs.
Note that you may need to wait about 30 seconds for a recently submitted job to show up.
- scancel - This command deletes a queued or running job
Throughout this tutorial, you will submit your jobs to the Undetermined queue with an account
number Undetermined. At various points in the practical exercises, we'll need to submit jobs to Taiwania-3's queueing system
using the sbatch command, and after doing so, we may often check on the status
of the job with the squeue command, monitoring the log files produced by the job
once we see that the job has begun to run. All exercises will run with one single node (56 cores, 172GB available memory)
except for 4DEnVar that uses 3 nodes with the analysis at 3 times.
You're now ready to begin with the practical exercises of this tutorial!
In this section, the goal is to obtain, compile, and test the MPAS-JEDI code
through cmake/ctest mechanism.
In order to build MPAS-JEDI and its dependencies, it is recommended to access source code from mpas-bundle.
We will use the 'release/2.0.0' branch in the mpas-bundle repository:
$ cd /work/$USER/mpas_jedi_tutorial
$ mkdir mpas_bundle_v2
$ cd mpas_bundle_v2
$ git clone -b release/2.0.0 https://github.com/JCSDA/mpas-bundle.git code
The git command will clone the mpas-bundle repository from github to a local directory 'code'
, then make the 'release/2.0.0' branch as the active branch. The output to the terminal should look as following:
Cloning into 'code'...
remote: Enumerating objects: 461, done.
remote: Counting objects: 100% (82/82), done.
remote: Compressing objects: 100% (44/44), done.
remote: Total 461 (delta 44), reused 71 (delta 38), pack-reused 379
Receiving objects: 100% (461/461), 144.69 KiB | 1.64 MiB/s, done.
Resolving deltas: 100% (273/273), done.
Note that the mpas-bundle repository does not contain actual source code.
Instead, the CMakeLists.txt file under code includes the github repositories' branch/tag
information needed to build MPAS-JEDI.
We now create a 'build' folder and load the spack-stack environment pre-built
using the GNU compiler on Taiwania-3:
$ mkdir build
$ cd build
$ source ../../gnu-openmpi-taiwania3.sh
$ module list
The output to the module list in the terminal should look like the following (90 modules!!!):
Currently Loaded Modules:
1) stack-gcc/11.3.0 19) git-lfs/2.10.0 37) eccodes/2.27.0 55) py-pycodestyle/2.8.0 73) py-numpy/1.22.3
2) stack-openmpi/4.1.5 20) netcdf-fortran/4.6.0 38) py-attrs/22.2.0 56) krb5/1.20.1 74) py-six/1.16.0
3) git/1.8.3.1 21) gsibec/1.1.2 39) py-pycparser/2.21 57) libtirpc/1.2.6 75) py-python-dateutil/2.8.2
4) nccmp/1.9.0.1 22) gsl-lite/0.37.0 40) py-cffi/1.15.1 58) hdf/4.2.15 76) py-pytz/2022.2.1
5) parallel-netcdf/1.12.2 23) jedi-cmake/1.4.0 41) py-findlibs/0.0.2 59) libjpeg/2.1.0 77) py-pandas/1.4.0
6) parallelio/2.5.9 24) libpng/1.6.37 42) py-eccodes/1.4.2 60) py-pyhdf/0.10.4 78) py-setuptools/59.4.0
7) py-pip/23.0 25) libxmu/1.1.2 43) py-f90nml/1.4.3 61) libyaml/0.2.5 79) tar/1.26
8) wget/1.14 26) libxpm/3.5.12 44) py-h5py/3.7.0 62) py-pyyaml/6.0 80) gettext/0.21.1
9) base-env/1.0.0 27) libxt/1.1.5 45) curl/7.29.0 63) py-pybind11/2.8.1 81) libxcrypt/4.4.33
10) bufr/12.0.0 28) libxaw/1.0.13 46) pkg-config/0.27.1 64) py-beniget/0.4.1 82) sqlite/3.40.1
11) cmake/3.23.1 29) ncview/2.1.8 47) hdf5/1.14.0 65) py-gast/0.5.3 83) util-linux-uuid/2.38.1
12) ecbuild/3.7.2 30) netcdf-cxx4/4.3.1 48) numactl/2.0.14 66) py-ply/3.11 84) zlib/1.2.13
13) boost/1.78.0 31) json/3.10.5 49) pmix/4.2.3 67) py-pythran/0.12.2 85) python/3.10.8
14) fiat/1.1.0 32) json-schema-validator/2.1.0 50) openmpi/4.1.5 68) py-scipy/1.9.3 86) py-xarray/2022.3.0
15) ectrans/1.2.0 33) eigen/3.4.0 51) zstd/1.5.2 69) py-bottleneck/1.3.5 87) sp/2.3.3
16) ecmwf-atlas/0.33.0 34) eckit/1.23.1 52) netcdf-c/4.9.2 70) py-packaging/23.0 88) udunits/2.2.28
17) fckit/0.10.1 35) odc/1.4.6 53) py-cftime/1.0.3.4 71) py-numexpr/2.8.3 89) jedi-base-env/1.0.0
18) fftw/3.3.10 36) openjpeg/2.3.1 54) py-netcdf4/1.5.3 72) openblas/0.3.19 90) jedi-mpas-env/skylab-dev
Now we are ready to 'download' actual source code through cmake.
MPAS-JEDI uses cmake to automatically download the code from various github repositories,
listed in CMakeLists.txt under ~code. Now type the command below under the build directory:
$ git lfs install
$ cmake ../code
After it completes, you will see the actual source code of various
reporitories (e.g., oops, saber, ufo, ioda, crtm, mpas-jedi, and MPAS) is now under the code directory.
Meanwhile, Makefile files to build executables are generated under the build directory.
Now it is ready to compile MPAS-JEDI under 'build' using the standard make.
However, it is NOT recommended to compile the code on the login node at the same time.
$ make -j14
This requests 20 min walltime of a login node with 14 cores.
Also note that the reason we do the cmake step (will clone various github repositories) on the login
node instead of a compute node is that
the Taiwania-3's compute nodes have no access to the internet (so cannot clone github repos).
WARNING:The compilation will take ~20 min to complete.
Once we reach 100% of the compilation, we can check mpas-jedi related executables:
$ ls bin/mpas*
bin/mpas_atmosphere bin/mpasjedi_convertstate.x bin/mpasjedi_error_covariance_training.x bin/mpasjedi_rtpp.x bin/mpas_parse_init_atmosphere
bin/mpas_atmosphere_build_tables bin/mpasjedi_dirac.x bin/mpasjedi_forecast.x bin/mpasjedi_staticbinit.x bin/mpas_streams_gen
bin/mpas_data_checker.py bin/mpasjedi_eda.x bin/mpasjedi_gen_ens_pert_B.x bin/mpasjedi_variational.x
bin/mpas_data_downloader.py bin/mpasjedi_enkf.x bin/mpasjedi_hofx3d.x bin/mpas_namelist_gen
bin/mpas_init_atmosphere bin/mpasjedi_enshofx.x bin/mpasjedi_hofx.x bin/mpas_parse_atmosphere
The last step is to ensure that the code was compiled properly by running the MPAS-JEDI ctests,
with two lines of simple command:
$ cd mpas-jedi
$ ctest
At the moment the tests are running (take ~8 min to finish), it indicates if it passes or fails. At the end, a summary is provided with a percentage of the tests that passed, failed and the processing times.
The output to the terminal should look as following:
Start 46: test_mpasjedi_3dvar_2pe
46/47 Test #46: test_mpasjedi_3dvar_2pe ...................... Passed 22.43 sec
Start 47: test_mpasjedi_3dhybrid_bumpcov_bumploc_2pe
47/47 Test #47: test_mpasjedi_3dhybrid_bumpcov_bumploc_2pe ... Passed 21.44 sec
100% tests passed, 0 tests failed out of 47
Label Time Summary:
executable = 27.56 secxproc (13 tests)
mpasjedi = 523.69 secxproc (47 tests)
mpi = 503.92 secxproc (43 tests)
script = 496.13 seckproc (34 tests)
Total Test time (real) = 523.74 sec
To determine if a test passes or fails, a comparison of the test log and reference file is done internally taking as reference a tolerance value. This tolerance is specified via YAML file. These files can be found under mpas-jedi/test/testoutput in the build folder:
3denvar_2stream_bumploc.ref 3denvar_dual_resolution.run.ref 4denvar_bumploc.ref dirac_bumploc.ref hofx3d.ref parameters_bumpcov.run
3denvar_2stream_bumploc.run 3dfgat.ref 4denvar_ID.ref dirac_bumploc.run hofx3d_ropp.ref parameters_bumpcov.run.ref
3denvar_2stream_bumploc.run.ref 3dhybrid_bumpcov_bumploc.ref convertstate_bumpinterp.ref dirac_bumploc.run.ref hofx3d_rttovcpp.ref parameters_bumploc.ref
3denvar_amsua_allsky.ref 3dvar_bumpcov.ref convertstate_bumpinterp.run dirac_noloc.ref hofx3d.run parameters_bumploc.run
3denvar_amsua_allsky.run 3dvar_bumpcov_ropp.ref convertstate_bumpinterp.run.ref dirac_noloc.run hofx3d.run.ref parameters_bumploc.run.ref
3denvar_amsua_bc.ref 3dvar_bumpcov_rttovcpp.ref convertstate_unsinterp.ref dirac_noloc.run.ref hofx.ref rtpp.ref
3denvar_bumploc.ref 3dvar_bumpcov.run convertstate_unsinterp.run eda_3dhybrid.ref hofx.run
3denvar_bumploc.run 3dvar_bumpcov.run.ref convertstate_unsinterp.run.ref forecast.ref hofx.run.ref
3denvar_bumploc.run.ref 3dvar.ref dirac_bumpcov.ref forecast.run letkf_3dloc.ref
3denvar_dual_resolution.ref 3dvar.run dirac_bumpcov.run forecast.run.ref lgetkf.ref
3denvar_dual_resolution.run 3dvar.run.ref dirac_bumpcov.run.ref gen_ens_pert_B.ref parameters_bumpcov.ref
The goal of this session is to create the observation input files needed for running
subsequent MPAS-JEDI test cases.
In mpas-jedi_tutorial directory,
we can create a directory for converting NCEP BUFR observation into IODA-HDF5 format.
Then, we can make a clone of the obs2ioda repository and compile the converter.
$ cd /work/${USER}/mpas_jedi_tutorial
$ git clone https://github.com/NCAR/obs2ioda
$ cd obs2ioda/obs2ioda-v2/src/
NCEP BUFR library (https://github.com/NOAA-EMC/NCEPLIBS-bufr) is required to compile
obs2ioda-v2.x.
Edit obs2ioda-v2/src/Makefile to set
proper BUFR_LIB before make.
Here we have pre-build BUFR_LIB on Taiwania-3 and
we put the pre-build BUFR_LIB path in Makefile.
Then compile the converter by using make command:
$ make
If the compilation is successful, the executable file obs2ioda-v2.x will be generated.
The NCEP bufr input files are under the existing ~obs_bufr folder. Now we need to convert
the bufr-format data into the ioda-hdf5 format under the ~obs_ioda folder.
$ mkdir /work/${USER}/mpas_jedi_tutorial/obs_ioda
$ cd /work/${USER}/mpas_jedi_tutorial/obs_ioda
$ mkdir 2018041500
$ cd 2018041500
Link the prepbufr/bufr files to the working directory:
$ ln -fs ../../obs_bufr/2018041500/prepbufr.gdas.20180415.t00z.nr.48h .
$ ln -fs ../../obs_bufr/2018041500/gdas.1bamua.t00z.20180415.bufr .
$ ln -fs ../../obs_bufr/2018041500/gdas.gpsro.t00z.20180415.bufr .
$ ln -fs ../../obs_bufr/2018041500/gdas.satwnd.t00z.20180415.bufr .
Link executables and files:
$ ln -fs ../../obs2ioda/obs2ioda-v2/src/obs2ioda-v2.x .
$ ln -fs ../../mpas_bundle_v2_SP/build/bin/ioda-upgrade-v1-to-v2.x .
$ ln -fs ../../mpas_bundle_v2_SP/build/bin/ioda-upgrade-v2-to-v3.x .
$ ln -fs ../../mpas_bundle_v2_SP/code/ioda/share/ioda/yaml/validation/ObsSpace.yaml .
Copy obs_errtable to work directory:
$ cp ../../obs_bufr/2018041500/obs_errtable .
Now we can run obs2ioda-v2.x.
The usage is:
$ ./obs2ioda-v2.x [-i input_dir] [-o output_dir] [bufr_filename(s)_to_convert]
If input_dir and output_dir are not specified, the default is the current working directory.
If bufr_filename(s)_to_convert is not specified, the code looks for file name, **prepbufr.bufr** (also **satwnd.bufr**, **gnssro.bufr**, **amsua.bufr**, **airs.bufr**, **mhs.bufr**, **iasi.bufr**, **cris.bufr**) in the input/working directory. If the file exists, do the conversion, otherwise skip it.
$ ./obs2ioda-v2.x prepbufr.gdas.20180415.t00z.nr.48h
$ ./obs2ioda-v2.x gdas.satwnd.t00z.20180415.bufr
$ ./obs2ioda-v2.x gdas.gpsro.t00z.20180415.bufr
$ ./obs2ioda-v2.x gdas.1bamua.t00z.20180415.bufr
Create the directory iodav2
and move amsua and gnssro h5 files to the directory iodav2.
$ mkdir iodav2
$ mv amsua_*obs* iodav2
$ mv gnssro_obs_*.h5 iodav2
Due to NetCDF-Fortran interface does not allow reading/writing NF90_STRING, so station_id and variable_names are still written out as
char station_id(nlocs, nstring)
char variable_names(nvars, nstring)
rather than
string station_id(nlocs)
string variable_names(nvars)
So, for aircraft/satwind/satwnd/sfc/sondes, we need run upgrade executable ioda-upgrade-v1-to-v2.x.
The output files for aircraft/satwind/satwnd/sfc/sondes are created in iodav2 directory.
The input *h5 files in currect directory can be deleted.
$ ./ioda-upgrade-v1-to-v2.x satwind_obs_2018041500.h5 iodav2/satwind_obs_2018041500.h5
$ ./ioda-upgrade-v1-to-v2.x satwnd_obs_2018041500.h5 iodav2/satwnd_obs_2018041500.h5
$ ./ioda-upgrade-v1-to-v2.x sfc_obs_2018041500.h5 iodav2/sfc_obs_2018041500.h5
$ ./ioda-upgrade-v1-to-v2.x aircraft_obs_2018041500.h5 iodav2/aircraft_obs_2018041500.h5
$ ./ioda-upgrade-v1-to-v2.x sondes_obs_2018041500.h5 iodav2/sondes_obs_2018041500.h5
$ rm *.h5
Now we can convert iodav2 to iodav3 by using the executable file ioda-upgrade-v2-to-v3.x
$ ./ioda-upgrade-v2-to-v3.x iodav2/aircraft_obs_2018041500.h5 ./aircraft_obs_2018041500.h5 ObsSpace.yaml
$ ./ioda-upgrade-v2-to-v3.x iodav2/gnssro_obs_2018041500.h5 ./gnssro_obs_2018041500.h5 ObsSpace.yaml
$ ./ioda-upgrade-v2-to-v3.x iodav2/satwind_obs_2018041500.h5 ./satwind_obs_2018041500.h5 ObsSpace.yaml
$ ./ioda-upgrade-v2-to-v3.x iodav2/satwnd_obs_2018041500.h5 ./satwnd_obs_2018041500.h5 ObsSpace.yaml
$ ./ioda-upgrade-v2-to-v3.x iodav2/sfc_obs_2018041500.h5 ./sfc_obs_2018041500.h5 ObsSpace.yaml
$ ./ioda-upgrade-v2-to-v3.x iodav2/sondes_obs_2018041500.h5 ./sondes_obs_2018041500.h5 ObsSpace.yaml
$ ./ioda-upgrade-v2-to-v3.x iodav2/amsua_n15_obs_2018041500.h5 ./amsua_n15_obs_2018041500.h5 ObsSpace.yaml
$ ./ioda-upgrade-v2-to-v3.x iodav2/amsua_n18_obs_2018041500.h5 ./amsua_n18_obs_2018041500.h5 ObsSpace.yaml
$ ./ioda-upgrade-v2-to-v3.x iodav2/amsua_n19_obs_2018041500.h5 ./amsua_n19_obs_2018041500.h5 ObsSpace.yaml
$ ./ioda-upgrade-v2-to-v3.x iodav2/amsua_metop-a_obs_2018041500.h5 ./amsua_metop-a_obs_2018041500.h5 ObsSpace.yaml
$ ./ioda-upgrade-v2-to-v3.x iodav2/amsua_metop-b_obs_2018041500.h5 ./amsua_metop-b_obs_2018041500.h5 ObsSpace.yaml
Now you have the proper format ioda files under ~obs_ioda/2018041500 for subsequent mpas-jedi tests:
aircraft_obs_2018041500.h5 amsua_metop-b_obs_2018041500.h5 amsua_n18_obs_2018041500.h5
amsua_metop-a_obs_2018041500.h5 amsua_n15_obs_2018041500.h5 amsua_n19_obs_2018041500.h5
gnssro_obs_2018041500.h5 satwnd_obs_2018041500.h5 satwind_obs_2018041500.h5
sfc_obs_2018041500.h5 sondes_obs_2018041500.h5
Copy graphics directory to current directory:
$ cp -r ../../graphics .
$ cd graphics/standalone
$ python plot_obs_loc_tut.py
Now, we can check the figures:
In this session, we will be running an application called 𝐻(𝐱).
Creating hofx directory in mpas_jedi_tutorial,
$ cd /work/${USER}/mpas_jedi_tutorial
$ mkdir hofx
$ cd hofx
Now, we need prepare input files to run hofx.
The following commands are for linking files, such as,
physics-related data and tables, MPAS graph, streams, and namelist files:
$ ln -fs ../MPAS_namelist_stream_physics_files/*TBL .
$ ln -fs ../MPAS_namelist_stream_physics_files/*DBL .
$ ln -fs ../MPAS_namelist_stream_physics_files/*DATA .
$ ln -fs ../MPAS_namelist_stream_physics_files/x1.10242.graph.info.part.36 .
$ ln -fs ../MPAS_namelist_stream_physics_files/namelist.atmosphere_240km .
$ ln -fs ../MPAS_namelist_stream_physics_files/streams.atmosphere_240km .
$ ln -fs ../MPAS_namelist_stream_physics_files/stream_list.atmosphere.analysis .
$ ln -fs ../MPAS_namelist_stream_physics_files/stream_list.atmosphere.background .
$ ln -fs ../MPAS_namelist_stream_physics_files/stream_list.atmosphere.control .
$ ln -fs ../MPAS_namelist_stream_physics_files/stream_list.atmosphere.ensemble .
Link yamls:
$ ln -fs ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/namelists/geovars.yaml .
$ ln -fs ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/namelists/keptvars.yaml .
$ ln -fs ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/obsop_name_map.yaml .
Link MPAS 6-h forecast background file into template file and an MPAS init file as static.nc for 2-stream I/O:
$ ln -fs ../background/2018041418/x1.10242.init.2018-04-14_18.00.00.nc static.10242.nc
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc templateFields.10242.nc
Link the background and observaion input files:
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./bg.2018-04-15_00.00.00.nc
$ ln -fs ../obs_ioda/2018041500/aircraft_obs_2018041500.h5 .
$ ln -fs ../obs_ioda/2018041500/gnssro_obs_2018041500.h5 .
$ ln -fs ../obs_ioda/2018041500/satwind_obs_2018041500.h5 .
$ ln -fs ../obs_ioda/2018041500/sfc_obs_2018041500.h5 .
$ ln -fs ../obs_ioda/2018041500/sondes_obs_2018041500.h5 .
Copy hofx3d.yaml and run_hofx3d.sh to the working directory:
$ cp ../MPAS_JEDI_yamls_scripts/hofx3d.yaml .
$ cp ../MPAS_JEDI_yamls_scripts/run_hofx3d.sh .
Submit a slurm job to run hofx,
$ sbatch run_hofx3d.sh
According to the setting in hofx3d.yaml, the output files are obsout_hofx_*.h5.
You may check the output file content using 'ncdump -h obsout_hofx_aircraft.h5'.
you can also check hofx by run plot_diag_hofx_tut.py.
The following steps are for copy graphics directory, and run plot_diag_hofx_tut.py:
$ cp -r ../graphics .
$ cd graphics/standalone
$ python plot_diag_hofx_tut.py
In this practice, we will generate the BUMP localization files to be used for spatial localization
of ensemble background error covariance. These files will be used in the following 3D/4DEnVar and hybrid-3DEnVar
data assimilation practice.
After changing to our /work/${USER}/mpas_jedi_tutorial directory,
create the working directory.
$ cd /work/$USER/mpas_jedi_tutorial
$ mkdir localization
$ cd localization
Link MPAS physics-related data and tables, MPAS graph, MPAS streams, MPAS namelist files.
$ ln -fs ../MPAS_namelist_stream_physics_files/*BL ./
$ ln -fs ../MPAS_namelist_stream_physics_files/*DATA ./
$ ln -fs ../MPAS_namelist_stream_physics_files/x1.10242.graph.info.part.36 ./
$ ln -fs ../MPAS_namelist_stream_physics_files/streams.atmosphere_240km ./
$ ln -fs ../MPAS_namelist_stream_physics_files/stream_list.atmosphere.analysis ./
$ ln -fs ../MPAS_namelist_stream_physics_files/stream_list.atmosphere.background ./
$ ln -fs ../MPAS_namelist_stream_physics_files/stream_list.atmosphere.control ./
$ ln -fs ../MPAS_namelist_stream_physics_files/stream_list.atmosphere.ensemble ./
$ ln -fs ../MPAS_namelist_stream_physics_files/namelist.atmosphere_240km ./
Link MPAS 2-stream files.
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./templateFields.10242.nc
$ ln -fs ../background/2018041418/x1.10242.init.2018-04-14_18.00.00.nc ./static.10242.nc
Link mpasjedi executable.
$ ln -fs ../mpas_bundle_v2_SP/build/bin/mpasjedi_error_covariance_training.x ./
Link or copy yaml file.
$ ln -fs ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/namelists/geovars.yaml ./
$ ln -fs ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/namelists/keptvars.yaml ./
$ cp ../MPAS_JEDI_yamls_scripts/bumploc.yaml ./
Copy slurm job script.
$ cp ../MPAS_JEDI_yamls_scripts/run_bumploc.sh ./
Submit the slurm job.
$ sbatch run_bumploc.sh
User may check the slurm job status by running squeue -u $USER,
or monitoring the log file jedi.log.
After finishing the slurm job, user can check the BUMP localization files generated. The "local" files are generated for each processor (36 processors in this practice).
$ ls bumploc_1200km6km_nicas_local_*.nc
bumploc_1200km6km_nicas_local_000036-000001.nc bumploc_1200km6km_nicas_local_000036-000019.nc
bumploc_1200km6km_nicas_local_000036-000002.nc bumploc_1200km6km_nicas_local_000036-000020.nc
bumploc_1200km6km_nicas_local_000036-000003.nc bumploc_1200km6km_nicas_local_000036-000021.nc
bumploc_1200km6km_nicas_local_000036-000004.nc bumploc_1200km6km_nicas_local_000036-000022.nc
bumploc_1200km6km_nicas_local_000036-000005.nc bumploc_1200km6km_nicas_local_000036-000023.nc
bumploc_1200km6km_nicas_local_000036-000006.nc bumploc_1200km6km_nicas_local_000036-000024.nc
bumploc_1200km6km_nicas_local_000036-000007.nc bumploc_1200km6km_nicas_local_000036-000025.nc
bumploc_1200km6km_nicas_local_000036-000008.nc bumploc_1200km6km_nicas_local_000036-000026.nc
bumploc_1200km6km_nicas_local_000036-000009.nc bumploc_1200km6km_nicas_local_000036-000027.nc
bumploc_1200km6km_nicas_local_000036-000010.nc bumploc_1200km6km_nicas_local_000036-000028.nc
bumploc_1200km6km_nicas_local_000036-000011.nc bumploc_1200km6km_nicas_local_000036-000029.nc
bumploc_1200km6km_nicas_local_000036-000012.nc bumploc_1200km6km_nicas_local_000036-000030.nc
bumploc_1200km6km_nicas_local_000036-000013.nc bumploc_1200km6km_nicas_local_000036-000031.nc
bumploc_1200km6km_nicas_local_000036-000014.nc bumploc_1200km6km_nicas_local_000036-000032.nc
bumploc_1200km6km_nicas_local_000036-000015.nc bumploc_1200km6km_nicas_local_000036-000033.nc
bumploc_1200km6km_nicas_local_000036-000016.nc bumploc_1200km6km_nicas_local_000036-000034.nc
bumploc_1200km6km_nicas_local_000036-000017.nc bumploc_1200km6km_nicas_local_000036-000035.nc
bumploc_1200km6km_nicas_local_000036-000018.nc bumploc_1200km6km_nicas_local_000036-000036.nc
For a Dirac function multiplied by a given localization function, user can make a plot with mpas.dirac_nicas.2018-04-15_00.00.00.nc .
Copy the plotting script for Dirac output.
$ cp ../MPAS_JEDI_yamls_scripts/plot_bumploc.py ./
$ python plot_bumploc.py
$ display figure_bumploc.png
User can play with the specified localization lengths.
Please modify the following YAML keys from bumploc.yaml, then re-submit the slurm job script.
[Note that mpas.dirac_nicas.2018-04-15_00.00.00.nc and
figure_bumploc.png will be overwritten, so it is a good idea
to rename the existing files to compare them to those from a new experiment with modified localization lengths.]
...
io:
files prefix: bumploc_1200km6km
...
nicas:
resolution: 8.0
explicit length-scales: true
horizontal length-scale:
- groups: [common]
value: 1200.0e3
vertical length-scale:
- groups: [common]
value: 6.0e3
...
Make a plot for Dirac result (with plot_bumploc.py) to see how the localization looks like.
We are now ready to set up and run our first data assimilation update using 3DEnVar. For this 3DEnVar update we will use a 240 km background forecast and an ensemble at the same resolution.
Let's create a new working directory to run our 3DEnVar update:
$ cd /work/${USER}/mpas_jedi_tutorial
$ mkdir 3denvar
$ cd 3denvar
We can then link the relevant MPAS stream, graph, namelist, and physics files:
$ ln -sf ../MPAS_namelist_stream_physics_files/*BL ./
$ ln -sf ../MPAS_namelist_stream_physics_files/*DATA ./
$ ln -sf ../MPAS_namelist_stream_physics_files/x1.10242.graph* ./
$ ln -sf ../MPAS_namelist_stream_physics_files/namelist.atmosphere_240km ./
$ ln -sf ../MPAS_namelist_stream_physics_files/streams.atmosphere_240km ./
$ ln -sf ../MPAS_namelist_stream_physics_files/stream_list* ./
Next, we need to link in a few yaml files that define MPAS variables for JEDI:
$ ln -sf ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/namelists/geovars.yaml ./
$ ln -sf ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/namelists/keptvars.yaml ./
$ ln -sf ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/obsop_name_map.yaml ./
Next, we need to link our MPAS 2-stream files corresponding with our MPAS grid. We need both a static.nc file and templateFields.nc file.
For this experiment we need to link the following files:
$ ln -fs ../background/2018041418/x1.10242.init.2018-04-14_18.00.00.nc ./static.10242.nc
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./templateFields.10242.nc
Next, we should link our background file mpasout.2018-04-15_00.00.00.nc into the working directory and make a copy of the background file as the analysis file to be overwritten.
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./bg.2018-04-15_00.00.00.nc
$ cp ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./an.2018-04-15_00.00.00.nc
Next, we need to link our bump localization files.
$ ln -fs ../localization ./BUMP_files
Next, we need to link our observation ioda files into the working directory.
$ ln -fs ../obs_ioda/2018041500/*.h5 ./
Next, we are ready to copy the prepared 3denvar.yaml file for 3DEnVar. 3denvar.yaml contains all of the settings we need for this 3DEnVar update.
$ cp ../MPAS_JEDI_yamls_scripts/3denvar.yaml ./
We are finally ready to link a script to submit a slurm job to run our update and then submit our job. Notice our job submission script runs the following executable mpasjedi_variational.x, which is used for all of our variational updates.
$ ln -fs ../MPAS_JEDI_yamls_scripts/run_3denvar.sh ./
$ sbatch run_3denvar.sh
Once your job is submitted you can check the status of the job using squeue -u ${USER}.
When the job is complete, you can check the log file by using vi mpasjedi_3denvar.log.
If you look in your working directory you should also now see the following file:
an.2018-04-15_00.00.00.nc
Now that our update was successful, we may be interested in plotting the cost function and gradient norm reduction. We can accompilsh this with a prepared python script using the following:
$ cp ../MPAS_JEDI_yamls_scripts/plot_cost_grad_tut.py .
$ python plot_cost_grad_tut.py
We can then view the two figures using the following:
$ display cost.png
$ display grad.png
We can also plot the analysis increment with the following.
$ cp ../MPAS_JEDI_yamls_scripts/plot_analysis_increment.py ./
$ python plot_analysis_increment.py
$ display figure_increment_uReconstructZonal_10.png
Now they we know how to run a 3DEnVar update, let's prepare a second test which uses a higher resolution background of 120 km but a coarser ensemble of 240 km.
Let's create a new working directory to run our 120-240km 3DEnVar update:
$ cd /work/${USER}/mpas_jedi_tutorial
$ mkdir 120km240km_3denvar
$ cd 120km240km_3denvar
We can then link the relevant MPAS stream, graph, namelist, and physics files. In addition to our files for the last update we also need files for our 120 km mesh.
$ ln -sf ../MPAS_namelist_stream_physics_files/*BL ./
$ ln -sf ../MPAS_namelist_stream_physics_files/*DATA ./
$ ln -sf ../MPAS_namelist_stream_physics_files/x1.10242.graph* ./
$ ln -sf ../MPAS_namelist_stream_physics_files/x1.40962.graph* ./
$ ln -sf ../MPAS_namelist_stream_physics_files/namelist.atmosphere_240km ./
$ ln -sf ../MPAS_namelist_stream_physics_files/namelist.atmosphere_120km ./
$ ln -sf ../MPAS_namelist_stream_physics_files/streams.atmosphere_240km ./
$ ln -sf ../MPAS_namelist_stream_physics_files/streams.atmosphere_120km ./
$ ln -sf ../MPAS_namelist_stream_physics_files/stream_list* ./
Next, we need to link in a few yaml files again that define MPAS variables for JEDI:
$ ln -sf ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/namelists/geovars.yaml ./
$ ln -sf ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/namelists/keptvars.yaml ./
$ ln -sf ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/obsop_name_map.yaml ./
Next, we need to link our MPAS 2-stream files corresponding with both our 240 km and 120 km MPAS grids. We need both a static.nc file and templateFields.nc file.
For this experiment we need to link the following files:
$ ln -fs ../background/2018041418/x1.10242.init.2018-04-14_18.00.00.nc ./static.10242.nc
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./templateFields.10242.nc
$ ln -fs ../background_120km/2018041418/x1.40962.init.2018-04-14_18.00.00.nc ./static.40962.nc
$ ln -fs ../background_120km/2018041418/mpasout.2018-04-15_00.00.00.nc ./templateFields.40962.nc
Next, we should link our 120 km background file into the working directory and make a copy of the background file as the analysis file to be overwritten.
$ ln -fs ../background_120km/2018041418/mpasout.2018-04-15_00.00.00.nc ./bg.2018-04-15_00.00.00.nc
$ cp ../background_120km/2018041418/mpasout.2018-04-15_00.00.00.nc ./an.2018-04-15_00.00.00.nc
Next, we need to link our bump localization files.
$ ln -fs ../localization ./BUMP_files
Next, we need to link our observation ioda files into the working directory.
$ ln -fs ../obs_ioda/2018041500/*.h5 ./
Next, we are ready to copy the prepared 120km240km_3denvar.yaml file for 3DEnVar. 3denvar.yaml contains all of the settings we need for this 3DEnVar update.
$ cp ../MPAS_JEDI_yamls_scripts/120km240km_3denvar.yaml ./3denvar.yaml
We are finally ready to link a script to submit a slurm job to run our update and then submit our job. Notice our job submission script runs the following executable mpasjedi_variational.x, which is used for all of our variational updates.
$ ln -fs ../MPAS_JEDI_yamls_scripts/run_3denvar.sh ./
$ sbatch run_3denvar.sh
Once your job is submitted you can check the status of the job using squeue -u ${USER}.
When the job is complete, you can check the log file by using vi mpasjedi_3denvar.log.
If you look in your working directory you should also now see the following file:
an.2018-04-15_00.00.00.nc
We can plot the analysis increment as before with one small change to our script.
$ cp ../MPAS_JEDI_yamls_scripts/plot_analysis_increment.py ./
We need to change STATIC = "static.10242.nc" to STATIC = "static.40962.nc" within plot_analysis_increment.py in order to now correspond with our 120 km mesh.
$ python plot_analysis_increment.py
$ display figure_increment_uReconstructZonal_10.png
Now they we know how to run a 3DEnVar update, let's prepare an update with 4DEnVar which uses 3 subwindows.
Let's create a new working directory to run our 240km 4DEnVar update:
$ cd /work/${USER}/mpas_jedi_tutorial
$ mkdir 4denvar
$ cd 4denvar
We can then link the relevant MPAS stream, graph, namelist, and physics files just as we did for 3DEnVar.
$ ln -sf ../MPAS_namelist_stream_physics_files/*BL ./
$ ln -sf ../MPAS_namelist_stream_physics_files/*DATA ./
$ ln -sf ../MPAS_namelist_stream_physics_files/x1.10242.graph* ./
$ ln -sf ../MPAS_namelist_stream_physics_files/namelist.atmosphere_240km ./
$ ln -sf ../MPAS_namelist_stream_physics_files/streams.atmosphere_240km ./
$ ln -sf ../MPAS_namelist_stream_physics_files/stream_list* ./
Next, we need to link in a few yaml files again that define MPAS variables for JEDI:
$ ln -sf ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/namelists/geovars.yaml ./
$ ln -sf ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/namelists/keptvars.yaml ./
$ ln -sf ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/obsop_name_map.yaml ./
Next, we need to link our MPAS 2-stream files corresponding with our 240 km MPAS grid. We need both a static.nc file and templateFields.nc file.
For this experiment we need to link the following files:
$ ln -fs ../background/2018041418/x1.10242.init.2018-04-14_18.00.00.nc ./static.10242.nc
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./templateFields.10242.nc
Next, we need to link our 240 km background files for each subwindow into the working directory and make a copy of each background file as the analysis file for each subwindow to be overwritten.
$ ln -fs ../background/2018041418/mpasout.2018-04-14_21.00.00.nc ./bg.2018-04-14_21.00.00.nc
$ cp ../background/2018041418/mpasout.2018-04-14_21.00.00.nc ./an.2018-04-14_21.00.00.nc
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./bg.2018-04-15_00.00.00.nc
$ cp ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./an.2018-04-15_00.00.00.nc
$ ln -fs ../background/2018041418/mpasout.2018-04-15_03.00.00.nc ./bg.2018-04-15_03.00.00.nc
$ cp ../background/2018041418/mpasout.2018-04-15_03.00.00.nc ./an.2018-04-15_03.00.00.nc
Next, we need to link our bump localization files.
$ ln -fs ../localization ./BUMP_files
Next, we need to link our observation ioda files.
$ ln -fs ../obs_ioda/2018041500/*.h5 ./
Next, we are ready to copy the prepared 4denvar.yaml file for 4DEnVar. 4denvar.yaml contains all of the settings we need for this 4DEnVar update.
$ cp ../MPAS_JEDI_yamls_scripts/4denvar.yaml ./4denvar.yaml
We are finally ready to link a script to submit a slurm job to run our update and then submit our job. Notice our job submission script runs the following executable mpasjedi_variational.x, which is used for all of our variational updates.
$ ln -fs ../MPAS_JEDI_yamls_scripts/run_4denvar.sh ./
$ sbatch run_4denvar.sh
Once your job is submitted you can check the status of the job using squeue -u ${USER}.
When the job is complete, you can check the log file by using vi mpasjedi_4denvar.log.
If you look in your working directory you should also now see the following files corresponding to each subwindow:
an.2018-04-14_21.00.00.nc,
an.2018-04-15_00.00.00.nc, and
an.2018-04-15_03.00.00.nc
We can also plot the analysis increment at one of the subwindows with the following.
$ cp ../MPAS_JEDI_yamls_scripts/plot_analysis_increment.py ./
$ python plot_analysis_increment.py
$ display figure_increment_uReconstructZonal_10.png
Feel free to modify the an and bg files in plot_analysis_increment.py to see how the increment changes across the subwindows.
In this practice, we will run MPAS-JEDI 3DVar at global 240 km mesh.
The overall procedure is similar to the 3D/4DEnVar practice. However, 3DVar
uses the static multivariate background error covariance (B), which contains the
climatological characteristics. Here, the pre-generated static B will be used.
We hope to provide a practice for static B training in future MPAS-JEDI tutorials.
After changing to our /work/${USER}/mpas_jedi_tutorial directory,
create the working directory for 3dvar.
$ cd /work/$USER/mpas_jedi_tutorial
$ mkdir 3dvar
$ cd 3dvar
Link MPAS physics-related data and tables, MPAS graph, MPAS streams, and MPAS namelist files.
$ ln -fs ../MPAS_namelist_stream_physics_files/*BL ./
$ ln -fs ../MPAS_namelist_stream_physics_files/*DATA ./
$ ln -fs ../MPAS_namelist_stream_physics_files/x1.10242.graph.info.part.36 ./
$ ln -fs ../MPAS_namelist_stream_physics_files/streams.atmosphere_240km ./
$ ln -fs ../MPAS_namelist_stream_physics_files/stream_list.atmosphere.analysis ./
$ ln -fs ../MPAS_namelist_stream_physics_files/stream_list.atmosphere.background ./
$ ln -fs ../MPAS_namelist_stream_physics_files/stream_list.atmosphere.control ./
$ ln -fs ../MPAS_namelist_stream_physics_files/stream_list.atmosphere.ensemble ./
$ ln -fs ../MPAS_namelist_stream_physics_files/namelist.atmosphere_240km ./
Link the background file and copy it as an analysis file.
The file with prefix mpasout contains the variables for
"da_state" stream of 2-stream I/O.
And the data assimilation will overwrite
several selected analyzed variables (defined in stream_list.atmosphere.analysis)
in the analysis file.
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./bg.2018-04-15_00.00.00.nc
$ cp ./bg.2018-04-15_00.00.00.nc ./an.2018-04-15_00.00.00.nc
Link MPAS 2-stream files.
$ ln -fs ../background/2018041418/x1.10242.init.2018-04-14_18.00.00.nc ./static.10242.nc
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./templateFields.10242.nc
Prepare the pre-generated static B files.
$ mkdir -p BUMP_files/
$ cd BUMP_files
$ ln -fs ../../B_Matrix/bump_vertical_balance ./
$ ln -fs ../../B_Matrix/bump_nicas ./
$ ln -fs ../../B_Matrix/stddev ./
$ cd ..
Create a directory for observation input, and link the observation files to be assimilated.
Along with the conventional observations, here AMSU-A radiance data from NOAA-18 is added.
$ ln -fs ../obs_ioda/2018041500/aircraft_obs_2018041500.h5 ./
$ ln -fs ../obs_ioda/2018041500/gnssro_obs_2018041500.h5 ./
$ ln -fs ../obs_ioda/2018041500/satwind_obs_2018041500.h5 ./
$ ln -fs ../obs_ioda/2018041500/sfc_obs_2018041500.h5 ./
$ ln -fs ../obs_ioda/2018041500/sondes_obs_2018041500.h5 ./
$ ln -fs ../obs_ioda/2018041500/amsua_n18_obs_2018041500.h5 ./
In case user has a trouble with converting the obs format from bufr to ioda,
please link them from obs_ioda_pregenerated .
Prepare the coefficient files for CRTM.
$ mkdir crtm_coeffs_v2
$ cd crtm_coeffs_v2
$ ln -fs ../../crtm_coeffs_v2.3/* ./
$ cd ..
Link the mpasjedi executable.
$ ln -fs ../mpas_bundle_v2_SP/build/bin/mpasjedi_variational.x ./
Link or Copy yaml file.
$ ln -fs ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/namelists/geovars.yaml ./
$ ln -fs ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/namelists/keptvars.yaml ./
$ ln -fs ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/obsop_name_map.yaml ./
$ cp ../MPAS_JEDI_yamls_scripts/3dvar.yaml ./
Copy slurm job script.
$ cp ../MPAS_JEDI_yamls_scripts/run_3dvar.sh ./
Submit the slurm job.
$ sbatch run_3dvar.sh
User may check the slurm job status by running squeue -u $USER,
or monitoring the log file mpasjedi_3dvar.log .
After finishing the slurm job, user may explore the following steps to check the result.
- Check the analysis file and obs feedback files created.
$ ls -l an.2018-04-15_00.00.00.nc obsout_da_*.h5
-rw-r--r-- 1 bjung ncar 53929600 Sep 18 13:58 an.2018-04-15_00.00.00.nc
-rw-r--r-- 1 bjung ncar 140490883 Sep 18 13:58 obsout_da_aircraft.h5
-rw-r--r-- 1 bjung ncar 95854693 Sep 18 13:58 obsout_da_amsua_n18.h5
-rw-r--r-- 1 bjung ncar 14843761 Sep 18 13:58 obsout_da_gnssro.h5
-rw-r--r-- 1 bjung ncar 59712043 Sep 18 13:58 obsout_da_satwind.h5
-rw-r--r-- 1 bjung ncar 48512797 Sep 18 13:58 obsout_da_sfc.h5
-rw-r--r-- 1 bjung ncar 34931556 Sep 18 13:58 obsout_da_sondes.h5
- Check the nonlinear Jo term (fit-to-obs statistics) from the log file.
$ grep 'Nonlinear Jo' mpasjedi_3dvar.log
CostJo : Nonlinear Jo(Aircraft) = 338059, nobs = 396763, Jo/n = 0.852042, err = 1.77066
CostJo : Nonlinear Jo(GnssroRefNCEP) = 20941.7, nobs = 18845, Jo/n = 1.11126, err = 1.08039
CostJo : Nonlinear Jo(Satwind) = 17292.9, nobs = 26858, Jo/n = 0.643864, err = 2.39044
CostJo : Nonlinear Jo(SfcPCorrected) = 42584.8, nobs = 36491, Jo/n = 1.167, err = 52.8086
CostJo : Nonlinear Jo(Radiosonde) = 99066.7, nobs = 113434, Jo/n = 0.873342, err = 1.65807
CostJo : Nonlinear Jo(amsua_n18) = 98556.9, nobs = 57870, Jo/n = 1.70307, err = 0.287548
CostJo : Nonlinear Jo = 616502
CostJo : Nonlinear Jo(Aircraft) = 267155, nobs = 393888, Jo/n = 0.678251, err = 1.76997
CostJo : Nonlinear Jo(GnssroRefNCEP) = 13803, nobs = 17656, Jo/n = 0.781775, err = 1.11408
CostJo : Nonlinear Jo(Satwind) = 16614.4, nobs = 26777, Jo/n = 0.620471, err = 2.39035
CostJo : Nonlinear Jo(SfcPCorrected) = 33234.9, nobs = 35945, Jo/n = 0.924603, err = 52.8168
CostJo : Nonlinear Jo(Radiosonde) = 92407.5, nobs = 112614, Jo/n = 0.820569, err = 1.65863
CostJo : Nonlinear Jo(amsua_n18) = 37441.5, nobs = 57307, Jo/n = 0.653349, err = 0.287535
CostJo : Nonlinear Jo = 460656
CostJo : Nonlinear Jo(Aircraft) = 256533, nobs = 393389, Jo/n = 0.652112, err = 1.76987
CostJo : Nonlinear Jo(GnssroRefNCEP) = 10744.8, nobs = 17471, Jo/n = 0.615006, err = 1.11985
CostJo : Nonlinear Jo(Satwind) = 16070.3, nobs = 26755, Jo/n = 0.600648, err = 2.39045
CostJo : Nonlinear Jo(SfcPCorrected) = 31846.6, nobs = 35819, Jo/n = 0.889099, err = 52.8222
CostJo : Nonlinear Jo(Radiosonde) = 89111.9, nobs = 112390, Jo/n = 0.792881, err = 1.65882
CostJo : Nonlinear Jo(amsua_n18) = 31802.5, nobs = 57194, Jo/n = 0.556046, err = 0.287525
CostJo : Nonlinear Jo = 436110
Here, the first 7 lines represent the nonlinear Jo terms from the background field.
The following two sections represent the nonlinear Jo terms after 1st- and 2nd outer loops.
- Check the convergence using plot_cost_grad_tut.py .
$ cp ../MPAS_JEDI_yamls_scripts/plot_cost_grad_tut.py ./
We need to modify the following line of plot_cost_grad_tut.py as below.
dalogfile = 'mpasjedi_3dvar.log'
Then, execute the python script and check the figures.
$ python plot_cost_grad_tut.py
$ display cost.png
$ display grad.png
- Make a plot for analysis increment using plot_analysis_increment.py .
$ cp ../MPAS_JEDI_yamls_scripts/plot_analysis_increment.py ./
$ python plot_analysis_increment.py
$ display figure_increment_uReconstructZonal_10.png
User may want to plot the other fields by modifying the following lines of plot_analysis_increment.py.
VARIABLE = "uReconstructZonal" # variable to plot
LEVEL = 10 # vertical level
- Make a plot for Omb/Oma using plot_diag_omboma_tut.py .
$ cp ../MPAS_JEDI_yamls_scripts/plot_diag_omboma_tut.py ./
$ python ./plot_diag_omboma_tut.py
$ display figure_obsout_simple.png
This script will plot the horizontal map of OMB for aircraft temperature observation.
User can plot OMA or OMB for different observed variables, or for different observation feedback files
by modifying the python function
def readdata():
in plot_diag_omboma_tut.py .
User may want to try addtional practice as below.
- Change the number of inner and/or outer loop.
The default 3dvar.yaml uses 15 inner loops for 2 outer loops. User can modify the number of
inner loops, or even add the 3rd outer loop.
variational:
minimizer:
algorithm: DRPCG
iterations:
- <<: *iterationConfig
diagnostics:
departures: ombg
ninner: 15
- <<: *iterationConfig
ninner: 15
Is there any change in the fit-to-obs statistics, convergence, or analysis increment fields?
- Assimilate additional satellite radiance observation.
To do this, additional obs file (for example, ../obs_ioda/2018041500/amsua_n19_obs_2018041500.h5)
should be linked to /work/$USER/mpas_jedi_tutorial/3dvar/.
Also, user need to add the corresponding obs space yaml key in 3dvar.yaml .
Is there any change in the fit-to-obs statistics, convergence, or analysis increment fields?
In this practice, we will run MPAS-JEDI Hybrid-3DEnVar at global 240 km mesh.
The hybrid covariance will be made from both static B and ensemble B.
After changing to our /work/${USER}/mpas_jedi_tutorial directory,
create the working directory for Hybrid-3DEnVar.
$ cd /work/$USER/mpas_jedi_tutorial
$ mkdir hybrid-3denvar
$ cd hybrid-3denvar
Link MPAS physics-related data and tables, MPAS graph, MPAS streams, and MPAS namelist files.
$ ln -fs ../MPAS_namelist_stream_physics_files/*BL ./
$ ln -fs ../MPAS_namelist_stream_physics_files/*DATA ./
$ ln -fs ../MPAS_namelist_stream_physics_files/x1.10242.graph.info.part.36 ./
$ ln -fs ../MPAS_namelist_stream_physics_files/streams.atmosphere_240km ./
$ ln -fs ../MPAS_namelist_stream_physics_files/stream_list.atmosphere.analysis ./
$ ln -fs ../MPAS_namelist_stream_physics_files/stream_list.atmosphere.background ./
$ ln -fs ../MPAS_namelist_stream_physics_files/stream_list.atmosphere.control ./
$ ln -fs ../MPAS_namelist_stream_physics_files/stream_list.atmosphere.ensemble ./
$ ln -fs ../MPAS_namelist_stream_physics_files/namelist.atmosphere_240km ./
Link the background file and copy it as an analysis file.
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./bg.2018-04-15_00.00.00.nc
$ cp ./bg.2018-04-15_00.00.00.nc ./an.2018-04-15_00.00.00.nc
Link MPAS 2-stream files.
$ ln -fs ../background/2018041418/x1.10242.init.2018-04-14_18.00.00.nc ./static.10242.nc
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./templateFields.10242.nc
Prepare the pre-generated static B files.
$ mkdir -p BUMP_files/
$ cd BUMP_files
$ ln -fs ../../B_Matrix/bump_vertical_balance ./
$ ln -fs ../../B_Matrix/bump_nicas ./
$ ln -fs ../../B_Matrix/stddev ./
$ cd ..
Prepare the ensemble and localization files.
$ mkdir -p ensemble
$ cd ensemble
$ ln -fs ../../ensemble/2018041418 ./
$ cd ../BUMP_files
$ mkdir -p localization
$ cd localization
$ ln -fs ../../../localization/bumploc_1200km6km_nicas_local_*.nc ./
$ cd ../..
In case user has a trouble with generating the BUMP localization files,
please link them from localization_pregenerated.
Link the observation files to be assimilated.
Along we the conventional observation, here a single AMSU-A radiance observation from NOAA-18 is added.
$ ln -fs ../obs_ioda/2018041500/aircraft_obs_2018041500.h5 ./
$ ln -fs ../obs_ioda/2018041500/gnssro_obs_2018041500.h5 ./
$ ln -fs ../obs_ioda/2018041500/satwind_obs_2018041500.h5 ./
$ ln -fs ../obs_ioda/2018041500/sfc_obs_2018041500.h5 ./
$ ln -fs ../obs_ioda/2018041500/sondes_obs_2018041500.h5 ./
$ ln -fs ../obs_ioda/2018041500/amsua_n18_obs_2018041500.h5 ./
In case user has a trouble with converting the obs format from bufr to ioda,
please link them from obs_ioda_pregenerated.
Prepare the coefficient files for CRTM.
$ mkdir crtm_coeffs_v2
$ cd crtm_coeffs_v2
$ ln -fs ../../crtm_coeffs_v2.3/* ./
$ cd ..
Link the mpasjedi executable.
$ ln -fs ../mpas_bundle_v2_SP/build/bin/mpasjedi_variational.x ./
Link or Copy yaml file.
$ ln -fs ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/namelists/geovars.yaml ./
$ ln -fs ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/namelists/keptvars.yaml ./
$ ln -fs ../mpas_bundle_v2_SP/code/mpas-jedi/test/testinput/obsop_name_map.yaml ./
$ cp ../MPAS_JEDI_yamls_scripts/3dhyb.yaml ./
Copy slurm job script.
$ cp ../MPAS_JEDI_yamls_scripts/run_3dhyb.sh ./
Submit the slurm job.
$ sbatch run_3dhyb.sh
User may check the slurm job status by running squeue -u $USER,
or monitoring the log file mpasjedi_3dhyb.log .
After finishing the slurm job, user may explore the following steps to check the result.
- Check the analysis file and obs feedback files created.
$ ls -l an.2018-04-15_00.00.00.nc obsout_da_*.h5
-rw-r--r-- 1 bjung ncar 53929600 Sep 18 14:19 an.2018-04-15_00.00.00.nc
-rw-r--r-- 1 bjung ncar 140490883 Sep 18 14:19 obsout_da_aircraft.h5
-rw-r--r-- 1 bjung ncar 95854693 Sep 18 14:19 obsout_da_amsua_n18.h5
-rw-r--r-- 1 bjung ncar 14843761 Sep 18 14:19 obsout_da_gnssro.h5
-rw-r--r-- 1 bjung ncar 59712043 Sep 18 14:19 obsout_da_satwind.h5
-rw-r--r-- 1 bjung ncar 48512797 Sep 18 14:19 obsout_da_sfc.h5
-rw-r--r-- 1 bjung ncar 34931556 Sep 18 14:19 obsout_da_sondes.h5
- Check the nonlinear Jo term (fit-to-obs statistics) from the log file.
$ grep 'Nonlinear Jo' mpasjedi_3dhyb.log
CostJo : Nonlinear Jo(Aircraft) = 338059, nobs = 396763, Jo/n = 0.852042, err = 1.77066
CostJo : Nonlinear Jo(GnssroRefNCEP) = 20941.7, nobs = 18845, Jo/n = 1.11126, err = 1.08039
CostJo : Nonlinear Jo(Satwind) = 17292.9, nobs = 26858, Jo/n = 0.643864, err = 2.39044
CostJo : Nonlinear Jo(SfcPCorrected) = 42584.8, nobs = 36491, Jo/n = 1.167, err = 52.8086
CostJo : Nonlinear Jo(Radiosonde) = 99066.7, nobs = 113434, Jo/n = 0.873342, err = 1.65807
CostJo : Nonlinear Jo(amsua_n18) = 98556.9, nobs = 57870, Jo/n = 1.70307, err = 0.287548
CostJo : Nonlinear Jo = 616502
CostJo : Nonlinear Jo(Aircraft) = 238095, nobs = 394058, Jo/n = 0.604212, err = 1.76997
CostJo : Nonlinear Jo(GnssroRefNCEP) = 14260.2, nobs = 17575, Jo/n = 0.811391, err = 1.11652
CostJo : Nonlinear Jo(Satwind) = 15894.1, nobs = 26815, Jo/n = 0.592731, err = 2.39039
CostJo : Nonlinear Jo(SfcPCorrected) = 28765.4, nobs = 35917, Jo/n = 0.800887, err = 52.8231
CostJo : Nonlinear Jo(Radiosonde) = 90401.4, nobs = 112784, Jo/n = 0.801544, err = 1.65907
CostJo : Nonlinear Jo(amsua_n18) = 36962.3, nobs = 57403, Jo/n = 0.643908, err = 0.287586
CostJo : Nonlinear Jo = 424378
CostJo : Nonlinear Jo(Aircraft) = 223128, nobs = 393571, Jo/n = 0.566933, err = 1.76989
CostJo : Nonlinear Jo(GnssroRefNCEP) = 10531.6, nobs = 17343, Jo/n = 0.607254, err = 1.12392
CostJo : Nonlinear Jo(Satwind) = 14661.5, nobs = 26801, Jo/n = 0.547051, err = 2.39033
CostJo : Nonlinear Jo(SfcPCorrected) = 27206.5, nobs = 35824, Jo/n = 0.75945, err = 52.8263
CostJo : Nonlinear Jo(Radiosonde) = 85232.9, nobs = 112620, Jo/n = 0.756819, err = 1.65929
CostJo : Nonlinear Jo(amsua_n18) = 30078, nobs = 57291, Jo/n = 0.525004, err = 0.287582
CostJo : Nonlinear Jo = 390839
How different is it compared to that from 3DVar or 3D/4DEnVar?
- Check the convergence using plot_cost_grad_tut.py .
$ cp ../MPAS_JEDI_yamls_scripts/plot_cost_grad_tut.py ./
We need to modify the following line of plot_cost_grad_tut.py as below.
dalogfile = 'mpasjedi_3dhyb.log'
Then, execute the python script and check the figures.
$ python plot_cost_grad_tut.py
$ display cost.png
$ display grad.png
- Make a plot for analysis increment using plot_analysis_increment.py .
$ cp ../MPAS_JEDI_yamls_scripts/plot_analysis_increment.py ./
$ python plot_analysis_increment.py
$ display figure_increment_uReconstructZonal_10.png
User may want to plot the other fields by modifying the following lines of plot_analysis_increment.py.
VARIABLE = "uReconstructZonal" # variable to plot
LEVEL = 10 # vertical level
- Make a plot for Omb/Oma using plot_diag_omboma_tut.py .
$ cp ../MPAS_JEDI_yamls_scripts/plot_diag_omboma_tut.py ./
$ python ./plot_diag_omboma_tut.py
$ display figure_obsout_simple.png
User may want to try addtional practice as below.
- Change the hybrid weights.
The default 3dhyb.yaml uses 0.5 and 0.5 weights for static B and ensemble B. User can modify these weights, for example, 0.2/0.8 , 0.8/0.2 , 1.0/0.0 , 0.0/1.0 , or even 2.0/1.5 .
background error:
covariance model: hybrid
components:
- weight:
value: 0.5
covariance:
covariance model: SABER
...
- weight:
value: 0.5
covariance:
covariance model: ensemble
...
Is there any change in the fit-to-obs statistics, convergence, or analysis increment fields?
In this section, we will make plots for the OmA and OmB diagnostic files from two 5 day
cycling experiments by using the graphics package. In order to do this,
we first generate observation-space statistics following the instructions in Section 6.1.
Once we verify that the statistics files were successfully generated, we proceed to create the plots as shown below in Section 6.2.
The statistics are generated for each cycle and all ObsSpace (e.g., sondes, aircraft, gnssro, etc.) found in the folder containing the data assimilation feedback files. The statistics are binned globally and by latitude bands as specified in config.py.
The binned statistics are stored in HDF5 files located in a folder called stats that is created for each cycle.
The statistics are generated by the script DiagnoseObsStatistics.py in graphics.
Let's begin by checking the directory called omboma_from2experiments in your mpas_jedi_tutorial directory:
$ cd /work/$USER/mpas_jedi_tutorial
$ cd omboma_from2experiments
This folder contains the data assimilation outputs and scripts to generate the statistics:
- Exp1 and Exp2 — subfolders containing the data assimilation diagnostics outputs, e.g.,:
Exp1/CyclingDA/YYYYMMDDHH
- drive_stats.csh — script that drives the statistics generation for each experiment Exp1 and Exp 2 ($experiment) and cycle date between 00Z 15 April - 00Z 20 April 2018 ($vDate), by submitting the script generate_stats_omaomb.csh, e.g.,:
sbatch --export=ALL,da="$vDate",exp="$experiment" generate_stats_omaomb.csh
- generate_stats_omaomb.csh — sets the input and output directories,
makes a local copy of the graphics package and executes the script DiagnoseObsStatistics.py with corresponding arguments:
python DiagnoseObsStatistics.py -n 36 -p ${daDiagDir} -o obsout -app variational -nout 1
where:
-n: number of tasks/processors for multiprocessing
-p: path to deterministic or mean state UFO feedback files (e.g., obsout_da_*.h5)
-o: prefix for ObsSpace files
-app: application
-nout: number of outer loops
daDiagDir is for example: Exp1/CyclingDA/2018041500
- advanceCYMDH.py — defines a python function to help advance the time for the next cycle
Execute the driving script:
$ ./drive_stats.csh
We can check the status of our slurm jobs by squeue -u $USER, then verify that the statistics were created:
$ ls -l Exp*/CyclingDA/*/stats/*.h5
We can also check the log file log.DiagnoseObsStatistics.py to confirm that they were successfully created, e.g.,:
$ cd Exp1/CyclingDA/2018041500/stats/
$ tail -f log.DiagnoseObsStatistics.py
2023-09-15_14:11:16 [INFO] : DiagnoseObsStatistics.diagnoseObsSpace(da_gnssrobndropp1d) : Writing statistics file
2023-09-15_14:11:16 [INFO] : DiagnoseObsStatistics.diagnoseObsSpace(da_gnssrobndropp1d) : Finished
2023-09-15_14:11:16 [INFO] : __main__ : Finished __main__ successfully
If our statistics files were created and we can confirm that in the log file, we can proceed to the next section, were we will finally make the plots!
In order to make the plots, we first need to make a local copy of the graphics package and then make the corresponding modifications to the top-level script analyze_config.py. The plots are made in parallel by executing the script SpawnAnalyzeStats.py with all necessary arguments.
In this practice, the script analyze_config.py in the folder graphics has already been modified but here we provide a list of changes.
Let's start by changing the directory to graphics:
$ cd /work/$USER/mpas_jedi_tutorial/graphics
Let's check the changes made in analyze_config.py:
- Set the first and last cycle date time:
dbConf['firstCycleDTime'] = dt.datetime(2018,4,15,0,0,0)
dbConf['lastCycleDTime'] = dt.datetime(2018,4,20,0,0,0)
- Set the time increment between valid date times to 6 hours:
dbConf['cyTimeInc'] = dt.timedelta(hours=6)
- Set the verification type to omb/oma:
VerificationType = 'omb/oma'
- Set the verification space to obs:
VerificationSpace = 'obs'
- Set the workflow type to MPAS-JEDI_tutorial:
workflowType = 'MPAS-JEDI_tutorial'
- Set the path to the experiment's folder:
dbConf['expDirectory'] = os.getenv('EXP_DIR','/work/'+user+'/mpas_jedi_tutorial/omboma_from2experiments')
- Set the the name of the experiment that will be used as reference or control:
dbConf['cntrlExpName'] = 'exp1'
- Specify the experiments to be analyzed. Make sure that the value given to cntrlExpName is one of the keys in the dictionary experiments:
experiments['exp1'] = 'Exp1' + deterministicVerifyDir
experiments['exp2'] = 'Exp2' + deterministicVerifyDir
- In this practice the focus will only be the oma and omb results and analysis types:
BinValAxisProfileDiffCI: vertical profiles of RMS relative differences from the RMS of OmB of Exp1 with confidence intervals (95%), as a function of the impact height for different latitude bands and binned by latitude at 10km
CYAxisExpLines: time series of the bias, count, RMS, and STD averaged globally and by latitude bands
For this, append analysis types CYAxisExpLines and BinValAxisProfileDiffCI to the variable analysisTypes in lines 421, and 424.
Now we can proceed to execute the script SpawnAnalyzeStats.py, which submits individual jobs to Cheyenne computing nodes to process plots for
each diagSpace and analysis type (i.e., binning method):
$ ./SpawnAnalyzeStats.py -app variational -nout 1 -d gnssrobndropp1d
It can take a number of optional arguments (similar as in previous section), such as:
-app: application (e.g., variational, hofx)
-nout: number of outer loops in the variational method
-d: diagSpaces (e.g., amsua_,sonde,airc,sfc,gnssrobndropp1d,satwind)
For command-line examples, please see the script analyze_config.py.
In order to check the status of our slurm jobs, we can issue squeue -u $USER. Once the jobs are running, we can check the results in the folder gnssrobndropp1d_analyses (i.e., a diagSpace_analyses folder is created for each ObsSpace) for each analysis or plot types as specified above. Inside gnssrobndropp1d_analyses we will find subfolders for each of the analysis types containing several subfolders with diagnostics. Here we will take a look at the results in the folder omm, which contains oma and omb statistics in the same plot. Given that plots were created successfully, the end of the output log file should look like, e.g., for CYAxisExpLines :
$ cd gnssrobndropp1d_analyses/CYAxisExpLines
$ tail -f an.log
2023-09-18_11:14:49 [INFO] : analysis.AnalysisBase.gnssrobndropp1d.CYAxisExpLines : sigmaob, QCflag, good
2023-09-18_11:14:49 [INFO] : analysis.AnalysisBase.gnssrobndropp1d.CYAxisExpLines : sigmaob, lat, PolarLatBands
2023-09-18_11:15:03 [INFO] : analysis.Analyses.gnssrobndropp1d : Exiting Analyses.analyze()
2023-09-18_11:15:03 [INFO] : __main__ : Finished main() successfully
For visualization of the results, we can use the command display to show the plots, e.g.,:
$ cd omm
$ display QCflag_good_TSeries_0min_gnssrobndropp1d_omm_RMS.png
Now that we have learned how to make the plots, how about trying to plot statistics for the diagSpace amsua_n18?
