
   !!======================================================================

  SUBROUTINE out_nud_createfile(nb)

   !&E---------------------------------------------------------------------
   !&E                 ***  ROUTINE out_nud_createfile  ***
   !&E
   !&E ** Purpose :
   !&E
   !&E ** Description :
   !&E
   !&E ** Called by :
   !&E
   !&E ** External calls :
   !&E
   !&E ** Used ij-arrays :
   !&E
   !&E ** Modified variables :
   !&E
   !&E ** Reference :
   !&E
   !&E ** History :
   !&E       !  2013-09  (V. Garnier) Original code
   !&E
   !&E---------------------------------------------------------------------
   !! * Modules used

#include "toolcpp.h"

   USE parameters,ONLY : limin,limax,ljmin,ljmax,kmax,lchain,riosh,   &
                         imin,imax,jmin,jmax,rg_valmanq_io,rsh
   USE comvars2d, ONLY : riog_valid_min_h0,riog_valid_max_h0,         &
                         riog_valid_min_u,riog_valid_max_u,           &
                         riog_valid_min_temp,riog_valid_max_temp,     &
                         riog_valid_min_sal,riog_valid_max_sal,       &
                         riog_valid_min_xe,riog_valid_max_xe,         &
                         riog_valid_min_bz,riog_valid_max_bz
   USE comvars2d, ONLY : name_out_h0,name_out_hx,name_out_hy,name_out_u, &
                         name_out_v,name_out_uz,name_out_vz,name_out_temp, &
                         name_out_sal,name_out_xe,name_out_bz
   USE comvars2d, ONLY : lon2d_g,lat2d_g,h0,hx,hy,h0fond,lon2d,lat2d
   USE comvars2d, ONLY : out_dir,suffix
   USE comvars3d, ONLY : rhoref
   USE ionc4,     ONLY : ionc4_createfile,ionc4_gatt,ionc4_createvar,&
                         ionc4_write, ionc4_write_sig,&
                         ionc4_write_xy, ionc4_write_xyt, ionc4_vatt_chain
   USE comvars3d, ONLY : sig
#if defined key_siggen || defined key_gencoord
   USE comsiggen, ONLY : hc_sig,Csu_sig,Csw_sig,b_sig,theta_sig
#endif
   USE outtools,  ONLY : out_grid
   USE outnetcdf, ONLY : out_netcdf_conv
   USE output,    ONLY : iextmin,iextmax,jextmin,jextmax,stepi,stepj, &
                         igextmin,igextmax,jgextmin,jgextmax,        &
                         kextmin,kextmax,stepk
   USE_MPI toolmpi, ONLY : IPROC,NPROCS,ADD_ALL_MPI_INT,MPI_COMM_MARS,ADD_TAG
#ifdef key_nudging
   USE nudging,   ONLY : nud_init,l_nud_ssh,l_nud_sal,l_nud_temp,l_nud_boz,l_nud_u3d,l_nud_ubt
#endif

#ifdef key_MPI_2D
   INCLUDE 'mpif.h'
#endif

   !! * Arguments
   INTEGER, INTENT(in) :: nb

#ifdef key_nudging
   !! * Local declarations
   LOGICAL               :: l_out_pack,l_out_nc4par
   CHARACTER(LEN=lchain) :: fileres
   INTEGER               :: communicator
   INTEGER               :: grp_world,grp_zoom,ncpu,ip,ierr_mpi
   INTEGER, DIMENSION(:), ALLOCATABLE  :: numrank,listrank
   INTEGER, DIMENSION(1) :: listrank_nb
   REAL(KIND=rsh),DIMENSION(kmax,limin:limax,ljmin:ljmax)   :: tmp3d
   REAL(KIND=riosh),DIMENSION(limin:limax,ljmin:ljmax)      :: tecrit
   REAL(KIND=riosh),DIMENSION(0:kmax+1) :: sigecrit
   REAL(KIND=riosh)      :: vecrit

   !!----------------------------------------------------------------------
   !! * Executable part

   PRINT_DBG*, 'enter OUT_NUD_CREATEFILE routine'

   ! Create the output file saving the nudging terms
   ! -----------------------------------------------
   IF_AGRIF(agrif_root()) THEN
     fileres=TRIM(out_dir)//'nudged_fields_'//TRIM(suffix)//'.nc'
   ELSE_AGRIF
     INST_AGRIF fileres=TRIM(out_dir)//'nudged_fields_zoom'//TRIM(Agrif_Cfixed())//'_'//TRIM(suffix)//'.nc'
   ENDIF_AGRIF

   communicator=0
