From 86ab37639737114e7bbc482a775de2ad17ef6095 Mon Sep 17 00:00:00 2001 From: Robert Hallberg Date: Sat, 20 Dec 2025 11:48:13 -0500 Subject: [PATCH] *Correct metadata or calculation of 19 diagnostics Corrected the longname, unit description, conversion factor or the calculation of 19 diagnostics in 4 files. In two cases, the calculations were revised to avoid array syntax calculation of whole arrays, including potentially uninitialized halo regions. In the case of two other diagnostics, necessary unit conversion factors were missing altogether. In yet other cases, the unit conversion factors were moved into the register_SIS_diagnostics conversion factors to facilitate comparison with the declared units. A total of 34 now unnecessary missing_value arguments to related calls to register_SIS_diag_field() were also removed. All solutions are bitwise identical, but the metadata associated with some diagnostics have been changed, and some diagnostics that previously were not passing unit testing are passing now. --- src/SIS_ctrl_types.F90 | 14 ++++++------ src/SIS_ice_diags.F90 | 47 ++++++++++++++++++++++++++++------------- src/SIS_slow_thermo.F90 | 19 ++++++++++------- src/SIS_transport.F90 | 10 ++++----- 4 files changed, 54 insertions(+), 36 deletions(-) diff --git a/src/SIS_ctrl_types.F90 b/src/SIS_ctrl_types.F90 index 041d073..1327e3d 100644 --- a/src/SIS_ctrl_types.F90 +++ b/src/SIS_ctrl_types.F90 @@ -351,13 +351,13 @@ subroutine ice_diagnostics_init(IOF, OSS, FIA, G, US, IG, diag, Time, Cgrid) ! following iceberg diagnostics should be offered. if (associated(IOF%ustar_berg)) & IOF%id_ustar_berg = register_SIS_diag_field('ice_model', 'USTAR_BERG', diag%axesT1, Time, & - 'iceberg ustar', 'm/s', missing_value=missing) + 'iceberg friction velocity', units='m s-1', conversion=1.0) if (associated(IOF%area_berg)) & IOF%id_area_berg = register_SIS_diag_field('ice_model', 'AREA_BERG', diag%axesT1, Time, & - 'icebergs area', 'm2/m2', missing_value=missing) + 'icebergs area per unit ocean area', units='nondim', conversion=1.0) if (associated(IOF%mass_berg)) & IOF%id_mass_berg = register_SIS_diag_field('ice_model', 'MASS_BERG', diag%axesT1, Time, & - 'icebergs mass', 'kg/m2', missing_value=missing) + 'icebergs mass per unit ocean area', units='kg m-2', conversion=1.0) ! Write out static fields. @@ -398,10 +398,10 @@ subroutine ice_diags_fast_init(Rad, G, IG, diag, Time, component) isc = G%isc ; iec = G%iec ; jsc = G%jsc ; jec = G%jec nLay = IG%NkIce - Rad%id_swdn = register_SIS_diag_field(trim(comp_name),'SWDN', diag%axesT1, Time, & - 'downward shortwave flux', 'W/m^2', missing_value=missing) - Rad%id_lwdn = register_SIS_diag_field(trim(comp_name),'LWDN', diag%axesT1, Time, & - 'downward longwave flux', 'W/m^2', missing_value=missing) + Rad%id_swdn = register_SIS_diag_field(trim(comp_name), 'SWDN', diag%axesT1, Time, & + 'downward shortwave flux', units='W m-2', conversion=1.0) + Rad%id_lwdn = register_SIS_diag_field(trim(comp_name), 'LWDN', diag%axesT1, Time, & + 'downward longwave flux', units='W m-2', conversion=1.0) Rad%id_alb = register_SIS_diag_field(trim(comp_name),'ALB',diag%axesT1, Time, & 'surface albedo','0-1', missing_value=missing ) diff --git a/src/SIS_ice_diags.F90 b/src/SIS_ice_diags.F90 index c549b04..2f185e3 100644 --- a/src/SIS_ice_diags.F90 +++ b/src/SIS_ice_diags.F90 @@ -86,6 +86,10 @@ subroutine post_ice_state_diagnostics(IDs, IST, OSS, IOF, dt_slow, Time, G, US, temp_snow ! A diagnostic array with the snow temperature [C ~> degC]. real, dimension(SZI_(G),SZJ_(G),IG%CatIce) :: & rdg_frac ! fraction of ridged ice per category [nondim] + real, dimension(SZI_(G),SZJ_(G),IG%CatIce) :: & + mass_by_cat ! Sea ice mass per unit ocean area by thickness category [R Z ~> kg m-2] + real, dimension(SZI_(G),SZJ_(G),IG%CatIce) :: & + thick_by_cat ! Sea ice thickness by thickness category [Z ~> m] real, dimension(SZI_(G),SZJ_(G)) :: diagVar ! A temporary array for diagnostics. real, dimension(IG%NkIce) :: S_col ! Specified thermodynamic salinity of each ! ice layer if spec_thermo_sal is true [S ~> gSalt kg-1] @@ -137,8 +141,20 @@ subroutine post_ice_state_diagnostics(IDs, IST, OSS, IOF, dt_slow, Time, G, US, ! Thermodynamic state diagnostics ! if (IDs%id_cn>0) call post_data(IDs%id_cn, IST%part_size(:,:,1:ncat), diag) - if (IDs%id_siitdthick>0) call post_data(IDs%id_siitdthick, IST%mH_ice * Spec_vol_ice, diag) - if (IDs%id_simass_n>0) call post_data(IDs%id_simass_n, IST%mH_ice * IST%part_size(:,:,1:ncat), diag) + if (IDs%id_siitdthick>0) then + thick_by_cat(:,:,:) = 0.0 + do k=1,ncat ; do j=jsc,jec ; do i=isc,iec + thick_by_cat(i,j,k) = IST%mH_ice(i,j,k) * Spec_vol_ice + enddo ; enddo ; enddo + call post_data(IDs%id_siitdthick, thick_by_cat, diag) + endif + if (IDs%id_simass_n>0) then + mass_by_cat(:,:,:) = 0.0 + do k=1,ncat ; do j=jsc,jec ; do i=isc,iec + mass_by_cat(i,j,k) = IST%mH_ice(i,j,k) * IST%part_size(i,j,k) + enddo ; enddo ; enddo + call post_data(IDs%id_simass_n, mass_by_cat, diag) + endif if ((IDs%id_siconc>0) .or. (IDs%id_siconc_CMOR>0)) then diagVar(:,:) = 0.0 do j=jsc,jec ; do i=isc,iec ; do k=1,ncat @@ -180,15 +196,15 @@ subroutine post_ice_state_diagnostics(IDs, IST, OSS, IOF, dt_slow, Time, G, US, call post_data(IDs%id_ext, diagVar, diag) endif if (IDs%id_hp>0) call post_avg(IDs%id_hp, IST%mH_pond, IST%part_size(:,:,1:), & ! mw/new - diag, G=G, scale=US%RZ_to_kg_m2/1e3, wtd=.true.) ! rho_water=1e3 + diag, G=G, scale=1.0/(1e3*US%kg_m3_to_R), wtd=.true.) ! rho_water=1e3 [kg m-3] if (IDs%id_hs>0) call post_avg(IDs%id_hs, IST%mH_snow, IST%part_size(:,:,1:), & - diag, G=G, scale=US%Z_to_m/Rho_snow, wtd=.true.) + diag, G=G, scale=1.0/Rho_snow, wtd=.true.) if (IDs%id_sisnthick>0) call post_avg(IDs%id_sisnthick, IST%mH_snow, IST%part_size(:,:,1:), & - diag, G=G, scale=US%Z_to_m/Rho_snow, wtd=.true.) + diag, G=G, scale=1.0/Rho_snow, wtd=.true.) if (IDs%id_hi>0) call post_avg(IDs%id_hi, IST%mH_ice, IST%part_size(:,:,1:), & - diag, G=G, scale=US%Z_to_m/Rho_ice, wtd=.true.) + diag, G=G, scale=1.0/Rho_ice, wtd=.true.) if (IDs%id_sithick>0) call post_avg(IDs%id_sithick, IST%mH_ice, IST%part_size(:,:,1:), & - diag, G=G, scale=US%Z_to_m/Rho_ice, wtd=.true.) + diag, G=G, scale=1.0/Rho_ice, wtd=.true.) if (IDs%id_tsn>0) call post_avg(IDs%id_tsn, temp_snow, IST%part_size(:,:,1:), & diag, G=G, wtd=.true.) if (IDs%id_sitimefrac>0) then @@ -324,17 +340,18 @@ subroutine register_ice_state_diagnostics(Time, IG, US, param_file, diag, IDs) ! Ice state diagnostics. IDs%id_ext = register_diag_field('ice_model', 'EXT', diag%axesT1, Time, & - 'ice modeled', '0 or 1', missing_value=missing) + 'ice extent, indicating cells with more than 15% sea ice cover', & + units='nondim') IDs%id_cn = register_diag_field('ice_model', 'CN', diag%axesTc, Time, & 'ice concentration', '0-1', missing_value=missing) IDs%id_hp = register_diag_field('ice_model', 'HP', diag%axesT1, Time, & - 'pond thickness', 'm-pond', missing_value=missing) ! mw/new + 'pond thickness', units='m', conversion=US%Z_to_m) IDs%id_hs = register_diag_field('ice_model', 'HS', diag%axesT1, Time, & - 'snow thickness', 'm-snow', missing_value=missing) + 'snow thickness', units='m', conversion=US%Z_to_m) IDs%id_tsn = register_diag_field('ice_model', 'TSN', diag%axesT1, Time, & 'snow layer temperature', 'C', conversion=US%C_to_degC, missing_value=missing) IDs%id_hi = register_diag_field('ice_model', 'HI', diag%axesT1, Time, & - 'ice thickness', 'm-ice', missing_value=missing) + 'ice thickness', units='m', conversion=US%Z_to_m) IDs%id_sitimefrac = register_diag_field('ice_model', 'sitimefrac', diag%axesT1, Time, & 'time fraction of ice cover', '0-1', missing_value=missing) IDs%id_siconc = register_diag_field('ice_model', 'siconc', diag%axesT1, Time, & @@ -343,16 +360,16 @@ subroutine register_ice_state_diagnostics(Time, IG, US, param_file, diag, IDs) 'Sea-Ice Area Percentage', '%', missing_value=missing, & standard_name="SeaIceAreaFraction", conversion=100.0) IDs%id_sithick = register_diag_field('ice_model', 'sithick', diag%axesT1, Time, & - 'ice thickness', 'm-ice', missing_value=missing) + 'ice thickness', units='m', conversion=US%Z_to_m) IDs%id_sivol = register_diag_field('ice_model', 'sivol', diag%axesT1, Time, & - 'ice volume', 'm-ice', missing_value=missing) + 'ice volume', units='m', conversion=US%Z_to_m) IDs%id_sisnconc = register_diag_field('ice_model', 'sisnconc', diag%axesT1, Time, & 'snow concentration', '0-1', missing_value=missing) IDs%id_sisnconc_CMOR = register_diag_field('ice_model', 'sisnconc_CMOR', diag%axesT1, Time, & 'Snow Area Percentage', '%', missing_value=missing, & standard_name="SurfaceSnowAreaFraction", conversion=100.0) IDs%id_sisnthick= register_diag_field('ice_model', 'sisnthick', diag%axesT1, Time, & - 'snow thickness', 'm-snow', missing_value=missing) + 'snow thickness', units='m', conversion=US%Z_to_m) IDs%id_t_iceav = register_diag_field('ice_model', 'T_bulkice', diag%axesT1, Time, & 'Volume-averaged ice temperature', 'C', conversion=US%C_to_degC, missing_value=missing) @@ -377,7 +394,7 @@ subroutine register_ice_state_diagnostics(Time, IG, US, param_file, diag, IDs) IDs%id_simass_n = register_diag_field('ice_model', 'simass_n', diag%axesTc, Time, & 'ice mass in categories', 'kg/m^2', conversion=US%RZ_to_kg_m2, missing_value=missing) IDs%id_siitdthick = register_diag_field('ice_model', 'siitdthick', diag%axesTc, Time, & - 'ice thickness in categories', 'm-ice', missing_value=missing) + 'ice thickness in categories', units='m', conversion=US%Z_to_m) IDs%id_sisnmass = register_diag_field('ice_model', 'sisnmass', diag%axesT1, Time, & 'snow mass', 'kg/m^2', conversion=US%RZ_to_kg_m2, missing_value=missing) IDs%id_mib = register_diag_field('ice_model', 'MIB', diag%axesT1, Time, & diff --git a/src/SIS_slow_thermo.F90 b/src/SIS_slow_thermo.F90 index a2d224e..cdbfe3e 100644 --- a/src/SIS_slow_thermo.F90 +++ b/src/SIS_slow_thermo.F90 @@ -580,7 +580,7 @@ subroutine SIS2_thermodynamics(IST, dt_slow, CS, OSS, FIA, IOF, G, US, IG) salt_change, & ! The change in integrated salinity [R Z S ~> gSalt m-2] h2o_change, & ! The change in water in the ice [R Z ~> kg m-2] bsnk, & ! The bottom melting mass sink [R Z T-1 ~> kg m-2 s-1] - tmp2d, & ! A temporary array for mass balance diagnostics [R Z yr-1 ~> kg m-2 yr-1] + tmp2d, & ! A temporary array for mass balance diagnostics [R Z s T-1 yr-1 ~> kg m-2 yr-1] qflx_lim_ice, & ! Ice limiting heat flux [Q R Z T-1 ~> W m-2] qflx_res_ice, & ! Ice restoring heat flux [Q R Z T-1 ~> W m-2] cool_nudge, & ! A heat flux out of the sea ice that @@ -632,8 +632,8 @@ subroutine SIS2_thermodynamics(IST, dt_slow, CS, OSS, FIA, IOF, G, US, IG) real :: rho_ice ! The nominal density of sea ice [R ~> kg m-3]. real :: Idt_slow ! The inverse of the thermodynamic step [T-1 ~> s-1]. - real :: yr_dtslow ! The ratio of 1 year to the thermodynamic time step times some scaling - ! factors, used to change the units of several diagnostics to rate yr-1. + real :: yr_dtslow ! The number of seconds in a year divided by the timestep [s yr-1 T-1 ~> yr-1], + ! used to change the units of several diagnostics to rates in yr-1. real :: heat_to_ocn ! The heat passed from the ice to the ocean [Q R Z ~> J m-2] real :: water_to_ocn ! The water passed to the ocean [R Z ~> kg m-2] real :: salt_to_ocn ! The salt passed to the ocean [R Z S ~> gSalt m-2] @@ -1357,7 +1357,7 @@ subroutine SIS2_thermodynamics(IST, dt_slow, CS, OSS, FIA, IOF, G, US, IG) ! output that has been requested. call enable_SIS_averaging(US%T_to_s*dt_slow, CS%Time, CS%diag) - yr_dtslow = US%RZ_T_to_kg_m2s*(864e2*365*Idt_slow) + yr_dtslow = (864e2*365*Idt_slow) if (CS%id_lsnk>0) then !$OMP parallel do default(shared) do j=jsc,jec ; do i=isc,iec @@ -1380,7 +1380,7 @@ subroutine SIS2_thermodynamics(IST, dt_slow, CS, OSS, FIA, IOF, G, US, IG) scale=Idt_slow, wtd=.true.) if (FIA%id_bheat>0) call post_data(FIA%id_bheat, OSS%bheat, CS%diag) if (CS%id_sn2ic>0) call post_avg(CS%id_sn2ic, snow_to_ice, IST%part_size(:,:,1:), CS%diag, G=G, & - scale=US%RZ_T_to_kg_m2s*Idt_slow) + scale=Idt_slow) if (CS%id_qflim>0) call post_data(CS%id_qflim, qflx_lim_ice, CS%diag) if (CS%id_qfres>0) call post_data(CS%id_qfres, qflx_res_ice, CS%diag) if (CS%id_net_melt>0) call post_data(CS%id_net_melt, net_melt, CS%diag) @@ -1587,15 +1587,18 @@ subroutine SIS_slow_thermo_init(Time, G, US, IG, param_file, diag, CS, tracer_fl endif CS%id_lsrc = register_diag_field('ice_model','LSRC', diag%axesT1, Time, & - 'frozen water local source', 'kg/(m^2*yr)', missing_value=missing) + 'frozen water local source', & + units='kg m-2 yr-1', conversion=US%RZ_T_to_kg_m2s) CS%id_lsnk = register_diag_field('ice_model','LSNK',diag%axesT1, Time, & - 'frozen water local sink', 'kg/(m^2*yr)', missing_value=missing) + 'frozen water local sink', & + units='kg m-2 yr-1', conversion=US%RZ_T_to_kg_m2s) CS%id_bsnk = register_diag_field('ice_model','BSNK',diag%axesT1, Time, & 'frozen water local bottom sink', & 'kg/(m^2*yr)', conversion= 864e2*365.*US%RZ_T_to_kg_m2s, & missing_value=missing) CS%id_sn2ic = register_diag_field('ice_model','SN2IC' ,diag%axesT1,Time, & - 'rate of snow to ice conversion', 'kg/(m^2*s)', missing_value=missing) + 'rate of snow to ice conversion', & + units='kg m-2 s-1', conversion=US%RZ_T_to_kg_m2s) CS%id_net_melt = register_diag_field('ice_model','net_melt' ,diag%axesT1, Time, & 'net mass flux from ice & snow to ocean due to melting & freezing', & 'kg m-2 s-1', conversion=US%RZ_T_to_kg_m2s, missing_value=missing) diff --git a/src/SIS_transport.F90 b/src/SIS_transport.F90 index d928b34..4a56c73 100644 --- a/src/SIS_transport.F90 +++ b/src/SIS_transport.F90 @@ -257,7 +257,7 @@ subroutine finish_ice_transport(CAS, IST, TrReg, G, US, IG, dt, CS, OSS, rdg_rat ! real, dimension(SZI_(G),SZJ_(G)) :: & ! rdg_open, & ! formation rate of open water due to ridging [T-1 ~> s-1] ! rdg_vosh ! rate of ice mass shifted from level to ridged ice [R Z T-1 ~> kg m-2 s-1] - real :: yr_dt ! Tne number of timesteps in a year [nondim]. + real :: yr_dt ! The number of seconds in a year divided by the timestep [s yr-1 T-1 ~> yr-1] real, dimension(SZI_(G),SZJ_(G)) :: trans_conv ! The convergence of frozen water transport [R Z ~> kg m-2]. real, dimension(SZI_(G),SZJ_(G)) :: ice_cover ! The summed fractional ice concentration [nondim]. type(EFP_type) :: tot_ice, tot_snow, enth_ice, enth_snow @@ -379,7 +379,7 @@ subroutine finish_ice_transport(CAS, IST, TrReg, G, US, IG, dt, CS, OSS, rdg_rat ! Calculate and send transport-related diagnostics. Idt = 0.0 ; if (CAS%dt_sum > 0.0) Idt = 1.0 / CAS%dt_sum if (CS%id_xprt>0) then - yr_dt = (8.64e4 * 365.0) * US%s_to_T * Idt + yr_dt = (8.64e4 * 365.0) * Idt call get_cell_mass(IST, G, IG, trans_conv) do j=jsc,jec ; do i=isc,iec trans_conv(i,j) = (trans_conv(i,j) - CAS%mass0(i,j)) * yr_dt @@ -1258,8 +1258,7 @@ subroutine SIS_transport_init(Time, G, IG, US, param_file, diag, CS, continuity_ 'y-direction ice transport', 'kg/s', conversion=US%RZ_T_to_kg_m2s*US%L_to_m**2, & missing_value=missing, interp_method='none') CS%id_xprt = register_diag_field('ice_model', 'XPRT', diag%axesT1, Time, & - 'frozen water transport convergence', 'kg/(m^2*yr)', conversion=US%RZ_to_kg_m2, & - missing_value=missing) + 'frozen water transport convergence', units='kg m-2 yr-1', conversion=US%RZ_T_to_kg_m2s) CS%id_rdgr = register_diag_field('ice_model', 'RDG_RATE', diag%axesT1, Time, & 'ice ridging rate', '1/sec', conversion=US%s_to_T, missing_value=missing) CS%id_rdgh = register_diag_field('ice_model', 'RDG_HEIGHT', diag%axesTc, Time, & @@ -1268,8 +1267,7 @@ subroutine SIS_transport_init(Time, G, IG, US, param_file, diag, CS, continuity_ ! CS%id_rdgo = register_diag_field('ice_model', 'RDG_OPEN', diag%axesT1, Time, & ! 'rate of opening due to ridging', '1/s', conversion=US%s_to_T, missing_value=missing) ! CS%id_rdgv = register_diag_field('ice_model', 'RDG_VOSH', diag%axesT1, Time, & -! 'volume shifted from level to ridged ice', 'm^3/s', conversion=US%RZ_T_to_kg_m2s*US%L_to_m**2, & -! missing_value=missing) +! 'Mass shifted from level to ridged ice', units='kg s-1', conversion=US%RZ_T_to_kg_m2s*US%L_to_m**2) end subroutine SIS_transport_init