Context Navigation

newmicro.F @ 735

Last change on this file since 735 was 704, checked in by Laurent Fairhead, 18 years ago
Inclusion des modifs de Y. Meurdesoif pour la version V3 LF
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 13.4 KB

Line
1	!
2	! $Header$
3	!
4	SUBROUTINE newmicro (paprs, pplay,ok_newmicro,
5	. t, pqlwp, pclc, pcltau, pclemi,
6	. pch, pcl, pcm, pct, pctlwp,
7	s xflwp, xfiwp, xflwc, xfiwc,
8	e ok_aie,
9	e sulfate, sulfate_pi,
10	e bl95_b0, bl95_b1,
11	s cldtaupi, re, fl)
12	USE dimphy
13	IMPLICIT none
14	c======================================================================
15	c Auteur(s): Z.X. Li (LMD/CNRS) date: 19930910
16	c Objet: Calculer epaisseur optique et emmissivite des nuages
17	c======================================================================
18	c Arguments:
19	c t-------input-R-temperature
20	c pqlwp---input-R-eau liquide nuageuse dans l'atmosphere (kg/kg)
21	c pclc----input-R-couverture nuageuse pour le rayonnement (0 a 1)
22	c
23	c ok_aie--input-L-apply aerosol indirect effect or not
24	c sulfate-input-R-sulfate aerosol mass concentration [um/m^3]
25	c sulfate_pi-input-R-dito, pre-industrial value
26	c bl95_b0-input-R-a parameter, may be varied for tests (s-sea, l-land)
27	c bl95_b1-input-R-a parameter, may be varied for tests ( -"- )
28	c
29	c cldtaupi-output-R-pre-industrial value of cloud optical thickness,
30	c needed for the diagnostics of the aerosol indirect
31	c radiative forcing (see radlwsw)
32	c re------output-R-Cloud droplet effective radius multiplied by fl [um]
33	c fl------output-R-Denominator to re, introduced to avoid problems in
34	c the averaging of the output. fl is the fraction of liquid
35	c water clouds within a grid cell
36	c pcltau--output-R-epaisseur optique des nuages
37	c pclemi--output-R-emissivite des nuages (0 a 1)
38	c======================================================================
39	C
40	#include "YOMCST.h"
41	c
42	cym#include "dimensions.h"
43	cym#include "dimphy.h"
44	#include "nuage.h"
45	cIM cf. CR: include pour NOVLP et ZEPSEC
46	#include "radepsi.h"
47	#include "radopt.h"
48	REAL paprs(klon,klev+1), pplay(klon,klev)
49	REAL t(klon,klev)
50	c
51	REAL pclc(klon,klev)
52	REAL pqlwp(klon,klev)
53	REAL pcltau(klon,klev), pclemi(klon,klev)
54	c
55	REAL pct(klon), pctlwp(klon), pch(klon), pcl(klon), pcm(klon)
56	c
57	LOGICAL lo
58	c
59	REAL cetahb, cetamb
60	PARAMETER (cetahb = 0.45, cetamb = 0.80)
61	C
62	INTEGER i, k
63	cIM: 091003 REAL zflwp, zradef, zfice, zmsac
64	REAL zflwp(klon), zradef, zfice, zmsac
65	cIM: 091003 rajout
66	REAL xflwp(klon), xfiwp(klon)
67	REAL xflwc(klon,klev), xfiwc(klon,klev)
68	c
69	REAL radius, rad_chaud
70	cc PARAMETER (rad_chau1=13.0, rad_chau2=9.0, rad_froid=35.0)
71	ccc PARAMETER (rad_chaud=15.0, rad_froid=35.0)
72	c sintex initial PARAMETER (rad_chaud=10.0, rad_froid=30.0)
73	REAL coef, coef_froi, coef_chau
74	PARAMETER (coef_chau=0.13, coef_froi=0.09)
75	REAL seuil_neb, t_glace
76	PARAMETER (seuil_neb=0.001, t_glace=273.0-15.0)
77	INTEGER nexpo ! exponentiel pour glace/eau
78	PARAMETER (nexpo=6)
79	ccc PARAMETER (nexpo=1)
80
81	c -- sb:
82	logical ok_newmicro
83	c parameter (ok_newmicro=.FALSE.)
84	cIM: 091003 real rel, tc, rei, zfiwp
85	real rel, tc, rei, zfiwp(klon)
86	real k_liq, k_ice0, k_ice, DF
87	parameter (k_liq=0.0903, k_ice0=0.005) ! units=m2/g
88	parameter (DF=1.66) ! diffusivity factor
89	c sb --
90	cjq for the aerosol indirect effect
91	cjq introduced by Johannes Quaas (quaas@lmd.jussieu.fr), 27/11/2003
92	cjq
93	LOGICAL ok_aie ! Apply AIE or not?
94	LOGICAL ok_a1lwpdep ! a1 LWP dependent?
95
96	REAL sulfate(klon, klev) ! sulfate aerosol mass concentration [ug m-3]
97	REAL cdnc(klon, klev) ! cloud droplet number concentration [m-3]
98	REAL re(klon, klev) ! cloud droplet effective radius [um]
99	REAL sulfate_pi(klon, klev) ! sulfate aerosol mass concentration [ug m-3] (pre-industrial value)
100	REAL cdnc_pi(klon, klev) ! cloud droplet number concentration [m-3] (pi value)
101	REAL re_pi(klon, klev) ! cloud droplet effective radius [um] (pi value)
102
103	REAL fl(klon, klev) ! xliq * rneb (denominator to re; fraction of liquid water clouds within the grid cell)
104
105	REAL bl95_b0, bl95_b1 ! Parameter in B&L 95-Formula
106
107	REAL cldtaupi(klon, klev) ! pre-industrial cloud opt thickness for diag
108	cjq-end
109	cIM cf. CR:parametres supplementaires
110	REAL zclear(klon)
111	REAL zcloud(klon)
112	c
113	c Calculer l'epaisseur optique et l'emmissivite des nuages
114	c
115	cIM inversion des DO
116	DO i = 1, klon
117	xflwp(i)=0.
118	xfiwp(i)=0.
119	!CDIR NOVECTOR
120	DO k = 1, klev
121	c
122	xflwc(i,k)=0.
123	xfiwc(i,k)=0.
124	c
125	rad_chaud = rad_chau1
126	IF (k.LE.3) rad_chaud = rad_chau2
127	pclc(i,k) = MAX(pclc(i,k), seuil_neb)
128	zflwp(i) = 1000.*pqlwp(i,k)/RG/pclc(i,k)
129	. *(paprs(i,k)-paprs(i,k+1))
130	zfice = 1.0 - (t(i,k)-t_glace) / (273.13-t_glace)
131	zfice = MIN(MAX(zfice,0.0),1.0)
132	zfice = zfice**nexpo
133	radius = rad_chaud * (1.-zfice) + rad_froid * zfice
134	coef = coef_chau * (1.-zfice) + coef_froi * zfice
135	pcltau(i,k) = 3.0/2.0 * zflwp(i) / radius
136	pclemi(i,k) = 1.0 - EXP( - coef * zflwp(i))
137
138	if (ok_newmicro) then
139
140	c -- liquid/ice cloud water paths:
141
142	zfice = 1.0 - (t(i,k)-t_glace) / (273.13-t_glace)
143	zfice = MIN(MAX(zfice,0.0),1.0)
144
145	zflwp(i) = 1000.(1.-zfice)pqlwp(i,k)/pclc(i,k)
146	: *(paprs(i,k)-paprs(i,k+1))/RG
147	zfiwp(i) = 1000.zficepqlwp(i,k)/pclc(i,k)
148	: *(paprs(i,k)-paprs(i,k+1))/RG
149
150	xflwp(i) = xflwp(i)+ (1.-zfice)*pqlwp(i,k)
151	: *(paprs(i,k)-paprs(i,k+1))/RG
152	xfiwp(i) = xfiwp(i)+ zfice*pqlwp(i,k)
153	: *(paprs(i,k)-paprs(i,k+1))/RG
154
155	cIM Total Liquid/Ice water content
156	xflwc(i,k) = xflwc(i,k)+(1.-zfice)*pqlwp(i,k)
157	xfiwc(i,k) = xfiwc(i,k)+zfice*pqlwp(i,k)
158	cIM In-Cloud Liquid/Ice water content
159	c xflwc(i,k) = xflwc(i,k)+(1.-zfice)*pqlwp(i,k)/pclc(i,k)
160	c xfiwc(i,k) = xfiwc(i,k)+zfice*pqlwp(i,k)/pclc(i,k)
161
162	c -- effective cloud droplet radius (microns):
163
164	c for liquid water clouds:
165	IF (ok_aie) THEN
166	! Formula "D" of Boucher and Lohmann, Tellus, 1995
167	!
168	cdnc(i,k) = 10.*(bl95_b0+bl95_b1
169	. log(MAX(sulfate(i,k),1.e-4))/log(10.))*1.e6 !-m-3
170	! Cloud droplet number concentration (CDNC) is restricted
171	! to be within [20, 1000 cm^3]
172	!
173	cdnc(i,k)=MIN(1000.e6,MAX(20.e6,cdnc(i,k)))
174	!
175	!
176	cdnc_pi(i,k) = 10.*(bl95_b0+bl95_b1
177	. log(MAX(sulfate_pi(i,k),1.e-4))/log(10.))*1.e6 !-m-3
178	cdnc_pi(i,k)=MIN(1000.e6,MAX(20.e6,cdnc_pi(i,k)))
179	!
180	!
181	! air density: pplay(i,k) / (RD * zT(i,k))
182	! factor 1.1: derive effective radius from volume-mean radius
183	! factor 1000 is the water density
184	! _chaud means that this is the CDR for liquid water clouds
185	!
186	rad_chaud =
187	. 1.1 * ( (pqlwp(i,k) * pplay(i,k) / (RD * T(i,k)) )
188	. / (4./3. * RPI * 1000. * cdnc(i,k)) )**(1./3.)
189	!
190	! Convert to um. CDR shall be at least 3 um.
191	!
192	c rad_chaud = MAX(rad_chaud*1.e6, 3.)
193	rad_chaud = MAX(rad_chaud*1.e6, 5.)
194
195	! Pre-industrial cloud opt thickness
196	!
197	! "radius" is calculated as rad_chaud above (plus the
198	! ice cloud contribution) but using cdnc_pi instead of
199	! cdnc.
200	radius =
201	. 1.1 * ( (pqlwp(i,k) * pplay(i,k) / (RD * T(i,k)) )
202	. / (4./3. * RPI * 1000. * cdnc_pi(i,k)) )**(1./3.)
203	radius = MAX(radius*1.e6, 5.)
204
205	tc = t(i,k)-273.15
206	rei = 0.71*tc + 61.29
207	if (tc.le.-81.4) rei = 3.5
208	if (zflwp(i).eq.0.) radius = 1.
209	if (zfiwp(i).eq.0. .or. rei.le.0.) rei = 1.
210	cldtaupi(i,k) = 3.0/2.0 * zflwp(i) / radius
211	. + zfiwp(i) * (3.448e-03 + 2.431/rei)
212	ENDIF ! ok_aie
213	! For output diagnostics
214	!
215	! Cloud droplet effective radius [um]
216	!
217	! we multiply here with f * xl (fraction of liquid water
218	! clouds in the grid cell) to avoid problems in the
219	! averaging of the output.
220	! In the output of IOIPSL, derive the real cloud droplet
221	! effective radius as re/fl
222	!
223	fl(i,k) = pclc(i,k)*(1.-zfice)
224	re(i,k) = rad_chaud*fl(i,k)
225
226	c-jq end
227
228	rel = rad_chaud
229	c for ice clouds: as a function of the ambiant temperature
230	c [formula used by Iacobellis and Somerville (2000), with an
231	c asymptotical value of 3.5 microns at T<-81.4 C added to be
232	c consistent with observations of Heymsfield et al. 1986]:
233	tc = t(i,k)-273.15
234	rei = 0.71*tc + 61.29
235	if (tc.le.-81.4) rei = 3.5
236
237	c -- cloud optical thickness :
238
239	c [for liquid clouds, traditional formula,
240	c for ice clouds, Ebert & Curry (1992)]
241
242	if (zflwp(i).eq.0.) rel = 1.
243	if (zfiwp(i).eq.0. .or. rei.le.0.) rei = 1.
244	pcltau(i,k) = 3.0/2.0 * ( zflwp(i)/rel )
245	. + zfiwp(i) * (3.448e-03 + 2.431/rei)
246
247	c -- cloud infrared emissivity:
248
249	c [the broadband infrared absorption coefficient is parameterized
250	c as a function of the effective cld droplet radius]
251
252	c Ebert and Curry (1992) formula as used by Kiehl & Zender (1995):
253	k_ice = k_ice0 + 1.0/rei
254
255	pclemi(i,k) = 1.0
256	. - EXP( - coef_chauzflwp(i) - DFk_ice*zfiwp(i) )
257
258	endif ! ok_newmicro
259
260	lo = (pclc(i,k) .LE. seuil_neb)
261	IF (lo) pclc(i,k) = 0.0
262	IF (lo) pcltau(i,k) = 0.0
263	IF (lo) pclemi(i,k) = 0.0
264
265	IF (lo) cldtaupi(i,k) = 0.0
266	IF (.NOT.ok_aie) cldtaupi(i,k)=pcltau(i,k)
267	ENDDO
268	ENDDO
269	ccc DO k = 1, klev
270	ccc DO i = 1, klon
271	ccc t(i,k) = t(i,k)
272	ccc pclc(i,k) = MAX( 1.e-5 , pclc(i,k) )
273	ccc lo = pclc(i,k) .GT. (2.*1.e-5)
274	ccc zflwp = pqlwp(i,k)1000.(paprs(i,k)-paprs(i,k+1))
275	ccc . /(rg*pclc(i,k))
276	ccc zradef = 10.0 + (1.-sigs(k))*45.0
277	ccc pcltau(i,k) = 1.5 * zflwp / zradef
278	ccc zfice=1.0-MIN(MAX((t(i,k)-263.)/(273.-263.),0.0),1.0)
279	ccc zmsac = 0.13(1.0-zfice) + 0.08zfice
280	ccc pclemi(i,k) = 1.-EXP(-zmsac*zflwp)
281	ccc if (.NOT.lo) pclc(i,k) = 0.0
282	ccc if (.NOT.lo) pcltau(i,k) = 0.0
283	ccc if (.NOT.lo) pclemi(i,k) = 0.0
284	ccc ENDDO
285	ccc ENDDO
286	cccccc print*, 'pas de nuage dans le rayonnement'
287	cccccc DO k = 1, klev
288	cccccc DO i = 1, klon
289	cccccc pclc(i,k) = 0.0
290	cccccc pcltau(i,k) = 0.0
291	cccccc pclemi(i,k) = 0.0
292	cccccc ENDDO
293	cccccc ENDDO
294	C
295	C COMPUTE CLOUD LIQUID PATH AND TOTAL CLOUDINESS
296	C
297	cIM cf. CR:test: calcul prenant ou non en compte le recouvrement
298	cinitialisations
299	DO i=1,klon
300	zclear(i)=1.
301	zcloud(i)=0.
302	pch(i)=1.0
303	pcm(i) = 1.0
304	pcl(i) = 1.0
305	pctlwp(i) = 0.0
306	ENDDO
307	C
308	cIM cf CR DO k=1,klev
309	DO k = klev, 1, -1
310	DO i = 1, klon
311	pctlwp(i) = pctlwp(i)
312	. + pqlwp(i,k)*(paprs(i,k)-paprs(i,k+1))/RG
313	cIM cf. CR
314	IF (NOVLP.EQ.1) THEN
315	zclear(i)=zclear(i)*(1.-MAX(pclc(i,k),zcloud(i)))
316	s /(1.-MIN(zcloud(i),1.-ZEPSEC))
317	pct(i)=1.-zclear(i)
318	if (pplay(i,k).LE.cetahb*paprs(i,1)) then
319	pch(i) = pch(i)*(1.-MAX(pclc(i,k),zcloud(i)))
320	s /(1.-MIN(zcloud(i),1.-ZEPSEC))
321	else if (pplay(i,k).GT.cetahb*paprs(i,1) .AND.
322	. pplay(i,k).LE.cetamb*paprs(i,1)) then
323	pcm(i) = pcm(i)*(1.-MAX(pclc(i,k),zcloud(i)))
324	s /(1.-MIN(zcloud(i),1.-ZEPSEC))
325	else if (pplay(i,k).GT.cetamb*paprs(i,1)) then
326	pcl(i) = pcl(i)*(1.-MAX(pclc(i,k),zcloud(i)))
327	s /(1.-MIN(zcloud(i),1.-ZEPSEC))
328	endif
329	zcloud(i)=pclc(i,k)
330	ELSE IF (NOVLP.EQ.2) THEN
331	zcloud(i)=MAX(pclc(i,k),zcloud(i))
332	pct(i)=zcloud(i)
333	if (pplay(i,k).LE.cetahb*paprs(i,1)) then
334	pch(i) = MIN(pclc(i,k),pch(i))
335	else if (pplay(i,k).GT.cetahb*paprs(i,1) .AND.
336	. pplay(i,k).LE.cetamb*paprs(i,1)) then
337	pcm(i) = MIN(pclc(i,k),pcm(i))
338	else if (pplay(i,k).GT.cetamb*paprs(i,1)) then
339	pcl(i) = MIN(pclc(i,k),pcl(i))
340	endif
341	ELSE IF (NOVLP.EQ.3) THEN
342	zclear(i)=zclear(i)*(1.-pclc(i,k))
343	pct(i)=1-zclear(i)
344	if (pplay(i,k).LE.cetahb*paprs(i,1)) then
345	pch(i) = pch(i)*(1.0-pclc(i,k))
346	else if (pplay(i,k).GT.cetahb*paprs(i,1) .AND.
347	. pplay(i,k).LE.cetamb*paprs(i,1)) then
348	pcm(i) = pcm(i)*(1.0-pclc(i,k))
349	else if (pplay(i,k).GT.cetamb*paprs(i,1)) then
350	pcl(i) = pcl(i)*(1.0-pclc(i,k))
351	endif
352	ENDIF
353	ENDDO
354	ENDDO
355	C
356	DO i = 1, klon
357	cIM cf. CR pct(i)=1.-pct(i)
358	pch(i)=1.-pch(i)
359	pcm(i)=1.-pcm(i)
360	pcl(i)=1.-pcl(i)
361	ENDDO
362	C
363	RETURN
364	END

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format