diff --git a/sorc/grid_tools.fd/filter_topo.fd/filter_topo.F90 b/sorc/grid_tools.fd/filter_topo.fd/filter_topo.F90
index ff21ffc3e..bbfc3946f 100644
--- a/sorc/grid_tools.fd/filter_topo.fd/filter_topo.F90
+++ b/sorc/grid_tools.fd/filter_topo.fd/filter_topo.F90
@@ -1,4 +1,11 @@
 !> @file
+!! @brief This program does ???
+!! @author R. J. Purser
+
+!> This program does ???
+!!
+!! @return 0 for success, error code otherwise.
+!! @author GFDL Programmer
 program filter_topo
 
 
@@ -78,7 +85,13 @@ program filter_topo
 
 contains
 
-  !#####################################################################
+  !> ???
+  !!
+  !! @param[in] q1
+  !! @param[in] q2
+  !! @param[in] radius
+  !! @return ???
+  !! @author GFDL Programmer
   real function great_circle_dist( q1, q2, radius )
     real, intent(IN)           :: q1(2), q2(2)
     real, intent(IN), optional :: radius
@@ -103,17 +116,25 @@ real function great_circle_dist( q1, q2, radius )
 
   end function great_circle_dist
 
-  !####################################################################
+  !> ???
+  !!
+  !! <pre>
+  !!           p3
+  !!         /
+  !!        /
+  !!       p1 ---> angle
+  !!         \
+  !!          \
+  !!           p2
+  !! </pre>
+  !!
+  !! @param[in] p1 ???
+  !! @param[in] p2 ???
+  !! @param[in] p3 ???
+  !! @return ???
+  !! @author GFDL Programmer
   real function spherical_angle(p1, p2, p3)
 
-    !           p3
-    !         /
-    !        /
-    !       p1 ---> angle
-    !         \
-    !          \
-    !           p2
-
     real p1(3), p2(3), p3(3)
 
     real (f_p):: e1(3), e2(3), e3(3)
@@ -164,7 +185,15 @@ real function spherical_angle(p1, p2, p3)
   end function spherical_angle
 
 
-  !####################################################################
+  !> ???
+  !!
+  !! @param[in] p1 ???
+  !! @param[in] p4 ???
+  !! @param[in] p2 ???
+  !! @param[in] p3 ???
+  !! @param[in] radius ???
+  !! @return area
+  !! @author GFDL Programmer
   real function get_area(p1, p4, p2, p3, radius)
     !-----------------------------------------------
     real, intent(in), dimension(2):: p1, p2, p3, p4
@@ -208,6 +237,16 @@ real function get_area(p1, p4, p2, p3, radius)
 
   end function get_area
 
+  !> ???
+  !!
+  !! @param[in] q ???
+  !! @param[in] ng ???
+  !! @param[in] npx ???
+  !! @param[in] npy ???
+  !! @param[in] isd ???
+  !! @param[in] jsd ???
+  !! @param[in] fill ???
+  !! @author GFDL Programmer
   subroutine fill_AGRID_scalar_corners(q, ng, npx, npy, isd, jsd, fill)
     integer, intent(in)  :: ng, npx, npy, isd, jsd
     integer, intent(in)  :: fill
@@ -248,8 +287,16 @@ subroutine fill_AGRID_scalar_corners(q, ng, npx, npy, isd, jsd, fill)
 
   end subroutine fill_AGRID_scalar_corners
 
-
-  !####################################################################
+  !> ???
+  !!
+  !! @param[in] q ???
+  !! @param[in] ng ???
+  !! @param[in] npx ???
+  !! @param[in] npy ???
+  !! @param[in] isd ???
+  !! @param[in] jsd ???
+  !! @param[in] fill ???
+  !! @author GFDL Programmer
   subroutine fill_BGRID_scalar_corners(q, ng, npx, npy, isd, jsd, fill)
     integer, intent(in)  :: ng, npx, npy, isd, jsd
     integer, intent(in)  :: fill
@@ -291,7 +338,16 @@ subroutine fill_BGRID_scalar_corners(q, ng, npx, npy, isd, jsd, fill)
 
   end subroutine fill_BGRID_scalar_corners
 
-  !####################################################################
+  !> ???
+  !!
+  !! @param[in] x ???
+  !! @param[in] y ???
+  !! @param[in] ng ???
+  !! @param[in] npx ???
+  !! @param[in] npy ???
+  !! @param[in] isd ???
+  !! @param[in] jsd ???
+  !! @author GFDL Programmer
   subroutine fill_AGRID_xy_corners(x, y, ng, npx, npy, isd, jsd)
     integer, intent(in)                      :: ng, npx, npy, isd, jsd
     real, DIMENSION(isd:,jsd:,:), intent(INOUT):: x