#ifdef key_MPI_2D
   communicator=MPI_COMM_MARS
   ! what a code ! Not a clever code probably but not clear in my head at the moment
   ! Count the number of processors active over the area to be saved
   ncpu=0
   IF ( jextmin(nb) <= jextmax(nb) .AND. iextmin(nb) <= iextmax(nb) ) ncpu=ncpu+1
   CALL ADD_ALL_MPI_INT(ncpu)
   ! Save the number of active processors
   ALLOCATE (numrank(ncpu))
   ALLOCATE (listrank(nprocs))
   numrank(:)=0   ! must be initialized if the proc 0 is active on the saved zoom
   listrank(:)=0
   ! which cpu is active ?
   IF ( jextmin(nb) <= jextmax(nb) .AND. iextmin(nb) <= iextmax(nb) ) listrank(1)=iproc
   ! all cpus get the number of active cpus plus sleeping cpus (set to 0)
   listrank_nb(1) = listrank(1)
   CALL MPI_ALLGATHER(listrank_nb(1),1,MPI_INTEGER,listrank,1,MPI_INTEGER,MPI_COMM_MARS,ierr_mpi)
   ! all cpus keep the active cpus only
   IF (ncpu/=nprocs) THEN
     ip=0
     DO i=1,nprocs
       IF_MPI ( listrank(i) /= 0 ) THEN
       ip=ip+1
       numrank(ip)=listrank(i)
       ENDIF
     END DO
     ! Group relative to the MPI_COMM_MARS
     CALL MPI_COMM_GROUP(MPI_COMM_MARS,grp_world,ierr_mpi)
     ! Create the group of active cpus
     CALL MPI_GROUP_INCL(grp_world,ncpu,numrank,grp_zoom,ierr_mpi)
     ! Create the group of active cpus
     CALL MPI_COMM_CREATE(MPI_COMM_MARS,grp_zoom,communicator,ierr_mpi)
   ENDIF
#endif

   l_out_nc4par=.FALSE.
#ifdef key_MPI_2D
   l_out_nc4par=.TRUE.    ! also works with l_out_nc4par=.FALSE. (must be changed line 448 too)
