GaCore Examples

From OpenGrADS
Jump to: navigation, search

In this document we will walk you trough a simple script to exercise the most basic functionality provided by package grads. The only requirement is a recent enough version of Python (2.3 or newer), and the standard python library. This example should also work with Jython 2.2.1. As for the GrADS application, this example works with GrADS v1.9.0-rc1 and GrADS v2.0.x and newer.

The code fragments below could be typed in a file called gacore_ex.py which could executed with

% python gacore_ex.py

or typed at the classic python command line or even using the improved python shell ipython.

Getting the sample datasets[edit]

For running the examples below you need to download the following files to your local disk

Front matters[edit]

In the very top of your script you start by loading the module grads

import grads
from sys import stdout     # will need this for nicer output to the screen

Starting the GrADS application[edit]

Your next step is to start the grads application, for example the gradsnc executable.

ga = grads.GrADS(Bin='gradsnc', Window=False)

The optional parameter Window=False tells gradsnc to start in batch mode. The object ga is callable, in the sense that you can use it to send commands to the GrADS application. For example,

 ga("help")

results in the following output:

For Complete Information See:  http://grads.iges.org/grads

Basic Commands:
 OPEN <descr>      opens a data file
 Query             shows current status
 Clear             clears graphics window
 SET <args>        sets options
 Display expr      displays a variable graphically
 QUIT              exits the program

Opening a file[edit]

Let's start by purposely trying to open a file that does not exist:

try:
    fh = ga.open("wrong_file.nc")  
    print ">>> NOT OK <<< file 'wrong_file.nc' should not exist"
except GrADSError:
    print ">>> OK <<< No such file, we meant that!\n"

This code fragment also illustrates how errors reported by GrADS can be handled. Every time an abnormal condition is encountered, functions in the grads module raises a GrADSError exception. The attribute ga.rc contains the actual return code from gradsnc. The expected result from this code fragment is

>>> OK <<< No such file, we meant that!

Next let's open a file that actually exists:

try:
    fh = ga.open("model.ctl")
    print fh.title
    print '              File Id: ', fh.fid
    print '         Dataset type: ', fh.type
    print '    No. of Time steps: ', fh.nt
    print '    No. of Longitudes: ', fh.nx
    print '    No. of  Latitudes: ', fh.ny
    print '    No. of     Levels: ', fh.nz
    print '      Variable  names: ', fh.vars
    print '      Variable levels: ', fh.var_levs
    print '      Variable titles: ', fh.var_titles
    print 
    print ">>> OK <<< open CTL file"
except:
    print ">>> NOT OK <<< cannot open CTL file"

The file handle fh contains important information about the file contents. This code fragment produces the following output:

