From c9699eb32980b843273a17fff005fc25da609dfe Mon Sep 17 00:00:00 2001
From: Philippe WAUTELET <philippe.wautelet@aero.obs-mip.fr>
Date: Fri, 12 May 2017 15:18:36 +0200
Subject: [PATCH] Philippe 12/05/2017: IO: use IO_WRITE_FIELD to write several
 fields

---
 src/LIB/SURCOUCHE/src/mode_field.f90 |   28 +
 src/MNH/write_lfifm1_for_diag.f90    | 2317 ++++++++++++++------------
 src/MNH/write_lfin.f90               |   25 +-
 3 files changed, 1314 insertions(+), 1056 deletions(-)

diff --git a/src/LIB/SURCOUCHE/src/mode_field.f90 b/src/LIB/SURCOUCHE/src/mode_field.f90
index d94d611b5..667862b5a 100644
--- a/src/LIB/SURCOUCHE/src/mode_field.f90
+++ b/src/LIB/SURCOUCHE/src/mode_field.f90
@@ -2640,6 +2640,34 @@ ALLOCATE(TFIELDLIST(IDX)%TFIELD_X2D(IMODEL))
 IDX = IDX+1
 !
 IF(IDX>MAXFIELDS) CALL ERR_INI_FIELD_LIST()
+TFIELDLIST(IDX)%CMNHNAME   = 'INPRT'
+TFIELDLIST(IDX)%CSTDNAME   = ''
+TFIELDLIST(IDX)%CLONGNAME  = 'MesoNH: INPRT'
+TFIELDLIST(IDX)%CUNITS     = 'm s-1'
+TFIELDLIST(IDX)%CDIR       = 'XY'
+TFIELDLIST(IDX)%CCOMMENT   = 'X_Y_Total INstantaneaous PRecipitation rate'
+TFIELDLIST(IDX)%NGRID      = 1
+TFIELDLIST(IDX)%NTYPE      = TYPEREAL
+TFIELDLIST(IDX)%NDIMS      = 2
+!No permanent variable associated to this field
+!ALLOCATE(TFIELDLIST(IDX)%TFIELD_X2D(IMODEL))
+IDX = IDX+1
+!
+IF(IDX>MAXFIELDS) CALL ERR_INI_FIELD_LIST()
+TFIELDLIST(IDX)%CMNHNAME   = 'ACPRT'
+TFIELDLIST(IDX)%CSTDNAME   = ''
+TFIELDLIST(IDX)%CLONGNAME  = 'MesoNH: ACPRT'
+TFIELDLIST(IDX)%CUNITS     = 'm'
+TFIELDLIST(IDX)%CDIR       = 'XY'
+TFIELDLIST(IDX)%CCOMMENT   = 'X_Y_Total ACcumulated PRecipitation rate'
+TFIELDLIST(IDX)%NGRID      = 1
+TFIELDLIST(IDX)%NTYPE      = TYPEREAL
+TFIELDLIST(IDX)%NDIMS      = 2
+!No permanent variable associated to this field
+!ALLOCATE(TFIELDLIST(IDX)%TFIELD_X2D(IMODEL))
+IDX = IDX+1
+!
+IF(IDX>MAXFIELDS) CALL ERR_INI_FIELD_LIST()
 TFIELDLIST(IDX)%CMNHNAME   = 'VT_FLX'
 TFIELDLIST(IDX)%CSTDNAME   = ''
 TFIELDLIST(IDX)%CLONGNAME  = 'MesoNH: VT_FLX'
diff --git a/src/MNH/write_lfifm1_for_diag.f90 b/src/MNH/write_lfifm1_for_diag.f90
index e596e458b..2e72cd5c7 100644
--- a/src/MNH/write_lfifm1_for_diag.f90
+++ b/src/MNH/write_lfifm1_for_diag.f90
@@ -251,13 +251,8 @@ CHARACTER(LEN=28), INTENT(IN) :: HDADFILE     ! corresponding FM-file name of
 !*       0.2   Declarations of local variables
 !
 INTEGER           :: IRESP          ! return-code for the file routines 
-INTEGER           :: IGRID          ! IGRID : grid indicator
-INTEGER           :: ILENCH         ! ILENCH : length of comment string 
 !
 CHARACTER(LEN=28) :: YFMFILE        ! Temporary variable to store FM-file name
-CHARACTER(LEN=28) :: YFMFILE2       ! Temporary variable to store FM-file name
-CHARACTER(LEN=16) :: YRECFM         ! Name of the article to be written
-CHARACTER(LEN=100):: YCOMMENT       ! Comment string
 !
 CHARACTER(LEN=3)  :: YFRC           ! to mark the time of the forcing
 CHARACTER(LEN=31) :: YFGRI          ! file name for GPS stations
@@ -334,8 +329,6 @@ REAL,DIMENSION(:,:,:,:), ALLOCATABLE :: ZTMP1, ZTMP2, ZTMP3, ZTMP4
 !
 !*       0.     ARRAYS BOUNDS INITIALIZATION
 !
-YFMFILE2 = TPFILE%CNAME
-!
 CALL GET_DIM_EXT_ll ('B',IIU,IJU)
 CALL GET_INDICE_ll (IIB,IJB,IIE,IJE)
 IKU=NKMAX+2*JPVEXT
@@ -349,7 +342,7 @@ CALL FMLOOK_ll(CLUOUT0,CLUOUT0,ILUOUT0,IRESP)
 !*       1.     WRITES IN THE LFI FILE
 !               ---------------------- 
 !
-!*       1.0    YFMFILE2 and HDADFILE :
+!*       1.0    TPFILE%CNAME and HDADFILE :
 !
 CALL IO_WRITE_FIELD(TPFILE,'MASDEV',   CLUOUT,IRESP,NMASDEV)
 CALL IO_WRITE_FIELD(TPFILE,'BUGFIX',   CLUOUT,IRESP,NBUGFIX)
@@ -558,32 +551,15 @@ ZTEMP(:,:,:)=XTHT(:,:,:)*(XPABST(:,:,:)/ XP00) **(XRD/XCPD)
 !
 IF (LVAR_TURB) THEN
   IF (CTURB /= 'NONE') THEN
-    YRECFM='TKET'
-    YCOMMENT='X_Y_Z_Turbulent Kinetic Energy (M**2/S**2)'
-    IGRID=1
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',XTKET,IGRID,ILENCH,YCOMMENT,IRESP)
+    CALL IO_WRITE_FIELD(TPFILE,'TKET',CLUOUT,IRESP,XTKET)
     !
     IF( NRR > 1 ) THEN
-      YRECFM='SRCT'
-      YCOMMENT='X_Y_Z_normalized 2nd_order moment s_r_c/2Sigma_s2 (KG/KG**2)'
-      IGRID=1
-      ILENCH=LEN(YCOMMENT)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',XSRCT,IGRID,ILENCH,YCOMMENT,IRESP)
-      !
-      YRECFM='SIGS'
-      YCOMMENT='X_Y_Z_Sigma_s from turbulence scheme (KG/KG**2)'
-      IGRID=1
-      ILENCH=LEN(YCOMMENT)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',XSIGS,IGRID,ILENCH,YCOMMENT,IRESP)
+      CALL IO_WRITE_FIELD(TPFILE,'SRCT',CLUOUT,IRESP,XSRCT)
+      CALL IO_WRITE_FIELD(TPFILE,'SIGS',CLUOUT,IRESP,XSIGS)
     END IF
     ! 
     IF(CTOM=='TM06') THEN
-      YRECFM='BL_DEPTH'
-      YCOMMENT='X_Y_BL_DEPTH (M)'
-      IGRID=1
-      ILENCH=LEN(YCOMMENT)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',XBL_DEPTH,IGRID,ILENCH,YCOMMENT,IRESP)
+      CALL IO_WRITE_FIELD(TPFILE,'BL_DEPTH',CLUOUT,IRESP,XBL_DEPTH)
     END IF
   END IF
 END IF
@@ -594,125 +570,86 @@ IF (LVAR_PR .AND. LUSERR .AND. SIZE(XINPRR)>0 ) THEN
   !
   ! explicit species
   !
-  ZWORK21(:,:) = XINPRR(:,:)*3.6E6
-  YRECFM      ='INPRR'
-  YCOMMENT    ='X_Y_INstantaneous PRecipitation Rate (MM/H)'
-  IGRID       =1
-  ILENCH      =LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
-  !
-  ZWORK31(:,:,:)  = XINPRR3D(:,:,:)
-  YRECFM      = 'INPRR3D'
-  YCOMMENT    = 'X_Y_INstantaneous 3D Rain Precipitation flux (M/S)'
-  IGRID       = 1
-  ILENCH      = LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH, &
-                                             YCOMMENT,IRESP) ! unit conversion
-!
-  ZWORK31(:,:,:)  = XEVAP3D(:,:,:)
-  YRECFM      = 'EVAP3D'
-  YCOMMENT    = 'X_Y_INstantaneous 3D Rain Evaporation flux (KG/KG/S)'
-  IGRID       = 1
-  ILENCH      = LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH, &
-                                             YCOMMENT,IRESP) ! unit conversion
-!
-  ZWORK21(:,:) = XACPRR(:,:)*1.E3
-  YRECFM      ='ACPRR'
-  YCOMMENT    ='X_Y_ACcumulated PRecipitation Rate (MM)'
-  IGRID       =1
-  ILENCH      =LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
+  CALL FIND_FIELD_ID_FROM_MNHNAME('INPRR',IID,IRESP)
+  TZFIELD = TFIELDLIST(IID)
+  TZFIELD%CUNITS = 'mm hour-1'
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XINPRR*3.6E6)
+  !
+  CALL IO_WRITE_FIELD(TPFILE,'INPRR3D',CLUOUT,IRESP,XINPRR3D)
+  CALL IO_WRITE_FIELD(TPFILE,'EVAP3D', CLUOUT,IRESP,XEVAP3D)
+  !
+  CALL FIND_FIELD_ID_FROM_MNHNAME('ACPRR',IID,IRESP)
+  TZFIELD = TFIELDLIST(IID)
+  TZFIELD%CUNITS = 'mm'
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XACPRR*1.0E3)
   !
   IF (CCLOUD(1:3) == 'ICE' .OR. CCLOUD == 'C2R2' .OR. CCLOUD == 'C3R5' .OR.&
       CCLOUD == 'KHKO' .OR. CCLOUD == 'LIMA') THEN 
     IF (SIZE(XINPRC) /= 0 ) THEN
-      ZWORK21(:,:) = XINPRC(:,:)*3.6E6
-      YRECFM      ='INPRC'
-      YCOMMENT    ='X_Y_INstantaneous Cloud PRecipitation Rate (MM/H)'
-      IGRID       =1
-      ILENCH      =LEN(YCOMMENT)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
-  !
-      ZWORK21(:,:) = XACPRC(:,:)*1.E3
-      YRECFM      ='ACPRC'
-      YCOMMENT    ='X_Y_ACcumulated Cloud PRecipitation Rate (MM)'
-      IGRID       =1
-      ILENCH      =LEN(YCOMMENT)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
+      CALL FIND_FIELD_ID_FROM_MNHNAME('INPRC',IID,IRESP)
+      TZFIELD = TFIELDLIST(IID)
+      TZFIELD%CUNITS = 'mm hour-1'
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XINPRC*3.6E6)
+      !
+      CALL FIND_FIELD_ID_FROM_MNHNAME('ACPRC',IID,IRESP)
+      TZFIELD = TFIELDLIST(IID)
+      TZFIELD%CUNITS = 'mm'
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XACPRC*1.0E3)
     END IF 
     IF (SIZE(XINDEP) /= 0 ) THEN
-      ZWORK21(:,:) = XINDEP(:,:)*3.6E6
-      YRECFM      ='INDEP'
-      YCOMMENT    ='X_Y_INstantaneous Cloud Deposition Rate (MM/H)'
-      IGRID       =1
-      ILENCH      =LEN(YCOMMENT)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
-  !
-      ZWORK21(:,:) = XACDEP(:,:)*1.E3
-      YRECFM      ='ACDEP'
-      YCOMMENT    ='X_Y_ACcumulated Cloud Deposition Rate (MM)'
-      IGRID       =1
-      ILENCH      =LEN(YCOMMENT)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
+      CALL FIND_FIELD_ID_FROM_MNHNAME('INDEP',IID,IRESP)
+      TZFIELD = TFIELDLIST(IID)
+      TZFIELD%CUNITS = 'mm hour-1'
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XINDEP*3.6E6)
+      !
+      CALL FIND_FIELD_ID_FROM_MNHNAME('ACDEP',IID,IRESP)
+      TZFIELD = TFIELDLIST(IID)
+      TZFIELD%CUNITS = 'mm'
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XACDEP*1.0E3)
     END IF 
   END IF 
   IF (CCLOUD(1:3) == 'ICE' .OR. CCLOUD == 'C3R5' .OR. CCLOUD == 'LIMA') THEN
-    ZWORK21(:,:)  = XINPRS(:,:)*3.6E6
-    YRECFM      = 'INPRS'
-    YCOMMENT    = 'X_Y_INstantaneaous PRecipitation rate for Snow (MM/H)'
-    IGRID       = 1
-    ILENCH      = LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
-  !
-    ZWORK21(:,:)  = XACPRS(:,:)*1.0E3
-    YRECFM      = 'ACPRS'
-    YCOMMENT    = 'X_Y_ACccumuated PRecipitation rate for Snow (MM)'
-    IGRID       = 1
-    ILENCH      = LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
-  !
-    ZWORK21(:,:)  = XINPRG(:,:)*3.6E6
-    YRECFM      = 'INPRG'
-    YCOMMENT    = 'X_Y_INstantaneaous PRecipitation rate for Graupel (MM/H)'
-    IGRID       = 1
-    ILENCH      = LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
-  !
-    ZWORK21(:,:)  = XACPRG(:,:)*1.0E3
-    YRECFM      = 'ACPRG'
-    YCOMMENT    = 'X_Y_ACccumuated PRecipitation rate for Graupel (MM)'
-    IGRID       = 1
-    ILENCH      = LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
+    CALL FIND_FIELD_ID_FROM_MNHNAME('INPRS',IID,IRESP)
+    TZFIELD = TFIELDLIST(IID)
+    TZFIELD%CUNITS = 'mm hour-1'
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XINPRS*3.6E6)
+    !
+    CALL FIND_FIELD_ID_FROM_MNHNAME('ACPRS',IID,IRESP)
+    TZFIELD = TFIELDLIST(IID)
+    TZFIELD%CUNITS = 'mm'
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XACPRS*1.0E3)
+    !
+    CALL FIND_FIELD_ID_FROM_MNHNAME('INPRG',IID,IRESP)
+    TZFIELD = TFIELDLIST(IID)
+    TZFIELD%CUNITS = 'mm hour-1'
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XINPRG*3.6E6)
+    !
+    CALL FIND_FIELD_ID_FROM_MNHNAME('ACPRG',IID,IRESP)
+    TZFIELD = TFIELDLIST(IID)
+    TZFIELD%CUNITS = 'mm'
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XACPRG*1.0E3)
   !
     IF (SIZE(XINPRH) /= 0 ) THEN
-      ZWORK21(:,:)  = XINPRH(:,:)*3.6E6
-      YRECFM      = 'INPRH'
-      YCOMMENT    = 'X_Y_INstantaneaous PRecipitation rate for Hail (MM/H)'
-      IGRID       = 1
-      ILENCH      = LEN(YCOMMENT)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
-    !
-      ZWORK21(:,:)  = XACPRH(:,:)*1.0E3
-      YRECFM      = 'ACPRH'
-      YCOMMENT    = 'X_Y_ACccumuated PRecipitation rate for Hail (MM)'
-      IGRID       = 1
-      ILENCH      = LEN(YCOMMENT)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
+      CALL FIND_FIELD_ID_FROM_MNHNAME('INPRH',IID,IRESP)
+      TZFIELD = TFIELDLIST(IID)
+      TZFIELD%CUNITS = 'mm hour-1'
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XINPRH*3.6E6)
+      !
+      CALL FIND_FIELD_ID_FROM_MNHNAME('ACPRH',IID,IRESP)
+      TZFIELD = TFIELDLIST(IID)
+      TZFIELD%CUNITS = 'mm'
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XACPRH*1.0E3)
     ENDIF
   !
-    ZWORK21(:,:) = (XINPRR(:,:) + XINPRS(:,:) + XINPRG(:,:))*3.6E6
+    ZWORK21(:,:) = XINPRR(:,:) + XINPRS(:,:) + XINPRG(:,:)
     IF (SIZE(XINPRC) /= 0 ) &     
-      ZWORK21(:,:) = ZWORK21(:,:) + XINPRC(:,:)*3.6E6          
+      ZWORK21(:,:) = ZWORK21(:,:) + XINPRC(:,:)
     IF (SIZE(XINPRH) /= 0 ) &       
-      ZWORK21(:,:) = ZWORK21(:,:) + XINPRH(:,:)*3.6E6
-    YRECFM      = 'INPRT'
-    YCOMMENT    = 'X_Y_Total INstantaneaous PRecipitation rate (MM/H)'
-    YCOMMENT    = 'X_Y_INPRT (MM/H)'
-    IGRID       = 1
-    ILENCH      = LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
+      ZWORK21(:,:) = ZWORK21(:,:) + XINPRH(:,:)
+    CALL FIND_FIELD_ID_FROM_MNHNAME('INPRT',IID,IRESP)
+    TZFIELD = TFIELDLIST(IID)
+    TZFIELD%CUNITS = 'mm hour-1'
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK21*3.6E6)
   !
     ZWORK21(:,:) = (XACPRR(:,:) + XACPRS(:,:) + XACPRG(:,:))*1.0E3
     IF (SIZE(XINPRC) /= 0 ) &      
@@ -720,38 +657,30 @@ IF (LVAR_PR .AND. LUSERR .AND. SIZE(XINPRR)>0 ) THEN
     IF (SIZE(XINPRH) /= 0 ) &        
       ZWORK21(:,:) = ZWORK21(:,:) + XACPRH(:,:)*1.0E3
   !
-    YRECFM      = 'ACPRT'
-    YCOMMENT    = 'X_Y_Total ACcumulated PRecipitation rate (MM)'
-    YCOMMENT    = 'X_Y_ACPRT (MM)'
-    IGRID       = 1
-    ILENCH      = LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
+    CALL FIND_FIELD_ID_FROM_MNHNAME('ACPRT',IID,IRESP)
+    TZFIELD = TFIELDLIST(IID)
+    TZFIELD%CUNITS = 'mm'
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK21*1.0E3)
   !
   END IF
   !
   !* Convective rain
   !
   IF (CDCONV /= 'NONE') THEN
