forrest@0
|
1 |
module clm_mod
|
forrest@0
|
2 |
|
forrest@0
|
3 |
use netcdf
|
forrest@0
|
4 |
use kind_mod
|
forrest@0
|
5 |
implicit none
|
forrest@0
|
6 |
|
forrest@0
|
7 |
integer, parameter :: namelen = 128 ! string length for paths
|
forrest@0
|
8 |
integer, parameter :: varnamelen = 32 ! string length for varnames
|
forrest@0
|
9 |
integer, parameter :: posnamelen = 4 ! string length for positive attribute
|
forrest@0
|
10 |
integer, parameter :: varattlen = 128 ! string length for variable attributes
|
forrest@0
|
11 |
character (len=*), parameter :: xname = 'lon'
|
forrest@0
|
12 |
character (len=*), parameter :: yname = 'lat'
|
forrest@0
|
13 |
character (len=*), parameter :: zsoi_name = 'levsoi'
|
forrest@0
|
14 |
character (len=*), parameter :: zlak_name = 'levlak'
|
forrest@0
|
15 |
character (len=*), parameter :: tname = 'time'
|
forrest@0
|
16 |
character (len=*), parameter :: hr_name = 'hour'
|
forrest@0
|
17 |
character (len=*), parameter :: nsamples_name = 'nsamples'
|
forrest@0
|
18 |
character (len=*), parameter :: hiname = 'hist_interval'
|
forrest@0
|
19 |
character (len=*), parameter :: tbname = 'time_bounds'
|
forrest@0
|
20 |
character (len=*), parameter :: aname = 'area'
|
forrest@0
|
21 |
character (len=*), parameter :: lname = 'landfrac'
|
forrest@0
|
22 |
character (len=*), parameter :: mname = 'landmask'
|
forrest@0
|
23 |
character (len=*), parameter :: oname = 'topo'
|
forrest@0
|
24 |
character (len=*), parameter :: long_name_name = 'long_name'
|
forrest@0
|
25 |
character (len=*), parameter :: missing_value_name = 'missing_value'
|
forrest@0
|
26 |
character (len=*), parameter :: units_name = 'units'
|
forrest@0
|
27 |
character (len=*), parameter :: case_id_name = 'case_id'
|
forrest@0
|
28 |
character (len=*), parameter :: inidat_name = 'Initial_conditions_dataset'
|
forrest@0
|
29 |
character (len=*), parameter :: surdat_name = 'Surface_dataset'
|
forrest@0
|
30 |
character (len=*), parameter :: pftdat_name = 'PFT_physiological_constants_dataset'
|
forrest@0
|
31 |
character (len=*), parameter :: rtmdat_name = 'RTM_input_datset'
|
forrest@0
|
32 |
|
forrest@0
|
33 |
type surface
|
forrest@0
|
34 |
integer :: xsize ! number of longitude points (lon)
|
forrest@0
|
35 |
integer :: ysize ! number of latitutde points (lat)
|
forrest@0
|
36 |
integer :: psize ! number of pfts
|
forrest@0
|
37 |
real(r8), pointer :: x(:) ! longitude values (lon)
|
forrest@0
|
38 |
real(r8), pointer :: y(:) ! latitude values (lat)
|
forrest@0
|
39 |
real(r8), pointer :: p(:) ! pft number
|
forrest@0
|
40 |
real(r4), pointer :: pft(:,:,:) ! percent of each pft
|
forrest@0
|
41 |
end type surface
|
forrest@0
|
42 |
|
forrest@0
|
43 |
type coordinates
|
forrest@0
|
44 |
integer :: xsize ! number of longitude points (lon)
|
forrest@0
|
45 |
integer :: ysize ! number of latitutde points (lat)
|
forrest@0
|
46 |
integer :: zsoi_size ! number of soil z-levels (levsoi)
|
forrest@0
|
47 |
integer :: zlak_size ! number of lake z-levels (levlak)
|
forrest@0
|
48 |
integer :: tsize ! number of time points (time)
|
forrest@0
|
49 |
integer :: hr_size ! number of hour points (hour = 24)
|
forrest@0
|
50 |
integer :: nsamples ! number of time samples in post-processing statistics
|
forrest@0
|
51 |
real(r8), pointer :: x(:) ! longitude values (lon)
|
forrest@0
|
52 |
real(r8), pointer :: y(:) ! latitude values (lat)
|
forrest@0
|
53 |
real(r8), pointer :: zsoi(:) ! soil z-level values (levsoi)
|
forrest@0
|
54 |
real(r8), pointer :: zlak(:) ! lake z-level values (levlak)
|
forrest@0
|
55 |
real(r8), pointer :: hr(:) ! hour-of-day values (hour)
|
forrest@0
|
56 |
real(r8), pointer :: area(:,:) ! gridcell areas
|
forrest@0
|
57 |
real(r8), pointer :: landfrac(:,:) ! gridcell land fractions
|
forrest@0
|
58 |
integer, pointer :: landmask(:,:) ! gridcell land/ocean mask
|
forrest@0
|
59 |
real(r8), pointer :: orog(:,:) ! gridcell orography/topography
|
forrest@0
|
60 |
character (len=namelen) :: case_id ! case name of run
|
forrest@0
|
61 |
character (len=namelen) :: inidat ! initial dataset filename
|
forrest@0
|
62 |
character (len=namelen) :: surdat ! surface dataset filename
|
forrest@0
|
63 |
character (len=namelen) :: pftdat ! pft dataset filename
|
forrest@0
|
64 |
character (len=namelen) :: rtmdat ! runoff dataset filename
|
forrest@0
|
65 |
logical :: hour_flag ! flag for hour-of-day statistical summary files
|
forrest@0
|
66 |
end type coordinates
|
forrest@0
|
67 |
|
forrest@0
|
68 |
type variable_data
|
forrest@0
|
69 |
character (len=varnamelen) :: varname
|
forrest@0
|
70 |
integer :: ndims
|
forrest@0
|
71 |
integer, pointer :: int2d(:,:)
|
forrest@0
|
72 |
real(r4), pointer :: var2d(:,:)
|
forrest@0
|
73 |
real(r4), pointer :: var3d(:,:,:)
|
forrest@0
|
74 |
real(r4), pointer :: var4d(:,:,:,:)
|
forrest@0
|
75 |
real(r8), pointer :: time(:)
|
forrest@0
|
76 |
real(r8), pointer :: time_bounds(:,:)
|
forrest@0
|
77 |
logical :: int_type ! set if type is integer
|
forrest@0
|
78 |
logical :: soil_layer_flag
|
forrest@0
|
79 |
logical :: hour_flag ! flag for hour-of-day statistical summary variable
|
forrest@0
|
80 |
character (len=varattlen) :: long_name
|
forrest@0
|
81 |
real(r4) :: missing_value
|
forrest@0
|
82 |
character (len=varattlen) :: units
|
forrest@0
|
83 |
character (len=posnamelen) :: positive
|
forrest@0
|
84 |
character (len=varnamelen) :: out_varname
|
forrest@0
|
85 |
end type variable_data
|
forrest@0
|
86 |
|
forrest@0
|
87 |
type(surface), target :: surf
|
forrest@0
|
88 |
type(coordinates), target :: coord
|
forrest@0
|
89 |
type(variable_data), target :: var_data
|
forrest@0
|
90 |
save
|
forrest@0
|
91 |
|
forrest@0
|
92 |
public clm_read_surf
|
forrest@0
|
93 |
public clm_free_surf
|
forrest@0
|
94 |
public clm_read_coord
|
forrest@0
|
95 |
public clm_copy_grid_data
|
forrest@0
|
96 |
public clm_read_data
|
forrest@0
|
97 |
public clm_convert_data
|
forrest@0
|
98 |
public clm_free_data
|
forrest@0
|
99 |
private nc_err
|
forrest@0
|
100 |
|
forrest@0
|
101 |
contains
|
forrest@0
|
102 |
|
forrest@0
|
103 |
!
|
forrest@0
|
104 |
!--------------------------------------------------------------------------
|
forrest@0
|
105 |
!
|
forrest@0
|
106 |
subroutine clm_read_surf(fsurdat)
|
forrest@0
|
107 |
|
forrest@0
|
108 |
implicit none
|
forrest@0
|
109 |
character (len=namelen), intent(in) :: fsurdat ! input surface dataset
|
forrest@0
|
110 |
integer :: i
|
forrest@0
|
111 |
integer :: ierr
|
forrest@0
|
112 |
integer :: ncid
|
forrest@0
|
113 |
integer :: xdimid, xvarid, xlen
|
forrest@0
|
114 |
integer :: ydimid, yvarid, ylen
|
forrest@0
|
115 |
integer :: pdimid, pvarid, plen
|
forrest@0
|
116 |
real(r8), allocatable :: longxy(:,:)
|
forrest@0
|
117 |
real(r8), allocatable :: latixy(:,:)
|
forrest@0
|
118 |
real(r4), allocatable :: pctpft(:,:,:)
|
forrest@0
|
119 |
|
forrest@0
|
120 |
! Read coordinated from surface dataset
|
forrest@0
|
121 |
print *, 'Opening ',trim(fsurdat)
|
forrest@0
|
122 |
call nc_err(nf90_open(path=trim(fsurdat), mode=nf90_nowrite, ncid=ncid))
|
forrest@0
|
123 |
call nc_err(nf90_inq_dimid(ncid, 'lsmlon', xdimid))
|
forrest@0
|
124 |
call nc_err(nf90_inq_dimid(ncid, 'lsmlat', ydimid))
|
forrest@0
|
125 |
call nc_err(nf90_inq_dimid(ncid, 'lsmpft', pdimid))
|
forrest@0
|
126 |
|
forrest@0
|
127 |
call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=xdimid, len=xlen))
|
forrest@0
|
128 |
call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=ydimid, len=ylen))
|
forrest@0
|
129 |
call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=pdimid, len=plen))
|
forrest@0
|
130 |
|
forrest@0
|
131 |
call nc_err(nf90_inq_varid(ncid, 'LONGXY', xvarid))
|
forrest@0
|
132 |
call nc_err(nf90_inq_varid(ncid, 'LATIXY', yvarid))
|
forrest@0
|
133 |
call nc_err(nf90_inq_varid(ncid, 'PCT_PFT', pvarid))
|
forrest@0
|
134 |
|
forrest@0
|
135 |
allocate(longxy(xlen, ylen), latixy(xlen, ylen), pctpft(xlen, ylen, plen), stat=ierr)
|
forrest@0
|
136 |
if (ierr /= 0) then
|
forrest@0
|
137 |
print *, 'clm_read_surf: Cannot allocate longxy, latixy, pctpft variables'
|
forrest@0
|
138 |
stop
|
forrest@0
|
139 |
end if
|
forrest@0
|
140 |
|
forrest@0
|
141 |
call nc_err(nf90_get_var(ncid, xvarid, longxy))
|
forrest@0
|
142 |
call nc_err(nf90_get_var(ncid, yvarid, latixy))
|
forrest@0
|
143 |
call nc_err(nf90_get_var(ncid, pvarid, pctpft))
|
forrest@0
|
144 |
|
forrest@0
|
145 |
print *, 'Closing ',trim(fsurdat)
|
forrest@0
|
146 |
call nc_err(nf90_close(ncid))
|
forrest@0
|
147 |
|
forrest@0
|
148 |
! Allocate and fill in the surface derived data type
|
forrest@0
|
149 |
surf%xsize = xlen
|
forrest@0
|
150 |
surf%ysize = ylen
|
forrest@0
|
151 |
surf%psize = plen
|
forrest@0
|
152 |
allocate(surf%x(xlen), surf%y(ylen), surf%p(plen), surf%pft(xlen, ylen, plen), stat=ierr)
|
forrest@0
|
153 |
if (ierr /= 0) then
|
forrest@0
|
154 |
print *, 'clm_read_surf: Cannot allocate surface type variables'
|
forrest@0
|
155 |
stop
|
forrest@0
|
156 |
end if
|
forrest@0
|
157 |
surf%x(:) = longxy(:,1)
|
forrest@0
|
158 |
surf%y(:) = latixy(1,:)
|
forrest@0
|
159 |
do i = 1, plen
|
forrest@0
|
160 |
surf%p(i) = real(i)
|
forrest@0
|
161 |
end do
|
forrest@0
|
162 |
surf%pft = pctpft
|
forrest@0
|
163 |
|
forrest@0
|
164 |
deallocate(longxy, latixy, pctpft)
|
forrest@0
|
165 |
|
forrest@0
|
166 |
end subroutine clm_read_surf
|
forrest@0
|
167 |
!
|
forrest@0
|
168 |
!--------------------------------------------------------------------------
|
forrest@0
|
169 |
!
|
forrest@0
|
170 |
subroutine clm_free_surf()
|
forrest@0
|
171 |
|
forrest@0
|
172 |
implicit none
|
forrest@0
|
173 |
|
forrest@0
|
174 |
if (allocated(surf%x)) deallocate(surf%x)
|
forrest@0
|
175 |
if (allocated(surf%y)) deallocate(surf%y)
|
forrest@0
|
176 |
if (allocated(surf%p)) deallocate(surf%p)
|
forrest@0
|
177 |
if (allocated(surf%pft)) deallocate(surf%pft)
|
forrest@0
|
178 |
|
forrest@0
|
179 |
end subroutine clm_free_surf
|
forrest@0
|
180 |
!
|
forrest@0
|
181 |
!--------------------------------------------------------------------------
|
forrest@0
|
182 |
!
|
forrest@0
|
183 |
subroutine clm_read_coord(numf, input_path, fnames, read_grid_all)
|
forrest@0
|
184 |
|
forrest@0
|
185 |
implicit none
|
forrest@0
|
186 |
integer, intent(in) :: numf ! number of filenames
|
forrest@0
|
187 |
character (len=namelen), intent(in) :: input_path ! input path
|
forrest@0
|
188 |
character (len=namelen), intent(in) :: fnames(numf) ! input filenames
|
forrest@0
|
189 |
logical, intent(in) :: read_grid_all ! flag to read area,...
|
forrest@0
|
190 |
integer :: ncid
|
forrest@0
|
191 |
integer :: xdimid, xvarid, xlen
|
forrest@0
|
192 |
integer :: ydimid, yvarid, ylen
|
forrest@0
|
193 |
integer :: zsoi_dimid, zsoi_varid, zsoi_len
|
forrest@0
|
194 |
integer :: zlak_dimid, zlak_varid, zlak_len
|
forrest@0
|
195 |
integer :: tdimid, tvarid, tlen
|
forrest@0
|
196 |
integer :: hr_dimid, hr_varid, hr_len
|
forrest@0
|
197 |
integer :: nsamples_dimid, nsamples_len
|
forrest@0
|
198 |
integer :: avarid ! area (when read_grid_all true)
|
forrest@0
|
199 |
integer :: lvarid ! landfrac (when read_grid_all true)
|
forrest@0
|
200 |
integer :: mvarid ! landmask (when read_grid_all true)
|
forrest@0
|
201 |
integer :: ovarid ! orog/topo (when read_grid_all true)
|
forrest@0
|
202 |
integer(i8) :: i
|
forrest@0
|
203 |
integer(i8) :: time_size
|
forrest@0
|
204 |
real(r4), allocatable :: xvals(:)
|
forrest@0
|
205 |
real(r4), allocatable :: yvals(:)
|
forrest@0
|
206 |
real(r4), allocatable :: zsoi_vals(:)
|
forrest@0
|
207 |
real(r4), allocatable :: zlak_vals(:)
|
forrest@0
|
208 |
real(r4), allocatable :: hr_vals(:)
|
forrest@0
|
209 |
real(r4), allocatable :: avals(:,:) ! area (when read_grid_all true)
|
forrest@0
|
210 |
real(r4), allocatable :: lvals(:,:) ! landfrac (when read_grid_all true)
|
forrest@0
|
211 |
integer, allocatable :: mvals(:,:) ! landmask (when read_grid_all true)
|
forrest@0
|
212 |
real(r4), allocatable :: ovals(:,:) ! orog/topo (when read_grid_all true)
|
forrest@0
|
213 |
character (len=namelen) :: case_id = ''
|
forrest@0
|
214 |
character (len=namelen) :: inidat = ''
|
forrest@0
|
215 |
character (len=namelen) :: surdat = ''
|
forrest@0
|
216 |
character (len=namelen) :: pftdat = ''
|
forrest@0
|
217 |
character (len=namelen) :: rtmdat = ''
|
forrest@0
|
218 |
logical :: hrflag = .false.
|
forrest@0
|
219 |
integer :: ierr
|
forrest@0
|
220 |
|
forrest@0
|
221 |
! Read coordinates from first file
|
forrest@0
|
222 |
print *, 'Opening ',trim(input_path)//'/'//trim(fnames(1))
|
forrest@0
|
223 |
call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(1)), mode=nf90_nowrite, ncid=ncid))
|
forrest@0
|
224 |
call nc_err(nf90_inq_dimid(ncid, xname, xdimid))
|
forrest@0
|
225 |
call nc_err(nf90_inq_dimid(ncid, yname, ydimid))
|
forrest@0
|
226 |
call nc_err(nf90_inq_dimid(ncid, zsoi_name, zsoi_dimid))
|
forrest@0
|
227 |
call nc_err(nf90_inq_dimid(ncid, zlak_name, zlak_dimid))
|
forrest@0
|
228 |
call nc_err(nf90_inq_dimid(ncid, tname, tdimid))
|
forrest@0
|
229 |
if (nf90_inq_dimid(ncid, hr_name, hr_dimid) == nf90_noerr) then
|
forrest@0
|
230 |
call nc_err(nf90_inq_dimid(ncid, nsamples_name, nsamples_dimid))
|
forrest@0
|
231 |
hrflag = .true.
|
forrest@0
|
232 |
end if
|
forrest@0
|
233 |
|
forrest@0
|
234 |
call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=xdimid, len=xlen))
|
forrest@0
|
235 |
call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=ydimid, len=ylen))
|
forrest@0
|
236 |
call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=zsoi_dimid, len=zsoi_len))
|
forrest@0
|
237 |
call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=zlak_dimid, len=zlak_len))
|
forrest@0
|
238 |
call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=tdimid, len=tlen))
|
forrest@0
|
239 |
if (hrflag) then
|
forrest@0
|
240 |
call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=hr_dimid, len=hr_len))
|
forrest@0
|
241 |
call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=nsamples_dimid, len=nsamples_len))
|
forrest@0
|
242 |
else
|
forrest@0
|
243 |
hr_len = 0
|
forrest@0
|
244 |
nsamples_len = 0
|
forrest@0
|
245 |
end if
|
forrest@0
|
246 |
print *, tlen,' values'
|
forrest@0
|
247 |
|
forrest@0
|
248 |
allocate(xvals(xlen), yvals(ylen), zsoi_vals(zsoi_len), zlak_vals(zlak_len), hr_vals(hr_len), stat=ierr)
|
forrest@0
|
249 |
if (ierr /= 0) then
|
forrest@0
|
250 |
print *, 'clm_read_coord: Cannot allocate grid variables'
|
forrest@0
|
251 |
stop
|
forrest@0
|
252 |
end if
|
forrest@0
|
253 |
if (read_grid_all) then
|
forrest@0
|
254 |
allocate(avals(xlen, ylen), lvals(xlen, ylen), mvals(xlen, ylen), ovals(xlen, ylen), stat=ierr)
|
forrest@0
|
255 |
if (ierr /= 0) then
|
forrest@0
|
256 |
print *, 'clm_read_coord: Cannot allocate extra grid variables'
|
forrest@0
|
257 |
stop
|
forrest@0
|
258 |
end if
|
forrest@0
|
259 |
end if
|
forrest@0
|
260 |
|
forrest@0
|
261 |
call nc_err(nf90_inq_varid(ncid, xname, xvarid))
|
forrest@0
|
262 |
call nc_err(nf90_inq_varid(ncid, yname, yvarid))
|
forrest@0
|
263 |
call nc_err(nf90_inq_varid(ncid, zsoi_name, zsoi_varid))
|
forrest@0
|
264 |
call nc_err(nf90_inq_varid(ncid, zlak_name, zlak_varid))
|
forrest@0
|
265 |
call nc_err(nf90_inq_varid(ncid, tname, tvarid))
|
forrest@0
|
266 |
if (hrflag) call nc_err(nf90_inq_varid(ncid, hr_name, hr_varid))
|
forrest@0
|
267 |
if (read_grid_all) then
|
forrest@0
|
268 |
call nc_err(nf90_inq_varid(ncid, aname, avarid))
|
forrest@0
|
269 |
call nc_err(nf90_inq_varid(ncid, lname, lvarid))
|
forrest@0
|
270 |
call nc_err(nf90_inq_varid(ncid, mname, mvarid))
|
forrest@0
|
271 |
call nc_err(nf90_inq_varid(ncid, oname, ovarid))
|
forrest@0
|
272 |
end if
|
forrest@0
|
273 |
|
forrest@0
|
274 |
call nc_err(nf90_get_var(ncid, xvarid, xvals))
|
forrest@0
|
275 |
call nc_err(nf90_get_var(ncid, yvarid, yvals))
|
forrest@0
|
276 |
call nc_err(nf90_get_var(ncid, zsoi_varid, zsoi_vals))
|
forrest@0
|
277 |
call nc_err(nf90_get_var(ncid, zlak_varid, zlak_vals))
|
forrest@0
|
278 |
if (hrflag) call nc_err(nf90_get_var(ncid, hr_varid, hr_vals))
|
forrest@0
|
279 |
if (read_grid_all) then
|
forrest@0
|
280 |
call nc_err(nf90_get_var(ncid, avarid, avals))
|
forrest@0
|
281 |
call nc_err(nf90_get_var(ncid, lvarid, lvals))
|
forrest@0
|
282 |
call nc_err(nf90_get_var(ncid, mvarid, mvals))
|
forrest@0
|
283 |
call nc_err(nf90_get_var(ncid, ovarid, ovals))
|
forrest@0
|
284 |
end if
|
forrest@0
|
285 |
|
forrest@0
|
286 |
! Read global attributes
|
forrest@0
|
287 |
call nc_err(nf90_get_att(ncid, nf90_global, case_id_name, case_id))
|
forrest@0
|
288 |
call nc_err(nf90_get_att(ncid, nf90_global, inidat_name, inidat))
|
forrest@0
|
289 |
call nc_err(nf90_get_att(ncid, nf90_global, surdat_name, surdat))
|
forrest@0
|
290 |
call nc_err(nf90_get_att(ncid, nf90_global, pftdat_name, pftdat))
|
forrest@0
|
291 |
call nc_err(nf90_get_att(ncid, nf90_global, rtmdat_name, rtmdat))
|
forrest@0
|
292 |
|
forrest@0
|
293 |
call nc_err(nf90_close(ncid))
|
forrest@0
|
294 |
|
forrest@0
|
295 |
time_size = tlen
|
forrest@0
|
296 |
! Loop over all remaining files, reading size of time dimension
|
forrest@0
|
297 |
do i = 2, numf
|
forrest@0
|
298 |
print *, 'Opening ',trim(input_path)//'/'//trim(fnames(i))
|
forrest@0
|
299 |
call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(i)), &
|
forrest@0
|
300 |
mode=nf90_nowrite, ncid=ncid))
|
forrest@0
|
301 |
call nc_err(nf90_inq_dimid(ncid, tname, tdimid))
|
forrest@0
|
302 |
call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=tdimid, len=tlen))
|
forrest@0
|
303 |
call nc_err(nf90_close(ncid))
|
forrest@0
|
304 |
print *, tlen,' values'
|
forrest@0
|
305 |
time_size = time_size + tlen
|
forrest@0
|
306 |
end do
|
forrest@0
|
307 |
|
forrest@0
|
308 |
! Allocate and fill the coordinates derived data type
|
forrest@0
|
309 |
coord%xsize = xlen
|
forrest@0
|
310 |
coord%ysize = ylen
|
forrest@0
|
311 |
coord%zsoi_size = zsoi_len
|
forrest@0
|
312 |
coord%zlak_size = zlak_len
|
forrest@0
|
313 |
coord%tsize = time_size
|
forrest@0
|
314 |
coord%hr_size = hr_len
|
forrest@0
|
315 |
coord%nsamples = nsamples_len
|
forrest@0
|
316 |
allocate(coord%x(xlen), coord%y(ylen), coord%zsoi(zsoi_len), coord%zlak(zlak_len), coord%hr(hr_len), stat=ierr)
|
forrest@0
|
317 |
if (ierr /= 0) then
|
forrest@0
|
318 |
print *, 'clm_read_coord: Cannot allocate grid variables'
|
forrest@0
|
319 |
stop
|
forrest@0
|
320 |
end if
|
forrest@0
|
321 |
if (read_grid_all) then
|
forrest@0
|
322 |
allocate(coord%area(xlen, ylen), coord%landfrac(xlen, ylen), coord%landmask(xlen, ylen), coord%orog(xlen, ylen), stat=ierr)
|
forrest@0
|
323 |
if (ierr /= 0) then
|
forrest@0
|
324 |
print *, 'clm_read_coord: Cannot allocate extra grid variables'
|
forrest@0
|
325 |
stop
|
forrest@0
|
326 |
end if
|
forrest@0
|
327 |
end if
|
forrest@0
|
328 |
coord%x = xvals
|
forrest@0
|
329 |
coord%y = yvals
|
forrest@0
|
330 |
coord%zsoi = zsoi_vals
|
forrest@0
|
331 |
coord%zlak = zlak_vals
|
forrest@0
|
332 |
if (hrflag) coord%hr = hr_vals
|
forrest@0
|
333 |
if (read_grid_all) then
|
forrest@0
|
334 |
coord%area = avals
|
forrest@0
|
335 |
coord%landfrac = lvals
|
forrest@0
|
336 |
coord%landmask = mvals
|
forrest@0
|
337 |
coord%orog = ovals
|
forrest@0
|
338 |
end if
|
forrest@0
|
339 |
deallocate(xvals, yvals, zsoi_vals, zlak_vals, hr_vals)
|
forrest@0
|
340 |
if (read_grid_all) then
|
forrest@0
|
341 |
deallocate(avals, lvals, mvals, ovals)
|
forrest@0
|
342 |
end if
|
forrest@0
|
343 |
coord%case_id = trim(case_id)
|
forrest@0
|
344 |
coord%inidat = trim(inidat)
|
forrest@0
|
345 |
coord%surdat = trim(surdat)
|
forrest@0
|
346 |
coord%pftdat = trim(pftdat)
|
forrest@0
|
347 |
coord%rtmdat = trim(rtmdat)
|
forrest@0
|
348 |
coord%hour_flag = hrflag
|
forrest@0
|
349 |
|
forrest@0
|
350 |
end subroutine clm_read_coord
|
forrest@0
|
351 |
!
|
forrest@0
|
352 |
!--------------------------------------------------------------------------
|
forrest@0
|
353 |
!
|
forrest@0
|
354 |
subroutine clm_copy_grid_data(varname)
|
forrest@0
|
355 |
|
forrest@0
|
356 |
implicit none
|
forrest@0
|
357 |
character (len=varnamelen), intent(in) :: varname ! variable name
|
forrest@0
|
358 |
integer :: ierr
|
forrest@0
|
359 |
|
forrest@0
|
360 |
var_data%varname = varname
|
forrest@0
|
361 |
var_data%ndims = 2
|
forrest@0
|
362 |
var_data%missing_value = 1.0e+36
|
forrest@0
|
363 |
if (varname == 'landmask') then
|
forrest@0
|
364 |
var_data%int_type = .true.
|
forrest@0
|
365 |
allocate(var_data%int2d(coord%xsize, coord%ysize), stat=ierr)
|
forrest@0
|
366 |
if (ierr /= 0) then
|
forrest@0
|
367 |
print *, 'clm_copy_grid_data: Cannot allocate 2d integer variable'
|
forrest@0
|
368 |
stop
|
forrest@0
|
369 |
end if
|
forrest@0
|
370 |
else
|
forrest@0
|
371 |
var_data%int_type = .false.
|
forrest@0
|
372 |
allocate(var_data%var2d(coord%xsize, coord%ysize), stat=ierr)
|
forrest@0
|
373 |
if (ierr /= 0) then
|
forrest@0
|
374 |
print *, 'clm_copy_grid_data: Cannot allocate 2d real variable'
|
forrest@0
|
375 |
stop
|
forrest@0
|
376 |
end if
|
forrest@0
|
377 |
end if
|
forrest@0
|
378 |
|
forrest@0
|
379 |
select case (varname)
|
forrest@0
|
380 |
case ('area')
|
forrest@0
|
381 |
var_data%var2d = coord%area
|
forrest@0
|
382 |
case ('landfrac')
|
forrest@0
|
383 |
var_data%var2d = coord%landfrac
|
forrest@0
|
384 |
case ('landmask')
|
forrest@0
|
385 |
var_data%int2d = coord%landmask
|
forrest@0
|
386 |
case ('topo')
|
forrest@0
|
387 |
var_data%var2d = coord%orog
|
forrest@0
|
388 |
case default
|
forrest@0
|
389 |
print *,'clm_copy_grid_data: ',varname,' is not a known grid variable'
|
forrest@0
|
390 |
stop
|
forrest@0
|
391 |
end select
|
forrest@0
|
392 |
|
forrest@0
|
393 |
end subroutine clm_copy_grid_data
|
forrest@0
|
394 |
!
|
forrest@0
|
395 |
!--------------------------------------------------------------------------
|
forrest@0
|
396 |
!
|
forrest@0
|
397 |
subroutine clm_read_data(numf, input_path, fnames, tshift, varname)
|
forrest@0
|
398 |
|
forrest@0
|
399 |
implicit none
|
forrest@0
|
400 |
integer, intent(in) :: numf ! number of filenames
|
forrest@0
|
401 |
character (len=namelen), intent(in) :: input_path ! input path
|
forrest@0
|
402 |
character (len=namelen), intent(in) :: fnames(numf) ! input filenames
|
forrest@0
|
403 |
integer, intent(in) :: tshift(numf) ! time shift by file
|
forrest@0
|
404 |
character (len=varnamelen), intent(in) :: varname ! variable name
|
forrest@0
|
405 |
integer :: ncid
|
forrest@0
|
406 |
integer :: tdimid, tvarid, tlen
|
forrest@0
|
407 |
integer :: hidimid, hilen
|
forrest@0
|
408 |
integer :: tbvarid
|
forrest@0
|
409 |
integer :: varid, dimids(nf90_max_var_dims), dimlen
|
forrest@0
|
410 |
integer :: zsize
|
forrest@0
|
411 |
integer(i8) :: i, j, k, l
|
forrest@0
|
412 |
integer(i8) :: time_size
|
forrest@0
|
413 |
real(r4), allocatable :: tvals(:) ! temporary time values
|
forrest@0
|
414 |
real(r4), allocatable :: tbvals(:,:) ! temporary time bounds values
|
forrest@0
|
415 |
real(r4), allocatable :: vals(:,:) ! temporary data values
|
forrest@0
|
416 |
character (len=varnamelen) :: dimname ! dimension name
|
forrest@0
|
417 |
integer :: ierr
|
forrest@0
|
418 |
|
forrest@0
|
419 |
! Figure out the dimensionality of the variable and variable attributes
|
forrest@0
|
420 |
! from the first file
|
forrest@0
|
421 |
print *, 'Opening ',trim(input_path)//'/'//trim(fnames(1)), &
|
forrest@0
|
422 |
' to check dimensionality of ',trim(varname)
|
forrest@0
|
423 |
call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(1)), mode=nf90_nowrite, ncid=ncid))
|
forrest@0
|
424 |
call nc_err(nf90_inq_varid(ncid, trim(varname), varid))
|
forrest@0
|
425 |
call nc_err(nf90_inquire_variable(ncid, varid, name=var_data%varname, ndims=var_data%ndims, dimids=dimids))
|
forrest@0
|
426 |
! Figure out if there is a Z-dimension and/or an hour-dimension
|
forrest@0
|
427 |
ierr = 0
|
forrest@0
|
428 |
var_data%int_type = .false.
|
forrest@0
|
429 |
var_data%hour_flag = .false.
|
forrest@0
|
430 |
var_data%soil_layer_flag = .false.
|
forrest@0
|
431 |
do i = 1, var_data%ndims
|
forrest@0
|
432 |
call nc_err(nf90_inquire_dimension(ncid, dimids(i), name=dimname, len=dimlen))
|
forrest@0
|
433 |
if (i == 1 .and. dimname /= xname) ierr = ierr + 1
|
forrest@0
|
434 |
if (i == 2 .and. dimname /= yname) ierr = ierr + 1
|
forrest@0
|
435 |
if (i > 2 .and. i < var_data%ndims) then
|
forrest@0
|
436 |
if (dimname == hr_name) then
|
forrest@0
|
437 |
var_data%hour_flag = .true.
|
forrest@0
|
438 |
else if (dimname == zsoi_name) then
|
forrest@0
|
439 |
var_data%soil_layer_flag = .true.
|
forrest@0
|
440 |
!else if (dimname == zlak_name) then
|
forrest@0
|
441 |
! var_data%zdim_name = dimname
|
forrest@0
|
442 |
else
|
forrest@0
|
443 |
ierr = ierr + 1
|
forrest@0
|
444 |
end if
|
forrest@0
|
445 |
end if
|
forrest@0
|
446 |
! This is not a good assumption (that time is always there), so skip it
|
forrest@0
|
447 |
!if (i == var_data%ndims .and. dimname /= tname) ierr = ierr + 1
|
forrest@0
|
448 |
print *,'Dimension ',i,' is ',dimname
|
forrest@0
|
449 |
end do
|
forrest@0
|
450 |
if (ierr /= 0) then
|
forrest@0
|
451 |
print *,'clm_read_data: Unexpected dimension combination'
|
forrest@0
|
452 |
stop
|
forrest@0
|
453 |
end if
|
forrest@0
|
454 |
!
|
forrest@0
|
455 |
call nc_err(nf90_get_att(ncid, varid, long_name_name, var_data%long_name))
|
forrest@0
|
456 |
call nc_err(nf90_get_att(ncid, varid, missing_value_name, var_data%missing_value))
|
forrest@0
|
457 |
call nc_err(nf90_get_att(ncid, varid, units_name, var_data%units))
|
forrest@0
|
458 |
call nc_err(nf90_close(ncid))
|
forrest@0
|
459 |
|
forrest@0
|
460 |
! Allocate space for the variable of the correct dimensionality
|
forrest@0
|
461 |
select case (var_data%ndims)
|
forrest@0
|
462 |
case (2)
|
forrest@0
|
463 |
allocate(var_data%var2d(coord%xsize, coord%ysize), stat=ierr)
|
forrest@0
|
464 |
if (ierr /= 0) then
|
forrest@0
|
465 |
print *, 'clm_read_data: Cannot allocate 2d variable'
|
forrest@0
|
466 |
stop
|
forrest@0
|
467 |
end if
|
forrest@0
|
468 |
case (3)
|
forrest@0
|
469 |
allocate(var_data%var3d(coord%xsize, coord%ysize, coord%tsize), stat=ierr)
|
forrest@0
|
470 |
if (ierr /= 0) then
|
forrest@0
|
471 |
print *, 'clm_read_data: Cannot allocate 3d variable'
|
forrest@0
|
472 |
stop
|
forrest@0
|
473 |
end if
|
forrest@0
|
474 |
case (4)
|
forrest@0
|
475 |
if (var_data%hour_flag) then
|
forrest@0
|
476 |
zsize = coord%hr_size
|
forrest@0
|
477 |
else if (var_data%soil_layer_flag) then
|
forrest@0
|
478 |
zsize = coord%zsoi_size
|
forrest@0
|
479 |
!else if (var_data%zdim_name == zlak_name) then
|
forrest@0
|
480 |
! zsize = coord%zlak_size
|
forrest@0
|
481 |
else
|
forrest@0
|
482 |
print *,'clm_read_data: Unexpected fourth dimension while allocating'
|
forrest@0
|
483 |
stop
|
forrest@0
|
484 |
end if
|
forrest@0
|
485 |
allocate(var_data%var4d(coord%xsize, coord%ysize, zsize, coord%tsize), stat=ierr)
|
forrest@0
|
486 |
if (ierr /= 0) then
|
forrest@0
|
487 |
print *, 'clm_read_data: Cannot allocate 4d variable'
|
forrest@0
|
488 |
stop
|
forrest@0
|
489 |
end if
|
forrest@0
|
490 |
case default
|
forrest@0
|
491 |
print *, 'clm_read_data: Cannot handle variables of ', var_data%ndims, ' dimensions'
|
forrest@0
|
492 |
stop
|
forrest@0
|
493 |
end select
|
forrest@0
|
494 |
|
forrest@0
|
495 |
! Allocate space for time stamps and bounds
|
forrest@0
|
496 |
allocate(var_data%time(coord%tsize), stat=ierr)
|
forrest@0
|
497 |
if (ierr /= 0) then
|
forrest@0
|
498 |
print *, 'clm_read_data: Cannot allocate time variable'
|
forrest@0
|
499 |
stop
|
forrest@0
|
500 |
end if
|
forrest@0
|
501 |
allocate(var_data%time_bounds(2,coord%tsize), stat=ierr)
|
forrest@0
|
502 |
if (ierr /= 0) then
|
forrest@0
|
503 |
print *, 'clm_read_data: Cannot allocate time bounds variable'
|
forrest@0
|
504 |
stop
|
forrest@0
|
505 |
end if
|
forrest@0
|
506 |
|
forrest@0
|
507 |
! Allocate temporary space for the variable level-slice/time-slice */
|
forrest@0
|
508 |
allocate(vals(coord%xsize, coord%ysize), stat=ierr)
|
forrest@0
|
509 |
if (ierr /= 0) then
|
forrest@0
|
510 |
print *, 'clm_read_data: Cannot allocate variable time-slice'
|
forrest@0
|
511 |
stop
|
forrest@0
|
512 |
end if
|
forrest@0
|
513 |
|
forrest@0
|
514 |
if (var_data%ndims == 2) then
|
forrest@0
|
515 |
! No time axis, so just read the first instance from the first file
|
forrest@0
|
516 |
print *, 'Opening ',trim(input_path)//'/'//trim(fnames(1))
|
forrest@0
|
517 |
call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(1)), mode=nf90_nowrite, ncid=ncid))
|
forrest@0
|
518 |
call nc_err(nf90_inq_varid(ncid, trim(varname), varid))
|
forrest@0
|
519 |
call nc_err(nf90_get_var(ncid, varid, vals, start = (/ 1, 1 /), count = (/ coord%xsize, coord%ysize /)))
|
forrest@0
|
520 |
var_data%var2d(:,:) = vals(:,:)
|
forrest@0
|
521 |
call nc_err(nf90_close(ncid))
|
forrest@0
|
522 |
else
|
forrest@0
|
523 |
! Otherwise, loop over all times
|
forrest@0
|
524 |
time_size = 0
|
forrest@0
|
525 |
! Loop over all files, reading data
|
forrest@0
|
526 |
do i = 1, numf
|
forrest@0
|
527 |
print *, 'Opening ',trim(input_path)//'/'//trim(fnames(i))
|
forrest@0
|
528 |
call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(i)), mode=nf90_nowrite, ncid=ncid))
|
forrest@0
|
529 |
call nc_err(nf90_inq_dimid(ncid, tname, tdimid))
|
forrest@0
|
530 |
call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=tdimid, len=tlen))
|
forrest@0
|
531 |
call nc_err(nf90_inq_dimid(ncid, hiname, hidimid))
|
forrest@0
|
532 |
call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=hidimid, len=hilen))
|
forrest@0
|
533 |
|
forrest@0
|
534 |
call nc_err(nf90_inq_varid(ncid, tname, tvarid))
|
forrest@0
|
535 |
call nc_err(nf90_inq_varid(ncid, tbname, tbvarid))
|
forrest@0
|
536 |
call nc_err(nf90_inq_varid(ncid, trim(varname), varid))
|
forrest@0
|
537 |
|
forrest@0
|
538 |
! Allocate temporary space for the time stamps and bounds
|
forrest@0
|
539 |
allocate(tvals(tlen), stat=ierr)
|
forrest@0
|
540 |
if (ierr /= 0) then
|
forrest@0
|
541 |
print *, 'clm_read_data: Cannot allocate temporary time variable'
|
forrest@0
|
542 |
stop
|
forrest@0
|
543 |
end if
|
forrest@0
|
544 |
allocate(tbvals(hilen,tlen), stat=ierr)
|
forrest@0
|
545 |
if (ierr /= 0) then
|
forrest@0
|
546 |
print *, 'clm_read_data: Cannot allocate temporary time bounds variable'
|
forrest@0
|
547 |
stop
|
forrest@0
|
548 |
end if
|
forrest@0
|
549 |
|
forrest@0
|
550 |
! Get the time stamps and bounds
|
forrest@0
|
551 |
call nc_err(nf90_get_var(ncid, tvarid, tvals))
|
forrest@0
|
552 |
call nc_err(nf90_get_var(ncid, tbvarid, tbvals))
|
forrest@0
|
553 |
|
forrest@0
|
554 |
do j = 1, tlen
|
forrest@0
|
555 |
! Get and store a slab of data
|
forrest@0
|
556 |
select case (var_data%ndims)
|
forrest@0
|
557 |
case (3)
|
forrest@0
|
558 |
!print *, 'Reading slice of data'
|
forrest@0
|
559 |
call nc_err(nf90_get_var(ncid, varid, vals, start = (/ 1, 1, j /), count = (/ coord%xsize, coord%ysize, 1 /)))
|
forrest@0
|
560 |
var_data%var3d(:,:,time_size + j) = vals(:,:)
|
forrest@0
|
561 |
case (4)
|
forrest@0
|
562 |
do k = 1, zsize
|
forrest@0
|
563 |
call nc_err(nf90_get_var(ncid, varid, vals, start = (/ 1, 1, k, j /), count = (/ coord%xsize, coord%ysize, 1, 1 /)))
|
forrest@0
|
564 |
var_data%var4d(:,:,k,time_size + j) = vals(:,:)
|
forrest@0
|
565 |
end do
|
forrest@0
|
566 |
case default
|
forrest@0
|
567 |
print *, 'clm_read_data: Cannot handle variables of ', var_data%ndims, ' dimensions'
|
forrest@0
|
568 |
stop
|
forrest@0
|
569 |
end select
|
forrest@0
|
570 |
! Store each time stamp and bounds, including the time offset
|
forrest@0
|
571 |
var_data%time(time_size + j) = tvals(j) + tshift(i)
|
forrest@0
|
572 |
var_data%time_bounds(:,time_size + j) = tbvals(:,j) + tshift(i)
|
forrest@0
|
573 |
|
forrest@0
|
574 |
end do
|
forrest@0
|
575 |
|
forrest@0
|
576 |
call nc_err(nf90_close(ncid))
|
forrest@0
|
577 |
|
forrest@0
|
578 |
deallocate(tvals)
|
forrest@0
|
579 |
deallocate(tbvals)
|
forrest@0
|
580 |
|
forrest@0
|
581 |
time_size = time_size + tlen
|
forrest@0
|
582 |
|
forrest@0
|
583 |
end do
|
forrest@0
|
584 |
end if
|
forrest@0
|
585 |
|
forrest@0
|
586 |
deallocate(vals)
|
forrest@0
|
587 |
|
forrest@0
|
588 |
end subroutine clm_read_data
|
forrest@0
|
589 |
!
|
forrest@0
|
590 |
!--------------------------------------------------------------------------
|
forrest@0
|
591 |
!
|
forrest@0
|
592 |
subroutine clm_convert_data(casa_flux_bug)
|
forrest@0
|
593 |
|
forrest@0
|
594 |
implicit none
|
forrest@0
|
595 |
logical, intent(in) :: casa_flux_bug
|
forrest@0
|
596 |
real(r8) :: sfactor
|
forrest@0
|
597 |
real(r8) :: offset
|
forrest@0
|
598 |
integer :: ierr
|
forrest@0
|
599 |
logical, allocatable :: mask2d(:,:)
|
forrest@0
|
600 |
logical, allocatable :: mask3d(:,:,:)
|
forrest@0
|
601 |
logical, allocatable :: mask4d(:,:,:,:)
|
forrest@0
|
602 |
|
forrest@0
|
603 |
! Default offset is zero
|
forrest@0
|
604 |
offset = 0.
|
forrest@0
|
605 |
|
forrest@0
|
606 |
! Default var_data structure values
|
forrest@0
|
607 |
var_data%positive = ''
|
forrest@0
|
608 |
|
forrest@0
|
609 |
! Set scale factor according to which variable is to be converted
|
forrest@0
|
610 |
select case (var_data%varname)
|
forrest@0
|
611 |
|
forrest@0
|
612 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
613 |
! Atmospheric forcing
|
forrest@0
|
614 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
615 |
! husf
|
forrest@0
|
616 |
case ('QBOT')
|
forrest@0
|
617 |
! QBOT: kg kg-1 -> kg kg-1
|
forrest@0
|
618 |
sfactor = 1.0
|
forrest@0
|
619 |
var_data%units = 'kg kg-1'
|
forrest@0
|
620 |
var_data%out_varname = 'husf'
|
forrest@0
|
621 |
! prra
|
forrest@0
|
622 |
case ('RAIN')
|
forrest@0
|
623 |
! RAIN: mm s-1 -> kg m-2 s-1
|
forrest@0
|
624 |
sfactor = 1.0
|
forrest@0
|
625 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
626 |
var_data%out_varname = 'prra'
|
forrest@0
|
627 |
! prsn
|
forrest@0
|
628 |
case ('SNOW')
|
forrest@0
|
629 |
! SNOW: mm s-1 -> kg m-2 s-1
|
forrest@0
|
630 |
sfactor = 1.0
|
forrest@0
|
631 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
632 |
var_data%out_varname = 'prsn'
|
forrest@0
|
633 |
|
forrest@0
|
634 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
635 |
! Biogeochemistry
|
forrest@0
|
636 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
637 |
! agbc
|
forrest@0
|
638 |
! aglbc
|
forrest@0
|
639 |
! agnpp
|
forrest@0
|
640 |
case ('AGNPP')
|
forrest@0
|
641 |
! AGNPP: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
642 |
sfactor = 1.0e-3
|
forrest@0
|
643 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
644 |
var_data%out_varname = 'agnpp'
|
forrest@0
|
645 |
! ar
|
forrest@0
|
646 |
case ('AR')
|
forrest@0
|
647 |
! AR: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
648 |
sfactor = 1.0e-3
|
forrest@0
|
649 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
650 |
var_data%out_varname = 'ar'
|
forrest@0
|
651 |
! bco
|
forrest@0
|
652 |
case ('BIOGENCO')
|
forrest@0
|
653 |
! BIOGENCO: ug m-2 h-1 -> kg m-2 s-1
|
forrest@0
|
654 |
sfactor = 1.0e-9 * (1.0 / 3600.)
|
forrest@0
|
655 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
656 |
var_data%out_varname = 'bco'
|
forrest@0
|
657 |
! bgbc
|
forrest@0
|
658 |
! bgnpp
|
forrest@0
|
659 |
case ('BGNPP')
|
forrest@0
|
660 |
! BGNPP: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
661 |
sfactor = 1.0e-3
|
forrest@0
|
662 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
663 |
var_data%out_varname = 'bgnpp'
|
forrest@0
|
664 |
! bipn
|
forrest@0
|
665 |
case ('ISOPRENE')
|
forrest@0
|
666 |
! ISOPRENE: ug m-2 h-1 -> kg m-2 s-1
|
forrest@0
|
667 |
sfactor = 1.0e-9 * (1.0 / 3600.)
|
forrest@0
|
668 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
669 |
var_data%out_varname = 'bipn'
|
forrest@0
|
670 |
! bmtpn
|
forrest@0
|
671 |
case ('MONOTERP')
|
forrest@0
|
672 |
! MONOTERP: ug m-2 h-1 -> kg m-2 s-1
|
forrest@0
|
673 |
sfactor = 1.0e-9 * (1.0 / 3600.)
|
forrest@0
|
674 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
675 |
var_data%out_varname = 'bmtpn'
|
forrest@0
|
676 |
! borvoc
|
forrest@0
|
677 |
case ('ORVOC')
|
forrest@0
|
678 |
! ORVOC: ug m-2 h-1 -> kg m-2 s-1
|
forrest@0
|
679 |
sfactor = 1.0e-9 * (1.0 / 3600.)
|
forrest@0
|
680 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
681 |
var_data%out_varname = 'borvoc'
|
forrest@0
|
682 |
! bovoc
|
forrest@0
|
683 |
case ('OVOC')
|
forrest@0
|
684 |
! OVOC: ug m-2 h-1 -> kg m-2 s-1
|
forrest@0
|
685 |
sfactor = 1.0e-9 * (1.0 / 3600.)
|
forrest@0
|
686 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
687 |
var_data%out_varname = 'bovoc'
|
forrest@0
|
688 |
! cwdc
|
forrest@0
|
689 |
case ('CWDC')
|
forrest@0
|
690 |
! CWDC: g m-2 -> kg m-2
|
forrest@0
|
691 |
sfactor = 1.0e-3
|
forrest@0
|
692 |
var_data%units = 'kg m-2'
|
forrest@0
|
693 |
var_data%out_varname = 'cwdc'
|
forrest@0
|
694 |
! cwdc_hr
|
forrest@0
|
695 |
case ('CWDC_HR')
|
forrest@0
|
696 |
! CWDC_HR: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
697 |
sfactor = 1.0e-3
|
forrest@0
|
698 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
699 |
var_data%out_varname = 'cwdc_hr'
|
forrest@0
|
700 |
! cwdc_loss
|
forrest@0
|
701 |
case ('CWDC_LOSS')
|
forrest@0
|
702 |
! CWDC_LOSS: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
703 |
sfactor = 1.0e-3
|
forrest@0
|
704 |
if (casa_flux_bug) then
|
forrest@0
|
705 |
print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
|
forrest@0
|
706 |
sfactor = sfactor * (1. / 1200.)
|
forrest@0
|
707 |
end if
|
forrest@0
|
708 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
709 |
var_data%out_varname = 'cwdc_loss'
|
forrest@0
|
710 |
! froot_litterc
|
forrest@0
|
711 |
! froot_litterc_hr
|
forrest@0
|
712 |
! froot_litterc_loss
|
forrest@0
|
713 |
! frootc
|
forrest@0
|
714 |
case ('FROOTC')
|
forrest@0
|
715 |
! FROOTC: g m-2 -> kg m-2
|
forrest@0
|
716 |
sfactor = 1.0e-3
|
forrest@0
|
717 |
var_data%units = 'kg m-2'
|
forrest@0
|
718 |
var_data%out_varname = 'frootc'
|
forrest@0
|
719 |
! frootc_alloc
|
forrest@0
|
720 |
case ('FROOTC_ALLOC')
|
forrest@0
|
721 |
! FROOTC_ALLOC: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
722 |
sfactor = 1.0e-3
|
forrest@0
|
723 |
if (casa_flux_bug) then
|
forrest@0
|
724 |
print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
|
forrest@0
|
725 |
sfactor = sfactor * (1. / 1200.)
|
forrest@0
|
726 |
end if
|
forrest@0
|
727 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
728 |
var_data%out_varname = 'frootc_alloc'
|
forrest@0
|
729 |
! frootc_loss
|
forrest@0
|
730 |
case ('FROOTC_LOSS')
|
forrest@0
|
731 |
! FROOTC_LOSS: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
732 |
sfactor = 1.0e-3
|
forrest@0
|
733 |
if (casa_flux_bug) then
|
forrest@0
|
734 |
print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
|
forrest@0
|
735 |
sfactor = sfactor * (1. / 1200.)
|
forrest@0
|
736 |
end if
|
forrest@0
|
737 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
738 |
var_data%out_varname = 'frootc_loss'
|
forrest@0
|
739 |
! gpp
|
forrest@0
|
740 |
case ('GPP')
|
forrest@0
|
741 |
! GPP: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
742 |
sfactor = 1.0e-3
|
forrest@0
|
743 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
744 |
var_data%out_varname = 'gpp'
|
forrest@0
|
745 |
! gnmin
|
forrest@0
|
746 |
case ('GROSS_NMIN')
|
forrest@0
|
747 |
! GROSS_NMIN: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
748 |
sfactor = 1.0e-3
|
forrest@0
|
749 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
750 |
var_data%out_varname = 'gnmin'
|
forrest@0
|
751 |
! hr
|
forrest@0
|
752 |
case ('HR')
|
forrest@0
|
753 |
! HR: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
754 |
sfactor = 1.0e-3
|
forrest@0
|
755 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
756 |
var_data%out_varname = 'hr'
|
forrest@0
|
757 |
! leaf_litterc
|
forrest@0
|
758 |
! leaf_litterc_hr
|
forrest@0
|
759 |
! leaf_litterc_loss
|
forrest@0
|
760 |
! leafc
|
forrest@0
|
761 |
case ('LEAFC')
|
forrest@0
|
762 |
! LEAFC: g m-2 -> kg m-2
|
forrest@0
|
763 |
sfactor = 1.0e-3
|
forrest@0
|
764 |
var_data%units = 'kg m-2'
|
forrest@0
|
765 |
var_data%out_varname = 'leafc'
|
forrest@0
|
766 |
! leafc_alloc
|
forrest@0
|
767 |
case ('LEAFC_ALLOC')
|
forrest@0
|
768 |
! LEAFC_ALLOC: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
769 |
sfactor = 1.0e-3
|
forrest@0
|
770 |
if (casa_flux_bug) then
|
forrest@0
|
771 |
print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
|
forrest@0
|
772 |
sfactor = sfactor * (1. / 1200.)
|
forrest@0
|
773 |
end if
|
forrest@0
|
774 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
775 |
var_data%out_varname = 'leafc_alloc'
|
forrest@0
|
776 |
! leafc_loss
|
forrest@0
|
777 |
case ('LEAFC_LOSS')
|
forrest@0
|
778 |
! LEAFC_LOSS: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
779 |
sfactor = 1.0e-3
|
forrest@0
|
780 |
if (casa_flux_bug) then
|
forrest@0
|
781 |
print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
|
forrest@0
|
782 |
sfactor = sfactor * (1. / 1200.)
|
forrest@0
|
783 |
end if
|
forrest@0
|
784 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
785 |
var_data%out_varname = 'leafc_loss'
|
forrest@0
|
786 |
! litterc
|
forrest@0
|
787 |
case ('LITTERC')
|
forrest@0
|
788 |
! LITTERC: g m-2 -> kg m-2
|
forrest@0
|
789 |
sfactor = 1.0e-3
|
forrest@0
|
790 |
var_data%units = 'kg m-2'
|
forrest@0
|
791 |
var_data%out_varname = 'litterc'
|
forrest@0
|
792 |
! litterc_hr
|
forrest@0
|
793 |
case ('LITTERC_HR')
|
forrest@0
|
794 |
! LITTERC_HR: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
795 |
sfactor = 1.0e-3
|
forrest@0
|
796 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
797 |
var_data%out_varname = 'litterc_hr'
|
forrest@0
|
798 |
! litterc_loss
|
forrest@0
|
799 |
case ('LITTERC_LOSS')
|
forrest@0
|
800 |
! LITTERC_LOSS: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
801 |
sfactor = 1.0e-3
|
forrest@0
|
802 |
if (casa_flux_bug) then
|
forrest@0
|
803 |
print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
|
forrest@0
|
804 |
sfactor = sfactor * (1. / 1200.)
|
forrest@0
|
805 |
end if
|
forrest@0
|
806 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
807 |
var_data%out_varname = 'litterc_loss'
|
forrest@0
|
808 |
! nee
|
forrest@0
|
809 |
case ('NEE')
|
forrest@0
|
810 |
! NEE: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
811 |
sfactor = 1.0e-3
|
forrest@0
|
812 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
813 |
var_data%out_varname = 'nee'
|
forrest@0
|
814 |
! nep
|
forrest@0
|
815 |
case ('NEP')
|
forrest@0
|
816 |
! NEP: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
817 |
sfactor = 1.0e-3
|
forrest@0
|
818 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
819 |
var_data%out_varname = 'nep'
|
forrest@0
|
820 |
! nnmin
|
forrest@0
|
821 |
case ('NET_NMIN')
|
forrest@0
|
822 |
! NEP: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
823 |
sfactor = 1.0e-3
|
forrest@0
|
824 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
825 |
var_data%out_varname = 'nnmin'
|
forrest@0
|
826 |
! npp
|
forrest@0
|
827 |
case ('NPP')
|
forrest@0
|
828 |
! NPP: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
829 |
sfactor = 1.0e-3
|
forrest@0
|
830 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
831 |
var_data%out_varname = 'npp'
|
forrest@0
|
832 |
! soilc
|
forrest@0
|
833 |
case ('SOILC')
|
forrest@0
|
834 |
! SOILC: g m-2 -> kg m-2
|
forrest@0
|
835 |
sfactor = 1.0e-3
|
forrest@0
|
836 |
var_data%units = 'kg m-2'
|
forrest@0
|
837 |
var_data%out_varname = 'soilc'
|
forrest@0
|
838 |
! soilc_hr
|
forrest@0
|
839 |
case ('SOILC_HR')
|
forrest@0
|
840 |
! SOILC_HR: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
841 |
sfactor = 1.0e-3
|
forrest@0
|
842 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
843 |
var_data%out_varname = 'soilc_hr'
|
forrest@0
|
844 |
! soilc_loss
|
forrest@0
|
845 |
case ('SOILC_LOSS')
|
forrest@0
|
846 |
! SOILC_LOSS: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
847 |
sfactor = 1.0e-3
|
forrest@0
|
848 |
if (casa_flux_bug) then
|
forrest@0
|
849 |
print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
|
forrest@0
|
850 |
sfactor = sfactor * (1. / 1200.)
|
forrest@0
|
851 |
end if
|
forrest@0
|
852 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
853 |
var_data%out_varname = 'soilc_loss'
|
forrest@0
|
854 |
! tbvoc
|
forrest@0
|
855 |
case ('VOCFLXT')
|
forrest@0
|
856 |
! VOCFLXT: ug m-2 h-1 -> kg m-2 s-1
|
forrest@0
|
857 |
sfactor = 1.0e-9 * (1.0 / 3600.)
|
forrest@0
|
858 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
859 |
var_data%out_varname = 'tbvoc'
|
forrest@0
|
860 |
! woodc
|
forrest@0
|
861 |
case ('WOODC')
|
forrest@0
|
862 |
! WOODC: g m-2 -> kg m-2
|
forrest@0
|
863 |
sfactor = 1.0e-3
|
forrest@0
|
864 |
var_data%units = 'kg m-2'
|
forrest@0
|
865 |
var_data%out_varname = 'woodc'
|
forrest@0
|
866 |
! woodc_alloc
|
forrest@0
|
867 |
case ('WOODC_ALLOC')
|
forrest@0
|
868 |
! WOODC_ALLOC: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
869 |
sfactor = 1.0e-3
|
forrest@0
|
870 |
if (casa_flux_bug) then
|
forrest@0
|
871 |
print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
|
forrest@0
|
872 |
sfactor = sfactor * (1. / 1200.)
|
forrest@0
|
873 |
end if
|
forrest@0
|
874 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
875 |
var_data%out_varname = 'woodc_alloc'
|
forrest@0
|
876 |
! woodc_loss
|
forrest@0
|
877 |
case ('WOODC_LOSS')
|
forrest@0
|
878 |
! WOODC_LOSS: g m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
879 |
sfactor = 1.0e-3
|
forrest@0
|
880 |
if (casa_flux_bug) then
|
forrest@0
|
881 |
print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
|
forrest@0
|
882 |
sfactor = sfactor * (1. / 1200.)
|
forrest@0
|
883 |
end if
|
forrest@0
|
884 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
885 |
var_data%out_varname = 'woodc_loss'
|
forrest@0
|
886 |
|
forrest@0
|
887 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
888 |
! Energy
|
forrest@0
|
889 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
890 |
! hfg
|
forrest@0
|
891 |
case ('FGR')
|
forrest@0
|
892 |
! FGR: W m-2 -> W m-2
|
forrest@0
|
893 |
sfactor = 1.0
|
forrest@0
|
894 |
var_data%units = 'W m-2'
|
forrest@0
|
895 |
var_data%positive = 'down'
|
forrest@0
|
896 |
var_data%out_varname = 'hfg'
|
forrest@0
|
897 |
! hfls
|
forrest@0
|
898 |
case ('LATENT')
|
forrest@0
|
899 |
! LATENT: W m-2 -> W m-2
|
forrest@0
|
900 |
sfactor = 1.0
|
forrest@0
|
901 |
var_data%units = 'W m-2'
|
forrest@0
|
902 |
var_data%positive = 'up'
|
forrest@0
|
903 |
var_data%out_varname = 'hfls'
|
forrest@0
|
904 |
! hfss
|
forrest@0
|
905 |
case ('FSH')
|
forrest@0
|
906 |
! FSH: W m-2 -> W m-2
|
forrest@0
|
907 |
sfactor = 1.0
|
forrest@0
|
908 |
var_data%units = 'W m-2'
|
forrest@0
|
909 |
var_data%positive = 'up'
|
forrest@0
|
910 |
var_data%out_varname = 'hfss'
|
forrest@0
|
911 |
! hfssg
|
forrest@0
|
912 |
case ('FSH_G')
|
forrest@0
|
913 |
! FSH_G: W m-2 -> W m-2
|
forrest@0
|
914 |
sfactor = 1.0
|
forrest@0
|
915 |
var_data%units = 'W m-2'
|
forrest@0
|
916 |
var_data%positive = 'up'
|
forrest@0
|
917 |
var_data%out_varname = 'hfssg'
|
forrest@0
|
918 |
! hfssv
|
forrest@0
|
919 |
case ('FSH_V')
|
forrest@0
|
920 |
! FSH_V: W m-2 -> W m-2
|
forrest@0
|
921 |
sfactor = 1.0
|
forrest@0
|
922 |
var_data%units = 'W m-2'
|
forrest@0
|
923 |
var_data%positive = 'up'
|
forrest@0
|
924 |
var_data%out_varname = 'hfssv'
|
forrest@0
|
925 |
|
forrest@0
|
926 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
927 |
! Grid
|
forrest@0
|
928 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
929 |
! gca
|
forrest@0
|
930 |
case ('area')
|
forrest@0
|
931 |
! area: km2 -> m2
|
forrest@0
|
932 |
sfactor = 1.0e+3
|
forrest@0
|
933 |
var_data%units = 'm2'
|
forrest@0
|
934 |
var_data%out_varname = 'gca'
|
forrest@0
|
935 |
! lbm
|
forrest@0
|
936 |
case ('landmask')
|
forrest@0
|
937 |
! landmask: 1 -> 1
|
forrest@0
|
938 |
var_data%units = '1'
|
forrest@0
|
939 |
var_data%out_varname = 'lbm'
|
forrest@0
|
940 |
! Do not convert integer landmask, just return
|
forrest@0
|
941 |
return
|
forrest@0
|
942 |
! lcf
|
forrest@0
|
943 |
! lcn
|
forrest@0
|
944 |
! orog
|
forrest@0
|
945 |
case ('topo')
|
forrest@0
|
946 |
! topo: m -> m
|
forrest@0
|
947 |
sfactor = 1.0
|
forrest@0
|
948 |
var_data%units = 'm'
|
forrest@0
|
949 |
var_data%out_varname = 'orog'
|
forrest@0
|
950 |
! sftlf
|
forrest@0
|
951 |
case ('landfrac')
|
forrest@0
|
952 |
! landfrac: 1 -> 1
|
forrest@0
|
953 |
sfactor = 1.0
|
forrest@0
|
954 |
var_data%units = '1'
|
forrest@0
|
955 |
var_data%out_varname = 'sftlf'
|
forrest@0
|
956 |
|
forrest@0
|
957 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
958 |
! Humidity
|
forrest@0
|
959 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
960 |
! hur2m
|
forrest@0
|
961 |
case ('RH2M')
|
forrest@0
|
962 |
! RH2M: % -> %
|
forrest@0
|
963 |
sfactor = 1.0
|
forrest@0
|
964 |
var_data%units = '%'
|
forrest@0
|
965 |
var_data%out_varname = 'hur2m'
|
forrest@0
|
966 |
! hus2m
|
forrest@0
|
967 |
case ('Q2M')
|
forrest@0
|
968 |
! Q2M: kg kg-1 -> kg kg-1
|
forrest@0
|
969 |
sfactor = 1.0
|
forrest@0
|
970 |
var_data%units = 'kg kg-1'
|
forrest@0
|
971 |
var_data%out_varname = 'hus2m'
|
forrest@0
|
972 |
|
forrest@0
|
973 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
974 |
! Hydrology
|
forrest@0
|
975 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
976 |
! btran
|
forrest@0
|
977 |
case ('BTRAN')
|
forrest@0
|
978 |
! BTRAN: 1 -> 1
|
forrest@0
|
979 |
sfactor = 1.0
|
forrest@0
|
980 |
var_data%units = '1'
|
forrest@0
|
981 |
var_data%out_varname = 'btran'
|
forrest@0
|
982 |
! et
|
forrest@0
|
983 |
case ('ET')
|
forrest@0
|
984 |
! ET: m s-1 -> kg m-2 s-1
|
forrest@0
|
985 |
sfactor = 1.0e+3
|
forrest@0
|
986 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
987 |
var_data%out_varname = 'et'
|
forrest@0
|
988 |
! evspsblsoi
|
forrest@0
|
989 |
case ('QSOIL')
|
forrest@0
|
990 |
! QSOIL: mm s-1 -> kg m-2 s-1
|
forrest@0
|
991 |
sfactor = 1.0
|
forrest@0
|
992 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
993 |
var_data%out_varname = 'evspsblsoi'
|
forrest@0
|
994 |
! evspsblveg
|
forrest@0
|
995 |
case ('QVEGE')
|
forrest@0
|
996 |
! QVEGE: mm s-1 -> kg m-2 s-1
|
forrest@0
|
997 |
sfactor = 1.0
|
forrest@0
|
998 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
999 |
var_data%out_varname = 'evspsblveg'
|
forrest@0
|
1000 |
! mrfsl (multi-level)
|
forrest@0
|
1001 |
case ('SOILICE')
|
forrest@0
|
1002 |
! SOILICE: kg m-2 -> kg m-2
|
forrest@0
|
1003 |
sfactor = 1.0
|
forrest@0
|
1004 |
var_data%units = 'kg m-2'
|
forrest@0
|
1005 |
var_data%out_varname = 'mrfsl'
|
forrest@0
|
1006 |
! mrfso
|
forrest@0
|
1007 |
case ('TOTAL_SOILICE')
|
forrest@0
|
1008 |
! TOTAL_SOILICE: kg m-2 -> kg m-2
|
forrest@0
|
1009 |
sfactor = 1.0
|
forrest@0
|
1010 |
var_data%units = 'kg m-2'
|
forrest@0
|
1011 |
var_data%out_varname = 'mrfso'
|
forrest@0
|
1012 |
! mrlsl (multi-level)
|
forrest@0
|
1013 |
case ('SOILLIQ')
|
forrest@0
|
1014 |
! SOILLIQ: kg m-2 -> kg m-2
|
forrest@0
|
1015 |
sfactor = 1.0
|
forrest@0
|
1016 |
var_data%units = 'kg m-2'
|
forrest@0
|
1017 |
var_data%out_varname = 'mrlsl'
|
forrest@0
|
1018 |
! mrlso
|
forrest@0
|
1019 |
case ('TOTAL_SOILLIQ')
|
forrest@0
|
1020 |
! TOTAL_SOILLIQ: kg m-2 -> kg m-2
|
forrest@0
|
1021 |
sfactor = 1.0
|
forrest@0
|
1022 |
var_data%units = 'kg m-2'
|
forrest@0
|
1023 |
var_data%out_varname = 'mrlso'
|
forrest@0
|
1024 |
! mrros
|
forrest@0
|
1025 |
case ('QOVER')
|
forrest@0
|
1026 |
! QOVER: mm s-1 -> kg m-2 s-1
|
forrest@0
|
1027 |
sfactor = 1.0
|
forrest@0
|
1028 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
1029 |
var_data%out_varname = 'mrros'
|
forrest@0
|
1030 |
! mrross
|
forrest@0
|
1031 |
case ('QDRAI')
|
forrest@0
|
1032 |
! QDRAI: mm s-1 -> kg m-2 s-1
|
forrest@0
|
1033 |
sfactor = 1.0
|
forrest@0
|
1034 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
1035 |
var_data%out_varname = 'mrross'
|
forrest@0
|
1036 |
! prveg
|
forrest@0
|
1037 |
case ('QINTR')
|
forrest@0
|
1038 |
! QINTR: mm s-1 -> kg m-2 s-1
|
forrest@0
|
1039 |
sfactor = 1.0
|
forrest@0
|
1040 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
1041 |
var_data%out_varname = 'prveg'
|
forrest@0
|
1042 |
! snd
|
forrest@0
|
1043 |
case ('H2OSNO')
|
forrest@0
|
1044 |
! H2OSNO: mm -> m
|
forrest@0
|
1045 |
sfactor = 1.0e-3
|
forrest@0
|
1046 |
var_data%units = 'm'
|
forrest@0
|
1047 |
var_data%out_varname = 'snd'
|
forrest@0
|
1048 |
! tran
|
forrest@0
|
1049 |
case ('QVEGT')
|
forrest@0
|
1050 |
! QVEGT: mm s-1 -> kg m-2 s-1
|
forrest@0
|
1051 |
sfactor = 1.0
|
forrest@0
|
1052 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
1053 |
var_data%out_varname = 'tran'
|
forrest@0
|
1054 |
! vfwsl (multi-level)
|
forrest@0
|
1055 |
case ('H2OSOI')
|
forrest@0
|
1056 |
! H2OSOI: mm3 mm-3 -> m3 m-3
|
forrest@0
|
1057 |
sfactor = 1.0
|
forrest@0
|
1058 |
var_data%units = 'm3 m-3'
|
forrest@0
|
1059 |
var_data%out_varname = 'vfwsl'
|
forrest@0
|
1060 |
! wpsl (multi-level)
|
forrest@0
|
1061 |
case ('SOILPSI')
|
forrest@1
|
1062 |
! SOILPSI: MPa -> Pa
|
forrest@1
|
1063 |
sfactor = 1.0e+6
|
forrest@0
|
1064 |
var_data%units = 'Pa'
|
forrest@0
|
1065 |
var_data%out_varname = 'wpsl'
|
forrest@0
|
1066 |
|
forrest@0
|
1067 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
1068 |
! Radiation
|
forrest@0
|
1069 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
1070 |
! alsas
|
forrest@0
|
1071 |
case ('ALL_SKY_ALBEDO')
|
forrest@0
|
1072 |
! ALL_SKY_ALBEDO: 1 -> 1
|
forrest@0
|
1073 |
sfactor = 1.0
|
forrest@0
|
1074 |
var_data%units = '1'
|
forrest@0
|
1075 |
var_data%out_varname = 'alsas'
|
forrest@0
|
1076 |
! alsbs
|
forrest@0
|
1077 |
case ('BLACK_SKY_ALBEDO')
|
forrest@0
|
1078 |
! BLACK_SKY_ALBEDO: 1 -> 1
|
forrest@0
|
1079 |
sfactor = 1.0
|
forrest@0
|
1080 |
var_data%units = '1'
|
forrest@0
|
1081 |
var_data%out_varname = 'alsbs'
|
forrest@0
|
1082 |
! laisha
|
forrest@0
|
1083 |
case ('LAISHA')
|
forrest@0
|
1084 |
! LAISHA: m2 m-2 -> m2 m-2
|
forrest@0
|
1085 |
sfactor = 1.0
|
forrest@0
|
1086 |
var_data%units = 'm2 m-2'
|
forrest@0
|
1087 |
var_data%out_varname = 'laisha'
|
forrest@0
|
1088 |
! laisun
|
forrest@0
|
1089 |
case ('LAISUN')
|
forrest@0
|
1090 |
! LAISUN: m2 m-2 -> m2 m-2
|
forrest@0
|
1091 |
sfactor = 1.0
|
forrest@0
|
1092 |
var_data%units = 'm2 m-2'
|
forrest@0
|
1093 |
var_data%out_varname = 'laisun'
|
forrest@0
|
1094 |
! rlds
|
forrest@0
|
1095 |
case ('FLDS')
|
forrest@0
|
1096 |
! FLDS: W m-2 -> W m-2
|
forrest@0
|
1097 |
sfactor = 1.0
|
forrest@0
|
1098 |
var_data%units = 'W m-2'
|
forrest@0
|
1099 |
var_data%positive = 'down'
|
forrest@0
|
1100 |
var_data%out_varname = 'rlds'
|
forrest@0
|
1101 |
! rls
|
forrest@0
|
1102 |
case ('FIRA')
|
forrest@0
|
1103 |
! FIRA: W m-2 -> W m-2
|
forrest@0
|
1104 |
sfactor = 1.0
|
forrest@0
|
1105 |
var_data%units = 'W m-2'
|
forrest@0
|
1106 |
var_data%positive = 'down'
|
forrest@0
|
1107 |
var_data%out_varname = 'rls'
|
forrest@0
|
1108 |
! rlus
|
forrest@0
|
1109 |
case ('FIRE')
|
forrest@0
|
1110 |
! FIRE: W m-2 -> W m-2
|
forrest@0
|
1111 |
sfactor = 1.0
|
forrest@0
|
1112 |
var_data%units = 'W m-2'
|
forrest@0
|
1113 |
var_data%positive = 'up'
|
forrest@0
|
1114 |
var_data%out_varname = 'rlus'
|
forrest@0
|
1115 |
! rs
|
forrest@0
|
1116 |
case ('NETRAD')
|
forrest@0
|
1117 |
! NETRAD: W m-2 -> W m-2
|
forrest@0
|
1118 |
sfactor = 1.0
|
forrest@0
|
1119 |
var_data%units = 'W m-2'
|
forrest@0
|
1120 |
var_data%positive = 'down'
|
forrest@0
|
1121 |
var_data%out_varname = 'rs'
|
forrest@0
|
1122 |
! rsas
|
forrest@0
|
1123 |
case ('FSA')
|
forrest@0
|
1124 |
! FSA: W m-2 -> W m-2
|
forrest@0
|
1125 |
sfactor = 1.0
|
forrest@0
|
1126 |
var_data%units = 'W m-2'
|
forrest@0
|
1127 |
var_data%positive = 'down'
|
forrest@0
|
1128 |
var_data%out_varname = 'rsas'
|
forrest@0
|
1129 |
! rsasg
|
forrest@0
|
1130 |
case ('SABG')
|
forrest@0
|
1131 |
! SABG: W m-2 -> W m-2
|
forrest@0
|
1132 |
sfactor = 1.0
|
forrest@0
|
1133 |
var_data%units = 'W m-2'
|
forrest@0
|
1134 |
var_data%positive = 'down'
|
forrest@0
|
1135 |
var_data%out_varname = 'rsasg'
|
forrest@0
|
1136 |
! rsasv
|
forrest@0
|
1137 |
case ('SABV')
|
forrest@0
|
1138 |
! SABV: W m-2 -> W m-2
|
forrest@0
|
1139 |
sfactor = 1.0
|
forrest@0
|
1140 |
var_data%units = 'W m-2'
|
forrest@0
|
1141 |
var_data%positive = 'down'
|
forrest@0
|
1142 |
var_data%out_varname = 'rsasv'
|
forrest@0
|
1143 |
! rsds
|
forrest@0
|
1144 |
case ('FSDS')
|
forrest@0
|
1145 |
! FSDS: W m-2 -> W m-2
|
forrest@0
|
1146 |
sfactor = 1.0
|
forrest@0
|
1147 |
var_data%units = 'W m-2'
|
forrest@0
|
1148 |
var_data%positive = 'down'
|
forrest@0
|
1149 |
var_data%out_varname = 'rsds'
|
forrest@0
|
1150 |
! rsus
|
forrest@0
|
1151 |
case ('FSR')
|
forrest@0
|
1152 |
! FSR: W m-2 -> W m-2
|
forrest@0
|
1153 |
sfactor = 1.0
|
forrest@0
|
1154 |
var_data%units = 'W m-2'
|
forrest@0
|
1155 |
var_data%positive = 'up'
|
forrest@0
|
1156 |
var_data%out_varname = 'rsus'
|
forrest@0
|
1157 |
! snc
|
forrest@0
|
1158 |
case ('FSNO')
|
forrest@0
|
1159 |
! FSNO: 1 -> 1
|
forrest@0
|
1160 |
sfactor = 1.0
|
forrest@0
|
1161 |
var_data%units = '1'
|
forrest@0
|
1162 |
var_data%out_varname = 'snc'
|
forrest@0
|
1163 |
|
forrest@0
|
1164 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
1165 |
! Temperature
|
forrest@0
|
1166 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
1167 |
! degday
|
forrest@0
|
1168 |
case ('DEGDAY')
|
forrest@0
|
1169 |
! DEGDAY: C -> K
|
forrest@0
|
1170 |
sfactor = 1.0
|
forrest@0
|
1171 |
offset = 273.15
|
forrest@0
|
1172 |
var_data%units = 'K'
|
forrest@0
|
1173 |
var_data%out_varname = 'degday'
|
forrest@0
|
1174 |
! ta2m
|
forrest@0
|
1175 |
case ('TSA')
|
forrest@0
|
1176 |
! TSA: K -> K
|
forrest@0
|
1177 |
sfactor = 1.0
|
forrest@0
|
1178 |
var_data%units = 'K'
|
forrest@0
|
1179 |
var_data%out_varname = 'ta2m'
|
forrest@0
|
1180 |
! ta2mdmax
|
forrest@0
|
1181 |
case ('TREFMXAV')
|
forrest@0
|
1182 |
! TREFMXAV: K -> K
|
forrest@0
|
1183 |
sfactor = 1.0
|
forrest@0
|
1184 |
var_data%units = 'K'
|
forrest@0
|
1185 |
var_data%out_varname = 'ta2mdmax'
|
forrest@0
|
1186 |
! ta2mdmin
|
forrest@0
|
1187 |
case ('TREFMNAV')
|
forrest@0
|
1188 |
! TREFMNAV: K -> K
|
forrest@0
|
1189 |
sfactor = 1.0
|
forrest@0
|
1190 |
var_data%units = 'K'
|
forrest@0
|
1191 |
var_data%out_varname = 'ta2mdmin'
|
forrest@0
|
1192 |
! tsl (multi-level)
|
forrest@0
|
1193 |
case ('TSOI')
|
forrest@0
|
1194 |
! TSOI: K -> K
|
forrest@0
|
1195 |
sfactor = 1.0
|
forrest@0
|
1196 |
var_data%units = 'K'
|
forrest@0
|
1197 |
var_data%out_varname = 'tsl'
|
forrest@0
|
1198 |
! tv
|
forrest@0
|
1199 |
case ('TV')
|
forrest@0
|
1200 |
! TV: K -> K
|
forrest@0
|
1201 |
sfactor = 1.0
|
forrest@0
|
1202 |
var_data%units = 'K'
|
forrest@0
|
1203 |
var_data%out_varname = 'tv'
|
forrest@0
|
1204 |
|
forrest@0
|
1205 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
1206 |
! Vegetation phenology
|
forrest@0
|
1207 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
1208 |
! laiexp
|
forrest@0
|
1209 |
case ('ELAI')
|
forrest@0
|
1210 |
! ELAI: m2 m-2 -> m2 m-2
|
forrest@0
|
1211 |
sfactor = 1.0
|
forrest@0
|
1212 |
var_data%units = 'm2 m-2'
|
forrest@0
|
1213 |
var_data%out_varname = 'laiexp'
|
forrest@0
|
1214 |
! laitot
|
forrest@0
|
1215 |
case ('TLAI')
|
forrest@0
|
1216 |
! TLAI: m2 m-2 -> m2 m-2
|
forrest@0
|
1217 |
sfactor = 1.0
|
forrest@0
|
1218 |
var_data%units = 'm2 m-2'
|
forrest@0
|
1219 |
var_data%out_varname = 'laitot'
|
forrest@0
|
1220 |
! saiexp
|
forrest@0
|
1221 |
case ('ESAI')
|
forrest@0
|
1222 |
! ESAI: m2 m-2 -> m2 m-2
|
forrest@0
|
1223 |
sfactor = 1.0
|
forrest@0
|
1224 |
var_data%units = 'm2 m-2'
|
forrest@0
|
1225 |
var_data%out_varname = 'saiexp'
|
forrest@0
|
1226 |
! saitot
|
forrest@0
|
1227 |
case ('TSAI')
|
forrest@0
|
1228 |
! TSAI: m2 m-2 -> m2 m-2
|
forrest@0
|
1229 |
sfactor = 1.0
|
forrest@0
|
1230 |
var_data%units = 'm2 m-2'
|
forrest@0
|
1231 |
var_data%out_varname = 'saitot'
|
forrest@0
|
1232 |
|
forrest@0
|
1233 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
1234 |
! Vegetation physiology
|
forrest@0
|
1235 |
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
forrest@0
|
1236 |
! psn
|
forrest@0
|
1237 |
case ('FPSN')
|
forrest@0
|
1238 |
! FPSN: umol m-2 s-1 -> kg m-2 s-1
|
forrest@0
|
1239 |
sfactor = 12.0 * 1.0e-6 * 1.0e-3
|
forrest@0
|
1240 |
var_data%units = 'kg m-2 s-1'
|
forrest@0
|
1241 |
var_data%out_varname = 'psn'
|
forrest@0
|
1242 |
! rssha
|
forrest@0
|
1243 |
case ('RSSHA')
|
forrest@0
|
1244 |
! RSSHA: s m-1 -> s m-1
|
forrest@0
|
1245 |
sfactor = 1.0
|
forrest@0
|
1246 |
var_data%units = 's m-1'
|
forrest@0
|
1247 |
var_data%out_varname = 'rssha'
|
forrest@0
|
1248 |
! rssun
|
forrest@0
|
1249 |
case ('RSSUN')
|
forrest@0
|
1250 |
! RSSUN: s m-1 -> s m-1
|
forrest@0
|
1251 |
sfactor = 1.0
|
forrest@0
|
1252 |
var_data%units = 's m-1'
|
forrest@0
|
1253 |
var_data%out_varname = 'rssun'
|
forrest@0
|
1254 |
|
forrest@0
|
1255 |
case default
|
forrest@0
|
1256 |
print *, 'clm_convert_data: Unknown variable name ', var_data%varname
|
forrest@0
|
1257 |
stop
|
forrest@0
|
1258 |
end select
|
forrest@0
|
1259 |
|
forrest@0
|
1260 |
! Convert to the desired units using the offset and scale factor (sfactor)
|
forrest@0
|
1261 |
select case (var_data%ndims)
|
forrest@0
|
1262 |
case (2)
|
forrest@0
|
1263 |
allocate(mask2d(coord%xsize, coord%ysize), stat=ierr)
|
forrest@0
|
1264 |
if (ierr /= 0) then
|
forrest@0
|
1265 |
print *, 'clm_convert_data: Cannot allocate 2d mask'
|
forrest@0
|
1266 |
stop
|
forrest@0
|
1267 |
end if
|
forrest@0
|
1268 |
mask2d = .false.
|
forrest@0
|
1269 |
if (var_data%int_type) then
|
forrest@0
|
1270 |
where (var_data%int2d /= var_data%missing_value)
|
forrest@0
|
1271 |
var_data%int2d = int(var_data%int2d * sfactor + offset)
|
forrest@0
|
1272 |
mask2d = .true.
|
forrest@0
|
1273 |
end where
|
forrest@0
|
1274 |
print *, var_data%out_varname, &
|
forrest@0
|
1275 |
' min:', minval(var_data%int2d, mask=mask2d), &
|
forrest@0
|
1276 |
' max:', maxval(var_data%int2d, mask=mask2d)
|
forrest@0
|
1277 |
else
|
forrest@0
|
1278 |
where (var_data%var2d /= var_data%missing_value)
|
forrest@0
|
1279 |
var_data%var2d = var_data%var2d * sfactor + offset
|
forrest@0
|
1280 |
mask2d = .true.
|
forrest@0
|
1281 |
end where
|
forrest@0
|
1282 |
print *, var_data%out_varname, &
|
forrest@0
|
1283 |
' min:', minval(var_data%var2d, mask=mask2d), &
|
forrest@0
|
1284 |
' max:', maxval(var_data%var2d, mask=mask2d)
|
forrest@0
|
1285 |
end if
|
forrest@0
|
1286 |
if (allocated(mask2d)) deallocate(mask2d)
|
forrest@0
|
1287 |
case (3)
|
forrest@0
|
1288 |
allocate(mask3d(coord%xsize, coord%ysize, coord%tsize), stat=ierr)
|
forrest@0
|
1289 |
if (ierr /= 0) then
|
forrest@0
|
1290 |
print *, 'clm_convert_data: Cannot allocate 3d mask'
|
forrest@0
|
1291 |
stop
|
forrest@0
|
1292 |
end if
|
forrest@0
|
1293 |
mask3d = .false.
|
forrest@0
|
1294 |
where (var_data%var3d /= var_data%missing_value)
|
forrest@0
|
1295 |
var_data%var3d = (var_data%var3d + offset) * sfactor
|
forrest@0
|
1296 |
mask3d = .true.
|
forrest@0
|
1297 |
end where
|
forrest@0
|
1298 |
print *, var_data%out_varname, &
|
forrest@0
|
1299 |
' min:', minval(var_data%var3d, mask=mask3d), &
|
forrest@0
|
1300 |
' max:', maxval(var_data%var3d, mask=mask3d)
|
forrest@0
|
1301 |
if (allocated(mask3d)) deallocate(mask3d)
|
forrest@0
|
1302 |
case (4)
|
forrest@0
|
1303 |
if (var_data%soil_layer_flag) then
|
forrest@0
|
1304 |
allocate(mask4d(coord%xsize, coord%ysize, coord%zsoi_size, coord%tsize), stat=ierr)
|
forrest@0
|
1305 |
if (ierr /= 0) then
|
forrest@0
|
1306 |
print *, 'clm_convert_data: Cannot allocate 4d mask'
|
forrest@0
|
1307 |
stop
|
forrest@0
|
1308 |
end if
|
forrest@0
|
1309 |
mask4d = .false.
|
forrest@0
|
1310 |
where (var_data%var4d /= var_data%missing_value)
|
forrest@0
|
1311 |
var_data%var4d = (var_data%var4d + offset) * sfactor
|
forrest@0
|
1312 |
mask4d = .true.
|
forrest@0
|
1313 |
end where
|
forrest@0
|
1314 |
print *, var_data%out_varname, &
|
forrest@0
|
1315 |
' min:', minval(var_data%var4d, mask=mask4d), &
|
forrest@0
|
1316 |
' max:', maxval(var_data%var4d, mask=mask4d)
|
forrest@0
|
1317 |
if (allocated(mask4d)) deallocate(mask4d)
|
forrest@0
|
1318 |
else
|
forrest@0
|
1319 |
print *,'clm_convert_data: Unknown data dimension'
|
forrest@0
|
1320 |
stop
|
forrest@0
|
1321 |
end if
|
forrest@0
|
1322 |
case default
|
forrest@0
|
1323 |
print *,'clm_convert_data: Dimension count not supported'
|
forrest@0
|
1324 |
stop
|
forrest@0
|
1325 |
end select
|
forrest@0
|
1326 |
|
forrest@0
|
1327 |
end subroutine clm_convert_data
|
forrest@0
|
1328 |
!
|
forrest@0
|
1329 |
!--------------------------------------------------------------------------
|
forrest@0
|
1330 |
!
|
forrest@0
|
1331 |
subroutine clm_free_data
|
forrest@0
|
1332 |
|
forrest@0
|
1333 |
implicit none
|
forrest@0
|
1334 |
|
forrest@0
|
1335 |
var_data%varname = ''
|
forrest@0
|
1336 |
var_data%ndims = 0
|
forrest@0
|
1337 |
if (allocated(var_data%int2d)) deallocate(var_data%int2d)
|
forrest@0
|
1338 |
if (allocated(var_data%var2d)) deallocate(var_data%var2d)
|
forrest@0
|
1339 |
if (allocated(var_data%var3d)) deallocate(var_data%var3d)
|
forrest@0
|
1340 |
if (allocated(var_data%var4d)) deallocate(var_data%var4d)
|
forrest@0
|
1341 |
if (allocated(var_data%time)) deallocate(var_data%time)
|
forrest@0
|
1342 |
if (allocated(var_data%time_bounds)) deallocate(var_data%time_bounds)
|
forrest@0
|
1343 |
var_data%int_type = .false.
|
forrest@0
|
1344 |
var_data%soil_layer_flag = .false.
|
forrest@0
|
1345 |
var_data%hour_flag = .false.
|
forrest@0
|
1346 |
var_data%long_name = ''
|
forrest@0
|
1347 |
var_data%missing_value = 0.
|
forrest@0
|
1348 |
var_data%units = ''
|
forrest@0
|
1349 |
var_data%positive = ''
|
forrest@0
|
1350 |
var_data%out_varname = ''
|
forrest@0
|
1351 |
|
forrest@0
|
1352 |
end subroutine clm_free_data
|
forrest@0
|
1353 |
!
|
forrest@0
|
1354 |
!--------------------------------------------------------------------------
|
forrest@0
|
1355 |
!
|
forrest@0
|
1356 |
subroutine nc_err(status)
|
forrest@0
|
1357 |
|
forrest@0
|
1358 |
implicit none
|
forrest@0
|
1359 |
integer :: status ! status code
|
forrest@0
|
1360 |
|
forrest@0
|
1361 |
if (status /= nf90_noerr) then
|
forrest@0
|
1362 |
print *, nf90_strerror(status)
|
forrest@0
|
1363 |
stop 'Program stopped after encountering netCDF error.'
|
forrest@0
|
1364 |
end if
|
forrest@0
|
1365 |
|
forrest@0
|
1366 |
end subroutine nc_err
|
forrest@0
|
1367 |
!
|
forrest@0
|
1368 |
!--------------------------------------------------------------------------
|
forrest@0
|
1369 |
!
|
forrest@0
|
1370 |
|
forrest@0
|
1371 |
end module clm_mod
|