From 3a4d02e6e22216811964bdf5d302399c3f860232 Mon Sep 17 00:00:00 2001 From: "David A. Bailey" Date: Fri, 16 Dec 2022 17:44:24 -0700 Subject: [PATCH] Fix to driver history (#421) --- configuration/driver/icedrv_history.F90 | 81 ++++++++++++++++++++++--- 1 file changed, 74 insertions(+), 7 deletions(-) diff --git a/configuration/driver/icedrv_history.F90 b/configuration/driver/icedrv_history.F90 index 6c30e5f91..e6abdae32 100644 --- a/configuration/driver/icedrv_history.F90 +++ b/configuration/driver/icedrv_history.F90 @@ -8,7 +8,7 @@ module icedrv_history use icedrv_kinds use icedrv_constants, only: nu_diag, nu_diag_out - use icedrv_domain_size, only: nx, ncat + use icedrv_domain_size, only: nx, ncat, nfsd use icedrv_diagnostics, only: nx_names use icepack_intfc, only: icepack_warnings_flush, icepack_warnings_aborted use icepack_intfc, only: icepack_query_parameters, icepack_query_tracer_sizes @@ -27,7 +27,7 @@ module icedrv_history character (len=char_len_long) :: hist_file ! hist file name integer (kind=int_kind) :: ncid ! cdf file id - integer (kind=int_kind) :: nxid, ncatid, ntrcrid, timid ! cdf dim ids + integer (kind=int_kind) :: nxid, ncatid, ntrcrid, nfsdid, timid ! cdf dim ids integer (kind=int_kind) :: timcnt ! time counter !======================================================================= @@ -50,6 +50,7 @@ subroutine history_write() use icedrv_flux, only: Tair, Qa, fsw, fcondtop use icedrv_flux, only: meltt, meltb, meltl, snoice use icedrv_flux, only: dsnow, congel, sst, sss, Tf, fhocn + use icedrv_arrays_column, only: d_afsd_newi, d_afsd_latg, d_afsd_latm, d_afsd_wave, d_afsd_weld #ifdef USE_NETCDF use netcdf #endif @@ -95,6 +96,14 @@ subroutine history_write() character(len=16), parameter :: fld_3d_ncat(num_3d_ncat) = & (/ 'aicen ', 'vicen ', 'vsnon ' /) + logical (kind=log_kind) :: & + tr_fsd ! flag for tracing fsd + + integer (kind=dbl_kind), parameter :: num_3d_nfsd = 5 + character(len=16), parameter :: fld_3d_nfsd(num_3d_nfsd) = & + (/ 'd_afsd_newi ', 'd_afsd_latg ', 'd_afsd_latm ', & + 'd_afsd_wave ', 'd_afsd_weld ' /) + integer (kind=dbl_kind), parameter :: num_3d_ntrcr = 1 character(len=16), parameter :: fld_3d_ntrcr(num_3d_ntrcr) = & (/ 'trcr ' /) @@ -109,6 +118,7 @@ subroutine history_write() #ifdef USE_NETCDF call icepack_query_tracer_sizes(ntrcr_out=ntrcr) + call icepack_query_tracer_flags(tr_fsd_out=tr_fsd) if (first_call) then timcnt = 0 write(hist_file,'(a,i8.8,a)') './history/icepack.h.',idate,'.nc' @@ -139,6 +149,14 @@ subroutine history_write() status = nf90_def_var(ncid,'ntrcr',NF90_INT,ntrcrid,varid) if (status /= nf90_noerr) call icedrv_system_abort(string=subname//' ERROR: def_var ntrcr') + if (tr_fsd) then + ! nfsd category dimension + status = nf90_def_dim(ncid,'nfsd',nfsd,nfsdid) + if (status /= nf90_noerr) call icedrv_system_abort(string=subname//' ERROR: def_dim nfsd') + status = nf90_def_var(ncid,'nfsd',NF90_INT,nfsdid,varid) + if (status /= nf90_noerr) call icedrv_system_abort(string=subname//' ERROR: def_var nfsd') + endif + ! time dimension status = nf90_def_dim(ncid,'time',NF90_UNLIMITED,timid) if (status /= nf90_noerr) call icedrv_system_abort(string=subname//' ERROR: def_dim time') @@ -189,6 +207,19 @@ subroutine history_write() if (status /= nf90_noerr) call icedrv_system_abort(string=subname//' ERROR in def_var '//trim(fld_3d_ncat(n))) enddo + if (tr_fsd) then + ! 3d nfsd fields + + dimid3(1) = nxid + dimid3(2) = nfsdid + dimid3(3) = timid + + do n = 1,num_3d_nfsd + status = nf90_def_var(ncid,trim(fld_3d_nfsd(n)),NF90_DOUBLE,dimid3,varid) + if (status /= nf90_noerr) call icedrv_system_abort(string=subname//' ERROR in def_var '//trim(fld_3d_nfsd(n))) + enddo + endif + ! 3d ntrcr fields dimid3(1) = nxid @@ -203,8 +234,8 @@ subroutine history_write() ! 4d ncat ntrcr fields dimid4(1) = nxid - dimid4(2) = ncatid - dimid4(3) = ntrcrid + dimid4(2) = ntrcrid + dimid4(3) = ncatid dimid4(4) = timid do n = 1,num_4d_ncat_ntrcr @@ -234,6 +265,13 @@ subroutine history_write() status = nf90_put_var(ncid,varid,(/(n,n=1,ntrcr)/)) if (status /= nf90_noerr) call icedrv_system_abort(string=subname//' ERROR: put_var '//'ntrcr') + if (tr_fsd) then + status = nf90_inq_varid(ncid,'nfsd',varid) + if (status /= nf90_noerr) call icedrv_system_abort(string=subname//' ERROR: inq_var '//'nfsd') + status = nf90_put_var(ncid,varid,(/(n,n=1,nfsd)/)) + if (status /= nf90_noerr) call icedrv_system_abort(string=subname//' ERROR: put_var '//'nfsd') + endif + endif first_call = .false. @@ -336,6 +374,35 @@ subroutine history_write() deallocate(value3) enddo + if (tr_fsd) then + ! 3d nfsd fields + + start3(1) = 1 + count3(1) = nx + start3(2) = 1 + count3(2) = nfsd + start3(3) = timcnt + count3(3) = 1 + + do n = 1,num_3d_nfsd + allocate(value3(count3(1),count3(2),1)) + + value3 = -9999._dbl_kind + if (trim(fld_3d_nfsd(n)) == 'd_afsd_newi') value3(1:count3(1),1:count3(2),1) = d_afsd_newi(1:count3(1),1:count3(2)) + if (trim(fld_3d_nfsd(n)) == 'd_afsd_latg') value3(1:count3(1),1:count3(2),1) = d_afsd_latg(1:count3(1),1:count3(2)) + if (trim(fld_3d_nfsd(n)) == 'd_afsd_latm') value3(1:count3(1),1:count3(2),1) = d_afsd_latm(1:count3(1),1:count3(2)) + if (trim(fld_3d_nfsd(n)) == 'd_afsd_wave') value3(1:count3(1),1:count3(2),1) = d_afsd_wave(1:count3(1),1:count3(2)) + if (trim(fld_3d_nfsd(n)) == 'd_afsd_weld') value3(1:count3(1),1:count3(2),1) = d_afsd_weld(1:count3(1),1:count3(2)) + + status = nf90_inq_varid(ncid,trim(fld_3d_nfsd(n)),varid) + if (status /= nf90_noerr) call icedrv_system_abort(string=subname//' ERROR: inq_var '//trim(fld_3d_nfsd(n))) + status = nf90_put_var(ncid,varid,value3,start=start3,count=count3) + if (status /= nf90_noerr) call icedrv_system_abort(string=subname//' ERROR: put_var '//trim(fld_3d_nfsd(n))) + + deallocate(value3) + enddo + endif + ! 3d ntrcr fields start3(1) = 1 @@ -364,9 +431,9 @@ subroutine history_write() start4(1) = 1 count4(1) = nx start4(2) = 1 - count4(2) = ncat + count4(2) = ntrcr start4(3) = 1 - count4(3) = ntrcr + count4(3) = ncat start4(4) = timcnt count4(4) = 1 @@ -374,7 +441,7 @@ subroutine history_write() allocate(value4(count4(1),count4(2),count4(3),1)) value4 = -9999._dbl_kind - if (trim(fld_4d_ncat_ntrcr(n)) == 'trcr') value4(1:count4(1),1:count4(2),1:count4(3),1) = trcrn(1:count4(1),1:count4(2),1:count4(3)) + if (trim(fld_4d_ncat_ntrcr(n)) == 'trcrn') value4(1:count4(1),1:count4(2),1:count4(3),1) = trcrn(1:count4(1),1:count4(2),1:count4(3)) status = nf90_inq_varid(ncid,trim(fld_4d_ncat_ntrcr(n)),varid) if (status /= nf90_noerr) call icedrv_system_abort(string=subname//' ERROR: inq_var '//trim(fld_4d_ncat_ntrcr(n)))