-    ZWORK21(:,:) = XPRCONV(:,:)*3.6E6 
-    YRECFM      ='PRCONV'
-    YCOMMENT    ='X_Y_CONVective instantaneous Precipitation Rate (MM/H)'
-    IGRID       =1
-    ILENCH      =LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
-  !
-    ZWORK21(:,:) = XPACCONV(:,:)*1.E3 
-    YRECFM      ='PACCONV'
-    YCOMMENT    ='X_Y_CONVective ACcumulated Precipitation rate (MM)'
-    IGRID       =1
-    ILENCH      =LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
-  !
-    ZWORK21(:,:) = XPRSCONV(:,:)*3.6E6 
-    YRECFM      ='PRSCONV'
-    YCOMMENT    ='X_Y_CONVective instantaneous Precipitation Rate for Snow (MM/H)'
-    IGRID       = 1
-    ILENCH      = LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
+    CALL FIND_FIELD_ID_FROM_MNHNAME('PRCONV',IID,IRESP)
+    TZFIELD = TFIELDLIST(IID)
+    TZFIELD%CUNITS = 'mm hour-1'
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XPRCONV*3.6E6)
+    !
+    CALL FIND_FIELD_ID_FROM_MNHNAME('PACCONV',IID,IRESP)
+    TZFIELD = TFIELDLIST(IID)
+    TZFIELD%CUNITS = 'mm'
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XPACCONV*1.0E3)
+    !
+    CALL FIND_FIELD_ID_FROM_MNHNAME('PRSCONV',IID,IRESP)
+    TZFIELD = TFIELDLIST(IID)
+    TZFIELD%CUNITS = 'mm hour-1'
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XPRSCONV*3.6E6)
   END IF
 END IF
 IF (LVAR_PR ) THEN
@@ -779,11 +708,16 @@ IF (LVAR_PR ) THEN
     ZWORK23(:,:) = 0.
   END DO
   !* Precipitable water in kg/m**2
-  YRECFM='PRECIP_WAT'
-  YCOMMENT='(kg/m²)'
-  IGRID=1
-  ILENCH=LEN(YCOMMENT)  
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK22,IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'PRECIP_WAT'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: PRECIP_WAT'
+  TZFIELD%CUNITS     = 'kg m-2'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = ''
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 2
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK22)
 ENDIF
 !
 !
@@ -852,84 +786,141 @@ IF (LHU_FLX) THEN
   ENDIF
   ! Ecriture
   !  composantes U et V du flux surfacique d'humidité  
-  YRECFM='UM90'
-  YCOMMENT='(kg / s / m²)'
-  IGRID=2
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'UM90'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: UM90'
+  TZFIELD%CUNITS     = 'kg s-1 m-2'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = ''
+  TZFIELD%NGRID      = 2
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
   !  
-  YRECFM='VM90'
-  YCOMMENT='(kg / s / m²)'
-  IGRID=3
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK32,IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'VM90'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: VM90'
+  TZFIELD%CUNITS     = 'kg s-1 m-2'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = ''
+  TZFIELD%NGRID      = 3
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK32)
   !  composantes U et V du flux d'humidité intégré sur 3000 metres
-  YRECFM='UM91'
-  YCOMMENT='(kg / s / m)'
-  IGRID=2
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
-  !
-  YRECFM='VM91'
-  YCOMMENT='(kg / s / m)'
-  IGRID=3
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK22,IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'UM91'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: UM91'
+  TZFIELD%CUNITS     = 'kg s-1 m-1'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = ''
+  TZFIELD%NGRID      = 2
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 2
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK21)
+  !
+  TZFIELD%CMNHNAME   = 'VM91'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: VM91'
+  TZFIELD%CUNITS     = 'kg s-1 m-1'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = ''
+  TZFIELD%NGRID      = 3
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 2
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK22)
   !
   !   Convergence d'humidité
-  YRECFM='HMCONV'
-  YCOMMENT='X_Y_Horizontal CONVergence of moisture flux (kg / s / m)'
-  YCOMMENT='(kg / s / m^3)'
-  IGRID=1
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK35*(-1),IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'HMCONV'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: HMCONV'
+  TZFIELD%CUNITS     = 'kg s-1 m-3'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = 'X_Y_Horizontal CONVergence of moisture flux'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,-ZWORK35)
   !
   !   Convergence d'humidité intégré sur 3000 mètres
-  YRECFM='HMCONV3000'
-  YCOMMENT='X_Y_Horizontal CONVergence of moisture flux (kg / s / m)'
-  YCOMMENT='(kg / s / m^3)'
-  IGRID=1
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK25*(-1),IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'HMCONV3000'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: HMCONV3000'
+  TZFIELD%CUNITS     = 'kg s-1 m-3'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = 'X_Y_Horizontal CONVergence of moisture flux'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 2
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,-ZWORK25)
   !
   IF  (CCLOUD(1:3) == 'ICE' .OR. CCLOUD == 'LIMA') THEN
     !  composantes U et V du flux surfacique d'hydrométéores  
-    YRECFM='UM92'
-    YCOMMENT='(kg / s / m²)'
-    IGRID=2
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK33,IGRID,ILENCH,YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'UM92'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: UM92'
+    TZFIELD%CUNITS     = 'kg s-1 m-2'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = ''
+    TZFIELD%NGRID      = 2
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 3
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK33)
     ! 
-    YRECFM='VM92'
-    YCOMMENT='(kg / s / m²)'
-    IGRID=3
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK34,IGRID,ILENCH,YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'VM92'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: VM92'
+    TZFIELD%CUNITS     = 'kg s-1 m-2'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = ''
+    TZFIELD%NGRID      = 3
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 3
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK34)
     !  composantes U et V du flux d'hydrométéores intégré sur 3000 metres
-    YRECFM='UM93'
-    YCOMMENT='(kg / s / m)'
-    IGRID=2
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK23,IGRID,ILENCH,YCOMMENT,IRESP)
-    YRECFM='VM93'
-    YCOMMENT='(kg / s / m)'
-    IGRID=3
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK24,IGRID,ILENCH,YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'UM93'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: UM93'
+    TZFIELD%CUNITS     = 'kg s-1 m-1'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = ''
+    TZFIELD%NGRID      = 2
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 2
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK23)
+    ! 
+    TZFIELD%CMNHNAME   = 'VM93'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: VM93'
+    TZFIELD%CUNITS     = 'kg s-1 m-1'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = ''
+    TZFIELD%NGRID      = 3
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 2
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK24)
     !   Convergence d'hydrométéores
-    YRECFM='HMCONV_TT'
-    YCOMMENT='X_Y_Horizontal CONVergence of hydrometeor flux (kg / s / m)'
-    YCOMMENT='(kg / s / m^3)'
-    IGRID=1
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK36*(-1),IGRID,ILENCH,YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'HMCONV_TT'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: HMCONV_TT'
+    TZFIELD%CUNITS     = 'kg s-1 m-3'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_Horizontal CONVergence of hydrometeor flux'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 3
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,-ZWORK36)
     !   Convergence d'hydrométéores intégré sur 3000 mètres
-    YRECFM='HMCONV3000_TT'
-    YCOMMENT='X_Y_Horizontal CONVergence of hydrometeor flux (kg / s / m)'
-    YCOMMENT='(kg / s / m^3)'
-    IGRID=1
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK26*(-1),IGRID,ILENCH,YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'HMCONV3000_TT'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: HMCONV3000_TT'
+    TZFIELD%CUNITS     = 'kg s-1 m-3'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_Horizontal CONVergence of hydrometeor flux'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 2
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,-ZWORK26)
   ENDIF
 ENDIF
 !
@@ -1014,129 +1005,160 @@ END IF
 !* Scalar Variables
 !
 ! User scalar variables
-IGRID=1                                     ! individually in the file
+! individually in the file
 IF (LVAR_MRSV) THEN
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CUNITS     = 'g kg-1'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  !
   DO JSV = 1,NSV_USER
-    WRITE(YRECFM,'(A4,I3.3)')'MRSV',JSV
-    WRITE(YCOMMENT,'(A6,A4,I3.3,A7)')'X_Y_Z_','MRSV',JSV,' (G/KG)'
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=XSVT(:,:,:,JSV)*1.E3
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-         YCOMMENT,IRESP)
+    WRITE(TZFIELD%CMNHNAME,'(A4,I3.3)')'MRSV',JSV
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    WRITE(TZFIELD%CCOMMENT,'(A6,A4,I3.3)')'X_Y_Z_','MRSV',JSV
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XSVT(:,:,:,JSV)*1.E3)
   END DO
 END IF
 ! microphysical C2R2 scheme scalar variables
 IF(LVAR_MRW) THEN
-  DO JSV = NSV_C2R2BEG,NSV_C2R2END
-    YRECFM=TRIM(C2R2NAMES(JSV-NSV_C2R2BEG+1))//'T'
-    IF (JSV < NSV_C2R2END) THEN
-      WRITE(YCOMMENT,'(A6,A4,I3.3,A7)')'X_Y_Z_','MRSV',JSV,' (/CM3)'
-    ZWORK31(:,:,:)=XSVT(:,:,:,JSV)*1.E-6
-    ELSE
-      WRITE(YCOMMENT,'(A6,A4,I3.3,A7)')'X_Y_Z_','MRSV',JSV,' (/L)'
-      ZWORK31(:,:,:)=XSVT(:,:,:,JSV)*1.E-3
-    ENDIF
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-         YCOMMENT,IRESP)
-  END DO
+  IF (NSV_C2R2END>=NSV_C2R2BEG) THEN
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 3
+    !
+    DO JSV = NSV_C2R2BEG,NSV_C2R2END
+      TZFIELD%CMNHNAME   = TRIM(C2R2NAMES(JSV-NSV_C2R2BEG+1))//'T'
+      TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+      IF (JSV < NSV_C2R2END) THEN
+        TZFIELD%CUNITS     = 'cm-3'
+        ZWORK31(:,:,:)=XSVT(:,:,:,JSV)*1.E-6
+      ELSE
+        TZFIELD%CUNITS     = 'l-1'
+        ZWORK31(:,:,:)=XSVT(:,:,:,JSV)*1.E-3
+      END IF
+      WRITE(TZFIELD%CCOMMENT,'(A6,A4,I3.3)')'X_Y_Z_','MRSV',JSV
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
+    END DO
+  END IF
   ! microphysical C3R5 scheme additional scalar variables
-  DO JSV = NSV_C1R3BEG,NSV_C1R3END
-    YRECFM=TRIM(C1R3NAMES(JSV-NSV_C1R3BEG+1))//'T'
-    WRITE(YCOMMENT,'(A6,A3,I3.3,A8)')'X_Y_Z_','SVT',JSV,' (/L)'
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=XSVT(:,:,:,JSV)*1.E-3
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-         YCOMMENT,IRESP)   
-  END DO
+  IF (NSV_C1R3END>=NSV_C1R3BEG) THEN
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CUNITS     = 'l-1'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 3
+    !
+    DO JSV = NSV_C1R3BEG,NSV_C1R3END
+      TZFIELD%CMNHNAME   = TRIM(C1R3NAMES(JSV-NSV_C1R3BEG+1))//'T'
+      TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+      WRITE(TZFIELD%CCOMMENT,'(A6,A3,I3.3)')'X_Y_Z_','SVT',JSV
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XSVT(:,:,:,JSV)*1.E-3)
+    END DO
+  END IF
 END IF
 !
 ! microphysical LIMA scheme scalar variables
 !
 IF (LLIMA_DIAG) THEN
-   DO JSV = NSV_LIMA_BEG,NSV_LIMA_END
+  IF (NSV_LIMA_END>=NSV_LIMA_BEG) THEN
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 3
+  END IF
+  !
+  DO JSV = NSV_LIMA_BEG,NSV_LIMA_END
+    !
+    TZFIELD%CUNITS     = 'cm-3'
+    WRITE(TZFIELD%CCOMMENT,'(A6,A3,I3.3)')'X_Y_Z_','SVT',JSV
+    !
 ! Nc
-      IF (JSV .EQ. NSV_LIMA_NC) THEN
-         YRECFM=TRIM(CLIMA_WARM_CONC(1))//'T'
-         WRITE(YCOMMENT,'(A6,A3,I3.3,A8)')'X_Y_Z_','SVT',JSV,' (/cm^3)'
-      END IF
+    IF (JSV .EQ. NSV_LIMA_NC) THEN
+      TZFIELD%CMNHNAME   = TRIM(CLIMA_WARM_CONC(1))//'T'
+    END IF
 ! Nr
-      IF (JSV .EQ. NSV_LIMA_NR) THEN
-         YRECFM=TRIM(CLIMA_WARM_CONC(2))//'T'
-         WRITE(YCOMMENT,'(A6,A3,I3.3,A8)')'X_Y_Z_','SVT',JSV,' (/cm^3)'
-      END IF
+    IF (JSV .EQ. NSV_LIMA_NR) THEN
+      TZFIELD%CMNHNAME   = TRIM(CLIMA_WARM_CONC(2))//'T'
+    END IF
 ! N CCN free
-      IF (JSV .GE. NSV_LIMA_CCN_FREE .AND. JSV .LT. NSV_LIMA_CCN_ACTI) THEN
-         WRITE(INDICE,'(I2.2)')(JSV - NSV_LIMA_CCN_FREE + 1)
-         YRECFM=TRIM(CLIMA_WARM_CONC(3))//INDICE//'T'
-         WRITE(YCOMMENT,'(A6,A3,I3.3,A8)')'X_Y_Z_','SVT',JSV,' (/cm^3)'
-      END IF
+    IF (JSV .GE. NSV_LIMA_CCN_FREE .AND. JSV .LT. NSV_LIMA_CCN_ACTI) THEN
+      WRITE(INDICE,'(I2.2)')(JSV - NSV_LIMA_CCN_FREE + 1)
+      TZFIELD%CMNHNAME   = TRIM(CLIMA_WARM_CONC(3))//INDICE//'T'
+    END IF
 ! N CCN acti
-      IF (JSV .GE. NSV_LIMA_CCN_ACTI .AND. JSV .LT. NSV_LIMA_CCN_ACTI + NMOD_CCN) THEN
-         WRITE(INDICE,'(I2.2)')(JSV - NSV_LIMA_CCN_ACTI + 1)
-         YRECFM=TRIM(CLIMA_WARM_CONC(4))//INDICE//'T'
-         WRITE(YCOMMENT,'(A6,A3,I3.3,A8)')'X_Y_Z_','SVT',JSV,' (/cm^3)'
-      END IF
+    IF (JSV .GE. NSV_LIMA_CCN_ACTI .AND. JSV .LT. NSV_LIMA_CCN_ACTI + NMOD_CCN) THEN
+      WRITE(INDICE,'(I2.2)')(JSV - NSV_LIMA_CCN_ACTI + 1)
+      TZFIELD%CMNHNAME   = TRIM(CLIMA_WARM_CONC(4))//INDICE//'T'
+    END IF
 ! Scavenging
-      IF (JSV .EQ. NSV_LIMA_SCAVMASS) THEN
-         YRECFM=TRIM(CAERO_MASS(1))//'T'
-         WRITE(YCOMMENT,'(A6,A3,I3.3,A8)')'X_Y_Z_','SVT',JSV,' (kg/cm^3)'
-      END IF
+    IF (JSV .EQ. NSV_LIMA_SCAVMASS) THEN
+      TZFIELD%CMNHNAME   = TRIM(CAERO_MASS(1))//'T'
+      TZFIELD%CUNITS     = 'kg cm-3'
+    END IF
 ! Ni
-      IF (JSV .EQ. NSV_LIMA_NI) THEN
-         YRECFM=TRIM(CLIMA_COLD_CONC(1))//'T'
-         WRITE(YCOMMENT,'(A6,A3,I3.3,A8)')'X_Y_Z_','SVT',JSV,' (/cm^3)'
-      END IF
+    IF (JSV .EQ. NSV_LIMA_NI) THEN
+      TZFIELD%CMNHNAME   = TRIM(CLIMA_COLD_CONC(1))//'T'
+    END IF
 ! N IFN free
-      IF (JSV .GE. NSV_LIMA_IFN_FREE .AND. JSV .LT. NSV_LIMA_IFN_NUCL) THEN
-         WRITE(INDICE,'(I2.2)')(JSV - NSV_LIMA_IFN_FREE + 1)
-         YRECFM=TRIM(CLIMA_COLD_CONC(2))//INDICE//'T'
-         WRITE(YCOMMENT,'(A6,A3,I3.3,A8)')'X_Y_Z_','SVT',JSV,' (/cm^3)'
-      END IF
+    IF (JSV .GE. NSV_LIMA_IFN_FREE .AND. JSV .LT. NSV_LIMA_IFN_NUCL) THEN
+      WRITE(INDICE,'(I2.2)')(JSV - NSV_LIMA_IFN_FREE + 1)
+      TZFIELD%CMNHNAME   = TRIM(CLIMA_COLD_CONC(2))//INDICE//'T'
+    END IF
 ! N IFN nucl
-      IF (JSV .GE. NSV_LIMA_IFN_NUCL .AND. JSV .LT. NSV_LIMA_IFN_NUCL + NMOD_IFN) THEN
-         WRITE(INDICE,'(I2.2)')(JSV - NSV_LIMA_IFN_NUCL + 1)
-         YRECFM=TRIM(CLIMA_COLD_CONC(3))//INDICE//'T'
-         WRITE(YCOMMENT,'(A6,A3,I3.3,A8)')'X_Y_Z_','SVT',JSV,' (/cm^3)'
-      END IF
+    IF (JSV .GE. NSV_LIMA_IFN_NUCL .AND. JSV .LT. NSV_LIMA_IFN_NUCL + NMOD_IFN) THEN
+      WRITE(INDICE,'(I2.2)')(JSV - NSV_LIMA_IFN_NUCL + 1)
+      TZFIELD%CMNHNAME   = TRIM(CLIMA_COLD_CONC(3))//INDICE//'T'
+    END IF
 ! N IMM nucl
