Return to Ferret FAQ


Outputting NetCDF?


Question:

How do I make my model write out NetCDF files that are suitable for Ferret?

Solution:

The NetCDF utilities ncgen and ncdump do most of the work for you. They actually generate the FORTRAN or C code that you will need. All you have to do is to "hack it" a little bit.

You can use Ferret commands to describe what your output files should look like -- everything about the output files except the actual data and time word values.

We proceed in 4 steps:

  1. describe the file using Ferret commands
  2. use ncdump to create a "CDL" language description of the file
  3. use ncgen to generate FORTRAN or C code
  4. hack that code and insert it into our model code

1. Describe the file using Ferret commands

For this example let us suppose that we will have variables U and TEMP and that they are on staggered grids.

We will define variables like this with Ferret commands. The **values** of these variables will never be used ... but the **shape** of the variables is important. Our example doesn't use a Z axis but it could be added in a straightforward way.

The encoding of time is discussed at some length in the Ferret Users Guide

! define the axes and grids
! temperature axes
define axis/x=130e:80w:2/unit="degrees" xt
define axis/y=20s:20n:2/unit="degrees" yt
! velocity axes (staggered relative to temperature)
define axis/x=131e:79w:2/unit="degrees" xu
define axis/y=19s:21n:2/unit="degrees" yu
! time axis
define axis/t="1-jan-1980:12:00":"1-jan-1990:12:00":1/unit="days" time
! define the grids
define grid/x=xt/y=yt/t=time gt
define grid/x=xu/y=yu/t=time gu

! define the (dummy) variables
let u = SIN((x[g=gu]+y[g=gu]+t[g=gu])/100)
let temp = SIN((x[g=gt]+y[g=gt]+t[g=gt])/100)
set variable/title="Temperature"/units="centigrade" temp
set variable/title="Zonal Velocity"/units="m/sec" u

! create the template NetCDF file
save/file=template.cdf/l=1 temp,u
2. Now create a CDL description of this file

! ... note -- you can edit template.cdl with your favorite editor

> ncdump template.cdf > template.cdl
3. now create FORTRAN code that would generate this NetCDF file

( use the "-c" qualifier instead of "-f" to get C code)

! ... note: you can "hack" this FORTRAN code to create the NetCDF output
!	parts of your model

> ncgen -f template.cdl > template.F
4. Hack the code from template.F into your model code

Usually you will

  1. create a subroutine to create (but not insert data into) the output file. It is executed only once.
  2. insert the NetCDF into your model code where you want to write out the time word and the values of your variables. It is executed once per output time step of your model

  1. The subroutine to create the output file is the FORTRAN code from the line
    	include 'netcdf.inc'
    
    through
       * leave define mode
    	call ncendf(ncid, iret)
    

    You will also want the blocks of code like store XT that write your models coordinate system into the output file.
    Note: For maximum accuracy you may want to write the coordinates directly from your model arrays instead of from the DATA statement. For this make sure you insert the matching data type (e.g. "NCDOUBLE") in the define variables section of the file creation code.

  2. The code to write the time word is hacked from store TIME. Note that the "corner" value will be "1" for time step 1, 2 for step 2, etc. As with the other coordinates make the data type of your model time word variable matches the data type specified during define variables

    The code will look something like

          corner(1) = istep		! the current output time step
          edges(1) = 1
          call ncvpt(ncid, TIMEid, corner, edges, time_word, iret)
    

    The code to write the model variables is similar:

    
    ! corner is the start point for writing -- normally "1" on each axis
    ! edges is the amount to write for each axis -- normally the full axis length
    
          corner(1) = 1
          corner(2) = 1
          corner(3) = istep
          edges(1) = 76
          edges(2) = 21
          edges(3) = 1
    
          call ncvpt(ncid, TEMPid, corner, edges, temp_var, iret)
    
    


Last modified: Dec 9, 1996