  MODULE diagdyn

!!======================================================================
!!                   ***  MODULE diagdyn  ***
!!
!! ** Description :
!!
!! ** Reference :
!!
!! ** History :  !  2011-04  (P. Lazure)  Original code
!!
!!======================================================================





























































!! * Modules used
   USE parameters, ONLY : lchain,rsh,rlg,riolg,riosh,kmax,imin,imax,jmin,jmax,     &
                          limin,limax,ljmin,ljmax,ishift_obc_west,ishift_obc_east, &
                          jshift_obc_south,jshift_obc_north


   IMPLICIT NONE


   INCLUDE 'mpif.h'


   PRIVATE

!! * Accessibility

   PUBLIC diagdyn_alloc,   &       ! called by alloc.F90
          diagdyn_dealloc, &       ! called by dealloc.F90
          diagdyn_init,   &       ! called by main.F90
          diagdyn_output          ! called by step.F90

!! * Shared module variables

   REAL(kind=rlg),ALLOCATABLE,DIMENSION(:,:,:),PUBLIC :: diagx_cor,diagx_nlt,diagx_ipg,diagx_hdif,diagx_epg, &
                                                         diagx_corcor,diagx_vdif,diagx_dudt,diag_boz,        &
                                                         diagx_nlt_ududx,diagx_nlt_vdudy,diagx_nlt_wdudz,    &
                                                         diagx_hdif_dnududx2,diagx_hdif_dnududy2
   REAL(kind=rlg),ALLOCATABLE,DIMENSION(:,:,:),PUBLIC :: diagy_cor,diagy_nlt,diagy_ipg,diagy_hdif,diagy_epg, &
                                                         diagy_corcor,diagy_vdif,diagy_dudt,                 &
                                                         diagy_nlt_udvdx,diagy_nlt_vdvdy,diagy_nlt_wdvdz,    &
                                                         diagy_hdif_dnudvdx2,diagy_hdif_dnudvdy2
   LOGICAL,ALLOCATABLE,DIMENSION(:,:),PUBLIC          :: mask_diagx,mask_diagy

!! * Private variables

   REAL(KIND=rlg)  :: tdebdiagdyn,tfindiagdyn,dtsordiagdyn,tnextout,tcumul,tprecmoy,dtmoy
   REAL(kind=rsh),ALLOCATABLE,DIMENSION(:,:,:) :: diagx_cor_save,diagx_nlt_save,diagx_ipg_save,diagx_hdif_save,diagx_epg_save,&
                                                  diagx_corcor_save,diagx_vdif_save,diagx_dudt_save,diagx_velnud_save, &
                                                  diagx_nlt_ududx_save,diagx_nlt_vdudy_save,diagx_nlt_wdudz_save, &
                                                  diag_boz_save,diagx_hdif_dnududx2_save,diagx_hdif_dnududy2_save
   REAL(kind=rsh),ALLOCATABLE,DIMENSION(:,:,:) :: diagy_cor_save,diagy_nlt_save,diagy_ipg_save,diagy_hdif_save,diagy_epg_save,&
                                                  diagy_corcor_save,diagy_vdif_save,diagy_dudt_save,diagy_velnud_save,           &
                                                  diagy_nlt_udvdx_save,diagy_nlt_vdvdy_save,diagy_nlt_wdvdz_save,             &
                                                  diagy_hdif_dnudvdx2_save,diagy_hdif_dnudvdy2_save
   REAL(kind=rsh),ALLOCATABLE,DIMENSION(:,:) :: diag_sshnud_save

   LOGICAL,ALLOCATABLE,DIMENSION(:,:)  :: mask_h0_diag,mask_hx_diag,mask_hy_diag

   LOGICAL     :: l_moy_save
   INTEGER     :: nbmoy
   CHARACTER(LEN=lchain) :: filediag
   INTEGER      :: iextmindiag,iextmaxdiag,stepidiag,jextmindiag,jextmaxdiag,stepjdiag,kextmindiag,kextmaxdiag,stepkdiag
   INTEGER      :: igextmindiag,igextmaxdiag,jgextmindiag,jgextmaxdiag
   LOGICAL      :: l_trend,l_nlt,l_cor,l_epgrad,l_ipg,l_hdif,l_epg,l_corcor,l_vdif,l_dudt,l_boz

   CHARACTER(len=lchain) :: name_outdiag_h0,name_outdiag_hx,name_outdiag_hy, &
                            name_outdiagx_cor,name_outdiagx_nlt,name_outdiagx_ipg,name_outdiagx_hdif, &
                            name_outdiagx_epg,name_outdiagx_corcor,name_outdiagx_vdif,name_outdiagx_dudt, &
                            name_outdiagx_velnud,name_outdiag_sshnud,name_outdiag_boz,name_outdiagx_nlt_ududx, &
                            name_outdiagx_nlt_vdudy,name_outdiagx_nlt_wdudz,name_outdiagx_hdif_dnududx2,  &
                            name_outdiagx_hdif_dnududy2

   CHARACTER(len=lchain) :: name_outdiagy_cor,name_outdiagy_nlt,name_outdiagy_ipg,name_outdiagy_hdif, &
                            name_outdiagy_epg,name_outdiagy_corcor,name_outdiagy_vdif,name_outdiagy_dudt,name_outdiagy_velnud, &
                            name_outdiagy_nlt_udvdx,name_outdiagy_nlt_vdvdy,name_outdiagy_nlt_wdvdz, &
                            name_outdiagy_hdif_dnudvdx2,name_outdiagy_hdif_dnudvdy2

   REAL(KIND=riosh)      :: riog_valid_min_diag,riog_valid_max_diag

  CONTAINS

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

  SUBROUTINE diagdyn_alloc

!&E---------------------------------------------------------------------
!&E                 ***  ROUTINE diagdyn_alloc  ***
!&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       !  2009-09 (S. Michel, HOCER) correct i dimension of the arrays
!&E
!&E---------------------------------------------------------------------
!! * Modules used

