Source code for metrics

# Written by Ori Adam Mar.21.2017
# Edited by Alison Ming Jul.4.2017
from __future__ import division
import numpy as np
from functions import *

[docs]def TropD_Metric_EDJ(U, lat, lev=np.array([1]), method='peak', n=0, n_fit=1): '''TropD Eddy Driven Jet (EDJ) metric Latitude of maximum of the zonal wind at the level closest to the 850 hPa level Args: U (lat,lev) or U (lat,): Zonal mean zonal wind. Also takes surface wind lat : latitude vector lev: vertical level vector in hPa units method (str, optional): 'peak' (default) | 'max' | 'fit' peak (Default): Latitude of the maximum of the zonal wind at the level closest to the 850 hPa level (smoothing parameter n=30) max: Latitude of the maximum of the zonal wind at the level closest to the 850 hPa level (smoothing parameter n=6) fit: Latitude of the maximum of the zonal wind at the level closest to the 850 hPa level using a quadratic polynomial fit of data from gridpoints surrounding the gridpoint of the maximum n (int, optional): If n is not set (0), n=6 (default) is used in TropD_Calculate_MaxLat. Rank of moment used to calculate the position of max value. n = 1,2,4,6,8,... Returns: tuple: PhiSH (ndarray), PhiNH (ndarray) Latitude of EDJ in SH and NH ''' try: assert (not hasattr(n, "__len__") and n >= 0) except AssertionError: print('TropD_Metric_EDJ: ERROR : the smoothing parameter n must be >= 0') try: assert(method in ['max','peak']) except AssertionError: print('TropD_Metric_EDJ: ERROR : unrecognized method ', method) eq_boundary = 15 polar_boundary = 70 if len(lev) > 1: u = U[:,find_nearest(lev, 850)] else: u = np.copy(U) if method == 'max': if n: PhiNH = TropD_Calculate_MaxLat(u[(lat > eq_boundary) & (lat < polar_boundary)],\ lat[(lat > eq_boundary) & (lat < polar_boundary)], n=n) PhiSH = TropD_Calculate_MaxLat(u[(lat > -polar_boundary) & (lat < -eq_boundary)],\ lat[(lat > -polar_boundary) & (lat < -eq_boundary)], n=n) else: #Default value of n=6 is used PhiNH = TropD_Calculate_MaxLat(u[(lat > eq_boundary) & (lat < polar_boundary)],\ lat[(lat > eq_boundary) & (lat < polar_boundary)]) PhiSH = TropD_Calculate_MaxLat(u[(lat > -polar_boundary) & (lat < -eq_boundary)],\ lat[(lat > -polar_boundary) & (lat < -eq_boundary)]) elif method == 'peak': if n: PhiNH = TropD_Calculate_MaxLat(u[(lat > eq_boundary) & (lat < polar_boundary)],\ lat[(lat > eq_boundary) & (lat < polar_boundary)], n=n) PhiSH = TropD_Calculate_MaxLat(u[(lat > -polar_boundary) & (lat < -eq_boundary)],\ lat[(lat > -polar_boundary) & (lat < -eq_boundary)], n=n) else: PhiNH = TropD_Calculate_MaxLat(u[(lat > eq_boundary) & (lat < polar_boundary)],\ lat[(lat > eq_boundary) & (lat < polar_boundary)],n=30) PhiSH = TropD_Calculate_MaxLat(u[(lat > -polar_boundary) & (lat < -eq_boundary)],\ lat[(lat > -polar_boundary) & (lat < -eq_boundary)],n=30) elif method == 'fit': Uh = u[(lat > eq_boundary) & (lat < polar_boundary)] Lat = lat[(lat > eq_boundary) & (lat < polar_boundary)] m = np.nanmax(Uh) Im = np.nanargmax(Uh) if (Im == 0 or Im == len(Uh)-1): PhiNH = Lat[Im] elif (n_fit > Im or n_fit > len(Uh)-Im+1): N = np.min(Im, len(Uh)-Im+1) p = np.polyfit(Lat[Im-N:Im+N+1], Uh[Im-N:Im+N+1],2) PhiNH = -p[1]/(2*p[0]) else: p = np.polyfit(Lat[Im-n_fit:Im+n_fit+1], Uh[Im-n_fit:Im+n_fit+1],2) PhiNH = -p[1]/(2*p[0]) Uh = u[(lat > -polar_boundary) & (lat < -eq_boundary)] Lat = lat[(lat > -polar_boundary) & (lat < -eq_boundary)] m = np.nanmax(Uh) Im = np.nanargmax(Uh) if (Im == 0 or Im == len(Uh)-1): PhiSH = Lat[Im] elif (n_fit > Im or n_fit > len(Uh)-Im+1): N = np.min(Im, len(Uh)-Im+1) p = np.polyfit(Lat[Im-N:Im+N+1], Uh[Im-N:Im+N+1],2) PhiSH = -p[1]/(2*p[0]) else: p = np.polyfit(Lat[Im-n_fit:Im+n_fit+1], Uh[Im-n_fit:Im+n_fit+1],2) PhiSH = -p[1]/(2*p[0]) else: print('TropD_Metric_EDJ: ERROR: unrecognized method ', method) return PhiSH, PhiNH
[docs]def TropD_Metric_OLR(olr, lat, method='250W', Cutoff=50, n=int(6)): """TropD Outgoing Longwave Radiation (OLR) metric Args: olr(lat,): zonal mean TOA olr (positive) lat: equally spaced latitude column vector method (str, optional): '250W'(Default): the first latitude poleward of the tropical OLR maximum in each hemisphere where OLR crosses 250W/m^2 '20W': the first latitude poleward of the tropical OLR maximum in each hemisphere where OLR crosses the tropical OLR max minus 20W/m^2 'cutoff': the first latitude poleward of the tropical OLR maximum in each hemisphere where OLR crosses a specified cutoff value '10Perc': the first latitude poleward of the tropical OLR maximum in each hemisphere where OLR is 10# smaller than the tropical OLR maximum 'max': the latitude of maximum of tropical olr in each hemisphere with the smoothing paramerer n=6 in TropD_Calculate_MaxLat 'peak': the latitude of maximum of tropical olr in each hemisphere with the smoothing parameter n=30 in TropD_Calculate_MaxLat Cutoff (float, optional): Scalar. For the method 'cutoff', Cutoff specifies the OLR cutoff value. n (int, optional): For the 'max' method, n is the smoothing parameter in TropD_Calculate_MaxLat Returns: tuple: PhiSH (ndarray), PhiNH (ndarray) Latitude of near equator OLR threshold crossing in SH and NH """ try: assert(isinstance(n, int)) except AssertionError: print('TropD_Metric_OLR: ERROR: the smoothing parameter n must be an integer') try: assert(n>=1) except AssertionError: print('TropD_Metric_OLR: ERROR: the smoothing parameter n must be >= 1') # make latitude vector monotonically increasing if lat[-1] < lat[0]: olr = np.flip(olr,0) lat = np.flip(lat,0) eq_boundary = 5 subpolar_boundary = 40 polar_boundary = 60 # NH olr_max_lat_NH = TropD_Calculate_MaxLat(olr[(lat > eq_boundary) & (lat < subpolar_boundary)],\ lat[(lat > eq_boundary) & (lat < subpolar_boundary)]) olr_max_NH = max(olr[(lat > eq_boundary) & (lat < subpolar_boundary)]) # SH olr_max_lat_SH = TropD_Calculate_MaxLat(olr[(lat > -subpolar_boundary) & (lat < -eq_boundary)],\ lat[(lat > -subpolar_boundary) & (lat < -eq_boundary)]) olr_max_SH = max(olr[(lat > -subpolar_boundary) & (lat < -eq_boundary)]) if method == '20W': PhiNH = TropD_Calculate_ZeroCrossing(olr[(lat > olr_max_lat_NH) & (lat < polar_boundary)] - olr_max_NH + 20,\ lat[(lat > olr_max_lat_NH) & (lat < polar_boundary)]) PhiSH = TropD_Calculate_ZeroCrossing(np.flip(olr[(lat < olr_max_lat_SH) & \ (lat > -polar_boundary)],0) - olr_max_SH + 20,\ np.flip(lat[(lat < olr_max_lat_SH) & (lat > -polar_boundary)],0)) elif method == '250W': PhiNH = TropD_Calculate_ZeroCrossing(olr[(lat > olr_max_lat_NH) & (lat < polar_boundary)] - 250,\ lat[(lat > olr_max_lat_NH) & (lat < polar_boundary)]) PhiSH = TropD_Calculate_ZeroCrossing(np.flip(olr[(lat < olr_max_lat_SH) & (lat > -polar_boundary)],0) - 250,\ np.flip(lat[(lat < olr_max_lat_SH) & (lat > -polar_boundary)],0)) elif method == 'cutoff': PhiNH = TropD_Calculate_ZeroCrossing(olr[(lat > olr_max_lat_NH) & (lat < polar_boundary)] - Cutoff,\ lat[(lat > olr_max_lat_NH) & (lat < polar_boundary)]) PhiSH = TropD_Calculate_ZeroCrossing(np.flip(olr[(lat < olr_max_lat_SH) & (lat > -polar_boundary)],0) - Cutoff,\ np.flip(lat[(lat < olr_max_lat_SH) & (lat > -polar_boundary)],0)) elif method == '10Perc': PhiNH = TropD_Calculate_ZeroCrossing(olr[(lat > olr_max_lat_NH) & (lat < polar_boundary)] / olr_max_NH - 0.9,\ lat[(lat > olr_max_lat_NH) & (lat < polar_boundary)]) PhiSH = TropD_Calculate_ZeroCrossing(np.flip(olr[(lat < olr_max_lat_SH) & (lat > -polar_boundary)],0) \ / olr_max_SH - 0.9, np.flip(lat[(lat < olr_max_lat_SH) & (lat > -polar_boundary)],0)) elif method == 'max': if Cutoff_is_set: PhiNH = TropD_Calculate_MaxLat(olr[(lat > eq_boundary) & (lat < subpolar_boundary)],\ lat[(lat > eq_boundary) & (lat < subpolar_boundary)], n=n) PhiSH = TropD_Calculate_MaxLat(olr[(lat > -subpolar_boundary) & (lat < -eq_boundary)],\ lat[(lat > -subpolar_boundary) & (lat < -eq_boundary)], n=n) else: PhiNH = np.copy(olr_max_lat_NH) PhiSH = np.copy(olr_max_lat_SH) elif method == 'peak': PhiNH = TropD_Calculate_MaxLat(olr[(lat > eq_boundary) & (lat < subpolar_boundary)],\ lat[(lat > eq_boundary) & (lat < subpolar_boundary)],n=30) PhiSH = TropD_Calculate_MaxLat(olr[(lat > -subpolar_boundary) & (lat < -eq_boundary)],\ lat[(lat > -subpolar_boundary) & (lat < -eq_boundary)], n=30) else: print('TropD_Metric_OLR: unrecognized method ', method) PhiNH = np.empty(0) PhiSH = np.empty(0) return PhiSH, PhiNH
[docs]def TropD_Metric_PE(pe,lat,method='zero_crossing',lat_uncertainty=0.0): ''' TropD Precipitation minus Evaporation (PE) metric Args: pe(lat,): zonal-mean precipitation minus evaporation lat: equally spaced latitude column vector method (str): 'zero_crossing': the first latitude poleward of the subtropical minimum where P-E changes from negative to positive values. Only one method so far. lat_uncertainty (float, optional): The minimal distance allowed between the first and second zero crossings along lat Returns: tuple: PhiSH (ndarray), PhiNH (ndarray) Latitude of first subtropical P-E zero crossing in SH and NH ''' try: assert(method in ['zero_crossing']) except AssertionError: print('TropD_Metric_PE: ERROR : unrecognized method ', method) # make latitude vector monotonically increasing if lat[-1] < lat[0]: pe = np.flip(pe) lat = np.flip(lat) # The gradient of PE is used to determine whether PE becomes positive at the zero crossing ped = np.interp(lat, (lat[:-1] + lat[1:])/2.0, np.diff(pe)) # define latitudes of boundaries certain regions eq_boundary = 5 subpolar_boundary = 50 polar_boundary = 60 # NH M1 = TropD_Calculate_MaxLat(-pe[(lat > eq_boundary) & (lat < subpolar_boundary)],\ lat[(lat > eq_boundary) & (lat < subpolar_boundary)], n=30) ZC1 = TropD_Calculate_ZeroCrossing(pe[(lat > M1) & (lat < polar_boundary)], \ lat[(lat > M1) & (lat < polar_boundary)], lat_uncertainty=lat_uncertainty) if np.interp(ZC1, lat, ped) > 0: PhiNH = ZC1 else: PhiNH = TropD_Calculate_ZeroCrossing(pe[(lat > ZC1) & (lat < polar_boundary)], \ lat[(lat > ZC1) & (lat < polar_boundary)], lat_uncertainty=lat_uncertainty) # SH # flip arrays to find the most equatorward zero crossing M1 = TropD_Calculate_MaxLat(np.flip(-pe[(lat < -eq_boundary) & (lat > -subpolar_boundary)],0),\ np.flip(lat[(lat < -eq_boundary) & (lat > -subpolar_boundary)],0), n=30) ZC1 = TropD_Calculate_ZeroCrossing(np.flip(pe[(lat < M1) & (lat > -polar_boundary)],0), \ np.flip(lat[(lat < M1) & (lat > -polar_boundary)],0), lat_uncertainty=lat_uncertainty) if np.interp(ZC1, lat, ped) < 0: PhiSH = ZC1 else: PhiSH = TropD_Calculate_ZeroCrossing(np.flip(pe[(lat < ZC1) & (lat > -polar_boundary)],0), \ np.flip(lat[(lat < ZC1) & (lat > -polar_boundary)],0), lat_uncertainty=lat_uncertainty) return PhiSH, PhiNH
[docs]def TropD_Metric_PSI(V, lat, lev, method='Psi_500', lat_uncertainty=0): ''' TropD Mass streamfunction (PSI) metric Latitude of the meridional mass streamfunction subtropical zero crossing Args: V(lat,lev): zonal-mean meridional wind lat: latitude vector lev: vertical level vector in hPa units method (str, optional): 'Psi_500'(default): Zero crossing of the stream function (Psi) at the 500hPa level 'Psi_500_10Perc': Crossing of 10# of the extremum value of Psi in each hemisphre at the 500hPa level 'Psi_300_700': Zero crossing of Psi vertically averaged between the 300hPa and 700 hPa levels 'Psi_500_Int': Zero crossing of the vertically-integrated Psi at the 500 hPa level 'Psi_Int' : Zero crossing of the column-averaged Psi lat_uncertainty (float, optional): The minimal distance allowed between the first and second zero crossings. For example, for lat_uncertainty = 10, the function will return a NaN value if a second zero crossings is found within 10 degrees of the most equatorward zero crossing. Returns: tuple: PhiSH (ndarray), PhiNH (ndarray) Latitude of Psi zero crossing in SH and NH ''' try: assert (lat_uncertainty >= 0) except AssertionError: print('TropD_Metric_PSI: ERROR : lat_uncertainty must be >= 0') try: assert(method in ['Psi_500','Psi_500_10Perc','Psi_300_700','Psi_500_Int','Psi_Int']) except AssertionError: print('TropD_Metric_PSI: ERROR : unrecognized method ', method) subpolar_boundary = 30 polar_boundary = 60 Psi = TropD_Calculate_StreamFunction(V, lat, lev) Psi[np.isnan(Psi)]=0 # make latitude vector monotonically increasing if lat[-1] < lat[0]: Psi = np.flip(Psi, 0) lat = np.flip(lat, 0) COS = np.repeat(np.cos(lat*np.pi/180), len(lev), axis=0).reshape(len(lat),len(lev)) if ( method == 'Psi_500' or method == 'Psi_500_10Perc'): # Use Psi at the level nearest to 500 hPa P = Psi[:,find_nearest(lev, 500)] elif method == 'Psi_300_700': # Use Psi averaged between the 300 and 700 hPa level P = np.trapz(Psi[:,(lev <= 700) & (lev >= 300)] * COS[:,(lev <= 700) & (lev >= 300)],\ lev[(lev <= 700) & (lev >= 300)]*100, axis=1) elif method == 'Psi_500_Int': # Use integrated Psi from p=0 to level mearest to 500 hPa PPsi_temp = sp.integrate.cumtrapz(Psi*COS, lev, axis=1) PPsi = np.zeros(np.shape(Psi)) PPsi[:,1:] = PPsi_temp P = PPsi[:,find_nearest(lev, 500)] elif method == 'Psi_Int': # Use vertical mean of Psi P = np.trapz(Psi*COS, lev, axis=1) else: print('TropD_Metric_PSI: ERROR : Unrecognized method ', method) # 1. Find latitude of maximal (minimal) tropical Psi in the NH (SH) # 2. Find latitude of minimal (maximal) subtropical Psi in the NH (SH) # 3. Find the zero crossing between the above latitudes # NH Lmax = TropD_Calculate_MaxLat(P[(lat > 0) & (lat < subpolar_boundary)],\ lat[(lat > 0) & (lat < subpolar_boundary)]) Lmin = TropD_Calculate_MaxLat(-P[(lat > Lmax) & (lat < polar_boundary)],\ lat[(lat > Lmax) & (lat < polar_boundary)]) if method == 'Psi_500_10Perc': Pmax = max(P[(lat > 0) & (lat < subpolar_boundary)]) PhiNH = TropD_Calculate_ZeroCrossing(P[(lat > Lmax) & (lat < Lmin)] - 0.1*Pmax,\ lat[(lat > Lmax) & (lat < Lmin)]) else: PhiNH = TropD_Calculate_ZeroCrossing(P[(lat > Lmax) & (lat < Lmin)],\ lat[(lat > Lmax) & (lat < Lmin)], lat_uncertainty=lat_uncertainty) # SH Lmax = TropD_Calculate_MaxLat(-P[(lat < 0) & (lat > -subpolar_boundary)],\ lat[(lat < 0) & (lat > -subpolar_boundary)]) Lmin = TropD_Calculate_MaxLat(P[(lat < Lmax) & (lat > -polar_boundary)],\ lat[(lat < Lmax) & (lat > -polar_boundary)]) if method == 'Psi_500_10Perc': Pmin = min(P[(lat < 0) & (lat > -subpolar_boundary)]) PhiSH = TropD_Calculate_ZeroCrossing(np.flip(P[(lat < Lmax) & (lat > Lmin)], 0) + 0.1*Pmin,\ np.flip(lat[(lat < Lmax) & (lat > Lmin)], 0)) else: PhiSH = TropD_Calculate_ZeroCrossing(np.flip(P[(lat < Lmax) & (lat > Lmin)], 0),\ np.flip(lat[(lat < Lmax) & (lat > Lmin)], 0), lat_uncertainty=lat_uncertainty) return PhiSH, PhiNH
[docs]def TropD_Metric_PSL(ps, lat, method='peak', n=0): ''' TropD Sea-level pressure (PSL) metric Latitude of maximum of the subtropical sea-level pressure Args: ps(lat,): sea-level pressure lat: equally spaced latitude column vector method (str, optional): 'peak' (default) | 'max' Returns: tuple: PhiSH (ndarray), PhiNH (ndarray) Latitude of subtropical sea-level pressure maximum SH and NH ''' try: assert(method in ['max','peak']) except AssertionError: print('TropD_Metric_PSL: ERROR : unrecognized method ', method) try: assert (not hasattr(n, "__len__") and n >= 0) except AssertionError: print 'TropD_Metric_PSL: ERROR : the smoothing parameter n must be >= 0' eq_boundary = 15 polar_boundary = 60 if method == 'max': if n: PhiNH = TropD_Calculate_MaxLat(ps[(lat > eq_boundary) & (lat < polar_boundary)],\ lat[(lat > eq_boundary) & (lat < polar_boundary)],n=n) PhiSH = TropD_Calculate_MaxLat(ps[(lat > -polar_boundary) & (lat < -eq_boundary)],\ lat[(lat > -polar_boundary) & (lat < -eq_boundary)], n=n) else: PhiNH = TropD_Calculate_MaxLat(ps[(lat > eq_boundary) & (lat < polar_boundary)],\ lat[(lat > eq_boundary) & (lat < polar_boundary)]) PhiSH = TropD_Calculate_MaxLat(ps[(lat > -polar_boundary) & (lat < -eq_boundary)],\ lat[(lat > -polar_boundary) & (lat < -eq_boundary)]) elif method == 'peak': if n: PhiNH = TropD_Calculate_MaxLat(ps[(lat > eq_boundary) & (lat < polar_boundary)],\ lat[(lat > eq_boundary) & (lat < polar_boundary)], n=n) PhiSH = TropD_Calculate_MaxLat(ps[(lat > -polar_boundary) & (lat < -eq_boundary)],\ lat[(lat > -polar_boundary) & (lat < -eq_boundary)], n=n) else: PhiNH = TropD_Calculate_MaxLat(ps[(lat > eq_boundary) & (lat < polar_boundary)],\ lat[(lat > eq_boundary) & (lat < polar_boundary)], n=30) PhiSH = TropD_Calculate_MaxLat(ps[(lat > -polar_boundary) & (lat < -eq_boundary)],\ lat[(lat > -polar_boundary) & (lat < -eq_boundary)], n=30) else: print('TropD_Metric_PSL: ERROR: unrecognized method ', method) return PhiSH, PhiNH
[docs]def TropD_Metric_STJ(U, lat, lev, method='adjusted_peak', n=0): ''' TropD Subtropical Jet (STJ) metric Args: U(lat,lev): zonal mean zonal wind lat: latitude vector lev: vertical level vector in hPa units method (str, optional): 'adjusted_peak': Latitude of maximum (smoothing parameter n=30) of the zonal wind averaged between the 100 and 400 hPa levels minus the zonal mean zonal wind at the level closes to the 850 hPa level, poleward of 10 degrees and equatorward of the Eddy Driven Jet latitude 'adjusted_max' : Latitude of maximum (smoothing parameter n=6) of the zonal wind averaged between the 100 and 400 hPa levels minus the zonal mean zonal wind at the level closes to the 850 hPa level, poleward of 10 degrees and equatorward of the Eddy Driven Jet latitude 'core_peak': Latitude of maximum of the zonal wind (smoothing parameter n=30) averaged between the 100 and 400 hPa levels, poleward of 10 degrees and equatorward of 70 degrees 'core_max': Latitude of maximum of the zonal wind (smoothing parameter n=6) averaged between the 100 and 400 hPa levels, poleward of 10 degrees and equatorward of 70 degrees Returns: tuple: PhiSH (ndarray), PhiNH (ndarray) Latitude of STJ SH and NH ''' try: assert (not hasattr(n, "__len__") and n >= 0) except AssertionError: print('TropD_Metric_STJ: ERROR : the smoothing parameter n must be >= 0') try: assert(method in ['adjusted_peak','core_peak','adjusted_max','core_max']) except AssertionError: print('TropD_Metric_STJ: ERROR : unrecognized method ', method) eq_boundary = 10 polar_boundary = 60 lev_int = lev[(lev >= 100) & (lev <= 400)] if (method == 'adjusted_peak' or method == 'adjusted_max'): idx_850 = find_nearest(lev, 850) # Pressure weighted vertical mean of U minus near surface U if len(lev_int) > 1: u = np.trapz(U[:, (lev >= 100) & (lev <= 400)], lev_int, axis=1) \ / (lev_int[-1] - lev_int[0]) - U[:,idx_850] else: u = np.mean(U[:,(lev >= 100) & (lev <= 400)], axis=1) - U[:,idx_850] elif (method == 'core_peak' or method == 'core_max'): # Pressure weighted vertical mean of U if len(lev_int) > 1: u = np.trapz(U[:, (lev >= 100) & (lev <= 400)], lev_int, axis=1) \ / (lev_int[-1] - lev_int[0]) else: u = np.mean(U[:, (lev >= 100) & (lev <= 400)], axis=1) else: print('TropD_Metric_STJ: unrecognized method ', method) print('TropD_Metric_STJ: optional methods are: adjusted_peak (default), adjusted_max, core_peak, core_max') if method == 'core_peak': if n: PhiNH = TropD_Calculate_MaxLat(u[(lat > eq_boundary) & (lat < polar_boundary)],\ lat[(lat > eq_boundary) & (lat < polar_boundary)], n=n) PhiSH = TropD_Calculate_MaxLat(u[(lat > -polar_boundary) & (lat < -eq_boundary)],\ lat[(lat > -polar_boundary) & (lat < -eq_boundary)], n=n) else: PhiNH = TropD_Calculate_MaxLat(u[(lat > eq_boundary) & (lat < polar_boundary)],\ lat[(lat > eq_boundary) & (lat < polar_boundary)],n=30) PhiSH = TropD_Calculate_MaxLat(u[(lat > -polar_boundary) & (lat < -eq_boundary)],\ lat[(lat > -polar_boundary) & (lat < -eq_boundary)], n=30) elif method == 'core_max': if n: PhiNH = TropD_Calculate_MaxLat(u[(lat > eq_boundary) & (lat < polar_boundary)],\ lat[(lat > eq_boundary) & (lat < polar_boundary)], n=n) PhiSH = TropD_Calculate_MaxLat(u[(lat > -polar_boundary) & (lat < - eq_boundary)],\ lat[(lat > -polar_boundary) & (lat < -eq_boundary)], n=n) else: PhiNH = TropD_Calculate_MaxLat(u[(lat > eq_boundary) & (lat < polar_boundary)],\ lat[(lat > eq_boundary) & (lat < polar_boundary)], n=6) PhiSH = TropD_Calculate_MaxLat(u[(lat > -polar_boundary) & (lat < -eq_boundary)],\ lat[(lat > -polar_boundary) & (lat < -eq_boundary)], n=6) elif method == 'adjusted_peak': PhiSH_EDJ, PhiNH_EDJ = TropD_Metric_EDJ(U,lat,lev) if n: PhiNH = TropD_Calculate_MaxLat(u[(lat > eq_boundary) & (lat < PhiNH_EDJ)],\ lat[(lat > eq_boundary) & (lat < PhiNH_EDJ)], n=n) PhiSH = TropD_Calculate_MaxLat(u[(lat > PhiSH_EDJ) & (lat < -eq_boundary)],\ lat[(lat > PhiSH_EDJ) & (lat < -eq_boundary)], n=n) else: PhiNH = TropD_Calculate_MaxLat(u[(lat > eq_boundary) & (lat < PhiNH_EDJ)],\ lat[(lat > eq_boundary) & (lat < PhiNH_EDJ)], n=30) PhiSH = TropD_Calculate_MaxLat(u[(lat > PhiSH_EDJ) & (lat < -eq_boundary)],\ lat[(lat > PhiSH_EDJ) & (lat < -eq_boundary)], n=30) elif method == 'adjusted_max': PhiSH_EDJ,PhiNH_EDJ = TropD_Metric_EDJ(U,lat,lev) if n: PhiNH = TropD_Calculate_MaxLat(u[(lat > eq_boundary) & (lat < PhiNH_EDJ)],\ lat[(lat > eq_boundary) & (lat < PhiNH_EDJ)], n=n) PhiSH = TropD_Calculate_MaxLat(u[(lat > PhiSH_EDJ) & (lat < -eq_boundary)],\ lat[(lat > PhiSH_EDJ) & (lat < -eq_boundary)], n=n) else: PhiNH = TropD_Calculate_MaxLat(u[(lat > eq_boundary) & (lat < PhiNH_EDJ)],\ lat[(lat > eq_boundary) & (lat < PhiNH_EDJ)], n=6) PhiSH = TropD_Calculate_MaxLat(u[(lat > PhiSH_EDJ) & (lat < -eq_boundary)],\ lat[(lat > PhiSH_EDJ) & (lat < -eq_boundary)], n=6) return PhiSH, PhiNH
[docs]def TropD_Metric_TPB(T, lat, lev, method='max_gradient', n=0, Z=None, Cutoff=15*1000): ''' TropD Tropopause break (TPB) metric Args: T(lat,lev): temperature (K) lat: latitude vector lev: pressure levels column vector in hPa method (str, optional): 'max_gradient' (default): The latitude of maximal poleward gradient of the tropopause height 'cutoff': The most equatorward latitude where the tropopause crosses a prescribed cutoff value 'max_potemp': The latitude of maximal difference between the potential temperature at the tropopause and at the surface Z(lat,lev) (optional): geopotential height (m) Cutoff (float, optional): geopotential height (m) cutoff that marks the location of the tropopause break Returns: tuple: PhiSH (ndarray), PhiNH (ndarray) Latitude of tropopause break SH and NH ''' Rd = 287.04 Cpd = 1005.7 k = Rd / Cpd try: assert (not hasattr(n, "__len__") and n >= 0) except AssertionError: print('TropD_Metric_TPB: ERROR : the smoothing parameter n must be >= 0') try: assert(method in ['max_gradient','max_potemp','cutoff']) except AssertionError: print('TropD_Metric_TPB: ERROR : unrecognized method ', method) polar_boundary = 60 if method == 'max_gradient': Pt = TropD_Calculate_TropopauseHeight(T,lev) Ptd = np.diff(Pt) / (lat[1] - lat[0]) lat2 = (lat[1:] + lat[:-1]) / 2 if (n >= 1): PhiNH = TropD_Calculate_MaxLat(Ptd[:,(lat2 > 0) & (lat2 < polar_boundary)],\ lat2[(lat2 > 0) & (lat2 < polar_boundary)], n=n) PhiSH = TropD_Calculate_MaxLat(-Ptd[:,(lat2 > -polar_boundary) & (lat2 < 0)],\ lat2[(lat2 > -polar_boundary) & (lat2 < 0)], n=n) else: PhiNH = TropD_Calculate_MaxLat(Ptd[:,(lat2 > 0) & (lat2 < polar_boundary)],\ lat2[(lat2 > 0) & (lat2 < polar_boundary)]) PhiSH = TropD_Calculate_MaxLat(-Ptd[:,(lat2 > -polar_boundary) & (lat2 < 0)],\ lat2[(lat2 > -polar_boundary) & (lat2 < 0)]) elif method == 'max_potemp': XF = np.tile((lev / 1000) ** k, (len(lat), 1)) PT = T / XF Pt, PTt = TropD_Calculate_TropopauseHeight(T, lev, Z=PT) PTdif = PTt - np.nanmin(PT, axis = 1) if (n >= 1): PhiNH = TropD_Calculate_MaxLat(PTdif[:,(lat > 0) & (lat < polar_boundary)],\ lat[(lat > 0) & (lat < polar_boundary)], n=n) PhiSH = TropD_Calculate_MaxLat(PTdif[:,(lat > - polar_boundary) & (lat < 0)],\ lat[(lat > -polar_boundary) & (lat < 0)], n=n) else: PhiNH = TropD_Calculate_MaxLat(PTdif[:,(lat > 0) & (lat < polar_boundary)],\ lat[(lat > 0) & (lat < polar_boundary)], n=30) PhiSH = TropD_Calculate_MaxLat(PTdif[:,(lat > - polar_boundary) & (lat < 0)],\ lat[(lat > -polar_boundary) & (lat < 0)], n=30) elif method == 'cutoff': Pt, Ht = TropD_Calculate_TropopauseHeight(T, lev, Z) # make latitude vector monotonically increasing if lat[-1] < lat[0]: Ht = np.flip(np.squeeze(Ht),0) lat = np.flip(lat,0) polar_boundary = 60 PhiNH = TropD_Calculate_ZeroCrossing(Ht[(lat > 0) & (lat < polar_boundary)] - Cutoff, lat[(lat > 0) & (lat < polar_boundary)]) PhiSH = TropD_Calculate_ZeroCrossing(np.flip(Ht[(lat < 0) & (lat > -polar_boundary)], 0) - Cutoff, np.flip(lat[(lat < 0) & (lat > -polar_boundary)], 0)) else: print('TropD_Metric_TPB: ERROR : Unrecognized method ', method) return PhiSH, PhiNH
[docs]def TropD_Metric_UAS(U, lat, lev=np.array([1]), method='zero_crossing', lat_uncertainty = 0): ''' TropD near-surface zonal wind (UAS) metric Args: U(lat,lev) or U (lat,)-- Zonal mean zonal wind. Also takes surface wind lat: latitude vector lev: vertical level vector in hPa units. lev=np.array([1]) for single-level input zonal wind U(lat,) method (str): 'zero_crossing': the first subtropical latitude where near-surface zonal wind changes from negative to positive lat_uncertainty (float, optional): the minimal distance allowed between the first and second zero crossings Returns: tuple: PhiSH (ndarray), PhiNH (ndarray) Latitude of first subtropical zero crossing of the near surface zonal wind in SH and NH ''' try: assert (lat_uncertainty >= 0) except AssertionError: print('TropD_Metric_PSI: ERROR : lat_uncertainty must be >= 0') try: assert(method in ['zero_crossing']) except AssertionError: print('TropD_Metric_PSI: ERROR : unrecognized method ', method) if len(lev) > 1: uas = U[:,find_nearest(lev, 850)] else: uas = np.copy(U) # make latitude vector monotonically increasing if lat[-1] < lat[0]: uas = np.flip(uas) lat = np.flip(lat) # define latitudes of boundaries certain regions eq_boundary = 5 subpolar_boundary = 30 polar_boundary = 60 # NH uas_min_lat_NH = TropD_Calculate_MaxLat(-uas[(lat > eq_boundary) & (lat < subpolar_boundary)],\ lat[(lat > eq_boundary) & (lat < subpolar_boundary)]) # SH uas_min_lat_SH = TropD_Calculate_MaxLat(-uas[(lat > -subpolar_boundary) & (lat < -eq_boundary)],\ lat[(lat > -subpolar_boundary) & (lat < -eq_boundary)]) try: assert(method == 'zero_crossing') PhiNH = TropD_Calculate_ZeroCrossing(uas[(lat > uas_min_lat_NH) & (lat < polar_boundary)],\ lat[(lat > uas_min_lat_NH) & (lat < polar_boundary)], lat_uncertainty=lat_uncertainty) # flip arrays to find the most equatorward zero crossing PhiSH = TropD_Calculate_ZeroCrossing(np.flip(uas[(lat < uas_min_lat_SH) & (lat > -polar_boundary)],0),\ np.flip(lat[(lat < uas_min_lat_SH) & (lat > -polar_boundary)],0), lat_uncertainty=lat_uncertainty) return PhiSH, PhiNH except AssertionError: print('TropD_Metric_UAS: ERROR : unrecognized method ', method)