diff --git a/LIBTOOLS/tools/lfi2cdf/src/mode_util.f90 b/LIBTOOLS/tools/lfi2cdf/src/mode_util.f90
index cc6146d7e9847bf21fe2b83984f49a98a44ba8cf..fd7bf509740c1311b6bd1ccf2baace47c1226775 100644
--- a/LIBTOOLS/tools/lfi2cdf/src/mode_util.f90
+++ b/LIBTOOLS/tools/lfi2cdf/src/mode_util.f90
@@ -12,11 +12,12 @@
 !  P. Wautelet 19/09/2019: add possibility to provide a fallback file if some information are not found in the input file
 !  P. Wautelet 21/10/2019: add OPTDIR option to set directory for writing outfiles
 !  P. Wautelet 21/10/2019: if DTMOD and DTCUR not found, try to read the time coordinate
+!  P. Wautelet 10/11/2020: new data structures for netCDF dimensions
 !-----------------------------------------------------------------
 MODULE mode_util
   use modd_field,      only: tfielddata, tfieldlist
   USE MODD_IO,         ONLY: TFILEDATA, TFILE_ELT
-  USE MODD_NETCDF,     ONLY: DIMCDF, CDFINT
+  USE MODD_NETCDF,     ONLY: CDFINT, tdimnc
   USE MODD_PARAMETERS, ONLY: NLFIMAXCOMMENTLENGTH, NMNHNAMELGTMAX
   use modd_precision,  only: LFIINT
 
@@ -56,7 +57,7 @@ MODULE mode_util
      INTEGER,DIMENSION(MAXRAW)             :: src    ! List of variables used to compute the variable (needed only if calc=.true.)
      INTEGER                               :: tgt    ! Target: id of the variable that use it (calc variable)
      TYPE(TFIELDDATA)                      :: TFIELD ! Metadata about the field
-     TYPE(DIMCDF),DIMENSION(:),ALLOCATABLE :: TDIMS  ! Dimensions of the field
+     TYPE(tdimnc),DIMENSION(:),ALLOCATABLE :: TDIMS  ! Dimensions of the field
   END TYPE workfield
 
   LOGICAL(KIND=LFIINT), PARAMETER :: ltrue  = .TRUE.
@@ -428,7 +429,7 @@ END DO
               CALL IO_Dimids_guess_nc4(outfiles(idx_out)%TFILE,tpreclist(ji)%TFIELD,&
                                       tpreclist(ji)%NSIZE,tpreclist(ji)%TDIMS,IRESP)
               !
-              IF (IRESP/=0 .OR. tpreclist(ji)%TDIMS(3)%LEN==1) THEN
+              IF (IRESP/=0 .OR. tpreclist(ji)%TDIMS(3)%nlen==1) THEN
                 CALL PRINT_MSG(NVERB_DEBUG,'IO','parse_infiles',tpreclist(ji)%TFIELD%CMNHNAME//': try 2D')
                 !Try again with 2D
                 tpreclist(ji)%TFIELD%NDIMS = 2
@@ -436,7 +437,7 @@ END DO
                                         tpreclist(ji)%NSIZE,tpreclist(ji)%TDIMS,IRESP)
               END IF
               !
