convert_soitex.f90

CCSM Research Tools: Community Atmosphere Model (CAM)

  name = 'percent clay'
  unit = 'unitless'
  call wrap_def_var (ncid ,'PCT_CLAY' ,nf_float, 2, dim2_id, pct_clay_id)
  call wrap_put_att_text (ncid, pct_clay_id, 'long_name', name)
  call wrap_put_att_text (ncid, pct_clay_id, 'units'    , unit)

! End of definition

  status = nf_enddef(ncid)

! Read in formatted surface data

  open (unit=ndata,file=trim(filei),status='old',form='formatted',iostat=status)
  if (status .ne. 0) then
     write (6,*)'failed to open ',trim(filei),' on unit ',ndata,' ierr=',status
     stop
  end if
  do j = 1, nlat
     if (lat(j) <= 84. .and. lat(j) >= -56.5) then
        read (ndata,*) (mapunit(i,j),i=1,nlon)
        do i = 1, nlon
           if (mapunit(i,j) ==    0. .or. mapunit(i,j) ==  794. .or. &
               mapunit(i,j) == 1972. .or. mapunit(i,j) == 3214. .or. &
               mapunit(i,j) == 6997. .or. mapunit(i,j) == 6998.) then
              landmask(i,j) = 0. !ocean, no soil data, lakes, glaciers
           else
              landmask(i,j) = 1.
           end if
        end do
     else
        landmask(i,j) = 0.
        mapunit(i,j) = 0.
     end if
  end do
  close(ndata)

  open (unit=11,file=trim(filei1),status='old')
  open (unit=12,file=trim(filei2),status='old')
  open (unit=13,file=trim(filei3),status='old')
  open (unit=14,file=trim(filei4),status='old')
  open (unit=15,file=trim(filei5),status='old')
  open (unit=16,file=trim(filei6),status='old')
  open (unit=17,file=trim(filei7),status='old')
  open (unit=18,file=trim(filei8),status='old')
  open (unit=19,file=trim(filei9),status='old')
  open (unit=20,file=trim(filei10),status='old')
  open (unit=21,file=trim(filei11),status='old')
  open (unit=22,file=trim(filei12),status='old')
  open (unit=23,file=trim(filei13),status='old')
  open (unit=24,file=trim(filei14),status='old')
  open (unit=25,file=trim(filei15),status='old')
  open (unit=26,file=trim(filei16),status='old')
  open (unit=27,file=trim(filei17),status='old')
  open (unit=28,file=trim(filei18),status='old')
  open (unit=29,file=trim(filei19),status='old')
  open (unit=30,file=trim(filei20),status='old')

! initialize first

  do j = 1, nlay
     do i = 1, mapunitmax
        pct_clay(i,j) = 0.
        pct_sand(i,j) = 0.
     end do
  end do

! first clay

  do j = 1, nlay
     read(10+j,*) ! clear the first line
     do i = 1, nmapunits
        read(10+j,*) mu, pct_clay(mu,j)
     end do
     close(10+j)
  end do

! then sand

  do j = 1, nlay
     read(20+j,*) ! clear the first line
     do i = 1, nmapunits
        read(20+j,*) mu, pct_sand(mu,j)
     end do
     close(20+j)
  end do

! make north to south back to south to north

  do j = 1, nlat
     do i = 1, nlon
        temp(i,j) = mapunit(i,nlat-j+1)
     end do
  end do
  do j = 1, nlat
     do i = 1, nlon
        mapunit(i,j) = temp(i,j)
     end do
  end do
  do j = 1, nlat
     do i = 1, nlon
        temp(i,j) = landmask(i,nlat-j+1)
     end do
  end do
  do j = 1, nlat
     do i = 1, nlon
        landmask(i,j) = temp(i,j)
     end do
  end do
  do j = 1, nlat
     do i = 1, nlon
        temp(i,j) = latixy(i,nlat-j+1)
     end do
  end do
  do j = 1, nlat
     do i = 1, nlon
        latixy(i,j) = temp(i,j)
     end do
  end do
  do j = 1, nlat
     do i = 1, nlon
        temp(i,j) = longxy(i,nlat-j+1)
     end do
  end do
  do j = 1, nlat
     do i = 1, nlon
        longxy(i,j) = temp(i,j)
     end do
  end do

  lat(:) = latixy(1,:)
  lon(:) = longxy(:,1)

