Browse Source

Merge branch 'feature/cape2d_1d' into develop

lon0
Bill Ladwig 7 years ago
parent
commit
35a792e62d
  1. 3
      doc/source/tutorial.rst
  2. 4
      doc/source/tutorials/wrf_workshop_2017.rst
  3. 410
      doc/source/tutorials/wrf_workshop_2018.rst
  4. 51
      src/wrf/computation.py
  5. 6
      src/wrf/g_cape.py
  6. 28
      src/wrf/metadecorators.py
  7. 45
      src/wrf/specialdec.py
  8. BIN
      test/ci_tests/ci_result_file.nc
  9. 65
      test/comp_utest.py

3
doc/source/tutorial.rst

@ -12,7 +12,7 @@ Upcoming Tutorials @@ -12,7 +12,7 @@ Upcoming Tutorials
.. toctree::
:maxdepth: 1
tutorials/tutorial_03_2018.rst
tutorials/wrf_workshop_2018.rst
Past Tutorials
@ -22,4 +22,5 @@ Past Tutorials @@ -22,4 +22,5 @@ Past Tutorials
:maxdepth: 1
tutorials/wrf_workshop_2017.rst
tutorials/tutorial_03_2018.rst

4
doc/source/tutorials/wrf_workshop_2017.rst

@ -1,5 +1,5 @@ @@ -1,5 +1,5 @@
WRF Workshop 2017
=====================
WRF Users' Workshop 2017
==========================
Welcome wrf-python tutorial attendees!

410
doc/source/tutorials/wrf_workshop_2018.rst