@@ -314,9 +370,16 @@ subroutine fill_AGRID_xy_corners(x, y, ng, npx, npy, isd, jsd)
 
   end subroutine fill_AGRID_xy_corners
 
-
-
-  !####################################################################
+  !> ???
+  !!
+  !! @param[in] x ???
+  !! @param[in] y ???
+  !! @param[in] ng ???
+  !! @param[in] npx ???
+  !! @param[in] npy ???
+  !! @param[in] isd ???
+  !! @param[in] jsd ???
+  !! @author GFDL Programmer
   subroutine fill_DGRID_xy_corners(x, y, ng, npx, npy, isd, jsd)
     integer, intent(in)                      :: ng, npx, npy, isd, jsd
     real, DIMENSION(isd:,jsd:,:), intent(INOUT):: x
@@ -338,7 +401,12 @@ subroutine fill_DGRID_xy_corners(x, y, ng, npx, npy, isd, jsd)
 
   end subroutine fill_DGRID_xy_corners
 
-  !###############################################################
+  !> ???
+  !!
+  !! @param[in] p1 ???
+  !! @param[in] p2 ???
+  !! @param[in] pm ???
+  !! @author GFDL Programmer
   subroutine mid_pt_sphere(p1, p2, pm)
     real, intent(IN)  :: p1(2), p2(2)
     real, intent(OUT) :: pm(2)
@@ -352,8 +420,12 @@ subroutine mid_pt_sphere(p1, p2, pm)
 
   end subroutine mid_pt_sphere
 
-
-  !#####################################################################
+  !> ???
+  !!
+  !! @param[in] p1 ???
+  !! @param[in] p2 ???
+  !! @param[in] e ???
+  !! @author GFDL Programmer
   subroutine mid_pt3_cart(p1, p2, e)
     real, intent(IN)  :: p1(3), p2(3)
     real, intent(OUT) :: e(3)
@@ -382,6 +454,11 @@ subroutine mid_pt3_cart(p1, p2, e)
 
   end subroutine mid_pt3_cart
 
+  !> ???
+  !!
+  !! @param[in] p ???
+  !! @param[in] e ???
+  !! @author GFDL Programmer
   subroutine latlon2xyz(p, e)
     !
     ! Routine to map (lon, lat) to (x,y,z)
@@ -409,8 +486,13 @@ subroutine latlon2xyz(p, e)
 
   end subroutine latlon2xyz
 
-
-
+  !> ???
+  !!
+  !! @param[in] np  ???
+  !! @param[in] q ???
+  !! @param[in] xs ???
+  !! @param[in] ys ???
+  !! @author GFDL Programmer
   subroutine cart_to_latlon(np, q, xs, ys)
     ! vector version of cart_to_latlon1
     integer, intent(in):: np
@@ -451,7 +533,13 @@ subroutine cart_to_latlon(np, q, xs, ys)
 
   end  subroutine cart_to_latlon
 
-  !#####################################################################
+  !> ???
+  !!
+  !! @param[in] p1 ???
+  !! @param[in] p2 ???
+  !! @param[in] p3 ???
+  !! @return ???
+  !! @author GFDL Programmer
   real function cos_angle(p1, p2, p3)
     ! As spherical_angle, but returns the cos(angle)
     !       p3
@@ -498,8 +586,14 @@ real function cos_angle(p1, p2, p3)
 
   end function cos_angle
 
-
-  !#####################################################################
+  !> ???
+  !!
+  !! @param[in] q1 ???
+  !! @param[in] q2 ???
+  !! @param[in] q3 ???
+  !! @param[in] q4 ???
+  !! @param[in] e2 ???
+  !! @author GFDL Programmer
   subroutine cell_center2(q1, q2, q3, q4, e2)
     real, intent(in ) :: q1(2), q2(2), q3(2), q4(2)
     real, intent(out) :: e2(2)
@@ -527,8 +621,10 @@ subroutine cell_center2(q1, q2, q3, q4, e2)
 
   end subroutine cell_center2
 
-
-  !#####################################################################
+  !> ???
+  !!
+  !! @param[in] regional ???
+  !! @author GFDL Programmer
   subroutine read_grid_file(regional)
 
     logical, intent(in) :: regional   ! Is this a regional run?
@@ -875,8 +971,10 @@ subroutine read_grid_file(regional)
 
   end subroutine read_grid_file
 
