 MODULE subreadflxatm

!!======================================================================
!!                   ***  MODULE  subreadflxatm  ***
!! Atmospherical fluxes of substances :  read of fluxes
!!
!!   History :
!!        !  2013-09 B.Thouvenin from sflxread.F90 and river.F90
!!
!!======================================================================





























































!! * Modules used
   USE parameters, ONLY : rsh,rlg,lchain,imin,imax,jmin,jmax,limin,limax, &
                          ljmin,ljmax,limaxp1,ljmaxp1,pi,valmanq
   USE ionc4, ONLY : ionc4_openr,ionc4_corres,ionc4_flimax,ionc4_fljmax, &
                     ionc4_read_lat,ionc4_read_lon,ionc4_read_dimt,ionc4_read_time,&
                     ionc4_read_subxyt
   USE comvars2d,  ONLY : h0,t,nbouc,iscreenlog,ierrorlog,iwarnlog

   USE comsubstance, ONLY : nv_adv,irk_fil,name_var

   USE toolmpi,  ONLY : MASTER


   IMPLICIT NONE
   PRIVATE

!! * Accessibility
   PUBLIC sub_sflx_read                           ! routine called by main and step

!! * Shared module variables
   INTEGER,PUBLIC                                     :: idim_flxatm,ilecf,ilecfmem,cdf_id_flxatm,nvlu
   INTEGER,DIMENSION(:),PUBLIC,ALLOCATABLE            :: ivlu_flxatm
   REAL(KIND=rsh),DIMENSION(:,:,:),PUBLIC,ALLOCATABLE :: sub_flx_atm_xyt,sub_flx_atm_xy1,sub_flx_atm_xy2
   REAL(kind=rlg),PUBLIC                              :: tnc1_flxatm,tnc2_flxatm

!! * Private variables

 CONTAINS

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

  SUBROUTINE sub_sflx_read


!&E---------------------------------------------------------------------
!&E                 ***  ROUTINE sub_sflx_read  ***
!&E
!&E ** Purpose : Read files generated by the user if he wants
!&E              to take into account temporal variations
!&E              of some substance concentration in rainwater
!&E              or of some atmospheric deposition of substances
!&E
!&E
!&E ** Description : If no spatial variation : one single value for all the grid
!&E                      read the following variables in 2 different files for each substance :
!&E                        concentration of one substance in rainwater  (mass unit/m3 of water)
!&E                        atmospheric deposition (mass unit/m2/s)
!&E                            These data are provided in a ASCII file
!&E                            One file per substance whose concentration in precipitation
!&E                            is nonzero in the data file variable.dat
!&E                            And One file per substance whose atmospheric deposition
!&E                            is nonzero in the data file variable.dat
!&E
!&E **            : If spatial variation (l_subflxatm_xyt=.true.) :
!&E                      read  atmospheric fluxes in one NETCDF file, interpolated on MARS grid
!&E
!&E                        These data are subsequently interpolated onto time.
!&E
!&E ** Called by : main , step
!&E
!&E ** External calls :
!&E
!&E ** Used ij-arrays :
!&E
!&E ** Modified variables : cv_rain, sub_flx_atm
!&E
!&E ** Reference :
!&E
!&E ** History :
!&E       !  2010-05  (B. Thouvenin)
!&E
!&E---------------------------------------------------------------------
!! * Modules used
   USE comvars2d, ONLY : nvalmax,tchrono,indchrono
   USE comsubstance, ONLY : l_cvrain_readfile,cv_rain_t, &
                            cv_rain,l_subflxatm_readfile,l_subflxatm_xyt,         &
                            sub_flx_atm_t,sub_flx_atm

   IMPLICIT NONE

!! * Arguments

!! * Local declarations
   LOGICAL                    :: ex
   CHARACTER(LEN=lchain)      :: rec,file_cvrain,file_flxatm
   CHARACTER(LEN=5)           :: comment
   INTEGER                    :: eof,i,iv,kk,IERR_MPI
   REAL(KIND=rsh)             :: val1,val2
   REAL(KIND=rlg)             :: tool_datosec,tint1,tint2,dt2,dt1,t1

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

! ****  Initialization  ****
!---------------------------

IF (nbouc == 0) THEN

!  Concentrations in rain water
!  ===============================
  IF (l_cvrain_readfile) THEN
     ALLOCATE(cv_rain_t(nvalmax,nv_adv))
    
