From 435b3bcc8a9574aa7283a0768cc8ae145490772b Mon Sep 17 00:00:00 2001 From: Ben Andre Date: Wed, 8 Feb 2017 15:06:16 -0700 Subject: [PATCH] Check fates parameter dimensions when reading from file. Automatically check the number of dimensions and their names when reading fates parameters from the file. Compare the data on file against what is expected by the code. Test suite: ed - yellowstone gnu, intel, pgi Test baseline: a651a4f Test namelist changes: add fates_paramfile Test answer changes: bit for bit Test summary: all tests pass --- .../src/ED/main/FatesParametersInterface.F90 | 4 +- components/clm/src/main/paramUtilMod.F90 | 160 ++++++++++++++++++ .../src/utils/clmfates_paraminterfaceMod.F90 | 5 +- 3 files changed, 166 insertions(+), 3 deletions(-) diff --git a/components/clm/src/ED/main/FatesParametersInterface.F90 b/components/clm/src/ED/main/FatesParametersInterface.F90 index cfebf64f7f..5da9bd7eff 100644 --- a/components/clm/src/ED/main/FatesParametersInterface.F90 +++ b/components/clm/src/ED/main/FatesParametersInterface.F90 @@ -315,7 +315,7 @@ subroutine SetDimensionSizes(this, is_host_file, num_used_dimensions, dimension_ end subroutine SetDimensionSizes !----------------------------------------------------------------------- - subroutine GetMetaData(this, index, name, dimension_shape, dimension_sizes, is_host_param) + subroutine GetMetaData(this, index, name, dimension_shape, dimension_sizes, dimension_names, is_host_param) implicit none @@ -324,11 +324,13 @@ subroutine GetMetaData(this, index, name, dimension_shape, dimension_sizes, is_h character(len=param_string_length), intent(out) :: name integer, intent(out) :: dimension_shape integer, intent(out) :: dimension_sizes(max_dimensions) + character(len=param_string_length), intent(out) :: dimension_names(max_dimensions) logical, intent(out) :: is_host_param name = this%parameters(index)%name dimension_shape = this%parameters(index)%dimension_shape dimension_sizes = this%parameters(index)%dimension_sizes + dimension_names = this%parameters(index)%dimension_names is_host_param = this%parameters(index)%sync_with_host end subroutine GetMetaData diff --git a/components/clm/src/main/paramUtilMod.F90 b/components/clm/src/main/paramUtilMod.F90 index 75a85e3e6c..96c95440e7 100644 --- a/components/clm/src/main/paramUtilMod.F90 +++ b/components/clm/src/main/paramUtilMod.F90 @@ -11,14 +11,22 @@ module paramUtilMod module procedure readNcdioScalar module procedure readNcdioArray1d module procedure readNcdioArray2d + module procedure readNcdioScalarCheckDimensions + module procedure readNcdioArray1dCheckDimensions + module procedure readNcdioArray2dCheckDimensions end interface public :: readNcdioScalar public :: readNcdioArray1d public :: readNcdioArray2d + public :: readNcdioScalarCheckDimensions + public :: readNcdioArray1dCheckDimensions + public :: readNcdioArray2dCheckDimensions public :: readNcdio + private :: checkDimensions + contains !----------------------------------------------------------------------- ! @@ -128,4 +136,156 @@ subroutine readNcdioArray2d(ncid, varName, callingName, retVal) end subroutine readNcdioArray2d !----------------------------------------------------------------------- + !----------------------------------------------------------------------- + ! + !----------------------------------------------------------------------- + subroutine readNcdioScalarCheckDimensions(ncid, varName, expected_numDims, expected_dimNames, & + callingName, retVal) + ! + ! read the netcdf file...generic, could be used for any parameter read + ! + use abortutils , only : endrun + use ncdio_pio , only : file_desc_t + + implicit none + + ! arguments + type(file_desc_t),intent(inout) :: ncid ! pio netCDF file id + character(len=*), intent(in) :: varName ! variable we are reading + integer, intent(in) :: expected_numDims + character(len=*), intent(in) :: expected_dimNames(:) ! expected dimension name + character(len=*), intent(in) :: callingName ! calling routine + real(r8), intent(inout) :: retVal + + ! local vars + character(len=32) :: subname = 'readNcdio::' + character(len=100) :: errCode = ' - Error reading. Var: ' + + ! + ! netcdf read here + ! + call checkDimensions(ncid, varName, expected_numDims, expected_dimNames, subname) + call readNcdio(ncid, varName, callingName, retVal) + + end subroutine readNcdioScalarCheckDimensions + !----------------------------------------------------------------------- + + !----------------------------------------------------------------------- + ! + !----------------------------------------------------------------------- + subroutine readNcdioArray1dCheckDimensions(ncid, varName, expected_numDims, expected_dimNames, & + callingName, retVal) + ! + ! read the netcdf file...generic, could be used for any parameter read + ! + use abortutils , only : endrun + use ncdio_pio , only : file_desc_t + + implicit none + + ! arguments + type(file_desc_t),intent(inout) :: ncid ! pio netCDF file id + character(len=*), intent(in) :: varName ! variable we are reading + integer, intent(in) :: expected_numDims + character(len=*), intent(in) :: expected_dimNames(:) ! expected dimension name + character(len=*), intent(in) :: callingName ! calling routine + real(r8), intent(inout) :: retVal( 1: ) + + ! local vars + character(len=32) :: subname = 'readNcdio::' + character(len=100) :: errCode = ' - Error reading. Var: ' + ! + ! netcdf read here + ! + call checkDimensions(ncid, varName, expected_numDims, expected_dimNames, subname) + call readNcdio(ncid, varName, callingName, retVal) + + end subroutine readNcdioArray1dCheckDimensions + !----------------------------------------------------------------------- + + !----------------------------------------------------------------------- + ! + !----------------------------------------------------------------------- + subroutine readNcdioArray2dCheckDimensions(ncid, varName, expected_numDims, expected_dimNames, & + callingName, retVal) + ! + ! read the netcdf file...generic, could be used for any parameter read + ! + use abortutils , only : endrun + use ncdio_pio , only : file_desc_t + + implicit none + + ! arguments + type(file_desc_t),intent(inout) :: ncid ! pio netCDF file id + character(len=*), intent(in) :: varName ! variable we are reading + integer, intent(in) :: expected_numDims + character(len=*), intent(in) :: expected_dimNames(:) ! expected dimension name + character(len=*), intent(in) :: callingName ! calling routine + real(r8), intent(inout) :: retVal(1:, : ) + + ! local vars + character(len=32) :: subname = 'readNcdio::' + character(len=100) :: errCode = ' - Error reading. Var: ' + ! + ! netcdf read here + ! + call checkDimensions(ncid, varName, expected_numDims, expected_dimNames, subname) + call readNcdio(ncid, varName, callingName, retVal) + + end subroutine readNcdioArray2dCheckDimensions + !----------------------------------------------------------------------- + + !----------------------------------------------------------------------- + ! + !----------------------------------------------------------------------- + subroutine checkDimensions(ncid, varName, expected_numDims, expected_dimNames, callingName) + ! + ! Assert that the expected number of dimensions and dimension + ! names for a variable match the actual names on the file. + ! + use abortutils , only : endrun + use ncdio_pio , only : file_desc_t, var_desc_t, check_var, ncd_inqvdname, ncd_inqvdims + + implicit none + + ! arguments + type(file_desc_t),intent(inout) :: ncid ! pio netCDF file id + character(len=*), intent(in) :: varName ! variable we are reading + integer, intent(in) :: expected_numDims ! number of expected dimensions on the variable + character(len=*), intent(in) :: expected_dimNames(:) ! expected dimension names + character(len=*), intent(in) :: callingName ! calling routine + integer :: error_num + + ! local vars + character(len=32) :: subname = 'checkDimensions::' + type(Var_desc_t) :: var_desc ! variable descriptor + logical :: readvar ! whether the variable was found + character(len=100) :: received_dimName + integer :: d, num_dims + character(len=256) :: msg + + call check_var(ncid, varName, var_desc, readvar) + if (readvar) then + call ncd_inqvdims(ncid, num_dims, var_desc) + if (num_dims /= expected_numDims) then + write(msg, *) trim(callingName)//trim(subname)//trim(varname)//":: expected number of dimensions = ", & + expected_numDims, " num dimensions received from file = ", num_dims + call endrun(msg) + end if + do d = 1, num_dims + received_dimName = '' + call ncd_inqvdname(ncid, varname=trim(varName), dimnum=d, dname=received_dimName, err_code=error_num) + if (trim(expected_dimNames(d)) /= trim(received_dimName)) then + write(msg, *) trim(callingName)//trim(subname)//trim(varname)//":: dimension ", d, & + " expected dimension name '"//trim(expected_dimNames(d))//& + "' dimension name received from file '"//trim(received_dimName)//"'." + call endrun(msg) + end if + end do + end if + + end subroutine checkDimensions + !----------------------------------------------------------------------- + end module paramUtilMod diff --git a/components/clm/src/utils/clmfates_paraminterfaceMod.F90 b/components/clm/src/utils/clmfates_paraminterfaceMod.F90 index 07fc107f06..0db35d1939 100644 --- a/components/clm/src/utils/clmfates_paraminterfaceMod.F90 +++ b/components/clm/src/utils/clmfates_paraminterfaceMod.F90 @@ -205,6 +205,7 @@ subroutine ParametersFromNetCDF(filename, is_host_file, fates_params) real(r8), allocatable :: data(:, :) character(len=param_string_length) :: name integer :: dimension_sizes(max_dimensions) + character(len=param_string_length) :: dimension_names(max_dimensions) integer :: size_dim_1, size_dim_2 logical :: is_host_param @@ -217,7 +218,7 @@ subroutine ParametersFromNetCDF(filename, is_host_file, fates_params) num_params = fates_params%num_params() do i = 1, num_params - call fates_params%GetMetaData(i, name, dimension_shape, dimension_sizes, is_host_param) + call fates_params%GetMetaData(i, name, dimension_shape, dimension_sizes, dimension_names, is_host_param) if (is_host_file .eqv. is_host_param) then select case(dimension_shape) case(dimension_shape_scalar) @@ -232,7 +233,7 @@ subroutine ParametersFromNetCDF(filename, is_host_file, fates_params) case default call endrun(msg='unsupported number of dimensions reading parameters.') end select - call readNcdio(ncid, name, subname, data(1:size_dim_1, 1:size_dim_2)) + call readNcdio(ncid, name, dimension_shape, dimension_names, subname, data(1:size_dim_1, 1:size_dim_2)) call fates_params%SetData(i, data(1:size_dim_1, 1:size_dim_2)) end if end do