Index: LMDZ5/branches/testing/libf/phylmd/rrtm/readaerosol_optic_rrtm.F90
===================================================================
--- LMDZ5/branches/testing/libf/phylmd/rrtm/readaerosol_optic_rrtm.F90	(revision 2056)
+++ LMDZ5/branches/testing/libf/phylmd/rrtm/readaerosol_optic_rrtm.F90	(revision 2160)
@@ -1,26 +1,33 @@
 ! $Id$
 !
-SUBROUTINE readaerosol_optic_rrtm(debut, new_aod, flag_aerosol, itap, rjourvrai, &
+SUBROUTINE readaerosol_optic_rrtm(debut, aerosol_couple,  &
+     new_aod, flag_aerosol, itap, rjourvrai, &
      pdtphys, pplay, paprs, t_seri, rhcl, presnivs, &
-     mass_solu_aero, mass_solu_aero_pi, &
+     tr_seri, mass_solu_aero, mass_solu_aero_pi, &
      tau_aero, piz_aero, cg_aero, &
      tausum_aero, tau3d_aero )
 
-! This routine will :
-! 1) recevie the aerosols(already read and interpolated) corresponding to flag_aerosol
-! 2) calculate the optical properties for the aerosols
-!
-  
+  ! This routine will :
+  ! 1) recevie the aerosols(already read and interpolated) corresponding to flag_aerosol
+  ! 2) calculate the optical properties for the aerosols
+  !
+
   USE dimphy
   USE aero_mod
-  USE phys_local_var_mod, only: sconcso4,sconcoa,sconcbc,sconcss,sconcdust, &
-      concso4,concoa,concbc,concss,concdust,loadso4,loadoa,loadbc,loadss,loaddust, &
-      load_tmp1,load_tmp2,load_tmp3,load_tmp4,load_tmp5,load_tmp6,load_tmp7
+  USE phys_local_var_mod, only: sconcso4,sconcno3,sconcoa,sconcbc,sconcss,sconcdust, &
+       concso4,concno3,concoa,concbc,concss,concdust,loadso4,loadoa,loadbc,loadss,loaddust, &
+       load_tmp1,load_tmp2,load_tmp3,load_tmp4,load_tmp5,load_tmp6,load_tmp7
+
+  USE infotrac
+
   IMPLICIT NONE
   include "clesphys.h"
-
-! Input arguments
-!****************************************************************************************
+  include "YOMCST.h"
+
+
+  ! Input arguments
+  !****************************************************************************************
   LOGICAL, INTENT(IN)                      :: debut
+  LOGICAL, INTENT(IN)                      :: aerosol_couple
   LOGICAL, INTENT(IN)                      :: new_aod
   INTEGER, INTENT(IN)                      :: flag_aerosol
@@ -33,7 +40,8 @@
   REAL, DIMENSION(klon,klev), INTENT(IN)   :: rhcl   ! humidite relative ciel clair
   REAL, DIMENSION(klev), INTENT(IN)        :: presnivs
-
-! Output arguments
-!****************************************************************************************
+  REAL, DIMENSION(klon,klev,nbtr), INTENT(IN) :: tr_seri ! concentration tracer
+
+  ! Output arguments
+  !****************************************************************************************
   REAL, DIMENSION(klon,klev), INTENT(OUT)     :: mass_solu_aero    ! Total mass for all soluble aerosols
   REAL, DIMENSION(klon,klev), INTENT(OUT)     :: mass_solu_aero_pi !     -"-     preindustrial values
@@ -41,23 +49,25 @@
   REAL, DIMENSION(klon,klev,2,NSW), INTENT(OUT) :: piz_aero    ! Single scattering albedo aerosol
   REAL, DIMENSION(klon,klev,2,NSW), INTENT(OUT) :: cg_aero     ! asymmetry parameter aerosol
-!  REAL, DIMENSION(klon,nwave,naero_spc), INTENT(OUT)       :: tausum_aero
-!  REAL, DIMENSION(klon,klev,nwave,naero_spc), INTENT(OUT)  :: tau3d_aero
-!--correction minibug OB
   REAL, DIMENSION(klon,nwave,naero_tot), INTENT(OUT)       :: tausum_aero
   REAL, DIMENSION(klon,klev,nwave,naero_tot), INTENT(OUT)  :: tau3d_aero
 