! Create output file

  call wrap_put_var_realx (ncid, lon_id     , lon)
  call wrap_put_var_realx (ncid, lat_id     , lat)
  call wrap_put_var_realx (ncid, longxy_id  , longxy)
  call wrap_put_var_realx (ncid, latixy_id  , latixy)
  call wrap_put_var_realx (ncid, landmask_id, landmask)
  call wrap_put_var_realx (ncid, edgen_id   , edge(1))
  call wrap_put_var_realx (ncid, edgee_id   , edge(2))
  call wrap_put_var_realx (ncid, edges_id   , edge(3))
  call wrap_put_var_realx (ncid, edgew_id   , edge(4))
  call wrap_put_var_realx (ncid, dzsoi_id   , dzsoi)
  call wrap_put_var_realx (ncid, zsoi_id    , zsoi)
  call wrap_put_var_realx (ncid, mapunit_id , mapunit)
  call wrap_put_var_realx (ncid, pct_sand_id, pct_sand)
  call wrap_put_var_realx (ncid, pct_clay_id, pct_clay)

  call wrap_close(ncid)

end program convert_soitex

!===============================================================================

subroutine wrap_create (path, cmode, ncid)
  implicit none
  include 'netcdf.inc'
  integer, parameter :: r8 = selected_real_kind(12)
  character(len=*) path
  integer cmode, ncid, ret
  ret = nf_create (path, cmode, ncid)
  if (ret.ne.NF_NOERR) call handle_error (ret)
end subroutine wrap_create

!===============================================================================

subroutine wrap_def_dim (nfid, dimname, len, dimid)
  implicit none
  include 'netcdf.inc'
  integer, parameter :: r8 = selected_real_kind(12)
  integer :: nfid, len, dimid
  character(len=*) :: dimname
  integer ret
  ret = nf_def_dim (nfid, dimname, len, dimid)
  if (ret.ne.NF_NOERR) call handle_error (ret)
end subroutine wrap_def_dim

!===============================================================================

subroutine wrap_def_var (nfid, name, xtype, nvdims, vdims, varid)
  implicit none
  include 'netcdf.inc'
  integer, parameter :: r8 = selected_real_kind(12)
  integer :: nfid, xtype, nvdims, varid
  integer :: vdims(nvdims)
  character(len=*) :: name
  integer ret
  ret = nf_def_var (nfid, name, xtype, nvdims, vdims, varid)
  if (ret.ne.NF_NOERR) call handle_error (ret)
end subroutine wrap_def_var

!===============================================================================

subroutine wrap_put_att_text (nfid, varid, attname, atttext)
  implicit none
  include 'netcdf.inc'
  integer, parameter :: r8 = selected_real_kind(12)
  integer :: nfid, varid
  character(len=*) :: attname, atttext
  integer :: ret, siz
  siz = len_trim(atttext)
  ret = nf_put_att_text (nfid, varid, attname, siz, atttext)
  if (ret.ne.NF_NOERR) call handle_error (ret)
end subroutine wrap_put_att_text

!===============================================================================

subroutine wrap_put_var_realx (nfid, varid, arr)
  implicit none
  include 'netcdf.inc'
  integer, parameter :: r8 = selected_real_kind(12)
  integer :: nfid, varid
  real(r8) :: arr(*)
  integer :: ret
#ifdef CRAY
  ret = nf_put_var_real (nfid, varid, arr)
#else
  ret = nf_put_var_double (nfid, varid, arr)
#endif
  if (ret.ne.NF_NOERR) call handle_error (ret)
end subroutine wrap_put_var_realx

!===============================================================================

subroutine wrap_put_var_int (nfid, varid, arr)
  implicit none
  include 'netcdf.inc'
  integer, parameter :: r8 = selected_real_kind(12)
  integer :: nfid, varid
  integer :: arr(*)
  integer :: ret
  ret = nf_put_var_int (nfid, varid, arr)
  if (ret.ne.NF_NOERR) call handle_error (ret)
end subroutine wrap_put_var_int
  
!===============================================================================

subroutine wrap_close (ncid)
  implicit none
  include 'netcdf.inc'
  integer, parameter :: r8 = selected_real_kind(12)
  integer :: ncid
  integer :: ret
  ret = nf_close (ncid)
  if (ret.ne.NF_NOERR) then
     write(6,*)'WRAP_CLOSE: nf_close failed for id ',ncid
     call handle_error (ret)
  end if
end subroutine wrap_close

!===============================================================================

subroutine handle_error(ret)
  implicit none
  include 'netcdf.inc'
  integer :: ret
  if (ret .ne. nf_noerr) then
     write(6,*) 'NCDERR: ERROR: ',nf_strerror(ret)
     call abort
  endif
end subroutine handle_error

!===============================================================================