offline.F90

Go to the documentation of this file.

!SFX_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
!SFX_LIC This is part of the SURFEX software governed by the CeCILL-C licence
!SFX_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt  
!SFX_LIC for details. version 1.
! -------------------------------------------------
PROGRAM offline
!
! -------------------------------------------------
! Driver structure
! ----------------
! 1. Initializations
! 2. Temporal loo1s
!   2.a Read forcing
!   2.b Interpolate forcing in time
!   2.c Run surface
!   2.d Write prognostics and diagnostics variables
!
! modifications 
! 09/2012 G. Pigeon: coherence between radiation and zenith angle because of
!                    trouble with radiation received by wall in TEB
! 03/2014 E. Martin change indices names in OMP module according to GMAP changes
! 05/2014 B. Decharme delete trip
!                     I5 format to print DAY
! 04/2013 P. Lemoigne Add XDELTA_OROG to fix the maximum difference allowed between
!                     forcing and surface file orographies if LSET_FORC_ZS=.F
! 12/2013 S.Senesi    Add call to Gelato diag files init and close
! 02/2016: replace DOUBLE PRECISION by REAL to handle problem for promotion of real with GMKPACK or IBM SP
! 06/2016 S.Senesi    Use XIOS for diags output
! -------------------------------------------------
!
USE modd_off_surfex_n
!
USE modd_type_date_surf, ONLY : date
!
USE modd_write_surf_atm, ONLY : lfirst_write, ncpt_write
!
USE modd_forc_atm,  ONLY: csv         ,&! name of all scalar variables
                            xdir_alb    ,&! direct albedo for each band
                            xsca_alb    ,&! diffuse albedo for each band
                            xemis       ,&! emissivity
                            xtsrad      ,&! radiative temperature
                            xtsun       ,&! solar time                    (s from midnight)
                            xzs         ,&! orography                             (m)
                            xzref       ,&! height of T,q forcing                 (m)
                            xuref       ,&! height of wind forcing                (m)
                            xta         ,&! air temperature forcing               (K)
                            xqa         ,&! air humidity forcing                  (kg/m3)
                            xsv         ,&! scalar variables
                            xu          ,&! zonal wind                            (m/s)
                            xv          ,&! meridian wind                         (m/s)
                            xdir_sw     ,&! direct  solar radiation (on horizontal surf.)
                            xsca_sw     ,&! diffuse solar radiation (on horizontal surf.)
                            xsw_bands   ,&! mean wavelength of each shortwave band (m)
                            xzenith     ,&! zenithal angle       (radian from the vertical)
                            xzenith2    ,&! zenithal angle       (radian from the vertical)
                            xazim       ,&! azimuthal angle      (radian from North, clockwise)
                            xlw         ,&! longwave radiation (on horizontal surf.)
                            xps         ,&! pressure at atmospheric model surface (Pa)
                            xpa         ,&! pressure at forcing level             (Pa)
                            xrhoa       ,&! density at forcing level              (kg/m3)
                            xco2        ,&! CO2 concentration in the air          (kg/m3)
                            xsnow       ,&! snow precipitation                    (kg/m2/s)
                            xrain       ,&! liquid precipitation                  (kg/m2/s)
                            xsfth       ,&! flux of heat                          (W/m2)
                            xsftq       ,&! flux of water vapor                   (kg/m2/s)
                            xsfu        ,&! zonal momentum flux                   (m/s)
                            xsfv        ,&! meridian momentum flux                (m/s)
                            xsfco2      ,&! flux of CO2                           (kg/m2/s)
                            xsfts       ,&! flux of scalar var.                   (kg/m2/s)
                            xpew_a_coef ,&! implicit coefficients
                            xpew_b_coef ,&! needed if HCOUPLING='I'
                            xpet_a_coef ,&
                            xpeq_a_coef ,&
                            xpet_b_coef ,&
                            xpeq_b_coef ,&
                            xtsurf      ,&! effective temperature                  (K)
                            xz0         ,&! surface roughness length for momentum  (m)
                            xz0h        ,&! surface roughness length for heat      (m)
                            xqsurf        ! specific humidity at surface           (kg/kg)
!
USE modd_surf_conf,  ONLY : cprogname, csoftware
USE modd_csts,       ONLY : xpi, xday, xrv, xrd, xg
USE modd_io_surf_asc,ONLY : cfilein,cfilein_save,cfileout,cfilepgd
USE modd_surf_par
USE modd_io_surf_fa, ONLY : cfilein_fa, cfilein_fa_save,       &
                            cfileout_fa, nunit_fa, cdnomc,     &
                            iverbfa, lfanocompact, cfilepgd_fa  
USE modd_io_surf_lfi,ONLY : cfilein_lfi, cfilein_lfi_save, cluout_lfi, cfileout_lfi, &
                            lmnh_compatible, cfilepgd_lfi  
USE modd_io_surf_nc, ONLY : cfilein_nc, cfilein_nc_save, cfileout_nc, cluout_nc, &
                            cfilepgd_nc, ldef_nc=>ldef
USE modd_io_surf_ol, ONLY : xstart, xcount, xstride, lpartw,    &
                              xstartw, xcountw, ltime_written,  &
                              nstep_output, ldef_ol=>ldef
USE modd_write_bin,  ONLY : nwrite
!
USE modd_surfex_mpi, ONLY : ncomm, nproc, nrank, npio, wlog_mpi, prep_log_mpi,   &
                            nindex, nsize_task, xtime_npio_read, xtime_npio_write, &
                            xtime_comm_read, xtime_comm_write, xtime_sea, &
                            xtime_nature, xtime_water, xtime_town, &
                            xtime_init_sea, xtime_init_water, &
                            xtime_init_nature, xtime_init_town, &
                            xtime_write, xtime_calc,  &
                            xtime_calc_write, idx_w, end_log_mpi
!
USE modd_surfex_omp, ONLY :  nblock, nblocktot
!
USE modd_coupling_topd, ONLY : nnb_topd, nnb_stp_restart, lbudget_topd, ltopd_step, &
                               lcoupl_topd, ntopd_step, nyear, nmonth, nday, nh, nm
USE modd_topodyn, ONLY : xtopd_step, nnb_topd_step, xqtot, xqb_run, xqb_dr
!
USE modd_slope_effect, ONLY: xzs_thread,xzs_xy_thread,xslopang_thread,&
                             xslopazi_thread,xsurf_triangle_thread
!
USE modd_sfx_oasis, ONLY : loasis, xruntime
!
USE modd_xios, ONLY : lxios, txios_context, lxios_def_closed, ladd_dim=>lallow_add_dim, ntimestep
!
USE mode_pos_surf
!
USE mode_crodebug
!
USE modn_io_offline
!
USE modi_get_luout
USE modi_open_namelist
USE modi_test_nam_var_surf
USE modi_close_namelist
USE modi_read_all_namelists
USE modi_open_close_bin_asc_forc
USE modi_open_filein_ol
USE modi_ol_read_atm_conf
USE modi_abor1_sfx
USE modi_ol_alloc_atm
USE modi_compare_orography
USE modi_sunpos
USE modi_init_index_mpi
USE modi_ol_read_atm
USE modi_io_buff_clean
USE modi_init_surf_atm_n
USE modi_init_surf_landuse_n
USE modi_ol_time_interp_atm
USE modi_coupling_surf_atm_n
USE modi_add_forecast_to_date_surf
USE modi_write_surf_atm_n
USE modi_write_header_mnh
USE modi_flag_update
USE modi_flag_diag_update
USE modi_diag_surf_atm_n
USE modi_write_diag_surf_atm_n
USE modi_get_surf_var_n
USE modi_gather_and_write_mpi
USE modi_close_filein_ol
USE modi_close_fileout_ol
USE modi_init_output_ol_n
USE modi_init_output_nc_n
!
USE modi_write_header_fa
USE modi_abor1_sfx
!
USE modi_write_discharge_file
USE modi_write_budget_coupl_rout
USE modi_prep_restart_coupl_topd
!
USE modi_init_slope_param
USE modi_slope_radiative_effect
!
USE modi_sfx_oasis_read_nam
USE modi_sfx_oasis_init
USE modi_sfx_oasis_def_ol
USE modi_sfx_oasis_recv_ol
USE modi_sfx_oasis_send_ol
USE modi_sfx_oasis_end
!RJ: missing modi
USE modi_local_slope_param
!
#ifdef WXIOS
USE xios, ONLY : xios_context_finalize, xios_close_context_definition, xios_update_calendar
#endif
USE modi_sfx_xios_readnam_ol
USE modi_sfx_xios_setup_ol
!
USE mode_glt_dia_lu
!
#ifdef SFX_MPI
#ifdef SFX_MPL
USE mpl_data_module, ONLY : lmplusercomm, mplusercomm
#endif
#endif
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!
#ifdef AIX64
!$ USE OMP_LIB
#endif
!
IMPLICIT NONE
!
#ifdef SFX_MPI
include 'mpif.h'
#endif
!
#ifndef AIX64
!$ INCLUDE 'omp_lib.h'
#endif
!
!*      0.    declarations of local variables
!
 CHARACTER(LEN=3), PARAMETER       :: YINIT     = 'ALL'
