diff --git a/src/LIB/SURCOUCHE/src/modd_errcodes.f90 b/src/LIB/SURCOUCHE/src/modd_errcodes.f90
deleted file mode 100644
index d580fbce6b436e8ae8d543e369ff1ee90f5fb329..0000000000000000000000000000000000000000
--- a/src/LIB/SURCOUCHE/src/modd_errcodes.f90
+++ /dev/null
@@ -1,28 +0,0 @@
-!MNH_LIC Copyright 1994-2014 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 for details. version 1.
-!-----------------------------------------------------------------
-!--------------- special set of characters for CVS information
-!-----------------------------------------------------------------
-! $Source$
-! $Name$ 
-! $Revision$ 
-! $Date$
-!-----------------------------------------------------------------
-!-----------------------------------------------------------------
-
-! USED ONLY BY IO
-! ---------------
-
-MODULE MODD_ERRCODES
-  IMPLICIT NONE 
-
-  !! Error codes
-  INTEGER, PARAMETER :: NOERROR    =  0
-  INTEGER, PARAMETER :: IOERROR    = -1
-  INTEGER, PARAMETER :: NOSLOTLEFT = -2
-  INTEGER, PARAMETER :: BADVALUE   = -3
-  INTEGER, PARAMETER :: UNDEFINED  = -999
-
-END MODULE MODD_ERRCODES
diff --git a/src/LIB/SURCOUCHE/src/modd_io.f90 b/src/LIB/SURCOUCHE/src/modd_io.f90
index da81705fd9f038825de6100ad5797140444387c7..21e82e1da640584501076fdca952c74ab22f228b 100644
--- a/src/LIB/SURCOUCHE/src/modd_io.f90
+++ b/src/LIB/SURCOUCHE/src/modd_io.f90
@@ -1,10 +1,11 @@
-!MNH_LIC Copyright 1994-2018 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2019 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
+!  Philippe Wautelet: 10/01/2019: use NEWUNIT argument of OPEN (removed ISTDOUT, ISTDERR, added NNULLUNIT, CNULLFILE)
 !-----------------------------------------------------------------
 
 MODULE MODD_IO_ll
@@ -17,8 +18,8 @@ IMPLICIT NONE
 !
 INTEGER, PARAMETER :: NVERB_NO=0, NVERB_FATAL=1, NVERB_ERROR=2, NVERB_WARNING=3, NVERB_INFO=4, NVERB_DEBUG=5
 
-
-INTEGER, SAVE :: ISTDOUT, ISTDERR
+INTEGER                     :: NNULLUNIT = -1  ! /dev/null fortran unit, value set in INITIO_ll
+CHARACTER(LEN=*), PARAMETER :: CNULLFILE = "/dev/null"
 
 INTEGER, SAVE :: ISIOP   !! IOproc number
 INTEGER, SAVE :: ISP     !! Actual proc number
diff --git a/src/LIB/SURCOUCHE/src/mode_fm.f90 b/src/LIB/SURCOUCHE/src/mode_fm.f90
index 4723ab0e7530bcd9a6f42e25bd6cddfc9261adb6..eea6d328d7a0a7b8509c4c600ab4ebe05871957f 100644
--- a/src/LIB/SURCOUCHE/src/mode_fm.f90
+++ b/src/LIB/SURCOUCHE/src/mode_fm.f90
@@ -1,4 +1,4 @@
-!MNH_LIC Copyright 1994-2018 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2019 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 for details. version 1.
@@ -9,10 +9,11 @@
 !  Philippe Wautelet: 05/2016-04/2018: new data structures and calls for I/O
 !  Philippe Wautelet: 29/10/2018: better detection of older MNH version numbers
 !  Philippe Wautelet: 13/12/2018: moved some operations to new mode_io_*_nc4 modules
+!  Philippe Wautelet: 10/01/2019: use NEWUNIT argument of OPEN + move management
+!                                 of NNCID and NLFIFLU to the nc4 and lfi subroutines
 !-----------------------------------------------------------------
 
 MODULE MODE_FM
-USE MODD_ERRCODES
 USE MODD_MPIF
 
 USE MODE_MSG
@@ -156,8 +157,8 @@ IF (PRESENT(KRESP)) KRESP = IRESP
 END SUBROUTINE IO_FILE_OPEN_ll
 
 SUBROUTINE FMOPEN_ll(TPFILE,KRESP,OPARALLELIO,HPROGRAM_ORIG)
-USE MODD_IO_ll,               ONLY: ISTDOUT,TFILEDATA
-USE MODE_IO_ll,               ONLY: OPEN_ll,GCONFIO,IOFREEFLU,IONEWFLU
+USE MODD_IO_ll,               ONLY: TFILEDATA
+USE MODE_IO_ll,               ONLY: OPEN_ll,GCONFIO
 !JUANZ
 USE MODD_CONFZ,ONLY  : NB_PROCIO_R,NB_PROCIO_W
 !JUANZ
@@ -234,28 +235,22 @@ IF (TPFILE%LMASTER) THEN
           CALL PRINT_MSG(NVERB_WARNING,'IO','FMOPEN_ll',TRIM(TPFILE%CNAME)// &
                          ': .nc file does not exist but .lfi exists -> forced to LFI')
           TPFILE%CFORMAT='LFI'
-          TPFILE%NLFIFLU = IONEWFLU()
         END IF
       CASE ('LFI')
         IF (.NOT.GEXIST_LFI .AND. GEXIST_NC4) THEN
           CALL PRINT_MSG(NVERB_WARNING,'IO','FMOPEN_ll',TRIM(TPFILE%CNAME)// &
                          ': .lfi file does not exist but .nc exists -> forced to NETCDF4')
           TPFILE%CFORMAT='NETCDF4'
-          CALL IOFREEFLU(INT(TPFILE%NLFIFLU))
-          TPFILE%NLFIFLU = -1
         END IF
       CASE ('LFICDF4')
         IF (GEXIST_NC4) THEN
           CALL PRINT_MSG(NVERB_WARNING,'IO','FMOPEN_ll',TRIM(TPFILE%CNAME)// &
                          ': LFICDF4 format is not allowed in READ mode -> forced to NETCDF4')
           TPFILE%CFORMAT='NETCDF4'
-          IF (TPFILE%NLFIFLU>0) CALL IOFREEFLU(INT(TPFILE%NLFIFLU))
-          TPFILE%NLFIFLU = -1
         ELSE IF (GEXIST_LFI) THEN
           CALL PRINT_MSG(NVERB_WARNING,'IO','FMOPEN_ll',TRIM(TPFILE%CNAME)// &
                          ': LFICDF4 format is not allowed in READ mode -> forced to LFI')
           TPFILE%CFORMAT='LFI'
-          TPFILE%NLFIFLU = IONEWFLU()
         END IF
       CASE DEFAULT
         IF (GEXIST_NC4) THEN
@@ -266,7 +261,6 @@ IF (TPFILE%LMASTER) THEN
           CALL PRINT_MSG(NVERB_ERROR,'IO','FMOPEN_ll',TRIM(TPFILE%CNAME)// &
                          ': invalid fileformat (-> forced to LFI if no abort)')
           TPFILE%CFORMAT='LFI'
-          TPFILE%NLFIFLU = IONEWFLU()
         END IF
     END SELECT
   END IF
@@ -421,9 +415,6 @@ SELECT CASE(TPFILE%CTYPE)
     !
     CALL FMCLOS_ll(TPFILE,KRESP=IRESP,OPARALLELIO=OPARALLELIO,HPROGRAM_ORIG=HPROGRAM_ORIG)
     !
-    TPFILE%NLFIFLU = -1
-    TPFILE%NNCID   = -1
-    !
     DO JI = 1,TPFILE%NSUBFILES_IOZ
       TZFILE_IOZ => TPFILE%TFILES_IOZ(JI)%TFILE
       IF (.NOT.TZFILE_IOZ%LOPENED) &
@@ -434,8 +425,6 @@ SELECT CASE(TPFILE%CTYPE)
       TZFILE_IOZ%LOPENED       = .FALSE.
       TZFILE_IOZ%NOPEN_CURRENT = 0
       TZFILE_IOZ%NCLOSE        = TZFILE_IOZ%NCLOSE + 1
-      TZFILE_IOZ%NLFIFLU       = -1
-      TZFILE_IOZ%NNCID         = -1
     END DO
 END SELECT
 !
@@ -465,7 +454,7 @@ USE MODI_SYSTEM_MNH
   use mode_io_file_nc4,  only: io_close_file_nc4
   use mode_io_write_nc4, only: io_write_coordvar_nc4
 #endif
-TYPE(TFILEDATA),      INTENT(IN) :: TPFILE ! File structure
+TYPE(TFILEDATA),      INTENT(INOUT)         :: TPFILE ! File structure
 INTEGER,              INTENT(OUT), OPTIONAL :: KRESP   ! return-code if problems araised
 LOGICAL,              INTENT(IN),  OPTIONAL :: OPARALLELIO
 CHARACTER(LEN=*),     INTENT(IN),  OPTIONAL :: HPROGRAM_ORIG !To emulate a file coming from this program
diff --git a/src/LIB/SURCOUCHE/src/mode_io.f90 b/src/LIB/SURCOUCHE/src/mode_io.f90
index 18060d752c455bdd9407802a95c0b9814e1ee788..1a087fd50519ee9bf16b8e004251bab140dbcc3c 100644
--- a/src/LIB/SURCOUCHE/src/mode_io.f90
+++ b/src/LIB/SURCOUCHE/src/mode_io.f90
@@ -19,10 +19,12 @@
 !     J. Pianezze 01/08/2016  add LOASIS flag
 !     Philippe Wautelet: 13/12/2018: moved some operations to new mode_io_*_nc4 modules
 !     Philippe Wautelet: 10/01/2019: bug correction: close correctly Z-split files
+!     Philippe Wautelet: 10/01/2019: use NEWUNIT argument of OPEN
+!                                    + move IOFREEFLU and IONEWFLU to mode_io_file_lfi.f90
+!                                    + move management of NNCID and NLFIFLU to the nc4 and lfi subroutines
 !
 MODULE MODE_IO_ll
 
-  USE MODD_ERRCODES
   USE MODD_MPIF
   USE MODD_VAR_ll, ONLY : NMNH_COMM_WORLD
 
@@ -32,59 +34,13 @@ MODULE MODE_IO_ll
 
   PRIVATE
 
-  INTEGER, PARAMETER :: JPFNULL = 9       !! /dev/null fortran unit
-  INTEGER, PARAMETER :: JPRESERVED_UNIT   = 11
-  INTEGER, PARAMETER :: JPMAX_UNIT_NUMBER = JPRESERVED_UNIT+300
-  ! 
-  LOGICAL,SAVE :: GALLOC(JPRESERVED_UNIT:JPMAX_UNIT_NUMBER) = .FALSE.
-  !
-  CHARACTER(LEN=*),PARAMETER      :: CFILENULL="/dev/null"
-  !
   LOGICAL,SAVE :: GCONFIO = .FALSE. ! Turn TRUE when SET_CONFIO_ll is called.
 
-  PUBLIC IOFREEFLU,IONEWFLU,UPCASE,INITIO_ll,OPEN_ll,CLOSE_ll
+  PUBLIC UPCASE,INITIO_ll,OPEN_ll,CLOSE_ll
   PUBLIC SET_CONFIO_ll,GCONFIO
 
 CONTAINS 
 
-  FUNCTION IONEWFLU()
-
-    INTEGER :: IONEWFLU
-
-    INTEGER :: JI
-    INTEGER :: IOS
-    LOGICAL :: GEXISTS, GOPENED, GFOUND
-
-    GFOUND = .FALSE.
-
-    DO JI=JPRESERVED_UNIT, JPMAX_UNIT_NUMBER
-       IF (GALLOC(JI)) CYCLE
-       INQUIRE(UNIT=JI, EXIST=GEXISTS, OPENED=GOPENED, IOSTAT=IOS)
-       IF (GEXISTS .AND. .NOT. GOPENED .AND. IOS == 0) THEN
-          IONEWFLU   = JI
-          GFOUND     = .TRUE.
-          GALLOC(JI) = .TRUE.
-          EXIT
-       END IF
-    END DO
-
-    IF (.NOT. GFOUND) IONEWFLU = NOSLOTLEFT
-
-  END FUNCTION IONEWFLU
-
-  SUBROUTINE IOFREEFLU(KOFLU)
-
-    INTEGER :: KOFLU
-
-    IF ((KOFLU .GE. JPRESERVED_UNIT) .AND. (KOFLU .LE. JPMAX_UNIT_NUMBER )) THEN 
-       GALLOC(KOFLU) = .FALSE.
-    ELSE
-       print*,"mode_io.f90: IOFREEFLU BAD IUNIT=",KOFLU
-       STOP "mode_io.f90: IOFREEFLU BAD IUNIT"
-    END IF
-
-  END SUBROUTINE IOFREEFLU
-
   FUNCTION UPCASE(HSTRING)
     CHARACTER(LEN=*)            :: HSTRING
     CHARACTER(LEN=LEN(HSTRING)) :: UPCASE
@@ -137,9 +93,9 @@ CONTAINS
     END IF
     
   END SUBROUTINE SET_CONFIO_INTERN_ll
-  
+
   SUBROUTINE INITIO_ll()
-    USE  MODE_MNH_WORLD , ONLY :  INIT_NMNH_COMM_WORLD
+    USE MODE_MNH_WORLD, ONLY: INIT_NMNH_COMM_WORLD
     USE MODD_IO_ll
     USE MODE_FIELD
     IMPLICIT NONE
@@ -148,10 +104,8 @@ CONTAINS
 
     CALL PRINT_MSG(NVERB_DEBUG,'IO','INITIO_ll','called')
 
-    ISTDERR = 0
-
     CALL INIT_NMNH_COMM_WORLD(IERR)
