Recipe: Fortran program as OPeNDAP client: Difference between revisions

From GEOS-5
Jump to navigation Jump to search
Pchakrab (talk | contribs)
m Pchakrab moved page Fortran to G5NR data access using Fortran
Move from http to https
 
(105 intermediate revisions by one other user not shown)
Line 1: Line 1:
Back to [[G5NR Data Access Guide]].


<nowiki>
== Problem ==
program g5nr_reader


   use netcdf          ! for reading the NR files                                                                              
By accessing the collection <code>inst01hr_3d_T_Cv</code> via the OPeNDAP server
https://opendap.nccs.nasa.gov/dods/OSSE/G5NR/Ganymed/7km
 
we want to read the surface temperature data inside the box bound by latitudes 25<sup>o</sup>N, 50<sup>o</sup>N and longitudes -130<sup>o</sup>W, -65<sup>o</sup>W 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 <code>nc-config</code> is a utility bundled with NetCDF-4 package.
 
The metadata for the collection <code>inst01hr_3d_T_Cv</code> is available at
https://opendap.nccs.nasa.gov/dods/OSSE/G5NR/Ganymed/7km/0.5000_deg/inst/inst01hr_3d_T_Cv.info
 
==== Code ====
This code accesses the collection <code>inst01hr_3d_T_Cv</code> 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:
# Instead of reading a downloaded NetCDF-4 file, we read an OPeNDAP URL.
# 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 <code>nf90_inq_varid</code>.
# Via the OPeNDAP URL, we now have access to all times for which data exists. The hourly <code>inst</code> files are available starting at 2005/15/15, 2200z. Our desired time, 2006/09/18, 0900z is then the 11772<sup>th</sup> step.
 