! save into simu.log
! ------------------
     IF (MASTER) THEN
     WRITE(iscreenlog,*) ' '
     WRITE(iscreenlog,*) ' '
     WRITE(iscreenlog,*) ' '
     WRITE(iscreenlog,*) '**************************************************'
     WRITE(iscreenlog,*) '********* sub_sflx_read.F90 : cv_rain************'
     WRITE(iscreenlog,*) '**************************************************'
     WRITE(iscreenlog,*) ' '
     ENDIF      
     DO iv=1,nv_adv
      IF (cv_rain(iv) == 0)then
        cv_rain_t(:,iv)=0.0_rsh
      ELSE
! reading values into an ASCII file
! only if the constant value given in variable.dat is not equal to zero
       file_cvrain='../inputs/cvrain_'//TRIM(ADJUSTL(ADJUSTR(name_var(irk_fil(iv)))))//'.dat'
       eof=0

       INQUIRE(file=file_cvrain,exist=ex)
       IF (ex) THEN
 
         IF (MASTER) THEN
         WRITE(iscreenlog,*) 'CV_rain FILE :',TRIM(file_cvrain)
         ENDIF      

         OPEN(29,file = file_cvrain,form='formatted')
         comment='debut'
         DO WHILE (comment /= '*****')
            READ(29,'(a)',iostat=eof) comment
         END DO
         i=1

! check the first available date in the cv_rain file
! ------------------------------------------------------
         READ(29,'(a)',END=999) rec
         kk=INDEX(rec,'/')-2
         READ(rec(kk+19:),*) val1
         tint1 = tool_datosec(rec(kk:kk+18))
         t1=tchrono(i)
         DO WHILE (t1 < tint1.AND.i<=nvalmax)
             cv_rain_t(i,iv)=val1
             i=i+1
             IF (i<=nvalmax) t1=tchrono(i)
         END DO

! If the simulation starts before the first data record,
! the cv_rain is set to the first available data.
         IF (i > 1) THEN
           WRITE(iwarnlog,*) ' '
           WRITE(iwarnlog,*) 'WARNING: The cv_rain data for the substance ',TRIM(name_var(irk_fil(iv)))
           WRITE(iwarnlog,*) 'is not available from the beginning of the simulation.'
           WRITE(iwarnlog,*) 'The first available value from the file : ',TRIM(file_cvrain)
           WRITE(iwarnlog,*) 'corresponds to the date : ',rec(kk:kk+18)
           WRITE(iwarnlog,*) 'It is equal to : ',val1,' mass unit/m3'
           WRITE(iwarnlog,*) 'and will be used at the beginning of the simulation.'
         END IF

! time interpolation of cv_rain data
         IF (i < nvalmax) THEN
           eof=0
           DO WHILE (eof==0)
             READ(29,'(a)',END=999) rec
             kk=INDEX(rec,'/')-2
             READ(rec(kk+19:),*) val2
             tint2 = tool_datosec(rec(kk:kk+18))
             DO WHILE ((tchrono(i)-tint1)*(tchrono(i)-tint2) <= 0.0_rlg)
               dt1=(tchrono(i)-tint1)/(tint2-tint1)
               dt2=(tint2-tchrono(i))/(tint2-tint1)
               cv_rain_t(i,iv)=val1*dt2+val2*dt1
               i=i+1
               IF (i > nvalmax) GOTO 999
             END DO 
             tint1=tint2
             val1=val2
           END DO
         END IF

 999     CONTINUE

! If the simulation ends after the last data record,
! set the cv_rain to the last available data.
         IF (i < nvalmax) THEN
            WRITE(iwarnlog,*) ' '
            WRITE(iwarnlog,*) 'WARNING: The cv_rain for the substance ',TRIM(name_var(irk_fil(iv)))
            WRITE(iwarnlog,*) 'is not available to the end of the simulation.'
            WRITE(iwarnlog,*) 'The first available value from the file : ',TRIM(file_cvrain)
            WRITE(iwarnlog,*) 'corresponds to the date : ',rec(kk:kk+18)
            WRITE(iwarnlog,*) 'It is equal to : ',val2,' mass unit/m3'
            WRITE(iwarnlog,*) 'and will be used until the end of the simulation.'
            cv_rain_t(i,iv)=val2
         END IF

         CLOSE(29)
    
       ELSE    
         WRITE(ierrorlog,*)'ERROR : the data file does not exist for cv_rain varying '
         WRITE(ierrorlog,*)'with time for the substance ',TRIM(name_var(irk_fil(iv)))
         WRITE(ierrorlog,*)'His name should be :',TRIM(file_cvrain)
         WRITE(ierrorlog,*)'Build or rename your data file'
         WRITE(ierrorlog,*)'or set to zero the concentration in rainwater for this'
         WRITE(ierrorlog,*)'substance into the data file variable.dat'
         WRITE(ierrorlog,*)
         WRITE(ierrorlog,*)'if this concentration is constant but nonzero'
         WRITE(ierrorlog,*)'you have to build a file with 2 identical values '
         WRITE(ierrorlog,*)'at the beginning and at the end of the simulation '
         CALL MPI_FINALIZE(IERR_MPI)
         STOP
       ENDIF    
      ENDIF    ! end test cv_rain cte=0
     END DO  
  ENDIF   ! end test l_cvrain_readfile

