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.
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 Derecho system.
Derecho 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 Derecho.
MPAS-JEDI code build for this tutorial is based upon spack-stack-1.6.0.
Logging onto Derecho from your laptop or desktop by
$ ssh -X your-username@derecho.hpc.ucar.edu
then using your password and DUO two-factor authentication to get onto Derecho. Note that
'-X' is necessary in order to use X11-forwarding for direct graphics display from Derecho.
It is recommended to at least login onto Derecho with two terminals for different tasks.
First, copying the pre-prepared tutorial test dataset directory to your own scratch disk space:
$ cd /glade/derecho/scratch/$USER
$ cp -r /glade/derecho/scratch/liuz/mpasjedi_tutorial202506NCAS ./mpas_jedi_tutorial
$ ls -l mpas_jedi_tutorial
total 208
drwxr-xr-x 3 liuz ncar 16384 Jun 18 15:08 background
drwxr-xr-x 3 liuz ncar 16384 Jun 18 15:07 background_120km
drwxr-xr-x 3 liuz ncar 16384 Jun 18 22:44 Bflow_global120km
drwxr-xr-x 2 liuz ncar 16384 Jun 18 21:55 Bflow_preprocessing
drwxr-xr-x 5 liuz ncar 16384 Jun 18 15:08 B_Matrix
drwxr-xr-x 4 liuz ncar 16384 Jun 18 21:04 conus15km
drwxr-xr-x 2 liuz ncar 16384 Jun 18 15:07 crtm_coeffs_v3
drwxr-xr-x 5 liuz ncar 16384 Jun 19 11:12 cyclingDA
drwxr-xr-x 3 liuz ncar 16384 Jun 18 15:07 ensemble
drwxr-xr-x 2 liuz ncar 16384 Jun 18 21:14 MPAS_JEDI_yamls_scripts
drwxr-xr-x 2 liuz ncar 16384 Jun 18 15:08 MPAS_namelist_stream_physics_files
drwxr-xr-x 3 liuz ncar 16384 Jun 18 15:07 obs_bufr
drwxr-xr-x 4 liuz ncar 16384 Jun 18 15:07 omboma_from2experiments
This copy will take some time as the size of the whole directory is ~15GB!
On derecho, the default shell is bash.
Derecho 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
Currently Loaded Modules:
1) ncarenv/24.12 (S) 3) intel/2024.2.1 5) libfabric/1.15.2.0 7) hdf5/1.12.3
2) craype/2.7.31 4) ncarcompilers/1.0.0 6) cray-mpich/8.1.29 8) netcdf/4.9.2
Where:
S: Module is Sticky, requires --force to unload or purge
Post-processing and graphics exercises will need Python and NCAR Command
Language (NCL).
On Derecho, these are available in
a Conda environment named 'npl'. The npl environment also provides
NCL. We can load the conda module and activate the npl environment with
the following commands:
$ module reset # this may or may not be needed
$ module load conda
$ conda activate npl
$ which ncl # check whether ncl exists
You will run all of the practical exercises in your own scratch space under
/glade/derecho/scratch/$USER/mpas_jedi_tutorial.
Running jobs on Derecho 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 Derecho, there are several key commands worth
noting:
- qsub job-script - This command submits a PBS job script, which describes a job
to be run on one or more batch nodes.
- qstat -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.
- qdel JobID - This command deletes a queued or running job
Throughout this tutorial, you will submit your jobs to the 'main' queue with the 'premium' priority, using an account number UMMM0012.
At various points in the practical exercises, we'll need to submit jobs to Derecho's queueing system
using the qsub command, and after doing so, we may check on the status
of the job with the qstat command, monitoring the log files produced by the job
once we see that the job has begun to run. Most exercises will run with 64 cores of one single node
though each of Derecho's nodes have 128 cores (235GB available memory).
Note that 4DEnVar uses 3 nodes with the analyses 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/3.0.2' branch from mpas-bundle-3.0.2:
$ cd /glade/derecho/scratch/$USER/mpas_jedi_tutorial
$ mkdir mpas_bundle_v3
$ cd mpas_bundle_v3
$ git clone -b release/3.0.2 https://github.com/JCSDA/mpas-bundle code
The git command will clone the mpas-bundle repository from github to a local directory 'code'
, then make the 'release/3.0.2' 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's 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:
$ source ./code/env-setup/gnu-derecho.sh
$ module list
The output to the module list in the terminal should look like the following (102 modules!!!).
These packages are from pre-installed spack-stack-1.6.0.
Currently Loaded Modules:
1) ecflow/5.8.4 35) base-env/1.0.0 69) json-schema-validator/2.1.0
2) mysql/8.0.33 36) boost/1.83.0 70) odc/1.4.6
3) ncarenv/23.09 (S) 37) openblas/0.3.24 71) py-attrs/21.4.0
4) gcc/12.2.0 38) py-setuptools/63.4.3 72) py-pycparser/2.21
5) stack-gcc/12.2.0 39) py-numpy/1.22.3 73) py-cffi/1.15.1
6) craype/2.7.20 40) bufr/12.0.1 74) py-findlibs/0.0.2
7) cray-mpich/8.1.25 41) ecbuild/3.7.2 75) py-eccodes/1.5.0
8) libfabric/1.15.2.0 42) libpng/1.6.37 76) py-f90nml/1.4.3
9) cray-pals/1.2.11 43) openjpeg/2.3.1 77) py-h5py/3.7.0
10) stack-cray-mpich/8.1.25 44) eccodes/2.32.0 78) py-cftime/1.0.3.4
11) tar/1.34 45) eigen/3.4.0 79) py-netcdf4/1.5.8
12) gettext/0.21.1 46) eckit/1.24.5 80) py-bottleneck/1.3.7
13) libxcrypt/4.4.35 47) fftw/3.3.10 81) py-numexpr/2.8.4
14) zlib/1.2.13 48) fckit/0.11.0 82) py-et-xmlfile/1.0.1
15) sqlite/3.43.2 49) fiat/1.2.0 83) py-openpyxl/3.1.2
16) util-linux-uuid/2.38.1 50) ectrans/1.2.0 84) py-six/1.16.0
17) python/3.10.13 51) qhull/2020.2 85) py-python-dateutil/2.8.2
18) stack-python/3.10.13 52) atlas/0.35.1 86) py-pytz/2023.3
19) nghttp2/1.57.0 53) git-lfs/3.3.0 87) py-pyxlsb/1.0.10
20) curl/8.4.0 54) gsibec/1.1.3 88) py-xlrd/2.0.1
21) cmake/3.23.1 55) gsl-lite/0.37.0 89) py-xlsxwriter/3.1.7
22) git/2.41.0 56) libjpeg/2.1.0 90) py-xlwt/1.3.0
23) pkg-config/0.29.2 57) krb5/1.19.2 91) py-pandas/1.5.3
24) hdf5/1.14.0 58) libtirpc/1.3.3 92) py-pybind11/2.11.0
25) snappy/1.1.10 59) hdf/4.2.15 93) py-pycodestyle/2.11.0
26) zstd/1.5.2 60) jedi-cmake/1.4.0 94) py-pyhdf/0.10.4
27) c-blosc/1.21.5 61) libxt/1.1.5 95) libyaml/0.2.5
28) netcdf-c/4.9.2 62) libxmu/1.1.4 96) py-pyyaml/6.0
29) nccmp/1.9.0.1 63) libxpm/3.5.12 97) py-scipy/1.11.3
30) netcdf-fortran/4.6.1 64) libxaw/1.0.13 98) py-packaging/23.1
31) parallel-netcdf/1.12.2 65) udunits/2.2.28 99) py-xarray/2023.7.0
32) parallelio/2.5.10 66) ncview/2.1.9 100) sp/2.5.0
33) py-pip/23.1.2 67) netcdf-cxx4/4.3.1 101) jedi-base-env/1.0.0
34) wget/1.20.3 68) json/3.10.5 102) jedi-mpas-env/1.0.0
Where:
S: Module is Sticky, requires --force to unload or purge
Now we are ready to 'fetch' actual source code from github repositories through cmake.
MPAS-JEDI uses cmake to automatically fetch the code from various github repositories,
listed in CMakeLists.txt under ~code. Now type the command below under the build directory:
$ git lfs install
$ mkdir build
$ cd build
$ cmake ../code
After it completes (may take 15-20 min depending on network connection), you will see the actual source code of various
repositories (e.g., oops, vader, 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. Instead, better compiling code
using a compute node with an interactive job by issuing the command below:
WARNING:This cmake step could take 15-20 min to complete. You may continue to read instructions while waiting.
$ qsub -A ummm0012 -N build-bundle -q main -l job_priority=premium -l walltime=03:00:00 -l select=1:ncpus=128:mem=235GB -I
This requests an interactive job to the 'main' queue, with a job name 'build-bundle', 3 hours walltime of
a compute node with 128 cores, the 'premium' priority and an account number 'UMMM0012'.
Note that the interactive job does not have the spack-stack environment loaded, so you have to reload
the spack-stack environment (under ~build), and then compiling code with the parallel make using 20 cores:
$ source ../code/env-setup/gnu-derecho.sh
$ make -j20
This could take ~25 min. Using more cores for build could speed up a bit, but will not help too much from our experience.
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 derecho's compute nodes have a much slower internet connection.
WARNING:The compilation could take ~25 min to complete. You may continue to read instructions while waiting.
Once we reach 100% of the compilation, many executables will be generated under ~build/bin.
mpas model and mpas-jedi DA related executables are:
$ ls bin/mpas*
bin/mpas_atmosphere bin/mpasjedi_enshofx.x bin/mpasjedi_rtpp.x
bin/mpas_atmosphere_build_tables bin/mpasjedi_ens_mean_variance.x bin/mpasjedi_saca.x
bin/mpas_init_atmosphere bin/mpasjedi_error_covariance_toolbox.x bin/mpasjedi_variational.x
bin/mpasjedi_convertstate.x bin/mpasjedi_forecast.x bin/mpas_namelist_gen
bin/mpasjedi_converttostructuredgrid.x bin/mpasjedi_gen_ens_pert_B.x bin/mpas_parse_atmosphere
bin/mpasjedi_eda.x bin/mpasjedi_hofx3d.x bin/mpas_parse_init_atmosphere
bin/mpasjedi_enkf.x bin/mpasjedi_hofx.x bin/mpas_streams_gen
The last step is to ensure that the code was compiled properly by running the MPAS-JEDI ctests,
with two lines of simple command:
$ export LD_LIBRARY_PATH=/glade/derecho/scratch/${USER}/mpas_jedi_tutorial/mpas_bundle_v3/build/lib:$LD_LIBRARY_PATH
$ cd mpas-jedi
$ ctest
At the moment the tests are running (take ~5 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:
......
100% tests passed, 0 tests failed out of 59
Label Time Summary:
executable = 18.76 sec*proc (13 tests)
mpasjedi = 237.40 sec*proc (59 tests)
mpi = 235.68 sec*proc (56 tests)
script = 218.64 sec*proc (44 tests)
Total Test time (real) = 237.49 sec
WARNING:You could run ctest just under the 'build' directory, but that will run a total of 2159 ctest cases for all
component packages (oops, vader, saber, ufo, ioda, crtm, mpas-jedi etc.) in mpas-bundle, which will take much longer time.
Maybe something you can play with after this tutorial.
You may use 'ctest -N, which will only list names of ctest cases, but not run them).
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 4denvar_ID.run ens_mean_variance.run.ref
3denvar_2stream_bumploc.run 4denvar_ID.run.ref forecast.ref
3denvar_2stream_bumploc.run.ref 4denvar_VarBC_nonpar.run forecast.run
3denvar_amsua_allsky.ref 4denvar_VarBC_nonpar.run.ref forecast.run.ref
3denvar_amsua_allsky.run 4denvar_VarBC.ref gen_ens_pert_B.ref
3denvar_amsua_allsky.run.ref 4denvar_VarBC.run gen_ens_pert_B.run
3denvar_amsua_bc.ref 4denvar_VarBC.run.ref gen_ens_pert_B.run.ref
3denvar_amsua_bc.run 4dfgat.ref hofx3d_nbam.ref
3denvar_amsua_bc.run.ref 4dfgat.run hofx3d_nbam.run
3denvar_bumploc.ref 4dfgat.run.ref hofx3d_nbam.run.ref
3denvar_bumploc.run 4dhybrid_bumpcov_bumploc.ref hofx3d.ref
3denvar_bumploc.run.ref 4dhybrid_bumpcov_bumploc.run hofx3d_ropp.ref
3denvar_dual_resolution.ref 4dhybrid_bumpcov_bumploc.run.ref hofx3d_rttovcpp.ref
3denvar_dual_resolution.run convertstate_bumpinterp.ref hofx3d.run
3denvar_dual_resolution.run.ref convertstate_bumpinterp.run hofx3d.run.ref
3dfgat_pseudo.ref convertstate_bumpinterp.run.ref hofx4d_pseudo.ref
3dfgat_pseudo.run convertstate_unsinterp.ref hofx4d_pseudo.run
3dfgat_pseudo.run.ref convertstate_unsinterp.run hofx4d_pseudo.run.ref
3dfgat.ref convertstate_unsinterp.run.ref hofx4d.ref
3dfgat.run converttostructuredgrid_latlon.ref hofx4d.run
3dfgat.run.ref converttostructuredgrid_latlon.run hofx4d.run.ref
3dhybrid_bumpcov_bumploc.ref converttostructuredgrid_latlon.run.ref letkf_3dloc.ref
3dhybrid_bumpcov_bumploc.run dirac_bumpcov.ref letkf_3dloc.run
3dhybrid_bumpcov_bumploc.run.ref dirac_bumpcov.run letkf_3dloc.run.ref
3dvar_bumpcov_nbam.ref dirac_bumpcov.run.ref lgetkf_height_vloc.ref
3dvar_bumpcov_nbam.run dirac_bumploc.ref lgetkf_height_vloc.run
3dvar_bumpcov_nbam.run.ref dirac_bumploc.run lgetkf_height_vloc.run.ref
3dvar_bumpcov.ref dirac_bumploc.run.ref lgetkf.ref
3dvar_bumpcov_ropp.ref dirac_noloc.ref lgetkf.run
3dvar_bumpcov_rttovcpp.ref dirac_noloc.run lgetkf.run.ref
3dvar_bumpcov.run dirac_noloc.run.ref parameters_bumpcov.ref
3dvar_bumpcov.run.ref dirac_spectral_1.ref parameters_bumpcov.run
3dvar.ref dirac_spectral_1.run parameters_bumpcov.run.ref
3dvar.run dirac_spectral_1.run.ref parameters_bumploc.ref
3dvar.run.ref eda_3dhybrid.ref parameters_bumploc.run
4denvar_bumploc.ref eda_3dhybrid.run parameters_bumploc.run.ref
4denvar_bumploc.run eda_3dhybrid.run.ref rtpp.ref
4denvar_bumploc.run.ref ens_mean_variance.ref rtpp.run
4denvar_ID.ref ens_mean_variance.run rtpp.run.ref
After running ctest, you can terminate your interactive job by simply typing 'exit' or by:
$ qstat -u $USER
$ qdel job-id-number
The qstat command will return your job ID in a form of 'job-id-number.desched*'.
The qdel command will kill your interactive job by only supplying 'job-id-number'.
Note that the tutorial test cases are designed with single-precision MPAS model test dataset.
Therefore, we pre-compiled mpas-bundle code with single-precision mode in
/glade/derecho/scratch/liuz/mpas_bundle_v3.0.2_public_gnuSP
Tutorial attendees should use this pre-build bundle for mpas-jedi practicals and should
set an environment variable 'bundle_dir' for convenience in the subsequent practices.
$ export bundle_dir=/glade/derecho/scratch/liuz/mpas_bundle_v3.0.2_public_gnuSP
The goal of this session is to create the observation input files needed for running
MPAS-JEDI test cases.
In mpas_jedi_tutorial directory,
create a directory for converting NCEP BUFR observation into IODA-HDF5 format.
Then, clone the obs2ioda repository and proceed for compiling the converter.
$ cd /glade/derecho/scratch/${USER}/mpas_jedi_tutorial
$ git clone -b release/1.0.0 https://github.com/NCAR/obs2ioda
NCEP BUFR library (https://github.com/NOAA-EMC/NCEPLIBS-bufr) along with NETCDF-C and NETCDF-Fortran libraries are required to compile obs2ioda.
For this practice, we will use a pre-built BUFR library using GNU compilers, but users are encouraged to see NCEP BUFR library for further instructions on how to install it.
obs2ioda has been updated to use CMake as its primary build system. This change enables the use of the same environment file used for compiling MPAS-JEDI, simplifying setup and ensuring consistency across builds. For convenience, a tailored environment file is provided specifically for obs2ioda, which loads only the modules and packages required for its compilation. GNU compilers are used by default—consistent with MPAS-JEDI—though Intel compilers are also supported.
Load the appropriate compiler and library modules, create a new build directory and navigate into it. Then, run CMake to configure the build—set the CMAKE_BUILD_TYPE option to Release and specify the path to the NCEP BUFR library. Finally, compile the converter using the make command and run the ctest to verify that the build works correctly and validate functionalities:
$ cd obs2ioda/env-setup
$ source gnu_derecho.sh
$ cd ..
$ mkdir build && cd build
$ cmake ../ -DNCEP_BUFR_LIB=/glade/campaign/mmm/parc/ivette/pandac/codeBuild/NCEPLIBS-bufr_gnu/build/src/libbufr_4.a -DCMAKE_BUILD_TYPE=Release
$ make
$ ctest
If the compilation is successful, the executable file obs2ioda_v3 will be generated under the bin directory.
Once the compilation is complete, deactivate the obs2ioda environment:
$ conda deactivate
The NCEP BUFR input files are located in the existing ~/obs_bufr directory. The next step is to convert this BUFR-format data into IODAv3-HDF5 format and store the output in the ~/obs_ioda directory.
$ mkdir /glade/derecho/scratch/${USER}/mpas_jedi_tutorial/obs_ioda
$ cd /glade/derecho/scratch/${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 the executable and obs_errtable files to the working directory:
$ ln -fs ../../obs2ioda/build/bin/obs2ioda_v3 .
$ ln -sf ../../obs_bufr/2018041500/obs_errtable .
Now we can run obs2ioda_v3.
The usage is:
$ ./obs2ioda_v3 [-i input_dir] [-o output_dir] [bufr_filename(s)_to_convert]
If input_dir and output_dir are not specified, the current working directory is used by default.
If no bufr_filename(s)_to_convert is 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. For each file found, the corresponding conversion is performed; any files not present are skipped. If no obs_errtable file is found in the running directory, the script uses the observation errors embedded in the prepbufr file. Otherwise, it uses the user-provided obs_errtable.
$ ./obs2ioda_v3 prepbufr.gdas.20180415.t00z.nr.48h
$ ./obs2ioda_v3 gdas.satwnd.t00z.20180415.bufr
$ ./obs2ioda_v3 gdas.gpsro.t00z.20180415.bufr
$ ./obs2ioda_v3 gdas.1bamua.t00z.20180415.bufr
You should see *.h5 files in the folder.
Setup Python environment for plotting.
Before running the Python script, make sure the appropriate Python environment is set up by loading the NPL Conda environment.
$ module load conda
$ conda activate npl
Copy the graphics directory to your current working directory, then run the Python script:
$ export bundle_dir=/glade/derecho/scratch/liuz/mpas_bundle_v3.0.2_public_gnuSP
$ cp -r ${bundle_dir}/code/mpas-jedi/graphics .
$ cd graphics/standalone
$ ln -sf ../../../../MPAS_JEDI_yamls_scripts/plot_obs_loc_tut.py .
$ python plot_obs_loc_tut.py
Now, one of the generated figures is shown here:
In this session, we will be running an application called 𝐻(𝐱).
Creating hofx directory in mpas_jedi_tutorial,
$ cd /glade/derecho/scratch/${USER}/mpas_jedi_tutorial
$ mkdir hofx
$ cd hofx
Set environment variable for the mpas-bundle directory,
$ export bundle_dir=/glade/derecho/scratch/liuz/mpas_bundle_v3.0.2_public_gnuSP
Now, we need to 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.64 .
$ 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 ${bundle_dir}/code/mpas-jedi/test/testinput/namelists/geovars.yaml .
$ ln -fs ${bundle_dir}/code/mpas-jedi/test/testinput/namelists/keptvars.yaml .
$ ln -fs ${bundle_dir}/code/mpas-jedi/test/testinput/obsop_name_map.yaml .
Link MPAS 6-h forecast background file into template file and an MPAS invariant file for 2-stream I/O:
$ ln -fs ../MPAS_namelist_stream_physics_files/x1.10242.invariant.nc .
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc templateFields.10242.nc
Link the background and observation 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 .
Link hofx3d.yaml and run_hofx3d.csh to the working directory:
$ ln -sf ../MPAS_JEDI_yamls_scripts/hofx3d.yaml .
$ ln -sf ../MPAS_JEDI_yamls_scripts/run_hofx3d.csh .
Submit a PBS job to run hofx,
$ qsub run_hofx3d.csh
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 running plot_diag.py.
The following steps are for set up python environment, copy graphics directory, and run plot_diag.py:
$ module reset
$ module load conda
$ conda activate npl
$ cp -r ${bundle_dir}/code/mpas-jedi/graphics .
$ cd graphics/standalone
$ python plot_diag.py
This will produce a number of figure files, display one of them
display distri_windEastward_hofx_aircraft_omb_allLevels.png
which will look like the figure below
To view more figures, it is more convenient to transfer them to your local computer via:
$ scp -r $USER@data-access.ucar.edu:/glade/derecho/scratch/$USER/mpas_jedi_tutorial/hofx/graphics/standalone .
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 /glade/derecho/scratch/${USER}/mpas_jedi_tutorial directory,
create the working directory.
$ cd /glade/derecho/scratch/$USER/mpas_jedi_tutorial
$ mkdir localization
$ cd localization
Set environment variable for the mpas-bundle directory,
$ export bundle_dir=/glade/derecho/scratch/liuz/mpas_bundle_v3.0.2_public_gnuSP
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.64 ./
$ 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 ../MPAS_namelist_stream_physics_files/x1.10242.invariant.nc .
Link yaml files.
$ ln -fs ${bundle_dir}/code/mpas-jedi/test/testinput/namelists/geovars.yaml ./
$ ln -fs ${bundle_dir}/code/mpas-jedi/test/testinput/namelists/keptvars.yaml ./
$ ln -fs ../MPAS_JEDI_yamls_scripts/bumploc.yaml ./
Link and submit pbs job script.
$ ln -sf ../MPAS_JEDI_yamls_scripts/run_bumploc.csh ./
$ qsub run_bumploc.csh
User may check the PBS job status by running qstat -u $USER,
or monitoring the log file jedi.log.
After finishing the PBS job, user can check the BUMP localization files generated. The "local" files are generated for each processor (64 processors in this practice).
$ ls bumploc_1200km6km_nicas_local_*.nc
bumploc_1200km6km_nicas_local_000064-000001.nc bumploc_1200km6km_nicas_local_000064-000033.nc
bumploc_1200km6km_nicas_local_000064-000002.nc bumploc_1200km6km_nicas_local_000064-000034.nc
bumploc_1200km6km_nicas_local_000064-000003.nc bumploc_1200km6km_nicas_local_000064-000035.nc
bumploc_1200km6km_nicas_local_000064-000004.nc bumploc_1200km6km_nicas_local_000064-000036.nc
bumploc_1200km6km_nicas_local_000064-000005.nc bumploc_1200km6km_nicas_local_000064-000037.nc
bumploc_1200km6km_nicas_local_000064-000006.nc bumploc_1200km6km_nicas_local_000064-000038.nc
bumploc_1200km6km_nicas_local_000064-000007.nc bumploc_1200km6km_nicas_local_000064-000039.nc
bumploc_1200km6km_nicas_local_000064-000008.nc bumploc_1200km6km_nicas_local_000064-000040.nc
bumploc_1200km6km_nicas_local_000064-000009.nc bumploc_1200km6km_nicas_local_000064-000041.nc
bumploc_1200km6km_nicas_local_000064-000010.nc bumploc_1200km6km_nicas_local_000064-000042.nc
bumploc_1200km6km_nicas_local_000064-000011.nc bumploc_1200km6km_nicas_local_000064-000043.nc
bumploc_1200km6km_nicas_local_000064-000012.nc bumploc_1200km6km_nicas_local_000064-000044.nc
bumploc_1200km6km_nicas_local_000064-000013.nc bumploc_1200km6km_nicas_local_000064-000045.nc
bumploc_1200km6km_nicas_local_000064-000014.nc bumploc_1200km6km_nicas_local_000064-000046.nc
bumploc_1200km6km_nicas_local_000064-000015.nc bumploc_1200km6km_nicas_local_000064-000047.nc
bumploc_1200km6km_nicas_local_000064-000016.nc bumploc_1200km6km_nicas_local_000064-000048.nc
bumploc_1200km6km_nicas_local_000064-000017.nc bumploc_1200km6km_nicas_local_000064-000049.nc
bumploc_1200km6km_nicas_local_000064-000018.nc bumploc_1200km6km_nicas_local_000064-000050.nc
bumploc_1200km6km_nicas_local_000064-000019.nc bumploc_1200km6km_nicas_local_000064-000051.nc
bumploc_1200km6km_nicas_local_000064-000020.nc bumploc_1200km6km_nicas_local_000064-000052.nc
bumploc_1200km6km_nicas_local_000064-000021.nc bumploc_1200km6km_nicas_local_000064-000053.nc
bumploc_1200km6km_nicas_local_000064-000022.nc bumploc_1200km6km_nicas_local_000064-000054.nc
bumploc_1200km6km_nicas_local_000064-000023.nc bumploc_1200km6km_nicas_local_000064-000055.nc
bumploc_1200km6km_nicas_local_000064-000024.nc bumploc_1200km6km_nicas_local_000064-000056.nc
bumploc_1200km6km_nicas_local_000064-000025.nc bumploc_1200km6km_nicas_local_000064-000057.nc
bumploc_1200km6km_nicas_local_000064-000026.nc bumploc_1200km6km_nicas_local_000064-000058.nc
bumploc_1200km6km_nicas_local_000064-000027.nc bumploc_1200km6km_nicas_local_000064-000059.nc
bumploc_1200km6km_nicas_local_000064-000028.nc bumploc_1200km6km_nicas_local_000064-000060.nc
bumploc_1200km6km_nicas_local_000064-000029.nc bumploc_1200km6km_nicas_local_000064-000061.nc
bumploc_1200km6km_nicas_local_000064-000030.nc bumploc_1200km6km_nicas_local_000064-000062.nc
bumploc_1200km6km_nicas_local_000064-000031.nc bumploc_1200km6km_nicas_local_000064-000063.nc
bumploc_1200km6km_nicas_local_000064-000032.nc bumploc_1200km6km_nicas_local_000064-000064.nc
For a Dirac function multiplied by a given localization function, user can make a plot with dirac_nicas.2018-04-15_00.00.00.nc .
Note that we need to load the python environment before executing the python script.
If you have already loaded the spack-stack module, please reset the module environment to properly load the NCAR Python Library (NPL).
$ module reset
$ module load conda
$ conda activate npl
Copy the plotting script for Dirac output.
$ cp ../MPAS_JEDI_yamls_scripts/plot_bumploc.py ./
$ python plot_bumploc.py
$ display figure_bumploc.png
User may play with different localization lengths (e.g., change the horizontal one from 1200km to 2500km) by
changing the following YAML keys from bumploc.yaml, then re-submit the PBS job script.
[Note that 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 the first analysis test using 3DEnVar. For this 3DEnVar, we will use a 240 km background forecast and ensemble input at the same resolution.
Let's create a new working directory for this test:
$ cd /glade/derecho/scratch/${USER}/mpas_jedi_tutorial
$ mkdir 3denvar
$ cd 3denvar
Set environment variable for the mpas-bundle directory,
$ export bundle_dir=/glade/derecho/scratch/liuz/mpas_bundle_v3.0.2_public_gnuSP
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.info.part.64 ./
$ 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 ${bundle_dir}/code/mpas-jedi/test/testinput/namelists/geovars.yaml ./
$ ln -sf ${bundle_dir}/code/mpas-jedi/test/testinput/namelists/keptvars.yaml ./
$ ln -sf ${bundle_dir}/code/mpas-jedi/test/testinput/obsop_name_map.yaml ./
Next, we need to link MPAS 2-stream files corresponding to the 240 km mesh.
For this experiment we need to link the following files:
$ ln -fs ../MPAS_namelist_stream_physics_files/x1.10242.invariant.nc .
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./templateFields.10242.nc
Next, we should link the 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 bump localization files.
$ ln -fs ../localization ./BUMP_files
Next, we need to link observation files into the working directory.
$ 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 .
Next, we are ready to link the pre-prepared 3denvar.yaml file for 3DEnVar. 3denvar.yaml contains all of the settings we need for this 3DEnVar update.
$ ln -sf ../MPAS_JEDI_yamls_scripts/3denvar.yaml ./
We are finally ready to link a PBS job script and then submit it. Notice the job 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.csh ./
$ qsub run_3denvar.csh
Once your job is submitted, you can check the status of the job using qstat -u ${USER}.
When the job is complete, you can check the log file by using vi mpasjedi_3denvar.log.
Now that the analysis was successful, you may plot the cost function and gradient norm reduction.
$ module load conda
$ conda activate npl
$ cp ${bundle_dir}/code/mpas-jedi/graphics/standalone/plot_cost_grad.py .
Changing 'jedi.log' to 'mpasjedi_3denvar.log' in plot_cost_grad.py, then do
$ python plot_cost_grad.py
You can then view the figure using the following:
$ display costgrad.png
You can also plot the analysis increment (i.e., an.*.nc - bg.*.nc) with the following.
$ cp ../MPAS_JEDI_yamls_scripts/plot_analysis_increment.py ./
$ python plot_analysis_increment.py
$ display figure_increment_uReconstructZonal_10.png
You will also see several so-called DA "feedback" files below:
obsout_da_aircraft.h5
obsout_da_gnssro.h5
obsout_da_satwind.h5
obsout_da_sfc.h5
obsout_da_sondes.h5
You can follow the same procedure as in section 2.2 to plot these 'feedback' files, e.g.,
$ module load conda
$ conda activate npl
$ cp -r ${bundle_dir}/code/mpas-jedi/graphics .
$ cd graphics/standalone
$ python plot_diag.py
This will produce a number of figure files, display one of them
display RMSofOMM_da_aircraft_airTemperature.png
which will look like the figure below
Now we know how to run a 3DEnVar analysis, 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 /glade/derecho/scratch/${USER}/mpas_jedi_tutorial
$ mkdir 3denvar_120km240km
$ cd 3denvar_120km240km
Set environment variable for the mpas-bundle directory,
$ export bundle_dir=/glade/derecho/scratch/liuz/mpas_bundle_v3.0.2_public_gnuSP
We can then link the relevant MPAS stream, graph, namelist, and physics files. In addition to files for 240km mesh, we also need files for 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.info.part.64 ./
$ ln -sf ../MPAS_namelist_stream_physics_files/x1.40962.graph.info.part.64 ./
$ 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 ${bundle_dir}/code/mpas-jedi/test/testinput/namelists/geovars.yaml ./
$ ln -sf ${bundle_dir}/code/mpas-jedi/test/testinput/namelists/keptvars.yaml ./
$ ln -sf ${bundle_dir}/code/mpas-jedi/test/testinput/obsop_name_map.yaml ./
Next, we need to link MPAS 2-stream files corresponding to both 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 ../MPAS_namelist_stream_physics_files/x1.10242.invariant.nc .
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./templateFields.10242.nc
$ ln -fs ../MPAS_namelist_stream_physics_files/x1.40962.invariant.nc .
$ ln -fs ../background_120km/2018041418/mpasout.2018-04-15_00.00.00.nc ./templateFields.40962.nc
Next, we should link 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 observation files into the working directory.
$ 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 .
Next, we are ready to copy the pre-prepared 3denvar_120km240km.yaml file for this dual-resolution 3DEnVar.
$ ln -sf ../MPAS_JEDI_yamls_scripts/3denvar_120km240km.yaml ./3denvar.yaml
We are finally ready to link a PBS job script and then submit the job.
$ ln -fs ../MPAS_JEDI_yamls_scripts/run_3denvar.csh ./
$ qsub run_3denvar.csh
Again, you can check the status of the job using qstat -u ${USER}.
When the job is complete, you can check the log file by using vi mpasjedi_3denvar.log.
Again, you may plot cost function and gradient norm reduction as well as the analysis increment.
You need to change "x1.10242.invariant.nc" to "x1.40962.invariant.nc" within plot_analysis_increment.py in order to plot increment at the 120 km mesh.
Now we know how to run a 3DEnVar analysis, let's try 4DEnVar with 3 subwindows.
Let's create a new working directory to run 240km 4DEnVar:
$ cd /glade/derecho/scratch/${USER}/mpas_jedi_tutorial
$ mkdir 4denvar
$ cd 4denvar
Set environment variable for the mpas-bundle directory,
$ export bundle_dir=/glade/derecho/scratch/liuz/mpas_bundle_v3.0.2_public_gnuSP
We can then link 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.info.part.64 ./
$ 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 a few yaml files again that define MPAS variables for JEDI:
$ ln -sf ${bundle_dir}/code/mpas-jedi/test/testinput/namelists/geovars.yaml ./
$ ln -sf ${bundle_dir}/code/mpas-jedi/test/testinput/namelists/keptvars.yaml ./
$ ln -sf ${bundle_dir}/code/mpas-jedi/test/testinput/obsop_name_map.yaml ./
Next, we need to link MPAS 2-stream files corresponding to 240km mesh.
For this experiment we need to link the following files:
$ ln -fs ../MPAS_namelist_stream_physics_files/x1.10242.invariant.nc .
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./templateFields.10242.nc
Next, we need to link 240km 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 bump localization files.
$ ln -fs ../localization ./BUMP_files
Next, we need to link observation files.
$ 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 .
Next, we are ready to link the pre-prepared 4denvar.yaml file.
$ ln -sf ../MPAS_JEDI_yamls_scripts/4denvar.yaml ./4denvar.yaml
We are finally ready to link a PBS job script and then submit a job. Note that this 4DEnVar job uses 3 nodes with 64 cores per node for 3 subwindows.
$ ln -fs ../MPAS_JEDI_yamls_scripts/run_4denvar.csh ./
$ qsub run_4denvar.csh
Again, you can check the status of the job using qstat -u ${USER}.
When the job is complete, you can check the log file by using vi mpasjedi_4denvar.log.
Again, you can plot cost function and gradient norm reduction as well as analysis increments at one of 3 subwindows.
You may 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 /glade/derecho/scratch/${USER}/mpas_jedi_tutorial directory,
create the working directory for 3dvar.
$ cd /glade/derecho/scratch/$USER/mpas_jedi_tutorial
$ mkdir 3dvar
$ cd 3dvar
Set environment variable for the mpas-bundle directory,
$ export bundle_dir=/glade/derecho/scratch/liuz/mpas_bundle_v3.0.2_public_gnuSP
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.64 ./
$ 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 ../MPAS_namelist_stream_physics_files/x1.10242.invariant.nc .
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./templateFields.10242.nc
Prepare the pre-generated static B files.
$ ln -fs ../B_Matrix ./
Next, link the observation files to be assimilated.
Along 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
For this 3DVar test, satellite radiance data from NOAA-18's AMSU-A sensor is assimilated, so need coefficient files for CRTM.
$ ln -fs ../crtm_coeffs_v3 ./
Link yaml files:
$ ln -fs ${bundle_dir}/code/mpas-jedi/test/testinput/namelists/geovars.yaml ./
$ ln -fs ${bundle_dir}/code/mpas-jedi/test/testinput/namelists/keptvars.yaml ./
$ ln -fs ${bundle_dir}/code/mpas-jedi/test/testinput/obsop_name_map.yaml ./
$ ln -sf ../MPAS_JEDI_yamls_scripts/3dvar.yaml ./
Link and submit pbs job script:
$ ln -sf ../MPAS_JEDI_yamls_scripts/run_3dvar.csh ./
$ qsub run_3dvar.csh
User may check the PBS job status by running qstat -u $USER,
or monitoring the log file tail -f mpasjedi_3dvar.log .
You may make figures following similar procedures in previous sections.
You may also try to assimilate additional satellite radiance data.
To do this, additional obs file (for example, ../obs_ioda/2018041500/amsua_n19_obs_2018041500.h5)
should be linked to /glade/derecho/scratch/$USER/mpas_jedi_tutorial/3dvar/.
Also, user need to add the corresponding obs space yaml key in 3dvar.yaml .
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 /glade/derecho/scratch/${USER}/mpas_jedi_tutorial directory,
create the working directory for Hybrid-3DEnVar.
$ cd /glade/derecho/scratch/$USER/mpas_jedi_tutorial
$ mkdir hybrid-3denvar
$ cd hybrid-3denvar
Set environment variable for the mpas-bundle directory,
$ export bundle_dir=/glade/derecho/scratch/liuz/mpas_bundle_v3.0.2_public_gnuSP
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.64 ./
$ 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 ../MPAS_namelist_stream_physics_files/x1.10242.invariant.nc .
$ ln -fs ../background/2018041418/mpasout.2018-04-15_00.00.00.nc ./templateFields.10242.nc
Prepare the pre-generated static B files.
$ ln -fs ../B_Matrix ./
Link the localization files.
$ ln -fs ../localization ./BUMP_files
Link the observation files to be assimilated.
Along 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
Prepare the coefficient files for CRTM.
$ ln -fs ../crtm_coeffs_v3 ./
Link yaml files.
$ ln -fs ${bundle_dir}/code/mpas-jedi/test/testinput/namelists/geovars.yaml ./
$ ln -fs ${bundle_dir}/code/mpas-jedi/test/testinput/namelists/keptvars.yaml ./
$ ln -fs ${bundle_dir}/code/mpas-jedi/test/testinput/obsop_name_map.yaml ./
$ ln -fs ../MPAS_JEDI_yamls_scripts/3dhyb.yaml ./
Link and submit pbs job script.
$ ln -sf ../MPAS_JEDI_yamls_scripts/run_3dhyb.csh ./
$ qsub run_3dhyb.csh
User may check the PBS job status by running qstat -u $USER,
or monitoring the log file tail -f mpasjedi_3dhyb.log .
User may want to try additional 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
...
You may see some changes in the fit-to-obs statistics, convergence, or analysis increment fields.
In this section, we will make plots for the statistics of OmB/OmA feedback files from two 2-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,
then create the plots from statistics files following the instructions 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 (i.e., obsout_da_*.h5). 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 named 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:
$ cd /glade/derecho/scratch/$USER/mpas_jedi_tutorial/omboma_from2experiments
This folder contains sample variational-based DA observation "feedback" files from two experiments Exp1 and Exp2, e.g.,
$ ls Exp1/CyclingDA/*
You will see obsout_da_amsua_n18.h5 and obsout_da_gnssrobndropp1d.h5 under each cycle time during a 2-day period (9 cycles).
We can generate the bias, RMSE, and standard deviation statistics of these OMB/OMA files using
$ ln -sf ../MPAS_JEDI_yamls_scripts/generate_stats_omboma.csh .
$ ln -sf ../MPAS_JEDI_yamls_scripts/advanceCYMDH.py .
$ qsub generate_stats_omboma.csh
The job will take ~5 minutes to complete, and you can take a look at the folder of the statistics files:
$ ls Exp1/CyclingDA/*/stats
You will see two files stats_da_amsua_n18.h5 and stats_da_gnssrobndropp1d.h5 in each folder.
You may check generate_stats_omboma.csh, in which a line below generated statistical files.
python DiagnoseObsStatistics.py -n ${NUMPROC} -p ${daDiagDir} -o obsout -app variational -nout 1
where
-n: number of tasks/processors for multiprocessing
-p: path to DA feedback files (e.g., Exp1/CyclingDA/2018041500)
-o: prefix for DA feedback files
-app: application (variational or hofx)
-nout: number of outer loops
Also note that advanceCYMDH.py is a python code to help advance the time.
If statistics files were created, we can proceed to make the plots in the next section!
Set environment variable for the mpas-bundle directory,
$ export bundle_dir=/glade/derecho/scratch/liuz/mpas_bundle_v3.0.2_public_gnuSP
Let us make a local copy of the graphics package first:
$ cd /glade/derecho/scratch/$USER/mpas_jedi_tutorial/omboma_from2experiments
$ cp -r ${bundle_dir}/code/mpas-jedi/graphics .
$ cd graphics
Then diff and replace analyze_config.py with a pre-prepared one:
$ diff analyze_config.py ../../MPAS_JEDI_yamls_scripts/analyze_config.py
$ mv analyze_config.py analyze_config.py_default
$ cp ../../MPAS_JEDI_yamls_scripts/analyze_config.py .
You can see the changes made in analyze_config.py by diff.
Below provides some explanations about settings 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,17,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','/glade/derecho/scratch/'+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
Last, setting python environment and do
$ module reset
$ module load conda
$ conda activate npl
$ ./SpawnAnalyzeStats.py -a UMMM0012 -app variational -nout 1 -d gnssrobndropp1d
$ ./SpawnAnalyzeStats.py -a UMMM0012 -app variational -nout 1 -d amsua_n18
Each SpawnAnalyzeStats.py command triggers two PBS jobs for making two types of
figures under gnssrobndropp1d_analyses and amsua_n18_analyses directories.
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
SpawnAnalyzeStats.py can take some arguments (similar as in previous section), such as:
-a: account
-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)
After jobs complete, you may check gnssrobndropp1d_analyses/CYAxisExpLines/omm, for example,
$ cd gnssrobndropp1d_analyses/CYAxisExpLines/omm
$ display QCflag_good_TSeries_0min_gnssrobndropp1d_omm_RMS.pdf
which will look like the figure below
To view more figures, it is more convenent to transfer them to your local computer via:
$ scp -r $USER@data-access.ucar.edu:/glade/derecho/scratch/$USER/mpas_jedi_tutorial/omboma_from2experiments/graphics/gnssrobndropp1d_analyses .
$ scp -r $USER@data-access.ucar.edu:/glade/derecho/scratch/$USER/mpas_jedi_tutorial/omboma_from2experiments/graphics/amsua_n18_analyses .
In this section, we will demonstrate how to run a cycling data assimilation experiment.
Let’s begin by setting up the folder structure needed for the cycling experiment. Within the mpas_jedi_tutorial directory, create a new directory named 3dvarCyclingDA to set up the data assimilation (DA) experiment for 0000 UTC on 15 April 2018 and run a forecast from the generated analysis. This experiment will use MPAS-JEDI with the 3DVar configuration on a global 240 km mesh. Additionally, create subdirectories for the DA and executing the 6-hour MPAS-A forecasts during each cycle.
$ export mpasjedi_tutorial=/glade/derecho/scratch/${USER}/mpas_jedi_tutorial
$ cd ${mpasjedi_tutorial}
$ mkdir 3dvarCyclingDA
$ cd 3dvarCyclingDA
$ mkdir CyclingFC
$ mkdir CyclingDA
For the analysis at 0000 UTC 15 April 2018, we need a 6-hour MPAS-A forecast initialized at 1800 UTC on 14 April 2018 to serve as the background. To generate this forecast, we'll apply the methods learned in the MPAS-A tutorial—using GFS initial conditions along with a pre-generated invariant file that contains all time-invariant fields.
Navigate to the CyclingFC directory and create a folder named 2018041418 for the first cycle date:
$ cd CyclingFC
$ mkdir 2018041418
$ cd 2018041418
For later substitutions in the configuration YAML file, we construct the corresponding date and time strings in the required formats.
$ export thisValidDate=2018041418
$ export thisMPASFileDate=2018-04-14_18.00.00
$ export thisMPASNamelistDate=2018-04-14_18:00:00
$ export nextMPASFileDate=2018-04-15_00.00.00
Link the invariant file and the GFS initial conditions on the MPAS-A mesh into the cycle directory and create a symbolic link to the appropriate init file that will be used to update sea surface temperature (SST) and sea-ice fractional area coverage (XICE). These init files, used for surface updates, have been pre-generated using methods covered in the MPAS-A tutorial. Specifically, the invariant file was placed in the "input" stream, and the "filename_template" for the "output" stream was set to produce the init files (details not covered here).
$ ln -sf ${mpasjedi_tutorial}/MPAS_namelist_stream_physics_files/x1.10242.invariant.nc .
$ ln -sf ${mpasjedi_tutorial}/cyclingDA/init/${thisValidDate}/x1.10242.init.${thisMPASFileDate}.nc ./mpasin.${thisMPASFileDate}.nc
$ ln -sf ${mpasjedi_tutorial}/cyclingDA/init/${thisValidDate}/x1.10242.init.${thisMPASFileDate}.nc ./x1.10242.sfc_update.nc
Link MPAS physics-related data and tables, MPAS-A graph, MPAS-A streams, and MPAS-A namelist files:
$ ln -sf ${mpasjedi_tutorial}/MPAS_namelist_stream_physics_files/*BL ./
$ ln -sf ${mpasjedi_tutorial}/MPAS_namelist_stream_physics_files/*DATA ./
$ ln -sf ${mpasjedi_tutorial}/MPAS_namelist_stream_physics_files/x1.10242.graph.info.part.64 ./
$ cp -v ${mpasjedi_tutorial}/cyclingDA/namelist.atmosphere_240km ./namelist.atmosphere
$ ln -sf ${mpasjedi_tutorial}/cyclingDA/streams.atmosphere_240km ./streams.atmosphere
$ ln -sf ${mpasjedi_tutorial}/cyclingDA/stream_list* ./
Note that we use the "invariant" and "da_state" streams in streams.atmosphere.
immutable_stream name="invariant"
type="input"
precision="single"
filename_template="x1.10242.invariant.nc"
io_type="pnetcdf,cdf5"
input_interval="initial_only"
immutable_stream name="da_state"
type="output"
precision="single"
clobber_mode="truncate"
filename_template="mpasout.$Y-$M-$D_$h.$m.$s.nc"
packages="jedi_da"
io_type="pnetcdf,cdf5"
output_interval="06:00:00"
Also, note that config_jedi_da is set to True in namelist.atmosphere. For cold-start model runs config_do_DAcycling needs to be set to False, set configDODACycling to False and make the variable substitution in the namelist. Also, substitute the start time in the namelist.
$ export configDODACycling=false
$ sed -i 's@{{configDODACycling}}@'${configDODACycling}'@g' namelist.atmosphere
$ sed -i 's@{{startTime}}@'${thisMPASNamelistDate}'@g' namelist.atmosphere
Now, run the model submitting a job with qsub, but first copy the prepared job script to the working directory.
$ ln -sf ${mpasjedi_tutorial}/cyclingDA/run_model.pbs .
$ qsub run_model.pbs
You can monitor the status of the job using qstat -u ${USER}, and follow the model's progress in real time with tail -f log.atmosphere.0000.out.
At this resolution, the model runs quite fast, so you should see it complete shortly. Once finished, the output directory will contain the expected mpasout files for the 0- and 6-hour forecasts:
$ ls -l mpasout.*
With the 6-hour background already generated, we're now ready to begin the DA cycle. Navigate to the CyclingDA directory and create a folder for the current cycle.
$ cd ${mpasjedi_tutorial}/3dvarCyclingDA/CyclingDA
$ mkdir 2018041500
$ cd 2018041500
For later substitutions in the configuration YAML file, we define several date related variables, including the current valid time, the previous valid time, and the initial time of the DA time window. We then construct the corresponding date and time strings in the required formats.
$ export prevValidDate=2018041418
$ export thisValidDate=2018041500
$ export thisFileDate=2018-04-15T00:00:00
$ export thisMPASFileDate=2018-04-15_00.00.00
$ export thisMPASNamelistDate=2018-04-15_00:00:00
$ export halfprevDate=2018-04-14T21:00:00
$ export nextMPASFileDate=2018-04-15_06.00.00
Now, follow similar steps as outlined in Section 5.1, but this time for the 3DVar configuration. Link the required files, including MPAS-A physics data and tables, the MPAS-A graph file, streams configuration, and namelist files.
$ ln -sf ${mpasjedi_tutorial}/MPAS_namelist_stream_physics_files/*BL ./
$ ln -sf ${mpasjedi_tutorial}/MPAS_namelist_stream_physics_files/*DATA ./
$ ln -sf ${mpasjedi_tutorial}/MPAS_namelist_stream_physics_files/x1.10242.graph.info.part.64 ./
$ cp -v ${mpasjedi_tutorial}/cyclingDA/namelist.atmosphere_240km ./
$ ln -sf ${mpasjedi_tutorial}/MPAS_namelist_stream_physics_files/streams.atmosphere_240km ./
$ ln -sf ${mpasjedi_tutorial}/MPAS_namelist_stream_physics_files/stream_list* ./
For cycling DA runs config_do_DAcycling needs to be set to true, thus, set the variable configDODACycling to True and make the variable substitution in the namelist. In addtion, substitute the start date entry in the namelist file to match the current cycle.
$ export configDODACycling=true
$ sed -i 's@{{configDODACycling}}@'${configDODACycling}'@g' namelist.atmosphere_240km
$ sed -i 's@{{startTime}}@'${thisMPASNamelistDate}'@g' namelist.atmosphere_240km
Link the background file mpasout.${thisMPASFileDate}.nc into the working directory, and make a copy of it to serve as the analysis file, which will be overwritten during the DA step.
$ ln -fs ../../CyclingFC/${prevValidDate}/mpasout.${thisMPASFileDate}.nc ./bg.${thisMPASFileDate}.nc
$ cp -v ../../CyclingFC/${prevValidDate}/mpasout.${thisMPASFileDate}.nc ./an.${thisMPASFileDate}.nc
Link MPAS-A 2-stream I/O files corresponding to the 240 km mesh.
$ ln -fs ${mpasjedi_tutorial}/MPAS_namelist_stream_physics_files/x1.10242.invariant.nc .
$ ln -fs ../../CyclingFC/${prevValidDate}/mpasout.${thisMPASFileDate}.nc ./templateFields.10242.nc
Note that the templateFields file is used as the "input" stream in streams.atmosphere.
immutable_stream name="input"
type="input"
precision="single"
filename_template="templateFields.10242.nc"
io_type="pnetcdf,cdf5"
input_interval="initial_only"
Link the pre-generated static B files.
$ ln -fs ${mpasjedi_tutorial}/B_Matrix ./
Link the observation files to be assimilated.
$ ln -fs ${mpasjedi_tutorial}/obs_ioda/${thisValidDate}/aircraft_obs_${thisValidDate}.h5 .
$ ln -fs ${mpasjedi_tutorial}/obs_ioda/${thisValidDate}/gnssro_obs_${thisValidDate}.h5 .
$ ln -fs ${mpasjedi_tutorial}/obs_ioda/${thisValidDate}/satwind_obs_${thisValidDate}.h5 .
$ ln -fs ${mpasjedi_tutorial}/obs_ioda/${thisValidDate}/sfc_obs_${thisValidDate}.h5 .
$ ln -fs ${mpasjedi_tutorial}/obs_ioda/${thisValidDate}/sondes_obs_${thisValidDate}.h5 .
$ ln -fs ${mpasjedi_tutorial}/obs_ioda/${thisValidDate}/amsua_n18_obs_${thisValidDate}.h5 .
Link the coefficient files for CRTM.
$ ln -fs ${mpasjedi_tutorial}/crtm_coeffs_v3 ./
Link YAML files that define MPAS-A variables for JEDI.
$ export bundle_dir=/glade/derecho/scratch/liuz/mpas_bundle_v3.0.2_public_gnuSP
$ ln -fs ${bundle_dir}/code/mpas-jedi/test/testinput/namelists/geovars.yaml ./
$ ln -fs ${bundle_dir}/code/mpas-jedi/test/testinput/namelists/keptvars.yaml ./
$ ln -fs ${bundle_dir}/code/mpas-jedi/test/testinput/obsop_name_map.yaml ./
Next, we are ready to copy the template 3dvar_cyclingDA.yaml file for 3DVar and make the following substitutions:
$ cp -v ${mpasjedi_tutorial}/cyclingDA/3dvar_cyclingDA.yaml ./
$ sed -i 's@{{prevDate}}@'${prevValidDate}'@g' 3dvar_cyclingDA.yaml
$ sed -i 's@{{thisValidDate}}@'${thisValidDate}'@g' 3dvar_cyclingDA.yaml
$ sed -i 's@{{thisFileDate}}@'${thisFileDate}'@g' 3dvar_cyclingDA.yaml
$ sed -i 's@{{thisMPASFileDate}}@'${thisMPASFileDate}'@g' 3dvar_cyclingDA.yaml
$ sed -i 's@{{begin}}@'${halfprevDate}'@g' 3dvar_cyclingDA.yaml
Link and submit pbs job script.
$ ln -sf ${mpasjedi_tutorial}/cyclingDA/run_3dvar.csh ./
$ qsub run_3dvar.csh
You can check the PBS job status using qstat -u $USER, or monitor its progress in real time by running tail -f mpasjedi_3dvar.log.
After the analysis step completes successfully, the updated analysis file is used as the initial condition to run MPAS-A and generate the next 6-hour forecast. To proceed, navigate to the CyclingFC directory and create a subdirectory for the current cycle, named 2018041500.
$ cd ${mpasjedi_tutorial}/3dvarCyclingDA/CyclingFC
$ mkdir ${thisValidDate}
$ cd ${thisValidDate}
Link the analysis file and rename it to mpasin.
$ ln -sf ../../CyclingDA/${thisValidDate}/an.${thisMPASFileDate}.nc mpasin.${thisMPASFileDate}.nc
Link MPAS-A physics-related data and tables, MPAS-A graph, MPAS-A streams, and MPAS-A namelist files.
$ ln -sf ${mpasjedi_tutorial}/MPAS_namelist_stream_physics_files/*BL ./
$ ln -sf ${mpasjedi_tutorial}/MPAS_namelist_stream_physics_files/*DATA ./
$ ln -sf ${mpasjedi_tutorial}/MPAS_namelist_stream_physics_files/x1.10242.graph.info.part.64 ./
$ cp -v ${mpasjedi_tutorial}/cyclingDA/namelist.atmosphere_240km ./namelist.atmosphere
$ ln -sf ${mpasjedi_tutorial}/cyclingDA/streams.atmosphere_240km ./streams.atmosphere
$ ln -sf ${mpasjedi_tutorial}/cyclingDA/stream_list* ./
For warm-start model runs set the variable configDODACycling to True and make the variable substitution in the namelist. Also, substitute the start time in the namelist.
$ export configDODACycling=true
$ sed -i 's@{{configDODACycling}}@'${configDODACycling}'@g' namelist.atmosphere
$ sed -i 's@{{startTime}}@'${thisMPASNamelistDate}'@g' namelist.atmosphere
Link the invariant file and GFS initial conditions on the MPAS-A mesh for this current cycle and create a link of the init file to be used for sea-surface temperature and fractional area coverage of sea-ice update:
$ ln -sf ${mpasjedi_tutorial}/MPAS_namelist_stream_physics_files/x1.10242.invariant.nc .
$ ln -sf ${mpasjedi_tutorial}/cyclingDA/init/${thisValidDate}/x1.10242.init.${thisMPASFileDate}.nc ./
$ ln -sf ${mpasjedi_tutorial}/cyclingDA/init/${thisValidDate}/x1.10242.init.${thisMPASFileDate}.nc ./x1.10242.sfc_update.nc
Run the model submitting a job with qsub, but first copy the prepared job script to the working directory.
$ ln -sf ${mpasjedi_tutorial}/cyclingDA/run_model.pbs .
$ qsub run_model.pbs
Again, monitor the status of the job by issuing qstat -u ${USER} and follow the progress of the model run with tail -f log.atmosphere.0000.out. Once the model has finished successfully, we should see the following mpasout files in the directory for the 0- and 6-h forecast:
$ ls -l mpasout.*
Once the 6-hour forecast is generated, subsequent cycles are similar to the first one.
And with that, we have completed a one-cycle data assimilation experiment!!!
In this final section, we will go through how to generate the multivariate background error covariance statistics from the training dataset.
MPAS-JEDI uses some of external data processing (i.e., NCL and NCO) out of JEDI framework on its training dataset.
So this section consists of two parts: one for preprocessing of training dataset and the other for actual training procedure using the MPAS-JEDI executable, mpasjedi_error_covariance_toolbox.x.
Let's reset some module environment and move to the working directory.
$ module reset
$ module load conda
$ conda activate npl
$ cd /glade/derecho/scratch/$USER/mpas_jedi_tutorial/Bflow_preprocessing
Users will see the five shell scripts, which do smaller steps.
The major step will be 3_convert_uv_to_psichi.bash, but other steps are also included here to simplify the following subsection.
$ bash 1_generate_ESMF_weights.bash
This script will generate the NetCDF files of ESMF interpolation weights between MPAS mesh (global quasi uniform 120 km in this hands-on) and regular latitude-longitude grids (1 degree in this hands-on) via NCL.
Two interpolation weight files will be generated ( both directions; from MPAS to lat/lon; from lat/lon to MPAS) as follow:
ESMF_weights/latlon_1p0_to_MPAS_x1.40962_bilinear.nc
ESMF_weights/MPAS_x1.40962_to_latlon_1p0_bilinear.nc
$ bash 2_generate_template_PTB.bash
This script will generate a template file only including “stream_function” and “velocity_potential” via NCO commands.
This is especially convenient when we work with a large size of file.
template_PTB.nc
$ bash 3_convert_uv_to_psichi.bash 24
$ bash 3_convert_uv_to_psichi.bash 48
This script will convert the zonal and meridional winds at MPAS cell center into stream function and velocity_potential via NCL.
The zonal and meridional winds on the MPAS mesh will be interpolated to lat/lon grid, then converted as “stream_function” and velocity potential with “uv2sfvpf“ function.
Then, the stream_function and velocity_potential on the 1 degree lat/lon grid will be interpolated back to MPAS mesh.
This will be performed across range of model initialization time (from 2018040100 to 2018043000 with 12 hour interval in this hands-on example) for two different forecast lead times (24 hr and 48 hr).
To use the computing resource more efficiently, multiple NCL executions are performed with a PBS job script.
Users can find the outputs as follow:
output/YYYYMMDDHH/FULL_f24.nc # 24 hour forecast lead time
output/YYYYMMDDHH/FULL_f48.nc # 48 hour forecast lead time
Here, YYYYMMDDHH represents the model valid time.
At this point, FULL_fHH.nc files only contain two full field variables, which are “steram_function” and “velocity_potential".
$ bash 4_add_variables.bash 24
$ bash 4_add_variables.bash 48
This step adds additional variables to FULL_fXX.nc files with NCO commands.
“Temperature” and “spechum" are converted from other variables that are available from MPAS initial files, while “surface_pressure” is directly added.
“uReconstructZonal,uReconstructMeridional” are added in case user may want generate the univariate B statistics.
Users may also add the cloud hydrometeor variables if they are available.
This will be performed across range of model initialization time (from 2018040100 to 2018043000 with 12 hour interval in this hands-on example) for two different forecast lead times (24 hr and 48 hr).
To use the computing resource more efficiently, multiple bash shell scripts are executed with a PBS job script.
Users can find the outputs as follow:
output/YYYYMMDDHH/FULL_f24.nc # 24 hour forecast lead time
output/YYYYMMDDHH/FULL_f48.nc # 48 hour forecast lead time
Here, YYYYMMDDHH represents the model valid time.
At this point, FULL_fHH.nc files contain the following full field variables: “steram_function”, “velocity_potential”, “temperature”, “spechum”, “surface_pressure”, and (optionally) uReconstructZonal,uReconstructMeridional, relhum, qc, qi, qr, qs, qg.
$ bash 5_ncdiff.bash
In this step, the purturbation files are generated by subtracting the full fields of 24 hour forecast lead time from the full fields of 48 hour forecast lead time with NCO commands.
These perturbation files will be used in the next section for further B training.
This will be performed across range of model **valid** times (from 2018040300 to 2018050100 with 12 hour interval in this hands-on example).
To use the computing resource more efficiently, multiple bash shell scripts are executed with a PBS job script.
Users can find the outputs as follow:
output/YYYYMMDDHH/PTB_f48mf24.nc
Here, YYYYMMDDHH represents the model valid time.
$ cd /glade/derecho/scratch/$USER/mpas_jedi_tutorial/Bflow_global120km
Users can find several shell scripts, which do smaller steps, and plot_diagnostics directory, which includes several diagnostic scripts.
$ bash 0_link_samples.bash
While the B training yaml files can handle the list of training sample files (this can be several hundreds of lines), it is convenient to link the training samples, generated from previous section, as an ordered file name as follows
samples/PTB_f48mf24_###.nc #where ### represent the given index, sequentially from 001 to total number of samples.
$ bash 1_run_vbal.bash
This step reads the training samples (under samples directory), calculates the vertical regression coefficients between the desirable variables, then writes out the “unbalanced” training samples (under samplesUnbalanced directory).
It uses the MPAS-JEDI executable, mpasjedi_error_covariance_toolbox.x.
Users will see the following major output files.
VBAL/mpas_vbal_local_000064-0000##.nc : local regression coefficient files for each processor
VBAL/mpas_sampling_local_000064-0000##.nc : local sampling files for each processor
VBAL/mpas_vbal.nc : global regression coefficient file
VBAL/mpas_sampling.nc : global sampling file
samplesUnbalanced/PTB_f48mf24_###.nc : "unbalanced" training samples for {stream_function, unbalanced velocity_potential, unbalanced temperature, specific humidity, unbalanced surface_pressure}
Users can make a plot for diagnosed vertical regression coefficients.
$ cd /glade/derecho/scratch/$USER/mpas_jedi_tutorial/Bflow_global120km/plot_diagnostics
$ ncl plot_vbal_regcoef.ncl
$ display -rotate 270 plot_vbal_regcoef_150.pdf
$ ncl plot_vbal_explained_variance.ncl
$ display -rotate 270 plot_vbal_explained_variance.pdf
$ cd /glade/derecho/scratch/$USER/mpas_jedi_tutorial/Bflow_global120km
$ bash 2a_run_hdiag_var.bash
This step reads the unbalanced training samples (under samplesUnbalanced directory) of {stream_function, unbalanced velocity_potential, unbalanced temperature, specific humidity, unbalanced surface_pressure}, calculates and writes out three B parameters, which are the error standard deviation, horizontal and vertical correlation lengths.
It uses the MPAS-JEDI executable, mpasjedi_error_covariance_toolbox.x.
Users will see the following major output files.
HDIAG_VAR/vargroup1/mpas.cor_rh.nc
HDIAG_VAR/vargroup1/mpas.cor_rv.nc
HDIAG_VAR/vargroup1/mpas.stddev.nc
Users will find the similar script, 2b_run_hdiag_var.bash, for cloud hydrometeor variables.
$ bash 2c_modify_diagnostics.bash
This step performs several manipulations to the diagnosed B parameter files.
- If there exist separate diagnostic files, it is convenient to merge those into a single file.
- If the diagnostics for cloud hydrometer variables contains missing values, it is replaced by “0.0” values for further processing (see etc_modify_missing.bash).
- If needed, the diagnosed B parameters are modified. In this hands-on, the horizontal correlation lengths for stream function and unbalanced velocity_potential is reduced by a factor of 2 (see etc_modify_cor.bash). Also, the standard deviations for all variables are reduced by a factor of 3 (see etc_modify_var.bash).
Users will find the results as follow:
HDIAG_VAR/merge/mpas.cor_rh.nc
HDIAG_VAR/merge/mpas.cor_rv.nc
HDIAG_VAR/merge/mpas.stddev.nc
Users can make a plot for diagnosed (and modified) B parameters.
$ cd /glade/derecho/scratch/$USER/mpas_jedi_tutorial/Bflow_global120km/plot_diagnostics
$ python plot_B_param_cor_rh.py
$ display plot_B_param_cor_rh.png
$ python plot_B_param_cor_rv.py
$ display plot_B_param_cor_rv.png
$ python plot_B_param_stddev.py
$ display plot_B_param_stddev.png
$ cd /glade/derecho/scratch/$USER/mpas_jedi_tutorial/Bflow_global120km
$ bash 3_run_nicas_split.bash
This step calculates the spatial correlation by reading the diagnosed horizontal and vertical correlation length files (i.e., mpas.cor_rh.nc and mpas.cor_rv.nc).
Because the executaion can take significantly long time in some cases, this step is separately done for each variable.
It uses the MPAS-JEDI executable, mpasjedi_error_covariance_toolbox.x.
Users will see the following major output files.
NICAS.split/{varname}/mpas_nicas_local_000064-0000##.nc : local NICAS files for each processor
NICAS.split/{varname}/mpas_nicas.nc : global NICAS files
NICAS.split/{varname}/mpas.dirac_nicas.nc : fields of BUMP-Dirac test
$ bash 4_merge_nicas.bash
This step merged the NICAS files generated in the previous section (for individual variable) into a single file with NCO commands.
Users can find the merged files as follow:
NICAS.split/merge/mpas_nicas_local_000064-0000##.nc : local NICAS files for each processor
NICAS.split/merge/mpas_nicas.nc : global NICAS files
NICAS.split/merge/mpas.dirac_nicas.nc : fields of BUMP-Dirac test
Users can make a plot for the fields of Dirac test.
Note that these Dirac test shows the univariate structure for each individual variable.
$ cd /glade/derecho/scratch/$USER/mpas_jedi_tutorial/Bflow_global120km/plot_diagnostics
$ ncl plot_dirac_nicas_group1.ncl
$ display -rotate 270 plot_dirac_nicas_group1.pdf
Now users have generated all necessary files for multivariate static background error covariance.
The essential files for actual 3DVar or Hybrid-3DEnVar are as follows.
VBAL/mpas_vbal_local_000064-0000##.nc : local regression coefficient files for each processor
VBAL/mpas_sampling_local_000064-0000##.nc : local sampling files for each processor
HDIAG_VAR/merge/mpas.stddev.nc
NICAS.split/merge/mpas_nicas_local_000064-0000##.nc : local NICAS files for each processor
Users can explore the effect of multivariate B with a single observation test.
Here, only two observations, one for temperature and the other for zonal wind, are assimilated with the executable mpasjedi_variational.x.
$ cd /glade/derecho/scratch/$USER/mpas_jedi_tutorial/Bflow_global120km
$ bash 5_SO.bash
Like the previous variational practices, users can find the output (i.e., analysis file) from SO/an.2018-04-15_00.00.00.nc.
Users can make a plot for the analysis increment fields, focused on the two observation locations.
$ cd /glade/derecho/scratch/$USER/mpas_jedi_tutorial/Bflow_global120km/plot_diagnostics
$ ncl plot_SO_T.ncl
$ display plot_SO_T.pdf
$ ncl plot_SO_U.ncl
$ display plot_SO_U.pdf
The observation can be configured directly from the yaml.
For example, users may test different locations ("lats", "lons", "vert coords") or different "obs values" and "obs errors".
observations:
observers:
- obs space:
name: SO_T
simulated variables: [airTemperature]
obsdatain:
engine:
type: GenList
lats: [30.3061]
lons: [130.085]
vert coord type: pressure
vert coords: [78775.95]
dateTimes: [0]
epoch: "seconds since 2018-04-15T00:00:00Z"
obs errors: [0.8]
obs values: [284.5912]
obsdataout:
engine:
type: H5File
obsfile: ./obsout_SO_T.h5
obs operator:
name: VertInterp
#obs filters:
#- filter: GOMsaver
# filename: ./geoval_SO_T.nc4
- obs space:
name: SO_U
simulated variables: [windEastward]
obsdatain:
engine:
type: GenList
lats: [57.7699]
lons: [357.713]
vert coord type: pressure
vert coords: [77693.09]
dateTimes: [0]
epoch: "seconds since 2018-04-15T00:00:00Z"
obs errors: [1.0]
obs values: [0.7250047]
obsdataout:
engine:
type: H5File
obsfile: ./obsout_SO_U.h5
obs operator:
name: VertInterp
#obs filters:
#- filter: GOMsaver
# filename: ./geoval_SO_U.nc4
