Fixed unit conversion for SOILPSI (wpsl) and regenerated files. Changed path to CMOR table in namelist files.
7 integer, parameter :: namelen = 128 ! string length for paths
8 integer, parameter :: varnamelen = 32 ! string length for varnames
9 integer, parameter :: posnamelen = 4 ! string length for positive attribute
10 integer, parameter :: varattlen = 128 ! string length for variable attributes
11 character (len=*), parameter :: xname = 'lon'
12 character (len=*), parameter :: yname = 'lat'
13 character (len=*), parameter :: zsoi_name = 'levsoi'
14 character (len=*), parameter :: zlak_name = 'levlak'
15 character (len=*), parameter :: tname = 'time'
16 character (len=*), parameter :: hr_name = 'hour'
17 character (len=*), parameter :: nsamples_name = 'nsamples'
18 character (len=*), parameter :: hiname = 'hist_interval'
19 character (len=*), parameter :: tbname = 'time_bounds'
20 character (len=*), parameter :: aname = 'area'
21 character (len=*), parameter :: lname = 'landfrac'
22 character (len=*), parameter :: mname = 'landmask'
23 character (len=*), parameter :: oname = 'topo'
24 character (len=*), parameter :: long_name_name = 'long_name'
25 character (len=*), parameter :: missing_value_name = 'missing_value'
26 character (len=*), parameter :: units_name = 'units'
27 character (len=*), parameter :: case_id_name = 'case_id'
28 character (len=*), parameter :: inidat_name = 'Initial_conditions_dataset'
29 character (len=*), parameter :: surdat_name = 'Surface_dataset'
30 character (len=*), parameter :: pftdat_name = 'PFT_physiological_constants_dataset'
31 character (len=*), parameter :: rtmdat_name = 'RTM_input_datset'
34 integer :: xsize ! number of longitude points (lon)
35 integer :: ysize ! number of latitutde points (lat)
36 integer :: psize ! number of pfts
37 real(r8), pointer :: x(:) ! longitude values (lon)
38 real(r8), pointer :: y(:) ! latitude values (lat)
39 real(r8), pointer :: p(:) ! pft number
40 real(r4), pointer :: pft(:,:,:) ! percent of each pft
44 integer :: xsize ! number of longitude points (lon)
45 integer :: ysize ! number of latitutde points (lat)
46 integer :: zsoi_size ! number of soil z-levels (levsoi)
47 integer :: zlak_size ! number of lake z-levels (levlak)
48 integer :: tsize ! number of time points (time)
49 integer :: hr_size ! number of hour points (hour = 24)
50 integer :: nsamples ! number of time samples in post-processing statistics
51 real(r8), pointer :: x(:) ! longitude values (lon)
52 real(r8), pointer :: y(:) ! latitude values (lat)
53 real(r8), pointer :: zsoi(:) ! soil z-level values (levsoi)
54 real(r8), pointer :: zlak(:) ! lake z-level values (levlak)
55 real(r8), pointer :: hr(:) ! hour-of-day values (hour)
56 real(r8), pointer :: area(:,:) ! gridcell areas
57 real(r8), pointer :: landfrac(:,:) ! gridcell land fractions
58 integer, pointer :: landmask(:,:) ! gridcell land/ocean mask
59 real(r8), pointer :: orog(:,:) ! gridcell orography/topography
60 character (len=namelen) :: case_id ! case name of run
61 character (len=namelen) :: inidat ! initial dataset filename
62 character (len=namelen) :: surdat ! surface dataset filename
63 character (len=namelen) :: pftdat ! pft dataset filename
64 character (len=namelen) :: rtmdat ! runoff dataset filename
65 logical :: hour_flag ! flag for hour-of-day statistical summary files
69 character (len=varnamelen) :: varname
71 integer, pointer :: int2d(:,:)
72 real(r4), pointer :: var2d(:,:)
73 real(r4), pointer :: var3d(:,:,:)
74 real(r4), pointer :: var4d(:,:,:,:)
75 real(r8), pointer :: time(:)
76 real(r8), pointer :: time_bounds(:,:)
77 logical :: int_type ! set if type is integer
78 logical :: soil_layer_flag
79 logical :: hour_flag ! flag for hour-of-day statistical summary variable
80 character (len=varattlen) :: long_name
81 real(r4) :: missing_value
82 character (len=varattlen) :: units
83 character (len=posnamelen) :: positive
84 character (len=varnamelen) :: out_varname
85 end type variable_data
87 type(surface), target :: surf
88 type(coordinates), target :: coord
89 type(variable_data), target :: var_data
95 public clm_copy_grid_data
97 public clm_convert_data
104 !--------------------------------------------------------------------------
106 subroutine clm_read_surf(fsurdat)
109 character (len=namelen), intent(in) :: fsurdat ! input surface dataset
113 integer :: xdimid, xvarid, xlen
114 integer :: ydimid, yvarid, ylen
115 integer :: pdimid, pvarid, plen
116 real(r8), allocatable :: longxy(:,:)
117 real(r8), allocatable :: latixy(:,:)
118 real(r4), allocatable :: pctpft(:,:,:)
120 ! Read coordinated from surface dataset
121 print *, 'Opening ',trim(fsurdat)
122 call nc_err(nf90_open(path=trim(fsurdat), mode=nf90_nowrite, ncid=ncid))
123 call nc_err(nf90_inq_dimid(ncid, 'lsmlon', xdimid))
124 call nc_err(nf90_inq_dimid(ncid, 'lsmlat', ydimid))
125 call nc_err(nf90_inq_dimid(ncid, 'lsmpft', pdimid))
127 call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=xdimid, len=xlen))
128 call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=ydimid, len=ylen))
129 call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=pdimid, len=plen))
131 call nc_err(nf90_inq_varid(ncid, 'LONGXY', xvarid))
132 call nc_err(nf90_inq_varid(ncid, 'LATIXY', yvarid))
133 call nc_err(nf90_inq_varid(ncid, 'PCT_PFT', pvarid))
135 allocate(longxy(xlen, ylen), latixy(xlen, ylen), pctpft(xlen, ylen, plen), stat=ierr)
137 print *, 'clm_read_surf: Cannot allocate longxy, latixy, pctpft variables'
141 call nc_err(nf90_get_var(ncid, xvarid, longxy))
142 call nc_err(nf90_get_var(ncid, yvarid, latixy))
143 call nc_err(nf90_get_var(ncid, pvarid, pctpft))
145 print *, 'Closing ',trim(fsurdat)
146 call nc_err(nf90_close(ncid))
148 ! Allocate and fill in the surface derived data type
152 allocate(surf%x(xlen), surf%y(ylen), surf%p(plen), surf%pft(xlen, ylen, plen), stat=ierr)
154 print *, 'clm_read_surf: Cannot allocate surface type variables'
157 surf%x(:) = longxy(:,1)
158 surf%y(:) = latixy(1,:)
164 deallocate(longxy, latixy, pctpft)
166 end subroutine clm_read_surf
168 !--------------------------------------------------------------------------
170 subroutine clm_free_surf()
174 if (allocated(surf%x)) deallocate(surf%x)
175 if (allocated(surf%y)) deallocate(surf%y)
176 if (allocated(surf%p)) deallocate(surf%p)
177 if (allocated(surf%pft)) deallocate(surf%pft)
179 end subroutine clm_free_surf
181 !--------------------------------------------------------------------------
183 subroutine clm_read_coord(numf, input_path, fnames, read_grid_all)
186 integer, intent(in) :: numf ! number of filenames
187 character (len=namelen), intent(in) :: input_path ! input path
188 character (len=namelen), intent(in) :: fnames(numf) ! input filenames
189 logical, intent(in) :: read_grid_all ! flag to read area,...
191 integer :: xdimid, xvarid, xlen
192 integer :: ydimid, yvarid, ylen
193 integer :: zsoi_dimid, zsoi_varid, zsoi_len
194 integer :: zlak_dimid, zlak_varid, zlak_len
195 integer :: tdimid, tvarid, tlen
196 integer :: hr_dimid, hr_varid, hr_len
197 integer :: nsamples_dimid, nsamples_len
198 integer :: avarid ! area (when read_grid_all true)
199 integer :: lvarid ! landfrac (when read_grid_all true)
200 integer :: mvarid ! landmask (when read_grid_all true)
201 integer :: ovarid ! orog/topo (when read_grid_all true)
203 integer(i8) :: time_size
204 real(r4), allocatable :: xvals(:)
205 real(r4), allocatable :: yvals(:)
206 real(r4), allocatable :: zsoi_vals(:)
207 real(r4), allocatable :: zlak_vals(:)
208 real(r4), allocatable :: hr_vals(:)
209 real(r4), allocatable :: avals(:,:) ! area (when read_grid_all true)
210 real(r4), allocatable :: lvals(:,:) ! landfrac (when read_grid_all true)
211 integer, allocatable :: mvals(:,:) ! landmask (when read_grid_all true)
212 real(r4), allocatable :: ovals(:,:) ! orog/topo (when read_grid_all true)
213 character (len=namelen) :: case_id = ''
214 character (len=namelen) :: inidat = ''
215 character (len=namelen) :: surdat = ''
216 character (len=namelen) :: pftdat = ''
217 character (len=namelen) :: rtmdat = ''
218 logical :: hrflag = .false.
221 ! Read coordinates from first file
222 print *, 'Opening ',trim(input_path)//'/'//trim(fnames(1))
223 call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(1)), mode=nf90_nowrite, ncid=ncid))
224 call nc_err(nf90_inq_dimid(ncid, xname, xdimid))
225 call nc_err(nf90_inq_dimid(ncid, yname, ydimid))
226 call nc_err(nf90_inq_dimid(ncid, zsoi_name, zsoi_dimid))
227 call nc_err(nf90_inq_dimid(ncid, zlak_name, zlak_dimid))
228 call nc_err(nf90_inq_dimid(ncid, tname, tdimid))
229 if (nf90_inq_dimid(ncid, hr_name, hr_dimid) == nf90_noerr) then
230 call nc_err(nf90_inq_dimid(ncid, nsamples_name, nsamples_dimid))
234 call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=xdimid, len=xlen))
235 call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=ydimid, len=ylen))
236 call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=zsoi_dimid, len=zsoi_len))
237 call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=zlak_dimid, len=zlak_len))
238 call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=tdimid, len=tlen))
240 call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=hr_dimid, len=hr_len))
241 call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=nsamples_dimid, len=nsamples_len))
246 print *, tlen,' values'
248 allocate(xvals(xlen), yvals(ylen), zsoi_vals(zsoi_len), zlak_vals(zlak_len), hr_vals(hr_len), stat=ierr)
250 print *, 'clm_read_coord: Cannot allocate grid variables'
253 if (read_grid_all) then
254 allocate(avals(xlen, ylen), lvals(xlen, ylen), mvals(xlen, ylen), ovals(xlen, ylen), stat=ierr)
256 print *, 'clm_read_coord: Cannot allocate extra grid variables'
261 call nc_err(nf90_inq_varid(ncid, xname, xvarid))
262 call nc_err(nf90_inq_varid(ncid, yname, yvarid))
263 call nc_err(nf90_inq_varid(ncid, zsoi_name, zsoi_varid))
264 call nc_err(nf90_inq_varid(ncid, zlak_name, zlak_varid))
265 call nc_err(nf90_inq_varid(ncid, tname, tvarid))
266 if (hrflag) call nc_err(nf90_inq_varid(ncid, hr_name, hr_varid))
267 if (read_grid_all) then
268 call nc_err(nf90_inq_varid(ncid, aname, avarid))
269 call nc_err(nf90_inq_varid(ncid, lname, lvarid))
270 call nc_err(nf90_inq_varid(ncid, mname, mvarid))
271 call nc_err(nf90_inq_varid(ncid, oname, ovarid))
274 call nc_err(nf90_get_var(ncid, xvarid, xvals))
275 call nc_err(nf90_get_var(ncid, yvarid, yvals))
276 call nc_err(nf90_get_var(ncid, zsoi_varid, zsoi_vals))
277 call nc_err(nf90_get_var(ncid, zlak_varid, zlak_vals))
278 if (hrflag) call nc_err(nf90_get_var(ncid, hr_varid, hr_vals))
279 if (read_grid_all) then
280 call nc_err(nf90_get_var(ncid, avarid, avals))
281 call nc_err(nf90_get_var(ncid, lvarid, lvals))
282 call nc_err(nf90_get_var(ncid, mvarid, mvals))
283 call nc_err(nf90_get_var(ncid, ovarid, ovals))
286 ! Read global attributes
287 call nc_err(nf90_get_att(ncid, nf90_global, case_id_name, case_id))
288 call nc_err(nf90_get_att(ncid, nf90_global, inidat_name, inidat))
289 call nc_err(nf90_get_att(ncid, nf90_global, surdat_name, surdat))
290 call nc_err(nf90_get_att(ncid, nf90_global, pftdat_name, pftdat))
291 call nc_err(nf90_get_att(ncid, nf90_global, rtmdat_name, rtmdat))
293 call nc_err(nf90_close(ncid))
296 ! Loop over all remaining files, reading size of time dimension
298 print *, 'Opening ',trim(input_path)//'/'//trim(fnames(i))
299 call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(i)), &
300 mode=nf90_nowrite, ncid=ncid))
301 call nc_err(nf90_inq_dimid(ncid, tname, tdimid))
302 call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=tdimid, len=tlen))
303 call nc_err(nf90_close(ncid))
304 print *, tlen,' values'
305 time_size = time_size + tlen
308 ! Allocate and fill the coordinates derived data type
311 coord%zsoi_size = zsoi_len
312 coord%zlak_size = zlak_len
313 coord%tsize = time_size
314 coord%hr_size = hr_len
315 coord%nsamples = nsamples_len
316 allocate(coord%x(xlen), coord%y(ylen), coord%zsoi(zsoi_len), coord%zlak(zlak_len), coord%hr(hr_len), stat=ierr)
318 print *, 'clm_read_coord: Cannot allocate grid variables'
321 if (read_grid_all) then
322 allocate(coord%area(xlen, ylen), coord%landfrac(xlen, ylen), coord%landmask(xlen, ylen), coord%orog(xlen, ylen), stat=ierr)
324 print *, 'clm_read_coord: Cannot allocate extra grid variables'
330 coord%zsoi = zsoi_vals
331 coord%zlak = zlak_vals
332 if (hrflag) coord%hr = hr_vals
333 if (read_grid_all) then
335 coord%landfrac = lvals
336 coord%landmask = mvals
339 deallocate(xvals, yvals, zsoi_vals, zlak_vals, hr_vals)
340 if (read_grid_all) then
341 deallocate(avals, lvals, mvals, ovals)
343 coord%case_id = trim(case_id)
344 coord%inidat = trim(inidat)
345 coord%surdat = trim(surdat)
346 coord%pftdat = trim(pftdat)
347 coord%rtmdat = trim(rtmdat)
348 coord%hour_flag = hrflag
350 end subroutine clm_read_coord
352 !--------------------------------------------------------------------------
354 subroutine clm_copy_grid_data(varname)
357 character (len=varnamelen), intent(in) :: varname ! variable name
360 var_data%varname = varname
362 var_data%missing_value = 1.0e+36
363 if (varname == 'landmask') then
364 var_data%int_type = .true.
365 allocate(var_data%int2d(coord%xsize, coord%ysize), stat=ierr)
367 print *, 'clm_copy_grid_data: Cannot allocate 2d integer variable'
371 var_data%int_type = .false.
372 allocate(var_data%var2d(coord%xsize, coord%ysize), stat=ierr)
374 print *, 'clm_copy_grid_data: Cannot allocate 2d real variable'
379 select case (varname)
381 var_data%var2d = coord%area
383 var_data%var2d = coord%landfrac
385 var_data%int2d = coord%landmask
387 var_data%var2d = coord%orog
389 print *,'clm_copy_grid_data: ',varname,' is not a known grid variable'
393 end subroutine clm_copy_grid_data
395 !--------------------------------------------------------------------------
397 subroutine clm_read_data(numf, input_path, fnames, tshift, varname)
400 integer, intent(in) :: numf ! number of filenames
401 character (len=namelen), intent(in) :: input_path ! input path
402 character (len=namelen), intent(in) :: fnames(numf) ! input filenames
403 integer, intent(in) :: tshift(numf) ! time shift by file
404 character (len=varnamelen), intent(in) :: varname ! variable name
406 integer :: tdimid, tvarid, tlen
407 integer :: hidimid, hilen
409 integer :: varid, dimids(nf90_max_var_dims), dimlen
411 integer(i8) :: i, j, k, l
412 integer(i8) :: time_size
413 real(r4), allocatable :: tvals(:) ! temporary time values
414 real(r4), allocatable :: tbvals(:,:) ! temporary time bounds values
415 real(r4), allocatable :: vals(:,:) ! temporary data values
416 character (len=varnamelen) :: dimname ! dimension name
419 ! Figure out the dimensionality of the variable and variable attributes
420 ! from the first file
421 print *, 'Opening ',trim(input_path)//'/'//trim(fnames(1)), &
422 ' to check dimensionality of ',trim(varname)
423 call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(1)), mode=nf90_nowrite, ncid=ncid))
424 call nc_err(nf90_inq_varid(ncid, trim(varname), varid))
425 call nc_err(nf90_inquire_variable(ncid, varid, name=var_data%varname, ndims=var_data%ndims, dimids=dimids))
426 ! Figure out if there is a Z-dimension and/or an hour-dimension
428 var_data%int_type = .false.
429 var_data%hour_flag = .false.
430 var_data%soil_layer_flag = .false.
431 do i = 1, var_data%ndims
432 call nc_err(nf90_inquire_dimension(ncid, dimids(i), name=dimname, len=dimlen))
433 if (i == 1 .and. dimname /= xname) ierr = ierr + 1
434 if (i == 2 .and. dimname /= yname) ierr = ierr + 1
435 if (i > 2 .and. i < var_data%ndims) then
436 if (dimname == hr_name) then
437 var_data%hour_flag = .true.
438 else if (dimname == zsoi_name) then
439 var_data%soil_layer_flag = .true.
440 !else if (dimname == zlak_name) then
441 ! var_data%zdim_name = dimname
446 ! This is not a good assumption (that time is always there), so skip it
447 !if (i == var_data%ndims .and. dimname /= tname) ierr = ierr + 1
448 print *,'Dimension ',i,' is ',dimname
451 print *,'clm_read_data: Unexpected dimension combination'
455 call nc_err(nf90_get_att(ncid, varid, long_name_name, var_data%long_name))
456 call nc_err(nf90_get_att(ncid, varid, missing_value_name, var_data%missing_value))
457 call nc_err(nf90_get_att(ncid, varid, units_name, var_data%units))
458 call nc_err(nf90_close(ncid))
460 ! Allocate space for the variable of the correct dimensionality
461 select case (var_data%ndims)
463 allocate(var_data%var2d(coord%xsize, coord%ysize), stat=ierr)
465 print *, 'clm_read_data: Cannot allocate 2d variable'
469 allocate(var_data%var3d(coord%xsize, coord%ysize, coord%tsize), stat=ierr)
471 print *, 'clm_read_data: Cannot allocate 3d variable'
475 if (var_data%hour_flag) then
476 zsize = coord%hr_size
477 else if (var_data%soil_layer_flag) then
478 zsize = coord%zsoi_size
479 !else if (var_data%zdim_name == zlak_name) then
480 ! zsize = coord%zlak_size
482 print *,'clm_read_data: Unexpected fourth dimension while allocating'
485 allocate(var_data%var4d(coord%xsize, coord%ysize, zsize, coord%tsize), stat=ierr)
487 print *, 'clm_read_data: Cannot allocate 4d variable'
491 print *, 'clm_read_data: Cannot handle variables of ', var_data%ndims, ' dimensions'
495 ! Allocate space for time stamps and bounds
496 allocate(var_data%time(coord%tsize), stat=ierr)
498 print *, 'clm_read_data: Cannot allocate time variable'
501 allocate(var_data%time_bounds(2,coord%tsize), stat=ierr)
503 print *, 'clm_read_data: Cannot allocate time bounds variable'
507 ! Allocate temporary space for the variable level-slice/time-slice */
508 allocate(vals(coord%xsize, coord%ysize), stat=ierr)
510 print *, 'clm_read_data: Cannot allocate variable time-slice'
514 if (var_data%ndims == 2) then
515 ! No time axis, so just read the first instance from the first file
516 print *, 'Opening ',trim(input_path)//'/'//trim(fnames(1))
517 call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(1)), mode=nf90_nowrite, ncid=ncid))
518 call nc_err(nf90_inq_varid(ncid, trim(varname), varid))
519 call nc_err(nf90_get_var(ncid, varid, vals, start = (/ 1, 1 /), count = (/ coord%xsize, coord%ysize /)))
520 var_data%var2d(:,:) = vals(:,:)
521 call nc_err(nf90_close(ncid))
523 ! Otherwise, loop over all times
525 ! Loop over all files, reading data
527 print *, 'Opening ',trim(input_path)//'/'//trim(fnames(i))
528 call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(i)), mode=nf90_nowrite, ncid=ncid))
529 call nc_err(nf90_inq_dimid(ncid, tname, tdimid))
530 call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=tdimid, len=tlen))
531 call nc_err(nf90_inq_dimid(ncid, hiname, hidimid))
532 call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=hidimid, len=hilen))
534 call nc_err(nf90_inq_varid(ncid, tname, tvarid))
535 call nc_err(nf90_inq_varid(ncid, tbname, tbvarid))
536 call nc_err(nf90_inq_varid(ncid, trim(varname), varid))
538 ! Allocate temporary space for the time stamps and bounds
539 allocate(tvals(tlen), stat=ierr)
541 print *, 'clm_read_data: Cannot allocate temporary time variable'
544 allocate(tbvals(hilen,tlen), stat=ierr)
546 print *, 'clm_read_data: Cannot allocate temporary time bounds variable'
550 ! Get the time stamps and bounds
551 call nc_err(nf90_get_var(ncid, tvarid, tvals))
552 call nc_err(nf90_get_var(ncid, tbvarid, tbvals))
555 ! Get and store a slab of data
556 select case (var_data%ndims)
558 !print *, 'Reading slice of data'
559 call nc_err(nf90_get_var(ncid, varid, vals, start = (/ 1, 1, j /), count = (/ coord%xsize, coord%ysize, 1 /)))
560 var_data%var3d(:,:,time_size + j) = vals(:,:)
563 call nc_err(nf90_get_var(ncid, varid, vals, start = (/ 1, 1, k, j /), count = (/ coord%xsize, coord%ysize, 1, 1 /)))
564 var_data%var4d(:,:,k,time_size + j) = vals(:,:)
567 print *, 'clm_read_data: Cannot handle variables of ', var_data%ndims, ' dimensions'
570 ! Store each time stamp and bounds, including the time offset
571 var_data%time(time_size + j) = tvals(j) + tshift(i)
572 var_data%time_bounds(:,time_size + j) = tbvals(:,j) + tshift(i)
576 call nc_err(nf90_close(ncid))
581 time_size = time_size + tlen
588 end subroutine clm_read_data
590 !--------------------------------------------------------------------------
592 subroutine clm_convert_data(casa_flux_bug)
595 logical, intent(in) :: casa_flux_bug
599 logical, allocatable :: mask2d(:,:)
600 logical, allocatable :: mask3d(:,:,:)
601 logical, allocatable :: mask4d(:,:,:,:)
603 ! Default offset is zero
606 ! Default var_data structure values
607 var_data%positive = ''
609 ! Set scale factor according to which variable is to be converted
610 select case (var_data%varname)
612 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
613 ! Atmospheric forcing
614 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
617 ! QBOT: kg kg-1 -> kg kg-1
619 var_data%units = 'kg kg-1'
620 var_data%out_varname = 'husf'
623 ! RAIN: mm s-1 -> kg m-2 s-1
625 var_data%units = 'kg m-2 s-1'
626 var_data%out_varname = 'prra'
629 ! SNOW: mm s-1 -> kg m-2 s-1
631 var_data%units = 'kg m-2 s-1'
632 var_data%out_varname = 'prsn'
634 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
636 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
641 ! AGNPP: g m-2 s-1 -> kg m-2 s-1
643 var_data%units = 'kg m-2 s-1'
644 var_data%out_varname = 'agnpp'
647 ! AR: g m-2 s-1 -> kg m-2 s-1
649 var_data%units = 'kg m-2 s-1'
650 var_data%out_varname = 'ar'
653 ! BIOGENCO: ug m-2 h-1 -> kg m-2 s-1
654 sfactor = 1.0e-9 * (1.0 / 3600.)
655 var_data%units = 'kg m-2 s-1'
656 var_data%out_varname = 'bco'
660 ! BGNPP: g m-2 s-1 -> kg m-2 s-1
662 var_data%units = 'kg m-2 s-1'
663 var_data%out_varname = 'bgnpp'
666 ! ISOPRENE: ug m-2 h-1 -> kg m-2 s-1
667 sfactor = 1.0e-9 * (1.0 / 3600.)
668 var_data%units = 'kg m-2 s-1'
669 var_data%out_varname = 'bipn'
672 ! MONOTERP: ug m-2 h-1 -> kg m-2 s-1
673 sfactor = 1.0e-9 * (1.0 / 3600.)
674 var_data%units = 'kg m-2 s-1'
675 var_data%out_varname = 'bmtpn'
678 ! ORVOC: ug m-2 h-1 -> kg m-2 s-1
679 sfactor = 1.0e-9 * (1.0 / 3600.)
680 var_data%units = 'kg m-2 s-1'
681 var_data%out_varname = 'borvoc'
684 ! OVOC: ug m-2 h-1 -> kg m-2 s-1
685 sfactor = 1.0e-9 * (1.0 / 3600.)
686 var_data%units = 'kg m-2 s-1'
687 var_data%out_varname = 'bovoc'
690 ! CWDC: g m-2 -> kg m-2
692 var_data%units = 'kg m-2'
693 var_data%out_varname = 'cwdc'
696 ! CWDC_HR: g m-2 s-1 -> kg m-2 s-1
698 var_data%units = 'kg m-2 s-1'
699 var_data%out_varname = 'cwdc_hr'
702 ! CWDC_LOSS: g m-2 s-1 -> kg m-2 s-1
704 if (casa_flux_bug) then
705 print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
706 sfactor = sfactor * (1. / 1200.)
708 var_data%units = 'kg m-2 s-1'
709 var_data%out_varname = 'cwdc_loss'
715 ! FROOTC: g m-2 -> kg m-2
717 var_data%units = 'kg m-2'
718 var_data%out_varname = 'frootc'
720 case ('FROOTC_ALLOC')
721 ! FROOTC_ALLOC: g m-2 s-1 -> kg m-2 s-1
723 if (casa_flux_bug) then
724 print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
725 sfactor = sfactor * (1. / 1200.)
727 var_data%units = 'kg m-2 s-1'
728 var_data%out_varname = 'frootc_alloc'
731 ! FROOTC_LOSS: g m-2 s-1 -> kg m-2 s-1
733 if (casa_flux_bug) then
734 print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
735 sfactor = sfactor * (1. / 1200.)
737 var_data%units = 'kg m-2 s-1'
738 var_data%out_varname = 'frootc_loss'
741 ! GPP: g m-2 s-1 -> kg m-2 s-1
743 var_data%units = 'kg m-2 s-1'
744 var_data%out_varname = 'gpp'
747 ! GROSS_NMIN: g m-2 s-1 -> kg m-2 s-1
749 var_data%units = 'kg m-2 s-1'
750 var_data%out_varname = 'gnmin'
753 ! HR: g m-2 s-1 -> kg m-2 s-1
755 var_data%units = 'kg m-2 s-1'
756 var_data%out_varname = 'hr'
762 ! LEAFC: g m-2 -> kg m-2
764 var_data%units = 'kg m-2'
765 var_data%out_varname = 'leafc'
768 ! LEAFC_ALLOC: g m-2 s-1 -> kg m-2 s-1
770 if (casa_flux_bug) then
771 print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
772 sfactor = sfactor * (1. / 1200.)
774 var_data%units = 'kg m-2 s-1'
775 var_data%out_varname = 'leafc_alloc'
778 ! LEAFC_LOSS: g m-2 s-1 -> kg m-2 s-1
780 if (casa_flux_bug) then
781 print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
782 sfactor = sfactor * (1. / 1200.)
784 var_data%units = 'kg m-2 s-1'
785 var_data%out_varname = 'leafc_loss'
788 ! LITTERC: g m-2 -> kg m-2
790 var_data%units = 'kg m-2'
791 var_data%out_varname = 'litterc'
794 ! LITTERC_HR: g m-2 s-1 -> kg m-2 s-1
796 var_data%units = 'kg m-2 s-1'
797 var_data%out_varname = 'litterc_hr'
799 case ('LITTERC_LOSS')
800 ! LITTERC_LOSS: g m-2 s-1 -> kg m-2 s-1
802 if (casa_flux_bug) then
803 print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
804 sfactor = sfactor * (1. / 1200.)
806 var_data%units = 'kg m-2 s-1'
807 var_data%out_varname = 'litterc_loss'
810 ! NEE: g m-2 s-1 -> kg m-2 s-1
812 var_data%units = 'kg m-2 s-1'
813 var_data%out_varname = 'nee'
816 ! NEP: g m-2 s-1 -> kg m-2 s-1
818 var_data%units = 'kg m-2 s-1'
819 var_data%out_varname = 'nep'
822 ! NEP: g m-2 s-1 -> kg m-2 s-1
824 var_data%units = 'kg m-2 s-1'
825 var_data%out_varname = 'nnmin'
828 ! NPP: g m-2 s-1 -> kg m-2 s-1
830 var_data%units = 'kg m-2 s-1'
831 var_data%out_varname = 'npp'
834 ! SOILC: g m-2 -> kg m-2
836 var_data%units = 'kg m-2'
837 var_data%out_varname = 'soilc'
840 ! SOILC_HR: g m-2 s-1 -> kg m-2 s-1
842 var_data%units = 'kg m-2 s-1'
843 var_data%out_varname = 'soilc_hr'
846 ! SOILC_LOSS: g m-2 s-1 -> kg m-2 s-1
848 if (casa_flux_bug) then
849 print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
850 sfactor = sfactor * (1. / 1200.)
852 var_data%units = 'kg m-2 s-1'
853 var_data%out_varname = 'soilc_loss'
856 ! VOCFLXT: ug m-2 h-1 -> kg m-2 s-1
857 sfactor = 1.0e-9 * (1.0 / 3600.)
858 var_data%units = 'kg m-2 s-1'
859 var_data%out_varname = 'tbvoc'
862 ! WOODC: g m-2 -> kg m-2
864 var_data%units = 'kg m-2'
865 var_data%out_varname = 'woodc'
868 ! WOODC_ALLOC: g m-2 s-1 -> kg m-2 s-1
870 if (casa_flux_bug) then
871 print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
872 sfactor = sfactor * (1. / 1200.)
874 var_data%units = 'kg m-2 s-1'
875 var_data%out_varname = 'woodc_alloc'
878 ! WOODC_LOSS: g m-2 s-1 -> kg m-2 s-1
880 if (casa_flux_bug) then
881 print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
882 sfactor = sfactor * (1. / 1200.)
884 var_data%units = 'kg m-2 s-1'
885 var_data%out_varname = 'woodc_loss'
887 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
889 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
892 ! FGR: W m-2 -> W m-2
894 var_data%units = 'W m-2'
895 var_data%positive = 'down'
896 var_data%out_varname = 'hfg'
899 ! LATENT: W m-2 -> W m-2
901 var_data%units = 'W m-2'
902 var_data%positive = 'up'
903 var_data%out_varname = 'hfls'
906 ! FSH: W m-2 -> W m-2
908 var_data%units = 'W m-2'
909 var_data%positive = 'up'
910 var_data%out_varname = 'hfss'
913 ! FSH_G: W m-2 -> W m-2
915 var_data%units = 'W m-2'
916 var_data%positive = 'up'
917 var_data%out_varname = 'hfssg'
920 ! FSH_V: W m-2 -> W m-2
922 var_data%units = 'W m-2'
923 var_data%positive = 'up'
924 var_data%out_varname = 'hfssv'
926 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
928 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
933 var_data%units = 'm2'
934 var_data%out_varname = 'gca'
939 var_data%out_varname = 'lbm'
940 ! Do not convert integer landmask, just return
949 var_data%out_varname = 'orog'
955 var_data%out_varname = 'sftlf'
957 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
959 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
965 var_data%out_varname = 'hur2m'
968 ! Q2M: kg kg-1 -> kg kg-1
970 var_data%units = 'kg kg-1'
971 var_data%out_varname = 'hus2m'
973 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
975 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
981 var_data%out_varname = 'btran'
984 ! ET: m s-1 -> kg m-2 s-1
986 var_data%units = 'kg m-2 s-1'
987 var_data%out_varname = 'et'
990 ! QSOIL: mm s-1 -> kg m-2 s-1
992 var_data%units = 'kg m-2 s-1'
993 var_data%out_varname = 'evspsblsoi'
996 ! QVEGE: mm s-1 -> kg m-2 s-1
998 var_data%units = 'kg m-2 s-1'
999 var_data%out_varname = 'evspsblveg'
1000 ! mrfsl (multi-level)
1002 ! SOILICE: kg m-2 -> kg m-2
1004 var_data%units = 'kg m-2'
1005 var_data%out_varname = 'mrfsl'
1007 case ('TOTAL_SOILICE')
1008 ! TOTAL_SOILICE: kg m-2 -> kg m-2
1010 var_data%units = 'kg m-2'
1011 var_data%out_varname = 'mrfso'
1012 ! mrlsl (multi-level)
1014 ! SOILLIQ: kg m-2 -> kg m-2
1016 var_data%units = 'kg m-2'
1017 var_data%out_varname = 'mrlsl'
1019 case ('TOTAL_SOILLIQ')
1020 ! TOTAL_SOILLIQ: kg m-2 -> kg m-2
1022 var_data%units = 'kg m-2'
1023 var_data%out_varname = 'mrlso'
1026 ! QOVER: mm s-1 -> kg m-2 s-1
1028 var_data%units = 'kg m-2 s-1'
1029 var_data%out_varname = 'mrros'
1032 ! QDRAI: mm s-1 -> kg m-2 s-1
1034 var_data%units = 'kg m-2 s-1'
1035 var_data%out_varname = 'mrross'
1038 ! QINTR: mm s-1 -> kg m-2 s-1
1040 var_data%units = 'kg m-2 s-1'
1041 var_data%out_varname = 'prveg'
1046 var_data%units = 'm'
1047 var_data%out_varname = 'snd'
1050 ! QVEGT: mm s-1 -> kg m-2 s-1
1052 var_data%units = 'kg m-2 s-1'
1053 var_data%out_varname = 'tran'
1054 ! vfwsl (multi-level)
1056 ! H2OSOI: mm3 mm-3 -> m3 m-3
1058 var_data%units = 'm3 m-3'
1059 var_data%out_varname = 'vfwsl'
1060 ! wpsl (multi-level)
1062 ! SOILPSI: MPa -> Pa
1064 var_data%units = 'Pa'
1065 var_data%out_varname = 'wpsl'
1067 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1069 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1071 case ('ALL_SKY_ALBEDO')
1072 ! ALL_SKY_ALBEDO: 1 -> 1
1074 var_data%units = '1'
1075 var_data%out_varname = 'alsas'
1077 case ('BLACK_SKY_ALBEDO')
1078 ! BLACK_SKY_ALBEDO: 1 -> 1
1080 var_data%units = '1'
1081 var_data%out_varname = 'alsbs'
1084 ! LAISHA: m2 m-2 -> m2 m-2
1086 var_data%units = 'm2 m-2'
1087 var_data%out_varname = 'laisha'
1090 ! LAISUN: m2 m-2 -> m2 m-2
1092 var_data%units = 'm2 m-2'
1093 var_data%out_varname = 'laisun'
1096 ! FLDS: W m-2 -> W m-2
1098 var_data%units = 'W m-2'
1099 var_data%positive = 'down'
1100 var_data%out_varname = 'rlds'
1103 ! FIRA: W m-2 -> W m-2
1105 var_data%units = 'W m-2'
1106 var_data%positive = 'down'
1107 var_data%out_varname = 'rls'
1110 ! FIRE: W m-2 -> W m-2
1112 var_data%units = 'W m-2'
1113 var_data%positive = 'up'
1114 var_data%out_varname = 'rlus'
1117 ! NETRAD: W m-2 -> W m-2
1119 var_data%units = 'W m-2'
1120 var_data%positive = 'down'
1121 var_data%out_varname = 'rs'
1124 ! FSA: W m-2 -> W m-2
1126 var_data%units = 'W m-2'
1127 var_data%positive = 'down'
1128 var_data%out_varname = 'rsas'
1131 ! SABG: W m-2 -> W m-2
1133 var_data%units = 'W m-2'
1134 var_data%positive = 'down'
1135 var_data%out_varname = 'rsasg'
1138 ! SABV: W m-2 -> W m-2
1140 var_data%units = 'W m-2'
1141 var_data%positive = 'down'
1142 var_data%out_varname = 'rsasv'
1145 ! FSDS: W m-2 -> W m-2
1147 var_data%units = 'W m-2'
1148 var_data%positive = 'down'
1149 var_data%out_varname = 'rsds'
1152 ! FSR: W m-2 -> W m-2
1154 var_data%units = 'W m-2'
1155 var_data%positive = 'up'
1156 var_data%out_varname = 'rsus'
1161 var_data%units = '1'
1162 var_data%out_varname = 'snc'
1164 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1166 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1172 var_data%units = 'K'
1173 var_data%out_varname = 'degday'
1178 var_data%units = 'K'
1179 var_data%out_varname = 'ta2m'
1184 var_data%units = 'K'
1185 var_data%out_varname = 'ta2mdmax'
1190 var_data%units = 'K'
1191 var_data%out_varname = 'ta2mdmin'
1196 var_data%units = 'K'
1197 var_data%out_varname = 'tsl'
1202 var_data%units = 'K'
1203 var_data%out_varname = 'tv'
1205 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1206 ! Vegetation phenology
1207 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1210 ! ELAI: m2 m-2 -> m2 m-2
1212 var_data%units = 'm2 m-2'
1213 var_data%out_varname = 'laiexp'
1216 ! TLAI: m2 m-2 -> m2 m-2
1218 var_data%units = 'm2 m-2'
1219 var_data%out_varname = 'laitot'
1222 ! ESAI: m2 m-2 -> m2 m-2
1224 var_data%units = 'm2 m-2'
1225 var_data%out_varname = 'saiexp'
1228 ! TSAI: m2 m-2 -> m2 m-2
1230 var_data%units = 'm2 m-2'
1231 var_data%out_varname = 'saitot'
1233 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1234 ! Vegetation physiology
1235 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1238 ! FPSN: umol m-2 s-1 -> kg m-2 s-1
1239 sfactor = 12.0 * 1.0e-6 * 1.0e-3
1240 var_data%units = 'kg m-2 s-1'
1241 var_data%out_varname = 'psn'
1244 ! RSSHA: s m-1 -> s m-1
1246 var_data%units = 's m-1'
1247 var_data%out_varname = 'rssha'
1250 ! RSSUN: s m-1 -> s m-1
1252 var_data%units = 's m-1'
1253 var_data%out_varname = 'rssun'
1256 print *, 'clm_convert_data: Unknown variable name ', var_data%varname
1260 ! Convert to the desired units using the offset and scale factor (sfactor)
1261 select case (var_data%ndims)
1263 allocate(mask2d(coord%xsize, coord%ysize), stat=ierr)
1265 print *, 'clm_convert_data: Cannot allocate 2d mask'
1269 if (var_data%int_type) then
1270 where (var_data%int2d /= var_data%missing_value)
1271 var_data%int2d = int(var_data%int2d * sfactor + offset)
1274 print *, var_data%out_varname, &
1275 ' min:', minval(var_data%int2d, mask=mask2d), &
1276 ' max:', maxval(var_data%int2d, mask=mask2d)
1278 where (var_data%var2d /= var_data%missing_value)
1279 var_data%var2d = var_data%var2d * sfactor + offset
1282 print *, var_data%out_varname, &
1283 ' min:', minval(var_data%var2d, mask=mask2d), &
1284 ' max:', maxval(var_data%var2d, mask=mask2d)
1286 if (allocated(mask2d)) deallocate(mask2d)
1288 allocate(mask3d(coord%xsize, coord%ysize, coord%tsize), stat=ierr)
1290 print *, 'clm_convert_data: Cannot allocate 3d mask'
1294 where (var_data%var3d /= var_data%missing_value)
1295 var_data%var3d = (var_data%var3d + offset) * sfactor
1298 print *, var_data%out_varname, &
1299 ' min:', minval(var_data%var3d, mask=mask3d), &
1300 ' max:', maxval(var_data%var3d, mask=mask3d)
1301 if (allocated(mask3d)) deallocate(mask3d)
1303 if (var_data%soil_layer_flag) then
1304 allocate(mask4d(coord%xsize, coord%ysize, coord%zsoi_size, coord%tsize), stat=ierr)
1306 print *, 'clm_convert_data: Cannot allocate 4d mask'
1310 where (var_data%var4d /= var_data%missing_value)
1311 var_data%var4d = (var_data%var4d + offset) * sfactor
1314 print *, var_data%out_varname, &
1315 ' min:', minval(var_data%var4d, mask=mask4d), &
1316 ' max:', maxval(var_data%var4d, mask=mask4d)
1317 if (allocated(mask4d)) deallocate(mask4d)
1319 print *,'clm_convert_data: Unknown data dimension'
1323 print *,'clm_convert_data: Dimension count not supported'
1327 end subroutine clm_convert_data
1329 !--------------------------------------------------------------------------
1331 subroutine clm_free_data
1335 var_data%varname = ''
1337 if (allocated(var_data%int2d)) deallocate(var_data%int2d)
1338 if (allocated(var_data%var2d)) deallocate(var_data%var2d)
1339 if (allocated(var_data%var3d)) deallocate(var_data%var3d)
1340 if (allocated(var_data%var4d)) deallocate(var_data%var4d)
1341 if (allocated(var_data%time)) deallocate(var_data%time)
1342 if (allocated(var_data%time_bounds)) deallocate(var_data%time_bounds)
1343 var_data%int_type = .false.
1344 var_data%soil_layer_flag = .false.
1345 var_data%hour_flag = .false.
1346 var_data%long_name = ''
1347 var_data%missing_value = 0.
1349 var_data%positive = ''
1350 var_data%out_varname = ''
1352 end subroutine clm_free_data
1354 !--------------------------------------------------------------------------
1356 subroutine nc_err(status)
1359 integer :: status ! status code
1361 if (status /= nf90_noerr) then
1362 print *, nf90_strerror(status)
1363 stop 'Program stopped after encountering netCDF error.'
1366 end subroutine nc_err
1368 !--------------------------------------------------------------------------