!  atmospheric fluxes
!  ===================

  IF (l_subflxatm_readfile) THEN
! variation in time and space
    IF (l_subflxatm_xyt) THEN
     CALL sub_sflx_readxyt
    ELSE
! on single value varying with time
     ALLOCATE(sub_flx_atm_t(nvalmax,nv_adv))
! save into simu.log
! ------------------

     IF (MASTER) THEN
     WRITE(iscreenlog,*) ' '
     WRITE(iscreenlog,*) ' '
     WRITE(iscreenlog,*) ' '
     WRITE(iscreenlog,*) '**************************************************'
     WRITE(iscreenlog,*) '*** sub_sflx_read.F90 : flx_atm (substance)**'
     WRITE(iscreenlog,*) '**************************************************'
     WRITE(iscreenlog,*) ' '
     ENDIF      
     
     DO iv=1,nv_adv
      IF (sub_flx_atm(iv) == 0)then
        sub_flx_atm_t(:,iv)=0.0_rlg
      ELSE
! reading values into an ASCII file
! only if the constant value given in variable.dat is not equal to zero
       file_flxatm='../inputs/flxatm_'//TRIM(ADJUSTL(ADJUSTR(name_var(irk_fil(iv)))))//'.dat'
       eof=0

       INQUIRE(file=file_flxatm,exist=ex)
       IF (ex) THEN
         
         IF (MASTER) THEN
         WRITE(iscreenlog,*) 'Flux_atm FILE :',TRIM(file_flxatm)
         ENDIF      

         OPEN(29,file = file_flxatm,form='formatted')
         comment='debut'
         DO WHILE (comment /= '*****')
           READ(29,'(a)',iostat=eof) comment
         END DO
         i=1

! check the first available date in the flxatm file
! ------------------------------------------------------
         READ(29,'(a)',END=998) rec
         kk=INDEX(rec,'/')-2
         READ(rec(kk+19:),*) val1
         tint1 = tool_datosec(rec(kk:kk+18))
         t1=tchrono(i)
         DO WHILE (t1 < tint1.AND.i<=nvalmax)
             sub_flx_atm_t(i,iv)=val1
             i=i+1
             IF (i<=nvalmax) t1=tchrono(i)
         END DO

! If the simulation starts before the first data record,
! the cv_rain is set to the first available data.
         IF (i > 1) THEN
           WRITE(iwarnlog,*) ' '
           WRITE(iwarnlog,*) 'WARNING: The sub_flx_atm data for the substance ',TRIM(name_var(irk_fil(iv)))
           WRITE(iwarnlog,*) 'is not available from the beginning of the simulation.'
           WRITE(iwarnlog,*) 'The first available value from the file : ',TRIM(file_flxatm)
           WRITE(iwarnlog,*) 'corresponds to the date : ',rec(kk:kk+18)
           WRITE(iwarnlog,*) 'It is equal to : ',val1,' mass unit/m2/s'
           WRITE(iwarnlog,*) 'and will be used at the beginning of the simulation.'
         END IF

! time interpolation of atm. fluxes
         IF (i < nvalmax) THEN
           eof=0
           DO WHILE (eof==0)
             READ(29,'(a)',END=998) rec
             kk=INDEX(rec,'/')-2
             READ(rec(kk+19:),*) val2
             tint2 = tool_datosec(rec(kk:kk+18))
             DO WHILE ((tchrono(i)-tint1)*(tchrono(i)-tint2) <= 0.0_rlg)
               dt1=(tchrono(i)-tint1)/(tint2-tint1)
               dt2=(tint2-tchrono(i))/(tint2-tint1)
               sub_flx_atm_t(i,iv)=val1*dt2+val2*dt1
               i=i+1
               IF (i > nvalmax) GOTO 998
             END DO 
             tint1=tint2
             val1=val2
           END DO
         END IF

 998     CONTINUE