#endif
   l_out_pack=.FALSE.
   sigecrit(:)=sig(:)

   ! to work with active cpus over a small area
   IF ( jextmin(nb) <= jextmax(nb) .AND. iextmin(nb) <= iextmax(nb) ) THEN

     IF (l_out_nc4par) THEN
       ! To allow the different cpus to write at the same time into a global
       ! file, set l_out_nc4par= .TRUE.
       ! The dimension of the file must be global (igextmin)
       CALL ionc4_createfile(fileres,lon2d_g(igextmin(nb):igextmax(nb):stepi(nb),jgextmin(nb):jgextmax(nb):stepj(nb)), &
                                     lat2d_g(igextmin(nb):igextmax(nb):stepi(nb),jgextmin(nb):jgextmax(nb):stepj(nb)), &
                                     igextmin(nb),igextmax(nb),stepi(nb), &
                                     jgextmin(nb),jgextmax(nb),stepj(nb), &
                                     kextmin(nb),kextmax(nb),stepk(nb),   &
                                     kmax,0,l_out_nc4par= .TRUE.,comm_active=communicator, &
                                     zaxis_var=sigecrit(1:kmax))
     ELSE
       ! If none l_out_nc4par argument, l_out_nc4par is considered as false
       ! Each cpu create its own output file inside the local domain (iextmin)
       IF_MPI (.NOT. l_out_nc4par) CALL ADD_TAG(fileres,LEN(fileres))
       CALL ionc4_createfile(fileres,lon2d(iextmin(nb):iextmax(nb):stepi(nb),jextmin(nb):jextmax(nb):stepj(nb)), &
                       lat2d(iextmin(nb):iextmax(nb):stepi(nb),jextmin(nb):jextmax(nb):stepj(nb)), &
                       iextmin(nb),iextmax(nb),stepi(nb),jextmin(nb),jextmax(nb),stepj(nb), &
                       kextmin(nb),kextmax(nb),stepk(nb),kmax,0,zaxis_var=sigecrit(1:kmax))
     ENDIF

   ! save attribute of the spatial extension of the domain
   ! -----------------------------------------------------
     CALL ionc4_gatt(fileres,'global_imin',igextmin(nb))
     CALL ionc4_gatt(fileres,'global_imax',igextmax(nb))
     CALL ionc4_gatt(fileres,'global_jmin',jgextmin(nb))
     CALL ionc4_gatt(fileres,'global_jmax',jgextmax(nb))

   ! save longitude, latitude fields and variables relative to the mesh grid
   ! -----------------------------------------------------------------------
     CALL out_grid(fileres)

     CALL ionc4_createvar(fileres,'SIG',' ','sigma variable',&
#if defined key_siggen || defined key_gencoord
                        standard_name= 'ocean_s_variable',&
#else
                        standard_name= 'ocean_sigma_variable',&
#endif                          
                        l_auxcoordinate=.false.,&
                        dims='z',valid_min= -1.0_riosh,valid_max= 0.0_riosh,fill_value= rg_valmanq_io)
     CALL ionc4_write_sig(fileres,'SIG',sigecrit,kmax)

#if defined key_siggen || defined key_gencoord

     sigecrit(:)=Csu_sig(:)
     CALL ionc4_createvar(fileres,'Csu_sig',' ','C(s) law at mid layer',dims='z', &
                        standard_name='ocean_s_coordinate_function_at_midlayer', &
                        valid_min=-2.0_riosh,valid_max=2.0_riosh,fill_value= rg_valmanq_io, &
                        l_auxcoordinate=.false.)
     CALL ionc4_write_sig(fileres,'Csu_sig',sigecrit,kmax)
     CALL ionc4_vatt_chain(fileres,'Csu_sig','formula_definition', &
          'C(k) = (1-b)*sinh(a*s(k))/sinh(a) + b*[tanh(a*(s(k)+0.5))-tanh(0.5*a)]/(2*tanh(0.5*a))')


     CALL ionc4_createvar(fileres,'hc','m', 'depth over which the resolution is increased', &
                        standard_name='ocean_s_coordinate_thick_of_surface_resolution',dims='xy', &
                        valid_min= -20.0_riosh,valid_max= 500.0_riosh,fill_value=-rg_valmanq_io,&
                        l_out_nc4par=l_out_nc4par,l_pack=l_out_pack )
     tecrit(limin:limax,ljmin:ljmax)=hc_sig(limin:limax,ljmin:ljmax)
     CALL ionc4_write_xy(fileres,'hc',tecrit,limin,limax,ljmin,ljmax,fill_value=-rg_valmanq_io,l_out_nc4par=l_out_nc4par)

     CALL ionc4_createvar(fileres,'b',' ','bottom control parameter',  &
                        standard_name='ocean_s_coordinate_bottom_control', &
                        valid_min=0.0_riosh,valid_max=1.0_riosh,fill_value=rg_valmanq_io)
     vecrit=b_sig
     CALL ionc4_write(fileres,'b',vecrit)

     CALL ionc4_createvar(fileres,'theta',' ','surface control parameter', &
                        standard_name ='ocean_s_coordinate_surface_control', &
                        valid_min=0.0_riosh,valid_max=20.0_riosh,fill_value=rg_valmanq_io)
     vecrit=theta_sig
     CALL ionc4_write(fileres,'theta',vecrit)