!
 CHARACTER(LEN=28)                 :: YLUOUT    = 'LISTING_OFFLINE             '
!
INTEGER                           :: IYEAR               ! current year (UTC)
INTEGER                           :: IMONTH              ! current month (UTC)
INTEGER                           :: IDAY                ! current day (UTC)
INTEGER                           :: IYEAR2              ! current year at end of timestep(UTC)
INTEGER                           :: IMONTH2             ! current month at end of timestep(UTC)
INTEGER                           :: IDAY2               ! current day at end of timestep(UTC)
REAL                              :: ZTIME               ! current time since start of the day (s)
REAL                              :: ZTIME2              ! current time since start of the day at end of timestep (s)
REAL                              :: ZTIMEC              ! current duration since start of the run (s)
!
INTEGER                           :: IYEAR_OUT           ! output year name
INTEGER                           :: IMONTH_OUT          ! output month name
INTEGER                           :: IDAY_OUT            ! output day name
REAL                              :: ZTIME_OUT           ! output time since start of the run (s)
!
INTEGER, DIMENSION(11)  :: IDATEF
!
 CHARACTER(LEN=28), PARAMETER      :: YATMFILE     = '                            '
 CHARACTER(LEN=6),  PARAMETER      :: YATMFILETYPE = '      '
 CHARACTER(LEN=2),  PARAMETER      :: YTEST        = 'OK'          ! must be equal to 'OK'
!
REAL, DIMENSION(:), POINTER       :: ZLAT                ! latitude                         (rad)
REAL, DIMENSION(:), POINTER       :: ZLON                ! longitude                        (rad)
REAL, DIMENSION(:), POINTER       :: ZZS_FORC            ! orography                        (m)  
REAL, DIMENSION(:), POINTER       :: ZZREF               ! Forcing level for T
REAL, DIMENSION(:), POINTER       :: ZUREF               ! Forcing level for U
!
REAL                              :: ZTSTEP              ! atmospheric time-step            (s)
!
INTEGER                           :: INI                 ! grid dimension
INTEGER                           :: JLOOP               ! loop counter
INTEGER                           :: IBANDS              ! Number of radiative bands 
INTEGER                           :: INB_STEP_ATM        ! Number of atmospheric time-steps
INTEGER                           :: INB_ATM             ! Number of Isba time-steps 
                                                         ! within a forcing time-step
INTEGER                           :: ID_FORC             ! indice of forcing in the file
INTEGER                           :: INB_LINES           ! nb of lines to read in the forcing file
INTEGER                           :: IDMAX               ! nb of lines to read in the forcing file at last 
INTEGER                           :: JFORC_STEP          ! atmospheric loop index
INTEGER                           :: JSURF_STEP          ! isba loop index
INTEGER                           :: ICOUNT              ! day counter 
INTEGER                           :: ITIMESTARTINDEX
REAL                              :: ZDURATION, ZDURATION2    ! duration of run                     (s)
!
REAL, DIMENSION(:,:), ALLOCATABLE :: ZTA                 ! air temperature forcing               (K)
REAL, DIMENSION(:,:), ALLOCATABLE :: ZQA                 ! air humidity forcing                  (kg/m3)
REAL, DIMENSION(:,:), ALLOCATABLE :: ZWIND               ! wind speed                            (m/s)
REAL, DIMENSION(:,:), ALLOCATABLE :: ZSCA_SW             ! diffuse solar radiation (on horizontal surf.)
REAL, DIMENSION(:,:), ALLOCATABLE :: ZDIR_SW             ! direct  solar radiation (on horizontal surf.)
REAL, DIMENSION(:,:), ALLOCATABLE :: ZLW                 ! longwave radiation (on horizontal surf.)
REAL, DIMENSION(:,:), ALLOCATABLE :: ZSNOW               ! snow precipitation                    (kg/m2/s)
REAL, DIMENSION(:,:), ALLOCATABLE :: ZRAIN               ! liquid precipitation                  (kg/m2/s)
REAL, DIMENSION(:,:), ALLOCATABLE :: ZPS                 ! pressure at forcing level             (Pa)
REAL, DIMENSION(:,:), ALLOCATABLE :: ZCO2                ! CO2 concentration in the air          (kg/m3)
REAL, DIMENSION(:,:), ALLOCATABLE :: ZDIR                ! wind direction
INTEGER                           :: ILUOUT              ! ascii output unit number
INTEGER                           :: ILUNAM              ! namelist unit number
INTEGER                           :: IRET                ! error return code
INTEGER                           :: INB 
INTEGER                           :: INW, JNW
 CHARACTER(LEN=14)                :: YTAG      
LOGICAL                           :: GFOUND              ! return logical when reading namelist
LOGICAL                           :: GSHADOWS    
REAL, DIMENSION(:),   ALLOCATABLE :: ZSW                 ! total solar radiation (on horizontal surf.)
REAL, DIMENSION(:),   ALLOCATABLE :: ZCOEF               ! coefficient for solar radiation interpolation near sunset/sunrise
!
TYPE(date) :: TDATE_END
! Flag diag :
!
INTEGER                           :: I2M, IBEQ, IDSTEQ
LOGICAL                           :: GFRAC, GDIAG_GRID, GSURF_BUDGET, GRAD_BUDGET, GCOEF,    &
                                     GSURF_VARS, GDIAG_OCEAN, GDIAG_SEAICE, GWATER_PROFILE,  &
                                     GINTERPOL_TS, GSURF_EVAP_BUDGET, GFLOOD, GPGD_ISBA,     &
                                     GCH_NO_FLUX_ISBA, GSURF_MISC_BUDGET_ISBA, GPGD_TEB,     &
                                     GSURF_MISC_BUDGET_TEB
!
! Inquiry mode arrays:
!
REAL, DIMENSION(:), ALLOCATABLE   :: ZSEA, ZWATER, ZNATURE, ZTOWN
REAL, DIMENSION(:), ALLOCATABLE   :: ZSEA_FULL, ZWATER_FULL, ZNATURE_FULL, ZTOWN_FULL
REAL, DIMENSION(:), ALLOCATABLE   :: ZT2M, ZQ2M
REAL, DIMENSION(:), ALLOCATABLE   :: ZZ0, ZZ0H, ZQS
REAL, DIMENSION(:), ALLOCATABLE   :: ZQS_SEA, ZQS_WATER, ZQS_NATURE, ZQS_TOWN
REAL, DIMENSION(:), ALLOCATABLE   :: ZPSNG, ZPSNV
REAL, DIMENSION(:), ALLOCATABLE   :: ZZ0EFF
REAL, DIMENSION(:), ALLOCATABLE   :: ZZS
REAL, DIMENSION(:), ALLOCATABLE   :: ZZ0_FULL, ZZ0EFF_FULL, ZZS_FULL
REAL, DIMENSION(:), ALLOCATABLE :: ZSUMZEN
INTEGER :: ISERIES, ISIZE
!
! MPI variables
!
 CHARACTER(LEN=100) :: YNAME
 CHARACTER(LEN=10)  :: YRANK
INTEGER :: ILEVEL, INFOMPI, J
REAL :: XTIME0, XTIME1, XTIME
!
! SFX - OASIS coupling variables
!
LOGICAL :: GSAVHOOK
INTEGER :: IBLOCKTOT, IBLOCK
!
REAL(KIND=JPRB) :: ZHOOK_HANDLE
!
! --------------------------------------------------------------------------------------
!
!*     0.1.   MPI, OASIS, XIOS and dr_hook initializations
!
csoftware='OFFLINE'
!
infompi=1
!
gsavhook = lhook
lhook = .false.
! There are issues with oasis if LHOOK=T during its init phase 
#ifdef WXIOS
 CALL sfx_xios_readnam_ol(cnamelist)
#else
lxios=.false.
#endif
lhook = gsavhook
!
!Must be call before DRHOOK !
 CALL sfx_oasis_init(cnamelist,ncomm)
#ifdef SFX_MPI
#ifdef SFX_MPL
IF (loasis.OR.lxios) THEN
  lmplusercomm = .true.
  mplusercomm = ncomm
ENDIF
#endif
#endif
!
#ifdef SFX_MPI
IF(.NOT.loasis.AND..NOT.lxios)THEN
 CALL mpi_init_thread(mpi_thread_multiple,ilevel,infompi)
 IF (infompi /= mpi_success) THEN 
    CALL abor1_sfx('OFFLINE: ERROR WHEN INITIALIZING MPI')
 ENDIF
 ncomm=mpi_comm_world
