!==================================================================
module radii_mod
!==================================================================
!  module to centralize the radii calculations for aerosols
!==================================================================

!     N2 cloud properties (initialized in inifis)
      real,save :: Nmix_n2 ! Number mixing ratio of N2 ice particles
!$OMP THREADPRIVATE(Nmix_n2)

      ! flag to specify if we assume a constant fixed radius for particles
      logical,save :: radfixed ! initialized in inifis
!$OMP THREADPRIVATE(radfixed)

contains


!==================================================================
   subroutine su_aer_radii(ngrid,nlayer,reffrad,nueffrad)
!==================================================================
!     Purpose
!     -------
!     Compute the effective radii of haze particles  (TB)
!
!
!==================================================================
      use mod_phys_lmdz_para, only : is_master
      use ioipsl_getin_p_mod, only: getin_p
      use radinc_h, only: naerkind
      use datafile_mod, only: hazerad_file
      use aerosol_mod, only: iaero_haze, i_haze, &
                              iaero_generic, i_rgcs_ice
      use callkeys_mod, only: nlayaero, aeronlay_size, &
                              aeronlay_nueff,aerogeneric
      use tracer_h, only: radius, nqtot, is_rgcs, nmono
      Implicit none

      integer,intent(in) :: ngrid
      integer,intent(in) :: nlayer

      real, intent(out) :: reffrad(ngrid,nlayer,naerkind)      !aerosols radii (K)
      real, intent(out) :: nueffrad(ngrid,nlayer,naerkind)     !variance

      integer :: iaer, ia , iq, i_rad

      do iaer=1,naerkind
         print*, 'TB22 iaer',iaer,iaero_haze
         if(iaer.eq.iaero_haze)then
           ! Equivalent sphere radius
           reffrad(1:ngrid,1:nlayer,iaer)=radius(i_haze)*nmono**(1./3.)
           !reffrad(1:ngrid,1:nlayer,iaer) = 2.e-6 ! haze
           nueffrad(1:ngrid,1:nlayer,iaer) = 0.02 ! haze
         endif
      enddo

   end subroutine su_aer_radii
!==================================================================


!==================================================================
   subroutine haze_reffrad_fix(ngrid,nlayer,zzlay,reffrad,nueffrad)
      !==================================================================
!     Purpose
!     -------
!     Compute the effective radii of haze particles
!     fixed profile of radius (TB)
!
!==================================================================
      use radinc_h, only: naerkind
      USE tracer_h, only:rho_n2,nmono
      use comcstfi_mod, only: pi
      use aerosol_mod, only: iaero_haze, i_haze
      use datafile_mod
      Implicit none

      integer,intent(in) :: ngrid,nlayer
      real,intent(in) :: zzlay(ngrid,nlayer)
      real, intent(out) :: reffrad(ngrid,nlayer,naerkind)      ! haze particles radii (m)
      real, intent(out) :: nueffrad(ngrid,nlayer,naerkind)     !

      real :: zrad
      real,external :: CBRT

      logical, save :: firstcall=.true.
      !$OMP THREADPRIVATE(firstcall)

!     Local variables
      integer :: iaer,l,ifine,ig
      real :: radvec(ngrid,nlayer)

      !!read altitudes and radius
      integer Nfine
      !parameter(Nfine=21)
      parameter(Nfine=701)
      character(len=100) :: file_path
      real,save :: levdat(Nfine),raddat(Nfine)

!---------------- INPUT ------------------------------------------------

      IF (firstcall) then
         firstcall=.false.
         file_path=trim(datadir)//'/haze_prop/'//hazerad_file
         open(223,file=file_path,form='formatted')
         do ifine=1,Nfine
            read(223,*) levdat(ifine), raddat(ifine)
         enddo
         close(223)
         print*, 'TB22 READ HAZERAD'
       ENDIF

       ! in radii mod levs has been put in km
       DO ig=1,ngrid
         CALL interp_line(levdat,raddat,Nfine,zzlay(ig,:)/1000,radvec(ig,:),nlayer)
       enddo

       do iaer=1,naerkind
             if(iaer.eq.iaero_haze)then
                  ! spherical radius or eq spherical radius
                  ! TB22: fractal has no impact on reffrad if haze_radproffix
                  do ig=1,ngrid
                     do l=1,nlayer
                        reffrad(ig,l,iaer)=radvec(ig,l)*1.e-9    !  nm => m
                     enddo
                  enddo
                  nueffrad(1:ngrid,1:nlayer,iaer) = 0.02 ! haze
             endif
       enddo

   end subroutine haze_reffrad_fix
!==================================================================


end module radii_mod
!==================================================================