"Sample Model Data for lats4d Tutorial"
              File Id:  1
         Dataset type:  Gridded
    No. of Time steps:  5
    No. of Longitudes:  72
    No. of  Latitudes:  46
    No. of     Levels:  7
      Variable  names:  ['ps', 'ua', 'va', 'zg', 'ta', 'hus', 'ts', 'pr']
      Variable levels:  [0, 7, 7, 7, 7, 7, 0, 0]
      Variable titles:  ['Surface pressure [hPa]', 'Eastward wind [m/s]', 'Northward wind [m/s]',
  'Geopotential height [m]', 'Air Temperature [K]', 'Specific humidity [kg/kg]', 'Surface (2m) air
  temperature [K]', 'Total precipitation rate [kg/(m^2*s)]']

The method open is also able to handle NetCDF or HDF files trough sdfopen:

try:
    fh = ga.open("model.nc")   # or "model.htf"
    print fh.title
    print '              File Id: ', fh.fid
    print '         Dataset type: ', fh.type
    print '    No. of Time steps: ', fh.nt
    print '    No. of Longitudes: ', fh.nx
    print '    No. of  Latitudes: ', fh.ny
    print '    No. of     Levels: ', fh.nz
    print '      Variable  names: ', fh.vars
    print '      Variable levels: ', fh.var_levs
    print '      Variable titles: ', fh.var_titles
    print 
    print ">>> OK <<< open NetCDF file"
except:
    print ">>> NOT OK <<< cannot open NetCDF file"

which again produces

              File Id:  2
         Dataset type:  Gridded
    No. of Time steps:  5
    No. of Longitudes:  72
    No. of  Latitudes:  46
    No. of     Levels:  7
      Variable  names:  ['ps', 'ts', 'pr', 'ua', 'va', 'zg', 'ta', 'hus']
      Variable levels:  [0, 0, 0, 7, 7, 7, 7, 7]
      Variable titles:  ['Surface pressure [hPa]', 'Surface (2m) air temperature [K]', 'Total
 precipitation rate [kg/(m^2*s)]', 'Eastward wind [m/s]', 'Northward wind [m/s]', 'Geopotential
 height [m]', 'Air Temperature [K]', 'Specific humidity [kg/kg]']

Displaying a variable and saving an image[edit]

You can simply execute a sequence of GrADS commands for accomplishing just that:

ga("set gxout shaded")
ga("d ps")
ga("printim ps.png")

Querying the GraDS state[edit]

The method query allows you to conveniently capture the output of a GrADS query as attributes of the a GrADSHandle object.


To query the GrADS dimension state call with the argument 'dims':

try:
    qh = ga.query('dims')
    print 'Current dimensional state: '
    print '   X is '+qh.x_state+'   Lon = ',qh.lon,'  X = ',qh.x
    print '   Y is '+qh.x_state+'   Lat = ',qh.lat,'  Y = ',qh.y
    print '   Z is '+qh.y_state+'   Lev = ',qh.lev,'  Z = ',qh.z
    print '   T is '+qh.z_state+'  Time = ',qh.time,'  T = ',qh.t
    print 
    print ">>> OK <<< query dimensions"
except:
    print ">>> NOT OK <<< cannot query dimensions"

This code fragment results in the following output:

Current dimensional state:
   X is varying   Lon =  (0.0, 360.0)   X =  (1.0, 73.0)
   Y is varying   Lat =  (-90.0, 90.0)   Y =  (1.0, 46.0)
   Z is varying   Lev =  (1000.0, 1000.0)   Z =  (1, 1)
   T is fixed  Time =  ('00Z01JAN1987', '00Z01JAN1987')   T =  (1, 1)


To query the file state call with the argument 'file' or 'file #':

try:
    qh = ga.query('file')
    print 'Current file state: '
    print '         Title: ', qh.title
    print '       File Id: ', qh.fid
    print '   Description: ', qh.desc
    print '        Binary: ', qh.bin
    print 
    print '   Variable  Num levels  Description'

    for (var, nlevels, desc) in qh.var_info:
        print '   '+var.rjust(8)+'  '+str(nlevels).rjust(10)+'  '+desc

    print 
    print ">>> OK <<< query file"
except:
    print ">>> NOT OK <<< cannot query file"

This code fragment results in the following output:

Current file state: 
         Title:  "Sample Model Data for lats4d Tutorial"
       File Id:  1
   Description:  ../data/model.ctl
        Binary:  ../data/model.grb

   Variable  Num levels  Description
         ps           0  Surface pressure [hPa]
         ua           7  Eastward wind [m/s]
         va           7  Northward wind [m/s]
         zg           7  Geopotential height [m]
         ta           7  Air Temperature [K]
        hus           7  Specific humidity [kg/kg]
         ts           0  Surface (2m) air temperature [K]
         pr           0  Total precipitation rate [kg/(m^2*s)]

Examining GrADS standard output: the traditional way[edit]

If you have used the GrADS scripting language you are likely familiar with functions sublin() and subwrd() that gives you access to lines and words within lines of GrADS output. Two similar methods are provided here: rline() which captures whole lines, and rword(i,j) which returns the jth word from the ith line.

The following code fragment demonstrate the line interface:

try:
    ga("q config")
    print "--------------------------------------------------------------"
    print "            Captured GrADS output: Line interface"
    print "--------------------------------------------------------------"
    for i in range(1,ga.nLines):
        print ga.rline(i)
    print "                          ---------"
except:
    print ">>> NOT OK <<< rline() failed"

The following output results:

--------------------------------------------------------------
            Captured GrADS output: Line interface
--------------------------------------------------------------
Config: v1.9.0-rc1 32-bit big-endian readline sdf/xdf hdf-sds netcdf lats athena printim

Grid Analysis and Display System (GrADS) Version 1.9.0-rc1
Copyright (c) 1988-2005 by Brian Doty and IGES
Center for Ocean-Land-Atmosphere Studies (COLA)
Institute for Global Environment and Society (IGES)
This program is distributed WITHOUT ANY WARRANTY
See file COPYRIGHT for more information.

Built Tue Jan  1 13:19:12 EST 2008 for powerpc-apple-darwin8.11.0

This version of GrADS has been configured with the following options:

  o This is a 32-bit BIG ENDIAN machine version.
  o Command line editing (readline) ENABLED.
  o CIRES/CDC (http://www.cdc.noaa.gov) SDF/XDF interface ENABLED.
    Use sdfopen or xdfopen to read both NetCDF and HDF-SDS files.
  o DTYPE netcdf and DTYPE hdfsds are ENABLED.
  o OPeNDAP (a.k.a. DODS) gridded data interface DISABLED.
  o OPeNDAP (a.k.a. DODS) station data interface DISABLED.
  o PCMDI (http://www-pcmdi.llnl.gov) LATS interface ENABLED.
    This version is configured to write GRIB and HDF-SDS files.
  o DAO (http://dao.gsfc.nasa.gov) Athena Widget GUI ENABLED.
  o NRL/DAO/PCMDI XA or ImageMagick Image Output DISABLED.
  o printim command for direct png/gif output ENABLED.
    (via the GD Library -- http://www.boutell.com/gd)

For additional information please consult http://grads.iges.org/grads/
                         ---------

Similarly, one could have used the word interface with similar results:

try:
    print ""
    print "--------------------------------------------------------------"
    print "            Captured GrADS output: Word interface"
    print "--------------------------------------------------------------"
    for i in range(1,ga.nLines):
        for j in range(1,20):     # 20 is an over estimate, but i is OK
            stdout.write(ga.rword(i,j)+' ')
        stdout.write('\n')
    print "                          ---------"
 except:
    print ">>> NOT OK <<< rword() failed"

Quiting the GrADS application[edit]

The simplest way is to destroy the ga object

del ga