[mpas-developers] attributes in grid generation code

Michael Duda duda at ucar.edu
Wed May 19 13:57:42 MDT 2010


Hi, Todd.

Thanks for the input -- we had concluded the same here.

I've implemented the reading and writing of the attributes
on_a_sphere and sphere_radius in a hard-coded way; attached are
the changes (trunk on the left, modified on the right) to
module_grid_types.F, module_io_input.F, and module_io_output.F.
If anyone would like, I can send the modified .F files for easier
testing. These changes should allow for backward compatibility
with grid.nc files on the unit sphere; the attribute on_a_sphere
is assumed to be 'YES' if not present, and the attribute
sphere_radius is assumed to be 1.0 if not present in the grid.nc
file. Again, I've updated all of the grid files at
http://www.mmm.ucar.edu/people/duda/files/mpas/, and any existing
files can be updated using the code from
http://www.mmm.ucar.edu/people/duda/files/mpas/add_sphere_atts.f90.

If anyone has any comments or suggestions for improving the
changes, please let me know.

Cheers,
Michael


On Wed, May 19, 2010 at 11:48:13AM -0600, Todd Ringler wrote:
> 
> Hi Michael,
> 
> I think that we will be talking about a relatively small amount of  
> code that grows slowly with regards to the number of attributes, so my  
> recommendation is that we simply hard-code the read (and write) of  
> these attributes from (and to) the netCDF files.
> 
> Cheers,
> Todd
> 
> On May 18, 2010, at 11:44 AM, Michael Duda wrote:
> 
> >Hi, Everyone.
> >
> >I've updated the grid generation programs to add the attributes
> >on_a_sphere and sphere_radius to grid.nc files. I've also updated
> >the existing grid.nc files available from
> >http://www.mmm.ucar.edu/people/duda/files/mpas/ and
> >http://www.mmm.ucar.edu/people/duda/files/mpas/single-tracer.
> >
> >If you'd like to add the two new attributes to any grid.nc files
> >that you might have already, I've placed a simple Fortran program
> >to do this at
> >http://www.mmm.ucar.edu/people/duda/files/mpas/add_sphere_atts.f90
> >
> >I think the next step will be to add code in the MPAS
> >infrastructure to read and write these attributes from and to
> >files. Ideally, we would have some general mechanism -- perhaps
> >through additional registry syntax -- to specify that arbitrary
> >attributes be read and written to MPAS files. However, this might
> >get a bit involved, and the most direction solution would be to
> >just add code to handle these two attributes specifically. How
> >often we expect to need to add attributes might determine which
> >approach we take. Do any of you have any thoughts on the need to
> >add and subtract attributes more than infrequently (e.g., as
> >frequently as you might add fields in the registry)?
> >
> >Cheers,
> >Michael
> >
> >
> >On Mon, May 17, 2010 at 07:28:34PM -0600, Michael Duda wrote:
> >>Hi, All.
> >>
> >>Toward the goal of being able to distinguish between grids in the
> >>plane and grids on the sphere in MPAS, I'm proposing that we add
> >>two attributes to the grid.nc file, and make changes to our grid
> >>generation codes (global_scvt, periodic_hex, and periodic_quad)
> >>that will add these attributes to future grid.nc files. The two
> >>attributes are
> >>
> >>  on_a_sphere -- a character string, either 'YES' or 'NO'
> >>  sphere_radius -- the radius of the sphere
> >>
> >>I'd prefer on_a_sphere to be a logical, but it appears that having
> >>a logical global attribute in netCDF is not possible (which isn't
> >>to say that global attributes must be illogical!).
> >>
> >>The code changes to module_write_netcdf.F in the three grid
> >>generation codes mentioned above would essentially look like
> >>
> >>95a96,97
> >>>     real (kind=8) :: sphere_radius
> >>>     character (len=16) :: on_a_sphere
> >>104a107,109
> >>>     on_a_sphere = 'YES             '
> >>>     sphere_radius = 1.0
> >>>
> >>110a116,117
> >>>     nferr = nf_put_att_text(wr_ncid, NF_GLOBAL, 'on_a_sphere',  
> >>>16, on_a_sphere)
> >>>     nferr = nf_put_att_double(wr_ncid, NF_GLOBAL,  
> >>>'sphere_radius', NF_DOUBLE, 1, sphere_radius)
> >>
> >>with YES being replaced by NO and 1.0 being replaced by 0.0 in the
> >>case of the doubly-periodic codes.
> >>
> >>Although 'YES' requires only three characters and 'NO' requires
> >>only two, I wanted both to be the same length, and I figured
> >>making the attribute be 16 characters rather than three would
> >>allow for some potential flexibility in future without adding an
> >>unreasonable amount of storage. For example, if not on the sphere,
> >>we might in future find it useful to be able to set on_a_sphere to
> >>something like 'NO (X-Y PLANE)'. Who knows.
> >>
> >>If there are no suggestions for improvements (or requests to see
> >>the changes as they appear in the module_write_netcdf.F files),
> >>I'll commit these changes. I can also update all of the existing
> >>grid.nc files available from
> >>http://www.mmm.ucar.edu/people/duda/files/mpas/ as well as provide
> >>a utility program to update any other grid.nc file with the
> >>new attributes.
> >>
> >>Cheers,
> >>Michael
> >_______________________________________________
> >mpas-developers mailing list
> >mpas-developers at mailman.ucar.edu
> >http://mailman.ucar.edu/mailman/listinfo/mpas-developers
> 
-------------- next part --------------
72a73,75
>       logical :: on_a_sphere
>       real (kind=RKIND) :: sphere_radius
> 
-------------- next part --------------
51a52,54
> 
>       character (len=16) :: c_on_a_sphere
>       real (kind=RKIND) :: r_sphere_radius
731a735,746
>       !
>       ! Read attributes
>       !
>       call io_input_get_att_text(input_obj, 'on_a_sphere', c_on_a_sphere)
>       call io_input_get_att_real(input_obj, 'sphere_radius', r_sphere_radius)
>       if (index(c_on_a_sphere, 'YES') /= 0) then
>          domain % blocklist % mesh % on_a_sphere = .true.
>       else
>          domain % blocklist % mesh % on_a_sphere = .false.
>       end if
>       domain % blocklist % mesh % sphere_radius = r_sphere_radius
> 
1102a1118,1169
>    
>    subroutine io_input_get_att_real(input_obj, attname, attvalue)
>       
>       implicit none
> 
>       type (io_input_object), intent(in) :: input_obj
>       character (len=*), intent(in) :: attname
>       real (kind=RKIND), intent(out) :: attvalue
> 
>       include 'netcdf.inc'
> 
>       integer :: nferr
> 
>       if (RKIND == 8) then
>          nferr = nf_get_att_double(input_obj % rd_ncid, NF_GLOBAL, attname, attvalue)
>       else
>          nferr = nf_get_att_real(input_obj % rd_ncid, NF_GLOBAL, attname, attvalue)
>       end if
>       if (nferr /= NF_NOERR) then
>          write(0,*) 'Warning: Attribute '//trim(attname)//' not found in '//trim(input_obj % filename)
>          if (index(attname, 'sphere_radius') /= 0) then
>             write(0,*) '   Setting '//trim(attname)//' to 1.0'
>             attvalue = 1.0
>          end if
>       end if
> 
>    end subroutine io_input_get_att_real
> 
>    
>    subroutine io_input_get_att_text(input_obj, attname, attvalue)
>       
>       implicit none
> 
>       type (io_input_object), intent(in) :: input_obj
>       character (len=*), intent(in) :: attname
>       character (len=*), intent(out) :: attvalue
> 
>       include 'netcdf.inc'
> 
>       integer :: nferr
> 
>       nferr = nf_get_att_text(input_obj % rd_ncid, NF_GLOBAL, attname, attvalue)
>       if (nferr /= NF_NOERR) then
>          write(0,*) 'Warning: Attribute '//trim(attname)//' not found in '//trim(input_obj % filename)
>          if (index(attname, 'on_a_sphere') /= 0) then
>             write(0,*) '   Setting '//trim(attname)//' to ''YES'''
>             attvalue = 'YES'
>          end if
>       end if
> 
>    end subroutine io_input_get_att_text
> 
-------------- next part --------------
92a93
>                           block_ptr % mesh, &
328a330
>                               mesh, &
337a340
>       type (grid_meta), intent(in) :: mesh
350a354,364
> 
>       if (mesh % on_a_sphere) then
>          nferr = nf_put_att_text(output_obj % wr_ncid, NF_GLOBAL, 'on_a_sphere', 16, 'YES             ')
>       else
>          nferr = nf_put_att_text(output_obj % wr_ncid, NF_GLOBAL, 'on_a_sphere', 16, 'NO              ')
>       end if
>       if (RKIND == 8) then
>          nferr = nf_put_att_double(output_obj % wr_ncid, NF_GLOBAL, 'sphere_radius', NF_DOUBLE, 1, mesh % sphere_radius)
>       else
>          nferr = nf_put_att_real(output_obj % wr_ncid, NF_GLOBAL, 'sphere_radius', NF_FLOAT, 1, mesh % sphere_radius)
>       end if


More information about the mpas-developers mailing list