-      I = 0
-      IF (JSV .GE. NSV_LIMA_IMM_NUCL .AND. JSV .LT. NSV_LIMA_IMM_NUCL + NMOD_IMM) THEN
-         I = I + 1
-         WRITE(INDICE,'(I2.2)')(NINDICE_CCN_IMM(I))
-         YRECFM=TRIM(CLIMA_COLD_CONC(4))//INDICE//'T'
-         WRITE(YCOMMENT,'(A6,A3,I3.3,A8)')'X_Y_Z_','SVT',JSV,' (/cm^3)'
-      END IF
+    I = 0
+    IF (JSV .GE. NSV_LIMA_IMM_NUCL .AND. JSV .LT. NSV_LIMA_IMM_NUCL + NMOD_IMM) THEN
+      I = I + 1
+      WRITE(INDICE,'(I2.2)')(NINDICE_CCN_IMM(I))
+      TZFIELD%CMNHNAME   = TRIM(CLIMA_COLD_CONC(4))//INDICE//'T'
+    END IF
 ! Hom. freez. of CCN
-      IF (JSV .EQ. NSV_LIMA_HOM_HAZE) THEN
-         YRECFM=TRIM(CLIMA_COLD_CONC(5))//'T'
-         WRITE(YCOMMENT,'(A6,A3,I3.3,A8)')'X_Y_Z_','SVT',JSV,' (/cm^3)'
-      END IF
-!
-!
-      ILENCH=LEN(YCOMMENT)
-      ZWORK31(:,:,:)=XSVT(:,:,:,JSV)*1.E-6*XRHODREF(:,:,:)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                  YCOMMENT,IRESP)
-   END DO
+    IF (JSV .EQ. NSV_LIMA_HOM_HAZE) THEN
+      TZFIELD%CMNHNAME   = TRIM(CLIMA_COLD_CONC(5))//'T'
+    END IF
+    !
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    ZWORK31(:,:,:)=XSVT(:,:,:,JSV)*1.E-6*XRHODREF(:,:,:)
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
+  END DO
 !
-   IF (LUSERC) THEN
-      YRECFM= 'LWC'
-      YCOMMENT='X_Y_Z_LWC (g/m^3)'
-      ILENCH=LEN(YCOMMENT)
-      ZWORK31(:,:,:)=XRT(:,:,:,2)*1.E3*XRHODREF(:,:,:)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,  &
-                  YCOMMENT,IRESP)
-   END IF
+  IF (LUSERC) THEN
+    TZFIELD%CMNHNAME   = 'LWC'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: LWC'
+    TZFIELD%CUNITS     = 'g m-3'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_Z_LWC'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 3
+    ZWORK31(:,:,:)=XRT(:,:,:,2)*1.E3*XRHODREF(:,:,:)
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
+  END IF
 !
-   IF (LUSERI) THEN
-      YRECFM= 'IWC'
-      YCOMMENT='X_Y_Z_MRI (g/m^3)'
-      ILENCH=LEN(YCOMMENT)
-      ZWORK31(:,:,:)=XRT(:,:,:,4)*1.E3*XRHODREF(:,:,:)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,  &
-                  YCOMMENT,IRESP)
-   END IF
+  IF (LUSERI) THEN
+    TZFIELD%CMNHNAME   = 'IWC'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: IWC'
+    TZFIELD%CUNITS     = 'g m-3'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_Z_MRI'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 3
+    ZWORK31(:,:,:)=XRT(:,:,:,4)*1.E3*XRHODREF(:,:,:)
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
+  END IF
 !
 END IF
 !
@@ -1156,45 +1178,51 @@ IF (LCHEMDIAG) THEN
   END DO
 END IF
 IF (LCHAQDIAG) THEN    !aqueous concentration in M
+  TZFIELD%CSTDNAME   = ''
+  !PW TODO: check units
+  TZFIELD%CUNITS     = ''
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  !
   ZWORK31(:,:,:)=0.
   DO JSV = NSV_CHACBEG, NSV_CHACBEG-1+NEQAQ/2   !cloud water
-    YRECFM=TRIM(CNAMES(JSV-NSV_CHACBEG+NSV_CHGS+1))//'M'
-    WRITE(YCOMMENT,'(A6,A4,I3.3,A4)')'X_Y_Z_','CHAQ',JSV,' (M)'
-    ILENCH=LEN(YCOMMENT)
+    TZFIELD%CMNHNAME   = TRIM(CNAMES(JSV-NSV_CHACBEG+NSV_CHGS+1))//'M'
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    WRITE(TZFIELD%CCOMMENT,'(A6,A4,I3.3,A4)')'X_Y_Z_','CHAQ',JSV,' (M)'
     WHERE(((XRT(:,:,:,2)*XRHODREF(:,:,:))/1.e3) .GE. XRTMIN_AQ)
       ZWORK31(:,:,:)=(XSVT(:,:,:,JSV)*1000.)/(XMD*1.E+3*XRT(:,:,:,2))
     ENDWHERE
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-         YCOMMENT,IRESP)
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
   END DO
+  !
   ZWORK31(:,:,:)=0.
   DO JSV = NSV_CHACBEG+NEQAQ/2, NSV_CHACEND    !rain water
-    YRECFM=TRIM(CNAMES(JSV-NSV_CHACBEG+NSV_CHGS+1))//'M'
-    WRITE(YCOMMENT,'(A6,A4,I3.3,A4)')'X_Y_Z_','CHAQ',JSV,' (M)'
-    ILENCH=LEN(YCOMMENT)
+    TZFIELD%CMNHNAME   = TRIM(CNAMES(JSV-NSV_CHACBEG+NSV_CHGS+1))//'M'
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    WRITE(TZFIELD%CCOMMENT,'(A6,A4,I3.3,A4)')'X_Y_Z_','CHAQ',JSV,' (M)'
     WHERE(((XRT(:,:,:,3)*XRHODREF(:,:,:))/1.e3) .GE. XRTMIN_AQ)
       ZWORK31(:,:,:)=(XSVT(:,:,:,JSV)*1000.)/(XMD*1.E+3*XRT(:,:,:,3))
     ENDWHERE
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-         YCOMMENT,IRESP)
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
   END DO
 !  ZWORK31(:,:,:)=0.
 !  DO JSV = NSV_CHICBEG,NSV_CHICEND   ! ice phase
-!    YRECFM=TRIM(CICNAMES(JSV-NSV_CHICBEG+1))//'T'
-!    WRITE(YCOMMENT,'(A6,A4,I3.3,A4)')'X_Y_Z_','CHIC',JSV,' (M)'
-!    ILENCH=LEN(YCOMMENT)
+!    TZFIELD%CMNHNAME   = TRIM(CICNAMES(JSV-NSV_CHICBEG+1))//'T'
+!    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+!    WRITE(TZFIELD%CCOMMENT,'(A6,A4,I3.3,A4)')'X_Y_Z_','CHIC',JSV,' (M)'
 !    WHERE(((XRT(:,:,:,3)*XRHODREF(:,:,:))/1.e3) .GE. XRTMIN_AQ)
 !      ZWORK31(:,:,:)=(XSVT(:,:,:,JSV)*1000.)/(XMD*1.E+3*XRT(:,:,:,3))
 !    ENDWHERE
-!    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-!         YCOMMENT,IRESP)
+!    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
 !  END DO
 END IF
 
 ! Passive polluant scalar variables
 IF (LPASPOL) THEN
   ALLOCATE(ZRHOT( SIZE(XTHT,1), SIZE(XTHT,2),SIZE(XTHT,3)))
-  ALLOCATE(ZTMP( SIZE(XTHT,1), SIZE(XTHT,2),SIZE(XTHT,3)))
+  ALLOCATE(ZTMP(  SIZE(XTHT,1), SIZE(XTHT,2),SIZE(XTHT,3)))
 !
 !*	Density                                          
 !
@@ -1203,40 +1231,57 @@ IF (LPASPOL) THEN
 !*	Conversion g/m3.
 !
   ZRHOT(:,:,:)=ZRHOT(:,:,:)*1000.0
-
+  !
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CUNITS     = 'g m-3'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  !
   DO JSV = 1,NSV_PP
     ZTMP(:,:,:)=ABS( XSVT(:,:,:,JSV+NSV_PPBEG-1)*ZRHOT(:,:,:) )
-    WRITE(YRECFM,'(A3,I3.3)')'PPT',JSV
-    WRITE(YCOMMENT,'(A6,A3,I3.3,A8)')'X_Y_Z_','PPT',JSV,' (G/M3) '
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZTMP,IGRID,ILENCH,       &
-     YCOMMENT,IRESP)
-    !
+    WRITE(TZFIELD%CMNHNAME,'(A3,I3.3)')'PPT',JSV
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CCOMMENT   = 'X_Y_Z_'//TRIM(TZFIELD%CMNHNAME)
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZTMP)
   END DO
   DEALLOCATE(ZTMP)
   DEALLOCATE(ZRHOT)
 END IF
 ! Conditional sampling variables
 IF (LCONDSAMP) THEN
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CUNITS     = '1'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  !
   DO JSV = NSV_CSBEG,NSV_CSEND
-    WRITE(YRECFM,'(A3,I3.3)')'CST',JSV
-    WRITE(YCOMMENT,'(A6,A3,I3.3,A8)')'X_Y_Z_','CST',JSV,' () '
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',XSVT(:,:,:,JSV),IGRID,ILENCH,       &
-     YCOMMENT,IRESP)
-    !
+    WRITE(TZFIELD%CMNHNAME,'(A3,I3.3)')'CST',JSV
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CCOMMENT   = 'X_Y_Z_'//TRIM(TZFIELD%CMNHNAME)
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XSVT(:,:,:,JSV))
   END DO
 END IF
 
 
 ! Lagrangian variables
 IF (LTRAJ) THEN
+  TZFIELD%CSTDNAME   = ''
+  !PW TODO: check units
+  TZFIELD%CUNITS     = ''
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  !
   DO JSV = NSV_LGBEG,NSV_LGEND
-    YRECFM=TRIM(CLGNAMES(JSV-NSV_LGBEG+1))//'T'
-    WRITE(YCOMMENT,'(A6,A20,I3.3,A7)')'X_Y_Z_','Lagrangian variable ',JSV,' (M)'
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',XSVT(:,:,:,JSV),IGRID,ILENCH,   &
-                YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = TRIM(CLGNAMES(JSV-NSV_LGBEG+1))//'T'
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    WRITE(TZFIELD%CCOMMENT,'(A6,A20,I3.3,A4)')'X_Y_Z_','Lagrangian variable ',JSV,' (M)'
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XSVT(:,:,:,JSV))
   END DO
   ! X coordinate
   DO JK=1,IKU
@@ -1247,11 +1292,10 @@ IF (LTRAJ) THEN
       ZWORK31(IIU,JJ,JK)=2.*ZWORK31(IIU-1,JJ,JK) - ZWORK31(IIU-2,JJ,JK)
     END DO
   END DO
-  WRITE(YRECFM,'(A1)')'X'
-  WRITE(YCOMMENT,'(A6,A13,A7)')'X_Y_Z_','X coordinate ',' (M)'
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-              YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'X'
+  TZFIELD%CLONGNAME  = 'MesoNH: X'
+  TZFIELD%CCOMMENT   = 'X_Y_Z_X coordinate'
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
   ! Y coordinate
   DO JK=1,IKU
     DO JI=1,IIU
@@ -1261,35 +1305,45 @@ IF (LTRAJ) THEN
       ZWORK31(JI,IJU,JK)=2.*ZWORK31(JI,IJU-1,JK) - ZWORK31(JI,IJU-2,JK)
     END DO
   END DO
-  WRITE(YRECFM,'(A1)')'Y'
-  WRITE(YCOMMENT,'(A6,A13,A7)')'X_Y_Z_','Y coordinate ',' (M)'
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-              YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'Y'
+  TZFIELD%CLONGNAME  = 'MesoNH: Y'
+  TZFIELD%CCOMMENT   = 'X_Y_Z_Y coordinate'
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
 END IF
 ! linox scalar variables
 IF (.NOT.(LUSECHEM .OR. LCHEMDIAG) .AND. LCH_CONV_LINOX) THEN
   DO JSV = NSV_LNOXBEG,NSV_LNOXEND
-    YRECFM='LINOXT'
-    WRITE(YCOMMENT,'(A6,A4,I3.3,A7)')'X_Y_Z_','LNOX',JSV,' (ppbv)'
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=XSVT(:,:,:,JSV)*1.E9
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-         YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'LINOXT'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'ppbv'
+    TZFIELD%CDIR       = 'XY'
+    WRITE(TZFIELD%CCOMMENT,'(A6,A4,I3.3)')'X_Y_Z_','LNOX',JSV
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 3
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XSVT(:,:,:,JSV)*1.E9)
   END DO
 END IF
 IF (LELECDIAG .AND. CELEC .NE. "NONE") THEN
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  !
   DO JSV = NSV_ELECBEG,NSV_ELECEND
-    YRECFM=TRIM(CELECNAMES(JSV-NSV_ELECBEG+1))//'T'
-    IF (JSV .GT. NSV_ELECBEG .AND. JSV .LT. NSV_ELECEND) THEN
-      WRITE(YCOMMENT,'(A6,A3,I3.3,A8)')'X_Y_Z_','SVT',JSV,' (C/m3)'
+    TZFIELD%CMNHNAME   = TRIM(CELECNAMES(JSV-NSV_ELECBEG+1))//'T'
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    IF (JSV .GT. NSV_ELECBEG .AND. JSV .LT. NSV_ELECEND) THEN 
+      TZFIELD%CUNITS     = 'C m-3'
+      WRITE(TZFIELD%CCOMMENT,'(A6,A3,I3.3,A8)')'X_Y_Z_','SVT',JSV,' (C/m3)'
     ELSE
-      WRITE(YCOMMENT,'(A6,A3,I3.3,A8)')'X_Y_Z_','SVT',JSV,' (nb ions/m3)'
+      TZFIELD%CUNITS     = 'm-3'
+      WRITE(TZFIELD%CCOMMENT,'(A6,A3,I3.3,A8)')'X_Y_Z_','SVT',JSV,' (nb ions/m3)'
     END IF
-    ILENCH=LEN(YCOMMENT)
     ZWORK31(:,:,:)=XSVT(:,:,:,JSV) * XRHODREF(:,:,:)  ! C/kg --> C/m3
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-         YCOMMENT,IRESP)
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
   END DO
 END IF
 ! Sea Salt variables
@@ -1301,54 +1355,63 @@ IF (LSALT) THEN
   IF(.NOT.ALLOCATED(ZN0_SLT))  &
     ALLOCATE(ZN0_SLT(SIZE(XSVT,1), SIZE(XSVT,2), SIZE(XSVT,3), NMODE_SLT))
   !
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CUNITS     = 'ppbv'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  !
   DO JSV = NSV_SLTBEG,NSV_SLTEND
-    YRECFM=TRIM(UPCASE(CSALTNAMES(JSV-NSV_SLTBEG+1)))//'T'
-    WRITE(YCOMMENT,'(A6,A4,I3.3,A7)')'X_Y_Z_','SALT',JSV,' (ppbv)'
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=XSVT(:,:,:,JSV)*1.E9
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = TRIM(UPCASE(CSALTNAMES(JSV-NSV_SLTBEG+1)))//'T'
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    WRITE(TZFIELD%CCOMMENT,'(A6,A4,I3.3)')'X_Y_Z_','SALT',JSV
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XSVT(:,:,:,JSV)*1.E9)
   END DO
   !
   CALL PPP2SALT(XSVT(:,:,:,NSV_SLTBEG:NSV_SLTEND),XRHODREF,&
                PSIG3D=ZSIG_SLT, PRG3D=ZRG_SLT, PN3D=ZN0_SLT)
+  !
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  !
   DO JJ=1,NMODE_SLT
-    WRITE(YRECFM,'(A6,I1)')'SLTRGA',JJ
-    WRITE(YCOMMENT,'(A18,I1,A5)')'RG (nb) SALT MODE ',JJ,' (um)'
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=ZRG_SLT(:,:,:,JJ)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'SLTRGA'
+    TZFIELD%CLONGNAME  = 'MesoNH: SLTRGA'
+    TZFIELD%CUNITS     = 'um'
+    WRITE(TZFIELD%CCOMMENT,'(A18,I1)')'RG (nb) SALT MODE ',JJ
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZRG_SLT(:,:,:,JJ))
     !
-    WRITE(YRECFM,'(A7,I1)')'SLTRGAM',JJ
-    WRITE(YCOMMENT,'(A17,I1,A5)')'RG (m) SALT MODE ',JJ,' (um)'
-    ILENCH=LEN(YCOMMENT)
+    TZFIELD%CMNHNAME   = 'SLTRGAM'
+    TZFIELD%CLONGNAME  = 'MesoNH: SLTRGAM'
+    TZFIELD%CUNITS     = 'um'
+    WRITE(TZFIELD%CCOMMENT,'(A17,I1)')'RG (m) SALT MODE ',JJ
     ZWORK31(:,:,:)=ZRG_SLT(:,:,:,JJ) / (EXP(-3.*(LOG(ZSIG_SLT(:,:,:,JJ)))**2))
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
     !
-    WRITE(YRECFM,'(A6,I1)')'SLTN0A',JJ
-    WRITE(YCOMMENT,'(A13,I1,A7)')'N0 SALT MODE ',JJ,' (1/m3)'
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=ZN0_SLT(:,:,:,JJ)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    WRITE(TZFIELD%CMNHNAME,'(A6,I1)')'SLTN0A',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'm-3'
+    WRITE(TZFIELD%CCOMMENT,'(A13,I1)')'N0 SALT MODE ',JJ
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZN0_SLT(:,:,:,JJ))
     !
-    WRITE(YRECFM,'(A7,I1)')'SLTSIGA',JJ
-    WRITE(YCOMMENT,'(A16,I1)')'SIGMA SALT MODE ',JJ
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=ZSIG_SLT(:,:,:,JJ)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    WRITE(TZFIELD%CMNHNAME,'(A7,I1)')'SLTSIGA',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = '1'
+    WRITE(TZFIELD%CCOMMENT,'(A16,I1)')'SIGMA SALT MODE ',JJ
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZSIG_SLT(:,:,:,JJ))
     !SALT MASS CONCENTRATION