#endif  

   ! Add attributes to follow OCO netCDF convention
     CALL out_netcdf_conv(fileres,nb,communicator,igextmin(nb),igextmax(nb),jgextmin(nb),jgextmax(nb))

     IF_MPI (l_out_nc4par .AND. communicator/=MPI_COMM_MARS) CALL MPI_COMM_FREE(communicator,ierr_mpi)


   ! Create all the variables in the file
   ! ------------------------------------

   ! Save the bathymetry
   ! -------------------

   ! h0
   ! --
     CALL ionc4_createvar(fileres,name_out_h0,units='m',&
                        long_name='bathymetry relative to the mean level',&
                        standard_name = 'model_sea_floor_depth_below_geoid',&
                        fill_value = -rg_valmanq_io,&
                        valid_min = riog_valid_min_h0,&
                        valid_max = riog_valid_max_h0,&
                        dims = 'xy',&
                        l_pack = l_out_pack,&
                        l_out_nc4par = l_out_nc4par)
     tecrit(:,:)=h0(limin:limax,ljmin:ljmax)
     CALL ionc4_write_xy(fileres,TRIM(name_out_h0),tecrit,limin,limax,ljmin,ljmax,&
                       fill_value=-rg_valmanq_io,l_out_nc4par = l_out_nc4par)


     IF (l_nud_ubt .OR. l_nud_u3d) THEN
       CALL ionc4_createvar(fileres,name_out_hx,units='m',location_in_cell='u',&
                           long_name = 'bathymetry at u-location relative to the mean level',&
                           standard_name = 'model_sea_floor_depth_below_geoid_at_u_location',&
                           fill_value = -REAL(rg_valmanq_io+h0fond,riosh),&
                           valid_min = riog_valid_min_h0,&
                           valid_max = riog_valid_max_h0,&
                           dims = 'xy',&
                           l_pack = l_out_pack,&
                           l_out_nc4par = l_out_nc4par)
       tecrit(:,:)=hx(limin:limax,ljmin:ljmax)
       CALL ionc4_write_xy(fileres,TRIM(name_out_hx),tecrit,limin,limax,ljmin,ljmax ,&
                         fill_value = -REAL(rg_valmanq_io+h0fond,riosh),l_out_nc4par = l_out_nc4par)
       CALL ionc4_createvar(fileres,name_out_hy,units='m',location_in_cell='v',&
                           long_name = 'bathymetry at v-location relative to the mean level',&
                           standard_name = 'model_sea_floor_depth_below_geoid_at_v_location',&
                           fill_value = -REAL(rg_valmanq_io+h0fond,riosh),&
                           valid_min = riog_valid_min_h0,&
                           valid_max = riog_valid_max_h0,&
                           dims = 'xy',&
                           l_pack = l_out_pack,&
                           l_out_nc4par = l_out_nc4par)
       tecrit(:,:)=hy(limin:limax,ljmin:ljmax)
       CALL ionc4_write_xy(fileres,TRIM(name_out_hy),tecrit,limin,limax,ljmin,ljmax,&
                         fill_value = -REAL(rg_valmanq_io+h0fond,riosh),l_out_nc4par = l_out_nc4par)
     END IF

   ! Create time-depending variables
   !--------------------------------

     IF (l_nud_ssh) THEN
       CALL ionc4_createvar(fileres,name_out_xe,units= 'm',&
                          long_name = 'sea surface height',&
                          standard_name = 'sea_surface_height_above_geoid',&
                          fill_value = rg_valmanq_io,&
                          valid_min = riog_valid_min_xe,&
                          valid_max = riog_valid_max_xe,&
                          dims= 'xyt',&
                          l_pack = l_out_pack,&
                          l_out_nc4par = l_out_nc4par)