! If the simulation ends after the last data record,
! set the sub_flx_atm to the last available data.
         IF (i < nvalmax) THEN
            WRITE(iwarnlog,*) ' '
            WRITE(iwarnlog,*) 'WARNING: The sub_flx_atm for the substance ',TRIM(name_var(irk_fil(iv)))
            WRITE(iwarnlog,*) 'is not available to the end of the simulation.'
            WRITE(iwarnlog,*) 'The first available value from the file : ',TRIM(file_flxatm)
            WRITE(iwarnlog,*) 'corresponds to the date : ',rec(kk:kk+18)
            WRITE(iwarnlog,*) 'It is equal to : ',val2,' mass unit/m2/s'
            WRITE(iwarnlog,*) 'and will be used until the end of the simulation.'
            sub_flx_atm_t(i,iv)=val2
         END IF

         CLOSE(29)
    
       ELSE    
         WRITE(ierrorlog,*)'ERROR : the data file does not exist for atmospheric deposition '
         WRITE(ierrorlog,*)'varying with time for the substance ',TRIM(name_var(irk_fil(iv)))
         WRITE(ierrorlog,*)'His name should be :',TRIM(file_flxatm)
         WRITE(ierrorlog,*)'Build or rename your data file'
         WRITE(ierrorlog,*)'or set to zero the atmospheric deposition for this'
         WRITE(ierrorlog,*)'substance into the data file variable.dat'
         WRITE(ierrorlog,*)
         WRITE(ierrorlog,*)'if this flux is constant but nonzero'
         WRITE(ierrorlog,*)'you have to build a file with 2 identical values '
         WRITE(ierrorlog,*)'at the beginning and at the end of the simulation '
         CALL MPI_FINALIZE(IERR_MPI)
         STOP
       ENDIF    
      ENDIF    ! end test sub_flx_atm cte=0
     END DO
   ENDIF   ! end test l_subflxatm_xyt
  ENDIF   ! end test l_subflxatm_readfile


END IF   ! ends test on nbouc

! **** End of initialization part ****
! ------------------------------------



! ****  Current progress  ****
! ----------------------------
  IF (l_cvrain_readfile) THEN
    DO iv=1,nv_adv
      dt1=(t-tchrono(indchrono))/(tchrono(indchrono+1)-tchrono(indchrono))
      cv_rain(iv)=dt1*cv_rain_t(indchrono+1,iv)+(1.d0-dt1)*cv_rain_t(indchrono,iv)
    END DO
  END IF
  IF (l_subflxatm_readfile) THEN
! variation in time and space
    IF (l_subflxatm_xyt) THEN
     IF(nbouc>0)CALL sub_sflx_readxyt
    ELSE
! on single value varying with time
      DO iv=1,nv_adv
        dt1=(t-tchrono(indchrono))/(tchrono(indchrono+1)-tchrono(indchrono))
        sub_flx_atm(iv)=dt1*sub_flx_atm_t(indchrono+1,iv)+(1.d0-dt1)*sub_flx_atm_t(indchrono,iv)
      END DO
    END IF
  END IF
  

  END SUBROUTINE sub_sflx_read

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

  SUBROUTINE sub_sflx_readxyt

!&E---------------------------------------------------------------------
!&E                 ***  ROUTINE sub_sflx_readxyt  ***
!&E
!&E ** Purpose : Read atmospherioc fluxes data for substances
!&E              These data are provided in a NetCDF file with an AO structure :
!&E              the grid is the same as MARS. interpolation has been made previously
!&E              from initial file (example : EMEP file)
!&E
!&E ** Description : Read the atmospheric fluxes for some substances
!&E                  in masse/m2/seconde (N: mmole/m2/s)
!&E
!&E ** Called by : sub_sflx_read
!&E
!&E ** External calls : tool_interp2d_sflx, ionc4_openr, ionc4_corres, ionc4_flimax
!&E                     ionc4_fljmax, ionc4_read_lat, ionc4_read_lon, ionc4_read_dimt
!&E                     ionc4_read_time, ionc4_read_subxyt
!&E
!&E ** Used ij-arrays :
!&E
!&E ** Modified variables : sub_flx_atm_xyt
!&E                         sub_flx_atm_xy1
!&E                in/out : sub_flx_atm_xy2,ilecfmem,tnc1_flxatm, tnc2_flxatm
!&E
!&E ** Reference :
!&E
!&E ** History :
!&E---------------------------------------------------------------------
!! * Modules used
   use netcdf
   USE parameters,ONLY : rsh,riosh,tref
   USE comsubstance, ONLY : file_flxatm_subs
   USE toolmpi, ONLY : MASTER,BCAST_MPI

   IMPLICIT NONE

