MODULE restart_mod
   USE prec, ONLY: rstd
   USE field_mod
   USE domain_mod
   USE xios_mod
   USE netcdf_mod
   USE omp_para, ONLY: is_omp_master
   USE mpipara, ONLY: is_mpi_master, is_mpi_root
   USE time_mod, ONLY: itaumax
   USE metric, ONLY: ncell_glo, cell_glo
   USE getin_mod, ONLY: getin
   USE spherical_geom_mod, ONLY: xyz2lonlat
   USE transfert_mod, ONLY: gather_field, scatter_field
   IMPLICIT NONE
   PRIVATE

   TYPE t_field_array
      TYPE(t_field), POINTER :: field(:)
   END TYPE t_field_array

   LOGICAL, SAVE :: is_write_start = .TRUE.
   LOGICAL, SAVE :: is_read_start = .FALSE.
   LOGICAL, SAVE :: is_read_interp = .FALSE.

   PUBLIC :: init_restart, restart_var_present, read_xios_start_field, &
             is_write_start, set_write_start, write_restart, write_restart_fields, &
             is_read_start, set_read_start, read_start, read_start_fields, &
             is_read_interp, set_read_interp

   INTERFACE write_restart
      MODULE PROCEDURE write_restart_fields, write_restart_r0, write_restart_i0
   END INTERFACE

   INTERFACE read_start
      MODULE PROCEDURE read_start_fields, read_start_r0, read_start_i0
   END INTERFACE

CONTAINS

   SUBROUTINE set_write_start(flag)
      LOGICAL :: flag
      is_write_start = flag
   END SUBROUTINE set_write_start

   SUBROUTINE set_read_start(flag)
      LOGICAL :: flag
      is_read_start = flag
   END SUBROUTINE set_read_start

   SUBROUTINE set_read_interp(flag)
      LOGICAL :: flag
      is_read_interp = flag
   END SUBROUTINE set_read_interp

   SUBROUTINE init_restart
      USE grid_param, ONLY: nqtot
      CHARACTER(LEN=255) :: start_file_name
      CHARACTER(LEN=255) :: restart_file_name

      IF (using_xios) THEN
         start_file_name = "start"
         CALL getin("start_file_name", start_file_name)
         restart_file_name = "restart"
         CALL getin("restart_file_name", restart_file_name)
         IF (is_omp_master) THEN
