print_update_clm_5.F90 Source File


Source Code

!-------------------------------------------------------------------------------------------
!Copyright (c) 2013-2016 by Wolfgang Kurtz and Guowei He (Forschungszentrum Juelich GmbH)
!
!This file is part of TSMP-PDAF
!
!TSMP-PDAF is free software: you can redistribute it and/or modify
!it under the terms of the GNU Lesser General Public License as published by
!the Free Software Foundation, either version 3 of the License, or
!(at your option) any later version.
!
!TSMP-PDAF is distributed in the hope that it will be useful,
!but WITHOUT ANY WARRANTY; without even the implied warranty of
!MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
!GNU LesserGeneral Public License for more details.
!
!You should have received a copy of the GNU Lesser General Public License
!along with TSMP-PDAF.  If not, see <http://www.gnu.org/licenses/>.
!-------------------------------------------------------------------------------------------
!
!
!-------------------------------------------------------------------------------------------
!print_update_clm_5.F90: Module for printing updated CLM5 ensemble
!-------------------------------------------------------------------------------------------

#if defined CLMSA
subroutine print_update_clm(ts,ttot) bind(C,name="print_update_clm")
    
    use iso_c_binding   
    use shr_kind_mod , only : r8 => shr_kind_r8
    use subgridavemod, only : p2g, c2g
    use domainMod    , only : ldomain 
    use clm_varpar   , only : nlevsoi
    use clm_varcon   , only : nameg 
    use decompmod    , only : get_proc_global, get_proc_bounds, ldecomp
    use spmdgathscatmod , only : gather_data_to_master
    use spmdmod      , only : masterproc
    use clm_time_manager        , only : get_nstep    
    use clm_instMod, only : soilstate_inst, waterstate_inst
    use netcdf
    use enkf_clm_mod, only : clmupdate_swc,clmupdate_texture,clmprint_swc

    implicit none

    integer(c_int), intent(in) :: ts,ttot

    ! *** local variables ***
    integer :: numg           ! total number of gridcells across all processors
    integer :: numl           ! total number of landunits across all processors
    integer :: numc           ! total number of columns across all processors
    integer :: nump           ! total number of pfts across all processors
    integer :: begg,endg      ! local beg/end gridcells gdc
    integer :: begl,endl      ! local beg/end landunits
    integer :: begc,endc      ! local beg/end columns
    integer :: begp,endp      ! local beg/end pfts
    
    integer ::   isec, info, jn, jj, ji, g1, jx    ! temporary integer
    real(r8), pointer :: swc(:,:)
    real(r8), pointer :: psand(:,:)
    real(r8), pointer :: pclay(:,:)
    real(r8), pointer :: porgm(:,:)
    ! real(r8), pointer :: clmstate_tmp_local(:)
    ! real(r8), pointer :: clmstate_tmp_global(:)
    real(r8), allocatable :: clmstate_out(:,:,:)
    integer ,dimension(4) :: dimids
    integer ,dimension(1) :: il_var_id
    integer :: il_file_id, ncvarid(4), status
    character(len = 300) :: update_filename
    integer :: nerror
    integer :: ndlon,ndlat


    call get_proc_global(ng=numg,nl=numl,nc=numc,np=nump)
    call get_proc_bounds(begg,endg,begl,endl,begc,endc,begp,endp)
    ! allocate(clmstate_tmp_local(nlevsoi*(-begc+endc)), stat=nerror)

    ndlon  = ldomain%ni
    ndlat  = ldomain%nj

    if (masterproc) then
      ! allocate(clmstate_tmp_global(nlevsoi*numg), stat=nerror)
      allocate(clmstate_out(ndlon,ndlat,nlevsoi), stat=nerror)
    end if

    if(masterproc) then
      call get_update_filename(update_filename)
      if(ts.eq.1) then
        status =  nf90_create(update_filename, NF90_CLOBBER, il_file_id)
        status =  nf90_def_dim(il_file_id, "x", ndlon, dimids(1))
        status =  nf90_def_dim(il_file_id, "y", ndlat, dimids(2))
        status =  nf90_def_dim(il_file_id, "z", nlevsoi, dimids(3))
        status =  nf90_def_dim(il_file_id, "t", ttot, dimids(4))

        if(clmprint_swc.eq.1) then
          status =  nf90_def_var(il_file_id, "swc", NF90_DOUBLE, dimids, ncvarid(1))
        endif
 
        if(clmupdate_texture.eq.1) then 
          status =  nf90_def_var(il_file_id, "sand", NF90_DOUBLE, dimids, ncvarid(2))
          status =  nf90_def_var(il_file_id, "clay", NF90_DOUBLE, dimids, ncvarid(3))
        endif

        ! write updates to sand, clay and organic matter
        if(clmupdate_texture.eq.2) then
          status =  nf90_def_var(il_file_id, "sand", NF90_DOUBLE, dimids, ncvarid(2))
          status =  nf90_def_var(il_file_id, "clay", NF90_DOUBLE, dimids, ncvarid(3))
          status =  nf90_def_var(il_file_id, "orgm", NF90_DOUBLE, dimids, ncvarid(4))
        endif
        status =  nf90_enddef(il_file_id)
      else
        status = nf90_open(update_filename,NF90_WRITE,il_file_id)
      endif
    endif
  
    
    if(clmprint_swc.eq.1) then
      swc  => waterstate_inst%h2osoi_vol_col
      ! swc