!! * Arguments

!! * Local declarations
   LOGICAL                                     :: l_pb
   CHARACTER(LEN=19)                           :: tool_sectodat,tdate
   CHARACTER(LEN=lchain)                       :: name_flux
   INTEGER                                     :: i,j,ncerror,IERR_MPI,ivv,iv,var_id
   INTEGER                                     :: iflxatmin,iflxatmax,jflxatmin,jflxatmax
   REAL(kind=riosh),ALLOCATABLE,DIMENSION(:,:) :: tab1,tab2
   REAL(kind=rlg)                              :: torigin,tmet,dt1,dt2
   REAL(kind=rlg)                              :: tool_datosec
   INTEGER                                     :: ixmet,jxmet
!!----------------------------------------------------------------------
!! * Executable part

   !PRINT*, 'enter sub_sflx_readxyt routine, nbouc=',nbouc

!$OMP PARALLEL DEFAULT(none) PRIVATE (i,j,dt1,dt2,IERR_MPI,l_pb,tdate, &
!$OMP& torigin,ncerror,name_flux,var_id)                               &
!$OMP& SHARED(t,nbouc,tmet,iscreenlog,ierrorlog,iwarnlog,irk_fil,      &

!$OMP& MASTER,                                                         &

!$OMP& tab1,tab2,ixmet,jxmet,cdf_id_flxatm,file_flxatm_subs,           &
!$OMP& idim_flxatm,ilecf,tnc1_flxatm,tnc2_flxatm,sub_flx_atm_xyt,      &
!$OMP& sub_flx_atm_xy2,sub_flx_atm_xy1,ilecfmem,nv_adv,                &
!$OMP& iflxatmin,iflxatmax,jflxatmin,jflxatmax,name_var,               &
!$OMP& limin,limax,ljmax,ljmin,nvlu,ivlu_flxatm)

   l_pb = .false. ! No problem

   tmet=t
   torigin=tref
   torigin=0.0_rlg
! include time lag if time origine is not 01/01/1900 in atmospheric file
!torigin=tool_datosec('01/01/1900 00:00:00')-tool_datosec('01/01/0000 00:00:00')

   IF(nbouc==0) THEN

!$OMP single
     CALL ionc4_openr(file_flxatm_subs,l_out_nc4par= .true.)
     CALL ionc4_corres(file_flxatm_subs,cdf_id_flxatm)

!  le fichier des flux atmospherique est deja un fichier netdcf, interpole sur la grille de MARS
! il doit donc avoir les memes dimensions que les fichiers de la configuration
     ixmet=ionc4_flimax(cdf_id_flxatm)
     jxmet=ionc4_fljmax(cdf_id_flxatm)

! dimensions of the atmospheric file usefull for each cpu
     iflxatmin=0               ! global dimensions of the atmospheric file by default
     iflxatmax=ixmet-1
     jflxatmin=0
     jflxatmax=jxmet-1
     WRITE(*,*) 'iflxatmin =',iflxatmin,' iflxatmax =',iflxatmax,' ixmet=',ixmet
     WRITE(*,*) 'jflxatmin =',jflxatmin,' jflxatmax =',jflxatmax,' jxmet=',jxmet


     IF (MASTER) THEN
      IF (iflxatmax /= imax .OR. jflxatmax /= jmax) THEN
       WRITE(ierrorlog,*) 'The size of the domain read from the flux atmospheric field is not the same as the configuration grid'
       WRITE(ierrorlog,*) 'iflxatmin =',iflxatmin,' iflxatmax =',iflxatmax
       WRITE(ierrorlog,*) 'jflxatmin =',jflxatmin,' jflxatmax =',jflxatmax
       WRITE(ierrorlog,*) 'Simulation stopped.'
       CALL MPI_FINALIZE(IERR_MPI)
       STOP
      END IF
     ENDIF

