Skip to content

UNST-9376: Accept both "waterlevel" and "waterlevelbnd" #369

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
adrimourits wants to merge 2 commits into
base: main
Choose a base branch
from
Open

UNST-9376: Accept both "waterlevel" and "waterlevelbnd" #369

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
12c3bfc
UNST-9376: Accept both "waterlevel" and "waterlevelbnd"
adrimourits Nov 7, 2025
7f11fd6
UNST-9376: Ugly temporary fix: ecNetCDFInit: skip reading layers
adrimourits Nov 13, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
70 changes: 46 additions & 24 deletions src/utils_lgpl/ec_module/packages/ec_module/src/ec_netcdf_timeseries.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -137,6 +137,7 @@ end function ecNetCDFFree1dArray
!! Open a netCDF file, store ncid, standard names and long names ....
!! Open only if this file is not already opened, so check the list of nc-objects first and return a pointer ....
function ecNetCDFInit(ncname, ncptr, iostat) result(success)
use m_alloc, only: realloc
logical :: success
logical :: station_id_found
logical :: time_found
Expand Down Expand Up @@ -190,12 +191,15 @@ function ecNetCDFInit(ncname, ncptr, iostat) result(success)
if (ierr /= 0) return
allocate (ncptr%offsets(nVars), stat=ierr)
if (ierr /= 0) return
allocate (ncptr%time_varying_var(nVars), stat=ierr)
if (ierr /= 0) return
ncptr%standard_names = ' '
ncptr%long_names = ' '
ncptr%variable_names = ' '
ncptr%fillvalues = -huge(dp)
ncptr%scales = 1.0_dp
ncptr%offsets = 0.0_dp
ncptr%time_varying_var = .false.
allocate (var_dimids(n_dims, nVars), stat=ierr) ! NOTE: n_dims is only an upper bound here!
if (ierr /= 0) return
allocate (var_ndims(nVars), stat=ierr)
Expand Down Expand Up @@ -263,31 +267,40 @@ function ecNetCDFInit(ncname, ncptr, iostat) result(success)
time_found = .true.
end if
end if

if (any(var_dimids(:, iVars) == ncptr%timedimid)) then
ncptr%time_varying_var(iVars) = .true.
end if

! Check for important var: was it vertical layering?
positive = ''
positive = '' ! axis = ''
zunits = ''
ierr = ncu_get_att(ncptr%ncid, iVars, 'positive', positive)
ierr = ncu_get_att(ncptr%ncid, iVars, 'positive', positive) ! ierr = ncu_get_att(ncptr%ncid, iVars, 'axis', axis)
if (len_trim(positive) > 0) then ! Identified a layercoord variable, by its positive:up/down attribute
! NOTE: officially, a vertical coord var may also be identified by a unit of pressure, but we don't support that here.
ncptr%layervarid = iVars
ncptr%layerdimid = var_dimids(1, iVars) ! For convenience also store the dimension ID explicitly
ncptr%nLayer = ncptr%dimlen(ncptr%layerdimid)
allocate (ncptr%vp(ncptr%nLayer), stat=ierr)
if (ierr /= 0) return
ierr = nf90_get_var(ncptr%ncid, ncptr%layervarid, ncptr%vp, (/1/), (/ncptr%nLayer/))
if (ierr /= NF90_NOERR) return
ierr = ncu_get_att(ncptr%ncid, iVars, 'units', zunits)
if (ierr /= NF90_NOERR) return
if (strcmpi(zunits, 'm')) then
if (strcmpi(positive, 'up')) ncptr%vptyp = BC_VPTYP_ZDATUM ! z upward from datum, unmodified z-values
if (strcmpi(positive, 'down')) ncptr%vptyp = BC_VPTYP_ZSURF ! z downward
else
if (strcmpi(positive, 'up')) ncptr%vptyp = BC_VPTYP_PERCBED ! sigma upward
if (strcmpi(positive, 'down')) ncptr%vptyp = BC_VPTYP_PERCSURF ! sigma downward
end if
if (ncptr%vptyp < 1) then
call setECMessage("ec_bcreader::ecNetCDFCreate: Unable to determine vertical coordinate system.")
! This line does not work: if (any((/ 'z ', 'zcoordinate ', 'zcoordinate_c'/) == ncptr%standard_names(iVars))) then
if (any((/ 'z ', 'zcoordinate '/) == ncptr%standard_names(iVars))) then
Copy link
Member

@arthurvd arthurvd Nov 14, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Maybe we need to do this differently, let's refine a proper 3D initialization under the other dedicated issue: UNST-8909 (so not do this here)

Copy link
Member

@arthurvd arthurvd Nov 14, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

p.s. Coding guidelines require the modern format for initialization arrays, namely: [ .. ] instead of (/ .. /).