@ -0,0 +1,410 @@ @@ -0,0 +1,410 @@
WRF Users' Workshop 2018
=========================
Welcome WRF-Python tutorial attendees!
The instructions below should be completed prior to arriving at the tutorial.
There will not be enough time to do this during the tutorial.
Prerequisites
---------------
This tutorial assumes that you have basic knowledge of how to type commands
in to a command terminal using your preferred operating system. You
should know some basic directory commands like *cd*, *mkdir*, *cp*, *mv*.
Regarding Python, to understand the examples in this tutorial, you
should have some experience with Python basics. Below is a list of some
Python concepts that you will see in the examples, but don't worry if you aren't
familiar with everything.
- Opening a Python interpreter and entering commands.
- Importing packages via the import statement.
- Familiarity with some of the basic Python types: str, list, tuple, dict, bool, float, int, None.
- Creating a list, tuple, or dict with "[ ]", "( )", "{ }" syntax (e.g. my_list = [1,2,3,4,5]).
- Accessing dict/list/tuple items with the "x[ ]" syntax (e.g. my_list_item = my_list[0]).
- Slicing str/list/tuple with the ":" syntax (e.g. my_slice = my_list[1:3]).
- Using object methods and attributes with the "x.y" syntax (e.g. my_list.append(6)).
- Calling functions (e.g. result = some_function(x, y))
- Familiarity with numpy would be helpful, as only a very brief introduction
is provided.
- Familiarity with matplotlib would be helpful, as only a very brief
introduction is provided.
If you are completely new to Python, that shouldn't be a problem, since
most of the examples consist of basic container types and function calls. It
would be helpful to look at some introductory material before arriving at the
tutorial. If you've programmed before, picking up Python isn't too difficult.
Here are some links:
https://www.learnpython.org/
https://developers.google.com/edu/python/
Step 1: Open a Command Terminal
--------------------------------
To begin, you will first need to know how to open a command line terminal for
your operating system.
For Windows:
.. code-block:: none
WINDOWS + r
type cmd in the run window
For Mac:
.. code-block:: none
Finder -> Applications -> Utilities -> Terminal
For Linux:
.. code-block:: none
Try one of the following:
CTRL + ALT + T
CTRL + ALT + F2
Step 2: Download Miniconda
----------------------------
For this tutorial, you will need to download and install Miniconda. We are
going to use Python 2.7, but it will also work with Python 3.5+.
Please use the appropriate link below to download Miniconda for your operating
system.
.. note::
64-bit OS recommended
`Win64 <https://repo.continuum.io/miniconda/Miniconda2-latest-Windows-x86_64.exe>`_
`Mac <https://repo.continuum.io/miniconda/Miniconda2-latest-MacOSX-x86_64.sh>`_
`Linux <https://repo.continuum.io/miniconda/Miniconda2-latest-Linux-x86_64.sh>`_
For more information, see: https://conda.io/miniconda.html
.. note::
**What is Miniconda?**
If you have used the Anaconda distribution for Python before, then you will be
familiar with Miniconda. The Anaconda Python distribution includes numerous
scientific packages out of the box, which can be difficult for users to build and
install. More importantly, Anaconda includes the conda package manager.
The conda package manager is a utility (similar to yum or apt-get) that installs
packages from a repository of pre-compiled Python packages. These repositories
are called channels. Conda makes it easy for Python users to install and
uninstall packages, and also can be used to create isolated Python environments
(more on that later).
Miniconda is a bare bones implementation of Anaconda and only includes the
conda package manager. Since we are going to use the conda-forge channel to
install our scientific packages, Miniconda avoids any complications between
packages provided by Anaconda and conda-forge.
Step 3: Install Miniconda
----------------------------
Windows:
1. Browse to the directory where you downloaded Miniconda2-latest-Windows-x86_64.exe.
2. Double click on Miniconda2-latest-Windows-x86_64.exe.
3. Follow the instructions.
4. Restart your command terminal.
Mac and Linux:
For Mac and Linux, the installer is a bash script.
1. Using a terminal, you need to execute the bash shell script that you downloaded by
doing::
bash /path/to/Miniconda2-latest-MacOSX-x86_64.sh [Mac]
bash /path/to/Miniconda2-latest-Linux-x86_64.sh [Linux]
2. Follow the instructions.
3. At the end of the installation, it will ask if you want to add the
miniconda2 path to your bash environment. If you are unsure what to do,
you should say "yes". If you say "no", we're going to assume you know
what you are doing.
If you said "yes", then once you restart your shell, the miniconda2 Python
will be found instead of the system Python when you type the "python"
command. If you want to undo this later, then you can edit
either ~/.bash_profile or ~/.bashrc (depending on OS used) and
comment out the line that looks similar to::
# added by Miniconda2 4.1.11 installer
export PATH="/path/to/miniconda2/bin:$PATH"
4. Restart your command terminal.
5. [Linux and Mac Users Only] Miniconda only works with bash. If bash is
not your default shell, then you need to activate the bash shell by typing
the following in to your command terminal::
bash
6. Verify that your system is using the correct Python interpreter by typing
the following in to your command terminal::
which python
You should see the path to your miniconda installation. If not, see the
note below.
.. note::
If you have already installed another Python distribution, like Enthought
Canopy, you will need to comment out any PATH entries for that distribution
in your .bashrc or .bash_profile. Otherwise, your shell environment may
pick to wrong Python installation.
If bash is not your default shell type, and the PATH variable has been
set in .bash_profile by the miniconda installer, try executing
"bash -l" instead of the "bash" command in step 5.
Step 4: Set Up the Conda Environment
--------------------------------------
If you are new to the conda package manager, one of the nice features of conda
is that you can create isolated Python environments that prevent package
incompatibilities. This is similar to the *virtualenv* package that some
Python users may be familiar with. However, conda is not compatible with
virtualenv, so only use conda environments when working with conda.
The name of our conda environment for this tutorial is: **tutorial_2018**.
Follow the instructions below to create the tutorial_2018 environment.
1. Open a command terminal if you haven't done so.
2. [Linux and Mac Users Only] The conda package manager only works with bash,
so if bash is not your current shell, type::
bash
3. Add the conda-forge channel to your conda package manager.
Type or copy the command below in to your command terminal. You should
run this command even if you have already done it in the past.
This will ensure that conda-forge is set as the highest priority channel.
::
conda config --add channels conda-forge
.. note::
Conda-forge is a community driven collection of packages that are
continually tested to ensure compatibility. We highly recommend using
conda-forge when working with conda. See https://conda-forge.github.io/
for more details on this excellent project.
4. Create the conda environment for the tutorial.
Type or copy this command in to your command terminal::
conda create -n tutorial_2018 python=2.7 matplotlib cartopy netcdf4 jupyter git ffmpeg wrf-python
Type "y" when prompted. It will take several minutes to install everything.
This command creates an isolated Python environment named *tutorial_2018*, and installs
the python interpreter, matplotlib, cartopy, netcdf4, jupyter, git, ffmpeg, and wrf-python
packages.
.. note::
When the installation completes, your command terminal might post a message similar to:
.. code-block:: none
If this is your first install of dbus, automatically load on login with:
mkdir -p ~/Library/LaunchAgents
cp /path/to/miniconda2/envs/tutorial_test/org.freedesktop.dbus-session.plist ~/Library/LaunchAgents/
launchctl load -w ~/Library/LaunchAgents/org.freedesktop.dbus-session.plist
This is indicating that the dbus package can be set up to automatically load on login. You
can either ignore this message or type in the commands as indicated on your command terminal.
The tutorial should work fine in either case.
5. Activate the conda environment.
To activate the tutorial_2018 Python environment, type the following
in to the command terminal:
For Linux and Mac (using bash)::
source activate tutorial_2018
For Windows::
activate tutorial_2018
You should see (tutorial_2018) on your command prompt.
To deactivate your conda environment, type the following in to the
command terminal:
For Linux and Mac::
source deactivate
For Windows::
deactivate tutorial_2018
Step 5: Download the Student Workbook
---------------------------------------
The student workbook for the tutorial is available on GitHub. The tutorial_2018
conda environment includes the git application needed to download the repository.
These instructions download the tutorial in to your home directory. If you want
to place the tutorial in to another directory, we're going to assume you know
how to do this yourself.
To download the student workbook, follow these instructions:
1. Activate the tutorial_2018 conda environment following the instructions
in the previous step (*source activate tutorial_2018* or
*activate tutorial_2018*).
2. Change your working directory to the home directory by typing the
following command in to the command terminal:
For Linux and Mac::
cd ~
For Windows::
cd %HOMEPATH%
3. Download the git repository for the tutorial by typing the following
in to the command terminal::
git clone https://github.com/NCAR/wrf_python_tutorial.git
4. There may be additional changes to the tutorial after you have downloaded
it. To pull down the latest changes, type the following in to the
command terminal:
For Linux and Mac::
source activate tutorial_2018
cd ~/wrf_python_tutorial/wrf_workshop_2018
git pull
For Windows::
activate tutorial_2018
cd %HOMEPATH%\wrf_python_tutorial\wrf_workshop_2018
git pull
.. note::
If you try the "git pull" command and it returns an error indicating
that you have made changes to the workbook, this is probably because
you ran the workbook and it contains the cell output. To fix this,
first do a checkout of the workbook, then do the pull.
.. code-block:: none
git checkout -- .
git pull
Step 6: Verify Your Environment
----------------------------------
Verifying that your environment is correct involves importing a few
packages and checking for errors (you may see some warnings for matplotlib
or xarray, but you can safely ignore these).
1. Activate the tutorial_2018 conda environment if it isn't already active
(see instructions above).
2. Open a python terminal by typing the following in to the command
terminal::
python
3. Now type the following in to the Python interpreter::
>>> import netCDF4
>>> import matplotlib
>>> import xarray
>>> import wrf
4. You can exit the Python interpreter using **CTRL + D**
Step 7: Obtain WRF Output Files
----------------------------------
For this tutorial, we strongly recommend that you use your own WRF output files.
The tutorial includes an easy way to point to your own data files. The WRF
output files should all be from the same WRF run and use the same domain.
If your files are located on another system (e.g. yellowstone), then copy 2 or
3 of these files to your local computer prior to the tutorial.
If you do not have any of your own WRF output files, then you can download the
instructor data files from a link that should have been provided to you in an
email prior to the tutorial.
If you are using the link provided in the email for your data, you can follow
the instructions below to place your data in the default location for your
workbook.
1. The link in the email should take you to a location on an Amazon cloud
drive.
2. If you hover your mouse over the wrf_tutorial_data.zip file, you'll see
an empty check box appear next to the file name. Click this check
box.
3. At the bottom of the screen, you'll see a Download button next to a
cloud icon. Click this button to start the download.
4. The download was most likely placed in to your ~/Downloads folder
[%HOMEPATH%\\Downloads for Windows]. Using your preferred method of choice
for unzipping files, unzip this file in to your home directory. Your data
should now be in ~/wrf_tutorial_data
[%HOMEPATH%\\wrf_tutorial_data for Windows].
5. Verify that you have three WRF output files in that directory.
Getting Help
----------------
If you experience problems during this installation, please send a question
to the :ref:`google-group` support mailing list.
We look forward to seeing you at the tutorial!