-    WRITE(YRECFM,'(A4,I1)')'SLTMSS',JJ
-    WRITE(YCOMMENT,'(A14,I1,A7)')'MASSCONC MODE ',JJ,'(ug/m3)'
-    ILENCH=LEN(YCOMMENT)
+    WRITE(TZFIELD%CMNHNAME,'(A4,I1)')'SLTMSS',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'ug m-3'
+    WRITE(TZFIELD%CCOMMENT,'(A14,I1)')'MASSCONC MODE ',JJ
     ZWORK31(:,:,:)= ZN0_SLT(:,:,:,JJ)*4./3.*3.14*2500.*1e9 & !kg-->ug
        * (ZRG_SLT(:,:,:,JJ)**3)*1.d-18 &  !um-->m
        * exp(4.5*log(ZSIG_SLT(:,:,:,JJ))*log(ZSIG_SLT(:,:,:,JJ)))
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
     !SALT BURDEN (g/m2)
     ZWORK21(:,:)=0.0
     DO JK=IKB,IKE
@@ -1362,106 +1425,110 @@ IF (LSALT) THEN
         ENDDO
       ENDDO
     ENDDO
-    WRITE(YRECFM,'(A7,I1)')'SLTBRDN',JJ
-    WRITE(YCOMMENT,'(A13,I1)')'BURDEN (g/m2)',JJ
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    WRITE(TZFIELD%CMNHNAME,'(A7,I1)')'SLTBRDN',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'g m-2'
+    WRITE(TZFIELD%CCOMMENT,'(A6,I1)')'BURDEN',JJ
+    TZFIELD%NDIMS      = 2
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK21)
+    !
+    TZFIELD%NDIMS      = 3
   ENDDO
 END IF
 IF (LSALT.AND.LDEPOS_SLT(IMI)) THEN
-    ! 
-     ZSSLTDEP=XSVT(:,:,:,NSV_SLTDEPBEG:NSV_SLTDEPEND)
-     DO JSV = 1,NSV_SLTDEP   
-      YRECFM=TRIM(UPCASE(CDESLTNAMES(JSV)))//'T'
-       WRITE(YCOMMENT,'(A6,A4,I3.3,A7)')'X_Y_Z_','SALTDEP',JSV,' (ppbv)'
-       ILENCH=LEN(YCOMMENT)
-       ZWORK31(:,:,:)=ZSSLTDEP(:,:,:,JSV)*1.E9
-       CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
-     END DO
-     DO JJ=1,NMODE_SLT           
-            ZWORK31(:,:,:)=0.0
-! FOR CLOUDS
-        WRITE(YRECFM,'(A9,I1)')'SLTDEPN0A',JJ
-        WRITE(YCOMMENT,'(A16,I1,A7)')'N0 DUSTDEP MODE ', &
-                                   JJ,' (1/m3)'            
-! CLOUD: CALCULATE MOMENT 3 FROM TOTAL AEROSOL MASS            
-            ZWORK31(:,:,:)=ZSSLTDEP(:,:,:,JJ)  &!==>molec_{aer}/molec_{air}
-            *(XMOLARWEIGHT_DUST/XMD)           &!==>kg_{aer}/kg_{air}
-            *XRHODREF(:,:,:)                   &!==>kg_{aer}/m3_{air}
-            *(1.d0/XDENSITY_DUST)              &!==>m3_{aer}/m3_{air}
-            *XM3TOUM3                          &!==>um3_{aer}/m3_{air}
-            /(XPI*4./3.)                        !==>um3_{aer}/m3_{air}
+  !
+  ZSSLTDEP=XSVT(:,:,:,NSV_SLTDEPBEG:NSV_SLTDEPEND)
+  !
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CUNITS     = 'ppbv'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  !
+  DO JSV = 1,NSV_SLTDEP
+    TZFIELD%CMNHNAME   = TRIM(UPCASE(CDESLTNAMES(JSV)))//'T'
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    WRITE(TZFIELD%CCOMMENT,'(A6,A4,I3.3)')'X_Y_Z_','SALTDEP',JSV
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZSSLTDEP(:,:,:,JSV)*1.E9)
+  END DO
+  !
+  DO JJ=1,NMODE_SLT
+    ! FOR CLOUDS
+    WRITE(TZFIELD%CMNHNAME,'(A9,I1)')'SLTDEPN0A',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    WRITE(TZFIELD%CCOMMENT,'(A16,I1)')'N0 DUSTDEP MODE ',JJ
+    TZFIELD%CUNITS     = 'm-3'
+    ! CLOUD: CALCULATE MOMENT 3 FROM TOTAL AEROSOL MASS
+    ZWORK31(:,:,:) = ZSSLTDEP(:,:,:,JJ)         &!==>molec_{aer}/molec_{air}
+                     *(XMOLARWEIGHT_DUST/XMD)   &!==>kg_{aer}/kg_{air}
+                     *XRHODREF(:,:,:)           &!==>kg_{aer}/m3_{air}
+                     /XDENSITY_DUST             &!==>m3_{aer}/m3_{air}
+                     *XM3TOUM3                  &!==>um3_{aer}/m3_{air}
+                     /(XPI*4./3.)                !==>um3_{aer}/m3_{air}
             !==>volume 3rd moment
-!CLOUD: CALCULATE MOMENT 0 FROM DISPERSION AND MEAN RADIUS            
-            ZWORK31(:,:,:)=  ZWORK31(:,:,:)/      &
-                    ((ZRG_SLT(:,:,:,JJ)**3)*      &
-                    EXP(4.5 * LOG(ZSIG_SLT(:,:,:,JJ))**2))
-!CLOUD: RETURN TO CONCENTRATION #/m3
-            ZWORK31(:,:,:)= ZWORK31(:,:,:) *   XMD/ &
-                     (XAVOGADRO*XRHODREF(:,:,:))
-!CLOUD:  Get number concentration (#/molec_{air}==>#/m3)  
-             ZWORK31(:,:,:)=                         &
-                    ZWORK31(:,:,:)                  & !#/molec_{air}
-                    * XAVOGADRO                     & !==>#/mole
-                    / XMD                           & !==>#/kg_{air}
-                    * XRHODREF(:,:,:)                 !==>#/m3  
-     ILENCH=LEN(YCOMMENT)               
-     CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
-! CLOUD:   DUST MASS CONCENTRATION
-       WRITE(YRECFM,'(A9,I1)')'SLTDEPMSS',JJ
-       WRITE(YCOMMENT,'(A17,I1,A7)')'DEPMASSCONC MODE ', &
-                                 JJ,'(ug/m3)'
-       ILENCH=LEN(YCOMMENT)
-       ZWORK31(:,:,:)= ZWORK31(:,:,:)*4./3.*3.14*2500.*1e9 & !kg-->ug
-          * (ZRG_SLT(:,:,:,JJ)**3)*1.d-18 &  !um-->m
-          * exp(4.5*log(ZSIG_SLT(:,:,:,JJ))*log(ZSIG_SLT(:,:,:,JJ)))
-       CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)                
-!   FOR RAIN DROPS                
-        WRITE(YRECFM,'(A9,I1)')'SLTDEPN0A',JJ+NMODE_SLT
-        WRITE(YCOMMENT,'(A16,I1,A7)')'N0 DUSTDEP MODE ', &
-                                   JJ+NMODE_SLT,' (1/m3)'
-        ILENCH=LEN(YCOMMENT) 
-        ZWORK31(:,:,:)=0.0
-! RAIN: CALCULATE MOMENT 3 FROM TOTAL AEROSOL MASS                
-        ZWORK31(:,:,:)=ZSSLTDEP(:,:,:,JJ+NMODE_SLT)  &!==>molec_{aer}/molec_{air}
-            *(XMOLARWEIGHT_DUST/XMD)           &!==>kg_{aer}/kg_{air}
-            *XRHODREF(:,:,:)                   &!==>kg_{aer}/m3_{air}
-            *(1.d0/XDENSITY_DUST)              &!==>m3_{aer}/m3_{air}
-            *XM3TOUM3                          &!==>um3_{aer}/m3_{air}
-            /(XPI*4./3.)                        !==>um3_{aer}/m3_{air}
+    !CLOUD: CALCULATE MOMENT 0 FROM DISPERSION AND MEAN RADIUS
+    ZWORK31(:,:,:) =  ZWORK31(:,:,:)/                        &
+                      ((ZRG_SLT(:,:,:,JJ)**3)*               &
+                      EXP(4.5 * LOG(ZSIG_SLT(:,:,:,JJ))**2))
+    !CLOUD: RETURN TO CONCENTRATION #/m3
+    ZWORK31(:,:,:)= ZWORK31(:,:,:) *   XMD/ &
+                    (XAVOGADRO*XRHODREF(:,:,:))
+    !CLOUD:  Get number concentration (#/molec_{air}==>#/m3)
+    ZWORK31(:,:,:)=                   &
+                    ZWORK31(:,:,:)    & !#/molec_{air}
+                    * XAVOGADRO       & !==>#/mole
+                    / XMD             & !==>#/kg_{air}
+                    * XRHODREF(:,:,:)   !==>#/m3
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
+    ! CLOUD:   DUST MASS CONCENTRATION
+    WRITE(TZFIELD%CMNHNAME,'(A9,I1)')'SLTDEPMSS',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    WRITE(TZFIELD%CCOMMENT,'(A17,I1)')'DEPMASSCONC MODE ',JJ
+    TZFIELD%CUNITS     = 'ug m-3'
+    ZWORK31(:,:,:)= ZWORK31(:,:,:)*4./3.*3.14*2500.*1e9 & !kg-->ug
+                    * (ZRG_SLT(:,:,:,JJ)**3)*1.d-18     &  !um-->m
+                    * exp(4.5*log(ZSIG_SLT(:,:,:,JJ))*log(ZSIG_SLT(:,:,:,JJ)))
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
+    !   FOR RAIN DROPS
+    WRITE(TZFIELD%CMNHNAME,'(A9,I1)')'SLTDEPN0A',JJ+NMODE_SLT
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    WRITE(TZFIELD%CCOMMENT,'(A16,I1)')'N0 DUSTDEP MODE ',JJ+NMODE_SLT
+    TZFIELD%CUNITS     = 'm-3'
+    ! RAIN: CALCULATE MOMENT 3 FROM TOTAL AEROSOL MASS
+    ZWORK31(:,:,:) = ZSSLTDEP(:,:,:,JJ+NMODE_SLT)  &!==>molec_{aer}/molec_{air}
+                     *(XMOLARWEIGHT_DUST/XMD)      &!==>kg_{aer}/kg_{air}
+                     *XRHODREF(:,:,:)              &!==>kg_{aer}/m3_{air}
+                     /XDENSITY_DUST                &!==>m3_{aer}/m3_{air}
+                     *XM3TOUM3                     &!==>um3_{aer}/m3_{air}
+                     /(XPI*4./3.)                   !==>um3_{aer}/m3_{air}
             !==>volume 3rd moment   
-!RAIN: CALCULATE MOMENT 0 FROM DISPERSION AND MEAN RADIUS      
-     ZWORK31(:,:,:)= ZWORK31(:,:,:)/ &
-             ((ZRG_SLT(:,:,:,JJ)**3)*       &
-              EXP(4.5 * LOG(ZSIG_SLT(:,:,:,JJ))**2))
-!RAIN: RETURN TO CONCENTRATION #/m3       
-      ZWORK31(:,:,:)= ZWORK31(:,:,:) *   XMD/ &
-                (XAVOGADRO*XRHODREF(:,:,:))
-!RAIN: Get number concentration (#/molec_{air}==>#/m3)  
-             ZWORK31(:,:,:)=                         &
-             ZWORK31(:,:,:)                  & !#/molec_{air}
-              * XAVOGADRO                           & !==>#/mole
-              / XMD                                 & !==>#/kg_{air}
-              * XRHODREF(:,:,:)                       !==>#/m3
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
-! RAIN:   DUST MASS CONCENTRATION
-       WRITE(YRECFM,'(A9,I1)')'SLTDEPMSS',JJ+NMODE_SLT
-       WRITE(YCOMMENT,'(A17,I1,A7)')'DEPMASSCONC MODE ', &
-                                 JJ+NMODE_SLT,'(ug/m3)'
-       ILENCH=LEN(YCOMMENT)
-       ZWORK31(:,:,:)= ZWORK31(:,:,:)*4./3.*3.14*2500.*1e9 & !kg-->ug
-          * (ZRG_SLT(:,:,:,JJ)**3)*1.d-18 &  !um-->m
-          * exp(4.5*log(ZSIG_SLT(:,:,:,JJ))*log(ZSIG_SLT(:,:,:,JJ)))
-       CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
-     END DO       
-!
- END IF
+    !RAIN: CALCULATE MOMENT 0 FROM DISPERSION AND MEAN RADIUS
+    ZWORK31(:,:,:)= ZWORK31(:,:,:)/                        &
+                    ((ZRG_SLT(:,:,:,JJ)**3)*               &
+                    EXP(4.5 * LOG(ZSIG_SLT(:,:,:,JJ))**2))
+    !RAIN: RETURN TO CONCENTRATION #/m3
+    ZWORK31(:,:,:)= ZWORK31(:,:,:) *   XMD/ &
+                    (XAVOGADRO*XRHODREF(:,:,:))
+    !RAIN: Get number concentration (#/molec_{air}==>#/m3)
+    ZWORK31(:,:,:)=                   &
+                    ZWORK31(:,:,:)    & !#/molec_{air}
+                    * XAVOGADRO       & !==>#/mole
+                    / XMD             & !==>#/kg_{air}
+                    * XRHODREF(:,:,:)   !==>#/m3
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
+    ! RAIN:   DUST MASS CONCENTRATION
+    WRITE(TZFIELD%CMNHNAME,'(A9,I1)')'SLTDEPMSS',JJ+NMODE_SLT
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    WRITE(TZFIELD%CCOMMENT,'(A17,I1)')'DEPMASSCONC MODE ',JJ+NMODE_SLT
+    TZFIELD%CUNITS     = 'ug m-3'
+    ZWORK31(:,:,:)= ZWORK31(:,:,:)*4./3.*3.14*2500.*1e9 & !kg-->ug
+                    * (ZRG_SLT(:,:,:,JJ)**3)*1.d-18     &  !um-->m
+                    * exp(4.5*log(ZSIG_SLT(:,:,:,JJ))*log(ZSIG_SLT(:,:,:,JJ)))
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
+  END DO
+!
+END IF
 ! Dust variables
 IF (LDUST) THEN
   IF(.NOT.ALLOCATED(ZSIG_DST)) &
@@ -1471,54 +1538,56 @@ IF (LDUST) THEN
   IF(.NOT.ALLOCATED(ZN0_DST))  &
     ALLOCATE(ZN0_DST(SIZE(XSVT,1), SIZE(XSVT,2), SIZE(XSVT,3), NMODE_DST))
   !
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CUNITS     = 'ppbv'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  !
   DO JSV = NSV_DSTBEG,NSV_DSTEND
-    YRECFM=TRIM(UPCASE(CDUSTNAMES(JSV-NSV_DSTBEG+1)))//'T'
-    WRITE(YCOMMENT,'(A6,A4,I3.3,A7)')'X_Y_Z_','DUST',JSV,' (ppbv)'
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=XSVT(:,:,:,JSV)*1.E9
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = TRIM(UPCASE(CDUSTNAMES(JSV-NSV_DSTBEG+1)))//'T'
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    WRITE(TZFIELD%CCOMMENT,'(A6,A4,I3.3)')'X_Y_Z_','DUST',JSV
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XSVT(:,:,:,JSV)*1.E9)
   END DO
   !
   CALL PPP2DUST(XSVT(:,:,:,NSV_DSTBEG:NSV_DSTEND),XRHODREF,&
                PSIG3D=ZSIG_DST, PRG3D=ZRG_DST, PN3D=ZN0_DST)
   DO JJ=1,NMODE_DST
-    WRITE(YRECFM,'(A6,I1)')'DSTRGA',JJ
-    WRITE(YCOMMENT,'(A18,I1,A5)')'RG (nb) DUST MODE ',JJ,' (um)'
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=ZRG_DST(:,:,:,JJ)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'DSTRGA'
+    TZFIELD%CLONGNAME  = 'MesoNH: DSTRGA'
+    TZFIELD%CUNITS     = 'um'
+    WRITE(TZFIELD%CCOMMENT,'(A18,I1)')'RG (nb) DUST MODE ',JJ
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZRG_DST(:,:,:,JJ))
     !
-    WRITE(YRECFM,'(A7,I1)')'DSTRGAM',JJ
-    WRITE(YCOMMENT,'(A17,I1,A5)')'RG (m) DUST MODE ',JJ,' (um)'
-    ILENCH=LEN(YCOMMENT)
+    TZFIELD%CMNHNAME   = 'DSTRGAM'
+    TZFIELD%CLONGNAME  = 'MesoNH: DSTRGAM'
+    TZFIELD%CUNITS     = 'um'
+    WRITE(TZFIELD%CCOMMENT,'(A17,I1)')'RG (m) DUST MODE ',JJ
     ZWORK31(:,:,:)=ZRG_DST(:,:,:,JJ) / (EXP(-3.*(LOG(ZSIG_DST(:,:,:,JJ)))**2))
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
     !
-    WRITE(YRECFM,'(A6,I1)')'DSTN0A',JJ
-    WRITE(YCOMMENT,'(A13,I1,A7)')'N0 DUST MODE ',JJ,' (1/m3)'
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=ZN0_DST(:,:,:,JJ)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    WRITE(TZFIELD%CMNHNAME,'(A6,I1)')'DSTN0A',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'm-3'
+    WRITE(TZFIELD%CCOMMENT,'(A13,I1)')'N0 DUST MODE ',JJ
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZN0_DST(:,:,:,JJ))
     !
-    WRITE(YRECFM,'(A7,I1)')'DSTSIGA',JJ
-    WRITE(YCOMMENT,'(A16,I1)')'SIGMA DUST MODE ',JJ
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=ZSIG_DST(:,:,:,JJ)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    WRITE(TZFIELD%CMNHNAME,'(A7,I1)')'DSTSIGA',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = '1'
+    WRITE(TZFIELD%CCOMMENT,'(A16,I1)')'SIGMA DUST MODE ',JJ
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZSIG_DST(:,:,:,JJ))
     !DUST MASS CONCENTRATION
