diff --git a/CMakeLists.txt b/CMakeLists.txt index 8dedcd26..bea151cd 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -142,6 +142,7 @@ endif() if(CICE_IO MATCHES "PIO" AND NOT TARGET PIO::PIO_Fortran) # Code has not been tested with versions older than 2.5.3, but probably still works fine + # For multiple timesteps per history file, needs latest from main (>=2.6.7) find_package(PIO 2.5.3 REQUIRED COMPONENTS Fortran) endif() if(CICE_IO MATCHES "^(NetCDF|PIO)$" AND NOT TARGET NetCDF::NetCDF_Fortran) diff --git a/io_netcdf/ice_history_write.F90 b/io_netcdf/ice_history_write.F90 index 14e27351..f36ee0fa 100644 --- a/io_netcdf/ice_history_write.F90 +++ b/io_netcdf/ice_history_write.F90 @@ -69,15 +69,6 @@ module ice_history_write contains -subroutine check(status, msg) - integer, intent (in) :: status - character(len=*), intent (in) :: msg - - if(status /= nf90_noerr) then - call abort_ice('ice: NetCDF error '//trim(nf90_strerror(status)//' '//trim(msg))) - end if -end subroutine check - !======================================================================= ! @@ -87,242 +78,367 @@ end subroutine check subroutine ice_write_hist (ns) - use ice_calendar, only: time, sec, idate, idate0, & - dayyr, days_per_year, use_leap_years - use ice_fileunits, only: nu_diag - use ice_restart_shared, only: runid - - integer (kind=int_kind), intent(in) :: ns - - ! local variables - - real (kind=dbl_kind), dimension(:,:), allocatable :: work_g1 - real (kind=real_kind), dimension(:,:), allocatable :: work_gr - real (kind=real_kind), dimension(:,:,:), allocatable :: work_gr3 - real (kind=dbl_kind), dimension(nx_block,ny_block,max_blocks) :: & - work1 - - integer (kind=int_kind) :: i,k,ic,n,nn, & - ncid,status,imtid,jmtid,kmtidi,kmtids,kmtidb, cmtid,timid,varid, & - nvertexid,ivertex - integer (kind=int_kind), dimension(3) :: dimid - integer (kind=int_kind), dimension(4) :: dimidz - integer (kind=int_kind), dimension(5) :: dimidcz - integer (kind=int_kind), dimension(3) :: dimid_nverts - integer (kind=int_kind), dimension(4) :: dimidex - real (kind=real_kind) :: ltime - character (char_len) :: title - character (char_len_long) :: ncfile(max_nstrm) - - integer (kind=int_kind) :: shuffle, deflate, deflate_level - - integer (kind=int_kind) :: ind,boundid - - character (char_len) :: start_time,current_date,current_time - character (len=8) :: cdate - - TYPE(req_attributes), dimension(nvar) :: var - TYPE(coord_attributes), dimension(ncoord) :: coord_var - TYPE(coord_attributes), dimension(nvar_verts) :: var_nverts - TYPE(coord_attributes), dimension(nvarz) :: var_nz - CHARACTER (char_len), dimension(ncoord) :: coord_bounds - - ! We leave shuffle at 0, this is only useful for integer data. - shuffle = 0 - - ! If history_deflate_level < 0 then don't do deflation, - ! otherwise it sets the deflate level - if (history_deflate_level < 0) then - deflate = 0 - deflate_level = 0 - else - deflate = 1 - deflate_level = history_deflate_level - endif - - if (my_task == master_task .or. history_parallel_io) then - - ltime=time/int(secday) + use ice_calendar, only: time, month, daymo + use ice_fileunits, only: nu_diag + + integer (kind=int_kind), intent(in) :: ns !history stream number + + ! local variables + + real (kind=real_kind) :: ltime !history timestamp in days + character (char_len_long) :: ncfile(max_nstrm) !filenames + character (char_len) :: time_string !model time for logging + logical :: file_exists + integer (kind=int_kind) :: & + ncid, & ! netcdf id + varid, & + i_time, & ! time index + timid ! time var id + + type(req_attributes), dimension(nvar) :: var + type(coord_attributes), dimension(ncoord) :: coord_var + type(coord_attributes), dimension(nvar_verts) :: var_nverts + type(coord_attributes), dimension(nvarz) :: var_nz + + if (my_task == master_task .or. history_parallel_io) then + ! set timestamp in middle of time interval + if (histfreq(ns) == 'm' .or. histfreq(ns) == 'M') then + if (month /= 1) then + ltime=time/int(secday)-real(daymo(month-1))/2.0 + else + ltime=time/int(secday)-real(daymo(12))/2.0 + endif + else if(histfreq(ns) == 'd' .or. histfreq(ns) == 'D') then + ltime=time/int(secday) - 0.5 + else + ltime=time/int(secday) + endif - call construct_filename(ncfile(ns),'nc',ns) + call construct_filename(ncfile(ns),'nc',ns,time_string) ! add local directory path name to ncfile if (write_ic) then - ncfile(ns) = trim(incond_dir)//ncfile(ns) + ncfile(ns) = trim(incond_dir)//ncfile(ns) else - ncfile(ns) = trim(history_dir)//ncfile(ns) + ncfile(ns) = trim(history_dir)//ncfile(ns) endif - ! create file - if (history_parallel_io) then - call check(nf90_create(ncfile(ns), ior(NF90_NETCDF4, NF90_MPIIO), ncid, & - comm=MPI_COMM_ICE, info=MPI_INFO_NULL), & - 'create history ncfile '//ncfile(ns)) - if (.not. equal_num_blocks_per_cpu) then - call abort_ice('ice: error history_parallel_io needs equal_num_blocks_per_cpu') - endif + inquire(file=trim(ncfile(ns)),exist=file_exists) + if (.not. file_exists) then + call ice_hist_create(ns, ncfile(ns), ncid, var, coord_var, var_nverts, var_nz) + write(nu_diag,*) 'Created:'//trim(ncfile(ns)) else - call check(nf90_create(ncfile(ns), ior(NF90_CLASSIC_MODEL, NF90_HDF5), ncid), & - 'create history ncfile '//ncfile(ns)) + if (history_parallel_io) then + call check(nf90_open(trim(ncfile(ns)), NF90_WRITE, ncid, & + comm=MPI_COMM_ICE, info=MPI_INFO_NULL), & + 'parallel open existing history file '//ncfile(ns)) + else + call check(nf90_open(trim(ncfile(ns)), NF90_WRITE, ncid), & + "opening existing history file "//ncfile(ns)) + endif endif - !----------------------------------------------------------------- - ! define dimensions - !----------------------------------------------------------------- + !----------------------------------------------------------------- + ! write time variable + !----------------------------------------------------------------- + call check(nf90_inq_dimid(ncid, 'time', timid), & + 'inq dimid time') + call check(nf90_inquire_dimension(ncid, timid, len=i_time), & + 'inquire dim time') + call check(nf90_inq_varid(ncid,'time',varid), & + 'inq varid time') + if (history_parallel_io) then + ! unlimited dimensions need to have collective access set + call check(nf90_var_par_access(ncid, varid, NF90_COLLECTIVE), & + 'parallel access time') + endif + i_time = i_time + 1 ! index of the current history time + call check(nf90_put_var(ncid,varid,ltime,start=(/i_time/)), & + 'put var time') + + !----------------------------------------------------------------- + ! write time_bounds info + !----------------------------------------------------------------- if (hist_avg) then - call check(nf90_def_dim(ncid,'d2',2,boundid), 'def dim d2') + call check(nf90_inq_varid(ncid,'time_bounds',varid), & + 'inq varid time_bounds') + if (history_parallel_io) then + call check(nf90_var_par_access(ncid, varid, NF90_COLLECTIVE), & + 'parallel access time_bounds') + endif + call check(nf90_put_var(ncid,varid,time_beg(ns),start=(/1,i_time/)), & + 'put var time_bounds beginning') + call check(nf90_put_var(ncid,varid,time_end(ns),start=(/2,i_time/)), & + 'put var time_bounds end') endif + endif ! master_task or history_parallel_io + + call broadcast_scalar(i_time, master_task) !we need this on every processor for parallel writes + if (i_time == 1) then + ! these variables are time-invariant, only write once per file + ! ice_hist_create is only run on master task, but these variables are distributed, so call on all tasks + !----------------------------------------------------------------- + ! write coordinate variables + !----------------------------------------------------------------- + if (history_parallel_io) then + call write_coordinate_variables_parallel(ncid, coord_var, var_nz) + else + call write_coordinate_variables(ncid, coord_var, var_nz) + endif + + !----------------------------------------------------------------- + ! write grid masks, area and rotation angle + !----------------------------------------------------------------- + if (history_parallel_io) then + call write_grid_variables_parallel(ncid, var, var_nverts) + else + call write_grid_variables(ncid, var, var_nverts) + endif - call check(nf90_def_dim(ncid, 'ni', nx_global, imtid), & - 'def dim ni') - call check(nf90_def_dim(ncid, 'nj', ny_global, jmtid), & - 'def dim nj') - call check(nf90_def_dim(ncid, 'nc', ncat_hist, cmtid), & - 'def dim nc') - call check(nf90_def_dim(ncid, 'nkice', nzilyr, kmtidi), & - 'def dim nkice') - call check(nf90_def_dim(ncid, 'nksnow', nzslyr, kmtids), & - 'def dim nksnow') - call check(nf90_def_dim(ncid, 'nkbio', nzblyr, kmtidb), & - 'def dim nkbio') - call check(nf90_def_dim(ncid, 'time', 1, timid), & - 'def dim time') - call check(nf90_def_dim(ncid, 'nvertices', nverts, nvertexid), & - 'def dim nverts') + endif + !----------------------------------------------------------------- + ! write variable data + !----------------------------------------------------------------- - !----------------------------------------------------------------- - ! define coordinate variables - !----------------------------------------------------------------- + if (history_parallel_io) then + call write_2d_variables_parallel(ns, ncid, i_time) + else + call write_2d_variables(ns, ncid, i_time) + endif + + if (history_parallel_io) then + call write_3d_and_4d_variables_parallel(ns, ncid, i_time) + else + call write_3d_and_4d_variables(ns, ncid, i_time) + endif - call check(nf90_def_var(ncid,'time',nf90_float,timid,varid), & - 'def var time') - call check(nf90_put_att(ncid,varid,'long_name','model time'), & - 'put_att long_name') + !----------------------------------------------------------------- + ! close output dataset + !----------------------------------------------------------------- - write(cdate,'(i8.8)') idate0 - write(title,'(a,a,a,a,a,a,a,a)') 'days since ', & - cdate(1:4),'-',cdate(5:6),'-',cdate(7:8),' 00:00:00' - call check(nf90_put_att(ncid,varid,'units',title), & - 'put_att time units') - - if (days_per_year == 360) then - call check(nf90_put_att(ncid,varid,'calendar','360_day'), & - 'att time calendar') - elseif (days_per_year == 365 .and. .not.use_leap_years ) then - call check(nf90_put_att(ncid,varid,'calendar','NoLeap'), & - 'att time calendar') - elseif (use_leap_years) then - call check(nf90_put_att(ncid,varid,'calendar','Gregorian'), & - 'att time calendar') - else - call abort_ice( 'ice Error: invalid calendar settings') - endif + if (my_task == master_task .or. history_parallel_io) then + call check(nf90_close(ncid), 'closing netCDF history file') + write(nu_diag,*) 'Wrote ',trim(ncfile(ns)),' at time ',trim(time_string) + endif - if (hist_avg) then - call check(nf90_put_att(ncid,varid,'bounds','time_bounds'), & - 'att time bounds') - endif +end subroutine ice_write_hist - !----------------------------------------------------------------- - ! Define attributes for time bounds if hist_avg is true - !----------------------------------------------------------------- +subroutine ice_hist_create(ns, ncfile, ncid, var, coord_var, var_nverts, var_nz) - if (hist_avg) then - dimid(1) = boundid - dimid(2) = timid - call check(nf90_def_var(ncid, 'time_bounds', & - nf90_float,dimid(1:2),varid), & - 'def var time_bounds') - - call check(nf90_put_att(ncid,varid,'long_name', & - 'boundaries for time-averaging interval'), & - 'att time_bounds long_name') - write(cdate,'(i8.8)') idate0 - write(title,'(a,a,a,a,a,a,a,a)') 'days since ', & - cdate(1:4),'-',cdate(5:6),'-',cdate(7:8),' 00:00:00' - call check(nf90_put_att(ncid,varid,'units',title), & - 'att time_bounds units') - endif + use ice_calendar, only: idate, idate0, & + dayyr, days_per_year, use_leap_years, histfreq_n + use ice_restart_shared, only: runid - !----------------------------------------------------------------- - ! define information for required time-invariant variables - !----------------------------------------------------------------- + integer (kind=int_kind), intent(in) :: ns + character (char_len_long), intent(in) :: ncfile + integer (kind=int_kind), intent(out) :: ncid + TYPE(req_attributes), dimension(nvar), intent(inout) :: var + TYPE(coord_attributes), dimension(ncoord), intent(inout) :: coord_var + TYPE(coord_attributes), dimension(nvar_verts), intent(inout) :: var_nverts + TYPE(coord_attributes), dimension(nvarz), intent(inout) :: var_nz - ind = 0 - ind = ind + 1 - coord_var(ind) = coord_attributes('TLON', & - 'T grid center longitude', 'degrees_east') - coord_bounds(ind) = 'lont_bounds' - ind = ind + 1 - coord_var(ind) = coord_attributes('TLAT', & - 'T grid center latitude', 'degrees_north') - coord_bounds(ind) = 'latt_bounds' - ind = ind + 1 - coord_var(ind) = coord_attributes('ULON', & - 'U grid center longitude', 'degrees_east') - coord_bounds(ind) = 'lonu_bounds' - ind = ind + 1 - coord_var(ind) = coord_attributes('ULAT', & - 'U grid center latitude', 'degrees_north') - coord_bounds(ind) = 'latu_bounds' - - var_nz(1) = coord_attributes('NCAT', 'category maximum thickness', 'm') - var_nz(2) = coord_attributes('VGRDi', 'vertical ice levels', '1') - var_nz(3) = coord_attributes('VGRDs', 'vertical snow levels', '1') - var_nz(4) = coord_attributes('VGRDb', 'vertical ice-bio levels', '1') + ! local variables - !----------------------------------------------------------------- - ! define information for optional time-invariant variables - !----------------------------------------------------------------- + integer (kind=int_kind) :: i, n, status, imtid, jmtid, kmtidi, kmtids, & + kmtidb, cmtid, timid, varid, nvertexid + integer (kind=int_kind), dimension(3) :: dimid, dimid_nverts + integer (kind=int_kind), dimension(4) :: dimidz, dimidex + integer (kind=int_kind), dimension(5) :: dimidcz - var(n_tarea)%req = coord_attributes('tarea', & - 'area of T grid cells', 'm^2') - var(n_tarea)%coordinates = 'TLON TLAT' - var(n_uarea)%req = coord_attributes('uarea', & - 'area of U grid cells', 'm^2') - var(n_uarea)%coordinates = 'ULON ULAT' - var(n_dxt)%req = coord_attributes('dxt', & - 'T cell width through middle', 'm') - var(n_dxt)%coordinates = 'TLON TLAT' - var(n_dyt)%req = coord_attributes('dyt', & - 'T cell height through middle', 'm') - var(n_dyt)%coordinates = 'TLON TLAT' - var(n_dxu)%req = coord_attributes('dxu', & - 'U cell width through middle', 'm') - var(n_dxu)%coordinates = 'ULON ULAT' - var(n_dyu)%req = coord_attributes('dyu', & - 'U cell height through middle', 'm') - var(n_dyu)%coordinates = 'ULON ULAT' - var(n_HTN)%req = coord_attributes('HTN', & - 'T cell width on North side','m') - var(n_HTN)%coordinates = 'TLON TLAT' - var(n_HTE)%req = coord_attributes('HTE', & - 'T cell width on East side', 'm') - var(n_HTE)%coordinates = 'TLON TLAT' - var(n_ANGLE)%req = coord_attributes('ANGLE', & - 'angle grid makes with latitude line on U grid', & - 'radians') - var(n_ANGLE)%coordinates = 'ULON ULAT' - var(n_ANGLET)%req = coord_attributes('ANGLET', & - 'angle grid makes with latitude line on T grid', & - 'radians') - var(n_ANGLET)%coordinates = 'TLON TLAT' - - ! These fields are required for CF compliance - ! dimensions (nx,ny,nverts) - var_nverts(n_lont_bnds) = coord_attributes('lont_bounds', & - 'longitude boundaries of T cells', 'degrees_east') - var_nverts(n_latt_bnds) = coord_attributes('latt_bounds', & - 'latitude boundaries of T cells', 'degrees_north') - var_nverts(n_lonu_bnds) = coord_attributes('lonu_bounds', & - 'longitude boundaries of U cells', 'degrees_east') - var_nverts(n_latu_bnds) = coord_attributes('latu_bounds', & - 'latitude boundaries of U cells', 'degrees_north') + integer (kind=int_kind) :: shuffle, deflate, deflate_level ! compression settings - !----------------------------------------------------------------- - ! define attributes for time-invariant variables - !----------------------------------------------------------------- + integer (kind=int_kind) :: ind,boundid + + character (char_len) :: title, start_time,current_date,current_time,time_period_freq + character (len=8) :: cdate + + + CHARACTER (char_len), dimension(ncoord) :: coord_bounds + + ! We leave shuffle at 0, this is only useful for integer data. + shuffle = 0 + + ! If history_deflate_level < 0 then don't do deflation, + ! otherwise it sets the deflate level + if (history_deflate_level < 0) then + deflate = 0 + deflate_level = 0 + else + deflate = 1 + deflate_level = history_deflate_level + endif + + ! create file + if (history_parallel_io) then + call check(nf90_create(ncfile, ior(NF90_NETCDF4, NF90_MPIIO), ncid, & + comm=MPI_COMM_ICE, info=MPI_INFO_NULL), & + 'create history ncfile '//ncfile) + if (.not. equal_num_blocks_per_cpu) then + call abort_ice('ice: error history_parallel_io needs equal_num_blocks_per_cpu') + endif + else + call check(nf90_create(ncfile, ior(NF90_CLASSIC_MODEL, NF90_HDF5), ncid), & + 'create history ncfile '//ncfile) + endif + + !----------------------------------------------------------------- + ! define dimensions + !----------------------------------------------------------------- + + if (hist_avg) then + call check(nf90_def_dim(ncid,'d2',2,boundid), 'def dim d2') + endif + + call check(nf90_def_dim(ncid, 'ni', nx_global, imtid), & + 'def dim ni') + call check(nf90_def_dim(ncid, 'nj', ny_global, jmtid), & + 'def dim nj') + call check(nf90_def_dim(ncid, 'nc', ncat_hist, cmtid), & + 'def dim nc') + call check(nf90_def_dim(ncid, 'nkice', nzilyr, kmtidi), & + 'def dim nkice') + call check(nf90_def_dim(ncid, 'nksnow', nzslyr, kmtids), & + 'def dim nksnow') + call check(nf90_def_dim(ncid, 'nkbio', nzblyr, kmtidb), & + 'def dim nkbio') + call check(nf90_def_dim(ncid, 'time', NF90_UNLIMITED, timid), & + 'def dim time') + call check(nf90_def_dim(ncid, 'nvertices', nverts, nvertexid), & + 'def dim nverts') + + !----------------------------------------------------------------- + ! define coordinate variables + !----------------------------------------------------------------- + + call check(nf90_def_var(ncid,'time',nf90_float,timid,varid), & + 'def var time') + call check(nf90_put_att(ncid,varid,'long_name','model time'), & + 'put_att long_name') + + write(cdate,'(i8.8)') idate0 + write(title,'(a,a,a,a,a,a,a,a)') 'days since ', & + cdate(1:4),'-',cdate(5:6),'-',cdate(7:8),' 00:00:00' + call check(nf90_put_att(ncid,varid,'units',title), & + 'put_att time units') + + if (days_per_year == 360) then + call check(nf90_put_att(ncid,varid,'calendar','360_day'), & + 'att time calendar') + elseif (days_per_year == 365 .and. .not.use_leap_years ) then + call check(nf90_put_att(ncid,varid,'calendar','NoLeap'), & + 'att time calendar') + elseif (use_leap_years) then + call check(nf90_put_att(ncid,varid,'calendar','proleptic_gregorian'), & + 'att time calendar') + else + call abort_ice( 'ice Error: invalid calendar settings') + endif + + if (hist_avg) then + call check(nf90_put_att(ncid,varid,'bounds','time_bounds'), & + 'att time bounds') + endif + + !----------------------------------------------------------------- + ! Define attributes for time bounds if hist_avg is true + !----------------------------------------------------------------- + + if (hist_avg) then + dimid(1) = boundid + dimid(2) = timid + call check(nf90_def_var(ncid, 'time_bounds', & + nf90_float,dimid(1:2),varid), & + 'def var time_bounds') + + call check(nf90_put_att(ncid,varid,'long_name', & + 'boundaries for time-averaging interval'), & + 'att time_bounds long_name') + write(cdate,'(i8.8)') idate0 + write(title,'(a,a,a,a,a,a,a,a)') 'days since ', & + cdate(1:4),'-',cdate(5:6),'-',cdate(7:8),' 00:00:00' + call check(nf90_put_att(ncid,varid,'units',title), & + 'att time_bounds units') + endif + + !----------------------------------------------------------------- + ! define information for required time-invariant variables + !----------------------------------------------------------------- + + ind = 0 + ind = ind + 1 + coord_var(ind) = coord_attributes('TLON', & + 'T grid center longitude', 'degrees_east') + coord_bounds(ind) = 'lont_bounds' + ind = ind + 1 + coord_var(ind) = coord_attributes('TLAT', & + 'T grid center latitude', 'degrees_north') + coord_bounds(ind) = 'latt_bounds' + ind = ind + 1 + coord_var(ind) = coord_attributes('ULON', & + 'U grid center longitude', 'degrees_east') + coord_bounds(ind) = 'lonu_bounds' + ind = ind + 1 + coord_var(ind) = coord_attributes('ULAT', & + 'U grid center latitude', 'degrees_north') + coord_bounds(ind) = 'latu_bounds' + + var_nz(1) = coord_attributes('NCAT', 'category maximum thickness', 'm') + var_nz(2) = coord_attributes('VGRDi', 'vertical ice levels', '1') + var_nz(3) = coord_attributes('VGRDs', 'vertical snow levels', '1') + var_nz(4) = coord_attributes('VGRDb', 'vertical ice-bio levels', '1') + + !----------------------------------------------------------------- + ! define information for optional time-invariant variables + !----------------------------------------------------------------- + + var(n_tarea)%req = coord_attributes('tarea', & + 'area of T grid cells', 'm^2') + var(n_tarea)%coordinates = 'TLON TLAT' + var(n_uarea)%req = coord_attributes('uarea', & + 'area of U grid cells', 'm^2') + var(n_uarea)%coordinates = 'ULON ULAT' + var(n_dxt)%req = coord_attributes('dxt', & + 'T cell width through middle', 'm') + var(n_dxt)%coordinates = 'TLON TLAT' + var(n_dyt)%req = coord_attributes('dyt', & + 'T cell height through middle', 'm') + var(n_dyt)%coordinates = 'TLON TLAT' + var(n_dxu)%req = coord_attributes('dxu', & + 'U cell width through middle', 'm') + var(n_dxu)%coordinates = 'ULON ULAT' + var(n_dyu)%req = coord_attributes('dyu', & + 'U cell height through middle', 'm') + var(n_dyu)%coordinates = 'ULON ULAT' + var(n_HTN)%req = coord_attributes('HTN', & + 'T cell width on North side','m') + var(n_HTN)%coordinates = 'TLON TLAT' + var(n_HTE)%req = coord_attributes('HTE', & + 'T cell width on East side', 'm') + var(n_HTE)%coordinates = 'TLON TLAT' + var(n_ANGLE)%req = coord_attributes('ANGLE', & + 'angle grid makes with latitude line on U grid', & + 'radians') + var(n_ANGLE)%coordinates = 'ULON ULAT' + var(n_ANGLET)%req = coord_attributes('ANGLET', & + 'angle grid makes with latitude line on T grid', & + 'radians') + var(n_ANGLET)%coordinates = 'TLON TLAT' + + ! These fields are required for CF compliance + ! dimensions (nx,ny,nverts) + var_nverts(n_lont_bnds) = coord_attributes('lont_bounds', & + 'longitude boundaries of T cells', 'degrees_east') + var_nverts(n_latt_bnds) = coord_attributes('latt_bounds', & + 'latitude boundaries of T cells', 'degrees_north') + var_nverts(n_lonu_bnds) = coord_attributes('lonu_bounds', & + 'longitude boundaries of U cells', 'degrees_east') + var_nverts(n_latu_bnds) = coord_attributes('latu_bounds', & + 'latitude boundaries of U cells', 'degrees_north') + + !----------------------------------------------------------------- + ! define attributes for time-invariant variables + !----------------------------------------------------------------- dimid(1) = imtid dimid(2) = jmtid @@ -382,9 +498,9 @@ subroutine ice_write_hist (ns) deflate_level), & 'deflate '//var_nz(i)%short_name) - call check(nf90_put_att(ncid,varid,'long_name',var_nz(i)%long_name), & + call check(nf90_put_att(ncid,varid,'long_name',var_nz(i)%long_name), & 'put att long_name '//var_nz(i)%short_name) - call check(nf90_put_att(ncid, varid, 'units', var_nz(i)%units), & + call check(nf90_put_att(ncid, varid, 'units', var_nz(i)%units), & 'for att units '//var_nz(i)%short_name) endif enddo @@ -530,6 +646,11 @@ subroutine ice_write_hist (ns) call check(nf90_put_att(ncid,varid,'cell_measures', & avail_hist_fields(n)%vcellmeas), & 'put att cell_measures '//avail_hist_fields(n)%vname) + if (avail_hist_fields(n)%vcomment /= "none") then + call check(nf90_put_att(ncid,varid,'comment', & + avail_hist_fields(n)%vcomment), & + 'put att comment '//avail_hist_fields(n)%vname) + endif call check(nf90_put_att(ncid,varid,'missing_value',spval), & 'put att missing_value '//avail_hist_fields(n)%vname) call check(nf90_put_att(ncid,varid,'_FillValue',spval), & @@ -541,7 +662,6 @@ subroutine ice_write_hist (ns) if (hist_avg) then if (TRIM(avail_hist_fields(n)%vname)/='sig1' .or. & TRIM(avail_hist_fields(n)%vname)/='sig2') then - call check(nf90_put_att(ncid,varid,'cell_methods','time: mean'), & 'put att cell methods time mean '//avail_hist_fields(n)%vname) endif @@ -601,6 +721,11 @@ subroutine ice_write_hist (ns) avail_hist_fields(n)%vcellmeas) if (status /= nf90_noerr) call abort_ice( & 'Error defining cell measures for '//avail_hist_fields(n)%vname) + if (avail_hist_fields(n)%vcomment /= "none") then + call check(nf90_put_att(ncid,varid,'comment', & + avail_hist_fields(n)%vcomment), & + 'put att comment '//avail_hist_fields(n)%vname) + endif status = nf90_put_att(ncid,varid,'missing_value',spval) if (status /= nf90_noerr) call abort_ice( & 'Error defining missing_value for '//avail_hist_fields(n)%vname) @@ -664,6 +789,11 @@ subroutine ice_write_hist (ns) avail_hist_fields(n)%vcellmeas) if (status /= nf90_noerr) call abort_ice( & 'Error defining cell measures for '//avail_hist_fields(n)%vname) + if (avail_hist_fields(n)%vcomment /= "none") then + call check(nf90_put_att(ncid,varid,'comment', & + avail_hist_fields(n)%vcomment), & + 'put att comment '//avail_hist_fields(n)%vname) + endif status = nf90_put_att(ncid,varid,'missing_value',spval) if (status /= nf90_noerr) call abort_ice( & 'Error defining missing_value for '//avail_hist_fields(n)%vname) @@ -713,6 +843,11 @@ subroutine ice_write_hist (ns) avail_hist_fields(n)%vcellmeas) if (status /= nf90_noerr) call abort_ice( & 'Error defining cell measures for '//avail_hist_fields(n)%vname) + if (avail_hist_fields(n)%vcomment /= "none") then + call check(nf90_put_att(ncid,varid,'comment', & + avail_hist_fields(n)%vcomment), & + 'put att comment '//avail_hist_fields(n)%vname) + endif status = nf90_put_att(ncid,varid,'missing_value',spval) if (status /= nf90_noerr) call abort_ice( & 'Error defining missing_value for '//avail_hist_fields(n)%vname) @@ -763,6 +898,11 @@ subroutine ice_write_hist (ns) avail_hist_fields(n)%vcellmeas) if (status /= nf90_noerr) call abort_ice( & 'Error defining cell measures for '//avail_hist_fields(n)%vname) + if (avail_hist_fields(n)%vcomment /= "none") then + call check(nf90_put_att(ncid,varid,'comment', & + avail_hist_fields(n)%vcomment), & + 'put att comment '//avail_hist_fields(n)%vname) + endif status = nf90_put_att(ncid,varid,'missing_value',spval) if (status /= nf90_noerr) call abort_ice( & 'Error defining missing_value for '//avail_hist_fields(n)%vname) @@ -828,6 +968,11 @@ subroutine ice_write_hist (ns) avail_hist_fields(n)%vcellmeas) if (status /= nf90_noerr) call abort_ice( & 'Error defining cell measures for '//avail_hist_fields(n)%vname) + if (avail_hist_fields(n)%vcomment /= "none") then + call check(nf90_put_att(ncid,varid,'comment', & + avail_hist_fields(n)%vcomment), & + 'put att comment '//avail_hist_fields(n)%vname) + endif status = nf90_put_att(ncid,varid,'missing_value',spval) if (status /= nf90_noerr) call abort_ice( & 'Error defining missing_value for '//avail_hist_fields(n)%vname) @@ -893,6 +1038,11 @@ subroutine ice_write_hist (ns) avail_hist_fields(n)%vcellmeas) if (status /= nf90_noerr) call abort_ice( & 'Error defining cell measures for '//avail_hist_fields(n)%vname) + if (avail_hist_fields(n)%vcomment /= "none") then + call check(nf90_put_att(ncid,varid,'comment', & + avail_hist_fields(n)%vcomment), & + 'put att comment '//avail_hist_fields(n)%vname) + endif status = nf90_put_att(ncid,varid,'missing_value',spval) if (status /= nf90_noerr) call abort_ice( & 'Error defining missing_value for '//avail_hist_fields(n)%vname) @@ -939,8 +1089,25 @@ subroutine ice_write_hist (ns) title = 'Los Alamos Sea Ice Model (CICE) Version 5' call check(nf90_put_att(ncid,nf90_global,'source',title), & 'global attribute source') - -#ifdef AusCOM + + select case (histfreq(ns)) + case ("y", "Y") + write(time_period_freq,'(a,i0)') 'year_',histfreq_n(ns) + case ("m", "M") + write(time_period_freq,'(a,i0)') 'month_',histfreq_n(ns) + case ("d", "D") + write(time_period_freq,'(a,i0)') 'day_',histfreq_n(ns) + case ("h", "H") + write(time_period_freq,'(a,i0)') 'hour_',histfreq_n(ns) + case ("1") + write(time_period_freq,'(a,i0)') 'step_',histfreq_n(ns) + end select + + call check(nf90_put_att(ncid,nf90_global,'time_period_freq',trim(time_period_freq)),& + 'global attribute time_period_freq') + + +#if defined(AUSCOM) && !defined(ACCESS) write(title,'(a,i3,a)') 'This Year Has ',int(dayyr),' days' #else if (use_leap_years) then @@ -952,17 +1119,14 @@ subroutine ice_write_hist (ns) call check(nf90_put_att(ncid,nf90_global,'comment',title), & 'global attribute comment') - write(title,'(a,i8.8)') 'File written on model date ',idate + write(title,'(a,i8.8)') 'File started on model date ',idate call check(nf90_put_att(ncid,nf90_global,'comment2',title), & - 'global attribute date1') - - write(title,'(a,i6)') 'seconds elapsed into model date: ',sec - call check(nf90_put_att(ncid,nf90_global,'comment3',title), & - 'global attribute date2') + 'global attribute comment2') - title = 'CF-1.0' - call check(nf90_put_att(ncid,nf90_global,'conventions',title), & - 'global attribute conventions') + ! TO-DO: Update output for CF compliance ! + ! title = 'CF-1.0' + ! call check(nf90_put_att(ncid,nf90_global,'conventions',title), & + ! 'global attribute conventions') call date_and_time(date=current_date, time=current_time) write(start_time,1000) current_date(1:4), current_date(5:6), & @@ -983,78 +1147,9 @@ subroutine ice_write_hist (ns) call check(nf90_enddef(ncid), 'enddef') - !----------------------------------------------------------------- - ! write time variable - !----------------------------------------------------------------- - - call check(nf90_inq_varid(ncid,'time',varid), & - 'inq varid time') - call check(nf90_put_var(ncid,varid,ltime), & - 'put var time') - - !----------------------------------------------------------------- - ! write time_bounds info - !----------------------------------------------------------------- - - if (hist_avg) then - call check(nf90_inq_varid(ncid,'time_bounds',varid), & - 'inq varid time_bounds') - call check(nf90_put_var(ncid,varid,time_beg(ns),start=(/1/)), & - 'put var time_bounds beginning') - call check(nf90_put_var(ncid,varid,time_end(ns),start=(/2/)), & - 'put var time_bounds end') - endif - endif ! master_task or history_parallel_io - - !----------------------------------------------------------------- - ! write coordinate variables - !----------------------------------------------------------------- - - if (history_parallel_io) then - call write_coordinate_variables_parallel(ncid, coord_var, var_nz) - else - call write_coordinate_variables(ncid, coord_var, var_nz) - endif - - !----------------------------------------------------------------- - ! write grid masks, area and rotation angle - !----------------------------------------------------------------- - - if (history_parallel_io) then - call write_grid_variables_parallel(ncid, var, var_nverts) - else - call write_grid_variables(ncid, var, var_nverts) - endif +end subroutine ice_hist_create - !----------------------------------------------------------------- - ! write 2d variable data - !----------------------------------------------------------------- - - if (history_parallel_io) then - call write_2d_variables_parallel(ns, ncid) - else - call write_2d_variables(ns, ncid) - endif - - if (history_parallel_io) then - call write_3d_and_4d_variables_parallel(ns, ncid) - else - call write_3d_and_4d_variables(ns, ncid) - endif - - !----------------------------------------------------------------- - ! close output dataset - !----------------------------------------------------------------- - - if (my_task == master_task .or. history_parallel_io) then - call check(nf90_close(ncid), 'closing netCDF history file') - write(nu_diag,*) ' ' - write(nu_diag,*) 'Finished writing ',trim(ncfile(ns)) - endif - -end subroutine ice_write_hist - subroutine write_coordinate_variables(ncid, coord_var, var_nz) integer, intent(in) :: ncid @@ -1080,7 +1175,7 @@ subroutine write_coordinate_variables(ncid, coord_var, var_nz) work_g1(:,:) = c0 do i = 1,ncoord - coord_var_name = coord_var(i)%short_name + if (my_task == master_task) coord_var_name = coord_var(i)%short_name call broadcast_scalar(coord_var_name, master_task) SELECT CASE (coord_var_name) CASE ('TLON') @@ -1110,28 +1205,27 @@ subroutine write_coordinate_variables(ncid, coord_var, var_nz) enddo ! Extra dimensions (NCAT, VGRD*) - - do i = 1, nvarz - if (igrdz(i)) then - call broadcast_scalar(var_nz(i)%short_name,master_task) - if (my_task == master_task) then - call check(nf90_inq_varid(ncid, var_nz(i)%short_name, varid), & - 'inq varid '//var_nz(i)%short_name) - SELECT CASE (var_nz(i)%short_name) - CASE ('NCAT') - status = nf90_put_var(ncid,varid,hin_max(1:ncat_hist)) - CASE ('VGRDi') ! index - needed for Met Office analysis code - status = nf90_put_var(ncid,varid,(/(k, k=1,nzilyr)/)) - CASE ('VGRDs') ! index - needed for Met Office analysis code - status = nf90_put_var(ncid,varid,(/(k, k=1,nzslyr)/)) - CASE ('VGRDb') - status = nf90_put_var(ncid,varid,(/(k, k=1,nzblyr)/)) - END SELECT - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error writing'//var_nz(i)%short_name) - endif - endif - enddo + if (my_task == master_task) then + do i = 1, nvarz + coord_var_name = var_nz(i)%short_name + if (igrdz(i)) then + call check(nf90_inq_varid(ncid, coord_var_name, varid), & + 'inq varid '//coord_var_name) + SELECT CASE (coord_var_name) + CASE ('NCAT') + status = nf90_put_var(ncid,varid,hin_max(1:ncat_hist)) + CASE ('VGRDi') ! index - needed for Met Office analysis code + status = nf90_put_var(ncid,varid,(/(k, k=1,nzilyr)/)) + CASE ('VGRDs') ! index - needed for Met Office analysis code + status = nf90_put_var(ncid,varid,(/(k, k=1,nzslyr)/)) + CASE ('VGRDb') + status = nf90_put_var(ncid,varid,(/(k, k=1,nzblyr)/)) + END SELECT + if (status /= nf90_noerr) call abort_ice( & + 'ice: Error writing'//coord_var_name) + endif + enddo + endif deallocate(work_g1) deallocate(work_gr) @@ -1139,7 +1233,6 @@ subroutine write_coordinate_variables(ncid, coord_var, var_nz) end subroutine write_coordinate_variables - subroutine write_grid_variables(ncid, var, var_nverts) integer, intent(in) :: ncid @@ -1194,8 +1287,7 @@ subroutine write_grid_variables(ncid, var, var_nverts) do i = 3, nvar ! note n_tmask=1, n_blkmask=2 if (igrd(i)) then - var_name = var(i)%req%short_name - + if (my_task == master_task) var_name = var(i)%req%short_name call broadcast_scalar(var_name,master_task) SELECT CASE (var_name) CASE ('tarea') @@ -1239,7 +1331,7 @@ subroutine write_grid_variables(ncid, var, var_nverts) work1 (:,:,:) = c0 do i = 1, nvar_verts - var_nverts_name = var_nverts(i)%short_name + if (my_task == master_task) var_nverts_name = var_nverts(i)%short_name call broadcast_scalar(var_nverts_name,master_task) SELECT CASE (var_nverts_name) CASE ('lont_bounds') @@ -1284,10 +1376,9 @@ subroutine write_grid_variables(ncid, var, var_nverts) end subroutine write_grid_variables -subroutine write_2d_variables(ns, ncid) +subroutine write_2d_variables(ns, ncid, i_time) - integer, intent(in) :: ns - integer, intent(in) :: ncid + integer, intent(in) :: ns, ncid, i_time real (kind=dbl_kind), dimension(:,:), allocatable :: work_g1 real (kind=real_kind), dimension(:,:), allocatable :: work_gr @@ -1315,6 +1406,7 @@ subroutine write_2d_variables(ns, ncid) call check(nf90_inq_varid(ncid,avail_hist_fields(n)%vname,varid), & 'inq varid '//avail_hist_fields(n)%vname) call check(nf90_put_var(ncid,varid,work_gr(:,:), & + start=(/1,1,i_time/), & count=(/nx_global,ny_global/)), & 'put var '//avail_hist_fields(n)%vname) endif @@ -1327,10 +1419,9 @@ subroutine write_2d_variables(ns, ncid) end subroutine write_2d_variables -subroutine write_3d_and_4d_variables(ns, ncid) +subroutine write_3d_and_4d_variables(ns, ncid, i_time) - integer, intent(in) :: ns - integer, intent(in) :: ncid + integer, intent(in) :: ns, ncid, i_time real (kind=dbl_kind), dimension(:,:), allocatable :: work_g1 real (kind=real_kind), dimension(:,:), allocatable :: work_gr @@ -1363,7 +1454,7 @@ subroutine write_3d_and_4d_variables(ns, ncid) if (my_task == master_task) then call check(nf90_put_var(ncid,varid,work_gr(:,:), & - start=(/ 1, 1, k/), & + start=(/1,1,k,i_time/), & count=(/nx_global,ny_global, 1/)), & 'put var '//avail_hist_fields(n)%vname) endif @@ -1388,7 +1479,7 @@ subroutine write_3d_and_4d_variables(ns, ncid) if (my_task == master_task) then call check(nf90_put_var(ncid,varid,work_gr(:,:), & - start=(/ 1, 1,k/), & + start=(/1,1,k,i_time/), & count=(/nx_global,ny_global,1/)), & 'put var '//avail_hist_fields(n)%vname) endif @@ -1413,7 +1504,7 @@ subroutine write_3d_and_4d_variables(ns, ncid) if (my_task == master_task) then call check(nf90_put_var(ncid,varid,work_gr(:,:), & - start=(/ 1, 1,k/), & + start=(/1,1,k,i_time/), & count=(/nx_global,ny_global,1/)), & 'put var '//avail_hist_fields(n)%vname) endif @@ -1438,7 +1529,7 @@ subroutine write_3d_and_4d_variables(ns, ncid) work_gr(:,:) = work_g1(:,:) if (my_task == master_task) then call check(nf90_put_var(ncid,varid,work_gr(:,:), & - start=(/ 1, 1,k,ic/), & + start=(/1,1,k,ic,i_time/), & count=(/nx_global,ny_global,1, 1/)), & 'put var '//avail_hist_fields(n)%vname) endif @@ -1464,7 +1555,7 @@ subroutine write_3d_and_4d_variables(ns, ncid) work_gr(:,:) = work_g1(:,:) if (my_task == master_task) then call check(nf90_put_var(ncid,varid,work_gr(:,:), & - start=(/ 1, 1,k,ic/), & + start=(/1,1,k,ic,i_time/), & count=(/nx_global,ny_global,1, 1/)), & 'put var '//avail_hist_fields(n)%vname) endif @@ -1490,7 +1581,7 @@ subroutine write_3d_and_4d_variables(ns, ncid) work_gr(:,:) = work_g1(:,:) if (my_task == master_task) then call check(nf90_put_var(ncid,varid,work_gr(:,:), & - start=(/ 1, 1,k,ic/), & + start=(/1,1,k,ic,i_time/), & count=(/nx_global,ny_global,1, 1/)), & 'put var '//avail_hist_fields(n)%vname) endif @@ -1530,7 +1621,7 @@ subroutine write_coordinate_variables_parallel(ncid, coord_var, var_nz) work1 = ULAT*rad_to_deg END SELECT - call put_2d_with_blocks(ncid, coord_var(i)%short_name, work1) + call put_2d_with_blocks(ncid, 1, coord_var(i)%short_name, work1) enddo ! Extra dimensions (NCAT, VGRD*) @@ -1538,6 +1629,8 @@ subroutine write_coordinate_variables_parallel(ncid, coord_var, var_nz) if (igrdz(i)) then call check(nf90_inq_varid(ncid, var_nz(i)%short_name, varid), & 'inq_varid '//var_nz(i)%short_name) + call check(nf90_var_par_access(ncid, varid, NF90_COLLECTIVE), & + 'parallel access '//var_nz(i)%short_name) SELECT CASE (var_nz(i)%short_name) CASE ('NCAT') call check(nf90_put_var(ncid, varid, hin_max(1:ncat_hist)), & @@ -1576,11 +1669,11 @@ subroutine write_grid_variables_parallel(ncid, var, var_nverts) integer :: varid if (igrd(n_tmask)) then - call put_2d_with_blocks(ncid, 'tmask', hm) + call put_2d_with_blocks(ncid, 1, 'tmask', hm) endif if (igrd(n_blkmask)) then - call put_2d_with_blocks(ncid, 'blkmask', bm) + call put_2d_with_blocks(ncid, 1, 'blkmask', bm) endif do i = 3, nvar ! note n_tmask=1, n_blkmask=2 @@ -1608,7 +1701,7 @@ subroutine write_grid_variables_parallel(ncid, var, var_nverts) work1 = ANGLET END SELECT - call put_2d_with_blocks(ncid, var(i)%req%short_name, work1) + call put_2d_with_blocks(ncid, 1, var(i)%req%short_name, work1) endif enddo @@ -1626,11 +1719,13 @@ subroutine write_grid_variables_parallel(ncid, var, var_nverts) CASE ('lonu_bounds') work2(:, :, :, :) = lonu_bounds(:, :, :, :) CASE ('latu_bounds') - work2(:, :, :, :) = lonu_bounds(:, :, :, :) + work2(:, :, :, :) = latu_bounds(:, :, :, :) END SELECT call check(nf90_inq_varid(ncid, var_nverts(i)%short_name, varid), & 'inq varid '//var_nverts(i)%short_name) + call check(nf90_var_par_access(ncid, varid, NF90_COLLECTIVE), & + 'parallel access '//var_nverts(i)%short_name) do iblk=1, nblocks the_block = get_block(blocks_ice(iblk), iblk) @@ -1657,17 +1752,16 @@ subroutine write_grid_variables_parallel(ncid, var, var_nverts) end subroutine write_grid_variables_parallel -subroutine write_2d_variables_parallel(ns, ncid) +subroutine write_2d_variables_parallel(ns, ncid, i_time) - integer, intent(in) :: ns - integer, intent(in) :: ncid + integer, intent(in) :: ns, ncid, i_time integer :: varid integer :: n do n=1, num_avail_hist_fields_2D if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then - call put_2d_with_blocks(ncid, avail_hist_fields(n)%vname, & + call put_2d_with_blocks(ncid, i_time, avail_hist_fields(n)%vname, & a2D(:, :, n, :)) endif enddo ! num_avail_hist_fields_2D @@ -1676,10 +1770,9 @@ end subroutine write_2d_variables_parallel -subroutine write_3d_and_4d_variables_parallel(ns, ncid) +subroutine write_3d_and_4d_variables_parallel(ns, ncid, i_time) - integer, intent(in) :: ns - integer, intent(in) :: ncid + integer, intent(in) :: ns, ncid, i_time integer :: varid integer :: n, nn, k, ic @@ -1687,7 +1780,7 @@ subroutine write_3d_and_4d_variables_parallel(ns, ncid) do n = n2D + 1, n3Dccum nn = n - n2D if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then - call put_3d_with_blocks(ncid, avail_hist_fields(n)%vname, & + call put_3d_with_blocks(ncid, i_time, avail_hist_fields(n)%vname, & ncat_hist, a3Dc(:, :, :, nn, :)) endif enddo ! num_avail_hist_fields_3Dc @@ -1696,7 +1789,7 @@ subroutine write_3d_and_4d_variables_parallel(ns, ncid) do n = n3Dccum+1, n3Dzcum nn = n - n3Dccum if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then - call put_3d_with_blocks(ncid, avail_hist_fields(n)%vname, & + call put_3d_with_blocks(ncid, i_time, avail_hist_fields(n)%vname, & nzilyr, a3Dz(:, :, :, nn, :)) endif enddo ! num_avail_hist_fields_3Dz @@ -1705,7 +1798,7 @@ subroutine write_3d_and_4d_variables_parallel(ns, ncid) do n = n3Dzcum+1, n3Dbcum nn = n - n3Dzcum if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then - call put_3d_with_blocks(ncid, avail_hist_fields(n)%vname, & + call put_3d_with_blocks(ncid, i_time, avail_hist_fields(n)%vname, & nzblyr, a3Db(:, :, :, nn, :)) endif enddo ! num_avail_hist_fields_3Db @@ -1714,7 +1807,7 @@ subroutine write_3d_and_4d_variables_parallel(ns, ncid) do n = n3Dbcum+1, n4Dicum nn = n - n3Dbcum if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then - call put_4d_with_blocks(ncid, avail_hist_fields(n)%vname, & + call put_4d_with_blocks(ncid, i_time, avail_hist_fields(n)%vname, & nzilyr, ncat_hist, a4Di(:, :, :, :, nn, :)) endif enddo ! num_avail_hist_fields_4Di @@ -1723,7 +1816,7 @@ subroutine write_3d_and_4d_variables_parallel(ns, ncid) do n = n4Dicum+1, n4Dscum nn = n - n4Dicum if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then - call put_4d_with_blocks(ncid, avail_hist_fields(n)%vname, & + call put_4d_with_blocks(ncid, i_time, avail_hist_fields(n)%vname, & nzslyr, ncat_hist, a4Ds(:, :, :, :, nn, :)) endif enddo ! num_avail_hist_fields_4Ds @@ -1731,7 +1824,7 @@ subroutine write_3d_and_4d_variables_parallel(ns, ncid) do n = n4Dscum+1, n4Dbcum nn = n - n4Dscum if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then - call put_4d_with_blocks(ncid, avail_hist_fields(n)%vname, & + call put_4d_with_blocks(ncid, i_time, avail_hist_fields(n)%vname, & nzblyr, ncat_hist, a4Db(:, :, :, :, nn, :)) endif enddo ! num_avail_hist_fields_4Db @@ -1739,20 +1832,26 @@ subroutine write_3d_and_4d_variables_parallel(ns, ncid) end subroutine write_3d_and_4d_variables_parallel -subroutine put_2d_with_blocks(ncid, var_name, data) +subroutine put_2d_with_blocks(ncid, i_start, var_name, data) - integer, intent(in) :: ncid + ! by convention only, 2d variables are actually 3d if you consider time + ! typically i_start is the current time index, but can be set to 1 to write + ! time-invarient (e.g. grid data) + + integer, intent(in) :: ncid, i_start character(len=*), intent(in) :: var_name real(kind=dbl_kind), dimension(nx_block, ny_block, max_blocks), intent(in) :: data integer :: varid integer :: iblk integer :: ilo, jlo, ihi, jhi, gilo, gjlo, gihi, gjhi - integer, dimension(2) :: start, count + integer, dimension(3) :: start, count type(block) :: the_block call check(nf90_inq_varid(ncid, var_name, varid), & 'inq varid for '//var_name) + call check(nf90_var_par_access(ncid, varid, NF90_COLLECTIVE), & + 'parallel access '//var_name) do iblk=1, nblocks the_block = get_block(blocks_ice(iblk), iblk) @@ -1766,8 +1865,8 @@ subroutine put_2d_with_blocks(ncid, var_name, data) gihi = the_block%i_glob(ihi) gjhi = the_block%j_glob(jhi) - start = (/ gilo, gjlo /) - count = (/ gihi - gilo + 1, gjhi - gjlo + 1 /) + start = (/ gilo, gjlo,i_start /) + count = (/ gihi - gilo + 1, gjhi - gjlo + 1 , 1/) call check(nf90_put_var(ncid, varid, real(data(ilo:ihi, jlo:jhi, iblk)), & start=start, count=count), & 'put_2d_with_blocks put '//trim(var_name)) @@ -1775,20 +1874,25 @@ subroutine put_2d_with_blocks(ncid, var_name, data) end subroutine put_2d_with_blocks -subroutine put_3d_with_blocks(ncid, var_name, len_3dim, data) +subroutine put_3d_with_blocks(ncid, i_time, var_name, len_3dim, data) + + ! by convention only, 3d variables are actually 4d if you consider time - integer, intent(in) :: ncid, len_3dim + + integer, intent(in) :: ncid, i_time, len_3dim character(len=*), intent(in) :: var_name real(kind=dbl_kind), dimension(nx_block, ny_block, len_3dim, max_blocks), intent(in) :: data integer :: varid integer :: iblk integer :: ilo, jlo, ihi, jhi, gilo, gjlo, gihi, gjhi - integer, dimension(3) :: start, count + integer, dimension(4) :: start, count type(block) :: the_block call check(nf90_inq_varid(ncid, var_name, varid), & 'inq varid for '//var_name) + call check(nf90_var_par_access(ncid, varid, NF90_COLLECTIVE), & + 'parallel access '//var_name) do iblk=1, nblocks the_block = get_block(blocks_ice(iblk), iblk) @@ -1802,8 +1906,8 @@ subroutine put_3d_with_blocks(ncid, var_name, len_3dim, data) gihi = the_block%i_glob(ihi) gjhi = the_block%j_glob(jhi) - start = (/ gilo, gjlo, 1 /) - count = (/ gihi - gilo + 1, gjhi - gjlo + 1, len_3dim /) + start = (/ gilo, gjlo, 1 , i_time/) + count = (/ gihi - gilo + 1, gjhi - gjlo + 1, len_3dim, 1 /) call check(nf90_put_var(ncid, varid, & real(data(ilo:ihi, jlo:jhi, 1:len_3dim, iblk)), & start=start, count=count), & @@ -1813,9 +1917,11 @@ subroutine put_3d_with_blocks(ncid, var_name, len_3dim, data) end subroutine put_3d_with_blocks -subroutine put_4d_with_blocks(ncid, var_name, len_3dim, len_4dim, data) +subroutine put_4d_with_blocks(ncid, i_time, var_name, len_3dim, len_4dim, data) + ! by convention only, 4d variables are actually 5d if you consider time - integer, intent(in) :: ncid, len_3dim, len_4dim + + integer, intent(in) :: ncid, i_time, len_3dim, len_4dim character(len=*), intent(in) :: var_name real(kind=dbl_kind), dimension(nx_block, ny_block, len_3dim, & len_4dim, max_blocks), intent(in) :: data @@ -1823,11 +1929,13 @@ subroutine put_4d_with_blocks(ncid, var_name, len_3dim, len_4dim, data) integer :: varid integer :: iblk integer :: ilo, jlo, ihi, jhi, gilo, gjlo, gihi, gjhi - integer, dimension(4) :: start, count + integer, dimension(5) :: start, count type(block) :: the_block call check(nf90_inq_varid(ncid, var_name, varid), & 'inq varid for '//var_name) + call check(nf90_var_par_access(ncid, varid, NF90_COLLECTIVE), & + 'parallel access '//var_name) do iblk=1, nblocks the_block = get_block(blocks_ice(iblk), iblk) @@ -1841,8 +1949,8 @@ subroutine put_4d_with_blocks(ncid, var_name, len_3dim, len_4dim, data) gihi = the_block%i_glob(ihi) gjhi = the_block%j_glob(jhi) - start = (/ gilo, gjlo, 1, 1 /) - count = (/ gihi - gilo + 1, gjhi - gjlo + 1, len_3dim, len_4dim /) + start = (/ gilo, gjlo, 1, 1 , i_time/) + count = (/ gihi - gilo + 1, gjhi - gjlo + 1, len_3dim, len_4dim, 1 /) call check(nf90_put_var(ncid, varid, & real(data(ilo:ihi, jlo:jhi, 1:len_3dim, 1:len_4dim, iblk)), & start=start, count=count), & diff --git a/io_pio/ice_history_write.F90 b/io_pio/ice_history_write.F90 index 2da823ad..8e233ec7 100644 --- a/io_pio/ice_history_write.F90 +++ b/io_pio/ice_history_write.F90 @@ -9,82 +9,46 @@ ! module ice_history_write - use ice_kinds_mod - use ice_pio, only: ice_pio_subsystem - use netcdf - - implicit none - private - public :: ice_write_hist - save - - logical (kind=log_kind), public :: lcdf64 - -!======================================================================= - - contains - -!======================================================================= -! -! write average ice quantities or snapshots -! -! author: Elizabeth C. Hunke, LANL - - subroutine ice_write_hist (ns) + use ice_pio + use pio + use netcdf, only: NF90_UNLIMITED, NF90_CHUNKED, nf90_global, nf90_noerr - use ice_blocks, only: nx_block, ny_block - use ice_broadcast, only: broadcast_scalar - use ice_calendar, only: time, sec, idate, idate0, write_ic, & - histfreq, histfreq_n, dayyr, days_per_year, use_leap_years - use ice_communicate, only: my_task, master_task - use ice_constants, only: c0, c360, secday, spval, spval_dbl, rad_to_deg - use ice_domain, only: distrb_info, nblocks - use ice_domain_size, only: nx_global, ny_global, max_blocks, max_nstrm - use ice_domain_size, only: block_size_x, block_size_y + use ice_kinds_mod + use ice_constants, only: c0, c360, secday, spval, rad_to_deg + use ice_blocks, only: nx_block, ny_block, block, get_block use ice_exit, only: abort_ice - use ice_fileunits, only: nu_diag + use ice_domain, only: distrb_info, nblocks, blocks_ice + use ice_domain, only: equal_num_blocks_per_cpu + use ice_communicate, only: my_task, master_task, MPI_COMM_ICE + use ice_broadcast, only: broadcast_scalar use ice_gather_scatter, only: gather_global + use ice_domain_size, only: nx_global, ny_global, max_nstrm, max_blocks use ice_grid, only: TLON, TLAT, ULON, ULAT, hm, bm, tarea, uarea, & - dxu, dxt, dyu, dyt, HTN, HTE, ANGLE, ANGLET, tmask, & - lont_bounds, latt_bounds, lonu_bounds, latu_bounds + dxu, dxt, dyu, dyt, HTN, HTE, ANGLE, ANGLET, & + lont_bounds, latt_bounds, lonu_bounds, latu_bounds use ice_history_shared use ice_itd, only: hin_max - use ice_restart_shared, only: runid - - use ice_pio - use pio - use pio_nf + use ice_calendar, only: write_ic, histfreq - integer (kind=int_kind), intent(in) :: ns + implicit none + private + public :: ice_write_hist + save - ! local variables + type coord_attributes ! netcdf coordinate attributes + character (len=11) :: short_name + character (len=45) :: long_name + character (len=20) :: units + end type coord_attributes - integer (kind=int_kind) :: i,j,k,ic,n,nn, & - ncid,status,imtid,jmtid,kmtidi,kmtids,kmtidb, cmtid,timid, & - length,nvertexid,ivertex - integer (kind=int_kind), dimension(2) :: dimid2 - integer (kind=int_kind), dimension(3) :: dimid3 - integer (kind=int_kind), dimension(4) :: dimidz - integer (kind=int_kind), dimension(5) :: dimidcz - integer (kind=int_kind), dimension(3) :: dimid_nverts - integer (kind=int_kind), dimension(4) :: dimidex - real (kind=real_kind) :: ltime - character (char_len) :: title - character (char_len) :: time_period_freq - character (char_len_long) :: ncfile(max_nstrm) - - integer (kind=int_kind) :: iyear, imonth, iday - integer (kind=int_kind) :: icategory,ind,i_aice,boundid - - character (char_len) :: start_time,current_date,current_time - character (len=16) :: c_aice - character (len=8) :: cdate + type req_attributes ! req'd netcdf attributes + type (coord_attributes) :: req + character (len=20) :: coordinates + end type req_attributes - type(file_desc_t) :: File type(io_desc_t) :: iodesc2d, & iodesc3dc, iodesc3dv, iodesc3di, iodesc3db, & iodesc4di, iodesc4ds - type(var_desc_t) :: varid ! 4 coordinate variables: TLON, TLAT, ULON, ULAT INTEGER (kind=int_kind), PARAMETER :: ncoord = 4 @@ -96,59 +60,54 @@ subroutine ice_write_hist (ns) ! lont_bounds, latt_bounds, lonu_bounds, latu_bounds INTEGER (kind=int_kind), PARAMETER :: nvar_verts = 4 - TYPE coord_attributes ! netcdf coordinate attributes - character (len=11) :: short_name - character (len=45) :: long_name - character (len=20) :: units - END TYPE coord_attributes - - TYPE req_attributes ! req'd netcdf attributes - type (coord_attributes) :: req - character (len=20) :: coordinates - END TYPE req_attributes - - TYPE(req_attributes), dimension(nvar) :: var - TYPE(coord_attributes), dimension(ncoord) :: coord_var - TYPE(coord_attributes), dimension(nvar_verts) :: var_nverts - TYPE(coord_attributes), dimension(nvarz) :: var_nz - CHARACTER (char_len), dimension(ncoord) :: coord_bounds + contains - real (kind=real_kind), allocatable :: workr2(:,:,:) - real (kind=real_kind), allocatable :: workr3(:,:,:,:) - real (kind=real_kind), allocatable :: workr4(:,:,:,:,:) - real (kind=real_kind), allocatable :: workr3v(:,:,:,:) +!======================================================================= - character(len=char_len_long) :: & - filename +! +! write average ice quantities or snapshots +! +! author: Elizabeth C. Hunke, LANL - integer (kind=int_kind), dimension(1) :: & - tim_start,tim_length ! dimension quantities for netCDF + subroutine ice_write_hist (ns) - integer (kind=int_kind), dimension(2) :: & - bnd_start,bnd_length ! dimension quantities for netCDF + use ice_calendar, only: time, month, daymo + use ice_fileunits, only: nu_diag - integer (kind=PIO_OFFSET_KIND) :: FRAME_1 = 1 - integer, dimension(2) :: chunksizes + integer (kind=int_kind), intent(in) :: ns !history stream number - integer (kind=int_kind) :: shuffle, deflate, deflate_level - character (len=9) :: ret_str - integer :: ierr + ! local variables - ! We leave shuffle at 0, this is only useful for integer data. - shuffle = 0 + real (kind=real_kind) :: ltime !history timestamp in days + character (char_len_long) :: ncfile(max_nstrm), filename !filenames + character (char_len) :: time_string !model time for logging + logical :: file_exists + integer (kind=int_kind) :: & + i_time, & ! time index + timid ! time var id + type(file_desc_t) :: File - ! If history_deflate_level < 0 then don't do deflation, - ! otherwise it sets the deflate level - if (history_deflate_level < 0) then - deflate = 0 - deflate_level = 0 - else - deflate = 1 - deflate_level = history_deflate_level - endif + type(var_desc_t) :: varid + TYPE(req_attributes), dimension(nvar) :: var + TYPE(coord_attributes), dimension(ncoord) :: coord_var + TYPE(coord_attributes), dimension(nvar_verts) :: var_nverts + TYPE(coord_attributes), dimension(nvarz) :: var_nz if (my_task == master_task) then - call construct_filename(ncfile(ns),'nc',ns) + ! set timestamp in middle of time interval + if (histfreq(ns) == 'm' .or. histfreq(ns) == 'M') then + if (month /= 1) then + ltime=time/int(secday)-real(daymo(month-1))/2.0 + else + ltime=time/int(secday)-real(daymo(12))/2.0 + endif + else if(histfreq(ns) == 'd' .or. histfreq(ns) == 'D') then + ltime=time/int(secday) - 0.5 + else + ltime=time/int(secday) + endif + + call construct_filename(ncfile(ns),'nc',ns,time_string) ! add local directory path name to ncfile if (write_ic) then @@ -156,109 +115,226 @@ subroutine ice_write_hist (ns) else ncfile(ns) = trim(history_dir)//ncfile(ns) endif - filename = ncfile(ns) - end if - call broadcast_scalar(filename, master_task) - ! create file + ! run inquire only on master_task, to avoid possible race condition with + ! multiple procs creating the file + inquire(file=trim(ncfile(ns)),exist=file_exists) + endif + + call broadcast_scalar(ncfile(ns), master_task) + call broadcast_scalar(file_exists, master_task) File%fh=-1 - call ice_pio_initfile(mode='write', filename=trim(filename), File=File, & - clobber=.true., cdf64=lcdf64) + call ice_pio_initfile(mode='write', filename=trim(ncfile(ns)), File=File) + + if (.not. file_exists) then + write(nu_diag,*) 'Writing dimensions and metadata: '//trim(ncfile(ns)) + call ice_hist_create(ns, ncfile(ns), File, var, coord_var, var_nverts, var_nz) + endif + + !----------------------------------------------------------------- + ! extent time coord by 1 and write time variable + !----------------------------------------------------------------- + call ice_pio_check(pio_inq_dimid(File, 'time', timid), & + 'inq dimid time') + call ice_pio_check(pio_inquire_dimension(File, timid, len=i_time), & + 'inquire dim time') + call ice_pio_check(pio_inq_varid(File,'time',varid), & + 'inq varid time') + i_time = i_time + 1 ! index of the current history time + call ice_pio_check(pio_put_var(File,varid,(/i_time/),ltime), & + 'put var time') + + !----------------------------------------------------------------- + ! write time_bounds info + !----------------------------------------------------------------- + + if (hist_avg) then + time_bounds=(/time_beg(ns),time_end(ns)/) + call ice_pio_check(pio_inq_varid(File,'time_bounds',varid), & + 'inq varid time_bounds') + call ice_pio_check(pio_put_var(File,varid,ival=(/time_beg(ns),time_end(ns)/), & + start=(/1,i_time/), count=(/2,1/)), & + 'put var time_bounds') + endif + + + ! these iodesc variables describe how to map from the local variable in blocks, + ! to the output variable call ice_pio_initdecomp(iodesc=iodesc2d) + call ice_pio_initdecomp(ndim3=nverts, iodesc=iodesc3dv, inner_dim=.true.) call ice_pio_initdecomp(ndim3=ncat_hist, iodesc=iodesc3dc) call ice_pio_initdecomp(ndim3=nzilyr, iodesc=iodesc3di) call ice_pio_initdecomp(ndim3=nzlyrb, iodesc=iodesc3db) - call ice_pio_initdecomp(ndim3=nverts, inner_dim=.true., iodesc=iodesc3dv) call ice_pio_initdecomp(ndim3=nzilyr, ndim4=ncat_hist, iodesc=iodesc4di) call ice_pio_initdecomp(ndim3=nzslyr, ndim4=ncat_hist, iodesc=iodesc4ds) -! ltime = time/int(secday) - ltime = real(time/int(secday),kind=real_kind) + if (i_time == 1) then + ! these variables are time-invariant, only write once per file + !----------------------------------------------------------------- + ! write coordinate variables + !----------------------------------------------------------------- + call write_coordinate_variables(File, coord_var, var_nz) + + !----------------------------------------------------------------- + ! write grid masks, area and rotation angle + !----------------------------------------------------------------- + call write_grid_variables(File, var, var_nverts) + + endif !----------------------------------------------------------------- - ! define dimensions + ! write variable data !----------------------------------------------------------------- - if (hist_avg .and. histfreq(ns) /= '1') then - status = pio_def_dim(File,'d2',2,boundid) - endif + call write_2d_variables(ns, File, i_time) + call write_3d_and_4d_variables(ns, File, i_time) + + !----------------------------------------------------------------- + ! close output dataset + !----------------------------------------------------------------- + + call pio_closefile(File) + if (my_task == master_task) then + write(nu_diag,*) 'Wrote ',trim(ncfile(ns)),' at time ',trim(time_string) + endif + + !----------------------------------------------------------------- + ! clean-up PIO descriptors + !----------------------------------------------------------------- + call pio_freedecomp(ice_pio_subsystem, iodesc2d) + call pio_freedecomp(ice_pio_subsystem, iodesc3dv) + call pio_freedecomp(ice_pio_subsystem, iodesc3dc) + call pio_freedecomp(ice_pio_subsystem, iodesc3di) + call pio_freedecomp(ice_pio_subsystem, iodesc3db) + call pio_freedecomp(ice_pio_subsystem, iodesc4di) + call pio_freedecomp(ice_pio_subsystem, iodesc4ds) - status = pio_def_dim(File,'ni',nx_global,imtid) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error defining dim ni') +end subroutine ice_write_hist - status = pio_def_dim(File,'nj',ny_global,jmtid) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error defining dim nj') - status = pio_def_dim(File,'nc',ncat_hist,cmtid) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error defining dim nc') +subroutine ice_hist_create(ns, ncfile, File, var, coord_var, var_nverts, var_nz) - status = pio_def_dim(File,'nkice',nzilyr,kmtidi) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error defining dim nkice') + use ice_calendar, only: idate, idate0, & + dayyr, days_per_year, use_leap_years, histfreq_n - status = pio_def_dim(File,'nksnow',nzslyr,kmtids) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error defining dim nksnow') + integer (kind=int_kind), intent(in) :: ns + character (char_len_long), intent(in) :: ncfile + type(file_desc_t), intent(inout) :: File + + ! local variables - status = pio_def_dim(File,'nkbio',nzblyr,kmtidb) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error defining dim nkbio') + integer (kind=int_kind) :: i, n, status, imtid, jmtid, kmtidi, kmtids, & + kmtidb, cmtid, timid, nvertexid + type(var_desc_t) :: varid + integer (kind=int_kind), dimension(3) :: dimid, dimid_nverts + integer (kind=int_kind), dimension(4) :: dimidz, dimidex + integer (kind=int_kind), dimension(5) :: dimidcz - status = pio_def_dim(File,'time',PIO_UNLIMITED,timid) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error defining dim time') + integer (kind=int_kind) :: deflate, deflate_level ! compression settings + ! We leave shuffle at 0, this is only useful for integer data. + integer (kind=int_kind), parameter :: shuffle = 0 - status = pio_def_dim(File,'nvertices',nverts,nvertexid) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error defining dim nvertices') + integer (kind=int_kind) :: ind,boundid + + character (char_len) :: title, start_time,current_date,current_time,time_period_freq + character (len=8) :: cdate + character (char_len), dimension(ncoord) :: coord_bounds + + type(req_attributes), dimension(nvar), intent(inout) :: var + type(coord_attributes), dimension(ncoord), intent(inout) :: coord_var + type(coord_attributes), dimension(nvar_verts), intent(inout) :: var_nverts + type(coord_attributes), dimension(nvarz), intent(inout) :: var_nz + + ! If history_deflate_level < 0 then don't do deflation, + ! otherwise it sets the deflate level + if (history_deflate_level < 0) then + deflate = 0 + deflate_level = 0 + else + deflate = 1 + deflate_level = history_deflate_level + endif !----------------------------------------------------------------- - ! define coordinate variables: time, time_bounds + ! define dimensions !----------------------------------------------------------------- + call ice_pio_check(pio_def_dim(File, 'time', NF90_UNLIMITED, timid), & + 'def dim time') - status = pio_def_var(File,'time',pio_real,(/timid/),varid) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error def var time') + if (hist_avg .and. histfreq(ns) /= '1') then + call ice_pio_check(pio_def_dim(File,'d2',2,boundid), 'def dim d2') + endif - status = pio_put_att(File,varid,'long_name','model time') + call ice_pio_check(pio_def_dim(File, 'ni', nx_global, imtid), & + 'def dim ni') + call ice_pio_check(pio_def_dim(File, 'nj', ny_global, jmtid), & + 'def dim nj') + call ice_pio_check(pio_def_dim(File, 'nc', ncat_hist, cmtid), & + 'def dim nc') + call ice_pio_check(pio_def_dim(File, 'nkice', nzilyr, kmtidi), & + 'def dim nkice') + call ice_pio_check(pio_def_dim(File, 'nksnow', nzslyr, kmtids), & + 'def dim nksnow') + call ice_pio_check(pio_def_dim(File, 'nkbio', nzblyr, kmtidb), & + 'def dim nkbio') + call ice_pio_check(pio_def_dim(File, 'nvertices', nverts, nvertexid), & + 'def dim nverts') - write(cdate,'(i8.8)') idate0 - write(title,'(a,a,a,a,a,a,a)') 'days since ', & - cdate(1:4),'-',cdate(5:6),'-',cdate(7:8),' 00:00:00' - status = pio_put_att(File,varid,'units',trim(title)) + !----------------------------------------------------------------- + ! define coordinate variables + !----------------------------------------------------------------- - if (days_per_year == 360) then - status = pio_put_att(File,varid,'calendar','360_day') - elseif (days_per_year == 365 .and. .not.use_leap_years ) then - status = pio_put_att(File,varid,'calendar','NoLeap') - elseif (use_leap_years) then - status = pio_put_att(File,varid,'calendar','Gregorian') - else - call abort_ice( 'ice Error: invalid calendar settings') - endif + call ice_pio_check(pio_def_var(File,'time',pio_real,(/timid/),varid), & + 'def var time') + call ice_pio_check(pio_put_att(File,varid,'long_name','model time'), & + 'put_att long_name') + + write(cdate,'(i8.8)') idate0 + write(title,'(a,a,a,a,a,a,a,a)') 'days since ', & + cdate(1:4),'-',cdate(5:6),'-',cdate(7:8),' 00:00:00' + call ice_pio_check(pio_put_att(File,varid,'units',title), & + 'put_att time units') + + if (days_per_year == 360) then + call ice_pio_check(pio_put_att(File,varid,'calendar','360_day'), & + 'att time calendar') + elseif (days_per_year == 365 .and. .not.use_leap_years ) then + call ice_pio_check(pio_put_att(File,varid,'calendar','NoLeap'), & + 'att time calendar') + elseif (use_leap_years) then + call ice_pio_check(pio_put_att(File,varid,'calendar','proleptic_gregorian'), & + 'att time calendar') + else + call abort_ice( 'ice Error: invalid calendar settings') + endif - if (hist_avg .and. histfreq(ns) /= '1') then - status = pio_put_att(File,varid,'bounds','time_bounds') - endif + if (hist_avg .and. histfreq(ns) /= '1' ) then + call ice_pio_check(pio_put_att(File,varid,'bounds','time_bounds'), & + 'att time bounds') + endif - ! Define attributes for time_bounds if hist_avg is true - if (hist_avg .and. histfreq(ns) /= '1') then - dimid2(1) = boundid - dimid2(2) = timid - status = pio_def_var(File,'time_bounds',pio_real,dimid2,varid) - - status = pio_put_att(File,varid,'long_name', & - 'boundaries for time-averaging interval') - write(cdate,'(i8.8)') idate0 - write(title,'(a,a,a,a,a,a,a,a)') 'days since ', & - cdate(1:4),'-',cdate(5:6),'-',cdate(7:8),' 00:00:00' - status = pio_put_att(File,varid,'units',trim(title)) - endif + !----------------------------------------------------------------- + ! Define attributes for time bounds if hist_avg is true + !----------------------------------------------------------------- + + if (hist_avg .and. histfreq(ns) /= '1') then + dimid(1) = boundid + dimid(2) = timid + call ice_pio_check(pio_def_var(File, 'time_bounds', & + pio_real,dimid(1:2),varid), & + 'def var time_bounds') + + call ice_pio_check(pio_put_att(File,varid,'long_name', & + 'boundaries for time-averaging interval'), & + 'att time_bounds long_name') + write(cdate,'(i8.8)') idate0 + write(title,'(a,a,a,a,a,a,a,a)') 'days since ', & + cdate(1:4),'-',cdate(5:6),'-',cdate(7:8),' 00:00:00' + call ice_pio_check(pio_put_att(File,varid,'units',title), & + 'att time_bounds units') + endif !----------------------------------------------------------------- ! define information for required time-invariant variables @@ -267,19 +343,19 @@ subroutine ice_write_hist (ns) ind = 0 ind = ind + 1 coord_var(ind) = coord_attributes('TLON', & - 'T grid center longitude', 'degrees_east') + 'T grid center longitude', 'degrees_east') coord_bounds(ind) = 'lont_bounds' ind = ind + 1 coord_var(ind) = coord_attributes('TLAT', & - 'T grid center latitude', 'degrees_north') + 'T grid center latitude', 'degrees_north') coord_bounds(ind) = 'latt_bounds' ind = ind + 1 coord_var(ind) = coord_attributes('ULON', & - 'U grid center longitude', 'degrees_east') + 'U grid center longitude', 'degrees_east') coord_bounds(ind) = 'lonu_bounds' ind = ind + 1 coord_var(ind) = coord_attributes('ULAT', & - 'U grid center latitude', 'degrees_north') + 'U grid center latitude', 'degrees_north') coord_bounds(ind) = 'latu_bounds' var_nz(1) = coord_attributes('NCAT', 'category maximum thickness', 'm') @@ -348,550 +424,679 @@ subroutine ice_write_hist (ns) ! define attributes for time-invariant variables !----------------------------------------------------------------- - dimid2(1) = imtid - dimid2(2) = jmtid - - do i = 1, ncoord - status = pio_def_var(File, trim(coord_var(i)%short_name), pio_real, & - dimid2, varid) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error def var '//trim(coord_var(i)%short_name)) - status = pio_def_var_deflate(File, varid, shuffle, deflate, & - deflate_level) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error deflating var '//trim(coord_var(i)%short_name)) - - if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then - status = pio_def_var_chunking(File, varid, NF90_CHUNKED, & - (/ history_chunksize_x, history_chunksize_y /)) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error chunking var '//trim(coord_var(i)%short_name)) - endif - - status = pio_put_att(File,varid,'long_name',trim(coord_var(i)%long_name)) - status = pio_put_att(File, varid, 'units', trim(coord_var(i)%units)) - status = pio_put_att(File, varid, 'missing_value', spval) - status = pio_put_att(File, varid,'_FillValue',spval) - if (coord_var(i)%short_name == 'ULAT') then - status = pio_put_att(File,varid,'comment', & - trim('Latitude of NE corner of T grid cell')) - endif - if (f_bounds) then - status = pio_put_att(File, varid, 'bounds', trim(coord_bounds(i))) - endif - enddo + dimid(1) = imtid + dimid(2) = jmtid + dimid(3) = timid - ! Extra dimensions (NCAT, NZILYR, NZSLYR, NZBLYR) - dimidex(1)=cmtid - dimidex(2)=kmtidi - dimidex(3)=kmtids - dimidex(4)=kmtidb + do i = 1, ncoord + call ice_pio_check(pio_def_var(File, coord_var(i)%short_name, pio_real, & + dimid(1:2), varid), & + 'def var '//coord_var(i)%short_name) - do i = 1, nvarz - if (igrdz(i)) then - status = pio_def_var(File, trim(var_nz(i)%short_name), pio_real, & - (/dimidex(i)/), varid) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error def var '//trim(var_nz(i)%short_name)) - - status = pio_put_att(File, varid, 'long_name', var_nz(i)%long_name) - status = pio_put_att(File, varid, 'units' , var_nz(i)%units) - endif - enddo + if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then + call ice_pio_check(pio_def_var_chunking(File, varid, NF90_CHUNKED, & + (/ history_chunksize_x, history_chunksize_y /)), & + 'def var chunking '//coord_var(i)%short_name) + endif - ! Attributes for tmask defined separately, since it has no units - if (igrd(n_tmask)) then - status = pio_def_var(File, 'tmask', pio_real, dimid2, varid) - if (status /= nf90_noerr) call abort_ice('ice: Error defining var tmask') - status = pio_def_var_deflate(File, varid, shuffle, deflate, & - deflate_level) - if (status /= nf90_noerr) call abort_ice('ice: Error deflating var tmask') - - if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then - status = pio_def_var_chunking(File, varid, NF90_CHUNKED, & - (/ history_chunksize_x, history_chunksize_y /)) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error def var chunking tmask') - endif + call ice_pio_check(pio_def_var_deflate(File, varid, shuffle, deflate, deflate_level), & + 'deflate '//coord_var(i)%short_name) - status = pio_put_att(File,varid, 'long_name', 'ocean grid mask') - if (status /= nf90_noerr) call abort_ice('ice: Error put_att') - status = pio_put_att(File, varid, 'coordinates', 'TLON TLAT') - if (status /= nf90_noerr) call abort_ice('ice: Error put_att') - status = pio_put_att(File, varid, 'missing_value', spval) - if (status /= nf90_noerr) call abort_ice('ice: Error put_att') - status = pio_put_att(File, varid,'_FillValue',spval) - if (status /= nf90_noerr) call abort_ice('ice: Error put_att') - status = pio_put_att(File,varid,'comment', '0 = land, 1 = ocean') - if (status /= nf90_noerr) call abort_ice('ice: Error put_att') + call ice_pio_check(pio_put_att(File, varid,'long_name',coord_var(i)%long_name), & + 'put att long_name '//coord_var(i)%short_name) + call ice_pio_check(pio_put_att(File, varid, 'units', coord_var(i)%units), & + 'put att units '//coord_var(i)%short_name) + call ice_pio_check(pio_put_att(File, varid,'missing_value',spval), & + 'put att missing_value '//coord_var(i)%short_name) - endif + call ice_pio_check(pio_put_att(File, varid, '_FillValue', spval), & + 'put att _FillValue '//coord_var(i)%short_name) - if (igrd(n_blkmask)) then - status = pio_def_var(File, 'blkmask', pio_real, dimid2, varid) - status = pio_def_var_deflate(File, varid, shuffle, deflate, & - deflate_level) - if (status /= nf90_noerr) call abort_ice('ice: Error deflating var blkmask') - if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then - status = pio_def_var_chunking(File, varid, NF90_CHUNKED, & - (/ history_chunksize_x, history_chunksize_y /)) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error def var chunking blkmask') - endif + if (coord_var(i)%short_name == 'ULAT') then + call ice_pio_check(pio_put_att(File,varid,'comment', & + 'Latitude of NE corner of T grid cell'), & + 'put att comment for '//coord_var(i)%short_name) + endif + if (f_bounds) then + call ice_pio_check(pio_put_att(File, varid, 'bounds', coord_bounds(i)), & + 'put att bounds '//coord_var(i)%short_name) + endif + enddo - status = pio_put_att(File,varid, 'long_name', 'ice grid block mask') - status = pio_put_att(File, varid, 'coordinates', 'TLON TLAT') - status = pio_put_att(File,varid,'comment', 'mytask + iblk/100') - status = pio_put_att(File, varid, 'missing_value', spval) - status = pio_put_att(File, varid,'_FillValue',spval) - endif + ! Extra dimensions (NCAT, NZILYR, NZSLYR, NZBLYR) + dimidex(1)=cmtid + dimidex(2)=kmtidi + dimidex(3)=kmtids + dimidex(4)=kmtidb + + do i = 1, nvarz + if (igrdz(i)) then + call ice_pio_check(pio_def_var(File, var_nz(i)%short_name, & + pio_real, (/dimidex(i)/), varid), & + 'def var '//var_nz(i)%short_name) + + call ice_pio_check(pio_def_var_deflate(File, varid, shuffle, deflate, & + deflate_level), & + 'deflate '//var_nz(i)%short_name) + + call ice_pio_check(pio_put_att(File, varid,'long_name',var_nz(i)%long_name), & + 'put att long_name '//var_nz(i)%short_name) + call ice_pio_check(pio_put_att(File, varid, 'units', var_nz(i)%units), & + 'for att units '//var_nz(i)%short_name) + endif + enddo - do i = 3, nvar ! note: n_tmask=1, n_blkmask=2 - if (igrd(i)) then - status = pio_def_var(File, trim(var(i)%req%short_name), & - pio_real, dimid2, varid) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error def var '//trim(var(i)%req%short_name)) - status = pio_def_var_deflate(File, varid, shuffle, deflate, & - deflate_level) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error deflating var '//trim(var(i)%req%short_name)) - - if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then - status = pio_def_var_chunking(File, varid, NF90_CHUNKED, & - (/ history_chunksize_x, history_chunksize_y /)) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error chunking var '//trim(var(i)%req%short_name)) - endif - - status = pio_put_att(File,varid, 'long_name', trim(var(i)%req%long_name)) - status = pio_put_att(File, varid, 'units', trim(var(i)%req%units)) - status = pio_put_att(File, varid, 'coordinates', trim(var(i)%coordinates)) - status = pio_put_att(File, varid, 'missing_value', spval) - status = pio_put_att(File, varid,'_FillValue',spval) - endif - enddo + ! Attributes for tmask, blkmask defined separately, since they have no units + if (igrd(n_tmask)) then + call ice_pio_check(pio_def_var(File, 'tmask', pio_real, dimid(1:2), varid), & + 'def var tmask') - ! Fields with dimensions (nverts,nx,ny) - dimid_nverts(1) = nvertexid - dimid_nverts(2) = imtid - dimid_nverts(3) = jmtid - do i = 1, nvar_verts - if (f_bounds) then - status = pio_def_var(File, trim(var_nverts(i)%short_name), & - pio_real,dimid_nverts, varid) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error def var '//trim(var_nverts(i)%short_name)) - status = pio_def_var_deflate(File, varid, shuffle, deflate, & - deflate_level) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error deflating var '//trim(var_nverts(i)%short_name)) - if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then - status = pio_def_var_chunking(File, varid, NF90_CHUNKED, & - (/ 1, history_chunksize_x, history_chunksize_y /)) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error chunking var '//trim(var_nverts(i)%short_name)) - endif - - status = & - pio_put_att(File,varid, 'long_name', trim(var_nverts(i)%long_name)) - status = & - pio_put_att(File, varid, 'units', trim(var_nverts(i)%units)) - status = pio_put_att(File, varid, 'missing_value', spval) - status = pio_put_att(File, varid,'_FillValue',spval) - endif - enddo - - !----------------------------------------------------------------- - ! define attributes for time-variant variables - !----------------------------------------------------------------- + if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then + call ice_pio_check(pio_def_var_chunking(File, varid, NF90_CHUNKED, & + (/ history_chunksize_x, history_chunksize_y /)), & + 'def var chunking tmask') + endif - !----------------------------------------------------------------- - ! 2D - !----------------------------------------------------------------- + call ice_pio_check(pio_def_var_deflate(File, varid, shuffle, deflate, & + deflate_level), 'deflating var tmask') + + call ice_pio_check(pio_put_att(File,varid, 'long_name', 'ocean grid mask'), & + 'put att tmask long_name') + call ice_pio_check(pio_put_att(File, varid, 'coordinates', 'TLON TLAT'), & + 'put att tmask units') + call ice_pio_check(pio_put_att(File,varid,'comment', '0 = land, 1 = ocean'), & + 'put att tmask comment') + call ice_pio_check(pio_put_att(File,varid,'missing_value',spval), & + 'put att missing_value for tmask') + call ice_pio_check(pio_put_att(File,varid,'_FillValue',spval), & + 'put att _FillValue for tmask') + endif + + if (igrd(n_blkmask)) then + call ice_pio_check(pio_def_var(File, 'blkmask', pio_real, dimid(1:2), varid), & + 'def var blkmask') - dimid3(1) = imtid - dimid3(2) = jmtid - dimid3(3) = timid - - do n=1,num_avail_hist_fields_2D - if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then - status = pio_def_var(File, trim(avail_hist_fields(n)%vname), & - pio_real, dimid3, varid) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error def var '//trim(avail_hist_fields(n)%vname)) - status = pio_def_var_deflate(File, varid, shuffle, deflate, & - deflate_level) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error deflating var '//trim(avail_hist_fields(n)%vname)) if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then - status = pio_def_var_chunking(File, varid, NF90_CHUNKED, & - (/ history_chunksize_x, history_chunksize_y, 1 /)) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error chunking var '//trim(avail_hist_fields(n)%vname)) + call ice_pio_check(pio_def_var_chunking(File, varid, NF90_CHUNKED, & + (/ history_chunksize_x, history_chunksize_y /)), & + 'def var chunking blkmask') endif - status = pio_put_att(File,varid,'units', & - trim(avail_hist_fields(n)%vunit)) - status = pio_put_att(File,varid, 'long_name', & - trim(avail_hist_fields(n)%vdesc)) - status = pio_put_att(File,varid,'coordinates', & - trim(avail_hist_fields(n)%vcoord)) - status = pio_put_att(File,varid,'cell_measures', & - trim(avail_hist_fields(n)%vcellmeas)) - status = pio_put_att(File,varid,'missing_value',spval) - status = pio_put_att(File,varid,'_FillValue',spval) - - ! Add cell_methods attribute to variables if averaged - if (hist_avg .and. histfreq(ns) /= '1') then - if (TRIM(avail_hist_fields(n)%vname)/='sig1' & - .or.TRIM(avail_hist_fields(n)%vname)/='sig2') then - status = pio_put_att(File,varid,'cell_methods','time: mean') - endif + call ice_pio_check(pio_def_var_deflate(File, varid, shuffle, deflate, & + deflate_level), & + 'deflating var blkmask') + + call ice_pio_check(pio_put_att(File,varid, 'long_name', 'ice grid block mask'), & + 'put att blkmask long_name') + call ice_pio_check(pio_put_att(File, varid, 'coordinates', 'TLON TLAT'), & + 'put att blkmask coordinates') + call ice_pio_check(pio_put_att(File,varid,'comment', 'mytask + iblk/100'), & + 'put att blkmask comment') + call ice_pio_check(pio_put_att(File,varid,'missing_value',spval), & + 'put att blkmask missing_value') + call ice_pio_check(pio_put_att(File,varid,'_FillValue',spval), & + 'put att blkmask _FillValue') + endif + + do i = 3, nvar ! note n_tmask=1, n_blkmask=2 + if (igrd(i)) then + call ice_pio_check(pio_def_var(File, var(i)%req%short_name, & + pio_real, dimid(1:2), varid), & + 'def variable '//var(i)%req%short_name) + + if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then + call ice_pio_check(pio_def_var_chunking(File, varid, NF90_CHUNKED, & + (/ history_chunksize_x, history_chunksize_y /)), & + 'def var chunking '//var(i)%req%short_name) + endif + + call ice_pio_check(pio_def_var_deflate(File, varid, shuffle, deflate, & + deflate_level), & + 'deflate var '//var(i)%req%short_name) + + call ice_pio_check(pio_put_att(File,varid, 'long_name', var(i)%req%long_name), & + 'put att long_name '//var(i)%req%short_name) + call ice_pio_check(pio_put_att(File, varid, 'units', var(i)%req%units), & + 'put att units '//var(i)%req%short_name) + call ice_pio_check(pio_put_att(File, varid, 'coordinates', var(i)%coordinates), & + 'put att coordinates '//var(i)%req%short_name) + call ice_pio_check(pio_put_att(File,varid,'missing_value',spval), & + 'put att missing_value '//var(i)%req%short_name) + call ice_pio_check(pio_put_att(File,varid,'_FillValue',spval), & + 'put att _FillValue '//var(i)%req%short_name) endif + enddo - if (histfreq(ns) == '1' .or. .not. hist_avg & - .or. n==n_divu(ns) .or. n==n_shear(ns) & ! snapshots - .or. n==n_sig1(ns) .or. n==n_sig2(ns) & - .or. n==n_trsig(ns) & - .or. n==n_mlt_onset(ns) .or. n==n_frz_onset(ns) & - .or. n==n_hisnap(ns) .or. n==n_aisnap(ns)) then - status = pio_put_att(File,varid,'time_rep','instantaneous') - else - status = pio_put_att(File,varid,'time_rep','averaged') + ! Fields with dimensions (nverts,nx,ny) + dimid_nverts(1) = nvertexid + dimid_nverts(2) = imtid + dimid_nverts(3) = jmtid + do i = 1, nvar_verts + if (f_bounds) then + call ice_pio_check(pio_def_var(File, var_nverts(i)%short_name, & + pio_real,dimid_nverts, varid), & + 'def var '//var_nverts(i)%short_name) + + if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then + call ice_pio_check(pio_def_var_chunking(File, varid, NF90_CHUNKED, & + (/ 1, history_chunksize_x, history_chunksize_y /)), & + 'def var chunking '//var_nverts(i)%short_name) + endif + + call ice_pio_check(pio_def_var_deflate(File, varid, shuffle, deflate, & + deflate_level), & + 'deflate var '//var_nverts(i)%short_name) + + call ice_pio_check(pio_put_att(File,varid, 'long_name', & + var_nverts(i)%long_name), & + 'put att long_name '//var_nverts(i)%short_name) + call ice_pio_check(pio_put_att(File, varid, 'units', var_nverts(i)%units), & + 'put att units '//var_nverts(i)%short_name) + call ice_pio_check(pio_put_att(File,varid,'missing_value',spval), & + 'put att missing_value '//var_nverts(i)%short_name) + call ice_pio_check(pio_put_att(File,varid,'_FillValue',spval), & + 'put att _FillValue '//var_nverts(i)%short_name) endif - endif - enddo ! num_avail_hist_fields_2D + enddo - !----------------------------------------------------------------- - ! 3D (category) - !----------------------------------------------------------------- + do n=1,num_avail_hist_fields_2D + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + call ice_pio_check(pio_def_var(File, avail_hist_fields(n)%vname, & + pio_real, dimid, varid), & + 'def var '//avail_hist_fields(n)%vname) + + if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then + call ice_pio_check(pio_def_var_chunking(File, varid, NF90_CHUNKED, & + (/ history_chunksize_x, history_chunksize_y, 1 /)), & + 'def var chunking '//avail_hist_fields(n)%vname) + endif + + call ice_pio_check(pio_def_var_deflate(File, varid, shuffle, deflate, & + deflate_level), & + 'deflate '//avail_hist_fields(n)%vname) + + call ice_pio_check(pio_put_att(File,varid,'units', & + avail_hist_fields(n)%vunit), & + 'put att units '//avail_hist_fields(n)%vname) + call ice_pio_check(pio_put_att(File,varid, 'long_name', & + avail_hist_fields(n)%vdesc), & + 'put att long_name '//avail_hist_fields(n)%vname) + call ice_pio_check(pio_put_att(File,varid,'coordinates', & + avail_hist_fields(n)%vcoord), & + 'put att coordinates '//avail_hist_fields(n)%vname) + call ice_pio_check(pio_put_att(File,varid,'cell_measures', & + avail_hist_fields(n)%vcellmeas), & + 'put att cell_measures '//avail_hist_fields(n)%vname) + if (avail_hist_fields(n)%vcomment /= "none") then + call ice_pio_check(pio_put_att(File,varid,'comment', & + avail_hist_fields(n)%vcomment), & + 'put att comment '//avail_hist_fields(n)%vname) + endif + call ice_pio_check(pio_put_att(File,varid,'missing_value',spval), & + 'put att missing_value '//avail_hist_fields(n)%vname) + call ice_pio_check(pio_put_att(File,varid,'_FillValue',spval), & + 'put att _FillValue '//avail_hist_fields(n)%vname) + + !--------------------------------------------------------------- + ! Add cell_methods attribute to variables if averaged + !--------------------------------------------------------------- + if (hist_avg) then + if (TRIM(avail_hist_fields(n)%vname)/='sig1' .or. & + TRIM(avail_hist_fields(n)%vname)/='sig2') then + call ice_pio_check(pio_put_att(File,varid,'cell_methods','time: mean'), & + 'put att cell methods time mean '//avail_hist_fields(n)%vname) + endif + endif + + if (histfreq(ns) == '1' .or. .not. hist_avg & + .or. n==n_divu(ns) .or. n==n_shear(ns) & ! snapshots + .or. n==n_sig1(ns) .or. n==n_sig2(ns) & + .or. n==n_trsig(ns) & + .or. n==n_mlt_onset(ns) .or. n==n_frz_onset(ns) & + .or. n==n_hisnap(ns) .or. n==n_aisnap(ns)) then + call ice_pio_check(pio_put_att(File,varid,'time_rep','instantaneous'), & + 'put att time_rep instantaneous') + else + call ice_pio_check(pio_put_att(File,varid,'time_rep','averaged'), & + 'put att time_rep averaged') + endif + endif + enddo ! num_avail_hist_fields_2D - dimidz(1) = imtid - dimidz(2) = jmtid - dimidz(3) = cmtid - dimidz(4) = timid - - do n = n2D + 1, n3Dccum - if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then - status = pio_def_var(File, trim(avail_hist_fields(n)%vname), & - pio_real, dimidz, varid) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error def var '//trim(avail_hist_fields(n)%vname)) - status = pio_def_var_deflate(File, varid, shuffle, deflate, & - deflate_level) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error deflating var '//trim(avail_hist_fields(n)%vname)) + ! 3D (category) - if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then - status = pio_def_var_chunking(File, varid, NF90_CHUNKED, & - (/ history_chunksize_x, history_chunksize_y, 1, 1 /)) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error chunking var '//trim(avail_hist_fields(n)%vname)) - endif + dimidz(1) = imtid + dimidz(2) = jmtid + dimidz(3) = cmtid + dimidz(4) = timid - status = pio_put_att(File,varid,'units', & - trim(avail_hist_fields(n)%vunit)) - status = pio_put_att(File,varid, 'long_name', & - trim(avail_hist_fields(n)%vdesc)) - status = pio_put_att(File,varid,'coordinates', & - trim(avail_hist_fields(n)%vcoord)) - status = pio_put_att(File,varid,'cell_measures', & - trim(avail_hist_fields(n)%vcellmeas)) - status = pio_put_att(File,varid,'missing_value',spval) - status = pio_put_att(File,varid,'_FillValue',spval) - - ! Add cell_methods attribute to variables if averaged - if (hist_avg .and. histfreq(ns) /= '1') then - status = pio_put_att(File,varid,'cell_methods','time: mean') - endif + do n = n2D + 1, n3Dccum + ! to-do: use check subroutine + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_def_var(File, avail_hist_fields(n)%vname, & + pio_real, dimidz, varid) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining variable '//avail_hist_fields(n)%vname) - if (histfreq(ns) == '1' .or. .not. hist_avg) then - status = pio_put_att(File,varid,'time_rep','instantaneous') - else - status = pio_put_att(File,varid,'time_rep','averaged') - endif - endif - enddo ! num_avail_hist_fields_3Dc + if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then + call ice_pio_check(pio_def_var_chunking(File, varid, NF90_CHUNKED, & + (/ history_chunksize_x, history_chunksize_y, 1, 1 /)), & + 'def var chunking '//avail_hist_fields(n)%vname) + endif - !----------------------------------------------------------------- - ! 3D (ice layers) - !----------------------------------------------------------------- + status = pio_def_var_deflate(File, varid, shuffle, deflate, & + deflate_level) + if (status /= nf90_noerr) call abort_ice( & + 'Error deflating variable '//avail_hist_fields(n)%vname) - dimidz(1) = imtid - dimidz(2) = jmtid - dimidz(3) = kmtidi - dimidz(4) = timid - - do n = n3Dccum + 1, n3Dzcum - if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then - status = pio_def_var(File, trim(avail_hist_fields(n)%vname), & - pio_real, dimidz, varid) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error def var '//trim(avail_hist_fields(n)%vname)) - status = pio_def_var_deflate(File, varid, shuffle, deflate, & - deflate_level) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error deflating var '//trim(avail_hist_fields(n)%vname)) - if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then - status = pio_def_var_chunking(File, varid, NF90_CHUNKED, & - (/ history_chunksize_x, history_chunksize_y, 1, 1 /)) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error chunking var '//trim(avail_hist_fields(n)%vname)) + status = pio_put_att(File,varid,'units', & + avail_hist_fields(n)%vunit) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining units for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid, 'long_name', & + avail_hist_fields(n)%vdesc) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining long_name for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid,'coordinates', & + avail_hist_fields(n)%vcoord) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining coordinates for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid,'cell_measures', & + avail_hist_fields(n)%vcellmeas) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining cell measures for '//avail_hist_fields(n)%vname) + if (avail_hist_fields(n)%vcomment /= "none") then + call ice_pio_check(pio_put_att(File,varid,'comment', & + avail_hist_fields(n)%vcomment), & + 'put att comment '//avail_hist_fields(n)%vname) + endif + status = pio_put_att(File,varid,'missing_value',spval) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining missing_value for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid,'_FillValue',spval) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining _FillValue for '//avail_hist_fields(n)%vname) + + !--------------------------------------------------------------- + ! Add cell_methods attribute to variables if averaged + !--------------------------------------------------------------- + if (hist_avg .and. histfreq(ns) /= '1') then + status = pio_put_att(File,varid,'cell_methods','time: mean') + if (status /= nf90_noerr) call abort_ice( & + 'Error defining cell methods for '//avail_hist_fields(n)%vname) + endif + + if (histfreq(ns) == '1' .or. .not. hist_avg) then + status = pio_put_att(File,varid,'time_rep','instantaneous') + else + status = pio_put_att(File,varid,'time_rep','averaged') + endif endif + enddo ! num_avail_hist_fields_3Dc - status = pio_put_att(File,varid,'units', & - trim(avail_hist_fields(n)%vunit)) - status = pio_put_att(File,varid, 'long_name', & - trim(avail_hist_fields(n)%vdesc)) - status = pio_put_att(File,varid,'coordinates', & - trim(avail_hist_fields(n)%vcoord)) - status = pio_put_att(File,varid,'cell_measures', & - trim(avail_hist_fields(n)%vcellmeas)) - status = pio_put_att(File,varid,'missing_value',spval) - status = pio_put_att(File,varid,'_FillValue',spval) - - ! Add cell_methods attribute to variables if averaged - if (hist_avg .and. histfreq(ns) /= '1') then - status = pio_put_att(File,varid,'cell_methods','time: mean') - endif + ! 3D (ice layers) - if (histfreq(ns) == '1' .or. .not. hist_avg) then - status = pio_put_att(File,varid,'time_rep','instantaneous') - else - status = pio_put_att(File,varid,'time_rep','averaged') - endif - endif - enddo ! num_avail_hist_fields_3Dz - - !----------------------------------------------------------------- - ! 3D (biology layers) - !----------------------------------------------------------------- + dimidz(1) = imtid + dimidz(2) = jmtid + dimidz(3) = kmtidi + dimidz(4) = timid - dimidz(1) = imtid - dimidz(2) = jmtid - dimidz(3) = kmtidb - dimidz(4) = timid - - do n = n3Dzcum + 1, n3Dbcum - if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then - status = pio_def_var(File, trim(avail_hist_fields(n)%vname), & - pio_real, dimidz, varid) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error def var '//trim(avail_hist_fields(n)%vname)) - status = pio_def_var_deflate(File, varid, shuffle, deflate, & - deflate_level) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error deflating var '//trim(avail_hist_fields(n)%vname)) + do n = n3Dccum + 1, n3Dzcum + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_def_var(File, avail_hist_fields(n)%vname, & + pio_real, dimidz, varid) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining variable '//avail_hist_fields(n)%vname) - if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then - status = pio_def_var_chunking(File, varid, NF90_CHUNKED, & - (/ history_chunksize_x, history_chunksize_y, 1, 1 /)) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error deflating var '//trim(avail_hist_fields(n)%vname)) - endif + if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then + call ice_pio_check(pio_def_var_chunking(File, varid, NF90_CHUNKED, & + (/ history_chunksize_x, history_chunksize_y, 1, 1 /)), & + 'def var chunking '//avail_hist_fields(n)%vname) + endif - status = pio_put_att(File,varid,'units', & - trim(avail_hist_fields(n)%vunit)) - status = pio_put_att(File,varid, 'long_name', & - trim(avail_hist_fields(n)%vdesc)) - status = pio_put_att(File,varid,'coordinates', & - trim(avail_hist_fields(n)%vcoord)) - status = pio_put_att(File,varid,'cell_measures', & - trim(avail_hist_fields(n)%vcellmeas)) - status = pio_put_att(File,varid,'missing_value',spval) - status = pio_put_att(File,varid,'_FillValue',spval) - - ! Add cell_methods attribute to variables if averaged - if (hist_avg .and. histfreq(ns) /= '1') then - status = pio_put_att(File,varid,'cell_methods','time: mean') - endif + status = pio_def_var_deflate(File, varid, shuffle, deflate, & + deflate_level) + if (status /= nf90_noerr) call abort_ice( & + 'Error deflating variable '//avail_hist_fields(n)%vname) + + status = pio_put_att(File,varid,'units', & + avail_hist_fields(n)%vunit) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining units for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid, 'long_name', & + avail_hist_fields(n)%vdesc) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining long_name for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid,'coordinates', & + avail_hist_fields(n)%vcoord) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining coordinates for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid,'cell_measures', & + avail_hist_fields(n)%vcellmeas) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining cell measures for '//avail_hist_fields(n)%vname) + if (avail_hist_fields(n)%vcomment /= "none") then + call ice_pio_check(pio_put_att(File,varid,'comment', & + avail_hist_fields(n)%vcomment), & + 'put att comment '//avail_hist_fields(n)%vname) + endif + status = pio_put_att(File,varid,'missing_value',spval) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining missing_value for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid,'_FillValue',spval) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining _FillValue for '//avail_hist_fields(n)%vname) + + ! Add cell_methods attribute to variables if averaged + if (hist_avg .and. histfreq(ns) /= '1') then + status = pio_put_att(File,varid,'cell_methods','time: mean') + endif + + if (histfreq(ns) == '1' .or. .not. hist_avg) then + status = pio_put_att(File,varid,'time_rep','instantaneous') + else + status = pio_put_att(File,varid,'time_rep','averaged') + endif - if (histfreq(ns) == '1' .or. .not. hist_avg) then - status = pio_put_att(File,varid,'time_rep','instantaneous') - else - status = pio_put_att(File,varid,'time_rep','averaged') endif - endif - enddo ! num_avail_hist_fields_3Db + enddo ! num_avail_hist_fields_3Dz - !----------------------------------------------------------------- - ! define attributes for 4D variables - ! time coordinate is dropped - !----------------------------------------------------------------- + ! 3D (biology layers) - !----------------------------------------------------------------- - ! 4D (ice categories and layers) - !----------------------------------------------------------------- + dimidz(1) = imtid + dimidz(2) = jmtid + dimidz(3) = kmtidb + dimidz(4) = timid - dimidcz(1) = imtid - dimidcz(2) = jmtid - dimidcz(3) = kmtidi - dimidcz(4) = cmtid - dimidcz(5) = timid - - do n = n3Dbcum + 1, n4Dicum - if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then - status = pio_def_var(File, trim(avail_hist_fields(n)%vname), & - pio_real, dimidcz, varid) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error def var '//trim(avail_hist_fields(n)%vname)) - status = pio_def_var_deflate(File, varid, shuffle, deflate, & - deflate_level) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error deflating var '//trim(avail_hist_fields(n)%vname)) + do n = n3Dzcum + 1, n3Dbcum + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_def_var(File, avail_hist_fields(n)%vname, & + pio_real, dimidz, varid) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining variable '//avail_hist_fields(n)%vname) - if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then - status = pio_def_var_chunking(File, varid, NF90_CHUNKED, & - (/ history_chunksize_x, history_chunksize_y, 1, 1, 1 /)) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error chunking var '//trim(avail_hist_fields(n)%vname)) - endif + if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then + call ice_pio_check(pio_def_var_chunking(File, varid, NF90_CHUNKED, & + (/ history_chunksize_x, history_chunksize_y, 1, 1 /)), & + 'def var chunking '//avail_hist_fields(n)%vname) + endif + + status = pio_def_var_deflate(File, varid, shuffle, deflate, & + deflate_level) + if (status /= nf90_noerr) call abort_ice( & + 'Error deflating variable '//avail_hist_fields(n)%vname) + + status = pio_put_att(File,varid,'units', & + avail_hist_fields(n)%vunit) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining units for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid, 'long_name', & + avail_hist_fields(n)%vdesc) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining long_name for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid,'coordinates', & + avail_hist_fields(n)%vcoord) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining coordinates for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid,'cell_measures', & + avail_hist_fields(n)%vcellmeas) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining cell measures for '//avail_hist_fields(n)%vname) + if (avail_hist_fields(n)%vcomment /= "none") then + call ice_pio_check(pio_put_att(File,varid,'comment', & + avail_hist_fields(n)%vcomment), & + 'put att comment '//avail_hist_fields(n)%vname) + endif + status = pio_put_att(File,varid,'missing_value',spval) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining missing_value for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid,'_FillValue',spval) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining _FillValue for '//avail_hist_fields(n)%vname) + + ! Add cell_methods attribute to variables if averaged + if (hist_avg .and. histfreq(ns) /= '1') then + status = pio_put_att(File,varid,'cell_methods','time: mean') + endif + + if (histfreq(ns) == '1' .or. .not. hist_avg) then + status = pio_put_att(File,varid,'time_rep','instantaneous') + else + status = pio_put_att(File,varid,'time_rep','averaged') + endif - status = pio_put_att(File,varid,'units', & - trim(avail_hist_fields(n)%vunit)) - status = pio_put_att(File,varid, 'long_name', & - trim(avail_hist_fields(n)%vdesc)) - status = pio_put_att(File,varid,'coordinates', & - trim(avail_hist_fields(n)%vcoord)) - status = pio_put_att(File,varid,'cell_measures', & - trim(avail_hist_fields(n)%vcellmeas)) - status = pio_put_att(File,varid,'missing_value',spval) - status = pio_put_att(File,varid,'_FillValue',spval) - - ! Add cell_methods attribute to variables if averaged - if (hist_avg .and. histfreq(ns) /= '1') then - status = pio_put_att(File,varid,'cell_methods','time: mean') endif + enddo ! num_avail_hist_fields_3Db - if (histfreq(ns) == '1' .or. .not. hist_avg) then - status = pio_put_att(File,varid,'time_rep','instantaneous') - else - status = pio_put_att(File,varid,'time_rep','averaged') + ! 4D (ice categories and layers) + + dimidcz(1) = imtid + dimidcz(2) = jmtid + dimidcz(3) = kmtidi + dimidcz(4) = cmtid + dimidcz(5) = timid + + do n = n3Dbcum + 1, n4Dicum + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_def_var(File, avail_hist_fields(n)%vname, & + pio_real, dimidcz, varid) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining variable '//avail_hist_fields(n)%vname) + + if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then + call ice_pio_check(pio_def_var_chunking(File, varid, NF90_CHUNKED, & + (/ history_chunksize_x, history_chunksize_y, 1, 1, 1 /)), & + 'def var chunking '//avail_hist_fields(n)%vname) + endif + + status = pio_def_var_deflate(File, varid, shuffle, deflate, & + deflate_level) + if (status /= nf90_noerr) call abort_ice( & + 'Error deflating variable '//avail_hist_fields(n)%vname) + + status = pio_put_att(File,varid,'units', & + avail_hist_fields(n)%vunit) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining units for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid, 'long_name', & + avail_hist_fields(n)%vdesc) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining long_name for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid,'coordinates', & + avail_hist_fields(n)%vcoord) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining coordinates for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid,'cell_measures', & + avail_hist_fields(n)%vcellmeas) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining cell measures for '//avail_hist_fields(n)%vname) + if (avail_hist_fields(n)%vcomment /= "none") then + call ice_pio_check(pio_put_att(File,varid,'comment', & + avail_hist_fields(n)%vcomment), & + 'put att comment '//avail_hist_fields(n)%vname) + endif + status = pio_put_att(File,varid,'missing_value',spval) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining missing_value for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid,'_FillValue',spval) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining _FillValue for '//avail_hist_fields(n)%vname) + + !--------------------------------------------------------------- + ! Add cell_methods attribute to variables if averaged + !--------------------------------------------------------------- + if (hist_avg .and. histfreq(ns) /= '1') then + status = pio_put_att(File,varid,'cell_methods','time: mean') + if (status /= nf90_noerr) call abort_ice( & + 'Error defining cell methods for '//avail_hist_fields(n)%vname) + endif + + if (histfreq(ns) == '1' .or. .not. hist_avg) then + status = pio_put_att(File,varid,'time_rep','instantaneous') + else + status = pio_put_att(File,varid,'time_rep','averaged') + endif endif - endif - enddo ! num_avail_hist_fields_4Di + enddo ! num_avail_hist_fields_4Di + - !----------------------------------------------------------------- ! 4D (ice categories and snow layers) - !----------------------------------------------------------------- - dimidcz(1) = imtid - dimidcz(2) = jmtid - dimidcz(3) = kmtids - dimidcz(4) = cmtid - dimidcz(5) = timid - - do n = n4Dicum + 1, n4Dscum - if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then - status = pio_def_var(File, trim(avail_hist_fields(n)%vname), & - pio_real, dimidcz, varid) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error def var '//trim(avail_hist_fields(n)%vname)) - status = pio_def_var_deflate(File, varid, shuffle, deflate, & - deflate_level) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error deflating var '//trim(avail_hist_fields(n)%vname)) - if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then - status = pio_def_var_chunking(File, varid, NF90_CHUNKED, & - (/ history_chunksize_x, history_chunksize_y, 1, 1, 1 /)) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error chunking var '//trim(avail_hist_fields(n)%vname)) - endif + dimidcz(1) = imtid + dimidcz(2) = jmtid + dimidcz(3) = kmtids + dimidcz(4) = cmtid + dimidcz(5) = timid - status = pio_put_att(File,varid,'units', & - trim(avail_hist_fields(n)%vunit)) - status = pio_put_att(File,varid, 'long_name', & - trim(avail_hist_fields(n)%vdesc)) - status = pio_put_att(File,varid,'coordinates', & - trim(avail_hist_fields(n)%vcoord)) - status = pio_put_att(File,varid,'cell_measures', & - trim(avail_hist_fields(n)%vcellmeas)) - status = pio_put_att(File,varid,'missing_value',spval) - status = pio_put_att(File,varid,'_FillValue',spval) - - ! Add cell_methods attribute to variables if averaged - if (hist_avg .and. histfreq(ns) /= '1') then - status = pio_put_att(File,varid,'cell_methods','time: mean') - endif + do n = n4Dicum + 1, n4Dscum + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_def_var(File, avail_hist_fields(n)%vname, & + pio_real, dimidcz, varid) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining variable '//avail_hist_fields(n)%vname) - if (histfreq(ns) == '1' .or. .not. hist_avg) then - status = pio_put_att(File,varid,'time_rep','instantaneous') - else - status = pio_put_att(File,varid,'time_rep','averaged') + if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then + call ice_pio_check(pio_def_var_chunking(File, varid, NF90_CHUNKED, & + (/ history_chunksize_x, history_chunksize_y, 1, 1 ,1 /)), & + 'def var chunking '//avail_hist_fields(n)%vname) + endif + + status = pio_def_var_deflate(File, varid, shuffle, deflate, & + deflate_level) + if (status /= nf90_noerr) call abort_ice( & + 'Error deflating variable '//avail_hist_fields(n)%vname) + + status = pio_put_att(File,varid,'units', & + avail_hist_fields(n)%vunit) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining units for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid, 'long_name', & + avail_hist_fields(n)%vdesc) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining long_name for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid,'coordinates', & + avail_hist_fields(n)%vcoord) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining coordinates for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid,'cell_measures', & + avail_hist_fields(n)%vcellmeas) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining cell measures for '//avail_hist_fields(n)%vname) + if (avail_hist_fields(n)%vcomment /= "none") then + call ice_pio_check(pio_put_att(File,varid,'comment', & + avail_hist_fields(n)%vcomment), & + 'put att comment '//avail_hist_fields(n)%vname) + endif + status = pio_put_att(File,varid,'missing_value',spval) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining missing_value for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid,'_FillValue',spval) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining _FillValue for '//avail_hist_fields(n)%vname) + + !--------------------------------------------------------------- + ! Add cell_methods attribute to variables if averaged + !--------------------------------------------------------------- + if (hist_avg .and. histfreq(ns) /= '1') then + status = pio_put_att(File,varid,'cell_methods','time: mean') + if (status /= nf90_noerr) call abort_ice( & + 'Error defining cell methods for '//avail_hist_fields(n)%vname) + endif + + if (histfreq(ns) == '1' .or. .not. hist_avg) then + status = pio_put_att(File,varid,'time_rep','instantaneous') + else + status = pio_put_att(File,varid,'time_rep','averaged') + endif endif - endif - enddo ! num_avail_hist_fields_4Ds + enddo ! num_avail_hist_fields_4Ds - !----------------------------------------------------------------- ! 4D (ice categories and biology layers) - !----------------------------------------------------------------- - dimidcz(1) = imtid - dimidcz(2) = jmtid - dimidcz(3) = kmtidb - dimidcz(4) = cmtid - dimidcz(5) = timid - - do n = n4Dscum + 1, n4Dbcum - if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then - status = pio_def_var(File, trim(avail_hist_fields(n)%vname), & - pio_real, dimidcz, varid) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error def var '//trim(avail_hist_fields(n)%vname)) - status = pio_def_var_deflate(File, varid, shuffle, deflate, & - deflate_level) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error deflating var '//trim(avail_hist_fields(n)%vname)) - if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then - status = pio_def_var_chunking(File, varid, NF90_CHUNKED, & - (/ history_chunksize_x, history_chunksize_y, 1, 1, 1 /)) - if (status /= nf90_noerr) call abort_ice( & - 'ice: Error chunking var '//trim(avail_hist_fields(n)%vname)) - endif + dimidcz(1) = imtid + dimidcz(2) = jmtid + dimidcz(3) = kmtidb + dimidcz(4) = cmtid + dimidcz(5) = timid - status = pio_put_att(File,varid,'units', & - trim(avail_hist_fields(n)%vunit)) - status = pio_put_att(File,varid, 'long_name', & - trim(avail_hist_fields(n)%vdesc)) - status = pio_put_att(File,varid,'coordinates', & - trim(avail_hist_fields(n)%vcoord)) - status = pio_put_att(File,varid,'cell_measures', & - trim(avail_hist_fields(n)%vcellmeas)) - status = pio_put_att(File,varid,'missing_value',spval) - status = pio_put_att(File,varid,'_FillValue',spval) - - ! Add cell_methods attribute to variables if averaged - if (hist_avg .and. histfreq(ns) /= '1') then - status = pio_put_att(File,varid,'cell_methods','time: mean') - endif + do n = n4Dscum + 1, n4Dbcum + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_def_var(File, avail_hist_fields(n)%vname, & + pio_real, dimidcz, varid) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining variable '//avail_hist_fields(n)%vname) - if (histfreq(ns) == '1' .or. .not. hist_avg) then - status = pio_put_att(File,varid,'time_rep','instantaneous') - else - status = pio_put_att(File,varid,'time_rep','averaged') + if (history_chunksize_x > 0 .and. history_chunksize_y > 0) then + call ice_pio_check(pio_def_var_chunking(File, varid, NF90_CHUNKED, & + (/ history_chunksize_x, history_chunksize_y, 1, 1 /)), & + 'def var chunking '//avail_hist_fields(n)%vname) + endif + + status = pio_def_var_deflate(File, varid, shuffle, deflate, & + deflate_level) + if (status /= nf90_noerr) call abort_ice( & + 'Error deflating variable '//avail_hist_fields(n)%vname) + + status = pio_put_att(File,varid,'units', & + avail_hist_fields(n)%vunit) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining units for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid, 'long_name', & + avail_hist_fields(n)%vdesc) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining long_name for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid,'coordinates', & + avail_hist_fields(n)%vcoord) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining coordinates for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid,'cell_measures', & + avail_hist_fields(n)%vcellmeas) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining cell measures for '//avail_hist_fields(n)%vname) + if (avail_hist_fields(n)%vcomment /= "none") then + call ice_pio_check(pio_put_att(File,varid,'comment', & + avail_hist_fields(n)%vcomment), & + 'put att comment '//avail_hist_fields(n)%vname) + endif + status = pio_put_att(File,varid,'missing_value',spval) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining missing_value for '//avail_hist_fields(n)%vname) + status = pio_put_att(File,varid,'_FillValue',spval) + if (status /= nf90_noerr) call abort_ice( & + 'Error defining _FillValue for '//avail_hist_fields(n)%vname) + + !--------------------------------------------------------------- + ! Add cell_methods attribute to variables if averaged + !--------------------------------------------------------------- + if (hist_avg .and. histfreq(ns) /= '1' ) then + status = pio_put_att(File,varid,'cell_methods','time: mean') + if (status /= nf90_noerr) call abort_ice( & + 'Error defining cell methods for '//avail_hist_fields(n)%vname) + endif + + if (histfreq(ns) == '1' .or. .not. hist_avg) then + call ice_pio_check(pio_put_att(File,varid,'time_rep','instantaneous'), & + 'put att time_rep instantaneous') + else + call ice_pio_check(pio_put_att(File,varid,'time_rep','averaged'), & + 'put att time_rep averaged') + endif endif - endif - enddo ! num_avail_hist_fields_4Db + enddo ! num_avail_hist_fields_4Db - !----------------------------------------------------------------- - ! global attributes - !----------------------------------------------------------------- - ! ... the user should change these to something useful ... - !----------------------------------------------------------------- -#ifdef CCSMCOUPLED - status = pio_put_att(File,pio_global,'title',runid) -#else - title = 'sea ice model output for CICE' - status = pio_put_att(File,pio_global,'title',trim(title)) -#endif - title = 'Diagnostic and Prognostic Variables' - status = pio_put_att(File,pio_global,'contents',trim(title)) + title = 'sea ice model output for CICE' + call ice_pio_check(pio_put_att(File,nf90_global,'title',title), & + 'global attribute title') + + title = 'Diagnostic and Prognostic Variables' + call ice_pio_check(pio_put_att(File,nf90_global,'contents',title), & + 'global attribute contents') - title = 'Los Alamos Sea Ice Model (CICE) Version 5' - status = pio_put_att(File,pio_global,'source',trim(title)) + title = 'Los Alamos Sea Ice Model (CICE) Version 5' + call ice_pio_check(pio_put_att(File,nf90_global,'source',title), & + 'global attribute source') - select case (histfreq(ns)) + select case (histfreq(ns)) case ("y", "Y") write(time_period_freq,'(a,i0)') 'year_',histfreq_n(ns) case ("m", "M") @@ -902,71 +1107,72 @@ subroutine ice_write_hist (ns) write(time_period_freq,'(a,i0)') 'hour_',histfreq_n(ns) case ("1") write(time_period_freq,'(a,i0)') 'step_',histfreq_n(ns) - end select - - status = pio_put_att(File,pio_global,'time_period_freq',trim(time_period_freq)) - - if (use_leap_years) then - write(title,'(a,i3,a)') 'This year has ',int(dayyr),' days' - else - write(title,'(a,i3,a)') 'All years have exactly ',int(dayyr),' days' - endif - status = pio_put_att(File,pio_global,'comment',trim(title)) - - write(title,'(a,i8.8)') 'File written on model date ',idate - status = pio_put_att(File,pio_global,'comment2',trim(title)) + end select - write(title,'(a,i6)') 'seconds elapsed into model date: ',sec - status = pio_put_att(File,pio_global,'comment3',trim(title)) + status = pio_put_att(File,pio_global,'time_period_freq',trim(time_period_freq)) - title = 'CF-1.0' - status = & - pio_put_att(File,pio_global,'conventions',trim(title)) - - call date_and_time(date=current_date, time=current_time) - write(start_time,1000) current_date(1:4), current_date(5:6), & - current_date(7:8), current_time(1:2), & - current_time(3:4) +#ifdef AusCOM + write(title,'(a,i3,a)') 'This Year Has ',int(dayyr),' days' +#else + if (use_leap_years) then + write(title,'(a,i3,a)') 'This year has ',int(dayyr),' days' + else + write(title,'(a,i3,a)') 'All years have exactly ',int(dayyr),' days' + endif +#endif + call ice_pio_check(pio_put_att(File,nf90_global,'comment',title), & + 'global attribute comment') + + write(title,'(a,i8.8)') 'File started on model date ',idate + call ice_pio_check(pio_put_att(File,nf90_global,'comment2',title), & + 'global attribute date1') + + ! TO-DO: Update output for CF compliance ! + ! title = 'CF-1.0' + ! call ice_pio_check(pio_put_att(File,nf90_global,'conventions',title), & + ! 'global attribute conventions') + + call date_and_time(date=current_date, time=current_time) + write(start_time,1000) current_date(1:4), current_date(5:6), & + current_date(7:8), current_time(1:2), & + current_time(3:4), current_time(5:8) 1000 format('This dataset was created on ', & - a,'-',a,'-',a,' at ',a,':',a) - status = pio_put_att(File,pio_global,'history',trim(start_time)) + a,'-',a,'-',a,' at ',a,':',a,':',a) + + call ice_pio_check(pio_put_att(File,nf90_global,'history',start_time), & + 'global attribute history') - status = pio_put_att(File,pio_global,'io_flavor','io_pio') + call ice_pio_check(pio_put_att(File,nf90_global,'io_flavor','io_pio'), & + 'global attribute io_flavor') !----------------------------------------------------------------- ! end define mode !----------------------------------------------------------------- - status = pio_enddef(File) + call ice_pio_check(pio_enddef(File), 'enddef') - !----------------------------------------------------------------- - ! write time variable - !----------------------------------------------------------------- +end subroutine ice_hist_create - status = pio_inq_varid(File,'time',varid) - status = pio_put_var(File,varid,(/1/),ltime) +subroutine write_coordinate_variables(File, coord_var, var_nz) - !----------------------------------------------------------------- - ! write time_bounds info - !----------------------------------------------------------------- + type(file_desc_t), intent(inout) :: File + type(coord_attributes), dimension(ncoord), intent(in) :: coord_var + type(coord_attributes), dimension(nvarz) :: var_nz - if (hist_avg .and. histfreq(ns) /= '1') then - status = pio_inq_varid(File,'time_bounds',varid) - time_bounds=(/time_beg(ns),time_end(ns)/) - bnd_start = (/1,1/) - bnd_length = (/2,1/) - status = pio_put_var(File,varid,ival=time_bounds, & - start=bnd_start(:),count=bnd_length(:)) - endif + real (kind=real_kind), allocatable :: workr2(:,:,:) + + integer :: i, k, status + type(var_desc_t) :: varid !----------------------------------------------------------------- ! write coordinate variables !----------------------------------------------------------------- - + call pio_seterrorhandling(File, PIO_RETURN_ERROR) allocate(workr2(nx_block,ny_block,nblocks)) do i = 1,ncoord - status = pio_inq_varid(File, coord_var(i)%short_name, varid) + call ice_pio_check(pio_inq_varid(File, coord_var(i)%short_name, varid), & + 'inquire varid '//coord_var(i)%short_name) SELECT CASE (coord_var(i)%short_name) CASE ('TLON') ! Convert T grid longitude from -180 -> 180 to 0 to 360 @@ -980,11 +1186,12 @@ subroutine ice_write_hist (ns) END SELECT call pio_write_darray(File, varid, iodesc2d, & workr2, status, fillval=spval) + call ice_pio_check(status, 'write darray '//coord_var(i)%short_name) enddo ! Extra dimensions (NCAT, VGRD*) - do i = 1, nvarz + do i = 1, nvarz if (igrdz(i)) then status = pio_inq_varid(File, var_nz(i)%short_name, varid) SELECT CASE (var_nz(i)%short_name) @@ -999,23 +1206,29 @@ subroutine ice_write_hist (ns) END SELECT endif enddo + call pio_seterrorhandling(File, PIO_INTERNAL_ERROR) - !----------------------------------------------------------------- - ! write grid masks, area and rotation angle - !----------------------------------------------------------------- + deallocate(workr2) -! if (igrd(n_tmask)) then -! status = pio_inq_varid(File, 'tmask', varid) -! call pio_write_darray(File, varid, iodesc2d, & -! hm(:,:,1:nblocks), status, fillval=spval) -! endif -! if (igrd(n_blkmask)) then -! status = pio_inq_varid(File, 'blkmask', varid) -! call pio_write_darray(File, varid, iodesc2d, & -! bm(:,:,1:nblocks), status, fillval=spval) -! endif - - do i = 1, nvar ! note: n_tmask=1, n_blkmask=2 +end subroutine write_coordinate_variables + + + +subroutine write_grid_variables(File, var, var_nverts) + + type(file_desc_t), intent(inout) :: File + type(req_attributes), dimension(nvar), intent(in) :: var + type(coord_attributes), dimension(nvar_verts), intent(in) :: var_nverts + + real (kind=real_kind), allocatable :: workr2(:,:,:) + real (kind=real_kind), allocatable :: workr3v(:,:,:,:) + + integer :: ivertex, i, status + type(var_desc_t) :: varid + + allocate(workr2(nx_block,ny_block,nblocks)) + + do i = 1, nvar ! note: n_tmask=1, n_blkmask=2 if (igrd(i)) then SELECT CASE (var(i)%req%short_name) CASE ('tmask') @@ -1083,12 +1296,26 @@ subroutine ice_write_hist (ns) deallocate(workr3v) endif ! f_bounds + deallocate(workr2) - !----------------------------------------------------------------- - ! write variable data - !----------------------------------------------------------------- +end subroutine write_grid_variables + + +subroutine write_2d_variables(ns, File, i_time) + + integer, intent(in) :: ns, i_time + type(file_desc_t), intent(inout) :: File + + real (kind=real_kind), allocatable :: workr2(:,:,:) + + integer :: n, status + type(var_desc_t) :: varid + integer (kind=PIO_OFFSET_KIND) :: FRAME_1 + + FRAME_1 = i_time + + allocate(workr2(nx_block,ny_block,nblocks)) - ! 2D do n=1,num_avail_hist_fields_2D if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then status = pio_inq_varid(File,avail_hist_fields(n)%vname,varid) @@ -1103,7 +1330,24 @@ subroutine ice_write_hist (ns) deallocate(workr2) - ! 3D (category) +end subroutine write_2d_variables + + +subroutine write_3d_and_4d_variables(ns, File, i_time) + + integer, intent(in) :: ns, i_time + type(file_desc_t), intent(inout) :: File + + real (kind=real_kind), allocatable :: workr3(:,:,:,:) + real (kind=real_kind), allocatable :: workr4(:,:,:,:,:) + + type(var_desc_t) :: varid + integer :: n, nn, i, j, k, status + integer (kind=PIO_OFFSET_KIND) :: FRAME_1 + + FRAME_1 = i_time + + ! 3D (category) allocate(workr3(nx_block,ny_block,nblocks,ncat_hist)) do n = n2D + 1, n3Dccum nn = n - n2D @@ -1207,30 +1451,8 @@ subroutine ice_write_hist (ns) enddo ! num_avail_hist_fields_4Ds deallocate(workr4) -! similarly for num_avail_hist_fields_4Db (define workr4b, iodesc4db) - - !----------------------------------------------------------------- - ! close output dataset - !----------------------------------------------------------------- - - call pio_closefile(File) - if (my_task == master_task) then - write(nu_diag,*) ' ' - write(nu_diag,*) 'Finished writing ',trim(ncfile(ns)) - endif - - !----------------------------------------------------------------- - ! clean-up PIO descriptors - !----------------------------------------------------------------- - call pio_freedecomp(ice_pio_subsystem, iodesc2d) - call pio_freedecomp(ice_pio_subsystem, iodesc3dv) - call pio_freedecomp(ice_pio_subsystem, iodesc3dc) - call pio_freedecomp(ice_pio_subsystem, iodesc3di) - call pio_freedecomp(ice_pio_subsystem, iodesc3db) - call pio_freedecomp(ice_pio_subsystem, iodesc4di) - +end subroutine write_3d_and_4d_variables - end subroutine ice_write_hist !======================================================================= diff --git a/io_pio/ice_pio.F90 b/io_pio/ice_pio.F90 index 62090241..626fd86a 100644 --- a/io_pio/ice_pio.F90 +++ b/io_pio/ice_pio.F90 @@ -19,7 +19,7 @@ module ice_pio implicit none private - public :: ice_pio_subsystem + public :: ice_pio_subsystem, ice_pio_check save interface ice_pio_initdecomp @@ -33,7 +33,7 @@ module ice_pio public ice_pio_initfile public ice_pio_initdecomp logical :: pio_initialized = .false. - integer :: pio_iotype + integer, parameter :: pio_iotype = pio_iotype_netcdf4p type(iosystem_desc_t) :: ice_pio_subsystem !=============================================================================== @@ -44,36 +44,32 @@ module ice_pio subroutine ice_pio_init() - integer :: num_iotasks, stride, ierr, num_agg + integer (kind=int_kind) :: num_iotasks, status + integer (kind=int_kind), parameter :: num_agg = 1 , stride = 1 character(*),parameter :: subName = '(ice_pio_init) ' if (pio_initialized) then return endif - - pio_iotype = pio_iotype_netcdf4p - - num_agg = 1 - stride = 1 num_iotasks = get_num_procs() / stride call pio_init(my_task, MPI_COMM_ICE, num_iotasks, num_agg, stride, PIO_rearr_box, ice_pio_subsystem) - ! PIO needs to be compiled with --enable-debug to use this - ierr = pio_set_log_level(0) + ! PIO needs to be compiled with -DPIO_ENABLE_LOGGING=ON (+logging in spack) to use this + status = pio_set_log_level(0) + call ice_pio_check(status, subname//' ERROR: Failed to set pio log level') pio_initialized = .true. end subroutine ice_pio_init -subroutine ice_pio_initfile(mode, filename, File, clobber, cdf64) +subroutine ice_pio_initfile(mode, filename, File, clobber) implicit none character(len=*) , intent(in), optional :: mode character(len=*) , intent(in), optional :: filename type(file_desc_t) , intent(inout), optional :: File logical , intent(in), optional :: clobber - logical , intent(in), optional :: cdf64 ! local variables @@ -85,6 +81,8 @@ subroutine ice_pio_initfile(mode, filename, File, clobber, cdf64) integer :: status character(*),parameter :: subName = '(ice_pio_wopen) ' + call pio_seterrorhandling(ice_pio_subsystem, PIO_RETURN_ERROR) + if (present(mode) .and. present(filename) .and. present(File)) then if (trim(mode) == 'write') then @@ -97,17 +95,20 @@ subroutine ice_pio_initfile(mode, filename, File, clobber, cdf64) if (exists) then if (lclobber) then status = pio_createfile(ice_pio_subsystem, File, pio_iotype, trim(filename), PIO_clobber) + call ice_pio_check(status, subname//' ERROR: Failed to overwrite file '//trim(filename)) if (my_task == master_task) then write(nu_diag,*) subname,' create file ',trim(filename) end if else status = pio_openfile(ice_pio_subsystem, File, pio_iotype, trim(filename), pio_write) + call ice_pio_check( status, subname//' ERROR: Failed to open file '//trim(filename)) if (my_task == master_task) then write(nu_diag,*) subname,' open file ',trim(filename) end if endif else status = pio_createfile(ice_pio_subsystem, File, pio_iotype, trim(filename), pio_noclobber) + call ice_pio_check( status, subname//' ERROR: Failed to create file '//trim(filename)) if (my_task == master_task) then write(nu_diag,*) subname,' create file ',trim(filename) end if @@ -115,22 +116,25 @@ subroutine ice_pio_initfile(mode, filename, File, clobber, cdf64) else ! filename is already open, just return endif - end if - - if (trim(mode) == 'read') then + else if (trim(mode) == 'read') then inquire(file=trim(filename),exist=exists) if (exists) then status = pio_openfile(ice_pio_subsystem, File, pio_iotype, trim(filename), pio_nowrite) + call ice_pio_check( status, subname//' ERROR: Failed to open file '//trim(filename)) else if(my_task==master_task) then write(nu_diag,*) 'ice_pio_ropen ERROR: file invalid ',trim(filename) end if call abort_ice('aborting in ice-pio_ropen with invalid file') endif + else + call abort_ice(subName//' : mode='//trim(mode)//' not recognised') end if end if + call pio_seterrorhandling(ice_pio_subsystem, PIO_INTERNAL_ERROR) + end subroutine ice_pio_initfile !================================================================================ @@ -142,7 +146,7 @@ subroutine ice_pio_initdecomp_2d(iodesc, use_double) logical :: luse_double integer (kind=int_kind) :: & - iblk,ilo,ihi,jlo,jhi,lon,lat,i,j,n,k + iblk,ilo,ihi,jlo,jhi,lon,lat,i,j,n type(block) :: this_block @@ -383,7 +387,32 @@ subroutine ice_pio_initdecomp_4d (ndim3, ndim4, iodesc) deallocate(dof4d) end subroutine ice_pio_initdecomp_4d - + + +!================================================================================ + + ! PIO Error handling + ! Author: Anton Steketee, ACCESS-NRI + + subroutine ice_pio_check(status, abort_msg) + use ice_exit, only: abort_ice + + integer(kind=int_kind), intent (in) :: status + character (len=*) , intent (in) :: abort_msg + + ! local variables + + character(len=pio_max_name) :: err_msg + integer(kind=int_kind) :: strerror_status + character(len=*), parameter :: subname = '(ice_pio_check)' + + if (status /= PIO_NOERR) then + strerror_status = pio_strerror(status, err_msg) + call abort_ice(subname//trim(err_msg)//', '//trim(abort_msg)) + endif + end subroutine ice_pio_check + + !================================================================================ end module ice_pio diff --git a/source/ice_calendar.F90 b/source/ice_calendar.F90 index a1bad8c5..9660d898 100644 --- a/source/ice_calendar.F90 +++ b/source/ice_calendar.F90 @@ -103,8 +103,9 @@ module ice_calendar write_history(max_nstrm) ! write history now character (len=1), public :: & - histfreq(max_nstrm), & ! history output frequency, 'y','m','d','h','1' - dumpfreq ! restart frequency, 'y','m','d' + histfreq(max_nstrm), & ! history output frequency, 'y','m','d','h','1' + hist_file_freq(max_nstrm), & ! history output file save frequency, 'y','m','d','h','1' + dumpfreq ! restart frequency, 'y','m','d' character (len=char_len),public :: calendar_type diff --git a/source/ice_history_shared.F90 b/source/ice_history_shared.F90 index d89f696b..da913c8b 100644 --- a/source/ice_history_shared.F90 +++ b/source/ice_history_shared.F90 @@ -25,15 +25,18 @@ module ice_history_shared use ice_kinds_mod - use ice_fileunits, only: nu_diag use ice_domain_size, only: ncat, nilyr, nslyr, nblyr, max_nstrm implicit none save private - public :: define_hist_field, accum_hist_field, icefields_nml, construct_filename - + public :: define_hist_field, accum_hist_field, icefields_nml, & +#ifdef ncdf + check, & +#endif + construct_filename + logical (kind=log_kind), public :: & hist_avg ! if true, write averaged data instead of snapshots @@ -465,20 +468,44 @@ module ice_history_shared !======================================================================= - subroutine construct_filename(ncfile,suffix,ns) +#ifdef ncdf + subroutine check(status, msg) + use netcdf, only: nf90_noerr, nf90_strerror + use ice_fileunits, only: nu_diag, ice_stderr, ice_stdout + use ice_exit, only: abort_ice + + integer, intent (in) :: status + character(len=*), intent (in) :: msg + + if(status /= nf90_noerr) then + !sometimes the netcdf error string is quite long, so print seperately to prevent overrun + write(nu_diag,*) trim(nf90_strerror(status)) + write (ice_stdout,*) trim(nf90_strerror(status)) + write (ice_stderr,*) trim(nf90_strerror(status)) + call abort_ice('ice: NetCDF error '//trim(msg)) + end if + end subroutine check +#endif + + subroutine construct_filename(ncfile,suffix,ns,time_string) + + ! construct filenames for history output + ! we follow cosima convention: + ! e.g. ice-1daily-mean_0001-01.nc for daily data in january 0001 use ice_calendar, only: time, sec, nyr, month, daymo, & - mday, write_ic, histfreq, histfreq_n, & + mday, write_ic, histfreq, hist_file_freq, histfreq_n, & year_init, new_year, new_month, new_day, & dt use ice_restart_shared, only: lenstr - character (len=*), intent(inout) :: ncfile + character (char_len_long), intent(inout) :: ncfile + character (char_len), intent(out), optional :: time_string + character (char_len) :: ldate_string character (len=2), intent(in) :: suffix integer (kind=int_kind), intent(in) :: ns integer (kind=int_kind) :: iyear, imonth, iday, isec - character (len=1) :: cstream iyear = nyr + year_init - 1 ! set year_init=1 in ice_in to get iyear=nyr imonth = month @@ -490,9 +517,9 @@ subroutine construct_filename(ncfile,suffix,ns) #endif ! construct filename if (write_ic) then - write(ncfile,'(a,a,i4.4,a,i2.2,a,i2.2,a,i5.5,a,a)') & - incond_file(1:lenstr(incond_file)),'.',iyear,'-', & - imonth,'-',iday,'-',isec,'.',suffix + ncfile=incond_file(1:lenstr(incond_file)) + write(ldate_string,'(i4.4,a,i2.2,a,i2.2,a,i5.5)') & + iyear,'-',imonth,'-',iday,'-',sec else if (hist_avg .and. histfreq(ns) /= '1') then @@ -510,41 +537,56 @@ subroutine construct_filename(ncfile,suffix,ns) endif endif - cstream = '' -!echmod ! this was implemented for CESM but it breaks post-processing software -!echmod ! of other groups (including RASM which uses CCSMCOUPLED) -!echmod if (ns > 1) write(cstream,'(i1.1)') ns-1 - - if (histfreq(ns) == '1') then ! instantaneous, write every dt - write(ncfile,'(a,a,i4.4,a,i2.2,a,i2.2,a,i5.5,a,a)') & - history_file(1:lenstr(history_file))//trim(cstream),'_inst.', & - iyear,'-',imonth,'-',iday,'-',sec,'.',suffix - - elseif (hist_avg) then ! write averaged data - - if (histfreq(ns) == 'd'.or.histfreq(ns) == 'D') then ! daily - write(ncfile,'(a,a,i4.4,a,i2.2,a,i2.2,a,a)') & - history_file(1:lenstr(history_file))//trim(cstream), & - '.',iyear,'-',imonth,'-',iday,'.',suffix - elseif (histfreq(ns) == 'h'.or.histfreq(ns) == 'H') then ! hourly - write(ncfile,'(a,a,i2.2,a,i4.4,a,i2.2,a,i2.2,a,i5.5,a,a)') & - history_file(1:lenstr(history_file))//trim(cstream),'_', & - histfreq_n(ns),'h.',iyear,'-',imonth,'-',iday,'-',sec,'.',suffix - elseif (histfreq(ns) == 'm'.or.histfreq(ns) == 'M') then ! monthly - write(ncfile,'(a,a,i4.4,a,i2.2,a,a)') & - history_file(1:lenstr(history_file))//trim(cstream),'.', & - iyear,'-',imonth,'.',suffix - elseif (histfreq(ns) == 'y'.or.histfreq(ns) == 'Y') then ! yearly - write(ncfile,'(a,a,i4.4,a,a)') & - history_file(1:lenstr(history_file))//trim(cstream),'.', & - iyear,'.',suffix - endif + ncfile=history_file(1:lenstr(history_file)) + + ! frequency of history ouput (typically 1) + if (histfreq_n(ns)>9) then + write(ldate_string,'(i2)') histfreq_n(ns) + else + write(ldate_string,'(i1)') histfreq_n(ns) + endif - else ! instantaneous with histfreq > dt - write(ncfile,'(a,a,i4.4,a,i2.2,a,i2.2,a,i5.5,a,a)') & - history_file(1:lenstr(history_file)),'_inst.', & - iyear,'-',imonth,'-',iday,'-',sec,'.',suffix + ncfile=ncfile(1:lenstr(ncfile))//'-'//trim(ldate_string) + + ! name file based on history frequency (e.g. "daily-mean") + if (histfreq(ns) == 'd'.or.histfreq(ns) == 'D') then ! daily + ncfile=ncfile(1:lenstr(ncfile))//'daily' + elseif (histfreq(ns) == 'h'.or.histfreq(ns) == 'H') then ! hourly + ncfile=ncfile(1:lenstr(ncfile))//'hourly' + elseif (histfreq(ns) == 'm'.or.histfreq(ns) == 'M') then ! monthly + ncfile=ncfile(1:lenstr(ncfile))//'monthly' + elseif (histfreq(ns) == 'y'.or.histfreq(ns) == 'Y') then ! yearly + ncfile=ncfile(1:lenstr(ncfile))//'yearly' endif + + if (hist_avg) then + ncfile=ncfile(1:lenstr(ncfile))//'-mean' + else + ncfile=ncfile(1:lenstr(ncfile))//'-snap' + endif + + ! date in filename is based on history file output frequency (e.g. "0001-01" for one file per month) + if (hist_file_freq(ns) == 'd'.or.hist_file_freq(ns) == 'D') then ! daily + write(ldate_string,'(i4.4,a,i2.2,a,i2.2)') iyear,'-',imonth,'-',iday + elseif (hist_file_freq(ns) == 'h'.or.hist_file_freq(ns) == 'H') then ! hourly + write(ldate_string,'(i4.4,a,i2.2,a,i2.2,a,i5.5)') iyear,'-',imonth,'-',iday,'-',sec + elseif (hist_file_freq(ns) == 'm'.or.hist_file_freq(ns) == 'M') then ! monthly + write(ldate_string,'(i4.4,a,i2.2)') iyear,'-',imonth + elseif (hist_file_freq(ns) == 'y'.or.hist_file_freq(ns) == 'Y') then ! yearly + write(ldate_string,'(i4.4)') iyear + else !instantaneous + write(ldate_string,'(i4.4,a,i2.2,a,i2.2,a,i5.5)') iyear,'-',imonth,'-',iday,'-',sec + endif + + endif + + ! join the pieces, typically iceh-1daily-mean_0001-01.nc + ncfile=ncfile(1:lenstr(ncfile))//'_'//ldate_string(1:lenstr(ldate_string))//'.'//suffix + + ! create a string of current time for debugging + if ( present(time_string) ) then + write(time_string,'(i4.4,a,i2.2,a,i2.2,a,i5.5)') & + iyear,'-',imonth,'-',iday,'-',sec endif end subroutine construct_filename diff --git a/source/ice_init.F90 b/source/ice_init.F90 index 951b4d7d..df5f0327 100644 --- a/source/ice_init.F90 +++ b/source/ice_init.F90 @@ -53,7 +53,7 @@ subroutine input_data(forcing_start_date, cur_exp_date, & use ice_calendar, only: year_init, istep0, histfreq, histfreq_n, & dumpfreq, dumpfreq_n, diagfreq, nstreams, & npt, dt, ndtd, days_per_year, use_leap_years, & - write_ic, dump_last + write_ic, dump_last, hist_file_freq use ice_restart_shared, only: & restart, restart_ext, input_dir, input_dir, restart_dir, restart_file, & pointer_file, runid, runtype, use_restart_time, restart_format @@ -142,6 +142,7 @@ subroutine input_data(forcing_start_date, cur_exp_date, & diagfreq, diag_type, diag_file, & print_global, print_points, latpnt, lonpnt, & dbug, histfreq, histfreq_n, hist_avg, & + hist_file_freq, & history_dir, history_file, history_deflate_level, & history_parallel_io, history_chunksize_x, history_chunksize_y, & write_ic, incond_dir, incond_file @@ -230,6 +231,7 @@ subroutine input_data(forcing_start_date, cur_exp_date, & histfreq(4) = 'm' ! output frequency option for different streams histfreq(5) = 'y' ! output frequency option for different streams histfreq_n(:) = 1 ! output frequency + hist_file_freq(:) = 'x' ! default to histfreq (below) hist_avg = .true. ! if true, write time-averages (not snapshots) history_dir = './' ! write to executable dir for default history_file = 'iceh' ! history file name prefix @@ -750,6 +752,12 @@ subroutine input_data(forcing_start_date, cur_exp_date, & fbot_xfer_type = 'constant' endif + !if hist_file_freq not set, default to histfreq + if (my_task == master_task) then + do n = 1, max_nstrm + if (hist_file_freq(n)=='x' .or. hist_file_freq(n) == 'X') hist_file_freq(n) = histfreq(n) + enddo + endif call broadcast_scalar(days_per_year, master_task) call broadcast_scalar(use_leap_years, master_task) @@ -765,7 +773,10 @@ subroutine input_data(forcing_start_date, cur_exp_date, & call broadcast_scalar(diag_file, master_task) do n = 1, max_nstrm call broadcast_scalar(histfreq(n), master_task) - enddo + enddo + do n = 1, max_nstrm + call broadcast_scalar(hist_file_freq(n), master_task) + enddo call broadcast_array(histfreq_n, master_task) call broadcast_scalar(hist_avg, master_task) call broadcast_scalar(history_dir, master_task) @@ -934,6 +945,7 @@ subroutine input_data(forcing_start_date, cur_exp_date, & write(nu_diag,1010) ' print_points = ', print_points write(nu_diag,1010) ' bfbflag = ', bfbflag write(nu_diag,1050) ' histfreq = ', histfreq(:) + write(nu_diag,1050) ' hist_file_freq = ', hist_file_freq(:) write(nu_diag,1040) ' histfreq_n = ', histfreq_n(:) write(nu_diag,1010) ' hist_avg = ', hist_avg if (.not. hist_avg) write (nu_diag,*) 'History data will be snapshots' @@ -941,8 +953,10 @@ subroutine input_data(forcing_start_date, cur_exp_date, & trim(history_dir) write(nu_diag,*) ' history_file = ', & trim(history_file) - write(nu_diag,*) ' history_deflate_level = ', & + write(nu_diag,1020) ' history_deflate_level = ', & history_deflate_level + write(nu_diag,1010) ' history_parallel_io = ', & + history_parallel_io if (write_ic) then write (nu_diag,*) 'Initial condition will be written in ', & trim(incond_dir)