CDAT GrADS Package
			  Mike Fiorino
			   PCMDI/LLNL
			  20 May, 2001

1.0  introduction
-----------------
This document assume you have a basic familiarity with cdat and
creating gridded data object using the cdms module.
The purpose of this add-on is to give users of the Grid Analysis
and Display System (GrADS) a more familiar way to query the state
of the cdat system regarding data files and variables and to
provide grads-style data slicing or subsetting.
2.0  startup
------------
The recommended way to import the grads package is:
import grads as g
from grads import q
The first command creates an object 'g' to invoke the grads
selectors and the second gives the grads 'Query' object as 'q'
3.O  query
----------
There a four basic methods (functions) in the grads query class:
1) 'files' to list all open cdms datasets; 2) 'file' to list
variables in an open dataset; 3) 'vars' to list variables that
have been defined; and 4) 'dim' to list the dimension state of
either a file or transient variable.
3.1  files
----------
Starting python and running the following:
import cdms 
import grads as g
from grads import q
#
#  data director
#
ddir='./'
f=cdms.open(ddir+'model.ctl')
#
#  list opened files
#
q.files
cdat output:
CDMS files/datasets open:
File/Dataset ID         Path
---------------         ----
f                       /home/fiorino/grads/tutorial/model.ctl
3.1.1 files -- grads 
--------------------
the comparable command and output from grads
'open 'ddir'model.ctl
'q files'
grads output:
ga-> q files
File 1 : "Sample Model Data for GrADS Tutorial"
  Descriptor: /home/fiorino/grads/tutorial/model.ctl
  Binary: /home/fiorino/grads/tutorial/model.grb
3.2 file
--------
Now list the variables in the file:
#
#  run the 'file' method of q on the file object 'f'
#
q.file(f)
>>> q.file(f)
output:
CDMS file: /home/fiorino/grads/tutorial/model.ctl contains the following variables:
Variable          Shape             Order   Description                                 units     
--------          -----             -----   -----------                                 -----     
ts              : (5, 46, 72)     : tyx   : Surface (2m) air temperature [K]         :            
ps              : (5, 46, 72)     : tyx   : Surface pressure [hPa]                   :            
z               : (5, 7, 46, 72)  : tzyx  : Geopotential height [m]                  :            
v               : (5, 7, 46, 72)  : tzyx  : Northward wind [m/s]                     :            
t               : (5, 7, 46, 72)  : tzyx  : Air Temperature [K]                      :            
u               : (5, 7, 46, 72)  : tzyx  : Eastward wind [m/s]                      :            
p               : (5, 46, 72)     : tyx   : Total precipitation rate [kg/(m^2*s)]    :            
q               : (5, 7, 46, 72)  : tzyx  : Specific humidity [kg/kg]                :        
    
the variable name, shape (size of each dimension or axis), axis
order, descriptions and units are returned, e.g.,
'ts' is a 5x46x72 array  where,
the first and slowest varying dimension is 't' or a time axis of size 5
the second  dimension is a latitude dimension 'y' with 46 points and
the third or fastest varying dimension is 'x' and is a  longitude dimension
the order of the dimension is important for slicing, but not when
using selectors.
>>> vf=f.getVariable ('v')
>>> q.vars
CDMS variables:
Variable (id)         Type         Shape                 Order
------------          ----         -----                 -----
vf                    file         (5, 7, 46, 72)        tzyx
>>> q.dim(vf)
Dimension State for: vf   which is a File      Variable: v
X is varying Lon =    0.0 to  355.0  X = 1 to 72
Y is varying Lat =  -90.0 to   90.0  Y = 1 to 46
Z is varying Lev =  1000.0 to  100.0  Z = 1 to 7
T is varying Time = 1987:1:2:0 to 1987:1:6:0  T = 1 to 5
3.2.2 file -- grads
-------------------
The comparable grads command and output is
'q file'
grads output:
ga-> q file
File 1 : "Sample Model Data for GrADS Tutorial"
  Descriptor: /home/fiorino/grads/tutorial/model.ctl
  Binary: /home/fiorino/grads/tutorial/model.grb
  Type = Gridded
  Xsize = 72  Ysize = 46  Zsize = 7  Tsize = 5
  Number of Variables = 8
    ps 0 1 Surface pressure [hPa]
    u 7 33 Eastward wind [m/s]
    v 7 34 Northward wind [m/s]
    z 7 7 Geopotential height [m]
    t 7 11 Air Temperature [K]
    q 7 51 Specific humidity [kg/kg]
    ts 0 11 Surface (2m) air temperature [K]
    p 0 59 Total precipitation rate [kg/(m^2*s)]