51
src/wrf/computation.py

@ -714,10 +714,10 @@ def smooth2d(field, passes, meta=True): @@ -714,10 +714,10 @@ def smooth2d(field, passes, meta=True):
@set_cape_alg_metadata(is2d=True, copyarg="pres_hpa")
def cape_2d(pres_hpa, tkel, qv, height, terrain, psfc_hpa, ter_follow,
missing=default_fill(np.float64), meta=True):
"""Return the two-dimensional CAPE, CIN, LCL, and LFC.
"""Return the two-dimensional MCAPE, MCIN, LCL, and LFC.
This function calculates the maximum convective available potential
energy (CAPE), maximum convective inhibition (CIN),
energy (MCAPE), maximum convective inhibition (MCIN),
lifted condensation level (LCL), and level of free convection (LFC). This
function uses the RIP [Read/Interpolate/plot] code to calculate
potential energy (CAPE) and convective inhibition
@ -725,18 +725,28 @@ def cape_2d(pres_hpa, tkel, qv, height, terrain, psfc_hpa, ter_follow, @@ -725,18 +725,28 @@ def cape_2d(pres_hpa, tkel, qv, height, terrain, psfc_hpa, ter_follow,
in the column (i.e. something akin to Colman's MCAPE). CAPE is defined as
the accumulated buoyant energy from the level of free convection (LFC) to
the equilibrium level (EL). CIN is defined as the accumulated negative
buoyant energy from the parcel starting point to the LFC. The word 'parcel'
here refers to a 500 meter deep parcel, with actual temperature and
moisture averaged over that depth.
buoyant energy from the parcel starting point to the LFC.
The cape_2d algorithm works by first finding the maximum theta-e height
level in the lowest 3000 m. A parcel with a depth of 500 m is then
calculated and centered over this maximum theta-e height level. The
parcel's moisture and temperature characteristics are calculated by
averaging over the depth of this 500 m parcel. This 'maximum' parcel
is then used to compute MCAPE, MCIN, LCL and LFC.
The leftmost dimension of the returned array represents four different
quantities:
- return_val[0,...] will contain CAPE [J kg-1]
- return_val[1,...] will contain CIN [J kg-1]
- return_val[0,...] will contain MCAPE [J kg-1]
- return_val[1,...] will contain MCIN [J kg-1]
- return_val[2,...] will contain LCL [m]
- return_val[3,...] will contain LFC [m]
This function also supports computing MCAPE along a single vertical
column. In this mode, the *pres_hpa*, *tkel*, *qv* and *height* arguments
must be one-dimensional vertical columns, and the *terrain* and
*psfc_hpa* arguments must be scalar values
(:obj:`float`, :class:`numpy.float32` or :class:`numpy.float64`).
This is the raw computational algorithm and does not extract any variables
from WRF output files. Use :meth:`wrf.getvar` to both extract and compute
@ -749,6 +759,9 @@ def cape_2d(pres_hpa, tkel, qv, height, terrain, psfc_hpa, ter_follow, @@ -749,6 +759,9 @@ def cape_2d(pres_hpa, tkel, qv, height, terrain, psfc_hpa, ter_follow,
least three dimensions. The rightmost dimensions can be
top_bottom x south_north x west_east or bottom_top x south_north x
west_east.
When operating on only a single column of values, the vertical
column can be bottom_top or top_bottom. In this case, *terrain*
and *psfc_hpa* must be scalars.
Note:
@ -774,12 +787,16 @@ def cape_2d(pres_hpa, tkel, qv, height, terrain, psfc_hpa, ter_follow, @@ -774,12 +787,16 @@ def cape_2d(pres_hpa, tkel, qv, height, terrain, psfc_hpa, ter_follow,
terrain (:class:`xarray.DataArray` or :class:`numpy.ndarray`):
Terrain height in [m]. This is at least a two-dimensional array
with the same dimensionality as *pres_hpa*, excluding the vertical
(bottom_top/top_bottom) dimension.
(bottom_top/top_bottom) dimension. When operating on a single
vertical column, this argument must be a scalar (:obj:`float`,
:class:`numpy.float32`, or :class:`numpy.float64`).
psfc_hpa (:class:`xarray.DataArray` or :class:`numpy.ndarray`):
The surface pressure in [hPa]. This is at least a two-dimensional
array with the same dimensionality as *pres_hpa*, excluding the
vertical (bottom_top/top_bottom) dimension.
vertical (bottom_top/top_bottom) dimension. When operating on a
singlevertical column, this argument must be a scalar
(:obj:`float`, :class:`numpy.float32`, or :class:`numpy.float64`).
Note:
@ -1111,12 +1128,16 @@ def ctt(pres_hpa, tkel, qv, qcld, height, terrain, qice=None, @@ -1111,12 +1128,16 @@ def ctt(pres_hpa, tkel, qv, qcld, height, terrain, qice=None,
True. Default is
:data:`wrf.default_fill(numpy.float64)`.
opt_thresh (:obj:`float`, optional): The required amount of optical
depth (looking from top down) required to trigger a cloud top
temperature calculation. Values less than this threshold will be
treated as no cloud areas. For thin cirrus, a value of .1 might be
appropriate. For large CB clouds, 1000.0 might be appropriate.
Default is 1.0.
opt_thresh (:obj:`float`, optional): The amount of optical
depth (integrated from top down) required to trigger a cloud top
temperature calculation. The cloud top temperature is calculated at
the vertical level where this threshold is met. Vertical columns
with less than this threshold will be treated as cloud free areas.
In general, the larger the value is for this
threshold, the lower the altitude will be for the cloud top
temperature calculation, and therefore higher cloud top
temperature values. Default is 1.0, which should be sufficient for
most users.
meta (:obj:`bool`): Set to False to disable metadata and return
:class:`numpy.ndarray` instead of

6
src/wrf/g_cape.py

@ -13,13 +13,13 @@ from .metadecorators import set_cape_metadata @@ -13,13 +13,13 @@ from .metadecorators import set_cape_metadata
@set_cape_metadata(is2d=True)
def get_2dcape(wrfin, timeidx=0, method="cat", squeeze=True, cache=None,
meta=True, _key=None, missing=default_fill(np.float64)):
"""Return the 2d fields of CAPE, CIN, LCL, and LFC.
"""Return the 2d fields of MCAPE, MCIN, LCL, and LFC.
The leftmost dimension of the returned array represents four different
quantities:
- return_val[0,...] will contain CAPE [J kg-1]
- return_val[1,...] will contain CIN [J kg-1]
- return_val[0,...] will contain MCAPE [J kg-1]
- return_val[1,...] will contain MCIN [J kg-1]
- return_val[2,...] will contain LCL [m]
- return_val[3,...] will contain LFC [m]

28
src/wrf/metadecorators.py

@ -1985,7 +1985,7 @@ def set_cape_alg_metadata(is2d, copyarg="pres_hpa"): @@ -1985,7 +1985,7 @@ def set_cape_alg_metadata(is2d, copyarg="pres_hpa"):
p = argvals[copyarg]
missing = argvals["missing"]
# Note: cape_3d supports using only a single column of data
# Note: 2D/3D cape supports using only a single column of data
is1d = p.ndim == 1
# Need to squeeze the right dimensions for 1D cape
@ -1997,31 +1997,34 @@ def set_cape_alg_metadata(is2d, copyarg="pres_hpa"): @@ -1997,31 +1997,34 @@ def set_cape_alg_metadata(is2d, copyarg="pres_hpa"):
outattrs = OrderedDict()
if is2d:
outname = "cape_2d"
if is1d:
outname = "cape_2d"
else:
outname = "cape_2d_column"
outattrs["description"] = "mcape ; mcin ; lcl ; lfc"
outattrs["units"] = "J kg-1 ; J kg-1 ; m ; m"
outattrs["MemoryOrder"] = "XY"
outattrs["MemoryOrder"] = ""
else:
if not is1d:
outname = "cape_3d"
outattrs["description"] = "cape; cin"
outattrs["units"] = "J kg-1 ; J kg-1"
outattrs["MemoryOrder"] = "XYZ"
else:
outname = "cape_column"
outattrs["description"] = "cape; cin"
outattrs["units"] = "J kg-1 ; J kg-1"
outname = "cape_3d_column"
outattrs["MemoryOrder"] = "Z"
outattrs["description"] = "cape; cin"
outattrs["units"] = "J kg-1 ; J kg-1"
if isinstance(p, DataArray):
if is2d:
# Right dims
outdims[-2:] = p.dims[-2:]
# Left dims
outdims[1:-2] = p.dims[0:-3]
if not is1d:
# Right dims
outdims[-2:] = p.dims[-2:]
# Left dims
outdims[1:-2] = p.dims[0:-3]
else:
if not is1d:
# Right dims
@ -2031,6 +2034,7 @@ def set_cape_alg_metadata(is2d, copyarg="pres_hpa"): @@ -2031,6 +2034,7 @@ def set_cape_alg_metadata(is2d, copyarg="pres_hpa"):
else:
outdims[1] = p.dims[0]
outcoords = {}
# Left-most is always cape_cin or cape_cin_lcl_lfc
if is2d:

45
src/wrf/specialdec.py

@ -229,7 +229,14 @@ def cape_left_iter(alg_dtype=np.float64): @@ -229,7 +229,14 @@ def cape_left_iter(alg_dtype=np.float64):
ter_follow = args[8]
is2d = i3dflag == 0
is1d = np.isscalar(sfp)
# Note: This should still work with DataArrays
is1d = np.isscalar(sfp) or np.size(sfp) == 1
# Make sure sfp and terrain are regular floats for 1D case
# This should also work with DataArrays
if is1d:
ter = float(ter)
sfp = float(sfp)
orig_dtype = p_hpa.dtype
@ -542,39 +549,27 @@ def check_cape_args(): @@ -542,39 +549,27 @@ def check_cape_args():
ter_follow = args[8]
is2d = False if i3dflag != 0 else True
is1d = ((np.isscalar(sfp) or np.size(sfp) == 1) or
(np.isscalar(ter) or np.size(ter) == 1))
if not (p_hpa.shape == tk.shape == qv.shape == ht.shape):
raise ValueError("arguments 0, 1, 2, 3 must be the same shape")
# 2D CAPE does not allow for scalars
if is2d:
if np.isscalar(ter) or np.isscalar(sfp):
raise ValueError("arguments 4 and 5 cannot be scalars with "
"cape_2d routine")
if not is1d:
if ter.ndim != p_hpa.ndim-1 or sfp.ndim != p_hpa.ndim-1:
raise ValueError("arguments 4 and 5 must have "
"{} dimensions".format(p_hpa.ndim-1))
# 3D cape is allowed to be just a vertical column with scalars
# for terrain and psfc_hpa.
else:
if np.isscalar(ter) and np.isscalar(sfp):
if p_hpa.ndim != 1:
raise ValueError("arguments 0-3 "
"must be 1-dimensional when "
"arguments 4 and 5 are scalars")
if is2d:
raise ValueError("argument 7 must be 0 when using 1D data")
else:
if ((np.isscalar(ter) and not np.isscalar(sfp)) or
(not np.isscalar(ter) and np.isscalar(sfp))):
raise ValueError("arguments 4 and 5 must both be scalars")
else:
if ter.ndim != p_hpa.ndim-1 or sfp.ndim != p_hpa.ndim-1:
raise ValueError("arguments 4 and 5 "
"must have {} dimensions".format(
p_hpa.ndim-1))
if np.size(ter) != np.size(sfp):
raise ValueError("arguments 4 and 5 must both be scalars or "
"both be arrays")
# Only need to test p_hpa since we assured args 0-3 have same ndim
if p_hpa.ndim != 1:
raise ValueError("arguments 0-3 "
"must be 1-dimensional when "
"arguments 4 and 5 are scalars")
return wrapped(*args, **kwargs)

BIN
test/ci_tests/ci_result_file.nc

Binary file not shown.

65
test/comp_utest.py

@ -599,7 +599,69 @@ def test_cape3d_1d(wrfnc): @@ -599,7 +599,69 @@ def test_cape3d_1d(wrfnc):
return func
def test_cape2d_1d(wrfnc):
def func(self):
varnames = ("T", "P", "PB", "QVAPOR", "PH", "PHB", "HGT", "PSFC")
ncvars = extract_vars(wrfnc, 0, varnames, method="cat", squeeze=True,
cache=None, meta=True)
t = ncvars["T"]
p = ncvars["P"]
pb = ncvars["PB"]
qv = ncvars["QVAPOR"]
ph = ncvars["PH"]
phb = ncvars["PHB"]
ter = ncvars["HGT"]
psfc = ncvars["PSFC"]
col_idxs = (slice(None), t.shape[-2]//2, t.shape[-1]//2)
t = t[col_idxs]
p = p[col_idxs]
pb = pb[col_idxs]
qv = qv[col_idxs]
ph = ph[col_idxs]
phb = phb[col_idxs]
ter = float(ter[col_idxs[1:]])
psfc = float(psfc[col_idxs[1:]])
full_t = t + Constants.T_BASE
full_p = p + pb
tkel = tk(full_p, full_t, meta=False)
geopt = ph + phb
geopt_unstag = destagger(geopt, -1)
z = geopt_unstag/Constants.G
# Convert pressure to hPa
p_hpa = ConversionFactors.PA_TO_HPA * full_p
psfc_hpa = ConversionFactors.PA_TO_HPA * psfc
i3dflag = 0
ter_follow = 1
result = cape_2d(p_hpa, tkel, qv, z, ter, psfc_hpa, ter_follow)
print ("RESULT", result)
ref = getvar(wrfnc, "cape_2d")
ref = ref[(slice(None),) + col_idxs[1:]]
print ("REF", ref)
nt.assert_allclose(to_np(result), to_np(ref))
return func
if __name__ == "__main__":
from wrf import (omp_set_num_threads, omp_set_schedule, omp_get_schedule,
omp_set_dynamic, omp_get_num_procs, OMP_SCHED_STATIC)
omp_set_num_threads(omp_get_num_procs()-1)
omp_set_schedule(OMP_SCHED_STATIC, 0)
omp_set_dynamic(False)
varnames = ["avo", "eth", "cape_2d", "cape_3d", "ctt", "dbz", "mdbz",
"geopt", "helicity", "lat", "lon", "omg", "p", "pressure",
"pvo", "pw", "rh2", "rh", "slp", "ter", "td2", "td", "tc",
@ -639,6 +701,9 @@ if __name__ == "__main__": @@ -639,6 +701,9 @@ if __name__ == "__main__":
func = test_cape3d_1d(wrfnc)
setattr(WRFVarsTest, "test_cape3d_1d", func)
func = test_cape2d_1d(wrfnc)
setattr(WRFVarsTest, "test_cape2d_1d", func)
ut.main()

Loading…
Cancel
Save