-    WRITE(YRECFM,'(A4,I1)')'DSTMSS',JJ
-    WRITE(YCOMMENT,'(A14,I1,A7)')'MASSCONC MODE ',JJ,'(ug/m3)'
-    ILENCH=LEN(YCOMMENT)
+    WRITE(TZFIELD%CMNHNAME,'(A4,I1)')'DSTMSS',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'ug m-3'
+    WRITE(TZFIELD%CCOMMENT,'(A14,I1)')'MASSCONC MODE ',JJ
     ZWORK31(:,:,:)= ZN0_DST(:,:,:,JJ)*4./3.*3.14*2500.*1e9 & !kg-->ug
        * (ZRG_DST(:,:,:,JJ)**3)*1.d-18 &  !um-->m
        * exp(4.5*log(ZSIG_DST(:,:,:,JJ))*log(ZSIG_DST(:,:,:,JJ)))
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
     !DUST BURDEN (g/m2)
     ZWORK21(:,:)=0.0
     DO JK=IKB,IKE
@@ -1532,115 +1601,124 @@ IF (LDUST) THEN
         ENDDO
       ENDDO
     ENDDO
-    WRITE(YRECFM,'(A7,I1)')'DSTBRDN',JJ
-    WRITE(YCOMMENT,'(A13,I1)')'BURDEN (g/m2)',JJ
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    WRITE(TZFIELD%CMNHNAME,'(A7,I1)')'DSTBRDN',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'g m-2'
+    WRITE(TZFIELD%CCOMMENT,'(A6,I1)')'BURDEN',JJ
+    TZFIELD%NDIMS      = 2
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK21)
+    !
+    TZFIELD%NDIMS      = 3
   ENDDO
 END IF
 IF (LDUST.AND.LDEPOS_DST(IMI)) THEN
-    ! 
-     ZSDSTDEP=XSVT(:,:,:,NSV_DSTDEPBEG:NSV_DSTDEPEND)
-     DO JSV = 1,NSV_DSTDEP   
-      YRECFM=TRIM(UPCASE(CDEDSTNAMES(JSV)))//'T'
-       WRITE(YCOMMENT,'(A6,A4,I3.3,A7)')'X_Y_Z_','DUSTDEP',JSV,' (ppbv)'
-       ILENCH=LEN(YCOMMENT)
-       ZWORK31(:,:,:)=ZSDSTDEP(:,:,:,JSV)*1.E9
-       CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
-     END DO
-     DO JJ=1,NMODE_DST           
-            ZWORK31(:,:,:)=0.0
-! FOR CLOUDS
-        WRITE(YRECFM,'(A9,I1)')'DSTDEPN0A',JJ
-        WRITE(YCOMMENT,'(A16,I1,A7)')'N0 DUSTDEP MODE ', &
-                                   JJ,' (1/m3)'            
-! CLOUD: CALCULATE MOMENT 3 FROM TOTAL AEROSOL MASS            
-            ZWORK31(:,:,:)=ZSDSTDEP(:,:,:,JJ)  &!==>molec_{aer}/molec_{air}
-            *(XMOLARWEIGHT_DUST/XMD)           &!==>kg_{aer}/kg_{air}
-            *XRHODREF(:,:,:)                   &!==>kg_{aer}/m3_{air}
-            *(1.d0/XDENSITY_DUST)              &!==>m3_{aer}/m3_{air}
-            *XM3TOUM3                          &!==>um3_{aer}/m3_{air}
-            /(XPI*4./3.)                        !==>um3_{aer}/m3_{air}
+  !
+  ZSDSTDEP=XSVT(:,:,:,NSV_DSTDEPBEG:NSV_DSTDEPEND)
+  !
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CUNITS     = 'ppbv'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  !
+  DO JSV = 1,NSV_DSTDEP
+    TZFIELD%CMNHNAME   = TRIM(UPCASE(CDEDSTNAMES(JSV)))//'T'
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    WRITE(TZFIELD%CCOMMENT,'(A6,A4,I3.3)')'X_Y_Z_','DUSTDEP',JSV
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZSDSTDEP(:,:,:,JSV)*1.E9)
+  END DO
+  !
+  DO JJ=1,NMODE_DST
+    ! FOR CLOUDS
+    WRITE(TZFIELD%CMNHNAME,'(A9,I1)')'DSTDEPN0A',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    WRITE(TZFIELD%CCOMMENT,'(A16,I1)')'N0 DUSTDEP MODE ',JJ
+    TZFIELD%CUNITS     = 'm-3'
+    ! CLOUD: CALCULATE MOMENT 3 FROM TOTAL AEROSOL MASS
+    ZWORK31(:,:,:) = ZSDSTDEP(:,:,:,JJ)         &!==>molec_{aer}/molec_{air}
+                     *(XMOLARWEIGHT_DUST/XMD)   &!==>kg_{aer}/kg_{air}
+                     *XRHODREF(:,:,:)           &!==>kg_{aer}/m3_{air}
+                     /XDENSITY_DUST             &!==>m3_{aer}/m3_{air}
+                     *XM3TOUM3                  &!==>um3_{aer}/m3_{air}
+                     /(XPI*4./3.)                !==>um3_{aer}/m3_{air}
             !==>volume 3rd moment
-!CLOUD: CALCULATE MOMENT 0 FROM DISPERSION AND MEAN RADIUS            
-            ZWORK31(:,:,:)=  ZWORK31(:,:,:)/      &
+    !CLOUD: CALCULATE MOMENT 0 FROM DISPERSION AND MEAN RADIUS
+    ZWORK31(:,:,:)=  ZWORK31(:,:,:)/      &
                     ((ZRG_DST(:,:,:,JJ)**3)*      &
                     EXP(4.5 * LOG(ZSIG_DST(:,:,:,JJ))**2))
-!CLOUD: RETURN TO CONCENTRATION #/m3
-            ZWORK31(:,:,:)= ZWORK31(:,:,:) *   XMD/ &
+    !CLOUD: RETURN TO CONCENTRATION #/m3
+    ZWORK31(:,:,:)= ZWORK31(:,:,:) *   XMD/ &
                      (XAVOGADRO*XRHODREF(:,:,:))
-!CLOUD:  Get number concentration (#/molec_{air}==>#/m3)  
-             ZWORK31(:,:,:)=                         &
+    !CLOUD:  Get number concentration (#/molec_{air}==>#/m3)
+    ZWORK31(:,:,:)=                         &
                     ZWORK31(:,:,:)                  & !#/molec_{air}
                     * XAVOGADRO                     & !==>#/mole
                     / XMD                           & !==>#/kg_{air}
                     * XRHODREF(:,:,:)                 !==>#/m3  
-     ILENCH=LEN(YCOMMENT)               
-     CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
-! CLOUD:   DUST MASS CONCENTRATION
-       WRITE(YRECFM,'(A9,I1)')'DSTDEPMSS',JJ
-       WRITE(YCOMMENT,'(A17,I1,A7)')'DEPMASSCONC MODE ', &
-                                 JJ,'(ug/m3)'
-       ILENCH=LEN(YCOMMENT)
-       ZWORK31(:,:,:)= ZWORK31(:,:,:)*4./3.*3.14*2500.*1e9 & !kg-->ug
-          * (ZRG_DST(:,:,:,JJ)**3)*1.d-18 &  !um-->m
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
+    ! CLOUD:   DUST MASS CONCENTRATION
+    WRITE(TZFIELD%CMNHNAME,'(A9,I1)')'DSTDEPMSS',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    WRITE(TZFIELD%CCOMMENT,'(A17,I1)')'DEPMASSCONC MODE ',JJ
+    TZFIELD%CUNITS     = 'ug m-3'
+    ZWORK31(:,:,:)= ZWORK31(:,:,:)*4./3.*3.14*2500.*1e9 & !kg-->ug
+          * (ZRG_DST(:,:,:,JJ)**3)*1.d-18               &  !um-->m
           * exp(4.5*log(ZSIG_DST(:,:,:,JJ))*log(ZSIG_DST(:,:,:,JJ)))
-       CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)                
-!   FOR RAIN DROPS                
-        WRITE(YRECFM,'(A9,I1)')'DSTDEPN0A',JJ+NMODE_DST
-        WRITE(YCOMMENT,'(A16,I1,A7)')'N0 DUSTDEP MODE ', &
-                                   JJ+NMODE_DST,' (1/m3)'
-        ILENCH=LEN(YCOMMENT) 
-        ZWORK31(:,:,:)=0.0
-! RAIN: CALCULATE MOMENT 3 FROM TOTAL AEROSOL MASS                
-        ZWORK31(:,:,:)=ZSDSTDEP(:,:,:,JJ+NMODE_DST)  &!==>molec_{aer}/molec_{air}
-            *(XMOLARWEIGHT_DUST/XMD)           &!==>kg_{aer}/kg_{air}
-            *XRHODREF(:,:,:)                   &!==>kg_{aer}/m3_{air}
-            *(1.d0/XDENSITY_DUST)              &!==>m3_{aer}/m3_{air}
-            *XM3TOUM3                          &!==>um3_{aer}/m3_{air}
-            /(XPI*4./3.)                        !==>um3_{aer}/m3_{air}
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
+    !   FOR RAIN DROPS
+    WRITE(TZFIELD%CMNHNAME,'(A9,I1)')'DSTDEPN0A',JJ+NMODE_DST
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    WRITE(TZFIELD%CCOMMENT,'(A16,I1)')'N0 DUSTDEP MODE ',JJ+NMODE_DST
+    TZFIELD%CUNITS     = 'm-3'
+    ! RAIN: CALCULATE MOMENT 3 FROM TOTAL AEROSOL MASS
+    ZWORK31(:,:,:)=ZSDSTDEP(:,:,:,JJ+NMODE_DST)  &!==>molec_{aer}/molec_{air}
+            *(XMOLARWEIGHT_DUST/XMD)             &!==>kg_{aer}/kg_{air}
+            *XRHODREF(:,:,:)                     &!==>kg_{aer}/m3_{air}
+            *(1.d0/XDENSITY_DUST)                &!==>m3_{aer}/m3_{air}
+            *XM3TOUM3                            &!==>um3_{aer}/m3_{air}
+            /(XPI*4./3.)                          !==>um3_{aer}/m3_{air}
             !==>volume 3rd moment   
-!RAIN: CALCULATE MOMENT 0 FROM DISPERSION AND MEAN RADIUS      
-     ZWORK31(:,:,:)= ZWORK31(:,:,:)/ &
-             ((ZRG_DST(:,:,:,JJ)**3)*       &
+    !RAIN: CALCULATE MOMENT 0 FROM DISPERSION AND MEAN RADIUS
+    ZWORK31(:,:,:)= ZWORK31(:,:,:)/                 &
+             ((ZRG_DST(:,:,:,JJ)**3)*               &
               EXP(4.5 * LOG(ZSIG_DST(:,:,:,JJ))**2))
-!RAIN: RETURN TO CONCENTRATION #/m3       
-      ZWORK31(:,:,:)= ZWORK31(:,:,:) *   XMD/ &
-                (XAVOGADRO*XRHODREF(:,:,:))
-!RAIN: Get number concentration (#/molec_{air}==>#/m3)  
-             ZWORK31(:,:,:)=                         &
-             ZWORK31(:,:,:)                  & !#/molec_{air}
-              * XAVOGADRO                           & !==>#/mole
-              / XMD                                 & !==>#/kg_{air}
-              * XRHODREF(:,:,:)                       !==>#/m3
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
-! RAIN:   DUST MASS CONCENTRATION
-       WRITE(YRECFM,'(A9,I1)')'DSTDEPMSS',JJ+NMODE_DST
-       WRITE(YCOMMENT,'(A17,I1,A7)')'DEPMASSCONC MODE ', &
-                                 JJ+NMODE_DST,'(ug/m3)'
-       ILENCH=LEN(YCOMMENT)
-       ZWORK31(:,:,:)= ZWORK31(:,:,:)*4./3.*3.14*2500.*1e9 & !kg-->ug
-          * (ZRG_DST(:,:,:,JJ)**3)*1.d-18 &  !um-->m
-          * exp(4.5*log(ZSIG_DST(:,:,:,JJ))*log(ZSIG_DST(:,:,:,JJ)))
-       CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
-     END DO       
+    !RAIN: RETURN TO CONCENTRATION #/m3
+    ZWORK31(:,:,:)= ZWORK31(:,:,:) *   XMD/ &
+                    (XAVOGADRO*XRHODREF(:,:,:))
+    !RAIN: Get number concentration (#/molec_{air}==>#/m3)
+    ZWORK31(:,:,:)=                   &
+                    ZWORK31(:,:,:)    & !#/molec_{air}
+                    * XAVOGADRO       & !==>#/mole
+                    / XMD             & !==>#/kg_{air}
+                    * XRHODREF(:,:,:)   !==>#/m3
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
+    ! RAIN:   DUST MASS CONCENTRATION
+    WRITE(TZFIELD%CMNHNAME,'(A9,I1)')'DSTDEPMSS',JJ+NMODE_DST
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    WRITE(TZFIELD%CCOMMENT,'(A17,I1)')'DEPMASSCONC MODE ',JJ+NMODE_DST
+    TZFIELD%CUNITS     = 'ug m-3'
+    ZWORK31(:,:,:)= ZWORK31(:,:,:)*4./3.*3.14*2500.*1e9 & !kg-->ug
+                    * (ZRG_DST(:,:,:,JJ)**3)*1.d-18     &  !um-->m
+                    * exp(4.5*log(ZSIG_DST(:,:,:,JJ))*log(ZSIG_DST(:,:,:,JJ)))
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
+  END DO
 !
- END IF
+END IF
 ! Aerosol
 IF ((LCHEMDIAG).AND.(LORILAM).AND.(LUSECHEM)) THEN
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CUNITS     = 'ppbv'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  !
   DO JSV = NSV_AERBEG,NSV_AEREND
-    YRECFM=TRIM(UPCASE(CAERONAMES(JSV-NSV_AERBEG+1)))//'T'
-    WRITE(YCOMMENT,'(A6,A4,I3.3,A7)')'X_Y_Z_','AERO',JSV,' (ppbv)'
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=XSVT(:,:,:,JSV)*1.E9
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-         YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = TRIM(UPCASE(CAERONAMES(JSV-NSV_AERBEG+1)))//'T'
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    WRITE(TZFIELD%CCOMMENT,'(A6,A4,I3.3)')'X_Y_Z_','AERO',JSV
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XSVT(:,:,:,JSV)*1.E9)
   END DO
   !
   IF (.NOT.(ASSOCIATED(XN3D)))   &
@@ -1653,164 +1731,136 @@ IF ((LCHEMDIAG).AND.(LORILAM).AND.(LUSECHEM)) THEN
   CALL  PPP2AERO(XSVT(:,:,:,NSV_AERBEG:NSV_AEREND), XRHODREF, &
                  PSIG3D=XSIG3D, PRG3D=XRG3D, PN3D=XN3D, PCTOTA=ZPTOTA) 
   DO JJ=1,JPMODE
-    WRITE(YRECFM,'(A3,I1)')'RGA',JJ
-    WRITE(YCOMMENT,'(A21,I1,A5)')'RG (nb) AEROSOL MODE ',JJ,' (um)'
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=XRG3D(:,:,:,JJ)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'RGA'
+    TZFIELD%CLONGNAME  = 'MesoNH: RGA'
+    TZFIELD%CUNITS     = 'um'
+    WRITE(TZFIELD%CCOMMENT,'(A21,I1)')'RG (nb) AEROSOL MODE ',JJ
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XRG3D(:,:,:,JJ))
     !
-    WRITE(YRECFM,'(A4,I1)')'RGAM',JJ
-    WRITE(YCOMMENT,'(A20,I1,A5)')'RG (m) AEROSOL MODE ',JJ,' (um)'
-    ILENCH=LEN(YCOMMENT)
+    TZFIELD%CMNHNAME   = 'RGAM'
+    TZFIELD%CLONGNAME  = 'MesoNH: RGAM'
+    TZFIELD%CUNITS     = 'um'
+    WRITE(TZFIELD%CCOMMENT,'(A20,I1)')'RG (m) AEROSOL MODE ',JJ
     ZWORK31(:,:,:)=XRG3D(:,:,:,JJ) / (EXP(-3.*(LOG(XSIG3D(:,:,:,JJ)))**2))
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
     !
-    WRITE(YRECFM,'(A3,I1)')'N0A',JJ
-    WRITE(YCOMMENT,'(A16,I1,A7)')'N0 AEROSOL MODE ',JJ,' (1/cc)'
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=XN3D(:,:,:,JJ)*1.E-6
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    WRITE(TZFIELD%CMNHNAME,'(A3,I1)')'N0A',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'cm-3'
+    WRITE(TZFIELD%CCOMMENT,'(A16,I1)')'N0 AEROSOL MODE ',JJ
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XN3D(:,:,:,JJ)*1.E-6)
     !
-    WRITE(YRECFM,'(A4,I1)')'SIGA',JJ
-    WRITE(YCOMMENT,'(A19,I1)')'SIGMA AEROSOL MODE ',JJ
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=XSIG3D(:,:,:,JJ)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    WRITE(TZFIELD%CMNHNAME,'(A4,I1)')'SIGA',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = '1'
+    WRITE(TZFIELD%CCOMMENT,'(A19,I1)')'SIGMA AEROSOL MODE ',JJ
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XSIG3D(:,:,:,JJ))
     !
-    WRITE(YRECFM,'(A4,I1)')'MSO4',JJ
-    WRITE(YCOMMENT,'(A22,I1,A5)')'MASS SO4 AEROSOL MODE ',JJ,'(ug/m3)'
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=ZPTOTA(:,:,:,JP_AER_SO4,JJ)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    WRITE(TZFIELD%CMNHNAME,'(A4,I1)')'MSO4',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'ug m-3'
+    WRITE(TZFIELD%CCOMMENT,'(A22,I1)')'MASS SO4 AEROSOL MODE ',JJ
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZPTOTA(:,:,:,JP_AER_SO4,JJ))
     !
-    WRITE(YRECFM,'(A4,I1)')'MNO3',JJ
-    WRITE(YCOMMENT,'(A22,I1,A5)')'MASS NO3 AEROSOL MODE ',JJ,'(ug/m3)'
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=ZPTOTA(:,:,:,JP_AER_NO3,JJ)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    WRITE(TZFIELD%CMNHNAME,'(A4,I1)')'MNO3',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'ug m-3'
+    WRITE(TZFIELD%CCOMMENT,'(A22,I1)')'MASS NO3 AEROSOL MODE ',JJ
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZPTOTA(:,:,:,JP_AER_NO3,JJ))
     !