-              IF (IRESP/=0 .OR. tpreclist(ji)%TDIMS(2)%LEN==1) THEN
+              IF (IRESP/=0 .OR. tpreclist(ji)%TDIMS(2)%nlen==1) THEN
                 CALL PRINT_MSG(NVERB_DEBUG,'IO','parse_infiles',tpreclist(ji)%TFIELD%CMNHNAME//': try 1D')
                 !Try again with 1D
                 tpreclist(ji)%TFIELD%NDIMS = 1
@@ -732,7 +733,7 @@ END DO
 
       SELECT CASE(tpreclist(ji)%TFIELD%NTYPE)
       CASE (TYPEINT)
-        IDIMLEN(1:IDIMS) = tpreclist(ji)%TDIMS(1:IDIMS)%LEN
+        IDIMLEN(1:IDIMS) = tpreclist(ji)%TDIMS(1:IDIMS)%nlen
 
         IF (.NOT.tpreclist(ji)%calc) THEN
           INSRC = 1
@@ -808,7 +809,7 @@ END DO
 
 
       CASE (TYPELOG)
-        IDIMLEN(1:IDIMS) = tpreclist(ji)%TDIMS(1:IDIMS)%LEN
+        IDIMLEN(1:IDIMS) = tpreclist(ji)%TDIMS(1:IDIMS)%nlen
 
         tpreclist(ji)%TFIELD%LTIMEDEP = gtimedep_in(ji)
         SELECT CASE(IDIMS)
@@ -832,7 +833,7 @@ END DO
 
 
       CASE (TYPEREAL)
-        IDIMLEN(1:IDIMS) = tpreclist(ji)%TDIMS(1:IDIMS)%LEN
+        IDIMLEN(1:IDIMS) = tpreclist(ji)%TDIMS(1:IDIMS)%nlen
 
         IF (.NOT.tpreclist(ji)%calc) THEN
           INSRC = 1
@@ -1593,14 +1594,14 @@ END DO
 
   SUBROUTINE IO_Dims_fill_nc4(TPFILE,TPREC,KRESP)
     USE MODD_IO,           ONLY: TFILEDATA
-    use mode_io_tools_nc4, only: IO_Dimcdf_get_nc4, IO_Dim_find_byname_nc4
+    use mode_io_tools_nc4, only: IO_Dim_find_create_nc4, IO_Dim_find_byname_nc4
 
     TYPE(TFILEDATA),INTENT(IN)    :: TPFILE
     TYPE(workfield),INTENT(INOUT) :: TPREC
     INTEGER,        INTENT(OUT)   :: KRESP
 
+    integer              :: iidx
     INTEGER              :: JJ
-    TYPE(DIMCDF),POINTER :: TZDIMPTR
 
     CALL PRINT_MSG(NVERB_DEBUG,'IO','IO_Dims_fill_nc4','called')
 
@@ -1622,12 +1623,12 @@ END DO
       CALL IO_Dim_find_byname_nc4(TPFILE,TPREC%CDIMNAMES_FILE(JJ),TPREC%TDIMS(JJ),KRESP)
       !If dimension not found => create it
       IF (KRESP/=0)  THEN
-        TZDIMPTR => IO_Dimcdf_get_nc4(TPFILE,TPREC%NDIMSIZES_FILE(JJ))
-        TPREC%TDIMS(JJ) = TZDIMPTR
+        call IO_Dim_find_create_nc4( tpfile, tprec%ndimsizes_file(jj), iidx )
+        tprec%tdims(jj) = tpfile%tncdims%tdims(iidx)
         KRESP = 0
       END IF
-      IF (TRIM(TPREC%TDIMS(JJ)%name)/='time' .AND. &
-        TPREC%TDIMS(JJ)%len /= TPREC%NDIMSIZES_FILE(JJ)) THEN
+      IF (TRIM(TPREC%TDIMS(JJ)%cname)/='time' .AND. &
+        TPREC%TDIMS(JJ)%nlen /= TPREC%NDIMSIZES_FILE(JJ)) THEN
         CALL PRINT_MSG(NVERB_WARNING,'IO','IO_FILL_DIMS_NC4','problem with dimensions for '//TPREC%TFIELD%CMNHNAME)
         KRESP = -3
         EXIT
diff --git a/src/LIB/SURCOUCHE/src/modd_field.f90 b/src/LIB/SURCOUCHE/src/modd_field.f90
index 62d113c9eee7fc0892e81b3c1c6f4bb786f4fab7..f2289401ea10dd101089e0c6ad7adc15acf425ee 100644
--- a/src/LIB/SURCOUCHE/src/modd_field.f90
+++ b/src/LIB/SURCOUCHE/src/modd_field.f90
@@ -11,6 +11,7 @@
 !  P. Wautelet 23/01/2020: split in modd_field.f90 and mode_field.f90
 !  P. Wautelet 27/01/2020: create the tfield_metadata_base abstract datatype
 !  P. Wautelet 14/09/2020: add ndimlist field to tfield_metadata_base
+!  P. Wautelet 10/11/2020: new data structures for netCDF dimensions
 !-----------------------------------------------------------------
 module modd_field
 
@@ -28,64 +29,83 @@ INTEGER,PARAMETER :: TYPEUNDEF = -1, TYPEINT = 1, TYPELOG = 2, TYPEREAL = 3, TYP
 !
 integer, parameter :: NMNHDIM_UNKNOWN             = -2
 
+!For efficient use of memory, it is better that all values for real dimensions be contiguous
+integer, parameter :: NMNHDIM_NI                  = 1
+integer, parameter :: NMNHDIM_NJ                  = 2
+integer, parameter :: NMNHDIM_NI_U                = 3
+integer, parameter :: NMNHDIM_NJ_U                = 4
+integer, parameter :: NMNHDIM_NI_V                = 5
+integer, parameter :: NMNHDIM_NJ_V                = 6
+integer, parameter :: NMNHDIM_LEVEL               = 7
+integer, parameter :: NMNHDIM_LEVEL_W             = 8
+integer, parameter :: NMNHDIM_TIME                = 9
+
 integer, parameter :: NMNHDIM_ONE                 = 10
+
+integer, parameter :: NMNHDIM_LASTDIM_NODIACHRO   = 10  ! Index of the last defined dimension for non-diachronic files
+
 integer, parameter :: NMNHDIM_COMPLEX             = 11
 
-integer, parameter :: NMNHDIM_NI                  = 21
-integer, parameter :: NMNHDIM_NJ                  = 22
-integer, parameter :: NMNHDIM_NI_U                = 23
-integer, parameter :: NMNHDIM_NJ_U                = 24
-integer, parameter :: NMNHDIM_NI_V                = 25
-integer, parameter :: NMNHDIM_NJ_V                = 26
-integer, parameter :: NMNHDIM_LEVEL               = 27
-integer, parameter :: NMNHDIM_LEVEL_W             = 28
-integer, parameter :: NMNHDIM_TIME                = 29
-
-integer, parameter :: NMNHDIM_BUDGET_CART_NI      = 30
-integer, parameter :: NMNHDIM_BUDGET_CART_NJ      = 31
-integer, parameter :: NMNHDIM_BUDGET_CART_NI_U    = 32
-integer, parameter :: NMNHDIM_BUDGET_CART_NJ_U    = 33
-integer, parameter :: NMNHDIM_BUDGET_CART_NI_V    = 34
-integer, parameter :: NMNHDIM_BUDGET_CART_NJ_V    = 35
-integer, parameter :: NMNHDIM_BUDGET_CART_LEVEL   = 36
-integer, parameter :: NMNHDIM_BUDGET_CART_LEVEL_W = 37
-
-integer, parameter :: NMNHDIM_BUDGET_MASK_LEVEL   = 40
-integer, parameter :: NMNHDIM_BUDGET_MASK_LEVEL_W = 41
-integer, parameter :: NMNHDIM_BUDGET_MASK_TIME    = 42
-integer, parameter :: NMNHDIM_BUDGET_MASK_NBUMASK = 43
-
-integer, parameter :: NMNHDIM_BUDGET_LES_TIME     = 50
-integer, parameter :: NMNHDIM_BUDGET_LES_AVG_TIME = 51
-integer, parameter :: NMNHDIM_BUDGET_LES_LEVEL    = 52
-integer, parameter :: NMNHDIM_BUDGET_LES_SV       = 53
-integer, parameter :: NMNHDIM_BUDGET_LES_MASK     = 54
-
-integer, parameter :: NMNHDIM_SPECTRA_2PTS_NI     = 60
-integer, parameter :: NMNHDIM_SPECTRA_2PTS_NJ     = 61
-integer, parameter :: NMNHDIM_SPECTRA_SPEC_NI     = 62
-integer, parameter :: NMNHDIM_SPECTRA_SPEC_NJ     = 63
-integer, parameter :: NMNHDIM_SPECTRA_LEVEL       = 64
-
-integer, parameter :: NMNHDIM_SERIES_LEVEL        = 70
-integer, parameter :: NMNHDIM_SERIES_LEVEL_W      = 71
-integer, parameter :: NMNHDIM_SERIES_TIME         = 72  ! Time dimension for time series
-
-integer, parameter :: NMNHDIM_FLYER_TIME          = 80  ! Time dimension for aircraft/balloon (dimension local to each flyer)
-integer, parameter :: NMNHDIM_PROFILER_TIME       = 81  ! Time dimension for profilers
-integer, parameter :: NMNHDIM_STATION_TIME        = 82  ! Time dimension for stations
-
-integer, parameter :: NMNHDIM_BUDGET_NGROUPS      = 100 ! This is not a true dimension
-integer, parameter :: NMNHDIM_FLYER_PROC          = 101 ! This is not a true dimension
-integer, parameter :: NMNHDIM_PROFILER_PROC       = 102 ! This is not a true dimension
-integer, parameter :: NMNHDIM_STATION_PROC        = 103 ! This is not a true dimension
-integer, parameter :: NMNHDIM_SERIES_PROC         = 104 ! This is not a true dimension
-integer, parameter :: NMNHDIM_BUDGET_TERM         = 105 ! This is not a true dimension
+integer, parameter :: NMNHDIM_BUDGET_CART_NI      = 12
+integer, parameter :: NMNHDIM_BUDGET_CART_NJ      = 13
+integer, parameter :: NMNHDIM_BUDGET_CART_NI_U    = 14
+integer, parameter :: NMNHDIM_BUDGET_CART_NJ_U    = 15
+integer, parameter :: NMNHDIM_BUDGET_CART_NI_V    = 16
+integer, parameter :: NMNHDIM_BUDGET_CART_NJ_V    = 17
+integer, parameter :: NMNHDIM_BUDGET_CART_LEVEL   = 18
+integer, parameter :: NMNHDIM_BUDGET_CART_LEVEL_W = 19
+
+integer, parameter :: NMNHDIM_BUDGET_MASK_LEVEL   = 20
+integer, parameter :: NMNHDIM_BUDGET_MASK_LEVEL_W = 21
+integer, parameter :: NMNHDIM_BUDGET_MASK_TIME    = 22
+integer, parameter :: NMNHDIM_BUDGET_MASK_NBUMASK = 23
+
+integer, parameter :: NMNHDIM_BUDGET_LES_TIME     = 24
+integer, parameter :: NMNHDIM_BUDGET_LES_AVG_TIME = 25
+integer, parameter :: NMNHDIM_BUDGET_LES_LEVEL    = 26
+integer, parameter :: NMNHDIM_BUDGET_LES_SV       = 27
+integer, parameter :: NMNHDIM_BUDGET_LES_MASK     = 100 ! This is not a true dimension
+
+integer, parameter :: NMNHDIM_SPECTRA_2PTS_NI     = 28
+integer, parameter :: NMNHDIM_SPECTRA_2PTS_NJ     = 29
+integer, parameter :: NMNHDIM_SPECTRA_SPEC_NI     = 30
+integer, parameter :: NMNHDIM_SPECTRA_SPEC_NJ     = 31
+integer, parameter :: NMNHDIM_SPECTRA_LEVEL       = 32
+
+integer, parameter :: NMNHDIM_SERIES_LEVEL        = 33
+integer, parameter :: NMNHDIM_SERIES_LEVEL_W      = 34
+integer, parameter :: NMNHDIM_SERIES_TIME         = 35  ! Time dimension for time series
+
+integer, parameter :: NMNHDIM_FLYER_TIME          = 36  ! Time dimension for aircraft/balloon (dimension local to each flyer)
+integer, parameter :: NMNHDIM_PROFILER_TIME       = 37  ! Time dimension for profilers
+integer, parameter :: NMNHDIM_STATION_TIME        = 38  ! Time dimension for stations
+
+integer, parameter :: NMNHDIM_LASTDIM_DIACHRO     = 38  ! Index of the last defined dimension for diachronic files
+
+integer, parameter :: NMNHDIM_BUDGET_NGROUPS      = 101 ! This is not a true dimension
+integer, parameter :: NMNHDIM_FLYER_PROC          = 102 ! This is not a true dimension
+integer, parameter :: NMNHDIM_PROFILER_PROC       = 103 ! This is not a true dimension
+integer, parameter :: NMNHDIM_STATION_PROC        = 104 ! This is not a true dimension
+integer, parameter :: NMNHDIM_SERIES_PROC         = 105 ! This is not a true dimension
+integer, parameter :: NMNHDIM_BUDGET_TERM         = 106 ! This is not a true dimension
 
 integer, parameter :: NMNHDIM_NOTLISTED           = 200 ! Special case for valid dimension but not in this parameter list
 
 integer, parameter :: NMNHDIM_UNUSED              = 300
-!
+
+!Array to allow easy identification of dimensions for Arakawa grid points
+integer, dimension(0:8,3), parameter :: NMNHDIM_ARAKAWA = reshape( [ &
+  NMNHDIM_UNKNOWN, NMNHDIM_UNKNOWN, NMNHDIM_UNKNOWN, & ! dummy point (to treat ngrid=0 without crash)
+  NMNHDIM_NI,      NMNHDIM_NJ,      NMNHDIM_LEVEL,   & ! mass point
+  NMNHDIM_NI_U,    NMNHDIM_NJ_U,    NMNHDIM_LEVEL,   & ! u point
+  NMNHDIM_NI_V,    NMNHDIM_NJ_V,    NMNHDIM_LEVEL,   & ! v point
+  NMNHDIM_NI,      NMNHDIM_NJ,      NMNHDIM_LEVEL_W, & ! w point
+  NMNHDIM_NI_U,    NMNHDIM_NJ_V,    NMNHDIM_LEVEL,   & ! xi vorticity point (=f point =uv point)
+  NMNHDIM_NI_U,    NMNHDIM_NJ_U,    NMNHDIM_LEVEL_W, & ! eta vorticity point (=uw point)
+  NMNHDIM_NI_V,    NMNHDIM_NJ_V,    NMNHDIM_LEVEL_W, & ! zeta vorticity point (=vw point)
+  NMNHDIM_NI_U,    NMNHDIM_NJ_V,    NMNHDIM_LEVEL_W] & ! fw point (=uvw point)
+  , shape = [ 9, 3 ], order = [ 2, 1 ] )
+
 TYPE TFIELDPTR_C0D
   CHARACTER(LEN=:),     POINTER :: DATA => NULL()
 END TYPE TFIELDPTR_C0D
diff --git a/src/LIB/SURCOUCHE/src/modd_io.f90 b/src/LIB/SURCOUCHE/src/modd_io.f90
index fe3bc68f3c7f9f271a76aec6c1785327aa0d5ddf..77c77b43787ca147b242fdfa67728b93de1ce56c 100644
--- a/src/LIB/SURCOUCHE/src/modd_io.f90
+++ b/src/LIB/SURCOUCHE/src/modd_io.f90
@@ -11,13 +11,14 @@
 !  P. Wautelet 12/03/2019: add TMAINFILE field in TFILEDATA
 !  P. Wautelet 17/01/2020: add 'BUD' category for Print_msg + corresponding namelist variables
 !  P. Wautelet 22/09/2020: add ldimreduced in tfiledata
+!  P. Wautelet 10/11/2020: new data structures for netCDF dimensions
 !-----------------------------------------------------------------
 
 #define MNH_REDUCE_DIMENSIONS_IN_FILES 1
 
 MODULE MODD_IO
 !
-USE MODD_NETCDF,     ONLY: IOCDF, TPTR2DIMCDF
+use modd_netcdf,     only: tdimsnc
 USE MODD_PARAMETERS, ONLY: NDIRNAMELGTMAX, NFILENAMELGTMAX
 use modd_precision,  only: CDFINT, LFIINT
 !
@@ -112,15 +113,16 @@ TYPE TFILEDATA
   INTEGER              :: NLFIVERB  = 1  !LFI verbosity level
   INTEGER(KIND=LFIINT) :: NLFIFLU   = -1 !File identifier
   !
+#ifdef MNH_IOCDF4
   ! Fields for netCDF files
-  INTEGER(KIND=CDFINT) :: NNCID = -1 !File identifier
-  INTEGER(KIND=CDFINT) :: NNCNAR = 0 !Number of articles of the netCDF file (only accurate if file opened in read mode)
-  LOGICAL              :: LNCREDUCE_FLOAT_PRECISION = .FALSE. ! Reduce the precision of floats to single precision
-                                                                  ! instead of double precision
-  LOGICAL              :: LNCCOMPRESS = .FALSE. ! Do compression on fields
-  INTEGER(KIND=CDFINT) :: NNCCOMPRESS_LEVEL = 0 ! Compression level
-  TYPE(IOCDF),POINTER  :: TNCDIMS => NULL()     ! Structure containing netCDF dimensions
-  TYPE(TPTR2DIMCDF),DIMENSION(:,:),ALLOCATABLE :: TNCCOORDS ! Structure pointing to the coordinates variables
+  INTEGER(KIND=CDFINT)   :: NNCID = -1 !File identifier (corresponding to the actual group)
+  INTEGER(KIND=CDFINT)   :: NNCNAR = 0 !Number of articles of the netCDF file (only accurate if file opened in read mode)
+  LOGICAL                :: LNCREDUCE_FLOAT_PRECISION = .FALSE. ! Reduce the precision of floats to single precision
+                                                                ! instead of double precision
+  LOGICAL                :: LNCCOMPRESS = .FALSE. ! Do compression on fields
+  INTEGER(KIND=CDFINT)   :: NNCCOMPRESS_LEVEL = 0 ! Compression level
+  type(tdimsnc), pointer :: tncdims => Null()     ! Dimensions of netCDF file
+#endif
   !
   !Fields for other files
   INTEGER :: NLU = -1                      !Logical unit number
@@ -150,6 +152,6 @@ TYPE(TFILEDATA),POINTER,SAVE :: TFILE_SURFEX  => NULL() !Pointer used to find th
 TYPE(TFILEDATA),POINTER,SAVE :: TFILE_OUTPUTLISTING  => NULL() !Pointer used to point to the file used when writing to OUTPUT_LISTINGn file
 
 !Non existing file which can be used as a dummy target
-TYPE(TFILEDATA),TARGET, SAVE :: TFILE_DUMMY = TFILEDATA(CNAME="dummy",CDIRNAME=NULL(),TFILES_IOZ=NULL(),TNCCOORDS=NULL())
+TYPE(TFILEDATA),TARGET, SAVE :: TFILE_DUMMY = TFILEDATA(CNAME="dummy",CDIRNAME=NULL(),TFILES_IOZ=NULL())
 
 END MODULE MODD_IO
diff --git a/src/LIB/SURCOUCHE/src/modd_netcdf.f90 b/src/LIB/SURCOUCHE/src/modd_netcdf.f90
index b73380b5636a81333b7eeafce263fc062f56bf5d..0d289c67a511f8b20d8399613ac766c906daccd7 100644
--- a/src/LIB/SURCOUCHE/src/modd_netcdf.f90
+++ b/src/LIB/SURCOUCHE/src/modd_netcdf.f90
@@ -1,42 +1,35 @@
-!MNH_LIC Copyright 1994-2018 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2020 CNRS, Meteo-France and Universite Paul Sabatier
 !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
-!MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt  
+!MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
 !MNH_LIC for details. version 1.
 !-----------------------------------------------------------------
 ! Modifications:
-!  Philippe Wautelet: 05/2016-04/2018: new data structures and calls for I/O
+!  P. Wautelet 05/2016-04/2018: new data structures and calls for I/O
+!  P. Wautelet 10/11/2020: new data structures for netCDF dimensions
 !-----------------------------------------------------------------
-MODULE MODD_NETCDF
+module modd_netcdf
 
 use modd_precision, only: CDFINT
 
-IMPLICIT NONE
+implicit none
 
-TYPE IOCDF
-   TYPE(DIMCDF), POINTER :: DIM_NI      => NULL()
-   TYPE(DIMCDF), POINTER :: DIM_NJ      => NULL()
-   TYPE(DIMCDF), POINTER :: DIM_LEVEL   => NULL()
-   TYPE(DIMCDF), POINTER :: DIM_NI_U    => NULL()
-   TYPE(DIMCDF), POINTER :: DIM_NJ_U    => NULL()
-   TYPE(DIMCDF), POINTER :: DIM_NI_V    => NULL()
-   TYPE(DIMCDF), POINTER :: DIM_NJ_V    => NULL()
-   TYPE(DIMCDF), POINTER :: DIM_LEVEL_W => NULL()
-   TYPE(DIMCDF), POINTER :: DIMTIME     => NULL()
-   TYPE(DIMCDF), POINTER :: DIMSTR      => NULL()
-   TYPE(DIMCDF), POINTER :: DIMLIST     => NULL()
-END TYPE IOCDF
+public
 
-TYPE DIMCDF
-   CHARACTER(LEN=32)     :: NAME = ''
-   INTEGER(KIND=CDFINT)  :: LEN  = 0
-   INTEGER(KIND=CDFINT)  :: ID   = -1
-   TYPE(DIMCDF), POINTER :: NEXT => NULL()
-END TYPE DIMCDF
+integer, parameter :: NMAXDIMNAMELGTNC4 = 16
 
-TYPE TPTR2DIMCDF
-  TYPE(DIMCDF),POINTER :: TDIM
-END TYPE TPTR2DIMCDF
+!Datatype to store metadata of 1 dimension
+type tdimnc
+  character(len=NMAXDIMNAMELGTNC4) :: cname = ''
+  integer(kind=CDFINT)             :: nlen  = -1
+  integer(kind=CDFINT)             :: nid   = -1
+end type tdimnc
 
-TYPE(DIMCDF),TARGET :: TDIM_DUMMY !Dummy dimension
+!Datatype to store dimension metadata of a netCDF file
+type tdimsnc
+  integer :: nmaxdims = 0
+  type(tdimnc), dimension(:), allocatable :: tdims
+  integer :: nmaxdims_str = 0                          ! For character strings
+  type(tdimnc), dimension(:), allocatable :: tdims_str ! For character strings
+end type tdimsnc
 
-END MODULE MODD_NETCDF
+end module modd_netcdf
diff --git a/src/LIB/SURCOUCHE/src/mode_io_tools_nc4.f90 b/src/LIB/SURCOUCHE/src/mode_io_tools_nc4.f90
index 2af6508f8ff0990d5b0133d86b5e8c514698ab3d..ecd39b6fd51ae887ff319b7963d4775ca79f8ed6 100644
--- a/src/LIB/SURCOUCHE/src/mode_io_tools_nc4.f90
+++ b/src/LIB/SURCOUCHE/src/mode_io_tools_nc4.f90
@@ -13,13 +13,14 @@
 !  P. Wautelet 18/09/2019: correct support of 64bit integers (MNH_INT=8)
 !  P. Wautelet 14/09/2020: IO_Knowndims_set_nc4: add new dimensions + remove 'time' dimension in diachronic files
 !  P. Wautelet 14/09/2020: IO_Vdims_fill_nc4: use ndimlist when provided to fill dimensions ids
+!  P. Wautelet 10/11/2020: new data structures for netCDF dimensions
 !-----------------------------------------------------------------
 #ifdef MNH_IOCDF4
 module mode_io_tools_nc4
 
 use modd_field,     only: tfielddata
 use modd_io,        only: tfiledata
-use modd_netcdf,    only: dimcdf, iocdf, tdim_dummy
+use modd_netcdf,    only: tdimnc, tdimsnc
 use modd_precision, only: CDFINT
 
 use mode_msg
@@ -33,59 +34,64 @@ private
 
 public :: IO_Dim_find_byname_nc4, IO_Dimids_guess_nc4, IO_Knowndims_set_nc4
 public :: IO_Iocdf_alloc_nc4, IO_Iocdf_dealloc_nc4, IO_Mnhname_clean
-public :: IO_Dimcdf_get_nc4, IO_Strdimid_get_nc4, IO_Vdims_fill_nc4, IO_Err_handle_nc4
+public :: IO_Dim_find_create_nc4, IO_Strdimid_get_nc4, IO_Vdims_fill_nc4, IO_Err_handle_nc4
 
 contains
 
-SUBROUTINE IO_Dim_find_byname_nc4(TPFILE, HDIMNAME, TPDIM, KRESP)
-TYPE(TFILEDATA),         INTENT(IN)  :: TPFILE
-CHARACTER(LEN=*),        INTENT(IN)  :: HDIMNAME
-TYPE(DIMCDF),            INTENT(OUT) :: TPDIM
-INTEGER,                 INTENT(OUT) :: KRESP
-!
-TYPE(DIMCDF), POINTER :: TMP
-!
-CALL PRINT_MSG(NVERB_DEBUG,'IO','IO_Dim_find_byname_nc4','called for dimension name '//TRIM(HDIMNAME))
-!
-KRESP = -2
-!
-IF(.NOT.ASSOCIATED(TPFILE%TNCDIMS%DIMLIST)) THEN
-  CALL PRINT_MSG(NVERB_WARNING,'IO','IO_Dim_find_byname_nc4','DIMLIST not associated for file  '//TRIM(TPFILE%CNAME))
-  KRESP = -1
-  RETURN
-END IF
-!
-TMP => TPFILE%TNCDIMS%DIMLIST
-!
-DO WHILE(ASSOCIATED(TMP))
-  IF (TRIM(TMP%NAME)==TRIM(HDIMNAME)) THEN
-    TPDIM = TMP
-    KRESP = 0
-    EXIT
-  END IF
-  TMP => TMP%NEXT
-END DO
-!
-END SUBROUTINE IO_Dim_find_byname_nc4
+subroutine IO_Dim_find_byname_nc4( tpfile, hdimname, tpdim, kresp )
+type(tfiledata),  intent(in)  :: tpfile
+character(len=*), intent(in)  :: hdimname
+type(tdimnc),     intent(out) :: tpdim
+integer,          intent(out) :: kresp
+
+integer :: idx
+integer :: ji
+
+call Print_msg( NVERB_DEBUG, 'IO', 'IO_Dim_find_byname_nc4', 'called for dimension name ' // Trim( hdimname ) )
+
+kresp = -2
+
+if ( .not.Associated( tpfile%tncdims ) ) then
+  call Print_msg( NVERB_ERROR, 'IO', 'IO_Dim_find_byname_nc4', 'tncdims not associated for file  '//trim(tpfile%cname) )
+  kresp = -3
+  return
+end if
+
+if ( .not.Allocated( tpfile%tncdims%tdims ) ) then
+  call Print_msg( NVERB_WARNING, 'IO', 'IO_Dim_find_byname_nc4', 'tdims not allocated for file  '//trim(tpfile%cname) )
+  kresp = -1
+  return
+end if
+
+idx = -1
+do ji = 1, tpfile%tncdims%nmaxdims
+  if ( Trim( hdimname ) == Trim( tpfile%tncdims%tdims(ji)%cname ) ) then
+    idx = ji
+    kresp = 0
+    exit
+  end if
+end do
+
+end subroutine IO_Dim_find_byname_nc4
 
 
 SUBROUTINE IO_Dimids_guess_nc4(TPFILE, TPFIELD, KLEN, TPDIMS, KRESP)
 !
-USE MODD_FIELD, ONLY: TYPECHAR
+USE MODD_FIELD, ONLY: NMNHDIM_ARAKAWA, TYPECHAR
 !
 !Used by LFI2CDF
 TYPE(TFILEDATA),                      INTENT(IN)  :: TPFILE
 TYPE(TFIELDDATA),                     INTENT(IN)  :: TPFIELD
 INTEGER,                              INTENT(IN)  :: KLEN
-TYPE(DIMCDF),DIMENSION(:),            INTENT(OUT) :: TPDIMS
+TYPE(tdimnc),DIMENSION(:),            INTENT(OUT) :: TPDIMS
 INTEGER,                              INTENT(OUT) :: KRESP
 !
 INTEGER               :: IGRID
+integer               :: iidx
 INTEGER(kind=CDFINT)  :: ILEN, ISIZE
 INTEGER               :: JI
 CHARACTER(LEN=32)     :: YINT
 CHARACTER(LEN=2)      :: YDIR
-TYPE(DIMCDF), POINTER :: PTDIM
 !
 CALL PRINT_MSG(NVERB_DEBUG,'IO','IO_Dimids_guess_nc4','called for '//TRIM(TPFIELD%CMNHNAME))
 !
@@ -94,7 +100,6 @@ YDIR   =  TPFIELD%CDIR
 !
 KRESP = 0
 ILEN = 0
-PTDIM => NULL()
 !
 IF(IGRID<0 .OR. IGRID>8) THEN
   WRITE(YINT,'( I0 )') IGRID
@@ -115,9 +120,9 @@ IF (IGRID==0) THEN
         ILEN = 1
       END IF
     CASE (1)
-      PTDIM => IO_Dimcdf_get_nc4(TPFILE, int( KLEN, kind=CDFINT ) )
-      TPDIMS(1) = PTDIM
-      ILEN      = PTDIM%LEN
+      call IO_Dim_find_create_nc4( tpfile, klen, iidx )
+      tpdims(1) = tpfile%tncdims%tdims(iidx)
+      ilen = tpdims(1)%nlen
     CASE DEFAULT
       CALL PRINT_MSG(NVERB_WARNING,'IO','IO_Dimids_guess_nc4','NGRID=0 and NDIMS>1 not yet supported (field '&
                      //TRIM(TPFIELD%CMNHNAME)//')')
@@ -129,29 +134,25 @@ ELSE IF (TPFIELD%CLBTYPE/='NONE') THEN
   END IF
   !
   IF (TPFIELD%CLBTYPE=='LBX' .OR. TPFIELD%CLBTYPE=='LBXU') THEN
-    PTDIM => TPFILE%TNCCOORDS(2,IGRID)%TDIM
-    TPDIMS(2) = PTDIM
-    PTDIM => TPFILE%TNCCOORDS(3,IGRID)%TDIM
-    TPDIMS(3) = PTDIM
-    ILEN = TPDIMS(2)%LEN * TPDIMS(3)%LEN
+    tpdims(2) = tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,2) )
+    tpdims(3) = tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,3) )
+    ilen = tpdims(2)%nlen * tpdims(3)%nlen
     ISIZE = KLEN/ILEN
     IF (MOD(KLEN,ILEN)/=0) CALL PRINT_MSG(NVERB_WARNING,'IO','IO_Dimids_guess_nc4', &
                                               'can not guess 1st dimension for field '//TRIM(TPFIELD%CMNHNAME))
-    PTDIM => IO_Dimcdf_get_nc4(TPFILE, ISIZE)
-    TPDIMS(1) = PTDIM
-    ILEN       = ILEN * PTDIM%LEN
+    call IO_Dim_find_create_nc4( tpfile, isize, iidx )
+    tpdims(1) = tpfile%tncdims%tdims(iidx)
+    ilen = ilen * tpdims(1)%nlen
   ELSE IF (TPFIELD%CLBTYPE=='LBY' .OR. TPFIELD%CLBTYPE=='LBYV') THEN
-    PTDIM => TPFILE%TNCCOORDS(1,IGRID)%TDIM
-    TPDIMS(1) = PTDIM
-    PTDIM => TPFILE%TNCCOORDS(3,IGRID)%TDIM
-    TPDIMS(3) = PTDIM
-    ILEN = TPDIMS(1)%LEN * TPDIMS(3)%LEN
+    tpdims(1) = tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,1) )
+    tpdims(3) = tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,3) )
+    ilen = tpdims(1)%nlen * tpdims(3)%nlen
     ISIZE = KLEN/ILEN
     IF (MOD(KLEN,ILEN)/=0) CALL PRINT_MSG(NVERB_WARNING,'IO','IO_Dimids_guess_nc4', &
                                               'can not guess 2nd dimension for field '//TRIM(TPFIELD%CMNHNAME))
-    PTDIM => IO_Dimcdf_get_nc4(TPFILE, ISIZE)
-    TPDIMS(2) = PTDIM
-    ILEN       = ILEN * PTDIM%LEN
+    call IO_Dim_find_create_nc4( tpfile, isize, iidx )
+    tpdims(2) = tpfile%tncdims%tdims(iidx)
+    ilen = ilen * tpdims(2)%nlen
   ELSE
     CALL PRINT_MSG(NVERB_WARNING,'IO','IO_Dimids_guess_nc4','invalid CLBTYPE ('//TPFIELD%CLBTYPE//') for field '&
                      //TRIM(TPFIELD%CMNHNAME))
@@ -162,19 +163,19 @@ ELSE
   DO JI=1,TPFIELD%NDIMS
     IF (JI == 1) THEN
       IF ( (YDIR == 'XX' .OR. YDIR == 'XY') ) THEN
-        PTDIM => TPFILE%TNCCOORDS(1,IGRID)%TDIM
+        tpdims(1) = tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,1) )
       ELSE IF ( YDIR == 'YY' ) THEN
-        PTDIM => TPFILE%TNCCOORDS(2,IGRID)%TDIM
+        tpdims(1) = tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,2) )
       ELSE IF ( YDIR == 'ZZ' ) THEN