!       CALL ionc4_vatt_axis(fileres,TRIM(name_out_xe),'Z')
!       CALL ionc4_vatt_positive(fileres,TRIM(name_out_xe),'up')
     END IF

     IF (l_nud_ubt) THEN
       CALL ionc4_createvar(fileres,name_out_u,units= 'm s-1',&
                            long_name = 'barotropic zonal velocity',&
                            standard_name = 'barotropic_sea_water_x_velocity_at_u_location',&
                            location_in_cell='u',&
                            fill_value = rg_valmanq_io,&
                            valid_min = riog_valid_min_u,&
                            valid_max = riog_valid_max_u,&
                            dims= 'xyt',&
                            l_pack = l_out_pack,&
                            l_out_nc4par = l_out_nc4par)
       CALL ionc4_createvar(fileres,name_out_v,units= 'm s-1',&
                            long_name = 'barotropic meridional velocity',&
                            standard_name = 'barotropic_sea_water_y_velocity_at_v_location',&
                            location_in_cell='v',&
                            fill_value = rg_valmanq_io,&
                            valid_min = riog_valid_min_u,&
                            valid_max = riog_valid_max_u,&
                            dims= 'xyt',&
                            l_pack = l_out_pack,&
                            l_out_nc4par = l_out_nc4par)                
     END IF

     IF (l_nud_u3d) THEN
       CALL ionc4_createvar(fileres,name_out_uz,units= 'm s-1',&
                            long_name = '3d zonal velocity',&
                            standard_name = 'sea_water_x_velocity_at_u_location',&
                            location_in_cell='u',&
                            fill_value = rg_valmanq_io,&
                            valid_min = riog_valid_min_u,&
                            valid_max = riog_valid_max_u,&
                            dims= 'zxyt',&
                            l_pack = l_out_pack,&
                            l_out_nc4par = l_out_nc4par)
       CALL ionc4_createvar(fileres,name_out_vz,units= 'm s-1',&
                            long_name = '3d meridional velocity',&
                            standard_name = 'sea_water_y_velocity_at_v_location',&
                            location_in_cell='v',&
                            fill_value = rg_valmanq_io,&
                            valid_min = riog_valid_min_u,&
                            valid_max = riog_valid_max_u,&
                            dims= 'zxyt',&
                            l_pack = l_out_pack,&
                            l_out_nc4par = l_out_nc4par)
     END IF

     IF (l_nud_sal) THEN
       CALL ionc4_createvar(fileres,name_out_sal,units= '1e-3',&
                            long_name = 'salinity',&
                            standard_name = 'sea_water_salinity',&
                            fill_value = rg_valmanq_io,&
                            valid_min = riog_valid_min_sal,&
                            valid_max = riog_valid_max_sal,&
                            dims= 'zxyt',&
                            l_pack = l_out_pack,&
                            l_out_nc4par = l_out_nc4par)
     END IF

     IF (l_nud_temp) THEN
       CALL ionc4_createvar(fileres,name_out_temp,units= 'degrees_Celsius',&
                            long_name = 'temperature',&
                            standard_name = 'sea_water_potential_temperature',&
                            fill_value = rg_valmanq_io,&
                            valid_min = riog_valid_min_temp,&
                            valid_max = riog_valid_max_temp,&
                            dims= 'zxyt',&
                            l_pack = l_out_pack,&
                            l_out_nc4par = l_out_nc4par)
     END IF

     IF (l_nud_boz) THEN
       CALL ionc4_createvar(fileres,name_out_bz,units= 'kg m-3',&
                            long_name = 'density of water at atmospheric pressure minus 1000 kg m-3',&
                            standard_name = 'sea_water_sigma_t',&
                            fill_value = rg_valmanq_io,&
                            valid_min = riog_valid_min_bz,&
                            valid_max = riog_valid_max_bz,&
                            dims= 'zxyt',&
                            l_pack = l_out_pack,&
                            l_out_nc4par = l_out_nc4par)
     END IF

   END IF
