module read_rsun_ecrad_m

  use radiation_config      
!  USE setup_gas_optics

  implicit none

  REAL, DIMENSION(14) :: SPSOLARSCAL 
!  type(config_type) :: config

contains

  SUBROUTINE read_rsun_ecrad(debut)

!****************************************************************************************
! This routine will read the solar constant fraction per band
!
! Olivier Boucher with inputs from Marion Marchand 
!****************************************************************************************

  USE netcdf95, ONLY: nf95_close, nf95_inq_varid, nf95_open, nf95_gw_var
  USE netcdf, ONLY: nf90_get_var, nf90_noerr, nf90_nowrite

  USE phys_cal_mod, ONLY : days_elapsed, year_len

  USE mod_phys_lmdz_mpi_data, ONLY: is_mpi_root
  USE mod_phys_lmdz_omp_data, ONLY: is_omp_root
  USE mod_phys_lmdz_para
  USE clesphys_mod_h

!  USE YOESW, ONLY : RSUN

  IMPLICIT NONE

  !INCLUDE "clesphys.h"

  ! Input arguments
  LOGICAL, INTENT(IN) :: debut
  
! Local variables
  INTEGER :: ncid, dimid, varid, ncerr, nbday
  REAL, ALLOCATABLE :: wlen(:), time(:)
  REAL, ALLOCATABLE, SAVE, DIMENSION(:,:) :: SSI_FRAC
!$OMP THREADPRIVATE(SSI_FRAC)
  REAL, ALLOCATABLE, SAVE, DIMENSION(:) :: TSI(:)
!$OMP THREADPRIVATE(TSI)

  INTEGER, SAVE :: day_pre=1
!$OMP THREADPRIVATE(day_pre)

  ! Output arguments
!  REAL, DIMENSION(14), INTENT(out) :: SPSOLARSCAL

!--only one processor reads
    IF (debut) THEN 

    ALLOCATE(SSI_FRAC(14,year_len))

    ALLOCATE(TSI(year_len))

    IF (is_mpi_root.AND.is_omp_root) THEN

       CALL nf95_open('solarforcing.nc', NF90_NOWRITE, ncid)

       CALL nf95_inq_varid(ncid, 'wlen', varid)
       CALL nf95_gw_var(ncid, varid, wlen)

       CALL nf95_inq_varid(ncid, 'time', varid)
       CALL nf95_gw_var(ncid, varid, time)

       IF (14.NE.size(wlen)) THEN 
         PRINT *,'read_rsun_ecrad 14 <> wlen = ',14, size(wlen)
         CALL abort_physic('read_rsun_ecrad','size of SSI is different from 14',1)
       ENDIF

!--test if time is different from year_len but allow a mismatch of 1 day
       IF (size(time).NE.year_len.AND.size(time).NE.year_len+1) THEN
         PRINT *,'read_rsun_ecrad time <> year_len = ', size(time), year_len
         CALL abort_physic('read_rsun_ecrad','time dim should be the number of days in year',1)
       ENDIF
!--warning only if forcing file has 366 days but year_len has only 365
       IF (size(time).EQ.year_len+1) THEN 
         PRINT *,'Warning read_rsun_ecrad uses a leap year rsun for a noleap year'
       ENDIF

       CALL nf95_inq_varid(ncid, 'ssi_frac', varid)
       ncerr = nf90_get_var(ncid, varid, SSI_FRAC)

       CALL nf95_inq_varid(ncid, 'tsi', varid)
       ncerr = nf90_get_var(ncid, varid, TSI)

       CALL nf95_close(ncid)

       DO nbday=1, year_len
! AI ATTENTION       
!         IF (ABS(SUM(SSI_FRAC(:,nbday))-1.).GT.1.e-6) THEN 
         IF (ABS(SUM(SSI_FRAC(:,nbday))-1.).GT.1.) THEN
           PRINT *,'somme SSI_FRAC=', SUM(SSI_FRAC(:,nbday))
           CALL abort_physic('read_rsun_ecrad','somme SSI_FRAC <> 1',1)
         ENDIF
       ENDDO
     
    ENDIF ! is_mpi_root .AND. is_omp_root

!$OMP BARRIER
    CALL bcast(SSI_FRAC)
    CALL bcast(TSI)

    ENDIF

!--only read at beginning of day
!--day in year is provided as days_elapsed since the beginning of the year +1
    IF (debut.OR.days_elapsed+1.NE.day_pre) THEN

!--keep memory of previous day
      day_pre=days_elapsed+1

!--copy 
      SPSOLARSCAL(1:14)=SSI_FRAC(:,days_elapsed+1)
      solaire=TSI(days_elapsed+1)

      print *,'READ_SPSOLARSCAL ECRAD day=', days_elapsed+1,' solaire=', &
                               solaire, ' SPSOLARSCAL=', SPSOLARSCAL(1:14)

    ENDIF !--fin allocation

  END SUBROUTINE read_rsun_ecrad

end module read_rsun_ecrad_m  