-        PTDIM => TPFILE%TNCCOORDS(3,IGRID)%TDIM
+        tpdims(1) = tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,3) )
       ELSE IF (JI==TPFIELD%NDIMS) THEN !Guess last dimension
-        PTDIM => IO_Dimcdf_get_nc4(TPFILE, int( KLEN, kind=CDFINT ) )
+        call IO_Dim_find_create_nc4( tpfile, klen, iidx )
+        tpdims(1) = tpfile%tncdims%tdims(iidx)
       END IF
-      ILEN       = PTDIM%LEN
-      TPDIMS(JI) = PTDIM
+      ilen = tpdims(1)%nlen
     ELSE IF (JI == 2) THEN
       IF ( YDIR == 'XY') THEN
-        PTDIM => TPFILE%TNCCOORDS(2,IGRID)%TDIM
+        tpdims(2) = tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,2) )
       ELSE IF (JI==TPFIELD%NDIMS) THEN !Guess last dimension
         ISIZE = KLEN/ILEN
         IF (MOD(KLEN,ILEN)/=0) THEN
@@ -182,21 +183,22 @@ ELSE
                                             'can not guess 2nd and last dimension for field '//TRIM(TPFIELD%CMNHNAME))
           EXIT
         END IF
-        PTDIM => IO_Dimcdf_get_nc4(TPFILE, ISIZE)
+        call IO_Dim_find_create_nc4( tpfile, isize, iidx )
+        tpdims(2) = tpfile%tncdims%tdims(iidx)
       ELSE
         CALL PRINT_MSG(NVERB_WARNING,'IO','IO_Dimids_guess_nc4','can not guess 2nd dimension for field '//TRIM(TPFIELD%CMNHNAME))
         EXIT
       END IF
-      ILEN       = ILEN * PTDIM%LEN
-      TPDIMS(JI) = PTDIM
+      ilen = ilen * tpdims(2)%nlen
     ELSE IF (JI == 3) THEN
       IF ( YDIR == 'XY' ) THEN
         IF (JI==TPFIELD%NDIMS .AND. KLEN/ILEN==1 .AND. MOD(KLEN,ILEN)==0) THEN
           !The last dimension is of size 1 => probably time dimension
           ISIZE = 1
-          PTDIM => IO_Dimcdf_get_nc4(TPFILE,ISIZE)
+          call IO_Dim_find_create_nc4( tpfile, isize, iidx )
+          tpdims(3) = tpfile%tncdims%tdims(iidx)
         ELSE
-          PTDIM => TPFILE%TNCCOORDS(3,IGRID)%TDIM
+          tpdims(3) = tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,3) )
         END IF
       ELSE IF (JI==TPFIELD%NDIMS) THEN !Guess last dimension
         ISIZE = KLEN/ILEN
@@ -205,13 +207,13 @@ ELSE
                                             'can not guess 3rd and last dimension for field '//TRIM(TPFIELD%CMNHNAME))
           EXIT
         END IF