ENDIF
 CALL mpi_comm_size(ncomm,nproc,infompi)
 CALL mpi_comm_rank(ncomm,nrank,infompi)
#endif
!
IF (lhook) CALL dr_hook('OFFLINE',0,zhook_handle)
!
!RJ: init modd_surefx_omp
!$OMP PARALLEL
!$ NBLOCKTOT = OMP_GET_NUM_THREADS()
!$ NBLOCK = OMP_GET_THREAD_NUM()
!$OMP END PARALLEL
!
iblocktot = 1
iblock = 0
!
 CALL prep_log_mpi
!
 CALL wlog_mpi(' ')
!
 CALL wlog_mpi('NBLOCKTOT ',klog=nblocktot)
!
#ifdef SFX_MPI
xtime0 = mpi_wtime()
#endif
!
!
!*      0.3.   Open ascii file for writing
!
WRITE(yrank,fmt='(I10)') nrank
yname=trim(yluout)//adjustl(yrank)
!
 cluout_lfi =  adjustl(adjustr(yname)//'.txt')
 cluout_nc  =  adjustl(adjustr(yname)//'.txt')
!
 CALL get_luout('ASCII ',iluout)
OPEN(unit=iluout,file=adjustl(adjustr(yname)//'.txt'),form='FORMATTED',action='WRITE')
!
!
IF ( nrank==npio ) THEN
  !
!RJ: be verbose just for openmp
  IF(nblocktot==1) THEN
!$  WRITE(*,*) "CAUTION: DID YOU THINK TO SET OMP_NUM_THREADS=1?"
!$  WRITE(*,*) "PLEASE VERIFY OMP_NUM_THREADS IS INITIALIZED : TYPE ECHO $OMP_NUM_THREADS IN A TERMINAL" 
  !
!$  WRITE(ILUOUT,*) "CAUTION: DID YOU THINK TO SET OMP_NUM_THREADS=1?"
!$  WRITE(ILUOUT,*) "PLEASE VERIFY OMP_NUM_THREADS IS INITIALIZED : TYPE ECHO $OMP_NUM_THREADS IN A TERMINAL" 
  ENDIF
  !
ENDIF
!
!*      0.4.   Reads namelists
!
 CALL open_namelist('ASCII ',ilunam,cnamelist)
!
 CALL posnam(ilunam,'NAM_IO_OFFLINE',gfound,iluout)
IF (gfound) READ (unit=ilunam,nml=nam_io_offline)
 CALL close_namelist('ASCII ',ilunam)
!
IF (nproc==1) THEN 
  xio_frac=1.
ELSE
  xio_frac = max(min(xio_frac,1.),0.)
ENDIF
!
 CALL test_nam_var_surf(iluout,'CSURF_FILETYPE',csurf_filetype,'ASCII ','LFI   ','FA    ','NC    ')
#ifdef WXIOS
 CALL test_nam_var_surf(iluout,'CTIMESERIES_FILETYPE',ctimeseries_filetype,'NETCDF','TEXTE ','BINARY',&
                                                                            'ASCII ','LFI   ','FA    ',&
                                                                            'NONE  ','OFFLIN','NC    '&
                                                                            ,'XIOS  ')
#else 
 CALL test_nam_var_surf(iluout,'CTIMESERIES_FILETYPE',ctimeseries_filetype,'NETCDF','TEXTE ','BINARY',&
                                                                            'ASCII ','LFI   ','FA    ',&
                                                                            'NONE  ','OFFLIN','NC    ')  
#endif                                                                    
 CALL test_nam_var_surf(iluout,'CFORCING_FILETYPE',cforcing_filetype,'NETCDF','ASCII ','BINARY')
!
IF (nscal>59) CALL abor1_sfx("OFFLINE: NSCAL MUST BE LOWER THAN OR EQUAL TO 59")
!
!
IF (ctimeseries_filetype=='NETCDF') ctimeseries_filetype='OFFLIN'
!
IF ((trim(ctimeseries_filetype) /= 'XIOS') .AND. ladd_dim) THEN
      CALL abor1_sfx('CANNOT YET SET LALLOW_ADD_DIM TO .TRUE. WITHOUT SETTING CTIMESERIES_FILETYPE to XIOS ')
ENDIF
!
!
 cfilepgd = adjustl(adjustr(cpgdfile)//'.txt')
 cfilein  = adjustl(adjustr(cprepfile)//'.txt')
 cfilein_save = cfilein
!
 cfilepgd_lfi = cpgdfile
 cfilein_lfi  = cprepfile
 cfilein_lfi_save = cfilein_lfi
!
 cfilepgd_fa = adjustl(adjustr(cpgdfile)//'.fa')
 cfilein_fa  = adjustl(adjustr(cprepfile)//'.fa')
 cfilein_fa_save  = cfilein_fa
!
 cfilepgd_nc = adjustl(adjustr(cpgdfile)//'.nc')
 cfilein_nc  = adjustl(adjustr(cprepfile)//'.nc')
 cfilein_nc_save  = cfilein_nc
!
!     Allocations of Surfex Types
 CALL surfex_alloc_list(1)
 ysc => ysurf_list(1)
!
!     Reading all namelist (also assimilation)
 CALL read_all_namelists(ysc, csurf_filetype,'ALL',.false.)
!
!
!*      0.5.   Reads SFX - OASIS coupling namelists
!
 CALL sfx_oasis_read_nam(csurf_filetype,xtstep_surf)
!
!*      0.6   Assume FA filetype consistency 
!
 cprogname = csurf_filetype
!
! --------------------------------------------------------------------------------------
!
!*      1.    Initializations
!
!       netcdf file handling
!
IF (nrank==npio) THEN
  !
  xstart            = nundef
  xstride           = nundef
  xcount            = nundef
  xstartw           = 0
  xcountw           = 1
  lpartw            = .true.
  !
ENDIF
!
#ifdef SFX_MPI
xtime = (mpi_wtime() - xtime0)
#endif
 CALL wlog_mpi('READ NAMELISTS ',plog=xtime)
#ifdef SFX_MPI
xtime0 = mpi_wtime()
#endif
!
!
!       splitting of the grid
!
gshadows = lshadows_slope .OR. lshadows_other
 CALL init_index_mpi(ysc%DTCO, ysc%U, ysc%UG, ysc%GCP, csurf_filetype, 'OFF', yalg_mpi, xio_frac, gshadows)
!
 CALL wlog_mpi(' ')
 CALL wlog_mpi('TIME_NPIO_READ init_index ',plog=xtime_npio_read)
 CALL wlog_mpi('TIME_COMM_READ init_index ',plog=xtime_comm_read)
xtime_npio_read = 0.
xtime_comm_read = 0.
!
#ifdef SFX_MPI
xtime = (mpi_wtime() - xtime0)
#endif
 CALL wlog_mpi(' ')
 CALL wlog_mpi('INIT_INDEX_MPI ',plog=xtime)
 CALL wlog_mpi(' ')
#ifdef SFX_MPI
xtime0 = mpi_wtime()
#endif
!
!       forcing file handling
!
IF (cforcing_filetype=='ASCII ' .OR. cforcing_filetype=='BINARY') CALL open_close_bin_asc_forc('CONF ',cforcing_filetype,'R')
IF (cforcing_filetype=='NETCDF') CALL open_filein_ol
!
!       configuration of run
!
 CALL ol_read_atm_conf(ysc%DTCO, ysc%U, ysc%UG%G%CGRID, csurf_filetype, cforcing_filetype,  &
                      ldelayedstart_nc, ndatestop, zduration, ztstep, ini,  &
                      iyear, imonth, iday, ztime, zlat, zlon, zzs_forc,     &
                      zzref, zuref, itimestartindex     )
!
tdate_end%YEAR = iyear
tdate_end%MONTH = imonth
tdate_end%DAY = iday
zduration2 = zduration
 CALL add_forecast_to_date_surf(tdate_end%YEAR, tdate_end%MONTH, tdate_end%DAY, zduration2)
!
 CALL wlog_mpi(' ')
 CALL wlog_mpi('TIME_NPIO_READ forc conf ',plog=xtime_npio_read)
 CALL wlog_mpi('TIME_COMM_READ forc conf ',plog=xtime_comm_read)
xtime_npio_read = 0.
xtime_comm_read = 0.
!
#ifdef SFX_MPI
xtime = (mpi_wtime() - xtime0)
#endif
 CALL wlog_mpi('OL_READ_ATM_CONF ',plog=xtime)
 CALL wlog_mpi(' ')
#ifdef SFX_MPI
xtime0 = mpi_wtime()
#endif
!
!*     time steps coherence check 
!
IF ( (mod(xtstep_output,ztstep)*mod(ztstep,xtstep_output) /= 0) .OR. (mod(ztstep,xtstep_surf) /= 0) ) THEN
   WRITE(iluout,*)' FORCING  AND OUTPUT/SURFACE TIME STEP SHOULD BE MULTIPLE', &
     nint(ztstep),nint(xtstep_output),nint(xtstep_surf)    
   CALL abor1_sfx('OFFLINE: FORCING  AND OUTPUT/SURFACE TIME STEP SHOULD BE MULTIPLE')
ENDIF
!
IF ( ztime /= 0. .AND. mod(ztime,xtstep_surf) /= 0  ) THEN
   WRITE(iluout,*)' INITIAL AND SURFACE TIME STEP SHOULD BE MULTIPLE', &
   nint(ztime),nint(xtstep_surf)  
   CALL abor1_sfx('OFFLINE: INITIAL AND SURFACE TIME STEP SHOULD BE MULTIPLE')
ENDIF
!
IF(loasis.AND.zduration/=xruntime)THEN
   WRITE(iluout,*)'Total simulated time given by Forcing field and OASIS namcouple are different'
   WRITE(iluout,*)'From Forcing (s) : ',zduration, 'From OASIS   (s) : ',xruntime
   CALL abor1_sfx('OFFLINE: TOTAL SIMULATED TIME DIFFERENT BETWEEN FORCING AND OASIS')
ENDIF
!
inb_step_atm  = int(zduration / ztstep)
inb_atm       = int(ztstep / xtstep_surf)
nstep_output  = int(zduration / xtstep_output)
!
xtopd_step = 0
nnb_topd_step = 0
ntopd_step = 0
IF ( lcoupl_topd ) THEN
  !
  xtopd_step = float(nnb_topd)* xtstep_surf
  nnb_topd_step = int( zduration / xtopd_step )
  !
  IF ( nnb_stp_restart==0 .AND. .NOT.lrestart ) nnb_stp_restart = -1
  !
  ntopd_step = 1
  !
ENDIF
!
!       allocation of variables
!
ibands = 1
!
 CALL ol_alloc_atm(ini,ibands,nscal)
!
xzs   = zzs_forc
xzref = zzref
xuref = zuref
!
!       compare orography
!
 CALL compare_orography(ysc%DTCO, ysc%U, csurf_filetype, lset_forc_zs, xdelta_orog)
!
!       miscellaneous initialization
!
icount = 0
ztimec = 0.
!
 CALL sunpos(iyear, imonth, iday, ztime, zlon, zlat, xtsun, xzenith, xazim)
!
!number of lines read in forcing files
inb_lines=1
IF (nb_read_forc.EQ.1) THEN
  inb_lines=inb_step_atm
ELSEIF (nb_read_forc.NE.0) THEN
  !to be sure the number of readings will be NB_READ_FORC as a maximum
  inb_lines=ceiling(1.*(inb_step_atm+1)/nb_read_forc)
ENDIF
!number of lines to be read effectively
idmax=inb_lines+1
!effective number of readings of the forcing files
nb_read_forc=ceiling(1.*(inb_step_atm+1)/inb_lines)
!
!     open Gelato specific diagnostic files (if requested by 
!     Gelato wizzard user)
!
#if ! defined in_arpege
 CALL opndia()
#endif
!
!       allocate local atmospheric variables
!
IF (.NOT.ALLOCATED(zta))    ALLOCATE(zta(ini,inb_lines+1)) 
IF (.NOT.ALLOCATED(zqa))    ALLOCATE(zqa(ini,inb_lines+1))
IF (.NOT.ALLOCATED(zwind))  ALLOCATE(zwind(ini,inb_lines+1))
IF (.NOT.ALLOCATED(zdir_sw))ALLOCATE(zdir_sw(ini,inb_lines+1))
IF (.NOT.ALLOCATED(zsca_sw))ALLOCATE(zsca_sw(ini,inb_lines+1))
IF (.NOT.ALLOCATED(zlw))    ALLOCATE(zlw(ini,inb_lines+1))
IF (.NOT.ALLOCATED(zsnow))  ALLOCATE(zsnow(ini,inb_lines+1))
IF (.NOT.ALLOCATED(zrain))  ALLOCATE(zrain(ini,inb_lines+1))
IF (.NOT.ALLOCATED(zps))    ALLOCATE(zps(ini,inb_lines+1))
IF (.NOT.ALLOCATED(zco2))   ALLOCATE(zco2(ini,inb_lines+1))
IF (.NOT.ALLOCATED(zdir))   ALLOCATE(zdir(ini,inb_lines+1))
IF (.NOT.ALLOCATED(zcoef))  ALLOCATE(zcoef(ini))
IF (.NOT.ALLOCATED(zsumzen))ALLOCATE(zsumzen(ini))
!
IF (.NOT.ALLOCATED(zsw))ALLOCATE(zsw(ini))
!
!      computes initial air co2 concentration and  density
!
#ifdef SFX_MPI
xtime = (mpi_wtime() - xtime0)
#endif
 CALL wlog_mpi('COMPARE_OROGRAPHY SUNPOS ',plog=xtime)
#ifdef SFX_MPI
xtime0 = mpi_wtime()
#endif
!
!* opens forcing files (if ASCII or BINARY)
!
IF (cforcing_filetype=='ASCII ' .OR. cforcing_filetype=='BINARY') &
        CALL open_close_bin_asc_forc('OPEN ',cforcing_filetype,'R')
!
 CALL ol_read_atm(csurf_filetype, cforcing_filetype, itimestartindex,&
                  zta,zqa,zwind,zdir_sw,zsca_sw,zlw,zsnow,zrain,zps,&
                  zco2,zdir,llimit_qair                           ) 
!
 CALL wlog_mpi(' ')
 CALL wlog_mpi('TIME_NPIO_READ forc ',plog=xtime_npio_read)
 CALL wlog_mpi('TIME_COMM_READ forc ',plog=xtime_comm_read)
xtime_npio_read = 0.
xtime_comm_read = 0.
!
#ifdef SFX_MPI
xtime = (mpi_wtime() - xtime0)
#endif
 CALL wlog_mpi(' ')
 CALL wlog_mpi('OL_READ_ATM0 ',plog=xtime)
 CALL wlog_mpi(' ')
#ifdef SFX_MPI
xtime0 = mpi_wtime()
#endif
!
xco2(:)  = zco2(:,1)
xrhoa(:) = zps(:,1) / (xrd * zta(:,1) * ( 1.+((xrv/xrd)-1.)*zqa(:,1) ) + xg * xzref )
!                 
!       surface Initialisation     
!
#ifdef SFX_MPI
xtime = (mpi_wtime() - xtime0)
#endif
 CALL wlog_mpi('CO2 RHOA ',plog=xtime)
!
 CALL io_buff_clean
!
 !CALL SURFEX_DEALLO_LIST
 !CALL SURFEX_ALLOC_LIST(IBLOCKTOT)
!
#ifdef SFX_MPI
xtime0 = mpi_wtime()
#endif
!
 CALL goto_model(1)
!
 CALL init_surf_atm_n(ysc, csurf_filetype, yinit, lland_use, ini, nscal, ibands,  &
                      csv,xco2(:),xrhoa(:),  xzenith(:),xazim(:),xsw_bands,  &
                      xdir_alb(:,:), xsca_alb(:,:), xemis(:), xtsrad(:),     &
                      xtsurf(:), iyear, imonth, iday, ztime, tdate_end,      &
                      yatmfile, yatmfiletype, ytest                          )
!
! initialization routines to compute shadows
IF (gshadows) THEN
  IF (iblock==0) THEN
    CALL init_slope_param(ysc%UG%G, ysc%UG%XGRID_FULL_PAR,zzs_forc,ini,zlat)
  END IF 
  CALL local_slope_param(1,ini)
END IF
!
#ifdef SFX_MPI
xtime = (mpi_wtime() - xtime0)
#endif
 CALL wlog_mpi(' ')
 CALL wlog_mpi('INIT_SURF_ATM ',plog=xtime)
 CALL wlog_mpi(' ')
!
 CALL wlog_mpi('TIME_NPIO_READ init ',plog=xtime_npio_read)
 CALL wlog_mpi('TIME_COMM_READ init ',plog=xtime_comm_read)
 CALL wlog_mpi(' ')
!
xtime_npio_read = 0.
xtime_comm_read = 0.
!
!   Land use or/and vegetation dynamic
!    
 CALL init_surf_landuse_n(ysc%DTCO, ysc%DUO%LREAD_BUDGETC, ysc%U, ysc%UG,   &
                         ysc%IM, ysc%SV, ysc%SLT, ysc%NDST, csurf_filetype,&
                         yinit, lland_use, ini, nscal, ibands, csv,        &
                         xco2, xrhoa, xzenith, xazim, xsw_bands, xdir_alb, &
                         xsca_alb, xemis, xtsrad, xtsurf, iyear, imonth,   &
                         iday, ztime, yatmfile, yatmfiletype, ytest )
!
#ifdef SFX_MPI
xtime0 = mpi_wtime()
#endif
!
 CALL init_crodebug(ysc%IM%NPE%AL(1)%TSNOW%SCHEME)
!
! * SURFEX - OASIS  grid, partitions and local field definitions
!
IF(loasis)THEN
  CALL sfx_oasis_def_ol(ysc%IM%O, ysc%U, csurf_filetype,yalg_mpi)
ENDIF
!
! --------------------------------------------------------------------------------------
!
 CALL sfx_xios_setup_ol(ysc,iluout,iyear,imonth,iday,ztime,xtstep_output)
!
nwrite = 0
!
#ifdef SFX_MPI
xtime = (mpi_wtime() - xtime0)
#endif
 CALL wlog_mpi('INIT FINISHED ',plog=xtime)
#ifdef SFX_MPI
xtime0 = mpi_wtime()
#endif
!*      2.    Temporal loops
!
xtime_calc(:) = 0.
xtime_write(:) = 0.
!
lfirst_write = .true.
!
ncpt_write = 0
!
DO jforc_step=1,inb_step_atm
  !
#ifdef SFX_MPI
  xtime1 = mpi_wtime()
#endif
  ! read Forcing
  !
  !indice of forcing line in forcing arrays
  id_forc=jforc_step-int(jforc_step/inb_lines)*inb_lines
  IF (id_forc==0) id_forc=inb_lines
  !new forcings to read
  IF (id_forc==1 .AND. jforc_step.NE.1) THEN
    !if last part of forcing, the last point has to be adjusted on the end of
    !files
    IF (jforc_step/inb_lines==nb_read_forc-1) THEN 
      idmax=inb_step_atm-jforc_step+1+1
      !for ascii and binary forcing files
      zta(:,idmax) = zta(:,SIZE(zta,2))
      zqa(:,idmax) = zqa(:,SIZE(zta,2))
      zwind(:,idmax) = zwind(:,SIZE(zta,2))
      zdir_sw(:,idmax) = zdir_sw(:,SIZE(zta,2))
      zsca_sw(:,idmax) = zsca_sw(:,SIZE(zta,2))
      zlw(:,idmax) = zlw(:,SIZE(zta,2))
      zsnow(:,idmax) = zsnow(:,SIZE(zta,2))
      zrain(:,idmax) = zrain(:,SIZE(zta,2))
      zps(:,idmax) = zps(:,SIZE(zta,2))
      zco2(:,idmax) = zco2(:,SIZE(zta,2))
      zdir(:,idmax) = zdir(:,SIZE(zta,2))
    ENDIF
    CALL ol_read_atm(csurf_filetype, cforcing_filetype, itimestartindex+jforc_step-1, &
                     zta(:,1:idmax),zqa(:,1:idmax),zwind(:,1:idmax),       &
                     zdir_sw(:,1:idmax),zsca_sw(:,1:idmax),zlw(:,1:idmax), &
                     zsnow(:,1:idmax),zrain(:,1:idmax),zps(:,1:idmax),     &
                     zco2(:,1:idmax),zdir(:,1:idmax),llimit_qair         )
  ENDIF

#ifdef SFX_MPI
  xtime_calc(1) = xtime_calc(1) + (mpi_wtime() - xtime1)
  xtime1 = mpi_wtime()
#endif
  !    
  !COMPUTE SUM ZENITH angle between 2 timestepA
  zsumzen(:)=0.0
  DO jsurf_step = 1,inb_atm
    iday2  = iday
    ztime2 = ztime + (jsurf_step-1.)*xtstep_surf
    IF (ztime2>86400.) THEN
      ztime2 = ztime2-86400
      iday2  = iday+1
    ENDIF
    CALL sunpos(iyear, imonth, iday2, ztime+(jsurf_step-1.)*xtstep_surf, &
                zlon, zlat, xtsun, xzenith, xazim)
    !
    zsumzen(:)= zsumzen(:) + max(cos(xzenith(:)+0.1),0.)/(inb_atm*1.0)
    !
  ENDDO
  WHERE ( zsumzen<0.01 ) zsumzen = 0.0
  !
  DO jsurf_step=1,inb_atm
    !
    ! time interpolation of the forcing
    !
#ifdef SFX_MPI
    xtime1 = mpi_wtime()
#endif
    !
    CALL sunpos(iyear, imonth, iday, ztime, zlon, zlat, xtsun, xzenith, xazim)
    iyear2 = iyear
    imonth2= imonth
    iday2  = iday
    ztime2 = ztime+xtstep_surf
    CALL add_forecast_to_date_surf(iyear2, imonth2, iday2, ztime2)
    CALL sunpos(iyear2, imonth2, iday2, ztime2, zlon, zlat, xtsun, xzenith2, xazim)
    !
#ifdef SFX_MPI
    xtime_calc(2) = xtime_calc(2) + (mpi_wtime() - xtime1)
    xtime1 = mpi_wtime()
#endif
    !interpolation between beginning and end of current forcing time step
    CALL ol_time_interp_atm(jsurf_step,inb_atm,                      &
                            zta(:,id_forc),zta(:,id_forc+1),         &
                            zqa(:,id_forc),zqa(:,id_forc+1),         &
                            zwind(:,id_forc),zwind(:,id_forc+1),     &
                            zdir_sw(:,id_forc),zdir_sw(:,id_forc+1), &
                            zsca_sw(:,id_forc),zsca_sw(:,id_forc+1), &
                            zlw(:,id_forc),zlw(:,id_forc+1),         &
                            zsnow(:,id_forc+1),zrain(:,id_forc+1),   &
                            zps(:,id_forc),zps(:,id_forc+1),         &
                            zco2(:,id_forc),zco2(:,id_forc+1),       &
                            zdir(:,id_forc),zdir(:,id_forc+1),       &
                            xzenith+0.1,zsumzen                      )
#ifdef SFX_MPI
    xtime_calc(3) = xtime_calc(3) + (mpi_wtime() - xtime1)
    xtime1 = mpi_wtime()
#endif
    !
    IF(ladapt_sw)THEN
      !
      ! coherence between solar zenithal angle and radiation
      ! when solar beam close to horizontal -> reduction of direct radiation to
      ! the benefit of scattered radiation
      ! when pi/2 - 0.1 < ZENITH < pi/2 - 0.05 => weight of direct to scattered radiation decreases linearly with zenith 
      ! when pi/2 - 0.05 < ZENITH => all the direct radiation is converted to scattered radiation
      ! coherence between solar zenithal angle and radiation
      !
      zcoef(:) = (xpi/2. - xzenith(:) - 0.05) / 0.05
      zcoef(:) = max(min(zcoef,1.),0.)
      DO jloop=1,SIZE(xdir_sw,2)
        xsca_sw(:,jloop) = xsca_sw(:,jloop) + xdir_sw(:,jloop) * (1 - zcoef)
        xdir_sw(:,jloop) = xdir_sw(:,jloop) * zcoef(:)
      ENDDO
      !
    ELSE
      !
      zsw(:) = 0.
      DO jloop=1,SIZE(xdir_sw,2)
        zsw(:) = zsw(:) + xdir_sw(:,jloop) + xsca_sw(:,jloop)
      END DO
      WHERE (zsw(:)>0.)
        xzenith  = min(xzenith ,xpi/2.-0.01)
        xzenith2 = min(xzenith2,xpi/2.-0.01)
      ELSEWHERE
        xzenith  = max(xzenith ,xpi/2.)
        xzenith2 = max(xzenith2,xpi/2.)
      END WHERE
      !
    ENDIF       
    !
    ! updates time
    ztimec= ztimec+xtstep_surf
    IF (lcoupl_topd) ltopd_step = ( mod((((jforc_step-1)*inb_atm)+jsurf_step),nnb_topd) == 0 )
    !
    ! run Surface
    !
#ifdef SFX_MPI
    xtime_calc(4) = xtime_calc(4) + (mpi_wtime() - xtime1)
#endif
    !
    CALL io_buff_clean
    !
    IF(loasis)THEN
     ! Receive fields to other models proc by proc
     CALL sfx_oasis_recv_ol(ysc%FM%F, ysc%IM, ysc%SM%S, ysc%U, ysc%WM%W,               &
                            csurf_filetype, ini, ibands, ztimec, xtstep_surf, xzenith, &
                            xsw_bands, xtsrad, xdir_alb, xsca_alb, xemis, xtsurf   )
    ENDIF
    !
#ifdef SFX_MPI
    xtime1 = mpi_wtime()
#endif
    !
    IF(gshadows) THEN 
      CALL slope_radiative_effect(xtstep_surf, xzenith, xazim, xps, xta, xrain, xdir_sw, xlw, &
                                  xzs_thread, xzs_xy_thread, xslopang_thread, xslopazi_thread,&
                                  xsurf_triangle_thread)
    END IF
    !
    CALL coupling_surf_atm_n(ysc, csurf_filetype, 'E', ztimec, xtstep_surf, iyear, imonth, iday, ztime, &
                             ini, nscal, ibands, xtsun, xzenith, xzenith2, xazim, xzref, xuref,         &
                             xzs, xu, xv, xqa, xta, xrhoa, xsv, xco2, csv, xrain, xsnow, xlw, xdir_sw,  &
                             xsca_sw, xsw_bands, xps, xpa, xsftq, xsfth, xsfts, xsfco2, xsfu, xsfv,     &
                             xtsrad, xdir_alb, xsca_alb, xemis, xtsurf, xz0, xz0h, xqsurf, xpew_a_coef, &
                             xpew_b_coef,xpet_a_coef,xpeq_a_coef,xpet_b_coef,xpeq_b_coef, ytest      )
    !
#ifdef SFX_MPI
    xtime_calc(5) = xtime_calc(5) + (mpi_wtime() - xtime1)
#endif
    !
#ifdef SFX_MPI
    xtime1 = mpi_wtime()
#endif
    !
    IF(loasis)THEN
     ! Send fields to other models proc by proc
     CALL sfx_oasis_send_ol(ysc%FM%F, ysc%IM, ysc%SM%S, ysc%U, ysc%WM%W, csurf_filetype,ini,ztimec,xtstep_surf)
    ENDIF
    !
    ztime = ztime + xtstep_surf
    CALL add_forecast_to_date_surf(iyear, imonth, iday, ztime)
#ifdef SFX_MPI
    xtime_calc(6) = xtime_calc(6) + (mpi_wtime() - xtime1)
    !
    xtime1 =  mpi_wtime()
#endif
    ! ecrit Surface
    !
    IF ( lcoupl_topd .AND. ltopd_step ) THEN
      !
      IF (.NOT.ALLOCATED(nyear))  ALLOCATE(nyear(nnb_topd_step))
      IF (.NOT.ALLOCATED(nmonth)) ALLOCATE(nmonth(nnb_topd_step))
      IF (.NOT.ALLOCATED(nday))   ALLOCATE(nday(nnb_topd_step))
      IF (.NOT.ALLOCATED(nh))     ALLOCATE(nh(nnb_topd_step))
      IF (.NOT.ALLOCATED(nm))     ALLOCATE(nm(nnb_topd_step))
      !
      nyear(ntopd_step) = iyear
      nmonth(ntopd_step) = imonth
      nday(ntopd_step) = iday
      nh(ntopd_step) = int(ztime/3600.)
      nm(ntopd_step) = int((ztime-nh(ntopd_step)*3600.)/60.)
      !
      IF ( nm(ntopd_step)==60 ) THEN
        !
        nm(ntopd_step) = 0
        nh(ntopd_step) = nh(ntopd_step)+1
        !
      ENDIF
      !
      IF ( nh(ntopd_step)==24 ) THEN
        !
        nh(ntopd_step) = 0
        nday(ntopd_step) = nday(ntopd_step)+1
        !
        !!AJOUT BEC 
        SELECT CASE (nmonth(ntopd_step))
          CASE(4,6,9,11)
            IF ( nday(ntopd_step)==31 ) THEN
              nmonth(ntopd_step) = nmonth(ntopd_step)+1
              nday(ntopd_step) = 1
            ENDIF
          CASE(1,3,5,7:8,10)
            IF ( nday(ntopd_step)==32 ) THEN
              nmonth(ntopd_step) = nmonth(ntopd_step)+1
              nday(ntopd_step) = 1
            ENDIF
          CASE(12)
            IF ( nday(ntopd_step)==32 ) THEN
              nyear(ntopd_step) = nyear(ntopd_step)+1
              nmonth(ntopd_step) = 1
              nday(ntopd_step) = 1
            ENDIF
          CASE(2)
            IF( mod(nyear(ntopd_step),4)==0 .AND. mod(nyear(ntopd_step),100)/=0 .OR. mod(nyear(ntopd_step),400)==0 ) THEN 
              IF (nday(ntopd_step)==30) THEN
                nmonth(ntopd_step) = nmonth(ntopd_step)+1
                nday(ntopd_step) = 1
              ENDIF
            ELSE
              IF (nday(ntopd_step)==29) THEN
                nmonth(ntopd_step) = nmonth(ntopd_step)+1
                nday(ntopd_step) = 1
              ENDIF
            ENDIF
        END SELECT
        !
      ENDIF
      !
      ! * 2. Stocking date of each time step
      !
      ntopd_step = ntopd_step + 1 
      !     
    ENDIF
    !
    IF (mod(ztimec,xtstep_output) == 0. .AND. ctimeseries_filetype/='NONE  ') THEN
      !
      IF (nrank==npio) THEN
        !
        !* name of the file
        IF (ctimeseries_filetype=="ASCII " .OR. &
            ctimeseries_filetype=="LFI   " .OR. &
            ctimeseries_filetype=="FA    " .OR. &
            ctimeseries_filetype=="NC    "    ) THEN  
          !
          !
          ztime_out  = ztime
          iday_out   = iday
          imonth_out = imonth
          iyear_out  = iyear
          !
          IF(lout_timename)THEN
            ! if true, change the name of output file at the end of a day
            ! (ex: 19860502_00h00 -> 19860501_24h00)                     
            IF(ztime==0.0)THEN
              ztime_out = 86400.
              iday_out   = iday-1
              IF(iday_out==0)THEN
                imonth_out = imonth - 1
                IF(imonth_out==0)THEN
                  imonth_out=12
                  iyear_out = iyear - 1
                ENDIF
                SELECT CASE (imonth_out)
                  CASE(4,6,9,11)
                    iday_out=30
                  CASE(1,3,5,7:8,10,12)
                    iday_out=31
                  CASE(2)
                    IF( ((mod(iyear_out,4)==0).AND.(mod(iyear_out,100)/=0)) .OR. (mod(iyear_out,400)==0))THEN 
                      iday_out=29
                    ELSE
                      iday_out=28
                    ENDIF
                END SELECT
              ENDIF
            ENDIF
            !
          ENDIF
          !
          WRITE(ytag,fmt='(I4.4,I2.2,I2.2,A1,I2.2,A1,I2.2)') iyear_out,imonth_out,iday_out,&
                '_',int(ztime_out/3600.),'h',nint(ztime_out)/60-60*int(ztime_out/3600.)  
          cfileout    = adjustl(adjustr(csurffile)//'.'//ytag//'.txt')
          cfileout_lfi= adjustl(adjustr(csurffile)//'.'//ytag)
          cfileout_fa = adjustl(adjustr(csurffile)//'.'//ytag//'.fa')
          cfileout_nc = adjustl(adjustr(csurffile)//'.'//ytag//'.nc')
          !
          IF (ctimeseries_filetype=='FA    ') THEN
#ifdef SFX_FA                  
            lfanocompact = ldiag_fa_nocompact
            idatef(1)= iyear!_OUT
            idatef(2)= imonth!_OUT
            idatef(3)= iday!_OUT
            !ZTIME instead of ZTIME_OUT (FA XRD39 do not like 24h)
            idatef(4)= floor(ztime/3600.)
            idatef(5)= floor(ztime/60.) - idatef(4) * 60 
            idatef(6)= nint(ztime) - idatef(4) * 3600 - idatef(5) * 60
            idatef(7:11) = 0
            nunit_fa = 19
            IF (csurf_filetype/='FA    ') THEN
              CALL write_header_fa(ysc%GCP, ysc%UG%G%CGRID, ysc%UG%XGRID_FULL_PAR, csurf_filetype,'ALL')
            ELSE
              CALL faitou(iret,nunit_fa,.true.,cfileout_fa,'UNKNOWN',.true.,.false.,iverbfa,0,inb,cdnomc)
            ENDIF
            CALL fandar(iret,nunit_fa,idatef)
#endif            
          END IF
          !
        END IF
        !
        xstartw = xstartw + 1
        nwrite  = nwrite  + 1
        ltime_written=.false.
        !
      ENDIF
      !
#ifdef SFX_MPI
      xtime_write(1) = xtime_write(1) + (mpi_wtime() - xtime1)
#endif
      !
      inw = 1
      !
      IF ( lxios .AND. .NOT. lxios_def_closed ) inw = 2      
      !
      ldef_nc = .false. 
      IF (ctimeseries_filetype=="NC    ") THEN
        ldef_nc = .true.       
        inw = 2
        CALL init_output_nc_n (ysc%TM%BDD, ysc%CHE, ysc%CHN, ysc%CHU, ysc%SM%DTS, &
                               ysc%TM%DTT, ysc%DTZ, ysc%IM, ysc%UG, ysc%U, ysc%DUO%CSELECT)                   
      ENDIF
      !
      ldef_ol = .false.
      IF (ctimeseries_filetype=="OFFLIN".AND.lfirst_write) THEN
        ldef_ol = .true.
        inw = 2
        IF (ctimeseries_filetype=="OFFLIN") CALL init_output_ol_n (ysc)
      ENDIF
      !
      idx_w = 0
      !
      DO jnw = 1,inw
        !
        CALL io_buff_clean
        !
#ifdef SFX_MPI
        xtime1 =  mpi_wtime()
#endif
        IF (lxios) THEN 
#ifdef WXIOS
          ntimestep=int(ztimec/xtstep_output + 1.)
#endif
        ENDIF
        CALL write_surf_atm_n(ysc, ctimeseries_filetype,'ALL',lland_use)
#ifdef SFX_MPI
        xtime_write(2) = xtime_write(2) + (mpi_wtime() - xtime1)
        xtime1 =  mpi_wtime()
#endif
        CALL diag_surf_atm_n(ysc, ctimeseries_filetype)
#ifdef SFX_MPI
        xtime_write(3) = xtime_write(3) + (mpi_wtime() - xtime1)
        xtime1 =  mpi_wtime()
#endif
        CALL write_diag_surf_atm_n(ysc, ctimeseries_filetype,'ALL')
        !
        IF (lxios) THEN 
#ifdef WXIOS
          IF (.NOT. lxios_def_closed) THEN 
            CALL xios_close_context_definition()
            lxios_def_closed=.true.
          ENDIF
          CALL xios_update_calendar(ntimestep)
#endif
        ENDIF
        !      
#ifdef SFX_MPI
        xtime_write(4) = xtime_write(4) + (mpi_wtime() - xtime1)
#endif
        !
        ldef_nc = .false.
        ldef_ol = .false.
        !
        ncpt_write = 0
        !
        lfirst_write = .false.
        !
      ENDDO
      !
      IF (lcoupl_topd .AND. ntopd_step > nnb_topd_step) THEN
        !
        ! Writing of file resulting of coupling with TOPMODEL or routing ****
        CALL write_discharge_file(csurf_filetype,'q_total.txt','FORMATTED',&
                                  nyear,nmonth,nday,nh,nm,xqtot)
        CALL write_discharge_file(csurf_filetype,'q_runoff.txt','FORMATTED',&
                                  nyear,nmonth,nday,nh,nm,xqb_run)
        CALL write_discharge_file(csurf_filetype,'q_drainage.txt','FORMATTED',&
                                  nyear,nmonth,nday,nh,nm,xqb_dr)
        ! Writing of budget files 
        IF (lbudget_topd) CALL write_budget_coupl_rout
        !
      ENDIF
      !
#ifdef SFX_MPI
      xtime1 =  mpi_wtime()
#endif
      !
      IF (nrank==npio) THEN
        IF (ctimeseries_filetype=='FA    ') THEN
#ifdef SFX_FA                
          CALL fairme(iret,nunit_fa,'UNKNOWN')
#endif          
        END IF
        !* add informations in the file
        IF (ctimeseries_filetype=='LFI   ' .AND. lmnh_compatible) CALL write_header_mnh
      ENDIF
#ifdef SFX_MPI
      xtime_write(5) = xtime_write(5) + (mpi_wtime() - xtime1)
#endif
      !
    ENDIF
    !
  END DO
  !
  IF (nrank==npio) THEN
    IF (lprint) THEN
      IF (mod(ztimec,xday) == 0.) THEN
        icount = icount + 1
        CALL wlog_mpi('SFX DAY :',klog=icount,klog2=int(zduration/xday))
        WRITE(*,'(A10,I5,A2,I5)')'SFX  DAY :',icount,' /',int(zduration/xday)
      ENDIF
    ENDIF
  ENDIF
  !
END DO
!
 CALL wlog_mpi(' ')
 CALL wlog_mpi('OL_READ_ATM ',plog=xtime_calc(1))
 CALL wlog_mpi('SUNPOS ',plog=xtime_calc(2))
 CALL wlog_mpi('OL_TIME_INTERP_ATM ',plog=xtime_calc(3))
 CALL wlog_mpi('')
 CALL wlog_mpi('ZENITH ',plog=xtime_calc(4))
 CALL wlog_mpi('')
 CALL wlog_mpi('COUPLING_SURF_ATM ',plog=xtime_calc(5))
 CALL wlog_mpi('')
 CALL wlog_mpi('ADD_FORECAST_TO_DATE_SURF ',plog=xtime_calc(6))
 CALL wlog_mpi('DEF_DATE ',plog=xtime_write(1))
 CALL wlog_mpi('')
 CALL wlog_mpi('WRITE_SURF_ATM ',plog=xtime_write(2))
 CALL wlog_mpi('DIAG_SURF_ATM ',plog=xtime_write(3))
 CALL wlog_mpi('WRITE_DIAG_SURF_ATM ',plog=xtime_write(4))
 CALL wlog_mpi('')
 CALL wlog_mpi('CLOSE FILES ',plog=xtime_write(5))
 CALL wlog_mpi('')
!
IF (cforcing_filetype=='ASCII ' .OR. cforcing_filetype=='BINARY') &
        CALL open_close_bin_asc_forc('CLOSE',cforcing_filetype,'R')
!
IF (cforcing_filetype=='NETCDF') CALL close_filein_ol
IF (ctimeseries_filetype=='OFFLIN') CALL close_fileout_ol
!
! --------------------------------------------------------------------------------------
!
!*    3.     write restart file
!            ------------------
!
IF ( lrestart ) THEN
  !
  lfirst_write = .true.
  !
  IF (nrank==npio) THEN
    !* name of the file
    cfileout    = adjustl(adjustr(csurffile)//'.txt')
    cfileout_lfi= csurffile
    cfileout_fa = adjustl(adjustr(csurffile)//'.fa')
    cfileout_nc = adjustl(adjustr(csurffile)//'.nc')

    !* opens the file
    IF (csurf_filetype=='FA    ') THEN
#ifdef SFX_FA            
      lfanocompact = .true.
      idatef(1)= iyear
      idatef(2)= imonth
      idatef(3)= iday
      idatef(4)= floor(ztime/3600.)
      idatef(5)= floor(ztime/60.) - idatef(4) * 60  
      idatef(6)= nint(ztime) - idatef(4) * 3600 - idatef(5) * 60
      idatef(7:11) = 0   
      nunit_fa = 19 
      CALL faitou(iret,nunit_fa,.true.,cfileout_fa,'UNKNOWN',.true.,.false.,iverbfa,0,inb,cdnomc)
      CALL fandar(iret,nunit_fa,idatef)
#endif      
    END IF
    !
  ENDIF
  !
  inw = 1
  IF (csurf_filetype=="NC    ") inw = 2
  !
  ldef_nc = .true.
  ldef_ol = .true.
  !
  IF (ASSOCIATED(ysc%DUO%CSELECT)) DEALLOCATE(ysc%DUO%CSELECT)
  ALLOCATE(ysc%DUO%CSELECT(0))
  !
  IF (csurf_filetype=="NC    ") THEN
    CALL init_output_nc_n (ysc%TM%BDD, ysc%CHE, ysc%CHN, ysc%CHU, ysc%SM%DTS, ysc%TM%DTT, &
                           ysc%DTZ, ysc%IM, ysc%UG, ysc%U, ysc%DUO%CSELECT)
  ENDIF
  !
  DO jnw = 1,inw
    !
    CALL io_buff_clean
    !
    CALL flag_update(ysc%IM%ID%O, ysc%DUO,.false.,.true.,.false.,.false.)
    !
    IF (lrestart_2m) THEN
      i2m       = 1
      gpgd_isba = .true.
    ELSE
      i2m       = 0
      gpgd_isba = .false.
    ENDIF
    gfrac                  = .true.
    gdiag_grid             = .true.
    gsurf_budget           = .false.
    grad_budget            = .false.
    gcoef                  = .false.
    gsurf_vars             = .false.
    ibeq                   = 0
    idsteq                 = 0
    gdiag_ocean            = .false.
    gdiag_seaice           = .false.
    gwater_profile         = .false.
    gsurf_evap_budget      = .false.
    gflood                 = .false.
    gch_no_flux_isba       = .false.
    gsurf_misc_budget_isba = .false.
    gpgd_teb               = .false.
    gsurf_misc_budget_teb  = .false.  
    !
    CALL flag_diag_update(ysc%FM, ysc%IM, ysc%SM, ysc%TM, ysc%WM, ysc%DUO, ysc%U, ysc%SV, &
                          gfrac, gdiag_grid, i2m, gsurf_budget, grad_budget, gcoef,  &
                          gsurf_vars, ibeq, idsteq, gdiag_ocean, gdiag_seaice,       &
                          gwater_profile, gsurf_evap_budget, gflood,  gpgd_isba,     &
                          gch_no_flux_isba, gsurf_misc_budget_isba, gpgd_teb,        &
                          gsurf_misc_budget_teb    )
    !
    ysc%DUO%LSNOWDIMNC = .false.
    !
    !* writes into the file
    CALL write_surf_atm_n(ysc, csurf_filetype,'ALL',lland_use)
    IF(csurf_filetype/='FA    ' .OR. lrestart_2m) THEN
       CALL write_diag_surf_atm_n(ysc, csurf_filetype,'ALL')
    ENDIF
    !
    ldef_nc = .false.
    ldef_ol = .false.
    !
    ncpt_write = 0
    !
    lfirst_write = .false.
    !
  ENDDO
  !
  !* closes the file
  IF (nrank==0 ) THEN
    IF (csurf_filetype=='FA    ') THEN
#ifdef SFX_FA            
      CALL fairme(iret,nunit_fa,'UNKNOWN')
#endif      
    END IF
    !* add informations in the file
    IF (csurf_filetype=='LFI   ' .AND. lmnh_compatible) CALL write_header_mnh
    !
  ENDIF
  !
  IF (lcoupl_topd .AND. ntopd_step > nnb_topd_step) &
          CALL prep_restart_coupl_topd(ysc%UG, ysc%U, csurf_filetype,ini)
  !
END IF
!
! --------------------------------------------------------------------------------------
!
!*    4.     inquiry mode
!            ------------
!
IF ( linquire ) THEN
  !
  ALLOCATE( zsea( ini ) )
  ALLOCATE( zwater( ini ) )
  ALLOCATE( znature( ini ) )
  ALLOCATE( ztown( ini ) )
  ALLOCATE( zt2m( ini ) )
  ALLOCATE( zq2m( ini ) )
  ALLOCATE( zz0( ini ) )
  ALLOCATE( zz0h( ini ) )
  ALLOCATE( zqs_sea( ini ) )
  ALLOCATE( zqs_water( ini ) )
  ALLOCATE( zqs_nature( ini ) )
  ALLOCATE( zqs_town( ini ) )
  ALLOCATE( zqs( ini ) )
  ALLOCATE( zpsng( ini ) )
  ALLOCATE( zpsnv( ini ) )
  ALLOCATE( zz0eff( ini ) )
  ALLOCATE( zzs( ini ) )
  !
  iseries = 0
  CALL get_surf_var_n(ysc%FM, ysc%IM, ysc%SM, ysc%TM, ysc%WM, ysc%DUO, ysc%DU,  ysc%UG, ysc%U, ysc%USS, &
                      csurf_filetype,ini,iseries,psea=zsea,pwater=zwater,pnature=znature,ptown=ztown, &
                      pt2m=zt2m,pq2m=zq2m,pqs=zqs,pz0=zz0,pz0h=zz0h,pz0eff=zz0eff,pqs_sea=zqs_sea,  &
                      pqs_water=zqs_water,pqs_nature=zqs_nature,pqs_town=zqs_town,                  &
                      ppsng=zpsng,ppsnv=zpsnv,pzs=zzs                                         )  
  !
  isize = SIZE(nindex)
  IF (nrank==npio) THEN
    ALLOCATE(zsea_full(isize))
    ALLOCATE(zwater_full(isize))
    ALLOCATE(znature_full(isize))
    ALLOCATE(ztown_full(isize))
    ALLOCATE(zz0_full(isize))
    ALLOCATE(zz0eff_full(isize))
    ALLOCATE(zzs_full(isize))
  ELSE
    ALLOCATE(zsea_full(0))
    ALLOCATE(zwater_full(0))
    ALLOCATE(znature_full(0))
    ALLOCATE(ztown_full(0))
    ALLOCATE(zz0_full(0))
    ALLOCATE(zz0eff_full(0))
    ALLOCATE(zzs_full(0))
  ENDIF
  CALL gather_and_write_mpi(zsea,zsea_full) 
  CALL gather_and_write_mpi(zwater,zwater_full)
  CALL gather_and_write_mpi(znature,znature_full)
  CALL gather_and_write_mpi(ztown,ztown_full)
  CALL gather_and_write_mpi(zz0,zz0_full)
  CALL gather_and_write_mpi(zz0eff,zz0eff_full)
  CALL gather_and_write_mpi(zzs,zzs_full)

  IF (nrank==npio) THEN
    WRITE(iluout,'(A32,I4,A3,I4)') ' GRID BOXES CONTAINING SEA    : ',count( zsea_full(:) > 0. ),' / ',isize
    WRITE(iluout,'(A32,I4,A3,I4)') ' GRID BOXES CONTAINING WATER  : ',count( zwater_full(:) > 0. ),' / ',isize
    WRITE(iluout,'(A32,I4,A3,I4)') ' GRID BOXES CONTAINING NATURE : ',count( znature_full(:) > 0. ),' / ',isize
    WRITE(iluout,'(A32,I4,A3,I4)') ' GRID BOXES CONTAINING TOWN   : ',count( ztown_full(:) > 0. ),' / ',isize
    WRITE(iluout,*)'ZZ0    = ',zz0_full
    WRITE(iluout,*)'ZZ0EFF = ',zz0eff_full
    WRITE(iluout,*)'ZZS = ',zzs_full
    WRITE(iluout,*)'MINVAL(ZZS) = ',minval(zzs_full),' MAXVAL(ZZS) = ',maxval(zzs_full)
  ENDIF
  !
  DEALLOCATE( zsea       )
  DEALLOCATE( zwater     )
  DEALLOCATE( znature    )
  DEALLOCATE( ztown      )
  DEALLOCATE( zt2m       )
  DEALLOCATE( zq2m       )
  DEALLOCATE( zz0        )
  DEALLOCATE( zz0h       )
  DEALLOCATE( zqs_sea    )
  DEALLOCATE( zqs_water  )
  DEALLOCATE( zqs_nature )
  DEALLOCATE( zqs_town   )
  DEALLOCATE( zqs        )
  DEALLOCATE( zpsng      )
  DEALLOCATE( zpsnv      )
  DEALLOCATE( zz0eff     )
  DEALLOCATE( zzs        )
  !
  IF (nrank==npio) THEN
    DEALLOCATE(zsea_full   )
    DEALLOCATE(zwater_full )
    DEALLOCATE(znature_full)
    DEALLOCATE(ztown_full  )
    DEALLOCATE(zz0_full    )
    DEALLOCATE(zz0eff_full )
    DEALLOCATE(zzs_full    )
  ENDIF
  !
ENDIF   
!
! --------------------------------------------------------------------------------------
!
!    4'    Close Gelato specific diagnostic 
#if ! defined in_arpege
 CALL clsdia()
#endif
!
!
!*    5.     Close parallelized I/O
!            ----------------------
!
IF (nrank==npio) THEN
  WRITE(iluout,*) ' '
  WRITE(iluout,*) '    --------------------------'
  WRITE(iluout,*) '    | OFFLINE ENDS CORRECTLY |'
  WRITE(iluout,*) '    --------------------------'
  WRITE(iluout,*) ' '
  CLOSE(iluout)
  WRITE(*,*) ' '
  WRITE(*,*) '    --------------------------'
  WRITE(*,*) '    | OFFLINE ENDS CORRECTLY |'
  WRITE(*,*) '    --------------------------'
  WRITE(*,*) ' '
ENDIF
!
 CALL surfex_deallo_list
!
IF (ALLOCATED(nindex)) DEALLOCATE(nindex)
IF (ALLOCATED(nsize_task)) DEALLOCATE(nsize_task)
!
 CALL end_log_mpi
!
IF (lhook) CALL dr_hook('OFFLINE',1,zhook_handle)
!
! * MPI and OASIS must be finalized after the last DR_HOOK call
!
IF (lxios) THEN 
#ifdef WXIOS
  CALL xios_context_finalize()
#endif
ENDIF
!
 CALL sfx_oasis_end
!
#ifdef SFX_MPI
IF(.NOT. loasis .AND. .NOT. lxios) THEN
  CALL mpi_finalize(infompi)
ENDIF
#endif
!
! --------------------------------------------------------------------------------------
!
END PROGRAM offline