-
-  !#####################################################################
+  !> ???
+  !!
+  !! @param[in] regional ???
+  !! @author GFDL Programmer
   subroutine read_topo_file(regional)
 
     logical,intent(in) :: regional   ! Is this a run with a regional domain?
@@ -961,8 +1059,16 @@ subroutine read_topo_file(regional)
 
   end subroutine read_topo_file
 
-  !##############################################################################
-  !--- replace the topo_field
+  !> Replace the topo_field.
+  !!
+  !! @param[in] is ???
+  !! @param[in] ie ???
+  !! @param[in] js ???
+  !! @param[in] je ???
+  !! @param[in] ntiles ???
+  !! @param[in] q ???
+  !! @param[in] regional ???
+  !! @author GFDL Programmer
   subroutine write_topo_file(is,ie,js,je,ntiles,q,regional)
      integer, intent(in) :: is,ie,js,je,ntiles                 
      real,    intent(in) :: q(is:ie,js:je,ntiles) 
@@ -1000,9 +1106,17 @@ subroutine write_topo_file(is,ie,js,je,ntiles,q,regional)
 
   end subroutine write_topo_file
 
-  !##############################################################################
-  ! this routine fill the halo points for the cubic grid. ioff and joff is used to distinguish
-  ! T, C, E, or N-cell
+  !> This routine fill the halo points for the cubic grid. ioff and
+  !! joff is used to distinguish T, C, E, or N-cell.
+  !!
+  !! @param[in] data  ???
+  !! @param[in] data2 ???
+  !! @param[in] halo ???
+  !! @param[in] ioff ???
+  !! @param[in] joff ???
+  !! @param[in] sign1 ???
+  !! @param[in] sign2 ???
+  !! @author GFDL Programmer
   subroutine fill_cubic_grid_halo(data, data2, halo, ioff, joff, sign1, sign2)
     integer, intent(in)                               :: halo
     real, dimension(1-halo:,1-halo:,:), intent(inout) :: data, data2
@@ -1044,7 +1158,34 @@ subroutine fill_cubic_grid_halo(data, data2, halo, ioff, joff, sign1, sign2)
 
   end subroutine fill_cubic_grid_halo
 
-  !#####################################################################
+  !> ???
+  !!
+  !! @param[in] is ???
+  !! @param[in] ie ???
+  !! @param[in] js ???
+  !! @param[in] je ???
+  !! @param[in] isd ???
+  !! @param[in] ied ???
+  !! @param[in] jsd ???
+  !! @param[in] jed ???
+  !! @param[in] npx ???
+  !! @param[in] npy ???
+  !! @param[in] npx_global ???
+  !! @param[in] ntiles ???
+  !! @param[in] grid_type ???
+  !! @param[in] stretch_fac ???
+  !! @param[in] nested ???
+  !! @param[in] area ???
+  !! @param[in] dxa ???
+  !! @param[in] dya ???
+  !! @param[in] dx ???
+  !! @param[in] dy ???
+  !! @param[in] dxc ???
+  !! @param[in] dyc ???
+  !! @param[in] sin_sg ???
+  !! @param[in] phis ???
+  !! @param[in] regional ???
+  !! @author GFDL Programmer
   subroutine FV3_zs_filter (is, ie, js, je, isd, ied, jsd, jed, npx, npy, npx_global, ntiles,  &
        grid_type, stretch_fac, nested, area, dxa, dya, dx, dy, dxc, dyc, &
        sin_sg,  phis, regional )
