source: LMDZ6/branches/IPSLCM6.0.15/libf/phylmd/macv2sp.F90 @ 5029

Last change on this file since 5029 was 3318, checked in by jghattas, 7 years ago

Integration des developpements fait sur la trunk par O. Boucher concernant le MACspV2 aerosol plume climatology(nouveau option flag_aersol=7). Les commits suivants fait sur le trunk sont ici merge : [3274], [3279], [3287], [3288], [3290], [3295], [3296], [3297].
Tout les modifications dans newmicro.f90 ne sont pas retenu mais les changemnets lie au flag_aerosol=7 sont prise.
File size: 7.1 KB
RevLine 
[3318]1SUBROUTINE MACv2SP(pphis,pplay,paprs,xlon,xlat,tau_allaer,piz_allaer,cg_allaer)
2  !
3  !--routine to read the MACv2SP plume and compute optical properties
4  !--requires flag_aerosol = 7
5  !--feeds into aerosol optical properties and newmicro cloud droplet size if ok_cdnc activated
6  !--for this one needs to feed natural (pre-industrial) aerosols twice for nat and 1980 files
7  !--pre-ind aerosols (index=1) are not changed, present-day aerosols (index=2) are incremented
8  !--uses model year so year_cur needs to be correct in the model simulation
9  !
10  !--aod_prof = AOD per layer
11  !--ssa_prof = SSA
12  !--asy_prof = asymetry parameter
13  !--dNovrN   = enhancement factor for CDNC
14  !
15  USE mo_simple_plumes, ONLY: sp_aop_profile
16  USE phys_cal_mod, ONLY : year_cur, day_cur, year_len
17  USE dimphy
18  USE aero_mod
19  USE phys_local_var_mod, ONLY: t_seri, od443aer, od550aer, od865aer, ec550aer, dryod550aer, od550lt1aer, dNovrN
20  !!USE YOMCST, ONLY : RD, RG
21  !
22  IMPLICIT NONE
23  !
24  include "YOMCST.h"
25  !
26  REAL,DIMENSION(klon),INTENT(IN)        :: pphis   ! Geopotentiel de surface
27  REAL,DIMENSION(klon,klev),INTENT(IN)   :: pplay   ! pression pour le mileu de chaque couche (en Pa)
28  REAL,DIMENSION(klon,klev+1),INTENT(IN) :: paprs   ! pression pour les interfaces de chaque couche (en Pa)
29  REAL,DIMENSION(klon),INTENT(IN)        :: xlat    ! latitudes pour chaque point
30  REAL,DIMENSION(klon),INTENT(IN)        :: xlon    ! longitudes pour chaque point
31  !
32  REAL, DIMENSION(klon,klev,2,nbands_sw_rrtm), INTENT(OUT) :: tau_allaer !  epaisseur optique aerosol
33  REAL, DIMENSION(klon,klev,2,nbands_sw_rrtm), INTENT(OUT) :: piz_allaer !  single scattering albedo aerosol
34  REAL, DIMENSION(klon,klev,2,nbands_sw_rrtm), INTENT(OUT) :: cg_allaer  !  asymmetry parameter aerosol
35  !
36  REAL,DIMENSION(klon,klev) :: aod_prof, ssa_prof, asy_prof
37  REAL,DIMENSION(klon,klev) :: z, dz
38  REAL,DIMENSION(klon)      :: oro, zrho, zt
39  !
40  INTEGER, PARAMETER :: nmon = 12
41  !
42  REAL, PARAMETER    :: l443 = 443.0, l550 = 550.0, l865 = 865.0 !--wavelengths in nm
43  !
44  INTEGER, PARAMETER :: Nwvmax=25
45  REAL, DIMENSION(0:Nwvmax), PARAMETER :: lambda=(/ 240.0, &  !--this one is for band 1
46                  280.0,  300.0,  330.0,  360.0,  400.0,   &  !--these are bounds of Streamer bands
47                  440.0,  480.0,  520.0,  570.0,  640.0,   &
48                  690.0,  750.0,  780.0,  870.0, 1000.0,   &
49                 1100.0, 1190.0, 1280.0, 1530.0, 1640.0,   &
50                 2130.0, 2380.0, 2910.0, 3420.0, 4000.0   /)
51  !
52  REAL, DIMENSION(1:Nwvmax-1), PARAMETER :: weight =(/    &   !--and the weights to be given to the bands
53                 0.01,  4.05,  9.51, 15.99, 26.07, 33.10, &   !--corresponding to a typical solar spectrum
54                33.07, 39.91, 52.67, 27.89, 43.60, 13.67, &
55                42.22, 40.12, 32.70, 14.44, 19.48, 14.23, &
56                13.43, 16.42,  8.33,  0.95,  0.65,  2.76  /)
57  !
58  REAL :: zlambda, zweight
59  REAL :: year_fr
60  !
61  INTEGER band, i, k, Nwv
62  !
63  ! define the height and dheight arrays
64  !
65  oro(:)  = pphis(:)/RG                             ! surface height in m
66  !
67  DO k = 1, klev
68    zrho(:) = pplay(:,k)/t_seri(:,k)/RD                         ! air density in kg/m3
69    dz(:,k) = (paprs(:,k)-paprs(:,k+1))/zrho(:)/RG              ! layer thickness in m
70    IF (k==1) THEN
71       z(:,1) = oro(:) + (paprs(:,1)-pplay(:,1))/zrho(:)/RG     ! altitude middle of first layer in m
72       zt(:)  = oro(:) + dz(:,1)                                ! altitude top of first layer in m
73    ELSE
74      z(:,k) = zt(:) + (paprs(:,k)-pplay(:,k))/zrho(:)/RG       ! altitude middle of layer k in m
75      zt(:)  = zt(:) + dz(:,k)                                  ! altitude top of layer k in m
76    ENDIF
77  ENDDO
78  !
79  !--fractional year
80  !
81  year_fr = FLOAT(year_cur) + (FLOAT(day_cur)-0.5) / FLOAT(year_len)
82  print *,'year_fr=',year_fr
83  !
84  !--call to sp routine -- 443 nm
85  !
86  CALL sp_aop_profile                                    ( &
87       klev     ,klon ,l443 ,oro    ,xlon     ,xlat      , &
88       year_fr  ,z    ,dz   ,dNovrN ,aod_prof ,ssa_prof  , &
89       asy_prof )
90  !
91  !--AOD calculations for diagnostics
92  od443aer(:)= od443aer(:)+SUM(aod_prof(:,:),dim=2)
93  !
94  !--call to sp routine -- 550 nm
95  !
96  CALL sp_aop_profile                                    ( &
97       klev     ,klon ,l550 ,oro    ,xlon     ,xlat      , &
98       year_fr  ,z    ,dz   ,dNovrN ,aod_prof ,ssa_prof  , &
99       asy_prof )
100  !
101  !--AOD calculations for diagnostics
102  od550aer(:)=od550aer(:)+SUM(aod_prof(:,:),dim=2)
103  !
104  !--dry AOD calculation for diagnostics
105  dryod550aer(:)=dryod550aer(:)+od550aer(:)
106  !
107  !--fine-mode AOD calculation for diagnostics
108  od550lt1aer(:)=od550lt1aer(:)+od550aer(:)
109  !
110  !--extinction coefficient for diagnostic
111  ec550aer(:,:)=ec550aer(:,:)+aod_prof(:,:)/dz(:,:)
112  !
113  !--call to sp routine -- 865 nm
114  !
115  CALL sp_aop_profile                                    ( &
116       klev     ,klon ,l865 ,oro    ,xlon     ,xlat      , &
117       year_fr  ,z    ,dz   ,dNovrN ,aod_prof ,ssa_prof  , &
118       asy_prof )
119  !
120  !--AOD calculations for diagnostics
121  od865aer(:)=od865aer(:)+SUM(aod_prof(:,:),dim=2)
122  !
123  !--re-weighting of piz and cg arrays before adding the anthropogenic aerosols
124  !--index 2 = all natural + anthropogenic aerosols
125  piz_allaer(:,:,2,:)=piz_allaer(:,:,2,:)*tau_allaer(:,:,2,:)
126  cg_allaer(:,:,2,:) =cg_allaer(:,:,2,:)*piz_allaer(:,:,2,:)
127  !
128  !--now computing the same at many wavelengths to fill the model bands
129  !
130  DO Nwv=0,Nwvmax-1
131
132    IF (Nwv.EQ.0) THEN          !--RRTM spectral band 1
133      zlambda=lambda(Nwv)
134      zweight=1.0
135      band=1
136    ELSEIF (Nwv.LE.5) THEN      !--RRTM spectral band 2
137      zlambda=0.5*(lambda(Nwv)+lambda(Nwv+1))
138      zweight=weight(Nwv)/SUM(weight(1:5))
139      band=2
140    ELSEIF (Nwv.LE.10) THEN     !--RRTM spectral band 3
141      zlambda=0.5*(lambda(Nwv)+lambda(Nwv+1))
142      zweight=weight(Nwv)/SUM(weight(6:10))
143      band=3
144    ELSEIF (Nwv.LE.16) THEN     !--RRTM spectral band 4
145      zlambda=0.5*(lambda(Nwv)+lambda(Nwv+1))
146      zweight=weight(Nwv)/SUM(weight(11:16))
147      band=4
148    ELSEIF (Nwv.LE.21) THEN     !--RRTM spectral band 5
149      zlambda=0.5*(lambda(Nwv)+lambda(Nwv+1))
150      zweight=weight(Nwv)/SUM(weight(17:21))
151      band=5
152    ELSE                        !--RRTM spectral band 6
153      zlambda=0.5*(lambda(Nwv)+lambda(Nwv+1))
154      zweight=weight(Nwv)/SUM(weight(22:Nwvmax-1))
155      band=6
156    ENDIF
157    !
158    CALL sp_aop_profile                                       ( &
159         klev     ,klon ,zlambda ,oro    ,xlon     ,xlat      , &
160         year_fr  ,z    ,dz      ,dNovrN ,aod_prof ,ssa_prof  , &
161         asy_prof )
162    !
163    !--adding up the quantities tau, piz*tau and cg*piz*tau
164    tau_allaer(:,:,2,band)=tau_allaer(:,:,2,band)+zweight*MAX(aod_prof(:,:),1.e-15)
165    piz_allaer(:,:,2,band)=piz_allaer(:,:,2,band)+zweight*MAX(aod_prof(:,:),1.e-15)*ssa_prof(:,:)
166    cg_allaer(:,:,2,band) =cg_allaer(:,:,2,band) +zweight*MAX(aod_prof(:,:),1.e-15)*ssa_prof(:,:)*asy_prof(:,:)
167    !
168  ENDDO
169  !
170  !--renpomalizing cg and piz now that MACv2SP increments have been added
171  cg_allaer(:,:,2,:) =cg_allaer(:,:,2,:) /piz_allaer(:,:,2,:)
172  piz_allaer(:,:,2,:)=piz_allaer(:,:,2,:)/tau_allaer(:,:,2,:)
173  !
174END SUBROUTINE MACv2SP
Note: See TracBrowser for help on using the repository browser.