<syntaxhighlight lang="fortran" line>
program g5nr_reader_dap
 
   use netcdf          ! for reading the NR files


   implicit none
   implicit none


   ! File name                                                                                                                
   ! File name
   ! ---------                                                                                                                
   ! ---------
   character(len=256) :: T_file
   character(len=256) :: T_file


   ! Global, 4D array: (lon,lat,lev,time)                                                                                    
   ! 4D array: (lon,lat,lev,time)
   ! ------------------------------------  
   ! ------------------------------------
   real, pointer :: T(:,:,:,:) => null()
   real, allocatable :: T(:,:,:,:)


   ! Miscellaneous                                                                                                            
   ! Miscellaneous
   ! -------------                                                                                                            
   ! -------------
   integer :: ierr
   integer :: ierr
   integer :: im, jm, lm
   integer :: im, jm, lm
   integer :: ncid, varid
   integer :: ncid, varid
   integer :: start(4), count(4)
   integer :: start(4), count(4)
  integer :: imin, imax, jmin, jmax
  real :: minlat, minlon, maxlat, maxlon


   ! For now hard code file name and dimensions                                                                               
   ! We replace filename by opendal url
   ! ------------------------------------------                                                                              
   ! ----------------------------------
   im = 720
  T_file = "https://opendap.nccs.nasa.gov:9090/dods/OSSE/G5NR/Ganymed/7km/0.5000_deg/inst/inst01hr_3d_T_Cv"
   jm = 361
 
  ! file name and (subset) dimensions
  ! bounding box:
  !  lons = -130:0.5:-65
  !  lats =  25:0.5:50
  ! ---------------------------------
  ! array indices from the
  ! bounding lat/lon values
  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
   lm = 72
   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                                                                     
   ! Open url and get var id
   ! ----------------------------------------------------                                                                     
   ! -----------------------
   allocate(T(im,jm,lm,1))
   call check(nf90_open(T_file,NF90_NOWRITE,ncid), "opening T file")
  call check(nf90_inq_varid(ncid,"t",varid), "getting T varid")


   ! Hypercube for reading one 1 array for a given time                                                                       
   ! Read temperature data
   ! -------------------------------------------------                                                                         
   ! ---------------------
   start = (/  1, 1, 1, 37 /)  ! time level 37                                                                               
  allocate(T(im,jm,lm,1))            ! global 4D array with 1 time level
   count = (/ im, jm, lm, 1 /)  ! 1 time level, 3D (lon,lat,lev) array                                                        
   start = [imin, jmin, 1, 11772]    ! time level 11772 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_get_var(ncid,varid,T,start=start,count=count), "reading T")
   write(*,*)'Reading T'
   write(*,*)'T: ', maxval(T),minval(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")


   ! Orint min/max of arrays                                                                                                   
   ! close file, release memory
   !  -----------------------                                                                                                   
   call check(nf90_close(ncid), "closing T file")
  write(*,*)'T: ', maxval(T),minval(T)
  deallocate(T)


   ! All done                                                                                                                  
   ! All done
   ! --------                                                                                                                  
   ! --------


contains
contains
Line 67: Line 105:
   end subroutine check
   end subroutine check


end program g5nr_reader
end program g5nr_reader_dap
</nowiki>
</syntaxhighlight>
 
==== Compile and link ====
 
We use the utility <code>nf-config</code> (included in the NetCDF-4 installation) to identify the linking rules
 
For a typical NetCDF-4 installation, the above code, <code>g5nr_reader_dap.f90</code> 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 <code>g5nr_reader_dap.x</code>.
 
NOTE:
# You can use your favorite Fortran compiler instead of <code>gfortran</code>.
# If the NetCDF-4 library was built with parallel I/O support, you will need to use <code>mpif90</code> 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 ==
 
# File Spec: [[File:G5NR-Ganymed-7km FileSpec-ON6-V1.0.pdf]]
# [[Recipe: Fortran program to read data from downloaded file]]
 
== No Warranty ==
 
== Copyright ==

Latest revision as of 10:24, 9 April 2019

Back to G5NR Data Access Guide.

Problem

By accessing the collection inst01hr_3d_T_Cv via the OPeNDAP server

https://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

https://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.
program g5nr_reader_dap

  use netcdf           ! for reading the NR files

  implicit none

  ! File name
  ! ---------
  character(len=256) :: T_file

  ! 4D array: (lon,lat,lev,time)
  ! ------------------------------------
  real, allocatable :: T(:,:,:,:)

  ! Miscellaneous
  ! -------------
  integer :: ierr
  integer :: im, jm, lm
  integer :: ncid, varid
  integer :: start(4), count(4)
  integer :: imin, imax, jmin, jmax
  real :: minlat, minlon, maxlat, maxlon

  ! We replace filename by opendal url
  ! ----------------------------------
  T_file = "https://opendap.nccs.nasa.gov:9090/dods/OSSE/G5NR/Ganymed/7km/0.5000_deg/inst/inst01hr_3d_T_Cv"

  ! file name and (subset) dimensions
  ! bounding box:
  !   lons = -130:0.5:-65
  !   lats =   25:0.5:50
  ! ---------------------------------
  ! array indices from the
  ! bounding lat/lon values
  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
  lm = 72

  ! Open url and get var id
  ! -----------------------
  call check(nf90_open(T_file,NF90_NOWRITE,ncid), "opening T file")
  call check(nf90_inq_varid(ncid,"t",varid), "getting T varid")

  ! Read temperature data
  ! ---------------------
  allocate(T(im,jm,lm,1))            ! global 4D array with 1 time level
  start = [imin, jmin, 1, 11772]     ! time level 11772 corresponds to 2006/09/18, 9z
  count = [im, jm, lm, 1]            ! 1 time level, 3D (lon,lat,lev) array
  write(*,*) 'Reading T...'
  call check(nf90_get_var(ncid,varid,T,start=start,count=count), "reading T")
  write(*,*)'T: ', maxval(T),minval(T)

  ! close file, release memory
  call check(nf90_close(ncid), "closing T file")
  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

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