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