-! Local variables
-!****************************************************************************************
-  REAL, DIMENSION(klon)        :: aerindex ! POLDER aerosol index 
-  REAL, DIMENSION(klon,klev)   :: sulfate  ! SO4 aerosol concentration [ug/m3]
-  REAL, DIMENSION(klon,klev)   :: bcsol    ! BC soluble concentration [ug/m3]
-  REAL, DIMENSION(klon,klev)   :: bcins    ! BC insoluble concentration [ug/m3]
-  REAL, DIMENSION(klon,klev)   :: pomsol   ! POM soluble concentration [ug/m3]
-  REAL, DIMENSION(klon,klev)   :: pomins   ! POM insoluble concentration [ug/m3]
-  REAL, DIMENSION(klon,klev)   :: cidust    ! DUST aerosol concentration  [ug/m3]
-  REAL, DIMENSION(klon,klev)   :: sscoarse  ! SS Coarse concentration [ug/m3]
-  REAL, DIMENSION(klon,klev)   :: sssupco   ! SS Super Coarse concentration [ug/m3]
-  REAL, DIMENSION(klon,klev)   :: ssacu     ! SS Acumulation concentration [ug/m3]
-  REAL, DIMENSION(klon,klev)   :: sulfate_pi
+  ! Local variables
+  !****************************************************************************************
+  REAL, DIMENSION(klon)        :: aerindex      ! POLDER aerosol index 
+  REAL, DIMENSION(klon,klev)   :: sulfacc       ! SO4 accumulation concentration [ug/m3]
+  REAL, DIMENSION(klon,klev)   :: sulfcoarse    ! SO4 coarse concentration [ug/m3]
+  REAL, DIMENSION(klon,klev)   :: bcsol         ! BC soluble concentration [ug/m3]
+  REAL, DIMENSION(klon,klev)   :: bcins         ! BC insoluble concentration [ug/m3]
+  REAL, DIMENSION(klon,klev)   :: pomsol        ! POM soluble concentration [ug/m3]
+  REAL, DIMENSION(klon,klev)   :: pomins        ! POM insoluble concentration [ug/m3]
+  REAL, DIMENSION(klon,klev)   :: cidust        ! DUST aerosol concentration  [ug/m3]
+  REAL, DIMENSION(klon,klev)   :: sscoarse      ! SS Coarse concentration [ug/m3]
+  REAL, DIMENSION(klon,klev)   :: sssupco       ! SS Super Coarse concentration [ug/m3]
+  REAL, DIMENSION(klon,klev)   :: ssacu         ! SS Acumulation concentration [ug/m3]
+  REAL, DIMENSION(klon,klev)   :: nitracc       ! nitrate accumulation concentration [ug/m3]
+  REAL, DIMENSION(klon,klev)   :: nitrcoarse    ! nitrate coarse concentration [ug/m3]
+  REAL, DIMENSION(klon,klev)   :: nitrinscoarse ! nitrate insoluble coarse concentration [ug/m3]
+  REAL, DIMENSION(klon,klev)   :: sulfacc_pi
+  REAL, DIMENSION(klon,klev)   :: sulfcoarse_pi
   REAL, DIMENSION(klon,klev)   :: bcsol_pi
   REAL, DIMENSION(klon,klev)   :: bcins_pi
@@ -68,118 +78,228 @@
   REAL, DIMENSION(klon,klev)   :: sssupco_pi
   REAL, DIMENSION(klon,klev)   :: ssacu_pi
-  REAL, DIMENSION(klon,klev)   :: pdel
-  REAL, DIMENSION(klon,klev,naero_spc) :: m_allaer
-  REAL, DIMENSION(klon,klev,naero_spc) :: m_allaer_pi !RAF  
-!  REAL, DIMENSION(klon,naero_tot)      :: fractnat_allaer !RAF delete??
-
+  REAL, DIMENSION(klon,klev)   :: nitracc_pi
+  REAL, DIMENSION(klon,klev)   :: nitrcoarse_pi
+  REAL, DIMENSION(klon,klev)   :: nitrinscoarse_pi
+  REAL, DIMENSION(klon,klev)   :: pdel, zrho
+  REAL, DIMENSION(klon,klev,naero_tot) :: m_allaer
+  REAL, DIMENSION(klon,klev,naero_tot) :: m_allaer_pi !RAF  
+  !  REAL, DIMENSION(klon,naero_tot)      :: fractnat_allaer !RAF delete??
+  character(len=8), dimension(nbtr) :: tracname
+  integer :: id_ASBCM, id_ASPOMM, id_ASSO4M, id_ASMSAM, id_CSSO4M, id_CSMSAM, id_SSSSM
+  integer :: id_CSSSM, id_ASSSM, id_CIDUSTM, id_AIBCM, id_AIPOMM, id_ASNO3M, id_CSNO3M, id_CINO3M
   INTEGER :: k, i