-        PTDIM => IO_Dimcdf_get_nc4(TPFILE, ISIZE)
+        call IO_Dim_find_create_nc4( tpfile, isize, iidx )
+        tpdims(3) = tpfile%tncdims%tdims(iidx)
       ELSE
         CALL PRINT_MSG(NVERB_WARNING,'IO','IO_Dimids_guess_nc4','can not guess 3rd dimension for field '//TRIM(TPFIELD%CMNHNAME))
         EXIT
       END IF
-      ILEN       = ILEN * PTDIM%LEN
-      TPDIMS(JI) = PTDIM
+      ilen = ilen * tpdims(3)%nlen
     ELSE IF (JI==4 .AND. JI==TPFIELD%NDIMS) THEN !Guess last dimension
       ISIZE = KLEN/ILEN
       IF (MOD(KLEN,ILEN)/=0) THEN
@@ -219,9 +221,9 @@ ELSE
                                           'can not guess 4th and last dimension for field '//TRIM(TPFIELD%CMNHNAME))
         EXIT
       END IF
-      PTDIM => IO_Dimcdf_get_nc4(TPFILE, ISIZE)
-      ILEN       = ILEN * PTDIM%LEN
-      TPDIMS(JI) = PTDIM
+      call IO_Dim_find_create_nc4( tpfile, isize, iidx )
+      tpdims(4) = tpfile%tncdims%tdims(iidx)
+      ilen = ilen * tpdims(4)%nlen
     ELSE
       CALL PRINT_MSG(NVERB_WARNING,'IO','IO_Dimids_guess_nc4','can not guess dimension above 4 for field '&
                      //TRIM(TPFIELD%CMNHNAME))
@@ -247,6 +249,23 @@ USE MODD_CONF_n,        ONLY: CSTORAGE_TYPE
 USE MODD_DIM_n,         ONLY: NIMAX_ll, NJMAX_ll, NKMAX
 use modd_dyn,           only: xseglen
 use modd_dyn_n,         only: xtstep
+use modd_field,         only: NMNHDIM_NI, NMNHDIM_NJ, NMNHDIM_NI_U, NMNHDIM_NJ_U, NMNHDIM_NI_V, NMNHDIM_NJ_V,   &
+                              NMNHDIM_LEVEL, NMNHDIM_LEVEL_W, NMNHDIM_TIME,                                     &
+                              NMNHDIM_ONE,  NMNHDIM_COMPLEX,                                                    &
+                              NMNHDIM_BUDGET_CART_NI, NMNHDIM_BUDGET_CART_NJ, NMNHDIM_BUDGET_CART_NI_U,         &
+                              NMNHDIM_BUDGET_CART_NJ_U, NMNHDIM_BUDGET_CART_NI_V, NMNHDIM_BUDGET_CART_NJ_V,     &
+                              NMNHDIM_BUDGET_CART_LEVEL, NMNHDIM_BUDGET_CART_LEVEL_W,                           &
+                              NMNHDIM_BUDGET_MASK_LEVEL, NMNHDIM_BUDGET_MASK_LEVEL_W, NMNHDIM_BUDGET_MASK_TIME, &
+                              NMNHDIM_BUDGET_MASK_NBUMASK,                                                      &
+                              NMNHDIM_BUDGET_LES_TIME, NMNHDIM_BUDGET_LES_AVG_TIME, NMNHDIM_BUDGET_LES_LEVEL,   &
+                              NMNHDIM_BUDGET_LES_SV,                                                            &
+                              NMNHDIM_SPECTRA_2PTS_NI, NMNHDIM_SPECTRA_2PTS_NJ,                                 &
+                              NMNHDIM_SPECTRA_SPEC_NI, NMNHDIM_SPECTRA_SPEC_NJ, NMNHDIM_SPECTRA_LEVEL,          &
+                              NMNHDIM_SERIES_LEVEL, NMNHDIM_SERIES_LEVEL_W,                                     &
+                              NMNHDIM_SERIES_TIME, NMNHDIM_PROFILER_TIME, NMNHDIM_STATION_TIME,                 &
+                              NMNHDIM_ARAKAWA,                                                                  &
+                              NMNHDIM_LASTDIM_NODIACHRO, NMNHDIM_LASTDIM_DIACHRO
+
 use modd_les,           only: nles_k, nspectra_k, xles_temp_mean_start, xles_temp_mean_step, xles_temp_mean_end
 use modd_les_n,         only: nles_times, nspectra_ni, nspectra_nj
 use modd_nsv,           only: nsv
@@ -263,13 +282,9 @@ CHARACTER(LEN=:),ALLOCATABLE :: YPROGRAM
 integer                      :: iavg, iprof, istation
 integer                      :: ispectra_ni, ispectra_nj
 INTEGER                      :: IIU_ll, IJU_ll, IKU
-TYPE(DIMCDF), POINTER        :: tzdimcdf
-TYPE(IOCDF),  POINTER        :: PIOCDF
 
 CALL PRINT_MSG(NVERB_DEBUG,'IO','IO_Knowndims_set_nc4','called for '//TRIM(TPFILE%CNAME))
 
-PIOCDF => TPFILE%TNCDIMS
-
 IF (PRESENT(HPROGRAM_ORIG)) THEN
   YPROGRAM = HPROGRAM_ORIG
 ELSE
@@ -280,75 +295,84 @@ IIU_ll = NIMAX_ll + 2*JPHEXT
 IJU_ll = NJMAX_ll + 2*JPHEXT
 IKU    = NKMAX    + 2*JPVEXT
 
-IF (.NOT. ASSOCIATED(PIOCDF%DIM_NI))      PIOCDF%DIM_NI      => IO_Dimcdf_get_nc4(TPFILE, int( IIU_ll, kind=CDFINT ), 'ni')
-IF (.NOT. ASSOCIATED(PIOCDF%DIM_NJ))      PIOCDF%DIM_NJ      => IO_Dimcdf_get_nc4(TPFILE, int( IJU_ll, kind=CDFINT ), 'nj')
-IF (.NOT. ASSOCIATED(PIOCDF%DIM_NI_U))    PIOCDF%DIM_NI_U    => IO_Dimcdf_get_nc4(TPFILE, int( IIU_ll, kind=CDFINT ), 'ni_u')
-IF (.NOT. ASSOCIATED(PIOCDF%DIM_NJ_U))    PIOCDF%DIM_NJ_U    => IO_Dimcdf_get_nc4(TPFILE, int( IJU_ll, kind=CDFINT ), 'nj_u')
-IF (.NOT. ASSOCIATED(PIOCDF%DIM_NI_V))    PIOCDF%DIM_NI_V    => IO_Dimcdf_get_nc4(TPFILE, int( IIU_ll, kind=CDFINT ), 'ni_v')
-IF (.NOT. ASSOCIATED(PIOCDF%DIM_NJ_V))    PIOCDF%DIM_NJ_V    => IO_Dimcdf_get_nc4(TPFILE, int( IJU_ll, kind=CDFINT ), 'nj_v')
-IF (TRIM(YPROGRAM)/='PGD' .AND. TRIM(YPROGRAM)/='NESPGD' .AND. TRIM(YPROGRAM)/='ZOOMPG' &
-    .AND. .NOT.(TRIM(YPROGRAM)=='REAL' .AND. CSTORAGE_TYPE=='SU') ) THEN !condition to detect PREP_SURFEX
-  IF (.NOT. ASSOCIATED(PIOCDF%DIM_LEVEL))   PIOCDF%DIM_LEVEL   => IO_Dimcdf_get_nc4(TPFILE, int( IKU, kind=CDFINT ), 'level')
-  IF (.NOT. ASSOCIATED(PIOCDF%DIM_LEVEL_W)) PIOCDF%DIM_LEVEL_W => IO_Dimcdf_get_nc4(TPFILE, int( IKU, kind=CDFINT ), 'level_w')
-  if ( tpfile%ctype /= 'MNHDIACHRONIC' ) then
-    IF (.NOT. ASSOCIATED(PIOCDF%DIMTIME)) PIOCDF%DIMTIME => IO_Dimcdf_get_nc4(TPFILE, NF90_UNLIMITED, 'time')
-  end if
-ELSE
-  !PGD and SURFEX files for MesoNH have no vertical levels or time scale
-  !These dimensions are allocated to default values
-  !(they need to be allocated when looking for dimensions of variables)
-  IF (.NOT. ASSOCIATED(PIOCDF%DIM_LEVEL))   ALLOCATE(PIOCDF%DIM_LEVEL)
-  IF (.NOT. ASSOCIATED(PIOCDF%DIM_LEVEL_W)) ALLOCATE(PIOCDF%DIM_LEVEL_W)
-  IF (.NOT. ASSOCIATED(PIOCDF%DIMTIME))     ALLOCATE(PIOCDF%DIMTIME)
-END IF
+if ( .not.Associated( tpfile%tncdims ) ) then
+  call Print_msg( NVERB_ERROR, 'IO', 'IO_Knowndims_set_nc4', 'tncdims not associated for ' // Trim( tpfile%cname ) )
+  Allocate( tpfile%tncdims )
+end if
+
+if ( Allocated( tpfile%tncdims%tdims ) ) then
+  call Print_msg( NVERB_ERROR, 'IO', 'IO_Knowndims_set_nc4', 'tdims already allocated for ' // Trim( tpfile%cname ) )
+  Deallocate( tpfile%tncdims%tdims )
+end if
 
+if ( tpfile%ctype /= 'MNHDIACHRONIC' ) then
+  tpfile%tncdims%nmaxdims = NMNHDIM_LASTDIM_NODIACHRO
+  Allocate( tpfile%tncdims%tdims(NMNHDIM_LASTDIM_NODIACHRO) )
+else
+  tpfile%tncdims%nmaxdims = NMNHDIM_LASTDIM_DIACHRO
+  Allocate( tpfile%tncdims%tdims(NMNHDIM_LASTDIM_DIACHRO) )
+end if
+
+call IO_Add_dim_nc4( tpfile, NMNHDIM_NI,   'ni',   IIU_ll )
+call IO_Add_dim_nc4( tpfile, NMNHDIM_NJ,   'nj',   IJU_ll )
+call IO_Add_dim_nc4( tpfile, NMNHDIM_NI_U, 'ni_u', IIU_ll )
+call IO_Add_dim_nc4( tpfile, NMNHDIM_NJ_U, 'nj_u', IJU_ll )
+call IO_Add_dim_nc4( tpfile, NMNHDIM_NI_V, 'ni_v', IIU_ll )
+call IO_Add_dim_nc4( tpfile, NMNHDIM_NJ_V, 'nj_v', IJU_ll )
+if ( Trim( yprogram ) /= 'PGD' .and. Trim( yprogram ) /= 'NESPGD' .and. Trim( yprogram ) /= 'ZOOMPG' &
+     .and. .not. ( Trim( yprogram ) == 'REAL' .and. cstorage_type == 'SU' ) ) then !condition to detect PREP_SURFEX
+  call IO_Add_dim_nc4( tpfile, NMNHDIM_LEVEL,   'level',   IKU )
+  call IO_Add_dim_nc4( tpfile, NMNHDIM_LEVEL_W, 'level_w', IKU )
+  if ( tpfile%ctype /= 'MNHDIACHRONIC' ) &
+    call IO_Add_dim_nc4( tpfile, NMNHDIM_TIME, 'time', NF90_UNLIMITED )
+end if
 
 if ( tpfile%ctype == 'MNHDIACHRONIC' .or. ( lpack .and. l2d ) ) then
   !Dimension of size 1 used for NMNHDIM_ONE
-  tzdimcdf => IO_Dimcdf_get_nc4( tpfile, 1_CDFINT, hdimname = 'one' )
+  call IO_Add_dim_nc4( tpfile, NMNHDIM_ONE, 'one', 1 )
 end if
 
 !Write dimensions used in diachronic files
 if ( tpfile%ctype == 'MNHDIACHRONIC' ) then
   !Dimension of size 2 used for NMNHDIM_COMPLEX
-  tzdimcdf => IO_Dimcdf_get_nc4( tpfile, 2_CDFINT, hdimname = 'real_imaginary' )
+  call IO_Add_dim_nc4( tpfile, NMNHDIM_COMPLEX, 'real_imaginary', 2 )
 
   !Dimensions for the budgets masks
   if ( cbutype == 'CART' .or. cbutype == 'SKIP' ) then
-    if ( .not. lbu_icp ) tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nbuimax_ll, kind = CDFINT ), hdimname = 'cart_ni'      )
-    if ( .not. lbu_jcp ) tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nbujmax_ll, kind = CDFINT ), hdimname = 'cart_nj'      )
-    if ( .not. lbu_icp ) tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nbuimax_ll, kind = CDFINT ), hdimname = 'cart_ni_u'    )
-    if ( .not. lbu_jcp ) tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nbujmax_ll, kind = CDFINT ), hdimname = 'cart_nj_u'    )
-    if ( .not. lbu_icp ) tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nbuimax_ll, kind = CDFINT ), hdimname = 'cart_ni_v'    )
-    if ( .not. lbu_jcp ) tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nbujmax_ll, kind = CDFINT ), hdimname = 'cart_nj_v'    )
-    if ( .not. lbu_kcp ) tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nbukmax,    kind = CDFINT ), hdimname = 'cart_level'   )
-    if ( .not. lbu_kcp ) tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nbukmax,    kind = CDFINT ), hdimname = 'cart_level_w' )
+    if ( .not. lbu_icp ) call IO_Add_dim_nc4( tpfile, NMNHDIM_BUDGET_CART_NI,      'cart_ni',      nbuimax_ll )
+    if ( .not. lbu_jcp ) call IO_Add_dim_nc4( tpfile, NMNHDIM_BUDGET_CART_NJ,      'cart_nj',      nbujmax_ll )
+    if ( .not. lbu_icp ) call IO_Add_dim_nc4( tpfile, NMNHDIM_BUDGET_CART_NI_U,    'cart_ni_u',    nbuimax_ll )
+    if ( .not. lbu_jcp ) call IO_Add_dim_nc4( tpfile, NMNHDIM_BUDGET_CART_NJ_U,    'cart_nj_u',    nbujmax_ll )
+    if ( .not. lbu_icp ) call IO_Add_dim_nc4( tpfile, NMNHDIM_BUDGET_CART_NI_V,    'cart_ni_v',    nbuimax_ll )
+    if ( .not. lbu_jcp ) call IO_Add_dim_nc4( tpfile, NMNHDIM_BUDGET_CART_NJ_V,    'cart_nj_v',    nbujmax_ll )
+    if ( .not. lbu_kcp ) call IO_Add_dim_nc4( tpfile, NMNHDIM_BUDGET_CART_LEVEL,   'cart_level',   nbukmax    )
+    if ( .not. lbu_kcp ) call IO_Add_dim_nc4( tpfile, NMNHDIM_BUDGET_CART_LEVEL_W, 'cart_level_w', nbukmax    )
   else if ( cbutype == 'MASK' ) then
-    if ( nbukmax > 0 ) tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nbukmax, kind = CDFINT ), hdimname = 'mask_level'   )
-    if ( nbukmax > 0 ) tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nbukmax, kind = CDFINT ), hdimname = 'mask_level_w' )
-    if ( nbuwrnb > 0 ) tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nbuwrnb, kind = CDFINT ), hdimname = 'time_mask'    )
-    if ( nbumask > 0 ) tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nbumask, kind = CDFINT ), hdimname = 'nbumask'      )
+    if ( nbukmax > 0 ) call IO_Add_dim_nc4( tpfile, NMNHDIM_BUDGET_MASK_LEVEL,   'mask_level',   nbukmax )
+    if ( nbukmax > 0 ) call IO_Add_dim_nc4( tpfile, NMNHDIM_BUDGET_MASK_LEVEL_W, 'mask_level_w', nbukmax )
+    if ( nbuwrnb > 0 ) call IO_Add_dim_nc4( tpfile, NMNHDIM_BUDGET_MASK_TIME,    'time_mask',    nbuwrnb )
+    if ( nbumask > 0 ) call IO_Add_dim_nc4( tpfile, NMNHDIM_BUDGET_MASK_NBUMASK, 'nbumask',      nbumask )
   end if
 
   !Dimension for the number of LES budget time samplings
