grid_tools  1.13.0
 All Data Structures Files Functions Variables Pages
global_equiv_resol.f90
Go to the documentation of this file.
1 
4 
16 program global_equiv_resol
17 
18  use netcdf
19 
20  implicit none
21 
22  integer, parameter :: dp = kind(1.0d0)
23  real(dp), parameter :: pi_geom = 4.0*atan(1.0), &
24  radius_earth = 6371200.0
25 
26  character(len=256) :: grid_fn
27  integer :: ncid, nxsg_dimid, nysg_dimid, dasg_varid, num_args
28  integer :: nxsg, nysg, nx, ny, res_equiv
29  real(dp) :: avg_cell_size, min_cell_size, max_cell_size
30  real(dp), dimension(:,:), allocatable :: &
31  quarter_da_ll, quarter_da_lr, quarter_da_ur, quarter_da_ul, &
32  dasg, da, sqrt_da
33 !
34 !=======================================================================
35 !
36 ! Read in the name of the file from the command line. The command-line
37 ! call to this program should have exactly one argument consisting of
38 ! the path to the NetCDF grid specification file to be read in. If this
39 ! is not the case, print out a usage message and exit.
40 !
41 !=======================================================================
42 !
43  num_args = command_argument_count()
44  if (num_args == 1) then
45  call get_command_argument(1, grid_fn)
46  else
47  WRITE(*,500)
48  WRITE(*,500) "Exactly one argument must be specified to program global_equiv_resol."
49  WRITE(*,500) "Usage:"
50  WRITE(*,500)
51  WRITE(*,500) " global_equiv_resol path_to_grid_file"
52  WRITE(*,500)
53  WRITE(*,500) "where path_to_grid_file is the path to the NetCDF grid file. Actual "
54  WRITE(*,500) "number of specified command line arguments is:"
55  WRITE(*,510) " num_args = ", num_args
56  WRITE(*,500) "Stopping."
57 500 FORMAT(a)
58 510 FORMAT(a, i3)
59  stop
60  end if
61 !
62 !=======================================================================
63 !
64 ! Open the grid file and read in the dimensions of the supergrid. The
65 ! supergrid is a grid that has twice the resolution of the actual/compu-
66 ! tational grid. In the file, the names of the supergrid dimensions are
67 ! nx and ny. Here, however, we reserve those names for the dimensions
68 ! of the actual grid (since in the FV3 code and in other data files, nx
69 ! and ny are used to denote the dimensions of the actual grid) and in-
70 ! stead use the variables nxSG and nySG to denote the dimensions of the
71 ! supergrid.
72 !
73 !=======================================================================
74 !
75  WRITE(*,500)
76  WRITE(*,500) "Opening NetCDF grid file for reading/writing:"
77  WRITE(*,500) " grid_fn = " // trim(grid_fn)
78 
79  call check( nf90_open(trim(grid_fn), nf90_write, ncid) )
80 
81  call check( nf90_inq_dimid(ncid, "nx", nxsg_dimid) )
82  call check( nf90_inquire_dimension(ncid, nxsg_dimid, len=nxsg) )
83 
84  call check( nf90_inq_dimid(ncid, "ny", nysg_dimid) )
85  call check( nf90_inquire_dimension(ncid, nysg_dimid, len=nysg) )
86 
87  WRITE(*,500)
88  WRITE(*,500) "Dimensions of supergrid are:"
89  WRITE(*,520) " nxSG = ", nxsg
90  WRITE(*,520) " nySG = ", nysg
91 520 FORMAT(a, i7)
92 !
93 !=======================================================================
94 !
95 ! Read in the cell areas on the supergrid. Then add the areas of the
96 ! four supergrid cells that make up one grid cell to obtain the cell
97 ! areas on the actual grid.
98 !
99 !=======================================================================
100 !
101  allocate(dasg(0:nxsg-1, 0:nysg-1))
102  call check( nf90_inq_varid(ncid, "area", dasg_varid) )
103  call check( nf90_get_var(ncid, dasg_varid, dasg) )
104 
105  nx = nxsg/2
106  ny = nysg/2
107 
108  WRITE(*,500)
109  WRITE(*,500) "Dimensions of (actual, i.e. computational) grid are:"
110  WRITE(*,520) " nx = ", nx
111  WRITE(*,520) " ny = ", ny
112 
113  allocate(quarter_da_ll(0:nx-1, 0:ny-1))
114  allocate(quarter_da_lr(0:nx-1, 0:ny-1))
115  allocate(quarter_da_ul(0:nx-1, 0:ny-1))
116  allocate(quarter_da_ur(0:nx-1, 0:ny-1))
117 
118  quarter_da_ll = dasg(0:nxsg-1:2, 0:nysg-1:2)
119  quarter_da_lr = dasg(0:nxsg-1:2, 1:nysg-1:2)
120  quarter_da_ur = dasg(1:nxsg-1:2, 1:nysg-1:2)
121  quarter_da_ul = dasg(1:nxsg-1:2, 0:nysg-1:2)
122 
123  allocate(da(0:nx-1, 0:ny-1))
124  allocate(sqrt_da(0:nx-1, 0:ny-1))
125 
126  da = quarter_da_ll + quarter_da_lr + quarter_da_ur + quarter_da_ul
127 !
128 !=======================================================================
129 !
130 ! Calculate a typical/representative cell size for each cell by taking
131 ! the square root of the area of the cell. Then calculate the minimum,
132 ! maximum, and average cell sizes over the whole grid.
133 !
134 !=======================================================================
135 !
136  sqrt_da = sqrt(da)
137  min_cell_size = minval(sqrt_da)
138  max_cell_size = maxval(sqrt_da)
139  avg_cell_size = sum(sqrt_da)/(nx*ny)
140 
141  WRITE(*,500)
142  WRITE(*,500) "Minimum, maximum, and average cell sizes are (based on square"
143  WRITE(*,500) "root of cell area):"
144  WRITE(*,530) " min_cell_size = ", min_cell_size
145  WRITE(*,530) " max_cell_size = ", max_cell_size
146  WRITE(*,530) " avg_cell_size = ", avg_cell_size
147 530 FORMAT(a, g11.4)
148 !
149 !=======================================================================
150 !
151 ! Use the average cell size to calculate an equivalent global uniform
152 ! cubed-sphere resolution (in units of number of cells) for the regional
153 ! grid. This is the RES that a global uniform (i.e. stretch factor of
154 ! 1) cubed-sphere grid would need to have in order to have the same no-
155 ! minal cell size as the average cell size of the regional grid.
156 !
157 !=======================================================================
158 !
159  res_equiv = nint( (2.0*pi_geom*radius_earth)/(4.0*avg_cell_size) )
160 
161  WRITE(*,500)
162  WRITE(*,500) "Equivalent global uniform cubed-sphere resolution is:"
163  WRITE(*,*) " RES_equiv = ", res_equiv
164 !
165 !=======================================================================
166 !
167 ! Write the average cell size and equivalent global resolution to the
168 ! grid file as a global attributes.
169 !
170 !=======================================================================
171 !
172  WRITE(*,500)
173  WRITE(*,500) "Writing avg_cell_size and RES_equiv to the grid specification"
174  WRITE(*,500) "file as global attributes..."
175 
176  call check( nf90_redef(ncid) )
177  call check( nf90_put_att(ncid, nf90_global, "avg_cell_size", avg_cell_size) )
178  call check( nf90_put_att(ncid, nf90_global, "RES_equiv", res_equiv) )
179  call check( nf90_enddef(ncid) )
180 
181  call check( nf90_close(ncid) )
182 
183  WRITE(*,500)
184  WRITE(*,500) "Done."
185 
186 end program global_equiv_resol
187 
192 subroutine check(status)
193  use netcdf
194  integer,intent(in) :: status
195 !
196  if(status /= nf90_noerr) then
197  write(0,*) ' check netcdf status = ', status
198  write(0,'("error ", a)') trim(nf90_strerror(status))
199  stop "Stopped"
200  endif
201 end subroutine check
global_equiv_resol
program global_equiv_resol
Compute the global equivalent resolution for regional grids using the average model grid cell size in...
Definition: global_equiv_resol.f90:16
check
subroutine check(status)
Check results of netCDF call.
Definition: filter_topo.F90:1958