Recipe: Fortran program to read data from downloaded file

From GEOS-5
Jump to: navigation, search

Back to G5NR Data Access Guide.

Problem

We want to read a downloaded data file using Fortran.

Solution

For the purpose of this example, we assume that we have already downloaded the file c1440_NR.inst01hr_3d_T_Cv.20060918_0900z.nc4 from the ftp server. For more information about file naming conventions, and how to download a file from the ftp server, please follow the links in the #See Also section.

Here we use the NetCDF-4 library to read the downloaded file (since NetCDF-4 uses HDF-5 as the underlying data format, HDF-5 library can also be used directly). First, we ensure that out NetCDF-4 library has been built with Fortran and HDF5 support. If both the queries

> nc-config --has-f90
> nc-config --has-hdf5

return "yes", we have a compatible NetCDF-4 library. Here nc-config is a utility bundled with NetCDF-4 package.

Code

This code reads the global temperature data from file and computes the maximum and minimum temperatures.

 1 program g5nr_reader_global
 2 
 3   use netcdf           ! for reading the NR files
 4 
 5   implicit none
 6 
 7   ! File name
 8   ! ---------
 9   character(len=256) :: T_file
10 
11   ! 4D array: (lon,lat,lev,time)
12   ! ----------------------------
13   real, allocatable :: T(:,:,:,:)    ! global data
14 
15   ! Miscellaneous
16   ! -------------
17   integer :: ierr
18   integer :: im, jm, lm
19   integer :: ncid, varid
20   integer :: start(4), count(4)
21 
22   ! file name and dimensions
23   ! ------------------------
24   T_file = "c1440_NR.inst01hr_3d_T_Cv.20060918_0900z.nc4"
25   im = 720
26   jm = 361
27   lm = 72
28 
29   ! Open file and get var id
30   ! ------------------------
31   call check(nf90_open(T_file,NF90_NOWRITE,ncid), "opening T file")
32   call check(nf90_inq_varid(ncid,"T",varid), "getting T varid")
33 
34   ! Read global temperature data
35   ! ----------------------------
36   allocate(T(im,jm,lm,1))            ! global 4D array with 1 time level
37   start = [1, 1, 1, 1]
38   count = [im, jm, lm, 1]            ! 1 time level, 3D (lon,lat,lev) array
39   write(*,*) 'Reading global T...'
40   call check(nf90_get_var(ncid,varid,T,start=start,count=count), "reading T")
41   write(*,*)'T: ', maxval(T),minval(T)
42 
43   ! close file, release memory
44   call check(nf90_close(ncid), "closing T file")
45   deallocate(T)
46 
47   ! All done
48   ! --------
49 
50 contains
51 
52   subroutine check(status, loc)
53 
54     integer, intent(in) :: status
55     character(len=*), intent(in) :: loc
56 
57     if(status /= NF90_NOERR) then
58        write (*,*) "Error at ", loc
59        write (*,*) NF90_STRERROR(status)
60     end if
61 
62   end subroutine check
63 
64 end program g5nr_reader_global
Compile and link

We use the utility nf-config (included in the NetCDF-4 installation) to identify the linking rules.

For a typical NetCDF-4 installation, the above code, g5nr_reader_global.f90 can be compiled and linked to the NetCDF-4 library via

gfortran -o g5nr_reader_global.x `nf-config --fflags` g5nr_reader_global.f90 `nf-config --flibs`

creating the executable g5nr_reader_global.x.

NOTE:

  1. You can use your favorite Fortran compiler instead of gfortran.
  2. If the NetCDF-4 library was built with parallel I/O support, you will need to use mpif90 to link, even if your code does not use the MPI library.
Run

Running the executable

> ./g5nr_reader_global.x

produces the output

Reading global T...
T:    315.6512       180.3667

Discussions

Modifications to read a subset of the data

The above Fortran code can be easily modified to read a subset of the temperature data instead. If we want to read temperature inside the box bounded by latitudes 25oN, 50oN and longitudes -130oW, -65oW, we need to modify the values of the dimensions im, jm and the argument start of the call to nf90_get_var.

Our new dimensions are

! array indices from boundary lat/lon
imin = nint((-130. + 180.)/0.5)
imax = nint((- 65. + 180.)/0.5)
jmin = nint((  25. +  90.)/0.5)
jmax = nint((  50. +  90.)/0.5)

! array sizes
im = imax-imin+1
jm = jmax-jmin+1

The new start argument to nf90_get_var defining the first indices of the slice of array to be read, is

start = [imin, jmin, 1, 1]

The modified version of the code is g5nr_reader_subset.txt.

Running this modified version

./g5nr_reader_subset.x

produces the output

Reading Tsub
Tsub:    305.6048       191.6956

See Also

  1. File Spec: File:G5NR-Ganymed-7km FileSpec-ON6-V1.0.pdf
  2. Recipe: Retrieve (global) data from FTP server
  3. Recipe: Fortran program as OPeNDAP client

No Warranty

Copyright