Skip to content

Commit

Permalink
io: allow disabling coordinates in history files (CICE-Consortium#935)
Browse files Browse the repository at this point in the history
* ice_history_shared: add namelist flags for coordinate variables

Users might not wish for all 8 coordinate variables to be written to
each history file, so add namelist flags allowing to disable each
coordinate variable independantly, following the approach used for the
other grid variables. Move the definition of 'ncoord' from each IO
backend to 'ice_history_shared'. Note that we already 'use' the whole of
ice_history_shared in ice_history_write, so we do not need to adjust the
use statement.

Note that the code that writes these variables in module
'ice_history_write' in each of the IO backends is not modified, it will
be adjusted in a following commit.

* io_{netcdf,pio2}/ice_history_write: define coordinates attributes in a loop

In ice_history_write::ice_write_hist, we initialize the 'var_coord' and
'coord_bounds' arrays by manually incrementing 'ind' and initializing
each elements with the corresponding information for each coordinate
variable. The rest of the code in that subroutine instead uses a loop on
'ncoord', which is a parameter holding the number of coordinates
variable, along with 'select' statements.

The latter approach is more elegant and also more flexible if we change
the number of coordinates variables. Refactor the code to use a loop and
leverage the integers n_{u,t,n,e}{lon,lat} added in the previous commit,
in both io_netcdf and io_pio2, which have identical code in this part of
the subroutine.

* io_{netcdf,pio2}/ice_history_write: use namelist flags for coordinate variables

In order to enable the namelist flags for each coordinate variables
added in the previous commit, adjust the code of the io_netcdf and
io_pio2 backends by checking the value of 'icoord(i)' for each loop on
'ncoord' that calls NetCDF or PIO functions.

Add the new flags to the reference namelist.
  • Loading branch information
phil-blain authored Feb 23, 2024
1 parent aca8357 commit 64177e3
Show file tree
Hide file tree
Showing 5 changed files with 221 additions and 163 deletions.
23 changes: 23 additions & 0 deletions cicecore/cicedyn/analysis/ice_history.F90
Original file line number Diff line number Diff line change
Expand Up @@ -447,6 +447,14 @@ subroutine init_hist (dt)
if (f_Tsnz (1:1) /= 'x') f_VGRDs = .true.
if (tr_fsd) f_NFSD = .true.

call broadcast_scalar (f_tlon, master_task)
call broadcast_scalar (f_tlat, master_task)
call broadcast_scalar (f_ulon, master_task)
call broadcast_scalar (f_ulat, master_task)
call broadcast_scalar (f_nlon, master_task)
call broadcast_scalar (f_nlat, master_task)
call broadcast_scalar (f_elon, master_task)
call broadcast_scalar (f_elat, master_task)
call broadcast_scalar (f_tmask, master_task)
call broadcast_scalar (f_umask, master_task)
call broadcast_scalar (f_nmask, master_task)
Expand Down Expand Up @@ -1973,6 +1981,21 @@ subroutine init_hist (dt)
! floe size distribution
call init_hist_fsd_4Df

!-----------------------------------------------------------------
! fill icoord array with namelist values
!-----------------------------------------------------------------

icoord=.true.

icoord(n_tlon ) = f_tlon
icoord(n_tlat ) = f_tlat
icoord(n_ulon ) = f_ulon
icoord(n_ulat ) = f_ulat
icoord(n_nlon ) = f_nlon
icoord(n_nlat ) = f_nlat
icoord(n_elon ) = f_elon
icoord(n_elat ) = f_elat

!-----------------------------------------------------------------
! fill igrd array with namelist values
!-----------------------------------------------------------------
Expand Down
19 changes: 19 additions & 0 deletions cicecore/cicedyn/analysis/ice_history_shared.F90
Original file line number Diff line number Diff line change
Expand Up @@ -131,6 +131,7 @@ module ice_history_shared
avail_hist_fields(max_avail_hist_fields)

integer (kind=int_kind), parameter, public :: &
ncoord = 8 , & ! number of coordinate variables: TLON, TLAT, ULON, ULAT, NLON, NLAT, ELON, ELAT
nvar_grd = 21 , & ! number of grid fields that can be written
! excluding grid vertices
nvar_grdz = 6 ! number of category/vertical grid fields written
Expand Down Expand Up @@ -165,6 +166,7 @@ module ice_history_shared
avgct(max_nstrm) ! average sample counter

logical (kind=log_kind), public :: &
icoord(ncoord) , & ! true if coord field is written to output file
igrd (nvar_grd), & ! true if grid field is written to output file
igrdz(nvar_grdz) ! true if category/vertical grid field is written

Expand Down Expand Up @@ -194,6 +196,10 @@ module ice_history_shared
!---------------------------------------------------------------

logical (kind=log_kind), public :: &
f_tlon = .true., f_tlat = .true., &
f_ulon = .true., f_ulat = .true., &
f_nlon = .true., f_nlat = .true., &
f_elon = .true., f_elat = .true., &
f_tmask = .true., f_umask = .true., &
f_nmask = .true., f_emask = .true., &
f_blkmask = .true., &
Expand Down Expand Up @@ -362,6 +368,10 @@ module ice_history_shared
!---------------------------------------------------------------

namelist / icefields_nml / &
f_tlon , f_tlat , &
f_ulon , f_ulat , &
f_nlon , f_nlat , &
f_elon , f_elat , &
f_tmask , f_umask , &
f_nmask , f_emask , &
f_blkmask , &
Expand Down Expand Up @@ -529,6 +539,15 @@ module ice_history_shared
!---------------------------------------------------------------

integer (kind=int_kind), parameter, public :: &
n_tlon = 1, &
n_tlat = 2, &
n_ulon = 3, &
n_ulat = 4, &
n_nlon = 5, &
n_nlat = 6, &
n_elon = 7, &
n_elat = 8, &

n_tmask = 1, &
n_umask = 2, &
n_nmask = 3, &
Expand Down
176 changes: 90 additions & 86 deletions cicecore/cicedyn/infrastructure/io/io_netcdf/ice_history_write.F90
Original file line number Diff line number Diff line change
Expand Up @@ -113,9 +113,6 @@ subroutine ice_write_hist (ns)
! time coord
TYPE(coord_attributes) :: time_coord

! 8 coordinate variables: TLON, TLAT, ULON, ULAT, NLON, NLAT, ELON, ELAT
INTEGER (kind=int_kind), PARAMETER :: ncoord = 8

! 4 vertices in each grid cell
INTEGER (kind=int_kind), PARAMETER :: nverts = 4

Expand Down Expand Up @@ -263,39 +260,42 @@ subroutine ice_write_hist (ns)
! define information for required time-invariant variables
!-----------------------------------------------------------------

ind = 0
ind = ind + 1
var_coord(ind) = coord_attributes('TLON', &
'T grid center longitude', 'degrees_east')
coord_bounds(ind) = 'lont_bounds'
ind = ind + 1
var_coord(ind) = coord_attributes('TLAT', &
'T grid center latitude', 'degrees_north')
coord_bounds(ind) = 'latt_bounds'
ind = ind + 1
var_coord(ind) = coord_attributes('ULON', &
'U grid center longitude', 'degrees_east')
coord_bounds(ind) = 'lonu_bounds'
ind = ind + 1
var_coord(ind) = coord_attributes('ULAT', &
'U grid center latitude', 'degrees_north')
coord_bounds(ind) = 'latu_bounds'
ind = ind + 1
var_coord(ind) = coord_attributes('NLON', &
'N grid center longitude', 'degrees_east')
coord_bounds(ind) = 'lonn_bounds'
ind = ind + 1
var_coord(ind) = coord_attributes('NLAT', &
'N grid center latitude', 'degrees_north')
coord_bounds(ind) = 'latn_bounds'
ind = ind + 1
var_coord(ind) = coord_attributes('ELON', &
'E grid center longitude', 'degrees_east')
coord_bounds(ind) = 'lone_bounds'
ind = ind + 1
var_coord(ind) = coord_attributes('ELAT', &
'E grid center latitude', 'degrees_north')
coord_bounds(ind) = 'late_bounds'
do ind = 1, ncoord
select case (ind)
case(n_tlon)
var_coord(ind) = coord_attributes('TLON', &
'T grid center longitude', 'degrees_east')
coord_bounds(ind) = 'lont_bounds'
case(n_tlat)
var_coord(ind) = coord_attributes('TLAT', &
'T grid center latitude', 'degrees_north')
coord_bounds(ind) = 'latt_bounds'
case(n_ulon)
var_coord(ind) = coord_attributes('ULON', &
'U grid center longitude', 'degrees_east')
coord_bounds(ind) = 'lonu_bounds'
case(n_ulat)
var_coord(ind) = coord_attributes('ULAT', &
'U grid center latitude', 'degrees_north')
coord_bounds(ind) = 'latu_bounds'
case(n_nlon)
var_coord(ind) = coord_attributes('NLON', &
'N grid center longitude', 'degrees_east')
coord_bounds(ind) = 'lonn_bounds'
case(n_nlat)
var_coord(ind) = coord_attributes('NLAT', &
'N grid center latitude', 'degrees_north')
coord_bounds(ind) = 'latn_bounds'
case(n_elon)
var_coord(ind) = coord_attributes('ELON', &
'E grid center longitude', 'degrees_east')
coord_bounds(ind) = 'lone_bounds'
case(n_elat)
var_coord(ind) = coord_attributes('ELAT', &
'E grid center latitude', 'degrees_north')
coord_bounds(ind) = 'late_bounds'
end select
end do

var_grdz(1) = coord_attributes('NCAT', 'category maximum thickness', 'm')
var_grdz(2) = coord_attributes('VGRDi', 'vertical ice levels', '1')
Expand Down Expand Up @@ -406,18 +406,20 @@ subroutine ice_write_hist (ns)
dimid(3) = timid

do i = 1, ncoord
call ice_hist_coord_def(ncid, var_coord(i), lprecision, dimid(1:2), varid)
call ice_write_hist_fill(ncid,varid,var_coord(i)%short_name,history_precision)
if (var_coord(i)%short_name == 'ULAT') then
status = nf90_put_att(ncid,varid,'comment', &
'Latitude of NE corner of T grid cell')
call ice_check_nc(status, subname// ' ERROR: defining comment for '//var_coord(i)%short_name, &
file=__FILE__, line=__LINE__)
endif
if (f_bounds) then
status = nf90_put_att(ncid, varid, 'bounds', coord_bounds(i))
call ice_check_nc(status, subname// ' ERROR: defining bounds for '//var_coord(i)%short_name, &
file=__FILE__, line=__LINE__)
if(icoord(i)) then
call ice_hist_coord_def(ncid, var_coord(i), lprecision, dimid(1:2), varid)
call ice_write_hist_fill(ncid,varid,var_coord(i)%short_name,history_precision)
if (var_coord(i)%short_name == 'ULAT') then
status = nf90_put_att(ncid,varid,'comment', &
'Latitude of NE corner of T grid cell')
call ice_check_nc(status, subname// ' ERROR: defining comment for '//var_coord(i)%short_name, &
file=__FILE__, line=__LINE__)
endif
if (f_bounds) then
status = nf90_put_att(ncid, varid, 'bounds', coord_bounds(i))
call ice_check_nc(status, subname// ' ERROR: defining bounds for '//var_coord(i)%short_name, &
file=__FILE__, line=__LINE__)
endif
endif
enddo

Expand Down Expand Up @@ -707,44 +709,46 @@ subroutine ice_write_hist (ns)
!-----------------------------------------------------------------

do i = 1,ncoord
call broadcast_scalar(var_coord(i)%short_name,master_task)
SELECT CASE (var_coord(i)%short_name)
CASE ('TLON')
! Convert T grid longitude from -180 -> 180 to 0 to 360
work1 = TLON*rad_to_deg + c360
where (work1 > c360) work1 = work1 - c360
where (work1 < c0 ) work1 = work1 + c360
call gather_global(work_g1,work1,master_task,distrb_info)
CASE ('TLAT')
work1 = TLAT*rad_to_deg
call gather_global(work_g1,work1,master_task,distrb_info)
CASE ('ULON')
work1 = ULON*rad_to_deg
call gather_global(work_g1,work1,master_task,distrb_info)
CASE ('ULAT')
work1 = ULAT*rad_to_deg
call gather_global(work_g1,work1,master_task,distrb_info)
CASE ('NLON')
work1 = NLON*rad_to_deg
call gather_global(work_g1,work1,master_task,distrb_info)
CASE ('NLAT')
work1 = NLAT*rad_to_deg
call gather_global(work_g1,work1,master_task,distrb_info)
CASE ('ELON')
work1 = ELON*rad_to_deg
call gather_global(work_g1,work1,master_task,distrb_info)
CASE ('ELAT')
work1 = ELAT*rad_to_deg
call gather_global(work_g1,work1,master_task,distrb_info)
END SELECT

if (my_task == master_task) then
status = nf90_inq_varid(ncid, var_coord(i)%short_name, varid)
call ice_check_nc(status, subname// ' ERROR: getting varid for '//var_coord(i)%short_name, &
file=__FILE__, line=__LINE__)
status = nf90_put_var(ncid,varid,work_g1)
call ice_check_nc(status, subname// ' ERROR: writing'//var_coord(i)%short_name, &
file=__FILE__, line=__LINE__)
if(icoord(i)) then
call broadcast_scalar(var_coord(i)%short_name,master_task)
SELECT CASE (var_coord(i)%short_name)
CASE ('TLON')
! Convert T grid longitude from -180 -> 180 to 0 to 360
work1 = TLON*rad_to_deg + c360
where (work1 > c360) work1 = work1 - c360
where (work1 < c0 ) work1 = work1 + c360
call gather_global(work_g1,work1,master_task,distrb_info)
CASE ('TLAT')
work1 = TLAT*rad_to_deg
call gather_global(work_g1,work1,master_task,distrb_info)
CASE ('ULON')
work1 = ULON*rad_to_deg
call gather_global(work_g1,work1,master_task,distrb_info)
CASE ('ULAT')
work1 = ULAT*rad_to_deg
call gather_global(work_g1,work1,master_task,distrb_info)
CASE ('NLON')
work1 = NLON*rad_to_deg
call gather_global(work_g1,work1,master_task,distrb_info)
CASE ('NLAT')
work1 = NLAT*rad_to_deg
call gather_global(work_g1,work1,master_task,distrb_info)
CASE ('ELON')
work1 = ELON*rad_to_deg
call gather_global(work_g1,work1,master_task,distrb_info)
CASE ('ELAT')
work1 = ELAT*rad_to_deg
call gather_global(work_g1,work1,master_task,distrb_info)
END SELECT

if (my_task == master_task) then
status = nf90_inq_varid(ncid, var_coord(i)%short_name, varid)
call ice_check_nc(status, subname// ' ERROR: getting varid for '//var_coord(i)%short_name, &
file=__FILE__, line=__LINE__)
status = nf90_put_var(ncid,varid,work_g1)
call ice_check_nc(status, subname// ' ERROR: writing'//var_coord(i)%short_name, &
file=__FILE__, line=__LINE__)
endif
endif
enddo

Expand Down
Loading

0 comments on commit 64177e3

Please sign in to comment.