-    IF (IERR .NE.0) CALL PRINT_MSG(NVERB_FATAL,'IO','SET_CONFIO_ll','problem with remapping of NMNH_COMM_WORLD')
+    IF (IERR .NE.0) CALL PRINT_MSG(NVERB_FATAL,'IO','INITIO_ll','problem with remapping of NMNH_COMM_WORLD')
 
     !! Now MPI is initialized for sure
 
@@ -168,22 +122,13 @@ CONTAINS
 
     !! Open /dev/null for GLOBAL mode
 #if defined(DEV_NULL)
-    OPEN(UNIT=JPFNULL,FILE=CFILENULL  ,ACTION='WRITE',IOSTAT=IOS)
+    OPEN(NEWUNIT=NNULLUNIT,FILE=CNULLFILE  ,ACTION='WRITE',IOSTAT=IOS)
 #else
-    OPEN(UNIT=JPFNULL,STATUS='SCRATCH',ACTION='WRITE',IOSTAT=IOS)
+    OPEN(NEWUNIT=NNULLUNIT,STATUS='SCRATCH',ACTION='WRITE',IOSTAT=IOS)
 #endif
     IF (IOS > 0) THEN
-       WRITE(ISTDERR,*) 'Error OPENING /dev/null...'
-       CALL MPI_ABORT(NMNH_COMM_WORLD, IOS, IERR)
-    END IF
-
-    !! Init STDOUT and PIPE
-    IF (ISP == ISIOP) THEN
-       ISTDOUT = 6
-    ELSE
-       ISTDOUT = JPFNULL
+       CALL PRINT_MSG(NVERB_FATAL,'IO','INITIO_ll','error opening /dev/null')
     END IF