-  if ( nles_times > 0 ) tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nles_times, kind = CDFINT ), hdimname = 'time_les' )
+  if ( nles_times > 0 ) call IO_Add_dim_nc4( tpfile, NMNHDIM_BUDGET_LES_TIME, 'time_les', nles_times )
 
   !Dimension for the number of LES budget time averages
   iavg = int( xles_temp_mean_end - 1.e-10 - xles_temp_mean_start ) / xles_temp_mean_step + 1
   !Condition also on nles_times to not create this dimension when not used (no time average if nles_times=0)
   if ( nles_times > 0 .and. iavg > 0 ) &
-    tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( iavg, kind = CDFINT ), hdimname = 'time_les_avg' )
+    call IO_Add_dim_nc4( tpfile, NMNHDIM_BUDGET_LES_AVG_TIME, 'time_les_avg', iavg )
 
   !Dimension for the number of vertical levels for local LES budgets
-  if ( nles_k > 0 ) tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nles_k, kind = CDFINT ), hdimname = 'level_les' )
+  if ( nles_k > 0 ) call IO_Add_dim_nc4( tpfile, NMNHDIM_BUDGET_LES_LEVEL, 'level_les', nles_k )
 
   !Dimension for the number of scalar variables
-  if ( nsv > 0 ) tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nsv, kind = CDFINT ), hdimname = 'nsv' )
+  if ( nsv > 0 ) call IO_Add_dim_nc4( tpfile, NMNHDIM_BUDGET_LES_SV, 'nsv', nsv )
 
   !Dimensions for the number of horizontal wavelengths for non-local LES budgets (2 points correlations)
-  if ( nspectra_ni > 0 ) tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nspectra_ni, kind = CDFINT ), hdimname = 'nspectra_2pts_ni' )
+  if ( nspectra_ni > 0 ) call IO_Add_dim_nc4( tpfile, NMNHDIM_SPECTRA_2PTS_NI, 'nspectra_2pts_ni', nspectra_ni )
   if ( nspectra_nj > 0 .and. .not. l2d ) &
-    tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nspectra_nj, kind = CDFINT ), hdimname = 'nspectra_2pts_nj' )
+    call IO_Add_dim_nc4( tpfile, NMNHDIM_SPECTRA_2PTS_NJ, 'nspectra_2pts_nj', nspectra_nj )
 
   !Dimensions for the number of horizontal wavelengths for LES spectra budgets
   if ( nspectra_ni > 0 ) then
@@ -358,7 +382,7 @@ if ( tpfile%ctype == 'MNHDIACHRONIC' ) then
       ispectra_ni =  nspectra_ni - 1
     end if
     if ( ispectra_ni > 0 ) &
-      tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( ispectra_ni, kind = CDFINT ), hdimname = 'nspectra_spec_ni' )
+      call IO_Add_dim_nc4( tpfile, NMNHDIM_SPECTRA_SPEC_NI, 'nspectra_spec_ni', ispectra_ni )
   end if
 
   if ( nspectra_nj > 0 .and. .not. l2d ) then
@@ -368,102 +392,77 @@ if ( tpfile%ctype == 'MNHDIACHRONIC' ) then
       ispectra_nj =  nspectra_nj - 1
     end if
     if ( ispectra_nj > 0 ) &
-      tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( ispectra_nj, kind = CDFINT ), hdimname = 'nspectra_spec_nj' )
+      call IO_Add_dim_nc4( tpfile, NMNHDIM_SPECTRA_SPEC_NJ, 'nspectra_spec_nj', ispectra_nj )
   end if
 
   !Dimension for the number of vertical levels for non-local LES budgets
-  if ( nspectra_k > 0 ) tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nspectra_k, kind = CDFINT ), hdimname = 'nspectra_level' )
+  if ( nspectra_k > 0 ) call IO_Add_dim_nc4( tpfile, NMNHDIM_SPECTRA_LEVEL, 'nspectra_level', nspectra_k )
 
   !Dimension for the number of profiler times
   if ( numbprofiler > 0 ) then
     iprof = Int ( ( xseglen - xtstep ) / tprofiler%step ) + 1
-    tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( iprof, kind = CDFINT ), hdimname = 'time_profiler' )
+    call IO_Add_dim_nc4( tpfile, NMNHDIM_PROFILER_TIME, 'time_profiler', iprof )
   end if
 
   !Dimension for the number of station times
   if ( numbstat > 0 ) then
     istation = Int ( ( xseglen - xtstep ) / tstation%step ) + 1
-    tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( istation, kind = CDFINT ), hdimname = 'time_station' )
+    call IO_Add_dim_nc4( tpfile, NMNHDIM_STATION_TIME, 'time_station', istation )
   end if
 
   !Dimension for the number of series times
   if ( lseries .and. nsnbstept > 0 ) then
-    tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nkmax, kind = CDFINT ), hdimname = 'series_level'   )
-    tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nkmax, kind = CDFINT ), hdimname = 'series_level_w' )
-    tzdimcdf => IO_Dimcdf_get_nc4( tpfile, int( nsnbstept, kind = CDFINT ), hdimname = 'time_series' )
+    call IO_Add_dim_nc4( tpfile, NMNHDIM_SERIES_LEVEL,   'series_level',   nkmax     )
+    call IO_Add_dim_nc4( tpfile, NMNHDIM_SERIES_LEVEL_W, 'series_level_w', nkmax     )
+    call IO_Add_dim_nc4( tpfile, NMNHDIM_SERIES_TIME,    'time_series',    nsnbstept )
   end if
 end if
 
-!Store X,Y,Z coordinates for the Arakawa points
-!0 2nd-dimension is to treat NGRID=0 case without crash
-IF (.NOT.ALLOCATED(TPFILE%TNCCOORDS)) ALLOCATE(TPFILE%TNCCOORDS(3,0:8))
-!Dummy point
-TPFILE%TNCCOORDS(1,0)%TDIM => TDIM_DUMMY
-TPFILE%TNCCOORDS(2,0)%TDIM => TDIM_DUMMY
-TPFILE%TNCCOORDS(3,0)%TDIM => TDIM_DUMMY
-! Mass point
-TPFILE%TNCCOORDS(1,1)%TDIM => PIOCDF%DIM_NI
-TPFILE%TNCCOORDS(2,1)%TDIM => PIOCDF%DIM_NJ
-TPFILE%TNCCOORDS(3,1)%TDIM => PIOCDF%DIM_LEVEL
-! u point
-TPFILE%TNCCOORDS(1,2)%TDIM => PIOCDF%DIM_NI_U
-TPFILE%TNCCOORDS(2,2)%TDIM => PIOCDF%DIM_NJ_U
-TPFILE%TNCCOORDS(3,2)%TDIM => PIOCDF%DIM_LEVEL
-! v point
-TPFILE%TNCCOORDS(1,3)%TDIM => PIOCDF%DIM_NI_V
-TPFILE%TNCCOORDS(2,3)%TDIM => PIOCDF%DIM_NJ_V
-TPFILE%TNCCOORDS(3,3)%TDIM => PIOCDF%DIM_LEVEL
-! w point
-TPFILE%TNCCOORDS(1,4)%TDIM => PIOCDF%DIM_NI
-TPFILE%TNCCOORDS(2,4)%TDIM => PIOCDF%DIM_NJ
-TPFILE%TNCCOORDS(3,4)%TDIM => PIOCDF%DIM_LEVEL_W
-! xi vorticity point (=f point =uv point)
-TPFILE%TNCCOORDS(1,5)%TDIM => PIOCDF%DIM_NI_U
-TPFILE%TNCCOORDS(2,5)%TDIM => PIOCDF%DIM_NJ_V
-TPFILE%TNCCOORDS(3,5)%TDIM => PIOCDF%DIM_LEVEL
-! eta vorticity point (=uw point)
-TPFILE%TNCCOORDS(1,6)%TDIM => PIOCDF%DIM_NI_U
-TPFILE%TNCCOORDS(2,6)%TDIM => PIOCDF%DIM_NJ_U
-TPFILE%TNCCOORDS(3,6)%TDIM => PIOCDF%DIM_LEVEL_W
-! zeta vorticity point (=vw point)
-TPFILE%TNCCOORDS(1,7)%TDIM => PIOCDF%DIM_NI_V
-TPFILE%TNCCOORDS(2,7)%TDIM => PIOCDF%DIM_NJ_V
-TPFILE%TNCCOORDS(3,7)%TDIM => PIOCDF%DIM_LEVEL_W
-! fw point (=uvw point)
-TPFILE%TNCCOORDS(1,8)%TDIM => PIOCDF%DIM_NI_U
-TPFILE%TNCCOORDS(2,8)%TDIM => PIOCDF%DIM_NJ_V
-TPFILE%TNCCOORDS(3,8)%TDIM => PIOCDF%DIM_LEVEL_W
+END SUBROUTINE IO_Knowndims_set_nc4
 
 
-END SUBROUTINE IO_Knowndims_set_nc4
+subroutine IO_Add_dim_nc4( tpfile, kidx, hdimname, klen )
+use NETCDF, only: NF90_DEF_DIM, NF90_NOERR
+
+type(tfiledata),  intent(inout) :: tpfile
+integer,          intent(in)    :: kidx     !Position of the dimension in the list
+character(len=*), intent(in)    :: hdimname !Name of the dimension
+integer,          intent(in)    :: klen     !Length of the dimension
 
+integer(kind=CDFINT)          :: istatus
 
-SUBROUTINE IO_Iocdf_dealloc_nc4(PIOCDF)
-TYPE(IOCDF),  POINTER :: PIOCDF
 
