| libbjt.gex - Ben-Jei Tsuang's Collection of GrADS v1.9 Extensions |
Minv: displaying the minimum value between EXP1 and EXP2Maxv: displaying the maximum value between EXP1 and EXP2If: selectively displaying values according to the relationship between two expressions.Which: label gridpoints according to specific conditionstfit: performing multiple linear regression on individual gridpointsfit: performing global multiple linear regressionmuadv: computing zonal advecionmvadv: computing meridional advecionmwadv: computing vertical advecion using an upwind schemecosz: computing the cosine of the solar zenith anglelt: computing the local timedayratio: computing the daytime ratiojd: computing the Julian day sice January 0001dew: computing the dew point temperaturesatvap: computing the saturated temperatureline: drawing linesvint2: modified vertical integrationzinterp: Vertical interpolation to a specified p level.zinterp: Vertical interpolation to a specified z level.
libbjt.gex - Ben-Jei Tsuang's Collection of GrADS v1.9 Extensions
display lt(EXPR) - Local time
display jd(EXPR) - Julian day
display cosz(EXPR, H|D|M) - Cosine of solar zenith angle
display dayratio(EXPR) - Daylight ratio
display ftest(RATIO,DF1,DF2) - F-Test
display ttest(DIFF,DF) - T-Test
display tfit(Y1,X1,X2,X3,..,XN,time=TBEG,time=TEND[,DT,FILENAME]) - Local multiple linear regression
display fit(Y1,X1,X2,X3,..,XN,time=TBEG,time=TEND[,DT,FILENAME]) - Global multiple linear regression
display tcorr2(X1,X2,time=TBEG,time=TEND [,DT]) - Time correlation (r)
display tregr2(X,Y,time=TBEG,time=TEND [,DT]) - Gives the expected value of y departure given a 1 unit departure in x. (Same as tregr with less contrained)
display tmave2(MASKEXPR,EXPR,time=TBEG,time=TEND [,DT]) - Time averaging with masking
display muadv(U,EXPR) - Changing rate [EXPR/s] due to zonal advection
display mvadv(V,EXPR) - Changing rate [EXPR/s] due to meridional advection
display mwadv(W,EXPR) - Changing rate [EXPR/s] due to vertical advection
display madvu(U,EXPR) - Zonal flux gradient: -d(U*EXPR)/dx
display madvv(V,EXPR) - Meridional flux gradient: -d(V*EXPR)/dy
display madvw(W,EXPR) - Vertical flux gradient: -d(W*EXPR)/dz
display satvap(T) - Saturated vapor pressure
display dew(PVAP) - Dew-point Temperature
display lw(PL,TA,WA,TAUCL,FCLD,PS,TG,EG,TB,WB,TAUCLB,FCLDB [,BINFILE,CTLFILE,-r|-c|-h|-l]) - Longwave radiative coolimg rate (K/d)
display lw2(PL,TA,WA,TAUCL,FCLD,PS,TG,EG,TB,WB,PM,TM,WM,TAUCLM,FCLDM [,BINFILE,CTLFILE,-r|-c|-h|-l]) - Longwave radiative coolimg rate (K/d) (same as lw but requires data at the top of mixed layer)
display pinterp(EXPR,PGRID,PLEV,[-l|-s|-p]) - Vertical interpolation within a 3-D grid to a specified p level.
display zinterp(field,zgrid,zlev,[-l|-s|-p]) - Vertical interpolation within a 3-D grid to a specified z level.
display vint2(PSEXPR,EXPR,TOP) - Modified vertical integration
display if(EXPR1,OP,EXPR2,TRUE_EXP,FALSE_EXP) - Compare expressions
display which(EXPR1,COND1,EXP1,COND2,EXPR2,...,ELSE_EXP) - Label gridpoints according to expressions
display maxv(EXPR1,EXPR2) - Maximum value
display minv(EXPR1,EXPR2) - Minimum value
display line(EXPR,LON1,LAT1,LON2,LAT2[,-l|-m|-r]) - Draw a Line
This library implements Ben-Jei Tsuang's collection GrADS extensions, including a variety of statistical, mathematical and otherwise utility functions.
These functions were initially built in an experimental version of GrADS v1.8 contributed by Ben-Jei Tsuang ca. 2003. With the introduction of dynamic linked extensions by the OpenGrADS project in 2006 these functions were among the first set of user defined functions to be implemented as such.
Usage notice for most of these functions can be obtained at the GrADS command line. For example,
ga-> d lw
produces
lw
Error from LW: 10-13 arguments expected
lw(pl,ta,wa,taucl,fcld,ps,tg,eg,tb,wb [,'binfile','ctlfile',-r-c-h-l])
usage:
Function returns: thermal infrared fluxes
Following the NASA Technical Memorandum (NASA/TM-2001-104606,
Vol. 19) of Chou, Suarez, Liang, and Yan (2001). This NASA TM
has been revised a few timessince. It computes thermal infrared
fluxes due to emission by water vapor, ozone, co2, o2, minor
trace gases, clouds, and aerosols and due to scattering by
clouds and aerosols.
Arguments:
pl = level pressure (hPa)
ta = layer temperature (K)
wa = layer specific humidity (kg/kg)
taucl= cloud optical thickness (dimensionless)
...
This function calculates local time (in hours); specify any valid GrADS expression on input.
This function computes the Julian day since January 0001; specify any valid GrADS expression on input.
This function calculates the cosine of solar zenith angle on various time scale. On input,
Any GrADS vaild expression, e.g., lat
Specify one of the following for the time scale:
Hourly value of the cosine of solar zenith angle.
Daily mean of the cosine of solar zenith angle.
XSMonthly mean of the cosine of solar zenith angle.
This function calculates daylight ratio of EXPR during the time interval of the variable, where 1 argument is required..
Any GrADS valid expression, e.g., lat.
This function conducts the F test. It calculates the probability for the ratio of two variances, where 3 arguments are expected.
The ration of variances: (variance 1)/(variance 2)
Degrees of freedom for variance 1
Degrees of freedom for variance 2
This function conducts the T test. It calculates the probability two means to be different, where 2 arguments are expected.
The normalized difference between 2 means:
DIFF = (m1-m2)/stdv
where m1 and m2 are the 2 means, and stdv is the standard deviation.
Degrees of freedom.
This function conducts point multiple linear regression for
Y1 = F(X1,X2,X2,..,XN)
for each grid point seperately during the period between TBEG and TEND, with time interval DT to skip. The correlation coefficients (r) at each grid are returned. If the option FILENAME is provided, many statictical summaries are written to the file.
Must be expressions.
Begining of the time, which can be any valid time variable with
formats like time=jan1990 or t=100.
End of the time period, which can be any valid time variable with
formats like time=dec2000, t=500 or t+100.
Time interval to skip, which can be any valid time variable
(e.g. 1mon, 1day).
Filename to which many statistical summaries are written to.
This function conducts global multiple linear regression for
Y1 = F(X1,X2,X2,..,XN)
for all the grid points simultaneously during the period between TBEG and TEND, with time interval DT to skip. The results are written to screen. If the option FILENAME is provided, many statictical summaries are written to the file.
Must be expressions.
Begining of the time, which can be any valid time variable with
formats like time=jan1990 or t=100.
End of the time period, which can be any valid time variable with
formats like time=dec2000, t=500 or t+100.
Time interval to skip, which can be any valid time variable
(e.g. 1mon, 1day).
Filename to which many statistical summaries are written to.
This function calculates m coefficient for the linear regression
Y = m * X + C
at each grid during the period between TBEG and TEND, with timeinterval DT to skip. (Same as tregr, but with DT option).
must be expressions.
Begining of the time, which can be any valid time variable with
formats like time=jan1990. or t=100.
End of the time, which can be any valid time variable with formats
like time=dec2000, t=500 or t+100.
Time interval to skip, which can be any valid time variable
(e.g. 1mon, 1day.).
This function does time averaging while applying a mask. (Same as tmave, but with DT option.)
The mask expression; when evaluated at a fixed time, it must give a single value.
The expression to be averaged
Begining of the time, which can be any valid time variable with
formats like time=jan1990 or t=100.
End of the time, which can be any valid time variable with formats
like time=dec2000, t=500 or t+100.
Time interval to skip, which can be any valid time variable
(e.g. 1mon, 1day).
This function works similarly to the ave function, except for the masking. Using tmave2 is much more efficient than using maskout with ave).
The function loops through the specified time steps, and evaluates MASKEXPR at each fixed time. MASKEXPR must yeild a single value. If this value is the undefined/missing data value, then expr for that time is not included in the average. If MASKEXPR is not the undefined data value, it is used as the weight for EXPR in the average. So if you define MASKEXPR accordingly, you can use the tmave2 function to do weighted time averaging.
The tricky aspect of using tmave2 is setting up MASKEXPR. If EXPR is a grid with X and/or Y and/or Z varying, then MASKEXPR must refer to either a defined variable or a file with only time varying. In general, you have to set up MASKEXPR in advance.
Say you want to average slp over some time range but only when sst
over some region is above some value. You can do this by:
ga-> set x 1 ga-> set y 1 ga-> set t 1 last ga-> define sstmask = aave(sst,lon=-180,lon=-90,lat=-20,lat=20) ga-> define sstmask = const(maskout(sstmask,sstmask-25.0),1)
Now SSTMASK is a time series where the value is 1 when the sst areal
average is above 25 and undefined when the value is below 25. maskout
set the values below 25 to missing; const set the non-missing values
to 1. We can now do our tmave2:
ga-> set lon -180 -90 ga-> set lat -20 20 ga-> set t 1 ga-> d tmave2(sstmask,slp,t=1,t=last)
The mask could also be written to a file with all dimensions nonvarying except for time. Here is what some of the records in the data descriptor file might look like:
DSET maskfilename XDEF 1 linear 1 1 YDEF 1 linear 1 1 ZDEF 1 linear 1 1 TDEF 100 linear
This function calculates
-d(u*EXPR)/dx
Zonal wind speed (u component), m/s.
An expression
This function calculates
-d(V*EXPR)/dy
Meridional wind speed (v component), m/s.
An expression
This function calculates
-d(W*EXPR)/dp
using an upwind scheme.
Vertical wind speed (w component) (unit of vertical coordinate/s, such as hPa/s for the vertical coordinate in hPa)
An expression.
Note that no-slip boundary conditions are assumed. i.e.,
w(lev=0)=0
C(lev=0)=C(lev=1)
concp, wp
|---cn
conc, w
|---cnm
concm, wm
This function computes the time rate of change EXPR at each grid point due to zonal advecion (plus increase) according to Bott's 4-level scheme.
Zonal wind speed (u component), m/s.
An expression
Bott, A., 1989a: A positive definite advection scheme obtained by nonlinear renormalization of the advection flux. Mon. Wea. Rev., 117, 1006-1015.
Bott, A., 1989b: Notes and correspondance. Mon. Wea. Rev., 117, 2633-2636.
This function computes the time rate of change EXPR at each grid point due to meridional advecion (plus increase) according to Bott's 4-level scheme.
Meridional wind speed (v component), m/s.
An expression
Bott, A., 1989a: A positive definite advection scheme obtained by nonlinear renormalization of the advection flux. Mon. Wea. Rev., 117, 1006-1015.
Bott, A., 1989b: Notes and correspondance. Mon. Wea. Rev., 117, 2633-2636.
This function computes the time rate of change EXPR at each grid point due to vertical advecion (plus increase) using an upwind scheme.
Vertical wind speed (w component) (unit of vertical coordinate/s, such as hPa/s for the vertical coordinate in hPa)
An expression
This function calculates the saturated vapor pressure (Pa)
Temperature (K)
This function calcuates the dew point temperature (K).
Vapor pressure (Pa)
This function computes the cooling rate [K/day] and thermal infrared fluxes and returns it in a file.
level pressure (hPa)
layer temperature (K)
layer specific humidity (kg/kg)
cloud optical thickness (dimensionless)
cloud fraction
surface pressure (hPa)
land or ocean surface temperature (K)
land or ocean surface emissivity (fraction)
surface air temperature (K)
surface air specific humidity (kg/kg)
surface cloud optical thickness (dimensionless)
surface cloud amount (fraction)
binary file name ('null' or blank for no output)
ctl file name ('null' or blank for no output)
using Roberts et al. (1976) water vaopr continuum
using Clough et al. (1989) water vaopr continuum (default)
using look-up table (high acc but slow) (default)
using k-dist method (low acc but fast)
Following the NASA Technical Memorandum (NASA/TM-2001-104606, Vol. 19) of Chou, Suarez, Liang, and Yan (2001). This NASA TM has been revised a few times since. It computes thermal infrared fluxes due to emission by water vapor, ozone, CO2, O2, minor trace gases, clouds, and aerosols and due to scattering by clouds and aerosols.
This function computes the cooling rate [K/day] and thermal infrared fluxes and returns it in a file. (Same as lw, but required data at the top of mixed layer)
level pressure (hPa)
layer temperature (K)
layer specific humidity (kg/kg)
cloud optical thickness (dimensionless)
cloud amount (fraction)
surface pressure (hPa)
land or ocean surface temperature (K)
land or ocean surface emissivity (fraction)
surface air temperature (K)
surface air specific humidity (kg/kg)
surface cloud optical thickness (dimensionless)
surface cloud amount (fraction)
pressure at the top of mixed layer (hPa)
temperature at the top of mixed layer (K)
specific humidity at the top of mixed layer (kg/kg)
cloud optical thickness at the top of mixed layer (dimensionless)
cloud amount at the top of mixed layer (fraction)
binary file name ('null' or blank for no output)
ctl file name ('null' or blank for no output)
using Roberts et al. (1976) water vaopr continuum
using Clough et al. (1989) water vaopr continuum (default)
using look-up table (high acc but slow) (default)
using k-dist method (low acc but fast)
Following the NASA Technical Memorandum (NASA/TM-2001-104606, Vol. 19) of Chou, Suarez, Liang, and Yan (2001). This NASA TM has been revised a few timessince. It computes thermal infrared fluxes due to emission by water vapor, ozone, co2, o2, minor trace gases, clouds, and aerosols and due to scattering by clouds and aerosols.
This function interpolates a 3-D field to a specified pressure level. Can also be used on non-pressure level data, such as sigma or eta-coordinate output where pressure is a function of time and grid level.
Defined grid vinterp holding interpolated values
Name of 3-D grid to interpolate
Name of 3-D field holding pressure values at each gridpoint. If you
are using regular pressure-level data, this should be set to the
builtin GrADS variable lev.
Pressure level at which to interpolate (can be 3-D, such as reserved
variable lev)
for piecewise linear interpolation (fast) (default)
for spline interpolation (slow & problematic while dx is small)
for polynomnial interpolation (fast but problematic with extrapolation)
This function interpolates a 3-D field to a specified z level. Can
also be used on non-z level data, such as sigma or eta-coordinate
output where z is a function of time and grid level.
Defined grid vinterp holding interpolated values
Name of 3-D grid to interpolate
Name of 3-D grid holding z values at each gridpoint
z level at which to interpolate
for piecewise linear interpolation (fast) (default)
for spline interpolation (slow & problematic while dx is small)
for polynomnial interpolation (fast but problematic with extrapolation)
This function performs a mass-weighted vertical integral in mb pressure coordinates.
The calculation is a sum of the mass-weighted layers:
f/g * sum(expr * Delta(level))
The bounds of the integration are the surface pressure and the
indicated top value. The scale factors are f=100 and g=9.8. The
summation is done for each layer present that is between the bounds.
The layers are determined by the Z levels of the default file. Each
layer is considered to be from the midpoints between the levels actually
present, and is assumed to have the same value throughout the layer, namely
the value of the gridpoint at the middle of the layer.
A GrADS expression for the surface pressure, in mb, which bounds the integral on the bottom.
A GrADS expression representing the quantity to be integrated.
The bounding top pressure, in mb. This value must be a constant and cannot be provided as an expression.
The summation is done using the Z levels from the default file, so it
is important that the default file have the same Z dimension
coordinates as EXPR. Data levels below and above the bounds of the
summation are ignored. The Z dimension in world-coordinate units is
assumed to be pressure values given in millibars (mb). The units of g
are such that when the expression integrated is specific humidity (q)
in units of g/g, the result is kg of water per square meter, or
precipitable water in mm. It is usually a good idea to make the top
pressure value to be at the top of a layer, which is midway between
grid points. For example, if the default file (and the data) have
pressure levels of ...,500,400,300,250,... then a good value for top
might be 275, the value at the top of the layer that extends from 350
to 275 mb. The vint2 function operates only in an X-Y varying dimension
environment.
A typical use of vint might be:
ga-> d vint(ps,q,275)
This expression will integrate specific humidity to obtain precipitable water, in mm.
This function compares EXPR1 and EXPR2, where 5 arguments are expected. if their relation is true, TRUE_EXP is returned, otherwise expression FALSE_EXPR is returned.
One of the following: ==, =, >, <, >=, <=, !=
These can be a expression, constant or -u.
This function labels gridpoints according to specified conditions, where an even number of arguments is expected. The number of arguments is limited by the maximum characters (about 128 charcters) being handled by the argument string. If EXPR equals to COND1, then EXPR1 is returned. Otherwise, if EXPR equals to COND2, then EXPR2 is returned. If none of the conditions are met, ELSE_EXPR is returned.
A valid GrADS expression
These can be a expression, constant or -u.
This function displays the maximum value between EXPR1 and EXPR2, where 2 arguments are expected.
Must be an expression
Can be an expression, constant or -u.
This function displays the minimum value between EXPR1 and EXPR2, where 2 arguments are expected.
Must be an expression
Can be an expression, constant or -u.
This function draws a line from (LON1,LAT1) to (LON2,LAT2). Set the line at the value VAL (default 1) using the template of EXPR.
Starting longitude/latitude
Ending longitude/latitude
return the length [km] in each grid
return the logics (1:pass, 0: non-passing)
return the length ratio of each segment (fraction)
Minv: displaying the minimum value between EXP1 and EXP2set clevs -12 -9 -6 -3 0 3 6 9 12 d minv(u,0) run cbarn.gs draw title d minv(u,0)
Maxv: displaying the maximum value between EXP1 and EXP2set clevs -12 -9 -6 -3 0 3 6 9 12 d maxv(u,0) run cbarn.gs draw title d maxv(u,0)
If: selectively displaying values according to the relationship between two expressions.d if(u,>,0,t,-u) run cbarn.gs draw title if(u,>,0,t,-u)
Which: label gridpoints according to specific conditionsset clevs 0 1 2 d which(minv(maxv(u,0),10),0,0,10,2,1) run cbarn.gs draw title u<=0: 0; 0<u<10: 1; u>=10: 2
tfit: performing multiple linear regression on individual gridpoints
set t 1
d tfit('var.0','var.3','var.4',t=1,t=100)
d tfit('var.0','var.3',t=1,t=100,1mo,"/silver/data/nchu/fit/dmdtE")
d tfit('var.0','var.3','var.4',t='mon',t+100,1mo,"/silver/data/nchu/fit/dmdtEP")
d tfit('var.0','var.3','var.4',time=jan1990,time=dec2002)
d tfit('var.0','var.3',time=jan1990,time=dec2002,1mo,"/silver/data/nchu/fit/dmdtE")
d tfit('var.0','var.3','var.4',time=jan1990,time=dec2002,1mo,"/silver/data/nchu/fit/dmdtEP")
fit: performing global multiple linear regression
set t 1
d fit('var.0','var.3','var.4',t=1,t=100)
d fit('var.0','var.3',t=1,t=100,1mo,"/silver/data/nchu/fit/dmdtE")
d fit('var.0','var.3','var.4',t='mon',t+100,1mo,"/silver/data/nchu/fit/dmdtEP")
d fit('var.0','var.3','var.4',time=jan1990,time=dec2002)
d fit('var.0','var.3',time=jan1990,time=dec2002,1mo,"/silver/data/nchu/fit/dmdtE")
d fit('var.0','var.3','var.4',time=jan1990,time=dec2002,1mo,"/silver/data/nchu/fit/dmdtEP")
muadv: computing zonal adveciond muadv(u,t) run cbarn.gs draw title -u*dT/dx (K s`a-1`n)
mvadv: computing meridional adveciond mvadv(v,t) run cbarn.gs draw title -v*dT/dy (K s`a-1`n)
mwadv: computing vertical advecion using an upwind schemed mwadv(w,t)/100 run cbarn.gs draw title -w*dT/dz (K s`a-1`n)
Notice that the 100 factor above converts the unit of w from Pa/s to
hPa/s if the vertical coordinate is in hPa.
cosz: computing the cosine of the solar zenith angled cosz(lat,h) run cbarn.gs draw title hourly cosine zenith angle
d cosz(lat,d) run cbarn.gs draw title daily mean cosine zenith angle
d cosz(lat,m) run cbarn.gs draw title monthly mean cosine zenith angle
lt: computing the local timed lt(lat) run cbarn.gs draw title local time (h)
dayratio: computing the daytime ratiod dayratio(lat) run cbarn.gs draw title daytime ratio (lat)
jd: computing the Julian day sice January 0001d jd(lat) run cbarn.gs draw title julian date from 1 January 0001 (jd=1)\d jd(lat)
dew: computing the dew point temperatured dew(q/0.622*lev*100) run cbarn.gs draw title dew point (K)
satvap: computing the saturated temperatured satvap(t) run cbarn.gs draw title saturated temperature (Pa)
line: drawing linesd line(u,30,60,200,-40,-r) run cbarn.gs draw title draw line\d line(u,lon0,lat0,lon1,lat1,-r)
vint2: modified vertical integration
d vint2(1000,dtdt,300)
run cbarn.gs
draw title Integrate[dTdt,{dp,1000,300}] (K/s*kg/m2)
d vint2(pressfc/100,dTdt,presML/100) run cbarn.gs draw title (K/s*kg/m2) d vint(ps,q,275) draw title Integrate q from surface to 275 hPa\ Total Precipitable Water (mm)
zinterp: Vertical interpolation to a specified p level.First, conduct 3-d vertical interpolation from file 3 p level to file 1 level.
define cbtzw=pinterp(cbtzw.3,lev.3,lev)
Then conduct vertical interpolation to 1000 hPa:
define sxdflmc=pinterp(sxdfl1,pl1,1000) define flxmc=pinterp(flx.'id',pl.'id',levmc.2) define tml=pinterp(tprs.2*pow(1000/lev.2,0.286),lev.2,levmh.2(z=1)) define tsfc=pinterp(tprs.2*pow(1000/lev.2,0.286),lev.2,pressfc.3(z=1)/100)
zinterp: Vertical interpolation to a specified z level.define tmp1=zinterp(tprs.2,zprs.2/9.8,h1) define tmp2=zinterp(tprs.2,zprs.2/9.8,h2) define pressfc=zinterp(lev.4,zprs.4/9.8,elv.1(t=1)) define levmh=zinterp(lev.4,zprs.4/9.8,elvmh)
http://opengrads.org/ - OpenGrADS Home Page
http://opengrads.org/wiki/index.php - OpenGrADS User Defined Extensions
http://www.iges.org/grads/ - Official GrADS Home Page
Ben-Jei Tsuang <btsuang@yahoo.com>. Arlindo da Silva <dasilva@opengrads.org> implemented these functions as dynamic User Defined Extensions.
Copyright (c) 2003-2008 by Ben-Jei Tsuang <btsuang@yahoo.com>. All Rights Reserved.
This is free software released under the GNU General Public License; see the source for copying conditions. There is NO warranty; not even for MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
| libbjt.gex - Ben-Jei Tsuang's Collection of GrADS v1.9 Extensions |