-  
-!****************************************************************************************
-! 1) Get aerosol mass
-!    
-!****************************************************************************************
-! Read and interpolate sulfate
-  IF ( flag_aerosol .EQ. 1 .OR. &
-       flag_aerosol .EQ. 6 ) THEN 
-
-     CALL readaerosol_interp(id_ASSO4M, itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, sulfate, sulfate_pi,loadso4)
-  ELSE
-     sulfate(:,:) = 0. ; sulfate_pi(:,:) = 0.
-     loadso4=0.
-  END IF
-
-! Read and interpolate bcsol and bcins
-  IF ( flag_aerosol .EQ. 2 .OR. &
-       flag_aerosol .EQ. 6 ) THEN 
-
-     ! Get bc aerosol distribution 
-     CALL readaerosol_interp(id_ASBCM, itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, bcsol, bcsol_pi, load_tmp1 )
-     CALL readaerosol_interp(id_AIBCM, itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, bcins, bcins_pi, load_tmp2 )
-     loadbc(:)=load_tmp1(:)+load_tmp2(:)
-  ELSE
-     bcsol(:,:) = 0. ; bcsol_pi(:,:) = 0.
-     bcins(:,:) = 0. ; bcins_pi(:,:) = 0.
-     loadbc=0.
-  END IF
-
-
-! Read and interpolate pomsol and pomins
-  IF ( flag_aerosol .EQ. 3 .OR. &
-       flag_aerosol .EQ. 6 ) THEN
-
-     CALL readaerosol_interp(id_ASPOMM, itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, pomsol, pomsol_pi, load_tmp3)
-     CALL readaerosol_interp(id_AIPOMM, itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, pomins, pomins_pi, load_tmp4)
-     loadoa(:)=load_tmp3(:)+load_tmp4(:)
-  ELSE
-     pomsol(:,:) = 0. ; pomsol_pi(:,:) = 0.
-     pomins(:,:) = 0. ; pomins_pi(:,:) = 0.
-     loadoa=0.
-  END IF
-
-
-! Read and interpolate csssm, ssssm, assssm
-  IF (flag_aerosol .EQ. 4 .OR. &
-      flag_aerosol .EQ. 6 ) THEN 
-
-      CALL readaerosol_interp(id_SSSSM ,itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, sssupco, sssupco_pi, load_tmp5) 
-      CALL readaerosol_interp(id_CSSSM ,itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, sscoarse,sscoarse_pi, load_tmp6) 
-      CALL readaerosol_interp(id_ASSSM ,itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, ssacu, ssacu_pi, load_tmp7) 
-     loadss(:)=load_tmp5(:)+load_tmp6(:)+load_tmp7(:)
-  ELSE
-     sscoarse(:,:) = 0. ; sscoarse_pi(:,:) = 0. 
-     ssacu(:,:)    = 0. ; ssacu_pi(:,:) = 0. 
-     sssupco(:,:)  = 0. ; sssupco_pi = 0. 
-     loadss=0.
-  ENDIF
-
-! Read and interpolate cidustm
-  IF (flag_aerosol .EQ. 5 .OR.  &
-      flag_aerosol .EQ. 6 ) THEN 
-
-      CALL readaerosol_interp(id_CIDUSTM, itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, cidust, cidust_pi, loaddust) 
-
-  ELSE
-      cidust(:,:) = 0. ; cidust_pi(:,:) = 0. 
-      loaddust=0.
-  ENDIF
-
-!
-! Store all aerosols in one variable
-!
-  m_allaer(:,:,id_ASBCM)  = bcsol(:,:)        ! ASBCM
-  m_allaer(:,:,id_ASPOMM) = pomsol(:,:)       ! ASPOMM
-  m_allaer(:,:,id_ASSO4M) = sulfate(:,:)      ! ASSO4M (= SO4) 
-  m_allaer(:,:,id_CSSO4M) = 0.                ! CSSO4M 
-  m_allaer(:,:,id_SSSSM)  = sssupco(:,:)      ! SSSSM
-  m_allaer(:,:,id_CSSSM)  = sscoarse(:,:)     ! CSSSM
-  m_allaer(:,:,id_ASSSM)  = ssacu(:,:)        ! ASSSM
-  m_allaer(:,:,id_CIDUSTM)= cidust(:,:)       ! CIDUSTM
-  m_allaer(:,:,id_AIBCM)  = bcins(:,:)        ! AIBCM
-  m_allaer(:,:,id_AIPOMM) = pomins(:,:)       ! AIPOMM
-
-!RAF
-  m_allaer_pi(:,:,1)  = bcsol_pi(:,:)        ! ASBCM pre-ind
-  m_allaer_pi(:,:,2)  = pomsol_pi(:,:)       ! ASPOMM pre-ind
-  m_allaer_pi(:,:,3)  = sulfate_pi(:,:)      ! ASSO4M (= SO4) pre-ind
-  m_allaer_pi(:,:,4)  = 0.                   ! CSSO4M pre-ind
-  m_allaer_pi(:,:,5)  = sssupco_pi(:,:)      ! SSSSM pre-ind
-  m_allaer_pi(:,:,6)  = sscoarse_pi(:,:)     ! CSSSM pre-ind
-  m_allaer_pi(:,:,7)  = ssacu_pi(:,:)        ! ASSSM pre-ind
-  m_allaer_pi(:,:,8)  = cidust_pi(:,:)       ! CIDUSTM pre-ind
-  m_allaer_pi(:,:,9)  = bcins_pi(:,:)        ! AIBCM pre-ind
-  m_allaer_pi(:,:,10) = pomins_pi(:,:)       ! AIPOMM pre-ind
-
-!
-! Calculate the total mass of all soluble aersosols
-!
-  mass_solu_aero(:,:)    = sulfate(:,:)    + bcsol(:,:)    + pomsol(:,:) !   + &
-!       sscoarse(:,:)    + ssacu(:,:)    + sssupco(:,:) 
-  mass_solu_aero_pi(:,:) = sulfate_pi(:,:) + bcsol_pi(:,:) + pomsol_pi(:,:) ! + &
-!       sscoarse_pi(:,:) + ssacu_pi(:,:) + sssupco_pi(:,:)
-
-!****************************************************************************************
-! 2) Calculate optical properties for the aerosols
-!
-!****************************************************************************************
+
+  !--air density
+  zrho(:,:)=pplay(:,:)/t_seri(:,:)/RD                     !--kg/m3
+
+  !****************************************************************************************
+  ! 1) Get aerosol mass
+  !    
+  !****************************************************************************************
+  !
+  !
+  IF (aerosol_couple) THEN   !--we get aerosols from tr_seri array from INCA
+     !
+     !--copy fields from INCA tr_seri 
+     !--convert to ug m-3 unit for consistency with offline fields
+     !
+#ifdef INCA
+     call tracinca_name(tracname)
+#endif
+
+     do i=1,nbtr
+        select case(trim(tracname(i)))
+           case ("ASBCM")
+              id_ASBCM = i
+           case ("ASPOMM") 
+              id_ASPOMM = i 
+           case ("ASSO4M")
+              id_ASSO4M = i 
+           case ("ASMSAM")
+              id_ASMSAM = i 
+           case ("CSSO4M")
+              id_CSSO4M = i 
+           case ("CSMSAM")
+              id_CSMSAM = i 
+           case ("SSSSM")
+              id_SSSSM = i 
+           case ("CSSSM")
+              id_CSSSM = i 
+           case ("ASSSM")
+              id_ASSSM = i 
+           case ("CIDUSTM")
+              id_CIDUSTM = i 
+           case ("AIBCM")
+              id_AIBCM = i 
+           case ("AIPOMM")
+              id_AIPOMM = i 
+           case ("ASNO3M")
+              id_ASNO3M = i 
+           case ("CSNO3M")
+              id_CSNO3M = i 
+           case ("CINO3M")
+              id_CINO3M = i 
+           end select
+     enddo
+
+
+     bcsol(:,:)        =   tr_seri(:,:,id_ASBCM)                         *zrho(:,:)*1.e9  ! ASBCM
+     pomsol(:,:)       =   tr_seri(:,:,id_ASPOMM)                        *zrho(:,:)*1.e9  ! ASPOMM
+     sulfacc(:,:)      =  (tr_seri(:,:,id_ASSO4M)+tr_seri(:,:,id_ASMSAM))*zrho(:,:)*1.e9  ! ASSO4M (=SO4) + ASMSAM (=MSA)
+     sulfcoarse(:,:)   =  (tr_seri(:,:,id_CSSO4M)+tr_seri(:,:,id_CSMSAM))*zrho(:,:)*1.e9  ! CSSO4M (=SO4) + CSMSAM (=MSA)
+     sssupco(:,:)      =   tr_seri(:,:,id_SSSSM)                         *zrho(:,:)*1.e9  ! SSSSM
+     sscoarse(:,:)     =   tr_seri(:,:,id_CSSSM)                         *zrho(:,:)*1.e9  ! CSSSM
+     ssacu(:,:)        =   tr_seri(:,:,id_ASSSM)                         *zrho(:,:)*1.e9  ! ASSSM
+     cidust(:,:)       =   tr_seri(:,:,id_CIDUSTM)                       *zrho(:,:)*1.e9  ! CIDUSTM
+     bcins(:,:)        =   tr_seri(:,:,id_AIBCM)                         *zrho(:,:)*1.e9  ! AIBCM
+     pomins(:,:)       =   tr_seri(:,:,id_AIPOMM)                        *zrho(:,:)*1.e9  ! AIPOMM
+     nitracc(:,:)      =   tr_seri(:,:,id_ASNO3M)                        *zrho(:,:)*1.e9  ! ASNO3M
+     nitrcoarse(:,:)   =   tr_seri(:,:,id_CSNO3M)                        *zrho(:,:)*1.e9  ! CSNO3M
+     nitrinscoarse(:,:)=   tr_seri(:,:,id_CINO3M)                        *zrho(:,:)*1.e9  ! CINO3M
+     !
+     bcsol_pi(:,:)        =   0.0 ! ASBCM pre-ind
+     pomsol_pi(:,:)       =   0.0 ! ASPOMM pre-ind
+     sulfacc_pi(:,:)      =   0.0 ! ASSO4M (=SO4) + ASMSAM (=MSA) pre-ind
+     sulfcoarse_pi(:,:)   =   0.0 ! CSSO4M (=SO4) + CSMSAM (=MSA) pre-ind
+     sssupco_pi(:,:)      =   0.0 ! SSSSM pre-ind
+     sscoarse_pi(:,:)     =   0.0 ! CSSSM pre-ind
+     ssacu_pi(:,:)        =   0.0 ! ASSSM pre-ind
+     cidust_pi(:,:)       =   0.0 ! CIDUSTM pre-ind
+     bcins_pi(:,:)        =   0.0 ! AIBCM pre-ind
+     pomins_pi(:,:)       =   0.0 ! AIPOMM pre-ind
+     nitracc_pi(:,:)      =   0.0 ! ASNO3M pre-ind
+     nitrcoarse_pi(:,:)   =   0.0 ! CSNO3M pre-ind
+     nitrinscoarse_pi(:,:)=   0.0 ! CINO3M
+     !
+  ELSE !--not aerosol_couple
+     !
+     ! Read and interpolate sulfate
+     IF ( flag_aerosol .EQ. 1 .OR. &
+          flag_aerosol .EQ. 6 ) THEN 
+
+        CALL readaerosol_interp(id_ASSO4M_phy, itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, sulfacc, sulfacc_pi,loadso4)
+     ELSE
+        sulfacc(:,:) = 0. ; sulfacc_pi(:,:) = 0.
+        loadso4=0.
+     END IF
+
+     ! Read and interpolate bcsol and bcins
+     IF ( flag_aerosol .EQ. 2 .OR. &
+          flag_aerosol .EQ. 6 ) THEN 
+
+        ! Get bc aerosol distribution 
+        CALL readaerosol_interp(id_ASBCM_phy, itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, bcsol, bcsol_pi, load_tmp1 )
+        CALL readaerosol_interp(id_AIBCM_phy, itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, bcins, bcins_pi, load_tmp2 )
+        loadbc(:)=load_tmp1(:)+load_tmp2(:)
+     ELSE
+        bcsol(:,:) = 0. ; bcsol_pi(:,:) = 0.
+        bcins(:,:) = 0. ; bcins_pi(:,:) = 0.
+        loadbc=0.
+     END IF
+
+
+     ! Read and interpolate pomsol and pomins
+     IF ( flag_aerosol .EQ. 3 .OR. &
+          flag_aerosol .EQ. 6 ) THEN
+
+        CALL readaerosol_interp(id_ASPOMM_phy, itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, pomsol, pomsol_pi, load_tmp3)
+        CALL readaerosol_interp(id_AIPOMM_phy, itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, pomins, pomins_pi, load_tmp4)
+        loadoa(:)=load_tmp3(:)+load_tmp4(:)
+     ELSE
+        pomsol(:,:) = 0. ; pomsol_pi(:,:) = 0.
+        pomins(:,:) = 0. ; pomins_pi(:,:) = 0.
+        loadoa=0.
+     END IF
+
+
+     ! Read and interpolate csssm, ssssm, assssm
+     IF (flag_aerosol .EQ. 4 .OR. &
+          flag_aerosol .EQ. 6 ) THEN 
+
+        CALL readaerosol_interp(id_SSSSM_phy ,itap, pdtphys,rjourvrai, &
+        debut, pplay, paprs, t_seri, sssupco, sssupco_pi, load_tmp5) 
+        CALL readaerosol_interp(id_CSSSM_phy ,itap, pdtphys,rjourvrai, &
+        debut, pplay, paprs, t_seri, sscoarse,sscoarse_pi, load_tmp6) 
+        CALL readaerosol_interp(id_ASSSM_phy ,itap, pdtphys, rjourvrai, &
+        debut, pplay, paprs, t_seri, ssacu, ssacu_pi, load_tmp7) 
+        loadss(:)=load_tmp5(:)+load_tmp6(:)+load_tmp7(:)
+     ELSE
+        sscoarse(:,:) = 0. ; sscoarse_pi(:,:) = 0. 
+        ssacu(:,:)    = 0. ; ssacu_pi(:,:) = 0. 
+        sssupco(:,:)  = 0. ; sssupco_pi = 0. 
+        loadss=0.
+     ENDIF
+
+     ! Read and interpolate cidustm
+     IF (flag_aerosol .EQ. 5 .OR.  &
+          flag_aerosol .EQ. 6 ) THEN 
+
+        CALL readaerosol_interp(id_CIDUSTM_phy, itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, cidust, cidust_pi, loaddust) 
+
+     ELSE
+        cidust(:,:) = 0. ; cidust_pi(:,:) = 0. 
+        loaddust=0.
+     ENDIF
+     !
+     sulfcoarse(:,:)      =   0.0 ! CSSO4M (=SO4) + CSMSAM (=MSA) 
+     sulfcoarse_pi(:,:)   =   0.0 ! CSSO4M (=SO4) + CSMSAM (=MSA) pre-ind
+     !
+     !--placeholder for offline nitrate   
+     !
+     nitracc(:,:)         =   0.0 
+     nitracc_pi(:,:)      =   0.0 
+     nitrcoarse(:,:)      =   0.0
+     nitrcoarse_pi(:,:)   =   0.0
+     nitrinscoarse(:,:)   =   0.0
+     nitrinscoarse_pi(:,:)=   0.0
+
+  ENDIF !--not aerosol_couple
+
+  !
+  ! Store all aerosols in one variable
+  !
+  m_allaer(:,:,id_ASBCM_phy)  = bcsol(:,:)        ! ASBCM
+  m_allaer(:,:,id_ASPOMM_phy) = pomsol(:,:)       ! ASPOMM
+  m_allaer(:,:,id_ASSO4M_phy) = sulfacc(:,:)      ! ASSO4M (= SO4) 
+  m_allaer(:,:,id_CSSO4M_phy) = sulfcoarse(:,:)   ! CSSO4M 
+  m_allaer(:,:,id_SSSSM_phy)  = sssupco(:,:)      ! SSSSM
+  m_allaer(:,:,id_CSSSM_phy)  = sscoarse(:,:)     ! CSSSM
+  m_allaer(:,:,id_ASSSM_phy)  = ssacu(:,:)        ! ASSSM
+  m_allaer(:,:,id_CIDUSTM_phy)= cidust(:,:)       ! CIDUSTM
+  m_allaer(:,:,id_AIBCM_phy)  = bcins(:,:)        ! AIBCM
+  m_allaer(:,:,id_ASNO3M_phy) = nitracc(:,:)      ! ASNO3M
+  m_allaer(:,:,id_CSNO3M_phy) = nitrcoarse(:,:)   ! CSNO3M
+  m_allaer(:,:,id_CINO3M_phy) = nitrinscoarse(:,:)! CINO3M
+  m_allaer(:,:,id_AIPOMM_phy) = pomins(:,:)       ! AIPOMM
+  m_allaer(:,:,id_STRAT_phy)  = 0.0
+
+  !RAF
+  m_allaer_pi(:,:,id_ASBCM_phy)  = bcsol_pi(:,:)        ! ASBCM pre-ind
+  m_allaer_pi(:,:,id_ASPOMM_phy) = pomsol_pi(:,:)       ! ASPOMM pre-ind
+  m_allaer_pi(:,:,id_ASSO4M_phy) = sulfacc_pi(:,:)      ! ASSO4M (= SO4) pre-ind
+  m_allaer_pi(:,:,id_CSSO4M_phy) = sulfcoarse_pi(:,:)   ! CSSO4M pre-ind
+  m_allaer_pi(:,:,id_SSSSM_phy)  = sssupco_pi(:,:)      ! SSSSM pre-ind
+  m_allaer_pi(:,:,id_CSSSM_phy)  = sscoarse_pi(:,:)     ! CSSSM pre-ind
+  m_allaer_pi(:,:,id_ASSSM_phy)  = ssacu_pi(:,:)        ! ASSSM pre-ind
+  m_allaer_pi(:,:,id_CIDUSTM_phy)= cidust_pi(:,:)       ! CIDUSTM pre-ind
+  m_allaer_pi(:,:,id_AIBCM_phy)  = bcins_pi(:,:)        ! AIBCM pre-ind
+  m_allaer_pi(:,:,id_ASNO3M_phy) = nitracc_pi(:,:)      ! ASNO3M pre-ind
+  m_allaer_pi(:,:,id_CSNO3M_phy) = nitrcoarse_pi(:,:)   ! CSNO3M pre-ind
+  m_allaer_pi(:,:,id_CINO3M_phy) = nitrinscoarse_pi(:,:)! CINO3M pre-ind
+  m_allaer_pi(:,:,id_AIPOMM_phy) = pomins_pi(:,:)       ! AIPOMM pre-ind
+  m_allaer_pi(:,:,id_STRAT_phy)  = 0.0
+
+  !
+  ! Calculate the total mass of all soluble aersosols
+  ! to be revisited for AR6
+  mass_solu_aero(:,:)    = sulfacc(:,:)    + bcsol(:,:)    + pomsol(:,:)   + nitracc(:,:)   ! + & 
+  !                           sscoarse(:,:)    + ssacu(:,:)    + sssupco(:,:) 
+  mass_solu_aero_pi(:,:) = sulfacc_pi(:,:) + bcsol_pi(:,:) + pomsol_pi(:,:) + nitracc_pi(:,:) ! + &
+  !                           sscoarse_pi(:,:) + ssacu_pi(:,:) + sssupco_pi(:,:)
+
+  !****************************************************************************************
+  ! 2) Calculate optical properties for the aerosols
+  !
+  !****************************************************************************************
   DO k = 1, klev
      DO i = 1, klon
