Convert a netCDF error code to a nice error message. Writes to error_unit
All the arguments except istatus
are optional and are used to
give more detailed information about the error
The abort
argument will abort the simulation
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
integer(kind=kind_nf), | intent(in) | :: | istatus |
The netCDF error code to convert to a message |
||
integer(kind=kind_nf), | intent(in), | optional | :: | ncid |
netCDF ID of the file |
|
integer(kind=kind_nf), | intent(in), | optional | :: | varid |
netCDF ID of the variable |
|
integer(kind=kind_nf), | intent(in), | optional | :: | dimid |
netCDF ID of the dimension |
|
character(len=*), | intent(in), | optional | :: | file |
Name of the file |
|
character(len=*), | intent(in), | optional | :: | dim |
Name of the dimension |
|
character(len=*), | intent(in), | optional | :: | var |
Name of the variable |
|
character(len=*), | intent(in), | optional | :: | att |
Name of the attribute |
|
character(len=*), | intent(in), | optional | :: | message |
Custom text to append to the error message |
|
logical, | intent(in), | optional | :: | abort |
Immediately abort the program if present and true |
subroutine netcdf_error &
(istatus, ncid, varid, dimid, file, dim, var, att, message, abort)
use mp, only: mp_abort
use file_utils, only: error_unit
use netcdf, only: NF90_GLOBAL
implicit none
!> The netCDF error code to convert to a message
integer (kind_nf), intent (in) :: istatus
!> netCDF ID of the file
integer (kind_nf), intent (in), optional :: ncid
!> netCDF ID of the variable
integer (kind_nf), intent (in), optional :: varid
!> netCDF ID of the dimension
integer (kind_nf), intent (in), optional :: dimid
!> Name of the file
character (*), intent (in), optional :: file
!> Name of the dimension
character (*), intent (in), optional :: dim
!> Name of the variable
character (*), intent (in), optional :: var
!> Name of the attribute
character (*), intent (in), optional :: att
!> Custom text to append to the error message
character (*), intent (in), optional :: message
!> Immediately abort the program if present and true
logical, intent (in), optional :: abort
integer (kind_nf) :: ist
integer :: ierr
character (20) :: varname, dimname
ierr = error_unit()
write (ierr, '(2a)', advance='no') 'ERROR: ', trim (nf90_strerror (istatus))
if (present(file)) &
write (ierr, '(2a)', advance='no') ' in file: ', trim (file)
if (present(dim)) &
write (ierr, '(2a)', advance='no') ' in dimension: ', trim (dim)
if (present(var)) &
write (ierr, '(2a)', advance='no') ' in variable: ', trim (var)
if (present(varid)) then
if (present(ncid)) then
if (present(att) ) then
if (varid == NF90_GLOBAL) then
write (ierr, '(2a)') ' in global attribute: ', trim(att)
else
write (ierr, '(2a)') ' with the attribute: ', trim(att)
end if
else
ist = nf90_inquire_variable (ncid, varid, varname)
if (ist == NF90_NOERR) then
write (ierr, '(a,i8,2a)', advance='no') ' in varid: ', varid, &
& ' variable name: ', trim (varname)
else
write (ierr, *) ''
write (ierr, '(3a,i8,a,i8)', advance='no') 'ERROR in netcdf_error: ', &
trim (nf90_strerror(ist)), ' in varid: ', varid, &
', ncid: ', ncid
end if
end if
else
write (ierr, *) ''
write (ierr, '(2a)', advance='no') 'ERROR in netcdf_error: ', &
& 'ncid missing while varid present in the argument'
end if
end if
if (present(dimid)) then
if (present(ncid)) then
ist = nf90_inquire_dimension (ncid, dimid, dimname)
if (ist == NF90_NOERR) then
write (ierr, '(a,i8,2a)', advance='no') ' in dimid: ', dimid, &
& ' dimension name: ', trim (dimname)
else
write (ierr, *) ''
write (ierr, '(3a,i8,a,i8)', advance='no') 'ERROR in netcdf_error: ', &
trim (nf90_strerror(ist)), ' in dimid: ', dimid, &
', ncid: ', ncid
end if
else
write (ierr, *) ''
write (ierr, '(2a)', advance='no') 'ERROR in netcdf_error: ', &
& 'ncid missing while dimid present in the argument'
end if
end if
if (present(message)) write (ierr, '(a)', advance='no') trim(message)
! append line-break
write(ierr,*)
! maybe this switching is not necessary.
! if error is detected, the program should abort immediately
if(present(abort)) then
if(abort) then
call mp_abort('Aborted by netcdf_error')
endif
endif
end subroutine netcdf_error