SUBROUTINE WRITE_RESTART()

  ! UW Land Surface Hydrology Group implementation of Catchment model
  ! modified from NLDAS implementation
  ! author: Kostas Andreadis, kostas@hydro.washington.edu
  ! updated by Ben Livneh to adhere to UW conventions netCDF

  ! WRITES RESTART FILE

  ! driverMod contains definitions of all global driver variables
  USE driverMod

  IMPLICIT NONE

  ! Define local variables
  INTEGER varid
  INTEGER J,I,K,L,start2d(2),count2d(2),start3d(3),count3d(3),start4d(4),count4d(4)
  REAL    TEMP(xlen,ylen)
  INTEGER TEMP_INT(xlen,ylen)
  INTEGER land_idx

  ! Write Ground Heat Content
  status = NF_INQ_VARID(RESTART_NCID,'ght',varid)
  DO L = 1, MAXNSOIL
      start3d(1) = 1
      start3d(2) = 1
      start3d(3) = L
      count3d(1) = xlen
      count3d(2) = ylen
      count3d(3) = 1
      land_idx = 1
      DO I = 1, ylen
        DO J = 1, xlen
          IF (LANDMASK(J,I) == 1) THEN
            TEMP(J,I) =  CAT_PROGN(land_idx)%ght(L)
            land_idx = land_idx + 1
          ELSE
            TEMP(J,I) = NODATA
          END IF
        END DO
      END DO
      status = NF_PUT_VARA_REAL(RESTART_NCID,varid,start3d,count3d,TEMP)
   END DO

  ! Write SWE
  status = NF_INQ_VARID(RESTART_NCID,'wesn',varid)
  DO L = 1, MAXNSNOW
      start3d(1) = 1
      start3d(2) = 1
      start3d(3) = L
      count3d(1) = xlen
      count3d(2) = ylen
      count3d(3) = 1
      land_idx = 1
      DO I = 1, ylen
        DO J = 1, xlen
          IF (LANDMASK(J,I) == 1) THEN
            TEMP(J,I) =  CAT_PROGN(land_idx)%wesn(L)
            land_idx = land_idx + 1
          ELSE
            TEMP(J,I) = NODATA
          END IF
        END DO
      END DO
      status = NF_PUT_VARA_REAL(RESTART_NCID,varid,start3d,count3d,TEMP)
   END DO

  ! Write Snow Heat Content
  status = NF_INQ_VARID(RESTART_NCID,'htsn',varid)
  DO L = 1, MAXNSNOW
      start3d(1) = 1
      start3d(2) = 1
      start3d(3) = L
      count3d(1) = xlen
      count3d(2) = ylen
      count3d(3) = 1
      land_idx = 1
      DO I = 1, ylen
        DO J = 1, xlen
          IF (LANDMASK(J,I) == 1) THEN
            TEMP(J,I) =  CAT_PROGN(land_idx)%htsn(L)
            land_idx = land_idx + 1
          ELSE
            TEMP(J,I) = NODATA
          END IF
        END DO
      END DO
      status = NF_PUT_VARA_REAL(RESTART_NCID,varid,start3d,count3d,TEMP)
   END DO


  ! Write Snow Depth
  status = NF_INQ_VARID(RESTART_NCID,'sndz',varid)
  DO L = 1, MAXNSNOW
      start3d(1) = 1
      start3d(2) = 1
      start3d(3) = L
      count3d(1) = xlen
      count3d(2) = ylen
      count3d(3) = 1
      land_idx = 1
      DO I = 1, ylen
        DO J = 1, xlen
          IF (LANDMASK(J,I) == 1) THEN
            TEMP(J,I) =  CAT_PROGN(land_idx)%sndz(L)
            land_idx = land_idx + 1
          ELSE
            TEMP(J,I) = NODATA
          END IF
        END DO
      END DO
      status = NF_PUT_VARA_REAL(RESTART_NCID,varid,start3d,count3d,TEMP)
   END DO


  ! Write surface/canopy temperature
  status = NF_INQ_VARID(RESTART_NCID,'tc1',varid)
  start2d(1) = 1
  start2d(2) = 1
  count2d(1) = xlen
  count2d(2) = ylen
  land_idx = 1
  DO I = 1, ylen
     DO J = 1, xlen
        IF (LANDMASK(J,I) == 1) THEN
           TEMP(J,I) =  CAT_PROGN(land_idx)%tc1
           land_idx = land_idx + 1
        ELSE
           TEMP(J,I) = NODATA
        END IF
     END DO
  END DO
  status = NF_PUT_VARA_REAL(RESTART_NCID,varid,start2d,count2d,TEMP)

  ! Write surface/canopy temperature
  status = NF_INQ_VARID(RESTART_NCID,'tc2',varid)
  start2d(1) = 1
  start2d(2) = 1
  count2d(1) = xlen
  count2d(2) = ylen
  land_idx = 1
  DO I = 1, ylen
     DO J = 1, xlen
        IF (LANDMASK(J,I) == 1) THEN
           TEMP(J,I) =  CAT_PROGN(land_idx)%tc2
           land_idx = land_idx + 1
        ELSE
           TEMP(J,I) = NODATA
        END IF
     END DO
  END DO
  status = NF_PUT_VARA_REAL(RESTART_NCID,varid,start2d,count2d,TEMP)

  ! Write surface/canopy temperature
  status = NF_INQ_VARID(RESTART_NCID,'tc4',varid)
  start2d(1) = 1
  start2d(2) = 1
  count2d(1) = xlen
  count2d(2) = ylen
  land_idx = 1
  DO I = 1, ylen
     DO J = 1, xlen
        IF (LANDMASK(J,I) == 1) THEN
           TEMP(J,I) =  CAT_PROGN(land_idx)%tc4
           land_idx = land_idx + 1
        ELSE
           TEMP(J,I) = NODATA
        END IF
     END DO
  END DO
  status = NF_PUT_VARA_REAL(RESTART_NCID,varid,start2d,count2d,TEMP)

  ! Write specific humidity in canopy air
  status = NF_INQ_VARID(RESTART_NCID,'qa1',varid)
  start2d(1) = 1
  start2d(2) = 1
  count2d(1) = xlen
  count2d(2) = ylen
  land_idx = 1
  DO I = 1, ylen
     DO J = 1, xlen
        IF (LANDMASK(J,I) == 1) THEN
           TEMP(J,I) =  CAT_PROGN(land_idx)%qa1
           land_idx = land_idx + 1
        ELSE
           TEMP(J,I) = NODATA
        END IF
     END DO
  END DO
  status = NF_PUT_VARA_REAL(RESTART_NCID,varid,start2d,count2d,TEMP)

  ! Write specific humidity in canopy air
  status = NF_INQ_VARID(RESTART_NCID,'qa2',varid)
  start2d(1) = 1
  start2d(2) = 1
  count2d(1) = xlen
  count2d(2) = ylen
  land_idx = 1
  DO I = 1, ylen
     DO J = 1, xlen
        IF (LANDMASK(J,I) == 1) THEN
           TEMP(J,I) =  CAT_PROGN(land_idx)%qa2
           land_idx = land_idx + 1
        ELSE
           TEMP(J,I) = NODATA
        END IF
     END DO
  END DO
  status = NF_PUT_VARA_REAL(RESTART_NCID,varid,start2d,count2d,TEMP)

  ! Write specific humidity in canopy air
  status = NF_INQ_VARID(RESTART_NCID,'qa4',varid)
  start2d(1) = 1
  start2d(2) = 1
  count2d(1) = xlen
  count2d(2) = ylen
  land_idx = 1
  DO I = 1, ylen
     DO J = 1, xlen
        IF (LANDMASK(J,I) == 1) THEN
           TEMP(J,I) =  CAT_PROGN(land_idx)%qa4
           land_idx = land_idx + 1
        ELSE
           TEMP(J,I) = NODATA
        END IF
     END DO
  END DO
  status = NF_PUT_VARA_REAL(RESTART_NCID,varid,start2d,count2d,TEMP)

  ! Write canopy interception water
  status = NF_INQ_VARID(RESTART_NCID,'capac',varid)
  start2d(1) = 1
  start2d(2) = 1
  count2d(1) = xlen
  count2d(2) = ylen
  land_idx = 1
  DO I = 1, ylen
     DO J = 1, xlen
        IF (LANDMASK(J,I) == 1) THEN
           TEMP(J,I) =  CAT_PROGN(land_idx)%capac
           land_idx = land_idx + 1
        ELSE
           TEMP(J,I) = NODATA
        END IF
     END DO
  END DO
  status = NF_PUT_VARA_REAL(RESTART_NCID,varid,start2d,count2d,TEMP)

  ! Write catchment defecit
  status = NF_INQ_VARID(RESTART_NCID,'catdef',varid)
  start2d(1) = 1
  start2d(2) = 1
  count2d(1) = xlen
  count2d(2) = ylen
  land_idx = 1
  DO I = 1, ylen
     DO J = 1, xlen
        IF (LANDMASK(J,I) == 1) THEN
           TEMP(J,I) =  CAT_PROGN(land_idx)%catdef
           land_idx = land_idx + 1
        ELSE
           TEMP(J,I) = NODATA
        END IF
     END DO
  END DO
  status = NF_PUT_VARA_REAL(RESTART_NCID,varid,start2d,count2d,TEMP)

  ! Write root zone excess
  status = NF_INQ_VARID(RESTART_NCID,'rzexc',varid)
  start2d(1) = 1
  start2d(2) = 1
  count2d(1) = xlen
  count2d(2) = ylen
  land_idx = 1
  DO I = 1, ylen
     DO J = 1, xlen
        IF (LANDMASK(J,I) == 1) THEN
           TEMP(J,I) =  CAT_PROGN(land_idx)%rzexc
           land_idx = land_idx + 1
        ELSE
           TEMP(J,I) = NODATA
        END IF
     END DO
  END DO
  status = NF_PUT_VARA_REAL(RESTART_NCID,varid,start2d,count2d,TEMP)

  ! Write surface excess
  status = NF_INQ_VARID(RESTART_NCID,'srfexc',varid)
  start2d(1) = 1
  start2d(2) = 1
  count2d(1) = xlen
  count2d(2) = ylen
  land_idx = 1
  DO I = 1, ylen
     DO J = 1, xlen
        IF (LANDMASK(J,I) == 1) THEN
           TEMP(J,I) =  CAT_PROGN(land_idx)%srfexc
           land_idx = land_idx + 1
        ELSE
           TEMP(J,I) = NODATA
        END IF
     END DO
  END DO
  status = NF_PUT_VARA_REAL(RESTART_NCID,varid,start2d,count2d,TEMP)

  ! Close file
  status = NF_CLOSE(RESTART_NCID)