#endif

   PRINT_DBG*, 'exit OUT_NUD_CREATEFILE routine'

  END SUBROUTINE out_nud_createfile

   !!======================================================================

  SUBROUTINE out_nud_save(nb)

   !&E---------------------------------------------------------------------
   !&E                 ***  ROUTINE out_nud_save  ***
   !&E
   !&E ** Purpose :
   !&E
   !&E ** Description :
   !&E
   !&E ** Called by :
   !&E
   !&E ** External calls :
   !&E
   !&E ** Used ij-arrays :
   !&E
   !&E ** Modified variables :
   !&E
   !&E ** Reference :
   !&E
   !&E ** History :
   !&E       !  2013-09  (V. Garnier) Original code
   !&E
   !&E---------------------------------------------------------------------
   !! * Modules used
   USE parameters, ONLY : limin,limax,ljmin,ljmax,kmax,valmanq,lchain,rg_valmanq_io,grav, &
                          rsh,rlg
   USE comvars2d,  ONLY : h0fond,h0,hx,hy,t
   USE comvars2d,  ONLY : name_out_u,name_out_v,name_out_uz,name_out_vz,name_out_temp, &
                          name_out_sal,name_out_xe,name_out_bz,out_dir,suffix
   USE comvarp2d,  ONLY : ssh
   USE comvars3d,  ONLY : rhoref
   USE ionc4,      ONLY : ionc4_write_time,ionc4_sync
#ifdef key_nudging
   USE nudging,    ONLY : nud_init,l_nud_ssh,l_nud_sal,l_nud_temp,l_nud_boz,l_nud_u3d,l_nud_ubt,l_nud_zonalmean
   USE nudging,    ONLY : ssh_clim_cum, u_clim_cum, v_clim_cum, &
                          temp_clim_cum, sal_clim_cum, boz_clim_cum, &
                          uz_clim_cum, vz_clim_cum, &
                          ssh_zonal_cum, u_zonal_cum, v_zonal_cum, &
                          uz_zonal_cum, vz_zonal_cum, &
                          temp_zonal_cum, sal_zonal_cum, boz_zonal_cum
   USE output,    ONLY : iextmin,iextmax,jextmin,jextmax
   USE outtools,  ONLY : out_write_zxytvar,out_write_xytvar,out_write_zytvar,out_write_ytvar
   USE_MPI toolmpi,ONLY : ADD_TAG
#endif

   !! * Arguments
   INTEGER, INTENT(in) :: nb

#ifdef key_nudging
   !! * Local declarations
   LOGICAL         :: l_out_nc4par
   LOGICAL, DIMENSION(limin:limax,ljmin:ljmax) :: mask_h0,mask_hx,mask_hy
   CHARACTER(LEN=lchain) :: fileres
   INTEGER         :: i,j
   REAL(KIND=rsh),DIMENSION(kmax,ljmin:ljmax)   :: tmp2d
   REAL(KIND=rsh),DIMENSION(kmax,limin:limax,ljmin:ljmax)   :: tmp3d

   !!----------------------------------------------------------------------
   !! * Executable part

   PRINT_DBG*, 'enter OUT_NUD_SAVE routine'

   ! Save the fields
   ! ---------------

   IF ( jextmin(nb) <= jextmax(nb) .AND. iextmin(nb) <= iextmax(nb) ) THEN

     l_out_nc4par=.FALSE.
#ifdef key_MPI_2D
   l_out_nc4par=.TRUE.    ! also works with l_out_nc4par=.FALSE. (must be changed line 132 too)
