Skip to content

Commit

Permalink
cicecore: correct initial condition metadata (CICE-Consortium#818)
Browse files Browse the repository at this point in the history
* ice_history_write: fix initial condition metadata under 'hist_avg'

When writing averaged history outputs (hist_avg=.true.), this setting
also affects the initial condition. Even if the actual data variables
written to the initial condition are not averaged (they are taken
more or less directly from the restart or the hard-coded defaults,
modulo aggregation over categories), their attributes ('cell_method' and
'time_rep') imply they are averaged, and the 'bound' attribute of the
'time' variable refers to the 'time_bounds' variable.

Make the metadata of the initial condition more correct by:
- not writing the 'time_bounds' variable (and the corresponding 'd2' dimension)
- not writing the 'bounds' attribute of the 'time' variable
- not writing the 'cell_method' attributes of each variable
- writing the 'time_rep' attribute of each variable as 'instantaneous'
instead of 'averaged'.

Do this by checking 'write_ic' at all places where we check for the
value of 'hist_avg' to write the above variables and attributes in each
of the 3 IO backends (binary, netcdf, pio2).

* drivers/{nemo_concepts,standalone}: write initial condition at initial time

In CICE_InitMod::cice_init, we call ice_calendar::advance_timestep
before writing the initial condition, such that the 'time' variable in
the initial condition is not zero; it has a value of 1*dt (the model
time step). The initial condition filename also reflects this, since
'msec' (model seconds) also has a value of 1*dt and is used in
ice_history_shared::construct_filename. This leads to the initial
condition filename not corresponding to the model initialization
date/time but rather 1*dt later.

Since we call 'accum_hist' after initializing the forcing, any forcing
field written to the initial condition has values corresponding to
msec=dt, whereas the ice state corresponds to msec=0, leading to an
inconsistency.

Fix that by calling 'accum_hist' to write the initial condition _before_
calling 'advance_timestep'. Since we now call 'accum_hist' before
initializing the forcing, any forcing field written to the initial
condition have its default, hard-coded value, instead of its value at
time=dt. An improvement would be to read the forcing at time=dt, write
the initial condition, advance the time step, and read the forcing
again, but let's not complicate things too much for now.
  • Loading branch information
phil-blain authored Mar 13, 2023
1 parent 0fcc140 commit 9a55ad9
Show file tree
Hide file tree
Showing 5 changed files with 27 additions and 24 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -158,6 +158,7 @@ subroutine ice_write_hist(ns)
trim(avail_hist_fields(n)%vcomment)

if (histfreq(ns) == '1' .or. .not. hist_avg &
.or. write_ic &
.or. n==n_divu(ns) .or. n==n_shear(ns) & ! snapshots
.or. n==n_sig1(ns) .or. n==n_sig2(ns) &
.or. n==n_sigP(ns) .or. n==n_trsig(ns) &
Expand Down Expand Up @@ -186,7 +187,7 @@ subroutine ice_write_hist(ns)
write (nu_hdr, 994) nrec,trim(avail_hist_fields(n)%vname), &
trim(avail_hist_fields(n)%vdesc),trim(avail_hist_fields(n)%vunit),nn

if (histfreq(ns) == '1' .or. .not. hist_avg) then
if (histfreq(ns) == '1' .or. .not. hist_avg .or. write_ic) then
write (nu_hdr, 996) nrec,trim(avail_hist_fields(n)%vname), &
'time_rep','instantaneous'
else
Expand All @@ -210,7 +211,7 @@ subroutine ice_write_hist(ns)
write (nu_hdr, 993) nrec,trim(avail_hist_fields(n)%vname), &
trim(avail_hist_fields(n)%vdesc),trim(avail_hist_fields(n)%vunit),k

if (histfreq(ns) == '1' .or. .not. hist_avg) then
if (histfreq(ns) == '1' .or. .not. hist_avg .or. write_ic) then
write (nu_hdr, 996) nrec,trim(avail_hist_fields(n)%vname), &
'time_rep','instantaneous'
else
Expand All @@ -234,7 +235,7 @@ subroutine ice_write_hist(ns)
write (nu_hdr, 993) nrec,trim(avail_hist_fields(n)%vname), &
trim(avail_hist_fields(n)%vdesc),trim(avail_hist_fields(n)%vunit),nn,k

if (histfreq(ns) == '1' .or. .not. hist_avg) then
if (histfreq(ns) == '1' .or. .not. hist_avg .or. write_ic) then
write (nu_hdr, 996) nrec,trim(avail_hist_fields(n)%vname), &
'time_rep','instantaneous'
else
Expand All @@ -258,7 +259,7 @@ subroutine ice_write_hist(ns)
write (nu_hdr, 993) nrec,trim(avail_hist_fields(n)%vname), &
trim(avail_hist_fields(n)%vdesc),trim(avail_hist_fields(n)%vunit),nn,k

if (histfreq(ns) == '1' .or. .not. hist_avg) then
if (histfreq(ns) == '1' .or. .not. hist_avg .or. write_ic) then
write (nu_hdr, 996) nrec,trim(avail_hist_fields(n)%vname), &
'time_rep','instantaneous'
else
Expand All @@ -282,7 +283,7 @@ subroutine ice_write_hist(ns)
write (nu_hdr, 993) nrec,trim(avail_hist_fields(n)%vname), &
trim(avail_hist_fields(n)%vdesc),trim(avail_hist_fields(n)%vunit),nn,k

if (histfreq(ns) == '1' .or. .not. hist_avg) then
if (histfreq(ns) == '1' .or. .not. hist_avg .or. write_ic) then
write (nu_hdr, 996) nrec,trim(avail_hist_fields(n)%vname), &
'time_rep','instantaneous'
else
Expand All @@ -307,7 +308,7 @@ subroutine ice_write_hist(ns)
write (nu_hdr, 993) nrec,trim(avail_hist_fields(n)%vname), &
trim(avail_hist_fields(n)%vdesc),trim(avail_hist_fields(n)%vunit),nn,k

if (histfreq(ns) == '1' .or. .not. hist_avg) then
if (histfreq(ns) == '1' .or. .not. hist_avg .or. write_ic) then
write (nu_hdr, 996) nrec,trim(avail_hist_fields(n)%vname), &
'time_rep','instantaneous'
else
Expand All @@ -333,7 +334,7 @@ subroutine ice_write_hist(ns)
write (nu_hdr, 993) nrec,trim(avail_hist_fields(n)%vname), &
trim(avail_hist_fields(n)%vdesc),trim(avail_hist_fields(n)%vunit),nn,k

if (histfreq(ns) == '1' .or. .not. hist_avg) then
if (histfreq(ns) == '1' .or. .not. hist_avg .or. write_ic) then
write (nu_hdr, 996) nrec,trim(avail_hist_fields(n)%vname), &
'time_rep','instantaneous'
else
Expand All @@ -359,7 +360,7 @@ subroutine ice_write_hist(ns)
write (nu_hdr, 993) nrec,trim(avail_hist_fields(n)%vname), &
trim(avail_hist_fields(n)%vdesc),trim(avail_hist_fields(n)%vunit),nn,k

if (histfreq(ns) == '1' .or. .not. hist_avg) then
if (histfreq(ns) == '1' .or. .not. hist_avg .or. write_ic) then
write (nu_hdr, 996) nrec,trim(avail_hist_fields(n)%vname), &
'time_rep','instantaneous'
else
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -159,7 +159,7 @@ subroutine ice_write_hist (ns)
! define dimensions
!-----------------------------------------------------------------

if (hist_avg) then
if (hist_avg .and. .not. write_ic) then
status = nf90_def_dim(ncid,'d2',2,boundid)
if (status /= nf90_noerr) call abort_ice(subname// &
'ERROR: defining dim d2')
Expand Down Expand Up @@ -241,7 +241,7 @@ subroutine ice_write_hist (ns)
call abort_ice(subname//'ERROR: invalid calendar settings')
endif

if (hist_avg) then
if (hist_avg .and. .not. write_ic) then
status = nf90_put_att(ncid,varid,'bounds','time_bounds')
if (status /= nf90_noerr) call abort_ice(subname// &
'ERROR: time bounds')
Expand All @@ -251,7 +251,7 @@ subroutine ice_write_hist (ns)
! Define attributes for time bounds if hist_avg is true
!-----------------------------------------------------------------

if (hist_avg) then
if (hist_avg .and. .not. write_ic) then
dimid(1) = boundid
dimid(2) = timid
status = nf90_def_var(ncid,'time_bounds',lprecision,dimid(1:2),varid)
Expand Down Expand Up @@ -745,7 +745,7 @@ subroutine ice_write_hist (ns)
! write time_bounds info
!-----------------------------------------------------------------

if (hist_avg) then
if (hist_avg .and. .not. write_ic) then
status = nf90_inq_varid(ncid,'time_bounds',varid)
if (status /= nf90_noerr) call abort_ice(subname// &
'ERROR: getting time_bounds id')
Expand Down Expand Up @@ -1236,7 +1236,7 @@ end subroutine ice_write_hist
subroutine ice_write_hist_attrs(ncid, varid, hfield, ns)

use ice_kinds_mod
use ice_calendar, only: histfreq, histfreq_n
use ice_calendar, only: histfreq, histfreq_n, write_ic
use ice_history_shared, only: ice_hist_field, history_precision, &
hist_avg
#ifdef USE_NETCDF
Expand Down Expand Up @@ -1279,7 +1279,7 @@ subroutine ice_write_hist_attrs(ncid, varid, hfield, ns)
call ice_write_hist_fill(ncid,varid,hfield%vname,history_precision)

! Add cell_methods attribute to variables if averaged
if (hist_avg) then
if (hist_avg .and. .not. write_ic) then
if (TRIM(hfield%vname(1:4))/='sig1' &
.and.TRIM(hfield%vname(1:4))/='sig2' &
.and.TRIM(hfield%vname(1:9))/='sistreave' &
Expand All @@ -1293,6 +1293,7 @@ subroutine ice_write_hist_attrs(ncid, varid, hfield, ns)

if ((histfreq(ns) == '1' .and. histfreq_n(ns) == 1) &
.or..not. hist_avg &
.or. write_ic &
.or.TRIM(hfield%vname(1:4))=='divu' &
.or.TRIM(hfield%vname(1:5))=='shear' &
.or.TRIM(hfield%vname(1:4))=='sig1' &
Expand Down
13 changes: 7 additions & 6 deletions cicecore/cicedyn/infrastructure/io/io_pio2/ice_history_write.F90
Original file line number Diff line number Diff line change
Expand Up @@ -195,7 +195,7 @@ subroutine ice_write_hist (ns)
! define dimensions
!-----------------------------------------------------------------

if (hist_avg) then
if (hist_avg .and. .not. write_ic) then
status = pio_def_dim(File,'d2',2,boundid)
endif

Expand Down Expand Up @@ -233,12 +233,12 @@ subroutine ice_write_hist (ns)
call abort_ice(subname//'ERROR: invalid calendar settings')
endif

if (hist_avg) then
if (hist_avg .and. .not. write_ic) then
status = pio_put_att(File,varid,'bounds','time_bounds')
endif

! Define attributes for time_bounds if hist_avg is true
if (hist_avg) then
if (hist_avg .and. .not. write_ic) then
dimid2(1) = boundid
dimid2(2) = timid
status = pio_def_var(File,'time_bounds',pio_double,dimid2,varid)
Expand Down Expand Up @@ -702,7 +702,7 @@ subroutine ice_write_hist (ns)
! write time_bounds info
!-----------------------------------------------------------------

if (hist_avg) then
if (hist_avg .and. .not. write_ic) then
status = pio_inq_varid(File,'time_bounds',varid)
time_bounds=(/time_beg(ns),time_end(ns)/)
bnd_start = (/1,1/)
Expand Down Expand Up @@ -1219,7 +1219,7 @@ end subroutine ice_write_hist
subroutine ice_write_hist_attrs(File, varid, hfield, ns)

use ice_kinds_mod
use ice_calendar, only: histfreq, histfreq_n
use ice_calendar, only: histfreq, histfreq_n, write_ic
use ice_history_shared, only: ice_hist_field, history_precision, &
hist_avg
use ice_pio
Expand Down Expand Up @@ -1250,7 +1250,7 @@ subroutine ice_write_hist_attrs(File, varid, hfield, ns)
call ice_write_hist_fill(File,varid,hfield%vname,history_precision)

! Add cell_methods attribute to variables if averaged
if (hist_avg) then
if (hist_avg .and. .not. write_ic) then
if (TRIM(hfield%vname(1:4))/='sig1' &
.and.TRIM(hfield%vname(1:4))/='sig2' &
.and.TRIM(hfield%vname(1:9))/='sistreave' &
Expand All @@ -1262,6 +1262,7 @@ subroutine ice_write_hist_attrs(File, varid, hfield, ns)

if ((histfreq(ns) == '1' .and. histfreq_n(ns) == 1) &
.or..not. hist_avg &
.or. write_ic &
.or.TRIM(hfield%vname(1:4))=='divu' &
.or.TRIM(hfield%vname(1:5))=='shear' &
.or.TRIM(hfield%vname(1:4))=='sig1' &
Expand Down
4 changes: 2 additions & 2 deletions cicecore/drivers/direct/nemo_concepts/CICE_InitMod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -185,6 +185,8 @@ subroutine cice_init
if (trim(runtype) == 'continue' .or. restart) &
call init_shortwave ! initialize radiative transfer

if (write_ic) call accum_hist(dt) ! write initial conditions

! determine the time and date at the end of the first timestep
call advance_timestep()

Expand Down Expand Up @@ -215,8 +217,6 @@ subroutine cice_init
call init_flux_atm ! initialize atmosphere fluxes sent to coupler
call init_flux_ocn ! initialize ocean fluxes sent to coupler

if (write_ic) call accum_hist(dt) ! write initial conditions

end subroutine cice_init

!=======================================================================
Expand Down
4 changes: 2 additions & 2 deletions cicecore/drivers/standalone/cice/CICE_InitMod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -199,6 +199,8 @@ subroutine cice_init
if (trim(runtype) == 'continue' .or. restart) &
call init_shortwave ! initialize radiative transfer

if (write_ic) call accum_hist(dt) ! write initial conditions

! tcraig, use advance_timestep here
! istep = istep + 1 ! update time step counters
! istep1 = istep1 + 1
Expand Down Expand Up @@ -243,8 +245,6 @@ subroutine cice_init
call init_flux_atm ! initialize atmosphere fluxes sent to coupler
call init_flux_ocn ! initialize ocean fluxes sent to coupler

if (write_ic) call accum_hist(dt) ! write initial conditions

if (my_task == master_task) then
call ice_memusage_print(nu_diag,subname//':end')
endif
Expand Down

0 comments on commit 9a55ad9

Please sign in to comment.