@@ -188,28 +308,31 @@
   END DO
 
+  ! aeropt_6bands for rrtm
   CALL aeropt_6bands_rrtm( &
-          pdel, m_allaer, pdtphys, rhcl, & 
-          tau_aero, piz_aero, cg_aero,   &
-          m_allaer_pi, flag_aerosol, &
-          pplay, t_seri, presnivs) 
-     
-     ! aeropt_5wv only for validation and diagnostics.
-     CALL aeropt_5wv_rrtm(                    &
-          pdel, m_allaer,                &
-          pdtphys, rhcl, aerindex,       & 
-          flag_aerosol, pplay, t_seri,   &
-          tausum_aero, tau3d_aero )
-
-! Diagnostics calculation for CMIP5 protocol
-  sconcso4(:)=m_allaer(:,1,id_ASSO4M)*1.e-9
-  sconcoa(:)=(m_allaer(:,1,id_ASPOMM)+m_allaer(:,1,id_AIPOMM))*1.e-9
-  sconcbc(:)=(m_allaer(:,1,id_ASBCM)+m_allaer(:,1,id_AIBCM))*1.e-9
-  sconcss(:)=(m_allaer(:,1,id_ASSSM)+m_allaer(:,1,id_CSSSM)+m_allaer(:,1,id_SSSSM))*1.e-9
-  sconcdust(:)=m_allaer(:,1,id_CIDUSTM)*1.e-9
-  concso4(:,:)=m_allaer(:,:,id_ASSO4M)*1.e-9
-  concoa(:,:)=(m_allaer(:,:,id_ASPOMM)+m_allaer(:,:,id_AIPOMM))*1.e-9
-  concbc(:,:)=(m_allaer(:,:,id_ASBCM)+m_allaer(:,:,id_AIBCM))*1.e-9
-  concss(:,:)=(m_allaer(:,:,id_ASSSM)+m_allaer(:,:,id_CSSSM)+m_allaer(:,:,id_SSSSM))*1.e-9
-  concdust(:,:)=m_allaer(:,:,id_CIDUSTM)*1.e-9
+       pdel, m_allaer, pdtphys, rhcl, & 
+       tau_aero, piz_aero, cg_aero,   &
+       m_allaer_pi, flag_aerosol, &
+       zrho ) 
+
+  ! aeropt_5wv only for validation and diagnostics
+  CALL aeropt_5wv_rrtm(                    &
+       pdel, m_allaer,                &
+       pdtphys, rhcl, aerindex,       & 
+       flag_aerosol, pplay, t_seri,   &
+       tausum_aero, tau3d_aero )
+
+  ! Diagnostics calculation for CMIP5 protocol
+  sconcso4(:)  =m_allaer(:,1,id_ASSO4M_phy)*1.e-9
+  sconcno3(:)  =(m_allaer(:,1,id_ASNO3M_phy)+m_allaer(:,1,id_CSNO3M_phy)+m_allaer(:,1,id_CINO3M_phy))*1.e-9
+  sconcoa(:)   =(m_allaer(:,1,id_ASPOMM_phy)+m_allaer(:,1,id_AIPOMM_phy))*1.e-9
+  sconcbc(:)   =(m_allaer(:,1,id_ASBCM_phy)+m_allaer(:,1,id_AIBCM_phy))*1.e-9
+  sconcss(:)   =(m_allaer(:,1,id_ASSSM_phy)+m_allaer(:,1,id_CSSSM_phy)+m_allaer(:,1,id_SSSSM_phy))*1.e-9
+  sconcdust(:) =m_allaer(:,1,id_CIDUSTM_phy)*1.e-9
+  concso4(:,:) =m_allaer(:,:,id_ASSO4M_phy)*1.e-9
+  concno3(:,:) =(m_allaer(:,:,id_ASNO3M_phy)+m_allaer(:,:,id_CSNO3M_phy)+m_allaer(:,:,id_CINO3M_phy))*1.e-9
+  concoa(:,:)  =(m_allaer(:,:,id_ASPOMM_phy)+m_allaer(:,:,id_AIPOMM_phy))*1.e-9
+  concbc(:,:)  =(m_allaer(:,:,id_ASBCM_phy)+m_allaer(:,:,id_AIBCM_phy))*1.e-9
+  concss(:,:)  =(m_allaer(:,:,id_ASSSM_phy)+m_allaer(:,:,id_CSSSM_phy)+m_allaer(:,:,id_SSSSM_phy))*1.e-9
+  concdust(:,:)=m_allaer(:,:,id_CIDUSTM_phy)*1.e-9
 
 END SUBROUTINE readaerosol_optic_rrtm