-INTEGER(KIND=CDFINT) :: IRESP
+if ( .not.Associated( tpfile%tncdims ) ) &
+  call Print_msg( NVERB_FATAL, 'IO', 'IO_Add_dim_nc4', 'tncdims not associated for ' // Trim( tpfile%cname ) )
 
-CALL PRINT_MSG(NVERB_DEBUG,'IO','IO_Iocdf_dealloc_nc4','called')
+if ( kidx < 1 .or. kidx > Size( tpfile%tncdims%tdims ) )                                                      &
+  call Print_msg( NVERB_FATAL, 'IO', 'IO_Add_dim_nc4', 'index out of range for dimension ' // Trim( hdimname ) // &
+                  ' of file ' //Trim( tpfile%cname ) )
+
+if ( tpfile%tncdims%tdims(kidx)%nlen /= -1 .or. tpfile%tncdims%tdims(kidx)%nid /= -1 ) &
+  call Print_msg( NVERB_WARNING, 'IO', 'IO_Add_dim_nc4', 'dimension ' // Trim( hdimname ) //   &
+                  ' already defined for file ' //Trim( tpfile%cname ) )
+
+tpfile%tncdims%tdims(kidx)%cname = hdimname
+tpfile%tncdims%tdims(kidx)%nlen  = Int( klen, kind = CDFINT )
+
+istatus = NF90_DEF_DIM( tpfile%nncid, Trim( hdimname ), Int( klen, kind = CDFINT ), tpfile%tncdims%tdims(kidx)%nid )
+if ( istatus /= NF90_NOERR ) &
+  call IO_Err_handle_nc4( istatus, 'IO_Add_dim_nc4', 'NF90_DEF_DIM', Trim( hdimname ) )
+
+end subroutine IO_Add_dim_nc4
 
-! Clean DIMLIST and DIMSTR
-CALL CLEANLIST(PIOCDF%DIMLIST)
-CALL CLEANLIST(PIOCDF%DIMSTR)
-! Then free iocdf
-DEALLOCATE(PIOCDF)
 
-CONTAINS
+SUBROUTINE IO_Iocdf_dealloc_nc4(tpdimsnc)
+TYPE(tdimsnc),  POINTER :: tpdimsnc
 
-SUBROUTINE CLEANLIST(PLIST)
-TYPE(DIMCDF), POINTER :: PLIST,TZDIMCUR, TZDIMNEXT
+CALL PRINT_MSG(NVERB_DEBUG,'IO','IO_Iocdf_dealloc_nc4','called')
 
-TZDIMCUR  => PLIST
-DO WHILE(ASSOCIATED(TZDIMCUR))
-   TZDIMNEXT => TZDIMCUR%NEXT
-   DEALLOCATE(TZDIMCUR)
-   TZDIMCUR => TZDIMNEXT
-END DO
+if ( Allocated( tpdimsnc%tdims     ) ) deallocate( tpdimsnc%tdims     )
+if ( Allocated( tpdimsnc%tdims_str ) ) deallocate( tpdimsnc%tdims_str )
 
-END SUBROUTINE CLEANLIST
+deallocate( tpdimsnc )
+tpdimsnc => Null()
 
 END SUBROUTINE IO_Iocdf_dealloc_nc4
 
@@ -487,8 +486,8 @@ use modd_field,  only: NMNHDIM_UNKNOWN, NMNHDIM_ONE, NMNHDIM_COMPLEX,
                        NMNHDIM_SPECTRA_LEVEL,                                                           &
                        NMNHDIM_SERIES_LEVEL,        NMNHDIM_SERIES_LEVEL_W, NMNHDIM_SERIES_TIME,        &
                        NMNHDIM_FLYER_TIME,          NMNHDIM_PROFILER_TIME,                              &
-                       NMNHDIM_STATION_TIME,                                                            &
-                       NMNHDIM_NOTLISTED, NMNHDIM_UNUSED
+                       NMNHDIM_STATION_TIME,        NMNHDIM_LASTDIM_DIACHRO,                            &
+                       NMNHDIM_NOTLISTED, NMNHDIM_UNUSED, NMNHDIM_ARAKAWA
 
 TYPE(TFILEDATA),                              INTENT(IN)  :: TPFILE
 TYPE(TFIELDDATA),                             INTENT(IN)  :: TPFIELD
@@ -499,12 +498,11 @@ CHARACTER(LEN=32)             :: YINT
 CHARACTER(LEN=2)              :: YDIR
 character(len=:), allocatable :: ydimname
 INTEGER                       :: IGRID
+integer                       :: iidx
 integer                       :: iresp
 INTEGER                       :: JI
 integer(kind=CDFINT)          :: ilen
 integer(kind=CDFINT)          :: istatus
-type(dimcdf)                  :: tzdim
-TYPE(DIMCDF), POINTER         :: PTDIM
 !
 CALL PRINT_MSG(NVERB_DEBUG,'IO','IO_Vdims_fill_nc4','called for '//TRIM(TPFIELD%CMNHNAME))
 !
@@ -522,7 +520,7 @@ END IF
 IF (TPFIELD%LTIMEDEP) THEN
   !Add time dimension
   ALLOCATE(KVDIMS(TPFIELD%NDIMS+1))
-  KVDIMS(TPFIELD%NDIMS+1) = TPFILE%TNCDIMS%DIMTIME%ID
+  KVDIMS(TPFIELD%NDIMS+1) = tpfile%tncdims%tdims(NMNHDIM_TIME)%nid
 ELSE
   ALLOCATE(KVDIMS(TPFIELD%NDIMS))
 END IF
@@ -543,7 +541,8 @@ if ( Any( tpfield%ndimlist(:) /= NMNHDIM_UNKNOWN ) ) then
       call Print_msg( NVERB_FATAL, 'IO', 'IO_Vdims_fill_nc4', 'ndimlist partially filled for field ' // Trim( tpfield%cmnhname ) )
 
     if ( tpfield%ndimlist(ji) == NMNHDIM_NOTLISTED ) then
-      ptdim => IO_Dimcdf_get_nc4( tpfile, kshape(ji) ); kvdims(ji) = ptdim%id
+      call IO_Dim_find_create_nc4( tpfile, kshape(ji), iidx )
+      kvdims(ji) = tpfile%tncdims%tdims(iidx)%nid
       cycle
     end if
 
@@ -563,290 +562,203 @@ if ( Any( tpfield%ndimlist(:) /= NMNHDIM_UNKNOWN ) ) then
       cycle
     end if
 
-    select case ( tpfield%ndimlist(ji) )
-      case ( NMNHDIM_ONE )
-        ydimname = 'one'
-
-      case ( NMNHDIM_COMPLEX )
-        ydimname = 'real_imaginary'
-
-      case ( NMNHDIM_NI )
-        ydimname = 'ni'
-
-      case ( NMNHDIM_NJ )
-        ydimname = 'nj'
-
-      case ( NMNHDIM_NI_U )
-        ydimname = 'ni_u'
-
-      case ( NMNHDIM_NJ_U )
-        ydimname = 'nj_u'
-
-      case ( NMNHDIM_NI_V )
-        ydimname = 'ni_v'
-
-      case ( NMNHDIM_NJ_V )
-        ydimname = 'nj_v'
-
-      case ( NMNHDIM_LEVEL )
-        ydimname = 'level'
-
-      case ( NMNHDIM_LEVEL_W )
-        ydimname = 'level_w'
-
-      case ( NMNHDIM_TIME )
-        ydimname = 'time'
-
-      case ( NMNHDIM_BUDGET_CART_NI )
-        ydimname = 'cart_ni'
-
-      case ( NMNHDIM_BUDGET_CART_NJ )
-        ydimname = 'cart_nj'
-
-      case ( NMNHDIM_BUDGET_CART_NI_U )
-        ydimname = 'cart_ni_u'
-
-      case ( NMNHDIM_BUDGET_CART_NJ_U )
-        ydimname = 'cart_nj_u'
-
-      case ( NMNHDIM_BUDGET_CART_NI_V )
-        ydimname = 'cart_ni_v'
-
-      case ( NMNHDIM_BUDGET_CART_NJ_V )
-        ydimname = 'cart_nj_v'
-
-      case ( NMNHDIM_BUDGET_CART_LEVEL )
-        ydimname = 'cart_level'
-
-      case ( NMNHDIM_BUDGET_CART_LEVEL_W )
-        ydimname = 'cart_level_w'
-
-      case ( NMNHDIM_BUDGET_MASK_LEVEL )
-        ydimname = 'mask_level'
-
-      case ( NMNHDIM_BUDGET_MASK_LEVEL_W )
-        ydimname = 'mask_level_w'
-
-      case ( NMNHDIM_BUDGET_MASK_TIME )
-        ydimname = 'time_mask'
-
-      case ( NMNHDIM_BUDGET_MASK_NBUMASK )
-        ydimname = 'nbumask'
+    if ( tpfield%ndimlist(ji) > 0 .and. tpfield%ndimlist(ji) <= tpfile%tncdims%nmaxdims ) then
+      kvdims(ji) = tpfile%tncdims%tdims(tpfield%ndimlist(ji))%nid
 
-      case ( NMNHDIM_BUDGET_LES_TIME )
-        ydimname = 'time_les'
-
-      case ( NMNHDIM_BUDGET_LES_AVG_TIME )
-        ydimname = 'time_les_avg'
-
-      case ( NMNHDIM_BUDGET_LES_LEVEL )
-        ydimname = 'level_les'
-
-      case ( NMNHDIM_BUDGET_LES_SV )
-        ydimname = 'nsv'
-
-      case ( NMNHDIM_SPECTRA_2PTS_NI )
-        ydimname = 'nspectra_2pts_ni'
-
-      case ( NMNHDIM_SPECTRA_2PTS_NJ )
-        ydimname = 'nspectra_2pts_nj'
-
-      case ( NMNHDIM_SPECTRA_SPEC_NI )
-        ydimname = 'nspectra_spec_ni'
-
-      case ( NMNHDIM_SPECTRA_SPEC_NJ )
-        ydimname = 'nspectra_spec_nj'
-
-      case ( NMNHDIM_SPECTRA_LEVEL )
-        ydimname = 'nspectra_level'
-
-      case ( NMNHDIM_PROFILER_TIME )
-        ydimname = 'time_profiler'
-
-      case ( NMNHDIM_STATION_TIME )
-        ydimname = 'time_station'
-
-      case ( NMNHDIM_SERIES_LEVEL )
-        ydimname = 'series_level'
-
-      case ( NMNHDIM_SERIES_LEVEL_W )
-        ydimname = 'series_level_w'
-
-      case ( NMNHDIM_SERIES_TIME )
-        ydimname = 'time_series'
-
-      case default
-        call Print_msg( NVERB_FATAL, 'IO', 'IO_Vdims_fill_nc4', &
-                        'ndimlist case not yet implemented for field ' // Trim( tpfield%cmnhname ) )
-    end select
-
-    call IO_Dim_find_byname_nc4( tpfile, ydimname, tzdim, iresp )
-    kvdims(ji) = tzdim%id
-
-    ! Check if dimension sizes are consistent with the declared dimensions ( skip ji>size(kshape), timedep dimension)
-    if ( ji <= Size( kshape ) ) then
-      if ( kshape(ji) /= tzdim%len ) then
-        call Print_msg( NVERB_FATAL, 'IO', 'IO_Vdims_fill_nc4', &
-                                     'wrong size for dimension '// Trim( tzdim%name ) // ' of field ' // Trim( tpfield%cmnhname ) )
+      ! Check if dimension sizes are consistent with the declared dimensions ( skip ji>size(kshape), timedep dimension)
+      if ( ji <= Size( kshape ) ) then
+        if ( kshape(ji) /= tpfile%tncdims%tdims(tpfield%ndimlist(ji))%nlen )              &
+          call Print_msg( NVERB_FATAL, 'IO', 'IO_Vdims_fill_nc4', 'wrong size for dimension ' &
+                          // Trim( tpfile%tncdims%tdims(tpfield%ndimlist(ji))%cname ) //  &
+                                       ' of field ' // Trim( tpfield%cmnhname ) )
       end if
+    else
+      call Print_msg( NVERB_FATAL, 'IO', 'IO_Vdims_fill_nc4', 'wrong dimension' )
     end if
-
   end do
 else !ndimlist not provided
-  DO JI=1,SIZE(KSHAPE)
-    IF (JI == 1) THEN
-      IF ( (YDIR == 'XX' .OR. YDIR == 'XY') .AND. KSHAPE(1)==TPFILE%TNCCOORDS(1,IGRID)%TDIM%LEN) THEN
-        KVDIMS(1) = TPFILE%TNCCOORDS(1,IGRID)%TDIM%ID
-      ELSE IF ( YDIR == 'YY'                .AND. KSHAPE(1)==TPFILE%TNCCOORDS(2,IGRID)%TDIM%LEN) THEN
-        KVDIMS(1) = TPFILE%TNCCOORDS(2,IGRID)%TDIM%ID
-      ELSE IF ( YDIR == 'ZZ'                .AND. KSHAPE(1)==TPFILE%TNCCOORDS(3,IGRID)%TDIM%LEN) THEN
-        KVDIMS(1) = TPFILE%TNCCOORDS(3,IGRID)%TDIM%ID
-      ELSE
-        PTDIM => IO_Dimcdf_get_nc4(TPFILE, KSHAPE(1)); KVDIMS(1) = PTDIM%ID
-      END IF
-    ELSE IF (JI == 2) THEN
-      IF ( YDIR == 'XY' .AND. KSHAPE(2)==TPFILE%TNCCOORDS(2,IGRID)%TDIM%LEN) THEN
-        KVDIMS(2) = TPFILE%TNCCOORDS(2,IGRID)%TDIM%ID
-      ELSE
-        PTDIM => IO_Dimcdf_get_nc4(TPFILE, KSHAPE(2)); KVDIMS(2) = PTDIM%ID
-      END IF
-    ELSE IF (JI == 3) THEN
-      IF ( YDIR == 'XY' .AND. KSHAPE(3)==TPFILE%TNCCOORDS(3,IGRID)%TDIM%LEN) THEN
-        KVDIMS(3) = TPFILE%TNCCOORDS(3,IGRID)%TDIM%ID
-      ELSE
-        PTDIM => IO_Dimcdf_get_nc4(TPFILE, KSHAPE(3)); KVDIMS(3) = PTDIM%ID
-      END IF
-    ELSE
-        PTDIM => IO_Dimcdf_get_nc4(TPFILE, KSHAPE(JI)); KVDIMS(JI) = PTDIM%ID
-    END IF
+  do ji = 1, Size( kshape )
+    if ( igrid > 0 .and. igrid < 9 ) then
+      if ( ji == 1 ) then
+        if ( ( ydir == 'XX' .or. ydir == 'XY') .and. kshape(1) == tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,1) )%nlen ) then
+          kvdims(1) = tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,1) )%nid
+        else if ( ydir == 'YY'                .and. kshape(1) == tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,2) )%nlen ) then
+          kvdims(1) = tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,2) )%nid
+        else if ( ydir == 'ZZ'                .and. kshape(1) == tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,3) )%nlen ) then
+          kvdims(1) = tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,3) )%nid
+        else
+          call IO_Dim_find_create_nc4( tpfile, kshape(1), iidx )
+          kvdims(1) = tpfile%tncdims%tdims(iidx)%nid
+        end if
+      else if ( ji == 2 ) then
+        if ( ydir == 'XY' .and. kshape(2) == tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,2) )%nlen ) then
+          kvdims(2) = tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,2) )%nid
+        else
+          call IO_Dim_find_create_nc4( tpfile, kshape(2), iidx )
+          kvdims(2) = tpfile%tncdims%tdims(iidx)%nid
+        end if
+      else if ( ji == 3 ) then
+        if ( ydir == 'XY' .and. kshape(3) == tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,3) )%nlen ) then
+          kvdims(3) = tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(igrid,3) )%nid
+        else
+          call IO_Dim_find_create_nc4( tpfile, kshape(3), iidx )
+          kvdims(3) = tpfile%tncdims%tdims(iidx)%nid
+        end if
+      else
+        call IO_Dim_find_create_nc4( tpfile, kshape(ji), iidx )
+        kvdims(ji) = tpfile%tncdims%tdims(iidx)%nid
+      end if
+    else
+      call IO_Dim_find_create_nc4( tpfile, kshape(ji), iidx )
+      kvdims(ji) = tpfile%tncdims%tdims(iidx)%nid
+    end if
   END DO
 end if
 
 END SUBROUTINE IO_Vdims_fill_nc4
 
 
-FUNCTION IO_Dimcdf_get_nc4(TPFILE, KLEN, HDIMNAME)
-TYPE(TFILEDATA),            INTENT(IN) :: TPFILE
-INTEGER(KIND=CDFINT),       INTENT(IN) :: KLEN
-CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: HDIMNAME ! When provided don't search but
-                                                   ! simply create with name HDIMNAME
-TYPE(DIMCDF), POINTER   :: IO_Dimcdf_get_nc4
+subroutine IO_Dim_find_create_nc4( tpfile, klen, kidx, hdimname)
+use modd_netcdf, only: tdimnc
 
-TYPE(DIMCDF), POINTER :: TMP
-INTEGER               :: COUNT
-CHARACTER(LEN=16)     :: YSUFFIX
-CHARACTER(LEN=20)     :: YDIMNAME
-INTEGER(KIND=CDFINT)  :: STATUS
-LOGICAL               :: GCHKLEN !Check if KLEN is valid
-TYPE(IOCDF), POINTER  :: PIOCDF
+type(tfiledata),            intent(in) :: tpfile
+integer,       intent(in) :: klen
+CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: HDIMNAME
+integer, intent(out) :: kidx !Position of the dimension in the dimension array
 
-CALL PRINT_MSG(NVERB_DEBUG,'IO','IO_Dimcdf_get_nc4','called')
+character(len=16)     :: ysuffix
+integer :: inewsize
+integer :: ji
+integer(kind=CDFINT)  :: istatus
+type(tdimnc), dimension(:), allocatable :: tzncdims
 
-PIOCDF => TPFILE%TNCDIMS
 
-GCHKLEN = .TRUE.
-!Do not check KLEN if 'time' (because NF90_UNLIMITED = 0)
-IF (PRESENT(HDIMNAME)) THEN
-  IF (TRIM(HDIMNAME)=='time') THEN
-    GCHKLEN = .FALSE.
-  END IF
-END IF
+kidx = -1
 