!  INTEGER :: n , k


!  OPEN(10, file=RESTFILE, form='unformatted', status='unknown', &
!       action='write')
!
!  WRITE (10) (CAT_PROGN(n)%tc1,    n=1,landlen)
!  WRITE (10) (CAT_PROGN(n)%tc2,    n=1,landlen)
!  WRITE (10) (CAT_PROGN(n)%tc4,    n=1,landlen)
!
!  WRITE (10) (CAT_PROGN(n)%qa1,    n=1,landlen)
!  WRITE (10) (CAT_PROGN(n)%qa2,    n=1,landlen)
!  WRITE (10) (CAT_PROGN(n)%qa4,    n=1,landlen)
!
!  WRITE (10) (CAT_PROGN(n)%capac,  n=1,landlen)
!
!  WRITE (10) (CAT_PROGN(n)%catdef, n=1,landlen)
!  WRITE (10) (CAT_PROGN(n)%rzexc,  n=1,landlen)
!  WRITE (10) (CAT_PROGN(n)%srfexc, n=1,landlen)
!
!  do k=1,N_gndtmp
!     WRITE (10) (CAT_PROGN(n)%ght(k),  n=1,landlen)
!  end do
!
!  do k=1,N_snow
!     WRITE (10) (CAT_PROGN(n)%wesn(k), n=1,landlen)
!  end do
!  do k=1,N_snow
!     WRITE (10) (CAT_PROGN(n)%htsn(k), n=1,landlen)
!  end do
!  do k=1,N_snow
!     WRITE (10) (CAT_PROGN(n)%sndz(k), n=1,landlen)
!  end do
!
!  ! Close file
!  CLOSE(10,status='keep')

END SUBROUTINE WRITE_RESTART