!! * Local declarations

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

   ALLOCATE(diagx_cor(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagx_ipg(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagx_epg(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagx_corcor(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagx_vdif(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagx_dudt(1:kmax,limin:limax,ljmin:ljmax) )

   ALLOCATE(diagx_nlt(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagx_hdif(1:kmax,limin:limax,ljmin:ljmax) )




   ALLOCATE(diagy_cor(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagy_ipg(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagy_epg(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagy_corcor(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagy_vdif(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagy_dudt(1:kmax,limin:limax,ljmin:ljmax) )

   ALLOCATE(diagy_nlt(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagy_hdif(1:kmax,limin:limax,ljmin:ljmax) )


   ALLOCATE(diagx_cor_save(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagx_ipg_save(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagx_epg_save(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagx_corcor_save(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagx_vdif_save(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagx_dudt_save(1:kmax,limin:limax,ljmin:ljmax) )

   ALLOCATE(diagx_nlt_save(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagx_hdif_save(1:kmax,limin:limax,ljmin:ljmax) )




   ALLOCATE(diagy_cor_save(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagy_ipg_save(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagy_epg_save(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagy_corcor_save(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagy_vdif_save(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagy_dudt_save(1:kmax,limin:limax,ljmin:ljmax) )

   ALLOCATE(diagy_nlt_save(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diagy_hdif_save(1:kmax,limin:limax,ljmin:ljmax) )




   ALLOCATE(diag_boz(1:kmax,limin:limax,ljmin:ljmax) )
   ALLOCATE(diag_boz_save(1:kmax,limin:limax,ljmin:ljmax) )

   ALLOCATE(mask_h0_diag(limin:limax,ljmin:ljmax) )
   ALLOCATE(mask_hx_diag(limin:limax,ljmin:ljmax) )
   ALLOCATE(mask_hy_diag(limin:limax,ljmin:ljmax) )
   ALLOCATE(mask_diagx(limin:limax,ljmin:ljmax) )
   ALLOCATE(mask_diagy(limin:limax,ljmin:ljmax) )

   diagx_cor(:,:,:)=0.0_rlg
   diagx_ipg(:,:,:)=0.0_rlg
   diagx_epg(:,:,:)=0.0_rlg
   diagx_corcor(:,:,:)=0.0_rlg
   diagx_vdif(:,:,:)=0.0_rlg
   diagx_dudt(:,:,:)=0.0_rlg

   diagx_nlt(:,:,:)=0.0_rlg
   diagx_hdif(:,:,:)=0.0_rlg




   diagy_cor(:,:,:)=0.0_rlg
   diagy_ipg(:,:,:)=0.0_rlg
   diagy_epg(:,:,:)=0.0_rlg
   diagy_corcor(:,:,:)=0.0_rlg
   diagy_vdif(:,:,:)=0.0_rlg
   diagy_dudt(:,:,:)=0.0_rlg

   diagy_nlt(:,:,:)=0.0_rlg
   diagy_hdif(:,:,:)=0.0_rlg




  END SUBROUTINE diagdyn_alloc

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

  SUBROUTINE diagdyn_dealloc

!&E---------------------------------------------------------------------
!&E                 ***  ROUTINE diagdyn_dealloc  ***
!&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       !
!&E
!&E---------------------------------------------------------------------
!! * Modules used

!! * Arguments

!! * Local declarations

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

   DEALLOCATE(diagx_cor,diagx_ipg,diagx_epg,diagx_corcor,diagx_vdif,diagx_dudt)


   DEALLOCATE(diagx_nlt,diagx_hdif,diagy_nlt,diagy_hdif)
   DEALLOCATE(diagx_nlt_save,diagx_hdif_save,diagy_nlt_save,diagy_hdif_save)


   DEALLOCATE(diagy_cor,diagy_ipg,diagy_epg,diagy_corcor,diagy_vdif,diagy_dudt)

   DEALLOCATE(diagx_cor_save,diagx_ipg_save,diagx_epg_save, &
              diagx_corcor_save,diagx_vdif_save,diagx_dudt_save)

   DEALLOCATE(diagy_cor_save,diagy_ipg_save,diagy_epg_save, &
              diagy_corcor_save,diagy_vdif_save,diagy_dudt_save)





   DEALLOCATE(diag_boz,diag_boz_save)

   DEALLOCATE(mask_h0_diag,mask_hx_diag,mask_hy_diag,mask_diagx,mask_diagy)

  END SUBROUTINE diagdyn_dealloc

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

  SUBROUTINE diagdyn_init

!&E---------------------------------------------------------------------
!&E                 ***  ROUTINE diagdyn_init  ***
!&E
!&E ** Purpose : Initialize diagnostics of dynamical terms
!&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 :  !  2011-04  (P. Lazure)  Original code
!&E               !  2014-09  (S. Petton)  Saving of a zoom only
!&E
!&E---------------------------------------------------------------------

!! * Modules used
   USE comvars2d,   ONLY : ierrorlog
   USE output,      ONLY : tdebsor,tfinsor,pasor

!! * Local declarations
   INTEGER :: i,j,IERR_MPI

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

   !PRINT*, 'enter DIAG_DYN_INIT routine'


!   tdebdiagdyn=tool_datosec('01/01/1900 00:00:00')
!   tdebdiagdyn=MAX(tdebdiagdyn, tdeb3d)
   tdebdiagdyn=tdebsor(1)
!   tfindiagdyn=tool_datosec('01/01/1904 01:00:00')
   tfindiagdyn=tfinsor(1)
!   dtsordiagdyn=24.0_rlg*3600.0_rlg
!   dtsordiagdyn=3.*3600.0_rlg
   dtsordiagdyn=pasor(1)
   l_moy_save=.TRUE.

   l_cor=.TRUE.
   l_ipg=.TRUE.
   l_nlt=.TRUE.
   l_hdif=.TRUE.
   l_epg=.TRUE.
   l_corcor=.TRUE.
   l_vdif=.TRUE.
   l_dudt=.TRUE.
   l_boz=.TRUE.

   igextmindiag=imin+ishift_obc_west
   igextmaxdiag=imax-ishift_obc_east
   iextmindiag=imin+ishift_obc_west
   iextmaxdiag=imax-ishift_obc_east
   stepidiag=1
   jextmindiag=jmin+jshift_obc_south
   jextmaxdiag=jmax-jshift_obc_north
   jgextmindiag=jmin+jshift_obc_south
   jgextmaxdiag=jmax-jshift_obc_north
   stepjdiag=1

! If we do not use the parallel writing, the spatial extension
! must be bounded by limin and limax
     DO WHILE ( iextmindiag < limin )
       iextmindiag=iextmindiag+stepidiag
     END DO
     iextmaxdiag=MIN(limax,iextmaxdiag)
     DO WHILE ( jextmindiag < ljmin )
       jextmindiag=jextmindiag+stepjdiag
     END DO
     jextmaxdiag=MIN(ljmax,jextmaxdiag)

   kextmindiag=1  ! should not be changed
   kextmaxdiag=kmax ! id
   stepkdiag=1 ! id


! Output names
! ------------

   riog_valid_min_diag=-10000.0_riosh
   riog_valid_max_diag=10000.0_riosh

   name_outdiag_h0='h0'
   name_outdiag_hx='hx'
   name_outdiag_hy='hy'


   name_outdiagx_cor='cor_x'
   name_outdiagx_ipg='ipg_x'
   name_outdiagx_nlt='nlt_x'
   name_outdiagx_nlt_ududx='nlt_x_ududx'
   name_outdiagx_nlt_vdudy='nlt_x_vdudy'
   name_outdiagx_nlt_wdudz='nlt_x_wdudz'
   name_outdiagx_hdif='hdif_x'
   name_outdiagx_hdif_dnududx2='hdif_x_dnududx2'
   name_outdiagx_hdif_dnududy2='hdif_x_dnududy2'
   name_outdiagx_epg='epg_x'
   name_outdiagx_corcor='cor2_x'
   name_outdiagx_vdif='vdif_x'
   name_outdiagx_dudt='dudt_x'

   name_outdiagy_cor='cor_y'
   name_outdiagy_ipg='ipg_y'
   name_outdiagy_nlt='nlt_y'
   name_outdiagy_nlt_udvdx='nlt_y_udvdx'
   name_outdiagy_nlt_vdvdy='nlt_y_vdvdy'
   name_outdiagy_nlt_wdvdz='nlt_y_wdvdz'
   name_outdiagy_hdif='hdif_y'
   name_outdiagy_hdif_dnudvdx2='hdif_y_dnudvdx2'
   name_outdiagy_hdif_dnudvdy2='hdif_y_dnudvdy2'
   name_outdiagy_epg='epg_y'
   name_outdiagy_corcor='cor2_y'
   name_outdiagy_vdif='vdif_y'
   name_outdiagy_dudt='dudt_y'

   name_outdiagx_velnud='nudgu3d'
   name_outdiagy_velnud='nudgv3d'
   name_outdiag_sshnud='nudssh'

   name_outdiag_boz='boz'


   tnextout=tdebdiagdyn+dtsordiagdyn


! Verifications
! -------------

   IF ( jgextmindiag > jgextmaxdiag .OR. igextmindiag > igextmaxdiag ) THEN
     WRITE(ierrorlog,*) ' '
     WRITE(ierrorlog,*) ' Routine DIAGDYN_INIT in DIAGDYN.F90'
     WRITE(ierrorlog,*) ' '
     WRITE(ierrorlog,*) 'Check your domain : because'
     WRITE(ierrorlog,*) ' jgextmindiag > jgextmaxdiag '
     WRITE(ierrorlog,*) '           or'
     WRITE(ierrorlog,*) ' igextmindiag > igextmaxdiag '
     WRITE(ierrorlog,*) ' '
     WRITE(ierrorlog,*) ' The simulation is stopped '
     WRITE(ierrorlog,*) ' '
     CALL MPI_FINALIZE(IERR_MPI)
     STOP
   END IF

! Masks to estimate the diagnoses
! -------------------------------

   mask_diagx(:,:)=.FALSE.
   DO j=jextmindiag,jextmaxdiag,stepjdiag
     DO i=iextmindiag,iextmaxdiag,stepidiag
       IF ( i >= limin .AND. i <= limax .AND. j >= ljmin .AND. j <= ljmax) THEN
         mask_diagx(i,j)=.TRUE.
       END IF
     END DO
   END DO

   mask_diagy(:,:)=.FALSE.
   DO j=jextmindiag,jextmaxdiag,stepjdiag
     DO i=iextmindiag,iextmaxdiag,stepidiag
       IF ( i >= limin .AND. i <= limax .AND. j >= ljmin .AND. j <= ljmax) THEN
         mask_diagy(i,j)=.TRUE.
       END IF
     END DO
   END DO

   !PRINT*, 'exit DIAG_DYN_INIT routine'

  END SUBROUTINE diagdyn_init

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

  SUBROUTINE diagdyn_output_create

!&E---------------------------------------------------------------------
!&E                 ***  ROUTINE diagdyn_output_create  ***
!&E
!&E ** Purpose : Create output file containing dynamical terms
!&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 :  !  2011-04  (P. Lazure)  Original code
!&E               !  2014-09  (S. Petton)  Saving of a zoom only
!&E
!&E---------------------------------------------------------------------

   USE comvars2d,   ONLY : lon2d_g,lat2d_g,lon2d,lat2d,out_dir,suffix,t
   USE parameters,  ONLY : rsh,rlg,riosh,kmax,imin,imax,limin,limax,ljmin,ljmax,  &
                           rg_valmanq_io,ljminp1,riolg,valmanq,lchain
   USE ionc4,       ONLY : ionc4_createfile,ionc4_gatt,ionc4_createvar,ionc4_open,&
                           ionc4_write,ionc4_write_sig,ionc4_write_sigw,          &
                           ionc4_write_xy,ionc4_write_time,ionc4_vatt_chain
   USE output,      ONLY : l_out_nc4par,l_out_pack,l_out_zaxis_w,fileresu, &
                           fileresu_verif,l_filebydate,nb_out
   USE outtools,    ONLY : out_grid
   USE outnetcdf,   ONLY : out_netcdf_conv
   USE comvars2d,   ONLY : hx,hy,h0,h0fond,riog_valid_min_h0,riog_valid_max_h0,   &
                           name_out_xe,name_out_h0
   USE comvars2d,   ONLY : ig,id,dx_u,dy_v,jb,jh,htx,hty,tdeb3d
   USE comvars3d,   ONLY : sig
   USE toolmpi, ONLY : ADD_TAG,IPROC,NPROCS,ADD_ALL_MPI_INT,MPI_COMM_MARS



!! * Local declarations

   REAL(kind=riosh),DIMENSION(limin:limax,ljmin:ljmax) :: tecrit
   REAL(KIND=riosh),DIMENSION(0:kmax+1)            :: sigecrit
   REAL(KIND=riosh),DIMENSION(limin:limax,ljmin:ljmax)      :: h0_local,hx_local,hy_local
   REAL(KIND=riosh),DIMENSION(kmax,limin:limax,ljmin:ljmax) :: tzecrit
   REAL(KIND=riolg),DIMENSION(imin:imax)                    :: lon
   REAL(kind=rlg)        :: tool_datosec
   CHARACTER (len=19)    :: tool_sectodat,cdateout
   INTEGER :: i,j,IERR_MPI
   INTEGER :: grp_world,grp_zoom,comm_zoom,ncpu,ip
   INTEGER, DIMENSION(:), ALLOCATABLE  :: numrank,listrank
   INTEGER, DIMENSION(1) :: listrank_nb
   REAL(KIND=riosh)                                              :: vecrit

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

   !PRINT*, 'enter DIAG_DYN_OUTPUT_CREATE routine'


! Files
! -----

   IF (l_filebydate(1)) THEN
     cdateout = tool_sectodat(t)
     filediag=TRIM(out_dir)//'diag_'//TRIM(suffix)

     filediag=TRIM(filediag)//'_'//cdateout(7:10)//cdateout(4:5)//cdateout(1:2)//'T'
     filediag=TRIM(filediag)//cdateout(12:13)//cdateout(15:16)//cdateout(18:19)//'Z'
     filediag=TRIM(filediag)//'.nc'
     IF (.NOT. l_out_nc4par) CALL ADD_TAG(filediag,LEN(filediag))
   ELSE
     !IF (.TRUE.) THEN
       filediag=TRIM(out_dir)//'diag_'//TRIM(suffix)//'.nc'
     !ELSE
       ! filediag=TRIM(out_dir)//'diag_zoom'//TRIM(Agrif_Cfixed())//'_'//TRIM(suffix)//'.nc'
     !ENDIF

   END IF



! Masks
! -----

   h0_local(limin:limax,ljmin:ljmax)=h0(limin:limax,ljmin:ljmax)
   hx_local(limin:limax,ljmin:ljmax)=hx(limin:limax,ljmin:ljmax)
   hy_local(limin:limax,ljmin:ljmax)=hy(limin:limax,ljmin:ljmax)

   mask_h0_diag(:,:)=.TRUE.
   DO j=ljmin,ljmax
     DO i=MAX(limin,ig(j)),MIN(limax,id(j))
       IF ( h0_local(i,j)> h0fond ) mask_h0_diag(i,j)=.FALSE.
     END DO
   END DO

! pascal
! mouais pas clair c t affaire
   mask_hx_diag(:,:)=.TRUE.
   DO j=ljminp1,ljmax
     DO i=MAX(limin,ig(j)-1),MIN(limax,id(j))
       IF ( hx_local(i,j) /= -valmanq -h0fond ) mask_hx_diag(i,j)=.FALSE.
     END DO
   END DO
!mask_hx_diag(:,:)=.FALSE.

   mask_hy_diag(:,:)=.TRUE.
   DO j=ljmin,ljmax
     DO i=MAX(limin,ig(j)),MIN(limax,id(j))
       IF ( (j >= jb(i)-1 .AND. j <= jh(i)) .AND. (hy_local(i,j) /= -valmanq -h0fond) ) mask_hy_diag(i,j)=.FALSE.
     END DO
   END DO


  IF ( TRIM(filediag) /= TRIM(fileresu_verif(nb_out)) .AND. TRIM(filediag) /= TRIM(fileresu_verif(nb_out+1)) ) THEN

   comm_zoom=0

     comm_zoom=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 ( jextmindiag <= jextmaxdiag .AND. iextmindiag <= iextmaxdiag ) 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 ( jextmindiag <= jextmaxdiag .AND. iextmindiag <= iextmaxdiag ) 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 ( 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,comm_zoom,ierr_mpi)
     ENDIF



! to work with active cpus over a small area
   IF ( jextmindiag <= jextmaxdiag .AND. iextmindiag <= iextmaxdiag ) THEN

   sigecrit(:)=sig(:)
   IF ( l_out_nc4par ) THEN
! Create a global file because of the parallel writing
     CALL ionc4_createfile(filediag,lon2d_g(igextmindiag:igextmaxdiag:stepidiag,jgextmindiag:jgextmaxdiag:stepjdiag),&
                                    lat2d_g(igextmindiag:igextmaxdiag:stepidiag,jgextmindiag:jgextmaxdiag:stepjdiag),&
                                    igextmindiag,igextmaxdiag,stepidiag, &
                                    jgextmindiag,jgextmaxdiag,stepjdiag, &
                                    kextmindiag,kextmaxdiag,stepkdiag,   &
                                    kmax,0,l_out_nc4par = .TRUE.,comm_active=comm_zoom, &
                                    zaxis_var=sigecrit(1:kmax))
   ELSE
! Create a local file bounded by limin,limax,ljmin,ljmax
     CALL ADD_TAG(filediag,LEN(filediag))
     CALL ionc4_createfile(filediag,lon2d(iextmindiag:iextmaxdiag:stepidiag,jextmindiag:jextmaxdiag:stepjdiag),&
                                    lat2d(iextmindiag:iextmaxdiag:stepidiag,jextmindiag:jextmaxdiag:stepjdiag),&
                                    iextmindiag,iextmaxdiag,stepidiag, &
                                    jextmindiag,jextmaxdiag,stepjdiag, &
                                    kextmindiag,kextmaxdiag,stepkdiag, &
                                    kmax,0,zaxis_var=sigecrit(1:kmax))
   END IF

! save attribute of the spatial extension of the domain
! -----------------------------------------------------
   CALL ionc4_gatt(filediag,'global_imin',igextmindiag)
   CALL ionc4_gatt(filediag,'global_imax',igextmaxdiag)
   CALL ionc4_gatt(filediag,'global_jmin',jgextmindiag)
   CALL ionc4_gatt(filediag,'global_jmax',jgextmaxdiag)

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

! Add the vertical axis
!CALL out_zaxis_w(filediag,kextmindiag,kextmaxdiag,stepkdiag,comm_zoom)

   CALL ionc4_createvar(filediag,'SIG',' ','sigma variable',&

                        standard_name= 'ocean_sigma_variable',&

                        dims='z',valid_min= -1.0_riosh,valid_max= 0.0_riosh,fill_value= rg_valmanq_io)
   CALL ionc4_write_sig(filediag,'SIG',sigecrit,kmax)



! Add attributes to follow OCO netCDF convention
   CALL out_netcdf_conv(filediag,1,comm_zoom,igextmindiag,igextmaxdiag,jgextmindiag,jgextmaxdiag)

   IF (l_out_nc4par .AND. comm_zoom/=MPI_COMM_MARS) CALL MPI_COMM_FREE(comm_zoom,ierr_mpi)


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

! h0
! --


! avec agrif on ne peut compiler car The type of the actual argument differs from the type of the dummy argument.   [ARRAY0]
!  &,  Agrif_tabvars(606) % var% array0 ,  Agrif_tabvars(605) &
! Pas de souci si on utilise riog_valid_min_diag a la place de riog_valid_min_h0 !!!!!!!!!!
! Cela m exaspere, et c est un pb recurant, voir mail du 3/3/2010 a Romain Ramel (PREVI_OMPMPI)
! Laurent peux-tu faire qq chose STP ou donner une explication rationnelle !!!!

   CALL ionc4_createvar(filediag,name_outdiag_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(limin:limax,ljmin:ljmax)=h0_local(limin:limax,ljmin:ljmax)
   CALL ionc4_write_xy(filediag,TRIM(name_outdiag_h0),tecrit,limin,limax,ljmin,ljmax,&
                       fill_value=-rg_valmanq_io,l_out_nc4par = l_out_nc4par )


! hx
! --
   CALL ionc4_createvar(filediag,name_outdiag_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_local(:,:)
   CALL ionc4_write_xy(filediag,TRIM(name_outdiag_hx),tecrit,limin,limax,ljmin,ljmax ,&
                       fill_value = -REAL(rg_valmanq_io+h0fond,riosh),l_out_nc4par = l_out_nc4par)

! hy
! --

   CALL ionc4_createvar(filediag,name_outdiag_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_local(:,:)
   CALL ionc4_write_xy(filediag,TRIM(name_outdiag_hy),tecrit,limin,limax,ljmin,ljmax,&
                       fill_value = -REAL(rg_valmanq_io+h0fond,riosh),l_out_nc4par = l_out_nc4par)



   END IF  ! test on active cpus
  END IF  ! test if file different from fileresu_verif

  IF ( jextmindiag <= jextmaxdiag .AND. iextmindiag <= iextmaxdiag .AND. &
       TRIM(filediag) /= TRIM(fileresu_verif(nb_out+1)) ) THEN

   IF ( TRIM(filediag) == TRIM(fileresu_verif(nb_out)) ) &
     CALL ionc4_open(fileresu_verif(nb_out),l_out_nc4par) 

   IF (l_cor) THEN
     CALL ionc4_createvar(filediag,name_outdiagx_cor,units= 'm s-2',&
                          long_name = 'Coriolis term at u location',&
                          standard_name = 'first_guess_of_x_coriolis_acceleration_at_u_location',&
                          fill_value = rg_valmanq_io,&
                          valid_min = riog_valid_min_diag,&
                          valid_max = riog_valid_max_diag,&
                          dims= 'zxyt',&
                          l_pack = l_out_pack,&
                          l_out_nc4par = l_out_nc4par)

     CALL ionc4_createvar(filediag,name_outdiagy_cor,units= 'm s-2',&
                          long_name = 'Coriolis term at v location',&
                          standard_name = 'first_guess_of_y_coriolis_acceleration_at_v_location',&
                          fill_value = rg_valmanq_io,&
                          valid_min = riog_valid_min_diag,&
                          valid_max = riog_valid_max_diag,&
                          dims= 'zxyt',&
                          l_pack = l_out_pack,&
                          l_out_nc4par = l_out_nc4par)
   END IF

   IF (l_ipg) THEN
     CALL ionc4_createvar(filediag,name_outdiagx_ipg,units= 'm s-2',&
                          long_name = '1/rho dPi/dx : internal pressure gradient at u location',&
                          standard_name = 'x_internal_pressure_gradient_at_u_location',&
                          fill_value = rg_valmanq_io,&
                          valid_min = riog_valid_min_diag,&
                          valid_max = riog_valid_max_diag,&
                          dims= 'zxyt',&
                          l_pack = l_out_pack,&
                          l_out_nc4par = l_out_nc4par)

     CALL ionc4_createvar(filediag,name_outdiagy_ipg,units= 'm s-2',&
                          long_name = '1/rho dPi/dy : internal pressure gradient at v location',&
                          standard_name = 'y_internal_pressure_gradient_at_v_location',&
                          fill_value = rg_valmanq_io,&
                          valid_min = riog_valid_min_diag,&
                          valid_max = riog_valid_max_diag,&
                          dims= 'zxyt',&
                          l_pack = l_out_pack,&
                          l_out_nc4par = l_out_nc4par)
   END IF

   IF (l_nlt) THEN

     CALL ionc4_createvar(filediag,name_outdiagx_nlt,units= 'm s-2',&
                          long_name = 'udu/dx+vdu/dy+wdu/dz : non linear terms at u location',&
                          standard_name = 'non_linear_advection_of_sea_water_x_velocity_at_u_location',&
                          fill_value = rg_valmanq_io,&
                          valid_min = riog_valid_min_diag,&
                          valid_max = riog_valid_max_diag,&
                          dims= 'zxyt',&
                          l_pack = l_out_pack,&
                          l_out_nc4par = l_out_nc4par)




     CALL ionc4_createvar(filediag,name_outdiagy_nlt,units= 'm s-2',&
                          long_name = 'udv/dx+vdv/dy+wdv/dz : non linear terms at v location',&
                          standard_name = 'non_linear_advection_of_sea_water_y_velocity_at_v_location',&
                          fill_value = rg_valmanq_io,&
                          valid_min = riog_valid_min_diag,&
                          valid_max = riog_valid_max_diag,&
                          dims= 'zxyt',&
                          l_pack = l_out_pack,&
                          l_out_nc4par = l_out_nc4par)


   END IF

   IF (l_hdif) THEN

     CALL ionc4_createvar(filediag,name_outdiagx_hdif,units= 'm s-2',&
                          long_name = 'd nu(du/dx) /dx + d nu(du/dy) /dy : horizontal diffusion at u location',&
                          standard_name = 'horizontal_diffusion_of_sea_water_x_velocity_at_u_location',&
                          fill_value = rg_valmanq_io,&
                          valid_min = riog_valid_min_diag,&
                          valid_max = riog_valid_max_diag,&
                          dims= 'zxyt',&
                          l_pack = l_out_pack,&
                          l_out_nc4par = l_out_nc4par)




     CALL ionc4_createvar(filediag,name_outdiagy_hdif,units= 'm s-2',&
                          long_name = 'd nu(dv/dx) /dx + d nu(dv/dy) /dy : horizontal diffusion at v location',&
                          standard_name = 'horizontal_diffusion_of_sea_water_y_velocity_at_v_location',&
                          fill_value = rg_valmanq_io,&
                          valid_min = riog_valid_min_diag,&
                          valid_max = riog_valid_max_diag,&
                          dims= 'zxyt',&
                          l_pack = l_out_pack,&
                          l_out_nc4par = l_out_nc4par)


   END IF

   IF (l_epg) THEN
     CALL ionc4_createvar(filediag,name_outdiagx_epg,units= 'm s-2',&
                          long_name = ' - g d ssh / dx - 1/rho d Pa /dx : external pressure gradient at u location',&
                          standard_name = 'x_external_pressure_gradient_at_u_location',&
                          fill_value = rg_valmanq_io,&
                          valid_min = riog_valid_min_diag,&
                          valid_max = riog_valid_max_diag,&
                          dims= 'zxyt',&
                          l_pack = l_out_pack,&
                          l_out_nc4par = l_out_nc4par)

     CALL ionc4_createvar(filediag,name_outdiagy_epg,units= 'm s-2',&
                          long_name = ' - g d ssh / dy - 1/rho d Pa /dy : external pressure gradient at v location',&
                          standard_name = 'y_external_pressure_gradient_at_v_location',&
                          fill_value = rg_valmanq_io,&
                          valid_min = riog_valid_min_diag,&
                          valid_max = riog_valid_max_diag,&
                          dims= 'zxyt',&
                          l_pack = l_out_pack,&
                          l_out_nc4par = l_out_nc4par)
   END IF

   IF (l_corcor) THEN
     CALL ionc4_createvar(filediag,name_outdiagx_corcor,units= 'm s-2',&
                          long_name = '-f*v : corrected Coriolis term at u location',&
                          standard_name = 'x_coriolis_acceleration_at_u_location',&
                          fill_value = rg_valmanq_io,&
                          valid_min = riog_valid_min_diag,&
                          valid_max = riog_valid_max_diag,&
                          dims= 'zxyt',&
                          l_pack = l_out_pack,&
                          l_out_nc4par = l_out_nc4par)

     CALL ionc4_createvar(filediag,name_outdiagy_corcor,units= 'm s-2',&
                          long_name = 'f*u : corrected Coriolis term at v location',&
                          standard_name = 'y_coriolis_acceleration_at_v_location',&
                          fill_value = rg_valmanq_io,&
                          valid_min = riog_valid_min_diag,&
                          valid_max = riog_valid_max_diag,&
                          dims= 'zxyt',&
                          l_pack = l_out_pack,&
                          l_out_nc4par = l_out_nc4par)
   END IF

   IF (l_vdif) THEN
     CALL ionc4_createvar(filediag,name_outdiagx_vdif,units= 'm s-2',&
                          long_name = 'dz(nz*du/dz)/dz : vertical diffusion at u location',&
                          standard_name = 'upward_diffusion_of_sea_water_x_velocity_at_u_location',&
                          fill_value = rg_valmanq_io,&
                          valid_min = riog_valid_min_diag,&
                          valid_max = riog_valid_max_diag,&
                          dims= 'zxyt',&
                          l_pack = l_out_pack,&
                          l_out_nc4par = l_out_nc4par)

     CALL ionc4_createvar(filediag,name_outdiagy_vdif,units= 'm s-2',&
                          long_name = 'dz(nz*dv/dz)/dz : vertical diffusion at v location',&
                          standard_name = 'upward_diffusion_of_sea_water_y_velocity_at_v_location',&
                          fill_value = rg_valmanq_io,&
                          valid_min = riog_valid_min_diag,&
                          valid_max = riog_valid_max_diag,&
                          dims= 'zxyt',&
                          l_pack = l_out_pack,&
                          l_out_nc4par = l_out_nc4par)
   END IF

   IF (l_dudt) THEN
     CALL ionc4_createvar(filediag,name_outdiagx_dudt,units= 'm s-2',&
                          long_name = 'du/dt : time trend at u location',&
                          standard_name = 'tendency_of_sea_water_x_velocity_at_u_location',&
                          fill_value = rg_valmanq_io,&
                          valid_min = riog_valid_min_diag,&
                          valid_max = riog_valid_max_diag,&
                          dims= 'zxyt',&
                          l_pack = l_out_pack,&
                          l_out_nc4par = l_out_nc4par)

     CALL ionc4_createvar(filediag,name_outdiagy_dudt,units= 'm s-2',&
                          long_name = 'dv/dt : time trend at v location',&
                          standard_name = 'tendency_of_sea_water_y_velocity_at_v_location',&
                          fill_value = rg_valmanq_io,&
                          valid_min = riog_valid_min_diag,&
                          valid_max = riog_valid_max_diag,&
                          dims= 'zxyt',&
                          l_pack = l_out_pack,&
                          l_out_nc4par = l_out_nc4par)

   END IF



   IF (l_boz) &
   CALL ionc4_createvar(filediag,name_outdiag_boz,units= 'm s-2',&
                        long_name = 'Buoyancy of sea water at atmospheric pressure',&
                        standard_name = 'sea_water_buoyancy',&
                        fill_value = rg_valmanq_io,&
                        valid_min = riog_valid_min_diag,&
                        valid_max = riog_valid_max_diag,&
                        dims= 'zxyt',&
                        l_pack = l_out_pack,&
                        l_out_nc4par = l_out_nc4par)

   END IF  ! ends test on iextmin...

   fileresu_verif(nb_out+1)=TRIM(filediag)

   !PRINT*, 'exit DIAG_DYN_OUTPUT_CREATE routine'

  END SUBROUTINE diagdyn_output_create

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

  SUBROUTINE diagdyn_output

!&E---------------------------------------------------------------------
!&E                 ***  ROUTINE outputdiagdyn  ***
!&E
!&E ** Purpose : Save diagnostics of dynamical terms
!&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 :  !  2011-04  (P. Lazure)  Original code
!&E               !  2014-09  (S. Petton)  Saving of a zoom only
!&E
!&E---------------------------------------------------------------------

!! * Modules used
   USE parameters,  ONLY : rsh,rlg,riosh,kmax,imin,imax,limin,limax,ljmin,ljmax,  &
                           rg_valmanq_io,ljminp1,riolg,valmanq
   USE outtools,    ONLY : out_write_xytvar,out_write_zxytvar
   USE output,      ONLY : fileresu_verif,nb_out
   USE ionc4,       ONLY : ionc4_write_time,ionc4_sync
   USE output,      ONLY : l_out_nc4par
   USE comvars2d,   ONLY : hx,hy,h0,h0fond,riog_valid_min_h0,riog_valid_max_h0
   USE comvars2d,   ONLY : ig,id,dx_u,dy_v,jb,jh,htx,hty
   USE comvars2d,   ONLY : t,h0fond
   USE comvars3d,   ONLY : sig
   USE toolmpi, ONLY : MASTER

!! * Local declarations
   REAL(kind=riosh),DIMENSION(imin:imax,jmin:jmax) :: tecrit
   REAL(KIND=riosh),DIMENSION(limin:limax,ljmin:ljmax)      :: h0_local,hx_local,hy_local
   REAL(KIND=riosh),DIMENSION(kmax,limin:limax,ljmin:ljmax) :: tzecrit
   REAL(KIND=riolg),DIMENSION(imin:imax)                    :: lon
   INTEGER :: i,j
   CHARACTER (len=19) :: tool_sectodat

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

   !PRINT*, 'enter OUTPUT_DIAG_DYN routine'

   IF ( t < tdebdiagdyn .OR. t > tfindiagdyn )  RETURN

!print*, ' prochaine ecriture:',tool_sectodat(tnextout)

   IF (l_moy_save) THEN
     nbmoy=nbmoy+1

     IF (nbmoy ==1) THEN
       tprecmoy=t
       tcumul=0
     ELSE
       dtmoy=t-tprecmoy
       tcumul=tcumul+dtmoy

       diagx_cor_save=diagx_cor_save+diagx_cor*dtmoy
       diagx_ipg_save=diagx_ipg_save+diagx_ipg*dtmoy
       diagx_epg_save=diagx_epg_save+diagx_epg*dtmoy
       diagx_corcor_save=diagx_corcor_save+diagx_corcor*dtmoy
       diagx_vdif_save=diagx_vdif_save+diagx_vdif*dtmoy
       diagx_dudt_save=diagx_dudt_save+diagx_dudt*dtmoy

       diagx_nlt_save=diagx_nlt_save+diagx_nlt*dtmoy
       diagx_hdif_save=diagx_hdif_save+diagx_hdif*dtmoy




       diagy_cor_save=diagy_cor_save+diagy_cor*dtmoy
       diagy_ipg_save=diagy_ipg_save+diagy_ipg*dtmoy
       diagy_epg_save=diagy_epg_save+diagy_epg*dtmoy
       diagy_corcor_save=diagy_corcor_save+diagy_corcor*dtmoy
       diagy_vdif_save=diagy_vdif_save+diagy_vdif*dtmoy
       diagy_dudt_save=diagy_dudt_save+diagy_dudt*dtmoy

       diagy_nlt_save=diagy_nlt_save+diagy_nlt*dtmoy
       diagy_hdif_save=diagy_hdif_save+diagy_hdif*dtmoy




       diag_boz_save=diag_boz_save+diag_boz*dtmoy
     END IF

   ELSE

     nbmoy=1
     tcumul=1._rlg
     diagx_cor_save=diagx_cor
     diagx_ipg_save=diagx_ipg
     diagx_epg_save=diagx_epg
     diagx_corcor_save=diagx_corcor
     diagx_vdif_save=diagx_vdif
     diagx_dudt_save=diagx_dudt

     diagx_nlt_save=diagx_nlt
     diagx_hdif_save=diagx_hdif




     diagy_cor_save=diagy_cor
     diagy_ipg_save=diagy_ipg
     diagy_epg_save=diagy_epg
     diagy_corcor_save=diagy_corcor
     diagy_vdif_save=diagy_vdif
     diagy_dudt_save=diagy_dudt

     diagy_nlt_save=diagy_nlt
     diagy_hdif_save=diagy_hdif




     diag_boz_save=diag_boz
   END IF


   IF ( t >= tnextout) THEN

! Create output file if required
     CALL diagdyn_output_create

     IF (MASTER) THEN
     print*,' OUTPUT_DIAG_DYN : on va ecrire dans fichier:',TRIM(filediag)
     print*,' date:',tool_sectodat(t)
     print*, 'date sortie prevue',tool_sectodat(tnextout)
     print*, 'longueur moyenne (en h)',tcumul/3600.
     ENDIF

     tnextout=tnextout+dtsordiagdyn

! Estimate output fields
     diagx_cor_save=diagx_cor_save/tcumul
     diagx_ipg_save=diagx_ipg_save/tcumul
     diagx_epg_save=diagx_epg_save/tcumul
     diagx_corcor_save=diagx_corcor_save/tcumul
     diagx_vdif_save=diagx_vdif_save/tcumul
     diagx_dudt_save=diagx_dudt_save/tcumul

     diagx_nlt_save=diagx_nlt_save/tcumul
     diagx_hdif_save=diagx_hdif_save/tcumul




     diagy_cor_save=diagy_cor_save/tcumul
     diagy_ipg_save=diagy_ipg_save/tcumul
     diagy_epg_save=diagy_epg_save/tcumul
     diagy_corcor_save=diagy_corcor_save/tcumul
     diagy_vdif_save=diagy_vdif_save/tcumul
     diagy_dudt_save=diagy_dudt_save/tcumul

     diagy_nlt_save=diagy_nlt_save/tcumul
     diagy_hdif_save=diagy_hdif_save/tcumul




     diag_boz_save=diag_boz_save/tcumul

   IF ( jextmindiag <= jextmaxdiag .AND. iextmindiag <= iextmaxdiag ) THEN

     IF ( TRIM(filediag) /= TRIM(fileresu_verif(nb_out)) ) &
     CALL ionc4_write_time(filediag,0,t)

     IF (l_cor) CALL out_write_zxytvar(filediag,name_outdiagx_cor,rg_valmanq_io,mask_hx_diag, &
                                       diagx_cor_save(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)
     IF (l_cor) CALL out_write_zxytvar(filediag,name_outdiagy_cor,rg_valmanq_io,mask_hy_diag, &
                                       diagy_cor_save(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)

     IF (l_ipg) CALL out_write_zxytvar(filediag,name_outdiagx_ipg,rg_valmanq_io,mask_hx_diag, &
                                       diagx_ipg_save(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)
     IF (l_ipg) CALL out_write_zxytvar(filediag,name_outdiagy_ipg,rg_valmanq_io,mask_hy_diag, &
                                       diagy_ipg_save(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)


     IF (l_nlt) CALL out_write_zxytvar(filediag,name_outdiagx_nlt,rg_valmanq_io,mask_hx_diag, &
                                       diagx_nlt_save(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)



     IF (l_nlt) CALL out_write_zxytvar(filediag,name_outdiagy_nlt,rg_valmanq_io,mask_hy_diag, &
                                       diagy_nlt_save(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)




     IF (l_hdif) CALL out_write_zxytvar(filediag,name_outdiagx_hdif,rg_valmanq_io,mask_hx_diag, &
                                        diagx_hdif_save(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)



     IF (l_hdif) CALL out_write_zxytvar(filediag,name_outdiagy_hdif,rg_valmanq_io,mask_hy_diag, &
                                        diagy_hdif_save(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)



     IF (l_epg) CALL out_write_zxytvar(filediag,name_outdiagx_epg,rg_valmanq_io,mask_hx_diag, &
                                       diagx_epg_save(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)
     IF (l_epg) CALL out_write_zxytvar(filediag,name_outdiagy_epg,rg_valmanq_io,mask_hy_diag, &
                                       diagy_epg_save(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)

     IF (l_corcor) CALL out_write_zxytvar(filediag,name_outdiagx_corcor,rg_valmanq_io,mask_hx_diag, &
                                          diagx_corcor_save(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)
     IF (l_corcor) CALL out_write_zxytvar(filediag,name_outdiagy_corcor,rg_valmanq_io,mask_hy_diag, &
                                          diagy_corcor_save(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)

     IF (l_vdif) CALL out_write_zxytvar(filediag,name_outdiagx_vdif,rg_valmanq_io,mask_hx_diag, &
                                        diagx_vdif_save(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)
     IF (l_vdif) CALL out_write_zxytvar(filediag,name_outdiagy_vdif,rg_valmanq_io,mask_hy_diag, &
                                        diagy_vdif_save(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)

     IF (l_dudt) CALL out_write_zxytvar(filediag,name_outdiagx_dudt,rg_valmanq_io,mask_hx_diag, &
                                        diagx_dudt_save(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)
     IF (l_dudt) CALL out_write_zxytvar(filediag,name_outdiagy_dudt,rg_valmanq_io,mask_hy_diag, &
                                        diagy_dudt_save(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)



     IF (l_boz) CALL out_write_zxytvar(filediag,name_outdiag_boz,rg_valmanq_io,mask_h0_diag, &
                                       diag_boz_save(:,limin:limax,ljmin:ljmax),0,1,l_out_nc4par = l_out_nc4par)

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

   END IF ! test on jextmin...

!print*, ' avant reinitialisation'

     nbmoy=0
     diagx_cor_save=0.0_rsh
     diagx_ipg_save=0.0_rsh
     diagx_epg_save=0.0_rsh
     diagx_corcor_save=0.0_rsh
     diagx_vdif_save=0.0_rsh
     diagx_dudt_save=0.0_rsh

     diagx_nlt_save=0.0_rsh
     diagx_hdif_save=0.0_rsh




     diagy_cor_save=0.0_rsh
     diagy_ipg_save=0.0_rsh
     diagy_epg_save=0.0_rsh
     diagy_corcor_save=0.0_rsh
     diagy_vdif_save=0.0_rsh
     diagy_dudt_save=0.0_rsh

     diagy_nlt_save=0.0_rsh
     diagy_hdif_save=0.0_rsh




     diag_boz_save=0.0_rsh

   END IF

   !PRINT*, 'exit OUTPUT_DIAG_DYN routine'

  END SUBROUTINE diagdyn_output

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

  END MODULE diagdyn