#endif

   ! create the mask at first

     mask_h0(:,:)=.TRUE.
     DO j=ljmin,ljmax
       DO i=limin,limax
       ! For heights smaller than h0fond, the cell is considered as wet
         IF ( h0(i,j)+ssh(i,j) > h0fond ) mask_h0(i,j)=.FALSE.
       END DO
     END DO

     mask_hx(:,:)=.TRUE.
     DO j=ljmin,ljmax
       DO i=limin,limax
         IF ( hx(i,j) /= -valmanq -h0fond ) mask_hx(i,j)=.FALSE.
       END DO
     END DO

     mask_hy(:,:)=.TRUE.
     DO j=ljmin,ljmax
       DO i=limin,limax
         IF (hy(i,j) /= -valmanq -h0fond) mask_hy(i,j)=.FALSE.
       END DO
     END DO

   ! save the fields

     fileres=TRIM(out_dir)//'nudged_fields_'//TRIM(suffix)//'.nc'
     IF_MPI (.NOT. l_out_nc4par) CALL ADD_TAG(fileres,LEN(fileres))
     CALL ionc4_write_time(fileres,0,t)

     IF (l_nud_zonalmean) THEN  ! Zonal fields

       IF (l_nud_ssh) CALL out_write_ytvar(fileres,name_out_xe,rg_valmanq_io,mask_h0(limin,:), &
                                            ssh_zonal_cum(ljmin:ljmax),0,l_out_nc4par = l_out_nc4par)

       IF (l_nud_ubt) THEN
         CALL out_write_ytvar(fileres,name_out_u,rg_valmanq_io,mask_hx(limin,:), &
                               u_zonal_cum(ljmin:ljmax),0,l_out_nc4par = l_out_nc4par)
         CALL out_write_ytvar(fileres,name_out_v,rg_valmanq_io,mask_hy(limin,:), &
                               v_zonal_cum(ljmin:ljmax),0,l_out_nc4par = l_out_nc4par)
       END IF

       IF (l_nud_u3d) THEN
         CALL out_write_zytvar(fileres,name_out_uz,rg_valmanq_io,mask_hx(limin,:), &
                                uz_zonal_cum(:,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)
         CALL out_write_zytvar(fileres,name_out_vz,rg_valmanq_io,mask_hy(limin,:), &
                                vz_zonal_cum(:,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)
       END IF

       IF (l_nud_sal) CALL out_write_zytvar(fileres,name_out_sal,rg_valmanq_io,mask_h0(limin,:), &
                              sal_zonal_cum(:,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)

       IF (l_nud_temp) CALL out_write_zytvar(fileres,name_out_temp,rg_valmanq_io,mask_h0(limin,:), &
                             temp_zonal_cum(:,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)

       IF (l_nud_boz) THEN
         tmp2d(:,:) = rhoref - rhoref*boz_zonal_cum(:,ljmin:ljmax)/grav
         CALL out_write_zytvar(fileres,name_out_bz,rg_valmanq_io,mask_h0(limin,:), &
                                tmp3d,0,1,l_out_nc4par = l_out_nc4par)
       END IF

     ELSE  ! 3d fields 

       IF (l_nud_ssh) CALL out_write_xytvar(fileres,name_out_xe,rg_valmanq_io,mask_h0, &
                                            ssh_clim_cum(limin:limax,ljmin:ljmax),0,l_out_nc4par = l_out_nc4par)

       IF (l_nud_ubt) THEN
         CALL out_write_xytvar(fileres,name_out_u,rg_valmanq_io,mask_hx, &
                               u_clim_cum(limin:limax,ljmin:ljmax),0,l_out_nc4par = l_out_nc4par)
         CALL out_write_xytvar(fileres,name_out_v,rg_valmanq_io,mask_hy, &
                               v_clim_cum(limin:limax,ljmin:ljmax),0,l_out_nc4par = l_out_nc4par)
       END IF
       IF (l_nud_u3d) THEN
         CALL out_write_zxytvar(fileres,name_out_uz,rg_valmanq_io,mask_hx, &
                                uz_clim_cum(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)
         CALL out_write_zxytvar(fileres,name_out_vz,rg_valmanq_io,mask_hy, &
                                vz_clim_cum(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)
       END IF

       IF (l_nud_sal) CALL out_write_zxytvar(fileres,name_out_sal,rg_valmanq_io,mask_h0, &
                              sal_clim_cum(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)

       IF (l_nud_temp) CALL out_write_zxytvar(fileres,name_out_temp,rg_valmanq_io,mask_h0, &
                             temp_clim_cum(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)

       IF (l_nud_boz) THEN
         tmp3d(:,:,:) = rhoref - rhoref*boz_clim_cum(:,limin:limax,ljmin:ljmax)/grav
         CALL out_write_zxytvar(fileres,name_out_bz,rg_valmanq_io,mask_h0, &
                                tmp3d,0,1,l_out_nc4par = l_out_nc4par)
       END IF

     END IF  ! ends test on 3D or zonal variables

     ! To write the data on the disk and not loose data in case of run crash
     CALL ionc4_sync(fileres)

   END IF
#endif

   PRINT_DBG*, 'exit OUT_NUD_SAVE routine'

  END SUBROUTINE out_nud_save