which contains similar information, but is not as complete as in cdat
3.3  vars and dim
-----------------
We will now create a 'file variable' to illustrate the vars and
dim methods (functions):
>>> vf=f.getVariable ('v')
>>> q.vars
CDMS variables:
Variable (id)         Type         Shape                 Order
------------          ----         -----                 -----
vf                    file         (5, 7, 46, 72)        tzyx
>>> q.dim(vf)
Dimension State for: vf   which is a File      Variable: v
X is varying Lon =    0.0 to  355.0  X = 1 to 72
Y is varying Lat =  -90.0 to   90.0  Y = 1 to 46
Z is varying Lev =  1000.0 to  100.0  Z = 1 to 7
T is varying Time = 1987:1:2:0 to 1987:1:6:0  T = 1 to 5
The vars method lists all the variables that have been created
and their type.  Since we have not done any slicing, the shape
and order information are the same as from the 'q.file(f)' command
The dim method applied to 'vf' gives information on the dimension
state of all dimensions of tv including the range in world
coordinates (e.g., latitude or Lat = -90.0 to 90.0) and grid
coordinates (e.g., X = 1 to 72 for longitude)
3.3.1  vars and dim -- grads
----------------------------
There is no comparable command to vars in GrADS, except for
defined variables where 'q define' list these.  Further, the
dimension state is global and linked to a file, not to individual
variables, but,
'q dims'
outputs:
Default file number is: 1 
X is varying   Lon = 0 to 360   X = 1 to 73
Y is varying   Lat = -90 to 90   Y = 1 to 46
Z is fixed     Lev = 1000  Z = 1
T is fixed     Time = 00Z02JAN1987  T = 1
which shows similar information, but that the global dimension
environment is different.  Here Z and T are fixed and X,Y vary to
cover a full cycle of the wrapped X or longitude dimension.
GrADS has other query functions and some may be implemented in
the future.
4.0  COORDINATE SELECTORS
-------------------------
cdat variables are 'data objects' that consist of a ' masked
array' (the data), metadata about the axis and coordinates and
methods (functions).
The most basic operation on a data object is to extract subsets
or 'slices' of the data array.  There are many options for
slicing, but the simplest for user applications is to use
'selectors.'  These selectors are actually objects which
calculate the grids points to extract based on some condition.  A
typical condition might be to specify a range in world
coordinates such as all grid points between latitude 30S and
30N.
The standard selector in cdms is the 'axis' selector which takes the form
axis=(axisvaluebeg,axisvalueend,axisintersection_option)
where,
axis			   name of the dimension or axis (e.g., 'longitude'), 
axisvaluebeg		   beginning value, assuming monotonic variation, e.g., 120W
axisvalueend		   ending value, e.g., 600W
axisintersection_option	   specify whether the endpoints of the
			   interval (axisvalubeg,axisvalueend)
			   are 'open' or 'closed' and how
			   extraction depends on the intersecion between
			   the world coordinate and grid point coordinate axis
