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 ifend 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  endifend subroutine handle_error!===============================================================================