Index: /LMDZ4/branches/LMDZ4-dev/libf/phylmd/readaerosol_interp.F90
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/phylmd/readaerosol_interp.F90	(revision 1215)
+++ /LMDZ4/branches/LMDZ4-dev/libf/phylmd/readaerosol_interp.F90	(revision 1216)
@@ -1,5 +1,5 @@
 ! $Id$
 !
-SUBROUTINE readaerosol_interp(id_aero, r_day, first, pplay, paprs, mass_out, pi_mass_out)
+SUBROUTINE readaerosol_interp(id_aero, r_day, first, pplay, paprs, t_seri, mass_out, pi_mass_out)
 !
 ! This routine will return the mass concentration at actual day(mass_out) and 
@@ -36,4 +36,5 @@
   REAL, DIMENSION(klon,klev), INTENT(IN) :: pplay  ! pression at model mid-layers
   REAL, DIMENSION(klon,klev+1),INTENT(IN):: paprs  ! pression between model layers
+  REAL, DIMENSION(klon,klev), INTENT(IN) :: t_seri ! air temperature
 !      
 ! Output:      
@@ -51,4 +52,6 @@
 !$OMP THREADPRIVATE(klev_src)
 
+  REAL                              :: zrho      ! Air density [kg/m3]
+  REAL                              :: volm      ! Volyme de melange [kg/kg]
   REAL, DIMENSION(klon)             :: psurf_day, pi_psurf_day
   REAL, DIMENSION(klon)             :: load_src, pi_load_src  ! Mass load at source grid
@@ -149,4 +152,6 @@
     
      IF (debug) THEN
+        CALL writefield_phy('var_year_jan',var_year(:,:,1,id_aero),klev_src)
+        CALL writefield_phy('var_year_dec',var_year(:,:,12,id_aero),klev_src)
         CALL writefield_phy('psurf_src',psurf_year(:,:,id_aero),1)
         CALL writefield_phy('pi_psurf_src',pi_psurf_year(:,:,id_aero),1)
@@ -300,5 +305,7 @@
         DO k= 1, klev
            DO i = 1, klon
-              load_tgt(i) = load_tgt(i) + 1/RG * var_day(i,k,id_aero)*delp(i,k)
+              zrho = pplay(i,k)/t_seri(i,k)/RD       ! [kg/m3]
+              volm = var_day(i,k,id_aero)*1.E-9/zrho ! [kg/kg]
+              load_tgt(i) = load_tgt(i) + 1/RG * volm *delp(i,k)
            END DO
         END DO
@@ -315,5 +322,7 @@
            DO k= 1, klev
               DO i = 1, klon
-                 load_tgt_test(i) = load_tgt_test(i) + 1/RG * var_day(i,k,id_aero)*delp(i,k)
+                 zrho = pplay(i,k)/t_seri(i,k)/RD       ! [kg/m3]
+                 volm = var_day(i,k,id_aero)*1.E-9/zrho ! [kg/kg]
+                 load_tgt_test(i) = load_tgt_test(i) + 1/RG * volm*delp(i,k)
               END DO
            END DO
@@ -353,5 +362,7 @@
         DO k = 1, klev
            DO i = 1, klon
-              load_tgt(i) = load_tgt(i) + 1/RG * pi_var_day(i,k,id_aero)*delp(i,k)
+              zrho = pplay(i,k)/t_seri(i,k)/RD          ! [kg/m3]
+              volm = pi_var_day(i,k,id_aero)*1.E-9/zrho ! [kg/kg]
+              load_tgt(i) = load_tgt(i) + 1/RG * volm * delp(i,k)
            END DO
         END DO
@@ -367,5 +378,7 @@
            DO k = 1, klev
               DO i = 1, klon
-                 load_tgt_test(i) = load_tgt_test(i) + 1/RG * pi_var_day(i,k,id_aero)*delp(i,k)
+                 zrho = pplay(i,k)/t_seri(i,k)/RD          ! [kg/m3]
+                 volm = pi_var_day(i,k,id_aero)*1.E-9/zrho ! [kg/kg]
+                 load_tgt_test(i) = load_tgt_test(i) + 1/RG * volm * delp(i,k)
               END DO
            END DO
Index: /LMDZ4/branches/LMDZ4-dev/libf/phylmd/readaerosol_optic.F90
===================================================================
--- /LMDZ4/branches/LMDZ4-dev/libf/phylmd/readaerosol_optic.F90	(revision 1215)
+++ /LMDZ4/branches/LMDZ4-dev/libf/phylmd/readaerosol_optic.F90	(revision 1216)
@@ -73,5 +73,5 @@
        flag_aerosol .EQ. 6 ) THEN 
 
-     CALL readaerosol_interp(5, rjourvrai, debut, pplay, paprs, sulfate, sulfate_pi)
+     CALL readaerosol_interp(5, rjourvrai, debut, pplay, paprs, t_seri, sulfate, sulfate_pi)
   ELSE
      sulfate(:,:) = 0. ; sulfate_pi(:,:) = 0.
@@ -83,6 +83,6 @@
 
      ! Get bc aerosol distribution 
-     CALL readaerosol_interp(3, rjourvrai, debut, pplay, paprs, bcsol, bcsol_pi )
-     CALL readaerosol_interp(7, rjourvrai, debut, pplay, paprs, bcins, bcins_pi )
+     CALL readaerosol_interp(3, rjourvrai, debut, pplay, paprs, t_seri, bcsol, bcsol_pi )
+     CALL readaerosol_interp(7, rjourvrai, debut, pplay, paprs, t_seri, bcins, bcins_pi )
   ELSE
      bcsol(:,:) = 0. ; bcsol_pi(:,:) = 0.
@@ -95,6 +95,6 @@
        flag_aerosol .EQ. 6 ) THEN
 
-     CALL readaerosol_interp(4, rjourvrai, debut, pplay, paprs, pomsol, pomsol_pi)
-     CALL readaerosol_interp(8, rjourvrai, debut, pplay, paprs, pomins, pomins_pi)
+     CALL readaerosol_interp(4, rjourvrai, debut, pplay, paprs, t_seri, pomsol, pomsol_pi)
+     CALL readaerosol_interp(8, rjourvrai, debut, pplay, paprs, t_seri, pomins, pomins_pi)
   ELSE
      pomsol(:,:) = 0. ; pomsol_pi(:,:) = 0.
@@ -107,6 +107,6 @@
       flag_aerosol .EQ. 6 ) THEN 
 
-      CALL readaerosol_interp(1 ,rjourvrai, debut, pplay, paprs, sscoarse,sscoarse_pi) 
-      CALL readaerosol_interp(2 ,rjourvrai, debut, pplay, paprs, ssacu, ssacu_pi) 
+      CALL readaerosol_interp(1 ,rjourvrai, debut, pplay, paprs, t_seri, sscoarse,sscoarse_pi) 
+      CALL readaerosol_interp(2 ,rjourvrai, debut, pplay, paprs, t_seri, ssacu, ssacu_pi) 
 
   ELSE
@@ -119,5 +119,5 @@
       flag_aerosol .EQ. 6 ) THEN 
 
-      CALL readaerosol_interp(6 , rjourvrai, debut, pplay, paprs, cidust, cidust_pi) 
+      CALL readaerosol_interp(6 , rjourvrai, debut, pplay, paprs, t_seri, cidust, cidust_pi) 
 
   ELSE