!! quelles sont les variables a lire disponible ?
     nvlu=0
     ALLOCATE(ivlu_flxatm(nv_adv))
     DO iv=1,nv_adv
        name_flux='Depo_'//TRIM(ADJUSTL(ADJUSTR(name_var(irk_fil(iv)))))
        ncerror = nf90_inq_varid(cdf_id_flxatm,TRIM(ADJUSTL(ADJUSTR(name_flux))),var_id)
        IF(ncerror == nf90_noerr) THEN
          nvlu=nvlu+1
          ivlu_flxatm(nvlu)=iv
        ENDIF
     ENDDO

     ALLOCATE(sub_flx_atm_xyt(nv_adv,limin:limax,ljmin:ljmax))
     ALLOCATE(sub_flx_atm_xy1(nvlu,imin:imax,jmin:jmax))
     ALLOCATE(sub_flx_atm_xy2(nvlu,imin:imax,jmin:jmax))
     ALLOCATE(tab1(imin:imax,jmin:jmax))
     ALLOCATE(tab2(imin:imax,jmin:jmax))
     sub_flx_atm_xyt(1:nv_adv,limin:limax,ljmin:ljmax)=0.0_rsh
     tab1(imin:imax,jmin:jmax)=0.0_rsh
     tab2(imin:imax,jmin:jmax)=0.0_rsh

!$OMP end single

 
! Time index in the atmospheric file
! -------------------------------------

     IF (MASTER) THEN

!$OMP single
     idim_flxatm=ionc4_read_dimt(file_flxatm_subs)
     ilecf=0 
     tnc1_flxatm=0.0_rlg
     tnc2_flxatm=0.0_rlg
     DO WHILE((tmet-(tnc1_flxatm-torigin))*(tmet-(tnc2_flxatm-torigin)) > 0.0_rlg .AND. ilecf < idim_flxatm)
       ilecf=ilecf+1 
       IF (ilecf == idim_flxatm) THEN
         WRITE(ierrorlog,*) ' '
         WRITE(ierrorlog,*) 'ERROR: routine  sub_sflx_readxyt in subreadflxatm.F90'
         WRITE(ierrorlog,*) 'All records are read, dates are not compatible with simulation time.'
         WRITE(ierrorlog,*) 'Make sure your atmospheric fluxes for substances file is correct: ',TRIM(file_flxatm_subs)
         WRITE(ierrorlog,*) 'Make sure the origin time is 01/01/1900 00:00:00 in the file.'
         WRITE(ierrorlog,*) 'If not, modify variable "torigin" in the routine subreadflxatm.F90.'
         WRITE(ierrorlog,*) ''
         tdate = tool_sectodat(tmet)
         WRITE(ierrorlog,*) 't in MARS model =',tdate
         tdate = tool_sectodat(tnc1_flxatm-torigin)
         WRITE(ierrorlog,*) 't1 in atmospherical fluxes file =',tdate
         tdate = tool_sectodat(tnc2_flxatm-torigin)
         WRITE(ierrorlog,*) 't2 in atmospherical fluxes file =',tdate
         WRITE(ierrorlog,*) 'Simulation stopped.'
         WRITE(ierrorlog,*) ''
         l_pb = .true.
       END IF
       IF (.NOT. l_pb) THEN
         CALL ionc4_read_time(file_flxatm_subs,ilecf,tnc1_flxatm)
         CALL ionc4_read_time(file_flxatm_subs,ilecf+1,tnc2_flxatm)
       END IF
     ENDDO
     ilecfmem=ilecf
!$OMP end single

     ENDIF

! the MASTER sends informations to every cpus
!$OMP barrier
!$OMP master
     CALL BCAST_MPI(l_pb)
     CALL BCAST_MPI(idim_flxatm)
     CALL BCAST_MPI(ilecf)
     CALL BCAST_MPI(ilecfmem)
     CALL BCAST_MPI(tnc1_flxatm)
     CALL BCAST_MPI(tnc2_flxatm)
!$OMP end master
!$OMP barrier
!$OMP flush(l_pb,idim_flxatm,ilecf,tnc1_flxatm,tnc2_flxatm)

     IF (l_pb) THEN
       CALL MPI_FINALIZE(IERR_MPI)
       STOP
     END IF

