atmqty Examples: Relative Vorticity

Description

This web page provides an example of calculating relative vorticity using the AtmQty class. The data used in this example comes from the ECMWF ERA-40 reanalysis. These plots operate as a "unit test" of the function.

A link to the full source code needed and the data file used in these examples is at the bottom of this page. For more information on the vort_rel method see the pydoc description for vort_rel.

Import Statements and Data

All examples below use the following import statements and data:

>>> import atmqty
>>> import cdms
>>> import MV
>>> import Numeric as N

The data is from a netCDF data file of fields from the ECMWF ERA-40 reanalysis for 15 Jan 2001 (0Z) on pressure levels. We read in the data using netCDF utilities from cdms and convert the arrays to Numeric arrays. Vectors lon and lat are longitude and latitude, respectively, in degrees, and lev is pressure levels, in hPa. We transpose the arrays so that the first index cycles through latitude, the second through longitude, and the third through levels. We also reverse the latitude order so that latitude increases with increasing index (in the dataset latitudes start from 90°N and progress south):

>>> f = cdms.open('era40_2001-01-15.nc')

>>> T = MV.filled(f('t', squeeze=1))
>>> u = MV.filled(f('u', squeeze=1))
>>> v = MV.filled(f('v', squeeze=1))
>>> q = MV.filled(f('q', squeeze=1))

>>> T_obj = f['t']
>>> lon = T_obj.getLongitude()[:]
>>> lat = T_obj.getLatitude()[::-1]
>>> lev = T_obj.getLevel()[:].astype(N.Float)

>>> T = N.transpose(T, axes=(1,2,0))[::-1,:,:]
>>> u = N.transpose(u, axes=(1,2,0))[::-1,:,:]
>>> v = N.transpose(v, axes=(1,2,0))[::-1,:,:]
>>> q = N.transpose(q, axes=(1,2,0))[::-1,:,:]

Temperature is in K, winds are in m/s, and specific humidity in kg/kg. The levels are at 250, 300, 400, 500, 600, 700, and 850 hPa. Longitude and latitude values are separated every 2.5 degrees.

Examples

Relative Vorticity Calculated By AtmQty Using ERA-40 Wind Etc. Data

We create an AtmQty object (x) for the above ECMWF domain and fill it with the above atmospheric quantities. Using those quantities, we compute relative vorticity:

>>> x = atmqty.AtmQty( lat, lon, lev, lev_type="pressure" )
>>> x.add_qty({'T':T, 'u':u, 'v':v, 'q':q})
>>> x.vort_rel()

The above calculations (which result in the plot below) are done on a platform with the sphere package installed, so the curl calculations use the spectral method. Note that using the method vort_rel in the fashion above calculates altitude automatically, using default settings. If the default settings (e.g. surface pressure) do not apply, you need to manually call method alt.

The contour intervals are irregular, to bring out the structure of small values of relative vorticity. Units are 1/s.

Relative Vorticity As Provided By ERA-40

The contour intervals are the same as in the AtmQty example above. Note that I do not know what the surface pressure and temperature values were used by ERA-40 to calculate altitude, but they probably were not the default values used by AtmQty. Nonetheless, the result is insensitive to that paramter and the ERA-40 relative vorticity value is nearly identical to the AtmQty value given above.

Full Code and Data File For Examples

• test_vort.py: Note that the code fragments in the examples on this page given above are slightly different from the full code in this file, since the code in test_vort.py generates the plots shown.
• era40_2001-01-15.nc (1.5 Mb): Please see here for the restrictions ECMWF has as to how this data can be used.