An alternative to the axis selector is to create standalone
selector objects (i.e., no '=' as in the axis selector) and in
the grads package we have implemented three kinds:
"world coordinate"  -	  lon,lat,lev,time
"grid coordinate"   -	  x,y,z,t
"index coordinate"  -	  ndx1,ndx2,ndx3,ndx4
the world and grid selectors operate in a similar to grads whereas the
index selector slices by order of the dimensions.
4.1  grid  and index coordinate selectors
-----------------------------------------
Since python is derived from the C programming language arrays
and indices follow the C convention where indices start with 0
and the first dimension is the slowest varying.  Further, python
slice intervals are closed on the right and opened on the left.
For example, to extract all points from a 20 point array A using
would be:
A[0:21]
In contrast, FORTRAN indices start with 1 and the left interval
is closed, e.g,
print*,(A(i),i=1,20)
Since GrADS was designed for the modeling community where FORTRAN
is still used, it slices following the FORTRAN convention.
Returning to the example in section 3 above where we opened the
tutorial file as f and then,
vf=f.getVariable ('v')
where vf is a file variable and when we
q.dim(vf) 
we find,
Dimension State for: fv   which is a File      Variable: v
X is varying Lon =    0.0 to  355.0  X = 1 to 72
Y is varying Lat =  -90.0 to   90.0  Y = 1 to 46
Z is varying Lev =  1000.0 to  100.0  Z = 1 to 7
T is varying Time = 1987:1:2:0 to 1987:1:6:0  T = 1 to 5
and from 
q.file(f)
v               : (5, 7, 46, 72)  : tzyx  : Northward wind [m/s]                     :      
or that vf is a 4-D array where t is the slowest varying (C
ordering) and x is the fastest.
Let's now extract the first time into a transient variable
(occupies memory) using a grads grid selector,
v1=vf(g.t(1))
v1 is the transient variable, g.t is the grads (g) selector t and
1 is the first time point.
q.dim(v1) shows:
Dimension State for: v1   which is a Transient Variable from File Variable: v
X is varying Lon =    0.0 to  355.0  X = 1 to 72
Y is varying Lat =  -90.0 to   90.0  Y = 1 to 46
Z is varying Lev =  1000.0 to  100.0  Z = 1 to 7
T is fixed   Time = 1987:1:2:0  T = 1
the 'T' dimension is fixed at the first point.  
To extract the 2nd and 3rd time:
v2=vf(g.t(2,3))
and q.dim(v2) shows:
Dimension State for: v2   which is a Transient Variable from File Variable: v
X is varying Lon =    0.0 to  355.0  X = 1 to 72
Y is varying Lat =  -90.0 to   90.0  Y = 1 to 46
Z is varying Lev =  1000.0 to  100.0  Z = 1 to 7
T is varying Time = 1987:1:3:0 to 1987:1:4:0  T = 1 to 2
note that the times are correct (daily data) but that t grid
coordinate ranges from '1 to 2' instead of '2 to 3'.  The reason
is that v2 is a full-fledged, standalone data object with its own
mapping between world and grid coordinates.  The variable in the
file ranges from T = 1 to 5 but now that we have an object, it
has its own coordinate.  
The grid coordinate selector also supports striding as in grads, e.g.,
v3=vf(g.t(1,5,2))
pulls out the first, third and fifth time, i.e, skipping by 2.
Note that the file variable vf has the order:
'tzyx'
the x,y,z,t indicators mean that axises have been associated
with the 'standard' world coordinates, 'longitude',
'latitude','level','time'.  If such associations where not found
then the order would have been:
'----'
or that slicing using the x,y,z,t, selectors would not work.
Instead use the 'index coordinate' selectors, ndx1, ndx2...
These are positional in that you have to ndx1 works on the first
dimension in array, e.g., the slice expression:
vf(g.t(1,5,2) 
is equivalent to
vf(g.ndx1(1,5,2)
One of the benefits of selectors is that there is no requirement
on ordering in the slice, e.g.,
v4=vf(g.z(1),g.x(1,30),g.y(1,10),g.t(1))
returns the same object as,
v4=vf(g.t(1),g.x(1,30),g.y(1,10),g.z(1))
4.2  world coordinate selectors
-------------------------------
The grads world coordinate selectors 'lon', 'lat','lev','time'
are included mostly for completeness since they operate
identically to axis selectors. i.e.,
vf(g.lon(0,180))
is equivalent to,
vf(longitude=(0,180))
and in,
vv=vf(g.lon(0,180),g.t(1))
where we mixed selectors,
q.dims(vv)
shows
Dimension State for: vv   which is a Transient Variable from File Variable: v
X is varying Lon =    0.0 to  180.0  X = 1 to 37
Y is varying Lat =  -90.0 to   90.0  Y = 1 to 46
Z is varying Lev =  1000.0 to  100.0  Z = 1 to 7
T is fixed   Time = 1987:1:2:0  T = 1