-    WRITE(YRECFM,'(A4,I1)')'MNH3',JJ
-    WRITE(YCOMMENT,'(A22,I1,A5)')'MASS NH3 AEROSOL MODE ',JJ,'(ug/m3)'
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=ZPTOTA(:,:,:,JP_AER_NH3,JJ)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    WRITE(TZFIELD%CMNHNAME,'(A4,I1)')'MNH3',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'ug m-3'
+    WRITE(TZFIELD%CCOMMENT,'(A22,I1)')'MASS NH3 AEROSOL MODE ',JJ
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZPTOTA(:,:,:,JP_AER_NH3,JJ))
     !
-    WRITE(YRECFM,'(A4,I1)')'MH2O',JJ
-    WRITE(YCOMMENT,'(A22,I1,A5)') 'MASS H2O AEROSOL MODE ',JJ,'(ug/m3)'
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=ZPTOTA(:,:,:,JP_AER_H2O,JJ)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    WRITE(TZFIELD%CMNHNAME,'(A4,I1)')'MH2O',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'ug m-3'
+    WRITE(TZFIELD%CCOMMENT,'(A22,I1)')'MASS H2O AEROSOL MODE ',JJ
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZPTOTA(:,:,:,JP_AER_H2O,JJ))
     !
     IF (NSOA .EQ. 10) THEN
-      WRITE(YRECFM,'(A5,I1)')'MSOA1',JJ
-      WRITE(YCOMMENT,'(A23,I1,A5)')'MASS SOA1 AEROSOL MODE ',JJ,'(ug/m3)'
-      ILENCH=LEN(YCOMMENT)
-      ZWORK31(:,:,:)=ZPTOTA(:,:,:,JP_AER_SOA1,JJ)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                  YCOMMENT,IRESP)
+      WRITE(TZFIELD%CMNHNAME,'(A5,I1)')'MSOA1',JJ
+      TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+      TZFIELD%CUNITS     = 'ug m-3'
+      WRITE(TZFIELD%CCOMMENT,'(A23,I1)')'MASS SOA1 AEROSOL MODE ',JJ
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZPTOTA(:,:,:,JP_AER_SOA1,JJ))
       !
-      WRITE(YRECFM,'(A5,I1)')'MSOA2',JJ
-      WRITE(YCOMMENT,'(A23,I1,A5)')'MASS SOA2 AEROSOL MODE ',JJ,'(ug/m3)'
-      ILENCH=LEN(YCOMMENT)
-      ZWORK31(:,:,:)=ZPTOTA(:,:,:,JP_AER_SOA2,JJ)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                  YCOMMENT,IRESP)
+      WRITE(TZFIELD%CMNHNAME,'(A5,I1)')'MSOA2',JJ
+      TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+      TZFIELD%CUNITS     = 'ug m-3'
+      WRITE(TZFIELD%CCOMMENT,'(A23,I1)')'MASS SOA2 AEROSOL MODE ',JJ
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZPTOTA(:,:,:,JP_AER_SOA2,JJ))
       !
-      WRITE(YRECFM,'(A5,I1)')'MSOA3',JJ
-      WRITE(YCOMMENT,'(A23,I1,A5)')'MASS SOA3 AEROSOL MODE ',JJ,'(ug/m3)'
-      ILENCH=LEN(YCOMMENT)
-      ZWORK31(:,:,:)=ZPTOTA(:,:,:,JP_AER_SOA3,JJ)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-           YCOMMENT,IRESP)
+      WRITE(TZFIELD%CMNHNAME,'(A5,I1)')'MSOA3',JJ
+      TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+      TZFIELD%CUNITS     = 'ug m-3'
+      WRITE(TZFIELD%CCOMMENT,'(A23,I1)')'MASS SOA3 AEROSOL MODE ',JJ
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZPTOTA(:,:,:,JP_AER_SOA3,JJ))
       !
-      WRITE(YRECFM,'(A5,I1)')'MSOA4',JJ
-      WRITE(YCOMMENT,'(A23,I1,A5)')'MASS SOA4 AEROSOL MODE ',JJ,'(ug/m3)'
-      ILENCH=LEN(YCOMMENT)
-      ZWORK31(:,:,:)=ZPTOTA(:,:,:,JP_AER_SOA4,JJ)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-           YCOMMENT,IRESP)
+      WRITE(TZFIELD%CMNHNAME,'(A5,I1)')'MSOA4',JJ
+      TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+      TZFIELD%CUNITS     = 'ug m-3'
+      WRITE(TZFIELD%CCOMMENT,'(A23,I1)')'MASS SOA4 AEROSOL MODE ',JJ
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZPTOTA(:,:,:,JP_AER_SOA4,JJ))
       !
-      WRITE(YRECFM,'(A5,I1)')'MSOA5',JJ
-      WRITE(YCOMMENT,'(A23,I1,A5)')'MASS SOA5 AEROSOL MODE ',JJ,'(ug/m3)'
-      ILENCH=LEN(YCOMMENT)
-      ZWORK31(:,:,:)=ZPTOTA(:,:,:,JP_AER_SOA5,JJ)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                  YCOMMENT,IRESP)
+      WRITE(TZFIELD%CMNHNAME,'(A5,I1)')'MSOA5',JJ
+      TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+      TZFIELD%CUNITS     = 'ug m-3'
+      WRITE(TZFIELD%CCOMMENT,'(A23,I1)')'MASS SOA5 AEROSOL MODE ',JJ
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZPTOTA(:,:,:,JP_AER_SOA5,JJ))
       !
-      WRITE(YRECFM,'(A5,I1)')'MSOA6',JJ
-      WRITE(YCOMMENT,'(A23,I1,A5)')'MASS SOA6 AEROSOL MODE ',JJ,'(ug/m3)'
-      ILENCH=LEN(YCOMMENT)
-      ZWORK31(:,:,:)=ZPTOTA(:,:,:,JP_AER_SOA6,JJ)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                  YCOMMENT,IRESP)
+      WRITE(TZFIELD%CMNHNAME,'(A5,I1)')'MSOA6',JJ
+      TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+      TZFIELD%CUNITS     = 'ug m-3'
+      WRITE(TZFIELD%CCOMMENT,'(A23,I1)')'MASS SOA6 AEROSOL MODE ',JJ
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZPTOTA(:,:,:,JP_AER_SOA6,JJ))
       !
-      WRITE(YRECFM,'(A5,I1)')'MSOA7',JJ
-      WRITE(YCOMMENT,'(A23,I1,A5)')'MASS SOA7 AEROSOL MODE ',JJ,'(ug/m3)'
-      ILENCH=LEN(YCOMMENT)
-      ZWORK31(:,:,:)=ZPTOTA(:,:,:,JP_AER_SOA7,JJ)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                  YCOMMENT,IRESP)
+      WRITE(TZFIELD%CMNHNAME,'(A5,I1)')'MSOA7',JJ
+      TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+      TZFIELD%CUNITS     = 'ug m-3'
+      WRITE(TZFIELD%CCOMMENT,'(A23,I1)')'MASS SOA7 AEROSOL MODE ',JJ
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZPTOTA(:,:,:,JP_AER_SOA7,JJ))
       !
-      WRITE(YRECFM,'(A5,I1)')'MSOA8',JJ
-      WRITE(YCOMMENT,'(A23,I1,A5)')'MASS SOA8 AEROSOL MODE ',JJ,'(ug/m3)'
-      ILENCH=LEN(YCOMMENT)
-      ZWORK31(:,:,:)=ZPTOTA(:,:,:,JP_AER_SOA8,JJ)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                  YCOMMENT,IRESP)
+      WRITE(TZFIELD%CMNHNAME,'(A5,I1)')'MSOA8',JJ
+      TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+      TZFIELD%CUNITS     = 'ug m-3'
+      WRITE(TZFIELD%CCOMMENT,'(A23,I1)')'MASS SOA8 AEROSOL MODE ',JJ
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZPTOTA(:,:,:,JP_AER_SOA8,JJ))
       !
-      WRITE(YRECFM,'(A5,I1)')'MSOA9',JJ
-      WRITE(YCOMMENT,'(A23,I1,A5)')'MASS SOA9 AEROSOL MODE ',JJ,'(ug/m3)'
-      ILENCH=LEN(YCOMMENT)
-      ZWORK31(:,:,:)=ZPTOTA(:,:,:,JP_AER_SOA9,JJ)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                  YCOMMENT,IRESP)
+      WRITE(TZFIELD%CMNHNAME,'(A5,I1)')'MSOA9',JJ
+      TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+      TZFIELD%CUNITS     = 'ug m-3'
+      WRITE(TZFIELD%CCOMMENT,'(A23,I1)')'MASS SOA9 AEROSOL MODE ',JJ
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZPTOTA(:,:,:,JP_AER_SOA9,JJ))
       !
-      WRITE(YRECFM,'(A6,I1)')'MSOA10',JJ
-      WRITE(YCOMMENT,'(A23,I1,A5)')'MASS SOA10 AEROSOL MODE ',JJ,'(ug/m3)'
-      ILENCH=LEN(YCOMMENT)
-      ZWORK31(:,:,:)=ZPTOTA(:,:,:,JP_AER_SOA10,JJ)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                  YCOMMENT,IRESP)
+      WRITE(TZFIELD%CMNHNAME,'(A6,I1)')'MSOA10',JJ
+      TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+      TZFIELD%CUNITS     = 'ug m-3'
+      WRITE(TZFIELD%CCOMMENT,'(A24,I1)')'MASS SOA10 AEROSOL MODE ',JJ
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZPTOTA(:,:,:,JP_AER_SOA10,JJ))
     END IF
     !
-    WRITE(YRECFM,'(A3,I1)')'MOC',JJ
-    WRITE(YCOMMENT,'(A21,I1,A5)')'MASS OC AEROSOL MODE ',JJ,'(ug/m3)'
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=ZPTOTA(:,:,:,JP_AER_OC,JJ)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    WRITE(TZFIELD%CMNHNAME,'(A3,I1)')'MOC',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'ug m-3'
+    WRITE(TZFIELD%CCOMMENT,'(A21,I1)')'MASS OC AEROSOL MODE ',JJ
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZPTOTA(:,:,:,JP_AER_OC,JJ))
     !
-    WRITE(YRECFM,'(A3,I1)')'MBC',JJ
-    WRITE(YCOMMENT,'(A21,I1,A5)')'MASS BC AEROSOL MODE ',JJ,'(ug/m3)'
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=ZPTOTA(:,:,:,JP_AER_BC,JJ)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,    &
-                YCOMMENT,IRESP)
+    WRITE(TZFIELD%CMNHNAME,'(A3,I1)')'MBC',JJ
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'ug m-3'
+    WRITE(TZFIELD%CCOMMENT,'(A21,I1)')'MASS BC AEROSOL MODE ',JJ
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZPTOTA(:,:,:,JP_AER_BC,JJ))
   ENDDO
 END IF
 !
 !* Large Scale variables
 !
 IF (LVAR_LS) THEN
-  YRECFM='LSUM'
-  YCOMMENT='X_Y_Z_Large Scale U component (M/S)'
-  IGRID=2
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',XLSUM,IGRID,ILENCH,YCOMMENT,IRESP)
-!
-  YRECFM='LSVM'
-  YCOMMENT='X_Y_Z_Large Scale V component (M/S)'
-  IGRID=3
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',XLSVM,IGRID,ILENCH,YCOMMENT,IRESP)
+  CALL IO_WRITE_FIELD(TPFILE,'LSUM', CLUOUT,IRESP,XLSUM)
+  CALL IO_WRITE_FIELD(TPFILE,'LSVM', CLUOUT,IRESP,XLSVM)
   !
   IF (LWIND_ZM) THEN
     TZFIELD2(1)%CMNHNAME   = 'LSUM_ZM'
@@ -1836,25 +1886,14 @@ IF (LVAR_LS) THEN
     CALL UV_TO_ZONAL_AND_MERID(XLSUM,XLSVM,23,TPFILE=TPFILE,TZFIELDS=TZFIELD2)
   ENDIF
   !
-  YRECFM='LSWM'
-  YCOMMENT='X_Y_Z_Large Scale vertical wind (M/S)'
-  IGRID=4
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',XLSWM,IGRID,ILENCH,YCOMMENT,IRESP)
-  !
-  YRECFM='LSTHM'
-  YCOMMENT='X_Y_Z_Large Scale potential Temperature (K)'
-  IGRID=1
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',XLSTHM,IGRID,ILENCH,YCOMMENT,IRESP)
+  CALL IO_WRITE_FIELD(TPFILE,'LSWM', CLUOUT,IRESP,XLSWM)
+  CALL IO_WRITE_FIELD(TPFILE,'LSTHM',CLUOUT,IRESP,XLSTHM)
 !
   IF (LUSERV) THEN
-    YRECFM='LSMRV'
-    YCOMMENT='X_Y_Z_Large Scale Vapor Mixing Ratio (G/KG)'
-    IGRID=1
-    ILENCH=LEN(YCOMMENT)
-    ZWORK31(:,:,:)=XLSRVM(:,:,:)*1.E3
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,YCOMMENT,IRESP)
+    CALL FIND_FIELD_ID_FROM_MNHNAME('LSRVM',IID,IRESP)
+    TZFIELD = TFIELDLIST(IID)
+    TZFIELD%CUNITS = 'g kg-1'
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XLSRVM(:,:,:)*1.E3)
   END IF
 END IF
 !
@@ -1865,75 +1904,115 @@ IF (LVAR_FRC .AND. LFORCING) THEN
   DO JT=1,NFRC
     WRITE (YFRC,'(I3.3)') JT
 !
-    YRECFM='UFRC'//YFRC
-    YCOMMENT=' '
-    IGRID=1
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'--',XUFRC(:,JT),IGRID,ILENCH,     &
-                                                            YCOMMENT,IRESP)
-!
-    YRECFM='VFRC'//YFRC
-    YCOMMENT=' '
-    IGRID=1
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'--',XVFRC(:,JT),IGRID,ILENCH,     &
-                                                            YCOMMENT,IRESP)
-!
-    YRECFM='WFRC'//YFRC
-    YCOMMENT=' '
-    IGRID=4
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'--',XWFRC(:,JT),IGRID,ILENCH,     &
-                                                            YCOMMENT,IRESP)
-!
-    YRECFM='THFRC'//YFRC
-    YCOMMENT=' '
-    IGRID=1
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'--',XTHFRC(:,JT),IGRID,ILENCH,    &
-                                                            YCOMMENT,IRESP)
-!
-    YRECFM='RVFRC'//YFRC
-    YCOMMENT=' '
-    IGRID=1
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'--',XRVFRC(:,JT),IGRID,ILENCH,    &
-                                                            YCOMMENT,IRESP)
-!
-    YRECFM='TENDTHFRC'//YFRC
-    YCOMMENT=' '
-    IGRID=1
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'--',XTENDTHFRC(:,JT),IGRID,ILENCH,  &
-                                                            YCOMMENT,IRESP)
-!
-    YRECFM='TENDRVFRC'//YFRC
-    YCOMMENT=' '
-    IGRID=1
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'--',XTENDRVFRC(:,JT),IGRID,ILENCH,  &
-                                                            YCOMMENT,IRESP)
-!
-    YRECFM='GXTHFRC'//YFRC
-    YCOMMENT=' '
-    IGRID=1
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'--',XGXTHFRC(:,JT),IGRID,ILENCH,  &
-                                                            YCOMMENT,IRESP)
-!
-    YRECFM='GYTHFRC'//YFRC
-    YCOMMENT=' '
-    IGRID=1
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'--',XGYTHFRC(:,JT),IGRID,ILENCH,  &
-                                                            YCOMMENT,IRESP)
-!
-    YRECFM='PGROUNDFRC'//YFRC
-    YCOMMENT=' '
-    IGRID=0
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'--',XPGROUNDFRC(JT),IGRID,ILENCH,  &
-                                                            YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'UFRC'//YFRC
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'm s-1'
+    TZFIELD%CDIR       = '--'
+    TZFIELD%CCOMMENT   = 'Zonal component of horizontal forcing wind'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 1
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XUFRC(:,JT))
+!
+    TZFIELD%CMNHNAME   = 'VFRC'//YFRC
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'm s-1'
+    TZFIELD%CDIR       = '--'
+    TZFIELD%CCOMMENT   = 'Meridian component of horizontal forcing wind'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 1
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XVFRC(:,JT))
+!
+    TZFIELD%CMNHNAME   = 'WFRC'//YFRC
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'm s-1'
+    TZFIELD%CDIR       = '--'
+    TZFIELD%CCOMMENT   = 'Vertical forcing wind'
+    TZFIELD%NGRID      = 4
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 1
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XWFRC(:,JT))
+!
+    TZFIELD%CMNHNAME   = 'THFRC'//YFRC
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'K'
+    TZFIELD%CDIR       = '--'
+    TZFIELD%CCOMMENT   = 'Forcing potential temperature'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 1
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XTHFRC(:,JT))
+!
+    TZFIELD%CMNHNAME   = 'RVFRC'//YFRC
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'kg kg-1'
+    TZFIELD%CDIR       = '--'
+    TZFIELD%CCOMMENT   = 'Forcing vapor mixing ratio'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 1
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XRVFRC(:,JT))
+!
+    TZFIELD%CMNHNAME   = 'TENDTHFRC'//YFRC
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'K s-1'
+    TZFIELD%CDIR       = '--'
+    TZFIELD%CCOMMENT   = 'Large-scale potential temperature tendency for forcing'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 1
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XTENDTHFRC(:,JT))
+!
+    TZFIELD%CMNHNAME   = 'TENDRVFRC'//YFRC
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'kg kg-1 s-1'
+    TZFIELD%CDIR       = '--'
+    TZFIELD%CCOMMENT   = 'Large-scale vapor mixing ratio tendency for forcing'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 1
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XTENDRVFRC(:,JT))
+!
+    TZFIELD%CMNHNAME   = 'GXTHFRC'//YFRC
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'K m-1'
+    TZFIELD%CDIR       = '--'
+    TZFIELD%CCOMMENT   = 'Large-scale potential temperature gradient for forcing'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 1
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XGXTHFRC(:,JT))
+!
+    TZFIELD%CMNHNAME   = 'GYTHFRC'//YFRC
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'K m-1'
+    TZFIELD%CDIR       = '--'
+    TZFIELD%CCOMMENT   = 'Large-scale potential temperature gradient for forcing'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 1
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XGYTHFRC(:,JT))
+!
+    TZFIELD%CMNHNAME   = 'PGROUNDFRC'//YFRC
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+    TZFIELD%CUNITS     = 'Pa'
+    TZFIELD%CDIR       = '--'
+    TZFIELD%CCOMMENT   = 'Forcing ground pressure'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 0
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XPGROUNDFRC(JT))
 !
   END DO
 END IF