-WRITE(YSUFFIX,'(I0)') KLEN
+do ji = 1, Size( tpfile%tncdims%tdims )
+  if ( tpfile%tncdims%tdims(ji)%nlen == klen ) then
+    if ( Present( hdimname ) ) then
+      if ( hdimname == Trim( tpfile%tncdims%tdims(ji)%cname ) ) then
+        kidx = ji
+        exit
+      end if
+    else
+      kidx = ji
+      exit
+    end if
+  end if
+end do
+
+if ( kidx == - 1 ) then
+  !Check if already exist with the provided name (if so => error)
+  if ( Present( hdimname ) ) then
+    do ji = 1, Size( tpfile%tncdims%tdims )
+      if ( hdimname == Trim( tpfile%tncdims%tdims(ji)%cname ) )                             &
+        call Print_msg( NVERB_ERROR, 'IO', 'IO_Dim_find_create_nc4', 'dimension '          &
+                        // Trim( hdimname ) // ' already exist but with a different size' )
+    end do
+  end if
 
-IF (GCHKLEN .AND. KLEN < 1) THEN
-  CALL PRINT_MSG(NVERB_FATAL,'IO','IO_Dimcdf_get_nc4','KLEN='//TRIM(YSUFFIX))
-END IF
+  !Create new dimension
+  inewsize = Size( tpfile%tncdims%tdims ) + 1
+  allocate( tzncdims(inewsize) )
+  tzncdims(1 : inewsize - 1) = tpfile%tncdims%tdims(:)
 
-IF (PRESENT(HDIMNAME)) THEN
-   NULLIFY(TMP)
-   YDIMNAME = TRIM(HDIMNAME)
-ELSE
-   ! Search dimension with KLEN length
-   COUNT = 1
-   TMP  => PIOCDF%DIMLIST
-   DO WHILE(ASSOCIATED(TMP))
-      IF (TMP%LEN == KLEN .AND. TMP%NAME(1:4) /= 'char') EXIT
-      TMP=>TMP%NEXT
-      COUNT = COUNT+1
-   END DO
-   YDIMNAME = 'size'//TRIM(YSUFFIX)
-END IF
+  if ( Present( hdimname ) ) then
+    tzncdims(inewsize)%cname = Trim( hdimname )
+  else
+    Write( ysuffix, '( i0 )' ) klen
+    tzncdims(inewsize)%cname = 'size' // Trim( ysuffix )
+  end if
+  tzncdims(inewsize)%nlen = Int( klen, kind = CDFINT )
 
-IF (.NOT. ASSOCIATED(TMP)) THEN
-   ! Not found then define new dimension
-   ALLOCATE(TMP)
-   TMP%NAME = YDIMNAME
-   TMP%LEN = KLEN
-   STATUS = NF90_DEF_DIM(TPFILE%NNCID, TMP%NAME, KLEN, TMP%ID)
-   IF (STATUS /= NF90_NOERR) CALL IO_Err_handle_nc4(status,'IO_Dimcdf_get_nc4','NF90_DEF_DIM',trim(TMP%NAME))
-   NULLIFY(TMP%NEXT)
-   TMP%NEXT       => PIOCDF%DIMLIST
-   PIOCDF%DIMLIST => TMP
-CALL PRINT_MSG(NVERB_DEBUG,'IO','IO_Dimcdf_get_nc4','new dimension: '//TRIM(TMP%NAME))
-END IF
+  istatus = NF90_DEF_DIM( tpfile%nncid, tzncdims(inewsize)%cname, Int( klen, kind = CDFINT ), tzncdims(inewsize)%nid )
+  if ( istatus /= NF90_NOERR ) &
+    call IO_Err_handle_nc4( istatus, 'IO_Dim_find_create_nc4', 'NF90_DEF_DIM', Trim( tzncdims(inewsize)%cname) )
 
-IO_Dimcdf_get_nc4 => TMP
+  call Move_alloc( from = tzncdims, to = tpfile%tncdims%tdims )
+  tpfile%tncdims%nmaxdims = inewsize
 
-END FUNCTION IO_Dimcdf_get_nc4
+  kidx = inewsize
+  call Print_msg( NVERB_DEBUG, 'IO', 'IO_Dim_find_create_nc4','new dimension: ' // Trim( tpfile%tncdims%tdims(inewsize)%cname ) )
+end if
+
+end subroutine IO_Dim_find_create_nc4
 
 
 FUNCTION IO_Strdimid_get_nc4(TPFILE,KLEN)
+use modd_netcdf, only: tdimnc
+
 TYPE(TFILEDATA),      INTENT(IN) :: TPFILE
 INTEGER(KIND=CDFINT), INTENT(IN) :: KLEN
 INTEGER(KIND=CDFINT)             :: IO_Strdimid_get_nc4
 
-TYPE(DIMCDF), POINTER :: TMP
-TYPE(IOCDF),  POINTER :: TZIOCDF
-CHARACTER(LEN=16)     :: YSUFFIX
-INTEGER(KIND=CDFINT)  :: STATUS
+CHARACTER(LEN=16)                       :: YSUFFIX
+integer                                 :: idx
+integer                                 :: inewsize
+integer                                 :: ji
+integer(kind=CDFINT)                    :: istatus
+type(tdimnc), dimension(:), allocatable :: tzncdims
 
 CALL PRINT_MSG(NVERB_DEBUG,'IO','IO_Strdimid_get_nc4','called')
 
-WRITE(YSUFFIX,'(I0)') KLEN
-
 IF (KLEN < 1) THEN
   CALL PRINT_MSG(NVERB_FATAL,'IO','IO_Strdimid_get_nc4','KLEN='//TRIM(YSUFFIX))
 END IF
 
 ! Search string dimension with KLEN length
-TMP  => TPFILE%TNCDIMS%DIMSTR
-DO WHILE(ASSOCIATED(TMP))
-   IF (TMP%LEN == KLEN) EXIT
-   TMP=>TMP%NEXT
-END DO
-
-IF (.NOT. ASSOCIATED(TMP)) THEN
-   ! Not found then define new dimension
-   ALLOCATE(TMP)
-   TMP%NAME = 'char'//TRIM(YSUFFIX)
-   TMP%LEN = KLEN
-   STATUS = NF90_DEF_DIM(TPFILE%NNCID, TMP%NAME, KLEN, TMP%ID)
-   IF (STATUS /= NF90_NOERR) CALL IO_Err_handle_nc4(status,'IO_Strdimid_get_nc4','NF90_DEF_DIM',trim(TMP%NAME))
-   NULLIFY(TMP%NEXT)
-   TMP%NEXT      => TPFILE%TNCDIMS%DIMSTR
-   TZIOCDF => TPFILE%TNCDIMS
-   TZIOCDF%DIMSTR => TMP
-END IF
+idx = -1
+if ( Allocated( tpfile%tncdims%tdims_str ) ) then
+  do ji = 1, Size( tpfile%tncdims%tdims_str )
+    if ( tpfile%tncdims%tdims_str(ji)%nlen == klen ) then
+      idx = ji
+      exit
+    end if
+  end do
+end if
+
+if ( idx == -1 ) then
+  !Create new dimension
+  inewsize = tpfile%tncdims%nmaxdims_str + 1
+  allocate( tzncdims(inewsize) )
+  if ( Allocated( tpfile%tncdims%tdims_str ) ) &
+    tzncdims(1 : inewsize - 1) = tpfile%tncdims%tdims_str(:)
+
+  Write( ysuffix, '( i0 )' ) klen
+  tzncdims(inewsize)%cname = 'char' // Trim( ysuffix )
+  tzncdims(inewsize)%nlen = Int( klen, kind = CDFINT )
+
+  istatus = NF90_DEF_DIM( tpfile%nncid, tzncdims(inewsize)%cname, Int( klen, kind = CDFINT ), tzncdims(inewsize)%nid )
+  if ( istatus /= NF90_NOERR ) &
+    call IO_Err_handle_nc4( istatus, 'IO_Strdimid_get_nc4', 'NF90_DEF_DIM', Trim( tzncdims(inewsize)%cname) )
+
+  call Move_alloc( from = tzncdims, to = tpfile%tncdims%tdims_str )
+  tpfile%tncdims%nmaxdims_str = inewsize
+
+  idx = inewsize
+
+  call Print_msg( NVERB_DEBUG, 'IO', 'IO_Dim_find_create_nc4','new dimension: ' &
+                  // Trim( tpfile%tncdims%tdims_str(inewsize)%cname ) )
+end if
 
-IO_Strdimid_get_nc4 = TMP%ID
+IO_Strdimid_get_nc4 = tpfile%tncdims%tdims_str(idx)%nid
 
 END FUNCTION IO_Strdimid_get_nc4
 
 
 FUNCTION IO_Iocdf_alloc_nc4()
-TYPE(IOCDF), POINTER :: IO_Iocdf_alloc_nc4
-TYPE(IOCDF), POINTER :: TZIOCDF
-INTEGER              :: IRESP
+type(tdimsnc), pointer :: IO_Iocdf_alloc_nc4
+type(tdimsnc), pointer :: tzdimsnc
+INTEGER                :: IRESP
 
 CALL PRINT_MSG(NVERB_DEBUG,'IO','IO_Iocdf_alloc_nc4','called')
 
-ALLOCATE(TZIOCDF, STAT=IRESP)
+ALLOCATE(tzdimsnc, STAT=IRESP)
 IF (IRESP > 0) THEN
   CALL PRINT_MSG(NVERB_FATAL,'IO','IO_Iocdf_alloc_nc4','memory allocation error')
 END IF
 
-IO_Iocdf_alloc_nc4=>TZIOCDF
+IO_Iocdf_alloc_nc4 => tzdimsnc
 
 END FUNCTION IO_Iocdf_alloc_nc4
 
diff --git a/src/LIB/SURCOUCHE/src/mode_io_write_nc4.f90 b/src/LIB/SURCOUCHE/src/mode_io_write_nc4.f90
index 7b49c02b5f327832ebef39e88f9a59c7815c55aa..808e451a2f6c45ba17b5bf6e18bfca54fc0ada8d 100644
--- a/src/LIB/SURCOUCHE/src/mode_io_write_nc4.f90
+++ b/src/LIB/SURCOUCHE/src/mode_io_write_nc4.f90
@@ -21,16 +21,16 @@
 !  P. Wautelet 25/06/2020: remove workaround for netCDF bug (see 19/09/2019)
 !  P. Wautelet 14/09/2020: IO_Coordvar_write_nc4: do not store 'time' coordinate in diachronic files
 !  P. Wautelet 22/09/2020: add ldimreduced field to allow reduction in the number of dimensions of fields (used by 2D simulations)
+!  P. Wautelet 10/11/2020: new data structures for netCDF dimensions
 !-----------------------------------------------------------------
 #ifdef MNH_IOCDF4
 module mode_io_write_nc4
 
 use modd_field,        only: tfielddata
 use modd_io,           only: gsmonoproc, tfiledata
-use modd_netcdf,       only: dimcdf, iocdf
 use modd_precision,    only: CDFINT, MNHINT_NF90, MNHREAL_MPI, MNHREAL_NF90
 
-use mode_io_tools_nc4, only: IO_Mnhname_clean, IO_Vdims_fill_nc4, IO_Dimcdf_get_nc4, IO_Strdimid_get_nc4, IO_Err_handle_nc4
+use mode_io_tools_nc4, only: IO_Mnhname_clean, IO_Vdims_fill_nc4, IO_Dim_find_create_nc4, IO_Strdimid_get_nc4, IO_Err_handle_nc4
 use mode_msg
 
 use NETCDF,            only: NF90_CHAR, NF90_FLOAT, NF90_INT1,                                    &
@@ -159,7 +159,7 @@ SUBROUTINE IO_Field_attr_write_nc4(TPFILE,TPFIELD,KVARID,OEXISTED,KSHAPE,HCALEND
 !
 USE MODD_CONF,   ONLY: CPROGRAM, LCARTESIAN
 USE MODD_CONF_n, ONLY: CSTORAGE_TYPE
-use modd_field,  only: TYPEINT, TYPEREAL
+use modd_field,  only: NMNHDIM_ARAKAWA, TYPEINT, TYPEREAL
 !
 TYPE(TFILEDATA),                              INTENT(IN) :: TPFILE
 TYPE(TFIELDDATA),                             INTENT(IN) :: TPFIELD
@@ -252,7 +252,8 @@ IF (.NOT.GISCOORD) THEN
   !1D: no direct correspondance with latitude(_x)/longitude(_x) 2D variables => nothing to do
   IF (.NOT.LCARTESIAN .AND. TPFIELD%NDIMS>1 .AND. TPFIELD%NGRID/=0) THEN
     IF (TPFIELD%CDIR=='XY') THEN
-      IF (KSHAPE(1)==TPFILE%TNCCOORDS(1,TPFIELD%NGRID)%TDIM%LEN .AND. KSHAPE(2)==TPFILE%TNCCOORDS(2,TPFIELD%NGRID)%TDIM%LEN ) THEN
+      if (       kshape(1) == tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(tpfield%ngrid,1) )%nlen &
+           .and. kshape(2) == tpfile%tncdims%tdims( NMNHDIM_ARAKAWA(tpfield%ngrid,2) )%nlen ) then
         SELECT CASE(TPFIELD%NGRID)
           CASE (0) !Not on Arakawa grid
             !Nothing to do
@@ -837,9 +838,11 @@ END SUBROUTINE IO_Field_write_nc4_X6
 SUBROUTINE IO_Field_write_nc4_N0(TPFILE,TPFIELD,KFIELD,KRESP)
 !
 #if 0
+use modd_field,          only: NMNHDIM_NI, NMNHDIM_NJ, NMNHDIM_LEVEL
 USE MODD_IO,             ONLY: LPACK,L1D,L2D
 USE MODD_PARAMETERS_ll,  ONLY: JPHEXT, JPVEXT
 #else
+use modd_field,          only: NMNHDIM_LEVEL
 USE MODD_PARAMETERS_ll,  ONLY: JPVEXT
 #endif
 !
@@ -848,6 +851,7 @@ TYPE(TFIELDDATA),      INTENT(IN) :: TPFIELD
 INTEGER,               INTENT(IN) :: KFIELD
 INTEGER,               INTENT(OUT):: KRESP
 !
+integer                                         :: iidx
 INTEGER(KIND=CDFINT)                            :: STATUS
 INTEGER(KIND=CDFINT)                            :: INCID
 CHARACTER(LEN=LEN(TPFIELD%CMNHNAME))            :: YVARNAME
@@ -855,7 +859,6 @@ INTEGER(KIND=CDFINT)                            :: IVARID
 INTEGER(KIND=CDFINT), DIMENSION(:), ALLOCATABLE :: IVDIMS
 INTEGER                                         :: IRESP
 LOGICAL                                         :: GEXISTED !True if variable was already defined
-TYPE(IOCDF), POINTER                            :: TZIOCDF
 !
 CALL PRINT_MSG(NVERB_DEBUG,'IO','IO_Field_write_nc4_N0',TRIM(TPFILE%CNAME)//': writing '//TRIM(TPFIELD%CMNHNAME))
 !
@@ -898,19 +901,20 @@ IF (status /= NF90_NOERR) CALL IO_Err_handle_nc4(status,'IO_Field_write_nc4_N0',
 ! /!\ Can only work if IMAX, JMAX or KMAX are written before any array
 !
 #if 0
-IF (YVARNAME == 'IMAX' .AND. .NOT. ASSOCIATED(TPFILE%TNCDIMS%DIM_NI)) TPFILE%TNCDIMS%DIM_NI=>IO_Dimcdf_get_nc4(TPFILE%TNCDIMS,KFIELD+2*JPHEXT,'X')
-IF (YVARNAME == 'JMAX' .AND. .NOT. ASSOCIATED(TPFILE%TNCDIMS%DIM_NJ)) THEN
-   IF (LPACK .AND. L2D) THEN
-      TPFILE%TNCDIMS%DIM_NJ=>IO_Dimcdf_get_nc4(TPFILE, 1,'Y')
-   ELSE
-      TPFILE%TNCDIMS%DIM_NJ=>IO_Dimcdf_get_nc4(TPFILE, KFIELD+2*JPHEXT, 'Y')
-   END IF
-END IF
+if ( yvarname == 'IMAX' .and. tpfile%tncdims%tdims(NMNHDIM_NI)%nid == -1 ) then
+  call IO_Dim_find_create_nc4( tpfile, kfield + 2 * jphext, iidx, 'X' )
+end if
+if ( yvarname == 'JMAX' .and. tpfile%tncdims%tdims(NMNHDIM_NJ)%nid == -1 ) then
+  if ( lpack .and. l2d ) then
+    call IO_Dim_find_create_nc4( tpfile, 1,                   iidx, 'Y' )
+  else
+    call IO_Dim_find_create_nc4( tpfile, kfield + 2 * jphext, iidx, 'Z' )
+  end if
+end if
 #endif
-IF (YVARNAME == 'KMAX' .AND. .NOT. ASSOCIATED(TPFILE%TNCDIMS%DIM_LEVEL)) THEN
-  TZIOCDF => TPFILE%TNCDIMS
-  TZIOCDF%DIM_LEVEL=>IO_Dimcdf_get_nc4(TPFILE,INT(KFIELD+2*JPVEXT,KIND=CDFINT),'Z')
-END IF
+if ( yvarname == 'KMAX' .and. tpfile%tncdims%tdims(NMNHDIM_LEVEL)%nid == -1 ) then
+  call IO_Dim_find_create_nc4( tpfile, kfield + 2 * JPVEXT, iidx, 'Z' )
+end if
 
 KRESP = IRESP
 END SUBROUTINE IO_Field_write_nc4_N0
@@ -1625,10 +1629,12 @@ END SUBROUTINE IO_Field_write_nc4_T1
 SUBROUTINE IO_Coordvar_write_nc4(TPFILE,HPROGRAM_ORIG)
 USE MODD_CONF,       ONLY: CPROGRAM, LCARTESIAN
 USE MODD_CONF_n,     ONLY: CSTORAGE_TYPE
-use modd_field,      only: tfieldlist
+use modd_field,      only: NMNHDIM_NI, NMNHDIM_NJ, NMNHDIM_NI_U, NMNHDIM_NJ_U, NMNHDIM_NI_V, NMNHDIM_NJ_V, &
+                           NMNHDIM_LEVEL, NMNHDIM_LEVEL_W, NMNHDIM_TIME,                                   &
+                           tfieldlist
 USE MODD_GRID,       ONLY: XLATORI, XLONORI
 USE MODD_GRID_n,     ONLY: LSLEVE, XXHAT, XYHAT, XZHAT
-use modd_netcdf,     only: dimcdf
+use modd_netcdf,     only: tdimnc
 USE MODD_PARAMETERS, ONLY: JPHEXT, JPVEXT
 
 use mode_field,      only: Find_field_id_from_mnhname
@@ -1650,8 +1656,7 @@ LOGICAL,POINTER                 :: GSLEVE
 REAL,DIMENSION(:),POINTER       :: ZXHAT, ZYHAT, ZZHAT
 REAL,DIMENSION(:),ALLOCATABLE   :: ZXHATM, ZYHATM,ZZHATM !Coordinates at mass points in the transformed space
 REAL,DIMENSION(:,:),POINTER     :: ZLAT, ZLON
-type(dimcdf), pointer           :: tzdim_ni, tzdim_nj, tzdim_ni_u, tzdim_nj_u, tzdim_ni_v, tzdim_nj_v
-TYPE(IOCDF),  POINTER           :: PIOCDF
+type(tdimnc), pointer           :: tzdim_ni, tzdim_nj, tzdim_ni_u, tzdim_nj_u, tzdim_ni_v, tzdim_nj_v
 
 !These variables are save: they are populated once for the master Z-split file and freed after the last file has been written
 real, dimension(:),   pointer, save :: zxhat_glob  => null(), zyhat_glob  => null()
@@ -1667,8 +1672,6 @@ ZXHAT => NULL()
 ZYHAT => NULL()
 ZZHAT => NULL()
 
-PIOCDF => TPFILE%TNCDIMS
-
 GCHANGEMODEL = .FALSE.
 
 IF (PRESENT(HPROGRAM_ORIG)) THEN
@@ -1719,23 +1722,22 @@ ELSE
   YSTDNAMEPREFIX = 'projection'
 ENDIF
 
-if(associated(piocdf)) then
-tzdim_ni   => piocdf%dim_ni
-tzdim_nj   => piocdf%dim_nj
-tzdim_ni_u => piocdf%dim_ni_u
-tzdim_nj_u => piocdf%dim_nj_u
-tzdim_ni_v => piocdf%dim_ni_v
-tzdim_nj_v => piocdf%dim_nj_v
+if ( Associated( tpfile%tncdims ) ) then
+  tzdim_ni   => tpfile%tncdims%tdims(NMNHDIM_NI)
+  tzdim_nj   => tpfile%tncdims%tdims(NMNHDIM_NJ)
+  tzdim_ni_u => tpfile%tncdims%tdims(NMNHDIM_NI_U)
+  tzdim_nj_u => tpfile%tncdims%tdims(NMNHDIM_NJ_U)
+  tzdim_ni_v => tpfile%tncdims%tdims(NMNHDIM_NI_V)
+  tzdim_nj_v => tpfile%tncdims%tdims(NMNHDIM_NJ_V)
 else
-tzdim_ni   => null()
-tzdim_nj   => null()
-tzdim_ni_u => null()
-tzdim_nj_u => null()
-tzdim_ni_v => null()
-tzdim_nj_v => null()
+  tzdim_ni   => Null()
+  tzdim_nj   => Null()
+  tzdim_ni_u => Null()
+  tzdim_nj_u => Null()
+  tzdim_ni_v => Null()
+  tzdim_nj_v => Null()
 end if
 
-
 !If the file is a Z-split subfile, coordinates are already collected
 if ( .not. associated( tpfile%tmainfile ) ) then
   call Gather_hor_coord1d( 'X', zxhat,  zxhat_glob  )
@@ -1809,10 +1811,10 @@ IF (TPFILE%LMASTER) THEN !vertical coordinates in the transformed space are the
     ZZHATM(1:IKU-1) = 0.5 * (ZZHAT(2:IKU)+ZZHAT(1:IKU-1))
     ZZHATM(IKU)     = 2.* ZZHAT(IKU) - ZZHATM(IKU-1)
     !
-    CALL WRITE_VER_COORD(PIOCDF%DIM_LEVEL,  'position z in the transformed space',              '', &
+    CALL WRITE_VER_COORD(tpfile%tncdims%tdims(NMNHDIM_LEVEL),  'position z in the transformed space',              '', &
                          'altitude',               0., JPVEXT,JPVEXT,ZZHATM)
     !
-    CALL WRITE_VER_COORD(PIOCDF%DIM_LEVEL_W,'position z in the transformed space at w location','', &
+    CALL WRITE_VER_COORD(tpfile%tncdims%tdims(NMNHDIM_LEVEL_W),'position z in the transformed space at w location','', &
                          'altitude_at_w_location',-0.5,JPVEXT,0,     ZZHAT)
     !
     DEALLOCATE(ZZHATM)
@@ -1824,7 +1826,7 @@ IF (TPFILE%LMASTER) THEN !Time scale is the same on all processes
   IF (TRIM(YPROGRAM)/='PGD' .AND. TRIM(YPROGRAM)/='NESPGD' .AND. TRIM(YPROGRAM)/='ZOOMPG' &
       .AND. .NOT.(TRIM(YPROGRAM)=='REAL' .AND. CSTORAGE_TYPE=='SU') ) THEN !condition to detect PREP_SURFEX
     if ( tpfile%ctype /= 'MNHDIACHRONIC' ) &
-      CALL WRITE_TIME_COORD(PIOCDF%DIMTIME)
+      CALL WRITE_TIME_COORD(tpfile%tncdims%tdims(NMNHDIM_TIME))
   END IF
 END IF
 
@@ -1931,7 +1933,7 @@ subroutine Write_hor_coord1d(TDIM,HLONGNAME,HSTDNAME,HAXIS,PSHIFT,KBOUNDLOW,KBOU
   USE MODE_ALLOCBUFFER_ll, ONLY: ALLOCBUFFER_ll
   USE MODE_GATHER_ll,      ONLY: GATHER_XXFIELD
 
-  TYPE(DIMCDF), POINTER,      INTENT(IN) :: TDIM
+  TYPE(tdimnc), POINTER,      INTENT(IN) :: TDIM
   CHARACTER(LEN=*),           INTENT(IN) :: HLONGNAME
   CHARACTER(LEN=*),           INTENT(IN) :: HSTDNAME
   CHARACTER(LEN=*),           INTENT(IN) :: HAXIS
@@ -1949,9 +1951,9 @@ subroutine Write_hor_coord1d(TDIM,HLONGNAME,HSTDNAME,HAXIS,PSHIFT,KBOUNDLOW,KBOU
   INTEGER(KIND=CDFINT)          :: ISTATUS
 
   IF (TPFILE%LMASTER) THEN
-    ISIZE = TDIM%LEN
-    YVARNAME = TRIM(TDIM%NAME)
-    IVDIM = TDIM%ID
+    isize    = tdim%nlen
+    yvarname = Trim( tdim%cname )
+    ivdim    = tdim%nid
 
     ISTATUS = NF90_INQ_VARID(INCID, YVARNAME, IVARID)
     IF (ISTATUS /= NF90_NOERR) THEN
@@ -2009,7 +2011,7 @@ end subroutine Write_hor_coord2d
 
 
 SUBROUTINE WRITE_VER_COORD(TDIM,HLONGNAME,HSTDNAME,HCOMPNAME,PSHIFT,KBOUNDLOW,KBOUNDHIGH,PCOORDS)
-  TYPE(DIMCDF), POINTER, INTENT(IN) :: TDIM
+  TYPE(tdimnc), POINTER, INTENT(IN) :: TDIM
   CHARACTER(LEN=*),      INTENT(IN) :: HLONGNAME
   CHARACTER(LEN=*),      INTENT(IN) :: HSTDNAME
   CHARACTER(LEN=*),      INTENT(IN) :: HCOMPNAME
@@ -2027,9 +2029,9 @@ SUBROUTINE WRITE_VER_COORD(TDIM,HLONGNAME,HSTDNAME,HCOMPNAME,PSHIFT,KBOUNDLOW,KB
   INTEGER(KIND=CDFINT)          :: IVDIM
   INTEGER(KIND=CDFINT)          :: STATUS
 
-  ISIZE = TDIM%LEN
-  YVARNAME = TRIM(TDIM%NAME)
-  IVDIM = TDIM%ID
+  isize    = tdim%nlen
+  yvarname = Trim( tdim%cname )
+  ivdim    = tdim%nid
 
   STATUS = NF90_INQ_VARID(INCID, YVARNAME, IVARID)
   IF (STATUS /= NF90_NOERR) THEN
@@ -2097,7 +2099,7 @@ SUBROUTINE WRITE_TIME_COORD(TDIM)
   use mode_field,      only: Find_field_id_from_mnhname
   USE MODE_GRIDPROJ
 
-  TYPE(DIMCDF), POINTER, INTENT(IN) :: TDIM
+  TYPE(tdimnc), POINTER, INTENT(IN) :: TDIM
 
   REAL                         :: ZDELTATIME
   CHARACTER(LEN=40)            :: YUNITS
@@ -2109,8 +2111,8 @@ SUBROUTINE WRITE_TIME_COORD(TDIM)
 
 
   IF (ASSOCIATED(TDTCUR) .AND. ASSOCIATED(TDTMOD)) THEN
-    YVARNAME = TRIM(TDIM%NAME)
-    IVDIM = TDIM%ID
+    yvarname = Trim( tdim%cname )
+    ivdim    = tdim%nid
 
     STATUS = NF90_INQ_VARID(INCID, YVARNAME, IVARID)
     IF (STATUS /= NF90_NOERR) THEN