-
   END SUBROUTINE INITIO_ll
 
   SUBROUTINE OPEN_ll(&
@@ -292,14 +237,14 @@ CONTAINS
     IF (YACTION /= "READ" .AND. YACTION /= "WRITE") THEN
        IOSTAT = 99
        TPFILE%NLU = -1
-       WRITE(ISTDERR,*) 'Erreur OPEN_ll : ACTION=',YACTION,' non supportee'
+       CALL PRINT_MSG(NVERB_ERROR,'IO','OPEN_ll','action='//TRIM(YACTION)//' not supported')
        RETURN
     END IF
 
     IF (.NOT. ANY(YMODE == (/'GLOBAL     ','SPECIFIC   ','DISTRIBUTED' , 'IO_ZSPLIT  '/))) THEN
        IOSTAT = 99
        TPFILE%NLU = -1
-       WRITE(ISTDERR,*) 'OPEN_ll error : MODE UNKNOWN'
+       CALL PRINT_MSG(NVERB_ERROR,'IO','OPEN_ll','ymode='//TRIM(YMODE)//' not supported')
        RETURN
     END IF
 
@@ -385,10 +330,8 @@ CONTAINS
 
        IF (TPFILE%LMASTER) THEN
           !! I/O processor case
-
-          TPFILE%NLU = IONEWFLU()
 #ifdef MNH_VPP
-          OPEN(UNIT=TPFILE%NLU,        &
+          OPEN(NEWUNIT=TPFILE%NLU,     &
                FILE=TRIM(YPREFILENAME),&
                STATUS=STATUS,          &
                ACCESS=ACCESS,          &
@@ -406,7 +349,7 @@ CONTAINS
 #if defined(MNH_SX5) || defined(MNH_SP4) || defined(NAGf95) || defined(MNH_LINUX)
           !JUAN : 31/03/2000 modif pour acces direct
           IF (YACCESS=='STREAM') THEN
-             OPEN(UNIT=TPFILE%NLU,        &
+             OPEN(NEWUNIT=TPFILE%NLU,     &
                   FILE=TRIM(YPREFILENAME),&
                   STATUS=YSTATUS,         &
                   ACCESS=YACCESS,         &
@@ -415,7 +358,7 @@ CONTAINS
                   FORM=YFORM,             &
                   ACTION=YACTION)
           ELSEIF (YACCESS=='DIRECT') THEN
-             OPEN(UNIT=TPFILE%NLU,        &
+             OPEN(NEWUNIT=TPFILE%NLU,     &
                   FILE=TRIM(YPREFILENAME),&
                   STATUS=YSTATUS,         &
                   ACCESS=YACCESS,         &
@@ -427,7 +370,7 @@ CONTAINS
           ELSE
              IF (YFORM=="FORMATTED") THEN
                IF (YACTION=='READ') THEN
-                OPEN(UNIT=TPFILE%NLU,        &
+                OPEN(NEWUNIT=TPFILE%NLU,     &
                      FILE=TRIM(YPREFILENAME),&
                      STATUS=YSTATUS,         &
                      ACCESS=YACCESS,         &
@@ -441,7 +384,7 @@ CONTAINS
                      !DELIM=YDELIM,          & !Philippe: commented because bug with GCC 5.X
                      PAD=YPAD)
                ELSE
-                OPEN(UNIT=TPFILE%NLU,        &
+                OPEN(NEWUNIT=TPFILE%NLU,     &
                      FILE=TRIM(YPREFILENAME),&
                      STATUS=YSTATUS,         &
                      ACCESS=YACCESS,         &
@@ -456,7 +399,7 @@ CONTAINS
                      PAD=YPAD)
                ENDIF
              ELSE
-                OPEN(UNIT=TPFILE%NLU,        &
+                OPEN(NEWUNIT=TPFILE%NLU,     &
                      FILE=TRIM(YPREFILENAME),&
                      STATUS=YSTATUS,         &
                      ACCESS=YACCESS,         &
@@ -471,7 +414,7 @@ CONTAINS
 
 
           !print*,' OPEN_ll'
-          !print*,' OPEN(UNIT=',TPFILE%NLU
+          !print*,' OPEN(NEWUNIT=',TPFILE%NLU
           !print*,' FILE=',TRIM(YPREFILENAME)
           !print*,' STATUS=',YSTATUS       
           !print*,' ACCESS=',YACCESS
@@ -484,7 +427,7 @@ CONTAINS
           !print*,' DELIM=',YDELIM
           !print*,' PAD=',YPAD
 #else
-          OPEN(UNIT=TPFILE%NLU,        &
+          OPEN(NEWUNIT=TPFILE%NLU,     &
                FILE=TRIM(YPREFILENAME),&
                STATUS=STATUS,          &
                ACCESS=ACCESS,          &
@@ -501,22 +444,21 @@ CONTAINS
 
 #endif
           IF (IOS/=0) CALL PRINT_MSG(NVERB_FATAL,'IO','OPEN_ll','Problem when opening '//TRIM(YPREFILENAME)//': '//TRIM(YIOERRMSG))
-       ELSE 
+       ELSE
           !! NON I/O processors case
           IOS = 0
-          TPFILE%NLU = JPFNULL
+          TPFILE%NLU = NNULLUNIT
        END IF
 
 
     CASE('SPECIFIC')
-       TPFILE%NLU = IONEWFLU()
        TPFILE%NMASTER_RANK  = -1
        TPFILE%LMASTER       = .TRUE. !Every process use the file
        TPFILE%LMULTIMASTERS = .TRUE.
        TPFILE%NSUBFILES_IOZ = 0
 
 #ifdef MNH_VPP
-       OPEN(UNIT=TPFILE%NLU,                       &
+       OPEN(NEWUNIT=TPFILE%NLU,                    &
             FILE=TRIM(YPREFILENAME)//SUFFIX(".P"), &
             STATUS=STATUS,                         &
             ACCESS=ACCESS,                         &
@@ -533,7 +475,7 @@ CONTAINS
 #else
 #if defined(MNH_SX5) || defined(MNH_SP4) || defined(NAGf95) || defined(MNH_LINUX)
        IF (ACCESS=='DIRECT') THEN
-          OPEN(UNIT=TPFILE%NLU,                       &
+          OPEN(NEWUNIT=TPFILE%NLU,                    &
                FILE=TRIM(YPREFILENAME)//SUFFIX(".P"), &
                STATUS=YSTATUS,                        &
                ACCESS=YACCESS,                        &
@@ -544,7 +486,7 @@ CONTAINS
                ACTION=YACTION)
        ELSE
         IF (YACTION=='READ') THEN
-          OPEN(UNIT=TPFILE%NLU,                        &
+          OPEN(NEWUNIT=TPFILE%NLU,                     &
                FILE=TRIM(YPREFILENAME)//SUFFIX(".P"),  &
                STATUS=YSTATUS,                         &
                ACCESS=YACCESS,                         &
@@ -558,7 +500,7 @@ CONTAINS
                !DELIM=YDELIM,         & !Philippe: commented because bug with GCC 5.X
                PAD=YPAD)
          ELSE
-          OPEN(UNIT=TPFILE%NLU,                        &
+          OPEN(NEWUNIT=TPFILE%NLU,                     &
                FILE=TRIM(YPREFILENAME)//SUFFIX(".P"),  &
                STATUS=YSTATUS,                         &
                ACCESS=YACCESS,                         &
@@ -574,7 +516,7 @@ CONTAINS
          ENDIF
        ENDIF
 #else
-       OPEN(UNIT=TPFILE%NLU,                       &
+       OPEN(NEWUNIT=TPFILE%NLU,                    &
             FILE=TRIM(YPREFILENAME)//SUFFIX(".P"), &
             STATUS=STATUS,                         &
             ACCESS=ACCESS,                         &
@@ -600,12 +542,9 @@ CONTAINS
        TPFILE%LMULTIMASTERS = .FALSE.
        TPFILE%NSUBFILES_IOZ = 0
 
-       IF (TPFILE%LMASTER) THEN
-          TPFILE%NLU = IONEWFLU()
-       ELSE 
+       IF (.NOT.TPFILE%LMASTER) THEN
           !! NON I/O processors case
           IOS = 0
-          TPFILE%NLU = -1
        END IF
 
 
@@ -626,11 +565,9 @@ CONTAINS
 #else
        IF (TPFILE%LMASTER) THEN
 #endif
-             TPFILE%NLFIFLU = IONEWFLU()
        ELSE 
           !! NON I/O processors OR netCDF read case
           IOS = 0
-          TPFILE%NLFIFLU = -1
        END IF
 
        IF (TPFILE%NSUBFILES_IOZ > 0) THEN
@@ -761,9 +698,8 @@ CONTAINS
     IGLOBALERR2 = 0
 
     IF (TPFILE%LMASTER) THEN
-      IF (TPFILE%NLU>0 .AND. TPFILE%NLU/=JPFNULL) THEN
+      IF (TPFILE%NLU/=-1 .AND. TPFILE%NLU/=NNULLUNIT) THEN
         CLOSE(UNIT=TPFILE%NLU, IOSTAT=IRESP,STATUS='KEEP')
-        CALL IOFREEFLU(TPFILE%NLU)
       END IF
     END IF
     !
@@ -781,9 +717,6 @@ CONTAINS
 #if defined(MNH_IOCDF4)
           if (tzfile%cformat == 'NETCDF4' .or. tzfile%cformat == 'LFICDF4') call io_close_file_nc4(tzfile,iresp2)
 #endif
-          IF (TZFILE%NLFIFLU > 0) THEN !if LFI
-            CALL IOFREEFLU(INT(TZFILE%NLFIFLU))
-          END IF
         END IF
       END DO
       !
diff --git a/src/LIB/SURCOUCHE/src/mode_io_file_lfi.f90 b/src/LIB/SURCOUCHE/src/mode_io_file_lfi.f90
index e477754859822771d91d2c6df108474ae1aa4b23..b79ff2ec1bb67c395249497bc4f61f403f30227b 100644
--- a/src/LIB/SURCOUCHE/src/mode_io_file_lfi.f90
+++ b/src/LIB/SURCOUCHE/src/mode_io_file_lfi.f90
@@ -1,4 +1,4 @@
-!MNH_LIC Copyright 2018-2018 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 2018-2019 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 for details. version 1.
@@ -9,6 +9,9 @@
 !           (was duplicated in the 2 files)
 !
 !  Modifications:
+!     Philippe Wautelet: 10/01/2019: use NEWUNIT argument of OPEN
+!                                    + move IOFREEFLU and IONEWFLU to mode_io_file_lfi.f90
+!                                    + move management of NNCID and NLFIFLU to the nc4 and lfi subroutines
 !
 !-----------------------------------------------------------------
 module mode_io_file_lfi
@@ -24,6 +27,11 @@ private
 
 public :: io_create_file_lfi, io_close_file_lfi, io_open_file_lfi
 
+integer, parameter :: JPRESERVED_UNIT   = 11
+integer, parameter :: JPMAX_UNIT_NUMBER = JPRESERVED_UNIT + 300
+
+logical,save :: galloc(JPRESERVED_UNIT:JPMAX_UNIT_NUMBER) = .false.
+
 contains
 
 subroutine io_create_file_lfi(tpfile, kstatus)
@@ -50,15 +58,17 @@ subroutine io_create_file_lfi(tpfile, kstatus)
     call io_construct_filename(tpfile, yfilem)
 
     iresou = 0
-    inumbr = tpfile%nlfiflu
+    if ( tpfile%nlfiflu /= -1 ) call print_msg(NVERB_ERROR,'IO', &
+                                               'io_create_file_lfi','file '//trim(yfilem)//'.lfi has already a unit number')
+    tpfile%nlfiflu = ionewflu()
     gnamfi = .true.
     yforstatus = 'REPLACE'
     gfater = .true.
 
     call io_prepare_verbosity_lfi(tpfile, imelev, gstats)
 
+    inumbr = tpfile%nlfiflu
     inprar = tpfile%nlfinprar
-
     call lfiouv(iresou, inumbr, gnamfi, trim(yfilem)//'.lfi', yforstatus, gfater, gstats, imelev, inprar, ininar)
 
     tpfile%nlfininar = ininar
@@ -68,7 +78,7 @@ subroutine io_create_file_lfi(tpfile, kstatus)
     !test if file is newly defined
     gnewfi = (ininar==0) .or. (imelev<2)
     if (.not.gnewfi) then
-      call print_msg(NVERB_INFO,'IO','file '//trim(yfilem)//'.lfi',' previously created with LFI')
+      call print_msg(NVERB_INFO,'IO','io_create_file_lfi','file '//trim(yfilem)//'.lfi previously created with LFI')
     endif
   end if
   call io_set_mnhversion(tpfile)
@@ -76,10 +86,8 @@ end subroutine io_create_file_lfi
 
 
 subroutine io_close_file_lfi(tpfile, kstatus)
-!   use mode_io_tools_nc4, only: cleaniocdf
-
-  type(tfiledata),   intent(in)  :: tpfile
-  integer, optional, intent(out) :: kstatus
+  type(tfiledata),   intent(inout)  :: tpfile
+  integer, optional, intent(out)    :: kstatus
 
   character(len=*), parameter :: YSTATUS = 'KEEP'
 
@@ -90,8 +98,10 @@ subroutine io_close_file_lfi(tpfile, kstatus)
   istatus = 0
 
   if (tpfile%lmaster) then
-    if ( tpfile%nlfiflu > 0 ) then
+    if ( tpfile%nlfiflu /= -1 ) then
       call lfifer(istatus, tpfile%nlfiflu, YSTATUS)
+      call iofreeflu(int(tpfile%nlfiflu))
+      tpfile%nlfiflu = -1
     else
       istatus = -1
       call print_msg(NVERB_WARNING, 'IO', 'io_close_file_lfi', 'file '//trim(tpfile%cname)//'.lfi is not opened')
@@ -127,13 +137,16 @@ subroutine io_open_file_lfi(tpfile, kstatus)
     call io_construct_filename(tpfile, yfilem)
 
     iresou = 0
-    inumbr = tpfile%nlfiflu
+    if ( tpfile%nlfiflu /= -1 ) call print_msg(NVERB_ERROR,'IO', &
+                                               'io_open_file_lfi','file '//trim(yfilem)//'.lfi has already a unit number')
+    tpfile%nlfiflu = ionewflu()
     gnamfi = .true.
     yforstatus = 'OLD'
     gfater = .true.
 
     call io_prepare_verbosity_lfi(tpfile, imelev, gstats)
 
+    inumbr = tpfile%nlfiflu
     inprar = tpfile%nlfinprar
 
     call lfiouv(iresou, inumbr, gnamfi, trim(yfilem)//'.lfi', yforstatus, gfater, gstats, imelev, inprar, ininar)
@@ -146,4 +159,44 @@ subroutine io_open_file_lfi(tpfile, kstatus)
 end subroutine io_open_file_lfi
 
 
+function ionewflu()
+  use modd_io_ll, only: nnullunit
+
+  integer :: ionewflu
+
+  integer :: ji
+  integer :: ios
+  logical :: gexists, gopened, gfound
+
+  gfound = .false.
+
+  do ji = JPRESERVED_UNIT, JPMAX_UNIT_NUMBER
+    if ( galloc(ji) ) cycle
+    inquire(unit=ji, exist=gexists, opened=gopened, iostat=ios)
+    if (gexists .and. .not. gopened .and. ios == 0) then
+      ionewflu   = ji
+      gfound     = .true.
+      galloc(ji) = .true.
+      exit
+    end if
+  end do
+
+  if (.not. gfound) then
+    call print_msg(NVERB_ERROR,'IO','ionewflu','wrong unit number')
+    ionewflu = nnullunit !/dev/null Fortran unit
+  end if
+end function ionewflu
+
+
+subroutine iofreeflu(koflu)
+  integer :: koflu
+
+  if ( (koflu >= JPRESERVED_UNIT) .and. (koflu <= JPMAX_UNIT_NUMBER) ) then
+    galloc(koflu) = .false.
+  else
+    call print_msg(NVERB_ERROR,'IO','iofreeflu','wrong unit number')
+  end if
+end subroutine iofreeflu
+
+
 end module mode_io_file_lfi
diff --git a/src/LIB/SURCOUCHE/src/mode_io_file_nc4.f90 b/src/LIB/SURCOUCHE/src/mode_io_file_nc4.f90
index 4738a2fc0b4e3537e7b498bba02b2a382d630615..8a6c45a711512944f68356fe00344cc7cf7cb66d 100644
--- a/src/LIB/SURCOUCHE/src/mode_io_file_nc4.f90
+++ b/src/LIB/SURCOUCHE/src/mode_io_file_nc4.f90
@@ -9,6 +9,9 @@
 !           (was duplicated in the 2 files)
 !
 !  Modifications:
+!     Philippe Wautelet: 10/01/2019: use NEWUNIT argument of OPEN
+!                                    + move IOFREEFLU and IONEWFLU to mode_io_file_lfi.f90
+!                                    + move management of NNCID and NLFIFLU to the nc4 and lfi subroutines
 !
 !-----------------------------------------------------------------
 #if defined(MNH_IOCDF4)
@@ -62,8 +65,8 @@ end subroutine io_create_file_nc4
 subroutine io_close_file_nc4(tpfile,kstatus)
   use mode_io_tools_nc4, only: cleaniocdf
 
-  type(tfiledata),                    intent(in)  :: tpfile
-  integer(kind=IDCDF_KIND), optional, intent(out) :: kstatus
+  type(tfiledata),                    intent(inout) :: tpfile
+  integer(kind=IDCDF_KIND), optional, intent(out)   :: kstatus
 
   integer(kind=IDCDF_KIND) :: istatus
 
@@ -80,6 +83,7 @@ subroutine io_close_file_nc4(tpfile,kstatus)
       if (istatus /= NF90_NOERR) then
         call print_msg(NVERB_WARNING, 'IO', 'io_close_file_nc4', 'NF90_CLOSE error: '//trim(NF90_STRERROR(istatus)))
       end if
+      tpfile%nncid = -1
       if (associated(tpfile%tncdims)) call cleaniocdf(tpfile%tncdims)
     end if
   end if
diff --git a/src/LIB/SURCOUCHE/src/mode_mnh_world.f90 b/src/LIB/SURCOUCHE/src/mode_mnh_world.f90
index 7873bdb8cbcbc730441646a98f6c16c7221991db..67ef014021be97382819f8df9f4d05f06ab1329a 100644
--- a/src/LIB/SURCOUCHE/src/mode_mnh_world.f90
+++ b/src/LIB/SURCOUCHE/src/mode_mnh_world.f90
@@ -1,6 +1,6 @@
-!MNH_LIC Copyright 1994-2018 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2019 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
@@ -8,6 +8,7 @@
 !!
 !!  J.Escobar 3/12/2014 : typo form -> from
 !!  Philippe 03/10/2017: set IP and NPROC in INIT_NMNH_COMM_WORLD
+!!  Philippe Wautelet: 10/01/2019: use NEWUNIT argument of OPEN
 !!
 MODULE MODE_MNH_WORLD
   IMPLICIT NONE
@@ -42,7 +43,8 @@ CONTAINS
 #ifdef MNH_GA
     INTEGER                         :: REQUIRED=MPI_THREAD_MULTIPLE,PROVIDED 
 #endif
-    !JUANZ    
+    !JUANZ
+    INTEGER :: ILU
 
     !
     KINFO_ll = 0
@@ -62,13 +64,13 @@ CONTAINS
        ! Read namelist config file
        !
        IF ( irank .EQ. 0 ) THEN
-          OPEN(unit=10,form="formatted",file=conf_mnh_world,STATUS='OLD',iostat=IERR)
+          OPEN(newunit=ILU,form="formatted",file=conf_mnh_world,STATUS='OLD',iostat=IERR)
           ! Read IO parameter
           IF (IERR.EQ.0) THEN
-             READ(10,NML=NAM_CONF_MNH_WORLD) 
+             READ(unit=ilu,NML=NAM_CONF_MNH_WORLD)
              WRITE(*,NAM_CONF_MNH_WORLD)
           ENDIF
-          CLOSE(10)
+          CLOSE(unit=ILU)
        ENDIF
        iroot = 0
        ! Brodcast mapping
@@ -85,12 +87,12 @@ CONTAINS
 !          myrank_key =  IY_COORD(IRANK+1) + IX_COORD(IRANK+1) * IX_DIM
 !          IF (IRANK .EQ. 0 ) THEN
 !             print *,"IX_DIM=",IX_DIM
-!             OPEN(unit=100,form="formatted",file="hilbert_2D")
+!             OPEN(newunit=ILU,form="formatted",file="hilbert_2D")
 !             DO  I=1,IPROC
-!                write(100,*) I, IX_COORD(I),IY_COORD(I) , &
+!                write(unit=ILU,fmt=*) I, IX_COORD(I),IY_COORD(I) , &
 !                     1+ IY_COORD(I) + IX_COORD(I) * IX_DIM 
 !             END DO
-!             CLOSE(100)
+!             CLOSE(unit=ILU)
 !          END IF
 !          color = 1
 !          call MPI_Comm_split( MPI_COMM_WORLD, color, myrank_key,NMNH_COMM_WORLD , KINFO_ll )
diff --git a/src/LIB/SURCOUCHE/src/mode_mppdb.f90 b/src/LIB/SURCOUCHE/src/mode_mppdb.f90
index c3181c2c73131f4998464b3e0e994f82335c21f7..d3693aa9c9224c780e95f0ec88edfb410e9e0864 100644
--- a/src/LIB/SURCOUCHE/src/mode_mppdb.f90
+++ b/src/LIB/SURCOUCHE/src/mode_mppdb.f90
@@ -1,6 +1,6 @@
-!MNH_LIC Copyright 1994-2018 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2019 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.
 !-----------------------------------------------------------------
 MODULE MODE_MPPDB
@@ -11,6 +11,7 @@ MODULE MODE_MPPDB
 !  J.Escobar : 15/09/2015 : WENO5 & JPHEXT <> 1 
 !  G.Delautier : 23/06/2016 : surfex v8
 !  Philippe Wautelet: 05/2016-04/2018: new data structures and calls for I/O
+!  Philippe Wautelet: 10/01/2019: use NEWUNIT argument of OPEN
 !
   IMPLICIT NONE
 
@@ -46,16 +47,15 @@ CONTAINS
     !pas de mpi_spawn sur IBM-SP ni MPI_ARGV_NULL etc ...
     RETURN           
 #else
-    !USE MPI
+    USE MODD_MPIF
     !JUANZ
     USE  MODE_MNH_WORLD , ONLY :  INIT_NMNH_COMM_WORLD
     USE  MODD_VAR_ll    , ONLY :  NMNH_COMM_WORLD
     !JUANZ
     IMPLICIT NONE
-    INCLUDE "mpif.h"
 
 
-    INTEGER                         :: IUNIT = 100
+    INTEGER                         :: IUNIT
     INTEGER                         :: IERR
     INTEGER                         :: IRANK_WORLD,IRANK_INTRA
     INTEGER                         :: NBPROC_WORLD,NBPROC_INTRA
@@ -104,7 +104,7 @@ CONTAINS
        !
        ! if no config file , inactive MPPDB routines
        !
-       OPEN(unit=IUNIT,file=MPPDB_CONF,STATUS='OLD',iostat=IERR)
+       OPEN(newunit=IUNIT,file=MPPDB_CONF,STATUS='OLD',iostat=IERR)
        IF (IERR.NE.0 ) THEN
           ! PRINT*,"IOSTAT=",IERR
           !
diff --git a/src/LIB/SURCOUCHE/src/system_mnh.f90 b/src/LIB/SURCOUCHE/src/system_mnh.f90
index d2355e28e9eb2bd489f06d464fe00c869880c672..ef39625e42eee1684e444aef8f0fc4562c7816c7 100644
--- a/src/LIB/SURCOUCHE/src/system_mnh.f90
+++ b/src/LIB/SURCOUCHE/src/system_mnh.f90
@@ -1,47 +1,47 @@
-!MNH_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2019 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.
-! $Source$
 !-----------------------------------------------------------------
-          SUBROUTINE SYSTEM_MNH(HCOMMAND)
+SUBROUTINE SYSTEM_MNH(HCOMMAND)
 !!
 !!    PURPOSE
 !!    -------
 !!    This subroutine writes the 1 command line HCOMMAND
 !!    in the file pipe_name and flushes the buffer.
-USE MODE_IO_ll
+!!    Modifications:
+!!      Philippe Wautelet: 10/01/2019: use NEWUNIT argument of OPEN
+!!
+  USE MODE_IO_ll
 !!
 !*       0.    DECLARATIONS
 !              ------------
 !
-IMPLICIT NONE
+  IMPLICIT NONE
 !
 !*       0.1   Declaration of arguments
 !              ------------------------
 !
-          CHARACTER(LEN=*)    :: HCOMMAND
+  CHARACTER(LEN=*)    :: HCOMMAND
 !
 !*       0.2   Declaration of local variables
 !              ------------------------------
 #if defined(MNH_LINUX) || defined(MNH_SP4)
-         CHARACTER(LEN=*),PARAMETER :: CFILE="file_for_xtransfer"
+  CHARACTER(LEN=*),PARAMETER :: CFILE="file_for_xtransfer"
 #else
 #if !defined(MNH_SX5)
-         CHARACTER(LEN=*),PARAMETER :: CFILE="file_for_fujitransfer"
+  CHARACTER(LEN=*),PARAMETER :: CFILE="file_for_fujitransfer"
 #else
-         CHARACTER(LEN=*),PARAMETER :: CFILE="file_for_nectransfer"
+  CHARACTER(LEN=*),PARAMETER :: CFILE="file_for_nectransfer"
 #endif
 #endif
-         INTEGER                    :: IUNIT
+  INTEGER                    :: IUNIT
 !
 !
 !
-          IUNIT=IONEWFLU()
-          OPEN(UNIT=IUNIT,FILE=CFILE,ACCESS="sequential" &
-&             ,FORM="formatted",POSITION="append")
-          WRITE(IUNIT,*) HCOMMAND
-          CLOSE(IUNIT)
-
-          END SUBROUTINE SYSTEM_MNH
+  IUNIT = -1
+  OPEN(NEWUNIT=IUNIT,FILE=CFILE,ACCESS="sequential",FORM="formatted",POSITION="append")
+  WRITE(IUNIT,*) HCOMMAND
+  CLOSE(IUNIT)
 
+END SUBROUTINE SYSTEM_MNH
diff --git a/src/MNH/compute_spectre.f90 b/src/MNH/compute_spectre.f90
index 3006a8268f89e5b9498a29d052964ccf4bcb9663..6a31130f07b6813137c6e781faf88fa1c42ee87a 100644
--- a/src/MNH/compute_spectre.f90
+++ b/src/MNH/compute_spectre.f90
@@ -1,6 +1,6 @@
-!MNH_LIC Copyright 1994-2018 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2019 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.
 !     ####################
       MODULE MODI_COMPUTE_SPECTRE
@@ -49,6 +49,7 @@ END MODULE MODI_COMPUTE_SPECTRE
 !	A. Mary, R. Legrand          **ENM**
 !       D. Ricard    **CNRM**            
 !!  Philippe Wautelet: 05/2016-04/2018: new data structures and calls for I/O
+!!  Philippe Wautelet: 10/01/2019: use NEWUNIT argument of OPEN
 !-------------------------------------------------------------------------------
 !
 !*       0.    DECLARATIONS
@@ -114,6 +115,7 @@ INTEGER :: IIX,IJX,IIY,IJY ! dimensions of the extended x or y slices subdomain
 REAL, DIMENSION(:,:,:), ALLOCATABLE :: ZBAND_YT  ! array in Y slices distribution transpose
 !
 ! local variables for section 5
+CHARACTER(len=9)                    :: YJNB       ! for error message
 INTEGER                             :: JNB        ! loop index : adimensionned wavenumbers
 INTEGER                             :: JERR       ! allocation or writing errors
 INTEGER                             :: IMIN       ! minimum domain dimension
@@ -142,6 +144,7 @@ REAL, DIMENSION(:,:), ALLOCATABLE   :: ZSP        ! results table
 !         ----!-------------------------------
 !         ... !        !         !
 !
+INTEGER                         :: ILU
 INTEGER                         :: ILUOUT    ! Logical unit number for
                                              ! output-listing   
 INTEGER                         :: IRESP     ! return code in FM routines         
@@ -469,29 +472,24 @@ ELSE
   WRITE(YFMT,FMT='(A,I1,A,I3,A)') "(I",JLOGIMIN,",F13.2,",IKU,"F25.16)"
 ENDIF
 YFMT=TRIM(ADJUSTL(YFMT))
-
-OPEN(UNIT=32, FILE=YOUTFILE, ACCESS='SEQUENTIAL', IOSTAT=JERR)
-IF (JERR /= 0) THEN 
-  print*,"error when open the  file ",YOUTFILE, " JERR=",JERR
-  CALL ABORT
-  STOP
+!
+OPEN(NEWUNIT=ILU, FILE=YOUTFILE, ACCESS='SEQUENTIAL', IOSTAT=JERR)
+IF (JERR /= 0) THEN
+  CALL PRINT_MSG(NVERB_FATAL,'IO','COMPUTE_SPECTRE','error when opening '//trim(YOUTFILE))
 ENDIF
+!
 DO JNB=1,IMIN-1
- WRITE(UNIT=32,IOSTAT=JERR,FMT=YFMT) INT(ZLGO(JNB,1)),ZLGO(JNB,2),ZSP(JNB,:)
- IF (JERR /= 0) THEN 
-   print*,"error when writing JNB=",JNB," JERR=",JERR
-   CALL ABORT
-   STOP
- ENDIF
-
+  WRITE(UNIT=ILU,IOSTAT=JERR,FMT=YFMT) INT(ZLGO(JNB,1)),ZLGO(JNB,2),ZSP(JNB,:)
+  IF (JERR /= 0) THEN
+    WRITE(YJNB,'( I9 )') JNB
+    CALL PRINT_MSG(NVERB_FATAL,'IO','COMPUTE_SPECTRE','error when writing JNB='//trim(YJNB)//' in file '//trim(YOUTFILE))
+  ENDIF
 END DO
-CLOSE(UNIT=32, IOSTAT=JERR)
+!
+CLOSE(UNIT=ILU, IOSTAT=JERR)
 IF (JERR /= 0) THEN 
-  print*,"error when closing the  file ",YOUTFILE, " JERR=",JERR
-  CALL ABORT
-  STOP
+  CALL PRINT_MSG(NVERB_ERROR,'IO','COMPUTE_SPECTRE','error when closing '//trim(YOUTFILE))
 ENDIF
-
 !
 !-------------------------------------------------------------------------------
 !
diff --git a/src/MNH/flash_geom_elec.f90 b/src/MNH/flash_geom_elec.f90
index bd800432fcc590afeedd516fde501d74df6b952c..868200a3186197b85fcc66567d9db984c5d555bd 100644
--- a/src/MNH/flash_geom_elec.f90
+++ b/src/MNH/flash_geom_elec.f90
@@ -1,17 +1,19 @@
-!MNH_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 2010-2019 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.
 !     #############################
       MODULE MODI_FLASH_GEOM_ELEC_n
 !     #############################
 !
 INTERFACE
-    SUBROUTINE FLASH_GEOM_ELEC_n (KTCOUNT, KMI, KRR, PTSTEP, OEXIT,    &
-                                  PRHODJ, PRHODREF,                    &
-                                  PRT, PCIT, PRSVS, PRS, PTHT, PPABST, &
-                                  PEFIELDU, PEFIELDV, PEFIELDW,        &
-                                  PZZ, PSVS_LINOX, PTOWN, PSEA         )
+    SUBROUTINE FLASH_GEOM_ELEC_n (KTCOUNT, KMI, KRR, PTSTEP, OEXIT,                      &
+                                  PRHODJ, PRHODREF, PRT, PCIT, PRSVS, PRS, PTHT, PPABST, &
+                                  PEFIELDU, PEFIELDV, PEFIELDW, PZZ, PSVS_LINOX,         &
+                                  TPFILE_FGEOM_DIAG, TPFILE_FGEOM_COORD, TPFILE_LMA,     &
+                                  PTOWN, PSEA                                            )
+!
+USE MODD_IO_ll, ONLY: TFILEDATA
 !
 INTEGER,                  INTENT(IN)    :: KTCOUNT  ! Temporal loop counter
 INTEGER,                  INTENT(IN)    :: KMI      ! current model index
@@ -32,6 +34,9 @@ REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PTHT     ! Theta (K) at time t
 REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PPABST   ! Absolute pressure at t
 REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PZZ      ! height
 REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PSVS_LINOX ! NOx source term
+TYPE(TFILEDATA),          INTENT(IN)    :: TPFILE_FGEOM_DIAG
+TYPE(TFILEDATA),          INTENT(IN)    :: TPFILE_FGEOM_COORD
+TYPE(TFILEDATA),          INTENT(IN)    :: TPFILE_LMA
 REAL, DIMENSION(:,:), OPTIONAL, INTENT(IN) :: PTOWN ! town fraction
 REAL, DIMENSION(:,:), OPTIONAL, INTENT(IN) :: PSEA  ! Land-sea mask
 !
@@ -40,13 +45,13 @@ END INTERFACE
 END MODULE MODI_FLASH_GEOM_ELEC_n
 !
 !
-!       ####################################################################
-        SUBROUTINE FLASH_GEOM_ELEC_n (KTCOUNT, KMI, KRR, PTSTEP, OEXIT,    &
-                                      PRHODJ, PRHODREF,                    &
-                                      PRT, PCIT, PRSVS, PRS, PTHT, PPABST, &
-                                      PEFIELDU, PEFIELDV, PEFIELDW,        &
-                                      PZZ, PSVS_LINOX, PTOWN, PSEA         )
-!       ####################################################################
+!   ######################################################################################
+    SUBROUTINE FLASH_GEOM_ELEC_n (KTCOUNT, KMI, KRR, PTSTEP, OEXIT,                      &
+                                  PRHODJ, PRHODREF, PRT, PCIT, PRSVS, PRS, PTHT, PPABST, &
+                                  PEFIELDU, PEFIELDV, PEFIELDW, PZZ, PSVS_LINOX,         &
+                                  TPFILE_FGEOM_DIAG, TPFILE_FGEOM_COORD, TPFILE_LMA,     &
+                                  PTOWN, PSEA                                            )
+!   ######################################################################################
 !
 !!****  * -
 !!
@@ -86,6 +91,7 @@ END MODULE MODI_FLASH_GEOM_ELEC_n
 !!      J.Escobar : 20/06/2018 : Correction of computation of global index I8VECT
 !!      J.Escobar : 10/12/2018 : // Correction , mpi_bcast CG & CG_POS parameter 
 !!                               & initialize INBLIGHT on all proc for filling/saving AREA* arrays
+!!      Philippe Wautelet: 10/01/2019: use NEWUNIT argument of OPEN
 !!
 !-------------------------------------------------------------------------------
 !
@@ -94,6 +100,7 @@ END MODULE MODI_FLASH_GEOM_ELEC_n
 !
 USE MODD_CST, ONLY : XAVOGADRO, XMD
 USE MODD_CONF, ONLY : CEXP, LCARTESIAN
+USE MODD_IO_ll, ONLY: TFILEDATA
 USE MODD_PARAMETERS, ONLY : JPHEXT, JPVEXT
 USE MODD_GRID, ONLY : XLATORI,XLONORI
 USE MODD_GRID_n, ONLY : XXHAT, XYHAT, XZHAT
@@ -112,8 +119,6 @@ USE MODD_RAIN_ICE_DESCR, ONLY : XLBR, XLBEXR, XLBS, XLBEXS, &
 USE MODD_NSV, ONLY : NSV_ELECBEG, NSV_ELECEND, NSV_ELEC
 USE MODD_VAR_ll, ONLY : NPROC,NMNH_COMM_WORLD
 USE MODD_ARGSLIST_ll, ONLY : LIST_ll
-USE MODD_PRINT_ELEC,  ONLY : NLU_fgeom_diag, NLU_fgeom_coord, &
-                             NIOSTAT_fgeom_diag, NIOSTAT_fgeom_coord
 USE MODD_SUB_ELEC_n
 USE MODD_TIME_n
 USE MODD_LMA_SIMULATOR
@@ -156,6 +161,9 @@ REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PTHT     ! Theta (K) at time t
 REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PPABST   ! Absolute pressure at t
 REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PZZ      ! height
 REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PSVS_LINOX ! NOx source term
+TYPE(TFILEDATA),          INTENT(IN)    :: TPFILE_FGEOM_DIAG
+TYPE(TFILEDATA),          INTENT(IN)    :: TPFILE_FGEOM_COORD
+TYPE(TFILEDATA),          INTENT(IN)    :: TPFILE_LMA
 REAL, DIMENSION(:,:), OPTIONAL, INTENT(IN) :: PTOWN ! town fraction
 REAL, DIMENSION(:,:), OPTIONAL, INTENT(IN) :: PSEA  ! Land-sea mask
 !
@@ -1375,15 +1383,9 @@ ENDIF
       PRSVS(:,:,:,1) = PRSVS(:,:,:,1) / XECHARGE
       PRSVS(:,:,:,NSV_ELEC) = - PRSVS(:,:,:,NSV_ELEC) / XECHARGE
 !
-      IF (PRESENT(PSEA)) THEN
-        CALL ION_ATTACH_ELEC(KTCOUNT, KRR, PTSTEP, PRHODREF,                   &
-                             PRHODJ, PRSVS, PRS, PTHT, PCIT, PPABST, PEFIELDU, &
-                             PEFIELDV, PEFIELDW, GATTACH, PTOWN, PSEA          )
-      ELSE
-        CALL ION_ATTACH_ELEC(KTCOUNT, KRR, PTSTEP, PRHODREF,                   &
-                             PRHODJ, PRSVS, PRS, PTHT, PCIT, PPABST, PEFIELDU, & 
-                             PEFIELDV, PEFIELDW, GATTACH                       )
-      ENDIF
+      CALL ION_ATTACH_ELEC(KTCOUNT, KRR, PTSTEP, PRHODREF,                   &
+                           PRHODJ, PRSVS, PRS, PTHT, PCIT, PPABST, PEFIELDU, &
+                           PEFIELDV, PEFIELDW, GATTACH, PTOWN, PSEA          )
 !
       PRSVS(:,:,:,1) = PRSVS(:,:,:,1) * XECHARGE
       PRSVS(:,:,:,NSV_ELEC) = - PRSVS(:,:,:,NSV_ELEC) * XECHARGE
@@ -2471,19 +2473,19 @@ END SUBROUTINE PT_DISCHARGE
 IMPLICIT NONE
 !
 INTEGER :: I1, I2
+INTEGER :: ILU ! unit number for IO
 !
 !
 !*      1.     FLASH PARAMETERS
 !              ----------------
 !
-OPEN (UNIT=NLU_fgeom_diag, FILE=CEXP//"_fgeom_diag.asc", ACTION="WRITE", &
-      STATUS="OLD", FORM="FORMATTED", POSITION="APPEND")
+ILU = TPFILE_FGEOM_DIAG%NLU
 !
 ! Ecriture ascii dans CEXP//'_fgeom_diag.asc" defini dans RESOLVED_ELEC
 !
 IF (LCARTESIAN) THEN
   DO I1 = 1, NNBLIGHT
-    WRITE (UNIT=NLU_fgeom_diag,FMT='(I8,F9.1,I4,I6,I4,I6,F9.3,F12.3,F12.3,F9.3,F8.2,F9.2,f9.4)') &
+    WRITE (UNIT=ILU,FMT='(I8,F9.1,I4,I6,I4,I6,F9.3,F12.3,F12.3,F9.3,F8.2,F9.2,f9.4)') &
           ISFLASH_NUMBER(I1),              &
           ISTCOUNT_NUMBER(I1) * PTSTEP,    &
           ISCELL_NUMBER(I1),               &
@@ -2504,7 +2506,7 @@ ELSE
                                    ZSCOORD_SEG(I1,1,2),&
                                    ZLAT,ZLON)
 !
-    WRITE (UNIT=NLU_fgeom_diag,FMT='(I8,F9.1,I4,I6,I4,I6,F9.3,F12.3,F12.3,F9.3,F8.2,F9.2,f9.4)') &
+    WRITE (UNIT=ILU,FMT='(I8,F9.1,I4,I6,I4,I6,F9.3,F12.3,F12.3,F9.3,F8.2,F9.2,f9.4)') &
           ISFLASH_NUMBER(I1),              &
           ISTCOUNT_NUMBER(I1) * PTSTEP,    &
           ISCELL_NUMBER(I1),               &
@@ -2520,7 +2522,7 @@ ELSE
   END DO
 END IF
 !
-CLOSE (UNIT=NLU_fgeom_diag)
+CALL FLUSH(UNIT=ILU)
 !
 !
 !*      2.     FLASH SEGMENT COORDINATES
@@ -2530,12 +2532,11 @@ IF (LSAVE_COORD) THEN
 !
 ! Ecriture ascii dans CEXP//'_fgeom_coord.asc" defini dans RESOLVED_ELEC
 !
-  OPEN (UNIT=NLU_fgeom_coord, FILE=CEXP//"_fgeom_coord.asc", ACTION="WRITE", &
-        STATUS="OLD", FORM="FORMATTED", POSITION="APPEND")
+  ILU = TPFILE_FGEOM_COORD%NLU
 !
   DO I1 = 1, NNBLIGHT
     DO I2 = 1, ISNBSEG(I1)
-      WRITE (NLU_fgeom_coord, FMT='(I4,F9.1,I4,F12.3,F12.3,F12.3)') &
+      WRITE (ILU, FMT='(I4,F9.1,I4,F12.3,F12.3,F12.3)') &
                  ISFLASH_NUMBER(I1),           & 
                  ISTCOUNT_NUMBER(I1) * PTSTEP, & 
                  ISTYPE(I1),                   &
@@ -2545,7 +2546,7 @@ IF (LSAVE_COORD) THEN
     END DO
   END DO
 !
-  CLOSE (UNIT=NLU_fgeom_coord)
+  CALL FLUSH(UNIT=ILU)
 END IF
 !
 END SUBROUTINE WRITE_OUT_ASCII
@@ -2564,6 +2565,7 @@ SUBROUTINE WRITE_OUT_LMA
 IMPLICIT NONE
 !
 INTEGER :: I1, I2
+INTEGER :: ILU ! unit number for IO
 !
 !
 !*      1.     LMA SIMULATOR
@@ -2571,9 +2573,12 @@ INTEGER :: I1, I2
 !
 CALL SM_LATLON(XLATORI,XLONORI,ZSCOORD_SEG(:,:,1),ZSCOORD_SEG(:,:,2), &
                                ZLMA_LAT(:,:),ZLMA_LON(:,:))
+!
+ILU = TPFILE_LMA%NLU
+!
 DO I1 = 1, NNBLIGHT
   DO I2 = 1, ISNBSEG(I1)
-    WRITE (ILMA_UNIT,FMT='(I6,F12.1,I6,2(F15.6),3(F15.3),3(I6),12(E15.4))') &
+    WRITE (UNIT=ILU,FMT='(I6,F12.1,I6,2(F15.6),3(F15.3),3(I6),12(E15.4))') &
                ISFLASH_NUMBER(I1),           &
                ISTCOUNT_NUMBER(I1) * PTSTEP, &
                ISTYPE(I1),                   &
@@ -2600,6 +2605,8 @@ DO I1 = 1, NNBLIGHT
   END DO
 END DO
 !
+CALL FLUSH(UNIT=ILU)
+!
 END SUBROUTINE WRITE_OUT_LMA
 !
 !-------------------------------------------------------------------------------
diff --git a/src/MNH/modd_lma_simulator.f90 b/src/MNH/modd_lma_simulator.f90
index 5cd2fe548f8c109e0edcfb4a56a504cf592cda8f..2946c6df964996d6e6060ed9844b77c9ee0693f7 100644
--- a/src/MNH/modd_lma_simulator.f90
+++ b/src/MNH/modd_lma_simulator.f90
@@ -1,6 +1,6 @@
-!MNH_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 2013-2019 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.
 !     ############################
       MODULE MODD_LMA_SIMULATOR
@@ -28,6 +28,7 @@
 !!    MODIFICATIONS
 !!    -------------
 !!      Original    15/02/13
+!!  Philippe Wautelet: 10/01/2019: use NEWUNIT argument of OPEN
 !!
 !-------------------------------------------------------------------------------
 !
@@ -47,8 +48,6 @@ LOGICAL, SAVE                          :: LLMA=.FALSE.! Flag to record LMA-like
 REAL, SAVE                             :: XDTLMA ! Time length of a LMA record
 TYPE (DATE_TIME), SAVE                 :: TDTLMA ! Date and Time of LMA file
 CHARACTER (LEN=31), SAVE               :: CLMA_FILE   ! File name
-INTEGER, SAVE                          :: ILMA_UNIT   ! File information
-INTEGER, SAVE                          :: ILMA_IOSTAT ! File information
 !
 !* storage monitoring
 !
diff --git a/src/MNH/modd_print_elec.f90 b/src/MNH/modd_print_elec.f90
deleted file mode 100644
index 070b835dfd8ddced8897e66bdbf7c2cccbb1b51f..0000000000000000000000000000000000000000
--- a/src/MNH/modd_print_elec.f90
+++ /dev/null
@@ -1,42 +0,0 @@
-!MNH_LIC Copyright 1994-2014 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 for details. version 1.
-!       #######################
-        MODULE  MODD_PRINT_ELEC
-!       #######################
-!
-!!****  *MODD_PRINT_ELEC* - declaration of the LU and IOSTAT for extra prints
-!!
-!!	PURPOSE
-!!	-------
-!
-!!**	IMPLICIT ARGUMENTS
-!!	------------------
-!!	  None
-!!
-!!	REFERENCE
-!!	---------
-!!
-!!	AUTHOR
-!!	------
-!!       Jean-Pierre Pinty    * Laboratoire d'Aerologie *
-!!
-!!	MODIFICATIONS
-!!	-------------
-!!	  Original	17/01/2012
-!!
-!-------------------------------------------------------------------------------
-!
-!*	0.	DECLARATIONS
-!		------------
-!
-IMPLICIT NONE
-!
-INTEGER :: NLU_series_cloud_elec, NIOSTAT_series_cloud_elec, &
-           NLU_fgeom_diag,        NIOSTAT_fgeom_diag,        &
-           NLU_fgeom_coord,       NIOSTAT_fgeom_coord,       &
-           NLU_light_diag,        NIOSTAT_light_diag,        &
-           NLU_light_coord,       NIOSTAT_light_coord
-!
-END MODULE MODD_PRINT_ELEC
diff --git a/src/MNH/mode_mnh_timing.f90 b/src/MNH/mode_mnh_timing.f90
index cb98e458f7c2c32c6c6db37cb909013091c133ed..35d7c559ad4193581e19a2e885659d9287214e37 100644
--- a/src/MNH/mode_mnh_timing.f90
+++ b/src/MNH/mode_mnh_timing.f90
@@ -1,6 +1,6 @@
-!MNH_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2019 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.
 MODULE MODE_MNH_TIMING
 !
@@ -8,6 +8,7 @@ MODULE MODE_MNH_TIMING
 ! J.ESCOBAR 13/11/2008 : change (2) in (:) for bug in IBM-SP6 compiler
 ! J.Escobar 1/09/2011  : reduce 'timing' format
 ! J.Escobar 12/02/2013 : tribulle to slow on large BG partition , inhib it by a early return in the code 
+! Philippe Wautelet: 10/01/2019: use NEWUNIT argument of OPEN
 !
 
 INTEGER     :: NLUOUT_TIMING
@@ -111,6 +112,7 @@ END SUBROUTINE SECOND_MNH2
 
   REAL*8,DIMENSION(2,NPROC)   :: ZSTAT_ALL 
   INTEGER, DIMENSION(NPROC) :: IND
+  INTEGER :: ILU
 !
 !-------------------------------------------------------------------------------
 !
@@ -164,14 +166,14 @@ INFO = -1
             END SELECT
 
          END DO
-         OPEN (100+NLUOUT_TIMING,file=FILE)
-         WRITE(100+NLUOUT_TIMING,FMT= "(10('#'),A30,10('#'))" ) HPRINT//VIDE
+         OPEN (newunit=ILU,file=FILE)
+         WRITE(unit=ILU,FMT= "(10('#'),A30,10('#'))" ) HPRINT//VIDE
          call  triabulle(ZSTAT_ALL(2,:),ind)
          DO JP=1,NPROC
-            WRITE(100+NLUOUT_TIMING,FMT= "(5(I8,' ',F15.3,' '))" ) &
+            WRITE(unit=ILU,FMT= "(5(I8,' ',F15.3,' '))" ) &
                  ind(JP),ZSTAT_ALL(2,ind(JP))
          END DO
-         CLOSE (100+NLUOUT_TIMING)
+         CLOSE (unit=ILU)
       END IF
    END IF
    !
diff --git a/src/MNH/read_asc_latpress.f90 b/src/MNH/read_asc_latpress.f90
index d34571053f25c318a522cf557b3620c29f61af40..61745f0d15b225622ff7c647f17b122796e25dfe 100644
--- a/src/MNH/read_asc_latpress.f90
+++ b/src/MNH/read_asc_latpress.f90
@@ -1,6 +1,6 @@
-!MNH_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2019 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.
 
 !     ########################
@@ -44,6 +44,7 @@ END MODULE MODI_READ_ASC_LATPRESS
 !!    MODIFICATION
 !!    ------------
 !!
+!!  Philippe Wautelet: 10/01/2019: use NEWUNIT argument of OPEN
 !!
 !----------------------------------------------------------------------------
 !
@@ -85,16 +86,14 @@ CALL GET_INDICE_ll(IIB,IJB,IIE,IJE)
 !
 !----------------------------------------------------------------------------
 !----------------------------------------------------------------------------
-KUNIT=221
-OPEN(KUNIT,FILE=HFILENAME)
+OPEN(NEWUNIT=KUNIT,FILE=HFILENAME)
 !
 !*    3.     Reading of a data point
 !            -----------------------
 !
 DO JI=IIB,IIE
   DO JK=1,KLEV
-  READ(KUNIT,*,END=99) PLAT(JI),PLEV(JK),                                 & 
-                      PTHFRC(JI,JK),PRVFRC(JI,JK)
+    READ(UNIT=KUNIT,FMT=*,END=99) PLAT(JI),PLEV(JK),PTHFRC(JI,JK),PRVFRC(JI,JK)
   END DO
 END DO
 !
@@ -104,7 +103,7 @@ END DO
 !*    8.    Closing of the data file
 !           ------------------------
 !
-99 CLOSE(KUNIT)
+99 CLOSE(UNIT=KUNIT)
 !
 !-------------------------------------------------------------------------------
 !
diff --git a/src/MNH/read_ascp.f90 b/src/MNH/read_ascp.f90
index 9af4d6802dfd8bef737bbf4cd29a0ceb10a9dd1b..e5490ed548be54325183e78431fd4b6afeef98ef 100644
--- a/src/MNH/read_ascp.f90
+++ b/src/MNH/read_ascp.f90
@@ -1,6 +1,6 @@
-!MNH_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2019 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.
 
 !     ########################
@@ -40,6 +40,7 @@ END MODULE MODI_READ_ASCP
 !!    MODIFICATION
 !!    ------------
 !!
+!!  Philippe Wautelet: 10/01/2019: use NEWUNIT argument of OPEN
 !!
 !----------------------------------------------------------------------------
 !
@@ -74,15 +75,14 @@ REAL  :: ZLEV
 !*    1.      Open the file
 !             -------------
 !
-KUNIT=222
-OPEN(KUNIT,FILE=HFILENAME)
+OPEN(NEWUNIT=KUNIT,FILE=HFILENAME)
 !
 !*    3.     Reading of a data point
 !            -----------------------
 !
-  DO JK=1,KLEV
-  READ(KUNIT,*) ZLEV, PTHDF(JK),PRVF(JK)
-  END DO
+DO JK=1,KLEV
+  READ(UNIT=KUNIT,FMT=*) ZLEV, PTHDF(JK),PRVF(JK)
+END DO
 !
 !----------------------------------------------------------------------------
 
@@ -90,8 +90,7 @@ OPEN(KUNIT,FILE=HFILENAME)
 !*    8.    Closing of the data file
 !           ------------------------
 !
-!99 CLOSE(KUNIT)
- CLOSE(KUNIT)
+ CLOSE(UNIT=KUNIT)
 !
 !-------------------------------------------------------------------------------
 !
diff --git a/src/MNH/resolved_elecn.f90 b/src/MNH/resolved_elecn.f90
index 69b069d27b314456a0c7a93d19b19cc50d3a45a2..5d6bd3e89391a8c7c5f23cb86ac4d5d90d4cb3a1 100644
--- a/src/MNH/resolved_elecn.f90
+++ b/src/MNH/resolved_elecn.f90
@@ -1,6 +1,6 @@
-!MNH_LIC Copyright 2009-2018 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 2009-2019 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.
 !     ###########################
       MODULE MODI_RESOLVED_ELEC_n
@@ -166,6 +166,7 @@ END MODULE MODI_RESOLVED_ELEC_n
 !!      M. Chong      31/07/14  Add explicit LiNOx
 !!      J.Escobar : 15/09/2015 : WENO5 & JPHEXT <> 1
 !!  Philippe Wautelet: 05/2016-04/2018: new data structures and calls for I/O
+!!  Philippe Wautelet: 10/01/2019: use NEWUNIT argument of OPEN
 !!
 !-------------------------------------------------------------------------------
 !
@@ -184,7 +185,7 @@ USE MODD_CST
 USE MODD_IO_ll,            ONLY: TFILEDATA
 USE MODD_PARAMETERS, ONLY : JPVEXT
 USE MODD_ELEC_DESCR
-USE MODD_ELEC_n          
+USE MODD_ELEC_n
 USE MODD_BUDGET
 USE MODD_NSV
 USE MODD_CH_MNHC_n,    ONLY: LUSECHEM,LCH_CONV_LINOX
@@ -299,6 +300,7 @@ INTEGER :: IKU
 INTEGER :: IINFO_ll      ! return code of parallel routine
 INTEGER :: IPROC         ! my proc number
 INTEGER :: IERR          ! error status
+INTEGER :: ILU           ! unit number for IO
 !
 REAL, DIMENSION(SIZE(PZZ,1),SIZE(PZZ,2),SIZE(PZZ,3)):: ZT,   &
                                                        ZEXN, &
@@ -342,10 +344,16 @@ CHARACTER (LEN=32) :: YASCFILE
 REAL               :: ZTEMP_DIST
 CHARACTER (LEN=18) :: YNAME
 LOGICAL            :: GLMA_FILE
-TYPE(TFILEDATA),POINTER :: TZFILE
+TYPE(TFILEDATA),POINTER :: TZFILE_FGEOM_COORD
+TYPE(TFILEDATA),POINTER :: TZFILE_FGEOM_DIAG
+TYPE(TFILEDATA),POINTER :: TZFILE_LMA
+TYPE(TFILEDATA),POINTER :: TZFILE_SERIES_CLOUD_ELEC
 !
 NULLIFY(TZFIELDS_ll)
-TZFILE => NULL()
+TZFILE_FGEOM_COORD       => NULL()
+TZFILE_FGEOM_DIAG        => NULL()
+TZFILE_LMA               => NULL()
+TZFILE_SERIES_CLOUD_ELEC => NULL()
 !
 !------------------------------------------------------------------------------
 !
@@ -835,93 +843,78 @@ ENDIF
 !*      9.      OPEN THE OUTPUT ASCII FILES
 !               ---------------------------
 !
-IF (KTCOUNT .EQ. 1) THEN
+IF (KTCOUNT==1 .AND. IPROC==0) THEN
   IF (LFLASH_GEOM) THEN
     YASCFILE = CEXP//"_fgeom_diag.asc"
-    CALL IO_FILE_ADD2LIST(TZFILE,YASCFILE,'TXT','WRITE')
-    CALL IO_FILE_OPEN_ll(TZFILE,HPOSITION='APPEND',HSTATUS='NEW',KRESP=NIOSTAT_fgeom_diag)
-    NLU_fgeom_diag = TZFILE%NLU
-    IF ( IPROC .EQ. 0) THEN
-      WRITE (NLU_fgeom_diag, FMT='(A)') '--------------------------------------------------------'
-      WRITE (NLU_fgeom_diag, FMT='(A)') '*FLASH CHARACTERISTICS FROM FLASH_GEOM_ELEC*'
-      WRITE (NLU_fgeom_diag, FMT='(A)') '-- Column 1 : total flash number          --'
-      WRITE (NLU_fgeom_diag, FMT='(A)') '-- Column 2 : time (s)                    --'
-      WRITE (NLU_fgeom_diag, FMT='(A)') '-- Column 3 : cell number                 --'
-      WRITE (NLU_fgeom_diag, FMT='(A)') '-- Column 4 : flash number/cell/time step --'
-      WRITE (NLU_fgeom_diag, FMT='(A)') '-- Column 5 : flash type 1=IC, 2=CGN, 3=CGP '
-      WRITE (NLU_fgeom_diag, FMT='(A)') '-- Column 6 : number of segments          --'
-      WRITE (NLU_fgeom_diag, FMT='(A)') '-- Column 7 : trig electric field (kV/m)  --'
-      WRITE (NLU_fgeom_diag, FMT='(A)') '-- Column 8 : x coord. trig. point        --'
-      WRITE (NLU_fgeom_diag, FMT='(A)') '-- Column 9 : y coord. trig. point        --'
-      WRITE (NLU_fgeom_diag, FMT='(A)') '--         --> x,y in km if lcartesian=t, --'
-      WRITE (NLU_fgeom_diag, FMT='(A)') '--                    deg otherwise       --'
-      WRITE (NLU_fgeom_diag, FMT='(A)') '-- Column 10 : z coord. trig. point (km)  --'
-      WRITE (NLU_fgeom_diag, FMT='(A)') '-- Column 11: neutr. positive charge (C)  --' 
-      WRITE (NLU_fgeom_diag, FMT='(A)') '-- Column 12: neutr. negative charge (C)  --'
-      WRITE (NLU_fgeom_diag, FMT='(A)') '--------------------------------------------'
-    END IF
-!  
-    CALL IO_FILE_CLOSE_ll(TZFILE)
-    TZFILE => NULL()
-    CALL MPI_BCAST (NLU_fgeom_diag,1, MPI_INTEGER, 0, NMNH_COMM_WORLD, IERR)
+    CALL IO_FILE_ADD2LIST(TZFILE_FGEOM_DIAG,YASCFILE,'TXT','WRITE')
+    CALL IO_FILE_OPEN_ll(TZFILE_FGEOM_DIAG,HPOSITION='APPEND',HSTATUS='NEW')
+    ILU = TZFILE_FGEOM_DIAG%NLU
+    WRITE (UNIT=ILU, FMT='(A)') '--------------------------------------------------------'
+    WRITE (UNIT=ILU, FMT='(A)') '*FLASH CHARACTERISTICS FROM FLASH_GEOM_ELEC*'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 1 : total flash number          --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 2 : time (s)                    --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 3 : cell number                 --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 4 : flash number/cell/time step --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 5 : flash type 1=IC, 2=CGN, 3=CGP '
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 6 : number of segments          --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 7 : trig electric field (kV/m)  --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 8 : x coord. trig. point        --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 9 : y coord. trig. point        --'
+    WRITE (UNIT=ILU, FMT='(A)') '--         --> x,y in km if lcartesian=t, --'
+    WRITE (UNIT=ILU, FMT='(A)') '--                    deg otherwise       --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 10 : z coord. trig. point (km)  --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 11: neutr. positive charge (C)  --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 12: neutr. negative charge (C)  --'
+    WRITE (UNIT=ILU, FMT='(A)') '--------------------------------------------'
+    CALL FLUSH(UNIT=ILU)
 !
     IF (LSAVE_COORD) THEN
       YASCFILE = CEXP//"_fgeom_coord.asc"
-      CALL IO_FILE_ADD2LIST(TZFILE,YASCFILE,'TXT','WRITE')
-      CALL IO_FILE_OPEN_ll(TZFILE,HPOSITION='APPEND',HSTATUS='NEW',KRESP=NIOSTAT_fgeom_coord)
-      NLU_fgeom_coord = TZFILE%NLU
-      IF ( IPROC .EQ. 0) THEN
-        WRITE (NLU_fgeom_coord,FMT='(A)') '------------------------------------------'
-        WRITE (NLU_fgeom_coord,FMT='(A)') '*****FLASH COORD. FROM FLASH_GEOM_ELEC****'
-        WRITE (NLU_fgeom_coord,FMT='(A)') '-- Column 1 : flash number             --'
-        WRITE (NLU_fgeom_coord,FMT='(A)') '-- Column 2 : time (s)                 --'
-        WRITE (NLU_fgeom_coord,FMT='(A)') '-- Column 3 : type                     --'
-        WRITE (NLU_fgeom_coord,FMT='(A)') '-- Column 4 : coordinate along X (km)  --'
-        WRITE (NLU_fgeom_coord,FMT='(A)') '-- Column 5 : coordinate along Y (km)  --'
-        WRITE (NLU_fgeom_coord,FMT='(A)') '-- Column 6 : coordinate along Z (km)  --'
-        WRITE (NLU_fgeom_coord,FMT='(A)') '------------------------------------------'
-      END IF
-!
-      CALL IO_FILE_CLOSE_ll(TZFILE)
-      TZFILE => NULL()
-      CALL MPI_BCAST (NLU_fgeom_coord,1, MPI_INTEGER, 0, NMNH_COMM_WORLD, IERR)
+      CALL IO_FILE_ADD2LIST(TZFILE_FGEOM_COORD,YASCFILE,'TXT','WRITE')
+      CALL IO_FILE_OPEN_ll(TZFILE_FGEOM_COORD,HPOSITION='APPEND',HSTATUS='NEW')
+      ILU = TZFILE_FGEOM_COORD%NLU
+      WRITE (UNIT=ILU,FMT='(A)') '------------------------------------------'
+      WRITE (UNIT=ILU,FMT='(A)') '*****FLASH COORD. FROM FLASH_GEOM_ELEC****'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 1 : flash number             --'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 2 : time (s)                 --'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 3 : type                     --'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 4 : coordinate along X (km)  --'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 5 : coordinate along Y (km)  --'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 6 : coordinate along Z (km)  --'
+      WRITE (UNIT=ILU,FMT='(A)') '------------------------------------------'
+      CALL FLUSH(UNIT=ILU)
     END IF
   END IF
 !
   IF (LSERIES_ELEC) THEN
-    YASCFILE = CEXP//"_series_cloud_elec.asc"                              
-    CALL IO_FILE_ADD2LIST(TZFILE,YASCFILE,'TXT','WRITE')
-    CALL IO_FILE_OPEN_ll(TZFILE,HPOSITION='APPEND',HSTATUS='NEW',KRESP=NIOSTAT_series_cloud_elec)
-    NLU_series_cloud_elec = TZFILE%NLU
-    IF ( IPROC .EQ. 0) THEN
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '----------------------------------------------------'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '********* RESULTS FROM of LSERIES_ELEC *************'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 1 : Time (s)                            --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 2 : Cloud top height / Z > 20 dBZ (m)   --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 3 : Cloud top height / m.r. > 1.e-4 (m) --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 4 : Maximum radar reflectivity (dBZ)    --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 5 : Maximum vertical velocity (m/s)     --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 6 : Updraft volume for W > 5 m/s        --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 7 : Updraft volume for W > 10 m/s       --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 8 : Cloud water mass (kg)               --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 9 : Rain water mass (kg)                --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 10 : Ice crystal mass (kg)              --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 11 : Snow mass (kg)                     --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 12 : Graupel mass (kg)                  --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 13 : Precipitation ice mass (kg)        --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 14 : Ice mass flux product (kg2 m2/s2)  --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 15 : Precip. ice mass flux (kg m/s)     --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 16 : Non-precip. ice mass flux (kg m/s) --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 17 : Ice water path (kg/m2)             --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 18 : Cloud volume (m3)                  --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 19 : Maximum rain inst. precip. (mm/H)  --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '-- Column 20 : Rain instant. precip. (mm/H)       --'
-      WRITE (NLU_series_cloud_elec, FMT='(A)') '----------------------------------------------------'
-    END IF
-!
-    CALL IO_FILE_CLOSE_ll(TZFILE)
-    TZFILE => NULL()
-    CALL MPI_BCAST (NLU_series_cloud_elec,1, MPI_INTEGER, 0, NMNH_COMM_WORLD, IERR)
+    YASCFILE = CEXP//"_series_cloud_elec.asc"
+    CALL IO_FILE_ADD2LIST(TZFILE_SERIES_CLOUD_ELEC,YASCFILE,'TXT','WRITE')
+    CALL IO_FILE_OPEN_ll(TZFILE_SERIES_CLOUD_ELEC,HPOSITION='APPEND',HSTATUS='NEW')
+    ILU = TZFILE_SERIES_CLOUD_ELEC%NLU
+    WRITE (UNIT=ILU, FMT='(A)') '----------------------------------------------------'
+    WRITE (UNIT=ILU, FMT='(A)') '********* RESULTS FROM of LSERIES_ELEC *************'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 1 : Time (s)                            --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 2 : Cloud top height / Z > 20 dBZ (m)   --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 3 : Cloud top height / m.r. > 1.e-4 (m) --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 4 : Maximum radar reflectivity (dBZ)    --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 5 : Maximum vertical velocity (m/s)     --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 6 : Updraft volume for W > 5 m/s        --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 7 : Updraft volume for W > 10 m/s       --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 8 : Cloud water mass (kg)               --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 9 : Rain water mass (kg)                --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 10 : Ice crystal mass (kg)              --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 11 : Snow mass (kg)                     --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 12 : Graupel mass (kg)                  --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 13 : Precipitation ice mass (kg)        --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 14 : Ice mass flux product (kg2 m2/s2)  --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 15 : Precip. ice mass flux (kg m/s)     --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 16 : Non-precip. ice mass flux (kg m/s) --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 17 : Ice water path (kg/m2)             --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 18 : Cloud volume (m3)                  --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 19 : Maximum rain inst. precip. (mm/H)  --'
+    WRITE (UNIT=ILU, FMT='(A)') '-- Column 20 : Rain instant. precip. (mm/H)       --'
+    WRITE (UNIT=ILU, FMT='(A)') '----------------------------------------------------'
+    CALL FLUSH(UNIT=ILU)
   END IF
 END IF
 !
@@ -942,9 +935,9 @@ IF (LFLASH_GEOM .AND. LLMA) THEN
 !
   IF (GLMA_FILE) THEN
     IF(CLMA_FILE(1:5) /= "BEGIN") THEN ! close previous file if exists
-      CALL IO_FILE_FIND_BYNAME(CLMA_FILE,TZFILE,IERR)
-      CALL IO_FILE_CLOSE_ll(TZFILE)
-      TZFILE => NULL()
+      CALL IO_FILE_FIND_BYNAME(CLMA_FILE,TZFILE_LMA,IERR)
+      CALL IO_FILE_CLOSE_ll(TZFILE_LMA)
+      TZFILE_LMA => NULL()
     ENDIF
 !
     TDTLMA%TIME = TDTLMA%TIME - XDTLMA
@@ -955,34 +948,32 @@ IF (LFLASH_GEOM .AND. LLMA) THEN
     TDTLMA%TIME = MOD(TDTLMA%TIME + XDTLMA,86400.)
     CLMA_FILE = CEXP//"_SIMLMA_"//YNAME//".dat"
 !
-    CALL IO_FILE_ADD2LIST(TZFILE,CLMA_FILE,'TXT','WRITE')
-    CALL IO_FILE_OPEN_ll(TZFILE,HPOSITION='APPEND',HSTATUS='NEW',KRESP=ILMA_IOSTAT)
-    ILMA_UNIT = TZFILE%NLU
     IF ( IPROC .EQ. 0 ) THEN
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '----------------------------------------'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '*** FLASH COORD. FROM LMA SIMULATOR ****'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '-- Column 1  : flash number           --'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '-- Column 2  : time (s)               --'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '-- Column 3  : type                   --'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '-- Column 4  : coordinate along X (km)--'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '-- Column 5  : coordinate along Y (km)--'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '-- Column 6  : coordinate along Z (km)--'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '-- Column 7  : cld drop. mixing ratio --'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '-- Column 8  : rain mixing ratio      --'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '-- Column 9  : ice cryst mixing ratio --'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '-- Column 10 : snow mixing ratio      --'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '-- Column 11 : graupel mixing ratio   --'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '-- Column 12 : rain charge neut       --'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '-- Column 13 : ice cryst. charge neut --'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '-- Column 14 : snow charge neut       --'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '-- Column 15 : graupel charge neut    --'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '-- Column 16 : positive ions neut     --'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '-- Column 17 : negative ions neut     --'
-      WRITE (UNIT=ILMA_UNIT,FMT='(A)') '----------------------------------------'
+      CALL IO_FILE_ADD2LIST(TZFILE_LMA,CLMA_FILE,'TXT','WRITE')
+      CALL IO_FILE_OPEN_ll(TZFILE_LMA,HPOSITION='APPEND',HSTATUS='NEW')
+      ILU = TZFILE_LMA%NLU
+      WRITE (UNIT=ILU,FMT='(A)') '----------------------------------------'
+      WRITE (UNIT=ILU,FMT='(A)') '*** FLASH COORD. FROM LMA SIMULATOR ****'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 1  : flash number           --'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 2  : time (s)               --'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 3  : type                   --'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 4  : coordinate along X (km)--'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 5  : coordinate along Y (km)--'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 6  : coordinate along Z (km)--'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 7  : cld drop. mixing ratio --'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 8  : rain mixing ratio      --'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 9  : ice cryst mixing ratio --'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 10 : snow mixing ratio      --'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 11 : graupel mixing ratio   --'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 12 : rain charge neut       --'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 13 : ice cryst. charge neut --'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 14 : snow charge neut       --'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 15 : graupel charge neut    --'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 16 : positive ions neut     --'
+      WRITE (UNIT=ILU,FMT='(A)') '-- Column 17 : negative ions neut     --'
+      WRITE (UNIT=ILU,FMT='(A)') '----------------------------------------'
+      CALL FLUSH(UNIT=ILU)
     END IF
-    CALL IO_FILE_CLOSE_ll(TZFILE)
-    TZFILE => NULL()
-    CALL MPI_BCAST (ILMA_UNIT,1, MPI_INTEGER, 0, NMNH_COMM_WORLD, IERR)
   END IF
 END IF
 !
@@ -1001,27 +992,14 @@ DO JSV = NSV_ELECBEG+1, NSV_ELECEND-1
 END DO
 !
 IF ((.NOT. LOCG) .AND. LELEC_FIELD .AND.  MAX_ll(ABS(ZQTOT),IINFO_ll)>0.) THEN
-  IF (PRESENT(PSEA)) THEN
-    IF (LFLASH_GEOM) THEN
-      CALL FLASH_GEOM_ELEC_n (KTCOUNT, KMI, KRR, PTSTEP, OEXIT,       &
-                              PRHODJ, PRHODREF,                       &
-                              PRT, PCIT,                              &
-                              PSVS(:,:,:,NSV_ELECBEG:NSV_ELECEND),    &
-                              PRS, PTHT, PPABST,                      & 
-                              XEFIELDU, XEFIELDV, XEFIELDW,           &
-                              PZZ, PSVS(:,:,:,NSV_LNOXBEG), PTOWN, PSEA)
-    END IF 
-  ELSE
-    IF (LFLASH_GEOM) THEN
-      CALL FLASH_GEOM_ELEC_n (KTCOUNT, KMI, KRR, PTSTEP, OEXIT,       &
-                              PRHODJ, PRHODREF,                       &
-                              PRT, PCIT,                              &
-                              PSVS(:,:,:,NSV_ELECBEG:NSV_ELECEND),    &
-                              PRS, PTHT, PPABST,                      & 
-                              XEFIELDU, XEFIELDV, XEFIELDW,           &
-                              PZZ, PSVS(:,:,:,NSV_LNOXBEG)            )
-    END IF
-  ENDIF
+  IF (LFLASH_GEOM) THEN
+    CALL FLASH_GEOM_ELEC_n (KTCOUNT, KMI, KRR, PTSTEP, OEXIT,                                 &
+                            PRHODJ, PRHODREF, PRT, PCIT, PSVS(:,:,:,NSV_ELECBEG:NSV_ELECEND), &
+                            PRS, PTHT, PPABST, XEFIELDU, XEFIELDV, XEFIELDW,                  &
+                            PZZ, PSVS(:,:,:,NSV_LNOXBEG),                                     &
+                            TZFILE_FGEOM_DIAG, TZFILE_FGEOM_COORD, TZFILE_LMA,                &
+                            PTOWN, PSEA)
+  END IF
 !
   PSVS(:,:,:,NSV_ELECBEG) = MAX(0., PSVS(:,:,:,NSV_ELECBEG))
   PSVS(:,:,:,NSV_ELECEND) = MAX(0., PSVS(:,:,:,NSV_ELECEND))
@@ -1038,18 +1016,20 @@ IF (LSERIES_ELEC) THEN
   CALL SERIES_CLOUD_ELEC (KTCOUNT, PTSTEP,                &
                           PZZ, PRHODJ, PRHODREF, PEXNREF, &
                           PRT, PRS, PSVT,                 &
-                          PTHT, PWT, PPABST, PCIT, PINPRR  )
+                          PTHT, PWT, PPABST, PCIT,        &
+                          TZFILE_SERIES_CLOUD_ELEC,       &
+                          PINPRR                          )
 END IF
 !
 !
 !-------------------------------------------------------------------------------
 !
 !   Close Ascii Files if KTCOUNT = NSTOP
-
-IF (OEXIT) THEN
-  IF (.NOT. LFLASH_GEOM) CLOSE (UNIT=NLU_light_diag)
-  IF (.NOT. LFLASH_GEOM .AND. LSAVE_COORD) CLOSE (UNIT=NLU_light_coord)
-  IF (LLMA) CLOSE (UNIT=ILMA_UNIT)
+IF (OEXIT .AND. IPROC==0) THEN
+  IF (LFLASH_GEOM)                  CALL IO_FILE_CLOSE_ll(TZFILE_FGEOM_DIAG)
+  IF(LFLASH_GEOM .AND. LSAVE_COORD) CALL IO_FILE_CLOSE_ll(TZFILE_FGEOM_COORD)
+  IF (LSERIES_ELEC)                 CALL IO_FILE_CLOSE_ll(TZFILE_SERIES_CLOUD_ELEC)
+  IF (LFLASH_GEOM .AND. LLMA)       CALL IO_FILE_CLOSE_ll(TZFILE_LMA)
 ENDIF
 !
 !
diff --git a/src/MNH/series_cloud_elec.f90 b/src/MNH/series_cloud_elec.f90
index ffe47a93a487e5ed7405530cd73f5ce39ea0f751..2b04e77e5962a14739ddf809509f2426809651f2 100644
--- a/src/MNH/series_cloud_elec.f90
+++ b/src/MNH/series_cloud_elec.f90
@@ -1,6 +1,6 @@
-!MNH_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 2010-2019 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.
 !     #############################
       MODULE MODI_SERIES_CLOUD_ELEC
@@ -10,40 +10,46 @@ INTERFACE
   SUBROUTINE SERIES_CLOUD_ELEC (KTCOUNT, PTSTEP,                &
                                 PZZ, PRHODJ, PRHODREF, PEXNREF, &
                                 PRT, PRS, PSVT,                 &
-                                PTHT, PWT, PPABST, PCIT, PINPRR )
+                                PTHT, PWT, PPABST, PCIT,        &
+                                TPFILE_SERIES_CLOUD_ELEC,       &
+                                PINPRR                          )
 !
+USE MODD_IO_ll, ONLY: TFILEDATA
 !
-INTEGER, INTENT(IN)  :: KTCOUNT  ! Temporal loop counter
+INTEGER,                  INTENT(IN)    :: KTCOUNT  ! Temporal loop counter
 !
-REAL, INTENT(IN)  :: PTSTEP   ! Double time step except for cold start
+REAL,                     INTENT(IN)    :: PTSTEP   ! Double time step except for cold start
 !
-REAL, DIMENSION(:,:,:), INTENT(IN) :: PZZ     ! Height (z)
-REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODJ  ! Dry density * Jacobian
-REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODREF! Reference density
-REAL, DIMENSION(:,:,:), INTENT(IN) :: PEXNREF ! Reference Exner function
+REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PZZ     ! Height (z)
+REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRHODJ  ! Dry density * Jacobian
+REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRHODREF! Reference density
+REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PEXNREF ! Reference Exner function
 !
-REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PRT     ! Moist variables at time t
-REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PRS     ! Moist  variable sources
-REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PSVT    ! Scalar variable at time t
+REAL, DIMENSION(:,:,:,:), INTENT(IN)    :: PRT     ! Moist variables at time t
+REAL, DIMENSION(:,:,:,:), INTENT(IN)    :: PRS     ! Moist  variable sources
+REAL, DIMENSION(:,:,:,:), INTENT(IN)    :: PSVT    ! Scalar variable at time t
 !
-REAL, DIMENSION(:,:,:), INTENT(IN) :: PTHT   ! Theta at time t
-REAL, DIMENSION(:,:,:), INTENT(IN) :: PWT    ! Vertical velocity at t-dt
-REAL, DIMENSION(:,:,:), INTENT(IN) :: PPABST ! ab. pressure at time t
-REAL, DIMENSION(:,:,:), INTENT(IN) :: PCIT   ! Pristine ice number
-                                                 ! concentration at time t
-REAL, DIMENSION(:,:), INTENT(INOUT) :: PINPRR  ! Rain instant precip
+REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PTHT   ! Theta at time t
+REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PWT    ! Vertical velocity at t-dt
+REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PPABST ! ab. pressure at time t
+REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PCIT   ! Pristine ice number
+                                                  ! concentration at time t
+TYPE(TFILEDATA),          INTENT(IN)    :: TPFILE_SERIES_CLOUD_ELEC
+REAL, DIMENSION(:,:),     INTENT(INOUT) :: PINPRR  ! Rain instant precip
 !
 END SUBROUTINE SERIES_CLOUD_ELEC
 END INTERFACE
 END MODULE MODI_SERIES_CLOUD_ELEC
 !
 !      
-!      ###############################################################
-       SUBROUTINE SERIES_CLOUD_ELEC (KTCOUNT, PTSTEP,                &
-                                     PZZ, PRHODJ, PRHODREF, PEXNREF, &
-                                     PRT, PRS, PSVT,                 &
-                                     PTHT, PWT, PPABST, PCIT, PINPRR )
-!      ###############################################################
+! ###############################################################
+  SUBROUTINE SERIES_CLOUD_ELEC (KTCOUNT, PTSTEP,                &
+                                PZZ, PRHODJ, PRHODREF, PEXNREF, &
+                                PRT, PRS, PSVT,                 &
+                                PTHT, PWT, PPABST, PCIT,        &
+                                TPFILE_SERIES_CLOUD_ELEC,       &
+                                PINPRR                          )
+! ###############################################################
 !
 !!****  * -
 !!
@@ -72,6 +78,7 @@ END MODULE MODI_SERIES_CLOUD_ELEC
 !!      Modifications:
 !!      C. Barthe  * LACy *  Dec. 2010    add some parameters
 !!      J.Escobar : 15/09/2015 : WENO5 & JPHEXT <> 1
+!!      Philippe Wautelet: 10/01/2019: use NEWUNIT argument of OPEN
 !!
 !-------------------------------------------------------------------------------
 !
@@ -80,6 +87,7 @@ END MODULE MODI_SERIES_CLOUD_ELEC
 !
 USE MODD_CONF, ONLY : CEXP
 USE MODD_CST
+USE MODD_IO_ll, ONLY: TFILEDATA
 USE MODD_REF
 USE MODD_PARAMETERS
 USE MODD_ELEC_DESCR
@@ -91,7 +99,6 @@ USE MODD_DYN_n, ONLY : XDXHATM, XDYHATM
 USE MODD_RAIN_ICE_DESCR
 USE MODD_RAIN_ICE_PARAM
 USE MODD_NSV, ONLY : NSV_ELECBEG, NSV_ELECEND
-USE MODD_PRINT_ELEC,  ONLY : NLU_series_cloud_elec, NIOSTAT_series_cloud_elec
 !
 USE MODI_MOMG
 USE MODI_RADAR_RAIN_ICE
@@ -103,25 +110,26 @@ IMPLICIT NONE
 !
 !*       0.1   Declarations of dummy arguments :
 !
-INTEGER, INTENT(IN) :: KTCOUNT  ! Temporal loop counter
+INTEGER,                  INTENT(IN)    :: KTCOUNT  ! Temporal loop counter
 !
-REAL, INTENT(IN) :: PTSTEP   ! Double time step except for
-                             ! cold start
+REAL,                     INTENT(IN)    :: PTSTEP   ! Double time step except for cold start
 !
-REAL, DIMENSION(:,:,:), INTENT(IN) :: PZZ      ! Height (z)
-REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODJ   ! Dry density * Jacobian
-REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODREF ! Reference density
-REAL, DIMENSION(:,:,:), INTENT(IN) :: PEXNREF  ! Reference Exner function
-REAL, DIMENSION(:,:,:), INTENT(IN) :: PTHT     ! Theta at time t
-REAL, DIMENSION(:,:,:), INTENT(IN) :: PWT      ! Vertical velocity at t
-REAL, DIMENSION(:,:,:), INTENT(IN) :: PPABST   ! abs. pressure at time t
-REAL, DIMENSION(:,:,:), INTENT(IN) :: PCIT     ! Pristine ice number
-                                               ! concentration at time t
-REAL, DIMENSION(:,:), INTENT(INOUT) :: PINPRR  ! Rain instant precip
+REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PZZ     ! Height (z)
+REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRHODJ  ! Dry density * Jacobian
+REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRHODREF! Reference density
+REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PEXNREF ! Reference Exner function
 !
-REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PRT    ! Moist variables at time t
-REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PRS    ! Moist  variable sources
-REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PSVT   ! Scalar variable at time t
+REAL, DIMENSION(:,:,:,:), INTENT(IN)    :: PRT     ! Moist variables at time t
+REAL, DIMENSION(:,:,:,:), INTENT(IN)    :: PRS     ! Moist  variable sources
+REAL, DIMENSION(:,:,:,:), INTENT(IN)    :: PSVT    ! Scalar variable at time t
+!
+REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PTHT   ! Theta at time t
+REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PWT    ! Vertical velocity at t-dt
+REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PPABST ! ab. pressure at time t
+REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PCIT   ! Pristine ice number
+                                                  ! concentration at time t
+TYPE(TFILEDATA),          INTENT(IN)    :: TPFILE_SERIES_CLOUD_ELEC
+REAL, DIMENSION(:,:),     INTENT(INOUT) :: PINPRR  ! Rain instant precip
 !
 !
 !*       0.2   Declarations of local variables :
@@ -135,6 +143,7 @@ INTEGER :: ICOUNT     ! counter for iwp computation
 INTEGER :: IPROC      ! my proc number
 INTEGER :: IPROC_MAX  ! proc that contains max value
 INTEGER :: IINFO_ll   ! return code of parallel routine
+INTEGER :: ILU        ! unit number for IO
 !
 INTEGER, SAVE :: JCOUNT
 !
@@ -532,9 +541,8 @@ IF (JCOUNT == JCOUNT_STOP) THEN
   CALL REDUCESUM_ll (ZINPRR,        IINFO_ll)
 !
   IF (IPROC == 0) THEN
-    OPEN (UNIT=NLU_series_cloud_elec, FILE=CEXP//"_series_cloud_elec.asc", ACTION="WRITE",   &
-                 STATUS="OLD", FORM="FORMATTED", POSITION="APPEND")
-    WRITE (NLU_series_cloud_elec, FMT='(I6,19(E12.4))') &
+    ILU = TPFILE_SERIES_CLOUD_ELEC%NLU
+    WRITE (ILU, FMT='(I6,19(E12.4))') &
              INT(KTCOUNT*PTSTEP),         & ! time
              ZCTH_REF/FLOAT(JCOUNT),      & ! cloud top height from Z
              ZCTH_MR/FLOAT(JCOUNT),       & ! cloud top height from m.r.
@@ -555,9 +563,9 @@ IF (JCOUNT == JCOUNT_STOP) THEN
              ZCLD_VOL/FLOAT(JCOUNT),      & ! cloud volume
              ZINPRR/FLOAT(JCOUNT),        & ! Rain instant precip
              ZMAX_INPRR/FLOAT(JCOUNT)       ! maximum rain instant. precip.
-    CLOSE (UNIT=NLU_series_cloud_elec)
+    CALL FLUSH(UNIT=ILU)
   END IF
-! 
+!
   JCOUNT = 0
   ZMASS_C     = 0.
   ZMASS_R     = 0.
diff --git a/src/MNH/spectre_arome.f90 b/src/MNH/spectre_arome.f90
index 0ebcff24209aa4787743477193686ece5c1a1106..532c3d9794e0859740b4b6115e3941c38de5d077 100644
--- a/src/MNH/spectre_arome.f90
+++ b/src/MNH/spectre_arome.f90
@@ -1,6 +1,6 @@
-!MNH_LIC Copyright 1994-2018 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2019 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.
 !     ####################
       MODULE MODI_SPECTRE_AROME
@@ -25,6 +25,7 @@ SUBROUTINE SPECTRE_AROME(HINIFILE,HOUTFILE,PDELTAX,PDELTAY,KI,KJ,KK)
 !     
 ! Modifications:
 !  Philippe Wautelet: 05/2016-04/2018: new data structures and calls for I/O
+!  Philippe Wautelet: 10/01/2019: use NEWUNIT argument of OPEN
 !
 USE MODD_CONF
 USE MODE_FM
@@ -46,6 +47,7 @@ REAL,DIMENSION(:,:,:),ALLOCATABLE:: ZWORK         ! work array
 REAL,DIMENSION(:,:,:),ALLOCATABLE:: ZWORK1         ! work array
 CHARACTER (LEN=32)                :: YFILE         ! file to open
 REAL :: ZLATB, ZLONB
+INTEGER :: ILU
 INTEGER :: JJJ,III,JERR
 !
 !
@@ -59,13 +61,13 @@ IF (LSPECTRE_U) THEN
   YFIELDSP="U"  
   ZWORK1(:,:,:)=999
   WRITE(YFILE, FMT='(A,A,A)') TRIM(ADJUSTL(HINIFILE)),"_",TRIM(ADJUSTL(YFIELDSP))
-  OPEN(UNIT=42,FILE=YFILE,ACCESS='SEQUENTIAL', IOSTAT=JERR)
+  OPEN(NEWUNIT=ILU,FILE=YFILE,ACCESS='SEQUENTIAL', IOSTAT=JERR)
   DO JJJ=1+JPHEXT,SIZE(ZWORK1,2)-JPHEXT
     DO III=1+JPHEXT,SIZE(ZWORK1,1)-JPHEXT
-       READ(UNIT=42,IOSTAT=JERR, FMT=*) ZLATB, ZLONB, ZWORK1(III,JJJ,1+JPVEXT:SIZE(ZWORK1,3)-JPVEXT)
+       READ(UNIT=ILU,IOSTAT=JERR, FMT=*) ZLATB, ZLONB, ZWORK1(III,JJJ,1+JPVEXT:SIZE(ZWORK1,3)-JPVEXT)
     END DO
   END DO
-  CLOSE(UNIT=42, IOSTAT=JERR)
+  CLOSE(UNIT=ILU, IOSTAT=JERR)
 !
    IF (LZOOM) THEN
      ALLOCATE(ZWORK(NITOT+2,NJTOT+2,SIZE(ZWORK1,3)))
@@ -83,13 +85,13 @@ IF (LSPECTRE_V) THEN
   YFIELDSP="V"  
   ZWORK1(:,:,:)=999.
   WRITE(YFILE, FMT='(A,A,A)') TRIM(ADJUSTL(HINIFILE)),"_",TRIM(ADJUSTL(YFIELDSP))
-  OPEN(UNIT=42,FILE=YFILE,ACCESS='SEQUENTIAL', IOSTAT=JERR)
+  OPEN(NEWUNIT=ILU,FILE=YFILE,ACCESS='SEQUENTIAL', IOSTAT=JERR)
   DO JJJ=1+JPHEXT,SIZE(ZWORK1,2)-JPHEXT
     DO III=1+JPHEXT,SIZE(ZWORK1,1)-JPHEXT
-       READ(UNIT=42,IOSTAT=JERR, FMT=*) ZLATB, ZLONB, ZWORK1(III,JJJ,1+JPVEXT:SIZE(ZWORK1,3)-JPVEXT)
+       READ(UNIT=ILU,IOSTAT=JERR, FMT=*) ZLATB, ZLONB, ZWORK1(III,JJJ,1+JPVEXT:SIZE(ZWORK1,3)-JPVEXT)
     END DO
   END DO
-  CLOSE(UNIT=42, IOSTAT=JERR)
+  CLOSE(UNIT=ILU, IOSTAT=JERR)
 !
    IF (LZOOM) THEN
      ALLOCATE(ZWORK(NITOT+2,NJTOT+2,SIZE(ZWORK1,3)))
@@ -107,13 +109,13 @@ IF (LSPECTRE_W) THEN
   YFIELDSP="W"  
   ZWORK1(:,:,:)=999.
   WRITE(YFILE, FMT='(A,A,A)') TRIM(ADJUSTL(HINIFILE)),"_",TRIM(ADJUSTL(YFIELDSP))
-  OPEN(UNIT=42,FILE=YFILE,ACCESS='SEQUENTIAL', IOSTAT=JERR)
+  OPEN(NEWUNIT=ILU,FILE=YFILE,ACCESS='SEQUENTIAL', IOSTAT=JERR)
   DO JJJ=1+JPHEXT,SIZE(ZWORK1,2)-JPHEXT
     DO III=1+JPHEXT,SIZE(ZWORK1,1)-JPHEXT
-       READ(UNIT=42,IOSTAT=JERR, FMT=*) ZLATB, ZLONB, ZWORK1(III,JJJ,1+JPVEXT:SIZE(ZWORK1,3)-JPVEXT)
+       READ(UNIT=ILU,IOSTAT=JERR, FMT=*) ZLATB, ZLONB, ZWORK1(III,JJJ,1+JPVEXT:SIZE(ZWORK1,3)-JPVEXT)
     END DO
   END DO
-  CLOSE(UNIT=42, IOSTAT=JERR)
+  CLOSE(UNIT=ILU, IOSTAT=JERR)
 !
    IF (LZOOM) THEN
      ALLOCATE(ZWORK(NITOT+2,NJTOT+2,SIZE(ZWORK1,3)))
@@ -131,13 +133,13 @@ IF (LSPECTRE_RV) THEN
   YFIELDSP="Rv"  
   ZWORK1(:,:,:)=999.
   WRITE(YFILE, FMT='(A,A,A)') TRIM(ADJUSTL(HINIFILE)),"_",TRIM(ADJUSTL(YFIELDSP))
-  OPEN(UNIT=42,FILE=YFILE,ACCESS='SEQUENTIAL', IOSTAT=JERR)
+  OPEN(NEWUNIT=ILU,FILE=YFILE,ACCESS='SEQUENTIAL', IOSTAT=JERR)
   DO JJJ=1+JPHEXT,SIZE(ZWORK1,2)-JPHEXT
     DO III=1+JPHEXT,SIZE(ZWORK1,1)-JPHEXT
-       READ(UNIT=42,IOSTAT=JERR, FMT=*) ZLATB, ZLONB, ZWORK1(III,JJJ,1+JPVEXT:SIZE(ZWORK1,3)-JPVEXT)
+       READ(UNIT=ILU,IOSTAT=JERR, FMT=*) ZLATB, ZLONB, ZWORK1(III,JJJ,1+JPVEXT:SIZE(ZWORK1,3)-JPVEXT)
     END DO
   END DO
-  CLOSE(UNIT=42, IOSTAT=JERR)
+  CLOSE(UNIT=ILU, IOSTAT=JERR)
 !
    IF (LZOOM) THEN
      ALLOCATE(ZWORK(NITOT+2,NJTOT+2,SIZE(ZWORK1,3)))
@@ -155,13 +157,13 @@ IF (LSPECTRE_TH) THEN
   YFIELDSP="Theta"  
   ZWORK1(:,:,:)=999.
   WRITE(YFILE, FMT='(A,A,A)') TRIM(ADJUSTL(HINIFILE)),"_",TRIM(ADJUSTL(YFIELDSP))
-  OPEN(UNIT=42,FILE=YFILE,ACCESS='SEQUENTIAL', IOSTAT=JERR)
+  OPEN(NEWUNIT=ILU,FILE=YFILE,ACCESS='SEQUENTIAL', IOSTAT=JERR)
   DO JJJ=1+JPHEXT,SIZE(ZWORK1,2)-JPHEXT
     DO III=1+JPHEXT,SIZE(ZWORK1,1)-JPHEXT
-       READ(UNIT=42,IOSTAT=JERR, FMT=*) ZLATB, ZLONB, ZWORK1(III,JJJ,1+JPVEXT:SIZE(ZWORK1,3)-JPVEXT)
+       READ(UNIT=ILU,IOSTAT=JERR, FMT=*) ZLATB, ZLONB, ZWORK1(III,JJJ,1+JPVEXT:SIZE(ZWORK1,3)-JPVEXT)
     END DO
   END DO
-  CLOSE(UNIT=42, IOSTAT=JERR)
+  CLOSE(UNIT=ILU, IOSTAT=JERR)
 !
    IF (LZOOM) THEN
      ALLOCATE(ZWORK(NITOT+2,NJTOT+2,SIZE(ZWORK1,3)))