@@ -2009,21 +2088,31 @@ IF (LTPZH .OR. LCOREF) THEN
       ZWORK33(:,:,:)=(77.6*( XPABST(:,:,:)*1E-2                &
                             +ZWORK33(:,:,:)*4810/ZTEMP(:,:,:)) &
                       -6*ZWORK33(:,:,:)                        )/ZTEMP(:,:,:)
-      YRECFM='COREF'
-      YCOMMENT='X_Y_Z_REFraction COindex (N-units)'
-      IGRID=1
-      ILENCH=LEN(YCOMMENT)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK33,IGRID,ILENCH,YCOMMENT,IRESP)
+      TZFIELD%CMNHNAME   = 'COREF'
+      TZFIELD%CSTDNAME   = ''
+      TZFIELD%CLONGNAME  = 'MesoNH: COREF'
+      TZFIELD%CUNITS     = '1'
+      TZFIELD%CDIR       = 'XY'
+      TZFIELD%CCOMMENT   = 'X_Y_Z_REFraction COindex (N-units)'
+      TZFIELD%NGRID      = 1
+      TZFIELD%NTYPE      = TYPEREAL
+      TZFIELD%NDIMS      = 3
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK33)
       !
       ZWORK33(:,:,:)=ZWORK33(:,:,:)+MZF(1,IKU,1,XZZ(:,:,:))*1E6/XRADIUS
-      YRECFM='MCOREF'
-      YCOMMENT='X_Y_Z_Modified REFraction COindex (M-units)'
-      IGRID=1
-      ILENCH=LEN(YCOMMENT)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK33,IGRID,ILENCH,YCOMMENT,IRESP)
+      TZFIELD%CMNHNAME   = 'MCOREF'
+      TZFIELD%CSTDNAME   = ''
+      TZFIELD%CLONGNAME  = 'MesoNH: MCOREF'
+      TZFIELD%CUNITS     = '1'
+      TZFIELD%CDIR       = 'XY'
+      TZFIELD%CCOMMENT   = 'X_Y_Z_Modified REFraction COindex (M-units)'
+      TZFIELD%NGRID      = 1
+      TZFIELD%NTYPE      = TYPEREAL
+      TZFIELD%NDIMS      = 3
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK33)
     END IF
   ELSE
-    PRINT*, 'NO WATER VAPOR IN ',YFMFILE2,' RELATIVE HUMIDITY IS NOT COMPUTED'
+    PRINT*, 'NO WATER VAPOR IN ',TPFILE%CNAME,' RELATIVE HUMIDITY IS NOT COMPUTED'
   END IF
 !
 END IF
@@ -2051,11 +2140,16 @@ IF ( LMOIST_V .OR. LMSLP .OR. LBLTOP ) THEN
 !
   IF (LMOIST_V .AND. NRR > 0) THEN
 ! Virtual potential temperature
-    YRECFM='THETAV'
-    YCOMMENT='X_Y_Z_Virtual potential temperature (K)'
-    IGRID=1
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZTHETAV,IGRID,ILENCH,YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'THETAV'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: THETAV'
+    TZFIELD%CUNITS     = 'K'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_Z_Virtual potential temperature'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 3
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZTHETAV)
   END IF
 !
 END IF
@@ -2078,11 +2172,16 @@ IF (LVISI) THEN
      ZVISIKUN(:,:,:) =0.027/(XRT(:,:,:,2)*XRHODREF(:,:,:))**0.88*1000.
     END WHERE
 ! Visibity Kunkel                     
-      YRECFM='VISIKUN'
-      YCOMMENT='X_Y_Z_Visibility Kunkel (m)'
-      IGRID=1
-      ILENCH=LEN(YCOMMENT)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZVISIKUN,IGRID,ILENCH,YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'VISIKUN'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: VISIKUN'
+    TZFIELD%CUNITS     = 'm'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_Z_Visibility Kunkel'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 3
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZVISIKUN)
 !
     IF ((CCLOUD == 'C2R2') .OR. (CCLOUD =='KHKO')) THEN
       ZVISIGUL(:,:,:) = 10000.
@@ -2092,17 +2191,27 @@ IF (LVISI) THEN
        ZVISIZHA(:,:,:) =0.187/(XRT(:,:,:,2)*XRHODREF(:,:,:)*XSVT(:,:,:,NSV_C2R2BEG+1))**0.34*1000.
       END WHERE
 ! Visibity Gultepe                    
-      YRECFM='VISIGUL'
-      YCOMMENT='X_Y_Z_Visibility Gultepe (m)'
-      IGRID=1
-      ILENCH=LEN(YCOMMENT)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZVISIGUL,IGRID,ILENCH,YCOMMENT,IRESP)
+      TZFIELD%CMNHNAME   = 'VISIGUL'
+      TZFIELD%CSTDNAME   = ''
+      TZFIELD%CLONGNAME  = 'MesoNH: VISIGUL'
+      TZFIELD%CUNITS     = 'm'
+      TZFIELD%CDIR       = 'XY'
+      TZFIELD%CCOMMENT   = 'X_Y_Z_Visibility Gultepe'
+      TZFIELD%NGRID      = 1
+      TZFIELD%NTYPE      = TYPEREAL
+      TZFIELD%NDIMS      = 3
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZVISIGUL)
 ! Visibity Zhang                      
-      YRECFM='VISIZHA'
-      YCOMMENT='X_Y_Z_Visibility Zhang (m)'
-      IGRID=1
-      ILENCH=LEN(YCOMMENT)
-      CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZVISIZHA,IGRID,ILENCH,YCOMMENT,IRESP)
+      TZFIELD%CMNHNAME   = 'VISIZHA'
+      TZFIELD%CSTDNAME   = ''
+      TZFIELD%CLONGNAME  = 'MesoNH: VISIZHA'
+      TZFIELD%CUNITS     = 'm'
+      TZFIELD%CDIR       = 'XY'
+      TZFIELD%CCOMMENT   = 'X_Y_Z_Visibility Zhang'
+      TZFIELD%NGRID      = 1
+      TZFIELD%NTYPE      = TYPEREAL
+      TZFIELD%NDIMS      = 3
+      CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZVISIZHA)
 !
       DEALLOCATE(ZVISIGUL,ZVISIZHA)
     END IF
@@ -2127,11 +2236,16 @@ IF (( LMOIST_E .OR. LBV_FR ) .AND. (NRR>0)) THEN
                  *ZWORK31(:,:,:) *(1. +0.81 *ZWORK31(:,:,:)) )
 !
   IF (LMOIST_E) THEN
-    YRECFM='THETAE'
-    YCOMMENT='X_Y_Z_Equivalent potential Temperature (K)'
-    IGRID=1
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZTHETAE,IGRID,ILENCH,YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'THETAE'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: THETAE'
+    TZFIELD%CUNITS     = 'K'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_Z_Equivalent potential temperature'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 3
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZTHETAE)
   END IF
 END IF
 !-------------------------------------------------------------------------------
@@ -2147,11 +2261,16 @@ IF (LMOIST_ES .AND. (NRR>0)) THEN
        -4.805   )    ) + 55.
   ZTHETAES(:,:,:)= XTHT(:,:,:) * EXP( (3376. / ZTHETAE(:,:,:) - 2.54)  &
                *ZWORK31(:,:,:) *(1. +0.81 *ZWORK31(:,:,:)) )
-  YRECFM='THETAES'
-  YCOMMENT='X_Y_Z_Equivalent Saturated potential Temperature(K)'
-  IGRID=1
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZTHETAES,IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'THETAES'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: THETAES'
+  TZFIELD%CUNITS     = 'K'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = 'X_Y_Z_Equivalent Saturated potential temperature'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZTHETAES)
 ENDIF
 !
 !-------------------------------------------------------------------------------
@@ -2161,19 +2280,29 @@ ENDIF
 IF (LVORT) THEN
 ! Vorticity x
   ZWORK31(:,:,:)=MYF(MZF(1,IKU,1,MXM(ZVOX(:,:,:))))
-  YRECFM='UM1'
-  YCOMMENT='X_Y_Z_x component of vorticity (/S)'
-  IGRID=2
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'UM1'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: UM1'
+  TZFIELD%CUNITS     = 's-1'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = 'X_Y_Z_x component of vorticity'
+  TZFIELD%NGRID      = 2
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
 !    
 ! Vorticity y
   ZWORK32(:,:,:)=MZF(1,IKU,1,MXF(MYM(ZVOY(:,:,:))))
-  YRECFM='VM1'
-  YCOMMENT='X_Y_Z_y component of vorticity (/S)'
-  IGRID=3
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK32,IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'VM1'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: VM1'
+  TZFIELD%CUNITS     = 's-1'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = 'X_Y_Z_y component of vorticity'
+  TZFIELD%NGRID      = 3
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK32)
   !
   IF (LWIND_ZM) THEN
     TZFIELD2(1)%CMNHNAME   = 'UM1_ZM'
@@ -2201,19 +2330,29 @@ IF (LVORT) THEN
 !    
 ! Vorticity z
   ZWORK31(:,:,:)=MXF(MYF(MZM(1,IKU,1,ZVOZ(:,:,:))))
-  YRECFM='WM1'
-  YCOMMENT='X_Y_Z_relative vorticity (/S)'
-  IGRID=4
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'WM1'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: WM1'
+  TZFIELD%CUNITS     = 's-1'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = 'X_Y_Z_relative vorticity'
+  TZFIELD%NGRID      = 4
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
 !
 ! Absolute Vorticity 
   ZWORK31(:,:,:)=MYF(MXF(ZVOZ(:,:,:))) + ZCORIOZ(:,:,:)
-  YRECFM='ABVOR'
-  YCOMMENT='X_Y_Z_z ABsolute VORticity (/S)'
-  IGRID=1
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'ABVOR'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: ABVOR'
+  TZFIELD%CUNITS     = 's-1'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = 'X_Y_Z_z ABsolute VORticity'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
 !
 END IF
 !    
@@ -2235,11 +2374,16 @@ IF ( LMEAN_POVO ) THEN
     END WHERE
   END DO
   WHERE (IWORK1(:,:)>0) ZWORK21(:,:)=ZWORK21(:,:)/REAL( IWORK1(:,:) )
-  YRECFM='MEAN_POVO'
-  YCOMMENT='X_Y_Z_MEAN of POtential VOrticity (PVU)'
-  IGRID=4
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'MEAN_POVO'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: MEAN_POVO'
+  TZFIELD%CUNITS     = 'PVU' ! 1 PVU = 1e-6 m^2 s^-1 K kg^-1
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = 'X_Y_Z_MEAN of POtential VOrticity'
+  TZFIELD%NGRID      = 4
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 2
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK21)
 END IF
 !
 ! Virtual Potential Vorticity in PV units
@@ -2251,11 +2395,16 @@ IF (LMOIST_V .AND. (NRR>0) ) THEN
                + ZWORK32(:,:,:)*MZF(1,IKU,1,MXF(ZVOY(:,:,:)))     &
                + ZWORK33(:,:,:)*(MYF(MXF(ZVOZ(:,:,:))) + ZCORIOZ(:,:,:))
   ZWORK34(:,:,:)=ZWORK34(:,:,:)*1E6/XRHODREF(:,:,:)
-  YRECFM='POVOV'
-  YCOMMENT='X_Y_Z_Virtual POtential VOrticity (PVU)'
-  IGRID=1
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK34,IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'POVOV'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: POVOV'
+  TZFIELD%CUNITS     = 'PVU' ! 1 PVU = 1e-6 m^2 s^-1 K kg^-1
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = 'X_Y_Z_Virtual POtential VOrticity'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK34)
 !
   IF (LMEAN_POVO) THEN
     IWORK1(:,:)=0
@@ -2267,11 +2416,16 @@ IF (LMOIST_V .AND. (NRR>0) ) THEN
       END WHERE
     END DO
     WHERE(IWORK1(:,:)>0) ZWORK21(:,:)=ZWORK21(:,:)/REAL( IWORK1(:,:) )
-    YRECFM='MEAN_POVOV'
-    YCOMMENT='X_Y_Z_MEAN of Virtual POtential VOrticity (PVU)'
-    IGRID=4
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'MEAN_POVOV'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: MEAN_POVOV'
+    TZFIELD%CUNITS     = 'PVU' ! 1 PVU = 1e-6 m^2 s^-1 K kg^-1
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_Z_MEAN of Virtual POtential VOrticity'
+    TZFIELD%NGRID      = 4
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 2
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK21)
   END IF
 END IF
 !
@@ -2285,11 +2439,16 @@ IF (LMOIST_E .AND. (NRR>0) ) THEN
                 + ZWORK32(:,:,:)*MZF(1,IKU,1,MXF(ZVOY(:,:,:)))     &
                 + ZWORK33(:,:,:)*(MYF(MXF(ZVOZ(:,:,:))) + ZCORIOZ(:,:,:))
   ZWORK34(:,:,:)=ZWORK34(:,:,:)*1E6/XRHODREF(:,:,:)
-  YRECFM='POVOE'
-  YCOMMENT='X_Y_Z_Equivalent POtential VOrticity (PVU)'
-  IGRID=1
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK34,IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'POVOE'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: POVOE'
+  TZFIELD%CUNITS     = 'PVU' ! 1 PVU = 1e-6 m^2 s^-1 K kg^-1
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = 'X_Y_Z_Equivalent POtential VOrticity'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK34)
 !
   IF (LMEAN_POVO) THEN
     IWORK1(:,:)=0
@@ -2301,11 +2460,16 @@ IF (LMOIST_E .AND. (NRR>0) ) THEN
       END WHERE
     END DO
     WHERE(IWORK1(:,:)>0) ZWORK21(:,:)=ZWORK21(:,:)/REAL( IWORK1(:,:) )
-    YRECFM='MEAN_POVOE'
-    YCOMMENT='X_Y_Z_MEAN of Equivalent POtential VOrticity (PVU)'
-    IGRID=4
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'MEAN_POVOE'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: MEAN_POVOE'
+    TZFIELD%CUNITS     = 'PVU' ! 1 PVU = 1e-6 m^2 s^-1 K kg^-1
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_Z_MEAN of Equivalent POtential VOrticity'
+    TZFIELD%NGRID      = 4
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 2
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK21)
     DEALLOCATE(IWORK1)
   END IF 
   !
@@ -2320,11 +2484,16 @@ IF (LMOIST_ES .AND. (NRR>0) ) THEN
                 + ZWORK32(:,:,:)*MZF(1,IKU,1,MXF(ZVOY(:,:,:)))     &
                 + ZWORK33(:,:,:)*(MYF(MXF(ZVOZ(:,:,:))) + ZCORIOZ(:,:,:))
   ZWORK34(:,:,:)=ZWORK34(:,:,:)*1E6/XRHODREF(:,:,:)
-  YRECFM='POVOES'
-  YCOMMENT='X_Y_Z_Equivalent Saturated POtential VOrticity (PVU)'
-  IGRID=1
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK34,IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'POVOES'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: POVOES'
+  TZFIELD%CUNITS     = 'PVU' ! 1 PVU = 1e-6 m^2 s^-1 K kg^-1
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = 'X_Y_Z_Equivalent Saturated POtential VOrticity'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK34)
 ENDIF
 !
 !
@@ -2335,23 +2504,32 @@ ENDIF
 IF (LDIV) THEN
 !
   ZWORK31=GX_U_M(1,IKU,1,XUT,XDXX,XDZZ,XDZX) + GY_V_M(1,IKU,1,XVT,XDYY,XDZZ,XDZY)
-  YRECFM='HDIV'
-  YCOMMENT='X_Y_Z_Horizontal DIVergence (/S)'
-  YCOMMENT='(1/s)'
-  IGRID=1
-  ILENCH=LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'HDIV'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: HDIV'
+  TZFIELD%CUNITS     = 's-1'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = 'X_Y_Z_Horizontal DIVergence'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
 !
   IF (LUSERV) THEN
-   YRECFM= 'HMDIV'
-   YCOMMENT='X_Y_Z_Horizontal Moisture DIVergence HMDIV (KG/M3/S)'
-   ILENCH=LEN(YCOMMENT)
-   ZWORK31=MXM(XRHODREF*XRT(:,:,:,1))*XUT
-   ZWORK32=MYM(XRHODREF*XRT(:,:,:,1))*XVT
-   ZWORK33=GX_U_M(1,IKU,1,ZWORK31,XDXX,XDZZ,XDZX) + GY_V_M(1,IKU,1,ZWORK32,XDYY,XDZZ,XDZY)
-   CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK33,IGRID,ILENCH,  &
-        YCOMMENT,IRESP)
- END IF
+    TZFIELD%CMNHNAME   = 'HMDIV'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: HMDIV'
+    TZFIELD%CUNITS     = 'kg m-3 s-1'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_Z_Horizontal Moisture DIVergence HMDIV'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 3
+    ZWORK31=MXM(XRHODREF*XRT(:,:,:,1))*XUT
+    ZWORK32=MYM(XRHODREF*XRT(:,:,:,1))*XVT
+    ZWORK33=GX_U_M(1,IKU,1,ZWORK31,XDXX,XDZZ,XDZX) + GY_V_M(1,IKU,1,ZWORK32,XDYY,XDZZ,XDZY)
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK33)
+  END IF
 !
 ENDIF
 !
@@ -2395,17 +2573,27 @@ IF (LGEO .OR. LAGEO) THEN
   DEALLOCATE(ZPHI)
 !
   IF (LGEO) THEN 
-    YRECFM='UM88'
-    YCOMMENT='X_Y_Z_U component of GEOstrophic wind (m/s)'
-    IGRID=2
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'UM88'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: UM88'
+    TZFIELD%CUNITS     = 'm s-1'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_Z_U component of GEOstrophic wind'
+    TZFIELD%NGRID      = 2
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 3
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
 ! 
