Skip to content

Commit

Permalink
Add check to 'read_global_mask' to check the i/j dimensions
Browse files Browse the repository at this point in the history
of the umd land mask file.

Fixes ufs-community#1000.
  • Loading branch information
George Gayno committed Jan 3, 2025
1 parent 31d5486 commit ddfbc61
Showing 1 changed file with 25 additions and 2 deletions.
27 changes: 25 additions & 2 deletions sorc/orog_mask_tools.fd/orog.fd/io_utils.F90
Original file line number Diff line number Diff line change
Expand Up @@ -597,16 +597,39 @@ subroutine read_global_mask(imn, jmn, mask)

integer(1), intent(out) :: mask(imn,jmn)

integer :: ncid, id_var, error
integer :: ncid, id_var, id_dim, error, idim, jdim

print*,"- OPEN AND READ ./landcover.umd.30s.nc"

error=nf90_open("./landcover.umd.30s.nc",nf90_nowrite,ncid)
call netcdf_err(error, 'Open file landcover.umd.30s.nc' )
call netcdf_err(error, 'Opening file landcover.umd.30s.nc' )

error=nf90_inq_dimid(ncid, 'idim', id_dim)
call netcdf_err(error, 'Inquire dimid of idim' )

error=nf90_inquire_dimension(ncid,id_dim,len=idim)
call netcdf_err(error, 'Reading idim' )

if (imn /= idim) then
print*,"FATAL ERROR: i-dimensions do not match."
endif

error=nf90_inq_dimid(ncid, 'jdim', id_dim)
call netcdf_err(error, 'Inquire dimid of jdim' )

error=nf90_inquire_dimension(ncid,id_dim,len=jdim)
call netcdf_err(error, 'Reading jdim' )

if (jmn /= jdim) then
print*,"FATAL ERROR: j-dimensions do not match."
endif

error=nf90_inq_varid(ncid, 'land_mask', id_var)
call netcdf_err(error, 'Inquire varid of land_mask')

error=nf90_get_var(ncid, id_var, mask)
call netcdf_err(error, 'Inquire data of land_mask')

error = nf90_close(ncid)

call transpose_mask(imn,jmn,mask)
Expand Down

0 comments on commit ddfbc61

Please sign in to comment.