!      clmstate_tmp_local = pack(swc,.true.)
!      clmstate_tmp_global = clmstate_tmp_local
!      call gather_data_to_master(clmstate_tmp_local,clmstate_tmp_global,clmlevel=nameg)

      if(masterproc) then
        do jn = 1, nlevsoi
          do g1 = 1, numg
             ji = mod(ldecomp%gdc2glo(g1)-1,ldomain%ni) + 1
             jj = (ldecomp%gdc2glo(g1) - 1)/ldomain%ni + 1
             clmstate_out(ji,jj,jn) = swc(g1, jn) !clmstate_tmp_global(jn+(g1-1)*nlevsoi)
          end do
        end do
        status = nf90_inq_varid(il_file_id, "swc" , ncvarid(1))
        status = nf90_put_var( il_file_id, ncvarid(1), clmstate_out(:,:,:), &
                 start = (/ 1, 1, 1, ts/), count = (/ ndlon, ndlat, nlevsoi, 1/) )
        !status = nf90_close(il_file_id)
      end if
    end if
    
    if((clmupdate_texture.eq.1) .or. (clmupdate_texture.eq.2)) then
      psand => soilstate_inst%cellsand_col
      pclay => soilstate_inst%cellclay_col

      ! sand
      !clmstate_tmp_local = pack(psand,.true.)
      !clmstate_tmp_global = clmstate_tmp_local
!      call gather_data_to_master(clmstate_tmp_local,clmstate_tmp_global, clmlevel=nameg)
      if(masterproc) then
        do jn = 1, nlevsoi
          do g1 = 1, numg
             ji = mod(ldecomp%gdc2glo(g1)-1,ldomain%ni) + 1
             jj = (ldecomp%gdc2glo(g1) - 1)/ldomain%ni + 1
             clmstate_out(ji,jj,jn) = psand(g1,jn)!clmstate_tmp_global(jn+(g1-1)*nlevsoi)
          end do
        end do
        status = nf90_inq_varid(il_file_id, "sand" , ncvarid(2))
        status = nf90_put_var( il_file_id, ncvarid(2), clmstate_out(:,:,:), &
                 start = (/ 1, 1, 1, ts/), count = (/ ndlon, ndlat, nlevsoi, 1/) )
        !status = nf90_close(il_file_id)
      end if

      ! clay
     ! clmstate_tmp_local = pack(pclay,.true.)
     ! clmstate_tmp_global = clmstate_tmp_local
!      call gather_data_to_master(clmstate_tmp_local,clmstate_tmp_global,clmlevel=nameg)

      if(masterproc) then
        do jn = 1, nlevsoi
          do g1 = 1, numg
             ji = mod(ldecomp%gdc2glo(g1)-1,ldomain%ni) + 1
             jj = (ldecomp%gdc2glo(g1) - 1)/ldomain%ni + 1
             clmstate_out(ji,jj,jn) = pclay(g1, jn) !clmstate_tmp_global(jn+(g1-1)*nlevsoi)
          end do
        end do
        status = nf90_inq_varid(il_file_id, "clay" , ncvarid(3))
        status = nf90_put_var( il_file_id, ncvarid(3), clmstate_out(:,:,:), &
                 start = (/ 1, 1, 1, ts/), count = (/ ndlon, ndlat, nlevsoi, 1/) )
        !status = nf90_close(il_file_id)
      end if
      
      ! organic matter
      if(clmupdate_texture.eq.2) then
        porgm => soilstate_inst%cellorg_col

        !clmstate_tmp_local = pack(porgm,.true.)
        !clmstate_tmp_global = clmstate_tmp_local
!        call gather_data_to_master(clmstate_tmp_local,clmstate_tmp_global, clmlevel=nameg)

        if(masterproc) then
          do jn = 1, nlevsoi
            do g1 = 1, numg
             ji = mod(ldecomp%gdc2glo(g1)-1,ldomain%ni) + 1
             jj = (ldecomp%gdc2glo(g1) - 1)/ldomain%ni + 1
             clmstate_out(ji,jj,jn) = porgm(g1, jn) !clmstate_tmp_global(jn+(g1-1)*nlevsoi)
            end do
          end do
          status = nf90_inq_varid(il_file_id, "orgm" , ncvarid(4))
          status = nf90_put_var( il_file_id, ncvarid(4), clmstate_out(:,:,:), &
                   start = (/ 1, 1, 1, ts/), count = (/ ndlon, ndlat, nlevsoi, 1/) )
          !status = nf90_close(il_file_id)
        end if
      endif

    end if 
    
    if(masterproc) then
      status = nf90_close(il_file_id)
      deallocate(clmstate_out)
      ! deallocate(clmstate_tmp_global)
    end if

    ! deallocate(clmstate_tmp_local)

end subroutine print_update_clm
#endif

subroutine get_update_filename (iofile)
    use clm_varctl, only : caseid
    use clm_time_manager, only : get_curr_date, get_prev_date
    ! !ARGUMENTS:
    implicit none
    character(len=300),intent(inout) :: iofile
    ! LOCAL VARIABLES:
    character(len=256) :: cdate       !date char string
    integer :: day                    !day (1 -> 31)
    integer :: mon                    !month (1 -> 12)
    integer :: yr                     !year (0 -> ...)
    integer :: sec                    !seconds into current day
    !-----------------------------------------------------------------------

    call get_prev_date (yr, mon, day, sec)
    write(cdate,'(i4.4,"-",i2.2)') yr,mon 
    call get_curr_date (yr, mon, day, sec)
    !write(cdate,'(i4.4,"-",i2.2,"-",i2.2,"-",i5.5)') yr,mon,day,sec
    write(cdate,'(i4.4)') yr
    iofile = trim(caseid)//".update."//trim(cdate)//".nc"
    !iofile = trim(caseid)//".update.nc"
end subroutine get_update_filename