ncptr%layervarid = iVars
ncptr%layerdimid = var_dimids(1, iVars) ! For convenience also store the dimension ID explicitly
ncptr%nLayer = ncptr%dimlen(ncptr%layerdimid)
if (ncptr%nLayer > size(ncptr%vp)) then
call realloc (ncptr%vp, ncptr%nLayer, stat=ierr, keepExisting=.true.)
end if
if (ierr /= 0) return
Comment on lines +286 to +289
Copy link
Member

@arthurvd arthurvd Nov 14, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If 3D work is postponed until UNST-8909, why not simply have an if .. cycle here?

Suggested change
if (ncptr%nLayer > size(ncptr%vp)) then
call realloc (ncptr%vp, ncptr%nLayer, stat=ierr, keepExisting=.true.)
end if
if (ierr /= 0) return
allocate (ncptr%vp(ncptr%nLayer), stat=ierr)
if (ierr /= 0) cycle

ierr = nf90_get_var(ncptr%ncid, ncptr%layervarid, ncptr%vp, (/1/), (/ncptr%nLayer/))
if (ierr /= NF90_NOERR) return
ierr = ncu_get_att(ncptr%ncid, iVars, 'units', zunits)
if (ierr /= NF90_NOERR) return
if (strcmpi(zunits, 'm')) then
if (strcmpi(positive, 'up')) ncptr%vptyp = BC_VPTYP_ZDATUM ! z upward from datum, unmodified z-values
if (strcmpi(positive, 'down')) ncptr%vptyp = BC_VPTYP_ZSURF ! z downward
else
if (strcmpi(positive, 'up')) ncptr%vptyp = BC_VPTYP_PERCBED ! sigma upward
if (strcmpi(positive, 'down')) ncptr%vptyp = BC_VPTYP_PERCSURF ! sigma downward
end if
if (ncptr%vptyp < 1) then
call setECMessage("ec_bcreader::ecNetCDFCreate: Unable to determine vertical coordinate system.")
end if
end if
end if
end do
Expand Down Expand Up @@ -320,26 +333,24 @@ recursive function ecNetCDFScan(ncptr, quantity, location, q_id, l_id, dimids, v
integer, dimension(:), allocatable :: dimids_check

success = .false.
ltl = len_trim_nobnd(quantity)

! search for standard_name
do ivar = 1, ncptr%nVars
ltl = len_trim(quantity)
if (strcmpi(trim(ncptr%standard_names(ivar)), trim(quantity))) exit
if (strcmpi(trim(ncptr%standard_names(ivar)), trim(quantity), ltl)) exit
Copy link
Member

@arthurvd arthurvd Nov 14, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I would hope we can make this work instead using the existing ncstdnames approach, by having something like this in ecProviderCreateNetCdfItems():

      case ('waterlevelbnd')
         ncvarnames(1) = 'waterlevel'
         ncstdnames(1) = 'sea_surface_height'

Let's refine this together under UNST-9376 for proper followup in next sprint.

end do
vmax = 1

! if standard_name not found, search for long_name
if (ivar > ncptr%nVars) then
do ivar = 1, ncptr%nVars
ltl = len_trim(quantity)
if (strcmpi(ncptr%long_names(ivar), quantity, ltl)) exit
end do
end if

! if also long_name not found, search for variable_name
if (ivar > ncptr%nVars .and. allocated(ncptr%variable_names)) then
do ivar = 1, ncptr%nVars
ltl = len_trim(quantity)
if (strcmpi(ncptr%variable_names(ivar), quantity, ltl)) exit
end do
end if
Expand Down Expand Up @@ -455,4 +466,15 @@ function ecNetCDFGetAttrib(ncptr, q_id, attribute_name, attribute_value) result(
success = .true.
end function ecNetCDFGetAttrib

function len_trim_nobnd(name) result(len_without_bnd)
integer :: len_without_bnd
character(len=*), intent(in) :: name

len_without_bnd = len_trim(name)
if (strcmpi(name(len_trim(name)-2:len_trim(name)), 'bnd')) then
len_without_bnd = len_without_bnd - 3
end if
end function len_trim_nobnd


end module m_ec_netcdf_timeseries
1 change: 1 addition & 0 deletions src/utils_lgpl/ec_module/packages/ec_module/src/ec_typedefs.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -201,6 +201,7 @@ module m_ec_typedefs
character(len=50) :: timeunit !< netcdf-convention time unit definition
integer :: vptyp = -1 !< vertical coordinate type
real(hp), allocatable, dimension(:) :: vp !< vertical coordinate (layers)
logical, allocatable, dimension(:) :: time_varying_var !< Indicates whether variable i is time varying. Needed for nesting.
end type tEcNetCDF

type tEcNetCDFPtr
Expand Down