!||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| module io_netcdf !BOP ! !MODULE: io_netcdf ! !DESCRIPTION: ! This module provides a generic input/output interface ! for writing arrays in netCDF format. ! ! !REVISION HISTORY: ! SVN:$Id: io_netcdf.F90 2337 2006-10-30 21:54:07Z njn01 $ ! !USES: use kinds_mod use domain_size use domain use constants use communicate use boundary use broadcast use gather_scatter use exit_mod use io_types use netcdf use shr_sys_mod implicit none private save ! !PUBLIC MEMBER FUNCTIONS: public :: open_read_netcdf, & open_netcdf, & close_netcdf, & define_field_netcdf, & read_field_netcdf, & write_field_netcdf, & define_nstd_netcdf, & write_nstd_netcdf !EOP !BOC !----------------------------------------------------------------------- ! ! module variables ! !----------------------------------------------------------------------- !EOC !*********************************************************************** contains !*********************************************************************** !BOP ! !IROUTINE: open_read_netcdf ! !INTERFACE: subroutine open_read_netcdf(data_file) ! !INPUT/OUTPUT PARAMETERS: type (datafile), intent (inout) :: data_file ! !DESCRIPTION: ! This routine opens a netcdf data file and extracts global file ! attributes. ! ! !REVISION HISTORY: ! same as module !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- character (char_len) :: & path ! filename to read character (80) :: & work_line, &! temporary to use for parsing file lines att_name ! temporary to use for attribute names integer (i4) :: & iostat, &! status flag ncid, &! netCDF file id nsize, &! size parameter returned by inquire function n, &! loop index itype, &! netCDF data type att_ival, &! netCDF data type num_atts ! number of global attributes logical (log_kind) :: & att_lval ! temp space for logical attribute real (r4) :: & att_rval ! temp space for real attribute real (r8) :: & att_dval ! temp space for double attribute logical (log_kind) :: & attrib_error ! error flag for reading attributes !----------------------------------------------------------------------- ! ! set the readonly flag in the data file descriptor ! !----------------------------------------------------------------------- data_file%readonly = .true. !----------------------------------------------------------------------- ! ! open the netCDF file ! !----------------------------------------------------------------------- iostat = nf90_noerr data_file%id = 0 if (my_task == master_task) then path = trim(data_file%full_name) iostat = nf90_open(path=trim(path), mode=nf90_nowrite, ncid=ncid) call check_status(iostat) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) & call exit_POP(sigAbort,'error opening netCDF file for reading') call broadcast_scalar(ncid, master_task) data_file%id(1) = ncid !----------------------------------------------------------------------- ! ! determine number of global file attributes ! !----------------------------------------------------------------------- if (my_task == master_task) then iostat = nf90_Inquire(ncid, nAttributes = num_atts) end if call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) & call exit_POP(sigAbort, & 'error getting number of netCDF global attributes') call broadcast_scalar(num_atts, master_task) !----------------------------------------------------------------------- ! ! now read each attribute and set attribute values ! !----------------------------------------------------------------------- do n=1,num_atts !*** !*** get attribute name !*** att_name = char_blank if (my_task == master_task) then iostat = nf90_inq_attname(ncid, NF90_GLOBAL, n, att_name) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) & call exit_POP(sigAbort, & 'error getting netCDF global attribute name') call broadcast_scalar(att_name, master_task) !*** !*** check to see if name matches any of the standard file !*** attributes !*** select case(trim(att_name)) case('title') data_file%title = char_blank if (my_task == master_task) then iostat = nf90_inquire_attribute(ncid, NF90_GLOBAL, & name='title', len=nsize) if (iostat == nf90_noerr) then iostat = nf90_get_att(ncid=ncid, varid=NF90_GLOBAL, & name='title',values=data_file%title(1:nsize)) call check_status(iostat) endif endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading title from netCDF file') endif call broadcast_scalar(data_file%title, master_task) case('history') data_file%history = char_blank if (my_task == master_task) then iostat = nf90_inquire_attribute(ncid, NF90_GLOBAL, & 'history',len=nsize) if (iostat == nf90_noerr) then iostat = nf90_get_att(ncid, NF90_GLOBAL, 'history', & data_file%history(1:nsize)) call check_status(iostat) endif endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading history from netCDF file') endif call broadcast_scalar(data_file%history, master_task) case('conventions') data_file%conventions = char_blank if (my_task == master_task) then iostat = nf90_inquire_attribute(ncid, NF90_GLOBAL, & 'conventions',len=nsize) if (iostat == nf90_noerr) then iostat = nf90_get_att(ncid, NF90_GLOBAL, 'conventions', & data_file%conventions(1:nsize)) call check_status(iostat) endif endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading conventions from netCDF file') endif call broadcast_scalar(data_file%conventions, master_task) case default !*** !*** if does not match any of the standard file attributes !*** add the attribute to the datafile !*** if (my_task == master_task) then iostat = nf90_Inquire_Attribute(ncid, NF90_GLOBAL, & trim(att_name), & xtype = itype, & len = nsize) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading netCDF file attribute') endif call broadcast_scalar(itype, master_task) select case (itype) case (NF90_CHAR) work_line = char_blank call broadcast_scalar(nsize, master_task) if (my_task == master_task) then iostat = nf90_get_att(ncid, NF90_GLOBAL, trim(att_name),& work_line(1:nsize)) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading netCDF file attribute') endif call broadcast_scalar(work_line, master_task) call add_attrib_file(data_file, trim(att_name), & trim(work_line)) case (NF90_INT) if (my_task == master_task) then iostat = nf90_get_att(ncid, NF90_GLOBAL, & trim(att_name), att_ival) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading netCDF file attribute') endif call broadcast_scalar(att_ival, master_task) if (att_name(1:4) == 'LOG_') then !*** attribute logical work_line = att_name work_line(1:4) = ' ' att_name = adjustl(work_line) if (att_ival == 1) then att_lval = .true. else att_lval = .false. endif call add_attrib_file(data_file, trim(att_name), att_lval) else call add_attrib_file(data_file, trim(att_name), att_ival) endif case (NF90_FLOAT) if (my_task == master_task) then iostat = nf90_get_att(ncid, NF90_GLOBAL, & trim(att_name), att_rval) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading netCDF file attribute') endif call broadcast_scalar(att_rval, master_task) call add_attrib_file(data_file, trim(att_name), att_rval) case (NF90_DOUBLE) if (my_task == master_task) then iostat = nf90_get_att(ncid, NF90_GLOBAL, & trim(att_name), att_dval) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading netCDF file attribute') endif call broadcast_scalar(att_dval, master_task) call add_attrib_file(data_file, trim(att_name), att_dval) end select end select end do ! num_atts !----------------------------------------------------------------------- !EOC end subroutine open_read_netcdf !*********************************************************************** !BOP ! !IROUTINE: open_netcdf ! !INTERFACE: subroutine open_netcdf(data_file) ! !INPUT/OUTPUT PARAMETERS: type (datafile), intent (inout) :: data_file ! !DESCRIPTION: ! This routine opens a data file for writing and ! writes global file attributes. ! ! !REVISION HISTORY: ! same as module !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- character (char_len) :: & path ! temp to use for filename character (255) :: & work_line ! temp to use for character manipulation integer (i4) :: & ncid, &! netCDF id for file iostat, &! status flag for netCDF function calls itmp, &! integer temp for equivalent logical attribute n, &! loop index ncvals, &! counter for number of character attributes nlvals, &! counter for number of logical attributes nivals, &! counter for number of integer attributes nrvals, &! counter for number of real attributes ndvals ! counter for number of double attributes logical (log_kind) :: & attrib_error ! error flag for reading attributes !----------------------------------------------------------------------- ! ! open the netCDF file ! !----------------------------------------------------------------------- iostat = nf90_noerr data_file%id = 0 if (my_task==master_task) then path = trim(data_file%full_name) iostat = nf90_create(path=trim(path), cmode=nf90_write, ncid=ncid) call check_status(iostat) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) call exit_POP(sigAbort, & 'Error opening file') call broadcast_scalar(ncid, master_task) data_file%id(1) = ncid data_file%ldefine = .true. ! file in netCDF define mode !----------------------------------------------------------------------- ! ! define global file attributes ! !----------------------------------------------------------------------- attrib_error = .false. if (my_task == master_task) then !*** standard attributes iostat = nf90_put_att(ncid, NF90_GLOBAL, 'title', & trim(data_file%title)) call check_status(iostat) if (iostat /= nf90_noerr) then write(stdout,*) 'Error writing TITLE to netCDF file' attrib_error = .true. endif iostat = nf90_put_att(ncid, NF90_GLOBAL, 'history', & trim(data_file%history)) call check_status(iostat) if (iostat /= nf90_noerr) then write(stdout,*) 'Error writing HISTORY to netCDF file' attrib_error = .true. endif iostat = nf90_put_att(ncid, NF90_GLOBAL, 'conventions', & trim(data_file%conventions)) call check_status(iostat) if (iostat /= nf90_noerr) then write(stdout,*) 'Error writing CONVENTIONS to netCDF file' attrib_error = .true. endif !*** additional attributes if (associated(data_file%add_attrib_cval)) then ncvals = size(data_file%add_attrib_cval) else ncvals = 0 endif if (associated(data_file%add_attrib_lval)) then nlvals = size(data_file%add_attrib_lval) else nlvals = 0 endif if (associated(data_file%add_attrib_ival)) then nivals = size(data_file%add_attrib_ival) else nivals = 0 endif if (associated(data_file%add_attrib_rval)) then nrvals = size(data_file%add_attrib_rval) else nrvals = 0 endif if (associated(data_file%add_attrib_dval)) then ndvals = size(data_file%add_attrib_dval) else ndvals = 0 endif do n=1,ncvals work_line = data_file%add_attrib_cname(n) iostat = nf90_put_att(ncid, NF90_GLOBAL, trim(work_line), & trim(data_file%add_attrib_cval(n))) call check_status(iostat) if (iostat /= nf90_noerr) then write(stdout,*) 'Error writing ',trim(work_line) attrib_error = .true. endif end do do n=1,nlvals work_line = 'LOG_'/& &/data_file%add_attrib_lname(n) if (data_file%add_attrib_lval(n)) then itmp = 1 else itmp = 0 endif iostat = nf90_put_att(ncid, NF90_GLOBAL, trim(work_line), & itmp) call check_status(iostat) if (iostat /= nf90_noerr) then write(stdout,*) 'Error writing ',trim(work_line) attrib_error = .true. endif end do do n=1,nivals work_line = data_file%add_attrib_iname(n) iostat = nf90_put_att(ncid, NF90_GLOBAL, trim(work_line), & data_file%add_attrib_ival(n)) call check_status(iostat) if (iostat /= nf90_noerr) then write(stdout,*) 'Error writing ',trim(work_line) attrib_error = .true. endif end do do n=1,nrvals work_line = data_file%add_attrib_rname(n) iostat = nf90_put_att(ncid, NF90_GLOBAL, trim(work_line), & data_file%add_attrib_rval(n)) call check_status(iostat) if (iostat /= nf90_noerr) then write(stdout,*) 'Error writing ',trim(work_line) attrib_error = .true. endif end do do n=1,ndvals work_line = data_file%add_attrib_dname(n) iostat = nf90_put_att(ncid, NF90_GLOBAL, trim(work_line), & data_file%add_attrib_dval(n)) call check_status(iostat) if (iostat /= nf90_noerr) then write(stdout,*) 'Error writing ',trim(work_line) attrib_error = .true. endif end do endif ! master task call broadcast_scalar(attrib_error, master_task) if (attrib_error) call exit_POP(sigAbort, & 'Error writing file attributes') !----------------------------------------------------------------------- !EOC end subroutine open_netcdf !*********************************************************************** !BOP ! !IROUTINE: close_netcdf ! !INTERFACE: subroutine close_netcdf(data_file) ! !INPUT/OUTPUT PARAMETERS: type (datafile), intent (inout) :: data_file ! !DESCRIPTION: ! This routine closes an open netcdf data file. ! ! !REVISION HISTORY: ! same as module !EOP !BOC !----------------------------------------------------------------------- ! ! close a data file ! !----------------------------------------------------------------------- if (my_task == master_task) then call check_status(nf90_close(data_file%id(1))) end if !----------------------------------------------------------------------- !EOC end subroutine close_netcdf !*********************************************************************** !BOP ! !IROUTINE: define_field_netcdf ! !INTERFACE: subroutine define_field_netcdf(data_file, io_field) ! !DESCRIPTION: ! This routine defines an io field for a netCDF file. ! When reading a file, the define routine will attempt to fill an ! io field structure with meta-data information from the netCDF file. ! When writing a file, it calls the appropriate netCDF routines ! to define all the field attributes and assign a field id. ! ! !REVISION HISTORY: ! same as module ! !INPUT/OUTPUT PARAMETERS: type (datafile), intent (inout) :: & data_file ! data file in which field contained type (io_field_desc), intent (inout) :: & io_field ! field descriptor for this field !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- character (80) :: & work_line, &! workspace for manipulating input string comp_line, &! comparison string att_name ! attribute name integer (i4) :: & iostat, &! status flag for netCDF calls ncid, &! file id for netcdf file varid, &! variable id for field ndims, &! number of dimensions dimid, &! dimension id n, &! loop index ncount, &! num additional attributes nsize, &! length of character strings itype, &! netCDF data type num_atts, &! number of variable attributes att_ival, &! temp for integer attribute ncvals, &! counter for number of character attributes nlvals, &! counter for number of logical attributes nivals, &! counter for number of integer attributes nrvals, &! counter for number of real attributes ndvals ! counter for number of double attributes logical (log_kind) :: & att_lval ! temp for logical attribute real (r4) :: & att_rval ! temp for real attribute real (r8) :: & att_dval ! temp for double attribute logical (log_kind) :: & define_error ! error flag define_error = .false. ncid = data_file%id(1) !----------------------------------------------------------------------- ! ! make sure file has been opened ! !----------------------------------------------------------------------- call check_file_open(data_file, 'define_field_netcdf') !----------------------------------------------------------------------- ! ! for input files, get the variable id and determine number of field ! attributes ! !----------------------------------------------------------------------- if (data_file%readonly) then if (my_task == master_task) then iostat = NF90_INQ_VARID(ncid, trim(io_field%short_name), & io_field%id) call check_status(iostat) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) & call exit_POP(sigAbort, & 'Error finding field in netCDF input file') call broadcast_scalar(io_field%id, master_task) if (my_task == master_task) then iostat = NF90_Inquire_Variable(ncid,io_field%id,nAtts=num_atts) call check_status(iostat) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) & call exit_POP(sigAbort, & 'Error getting attribute count for netCDF field') call broadcast_scalar(num_atts, master_task) !*** !*** for each attribute, define standard attributes or add !*** attribute to io_field !*** do n=1,num_atts !*** !*** get attribute name !*** att_name = char_blank if (my_task == master_task) then iostat = nf90_inq_attname(ncid, io_field%id, n, att_name) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) & call exit_POP(sigAbort, & 'error getting netCDF field attribute name') call broadcast_scalar(att_name, master_task) !*** !*** check to see if name matches any of the standard field !*** attributes !*** select case(trim(att_name)) case('long_name') io_field%long_name = char_blank if (my_task == master_task) then iostat = nf90_inquire_attribute(ncid, io_field%id, & 'long_name', len=nsize) if (iostat == nf90_noerr) then iostat = nf90_get_att(ncid, io_field%id, 'long_name',& io_field%long_name(1:nsize)) call check_status(iostat) endif endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading long_name from netCDF file') endif call broadcast_scalar(io_field%long_name, master_task) case('units') io_field%units = char_blank if (my_task == master_task) then iostat = nf90_inquire_attribute(ncid, io_field%id, & 'units', len=nsize) if (iostat == nf90_noerr) then iostat = nf90_get_att(ncid, io_field%id, 'units', & io_field%units(1:nsize)) call check_status(iostat) endif endif call broadcast_scalar(io_field%units, master_task) case('coordinates') io_field%coordinates = char_blank if (my_task == master_task) then iostat = nf90_get_att(ncid, io_field%id, 'coordinates', & io_field%coordinates) call check_status(iostat) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading coordinates from netCDF file') endif call broadcast_scalar(io_field%coordinates, master_task) case('grid_loc') io_field%grid_loc = ' ' if (my_task == master_task) then iostat = nf90_get_att(ncid, io_field%id, 'grid_loc', & io_field%grid_loc) call check_status(iostat) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading grid_loc from netCDF file') endif call broadcast_scalar(io_field%grid_loc, master_task) case('missing_value') if (my_task == master_task) then iostat = nf90_get_att(ncid, io_field%id, & 'missing_value', & io_field%missing_value) call check_status(iostat) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading missing_value from netCDF file') endif call broadcast_scalar(io_field%missing_value, master_task) case('missing_value_i') if (my_task == master_task) then iostat = nf90_get_att(ncid, io_field%id, & 'missing_value_i', & io_field%missing_value_i) call check_status(iostat) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading missing_value_i from netCDF file') endif call broadcast_scalar(io_field%missing_value_i, master_task) case('valid_range') if (my_task == master_task) then iostat = nf90_get_att(ncid, io_field%id, & 'valid_range', & io_field%valid_range) call check_status(iostat) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading valid_range from netCDF file') endif call broadcast_array(io_field%valid_range, master_task) case default !*** !*** if does not match any of the standard file attributes !*** add the attribute to the datafile !*** if (my_task == master_task) then iostat = nf90_Inquire_Attribute(ncid, io_field%id, & trim(att_name), & xtype = itype, & len = nsize) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading netCDF file attribute') endif call broadcast_scalar(itype, master_task) select case (itype) case (NF90_CHAR) work_line = char_blank call broadcast_scalar(nsize, master_task) if (my_task == master_task) then iostat = nf90_get_att(ncid, io_field%id, & trim(att_name), & work_line(1:nsize)) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading netCDF file attribute') endif call broadcast_scalar(work_line, master_task) call add_attrib_io_field(io_field, trim(att_name), & trim(work_line)) case (NF90_INT) !*** both integer and logical attributes if (my_task == master_task) then iostat = nf90_get_att(ncid, io_field%id, & trim(att_name), att_ival) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading netCDF file attribute') endif call broadcast_scalar(att_ival, master_task) if (att_name(1:4) == 'LOG_') then !*** attribute logical work_line = att_name work_line(1:4) = ' ' att_name = adjustl(work_line) if (att_ival == 1) then att_lval = .true. else att_lval = .false. endif call add_attrib_file(data_file, trim(att_name), & att_lval) else call add_attrib_file(data_file, trim(att_name), & att_ival) endif case (NF90_FLOAT) if (my_task == master_task) then iostat = nf90_get_att(ncid, io_field%id, & trim(att_name), att_rval) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading netCDF file attribute') endif call broadcast_scalar(att_rval, master_task) call add_attrib_io_field(io_field, trim(att_name), & att_rval) case (NF90_DOUBLE) if (my_task == master_task) then iostat = nf90_get_att(ncid, io_field%id, & trim(att_name), att_dval) endif call broadcast_scalar(iostat, master_task) if (iostat /= nf90_noerr) then call exit_POP(sigAbort, & 'Error reading netCDF file attribute') endif call broadcast_scalar(att_dval, master_task) call add_attrib_io_field(io_field, trim(att_name), & att_dval) end select end select end do ! num_atts !----------------------------------------------------------------------- ! ! for output files, need to define everything ! make sure file is in define mode ! !----------------------------------------------------------------------- else ! output file if (.not. data_file%ldefine) & call exit_POP(sigAbort, & 'attempt to define field but not in define mode') !----------------------------------------------------------------------- ! ! define the dimensions ! !----------------------------------------------------------------------- ndims = io_field%nfield_dims if (my_task == master_task) then do n = 1,ndims dimid = 0 !*** check to see whether already defined iostat = NF90_INQ_DIMID(ncid=ncid, & name=trim(io_field%field_dim(n)%name),& dimid=dimid) if (iostat /= NF90_NOERR) then ! dimension not yet defined iostat = NF90_DEF_DIM (ncid=ncid, & name=trim(io_field%field_dim(n)%name), & len=io_field%field_dim(n)%length, & dimid=io_field%field_dim(n)%id) else io_field%field_dim(n)%id = dimid end if end do !----------------------------------------------------------------------- ! ! now define the field ! !----------------------------------------------------------------------- !*** check to see whether field of this name already defined. iostat = NF90_INQ_VARID(ncid, trim(io_field%short_name), varid) if (iostat /= NF90_NOERR) then ! variable was not yet defined if (associated (io_field%field_r_1d).or. & associated (io_field%field_r_2d).or. & associated (io_field%field_r_3d)) then iostat = NF90_DEF_VAR (ncid=ncid, & name=trim(io_field%short_name), & xtype=NF90_FLOAT, & dimids=(/ (io_field%field_dim(n)%id, n=1,ndims) /),& varid=io_field%id) else if ( io_field%nfield_dims == c0) then ! do not supply optional dimids for scalars iostat = NF90_DEF_VAR (ncid=ncid, & name=trim(io_field%short_name), & xtype=NF90_DOUBLE, & varid=io_field%id) else if (associated (io_field%field_d_1d).or. & associated (io_field%field_d_2d).or. & associated (io_field%field_d_3d)) then iostat = NF90_DEF_VAR (ncid=ncid, & name=trim(io_field%short_name), & xtype=NF90_DOUBLE, & dimids=(/ (io_field%field_dim(n)%id, n=1,ndims) /),& varid=io_field%id) else if (associated (io_field%field_i_1d).or. & associated (io_field%field_i_2d).or. & associated (io_field%field_i_3d)) then iostat = NF90_DEF_VAR (ncid=ncid, & name=trim(io_field%short_name), & xtype=NF90_INT, & dimids=(/ (io_field%field_dim(n)%id, n=1,ndims) /),& varid=io_field%id) else define_error = .true. end if call check_status(iostat) if (iostat /= nf90_noerr) define_error = .true. varid = io_field%id else ! Variable was previously defined, OK to use it io_field%id = varid end if end if ! master task call broadcast_scalar(define_error, master_task) if (define_error) then write(stdout,*) '(define_field_netcdf) ', trim(io_field%short_name) call exit_POP(sigAbort, 'Error defining netCDF field') endif !----------------------------------------------------------------------- ! ! Now define the field attributes ! !----------------------------------------------------------------------- if (my_task == master_task) then !*** long_name if (io_field%long_name /= char_blank) then iostat = NF90_INQUIRE_ATTRIBUTE(ncid=NCID, varid=varid, & name='long_name') if (iostat /= NF90_NOERR) then ! attrib probably not defined iostat = NF90_PUT_ATT(ncid=NCID, varid=varid, & name='long_name', & values=trim(io_field%long_name)) call check_status(iostat) if (iostat /= NF90_NOERR) define_error = .true. end if endif !*** units if (io_field%units /= char_blank) then iostat = NF90_INQUIRE_ATTRIBUTE(ncid=NCID, varid=varid, & name='units') if (iostat /= NF90_NOERR) then ! attrib probably not defined iostat = NF90_PUT_ATT(ncid=NCID, varid=varid, & name='units', & values=trim(io_field%units)) call check_status(iostat) if (iostat /= NF90_NOERR) define_error = .true. end if endif !*** coordinates if (io_field%coordinates /= char_blank) then iostat = NF90_INQUIRE_ATTRIBUTE(ncid=NCID, varid=varid, & name='coordinates') if (iostat /= NF90_NOERR) then ! attrib probably not defined iostat = NF90_PUT_ATT(ncid=NCID, varid=varid, & name='coordinates', & values=trim(io_field%coordinates)) call check_status(iostat) if (iostat /= NF90_NOERR) define_error = .true. end if endif !*** grid_loc if (io_field%grid_loc /= ' ') then iostat = NF90_INQUIRE_ATTRIBUTE(ncid=NCID, varid=varid, & name='grid_loc') if (iostat /= NF90_NOERR) then ! attrib probably not defined iostat = NF90_PUT_ATT(ncid=NCID, varid=varid, & name='grid_loc', & values=io_field%grid_loc) call check_status(iostat) if (iostat /= NF90_NOERR) define_error = .true. end if endif !*** missing_value if (io_field%missing_value /= undefined) then iostat = NF90_INQUIRE_ATTRIBUTE(ncid=NCID, varid=varid, & name='missing_value') if (iostat /= NF90_NOERR) then ! attrib probably not defined iostat = NF90_PUT_ATT(ncid=NCID, varid=varid, & name='missing_value', & values=io_field%missing_value) call check_status(iostat) if (iostat /= NF90_NOERR) define_error = .true. end if endif !*** missing_value_i if (io_field%missing_value_i == undefined_nf_int) then iostat = NF90_INQUIRE_ATTRIBUTE(ncid=NCID, varid=varid, & name='missing_value') if (iostat /= NF90_NOERR) then ! attrib probably not defined iostat = NF90_PUT_ATT(ncid=NCID, varid=varid, & name='missing_value', & values=io_field%missing_value_i) call check_status(iostat) if (iostat /= NF90_NOERR) define_error = .true. end if endif !*** valid_range(1:2) if (any(io_field%valid_range /= undefined)) then iostat = NF90_INQUIRE_ATTRIBUTE(ncid=NCID, varid=varid, & name='valid_range') if (iostat /= NF90_NOERR) then ! attrib probably not yet defined iostat = NF90_PUT_ATT(ncid=NCID, varid=varid, & name='valid_range', & values=io_field%valid_range(:)) call check_status(iostat) if (iostat /= NF90_NOERR) define_error = .true. end if endif !*** additional attributes if defined ncvals = 0 nlvals = 0 nivals = 0 nrvals = 0 ndvals = 0 if (associated(io_field%add_attrib_cval)) & ncvals = size(io_field%add_attrib_cval) if (associated(io_field%add_attrib_lval)) & nlvals = size(io_field%add_attrib_lval) if (associated(io_field%add_attrib_ival)) & nivals = size(io_field%add_attrib_ival) if (associated(io_field%add_attrib_rval)) & nrvals = size(io_field%add_attrib_rval) if (associated(io_field%add_attrib_dval)) & ndvals = size(io_field%add_attrib_dval) do n=1,ncvals iostat = NF90_PUT_ATT(ncid=NCID, varid=varid, & name=trim(io_field%add_attrib_cname(n)), & values=trim(io_field%add_attrib_cval(n))) call check_status(iostat) if (iostat /= NF90_NOERR) define_error = .true. end do do n=1,nlvals work_line = 'LOG_'/& &/trim(io_field%add_attrib_lname(n)) iostat = NF90_PUT_ATT(ncid=NCID, varid=varid, & name=trim(work_line), & values=io_field%add_attrib_ival(n)) call check_status(iostat) if (iostat /= NF90_NOERR) define_error = .true. end do do n=1,nivals iostat = NF90_PUT_ATT(ncid=NCID, varid=varid, & name=trim(io_field%add_attrib_iname(n)), & values=io_field%add_attrib_ival(n)) call check_status(iostat) if (iostat /= NF90_NOERR) define_error = .true. end do do n=1,nrvals iostat = NF90_PUT_ATT(ncid=NCID, varid=varid, & name=trim(io_field%add_attrib_rname(n)), & values=io_field%add_attrib_rval(n)) call check_status(iostat) if (iostat /= NF90_NOERR) define_error = .true. end do do n=1,ndvals iostat = NF90_PUT_ATT(ncid=NCID, varid=varid, & name=trim(io_field%add_attrib_dname(n)), & values=io_field%add_attrib_dval(n)) call check_status(iostat) if (iostat /= NF90_NOERR) define_error = .true. end do endif ! master_task call broadcast_scalar(define_error, master_task) if (define_error) call exit_POP(sigAbort, & 'Error adding attributes to field') endif ! input/output file !----------------------------------------------------------------------- !EOC end subroutine define_field_netcdf !*********************************************************************** !BOP ! !IROUTINE: write_field_netcdf ! !INTERFACE: subroutine write_field_netcdf(data_file, io_field) ! !DESCRIPTION: ! This routine writes a field to a netCDF data file. ! ! !REVISION HISTORY: ! same as module ! !INPUT/OUTPUT PARAMETERS: type (datafile), intent (inout) :: & data_file ! file to which field will be written type (io_field_desc), intent (inout) :: & io_field ! field to write to file !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- integer (i4), dimension(:,:), allocatable :: global_i_2d real (r4), dimension(:,:), allocatable :: global_r_2d real (r8), dimension(:,:), allocatable :: global_d_2d integer (i4), dimension(:), allocatable :: & start,length ! dimension quantities for netCDF integer (i4) :: & iostat, &! netCDF status flag k,n ! loop counters logical (log_kind) :: & write_error ! error flag !----------------------------------------------------------------------- ! ! exit define mode if necessary ! !----------------------------------------------------------------------- write_error = .false. if (my_task == master_task) then if (data_file%ldefine) then iostat = nf90_enddef(data_file%id(1)) data_file%ldefine = .false. call check_status(iostat) if (iostat /= nf90_noerr) write_error = .true. endif endif call broadcast_scalar(write_error, master_task) if (write_error) then write(stdout,*) '(write_field_netcdf) filename = ', & trim(data_file%full_name) call exit_POP(sigAbort, & 'Error exiting define mode in netCDF write') endif !----------------------------------------------------------------------- ! ! make sure field has been defined ! !----------------------------------------------------------------------- if (my_task == master_task) then if (io_field%id == 0) write_error = .true. endif call broadcast_scalar(write_error, master_task) if (write_error) & call exit_POP(sigAbort, & 'Attempt to write undefined field in netCDF write') !----------------------------------------------------------------------- ! ! allocate dimension start,stop quantities ! !----------------------------------------------------------------------- if (my_task == master_task) then if (associated(io_field%field_r_3d) .or. & associated(io_field%field_d_3d) .or. & associated(io_field%field_i_3d) ) then allocate(start(3), length(3) ) endif !----------------------------------------------------------------------- ! ! allocate global arrays - these are for 2-d slices of data which ! are gathered to the master task ! !----------------------------------------------------------------------- if (associated(io_field%field_r_3d) .or. & associated(io_field%field_r_2d)) then allocate(global_r_2d(nx_global,ny_global)) else if (associated(io_field%field_d_3d) .or. & associated(io_field%field_d_2d)) then allocate(global_d_2d(nx_global,ny_global)) else if (associated(io_field%field_i_3d) .or. & associated(io_field%field_i_2d)) then allocate(global_i_2d(nx_global,ny_global)) endif endif ! master task !----------------------------------------------------------------------- ! ! write data based on type ! !----------------------------------------------------------------------- !----------------------------------------------------------------------- ! ! real 3d array ! !----------------------------------------------------------------------- if (associated(io_field%field_r_3d)) then do k = 1,size(io_field%field_r_3d,dim=3) call gather_global(global_r_2d, io_field%field_r_3d(:,:,k,:), & master_task, distrb_clinic) if (my_task == master_task) then !*** tell netCDF to only write slice n io_field%field_dim(3)%start = k io_field%field_dim(3)%stop = k do n=1,3 start (n) = io_field%field_dim(n)%start length(n) = io_field%field_dim(n)%stop - start(n) + 1 end do iostat = NF90_PUT_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=global_r_2d, & start=start(:), count=length(:)) if (iostat /= nf90_noerr) then call check_status(iostat) write_error = .true. endif endif ! master task end do ! slice loop !----------------------------------------------------------------------- ! ! real 2d array ! !----------------------------------------------------------------------- else if (associated(io_field%field_r_2d)) then call gather_global(global_r_2d, io_field%field_r_2d, & master_task, distrb_clinic) if (my_task == master_task) then iostat = NF90_PUT_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=global_r_2d) if (iostat /= nf90_noerr) then call check_status(iostat) write_error = .true. endif endif ! master task !----------------------------------------------------------------------- ! ! real 1d array ! !----------------------------------------------------------------------- else if (associated(io_field%field_r_1d)) then if (my_task == master_task) then ! 1d vectors are not distributed to blocks; no need for gather_global iostat = NF90_PUT_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=io_field%field_r_1d) if (iostat /= nf90_noerr) then call check_status(iostat) write_error = .true. endif endif ! master task !----------------------------------------------------------------------- ! ! real 0d array ! !----------------------------------------------------------------------- ! deferred !----------------------------------------------------------------------- ! ! double 3d array ! !----------------------------------------------------------------------- else if (associated(io_field%field_d_3d)) then do k = 1,size(io_field%field_d_3d,dim=3) call gather_global(global_d_2d, io_field%field_d_3d(:,:,k,:), & master_task, distrb_clinic) if (my_task == master_task) then !*** tell netCDF to only write slice n io_field%field_dim(3)%start = k io_field%field_dim(3)%stop = k do n=1,io_field%nfield_dims start (n) = io_field%field_dim(n)%start length(n) = io_field%field_dim(n)%stop - start(n) + 1 end do iostat = NF90_PUT_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=global_d_2d, & start=start(:), count=length(:)) if (iostat /= nf90_noerr) then call check_status(iostat) write_error = .true. endif endif ! master task end do ! slice loop !----------------------------------------------------------------------- ! ! double 2d array ! !----------------------------------------------------------------------- else if (associated(io_field%field_d_2d)) then call gather_global(global_d_2d, io_field%field_d_2d, & master_task, distrb_clinic) if (my_task == master_task) then iostat = NF90_PUT_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=global_d_2d) if (iostat /= nf90_noerr) then call check_status(iostat) write_error = .true. endif endif ! master task !----------------------------------------------------------------------- ! ! double 1d array ! !----------------------------------------------------------------------- else if (associated(io_field%field_d_1d)) then if (my_task == master_task) then ! 1d vectors are not distributed to blocks; no need for gather_global iostat = NF90_PUT_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=io_field%field_d_1d) if (iostat /= nf90_noerr) then call check_status(iostat) write_error = .true. endif endif ! master task !----------------------------------------------------------------------- ! ! double 0d array ! !----------------------------------------------------------------------- else if ( io_field%nfield_dims == c0) then if (my_task == master_task) then ! scalars are not distributed to blocks; no need for gather_global ! for now, all scalars are r8 and are not pointers or targets iostat = NF90_PUT_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=io_field%field_d_0d) if (iostat /= nf90_noerr) then call check_status(iostat) write_error = .true. endif endif ! master task !----------------------------------------------------------------------- ! ! integer 3d array ! !----------------------------------------------------------------------- else if (associated(io_field%field_i_3d)) then do k = 1,size(io_field%field_i_3d,dim=3) call gather_global(global_i_2d, io_field%field_i_3d(:,:,k,:), & master_task, distrb_clinic) if (my_task == master_task) then !*** tell netCDF to only write slice n io_field%field_dim(3)%start = k io_field%field_dim(3)%stop = k do n=1,io_field%nfield_dims start (n) = io_field%field_dim(n)%start length(n) = io_field%field_dim(n)%stop - start(n) + 1 end do iostat = NF90_PUT_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=global_i_2d, & start=start(:), count=length(:)) if (iostat /= nf90_noerr) then call check_status(iostat) write_error = .true. endif endif ! master task end do ! slice loop !----------------------------------------------------------------------- ! ! integer 2d array ! !----------------------------------------------------------------------- else if (associated(io_field%field_i_2d)) then call gather_global(global_i_2d, io_field%field_i_2d, & master_task, distrb_clinic) if (my_task == master_task) then iostat = NF90_PUT_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=global_i_2d) if (iostat /= nf90_noerr) then call check_status(iostat) write_error = .true. endif endif ! master task !----------------------------------------------------------------------- ! ! integer 1d array ! !----------------------------------------------------------------------- else if (associated(io_field%field_i_1d)) then if (my_task == master_task) then ! 1d vectors are not distributed to blocks; no need for gather_global iostat = NF90_PUT_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=io_field%field_i_1d) if (iostat /= nf90_noerr) then call check_status(iostat) write_error = .true. endif endif ! master task !----------------------------------------------------------------------- ! ! check for write errors ! !----------------------------------------------------------------------- else call exit_POP(sigAbort, & 'No field associated for writing to netCDF') end if call broadcast_scalar(write_error, master_task) if (write_error) call exit_POP(sigAbort, & 'Error writing field to netCDF file') !----------------------------------------------------------------------- ! ! deallocate quantities ! !----------------------------------------------------------------------- if (my_task == master_task) then if (allocated(start)) deallocate(start) if (allocated(length)) deallocate(length) if (allocated(global_r_2d)) deallocate(global_r_2d) if (allocated(global_d_2d)) deallocate(global_d_2d) if (allocated(global_i_2d)) deallocate(global_i_2d) endif !----------------------------------------------------------------------- !EOC end subroutine write_field_netcdf !*********************************************************************** !BOP ! !IROUTINE: read_field_netcdf ! !INTERFACE: subroutine read_field_netcdf(data_file, io_field) ! !DESCRIPTION: ! This routine reads a field from a netcdf input file. ! ! !REVISION HISTORY: ! same as module ! !INPUT/OUTPUT PARAMETERS: type (datafile), intent (inout) :: & data_file ! file from which to read field type (io_field_desc), intent (inout) :: & io_field ! field to be read !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- integer (i4), dimension(:,:), allocatable :: global_i_2d real (r4), dimension(:,:), allocatable :: global_r_2d real (r8), dimension(:,:), allocatable :: global_d_2d integer (i4), dimension(:), allocatable :: & start,length ! dimension quantities for netCDF integer (i4) :: & iostat, &! netCDF status flag k,n ! loop counters logical (log_kind) :: & read_error ! error flag !----------------------------------------------------------------------- ! ! make sure field has been defined ! !----------------------------------------------------------------------- read_error = .false. if (my_task == master_task) then if (io_field%id == 0) read_error = .true. endif call broadcast_scalar(read_error, master_task) if (read_error) & call exit_POP(sigAbort, & 'Attempt to read undefined field in netCDF read') !----------------------------------------------------------------------- ! ! if no boundary update type defined, assume center location scalar ! !----------------------------------------------------------------------- if (io_field%field_loc == field_loc_unknown) then io_field%field_loc = field_loc_center io_field%field_type = field_type_scalar endif !----------------------------------------------------------------------- ! ! allocate dimension start,stop quantities ! !----------------------------------------------------------------------- if (my_task == master_task) then if (associated(io_field%field_r_3d) .or. & associated(io_field%field_d_3d) .or. & associated(io_field%field_i_3d) ) then allocate(start(3), length(3) ) endif !----------------------------------------------------------------------- ! ! allocate global arrays - these are for 2-d slices of data which ! are gathered to the master task ! !----------------------------------------------------------------------- if (associated(io_field%field_r_3d) .or. & associated(io_field%field_r_2d)) then allocate(global_r_2d(nx_global,ny_global)) else if (associated(io_field%field_d_3d) .or. & associated(io_field%field_d_2d)) then allocate(global_d_2d(nx_global,ny_global)) else if (associated(io_field%field_i_3d) .or. & associated(io_field%field_i_2d)) then allocate(global_i_2d(nx_global,ny_global)) endif endif ! master task !----------------------------------------------------------------------- ! ! read data based on type ! !----------------------------------------------------------------------- !----------------------------------------------------------------------- ! ! real 3d array ! !----------------------------------------------------------------------- if (associated(io_field%field_r_3d)) then do k = 1,size(io_field%field_r_3d,dim=3) if (my_task == master_task) then !*** tell netCDF to only read slice n io_field%field_dim(3)%start = k io_field%field_dim(3)%stop = k do n=1,io_field%nfield_dims start (n) = io_field%field_dim(n)%start length(n) = io_field%field_dim(n)%stop - start(n) + 1 end do iostat = NF90_GET_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=global_r_2d, & start=start(:), count=length(:)) if (iostat /= nf90_noerr) then call check_status(iostat) read_error = .true. endif endif ! master task call broadcast_scalar(read_error, master_task) if (.not. read_error) & call scatter_global(io_field%field_r_3d(:,:,k,:), & global_r_2d, master_task, distrb_clinic, & io_field%field_loc, io_field%field_type) end do ! slice loop !----------------------------------------------------------------------- ! ! real 2d array ! !----------------------------------------------------------------------- else if (associated(io_field%field_r_2d)) then if (my_task == master_task) then iostat = NF90_GET_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=global_r_2d) if (iostat /= nf90_noerr) then call check_status(iostat) read_error = .true. endif endif ! master task call broadcast_scalar(read_error, master_task) if (.not. read_error) then call scatter_global(io_field%field_r_2d, & global_r_2d, master_task, distrb_clinic, & io_field%field_loc, io_field%field_type) endif !----------------------------------------------------------------------- ! ! real 1d array ! !----------------------------------------------------------------------- else if (associated(io_field%field_r_1d)) then ! 1d vectors are not distributed to blocks; therefore, no scatter_global needed if (my_task == master_task) then iostat = NF90_GET_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=io_field%field_r_1d) if (iostat /= nf90_noerr) then call check_status(iostat) read_error = .true. endif endif ! master task call broadcast_scalar(read_error, master_task) !----------------------------------------------------------------------- ! ! real 0d array (scalar) ! !----------------------------------------------------------------------- else if (associated(io_field%field_r_1d)) then ! scalars are not distributed to blocks; therefore, no scatter_global needed if (my_task == master_task) then iostat = NF90_GET_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=io_field%field_r_0d) if (iostat /= nf90_noerr) then call check_status(iostat) read_error = .true. endif endif ! master task call broadcast_scalar(read_error, master_task) !----------------------------------------------------------------------- ! ! double 3d array ! !----------------------------------------------------------------------- else if (associated(io_field%field_d_3d)) then do k = 1,size(io_field%field_d_3d,dim=3) if (my_task == master_task) then !*** tell netCDF to only read slice n io_field%field_dim(3)%start = k io_field%field_dim(3)%stop = k do n=1,io_field%nfield_dims start (n) = io_field%field_dim(n)%start length(n) = io_field%field_dim(n)%stop - start(n) + 1 end do iostat = NF90_GET_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=global_d_2d, & start=start(:), count=length(:)) if (iostat /= nf90_noerr) then call check_status(iostat) read_error = .true. endif endif ! master task call broadcast_scalar(read_error, master_task) if (.not. read_error) & call scatter_global(io_field%field_d_3d(:,:,k,:), & global_d_2d, master_task, distrb_clinic, & io_field%field_loc, io_field%field_type) end do ! slice loop !----------------------------------------------------------------------- ! ! double 2d array ! !----------------------------------------------------------------------- else if (associated(io_field%field_d_2d)) then if (my_task == master_task) then iostat = NF90_GET_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=global_d_2d) if (iostat /= nf90_noerr) then call check_status(iostat) read_error = .true. endif endif ! master task call broadcast_scalar(read_error, master_task) if (.not. read_error) then call scatter_global(io_field%field_d_2d, & global_d_2d, master_task, distrb_clinic, & io_field%field_loc, io_field%field_type) endif !----------------------------------------------------------------------- ! ! double 1d array ! !----------------------------------------------------------------------- else if (associated(io_field%field_d_1d)) then ! 1d vectors are not distributed to blocks; therefore, no scatter_global needed if (my_task == master_task) then iostat = NF90_GET_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=io_field%field_d_1d) if (iostat /= nf90_noerr) then call check_status(iostat) read_error = .true. endif endif ! master task call broadcast_scalar(read_error, master_task) !----------------------------------------------------------------------- ! ! double 0d array (scalar) ! !----------------------------------------------------------------------- else if (associated(io_field%field_d_1d)) then ! scalars are not distributed to blocks; therefore, no scatter_global needed if (my_task == master_task) then iostat = NF90_GET_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=io_field%field_d_0d) if (iostat /= nf90_noerr) then call check_status(iostat) read_error = .true. endif endif ! master task call broadcast_scalar(read_error, master_task) !----------------------------------------------------------------------- ! ! integer 3d array ! !----------------------------------------------------------------------- else if (associated(io_field%field_i_3d)) then do k = 1,size(io_field%field_i_3d,dim=3) if (my_task == master_task) then !*** tell netCDF to only read slice n io_field%field_dim(3)%start = k io_field%field_dim(3)%stop = k do n=1,io_field%nfield_dims start (n) = io_field%field_dim(n)%start length(n) = io_field%field_dim(n)%stop - start(n) + 1 end do iostat = NF90_GET_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=global_i_2d, & start=start(:), count=length(:)) if (iostat /= nf90_noerr) then call check_status(iostat) read_error = .true. endif endif ! master task call broadcast_scalar(read_error, master_task) if (.not. read_error) & call scatter_global(io_field%field_i_3d(:,:,k,:), & global_i_2d, master_task, distrb_clinic, & io_field%field_loc, io_field%field_type) end do ! slice loop !----------------------------------------------------------------------- ! ! integer 2d array ! !----------------------------------------------------------------------- else if (associated(io_field%field_i_2d)) then if (my_task == master_task) then iostat = NF90_GET_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=global_i_2d) if (iostat /= nf90_noerr) then call check_status(iostat) read_error = .true. endif endif ! master task call broadcast_scalar(read_error, master_task) if (.not. read_error) then call scatter_global(io_field%field_i_2d, & global_i_2d, master_task, distrb_clinic, & io_field%field_loc, io_field%field_type) endif !----------------------------------------------------------------------- ! ! integer 1d array ! !----------------------------------------------------------------------- else if (associated(io_field%field_i_1d)) then ! 1d vectors are not distributed to blocks; therefore, no scatter_global needed if (my_task == master_task) then iostat = NF90_GET_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=io_field%field_i_1d) if (iostat /= nf90_noerr) then call check_status(iostat) read_error = .true. endif endif ! master task call broadcast_scalar(read_error, master_task) !----------------------------------------------------------------------- ! ! integer 0d array (scalar) ! !----------------------------------------------------------------------- else if (associated(io_field%field_i_1d)) then ! scalars are not distributed to blocks; therefore, no scatter_global needed if (my_task == master_task) then iostat = NF90_GET_VAR (ncid=data_file%id(1), & varid=io_field%id, & values=io_field%field_i_0d) if (iostat /= nf90_noerr) then call check_status(iostat) read_error = .true. endif endif ! master task call broadcast_scalar(read_error, master_task) !----------------------------------------------------------------------- ! ! check for read errors ! !----------------------------------------------------------------------- else call exit_POP(sigAbort, & 'No field associated for reading from netCDF') end if call broadcast_scalar(read_error, master_task) if (read_error) & call exit_POP(sigAbort,'Error reading field from netCDF file') !----------------------------------------------------------------------- ! ! deallocate quantities ! !----------------------------------------------------------------------- if (my_task == master_task) then if (allocated(start)) deallocate(start) if (allocated(length)) deallocate(length) if (allocated(global_r_2d)) deallocate(global_r_2d) if (allocated(global_d_2d)) deallocate(global_d_2d) if (allocated(global_i_2d)) deallocate(global_i_2d) endif !----------------------------------------------------------------------- !EOC end subroutine read_field_netcdf !*********************************************************************** !BOP ! !IROUTINE: check_status ! !INTERFACE: subroutine check_status(status) ! !DESCRIPTION: ! This exception handler subroutine can be used to check error status ! after a netcdf call. It prints out a text message assigned to ! an error code but does not exit because this routine is typically ! only called from a single process. ! ! !REVISION HISTORY: ! same as module ! !INPUT PARAMETERS: integer (i4), intent (in) :: & status ! status returned by netCDF call !EOP !BOC !----------------------------------------------------------------------- ! ! call netCDF routine to return error message ! !----------------------------------------------------------------------- if (status /= nf90_noerr) then write(stdout,*) trim(nf90_strerror(status)) call shr_sys_flush(stdout) end if !----------------------------------------------------------------------- !EOC end subroutine check_status !*********************************************************************** !BOP ! !IROUTINE: define_nstd_netcdf ! !INTERFACE: subroutine define_nstd_netcdf(data_file,ndims,io_dims,field_id, & short_name,long_name,units,coordinates, & fill_value,missing_value,nftype) ! !DESCRIPTION: ! This routine defines the nonstandard CCSM time-averaged diagnostic fields ! on nonstandard grids: MOC, N_HEAT, and N_SALT ! This routine is totally CCSM-specific ! ! ! !REVISION HISTORY: ! same as module ! !INPUT PARAMETERS: type (datafile), intent (in) :: & data_file ! data file in which field contained real (r4), intent (in) :: & fill_value, & missing_value integer (int_kind), intent(in) :: & ndims ! number of dimensions for nonstandard field character (*), intent (in) :: & short_name, & long_name, & units, & coordinates, & nftype ! !INPUT/OUTPUT PARAMETERS: type (io_dim), dimension(:), intent (inout) :: & io_dims integer (i4), intent (inout) :: & field_id ! variable id optional :: coordinates,fill_value,missing_value,nftype !EOP !BOP !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- integer (i4) :: & iostat, &! status flag for netCDF calls ncid, &! file id for netcdf file n, &! loop index xtype logical (log_kind) :: & define_error ! error flag define_error = .false. ncid = data_file%id(1) !----------------------------------------------------------------------- ! ! make sure file has been opened and is in define mode ! !----------------------------------------------------------------------- call check_file_open (data_file, 'define_nstd_netcdf') call check_definemode (data_file, 'define_nstd_netcdf') !----------------------------------------------------------------------- ! ! define the dimensions ! !----------------------------------------------------------------------- call define_dimensions(data_file,ndims,io_dims) !----------------------------------------------------------------------- ! ! define the field ! !----------------------------------------------------------------------- if (present(nftype)) then select case (trim(nftype)) case ('float','FLOAT') xtype = NF90_FLOAT case ('double','DOUBLE') xtype = NF90_DOUBLE case ('integer','INTEGER') xtype = NF90_INT case ('char','CHAR','character', 'CHARACTER') xtype = NF90_CHAR case default call exit_POP(sigAbort,'unknown nftype') end select else xtype = NF90_FLOAT endif call define_var (data_file,trim(short_name),ndims,io_dims, & xtype,field_id) !----------------------------------------------------------------------- ! ! Now define the field attributes ! !----------------------------------------------------------------------- if (my_task == master_task) then !*** long_name iostat = NF90_INQUIRE_ATTRIBUTE(ncid=NCID, varid=field_id, & name='long_name') if (iostat /= NF90_NOERR) then ! attrib probably not defined iostat = NF90_PUT_ATT(ncid=NCID, varid=field_id, & name='long_name', & values=trim(long_name)) call check_status(iostat) if (iostat /= NF90_NOERR) define_error = .true. end if !*** units iostat = NF90_INQUIRE_ATTRIBUTE(ncid=NCID, varid=field_id, & name='units') if (iostat /= NF90_NOERR) then ! attrib probably not defined iostat = NF90_PUT_ATT(ncid=NCID, varid=field_id, & name='units', & values=trim(units)) call check_status(iostat) if (iostat /= NF90_NOERR) define_error = .true. end if !*** coordinates if (present(coordinates)) then iostat = NF90_INQUIRE_ATTRIBUTE(ncid=NCID, varid=field_id, & name='coordinates') if (iostat /= NF90_NOERR) then ! attrib probably not defined iostat = NF90_PUT_ATT(ncid=NCID, varid=field_id, & name='coordinates', & values=trim(coordinates)) call check_status(iostat) if (iostat /= NF90_NOERR) define_error = .true. end if endif !*** missing_value if (present(missing_value)) then iostat = NF90_INQUIRE_ATTRIBUTE(ncid=NCID, varid=field_id, & name='missing_value') if (iostat /= NF90_NOERR) then ! attrib probably not defined iostat = NF90_PUT_ATT(ncid=NCID, varid=field_id, & name='missing_value', & values=missing_value) call check_status(iostat) if (iostat /= NF90_NOERR) define_error = .true. end if endif endif ! master_task call broadcast_scalar(define_error, master_task) if (define_error) call exit_POP(sigAbort, & '(define_nstd_netcdf) Error adding attributes to field') !----------------------------------------------------------------------- !EOC end subroutine define_nstd_netcdf !*********************************************************************** !BOP ! !IROUTINE: write_nstd_netcdf ! !INTERFACE: subroutine write_nstd_netcdf(data_file,field_id,num_writes, & ndims, io_dims, & nftype, & implied_time_dim, & indata_1d_r8, & indata_2d_r8, & indata_2d_r4, & indata_3d_r4 , & indata_4d_r4, & indata_1d_ch, & indata_2d_ch ) ! !DESCRIPTION: ! This is a specialized, CCSM-speicific routine to write any desired ! output field that cannot presently be defined through construct_io_field ! to the CCSM version of the netCDF time-averaged history output files ! ! !REVISION HISTORY: ! same as module ! !INPUT PARAMETERS: character (*), intent (in) :: & nftype integer (i4), intent (in) :: & field_id ! netCDF id for the nonstandard variables integer (int_kind), intent (in) :: & num_writes, & ndims type (io_dim), dimension(:), intent (in) :: & io_dims real (r8), dimension(:,:),intent (in) :: & indata_2d_r8 real (r8), dimension(:), intent (in) :: & indata_1d_r8 real (r4), dimension(:,:,:,:), intent (in) :: & indata_4d_r4 real (r4), dimension(:,:,:), intent (in) :: & indata_3d_r4 real (r4), dimension(:,:), intent (in) :: & indata_2d_r4 character (*), dimension(:,:), intent (in) :: & indata_2d_ch character (*), dimension(:), intent (in) :: & indata_1d_ch ! !INPUT/OUTPUT PARAMETERS: type (datafile), intent (inout) :: & data_file ! file to which field will be written logical (log_kind), intent(inout) :: & implied_time_dim optional :: & implied_time_dim, & indata_1d_r8, & indata_2d_r8, & indata_2d_r4, & indata_3d_r4, & indata_4d_r4, & indata_1d_ch, & indata_2d_ch !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- integer , dimension(2) :: & start,count ! dimension quantities for netCDF integer :: & iostat, &! netCDF status flag n ! index integer :: ncid, nout(5) logical (log_kind) :: & write_error, &! error flag supported real (r4), allocatable, dimension (:,:,:,:,:) :: & outdata_5d_r4 real (r4), allocatable, dimension (:,:,:,:) :: & outdata_4d_r4 real (r4), allocatable, dimension (:,:,:) :: & outdata_3d_r4 real (r4), allocatable, dimension (:,:) :: & outdata_2d_r4 real (r8), allocatable, dimension (:) :: & outdata_1d_r8 real (r8), allocatable, dimension (:,:) :: & outdata_2d_r8 character(char_len), allocatable, dimension (:,:) :: & outdata_2d_ch !----------------------------------------------------------------------- ! ! exit define mode if necessary ! !----------------------------------------------------------------------- write_error = .false. if (my_task == master_task) then if (data_file%ldefine) then iostat = nf90_enddef(data_file%id(1)) data_file%ldefine = .false. call check_status(iostat) if (iostat /= nf90_noerr) write_error = .true. endif endif call broadcast_scalar(write_error, master_task) if (write_error) then write(stdout,*) '(write_nstd_netcdf) filename = ', & trim(data_file%full_name) call exit_POP(sigAbort, & '(write_nstd_netcdf) Error exiting define mode in netCDF write') endif !----------------------------------------------------------------------- ! ! make sure field has been defined ! !----------------------------------------------------------------------- if (my_task == master_task) then if (field_id == 0) write_error = .true. endif call broadcast_scalar(write_error, master_task) if (write_error) & call exit_POP(sigAbort, & '(write_nstd_netcdf) Attempt to write undefined field in netCDF write') !----------------------------------------------------------------------- ! ! determine if the variable has the unlimited time dimension ! as an implicit dimension; if so, it must be made explicit ! in the outdata array ! !----------------------------------------------------------------------- if (.not. present(implied_time_dim)) then implied_time_dim = .false. endif !----------------------------------------------------------------------- ! NOTE: this version does not yet support multiple writes to the same ! netCDF file, but neither does basic pop2... !----------------------------------------------------------------------- supported = .true. if (my_task == master_task) then ncid = data_file%id(1) select case (trim(nftype)) case('double','DOUBLE') select case (implied_time_dim) case (.true.) select case (ndims) case(2) nout(1) = size(indata_1d_r8,DIM=1) allocate (outdata_2d_r8(nout(1),1)) outdata_2d_r8(:,1) = indata_1d_r8(:) iostat = NF90_PUT_VAR (ncid, field_id, outdata_2d_r8 ) deallocate (outdata_2d_r8) case default supported = .false. end select ! ndims case (.false.) select case (ndims) case(1) iostat = NF90_PUT_VAR (ncid, field_id, indata_1d_r8 ) case(2) iostat = NF90_PUT_VAR (ncid, field_id, indata_2d_r8 ) case default supported = .false. end select ! ndims end select ! implied_time_dim case('float','FLOAT') select case (implied_time_dim) case (.true.) select case (ndims) case(1) supported = .false. case(2) supported = .false. case(3) nout(1) = size(indata_3d_r4,DIM=1) nout(2) = size(indata_3d_r4,DIM=2) allocate (outdata_3d_r4(nout(1),nout(2),1)) outdata_3d_r4(:,:,1) = indata_2d_r4(:,:) iostat = NF90_PUT_VAR (ncid, field_id, outdata_3d_r4 ) deallocate (outdata_3d_r4) case(4) nout(1) = size(indata_3d_r4,DIM=1) nout(2) = size(indata_3d_r4,DIM=2) nout(3) = size(indata_3d_r4,DIM=3) allocate (outdata_4d_r4(nout(1),nout(2),nout(3),1)) outdata_4d_r4(:,:,:,1) = indata_3d_r4(:,:,:) iostat = NF90_PUT_VAR (ncid, field_id, outdata_4d_r4 ) deallocate (outdata_4d_r4) case(5) nout(1) = size(indata_4d_r4,DIM=1) nout(2) = size(indata_4d_r4,DIM=2) nout(3) = size(indata_4d_r4,DIM=3) nout(4) = size(indata_4d_r4,DIM=4) allocate (outdata_5d_r4(nout(1),nout(2),nout(3),nout(4),1)) outdata_5d_r4(:,:,:,:,1) = indata_4d_r4(:,:,:,:) iostat = NF90_PUT_VAR (ncid, field_id, outdata_5d_r4 ) deallocate (outdata_5d_r4) case default supported = .false. end select ! ndims case (.false.) select case (ndims) case(1) supported = .false. case(2) iostat = NF90_PUT_VAR (ncid, field_id, indata_2d_r4 ) case(3) iostat = NF90_PUT_VAR (ncid, field_id, indata_3d_r4 ) case(4) iostat = NF90_PUT_VAR (ncid, field_id, indata_4d_r4 ) case default supported = .false. end select ! ndims end select ! implied_time_dim case('char','character','CHAR','CHARACTER') select case (implied_time_dim) case (.true.) select case (ndims) case default supported = .false. end select ! ndims case (.false.) select case (ndims) case(2) do n=1,io_dims(2)%length start(1) = 1 start(2) = n count(1)=len_trim(indata_1d_ch(n)) count(2)=1 iostat = NF90_PUT_VAR (ncid, field_id, & trim(indata_1d_ch(n)), & start=start,count=count) enddo case default supported = .false. end select ! ndims end select ! implied_time_dim case default end select ! nftype if (iostat /= nf90_noerr) then call check_status(iostat) write_error = .true. endif endif ! master task !----------------------------------------------------------------------- ! ! check for errors ! !----------------------------------------------------------------------- call broadcast_scalar(write_error, master_task) if (write_error) call exit_POP(sigAbort, & '(write_nstd_netcdf) Error writing field to netCDF file') call broadcast_scalar(supported, master_task) if (.not. supported) call exit_POP(sigAbort, & '(write_nstd_netcdf) option not supported') !----------------------------------------------------------------------- !EOC end subroutine write_nstd_netcdf !*********************************************************************** !BOP ! !IROUTINE: define_dimensions ! !INTERFACE: subroutine define_dimensions(data_file,ndims,io_dims) ! !DESCRIPTION: ! This routine defines dimensions, if they are not already defined ! ! !REVISION HISTORY: ! same as module ! !INPUT PARAMETERS: type (datafile), intent (in) :: & data_file integer (int_kind), intent (in) :: & ndims ! !INPUT/OUTPUT PARAMETERS: type (io_dim), dimension(ndims), intent(inout) :: & io_dims !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- integer :: & iostat, & ! netCDF status flag dimid, & ncid, & n ncid = data_file%id(1) !----------------------------------------------------------------------- ! ! define the dimensions ! !----------------------------------------------------------------------- if (my_task == master_task) then do n = 1,ndims dimid = 0 !*** check to see whether dimension is already defined iostat = NF90_INQ_DIMID(ncid=ncid, name=trim(io_dims(n)%name),& dimid=dimid) if (iostat /= NF90_NOERR) then ! dimension not yet defined iostat = NF90_DEF_DIM (ncid=ncid, name=trim(io_dims(n)%name), & len=io_dims(n)%length, dimid=io_dims(n)%id) else io_dims(n)%id = dimid end if end do endif ! master_task !----------------------------------------------------------------------- !EOC end subroutine define_dimensions !*********************************************************************** !BOP ! !IROUTINE: define_var ! !INTERFACE: subroutine define_var (data_file,short_name,ndims,io_dims, & xtype,field_id) ! !DESCRIPTION: ! This routine defines a netCDF variable ! ! !REVISION HISTORY: ! same as module ! !INPUT PARAMETERS: type (datafile), intent (in) :: & data_file character(*), intent (in) :: & short_name integer (int_kind), intent (in) :: & ndims, & xtype ! !INPUT/OUTPUT PARAMETERS: type (io_dim), dimension(ndims) :: & io_dims ! !OUTPUT PARAMETERS: integer (i4), intent(out) :: & field_id !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- integer :: & iostat ! netCDF status flag integer :: & ncid, & dimid, & n logical (log_kind) :: & define_error ! error flag !----------------------------------------------------------------------- ! ! define the field ! !----------------------------------------------------------------------- define_error = .false. ncid = data_file%id(1) if (my_task == master_task) then !*** check to see whether field of this name already defined. iostat = NF90_INQ_VARID(ncid, trim(short_name), field_id) if (iostat /= NF90_NOERR) then ! variable was not yet defined iostat = NF90_DEF_VAR (ncid=ncid,name=trim(short_name), & xtype=xtype, & dimids=(/ (io_dims(n)%id, n=1,ndims) /),& varid=field_id) call check_status(iostat) if (iostat /= nf90_noerr) define_error = .true. end if end if ! master task call broadcast_scalar(define_error, master_task) if (define_error) then write(stdout,*) '(define_var) Error for field = ', trim(short_name) call exit_POP(sigAbort, 'Error defining nonstandard CCSM netCDF field') endif end subroutine define_var !*********************************************************************** !BOP ! !IROUTINE: check_file_open ! !INTERFACE: subroutine check_file_open(data_file, name) ! !DESCRIPTION: ! This utility routine checks if the data file has been opened ! ! !REVISION HISTORY: ! same as module ! !INPUT/OUTPUT PARAMETERS: type (datafile), intent (in) :: & data_file character(*),intent (in) :: name !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- integer :: & iostat ! netCDF status flag logical (log_kind) :: & define_error ! error flag character (char_len) :: string !----------------------------------------------------------------------- ! ! make sure file has been opened ! !----------------------------------------------------------------------- define_error = .false. if (data_file%id(1) <= 0) then define_error = .true. endif call broadcast_scalar(define_error, master_task) if (define_error) & call exit_POP(sigAbort, & '('//trim(name)//') attempt to define field without opening file first') !----------------------------------------------------------------------- !EOC end subroutine check_file_open !*********************************************************************** !BOP ! !IROUTINE: check_definemode ! !INTERFACE: subroutine check_definemode (data_file, name) ! !DESCRIPTION: ! This utility routine checks if the data file is in define mode ! ! !REVISION HISTORY: ! same as module ! !INPUT/OUTPUT PARAMETERS: type (datafile), intent (in) :: & data_file character(*),intent (in):: name !EOP !BOC !----------------------------------------------------------------------- ! ! local variables ! !----------------------------------------------------------------------- integer :: & iostat ! netCDF status flag logical (log_kind) :: & write_error ! error flag character (char_len) :: string !----------------------------------------------------------------------- ! ! make sure file is in define mode ! !----------------------------------------------------------------------- if (.not. data_file%ldefine) & call exit_POP(sigAbort, & '('//trim(name)//') attempt to define field but not in define mode') !----------------------------------------------------------------------- !EOC end subroutine check_definemode !*********************************************************************** end module io_netcdf !|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||