My Project
 All Classes Files Functions Variables Macros
aeropt.F
Go to the documentation of this file.
1 !
2 ! $Id: aeropt.F 1403 2010-07-01 09:02:53Z fairhead $
3 !
4  SUBROUTINE aeropt(pplay, paprs, t_seri, msulfate, RHcl,
5  . tau_ae, piz_ae, cg_ae, ai )
6 c
7  USE dimphy
8  IMPLICIT none
9 c
10 c
11 c
12 cym#include "dimensions.h"
13 cym#include "dimphy.h"
14 #include "YOMCST.h"
15 c
16 c Arguments:
17 c
18  REAL, INTENT(in) :: paprs(klon,klev+1)
19  REAL, INTENT(in) :: pplay(klon,klev), t_seri(klon,klev)
20  REAL, INTENT(in) :: msulfate(klon,klev) ! masse sulfate ug SO4/m3 [ug/m^3]
21  REAL, INTENT(in) :: rhcl(klon,klev) ! humidite relative ciel clair
22  REAL, INTENT(out) :: tau_ae(klon,klev,2) ! epaisseur optique aerosol
23  REAL, INTENT(out) :: piz_ae(klon,klev,2) ! single scattering albedo aerosol
24  REAL, INTENT(out) :: cg_ae(klon,klev,2) ! asymmetry parameter aerosol
25  REAL, INTENT(out) :: ai(klon) ! POLDER aerosol index
26 c
27 c Local
28 c
29  INTEGER i, k, inu
30  INTEGER rh_num, nbre_rh
31  parameter(nbre_rh=12)
32  REAL rh_tab(nbre_rh)
33  REAL rh_max, delta, rh
34  parameter(rh_max=95.)
35  DATA rh_tab/0.,10.,20.,30.,40.,50.,60.,70.,80.,85.,90.,95./
36  REAL zrho, zdz
37  REAL taue670(klon) ! epaisseur optique aerosol absorption 550 nm
38  REAL taue865(klon) ! epaisseur optique aerosol extinction 865 nm
39  REAL alpha_aer_sulfate(nbre_rh,5) !--unit m2/g SO4
40  REAL alphasulfate
41 
42  CHARACTER (LEN=20) :: modname='aeropt'
43  CHARACTER (LEN=80) :: abort_message
44 
45 c
46 c Proprietes optiques
47 c
48  REAL alpha_aer(nbre_rh,2) !--unit m2/g SO4
49  REAL cg_aer(nbre_rh,2)
50  DATA alpha_aer/.500130e+01, .500130e+01, .500130e+01,
51  . .500130e+01, .500130e+01, .616710e+01,
52  . .826850e+01, .107687e+02, .136976e+02,
53  . .162972e+02, .211690e+02, .354833e+02,
54  . .139460e+01, .139460e+01, .139460e+01,
55  . .139460e+01, .139460e+01, .173910e+01,
56  . .244380e+01, .332320e+01, .440120e+01,
57  . .539570e+01, .734580e+01, .136038e+02 /
58  DATA cg_aer/.619800e+00, .619800e+00, .619800e+00,
59  . .619800e+00, .619800e+00, .662700e+00,
60  . .682100e+00, .698500e+00, .712500e+00,
61  . .721800e+00, .734600e+00, .755800e+00,
62  . .545600e+00, .545600e+00, .545600e+00,
63  . .545600e+00, .545600e+00, .583700e+00,
64  . .607100e+00, .627700e+00, .645800e+00,
65  . .658400e+00, .676500e+00, .708500e+00 /
66  DATA alpha_aer_sulfate/
67  . 4.910,4.910,4.910,4.910,6.547,7.373,
68  . 8.373,9.788,12.167,14.256,17.924,28.433,
69  . 1.453,1.453,1.453,1.453,2.003,2.321,
70  . 2.711,3.282,4.287,5.210,6.914,12.305,
71  . 4.308,4.308,4.308,4.308,5.753,6.521,
72  . 7.449,8.772,11.014,12.999,16.518,26.772,
73  . 3.265,3.265,3.265,3.265,4.388,5.016,
74  . 5.775,6.868,8.745,10.429,13.457,22.538,
75  . 2.116,2.116,2.116,2.116,2.882,3.330,
76  . 3.876,4.670,6.059,7.327,9.650,16.883/
77 c
78  DO i=1, klon
79  taue670(i)=0.0
80  taue865(i)=0.0
81  ENDDO
82 c
83  DO k=1, klev
84  DO i=1, klon
85  if (t_seri(i,k).eq.0) write (*,*) 'aeropt T ',i,k,t_seri(i,k)
86  if (pplay(i,k).eq.0) write (*,*) 'aeropt p ',i,k,pplay(i,k)
87  zrho=pplay(i,k)/t_seri(i,k)/rd ! kg/m3
88  zdz=(paprs(i,k)-paprs(i,k+1))/zrho/rg ! m
89  rh=min(rhcl(i,k)*100.,rh_max)
90  rh_num = int( rh/10. + 1.)
91  IF (rh.LT.0.) THEN
92  abort_message = 'aeropt: RH < 0 not possible'
93  CALL abort_gcm(modname,abort_message,1)
94  ENDIF
95  IF (rh.gt.85.) rh_num=10
96  IF (rh.gt.90.) rh_num=11
97  delta=(rh-rh_tab(rh_num))/(rh_tab(rh_num+1)-rh_tab(rh_num))
98 c
99  inu=1
100  tau_ae(i,k,inu)=alpha_aer(rh_num,inu) +
101  . delta*(alpha_aer(rh_num+1,inu)-alpha_aer(rh_num,inu))
102  tau_ae(i,k,inu)=tau_ae(i,k,inu)*msulfate(i,k)*zdz*1.e-6
103  piz_ae(i,k,inu)=1.0
104  cg_ae(i,k,inu)=cg_aer(rh_num,inu) +
105  . delta*(cg_aer(rh_num+1,inu)-cg_aer(rh_num,inu))
106 c
107  inu=2
108  tau_ae(i,k,inu)=alpha_aer(rh_num,inu) +
109  . delta*(alpha_aer(rh_num+1,inu)-alpha_aer(rh_num,inu))
110  tau_ae(i,k,inu)=tau_ae(i,k,inu)*msulfate(i,k)*zdz*1.e-6
111  piz_ae(i,k,inu)=1.0
112  cg_ae(i,k,inu)=cg_aer(rh_num,inu) +
113  . delta*(cg_aer(rh_num+1,inu)-cg_aer(rh_num,inu))
114 cjq
115 cjq for aerosol index
116 c
117  alphasulfate=alpha_aer_sulfate(rh_num,4) +
118  . delta*(alpha_aer_sulfate(rh_num+1,4)-
119  . alpha_aer_sulfate(rh_num,4)) !--m2/g
120 c
121  taue670(i)=taue670(i)+
122  . alphasulfate*msulfate(i,k)*zdz*1.e-6
123 c
124  alphasulfate=alpha_aer_sulfate(rh_num,5) +
125  . delta*(alpha_aer_sulfate(rh_num+1,5)-
126  . alpha_aer_sulfate(rh_num,5)) !--m2/g
127 c
128  taue865(i)=taue865(i)+
129  . alphasulfate*msulfate(i,k)*zdz*1.e-6
130 
131  ENDDO
132  ENDDO
133 c
134  DO i=1, klon
135  ai(i)=(-log(max(taue670(i),0.0001)/
136  . max(taue865(i),0.0001))/log(670./865.)) *
137  . taue865(i)
138  ENDDO
139 c
140  RETURN
141  END