forrest@0: module clm_mod forrest@0: forrest@0: use netcdf forrest@0: use kind_mod forrest@0: implicit none forrest@0: forrest@0: integer, parameter :: namelen = 128 ! string length for paths forrest@0: integer, parameter :: varnamelen = 32 ! string length for varnames forrest@0: integer, parameter :: posnamelen = 4 ! string length for positive attribute forrest@0: integer, parameter :: varattlen = 128 ! string length for variable attributes forrest@0: character (len=*), parameter :: xname = 'lon' forrest@0: character (len=*), parameter :: yname = 'lat' forrest@0: character (len=*), parameter :: zsoi_name = 'levsoi' forrest@0: character (len=*), parameter :: zlak_name = 'levlak' forrest@0: character (len=*), parameter :: tname = 'time' forrest@0: character (len=*), parameter :: hr_name = 'hour' forrest@0: character (len=*), parameter :: nsamples_name = 'nsamples' forrest@0: character (len=*), parameter :: hiname = 'hist_interval' forrest@0: character (len=*), parameter :: tbname = 'time_bounds' forrest@0: character (len=*), parameter :: aname = 'area' forrest@0: character (len=*), parameter :: lname = 'landfrac' forrest@0: character (len=*), parameter :: mname = 'landmask' forrest@0: character (len=*), parameter :: oname = 'topo' forrest@0: character (len=*), parameter :: long_name_name = 'long_name' forrest@0: character (len=*), parameter :: missing_value_name = 'missing_value' forrest@0: character (len=*), parameter :: units_name = 'units' forrest@0: character (len=*), parameter :: case_id_name = 'case_id' forrest@0: character (len=*), parameter :: inidat_name = 'Initial_conditions_dataset' forrest@0: character (len=*), parameter :: surdat_name = 'Surface_dataset' forrest@0: character (len=*), parameter :: pftdat_name = 'PFT_physiological_constants_dataset' forrest@0: character (len=*), parameter :: rtmdat_name = 'RTM_input_datset' forrest@0: forrest@0: type surface forrest@0: integer :: xsize ! number of longitude points (lon) forrest@0: integer :: ysize ! number of latitutde points (lat) forrest@0: integer :: psize ! number of pfts forrest@0: real(r8), pointer :: x(:) ! longitude values (lon) forrest@0: real(r8), pointer :: y(:) ! latitude values (lat) forrest@0: real(r8), pointer :: p(:) ! pft number forrest@0: real(r4), pointer :: pft(:,:,:) ! percent of each pft forrest@0: end type surface forrest@0: forrest@0: type coordinates forrest@0: integer :: xsize ! number of longitude points (lon) forrest@0: integer :: ysize ! number of latitutde points (lat) forrest@0: integer :: zsoi_size ! number of soil z-levels (levsoi) forrest@0: integer :: zlak_size ! number of lake z-levels (levlak) forrest@0: integer :: tsize ! number of time points (time) forrest@0: integer :: hr_size ! number of hour points (hour = 24) forrest@0: integer :: nsamples ! number of time samples in post-processing statistics forrest@0: real(r8), pointer :: x(:) ! longitude values (lon) forrest@0: real(r8), pointer :: y(:) ! latitude values (lat) forrest@0: real(r8), pointer :: zsoi(:) ! soil z-level values (levsoi) forrest@0: real(r8), pointer :: zlak(:) ! lake z-level values (levlak) forrest@0: real(r8), pointer :: hr(:) ! hour-of-day values (hour) forrest@0: real(r8), pointer :: area(:,:) ! gridcell areas forrest@0: real(r8), pointer :: landfrac(:,:) ! gridcell land fractions forrest@0: integer, pointer :: landmask(:,:) ! gridcell land/ocean mask forrest@0: real(r8), pointer :: orog(:,:) ! gridcell orography/topography forrest@0: character (len=namelen) :: case_id ! case name of run forrest@0: character (len=namelen) :: inidat ! initial dataset filename forrest@0: character (len=namelen) :: surdat ! surface dataset filename forrest@0: character (len=namelen) :: pftdat ! pft dataset filename forrest@0: character (len=namelen) :: rtmdat ! runoff dataset filename forrest@0: logical :: hour_flag ! flag for hour-of-day statistical summary files forrest@0: end type coordinates forrest@0: forrest@0: type variable_data forrest@0: character (len=varnamelen) :: varname forrest@0: integer :: ndims forrest@0: integer, pointer :: int2d(:,:) forrest@0: real(r4), pointer :: var2d(:,:) forrest@0: real(r4), pointer :: var3d(:,:,:) forrest@0: real(r4), pointer :: var4d(:,:,:,:) forrest@0: real(r8), pointer :: time(:) forrest@0: real(r8), pointer :: time_bounds(:,:) forrest@0: logical :: int_type ! set if type is integer forrest@0: logical :: soil_layer_flag forrest@0: logical :: hour_flag ! flag for hour-of-day statistical summary variable forrest@0: character (len=varattlen) :: long_name forrest@0: real(r4) :: missing_value forrest@0: character (len=varattlen) :: units forrest@0: character (len=posnamelen) :: positive forrest@0: character (len=varnamelen) :: out_varname forrest@0: end type variable_data forrest@0: forrest@0: type(surface), target :: surf forrest@0: type(coordinates), target :: coord forrest@0: type(variable_data), target :: var_data forrest@0: save forrest@0: forrest@0: public clm_read_surf forrest@0: public clm_free_surf forrest@0: public clm_read_coord forrest@0: public clm_copy_grid_data forrest@0: public clm_read_data forrest@0: public clm_convert_data forrest@0: public clm_free_data forrest@0: private nc_err forrest@0: forrest@0: contains forrest@0: forrest@0: ! forrest@0: !-------------------------------------------------------------------------- forrest@0: ! forrest@0: subroutine clm_read_surf(fsurdat) forrest@0: forrest@0: implicit none forrest@0: character (len=namelen), intent(in) :: fsurdat ! input surface dataset forrest@0: integer :: i forrest@0: integer :: ierr forrest@0: integer :: ncid forrest@0: integer :: xdimid, xvarid, xlen forrest@0: integer :: ydimid, yvarid, ylen forrest@0: integer :: pdimid, pvarid, plen forrest@0: real(r8), allocatable :: longxy(:,:) forrest@0: real(r8), allocatable :: latixy(:,:) forrest@0: real(r4), allocatable :: pctpft(:,:,:) forrest@0: forrest@0: ! Read coordinated from surface dataset forrest@0: print *, 'Opening ',trim(fsurdat) forrest@0: call nc_err(nf90_open(path=trim(fsurdat), mode=nf90_nowrite, ncid=ncid)) forrest@0: call nc_err(nf90_inq_dimid(ncid, 'lsmlon', xdimid)) forrest@0: call nc_err(nf90_inq_dimid(ncid, 'lsmlat', ydimid)) forrest@0: call nc_err(nf90_inq_dimid(ncid, 'lsmpft', pdimid)) forrest@0: forrest@0: call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=xdimid, len=xlen)) forrest@0: call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=ydimid, len=ylen)) forrest@0: call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=pdimid, len=plen)) forrest@0: forrest@0: call nc_err(nf90_inq_varid(ncid, 'LONGXY', xvarid)) forrest@0: call nc_err(nf90_inq_varid(ncid, 'LATIXY', yvarid)) forrest@0: call nc_err(nf90_inq_varid(ncid, 'PCT_PFT', pvarid)) forrest@0: forrest@0: allocate(longxy(xlen, ylen), latixy(xlen, ylen), pctpft(xlen, ylen, plen), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_read_surf: Cannot allocate longxy, latixy, pctpft variables' forrest@0: stop forrest@0: end if forrest@0: forrest@0: call nc_err(nf90_get_var(ncid, xvarid, longxy)) forrest@0: call nc_err(nf90_get_var(ncid, yvarid, latixy)) forrest@0: call nc_err(nf90_get_var(ncid, pvarid, pctpft)) forrest@0: forrest@0: print *, 'Closing ',trim(fsurdat) forrest@0: call nc_err(nf90_close(ncid)) forrest@0: forrest@0: ! Allocate and fill in the surface derived data type forrest@0: surf%xsize = xlen forrest@0: surf%ysize = ylen forrest@0: surf%psize = plen forrest@0: allocate(surf%x(xlen), surf%y(ylen), surf%p(plen), surf%pft(xlen, ylen, plen), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_read_surf: Cannot allocate surface type variables' forrest@0: stop forrest@0: end if forrest@0: surf%x(:) = longxy(:,1) forrest@0: surf%y(:) = latixy(1,:) forrest@0: do i = 1, plen forrest@0: surf%p(i) = real(i) forrest@0: end do forrest@0: surf%pft = pctpft forrest@0: forrest@0: deallocate(longxy, latixy, pctpft) forrest@0: forrest@0: end subroutine clm_read_surf forrest@0: ! forrest@0: !-------------------------------------------------------------------------- forrest@0: ! forrest@0: subroutine clm_free_surf() forrest@0: forrest@0: implicit none forrest@0: forrest@0: if (allocated(surf%x)) deallocate(surf%x) forrest@0: if (allocated(surf%y)) deallocate(surf%y) forrest@0: if (allocated(surf%p)) deallocate(surf%p) forrest@0: if (allocated(surf%pft)) deallocate(surf%pft) forrest@0: forrest@0: end subroutine clm_free_surf forrest@0: ! forrest@0: !-------------------------------------------------------------------------- forrest@0: ! forrest@0: subroutine clm_read_coord(numf, input_path, fnames, read_grid_all) forrest@0: forrest@0: implicit none forrest@0: integer, intent(in) :: numf ! number of filenames forrest@0: character (len=namelen), intent(in) :: input_path ! input path forrest@0: character (len=namelen), intent(in) :: fnames(numf) ! input filenames forrest@0: logical, intent(in) :: read_grid_all ! flag to read area,... forrest@0: integer :: ncid forrest@0: integer :: xdimid, xvarid, xlen forrest@0: integer :: ydimid, yvarid, ylen forrest@0: integer :: zsoi_dimid, zsoi_varid, zsoi_len forrest@0: integer :: zlak_dimid, zlak_varid, zlak_len forrest@0: integer :: tdimid, tvarid, tlen forrest@0: integer :: hr_dimid, hr_varid, hr_len forrest@0: integer :: nsamples_dimid, nsamples_len forrest@0: integer :: avarid ! area (when read_grid_all true) forrest@0: integer :: lvarid ! landfrac (when read_grid_all true) forrest@0: integer :: mvarid ! landmask (when read_grid_all true) forrest@0: integer :: ovarid ! orog/topo (when read_grid_all true) forrest@0: integer(i8) :: i forrest@0: integer(i8) :: time_size forrest@0: real(r4), allocatable :: xvals(:) forrest@0: real(r4), allocatable :: yvals(:) forrest@0: real(r4), allocatable :: zsoi_vals(:) forrest@0: real(r4), allocatable :: zlak_vals(:) forrest@0: real(r4), allocatable :: hr_vals(:) forrest@0: real(r4), allocatable :: avals(:,:) ! area (when read_grid_all true) forrest@0: real(r4), allocatable :: lvals(:,:) ! landfrac (when read_grid_all true) forrest@0: integer, allocatable :: mvals(:,:) ! landmask (when read_grid_all true) forrest@0: real(r4), allocatable :: ovals(:,:) ! orog/topo (when read_grid_all true) forrest@0: character (len=namelen) :: case_id = '' forrest@0: character (len=namelen) :: inidat = '' forrest@0: character (len=namelen) :: surdat = '' forrest@0: character (len=namelen) :: pftdat = '' forrest@0: character (len=namelen) :: rtmdat = '' forrest@0: logical :: hrflag = .false. forrest@0: integer :: ierr forrest@0: forrest@0: ! Read coordinates from first file forrest@0: print *, 'Opening ',trim(input_path)//'/'//trim(fnames(1)) forrest@0: call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(1)), mode=nf90_nowrite, ncid=ncid)) forrest@0: call nc_err(nf90_inq_dimid(ncid, xname, xdimid)) forrest@0: call nc_err(nf90_inq_dimid(ncid, yname, ydimid)) forrest@0: call nc_err(nf90_inq_dimid(ncid, zsoi_name, zsoi_dimid)) forrest@0: call nc_err(nf90_inq_dimid(ncid, zlak_name, zlak_dimid)) forrest@0: call nc_err(nf90_inq_dimid(ncid, tname, tdimid)) forrest@0: if (nf90_inq_dimid(ncid, hr_name, hr_dimid) == nf90_noerr) then forrest@0: call nc_err(nf90_inq_dimid(ncid, nsamples_name, nsamples_dimid)) forrest@0: hrflag = .true. forrest@0: end if forrest@0: forrest@0: call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=xdimid, len=xlen)) forrest@0: call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=ydimid, len=ylen)) forrest@0: call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=zsoi_dimid, len=zsoi_len)) forrest@0: call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=zlak_dimid, len=zlak_len)) forrest@0: call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=tdimid, len=tlen)) forrest@0: if (hrflag) then forrest@0: call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=hr_dimid, len=hr_len)) forrest@0: call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=nsamples_dimid, len=nsamples_len)) forrest@0: else forrest@0: hr_len = 0 forrest@0: nsamples_len = 0 forrest@0: end if forrest@0: print *, tlen,' values' forrest@0: forrest@0: allocate(xvals(xlen), yvals(ylen), zsoi_vals(zsoi_len), zlak_vals(zlak_len), hr_vals(hr_len), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_read_coord: Cannot allocate grid variables' forrest@0: stop forrest@0: end if forrest@0: if (read_grid_all) then forrest@0: allocate(avals(xlen, ylen), lvals(xlen, ylen), mvals(xlen, ylen), ovals(xlen, ylen), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_read_coord: Cannot allocate extra grid variables' forrest@0: stop forrest@0: end if forrest@0: end if forrest@0: forrest@0: call nc_err(nf90_inq_varid(ncid, xname, xvarid)) forrest@0: call nc_err(nf90_inq_varid(ncid, yname, yvarid)) forrest@0: call nc_err(nf90_inq_varid(ncid, zsoi_name, zsoi_varid)) forrest@0: call nc_err(nf90_inq_varid(ncid, zlak_name, zlak_varid)) forrest@0: call nc_err(nf90_inq_varid(ncid, tname, tvarid)) forrest@0: if (hrflag) call nc_err(nf90_inq_varid(ncid, hr_name, hr_varid)) forrest@0: if (read_grid_all) then forrest@0: call nc_err(nf90_inq_varid(ncid, aname, avarid)) forrest@0: call nc_err(nf90_inq_varid(ncid, lname, lvarid)) forrest@0: call nc_err(nf90_inq_varid(ncid, mname, mvarid)) forrest@0: call nc_err(nf90_inq_varid(ncid, oname, ovarid)) forrest@0: end if forrest@0: forrest@0: call nc_err(nf90_get_var(ncid, xvarid, xvals)) forrest@0: call nc_err(nf90_get_var(ncid, yvarid, yvals)) forrest@0: call nc_err(nf90_get_var(ncid, zsoi_varid, zsoi_vals)) forrest@0: call nc_err(nf90_get_var(ncid, zlak_varid, zlak_vals)) forrest@0: if (hrflag) call nc_err(nf90_get_var(ncid, hr_varid, hr_vals)) forrest@0: if (read_grid_all) then forrest@0: call nc_err(nf90_get_var(ncid, avarid, avals)) forrest@0: call nc_err(nf90_get_var(ncid, lvarid, lvals)) forrest@0: call nc_err(nf90_get_var(ncid, mvarid, mvals)) forrest@0: call nc_err(nf90_get_var(ncid, ovarid, ovals)) forrest@0: end if forrest@0: forrest@0: ! Read global attributes forrest@0: call nc_err(nf90_get_att(ncid, nf90_global, case_id_name, case_id)) forrest@0: call nc_err(nf90_get_att(ncid, nf90_global, inidat_name, inidat)) forrest@0: call nc_err(nf90_get_att(ncid, nf90_global, surdat_name, surdat)) forrest@0: call nc_err(nf90_get_att(ncid, nf90_global, pftdat_name, pftdat)) forrest@0: call nc_err(nf90_get_att(ncid, nf90_global, rtmdat_name, rtmdat)) forrest@0: forrest@0: call nc_err(nf90_close(ncid)) forrest@0: forrest@0: time_size = tlen forrest@0: ! Loop over all remaining files, reading size of time dimension forrest@0: do i = 2, numf forrest@0: print *, 'Opening ',trim(input_path)//'/'//trim(fnames(i)) forrest@0: call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(i)), & forrest@0: mode=nf90_nowrite, ncid=ncid)) forrest@0: call nc_err(nf90_inq_dimid(ncid, tname, tdimid)) forrest@0: call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=tdimid, len=tlen)) forrest@0: call nc_err(nf90_close(ncid)) forrest@0: print *, tlen,' values' forrest@0: time_size = time_size + tlen forrest@0: end do forrest@0: forrest@0: ! Allocate and fill the coordinates derived data type forrest@0: coord%xsize = xlen forrest@0: coord%ysize = ylen forrest@0: coord%zsoi_size = zsoi_len forrest@0: coord%zlak_size = zlak_len forrest@0: coord%tsize = time_size forrest@0: coord%hr_size = hr_len forrest@0: coord%nsamples = nsamples_len forrest@0: allocate(coord%x(xlen), coord%y(ylen), coord%zsoi(zsoi_len), coord%zlak(zlak_len), coord%hr(hr_len), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_read_coord: Cannot allocate grid variables' forrest@0: stop forrest@0: end if forrest@0: if (read_grid_all) then forrest@0: allocate(coord%area(xlen, ylen), coord%landfrac(xlen, ylen), coord%landmask(xlen, ylen), coord%orog(xlen, ylen), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_read_coord: Cannot allocate extra grid variables' forrest@0: stop forrest@0: end if forrest@0: end if forrest@0: coord%x = xvals forrest@0: coord%y = yvals forrest@0: coord%zsoi = zsoi_vals forrest@0: coord%zlak = zlak_vals forrest@0: if (hrflag) coord%hr = hr_vals forrest@0: if (read_grid_all) then forrest@0: coord%area = avals forrest@0: coord%landfrac = lvals forrest@0: coord%landmask = mvals forrest@0: coord%orog = ovals forrest@0: end if forrest@0: deallocate(xvals, yvals, zsoi_vals, zlak_vals, hr_vals) forrest@0: if (read_grid_all) then forrest@0: deallocate(avals, lvals, mvals, ovals) forrest@0: end if forrest@0: coord%case_id = trim(case_id) forrest@0: coord%inidat = trim(inidat) forrest@0: coord%surdat = trim(surdat) forrest@0: coord%pftdat = trim(pftdat) forrest@0: coord%rtmdat = trim(rtmdat) forrest@0: coord%hour_flag = hrflag forrest@0: forrest@0: end subroutine clm_read_coord forrest@0: ! forrest@0: !-------------------------------------------------------------------------- forrest@0: ! forrest@0: subroutine clm_copy_grid_data(varname) forrest@0: forrest@0: implicit none forrest@0: character (len=varnamelen), intent(in) :: varname ! variable name forrest@0: integer :: ierr forrest@0: forrest@0: var_data%varname = varname forrest@0: var_data%ndims = 2 forrest@0: var_data%missing_value = 1.0e+36 forrest@0: if (varname == 'landmask') then forrest@0: var_data%int_type = .true. forrest@0: allocate(var_data%int2d(coord%xsize, coord%ysize), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_copy_grid_data: Cannot allocate 2d integer variable' forrest@0: stop forrest@0: end if forrest@0: else forrest@0: var_data%int_type = .false. forrest@0: allocate(var_data%var2d(coord%xsize, coord%ysize), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_copy_grid_data: Cannot allocate 2d real variable' forrest@0: stop forrest@0: end if forrest@0: end if forrest@0: forrest@0: select case (varname) forrest@0: case ('area') forrest@0: var_data%var2d = coord%area forrest@0: case ('landfrac') forrest@0: var_data%var2d = coord%landfrac forrest@0: case ('landmask') forrest@0: var_data%int2d = coord%landmask forrest@0: case ('topo') forrest@0: var_data%var2d = coord%orog forrest@0: case default forrest@0: print *,'clm_copy_grid_data: ',varname,' is not a known grid variable' forrest@0: stop forrest@0: end select forrest@0: forrest@0: end subroutine clm_copy_grid_data forrest@0: ! forrest@0: !-------------------------------------------------------------------------- forrest@0: ! forrest@0: subroutine clm_read_data(numf, input_path, fnames, tshift, varname) forrest@0: forrest@0: implicit none forrest@0: integer, intent(in) :: numf ! number of filenames forrest@0: character (len=namelen), intent(in) :: input_path ! input path forrest@0: character (len=namelen), intent(in) :: fnames(numf) ! input filenames forrest@0: integer, intent(in) :: tshift(numf) ! time shift by file forrest@0: character (len=varnamelen), intent(in) :: varname ! variable name forrest@0: integer :: ncid forrest@0: integer :: tdimid, tvarid, tlen forrest@0: integer :: hidimid, hilen forrest@0: integer :: tbvarid forrest@0: integer :: varid, dimids(nf90_max_var_dims), dimlen forrest@0: integer :: zsize forrest@0: integer(i8) :: i, j, k, l forrest@0: integer(i8) :: time_size forrest@0: real(r4), allocatable :: tvals(:) ! temporary time values forrest@0: real(r4), allocatable :: tbvals(:,:) ! temporary time bounds values forrest@0: real(r4), allocatable :: vals(:,:) ! temporary data values forrest@0: character (len=varnamelen) :: dimname ! dimension name forrest@0: integer :: ierr forrest@0: forrest@0: ! Figure out the dimensionality of the variable and variable attributes forrest@0: ! from the first file forrest@0: print *, 'Opening ',trim(input_path)//'/'//trim(fnames(1)), & forrest@0: ' to check dimensionality of ',trim(varname) forrest@0: call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(1)), mode=nf90_nowrite, ncid=ncid)) forrest@0: call nc_err(nf90_inq_varid(ncid, trim(varname), varid)) forrest@0: call nc_err(nf90_inquire_variable(ncid, varid, name=var_data%varname, ndims=var_data%ndims, dimids=dimids)) forrest@0: ! Figure out if there is a Z-dimension and/or an hour-dimension forrest@0: ierr = 0 forrest@0: var_data%int_type = .false. forrest@0: var_data%hour_flag = .false. forrest@0: var_data%soil_layer_flag = .false. forrest@0: do i = 1, var_data%ndims forrest@0: call nc_err(nf90_inquire_dimension(ncid, dimids(i), name=dimname, len=dimlen)) forrest@0: if (i == 1 .and. dimname /= xname) ierr = ierr + 1 forrest@0: if (i == 2 .and. dimname /= yname) ierr = ierr + 1 forrest@0: if (i > 2 .and. i < var_data%ndims) then forrest@0: if (dimname == hr_name) then forrest@0: var_data%hour_flag = .true. forrest@0: else if (dimname == zsoi_name) then forrest@0: var_data%soil_layer_flag = .true. forrest@0: !else if (dimname == zlak_name) then forrest@0: ! var_data%zdim_name = dimname forrest@0: else forrest@0: ierr = ierr + 1 forrest@0: end if forrest@0: end if forrest@0: ! This is not a good assumption (that time is always there), so skip it forrest@0: !if (i == var_data%ndims .and. dimname /= tname) ierr = ierr + 1 forrest@0: print *,'Dimension ',i,' is ',dimname forrest@0: end do forrest@0: if (ierr /= 0) then forrest@0: print *,'clm_read_data: Unexpected dimension combination' forrest@0: stop forrest@0: end if forrest@0: ! forrest@0: call nc_err(nf90_get_att(ncid, varid, long_name_name, var_data%long_name)) forrest@0: call nc_err(nf90_get_att(ncid, varid, missing_value_name, var_data%missing_value)) forrest@0: call nc_err(nf90_get_att(ncid, varid, units_name, var_data%units)) forrest@0: call nc_err(nf90_close(ncid)) forrest@0: forrest@0: ! Allocate space for the variable of the correct dimensionality forrest@0: select case (var_data%ndims) forrest@0: case (2) forrest@0: allocate(var_data%var2d(coord%xsize, coord%ysize), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_read_data: Cannot allocate 2d variable' forrest@0: stop forrest@0: end if forrest@0: case (3) forrest@0: allocate(var_data%var3d(coord%xsize, coord%ysize, coord%tsize), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_read_data: Cannot allocate 3d variable' forrest@0: stop forrest@0: end if forrest@0: case (4) forrest@0: if (var_data%hour_flag) then forrest@0: zsize = coord%hr_size forrest@0: else if (var_data%soil_layer_flag) then forrest@0: zsize = coord%zsoi_size forrest@0: !else if (var_data%zdim_name == zlak_name) then forrest@0: ! zsize = coord%zlak_size forrest@0: else forrest@0: print *,'clm_read_data: Unexpected fourth dimension while allocating' forrest@0: stop forrest@0: end if forrest@0: allocate(var_data%var4d(coord%xsize, coord%ysize, zsize, coord%tsize), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_read_data: Cannot allocate 4d variable' forrest@0: stop forrest@0: end if forrest@0: case default forrest@0: print *, 'clm_read_data: Cannot handle variables of ', var_data%ndims, ' dimensions' forrest@0: stop forrest@0: end select forrest@0: forrest@0: ! Allocate space for time stamps and bounds forrest@0: allocate(var_data%time(coord%tsize), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_read_data: Cannot allocate time variable' forrest@0: stop forrest@0: end if forrest@0: allocate(var_data%time_bounds(2,coord%tsize), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_read_data: Cannot allocate time bounds variable' forrest@0: stop forrest@0: end if forrest@0: forrest@0: ! Allocate temporary space for the variable level-slice/time-slice */ forrest@0: allocate(vals(coord%xsize, coord%ysize), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_read_data: Cannot allocate variable time-slice' forrest@0: stop forrest@0: end if forrest@0: forrest@0: if (var_data%ndims == 2) then forrest@0: ! No time axis, so just read the first instance from the first file forrest@0: print *, 'Opening ',trim(input_path)//'/'//trim(fnames(1)) forrest@0: call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(1)), mode=nf90_nowrite, ncid=ncid)) forrest@0: call nc_err(nf90_inq_varid(ncid, trim(varname), varid)) forrest@0: call nc_err(nf90_get_var(ncid, varid, vals, start = (/ 1, 1 /), count = (/ coord%xsize, coord%ysize /))) forrest@0: var_data%var2d(:,:) = vals(:,:) forrest@0: call nc_err(nf90_close(ncid)) forrest@0: else forrest@0: ! Otherwise, loop over all times forrest@0: time_size = 0 forrest@0: ! Loop over all files, reading data forrest@0: do i = 1, numf forrest@0: print *, 'Opening ',trim(input_path)//'/'//trim(fnames(i)) forrest@0: call nc_err(nf90_open(path=trim(input_path)//'/'//trim(fnames(i)), mode=nf90_nowrite, ncid=ncid)) forrest@0: call nc_err(nf90_inq_dimid(ncid, tname, tdimid)) forrest@0: call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=tdimid, len=tlen)) forrest@0: call nc_err(nf90_inq_dimid(ncid, hiname, hidimid)) forrest@0: call nc_err(nf90_inquire_dimension(ncid=ncid, dimid=hidimid, len=hilen)) forrest@0: forrest@0: call nc_err(nf90_inq_varid(ncid, tname, tvarid)) forrest@0: call nc_err(nf90_inq_varid(ncid, tbname, tbvarid)) forrest@0: call nc_err(nf90_inq_varid(ncid, trim(varname), varid)) forrest@0: forrest@0: ! Allocate temporary space for the time stamps and bounds forrest@0: allocate(tvals(tlen), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_read_data: Cannot allocate temporary time variable' forrest@0: stop forrest@0: end if forrest@0: allocate(tbvals(hilen,tlen), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_read_data: Cannot allocate temporary time bounds variable' forrest@0: stop forrest@0: end if forrest@0: forrest@0: ! Get the time stamps and bounds forrest@0: call nc_err(nf90_get_var(ncid, tvarid, tvals)) forrest@0: call nc_err(nf90_get_var(ncid, tbvarid, tbvals)) forrest@0: forrest@0: do j = 1, tlen forrest@0: ! Get and store a slab of data forrest@0: select case (var_data%ndims) forrest@0: case (3) forrest@0: !print *, 'Reading slice of data' forrest@0: call nc_err(nf90_get_var(ncid, varid, vals, start = (/ 1, 1, j /), count = (/ coord%xsize, coord%ysize, 1 /))) forrest@0: var_data%var3d(:,:,time_size + j) = vals(:,:) forrest@0: case (4) forrest@0: do k = 1, zsize forrest@0: call nc_err(nf90_get_var(ncid, varid, vals, start = (/ 1, 1, k, j /), count = (/ coord%xsize, coord%ysize, 1, 1 /))) forrest@0: var_data%var4d(:,:,k,time_size + j) = vals(:,:) forrest@0: end do forrest@0: case default forrest@0: print *, 'clm_read_data: Cannot handle variables of ', var_data%ndims, ' dimensions' forrest@0: stop forrest@0: end select forrest@0: ! Store each time stamp and bounds, including the time offset forrest@0: var_data%time(time_size + j) = tvals(j) + tshift(i) forrest@0: var_data%time_bounds(:,time_size + j) = tbvals(:,j) + tshift(i) forrest@0: forrest@0: end do forrest@0: forrest@0: call nc_err(nf90_close(ncid)) forrest@0: forrest@0: deallocate(tvals) forrest@0: deallocate(tbvals) forrest@0: forrest@0: time_size = time_size + tlen forrest@0: forrest@0: end do forrest@0: end if forrest@0: forrest@0: deallocate(vals) forrest@0: forrest@0: end subroutine clm_read_data forrest@0: ! forrest@0: !-------------------------------------------------------------------------- forrest@0: ! forrest@0: subroutine clm_convert_data(casa_flux_bug) forrest@0: forrest@0: implicit none forrest@0: logical, intent(in) :: casa_flux_bug forrest@0: real(r8) :: sfactor forrest@0: real(r8) :: offset forrest@0: integer :: ierr forrest@0: logical, allocatable :: mask2d(:,:) forrest@0: logical, allocatable :: mask3d(:,:,:) forrest@0: logical, allocatable :: mask4d(:,:,:,:) forrest@0: forrest@0: ! Default offset is zero forrest@0: offset = 0. forrest@0: forrest@0: ! Default var_data structure values forrest@0: var_data%positive = '' forrest@0: forrest@0: ! Set scale factor according to which variable is to be converted forrest@0: select case (var_data%varname) forrest@0: forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! Atmospheric forcing forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! husf forrest@0: case ('QBOT') forrest@0: ! QBOT: kg kg-1 -> kg kg-1 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'kg kg-1' forrest@0: var_data%out_varname = 'husf' forrest@0: ! prra forrest@0: case ('RAIN') forrest@0: ! RAIN: mm s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'prra' forrest@0: ! prsn forrest@0: case ('SNOW') forrest@0: ! SNOW: mm s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'prsn' forrest@0: forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! Biogeochemistry forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! agbc forrest@0: ! aglbc forrest@0: ! agnpp forrest@0: case ('AGNPP') forrest@0: ! AGNPP: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'agnpp' forrest@0: ! ar forrest@0: case ('AR') forrest@0: ! AR: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'ar' forrest@0: ! bco forrest@0: case ('BIOGENCO') forrest@0: ! BIOGENCO: ug m-2 h-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-9 * (1.0 / 3600.) forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'bco' forrest@0: ! bgbc forrest@0: ! bgnpp forrest@0: case ('BGNPP') forrest@0: ! BGNPP: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'bgnpp' forrest@0: ! bipn forrest@0: case ('ISOPRENE') forrest@0: ! ISOPRENE: ug m-2 h-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-9 * (1.0 / 3600.) forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'bipn' forrest@0: ! bmtpn forrest@0: case ('MONOTERP') forrest@0: ! MONOTERP: ug m-2 h-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-9 * (1.0 / 3600.) forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'bmtpn' forrest@0: ! borvoc forrest@0: case ('ORVOC') forrest@0: ! ORVOC: ug m-2 h-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-9 * (1.0 / 3600.) forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'borvoc' forrest@0: ! bovoc forrest@0: case ('OVOC') forrest@0: ! OVOC: ug m-2 h-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-9 * (1.0 / 3600.) forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'bovoc' forrest@0: ! cwdc forrest@0: case ('CWDC') forrest@0: ! CWDC: g m-2 -> kg m-2 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2' forrest@0: var_data%out_varname = 'cwdc' forrest@0: ! cwdc_hr forrest@0: case ('CWDC_HR') forrest@0: ! CWDC_HR: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'cwdc_hr' forrest@0: ! cwdc_loss forrest@0: case ('CWDC_LOSS') forrest@0: ! CWDC_LOSS: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: if (casa_flux_bug) then forrest@0: print *,'WARNING: casa_flux_bug set to true; dividing data by 1200' forrest@0: sfactor = sfactor * (1. / 1200.) forrest@0: end if forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'cwdc_loss' forrest@0: ! froot_litterc forrest@0: ! froot_litterc_hr forrest@0: ! froot_litterc_loss forrest@0: ! frootc forrest@0: case ('FROOTC') forrest@0: ! FROOTC: g m-2 -> kg m-2 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2' forrest@0: var_data%out_varname = 'frootc' forrest@0: ! frootc_alloc forrest@0: case ('FROOTC_ALLOC') forrest@0: ! FROOTC_ALLOC: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: if (casa_flux_bug) then forrest@0: print *,'WARNING: casa_flux_bug set to true; dividing data by 1200' forrest@0: sfactor = sfactor * (1. / 1200.) forrest@0: end if forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'frootc_alloc' forrest@0: ! frootc_loss forrest@0: case ('FROOTC_LOSS') forrest@0: ! FROOTC_LOSS: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: if (casa_flux_bug) then forrest@0: print *,'WARNING: casa_flux_bug set to true; dividing data by 1200' forrest@0: sfactor = sfactor * (1. / 1200.) forrest@0: end if forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'frootc_loss' forrest@0: ! gpp forrest@0: case ('GPP') forrest@0: ! GPP: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'gpp' forrest@0: ! gnmin forrest@0: case ('GROSS_NMIN') forrest@0: ! GROSS_NMIN: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'gnmin' forrest@0: ! hr forrest@0: case ('HR') forrest@0: ! HR: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'hr' forrest@0: ! leaf_litterc forrest@0: ! leaf_litterc_hr forrest@0: ! leaf_litterc_loss forrest@0: ! leafc forrest@0: case ('LEAFC') forrest@0: ! LEAFC: g m-2 -> kg m-2 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2' forrest@0: var_data%out_varname = 'leafc' forrest@0: ! leafc_alloc forrest@0: case ('LEAFC_ALLOC') forrest@0: ! LEAFC_ALLOC: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: if (casa_flux_bug) then forrest@0: print *,'WARNING: casa_flux_bug set to true; dividing data by 1200' forrest@0: sfactor = sfactor * (1. / 1200.) forrest@0: end if forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'leafc_alloc' forrest@0: ! leafc_loss forrest@0: case ('LEAFC_LOSS') forrest@0: ! LEAFC_LOSS: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: if (casa_flux_bug) then forrest@0: print *,'WARNING: casa_flux_bug set to true; dividing data by 1200' forrest@0: sfactor = sfactor * (1. / 1200.) forrest@0: end if forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'leafc_loss' forrest@0: ! litterc forrest@0: case ('LITTERC') forrest@0: ! LITTERC: g m-2 -> kg m-2 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2' forrest@0: var_data%out_varname = 'litterc' forrest@0: ! litterc_hr forrest@0: case ('LITTERC_HR') forrest@0: ! LITTERC_HR: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'litterc_hr' forrest@0: ! litterc_loss forrest@0: case ('LITTERC_LOSS') forrest@0: ! LITTERC_LOSS: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: if (casa_flux_bug) then forrest@0: print *,'WARNING: casa_flux_bug set to true; dividing data by 1200' forrest@0: sfactor = sfactor * (1. / 1200.) forrest@0: end if forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'litterc_loss' forrest@0: ! nee forrest@0: case ('NEE') forrest@0: ! NEE: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'nee' forrest@0: ! nep forrest@0: case ('NEP') forrest@0: ! NEP: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'nep' forrest@0: ! nnmin forrest@0: case ('NET_NMIN') forrest@0: ! NEP: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'nnmin' forrest@0: ! npp forrest@0: case ('NPP') forrest@0: ! NPP: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'npp' forrest@0: ! soilc forrest@0: case ('SOILC') forrest@0: ! SOILC: g m-2 -> kg m-2 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2' forrest@0: var_data%out_varname = 'soilc' forrest@0: ! soilc_hr forrest@0: case ('SOILC_HR') forrest@0: ! SOILC_HR: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'soilc_hr' forrest@0: ! soilc_loss forrest@0: case ('SOILC_LOSS') forrest@0: ! SOILC_LOSS: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: if (casa_flux_bug) then forrest@0: print *,'WARNING: casa_flux_bug set to true; dividing data by 1200' forrest@0: sfactor = sfactor * (1. / 1200.) forrest@0: end if forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'soilc_loss' forrest@0: ! tbvoc forrest@0: case ('VOCFLXT') forrest@0: ! VOCFLXT: ug m-2 h-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-9 * (1.0 / 3600.) forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'tbvoc' forrest@0: ! woodc forrest@0: case ('WOODC') forrest@0: ! WOODC: g m-2 -> kg m-2 forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'kg m-2' forrest@0: var_data%out_varname = 'woodc' forrest@0: ! woodc_alloc forrest@0: case ('WOODC_ALLOC') forrest@0: ! WOODC_ALLOC: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: if (casa_flux_bug) then forrest@0: print *,'WARNING: casa_flux_bug set to true; dividing data by 1200' forrest@0: sfactor = sfactor * (1. / 1200.) forrest@0: end if forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'woodc_alloc' forrest@0: ! woodc_loss forrest@0: case ('WOODC_LOSS') forrest@0: ! WOODC_LOSS: g m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e-3 forrest@0: if (casa_flux_bug) then forrest@0: print *,'WARNING: casa_flux_bug set to true; dividing data by 1200' forrest@0: sfactor = sfactor * (1. / 1200.) forrest@0: end if forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'woodc_loss' forrest@0: forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! Energy forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! hfg forrest@0: case ('FGR') forrest@0: ! FGR: W m-2 -> W m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'W m-2' forrest@0: var_data%positive = 'down' forrest@0: var_data%out_varname = 'hfg' forrest@0: ! hfls forrest@0: case ('LATENT') forrest@0: ! LATENT: W m-2 -> W m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'W m-2' forrest@0: var_data%positive = 'up' forrest@0: var_data%out_varname = 'hfls' forrest@0: ! hfss forrest@0: case ('FSH') forrest@0: ! FSH: W m-2 -> W m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'W m-2' forrest@0: var_data%positive = 'up' forrest@0: var_data%out_varname = 'hfss' forrest@0: ! hfssg forrest@0: case ('FSH_G') forrest@0: ! FSH_G: W m-2 -> W m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'W m-2' forrest@0: var_data%positive = 'up' forrest@0: var_data%out_varname = 'hfssg' forrest@0: ! hfssv forrest@0: case ('FSH_V') forrest@0: ! FSH_V: W m-2 -> W m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'W m-2' forrest@0: var_data%positive = 'up' forrest@0: var_data%out_varname = 'hfssv' forrest@0: forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! Grid forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! gca forrest@0: case ('area') forrest@0: ! area: km2 -> m2 forrest@0: sfactor = 1.0e+3 forrest@0: var_data%units = 'm2' forrest@0: var_data%out_varname = 'gca' forrest@0: ! lbm forrest@0: case ('landmask') forrest@0: ! landmask: 1 -> 1 forrest@0: var_data%units = '1' forrest@0: var_data%out_varname = 'lbm' forrest@0: ! Do not convert integer landmask, just return forrest@0: return forrest@0: ! lcf forrest@0: ! lcn forrest@0: ! orog forrest@0: case ('topo') forrest@0: ! topo: m -> m forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'm' forrest@0: var_data%out_varname = 'orog' forrest@0: ! sftlf forrest@0: case ('landfrac') forrest@0: ! landfrac: 1 -> 1 forrest@0: sfactor = 1.0 forrest@0: var_data%units = '1' forrest@0: var_data%out_varname = 'sftlf' forrest@0: forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! Humidity forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! hur2m forrest@0: case ('RH2M') forrest@0: ! RH2M: % -> % forrest@0: sfactor = 1.0 forrest@0: var_data%units = '%' forrest@0: var_data%out_varname = 'hur2m' forrest@0: ! hus2m forrest@0: case ('Q2M') forrest@0: ! Q2M: kg kg-1 -> kg kg-1 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'kg kg-1' forrest@0: var_data%out_varname = 'hus2m' forrest@0: forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! Hydrology forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! btran forrest@0: case ('BTRAN') forrest@0: ! BTRAN: 1 -> 1 forrest@0: sfactor = 1.0 forrest@0: var_data%units = '1' forrest@0: var_data%out_varname = 'btran' forrest@0: ! et forrest@0: case ('ET') forrest@0: ! ET: m s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0e+3 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'et' forrest@0: ! evspsblsoi forrest@0: case ('QSOIL') forrest@0: ! QSOIL: mm s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'evspsblsoi' forrest@0: ! evspsblveg forrest@0: case ('QVEGE') forrest@0: ! QVEGE: mm s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'evspsblveg' forrest@0: ! mrfsl (multi-level) forrest@0: case ('SOILICE') forrest@0: ! SOILICE: kg m-2 -> kg m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'kg m-2' forrest@0: var_data%out_varname = 'mrfsl' forrest@0: ! mrfso forrest@0: case ('TOTAL_SOILICE') forrest@0: ! TOTAL_SOILICE: kg m-2 -> kg m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'kg m-2' forrest@0: var_data%out_varname = 'mrfso' forrest@0: ! mrlsl (multi-level) forrest@0: case ('SOILLIQ') forrest@0: ! SOILLIQ: kg m-2 -> kg m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'kg m-2' forrest@0: var_data%out_varname = 'mrlsl' forrest@0: ! mrlso forrest@0: case ('TOTAL_SOILLIQ') forrest@0: ! TOTAL_SOILLIQ: kg m-2 -> kg m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'kg m-2' forrest@0: var_data%out_varname = 'mrlso' forrest@0: ! mrros forrest@0: case ('QOVER') forrest@0: ! QOVER: mm s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'mrros' forrest@0: ! mrross forrest@0: case ('QDRAI') forrest@0: ! QDRAI: mm s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'mrross' forrest@0: ! prveg forrest@0: case ('QINTR') forrest@0: ! QINTR: mm s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'prveg' forrest@0: ! snd forrest@0: case ('H2OSNO') forrest@0: ! H2OSNO: mm -> m forrest@0: sfactor = 1.0e-3 forrest@0: var_data%units = 'm' forrest@0: var_data%out_varname = 'snd' forrest@0: ! tran forrest@0: case ('QVEGT') forrest@0: ! QVEGT: mm s-1 -> kg m-2 s-1 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'tran' forrest@0: ! vfwsl (multi-level) forrest@0: case ('H2OSOI') forrest@0: ! H2OSOI: mm3 mm-3 -> m3 m-3 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'm3 m-3' forrest@0: var_data%out_varname = 'vfwsl' forrest@0: ! wpsl (multi-level) forrest@0: case ('SOILPSI') forrest@0: ! SOILPSI: Pa -> Pa forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'Pa' forrest@0: var_data%out_varname = 'wpsl' forrest@0: forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! Radiation forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! alsas forrest@0: case ('ALL_SKY_ALBEDO') forrest@0: ! ALL_SKY_ALBEDO: 1 -> 1 forrest@0: sfactor = 1.0 forrest@0: var_data%units = '1' forrest@0: var_data%out_varname = 'alsas' forrest@0: ! alsbs forrest@0: case ('BLACK_SKY_ALBEDO') forrest@0: ! BLACK_SKY_ALBEDO: 1 -> 1 forrest@0: sfactor = 1.0 forrest@0: var_data%units = '1' forrest@0: var_data%out_varname = 'alsbs' forrest@0: ! laisha forrest@0: case ('LAISHA') forrest@0: ! LAISHA: m2 m-2 -> m2 m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'm2 m-2' forrest@0: var_data%out_varname = 'laisha' forrest@0: ! laisun forrest@0: case ('LAISUN') forrest@0: ! LAISUN: m2 m-2 -> m2 m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'm2 m-2' forrest@0: var_data%out_varname = 'laisun' forrest@0: ! rlds forrest@0: case ('FLDS') forrest@0: ! FLDS: W m-2 -> W m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'W m-2' forrest@0: var_data%positive = 'down' forrest@0: var_data%out_varname = 'rlds' forrest@0: ! rls forrest@0: case ('FIRA') forrest@0: ! FIRA: W m-2 -> W m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'W m-2' forrest@0: var_data%positive = 'down' forrest@0: var_data%out_varname = 'rls' forrest@0: ! rlus forrest@0: case ('FIRE') forrest@0: ! FIRE: W m-2 -> W m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'W m-2' forrest@0: var_data%positive = 'up' forrest@0: var_data%out_varname = 'rlus' forrest@0: ! rs forrest@0: case ('NETRAD') forrest@0: ! NETRAD: W m-2 -> W m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'W m-2' forrest@0: var_data%positive = 'down' forrest@0: var_data%out_varname = 'rs' forrest@0: ! rsas forrest@0: case ('FSA') forrest@0: ! FSA: W m-2 -> W m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'W m-2' forrest@0: var_data%positive = 'down' forrest@0: var_data%out_varname = 'rsas' forrest@0: ! rsasg forrest@0: case ('SABG') forrest@0: ! SABG: W m-2 -> W m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'W m-2' forrest@0: var_data%positive = 'down' forrest@0: var_data%out_varname = 'rsasg' forrest@0: ! rsasv forrest@0: case ('SABV') forrest@0: ! SABV: W m-2 -> W m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'W m-2' forrest@0: var_data%positive = 'down' forrest@0: var_data%out_varname = 'rsasv' forrest@0: ! rsds forrest@0: case ('FSDS') forrest@0: ! FSDS: W m-2 -> W m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'W m-2' forrest@0: var_data%positive = 'down' forrest@0: var_data%out_varname = 'rsds' forrest@0: ! rsus forrest@0: case ('FSR') forrest@0: ! FSR: W m-2 -> W m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'W m-2' forrest@0: var_data%positive = 'up' forrest@0: var_data%out_varname = 'rsus' forrest@0: ! snc forrest@0: case ('FSNO') forrest@0: ! FSNO: 1 -> 1 forrest@0: sfactor = 1.0 forrest@0: var_data%units = '1' forrest@0: var_data%out_varname = 'snc' forrest@0: forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! Temperature forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! degday forrest@0: case ('DEGDAY') forrest@0: ! DEGDAY: C -> K forrest@0: sfactor = 1.0 forrest@0: offset = 273.15 forrest@0: var_data%units = 'K' forrest@0: var_data%out_varname = 'degday' forrest@0: ! ta2m forrest@0: case ('TSA') forrest@0: ! TSA: K -> K forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'K' forrest@0: var_data%out_varname = 'ta2m' forrest@0: ! ta2mdmax forrest@0: case ('TREFMXAV') forrest@0: ! TREFMXAV: K -> K forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'K' forrest@0: var_data%out_varname = 'ta2mdmax' forrest@0: ! ta2mdmin forrest@0: case ('TREFMNAV') forrest@0: ! TREFMNAV: K -> K forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'K' forrest@0: var_data%out_varname = 'ta2mdmin' forrest@0: ! tsl (multi-level) forrest@0: case ('TSOI') forrest@0: ! TSOI: K -> K forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'K' forrest@0: var_data%out_varname = 'tsl' forrest@0: ! tv forrest@0: case ('TV') forrest@0: ! TV: K -> K forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'K' forrest@0: var_data%out_varname = 'tv' forrest@0: forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! Vegetation phenology forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! laiexp forrest@0: case ('ELAI') forrest@0: ! ELAI: m2 m-2 -> m2 m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'm2 m-2' forrest@0: var_data%out_varname = 'laiexp' forrest@0: ! laitot forrest@0: case ('TLAI') forrest@0: ! TLAI: m2 m-2 -> m2 m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'm2 m-2' forrest@0: var_data%out_varname = 'laitot' forrest@0: ! saiexp forrest@0: case ('ESAI') forrest@0: ! ESAI: m2 m-2 -> m2 m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'm2 m-2' forrest@0: var_data%out_varname = 'saiexp' forrest@0: ! saitot forrest@0: case ('TSAI') forrest@0: ! TSAI: m2 m-2 -> m2 m-2 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 'm2 m-2' forrest@0: var_data%out_varname = 'saitot' forrest@0: forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! Vegetation physiology forrest@0: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! forrest@0: ! psn forrest@0: case ('FPSN') forrest@0: ! FPSN: umol m-2 s-1 -> kg m-2 s-1 forrest@0: sfactor = 12.0 * 1.0e-6 * 1.0e-3 forrest@0: var_data%units = 'kg m-2 s-1' forrest@0: var_data%out_varname = 'psn' forrest@0: ! rssha forrest@0: case ('RSSHA') forrest@0: ! RSSHA: s m-1 -> s m-1 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 's m-1' forrest@0: var_data%out_varname = 'rssha' forrest@0: ! rssun forrest@0: case ('RSSUN') forrest@0: ! RSSUN: s m-1 -> s m-1 forrest@0: sfactor = 1.0 forrest@0: var_data%units = 's m-1' forrest@0: var_data%out_varname = 'rssun' forrest@0: forrest@0: case default forrest@0: print *, 'clm_convert_data: Unknown variable name ', var_data%varname forrest@0: stop forrest@0: end select forrest@0: forrest@0: ! Convert to the desired units using the offset and scale factor (sfactor) forrest@0: select case (var_data%ndims) forrest@0: case (2) forrest@0: allocate(mask2d(coord%xsize, coord%ysize), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_convert_data: Cannot allocate 2d mask' forrest@0: stop forrest@0: end if forrest@0: mask2d = .false. forrest@0: if (var_data%int_type) then forrest@0: where (var_data%int2d /= var_data%missing_value) forrest@0: var_data%int2d = int(var_data%int2d * sfactor + offset) forrest@0: mask2d = .true. forrest@0: end where forrest@0: print *, var_data%out_varname, & forrest@0: ' min:', minval(var_data%int2d, mask=mask2d), & forrest@0: ' max:', maxval(var_data%int2d, mask=mask2d) forrest@0: else forrest@0: where (var_data%var2d /= var_data%missing_value) forrest@0: var_data%var2d = var_data%var2d * sfactor + offset forrest@0: mask2d = .true. forrest@0: end where forrest@0: print *, var_data%out_varname, & forrest@0: ' min:', minval(var_data%var2d, mask=mask2d), & forrest@0: ' max:', maxval(var_data%var2d, mask=mask2d) forrest@0: end if forrest@0: if (allocated(mask2d)) deallocate(mask2d) forrest@0: case (3) forrest@0: allocate(mask3d(coord%xsize, coord%ysize, coord%tsize), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_convert_data: Cannot allocate 3d mask' forrest@0: stop forrest@0: end if forrest@0: mask3d = .false. forrest@0: where (var_data%var3d /= var_data%missing_value) forrest@0: var_data%var3d = (var_data%var3d + offset) * sfactor forrest@0: mask3d = .true. forrest@0: end where forrest@0: print *, var_data%out_varname, & forrest@0: ' min:', minval(var_data%var3d, mask=mask3d), & forrest@0: ' max:', maxval(var_data%var3d, mask=mask3d) forrest@0: if (allocated(mask3d)) deallocate(mask3d) forrest@0: case (4) forrest@0: if (var_data%soil_layer_flag) then forrest@0: allocate(mask4d(coord%xsize, coord%ysize, coord%zsoi_size, coord%tsize), stat=ierr) forrest@0: if (ierr /= 0) then forrest@0: print *, 'clm_convert_data: Cannot allocate 4d mask' forrest@0: stop forrest@0: end if forrest@0: mask4d = .false. forrest@0: where (var_data%var4d /= var_data%missing_value) forrest@0: var_data%var4d = (var_data%var4d + offset) * sfactor forrest@0: mask4d = .true. forrest@0: end where forrest@0: print *, var_data%out_varname, & forrest@0: ' min:', minval(var_data%var4d, mask=mask4d), & forrest@0: ' max:', maxval(var_data%var4d, mask=mask4d) forrest@0: if (allocated(mask4d)) deallocate(mask4d) forrest@0: else forrest@0: print *,'clm_convert_data: Unknown data dimension' forrest@0: stop forrest@0: end if forrest@0: case default forrest@0: print *,'clm_convert_data: Dimension count not supported' forrest@0: stop forrest@0: end select forrest@0: forrest@0: end subroutine clm_convert_data forrest@0: ! forrest@0: !-------------------------------------------------------------------------- forrest@0: ! forrest@0: subroutine clm_free_data forrest@0: forrest@0: implicit none forrest@0: forrest@0: var_data%varname = '' forrest@0: var_data%ndims = 0 forrest@0: if (allocated(var_data%int2d)) deallocate(var_data%int2d) forrest@0: if (allocated(var_data%var2d)) deallocate(var_data%var2d) forrest@0: if (allocated(var_data%var3d)) deallocate(var_data%var3d) forrest@0: if (allocated(var_data%var4d)) deallocate(var_data%var4d) forrest@0: if (allocated(var_data%time)) deallocate(var_data%time) forrest@0: if (allocated(var_data%time_bounds)) deallocate(var_data%time_bounds) forrest@0: var_data%int_type = .false. forrest@0: var_data%soil_layer_flag = .false. forrest@0: var_data%hour_flag = .false. forrest@0: var_data%long_name = '' forrest@0: var_data%missing_value = 0. forrest@0: var_data%units = '' forrest@0: var_data%positive = '' forrest@0: var_data%out_varname = '' forrest@0: forrest@0: end subroutine clm_free_data forrest@0: ! forrest@0: !-------------------------------------------------------------------------- forrest@0: ! forrest@0: subroutine nc_err(status) forrest@0: forrest@0: implicit none forrest@0: integer :: status ! status code forrest@0: forrest@0: if (status /= nf90_noerr) then forrest@0: print *, nf90_strerror(status) forrest@0: stop 'Program stopped after encountering netCDF error.' forrest@0: end if forrest@0: forrest@0: end subroutine nc_err forrest@0: ! forrest@0: !-------------------------------------------------------------------------- forrest@0: ! forrest@0: forrest@0: end module clm_mod