1.1 --- /dev/null Thu Jan 01 00:00:00 1970 +0000
1.2 +++ b/clm_mod.f90 Sun Sep 21 21:59:01 2008 -0400
1.3 @@ -0,0 +1,1371 @@
1.4 +module clm_mod
1.5 +
1.6 + use netcdf
1.7 + use kind_mod
1.8 + implicit none
1.9 +
1.10 + integer, parameter :: namelen = 128 ! string length for paths
1.11 + integer, parameter :: varnamelen = 32 ! string length for varnames
1.12 + integer, parameter :: posnamelen = 4 ! string length for positive attribute
1.13 + integer, parameter :: varattlen = 128 ! string length for variable attributes
1.14 + character (len=*), parameter :: xname = 'lon'
1.15 + character (len=*), parameter :: yname = 'lat'
1.16 + character (len=*), parameter :: zsoi_name = 'levsoi'
1.17 + character (len=*), parameter :: zlak_name = 'levlak'
1.18 + character (len=*), parameter :: tname = 'time'
1.19 + character (len=*), parameter :: hr_name = 'hour'
1.20 + character (len=*), parameter :: nsamples_name = 'nsamples'
1.21 + character (len=*), parameter :: hiname = 'hist_interval'
1.22 + character (len=*), parameter :: tbname = 'time_bounds'
1.23 + character (len=*), parameter :: aname = 'area'
1.24 + character (len=*), parameter :: lname = 'landfrac'
1.25 + character (len=*), parameter :: mname = 'landmask'
1.26 + character (len=*), parameter :: oname = 'topo'
1.27 + character (len=*), parameter :: long_name_name = 'long_name'
1.28 + character (len=*), parameter :: missing_value_name = 'missing_value'
1.29 + character (len=*), parameter :: units_name = 'units'
1.30 + character (len=*), parameter :: case_id_name = 'case_id'
1.31 + character (len=*), parameter :: inidat_name = 'Initial_conditions_dataset'
1.32 + character (len=*), parameter :: surdat_name = 'Surface_dataset'
1.33 + character (len=*), parameter :: pftdat_name = 'PFT_physiological_constants_dataset'
1.34 + character (len=*), parameter :: rtmdat_name = 'RTM_input_datset'
1.35 +
1.36 + type surface
1.37 + integer :: xsize ! number of longitude points (lon)
1.38 + integer :: ysize ! number of latitutde points (lat)
1.39 + integer :: psize ! number of pfts
1.40 + real(r8), pointer :: x(:) ! longitude values (lon)
1.41 + real(r8), pointer :: y(:) ! latitude values (lat)
1.42 + real(r8), pointer :: p(:) ! pft number
1.43 + real(r4), pointer :: pft(:,:,:) ! percent of each pft
1.44 + end type surface
1.45 +
1.46 + type coordinates
1.47 + integer :: xsize ! number of longitude points (lon)
1.48 + integer :: ysize ! number of latitutde points (lat)
1.49 + integer :: zsoi_size ! number of soil z-levels (levsoi)
1.50 + integer :: zlak_size ! number of lake z-levels (levlak)
1.51 + integer :: tsize ! number of time points (time)
1.52 + integer :: hr_size ! number of hour points (hour = 24)
1.53 + integer :: nsamples ! number of time samples in post-processing statistics
1.54 + real(r8), pointer :: x(:) ! longitude values (lon)
1.55 + real(r8), pointer :: y(:) ! latitude values (lat)
1.56 + real(r8), pointer :: zsoi(:) ! soil z-level values (levsoi)
1.57 + real(r8), pointer :: zlak(:) ! lake z-level values (levlak)
1.58 + real(r8), pointer :: hr(:) ! hour-of-day values (hour)
1.59 + real(r8), pointer :: area(:,:) ! gridcell areas
1.60 + real(r8), pointer :: landfrac(:,:) ! gridcell land fractions
1.61 + integer, pointer :: landmask(:,:) ! gridcell land/ocean mask
1.62 + real(r8), pointer :: orog(:,:) ! gridcell orography/topography
1.63 + character (len=namelen) :: case_id ! case name of run
1.64 + character (len=namelen) :: inidat ! initial dataset filename
1.65 + character (len=namelen) :: surdat ! surface dataset filename
1.66 + character (len=namelen) :: pftdat ! pft dataset filename
1.67 + character (len=namelen) :: rtmdat ! runoff dataset filename
1.68 + logical :: hour_flag ! flag for hour-of-day statistical summary files
1.69 + end type coordinates
1.70 +
1.71 + type variable_data
1.72 + character (len=varnamelen) :: varname
1.73 + integer :: ndims
1.74 + integer, pointer :: int2d(:,:)
1.75 + real(r4), pointer :: var2d(:,:)
1.76 + real(r4), pointer :: var3d(:,:,:)
1.77 + real(r4), pointer :: var4d(:,:,:,:)
1.78 + real(r8), pointer :: time(:)
1.79 + real(r8), pointer :: time_bounds(:,:)
1.80 + logical :: int_type ! set if type is integer
1.81 + logical :: soil_layer_flag
1.82 + logical :: hour_flag ! flag for hour-of-day statistical summary variable
1.83 + character (len=varattlen) :: long_name
1.84 + real(r4) :: missing_value
1.85 + character (len=varattlen) :: units
1.86 + character (len=posnamelen) :: positive
1.87 + character (len=varnamelen) :: out_varname
1.88 + end type variable_data
1.89 +
1.90 + type(surface), target :: surf
1.91 + type(coordinates), target :: coord
1.92 + type(variable_data), target :: var_data
1.93 + save
1.94 +
1.95 + public clm_read_surf
1.96 + public clm_free_surf
1.97 + public clm_read_coord
1.98 + public clm_copy_grid_data
1.99 + public clm_read_data
1.100 + public clm_convert_data
1.101 + public clm_free_data
1.102 + private nc_err
1.103 +
1.104 + contains
1.105 +
1.106 + !
1.107 + !--------------------------------------------------------------------------
1.108 + !
1.109 + subroutine clm_read_surf(fsurdat)
1.110 +
1.111 + implicit none
1.112 + character (len=namelen), intent(in) :: fsurdat ! input surface dataset
1.113 + integer :: i
1.114 + integer :: ierr
1.115 + integer :: ncid
1.116 + integer :: xdimid, xvarid, xlen
1.117 + integer :: ydimid, yvarid, ylen
1.118 + integer :: pdimid, pvarid, plen
1.119 + real(r8), allocatable :: longxy(:,:)
1.120 + real(r8), allocatable :: latixy(:,:)
1.121 + real(r4), allocatable :: pctpft(:,:,:)
1.122 +
1.123 + ! Read coordinated from surface dataset
1.124 + print *, 'Opening ',trim(fsurdat)
1.125 + call nc_err(nf90_open(path=trim(fsurdat), mode=nf90_nowrite, ncid=ncid))
1.126 + call nc_err(nf90_inq_dimid(ncid, 'lsmlon', xdimid))
1.127 + call nc_err(nf90_inq_dimid(ncid, 'lsmlat', ydimid))
1.128 + call nc_err(nf90_inq_dimid(ncid, 'lsmpft', pdimid))
1.129 +
1.130 + call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=xdimid, len=xlen))
1.131 + call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=ydimid, len=ylen))
1.132 + call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=pdimid, len=plen))
1.133 +
1.134 + call nc_err(nf90_inq_varid(ncid, 'LONGXY', xvarid))
1.135 + call nc_err(nf90_inq_varid(ncid, 'LATIXY', yvarid))
1.136 + call nc_err(nf90_inq_varid(ncid, 'PCT_PFT', pvarid))
1.137 +
1.138 + allocate(longxy(xlen, ylen), latixy(xlen, ylen), pctpft(xlen, ylen, plen), stat=ierr)
1.139 + if (ierr /= 0) then
1.140 + print *, 'clm_read_surf: Cannot allocate longxy, latixy, pctpft variables'
1.141 + stop
1.142 + end if
1.143 +
1.144 + call nc_err(nf90_get_var(ncid, xvarid, longxy))
1.145 + call nc_err(nf90_get_var(ncid, yvarid, latixy))
1.146 + call nc_err(nf90_get_var(ncid, pvarid, pctpft))
1.147 +
1.148 + print *, 'Closing ',trim(fsurdat)
1.149 + call nc_err(nf90_close(ncid))
1.150 +
1.151 + ! Allocate and fill in the surface derived data type
1.152 + surf%xsize = xlen
1.153 + surf%ysize = ylen
1.154 + surf%psize = plen
1.155 + allocate(surf%x(xlen), surf%y(ylen), surf%p(plen), surf%pft(xlen, ylen, plen), stat=ierr)
1.156 + if (ierr /= 0) then
1.157 + print *, 'clm_read_surf: Cannot allocate surface type variables'
1.158 + stop
1.159 + end if
1.160 + surf%x(:) = longxy(:,1)
1.161 + surf%y(:) = latixy(1,:)
1.162 + do i = 1, plen
1.163 + surf%p(i) = real(i)
1.164 + end do
1.165 + surf%pft = pctpft
1.166 +
1.167 + deallocate(longxy, latixy, pctpft)
1.168 +
1.169 + end subroutine clm_read_surf
1.170 + !
1.171 + !--------------------------------------------------------------------------
1.172 + !
1.173 + subroutine clm_free_surf()
1.174 +
1.175 + implicit none
1.176 +
1.177 + if (allocated(surf%x)) deallocate(surf%x)
1.178 + if (allocated(surf%y)) deallocate(surf%y)
1.179 + if (allocated(surf%p)) deallocate(surf%p)
1.180 + if (allocated(surf%pft)) deallocate(surf%pft)
1.181 +
1.182 + end subroutine clm_free_surf
1.183 + !
1.184 + !--------------------------------------------------------------------------
1.185 + !
1.186 + subroutine clm_read_coord(numf, input_path, fnames, read_grid_all)
1.187 +
1.188 + implicit none
1.189 + integer, intent(in) :: numf ! number of filenames
1.190 + character (len=namelen), intent(in) :: input_path ! input path
1.191 + character (len=namelen), intent(in) :: fnames(numf) ! input filenames
1.192 + logical, intent(in) :: read_grid_all ! flag to read area,...
1.193 + integer :: ncid
1.194 + integer :: xdimid, xvarid, xlen
1.195 + integer :: ydimid, yvarid, ylen
1.196 + integer :: zsoi_dimid, zsoi_varid, zsoi_len
1.197 + integer :: zlak_dimid, zlak_varid, zlak_len
1.198 + integer :: tdimid, tvarid, tlen
1.199 + integer :: hr_dimid, hr_varid, hr_len
1.200 + integer :: nsamples_dimid, nsamples_len
1.201 + integer :: avarid ! area (when read_grid_all true)
1.202 + integer :: lvarid ! landfrac (when read_grid_all true)
1.203 + integer :: mvarid ! landmask (when read_grid_all true)
1.204 + integer :: ovarid ! orog/topo (when read_grid_all true)
1.205 + integer(i8) :: i
1.206 + integer(i8) :: time_size
1.207 + real(r4), allocatable :: xvals(:)
1.208 + real(r4), allocatable :: yvals(:)
1.209 + real(r4), allocatable :: zsoi_vals(:)
1.210 + real(r4), allocatable :: zlak_vals(:)
1.211 + real(r4), allocatable :: hr_vals(:)
1.212 + real(r4), allocatable :: avals(:,:) ! area (when read_grid_all true)
1.213 + real(r4), allocatable :: lvals(:,:) ! landfrac (when read_grid_all true)
1.214 + integer, allocatable :: mvals(:,:) ! landmask (when read_grid_all true)
1.215 + real(r4), allocatable :: ovals(:,:) ! orog/topo (when read_grid_all true)
1.216 + character (len=namelen) :: case_id = ''
1.217 + character (len=namelen) :: inidat = ''
1.218 + character (len=namelen) :: surdat = ''
1.219 + character (len=namelen) :: pftdat = ''
1.220 + character (len=namelen) :: rtmdat = ''
1.221 + logical :: hrflag = .false.
1.222 + integer :: ierr
1.223 +
1.224 + ! Read coordinates from first file
1.225 + print *, 'Opening ',trim(input_path)//'/'//trim(fnames(1))
1.226 + call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(1)), mode=nf90_nowrite, ncid=ncid))
1.227 + call nc_err(nf90_inq_dimid(ncid, xname, xdimid))
1.228 + call nc_err(nf90_inq_dimid(ncid, yname, ydimid))
1.229 + call nc_err(nf90_inq_dimid(ncid, zsoi_name, zsoi_dimid))
1.230 + call nc_err(nf90_inq_dimid(ncid, zlak_name, zlak_dimid))
1.231 + call nc_err(nf90_inq_dimid(ncid, tname, tdimid))
1.232 + if (nf90_inq_dimid(ncid, hr_name, hr_dimid) == nf90_noerr) then
1.233 + call nc_err(nf90_inq_dimid(ncid, nsamples_name, nsamples_dimid))
1.234 + hrflag = .true.
1.235 + end if
1.236 +
1.237 + call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=xdimid, len=xlen))
1.238 + call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=ydimid, len=ylen))
1.239 + call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=zsoi_dimid, len=zsoi_len))
1.240 + call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=zlak_dimid, len=zlak_len))
1.241 + call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=tdimid, len=tlen))
1.242 + if (hrflag) then
1.243 + call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=hr_dimid, len=hr_len))
1.244 + call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=nsamples_dimid, len=nsamples_len))
1.245 + else
1.246 + hr_len = 0
1.247 + nsamples_len = 0
1.248 + end if
1.249 + print *, tlen,' values'
1.250 +
1.251 + allocate(xvals(xlen), yvals(ylen), zsoi_vals(zsoi_len), zlak_vals(zlak_len), hr_vals(hr_len), stat=ierr)
1.252 + if (ierr /= 0) then
1.253 + print *, 'clm_read_coord: Cannot allocate grid variables'
1.254 + stop
1.255 + end if
1.256 + if (read_grid_all) then
1.257 + allocate(avals(xlen, ylen), lvals(xlen, ylen), mvals(xlen, ylen), ovals(xlen, ylen), stat=ierr)
1.258 + if (ierr /= 0) then
1.259 + print *, 'clm_read_coord: Cannot allocate extra grid variables'
1.260 + stop
1.261 + end if
1.262 + end if
1.263 +
1.264 + call nc_err(nf90_inq_varid(ncid, xname, xvarid))
1.265 + call nc_err(nf90_inq_varid(ncid, yname, yvarid))
1.266 + call nc_err(nf90_inq_varid(ncid, zsoi_name, zsoi_varid))
1.267 + call nc_err(nf90_inq_varid(ncid, zlak_name, zlak_varid))
1.268 + call nc_err(nf90_inq_varid(ncid, tname, tvarid))
1.269 + if (hrflag) call nc_err(nf90_inq_varid(ncid, hr_name, hr_varid))
1.270 + if (read_grid_all) then
1.271 + call nc_err(nf90_inq_varid(ncid, aname, avarid))
1.272 + call nc_err(nf90_inq_varid(ncid, lname, lvarid))
1.273 + call nc_err(nf90_inq_varid(ncid, mname, mvarid))
1.274 + call nc_err(nf90_inq_varid(ncid, oname, ovarid))
1.275 + end if
1.276 +
1.277 + call nc_err(nf90_get_var(ncid, xvarid, xvals))
1.278 + call nc_err(nf90_get_var(ncid, yvarid, yvals))
1.279 + call nc_err(nf90_get_var(ncid, zsoi_varid, zsoi_vals))
1.280 + call nc_err(nf90_get_var(ncid, zlak_varid, zlak_vals))
1.281 + if (hrflag) call nc_err(nf90_get_var(ncid, hr_varid, hr_vals))
1.282 + if (read_grid_all) then
1.283 + call nc_err(nf90_get_var(ncid, avarid, avals))
1.284 + call nc_err(nf90_get_var(ncid, lvarid, lvals))
1.285 + call nc_err(nf90_get_var(ncid, mvarid, mvals))
1.286 + call nc_err(nf90_get_var(ncid, ovarid, ovals))
1.287 + end if
1.288 +
1.289 + ! Read global attributes
1.290 + call nc_err(nf90_get_att(ncid, nf90_global, case_id_name, case_id))
1.291 + call nc_err(nf90_get_att(ncid, nf90_global, inidat_name, inidat))
1.292 + call nc_err(nf90_get_att(ncid, nf90_global, surdat_name, surdat))
1.293 + call nc_err(nf90_get_att(ncid, nf90_global, pftdat_name, pftdat))
1.294 + call nc_err(nf90_get_att(ncid, nf90_global, rtmdat_name, rtmdat))
1.295 +
1.296 + call nc_err(nf90_close(ncid))
1.297 +
1.298 + time_size = tlen
1.299 + ! Loop over all remaining files, reading size of time dimension
1.300 + do i = 2, numf
1.301 + print *, 'Opening ',trim(input_path)//'/'//trim(fnames(i))
1.302 + call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(i)), &
1.303 + mode=nf90_nowrite, ncid=ncid))
1.304 + call nc_err(nf90_inq_dimid(ncid, tname, tdimid))
1.305 + call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=tdimid, len=tlen))
1.306 + call nc_err(nf90_close(ncid))
1.307 + print *, tlen,' values'
1.308 + time_size = time_size + tlen
1.309 + end do
1.310 +
1.311 + ! Allocate and fill the coordinates derived data type
1.312 + coord%xsize = xlen
1.313 + coord%ysize = ylen
1.314 + coord%zsoi_size = zsoi_len
1.315 + coord%zlak_size = zlak_len
1.316 + coord%tsize = time_size
1.317 + coord%hr_size = hr_len
1.318 + coord%nsamples = nsamples_len
1.319 + allocate(coord%x(xlen), coord%y(ylen), coord%zsoi(zsoi_len), coord%zlak(zlak_len), coord%hr(hr_len), stat=ierr)
1.320 + if (ierr /= 0) then
1.321 + print *, 'clm_read_coord: Cannot allocate grid variables'
1.322 + stop
1.323 + end if
1.324 + if (read_grid_all) then
1.325 + allocate(coord%area(xlen, ylen), coord%landfrac(xlen, ylen), coord%landmask(xlen, ylen), coord%orog(xlen, ylen), stat=ierr)
1.326 + if (ierr /= 0) then
1.327 + print *, 'clm_read_coord: Cannot allocate extra grid variables'
1.328 + stop
1.329 + end if
1.330 + end if
1.331 + coord%x = xvals
1.332 + coord%y = yvals
1.333 + coord%zsoi = zsoi_vals
1.334 + coord%zlak = zlak_vals
1.335 + if (hrflag) coord%hr = hr_vals
1.336 + if (read_grid_all) then
1.337 + coord%area = avals
1.338 + coord%landfrac = lvals
1.339 + coord%landmask = mvals
1.340 + coord%orog = ovals
1.341 + end if
1.342 + deallocate(xvals, yvals, zsoi_vals, zlak_vals, hr_vals)
1.343 + if (read_grid_all) then
1.344 + deallocate(avals, lvals, mvals, ovals)
1.345 + end if
1.346 + coord%case_id = trim(case_id)
1.347 + coord%inidat = trim(inidat)
1.348 + coord%surdat = trim(surdat)
1.349 + coord%pftdat = trim(pftdat)
1.350 + coord%rtmdat = trim(rtmdat)
1.351 + coord%hour_flag = hrflag
1.352 +
1.353 + end subroutine clm_read_coord
1.354 + !
1.355 + !--------------------------------------------------------------------------
1.356 + !
1.357 + subroutine clm_copy_grid_data(varname)
1.358 +
1.359 + implicit none
1.360 + character (len=varnamelen), intent(in) :: varname ! variable name
1.361 + integer :: ierr
1.362 +
1.363 + var_data%varname = varname
1.364 + var_data%ndims = 2
1.365 + var_data%missing_value = 1.0e+36
1.366 + if (varname == 'landmask') then
1.367 + var_data%int_type = .true.
1.368 + allocate(var_data%int2d(coord%xsize, coord%ysize), stat=ierr)
1.369 + if (ierr /= 0) then
1.370 + print *, 'clm_copy_grid_data: Cannot allocate 2d integer variable'
1.371 + stop
1.372 + end if
1.373 + else
1.374 + var_data%int_type = .false.
1.375 + allocate(var_data%var2d(coord%xsize, coord%ysize), stat=ierr)
1.376 + if (ierr /= 0) then
1.377 + print *, 'clm_copy_grid_data: Cannot allocate 2d real variable'
1.378 + stop
1.379 + end if
1.380 + end if
1.381 +
1.382 + select case (varname)
1.383 + case ('area')
1.384 + var_data%var2d = coord%area
1.385 + case ('landfrac')
1.386 + var_data%var2d = coord%landfrac
1.387 + case ('landmask')
1.388 + var_data%int2d = coord%landmask
1.389 + case ('topo')
1.390 + var_data%var2d = coord%orog
1.391 + case default
1.392 + print *,'clm_copy_grid_data: ',varname,' is not a known grid variable'
1.393 + stop
1.394 + end select
1.395 +
1.396 + end subroutine clm_copy_grid_data
1.397 + !
1.398 + !--------------------------------------------------------------------------
1.399 + !
1.400 + subroutine clm_read_data(numf, input_path, fnames, tshift, varname)
1.401 +
1.402 + implicit none
1.403 + integer, intent(in) :: numf ! number of filenames
1.404 + character (len=namelen), intent(in) :: input_path ! input path
1.405 + character (len=namelen), intent(in) :: fnames(numf) ! input filenames
1.406 + integer, intent(in) :: tshift(numf) ! time shift by file
1.407 + character (len=varnamelen), intent(in) :: varname ! variable name
1.408 + integer :: ncid
1.409 + integer :: tdimid, tvarid, tlen
1.410 + integer :: hidimid, hilen
1.411 + integer :: tbvarid
1.412 + integer :: varid, dimids(nf90_max_var_dims), dimlen
1.413 + integer :: zsize
1.414 + integer(i8) :: i, j, k, l
1.415 + integer(i8) :: time_size
1.416 + real(r4), allocatable :: tvals(:) ! temporary time values
1.417 + real(r4), allocatable :: tbvals(:,:) ! temporary time bounds values
1.418 + real(r4), allocatable :: vals(:,:) ! temporary data values
1.419 + character (len=varnamelen) :: dimname ! dimension name
1.420 + integer :: ierr
1.421 +
1.422 + ! Figure out the dimensionality of the variable and variable attributes
1.423 + ! from the first file
1.424 + print *, 'Opening ',trim(input_path)//'/'//trim(fnames(1)), &
1.425 + ' to check dimensionality of ',trim(varname)
1.426 + call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(1)), mode=nf90_nowrite, ncid=ncid))
1.427 + call nc_err(nf90_inq_varid(ncid, trim(varname), varid))
1.428 + call nc_err(nf90_inquire_variable(ncid, varid, name=var_data%varname, ndims=var_data%ndims, dimids=dimids))
1.429 + ! Figure out if there is a Z-dimension and/or an hour-dimension
1.430 + ierr = 0
1.431 + var_data%int_type = .false.
1.432 + var_data%hour_flag = .false.
1.433 + var_data%soil_layer_flag = .false.
1.434 + do i = 1, var_data%ndims
1.435 + call nc_err(nf90_inquire_dimension(ncid, dimids(i), name=dimname, len=dimlen))
1.436 + if (i == 1 .and. dimname /= xname) ierr = ierr + 1
1.437 + if (i == 2 .and. dimname /= yname) ierr = ierr + 1
1.438 + if (i > 2 .and. i < var_data%ndims) then
1.439 + if (dimname == hr_name) then
1.440 + var_data%hour_flag = .true.
1.441 + else if (dimname == zsoi_name) then
1.442 + var_data%soil_layer_flag = .true.
1.443 + !else if (dimname == zlak_name) then
1.444 + ! var_data%zdim_name = dimname
1.445 + else
1.446 + ierr = ierr + 1
1.447 + end if
1.448 + end if
1.449 + ! This is not a good assumption (that time is always there), so skip it
1.450 + !if (i == var_data%ndims .and. dimname /= tname) ierr = ierr + 1
1.451 + print *,'Dimension ',i,' is ',dimname
1.452 + end do
1.453 + if (ierr /= 0) then
1.454 + print *,'clm_read_data: Unexpected dimension combination'
1.455 + stop
1.456 + end if
1.457 + !
1.458 + call nc_err(nf90_get_att(ncid, varid, long_name_name, var_data%long_name))
1.459 + call nc_err(nf90_get_att(ncid, varid, missing_value_name, var_data%missing_value))
1.460 + call nc_err(nf90_get_att(ncid, varid, units_name, var_data%units))
1.461 + call nc_err(nf90_close(ncid))
1.462 +
1.463 + ! Allocate space for the variable of the correct dimensionality
1.464 + select case (var_data%ndims)
1.465 + case (2)
1.466 + allocate(var_data%var2d(coord%xsize, coord%ysize), stat=ierr)
1.467 + if (ierr /= 0) then
1.468 + print *, 'clm_read_data: Cannot allocate 2d variable'
1.469 + stop
1.470 + end if
1.471 + case (3)
1.472 + allocate(var_data%var3d(coord%xsize, coord%ysize, coord%tsize), stat=ierr)
1.473 + if (ierr /= 0) then
1.474 + print *, 'clm_read_data: Cannot allocate 3d variable'
1.475 + stop
1.476 + end if
1.477 + case (4)
1.478 + if (var_data%hour_flag) then
1.479 + zsize = coord%hr_size
1.480 + else if (var_data%soil_layer_flag) then
1.481 + zsize = coord%zsoi_size
1.482 + !else if (var_data%zdim_name == zlak_name) then
1.483 + ! zsize = coord%zlak_size
1.484 + else
1.485 + print *,'clm_read_data: Unexpected fourth dimension while allocating'
1.486 + stop
1.487 + end if
1.488 + allocate(var_data%var4d(coord%xsize, coord%ysize, zsize, coord%tsize), stat=ierr)
1.489 + if (ierr /= 0) then
1.490 + print *, 'clm_read_data: Cannot allocate 4d variable'
1.491 + stop
1.492 + end if
1.493 + case default
1.494 + print *, 'clm_read_data: Cannot handle variables of ', var_data%ndims, ' dimensions'
1.495 + stop
1.496 + end select
1.497 +
1.498 + ! Allocate space for time stamps and bounds
1.499 + allocate(var_data%time(coord%tsize), stat=ierr)
1.500 + if (ierr /= 0) then
1.501 + print *, 'clm_read_data: Cannot allocate time variable'
1.502 + stop
1.503 + end if
1.504 + allocate(var_data%time_bounds(2,coord%tsize), stat=ierr)
1.505 + if (ierr /= 0) then
1.506 + print *, 'clm_read_data: Cannot allocate time bounds variable'
1.507 + stop
1.508 + end if
1.509 +
1.510 + ! Allocate temporary space for the variable level-slice/time-slice */
1.511 + allocate(vals(coord%xsize, coord%ysize), stat=ierr)
1.512 + if (ierr /= 0) then
1.513 + print *, 'clm_read_data: Cannot allocate variable time-slice'
1.514 + stop
1.515 + end if
1.516 +
1.517 + if (var_data%ndims == 2) then
1.518 + ! No time axis, so just read the first instance from the first file
1.519 + print *, 'Opening ',trim(input_path)//'/'//trim(fnames(1))
1.520 + call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(1)), mode=nf90_nowrite, ncid=ncid))
1.521 + call nc_err(nf90_inq_varid(ncid, trim(varname), varid))
1.522 + call nc_err(nf90_get_var(ncid, varid, vals, start = (/ 1, 1 /), count = (/ coord%xsize, coord%ysize /)))
1.523 + var_data%var2d(:,:) = vals(:,:)
1.524 + call nc_err(nf90_close(ncid))
1.525 + else
1.526 + ! Otherwise, loop over all times
1.527 + time_size = 0
1.528 + ! Loop over all files, reading data
1.529 + do i = 1, numf
1.530 + print *, 'Opening ',trim(input_path)//'/'//trim(fnames(i))
1.531 + call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(i)), mode=nf90_nowrite, ncid=ncid))
1.532 + call nc_err(nf90_inq_dimid(ncid, tname, tdimid))
1.533 + call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=tdimid, len=tlen))
1.534 + call nc_err(nf90_inq_dimid(ncid, hiname, hidimid))
1.535 + call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=hidimid, len=hilen))
1.536 +
1.537 + call nc_err(nf90_inq_varid(ncid, tname, tvarid))
1.538 + call nc_err(nf90_inq_varid(ncid, tbname, tbvarid))
1.539 + call nc_err(nf90_inq_varid(ncid, trim(varname), varid))
1.540 +
1.541 + ! Allocate temporary space for the time stamps and bounds
1.542 + allocate(tvals(tlen), stat=ierr)
1.543 + if (ierr /= 0) then
1.544 + print *, 'clm_read_data: Cannot allocate temporary time variable'
1.545 + stop
1.546 + end if
1.547 + allocate(tbvals(hilen,tlen), stat=ierr)
1.548 + if (ierr /= 0) then
1.549 + print *, 'clm_read_data: Cannot allocate temporary time bounds variable'
1.550 + stop
1.551 + end if
1.552 +
1.553 + ! Get the time stamps and bounds
1.554 + call nc_err(nf90_get_var(ncid, tvarid, tvals))
1.555 + call nc_err(nf90_get_var(ncid, tbvarid, tbvals))
1.556 +
1.557 + do j = 1, tlen
1.558 + ! Get and store a slab of data
1.559 + select case (var_data%ndims)
1.560 + case (3)
1.561 + !print *, 'Reading slice of data'
1.562 + call nc_err(nf90_get_var(ncid, varid, vals, start = (/ 1, 1, j /), count = (/ coord%xsize, coord%ysize, 1 /)))
1.563 + var_data%var3d(:,:,time_size + j) = vals(:,:)
1.564 + case (4)
1.565 + do k = 1, zsize
1.566 + call nc_err(nf90_get_var(ncid, varid, vals, start = (/ 1, 1, k, j /), count = (/ coord%xsize, coord%ysize, 1, 1 /)))
1.567 + var_data%var4d(:,:,k,time_size + j) = vals(:,:)
1.568 + end do
1.569 + case default
1.570 + print *, 'clm_read_data: Cannot handle variables of ', var_data%ndims, ' dimensions'
1.571 + stop
1.572 + end select
1.573 + ! Store each time stamp and bounds, including the time offset
1.574 + var_data%time(time_size + j) = tvals(j) + tshift(i)
1.575 + var_data%time_bounds(:,time_size + j) = tbvals(:,j) + tshift(i)
1.576 +
1.577 + end do
1.578 +
1.579 + call nc_err(nf90_close(ncid))
1.580 +
1.581 + deallocate(tvals)
1.582 + deallocate(tbvals)
1.583 +
1.584 + time_size = time_size + tlen
1.585 +
1.586 + end do
1.587 + end if
1.588 +
1.589 + deallocate(vals)
1.590 +
1.591 + end subroutine clm_read_data
1.592 + !
1.593 + !--------------------------------------------------------------------------
1.594 + !
1.595 + subroutine clm_convert_data(casa_flux_bug)
1.596 +
1.597 + implicit none
1.598 + logical, intent(in) :: casa_flux_bug
1.599 + real(r8) :: sfactor
1.600 + real(r8) :: offset
1.601 + integer :: ierr
1.602 + logical, allocatable :: mask2d(:,:)
1.603 + logical, allocatable :: mask3d(:,:,:)
1.604 + logical, allocatable :: mask4d(:,:,:,:)
1.605 +
1.606 + ! Default offset is zero
1.607 + offset = 0.
1.608 +
1.609 + ! Default var_data structure values
1.610 + var_data%positive = ''
1.611 +
1.612 + ! Set scale factor according to which variable is to be converted
1.613 + select case (var_data%varname)
1.614 +
1.615 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.616 + ! Atmospheric forcing
1.617 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.618 + ! husf
1.619 + case ('QBOT')
1.620 + ! QBOT: kg kg-1 -> kg kg-1
1.621 + sfactor = 1.0
1.622 + var_data%units = 'kg kg-1'
1.623 + var_data%out_varname = 'husf'
1.624 + ! prra
1.625 + case ('RAIN')
1.626 + ! RAIN: mm s-1 -> kg m-2 s-1
1.627 + sfactor = 1.0
1.628 + var_data%units = 'kg m-2 s-1'
1.629 + var_data%out_varname = 'prra'
1.630 + ! prsn
1.631 + case ('SNOW')
1.632 + ! SNOW: mm s-1 -> kg m-2 s-1
1.633 + sfactor = 1.0
1.634 + var_data%units = 'kg m-2 s-1'
1.635 + var_data%out_varname = 'prsn'
1.636 +
1.637 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.638 + ! Biogeochemistry
1.639 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.640 + ! agbc
1.641 + ! aglbc
1.642 + ! agnpp
1.643 + case ('AGNPP')
1.644 + ! AGNPP: g m-2 s-1 -> kg m-2 s-1
1.645 + sfactor = 1.0e-3
1.646 + var_data%units = 'kg m-2 s-1'
1.647 + var_data%out_varname = 'agnpp'
1.648 + ! ar
1.649 + case ('AR')
1.650 + ! AR: g m-2 s-1 -> kg m-2 s-1
1.651 + sfactor = 1.0e-3
1.652 + var_data%units = 'kg m-2 s-1'
1.653 + var_data%out_varname = 'ar'
1.654 + ! bco
1.655 + case ('BIOGENCO')
1.656 + ! BIOGENCO: ug m-2 h-1 -> kg m-2 s-1
1.657 + sfactor = 1.0e-9 * (1.0 / 3600.)
1.658 + var_data%units = 'kg m-2 s-1'
1.659 + var_data%out_varname = 'bco'
1.660 + ! bgbc
1.661 + ! bgnpp
1.662 + case ('BGNPP')
1.663 + ! BGNPP: g m-2 s-1 -> kg m-2 s-1
1.664 + sfactor = 1.0e-3
1.665 + var_data%units = 'kg m-2 s-1'
1.666 + var_data%out_varname = 'bgnpp'
1.667 + ! bipn
1.668 + case ('ISOPRENE')
1.669 + ! ISOPRENE: ug m-2 h-1 -> kg m-2 s-1
1.670 + sfactor = 1.0e-9 * (1.0 / 3600.)
1.671 + var_data%units = 'kg m-2 s-1'
1.672 + var_data%out_varname = 'bipn'
1.673 + ! bmtpn
1.674 + case ('MONOTERP')
1.675 + ! MONOTERP: ug m-2 h-1 -> kg m-2 s-1
1.676 + sfactor = 1.0e-9 * (1.0 / 3600.)
1.677 + var_data%units = 'kg m-2 s-1'
1.678 + var_data%out_varname = 'bmtpn'
1.679 + ! borvoc
1.680 + case ('ORVOC')
1.681 + ! ORVOC: ug m-2 h-1 -> kg m-2 s-1
1.682 + sfactor = 1.0e-9 * (1.0 / 3600.)
1.683 + var_data%units = 'kg m-2 s-1'
1.684 + var_data%out_varname = 'borvoc'
1.685 + ! bovoc
1.686 + case ('OVOC')
1.687 + ! OVOC: ug m-2 h-1 -> kg m-2 s-1
1.688 + sfactor = 1.0e-9 * (1.0 / 3600.)
1.689 + var_data%units = 'kg m-2 s-1'
1.690 + var_data%out_varname = 'bovoc'
1.691 + ! cwdc
1.692 + case ('CWDC')
1.693 + ! CWDC: g m-2 -> kg m-2
1.694 + sfactor = 1.0e-3
1.695 + var_data%units = 'kg m-2'
1.696 + var_data%out_varname = 'cwdc'
1.697 + ! cwdc_hr
1.698 + case ('CWDC_HR')
1.699 + ! CWDC_HR: g m-2 s-1 -> kg m-2 s-1
1.700 + sfactor = 1.0e-3
1.701 + var_data%units = 'kg m-2 s-1'
1.702 + var_data%out_varname = 'cwdc_hr'
1.703 + ! cwdc_loss
1.704 + case ('CWDC_LOSS')
1.705 + ! CWDC_LOSS: g m-2 s-1 -> kg m-2 s-1
1.706 + sfactor = 1.0e-3
1.707 + if (casa_flux_bug) then
1.708 + print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
1.709 + sfactor = sfactor * (1. / 1200.)
1.710 + end if
1.711 + var_data%units = 'kg m-2 s-1'
1.712 + var_data%out_varname = 'cwdc_loss'
1.713 + ! froot_litterc
1.714 + ! froot_litterc_hr
1.715 + ! froot_litterc_loss
1.716 + ! frootc
1.717 + case ('FROOTC')
1.718 + ! FROOTC: g m-2 -> kg m-2
1.719 + sfactor = 1.0e-3
1.720 + var_data%units = 'kg m-2'
1.721 + var_data%out_varname = 'frootc'
1.722 + ! frootc_alloc
1.723 + case ('FROOTC_ALLOC')
1.724 + ! FROOTC_ALLOC: g m-2 s-1 -> kg m-2 s-1
1.725 + sfactor = 1.0e-3
1.726 + if (casa_flux_bug) then
1.727 + print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
1.728 + sfactor = sfactor * (1. / 1200.)
1.729 + end if
1.730 + var_data%units = 'kg m-2 s-1'
1.731 + var_data%out_varname = 'frootc_alloc'
1.732 + ! frootc_loss
1.733 + case ('FROOTC_LOSS')
1.734 + ! FROOTC_LOSS: g m-2 s-1 -> kg m-2 s-1
1.735 + sfactor = 1.0e-3
1.736 + if (casa_flux_bug) then
1.737 + print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
1.738 + sfactor = sfactor * (1. / 1200.)
1.739 + end if
1.740 + var_data%units = 'kg m-2 s-1'
1.741 + var_data%out_varname = 'frootc_loss'
1.742 + ! gpp
1.743 + case ('GPP')
1.744 + ! GPP: g m-2 s-1 -> kg m-2 s-1
1.745 + sfactor = 1.0e-3
1.746 + var_data%units = 'kg m-2 s-1'
1.747 + var_data%out_varname = 'gpp'
1.748 + ! gnmin
1.749 + case ('GROSS_NMIN')
1.750 + ! GROSS_NMIN: g m-2 s-1 -> kg m-2 s-1
1.751 + sfactor = 1.0e-3
1.752 + var_data%units = 'kg m-2 s-1'
1.753 + var_data%out_varname = 'gnmin'
1.754 + ! hr
1.755 + case ('HR')
1.756 + ! HR: g m-2 s-1 -> kg m-2 s-1
1.757 + sfactor = 1.0e-3
1.758 + var_data%units = 'kg m-2 s-1'
1.759 + var_data%out_varname = 'hr'
1.760 + ! leaf_litterc
1.761 + ! leaf_litterc_hr
1.762 + ! leaf_litterc_loss
1.763 + ! leafc
1.764 + case ('LEAFC')
1.765 + ! LEAFC: g m-2 -> kg m-2
1.766 + sfactor = 1.0e-3
1.767 + var_data%units = 'kg m-2'
1.768 + var_data%out_varname = 'leafc'
1.769 + ! leafc_alloc
1.770 + case ('LEAFC_ALLOC')
1.771 + ! LEAFC_ALLOC: g m-2 s-1 -> kg m-2 s-1
1.772 + sfactor = 1.0e-3
1.773 + if (casa_flux_bug) then
1.774 + print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
1.775 + sfactor = sfactor * (1. / 1200.)
1.776 + end if
1.777 + var_data%units = 'kg m-2 s-1'
1.778 + var_data%out_varname = 'leafc_alloc'
1.779 + ! leafc_loss
1.780 + case ('LEAFC_LOSS')
1.781 + ! LEAFC_LOSS: g m-2 s-1 -> kg m-2 s-1
1.782 + sfactor = 1.0e-3
1.783 + if (casa_flux_bug) then
1.784 + print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
1.785 + sfactor = sfactor * (1. / 1200.)
1.786 + end if
1.787 + var_data%units = 'kg m-2 s-1'
1.788 + var_data%out_varname = 'leafc_loss'
1.789 + ! litterc
1.790 + case ('LITTERC')
1.791 + ! LITTERC: g m-2 -> kg m-2
1.792 + sfactor = 1.0e-3
1.793 + var_data%units = 'kg m-2'
1.794 + var_data%out_varname = 'litterc'
1.795 + ! litterc_hr
1.796 + case ('LITTERC_HR')
1.797 + ! LITTERC_HR: g m-2 s-1 -> kg m-2 s-1
1.798 + sfactor = 1.0e-3
1.799 + var_data%units = 'kg m-2 s-1'
1.800 + var_data%out_varname = 'litterc_hr'
1.801 + ! litterc_loss
1.802 + case ('LITTERC_LOSS')
1.803 + ! LITTERC_LOSS: g m-2 s-1 -> kg m-2 s-1
1.804 + sfactor = 1.0e-3
1.805 + if (casa_flux_bug) then
1.806 + print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
1.807 + sfactor = sfactor * (1. / 1200.)
1.808 + end if
1.809 + var_data%units = 'kg m-2 s-1'
1.810 + var_data%out_varname = 'litterc_loss'
1.811 + ! nee
1.812 + case ('NEE')
1.813 + ! NEE: g m-2 s-1 -> kg m-2 s-1
1.814 + sfactor = 1.0e-3
1.815 + var_data%units = 'kg m-2 s-1'
1.816 + var_data%out_varname = 'nee'
1.817 + ! nep
1.818 + case ('NEP')
1.819 + ! NEP: g m-2 s-1 -> kg m-2 s-1
1.820 + sfactor = 1.0e-3
1.821 + var_data%units = 'kg m-2 s-1'
1.822 + var_data%out_varname = 'nep'
1.823 + ! nnmin
1.824 + case ('NET_NMIN')
1.825 + ! NEP: g m-2 s-1 -> kg m-2 s-1
1.826 + sfactor = 1.0e-3
1.827 + var_data%units = 'kg m-2 s-1'
1.828 + var_data%out_varname = 'nnmin'
1.829 + ! npp
1.830 + case ('NPP')
1.831 + ! NPP: g m-2 s-1 -> kg m-2 s-1
1.832 + sfactor = 1.0e-3
1.833 + var_data%units = 'kg m-2 s-1'
1.834 + var_data%out_varname = 'npp'
1.835 + ! soilc
1.836 + case ('SOILC')
1.837 + ! SOILC: g m-2 -> kg m-2
1.838 + sfactor = 1.0e-3
1.839 + var_data%units = 'kg m-2'
1.840 + var_data%out_varname = 'soilc'
1.841 + ! soilc_hr
1.842 + case ('SOILC_HR')
1.843 + ! SOILC_HR: g m-2 s-1 -> kg m-2 s-1
1.844 + sfactor = 1.0e-3
1.845 + var_data%units = 'kg m-2 s-1'
1.846 + var_data%out_varname = 'soilc_hr'
1.847 + ! soilc_loss
1.848 + case ('SOILC_LOSS')
1.849 + ! SOILC_LOSS: g m-2 s-1 -> kg m-2 s-1
1.850 + sfactor = 1.0e-3
1.851 + if (casa_flux_bug) then
1.852 + print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
1.853 + sfactor = sfactor * (1. / 1200.)
1.854 + end if
1.855 + var_data%units = 'kg m-2 s-1'
1.856 + var_data%out_varname = 'soilc_loss'
1.857 + ! tbvoc
1.858 + case ('VOCFLXT')
1.859 + ! VOCFLXT: ug m-2 h-1 -> kg m-2 s-1
1.860 + sfactor = 1.0e-9 * (1.0 / 3600.)
1.861 + var_data%units = 'kg m-2 s-1'
1.862 + var_data%out_varname = 'tbvoc'
1.863 + ! woodc
1.864 + case ('WOODC')
1.865 + ! WOODC: g m-2 -> kg m-2
1.866 + sfactor = 1.0e-3
1.867 + var_data%units = 'kg m-2'
1.868 + var_data%out_varname = 'woodc'
1.869 + ! woodc_alloc
1.870 + case ('WOODC_ALLOC')
1.871 + ! WOODC_ALLOC: g m-2 s-1 -> kg m-2 s-1
1.872 + sfactor = 1.0e-3
1.873 + if (casa_flux_bug) then
1.874 + print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
1.875 + sfactor = sfactor * (1. / 1200.)
1.876 + end if
1.877 + var_data%units = 'kg m-2 s-1'
1.878 + var_data%out_varname = 'woodc_alloc'
1.879 + ! woodc_loss
1.880 + case ('WOODC_LOSS')
1.881 + ! WOODC_LOSS: g m-2 s-1 -> kg m-2 s-1
1.882 + sfactor = 1.0e-3
1.883 + if (casa_flux_bug) then
1.884 + print *,'WARNING: casa_flux_bug set to true; dividing data by 1200'
1.885 + sfactor = sfactor * (1. / 1200.)
1.886 + end if
1.887 + var_data%units = 'kg m-2 s-1'
1.888 + var_data%out_varname = 'woodc_loss'
1.889 +
1.890 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.891 + ! Energy
1.892 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.893 + ! hfg
1.894 + case ('FGR')
1.895 + ! FGR: W m-2 -> W m-2
1.896 + sfactor = 1.0
1.897 + var_data%units = 'W m-2'
1.898 + var_data%positive = 'down'
1.899 + var_data%out_varname = 'hfg'
1.900 + ! hfls
1.901 + case ('LATENT')
1.902 + ! LATENT: W m-2 -> W m-2
1.903 + sfactor = 1.0
1.904 + var_data%units = 'W m-2'
1.905 + var_data%positive = 'up'
1.906 + var_data%out_varname = 'hfls'
1.907 + ! hfss
1.908 + case ('FSH')
1.909 + ! FSH: W m-2 -> W m-2
1.910 + sfactor = 1.0
1.911 + var_data%units = 'W m-2'
1.912 + var_data%positive = 'up'
1.913 + var_data%out_varname = 'hfss'
1.914 + ! hfssg
1.915 + case ('FSH_G')
1.916 + ! FSH_G: W m-2 -> W m-2
1.917 + sfactor = 1.0
1.918 + var_data%units = 'W m-2'
1.919 + var_data%positive = 'up'
1.920 + var_data%out_varname = 'hfssg'
1.921 + ! hfssv
1.922 + case ('FSH_V')
1.923 + ! FSH_V: W m-2 -> W m-2
1.924 + sfactor = 1.0
1.925 + var_data%units = 'W m-2'
1.926 + var_data%positive = 'up'
1.927 + var_data%out_varname = 'hfssv'
1.928 +
1.929 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.930 + ! Grid
1.931 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.932 + ! gca
1.933 + case ('area')
1.934 + ! area: km2 -> m2
1.935 + sfactor = 1.0e+3
1.936 + var_data%units = 'm2'
1.937 + var_data%out_varname = 'gca'
1.938 + ! lbm
1.939 + case ('landmask')
1.940 + ! landmask: 1 -> 1
1.941 + var_data%units = '1'
1.942 + var_data%out_varname = 'lbm'
1.943 + ! Do not convert integer landmask, just return
1.944 + return
1.945 + ! lcf
1.946 + ! lcn
1.947 + ! orog
1.948 + case ('topo')
1.949 + ! topo: m -> m
1.950 + sfactor = 1.0
1.951 + var_data%units = 'm'
1.952 + var_data%out_varname = 'orog'
1.953 + ! sftlf
1.954 + case ('landfrac')
1.955 + ! landfrac: 1 -> 1
1.956 + sfactor = 1.0
1.957 + var_data%units = '1'
1.958 + var_data%out_varname = 'sftlf'
1.959 +
1.960 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.961 + ! Humidity
1.962 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.963 + ! hur2m
1.964 + case ('RH2M')
1.965 + ! RH2M: % -> %
1.966 + sfactor = 1.0
1.967 + var_data%units = '%'
1.968 + var_data%out_varname = 'hur2m'
1.969 + ! hus2m
1.970 + case ('Q2M')
1.971 + ! Q2M: kg kg-1 -> kg kg-1
1.972 + sfactor = 1.0
1.973 + var_data%units = 'kg kg-1'
1.974 + var_data%out_varname = 'hus2m'
1.975 +
1.976 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.977 + ! Hydrology
1.978 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.979 + ! btran
1.980 + case ('BTRAN')
1.981 + ! BTRAN: 1 -> 1
1.982 + sfactor = 1.0
1.983 + var_data%units = '1'
1.984 + var_data%out_varname = 'btran'
1.985 + ! et
1.986 + case ('ET')
1.987 + ! ET: m s-1 -> kg m-2 s-1
1.988 + sfactor = 1.0e+3
1.989 + var_data%units = 'kg m-2 s-1'
1.990 + var_data%out_varname = 'et'
1.991 + ! evspsblsoi
1.992 + case ('QSOIL')
1.993 + ! QSOIL: mm s-1 -> kg m-2 s-1
1.994 + sfactor = 1.0
1.995 + var_data%units = 'kg m-2 s-1'
1.996 + var_data%out_varname = 'evspsblsoi'
1.997 + ! evspsblveg
1.998 + case ('QVEGE')
1.999 + ! QVEGE: mm s-1 -> kg m-2 s-1
1.1000 + sfactor = 1.0
1.1001 + var_data%units = 'kg m-2 s-1'
1.1002 + var_data%out_varname = 'evspsblveg'
1.1003 + ! mrfsl (multi-level)
1.1004 + case ('SOILICE')
1.1005 + ! SOILICE: kg m-2 -> kg m-2
1.1006 + sfactor = 1.0
1.1007 + var_data%units = 'kg m-2'
1.1008 + var_data%out_varname = 'mrfsl'
1.1009 + ! mrfso
1.1010 + case ('TOTAL_SOILICE')
1.1011 + ! TOTAL_SOILICE: kg m-2 -> kg m-2
1.1012 + sfactor = 1.0
1.1013 + var_data%units = 'kg m-2'
1.1014 + var_data%out_varname = 'mrfso'
1.1015 + ! mrlsl (multi-level)
1.1016 + case ('SOILLIQ')
1.1017 + ! SOILLIQ: kg m-2 -> kg m-2
1.1018 + sfactor = 1.0
1.1019 + var_data%units = 'kg m-2'
1.1020 + var_data%out_varname = 'mrlsl'
1.1021 + ! mrlso
1.1022 + case ('TOTAL_SOILLIQ')
1.1023 + ! TOTAL_SOILLIQ: kg m-2 -> kg m-2
1.1024 + sfactor = 1.0
1.1025 + var_data%units = 'kg m-2'
1.1026 + var_data%out_varname = 'mrlso'
1.1027 + ! mrros
1.1028 + case ('QOVER')
1.1029 + ! QOVER: mm s-1 -> kg m-2 s-1
1.1030 + sfactor = 1.0
1.1031 + var_data%units = 'kg m-2 s-1'
1.1032 + var_data%out_varname = 'mrros'
1.1033 + ! mrross
1.1034 + case ('QDRAI')
1.1035 + ! QDRAI: mm s-1 -> kg m-2 s-1
1.1036 + sfactor = 1.0
1.1037 + var_data%units = 'kg m-2 s-1'
1.1038 + var_data%out_varname = 'mrross'
1.1039 + ! prveg
1.1040 + case ('QINTR')
1.1041 + ! QINTR: mm s-1 -> kg m-2 s-1
1.1042 + sfactor = 1.0
1.1043 + var_data%units = 'kg m-2 s-1'
1.1044 + var_data%out_varname = 'prveg'
1.1045 + ! snd
1.1046 + case ('H2OSNO')
1.1047 + ! H2OSNO: mm -> m
1.1048 + sfactor = 1.0e-3
1.1049 + var_data%units = 'm'
1.1050 + var_data%out_varname = 'snd'
1.1051 + ! tran
1.1052 + case ('QVEGT')
1.1053 + ! QVEGT: mm s-1 -> kg m-2 s-1
1.1054 + sfactor = 1.0
1.1055 + var_data%units = 'kg m-2 s-1'
1.1056 + var_data%out_varname = 'tran'
1.1057 + ! vfwsl (multi-level)
1.1058 + case ('H2OSOI')
1.1059 + ! H2OSOI: mm3 mm-3 -> m3 m-3
1.1060 + sfactor = 1.0
1.1061 + var_data%units = 'm3 m-3'
1.1062 + var_data%out_varname = 'vfwsl'
1.1063 + ! wpsl (multi-level)
1.1064 + case ('SOILPSI')
1.1065 + ! SOILPSI: Pa -> Pa
1.1066 + sfactor = 1.0
1.1067 + var_data%units = 'Pa'
1.1068 + var_data%out_varname = 'wpsl'
1.1069 +
1.1070 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.1071 + ! Radiation
1.1072 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.1073 + ! alsas
1.1074 + case ('ALL_SKY_ALBEDO')
1.1075 + ! ALL_SKY_ALBEDO: 1 -> 1
1.1076 + sfactor = 1.0
1.1077 + var_data%units = '1'
1.1078 + var_data%out_varname = 'alsas'
1.1079 + ! alsbs
1.1080 + case ('BLACK_SKY_ALBEDO')
1.1081 + ! BLACK_SKY_ALBEDO: 1 -> 1
1.1082 + sfactor = 1.0
1.1083 + var_data%units = '1'
1.1084 + var_data%out_varname = 'alsbs'
1.1085 + ! laisha
1.1086 + case ('LAISHA')
1.1087 + ! LAISHA: m2 m-2 -> m2 m-2
1.1088 + sfactor = 1.0
1.1089 + var_data%units = 'm2 m-2'
1.1090 + var_data%out_varname = 'laisha'
1.1091 + ! laisun
1.1092 + case ('LAISUN')
1.1093 + ! LAISUN: m2 m-2 -> m2 m-2
1.1094 + sfactor = 1.0
1.1095 + var_data%units = 'm2 m-2'
1.1096 + var_data%out_varname = 'laisun'
1.1097 + ! rlds
1.1098 + case ('FLDS')
1.1099 + ! FLDS: W m-2 -> W m-2
1.1100 + sfactor = 1.0
1.1101 + var_data%units = 'W m-2'
1.1102 + var_data%positive = 'down'
1.1103 + var_data%out_varname = 'rlds'
1.1104 + ! rls
1.1105 + case ('FIRA')
1.1106 + ! FIRA: W m-2 -> W m-2
1.1107 + sfactor = 1.0
1.1108 + var_data%units = 'W m-2'
1.1109 + var_data%positive = 'down'
1.1110 + var_data%out_varname = 'rls'
1.1111 + ! rlus
1.1112 + case ('FIRE')
1.1113 + ! FIRE: W m-2 -> W m-2
1.1114 + sfactor = 1.0
1.1115 + var_data%units = 'W m-2'
1.1116 + var_data%positive = 'up'
1.1117 + var_data%out_varname = 'rlus'
1.1118 + ! rs
1.1119 + case ('NETRAD')
1.1120 + ! NETRAD: W m-2 -> W m-2
1.1121 + sfactor = 1.0
1.1122 + var_data%units = 'W m-2'
1.1123 + var_data%positive = 'down'
1.1124 + var_data%out_varname = 'rs'
1.1125 + ! rsas
1.1126 + case ('FSA')
1.1127 + ! FSA: W m-2 -> W m-2
1.1128 + sfactor = 1.0
1.1129 + var_data%units = 'W m-2'
1.1130 + var_data%positive = 'down'
1.1131 + var_data%out_varname = 'rsas'
1.1132 + ! rsasg
1.1133 + case ('SABG')
1.1134 + ! SABG: W m-2 -> W m-2
1.1135 + sfactor = 1.0
1.1136 + var_data%units = 'W m-2'
1.1137 + var_data%positive = 'down'
1.1138 + var_data%out_varname = 'rsasg'
1.1139 + ! rsasv
1.1140 + case ('SABV')
1.1141 + ! SABV: W m-2 -> W m-2
1.1142 + sfactor = 1.0
1.1143 + var_data%units = 'W m-2'
1.1144 + var_data%positive = 'down'
1.1145 + var_data%out_varname = 'rsasv'
1.1146 + ! rsds
1.1147 + case ('FSDS')
1.1148 + ! FSDS: W m-2 -> W m-2
1.1149 + sfactor = 1.0
1.1150 + var_data%units = 'W m-2'
1.1151 + var_data%positive = 'down'
1.1152 + var_data%out_varname = 'rsds'
1.1153 + ! rsus
1.1154 + case ('FSR')
1.1155 + ! FSR: W m-2 -> W m-2
1.1156 + sfactor = 1.0
1.1157 + var_data%units = 'W m-2'
1.1158 + var_data%positive = 'up'
1.1159 + var_data%out_varname = 'rsus'
1.1160 + ! snc
1.1161 + case ('FSNO')
1.1162 + ! FSNO: 1 -> 1
1.1163 + sfactor = 1.0
1.1164 + var_data%units = '1'
1.1165 + var_data%out_varname = 'snc'
1.1166 +
1.1167 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.1168 + ! Temperature
1.1169 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.1170 + ! degday
1.1171 + case ('DEGDAY')
1.1172 + ! DEGDAY: C -> K
1.1173 + sfactor = 1.0
1.1174 + offset = 273.15
1.1175 + var_data%units = 'K'
1.1176 + var_data%out_varname = 'degday'
1.1177 + ! ta2m
1.1178 + case ('TSA')
1.1179 + ! TSA: K -> K
1.1180 + sfactor = 1.0
1.1181 + var_data%units = 'K'
1.1182 + var_data%out_varname = 'ta2m'
1.1183 + ! ta2mdmax
1.1184 + case ('TREFMXAV')
1.1185 + ! TREFMXAV: K -> K
1.1186 + sfactor = 1.0
1.1187 + var_data%units = 'K'
1.1188 + var_data%out_varname = 'ta2mdmax'
1.1189 + ! ta2mdmin
1.1190 + case ('TREFMNAV')
1.1191 + ! TREFMNAV: K -> K
1.1192 + sfactor = 1.0
1.1193 + var_data%units = 'K'
1.1194 + var_data%out_varname = 'ta2mdmin'
1.1195 + ! tsl (multi-level)
1.1196 + case ('TSOI')
1.1197 + ! TSOI: K -> K
1.1198 + sfactor = 1.0
1.1199 + var_data%units = 'K'
1.1200 + var_data%out_varname = 'tsl'
1.1201 + ! tv
1.1202 + case ('TV')
1.1203 + ! TV: K -> K
1.1204 + sfactor = 1.0
1.1205 + var_data%units = 'K'
1.1206 + var_data%out_varname = 'tv'
1.1207 +
1.1208 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.1209 + ! Vegetation phenology
1.1210 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.1211 + ! laiexp
1.1212 + case ('ELAI')
1.1213 + ! ELAI: m2 m-2 -> m2 m-2
1.1214 + sfactor = 1.0
1.1215 + var_data%units = 'm2 m-2'
1.1216 + var_data%out_varname = 'laiexp'
1.1217 + ! laitot
1.1218 + case ('TLAI')
1.1219 + ! TLAI: m2 m-2 -> m2 m-2
1.1220 + sfactor = 1.0
1.1221 + var_data%units = 'm2 m-2'
1.1222 + var_data%out_varname = 'laitot'
1.1223 + ! saiexp
1.1224 + case ('ESAI')
1.1225 + ! ESAI: m2 m-2 -> m2 m-2
1.1226 + sfactor = 1.0
1.1227 + var_data%units = 'm2 m-2'
1.1228 + var_data%out_varname = 'saiexp'
1.1229 + ! saitot
1.1230 + case ('TSAI')
1.1231 + ! TSAI: m2 m-2 -> m2 m-2
1.1232 + sfactor = 1.0
1.1233 + var_data%units = 'm2 m-2'
1.1234 + var_data%out_varname = 'saitot'
1.1235 +
1.1236 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.1237 + ! Vegetation physiology
1.1238 + !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1.1239 + ! psn
1.1240 + case ('FPSN')
1.1241 + ! FPSN: umol m-2 s-1 -> kg m-2 s-1
1.1242 + sfactor = 12.0 * 1.0e-6 * 1.0e-3
1.1243 + var_data%units = 'kg m-2 s-1'
1.1244 + var_data%out_varname = 'psn'
1.1245 + ! rssha
1.1246 + case ('RSSHA')
1.1247 + ! RSSHA: s m-1 -> s m-1
1.1248 + sfactor = 1.0
1.1249 + var_data%units = 's m-1'
1.1250 + var_data%out_varname = 'rssha'
1.1251 + ! rssun
1.1252 + case ('RSSUN')
1.1253 + ! RSSUN: s m-1 -> s m-1
1.1254 + sfactor = 1.0
1.1255 + var_data%units = 's m-1'
1.1256 + var_data%out_varname = 'rssun'
1.1257 +
1.1258 + case default
1.1259 + print *, 'clm_convert_data: Unknown variable name ', var_data%varname
1.1260 + stop
1.1261 + end select
1.1262 +
1.1263 + ! Convert to the desired units using the offset and scale factor (sfactor)
1.1264 + select case (var_data%ndims)
1.1265 + case (2)
1.1266 + allocate(mask2d(coord%xsize, coord%ysize), stat=ierr)
1.1267 + if (ierr /= 0) then
1.1268 + print *, 'clm_convert_data: Cannot allocate 2d mask'
1.1269 + stop
1.1270 + end if
1.1271 + mask2d = .false.
1.1272 + if (var_data%int_type) then
1.1273 + where (var_data%int2d /= var_data%missing_value)
1.1274 + var_data%int2d = int(var_data%int2d * sfactor + offset)
1.1275 + mask2d = .true.
1.1276 + end where
1.1277 + print *, var_data%out_varname, &
1.1278 + ' min:', minval(var_data%int2d, mask=mask2d), &
1.1279 + ' max:', maxval(var_data%int2d, mask=mask2d)
1.1280 + else
1.1281 + where (var_data%var2d /= var_data%missing_value)
1.1282 + var_data%var2d = var_data%var2d * sfactor + offset
1.1283 + mask2d = .true.
1.1284 + end where
1.1285 + print *, var_data%out_varname, &
1.1286 + ' min:', minval(var_data%var2d, mask=mask2d), &
1.1287 + ' max:', maxval(var_data%var2d, mask=mask2d)
1.1288 + end if
1.1289 + if (allocated(mask2d)) deallocate(mask2d)
1.1290 + case (3)
1.1291 + allocate(mask3d(coord%xsize, coord%ysize, coord%tsize), stat=ierr)
1.1292 + if (ierr /= 0) then
1.1293 + print *, 'clm_convert_data: Cannot allocate 3d mask'
1.1294 + stop
1.1295 + end if
1.1296 + mask3d = .false.
1.1297 + where (var_data%var3d /= var_data%missing_value)
1.1298 + var_data%var3d = (var_data%var3d + offset) * sfactor
1.1299 + mask3d = .true.
1.1300 + end where
1.1301 + print *, var_data%out_varname, &
1.1302 + ' min:', minval(var_data%var3d, mask=mask3d), &
1.1303 + ' max:', maxval(var_data%var3d, mask=mask3d)
1.1304 + if (allocated(mask3d)) deallocate(mask3d)
1.1305 + case (4)
1.1306 + if (var_data%soil_layer_flag) then
1.1307 + allocate(mask4d(coord%xsize, coord%ysize, coord%zsoi_size, coord%tsize), stat=ierr)
1.1308 + if (ierr /= 0) then
1.1309 + print *, 'clm_convert_data: Cannot allocate 4d mask'
1.1310 + stop
1.1311 + end if
1.1312 + mask4d = .false.
1.1313 + where (var_data%var4d /= var_data%missing_value)
1.1314 + var_data%var4d = (var_data%var4d + offset) * sfactor
1.1315 + mask4d = .true.
1.1316 + end where
1.1317 + print *, var_data%out_varname, &
1.1318 + ' min:', minval(var_data%var4d, mask=mask4d), &
1.1319 + ' max:', maxval(var_data%var4d, mask=mask4d)
1.1320 + if (allocated(mask4d)) deallocate(mask4d)
1.1321 + else
1.1322 + print *,'clm_convert_data: Unknown data dimension'
1.1323 + stop
1.1324 + end if
1.1325 + case default
1.1326 + print *,'clm_convert_data: Dimension count not supported'
1.1327 + stop
1.1328 + end select
1.1329 +
1.1330 + end subroutine clm_convert_data
1.1331 + !
1.1332 + !--------------------------------------------------------------------------
1.1333 + !
1.1334 + subroutine clm_free_data
1.1335 +
1.1336 + implicit none
1.1337 +
1.1338 + var_data%varname = ''
1.1339 + var_data%ndims = 0
1.1340 + if (allocated(var_data%int2d)) deallocate(var_data%int2d)
1.1341 + if (allocated(var_data%var2d)) deallocate(var_data%var2d)
1.1342 + if (allocated(var_data%var3d)) deallocate(var_data%var3d)
1.1343 + if (allocated(var_data%var4d)) deallocate(var_data%var4d)
1.1344 + if (allocated(var_data%time)) deallocate(var_data%time)
1.1345 + if (allocated(var_data%time_bounds)) deallocate(var_data%time_bounds)
1.1346 + var_data%int_type = .false.
1.1347 + var_data%soil_layer_flag = .false.
1.1348 + var_data%hour_flag = .false.
1.1349 + var_data%long_name = ''
1.1350 + var_data%missing_value = 0.
1.1351 + var_data%units = ''
1.1352 + var_data%positive = ''
1.1353 + var_data%out_varname = ''
1.1354 +
1.1355 + end subroutine clm_free_data
1.1356 + !
1.1357 + !--------------------------------------------------------------------------
1.1358 + !
1.1359 + subroutine nc_err(status)
1.1360 +
1.1361 + implicit none
1.1362 + integer :: status ! status code
1.1363 +
1.1364 + if (status /= nf90_noerr) then
1.1365 + print *, nf90_strerror(status)
1.1366 + stop 'Program stopped after encountering netCDF error.'
1.1367 + end if
1.1368 +
1.1369 + end subroutine nc_err
1.1370 + !
1.1371 + !--------------------------------------------------------------------------
1.1372 + !
1.1373 +
1.1374 +end module clm_mod