Recipe: Fortran program as OPeNDAP client
Back to G5NR Data Access Guide.
First, we ensure that out NetCDF-4 library has been built with Fortran, HDF5 and OPeNDAP support. All three queries
> nc-config --has-f90 > nc-config --has-hdf5 > nc-config --has-dap
should return yes.
Read downloaded nc4 file
Download the collection inst01hr_3d_T_Cv for 2006-Sep-18/0900z as described here. The downloaded file is c1440_NR.inst01hr_3d_T_Cv.20060918_0900z.nc4
. The program to read this file is
Code
program g5nr_reader
use netcdf ! for reading the NR files
implicit none
! File name
! ---------
character(len=256) :: T_file
! Global, 4D array: (lon,lat,lev,time)
! ------------------------------------
real, allocatable :: T(:,:,:,:)
! Miscellaneous
! -------------
integer :: ierr
integer :: im, jm, lm
integer :: ncid, varid
integer :: start(4), count(4)
! file name and dimensions
! ------------------------
im = 720
jm = 361
lm = 72
T_file = "c1440_NR.inst01hr_3d_T_Cv.20060918_0900z.nc4"
! Allocate the Global 4-D array with only 1 time level
! ----------------------------------------------------
allocate(T(im,jm,lm,1))
! Hypercube for reading one 1 array
! ---------------------------------
start = (/ 1, 1, 1, 1 /)
count = (/ im, jm, lm, 1 /) ! 1 time level, 3D (lon,lat,lev) array
! Read the data file
! ------------------
write(*,*)'Reading T'
call check(nf90_open(T_file,NF90_NOWRITE,ncid), "opening T file")
call check(nf90_inq_varid(ncid,"T",varid), "getting T varid")
call check(nf90_get_var(ncid,varid,T,start=start,count=count), "reading T")
call check(nf90_close(ncid), "closing T file")
! Print min/max of arrays
! -----------------------
write(*,*)'T: ', maxval(T),minval(T)
deallocate(T)
! All done
! --------
contains
subroutine check(status, loc)
integer, intent(in) :: status
character(len=*), intent(in) :: loc
if(status /= NF90_NOERR) then
write (*,*) "Error at ", loc
write (*,*) NF90_STRERROR(status)
end if
end subroutine check
end program g5nr_reader
Compile and link
We use the utility nf-config
(included in the NetCDF-4 installation) to identify the linker rules
For a typical NetCDF-4 installation, the above code, g5nr_reader.f90
can be compiled and linked to the NetCDF-4 library via
gfortran -o g5nr_reader.x `nf-config --fflags` g5nr_reader.f90 `nf-config --flibs`
creating the executable g5nr_reader.x
.
NOTE:
- You can use your favorite Fortran compiler instead of
gfortran
. - 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.x
produces the output
Reading T T: 315.6512 180.3667
Access data via OPeNDAP server
Code
The code to read data from the OPeNDAP server is a minor modification (3 lines) of the above program and is reproduced below
program g5nr_reader_dap
use netcdf ! for reading the NR files
implicit none
! File name
! ---------
character(len=256) :: T_file
! Global, 4D array: (lon,lat,lev,time)
! ------------------------------------
real, allocatable :: T(:,:,:,:)
! Miscellaneous
! -------------
integer :: ierr
integer :: im, jm, lm
integer :: ncid, varid
integer :: start(4), count(4)
! file name and dimensions
! ------------------------
im = 720
jm = 361
lm = 72
T_file = "http://opendap.nccs.nasa.gov:9090/dods/OSSE/GEOS-5.12/BETA9/0.5000_deg/inst/inst01hr_3d_T_Cv"
! Allocate the Global 4-D array with only 1 time level
! ----------------------------------------------------
allocate(T(im,jm,lm,1))
! Hypercube for reading one 1 array
! ---------------------------------
start = (/ 1, 1, 1, 11772 /) ! time level 11748 corresponds to 2006/09/18, 9z
count = (/ im, jm, lm, 1 /) ! 1 time level, 3D (lon,lat,lev) array
! Read the data file
! ------------------
write(*,*)'Reading T'
call check(nf90_open(T_file,NF90_NOWRITE,ncid), "opening T file")
call check(nf90_inq_varid(ncid,"t",varid), "getting T varid")
call check(nf90_get_var(ncid,varid,T,start=start,count=count), "reading T")
call check(nf90_close(ncid), "closing T file")
! Print min/max of arrays
! -----------------------
write(*,*)'T: ', maxval(T),minval(T)
deallocate(T)
! All done
! --------
contains
subroutine check(status, loc)
integer, intent(in) :: status
character(len=*), intent(in) :: loc
if(status /= NF90_NOERR) then
write (*,*) "Error at ", loc
write (*,*) NF90_STRERROR(status)
end if
end subroutine check
end program g5nr_reader_dap