subroutine GridReadCoords (weightFile, src_lat, src_lon, src_area, &
src_mask, src_frac, &
dst_lat, dst_lon, dst_area, dst_mask, &
dst_frac, rc)
character(len=*), intent(in) :: weightFile
real(ESMF_KIND_R8), pointer :: src_lat(:)
real(ESMF_KIND_R8), pointer :: src_lon(:)
real(ESMF_KIND_R8), pointer :: src_area(:)
real(ESMF_KIND_R8), pointer :: src_mask(:)
real(ESMF_KIND_R8), pointer :: src_frac(:)
real(ESMF_KIND_R8), pointer :: dst_lat(:)
real(ESMF_KIND_R8), pointer :: dst_lon(:)
real(ESMF_KIND_R8), pointer :: dst_area(:)
real(ESMF_KIND_R8), pointer :: dst_mask(:)
real(ESMF_KIND_R8), pointer :: dst_frac(:)
integer, intent(out), optional :: rc
integer :: ncstat, nc_file_id
integer :: nc_srcgridlat_id, nc_srcgridlon_id, &
nc_dstgridlat_id, nc_dstgridlon_id, &
nc_srcarea_id, nc_dstarea_id, &
nc_srcmask_id, nc_dstmask_id, &
nc_srcfrac_id, nc_dstfrac_id
integer :: unitsLen
character(ESMF_MAXPATHLEN) :: units, buffer
character(ESMF_MAXPATHLEN) :: msg
logical :: cart
real(ESMF_KIND_R8), parameter :: d2r = 3.141592653589793238/180
#ifdef ESMF_NETCDF
!-----------------------------------------------------------------
! open netcdf file
!-----------------------------------------------------------------
ncstat = nf90_open(weightFile, NF90_NOWRITE, nc_file_id)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_open: '//weightFile, &
ESMF_CONTEXT, rcToReturn=rc)) return
!-----------------------------------------------------------------
! get the grid coordinates
!-----------------------------------------------------------------
cart = .false.
ncstat = nf90_inq_varid(nc_file_id, 'yc_a', nc_srcgridlat_id)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_inq_varid', &
ESMF_CONTEXT, rcToReturn=rc)) return
ncstat = nf90_get_var(ncid=nc_file_id, varid=nc_srcgridlat_id, &
values=src_lat)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_get_var', &
ESMF_CONTEXT, rcToReturn=rc)) return
! get the units of the grid coordinates
ncstat = nf90_inquire_attribute(nc_file_id, nc_srcgridlat_id, 'units', &
len=unitsLen)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_inquire_attribute', &
ESMF_CONTEXT, rcToReturn=rc)) return
if(len(units) < unitsLen) then
print *, "Not enough space to get units."
return
endif
ncstat = nf90_get_att(nc_file_id, nc_srcgridlat_id, "units", buffer)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_get_att', &
ESMF_CONTEXT, rcToReturn=rc)) return
units = buffer(1:unitsLen)
! convert to radians if coordinates are in degrees
if (trim(units)==trim("degrees")) then
src_lat = src_lat*d2r
else if (units(1:1) .eq. 'm' .or. units(1:1) .eq. 'k') then
! Cartesian coordinate, convert into spherical
cart = .true.
else if (trim(units) /= "radians") then
write (msg, '(a,i4)') "- units are not 'degrees' or 'radians'"
call ESMF_LogSetError(ESMF_RC_OBJ_BAD, msg=msg, &
line=__LINE__, file=__FILE__ , rcToReturn=rc)
return
endif
ncstat = nf90_inq_varid(nc_file_id, 'xc_a', nc_srcgridlon_id)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_inq_varid', &
ESMF_CONTEXT, rcToReturn=rc)) return
ncstat = nf90_get_var(ncid=nc_file_id, varid=nc_srcgridlon_id, &
values=src_lon)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_get_var', &
ESMF_CONTEXT, rcToReturn=rc)) return
! get the units of the grid coordinates
ncstat = nf90_inquire_attribute(nc_file_id, nc_srcgridlon_id, 'units', &
len=unitsLen)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_inquire_attribute', &
ESMF_CONTEXT, rcToReturn=rc)) return
if(len(units) < unitsLen) then
print *, "Not enough space to get units."
return
endif
ncstat = nf90_get_att(nc_file_id, nc_srcgridlon_id, "units", buffer)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_get_att', &
ESMF_CONTEXT, rcToReturn=rc)) return
units = buffer(1:unitsLen)
! convert to radians if coordinates are in degrees
if ((trim(units) == "degrees" .or. trim(units) == "radians") .and. &
cart) then
write (msg, '(a,i4)') "- source grid lat/lon units are inconsistent"
call ESMF_LogSetError(ESMF_RC_OBJ_BAD, msg=msg, &
line=__LINE__, file=__FILE__ , rcToReturn=rc)
return
endif
if (trim(units)==trim("degrees")) then
src_lon = src_lon*d2r
else if (units(1:1) .eq. 'm' .or. units(1:1) .eq. 'k') then
! Cartesian coordinate, convert into spherical
if (.not. cart) then
write (msg, '(a,i4)') "- source grid lat/lon units are inconsistent"
call ESMF_LogSetError(ESMF_RC_OBJ_BAD, msg=msg, &
line=__LINE__, file=__FILE__ , rcToReturn=rc)
return
endif
else if (trim(units)/=trim("radians")) then
write (msg, '(a,i4)') "- units are not 'degrees' or 'radians'"
call ESMF_LogSetError(ESMF_RC_OBJ_BAD, msg=msg, &
line=__LINE__, file=__FILE__ , rcToReturn=rc)
return
endif
! Convert Cartesian to Spherical
if (cart) then
call convertCart2Sph(src_lon, src_lat,units)
endif
cart = .false.
!-----------------------------------------------------------------
! get the destination grid coordinates
!-----------------------------------------------------------------
ncstat = nf90_inq_varid(nc_file_id, 'yc_b', nc_dstgridlat_id)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_inq_varid', &
ESMF_CONTEXT, rcToReturn=rc)) return
ncstat = nf90_get_var(ncid=nc_file_id, varid=nc_dstgridlat_id, &
values=dst_lat)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_get_var', &
ESMF_CONTEXT, rcToReturn=rc)) return
! get the units of the grid coordinates
ncstat = nf90_inquire_attribute(nc_file_id, nc_dstgridlat_id, 'units', &
len=unitsLen)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_inquire_attribute', &
ESMF_CONTEXT, rcToReturn=rc)) return
if(len(units) < unitsLen) then
print *, "Not enough space to get units."
return
endif
ncstat = nf90_get_att(nc_file_id, nc_dstgridlat_id, "units", buffer)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_get_att', &
ESMF_CONTEXT, rcToReturn=rc)) return
units = buffer(1:unitsLen)
! convert to radians if coordinates are in degrees
if (trim(units)==trim("degrees")) then
dst_lat = dst_lat*d2r
else if (units(1:1) .eq. 'm' .or. units(1:1) .eq. 'k') then
! Cartesian coordinate, convert into spherical
cart = .true.
else if (trim(units)/=trim("radians")) then
write (msg, '(a,i4)') "- units are not 'degrees' or 'radians'"
call ESMF_LogSetError(ESMF_RC_OBJ_BAD, msg=msg, &
line=__LINE__, file=__FILE__ , rcToReturn=rc)
return
endif
ncstat = nf90_inq_varid(nc_file_id, 'xc_b', nc_dstgridlon_id)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_inq_varid', &
ESMF_CONTEXT, rcToReturn=rc)) return
ncstat = nf90_get_var(ncid=nc_file_id, varid=nc_dstgridlon_id, &
values=dst_lon)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_get_var', &
ESMF_CONTEXT, rcToReturn=rc)) return
! get the units of the grid coordinates
ncstat = nf90_inquire_attribute(nc_file_id, nc_dstgridlon_id, 'units', &
len=unitsLen)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_inquire_attribute', &
ESMF_CONTEXT, rcToReturn=rc)) return
if(len(units) < unitsLen) then
print *, "Not enough space to get units."
return
endif
ncstat = nf90_get_att(nc_file_id, nc_dstgridlon_id, "units", buffer)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_get_att', &
ESMF_CONTEXT, rcToReturn=rc)) return
units = buffer(1:unitsLen)
if ((trim(units) == "degrees" .or. trim(units) == "radians") .and. &
cart) then
write (msg, '(a,i4)') "- source grid lat/lon units are inconsistent"
call ESMF_LogSetError(ESMF_RC_OBJ_BAD, msg=msg, &
line=__LINE__, file=__FILE__ , rcToReturn=rc)
return
endif
! convert to radians if coordinates are in degrees
if (trim(units)==trim("degrees")) then
dst_lon = dst_lon*d2r
else if (units(1:1) .eq. 'm' .or. units(1:1) .eq. 'k') then
! Cartesian coordinate, convert into spherical
if (.not. cart) then
write (msg, '(a,i4)') "- source grid lat/lon units are inconsistent"
call ESMF_LogSetError(ESMF_RC_OBJ_BAD, msg=msg, &
line=__LINE__, file=__FILE__ , rcToReturn=rc)
return
endif
else if (trim(units)/=trim("radians")) then
write (msg, '(a,i4)') "- units are not 'degrees' or 'radians'"
call ESMF_LogSetError(ESMF_RC_OBJ_BAD, msg=msg, &
line=__LINE__, file=__FILE__ , rcToReturn=rc)
return
endif
! Convert Cartesian to Spherical
if (cart) then
call convertCart2Sph(dst_lon, dst_lat, units)
endif
!-----------------------------------------------------------------
! get the grid areas
!-----------------------------------------------------------------
ncstat = nf90_inq_varid(nc_file_id, 'area_a', nc_srcarea_id)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_inq_varid', &
ESMF_CONTEXT, rcToReturn=rc)) return
ncstat = nf90_get_var(ncid=nc_file_id, varid=nc_srcarea_id, &
values=src_area)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_get_var', &
ESMF_CONTEXT, rcToReturn=rc)) return
ncstat = nf90_inq_varid(nc_file_id, 'area_b', nc_dstarea_id)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_inq_varid error', &
ESMF_CONTEXT, rcToReturn=rc)) return
ncstat = nf90_get_var(ncid=nc_file_id, varid=nc_dstarea_id, &
values=dst_area)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_get_var', &
ESMF_CONTEXT, rcToReturn=rc)) return
!-----------------------------------------------------------------
! get the grid masks
!-----------------------------------------------------------------
ncstat = nf90_inq_varid(nc_file_id, 'mask_a', nc_srcmask_id)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_inq_varid', &
ESMF_CONTEXT, rcToReturn=rc)) return
ncstat = nf90_get_var(ncid=nc_file_id, varid=nc_srcmask_id, &
values=src_mask)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_get_var', &
ESMF_CONTEXT, rcToReturn=rc)) return
ncstat = nf90_inq_varid(nc_file_id, 'mask_b', nc_dstmask_id)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_inq_varid', &
ESMF_CONTEXT, rcToReturn=rc)) return
ncstat = nf90_get_var(ncid=nc_file_id, varid=nc_dstmask_id, &
values=dst_mask)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_get_var', &
ESMF_CONTEXT, rcToReturn=rc)) return
!-----------------------------------------------------------------
! get the grid fracs
!-----------------------------------------------------------------
ncstat = nf90_inq_varid(nc_file_id, 'frac_a', nc_srcfrac_id)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_inq_varid', &
ESMF_CONTEXT, rcToReturn=rc)) return
ncstat = nf90_get_var(ncid=nc_file_id, varid=nc_srcfrac_id, &
values=src_frac)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_get_var', &
ESMF_CONTEXT, rcToReturn=rc)) return
ncstat = nf90_inq_varid(nc_file_id, 'frac_b', nc_dstfrac_id)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_inq_varid', &
ESMF_CONTEXT, rcToReturn=rc)) return
ncstat = nf90_get_var(ncid=nc_file_id, varid=nc_dstfrac_id, &
values=dst_frac)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_get_var', &
ESMF_CONTEXT, rcToReturn=rc)) return
!------------------------------------------------------------------------
! close input file
!------------------------------------------------------------------------
ncstat = nf90_close(nc_file_id)
if (ESMF_LogFoundNetCDFError (ncstat, msg='nf90_close', &
ESMF_CONTEXT, rcToReturn=rc)) return
if(present(rc)) rc = ESMF_SUCCESS
#else
call ESMF_LogSetError(rcToCheck=ESMF_RC_LIB_NOT_PRESENT, &
msg="- ESMF_NETCDF not defined when lib was compiled", &
ESMF_CONTEXT, rcToReturn=rc)
return
#endif
end subroutine GridReadCoords