Recipe: Fortran program as OPeNDAP client

From GEOS-5
Jump to: navigation, search

Back to G5NR Data Access Guide.

Problem

By accessing the collection inst01hr_3d_T_Cv via the OPeNDAP server

http://opendap.nccs.nasa.gov/dods/OSSE/G5NR/Ganymed/7km

we want to read the surface temperature data inside the box bound by latitudes 25oN, 50oN and longitudes -130oW, -65oW for 2006/Sep/18, 9z and compute its min/max.

Solution

First, we ensure that our NetCDF-4 library has been built with Fortran and OPeNDAP support. If both the queries

> nc-config --has-f90
> nc-config --has-dap

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

The metadata for the collection inst01hr_3d_T_Cv is available at

http://opendap.nccs.nasa.gov/dods/OSSE/G5NR/Ganymed/7km/0.5000_deg/inst/inst01hr_3d_T_Cv.info

Code

This code accesses the collection inst01hr_3d_T_Cv from the OPeNDAP server, reads a subset of the temperature data (all levels inside the bounding box specified above) and computes its max/min.

NOTE:

  1. Instead of reading a downloaded NetCDF-4 file, we read an OPeNDAP URL.
  2. While in the downloaded file, the temperature variable appears in the uppercase (T), on the OPeNDAP server, this variable is in lowercase. This is reflected in the call to nf90_inq_varid.
  3. Via the OPeNDAP URL, we now have access to all times for which data exists. The hourly inst files are available starting at 2005/15/15, 2200z. Our desired time, 2006/09/18, 0900z is then the 11772th step.
 1 program g5nr_reader_dap
 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(:,:,:,:)
14 
15   ! Miscellaneous
16   ! -------------
17   integer :: ierr
18   integer :: im, jm, lm
19   integer :: ncid, varid
20   integer :: start(4), count(4)
21   integer :: imin, imax, jmin, jmax
22   real :: minlat, minlon, maxlat, maxlon
23 
24   ! We replace filename by opendal url
25   ! ----------------------------------
26   T_file = "http://opendap.nccs.nasa.gov:9090/dods/OSSE/G5NR/Ganymed/7km/0.5000_deg/inst/inst01hr_3d_T_Cv"
27 
28   ! file name and (subset) dimensions
29   ! bounding box:
30   !   lons = -130:0.5:-65
31   !   lats =   25:0.5:50
32   ! ---------------------------------
33   ! array indices from the
34   ! bounding lat/lon values
35   imin = nint((-130. + 180.)/0.5)
36   imax = nint(( -65. + 180.)/0.5)
37   jmin = nint((  25 +  90.)/0.5)
38   jmax = nint((  50 +  90.)/0.5)
39   ! array sizes
40   im = imax-imin+1
41   jm = jmax-jmin+1
42   lm = 72
43 
44   ! Open url and get var id
45   ! -----------------------
46   call check(nf90_open(T_file,NF90_NOWRITE,ncid), "opening T file")
47   call check(nf90_inq_varid(ncid,"t",varid), "getting T varid")
48 
49   ! Read temperature data
50   ! ---------------------
51   allocate(T(im,jm,lm,1))            ! global 4D array with 1 time level
52   start = [imin, jmin, 1, 11772]     ! time level 11772 corresponds to 2006/09/18, 9z
53   count = [im, jm, lm, 1]            ! 1 time level, 3D (lon,lat,lev) array
54   write(*,*) 'Reading T...'
55   call check(nf90_get_var(ncid,varid,T,start=start,count=count), "reading T")
56   write(*,*)'T: ', maxval(T),minval(T)
57 
58   ! close file, release memory
59   call check(nf90_close(ncid), "closing T file")
60   deallocate(T)
61 
62   ! All done
63   ! --------
64 
65 contains
66 
67   subroutine check(status, loc)
68 
69     integer, intent(in) :: status
70     character(len=*), intent(in) :: loc
71 
72     if(status /= NF90_NOERR) then
73        write (*,*) "Error at ", loc
74        write (*,*) NF90_STRERROR(status)
75     end if
76 
77   end subroutine check
78 
79 end program g5nr_reader_dap

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_dap.f90 can be compiled and linked to the NetCDF-4 library via

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

creating the executable g5nr_reader_dap.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_dap.x

produces the output

Reading T (subset)...
T (subset):    305.6048       191.6956

Discussion

See Also

  1. File Spec: File:G5NR-Ganymed-7km FileSpec-ON6-V1.0.pdf
  2. Recipe: Fortran program to read data from downloaded file

No Warranty

Copyright