-    YRECFM='VM88'
-    YCOMMENT='X_Y_Z_V component of GEOstrophic wind (m/s)'
-    IGRID=3
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK32,IGRID,ILENCH,YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'VM88'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: VM88'
+    TZFIELD%CUNITS     = 'm s-1'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_Z_V component of GEOstrophic wind'
+    TZFIELD%NGRID      = 3
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 3
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK32)
     !
     IF (LWIND_ZM) THEN
       TZFIELD2(1)%CMNHNAME   = 'UM88_ZM'
@@ -2432,12 +2620,11 @@ IF (LGEO .OR. LAGEO) THEN
     ENDIF
 !
 ! wm necessary to plot vertical cross sections of wind vectors
-    YRECFM='WM88'
-    YCOMMENT='X_Y_Z_vertical wind (m/s)'
-    YCOMMENT='(M/S)'
-    IGRID=4
-    ILENCH=LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',XWT,IGRID,ILENCH,YCOMMENT,IRESP)
+    CALL FIND_FIELD_ID_FROM_MNHNAME('WT',IID,IRESP)
+    TZFIELD = TFIELDLIST(IID)
+    TZFIELD%CMNHNAME  = 'WM88'
+    TZFIELD%CLONGNAME = 'MesoNH: WM88'
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,XWT)
   END IF
 !
   IF (LAGEO) THEN
@@ -2676,10 +2863,6 @@ END IF
 !* Accumulated and instantaneous total precip rates  in mm and mm/h
 !
 IF (LTOTAL_PR .AND. SIZE (XACPRR)>0 ) THEN
-  YRECFM      ='ACTOPR'
-  YCOMMENT    ='X_Y_ACccumulated TOtal Precipitation Rate (MM)'
-  IGRID       =1
-  ILENCH      =LEN(YCOMMENT)
   ZWORK21(:,:) = 0.
   !
   IF (LUSERR) THEN
@@ -2702,7 +2885,16 @@ IF (LTOTAL_PR .AND. SIZE (XACPRR)>0 ) THEN
   END IF
   IF (LUSERR .OR. CCLOUD(1:3) == 'ICE' .OR. CCLOUD == 'C3R5' .OR. &
                   CCLOUD == 'LIMA' .OR. CDCONV /= 'NONE') THEN
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'ACTOPR'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: ACTOPR'
+    TZFIELD%CUNITS     = 'mm'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_ACccumulated TOtal Precipitation Rate'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 2
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK21)
   ELSE
     PRINT * ,'YOU WANT TO COMPUTE THE ACCUMULATED RAIN'
     PRINT * ,'BUT NO RAIN IS PRESENT IN THE MODEL' 
@@ -2711,27 +2903,29 @@ IF (LTOTAL_PR .AND. SIZE (XACPRR)>0 ) THEN
   ! calculation of the mean accumulated precipitations in the mesh-grid of a 
   !large-scale model
   IF (LMEAN_PR .AND. LUSERR) THEN
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CUNITS     = 'mm'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_Large Scale ACccumulated TOtal Precipitation Rate'
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 2
+    !
     DO JK=1,SIZE(XMEAN_PR),2      
       IF (XMEAN_PR(JK) .NE. XUNDEF .AND. XMEAN_PR(JK+1) .NE. XUNDEF) THEN
         PRINT * ,'MEAN accumulated RAIN: GRID ', XMEAN_PR(JK), XMEAN_PR(JK+1)
-        CALL COMPUTE_MEAN_PRECIP(ZWORK21,XMEAN_PR(JK:JK+1),ZWORK22,IGRID)
+        CALL COMPUTE_MEAN_PRECIP(ZWORK21,XMEAN_PR(JK:JK+1),ZWORK22,TZFIELD%NGRID)
         !
         JI=INT(XMEAN_PR(JK))
         JJ=INT(XMEAN_PR(JK+1))
-        WRITE(YRECFM,'(A9,2I2.2)')'LS_ACTOPR',JI,JJ
-        YCOMMENT    ='X_Y_Large Scale ACccumulated TOtal Precipitation Rate (MM)'
-        ILENCH      =LEN(YCOMMENT)
-        CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK22,IGRID,ILENCH,YCOMMENT,IRESP)
+        WRITE(TZFIELD%CMNHNAME,'(A9,2I2.2)')'LS_ACTOPR',JI,JJ
+        TZFIELD%CLONGNAME  = 'MesoNH: '//TRIM(TZFIELD%CMNHNAME)
+        CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK22)
       END IF
     END DO
     !
   END IF
   !
   !
-  YRECFM      ='INTOPR'
-  YCOMMENT    ='X_Y_INstantaneous TOtal Precipitation Rate (MM/H)'
-  IGRID       =1
-  ILENCH      =LEN(YCOMMENT)
   ZWORK21(:,:) = 0.
   !    
   IF (LUSERR) THEN
@@ -2754,7 +2948,16 @@ IF (LTOTAL_PR .AND. SIZE (XACPRR)>0 ) THEN
   END IF
   IF (LUSERR .OR. CCLOUD(1:3) == 'ICE' .OR. CCLOUD == 'C3R5' .OR. &
                   CCLOUD == 'LIMA' .OR. CDCONV /= 'NONE') THEN
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK21,IGRID,ILENCH,YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'INTOPR'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: INTOPR'
+    TZFIELD%CUNITS     = 'mm hour-1'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_INstantaneous TOtal Precipitation Rate'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 2
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK21)
   ELSE
     PRINT * ,'YOU WANT TO COMPUTE THE RAIN RATE'
     PRINT * ,'BUT NO RAIN IS PRESENT IN THE MODEL' 
@@ -2763,12 +2966,17 @@ IF (LTOTAL_PR .AND. SIZE (XACPRR)>0 ) THEN
   ! calculation of the mean instantaneous precipitations in the mesh-grid of a 
   ! large-scale model
   IF (LMEAN_PR .AND. LUSERR) THEN
-    CALL COMPUTE_MEAN_PRECIP(ZWORK21,XMEAN_PR,ZWORK22,IGRID)
+    CALL COMPUTE_MEAN_PRECIP(ZWORK21,XMEAN_PR,ZWORK22,TZFIELD%NGRID)
 !
-    YRECFM      ='LS_INTOPR'
-    YCOMMENT    ='X_Y_Large Scale INstantaneous TOtal Precipitation Rate (MM/H)'
-    ILENCH      =LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK22,IGRID,ILENCH,YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'LS_INTOPR'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: LS_INTOPR'
+    TZFIELD%CUNITS     = 'mm hour-1'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_Large Scale INstantaneous TOtal Precipitation Rate'
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 2
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK22)
   END IF
 !
 END IF
@@ -2880,11 +3088,16 @@ IF (LBV_FR) THEN
    ENDDO
   ENDDO
   !
-  YRECFM      ='BV'   
-  YCOMMENT    ='X_Y_Z_Brunt-Vaissala frequency (/S)'
-  IGRID       =4
-  ILENCH      =LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'BV'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: BV'
+  TZFIELD%CUNITS     = 's-1'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = 'X_Y_Z_Brunt-Vaissala frequency'
+  TZFIELD%NGRID      = 4
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
 !  
   IF (NRR > 0) THEN
     ZWORK32(:,:,:)=DZM(1,IKU,1,ZTHETAE(:,:,:))/ MZM(1,IKU,1,ZTHETAE(:,:,:))
@@ -2900,11 +3113,16 @@ IF (LBV_FR) THEN
      ENDDO
     ENDDO
 !
-    YRECFM      ='BVE'   
-    YCOMMENT    ='X_Y_Z_Equivalent Brunt-Vaissala frequency (/S)'
-    IGRID       =4
-    ILENCH      =LEN(YCOMMENT)
-    CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'BVE'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: BVE'
+    TZFIELD%CUNITS     = 's-1'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_Z_Equivalent Brunt-Vaissala frequency'
+    TZFIELD%NGRID      = 4
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 3
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
   END IF
 END IF
 !
@@ -2919,27 +3137,42 @@ IF ( NGPS>=0 ) THEN
   ZGAMREF=-6.5E-3
   ZWORK21(:,:) = ZTEMP(:,:,IKB) - ZGAMREF*((XZZ(:,:,IKB)+XZZ(:,:,IKB+1))/2.-XZS(:,:))
   !
-  YFGRI=ADJUSTL(ADJUSTR(YFMFILE2)//'GPS')
+  YFGRI=ADJUSTL(ADJUSTR(TPFILE%CNAME)//'GPS')
   CALL GPS_ZENITH (YFGRI,XRT(:,:,:,1),ZTEMP,XPABST,ZWORK21,ZWORK22,ZWORK23,ZWORK24)    
   !
-  YRECFM       ='ZTD'
-  YCOMMENT     ='X_Y_Z_Zenithal Total Delay (m)'
-  IGRID        =1
-  ILENCH       =LEN(YCOMMENT)
-  CALL FMWRIT (YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK22,IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'ZTD'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: ZTD'
+  TZFIELD%CUNITS     = 'm'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = 'X_Y_Z_Zenithal Total Delay'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 2
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK22)
   !
   IF (NGPS>=1) THEN
-    YRECFM       ='ZHD'
-    YCOMMENT     ='X_Y_Z_Zenithal Hydrostatic Delay (m)'
-    IGRID        =1
-    ILENCH       =LEN(YCOMMENT)
-    CALL FMWRIT (YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK23,IGRID,ILENCH,YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'ZHD'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: ZHD'
+    TZFIELD%CUNITS     = 'm'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_Z_Zenithal Hydrostatic Delay'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 2
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK23)
     !
-    YRECFM       ='ZWD'
-    YCOMMENT     ='X_Y_Z_Zenithal Wet Delay (m)'
-    IGRID        =1
-    ILENCH       =LEN(YCOMMENT)
-    CALL FMWRIT (YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK24,IGRID,ILENCH,YCOMMENT,IRESP)
+    TZFIELD%CMNHNAME   = 'ZWD'
+    TZFIELD%CSTDNAME   = ''
+    TZFIELD%CLONGNAME  = 'MesoNH: ZWD'
+    TZFIELD%CUNITS     = 'm'
+    TZFIELD%CDIR       = 'XY'
+    TZFIELD%CCOMMENT   = 'X_Y_Z_Zenithal Wet Delay'
+    TZFIELD%NGRID      = 1
+    TZFIELD%NTYPE      = TYPEREAL
+    TZFIELD%NDIMS      = 2
+    CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK24)
     !
   END IF
   !
@@ -3141,7 +3374,7 @@ IF(LRADAR .AND. LUSERR) THEN
                  INT(ANINT(10.*XELEV(JI,JEL))-10*INT(XELEV(JI,JEL)))
           WRITE(YGRID_SIZE,'(I3.3)') 2*NMAX
           DO JJ=1,SIZE(ZWORK42(:,:,:,:,:),5)
-            YRS=YRAD(JJ)//CNAME_RAD(JI)(1:3)//YELEV//YGRID_SIZE//YFMFILE2
+            YRS=YRAD(JJ)//CNAME_RAD(JI)(1:3)//YELEV//YGRID_SIZE//TRIM(TPFILE%CNAME)
             CALL OPEN_ll(UNIT=ILURS,FILE=YRS,IOSTAT=IRESP,STATUS="NEW",ACTION='WRITE', &
                          FORM="FORMATTED",RECL=8192)
             WRITE(ILURS,'(A,4F12.6,2I5)') '**domaine LATLON ',ZWORK43(JI,1,1),ZWORK43(JI,4*NMAX-1,2*NMAX), &
@@ -3169,7 +3402,7 @@ IF(LRADAR .AND. LUSERR) THEN
           WRITE(YELEV,'(I2.2,A1,I1.1)') FLOOR(XELEV(JI,JEL)),'.',&
                 INT(ANINT(10.*XELEV(JI,JEL))-10*INT(XELEV(JI,JEL)))
           DO JJ=1,SIZE(ZWORK42(:,:,:,:,:),5)
-            YRS="P"//YRAD(JJ)//CNAME_RAD(JI)(1:3)//YELEV//YFMFILE2
+            YRS="P"//YRAD(JJ)//CNAME_RAD(JI)(1:3)//YELEV//TRIM(TPFILE%CNAME)
             CALL OPEN_ll(UNIT=ILURS,FILE=YRS,IOSTAT=IRESP,ACTION='WRITE',MODE=GLOBAL)
             DO JH=1,NBAZIM
               DO JV=1,NBSTEPMAX+1
@@ -3186,7 +3419,7 @@ IF(LRADAR .AND. LUSERR) THEN
 END IF
 !
 IF (LLIDAR) THEN
-  PRINT *,'CALL LIDAR/RADAR with YFMFILE2 =',YFMFILE2
+  PRINT *,'CALL LIDAR/RADAR with TPFILE%CNAME =',TPFILE%CNAME
   YVIEW='     '
   YVIEW=TRIM(CVIEW_LIDAR)
   PRINT *,'CVIEW_LIDAR REQUESTED ',YVIEW
@@ -3278,17 +3511,27 @@ IF (LLIDAR) THEN
   IF( ALLOCATED(ZTMP3) ) DEALLOCATE(ZTMP3)
   IF( ALLOCATED(ZTMP4) ) DEALLOCATE(ZTMP4)
 !
-  YRECFM       ='LIDAR'
-  YCOMMENT     ='X_Y_Z_Normalized_Lidar_Profile (1/m/sr)'
-  IGRID        =1
-  ILENCH       =LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK31,IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'LIDAR'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: LIDAR'
+  TZFIELD%CUNITS     = 'm-1 sr-1'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = 'X_Y_Z_Normalized_Lidar_Profile'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK31)
 !
-  YRECFM       ='LIPAR'
-  YCOMMENT     ='X_Y_Z_Particle_Lidar_Profile (1/m/sr)'
-  IGRID        =1
-  ILENCH       =LEN(YCOMMENT)
-  CALL FMWRIT(YFMFILE2,YRECFM,CLUOUT,'XY',ZWORK32,IGRID,ILENCH,YCOMMENT,IRESP)
+  TZFIELD%CMNHNAME   = 'LIPAR'
+  TZFIELD%CSTDNAME   = ''
+  TZFIELD%CLONGNAME  = 'MesoNH: LIPAR'
+  TZFIELD%CUNITS     = 'm-1 sr-1'
+  TZFIELD%CDIR       = 'XY'
+  TZFIELD%CCOMMENT   = 'X_Y_Z_Particle_Lidar_Profile'
+  TZFIELD%NGRID      = 1
+  TZFIELD%NTYPE      = TYPEREAL
+  TZFIELD%NDIMS      = 3
+  CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK32)
 !
 END IF
 !
@@ -3301,7 +3544,7 @@ IF (LBLTOP) THEN
   ZWORK31(:,:,1:IKU-1)=0.5*(XZZ(:,:,1:IKU-1)+XZZ(:,:,2:IKU))
   ZWORK31(:,:,IKU)=2.*ZWORK31(:,:,IKU-1)-ZWORK31(:,:,IKU-2)
   YFMFILE=CINIFILE
-  CINIFILE=YFMFILE2
+  CINIFILE=TPFILE%CNAME
   CALL FREE_ATM_PROFILE(TPFILE,ZTHETAV,ZWORK31,XZS,XZSMT,ZGAMREF,ZWORK32,ZWORK33)
   CINIFILE=YFMFILE
 END IF
diff --git a/src/MNH/write_lfin.f90 b/src/MNH/write_lfin.f90
index 5e66d2934..b462e2647 100644
--- a/src/MNH/write_lfin.f90
+++ b/src/MNH/write_lfin.f90
@@ -717,7 +717,6 @@ IF (NSV >=1) THEN
 !
   IF (NSV_LIMA_END>=NSV_LIMA_BEG) THEN
     TZFIELD%CSTDNAME   = ''
-    TZFIELD%CUNITS     = 'kg-1'
     TZFIELD%CDIR       = 'XY'
     TZFIELD%NGRID      = 1
     TZFIELD%NTYPE      = TYPEREAL
@@ -1683,30 +1682,18 @@ IF (CPROGRAM /= 'IDEAL') THEN
     IF (SIZE(XINPRG) /= 0 ) ZWORK2D = ZWORK2D + XINPRG
     IF (SIZE(XINPRH) /= 0 ) ZWORK2D = ZWORK2D + XINPRH
     IF (SIZE(XINPRC) /= 0 ) ZWORK2D = ZWORK2D + XINPRC
-    TZFIELD%CMNHNAME   = 'INPRT'
-    TZFIELD%CSTDNAME   = ''
-    TZFIELD%CLONGNAME  = 'MesoNH: INPRT'
-    TZFIELD%CUNITS     = 'mm hour-1'
-    TZFIELD%CDIR       = 'XY'
-    TZFIELD%CCOMMENT   = 'X_Y_Total INstantaneaous PRecipitation rate'
-    TZFIELD%NGRID      = 1
-    TZFIELD%NTYPE      = TYPEREAL
-    TZFIELD%NDIMS      = 2
+    CALL FIND_FIELD_ID_FROM_MNHNAME('INPRT',IID,IRESP)
+    TZFIELD = TFIELDLIST(IID)
+    TZFIELD%CUNITS = 'mm hour-1'
     CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK2D*3.6E6)
 !
     ZWORK2D = XACPRR + XACPRS
     IF (SIZE(XINPRG) /= 0 ) ZWORK2D = ZWORK2D + XACPRG
     IF (SIZE(XINPRH) /= 0 ) ZWORK2D = ZWORK2D + XACPRH
     IF (SIZE(XINPRC) /= 0 ) ZWORK2D = ZWORK2D + XACPRC
-    TZFIELD%CMNHNAME   = 'ACPRT'
-    TZFIELD%CSTDNAME   = ''
-    TZFIELD%CLONGNAME  = 'MesoNH: ACPRT'
-    TZFIELD%CUNITS     = 'mm'
-    TZFIELD%CDIR       = 'XY'
-    TZFIELD%CCOMMENT   = 'X_Y_Total ACcumulated PRecipitation rate'
-    TZFIELD%NGRID      = 1
-    TZFIELD%NTYPE      = TYPEREAL
-    TZFIELD%NDIMS      = 2
+    CALL FIND_FIELD_ID_FROM_MNHNAME('ACPRT',IID,IRESP)
+    TZFIELD = TFIELDLIST(IID)
+    TZFIELD%CUNITS = 'mm'
     CALL IO_WRITE_FIELD(TPFILE,TZFIELD,CLUOUT,IRESP,ZWORK2D*1.0E3)
   END IF
   END IF
-- 
GitLab