! -3- lecture des tableaux deux indices qui encadrent le pas de temps courant

     DO ivv=1,nvlu
       iv=ivlu_flxatm(ivv)
       name_flux='Depo_'//TRIM(ADJUSTL(ADJUSTR(name_var(irk_fil(iv)))))
       write(*,*)iv,ilecf
!$OMP single
       CALL ionc4_read_subxyt(file_flxatm_subs,TRIM(ADJUSTL(ADJUSTR(name_flux))),tab1,imin,imax,jmin,jmax,ilecf,0,0)
       CALL ionc4_read_subxyt(file_flxatm_subs,TRIM(ADJUSTL(ADJUSTR(name_flux))),tab2,imin,imax,jmin,jmax,ilecf+1,0,0)
!$OMP end single
       sub_flx_atm_xy1(ivv,imin:imax,jmin:jmax)=tab1(imin:imax,jmin:jmax)
       sub_flx_atm_xy2(ivv,imin:imax,jmin:jmax)=tab2(imin:imax,jmin:jmax)
     ENDDO

!$OMP single
     DEALLOCATE(tab1)
     DEALLOCATE(tab2)
!$OMP end single

   END IF !endif for nbouc=1?

! deroulement normal

!$OMP single
   DO WHILE((tmet-(tnc1_flxatm-torigin))*(tmet-(tnc2_flxatm-torigin)) > 0.0_rlg .AND. ilecf < idim_flxatm)
     ilecf=ilecf+1 
     IF (ilecf == idim_flxatm) THEN
       WRITE(ierrorlog,*) ' '
       WRITE(ierrorlog,*) 'ERROR: routine subreadflxatm.F90'
       WRITE(ierrorlog,*) 'All records are already read.'
       WRITE(ierrorlog,*) 'The atmospheric file is not long enough: ',TRIM(file_flxatm_subs)
       WRITE(ierrorlog,*) 'Simulation stopped.'
       CALL MPI_FINALIZE(IERR_MPI)
       STOP
     END IF
     CALL ionc4_read_time(file_flxatm_subs,ilecf,tnc1_flxatm)
     CALL ionc4_read_time(file_flxatm_subs,ilecf+1,tnc2_flxatm)
   ENDDO
!$OMP end single

   IF(ilecfmem/=ilecf) THEN
     IF(nbouc/=1) THEN
!$OMP single
      ALLOCATE(  tab1(imin:imax,jmin:jmax))
      tab1(imin:imax,jmin:jmax)=0.0_rsh
!$OMP end single
     END IF

     ilecfmem=ilecf

! -1- on transfert les tableaux 2 en 1

!$OMP do schedule(static)
     DO j=ljmin,ljmax
       sub_flx_atm_xy1(:,:,j)=sub_flx_atm_xy2(:,:,j)
     ENDDO
!$OMP end do

     DO ivv=1,nvlu
       iv=ivlu_flxatm(ivv)
       name_flux='Depo_'//TRIM(ADJUSTL(ADJUSTR(name_var(irk_fil(iv)))))
!$OMP single
       CALL ionc4_read_subxyt(file_flxatm_subs,TRIM(ADJUSTL(ADJUSTR(name_flux))),tab1,imin,imax,jmin,jmax,ilecf+1,0,0)

!$OMP end single
       sub_flx_atm_xy2(ivv,imin:imax,jmin:jmax)=tab1(imin:imax,jmin:jmax)

     ENDDO

!$OMP single
     DEALLOCATE(tab1)
!$OMP end single

   END IF

   dt1=(tnc2_flxatm-torigin-tmet)/(tnc2_flxatm-tnc1_flxatm)
   dt2=1.0_rlg-dt1

!$OMP do schedule(static)
   DO ivv=1,nvlu
     iv=ivlu_flxatm(ivv)

     DO j=ljmin,ljmax
        sub_flx_atm_xyt(iv,limin:limax,j) = (sub_flx_atm_xy1(ivv,limin:limax,j)*dt1 + sub_flx_atm_xy2(ivv,limin:limax,j)*dt2)
!IF(j==100.and.ivv==1)write(*,*)dt1,dt2,sub_flx_atm_xy1(ivv,limin,j),sub_flx_atm_xy2(ivv,limin,j),sub_flx_atm_xyt(iv,limin,j)
     ENDDO

   ENDDO
!$OMP end do
!$OMP end parallel

   !PRINT*, 'exit sub_sflx_readxyt routine'

  END SUBROUTINE sub_sflx_readxyt

 END MODULE
