Utility Functions
-
utilities.
compute_qgpv_givenvort
(omega, nlat, nlon, kmax, unih, ylat, avort, potential_temp, t0_cn, t0_cs, stat_cn, stat_cs, nlat_s=None, scale_height=7000.0)[source] The function “compute_qgpv_givenvort” computes the quasi-geostrophic potential vorticity based on the absolute vorticity, potential temperature and static stability given in a hemispheric domain.
Please make inquiries and report issues via Github: https://github.com/csyhuang/hn2016_falwa/issues
- Parameters
omega (float, optional) – Rotation rate of the planet.
nlat (int) – Latitudinal dimension of the latitude grid.
nlon (int) – Longitudinal dimension of the longitude grid.
kmax (int) – Vertical dimension of the height grid.
unih (numpy.array) – Numpy array of height in [meters]; dimension = (kmax)
ylat (numpy.array) – Numpy array of latitudes in [degrees]; dimension = (nlat)
avort (numpy.ndarray) – Three-dimension numpy array of absolute vorticity (i.e. relative vorticity + 2*Omega*sin(lat)) in [1/s]; dimension = (kmax x nlat x nlon)
potential_temp (numpy.ndarray) – Three-dimension numpy array of potential temperature in [K]; dimension = (kmax x nlat x nlon)
t0_cn (numpy.array) – Area-weighted average of potential temperature ( ilde{ heta} in HN16) in the Northern hemispheric domain with dimension = (kmax)
t0_cs (numpy.array) – Area-weighted average of potential temperature ( ilde{ heta} in HN16) in the Southern hemispheric domain with dimension = (kmax)
stat_cn (numpy.array) – Static stability (d ilde{ heta}/dz in HN16) in the Northern hemispheric domain with dimension = (kmax)
stat_cs (numpy.array) – Static stability (d ilde{ heta}/dz in HN16) in the Southern hemispheric domain with dimension = (kmax)
scale_height (float) – Scale height of the atmosphere in [m] with default value 7000.
- Returns
QGPV (numpy.ndarray) – Three-dimension numpy array of quasi-geostrophic potential vorticity; dimension = (kmax x nlat x nlon)
dzdiv (numpy.ndarray) – Three-dimension numpy array of the stretching term in QGPV; dimension = (kmax x nlat x nlon)
-
utilities.
curl_2d
(ufield, vfield, clat, dlambda, dphi, planet_radius=6378000.0)[source] Assuming regular latitude and longitude [in degree] grid, compute the curl of velocity on a pressure level in spherical coordinates.
-
utilities.
static_stability
(height, area, theta, s_et=None, n_et=None)[source] The function “static_stability” computes the vertical gradient (z-derivative) of hemispheric-averaged potential temperature, i.e. d ilde{theta}/dz in the def- inition of QGPV in eq.(3) of Huang and Nakamura (2016), by central differencing. At the boundary, the static stability is estimated by forward/backward differen- cing involving two adjacent z-grid points:
- i.e. stat_n[0] = (t0_n[1] - t0_n[0]) / (height[1] - height[0])
stat_n[-1] = (t0_n[-2] - t0_n[-1]) / (height[-2] - height[-1])
Please make inquiries and report issues via Github: https://github.com/csyhuang/hn2016_falwa/issues
- Parameters
height (numpy.array) – Array of z-coordinate [in meters] with dimension = (kmax), equally spaced
area (numpy.ndarray) – Two-dimension numpy array specifying differential areal element of each grid point; dimension = (nlat, nlon).
theta (numpy.ndarray) – Matrix of potential temperature [K] with dimension (kmax,nlat,nlon) or (kmax,nlat)
s_et (int, optional) – Index of the latitude that defines the boundary of the Southern hemispheric domain; initialized as nlat/2 if not input
n_et (int, optional) – Index of the latitude that defines the boundary of the Southern hemispheric domain; initialized as nlat/2 if not input
- Returns
t0_n (numpy.array) – Area-weighted average of potential temperature ( ilde{ heta} in HN16) in the Northern hemispheric domain with dimension = (kmax)
t0_s (numpy.array) – Area-weighted average of potential temperature ( ilde{ heta} in HN16) in the Southern hemispheric domain with dimension = (kmax)
stat_n (numpy.array) – Static stability (d ilde{ heta}/dz in HN16) in the Northern hemispheric domain with dimension = (kmax)
stat_s (numpy.array) – Static stability (d ilde{ heta}/dz in HN16) in the Southern hemispheric domain with dimension = (kmax)
-
utilities.
zonal_convergence
(field, clat, dlambda, planet_radius=6378000.0, tol=1e-05)[source] The function “zonal_convergence” computes the zonal component of the convergence operator of an arbitrary field f(lat, lon), i.e. it computes on the spherical surface: -1/(planet_radius * cos(lat)) * partial d(f(lat, lon))/partial d(lon)
Please make inquiries and report issues via Github: https://github.com/csyhuang/hn2016_falwa/issues
- Parameters
field (numpy.ndarray) – An arbitrary field that one needs to compute zonal divergence with dimension [nlat, nlon]
clat (numpy.array) – Numpy array of cosine latitude; dimension [nlat]
dlambda (float) – Differential element of longitude
planet_radius (float, optional) – Radius of the planet in meters. Default = 6.378e+6 (Earth’s radius)
tol (float, optional) – Tolerance below which clat is considered infinitely small that the corresponding grid points will have returned result = 0 (to avoid division by zero). Default = 1.e-5
- Returns
ans – Zonal convergence of field with the dimension same as field, i.e. [nlat, nlon]
- Return type
numpy.ndarray