!        CALL xios_set_file_attr("start",name=TRIM(ADJUSTL(start_file_name)),output_freq=(itaumax+1)*xios_timestep)
            CALL xios_set_file_attr("start", name=TRIM(ADJUSTL(start_file_name)), output_freq=1*xios_timestep)
            CALL xios_set_file_attr("restart", name=TRIM(ADJUSTL(restart_file_name)), output_freq=itaumax*xios_timestep)
            CALL xios_set_fieldgroup_attr("group_restart", freq_op=itaumax*xios_timestep)
            IF (nqtot == 0) CALL xios_set_field_attr("q_restart", enabled=.FALSE.)
         END IF
      END IF

   END SUBROUTINE init_restart

   FUNCTION restart_var_present(name)
      USE netcdf_mod
      USE omp_para
      USE mpipara
      USE icosa
      IMPLICIT NONE
      LOGICAL :: restart_var_present
      CHARACTER(LEN=*), INTENT(IN) :: name

      CHARACTER(LEN=255) :: start_file_name
      INTEGER :: status, varid, ncid

      restart_var_present = .FALSE.

      IF (no_io) RETURN
      start_file_name = "start"
      CALL getin("start_file_name", start_file_name)

      IF (is_omp_master) THEN
         IF (is_read_start) THEN
            IF (is_mpi_master) THEN
               status = NF90_OPEN(TRIM(start_file_name)//'.nc', NF90_NOWRITE, ncid)
               IF (status == NF90_NOERR) THEN
                  status = NF90_INQ_VARID(ncid, name, varid)
                  IF (status == NF90_NOERR) restart_var_present = .TRUE.
                  status = NF90_CLOSE(ncid)
               END IF
            END IF
         END IF
      END IF
      CALL bcast(restart_var_present)

   END FUNCTION restart_var_present

   SUBROUTINE write_restart_fields(it, field0, field1, field2, field3, field4, field5, field6, field7, field8, field9, &
                            field10, field11, field12, field13, field14, field15, field16, field17, field18, field19)
      USE math_const, ONLY: pi
      USE grid_param, ONLY: llm, nqtot
      USE prec, ONLY: r8
      USE disvert_mod, ONLY: presnivs

      INTEGER, INTENT(IN)     :: it

    TYPE(t_field), POINTER, OPTIONAL, DIMENSION(:) :: field0, field1, field2, field3, field4, field5, field6, field7, field8, field9
      TYPE(t_field),POINTER,OPTIONAL,DIMENSION(:) :: field10,field11,field12,field13,field14,field15,field16,field17,field18,field19

      TYPE(t_field_array) :: field_array(20)
      INTEGER             :: nfield
      INTEGER             :: fieldId(20)

      TYPE(t_domain), POINTER :: d
      TYPE(t_field), POINTER :: field_glo(:)
      TYPE(t_field), POINTER :: field(:)

      CHARACTER(LEN=255) :: file_name
      CHARACTER(LEN=255) :: suffix
      INTEGER, PARAMETER  :: nvert = 6
      INTEGER    ::  ncid, cellId, levId, edgeId, vertid, lonId, latId, bounds_lonId, bounds_latId, nqId, levAxisId
      INTEGER    :: ind, ind_glo, i, j, k, nf
      INTEGER    :: status
      REAL(rstd), ALLOCATABLE :: lon(:), lat(:), bounds_lon(:, :), bounds_lat(:, :)

      IF (no_io) RETURN

      IF (it == 0) THEN
         file_name = "start"
         CALL getin("start_file_name", file_name)
         suffix = "_start"
      ELSE
         file_name = "restart"
         CALL getin("restart_file_name", file_name)
         suffix = "_restart"
      END IF

      IF (using_xios) THEN
         IF (PRESENT(field0)) THEN; CALL xios_write_field(TRIM(field0(1)%name)//TRIM(suffix), field0); END IF
         IF (PRESENT(field1)) THEN; CALL xios_write_field(TRIM(field1(1)%name)//TRIM(suffix), field1); END IF
         IF (PRESENT(field2)) THEN; CALL xios_write_field(TRIM(field2(1)%name)//TRIM(suffix), field2); END IF
         IF (PRESENT(field3)) THEN; CALL xios_write_field(TRIM(field3(1)%name)//TRIM(suffix), field3); END IF
         IF (PRESENT(field4)) THEN; CALL xios_write_field(TRIM(field4(1)%name)//TRIM(suffix), field4); END IF
         IF (PRESENT(field5)) THEN; CALL xios_write_field(TRIM(field5(1)%name)//TRIM(suffix), field5); END IF
         IF (PRESENT(field6)) THEN; CALL xios_write_field(TRIM(field6(1)%name)//TRIM(suffix), field6); END IF
         IF (PRESENT(field7)) THEN; CALL xios_write_field(TRIM(field7(1)%name)//TRIM(suffix), field7); END IF
         IF (PRESENT(field8)) THEN; CALL xios_write_field(TRIM(field8(1)%name)//TRIM(suffix), field8); END IF
         IF (PRESENT(field9)) THEN; CALL xios_write_field(TRIM(field9(1)%name)//TRIM(suffix), field9); END IF
         IF (PRESENT(field10)) THEN; CALL xios_write_field(TRIM(field10(1)%name)//TRIM(suffix), field10); END IF
         IF (PRESENT(field11)) THEN; CALL xios_write_field(TRIM(field11(1)%name)//TRIM(suffix), field11); END IF
         IF (PRESENT(field12)) THEN; CALL xios_write_field(TRIM(field12(1)%name)//TRIM(suffix), field12); END IF
         IF (PRESENT(field13)) THEN; CALL xios_write_field(TRIM(field13(1)%name)//TRIM(suffix), field13); END IF
         IF (PRESENT(field14)) THEN; CALL xios_write_field(TRIM(field14(1)%name)//TRIM(suffix), field14); END IF
         IF (PRESENT(field15)) THEN; CALL xios_write_field(TRIM(field15(1)%name)//TRIM(suffix), field15); END IF
         IF (PRESENT(field16)) THEN; CALL xios_write_field(TRIM(field16(1)%name)//TRIM(suffix), field16); END IF
         IF (PRESENT(field17)) THEN; CALL xios_write_field(TRIM(field17(1)%name)//TRIM(suffix), field17); END IF
         IF (PRESENT(field18)) THEN; CALL xios_write_field(TRIM(field18(1)%name)//TRIM(suffix), field18); END IF
         IF (PRESENT(field19)) THEN; CALL xios_write_field(TRIM(field19(1)%name)//TRIM(suffix), field19); END IF
         IF (is_omp_master) CALL xios_send_field('it'//TRIM(suffix), it*1.0)

      ELSE

         !$OMP MASTER

         nfield = 0
         IF (PRESENT(field0)) THEN; nfield = nfield + 1; field_array(nfield)%field => field0; END IF
         IF (PRESENT(field1)) THEN; nfield = nfield + 1; field_array(nfield)%field => field1; END IF
         IF (PRESENT(field2)) THEN; nfield = nfield + 1; field_array(nfield)%field => field2; END IF
         IF (PRESENT(field3)) THEN; nfield = nfield + 1; field_array(nfield)%field => field3; END IF
         IF (PRESENT(field4)) THEN; nfield = nfield + 1; field_array(nfield)%field => field4; END IF
         IF (PRESENT(field5)) THEN; nfield = nfield + 1; field_array(nfield)%field => field5; END IF
         IF (PRESENT(field6)) THEN; nfield = nfield + 1; field_array(nfield)%field => field6; END IF
         IF (PRESENT(field7)) THEN; nfield = nfield + 1; field_array(nfield)%field => field7; END IF
         IF (PRESENT(field8)) THEN; nfield = nfield + 1; field_array(nfield)%field => field8; END IF
         IF (PRESENT(field9)) THEN; nfield = nfield + 1; field_array(nfield)%field => field9; END IF
         IF (PRESENT(field10)) THEN; nfield = nfield + 1; field_array(nfield)%field => field10; END IF
         IF (PRESENT(field11)) THEN; nfield = nfield + 1; field_array(nfield)%field => field11; END IF
         IF (PRESENT(field12)) THEN; nfield = nfield + 1; field_array(nfield)%field => field12; END IF
         IF (PRESENT(field13)) THEN; nfield = nfield + 1; field_array(nfield)%field => field13; END IF
         IF (PRESENT(field14)) THEN; nfield = nfield + 1; field_array(nfield)%field => field14; END IF
         IF (PRESENT(field15)) THEN; nfield = nfield + 1; field_array(nfield)%field => field15; END IF
         IF (PRESENT(field16)) THEN; nfield = nfield + 1; field_array(nfield)%field => field16; END IF
         IF (PRESENT(field17)) THEN; nfield = nfield + 1; field_array(nfield)%field => field17; END IF
         IF (PRESENT(field18)) THEN; nfield = nfield + 1; field_array(nfield)%field => field18; END IF
         IF (PRESENT(field19)) THEN; nfield = nfield + 1; field_array(nfield)%field => field19; END IF

         IF (is_mpi_root) THEN
            status = NF90_CREATE(TRIM(ADJUSTL(file_name))//'.nc', NF90_CLOBBER, ncid)
            status = NF90_DEF_DIM(ncid, 'cell', ncell_glo, cellId)
            status = NF90_DEF_DIM(ncid, 'edge', 3*ncell_glo, edgeId)
            status = NF90_DEF_DIM(ncid, 'lev', llm, levId)
            status = NF90_DEF_DIM(ncid, 'nvert', nvert, vertId)
            status = NF90_DEF_DIM(ncid, 'nq', nqtot, nqId)
            status = NF90_PUT_ATT(ncid, NF90_GLOBAL, "iteration", it)

            status = NF90_DEF_VAR(ncid, 'lon', NF90_DOUBLE, (/cellId/), lonId)
            status = NF90_PUT_ATT(ncid, lonId, "long_name", "longitude")
            status = NF90_PUT_ATT(ncid, lonId, "units", "degrees_east")
            status = NF90_PUT_ATT(ncid, lonId, "bounds", "bounds_lon")
            status = NF90_DEF_VAR(ncid, 'lat', NF90_DOUBLE, (/cellId/), latId)
            status = NF90_PUT_ATT(ncid, latId, "long_name", "latitude")
            status = NF90_PUT_ATT(ncid, latId, "units", "degrees_north")
            status = NF90_PUT_ATT(ncid, latId, "bounds", "bounds_lat")
            status = NF90_DEF_VAR(ncid, 'bounds_lon', NF90_DOUBLE, (/vertId, cellId/), bounds_lonId)
            status = NF90_DEF_VAR(ncid, 'bounds_lat', NF90_DOUBLE, (/vertId, cellId/), bounds_latId)
            status = NF90_DEF_VAR(ncid, 'lev', NF90_DOUBLE, (/levId/), levAxisId)
            status = NF90_PUT_ATT(ncid, levAxisId, "axis", "Z")
            status = NF90_PUT_ATT(ncid, levAxisId, "units", "Pa")
            status = NF90_PUT_ATT(ncid, levAxisId, "positive", "down")

            DO nf = 1, nfield
               field => field_array(nf)%field
               IF (field(1)%field_type == field_T) THEN
                  IF (field(1)%ndim == 2) THEN
                     status = NF90_DEF_VAR(ncid, TRIM(ADJUSTL(field(1)%name)), NF90_DOUBLE, (/cellId/), fieldId(nf))
                     status = NF90_PUT_ATT(ncid, FieldId(nf), "coordinates", "lon lat")
                  ELSE IF (field(1)%ndim == 3) THEN
                     status = NF90_DEF_VAR(ncid, TRIM(ADJUSTL(field(1)%name)), NF90_DOUBLE, (/cellId, levId/), fieldId(nf))
                     status = NF90_PUT_ATT(ncid, FieldId(nf), "coordinates", "lev lon lat")
                  ELSE IF (field(1)%ndim == 4) THEN
                     status = NF90_DEF_VAR(ncid, TRIM(ADJUSTL(field(1)%name)), NF90_DOUBLE, (/cellId, levId, nqId/), fieldId(nf))
                     status = NF90_PUT_ATT(ncid, FieldId(nf), "coordinates", "nq lev lon lat")
                  END IF
               ELSE IF (field(1)%field_type == field_U) THEN
                  IF (field(1)%ndim == 2) THEN
                     status = NF90_DEF_VAR(ncid, TRIM(ADJUSTL(field(1)%name)), NF90_DOUBLE, (/edgeId/), fieldId(nf))
                  ELSE IF (field(1)%ndim == 3) THEN
                     status = NF90_DEF_VAR(ncid, TRIM(ADJUSTL(field(1)%name)), NF90_DOUBLE, (/edgeId, levId/), fieldId(nf))
                  ELSE IF (field(1)%ndim == 4) THEN
                     status = NF90_DEF_VAR(ncid, TRIM(ADJUSTL(field(1)%name)), NF90_DOUBLE, (/edgeId, levId, nqId/), fieldId(nf))
                  END IF
               END IF
            END DO

            status = NF90_ENDDEF(ncid)

            ALLOCATE (lon(ncell_glo), lat(ncell_glo), bounds_lon(0:nvert - 1, ncell_glo), bounds_lat(0:nvert - 1, ncell_glo))
            DO ind = 1, ndomain_glo
               d => domain_glo(ind)
               DO j = d%jj_begin, d%jj_end
                  DO i = d%ii_begin, d%ii_end
                     ind_glo = d%assign_cell_glo(i, j)
                     CALL xyz2lonlat(d%xyz(:, i, j), lon(ind_glo), lat(ind_glo))
                     lon(ind_glo) = lon(ind_glo)*180/Pi
                     lat(ind_glo) = lat(ind_glo)*180/Pi
                     DO k = 0, 5
                        CALL xyz2lonlat(d%vertex(:, k, i, j), bounds_lon(k, ind_glo), bounds_lat(k, ind_glo))
                        bounds_lat(k, ind_glo) = bounds_lat(k, ind_glo)*180/Pi
                        bounds_lon(k, ind_glo) = bounds_lon(k, ind_glo)*180/Pi
                     END DO
                  END DO
               END DO
            END DO

            status = NF90_PUT_VAR(ncid, lonid, REAL(lon, r8), start=(/1/), count=(/ncell_glo/))
            status = NF90_PUT_VAR(ncid, latid, REAL(lat, r8), start=(/1/), count=(/ncell_glo/))
            status = NF90_PUT_VAR(ncid, bounds_lonId, REAL(bounds_lon, r8), start=(/1, 1/), count=(/nvert, ncell_glo/))
            status = NF90_PUT_VAR(ncid, bounds_latId, REAL(bounds_lat, r8), start=(/1, 1/), count=(/nvert, ncell_glo/))
            status = NF90_PUT_VAR(ncid, levAxisId, REAL(presnivs, r8), start=(/1/), count=(/llm/))
         END IF

         DO nf = 1, nfield
            field => field_array(nf)%field
            CALL write_restart_field(field, fieldId(nf), ncid)
         END DO

         !          status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(field(1)%name),NF90_DOUBLE,(/ ncellId /),fieldId(nf))
         !        ELSE IF (field(1)%ndim==3) THEN
         !          status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(field(1)%name),NF90_DOUBLE,(/ ncellId, levId /),fieldId(nf))
         !        ENDIF
         !      ENDDO

         IF (is_mpi_root) THEN
            status = NF90_CLOSE(ncid)
         END IF

         !$OMP END MASTER

      END IF

   END SUBROUTINE write_restart_fields

   SUBROUTINE write_restart_i0(it, name, var)
      USE xios_mod
      USE getin_mod
      USE time_mod
      USE netcdf_mod
      USE omp_para
      IMPLICIT NONE
      INTEGER, INTENT(IN)               :: it
      CHARACTER(LEN=*), INTENT(IN)     :: name
      INTEGER, INTENT(IN)                :: var
      CHARACTER(LEN=255) :: file_name
      CHARACTER(LEN=255) :: suffix

      IF (no_io) RETURN

      IF (it == itau0 .AND. is_write_start) THEN
         file_name = "start"
         CALL getin("start_file_name", file_name)
         suffix = "_start0"
      ELSE
         file_name = "restart"
         CALL getin("restart_file_name", file_name)
         suffix = "_restart"
      END IF

      IF (is_omp_master) CALL xios_send_field(TRIM(name)//TRIM(suffix), var*1.0)

   END SUBROUTINE write_restart_i0

   SUBROUTINE write_restart_r0(it, name, var)
      USE xios_mod
      USE getin_mod
      USE time_mod
      USE netcdf_mod
      USE omp_para
      IMPLICIT NONE
      INTEGER, INTENT(IN)               :: it
      CHARACTER(LEN=*), INTENT(IN)     :: name
      REAL, INTENT(IN)                   :: var

      CHARACTER(LEN=255) :: file_name
      CHARACTER(LEN=255) :: suffix

      IF (no_io) RETURN

      IF (it == itau0 .AND. is_write_start) THEN
         file_name = "start"
         CALL getin("start_file_name", file_name)
         suffix = "_start0"
      ELSE
         file_name = "restart"
         CALL getin("restart_file_name", file_name)
         suffix = "_restart"
      END IF

      IF (is_omp_master) CALL xios_send_field(TRIM(name)//TRIM(suffix), var)

   END SUBROUTINE write_restart_r0

   SUBROUTINE write_xios_restart_field(fieldId, suffix, f_field)
      USE field_mod
      USE grid_param
      USE domain_mod
      USE xios_mod
      USE omp_para
      USE tracer_icosa_mod
      IMPLICIT NONE
      CHARACTER(LEN=*), INTENT(IN)        :: fieldId
      CHARACTER(LEN=*), INTENT(IN)        :: suffix
      TYPE(t_field), POINTER, SAVE         :: f_tracer(:)
      TYPE(t_field), POINTER, DIMENSION(:) :: f_field
      INTEGER                            :: iq, ind, l

      IF (TRIM(fieldId) == "q") THEN! hook for tracers
         CALL allocate_field(f_tracer, field_t, type_real, llm)
         DO iq = 1, nqtot
            DO ind = 1, ndomain
               IF (.NOT. assigned_domain(ind)) CYCLE
               DO l = ll_begin, ll_end
                  f_tracer(ind)%rval3d(:, l) = f_field(ind)%rval4d(:, l, iq)
               END DO
            END DO
            CALL xios_write_field(TRIM(tracers(iq)%name)//TRIM(suffix), f_tracer)
         END DO
         CALL deallocate_field(f_tracer)
      ELSE
         CALL xios_write_field(TRIM(fieldId)//TRIM(suffix), f_field)
      END IF

   END SUBROUTINE write_xios_restart_field

   SUBROUTINE write_restart_field(field, fieldId, ncid)
      USE prec, ONLY: r8
      USE grid_param, ONLY: llm, nqtot
      TYPE(t_field), POINTER :: field(:)
      INTEGER, INTENT(IN)     :: fieldId
      INTEGER, INTENT(IN)     :: ncid

      TYPE(t_domain), POINTER :: d
      TYPE(t_field), POINTER :: field_glo(:)
      REAL(rstd), ALLOCATABLE :: global_field2d(:)
      REAL(rstd), ALLOCATABLE :: global_field3d(:, :)
      REAL(rstd), ALLOCATABLE :: global_field4d(:, :, :)
      INTEGER :: i, j, ij, k, e, ind, ind_glo
      INTEGER :: ndim, field_type
      INTEGER :: status

      ndim = field(1)%ndim
      field_Type = field(1)%field_type

      IF (is_mpi_root) THEN

         IF (ndim == 2) THEN
            CALL allocate_field_glo(field_glo, field(1)%field_type, field(1)%data_type, name=field(1)%name)
         ELSE IF (ndim == 3) THEN
            CALL allocate_field_glo(field_glo, field(1)%field_type, field(1)%data_type, field(1)%dim3, name=field(1)%name)
         ELSE IF (ndim == 4) THEN
       CALL allocate_field_glo(field_glo, field(1)%field_type, field(1)%data_type, field(1)%dim3, field(1)%dim4, name=field(1)%name)
         END IF
      END IF

      CALL gather_field(field, field_glo)

      IF (is_mpi_root) THEN

         IF (field_type == field_T) THEN
            IF (ndim == 2) THEN
               ALLOCATE (global_field2d(ncell_glo))
               DO ind = 1, ndomain_glo
                  d => domain_glo(ind)
                  DO j = d%jj_begin, d%jj_end
                     DO i = d%ii_begin, d%ii_end
                        IF (d%own(i, j)) THEN
                           ij = (j - 1)*d%iim + i
                           ind_glo = d%assign_cell_glo(i, j)
                           global_field2d(ind_glo) = field_glo(ind)%rval2d(ij)
                        END IF
                     END DO
                  END DO
               END DO
               status = NF90_PUT_VAR(ncid, fieldid, global_field2d, start=(/1/), count=(/ncell_glo/))

            ELSE IF (ndim == 3) THEN
               ALLOCATE (global_field3d(ncell_glo, llm))
               DO ind = 1, ndomain_glo
                  d => domain_glo(ind)
                  DO j = d%jj_begin, d%jj_end
                     DO i = d%ii_begin, d%ii_end
                        IF (d%own(i, j)) THEN
                           ij = (j - 1)*d%iim + i
                           ind_glo = d%assign_cell_glo(i, j)
                           global_field3d(ind_glo, :) = field_glo(ind)%rval3d(ij, :)
                        END IF
                     END DO
                  END DO
               END DO
               status = NF90_PUT_VAR(ncid, fieldid, global_field3d, start=(/1, 1/), count=(/ncell_glo, llm/))
            ELSE IF (ndim == 4) THEN
               ALLOCATE (global_field4d(ncell_glo, llm, nqtot))
               DO ind = 1, ndomain_glo
                  d => domain_glo(ind)
                  DO j = d%jj_begin, d%jj_end
                     DO i = d%ii_begin, d%ii_end
                        IF (d%own(i, j)) THEN
                           ij = (j - 1)*d%iim + i
                           ind_glo = d%assign_cell_glo(i, j)
                           global_field4d(ind_glo, :, :) = field_glo(ind)%rval4d(ij, :, :)
                        END IF
                     END DO
                  END DO
               END DO
               status = NF90_PUT_VAR(ncid, fieldid, global_field4d, start=(/1, 1, 1/), count=(/ncell_glo, llm, nqtot/))
            END IF

         ELSE IF (field_type == field_U) THEN

            IF (ndim == 2) THEN
               ALLOCATE (global_field2d(3*ncell_glo))
               DO ind = 1, ndomain_glo
                  d => domain_glo(ind)
                  DO j = d%jj_begin, d%jj_end
                     DO i = d%ii_begin, d%ii_end
                        DO k = 0, 5
                  IF (d%edge_assign_domain(k, i, j) == ind .AND. d%edge_assign_i(k, i, j) == i .AND. d%edge_assign_j(k, i, j) == j &
                            .AND. d%edge_assign_pos(k, i, j) == k) THEN
                           ij = (j - 1)*d%iim + i
                           ind_glo = d%assign_cell_glo(i, j)
                           e = cell_glo(ind_glo)%edge(MOD(k + d%delta(i, j) + 6, 6))
                           global_field2d(e) = d%edge_assign_sign(k, i, j)*field_glo(ind)%rval2d(ij + d%u_pos(k))
                        END IF
                        END DO
                     END DO
                  END DO
               END DO
               status = NF90_PUT_VAR(ncid, fieldid, REAL(global_field2d, r8), start=(/1/), count=(/3*ncell_glo/))
            ELSE IF (ndim == 3) THEN
               ALLOCATE (global_field3d(3*ncell_glo, llm))
               DO ind = 1, ndomain_glo
                  d => domain_glo(ind)
                  DO j = d%jj_begin, d%jj_end
                     DO i = d%ii_begin, d%ii_end
                        DO k = 0, 5
                  IF (d%edge_assign_domain(k, i, j) == ind .AND. d%edge_assign_i(k, i, j) == i .AND. d%edge_assign_j(k, i, j) == j &
                            .AND. d%edge_assign_pos(k, i, j) == k) THEN
                           ij = (j - 1)*d%iim + i
                           ind_glo = d%assign_cell_glo(i, j)
                           e = cell_glo(ind_glo)%edge(MOD(k + d%delta(i, j) + 6, 6))
                           global_field3d(e, :) = d%edge_assign_sign(k, i, j)*field_glo(ind)%rval3d(ij + d%u_pos(k + 1), :)
                        END IF
                        END DO
                     END DO
                  END DO
               END DO
               status = NF90_PUT_VAR(ncid, fieldid, global_field3d, start=(/1, 1/), count=(/3*ncell_glo, llm/))
            ELSE IF (ndim == 4) THEN
               ALLOCATE (global_field4d(3*ncell_glo, llm, nqtot))
               DO ind = 1, ndomain_glo
                  d => domain_glo(ind)
                  DO j = d%jj_begin, d%jj_end
                     DO i = d%ii_begin, d%ii_end
                        DO k = 0, 5
                  IF (d%edge_assign_domain(k, i, j) == ind .AND. d%edge_assign_i(k, i, j) == i .AND. d%edge_assign_j(k, i, j) == j &
                            .AND. d%edge_assign_pos(k, i, j) == k) THEN
                           ij = (j - 1)*d%iim + i
                           ind_glo = d%assign_cell_glo(i, j)
                           e = cell_glo(ind_glo)%edge(MOD(k + d%delta(i, j) + 6, 6))
                           global_field4d(e, :, :) = d%edge_assign_sign(k, i, j)*field_glo(ind)%rval4d(ij + d%u_pos(k + 1), :, :)
                        END IF
                        END DO
                     END DO
                  END DO
               END DO
               status = NF90_PUT_VAR(ncid, fieldid, global_field3d, start=(/1, 1, 1/), count=(/3*ncell_glo, llm, nqtot/))
            END IF

         END IF

         CALL deallocate_field_glo(field_glo)

      END IF

   END SUBROUTINE write_restart_field

   SUBROUTINE read_start_fields(it, field0, field1, field2, field3, field4, field5, field6, field7, field8, field9, &
                         field10, field11, field12, field13, field14, field15, field16, field17, field18, field19)

      INTEGER, INTENT(OUT)  :: it
    TYPE(t_field), POINTER, OPTIONAL, DIMENSION(:) :: field0, field1, field2, field3, field4, field5, field6, field7, field8, field9
      TYPE(t_field),POINTER,OPTIONAL,DIMENSION(:) :: field10,field11,field12,field13,field14,field15,field16,field17,field18,field19

      TYPE(t_field_array) :: field_array(20)
      INTEGER             :: nfield
      INTEGER             :: fieldId(20)

      TYPE(t_domain), POINTER :: d
      TYPE(t_field), POINTER :: field_glo(:)
      TYPE(t_field), POINTER :: field(:)

      CHARACTER(LEN=255) :: start_file_name
      INTEGER, PARAMETER  :: nvert = 6
      INTEGER    ::  ncid, cellId, levId, edgeId, vertid, lonId, latId, bounds_lonId, bounds_latId
      INTEGER    :: ind, ind_glo, i, j, k, nf
      INTEGER    :: status
      REAL(rstd), ALLOCATABLE :: lon(:), lat(:), bounds_lon(:, :), bounds_lat(:, :)
      REAL(rstd) :: it_real

      IF (no_io) RETURN

      start_file_name = "start"
      CALL getin("start_file_name", start_file_name)

      IF (using_xios) THEN
         IF (PRESENT(field0)) THEN; CALL xios_read_field(TRIM(field0(1)%name)//'_start', field0); END IF
         IF (PRESENT(field1)) THEN; CALL xios_read_field(TRIM(field1(1)%name)//'_start', field1); END IF
         IF (PRESENT(field2)) THEN; CALL xios_read_field(TRIM(field2(1)%name)//'_start', field2); END IF
         IF (PRESENT(field3)) THEN; CALL xios_read_field(TRIM(field3(1)%name)//'_start', field3); END IF
         IF (PRESENT(field4)) THEN; CALL xios_read_field(TRIM(field4(1)%name)//'_start', field4); END IF
         IF (PRESENT(field5)) THEN; CALL xios_read_field(TRIM(field5(1)%name)//'_start', field5); END IF
         IF (PRESENT(field6)) THEN; CALL xios_read_field(TRIM(field6(1)%name)//'_start', field6); END IF
         IF (PRESENT(field7)) THEN; CALL xios_read_field(TRIM(field7(1)%name)//'_start', field7); END IF
         IF (PRESENT(field8)) THEN; CALL xios_read_field(TRIM(field8(1)%name)//'_start', field8); END IF
         IF (PRESENT(field9)) THEN; CALL xios_read_field(TRIM(field9(1)%name)//'_start', field9); END IF
         IF (PRESENT(field10)) THEN; CALL xios_read_field(TRIM(field10(1)%name)//'_start', field10); END IF
         IF (PRESENT(field11)) THEN; CALL xios_read_field(TRIM(field11(1)%name)//'_start', field11); END IF
         IF (PRESENT(field12)) THEN; CALL xios_read_field(TRIM(field12(1)%name)//'_start', field12); END IF
         IF (PRESENT(field13)) THEN; CALL xios_read_field(TRIM(field13(1)%name)//'_start', field13); END IF
         IF (PRESENT(field14)) THEN; CALL xios_read_field(TRIM(field14(1)%name)//'_start', field14); END IF
         IF (PRESENT(field15)) THEN; CALL xios_read_field(TRIM(field15(1)%name)//'_start', field15); END IF
         IF (PRESENT(field16)) THEN; CALL xios_read_field(TRIM(field16(1)%name)//'_start', field16); END IF
         IF (PRESENT(field17)) THEN; CALL xios_read_field(TRIM(field17(1)%name)//'_start', field17); END IF
         IF (PRESENT(field18)) THEN; CALL xios_read_field(TRIM(field18(1)%name)//'_start', field18); END IF
         IF (PRESENT(field19)) THEN; CALL xios_read_field(TRIM(field19(1)%name)//'_start', field19); END IF

!      CALL xios_recv_field("it_start",it_real)
         CALL xios_read_var("it_start", it_real)
         it = it_real
      ELSE

         !$OMP MASTER

         nfield = 0
         IF (PRESENT(field0)) THEN; nfield = nfield + 1; field_array(nfield)%field => field0; END IF
         IF (PRESENT(field1)) THEN; nfield = nfield + 1; field_array(nfield)%field => field1; END IF
         IF (PRESENT(field2)) THEN; nfield = nfield + 1; field_array(nfield)%field => field2; END IF
         IF (PRESENT(field3)) THEN; nfield = nfield + 1; field_array(nfield)%field => field3; END IF
         IF (PRESENT(field4)) THEN; nfield = nfield + 1; field_array(nfield)%field => field4; END IF
         IF (PRESENT(field5)) THEN; nfield = nfield + 1; field_array(nfield)%field => field5; END IF
         IF (PRESENT(field6)) THEN; nfield = nfield + 1; field_array(nfield)%field => field6; END IF
         IF (PRESENT(field7)) THEN; nfield = nfield + 1; field_array(nfield)%field => field7; END IF
         IF (PRESENT(field8)) THEN; nfield = nfield + 1; field_array(nfield)%field => field8; END IF
         IF (PRESENT(field9)) THEN; nfield = nfield + 1; field_array(nfield)%field => field9; END IF
         IF (PRESENT(field10)) THEN; nfield = nfield + 1; field_array(nfield)%field => field10; END IF
         IF (PRESENT(field11)) THEN; nfield = nfield + 1; field_array(nfield)%field => field11; END IF
         IF (PRESENT(field12)) THEN; nfield = nfield + 1; field_array(nfield)%field => field12; END IF
         IF (PRESENT(field13)) THEN; nfield = nfield + 1; field_array(nfield)%field => field13; END IF
         IF (PRESENT(field14)) THEN; nfield = nfield + 1; field_array(nfield)%field => field14; END IF
         IF (PRESENT(field15)) THEN; nfield = nfield + 1; field_array(nfield)%field => field15; END IF
         IF (PRESENT(field16)) THEN; nfield = nfield + 1; field_array(nfield)%field => field16; END IF
         IF (PRESENT(field17)) THEN; nfield = nfield + 1; field_array(nfield)%field => field17; END IF
         IF (PRESENT(field18)) THEN; nfield = nfield + 1; field_array(nfield)%field => field18; END IF
         IF (PRESENT(field19)) THEN; nfield = nfield + 1; field_array(nfield)%field => field19; END IF

         IF (is_mpi_root) THEN
            status = NF90_OPEN(TRIM(ADJUSTL(start_file_name))//'.nc', NF90_NOWRITE, ncid)
         END IF

         DO nf = 1, nfield
            field => field_array(nf)%field
            status = nf90_inq_varid(ncid, TRIM(ADJUSTL(field(1)%name)), fieldId(nf))
            status = nf90_get_att(ncid, NF90_GLOBAL, "iteration", it)
            CALL read_start_field(field, fieldId(nf), ncid)
         END DO

         IF (is_mpi_root) THEN
            status = NF90_CLOSE(ncid)
         END IF

         !$OMP END MASTER

      END IF

   END SUBROUTINE read_start_fields

   SUBROUTINE read_start_i0(name, var)
      USE prec
      USE metric
      USE grid_param
      USE field_mod
      USE domain_mod
      USE netcdf_mod
      USE mpipara
      USE getin_mod
      USE spherical_geom_mod
      USE transfert_mod
      USE xios_mod
      USE omp_para

      IMPLICIT NONE
      CHARACTER(LEN=*), INTENT(IN)  :: name
      INTEGER, INTENT(OUT)           :: var
      CHARACTER(LEN=255) :: start_file_name
      REAL(rstd) :: var_real
      CHARACTER(LEN=255) :: suffix
      IF (no_io) RETURN

      start_file_name = "start"
      CALL getin("start_file_name", start_file_name)

      IF (is_read_interp) THEN
         suffix = "_start_interp"
      ELSE
         suffix = "_start"
      END IF

      CALL xios_read_var(TRIM(name)//TRIM(suffix), var_real)
      var = INT(var_real)
   END SUBROUTINE read_start_i0

   SUBROUTINE read_start_r0(name, var)
      USE prec
      USE metric
      USE grid_param
      USE field_mod
      USE domain_mod
      USE netcdf_mod
      USE mpipara
      USE getin_mod
      USE spherical_geom_mod
      USE transfert_mod
      USE xios_mod
      USE omp_para

      IMPLICIT NONE
      CHARACTER(LEN=*), INTENT(IN)  :: name
      REAL, INTENT(OUT)           :: var
      CHARACTER(LEN=255) :: start_file_name
      CHARACTER(LEN=255) :: suffix
      IF (no_io) RETURN

      start_file_name = "start"
      CALL getin("start_file_name", start_file_name)

      IF (is_read_interp) THEN
         suffix = "_start_interp"
      ELSE
         suffix = "_start"
      END IF

      CALL xios_read_var(TRIM(name)//TRIM(suffix), var)
   END SUBROUTINE read_start_r0

   SUBROUTINE read_xios_start_field(fieldId, suffix, f_field)
      USE field_mod
      USE grid_param
      USE domain_mod
      USE xios_mod
      USE omp_para
      USE tracer_icosa_mod
      IMPLICIT NONE
      CHARACTER(LEN=*), INTENT(IN)        :: fieldId
      CHARACTER(LEN=*), INTENT(IN)        :: suffix
      TYPE(t_field), POINTER, SAVE         :: f_tracer(:)
      TYPE(t_field), POINTER, DIMENSION(:) :: f_field
      INTEGER                            :: iq, ind, l

      IF (TRIM(fieldId) == "q") THEN! hook for tracers
         CALL allocate_field(f_tracer, field_t, type_real, llm)

         DO iq = 1, nqtot
            IF (xios_is_valid_field(TRIM(tracers(iq)%name)//TRIM(suffix))) THEN
               CALL xios_read_field(TRIM(tracers(iq)%name)//TRIM(suffix), f_tracer)
               DO ind = 1, ndomain
                  IF (.NOT. assigned_domain(ind)) CYCLE
                  DO l = ll_begin, ll_end
                     f_field(ind)%rval4d(:, l, iq) = f_tracer(ind)%rval3d(:, l)
                  END DO
               END DO
               tracers(iq)%already_initialized = .TRUE.
            ELSE
               IF (.NOT. tracers(iq)%already_initialized .AND. tracers(iq)%has_default_init_value) THEN
                  DO ind = 1, ndomain
                     IF (.NOT. assigned_domain(ind)) CYCLE
                     DO l = ll_begin, ll_end
                        f_field(ind)%rval4d(:, l, iq) = tracers(iq)%default_init_value
                     END DO
                  END DO
                  tracers(iq)%already_initialized = .TRUE.
               END IF
            END IF
         END DO
         CALL deallocate_field(f_tracer)
      ELSE
         CALL xios_read_field(TRIM(fieldId)//TRIM(suffix), f_field)
      END IF

   END SUBROUTINE read_xios_start_field

   SUBROUTINE read_start_field(field, fieldId, ncid)
      USE grid_param, ONLY: llm, nqtot
      TYPE(t_field), POINTER :: field(:)
      INTEGER, INTENT(IN)     :: fieldId
      INTEGER, INTENT(IN)     :: ncid

      TYPE(t_domain), POINTER :: d
      TYPE(t_field), POINTER :: field_glo(:)
      REAL(rstd), ALLOCATABLE :: global_field2d(:)
      REAL(rstd), ALLOCATABLE :: global_field3d(:, :)
      REAL(rstd), ALLOCATABLE :: global_field4d(:, :, :)
      INTEGER :: i, j, ij, k, e, ind, ind_glo
      INTEGER :: ndim, field_type
      INTEGER :: status

      ndim = field(1)%ndim
      field_Type = field(1)%field_type

      IF (is_mpi_root) THEN

         IF (ndim == 2) THEN
            CALL allocate_field_glo(field_glo, field(1)%field_type, field(1)%data_type, name=field(1)%name)
         ELSE IF (ndim == 3) THEN
            CALL allocate_field_glo(field_glo, field(1)%field_type, field(1)%data_type, field(1)%dim3, name=field(1)%name)
         ELSE IF (ndim == 4) THEN
       CALL allocate_field_glo(field_glo, field(1)%field_type, field(1)%data_type, field(1)%dim3, field(1)%dim4, name=field(1)%name)
         END IF
      END IF

      IF (is_mpi_root) THEN

         IF (field_type == field_T) THEN
            IF (ndim == 2) THEN
               ALLOCATE (global_field2d(ncell_glo))
               status = NF90_GET_VAR(ncid, fieldid, global_field2d, start=(/1/), count=(/ncell_glo/))
               DO ind = 1, ndomain_glo
                  d => domain_glo(ind)
                  DO j = d%jj_begin, d%jj_end
                     DO i = d%ii_begin, d%ii_end
                        ij = (j - 1)*d%iim + i
                        ind_glo = d%assign_cell_glo(i, j)
                        field_glo(ind)%rval2d(ij) = global_field2d(ind_glo)
                     END DO
                  END DO
               END DO

            ELSE IF (ndim == 3) THEN
               ALLOCATE (global_field3d(ncell_glo, llm))
               status = NF90_GET_VAR(ncid, fieldid, global_field3d, start=(/1, 1/), count=(/ncell_glo, llm/))
               DO ind = 1, ndomain_glo
                  d => domain_glo(ind)
                  DO j = d%jj_begin, d%jj_end
                     DO i = d%ii_begin, d%ii_end
                        ij = (j - 1)*d%iim + i
                        ind_glo = d%assign_cell_glo(i, j)
                        field_glo(ind)%rval3d(ij, :) = global_field3d(ind_glo, :)
                     END DO
                  END DO
               END DO
            ELSE IF (ndim == 4) THEN
               ALLOCATE (global_field4d(ncell_glo, llm, nqtot))
               status = NF90_GET_VAR(ncid, fieldid, global_field4d, start=(/1, 1, 1/), count=(/ncell_glo, llm, nqtot/))
               DO ind = 1, ndomain_glo
                  d => domain_glo(ind)
                  DO j = d%jj_begin, d%jj_end
                     DO i = d%ii_begin, d%ii_end
                        ij = (j - 1)*d%iim + i
                        ind_glo = d%assign_cell_glo(i, j)
                        field_glo(ind)%rval4d(ij, :, :) = global_field4d(ind_glo, :, :)
                     END DO
                  END DO
               END DO
            END IF

         ELSE IF (field_type == field_U) THEN

            IF (ndim == 2) THEN
               ALLOCATE (global_field2d(3*ncell_glo))
               status = NF90_GET_VAR(ncid, fieldid, global_field2d, start=(/1/), count=(/3*ncell_glo/))

               DO ind = 1, ndomain_glo
                  d => domain_glo(ind)
                  DO j = d%jj_begin, d%jj_end
                     DO i = d%ii_begin, d%ii_end
                        DO k = 0, 5
                           ij = (j - 1)*d%iim + i
                           ind_glo = d%assign_cell_glo(i, j)
                           e = cell_glo(ind_glo)%edge(MOD(k + d%delta(i, j) + 6, 6))
                           field_glo(ind)%rval2d(ij + d%u_pos(k)) = global_field2d(ind_glo)*d%edge_assign_sign(k, i, j)
                        END DO
                     END DO
                  END DO
               END DO
            ELSE IF (ndim == 3) THEN
               ALLOCATE (global_field3d(3*ncell_glo, llm))
               status = NF90_GET_VAR(ncid, fieldid, global_field3d, start=(/1, 1/), count=(/3*ncell_glo, llm/))
               DO ind = 1, ndomain_glo
                  d => domain_glo(ind)
                  DO j = d%jj_begin, d%jj_end
                     DO i = d%ii_begin, d%ii_end
                        DO k = 0, 5
                           ij = (j - 1)*d%iim + i
                           ind_glo = d%assign_cell_glo(i, j)
                           e = cell_glo(ind_glo)%edge(MOD(k + d%delta(i, j) + 6, 6))
                           field_glo(ind)%rval3d(ij + d%u_pos(k + 1), :) = global_field3d(e, :)*d%edge_assign_sign(k, i, j)
                        END DO
                     END DO
                  END DO
               END DO
            ELSE IF (ndim == 4) THEN
               ALLOCATE (global_field4d(3*ncell_glo, llm, nqtot))
               status = NF90_GET_VAR(ncid, fieldid, global_field4d, start=(/1, 1, 1/), count=(/3*ncell_glo, llm, nqtot/))
               DO ind = 1, ndomain_glo
                  d => domain_glo(ind)
                  DO j = d%jj_begin, d%jj_end
                     DO i = d%ii_begin, d%ii_end
                        DO k = 0, 5
                           ij = (j - 1)*d%iim + i
                           ind_glo = d%assign_cell_glo(i, j)
                           e = cell_glo(ind_glo)%edge(MOD(k + d%delta(i, j) + 6, 6))
                           field_glo(ind)%rval4d(ij + d%u_pos(k + 1), :, :) = global_field4d(e, :, :)*d%edge_assign_sign(k, i, j)
                        END DO
                     END DO
                  END DO
               END DO
            END IF

         END IF
      END IF

      CALL scatter_field(field_glo, field)

      IF (is_mpi_root) THEN
         CALL deallocate_field_glo(field_glo)
      END IF

   END SUBROUTINE read_start_field

END MODULE restart_mod

