MODULE ReadNETCDF IMPLICIT NONE ! This is part of the netCDF package. ! Copyright 2006 University Corporation for Atmospheric Research/Unidata. ! See COPYRIGHT file for conditions of use. ! This is an example which reads some 4D pressure and ! temperatures. The data file read by this program is produced by ! the companion program pres_temp_4D_wr.f90. It is intended to ! illustrate the use of the netCDF Fortran 90 API. ! This program is part of the netCDF tutorial: ! http://www.unidata.ucar.edu/software/netcdf/docs/netcdf-tutorial ! Full documentation of the netCDF Fortran 90 API can be found at: ! http://www.unidata.ucar.edu/software/netcdf/docs/netcdf-f90 ! $Id: pres_temp_4D_rd.f90,v 1.9 2010/04/06 19:32:09 ed Exp $ PRIVATE ! We are reading 4D data, a 12 x 6 x 2 lon-lat-lvl grid, with 2 ! timesteps of data. integer, parameter :: nMaxFile = 6 integer, parameter :: nMaxFile2 = 7 integer, parameter :: nvar = 4 integer, parameter :: NDIMS = 3 integer,public, parameter :: NTIME = 73700 integer, parameter :: NLATS = 96 integer, parameter :: NLONS = 192 integer :: lon_dimid integer :: lat_dimid integer :: tim_dimid !global attribute CHARACTER (len = *), parameter :: history = "history" CHARACTER (len = *), parameter :: institution = "CPTEC" CHARACTER (len = *), parameter :: source = "source" CHARACTER (len = *), parameter :: Conventions = "Conventions" INTEGER, PARAMETER :: MAX_ATT_LEN = 800 CHARACTER (len = MAX_ATT_LEN) :: gbl_history != "NASA LaRC/GEWEX SRB" ; CHARACTER (len = MAX_ATT_LEN) :: gbl_Conventions != "CF-1.1" ; ! In addition to the latitude and longitude dimensions, we will also ! create latitude and longitude variables which will hold the actual ! latitudes and longitudes. Since they hold data about the ! coordinate system, the netCDF term for these is: "coordinate ! variables." REAL(KIND=8) :: lats ( NLATS) REAL(KIND=8) :: lons ( NLONS) REAL(KIND=8) :: timer ( NTIME) ! We recommend that each variable carry a "units" attribute. TYPE :: attribute character (len = 5) :: UNITS ="units" character (len = 10) :: FillValue ="_FillValue" character (len = 10) :: add_offset ="add_offset" character (len = 12) :: scale_factor ="scale_factor" character (len = 13) :: standard_name ="standard_name" character (len = 9) :: long_name ="long_name" character (len = 12) :: actual_range ="actual_range" character (len = 7) :: delta_t ="delta_t" character(LEN=400) :: var_units != "W m-2" ; INTEGER(KIND=2) :: var_FillValue != -100s ; REAL :: lat_actual_range(2)! = -89.5f, 89.5f ; REAL :: lon_actual_range(2)! = -89.5f, 89.5f ; REAL :: var_add_offset != 0.f ; REAL :: var_scale_factor != 0.1f ; character(LEN=400) :: var_standard_name != "toa_incoming_shortwave_flux" ; character(LEN=400) :: var_long_name != "All-Sky Shortwave TOA Downward Flux" ; character(LEN=400) :: time_delta_t != "CF-1.1" ; REAL(KIND=8) :: time_actual_range(NTIME) ! = -89.5f, 89.5f ; REAL(KIND=8) :: time_actual_range_out(3*NTIME) ! = -89.5f, 89.5f ; END TYPE attribute TYPE(attribute), ALLOCATABLE :: var(:,:) ! The start and count arrays will tell the netCDF library where to ! read our data. integer :: start (NDIMS) integer :: start2 (NDIMS) integer :: start_out (NDIMS) INTEGER :: count (NDIMS) integer :: dimids (NDIMS) REAL(KIND=8), DIMENSION(:,:), ALLOCATABLE :: var_c REAL(KIND=8), DIMENSION(:,:), ALLOCATABLE :: var_p ! We will read surface temperature and pressure fields. In netCDF ! terminology these are called "variables." CHARACTER (len = 20),PUBLIC, PARAMETER :: VAR_NAME(1:nvar,1:nMaxFile)=RESHAPE((/& "lon ","lat ","time ","precip ",& "lon ","lat ","time ","temp ",& "lon ","lat ","time ","rhum ",& "lon ","lat ","time ","windspeed ",& "lon ","lat ","time ","swrad ",& "lon ","lat ","time ","lwrad "/),(/nvar,nMaxFile/)) CHARACTER (len = 20),PUBLIC, PARAMETER :: VAR_NAME2(1:nvar,1:nMaxFile2)=RESHAPE((/& "lon ","lat ","time ","pres ",& "lon ","lat ","time ","tems ",& "lon ","lat ","time ","umes ",& "lon ","lat ","time ","wind ",& "lon ","lat ","time ","prec ",& "lon ","lat ","time ","ocis ",& "lon ","lat ","time ","olis "/),(/nvar,nMaxFile2/)) INTEGER,PUBLIC :: ncid_pres INTEGER,PUBLIC :: ncid_prec INTEGER,PUBLIC :: ncid_temp INTEGER,PUBLIC :: ncid_rhum INTEGER,PUBLIC :: ncid_umes INTEGER,PUBLIC :: ncid_wind INTEGER,PUBLIC :: ncid_swrd INTEGER,PUBLIC :: ncid_lwrd integer,PUBLIC, ALLOCATABLE :: var_varid (:,:) LOGICAL :: ALLOC_INIT=.TRUE. LOGICAL :: DEALLOC_FINALIZE=.TRUE. PUBLIC :: InitReadNETCDF PUBLIC :: FinalizeReadNETCDF PUBLIC :: RunReadNETCDF CONTAINS SUBROUTINE InitReadNETCDF(FILE_NAME,icn_data,nFile) use netcdf IMPLICIT NONE ! This is the name of the data file we will read. CHARACTER (len = *), INTENT(IN ) :: FILE_NAME INTEGER , INTENT(IN ) :: nFile CHARACTER (len = *), INTENT(IN ) :: icn_data ! We will learn about the data file and store results in these ! program variables. INTEGER :: ndims_in INTEGER :: nvars_in INTEGER :: ngatts_in INTEGER :: unlimdimid_in INTEGER :: ivar INTEGER :: ncid ! Allocate memory. IF(ALLOC_INIT)THEN IF(TRIM(icn_data) == 'ERA5')THEN ALLOCATE(var (1:nvar,1:nMaxFile)) ALLOCATE(var_varid (1:nvar,1:nMaxFile)) ELSE IF(TRIM(icn_data) == 'GLDAS')THEN ALLOCATE(var (1:nvar,1:nMaxFile2)) ALLOCATE(var_varid (1:nvar,1:nMaxFile2)) ELSE PRINT*,'ERROR',TRIM(icn_data),'diff','ERA5','GLDAS' STOP END IF allocate(var_c (NLONS, NLATS)) allocate(var_p (NLONS, NLATS)) ALLOC_INIT=.FALSE. END IF !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Open the input file. CALL check( nf90_open (TRIM(FILE_NAME) , nf90_nowrite, ncid) ) !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! ! ! There are a number of inquiry functions in netCDF which can be ! used to learn about an unknown netCDF file. NF90_INQ tells how many ! netCDF variables, dimensions, and global attributes are in the ! file; also the dimension id of the unlimited dimension, if there ! is one. call check( nf90_inquire(ncid, ndims_in, nvars_in, ngatts_in, unlimdimid_in)) ! !call check( nf90_get_att(ncid, nf90_global, TRIM(history) , gbl_history)) !call check( nf90_get_att(ncid, nf90_global, TRIM(Conventions), gbl_Conventions)) ! Get the varids of the latitude and longitude coordinate variables. IF(TRIM(icn_data) == 'ERA5')THEN call check( nf90_inq_varid(ncid, VAR_NAME(1,nFile) , var_varid(1,nFile)) ) call check( nf90_inq_varid(ncid, VAR_NAME(2,nFile) , var_varid(2,nFile)) ) call check( nf90_inq_varid(ncid, VAR_NAME(3,nFile) , var_varid(3,nFile)) ) ELSE IF(TRIM(icn_data) == 'GLDAS')THEN call check( nf90_inq_varid(ncid, VAR_NAME2(1,nFile) , var_varid(1,nFile)) ) call check( nf90_inq_varid(ncid, VAR_NAME2(2,nFile) , var_varid(2,nFile)) ) call check( nf90_inq_varid(ncid, VAR_NAME2(3,nFile) , var_varid(3,nFile)) ) ELSE PRINT*,'ERROR',TRIM(icn_data),'diff','ERA5','GLDAS' STOP END IF PRINT*,'pkubota ',TRIM(FILE_NAME),TRIM(icn_data) ! Define the dimensions. The record dimension is defined to have ! unlimited length - it can grow as needed. In this example it is ! the time dimension. ! Read the timer, latitude and longitude data. IF(TRIM(icn_data) == 'ERA5')THEN call check( nf90_get_var(ncid, var_varid(1,nFile), lons) ) call check( nf90_get_var(ncid, var_varid(2,nFile), lats) ) ! call check( nf90_get_var(ncid, var_varid(3,nFile),timer) ) ELSE IF(TRIM(icn_data) == 'GLDAS')THEN call check( nf90_get_var(ncid, var_varid(1,nFile), lons) ) call check( nf90_get_var(ncid, var_varid(2,nFile), lats) ) ! call check( nf90_get_var(ncid, var_varid(3,nFile),timer) ) ELSE PRINT*,'ERROR',TRIM(icn_data),'diff','ERA5','GLDAS' STOP END IF PRINT*,'pkubota0 ',TRIM(FILE_NAME),TRIM(icn_data) ! Define the coordinate variables. We will only define coordinate ! variables for lat and lon. Ordinarily we would need to provide ! an array of dimension IDs for each variable's dimensions, but ! since coordinate variables only have one dimension, we can ! simply provide the address of that dimension ID (lat_dimid) and ! similarly for (lon_dimid). ! Assign units attributes to the netCDF variables. IF(TRIM(icn_data) == 'ERA5')THEN call check( nf90_get_att(ncid ,var_varid(1,nFile), var(1,nFile)%UNITS , var(1,nFile)%var_units ) ) call check( nf90_get_att(ncid ,var_varid(1,nFile), var(1,nFile)%long_name , var(1,nFile)%var_long_name ) ) call check( nf90_get_att(ncid ,var_varid(2,nFile), var(2,nFile)%UNITS , var(2,nFile)%var_units ) ) call check( nf90_get_att(ncid ,var_varid(2,nFile), var(2,nFile)%long_name , var(2,nFile)%var_long_name ) ) call check( nf90_get_att(ncid ,var_varid(3,nFile), var(3,nFile)%UNITS , var(3,nFile)%var_units ) ) call check( nf90_get_att(ncid ,var_varid(3,nFile), var(3,nFile)%long_name , var(3,nFile)%var_long_name ) ) ELSE IF(TRIM(icn_data) == 'GLDAS')THEN call check( nf90_get_att(ncid ,var_varid(1,nFile), var(1,nFile)%UNITS , var(1,nFile)%var_units ) ) call check( nf90_get_att(ncid ,var_varid(1,nFile), var(1,nFile)%long_name , var(1,nFile)%var_long_name ) ) call check( nf90_get_att(ncid ,var_varid(2,nFile), var(2,nFile)%UNITS , var(2,nFile)%var_units ) ) call check( nf90_get_att(ncid ,var_varid(2,nFile), var(2,nFile)%long_name , var(2,nFile)%var_long_name ) ) call check( nf90_get_att(ncid ,var_varid(3,nFile), var(3,nFile)%UNITS , var(3,nFile)%var_units ) ) call check( nf90_get_att(ncid ,var_varid(3,nFile), var(3,nFile)%long_name , var(3,nFile)%var_long_name ) ) ELSE PRINT*,'ERROR',TRIM(icn_data),'diff','ERA5','GLDAS' STOP END IF PRINT*,'pkubota1 ',TRIM(FILE_NAME),TRIM(icn_data) ! The dimids array is used to pass the dimids of the dimensions of ! the netCDF variables. Both of the netCDF variables we are creating ! share the same four dimensions. In Fortran, the unlimited ! dimension must come last on the list of dimids. dimids = (/ lon_dimid, lat_dimid, tim_dimid /) DO ivar=4,nvar ! ! Get the varids of the pressure and temperature netCDF variables. ! IF(TRIM(icn_data) == 'ERA5')THEN call check( nf90_inq_varid(ncid, VAR_NAME(ivar,nFile) , var_varid(ivar,nFile)) ) ELSE IF(TRIM(icn_data) == 'GLDAS')THEN call check( nf90_inq_varid(ncid, VAR_NAME2(ivar,nFile) , var_varid(ivar,nFile)) ) ELSE PRINT*,'ERROR',TRIM(icn_data),'diff','ERA5','GLDAS' STOP END IF ! ! Each of the netCDF variables has a "units" attribute. Let's read ! them and check them. ! ! call check( nf90_get_att(ncid, var_varid(ivar,nFile), var(ivar,nFile)%UNITS , var(ivar,nFile)%var_units ) ) ! call check( nf90_get_att(ncid, var_varid(ivar,nFile), var(ivar,nFile)%FillValue , var(ivar,nFile)%var_FillValue ) ) !call check( nf90_get_att(ncid, var_varid(ivar,nFile), var(ivar,nFile)%add_offset , var(ivar,nFile)%var_add_offset ) ) !call check( nf90_get_att(ncid, var_varid(ivar,nFile), var(ivar,nFile)%scale_factor , var(ivar,nFile)%var_scale_factor) ) ! call check( nf90_get_att(ncid, var_varid(ivar,nFile), var(ivar,nFile)%long_name , var(ivar,nFile)%var_long_name ) ) END DO PRINT*,'pkubota2 ',TRIM(FILE_NAME),TRIM(icn_data) IF(TRIM(icn_data) == 'ERA5')THEN IF (nFile ==1) THEN ncid_prec=ncid ELSE IF(nFile ==2)THEN ncid_temp=ncid ELSE IF(nFile ==3) THEN ncid_rhum=ncid ELSE IF(nFile ==4)THEN ncid_wind=ncid ELSE IF(nFile ==5)THEN ncid_swrd=ncid ELSE IF(nFile ==6)THEN ncid_lwrd=ncid END IF ELSE IF(TRIM(icn_data) == 'GLDAS')THEN IF (nFile ==1) THEN ncid_pres=ncid ELSE IF(nFile ==2)THEN ncid_temp=ncid ELSE IF(nFile ==3) THEN ncid_umes=ncid ELSE IF(nFile ==4)THEN ncid_wind=ncid ELSE IF(nFile ==5)THEN ncid_prec=ncid ELSE IF(nFile ==6)THEN ncid_swrd=ncid ELSE IF(nFile ==7)THEN ncid_lwrd=ncid END IF ELSE PRINT*,'ERROR',TRIM(icn_data),'diff','ERA5','GLDAS' STOP END IF END SUBROUTINE InitReadNETCDF SUBROUTINE RunReadNETCDF(rec,ncid,latic,latip,longc,longp,var_stepc,var_stepp,nFile) use netcdf IMPLICIT NONE INTEGER , INTENT(IN ) :: rec INTEGER , INTENT(IN ) :: nFile INTEGER , INTENT(IN ) :: ncid REAL(KIND=4) , INTENT(OUT ) :: latic( NLONS, NLATS) REAL(KIND=4) , INTENT(OUT ) :: latip( NLONS, NLATS) REAL(KIND=4) , INTENT(OUT ) :: longc( NLONS, NLATS) REAL(KIND=4) , INTENT(OUT ) :: longp( NLONS, NLATS) REAL(KIND=8) , INTENT(OUT ) :: var_stepc( NLONS, NLATS) REAL(KIND=8) , INTENT(OUT ) :: var_stepp( NLONS, NLATS) INTEGER ::ivar,i,j ! ! Read 1 record of NLONS*NLATS values, starting at the beginning ! of the record (the (1, 1, rec) element in the netCDF file). ! count = (/ NLONS, NLATS, 1 /) start = (/ 1, 1, 1 /) start2 = (/ 1, 1, 1 /) start_out = (/ 1, 1, 1 /) ! ! Read the surface pressure and temperature data from the file, one ! record at a time. ! start(3) = rec IF(rec==NTIME)THEN start2(3) = rec ELSE start2(3) = rec+1 END IF DO j=1,NLATS DO i=1,NLONS longp(i,j) =lons(i) latip(i,j) =lats(j) longc(i,j) = longp(i,j) latic(i,j) = latip(i,j) END DO END DO DO ivar=4,nvar call check( nf90_get_var(ncid, var_varid(ivar,nFile), var_c(1:NLONS, 1:NLATS), start = start , & count = count) ) call check( nf90_get_var(ncid, var_varid(ivar,nFile), var_p(1:NLONS, 1:NLATS), start = start2, & count = count) ) DO j=1,NLATS DO i=1,NLONS var_stepc(i,j) =var_c(i,j) var_stepp(i,j) =var_p(i,j) END DO END DO END DO END SUBROUTINE RunReadNETCDF SUBROUTINE FinalizeReadNETCDF(ncid,nFile) use netcdf IMPLICIT NONE INTEGER , INTENT(IN ) :: nFile INTEGER , INTENT(IN ) :: ncid IF(DEALLOC_FINALIZE)THEN DEALLOCATE(var) DEALLOCATE(var_varid) DEALLOCATE(var_c ) DEALLOCATE(var_p ) DEALLOC_FINALIZE=.FALSE. END IF ! Close the file. This frees up any internal netCDF resources ! associated with the file. IF(nFile==1) call check( nf90_close(ncid) ) END SUBROUTINE FinalizeReadNETCDF subroutine check(status) use netcdf integer, intent ( in) :: status if(status /= nf90_noerr) then stop 2 end if end subroutine check END MODULE ReadNETCDF