@@ -1096,7 +1237,38 @@ subroutine FV3_zs_filter (is, ie, js, je, isd, ied, jsd, jed, npx, npy, npx_glob
 
   end subroutine FV3_zs_filter
 
-  !#####################################################################
+  !> ???
+  !!
+  !! @param[in] is ???
+  !! @param[in] ie ???
+  !! @param[in] js ???
+  !! @param[in] je ???
+  !! @param[in] isd ???
+  !! @param[in] ied ???
+  !! @param[in] jsd ???
+  !! @param[in] jed ???
+  !! @param[in] npx ???
+  !! @param[in] npy ???
+  !! @param[in] ntiles ???
+  !! @param[in] q ???
+  !! @param[in] area ???
+  !! @param[in] dx ???
+  !! @param[in] dy ???
+  !! @param[in] dxa ???
+  !! @param[in] dya ???
+  !! @param[in] dxc ???
+  !! @param[in] dyc ???
+  !! @param[in] sin_sg ???
+  !! @param[in] cd ???
+  !! @param[in] zero_ocean ???
+  !! @param[in] check_slope ???
+  !! @param[in] filter_type ???
+  !! @param[in] grid_type ???
+  !! @param[in] mask ???
+  !! @param[in] nested ???
+  !! @param[in] ntmax ???
+  !! @param[in] regional ???
+  !! @author GFDL Programmer
   subroutine two_delta_filter(is, ie, js, je, isd, ied, jsd, jed, npx, npy, ntiles, &
        q, area, dx, dy, dxa, dya, dxc, dyc, sin_sg, cd, zero_ocean,  &
         check_slope, filter_type, grid_type, mask, nested, ntmax, regional)
@@ -1376,7 +1548,33 @@ subroutine two_delta_filter(is, ie, js, je, isd, ied, jsd, jed, npx, npy, ntiles
 
   end subroutine two_delta_filter
 
-  !#####################################################################
+  !> ???
+  !!
+  !! @param[in] is ???
+  !! @param[in] ie ???
+  !! @param[in] js ???
+  !! @param[in] je ???
+  !! @param[in] isd ???
+  !! @param[in] ied ???
+  !! @param[in] jsd ???
+  !! @param[in] jed ???
+  !! @param[in] npx ???
+  !! @param[in] npy ???
+  !! @param[in] ntiles ???
+  !! @param[in] q ???
+  !! @param[in] area ???
+  !! @param[in] dx ???
+  !! @param[in] dy ???
+  !! @param[in] dxc ???
+  !! @param[in] dyc ???
+  !! @param[in] sin_sg ???
+  !! @param[in] nmax ???
+  !! @param[in] cd ???
+  !! @param[in] zero_ocean ???
+  !! @param[in] mask ???
+  !! @param[in] nested ???
+  !! @param[in] regional ???
+  !! @author GFDL Programmer
   subroutine del2_cubed_sphere(is, ie, js, je, isd, ied, jsd, jed, npx, npy, ntiles,&
           q, area, dx, dy, dxc, dyc, sin_sg, nmax, cd, zero_ocean, mask, nested, regional)
     integer, intent(in) :: is,  ie,  js,  je
@@ -1474,7 +1672,32 @@ subroutine del2_cubed_sphere(is, ie, js, je, isd, ied, jsd, jed, npx, npy, ntile
 
   end subroutine del2_cubed_sphere
 
-  !#####################################################################
+  !> ???
+  !!
+  !! @param[in] is ???
+  !! @param[in] ie ???
+  !! @param[in] js ???
+  !! @param[in] je ???
+  !! @param[in] isd ???
+  !! @param[in] ied ???
+  !! @param[in] jsd ???
+  !! @param[in] jed ???
+  !! @param[in] npx ???
+  !! @param[in] npy ???
+  !! @param[in] ntiles ???
+  !! @param[in] q ???
+  !! @param[in] area ???
+  !! @param[in] dx ???
+  !! @param[in] dy ???
+  !! @param[in] dxc ???
+  !! @param[in] dyc ???
+  !! @param[in] sin_sg ???
+  !! @param[in] nmax ???
+  !! @param[in] zero_ocean ???
+  !! @param[in] mask ???
+  !! @param[in] nested ???
+  !! @param[in] regional ???
+  !! @author GFDL Programmer
   subroutine del4_cubed_sphere(is, ie, js, je, isd, ied, jsd, jed, npx, npy, ntiles, &
        q, area, dx, dy, dxc, dyc, sin_sg, nmax, zero_ocean, mask, nested, regional)
     integer, intent(in) :: is,  ie,  js,  je
@@ -1709,7 +1932,11 @@ subroutine del4_cubed_sphere(is, ie, js, je, isd, ied, jsd, jed, npx, npy, ntile
 
   end subroutine del4_cubed_sphere
 
-  !#####################################################################
+  !> ???
+  !!
+  !! @param[in] status ???
+  !! @param[in] string ???
+  !! @author GFDL Programmer
   subroutine handle_err(status, string)
     integer,          intent(in) :: status
     character(len=*), intent(in) :: string
@@ -1724,10 +1951,10 @@ subroutine handle_err(status, string)
 
   end subroutine  handle_err
 
-
-  !#######################################################################
-  ! reads the namelist file, write namelist to log file,
-
+  !> Reads the namelist file, write namelist to log file.
+  !!
+  !! 
+  !! @author GFDL Programmer
   subroutine read_namelist
 
     !  read namelist
@@ -1750,6 +1977,10 @@ subroutine read_namelist
 
   end subroutine read_namelist
 
+  !> Check results of netCDF call.
+  !!
+  !! @param[in] status return code to check
+  !! @author GFDL Programmer
   subroutine check(status)
   use netcdf
   integer,intent(in) :: status
@@ -1762,9 +1993,10 @@ subroutine check(status)
   endif
   end subroutine check
 
-  !#######################################################################
-  ! compute resolution-dependent values for the filtering.
-
+  !> Compute resolution-dependent values for the filtering.
+  !!
+  !! 
+  !! @author GFDL Programmer
   subroutine compute_filter_constants
 
   ! set the given values for various cube resolutions (c48, c96, c192, c384, c768, c1152, c3072)
diff --git a/sorc/grid_tools.fd/global_equiv_resol.fd/global_equiv_resol.f90 b/sorc/grid_tools.fd/global_equiv_resol.fd/global_equiv_resol.f90
index 08ad8a1d4..6a94e852d 100644
--- a/sorc/grid_tools.fd/global_equiv_resol.fd/global_equiv_resol.f90
+++ b/sorc/grid_tools.fd/global_equiv_resol.fd/global_equiv_resol.f90
@@ -1,7 +1,19 @@
 !> @file
-!=======================================================================
+!! @brief Compute the global equivalent resolution for regional grids.
+!! @author Gerard Ketefian NOAA/GSL
+
+!> Compute the global equivalent resolution for regional grids
+!! using the average model grid cell size in meters squared.
+!!
+!! - RES(in km) = (360 degrees / 4*CRES) * 111 km
+!!
+!! Using the average cell size the equivalent global resolution is:
+!!
+!! - CRES = NINT((2*PI*RADIUS_OF_EARTH)/(4*AVG_CELL_SIZE)
+!!
+!! @author Gerard Ketefian NOAA/GSL
+!! @return 0 for success, error code otherwise
 program global_equiv_resol
-!=======================================================================
 
   use netcdf
 
@@ -173,7 +185,10 @@ program global_equiv_resol
 
 end program global_equiv_resol
 
-
+!> Check results of netCDF call.
+!!
+!! @param[in] status return code to check
+!! @author Gerard Ketefian NOAA/GSL
 subroutine check(status)
   use netcdf
   integer,intent(in) :: status
diff --git a/sorc/grid_tools.fd/regional_esg_grid.fd/pmat4.f90 b/sorc/grid_tools.fd/regional_esg_grid.fd/pmat4.f90
index 3cd11cb43..b444f9455 100644
--- a/sorc/grid_tools.fd/regional_esg_grid.fd/pmat4.f90
+++ b/sorc/grid_tools.fd/regional_esg_grid.fd/pmat4.f90
@@ -1657,8 +1657,8 @@ end subroutine setem
 !!
 !! @param[in] a real quaternion
 !! @param[in] b real quaternion
+!! @return c result
 !! @author R. J. Purser
-!! @return c
 function mulqq(a,b)result(c)!                                         [mulqq]
 implicit none
 real(dp),dimension(0:3),intent(IN ):: a,b
diff --git a/sorc/grid_tools.fd/regional_esg_grid.fd/regional_esg_grid.f90 b/sorc/grid_tools.fd/regional_esg_grid.fd/regional_esg_grid.f90
index b66fe5c1b..264193d5c 100644
--- a/sorc/grid_tools.fd/regional_esg_grid.fd/regional_esg_grid.f90
+++ b/sorc/grid_tools.fd/regional_esg_grid.fd/regional_esg_grid.f90
@@ -1,9 +1,18 @@
 !> @file
-!! @brief ???
+!! @brief Compute geo-referencing parameters for the Extended
+!! Schmidt Gnomonic (ESG) regional grid.
 !! @author R. J. Purser
 
-!! Program which does ???
-!!
+!> Driver routine to compute geo-referencing parameters for
+!! the Extended Schmidt Gnomonic (ESG) regional grid.
+!! The parameters are:
+!! - Geographic longitude (degrees)
+!! - Geographic latitude (degrees)
+!! - Grid edge 'x' distance (meters)
+!! - Grid edge 'y' distance (meters)
+!! - Area (meters squared)
+!! - Grid vertex 'x' angle with respect to geographic east (degrees)
+!! - Grid vertex 'y' angle with respect to geographic north (degrees)
 !! @author R. J. Purser
 !! @return 0 for success, error code otherwise.
 program regional_grid
@@ -174,6 +183,10 @@ program regional_grid
 
 end program regional_grid
 
+!> Check results of netCDF call.
+!!
+!! @param[in] status return code to check
+!! @author R. J. Purser
 subroutine check(status)
 use netcdf
 integer,intent(in) :: status