source: LMDZ5/trunk/libf/phylmd/newmicro.F90 @ 2079

Last change on this file since 2079 was 2077, checked in by fhourdin, 10 years ago

Modification relative à la phase mixte des nuages :

  1. on fait tendre t_glace_min vers t_glace_max (en principe 0°C) linéairerment entre p/ps = 0.8 et p/ps=1.
  2. Passage de tableau à la routine icefrac_lsc en remplacement d'une fonction scalaire.
  3. Changement de nom pour le module (le module icefrac_lsc_mod.F90 contenant maintenant la routine icefrac_lsc).

Mofication concerning the mixte phase of clouds

  1. t_glace_min -> t_glace_max when p/ps : 0.8 -> 1
  2. passing arrays instead of scalars to icefrac_lsc (which becomes a routine rather than a function).
  3. Changing the name microphys_mod.F90 -> icerfrac_lsc_mod.F90
  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 20.2 KB
Line 
1! $Id: newmicro.F90 2077 2014-07-03 13:38:05Z lguez $
2
3
4
5SUBROUTINE newmicro(ok_cdnc, bl95_b0, bl95_b1, paprs, pplay, t, pqlwp, pclc, &
6    pcltau, pclemi, pch, pcl, pcm, pct, pctlwp, xflwp, xfiwp, xflwc, xfiwc, &
7    mass_solu_aero, mass_solu_aero_pi, pcldtaupi, re, fl, reliq, reice, &
8    reliq_pi, reice_pi)
9
10  USE dimphy
11  USE phys_local_var_mod, ONLY: scdnc, cldncl, reffclwtop, lcc, reffclws, &
12    reffclwc, cldnvi, lcc3d, lcc3dcon, lcc3dstra
13  USE phys_state_var_mod, ONLY: rnebcon, clwcon
14  USE icefrac_lsc_mod ! cloud microphysics (JBM 3/14)
15  IMPLICIT NONE
16  ! ======================================================================
17  ! Auteur(s): Z.X. Li (LMD/CNRS) date: 19930910
18  ! O.   Boucher (LMD/CNRS) mise a jour en 201212
19  ! Objet: Calculer epaisseur optique et emmissivite des nuages
20  ! ======================================================================
21  ! Arguments:
22  ! ok_cdnc-input-L-flag pour calculer les rayons a partir des aerosols
23
24  ! t-------input-R-temperature
25  ! pqlwp---input-R-eau liquide nuageuse dans l'atmosphere dans la partie
26  ! nuageuse (kg/kg)
27  ! pclc----input-R-couverture nuageuse pour le rayonnement (0 a 1)
28  ! mass_solu_aero-----input-R-total mass concentration for all soluble
29  ! aerosols[ug/m^3]
30  ! mass_solu_aero_pi--input-R-ditto, pre-industrial value
31
32  ! bl95_b0-input-R-a PARAMETER, may be varied for tests (s-sea, l-land)
33  ! bl95_b1-input-R-a PARAMETER, may be varied for tests (    -"-      )
34
35  ! re------output-R-Cloud droplet effective radius multiplied by fl [um]
36  ! fl------output-R-Denominator to re, introduced to avoid problems in
37  ! the averaging of the output. fl is the fraction of liquid
38  ! water clouds within a grid cell
39
40  ! pcltau--output-R-epaisseur optique des nuages
41  ! pclemi--output-R-emissivite des nuages (0 a 1)
42  ! pcldtaupi-output-R-pre-industrial value of cloud optical thickness,
43
44  ! pcl-output-R-2D low-level cloud cover
45  ! pcm-output-R-2D mid-level cloud cover
46  ! pch-output-R-2D high-level cloud cover
47  ! pct-output-R-2D total cloud cover
48  ! ======================================================================
49
50  include "YOMCST.h"
51  include "nuage.h"
52  include "radepsi.h"
53  include "radopt.h"
54
55  ! choix de l'hypothese de recouvrememnt nuageuse
56  LOGICAL random, maximum_random, maximum
57  PARAMETER (random=.FALSE., maximum_random=.TRUE., maximum=.FALSE.)
58
59  LOGICAL, SAVE :: first = .TRUE.
60  !$OMP THREADPRIVATE(FIRST)
61  INTEGER flag_max
62
63  ! threshold PARAMETERs
64  REAL thres_tau, thres_neb
65  PARAMETER (thres_tau=0.3, thres_neb=0.001)
66
67  REAL phase3d(klon, klev)
68  REAL tcc(klon), ftmp(klon), lcc_integrat(klon), height(klon)
69
70  REAL paprs(klon, klev+1)
71  REAL pplay(klon, klev)
72  REAL t(klon, klev)
73  REAL pclc(klon, klev)
74  REAL pqlwp(klon, klev)
75  REAL pcltau(klon, klev)
76  REAL pclemi(klon, klev)
77  REAL pcldtaupi(klon, klev)
78
79  REAL pct(klon)
80  REAL pcl(klon)
81  REAL pcm(klon)
82  REAL pch(klon)
83  REAL pctlwp(klon)
84
85  LOGICAL lo
86
87  ! !Abderr modif JL mail du 19.01.2011 18:31
88  ! REAL cetahb, cetamb
89  ! PARAMETER (cetahb = 0.45, cetamb = 0.80)
90  ! Remplacer
91  ! cetahb*paprs(i,1) par  prmhc
92  ! cetamb*paprs(i,1) par  prlmc
93  REAL prmhc ! Pressure between medium and high level cloud in Pa
94  REAL prlmc ! Pressure between low and medium level cloud in Pa
95  PARAMETER (prmhc=440.*100., prlmc=680.*100.)
96
97  INTEGER i, k
98  REAL xflwp(klon), xfiwp(klon)
99  REAL xflwc(klon, klev), xfiwc(klon, klev)
100
101  REAL radius
102
103  REAL coef_froi, coef_chau
104  PARAMETER (coef_chau=0.13, coef_froi=0.09)
105
106  REAL seuil_neb
107  PARAMETER (seuil_neb=0.001)
108
109! JBM (3/14) nexpo is replaced by exposant_glace
110! INTEGER nexpo ! exponentiel pour glace/eau
111! PARAMETER (nexpo=6)
112! PARAMETER (nexpo=1)
113! if iflag_t_glace=0, the old values are used:
114  REAL, PARAMETER :: t_glace_min_old = 258.
115  REAL, PARAMETER :: t_glace_max_old = 273.13
116
117  REAL rel, tc, rei
118  REAL k_ice0, k_ice, df
119  PARAMETER (k_ice0=0.005) ! units=m2/g
120  PARAMETER (df=1.66) ! diffusivity factor
121
122  ! jq for the aerosol indirect effect
123  ! jq introduced by Johannes Quaas (quaas@lmd.jussieu.fr), 27/11/2003
124  ! jq
125  REAL mass_solu_aero(klon, klev) ! total mass concentration for all soluble aerosols [ug m-3]
126  REAL mass_solu_aero_pi(klon, klev) ! - " - (pre-industrial value)
127  REAL cdnc(klon, klev) ! cloud droplet number concentration [m-3]
128  REAL re(klon, klev) ! cloud droplet effective radius [um]
129  REAL cdnc_pi(klon, klev) ! cloud droplet number concentration [m-3] (pi value)
130  REAL re_pi(klon, klev) ! cloud droplet effective radius [um] (pi value)
131
132  REAL fl(klon, klev) ! xliq * rneb (denominator to re; fraction of liquid water clouds
133  ! within the grid cell)
134
135  LOGICAL ok_cdnc
136  REAL bl95_b0, bl95_b1 ! Parameter in B&L 95-Formula
137
138  ! jq-end
139  ! IM cf. CR:parametres supplementaires
140  REAL zclear(klon)
141  REAL zcloud(klon)
142  REAL zcloudh(klon)
143  REAL zcloudm(klon)
144  REAL zcloudl(klon)
145  REAL rhodz(klon, klev) !--rho*dz pour la couche
146  REAL zrho(klon, klev) !--rho pour la couche
147  REAL dh(klon, klev) !--dz pour la couche
148  REAL zfice(klon, klev)
149  REAL rad_chaud(klon, klev) !--rayon pour les nuages chauds
150  REAL rad_chaud_pi(klon, klev) !--rayon pour les nuages chauds pre-industriels
151  REAL zflwp_var, zfiwp_var
152  REAL d_rei_dt
153
154  ! Abderrahmane oct 2009
155  REAL reliq(klon, klev), reice(klon, klev)
156  REAL reliq_pi(klon, klev), reice_pi(klon, klev)
157
158  ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
159  ! FH : 2011/05/24
160
161  ! rei = ( rei_max - rei_min ) * T(°C) / 81.4 + rei_max
162  ! to be used for a temperature in celcius T(°C) < 0
163  ! rei=rei_min for T(°C) < -81.4
164
165  ! Calcul de la pente de la relation entre rayon effective des cristaux
166  ! et la température.
167  ! Pour retrouver les résultats numériques de la version d'origine,
168  ! on impose 0.71 quand on est proche de 0.71
169
170  d_rei_dt = (rei_max-rei_min)/81.4
171  IF (abs(d_rei_dt-0.71)<1.E-4) d_rei_dt = 0.71
172  ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
173
174  ! Calculer l'epaisseur optique et l'emmissivite des nuages
175  ! IM inversion des DO
176
177  xflwp = 0.D0
178  xfiwp = 0.D0
179  xflwc = 0.D0
180  xfiwc = 0.D0
181
182  reliq = 0.
183  reice = 0.
184  reliq_pi = 0.
185  reice_pi = 0.
186
187  IF (iflag_t_glace.EQ.0) THEN
188    DO k = 1, klev
189      DO i = 1, klon
190        ! -layer calculation
191        rhodz(i, k) = (paprs(i,k)-paprs(i,k+1))/rg ! kg/m2
192        zrho(i, k) = pplay(i, k)/t(i, k)/rd ! kg/m3
193        dh(i, k) = rhodz(i, k)/zrho(i, k) ! m
194        ! -Fraction of ice in cloud using a linear transition
195        zfice(i, k) = 1.0 - (t(i,k)-t_glace_min_old)/(t_glace_max_old-t_glace_min_old)
196        zfice(i, k) = min(max(zfice(i,k),0.0), 1.0)
197        ! -IM Total Liquid/Ice water content
198        xflwc(i, k) = (1.-zfice(i,k))*pqlwp(i, k)
199        xfiwc(i, k) = zfice(i, k)*pqlwp(i, k)
200      END DO
201    END DO
202  ELSE ! of IF (iflag_t_glace.EQ.0)
203    DO k = 1, klev
204        CALL icefrac_lsc(klon,t(:,k),pplay(:,k)/paprs(:,1), &
205&               t_glace_min,t_glace_max,exposant_glace,zfice(:,k))
206 
207
208        ! JBM: icefrac_lsc is now a function contained in microphys_mod
209!       zfice(i, k) = icefrac_lsc(t(i,k), t_glace_min, &
210!                                 t_glace_max, exposant_glace)
211      DO i = 1, klon
212        ! -layer calculation
213        rhodz(i, k) = (paprs(i,k)-paprs(i,k+1))/rg ! kg/m2
214        zrho(i, k) = pplay(i, k)/t(i, k)/rd ! kg/m3
215        dh(i, k) = rhodz(i, k)/zrho(i, k) ! m
216        ! -IM Total Liquid/Ice water content
217        xflwc(i, k) = (1.-zfice(i,k))*pqlwp(i, k)
218        xfiwc(i, k) = zfice(i, k)*pqlwp(i, k)
219      END DO
220    END DO
221  ENDIF
222
223  IF (ok_cdnc) THEN
224
225    ! --we compute cloud properties as a function of the aerosol load
226
227    DO k = 1, klev
228      DO i = 1, klon
229
230        ! Formula "D" of Boucher and Lohmann, Tellus, 1995
231        ! Cloud droplet number concentration (CDNC) is restricted
232        ! to be within [20, 1000 cm^3]
233
234        ! --present-day case
235        cdnc(i, k) = 10.**(bl95_b0+bl95_b1*log(max(mass_solu_aero(i,k), &
236          1.E-4))/log(10.))*1.E6 !-m-3
237        cdnc(i, k) = min(1000.E6, max(20.E6,cdnc(i,k)))
238
239        ! --pre-industrial case
240        cdnc_pi(i, k) = 10.**(bl95_b0+bl95_b1*log(max(mass_solu_aero_pi(i,k), &
241          1.E-4))/log(10.))*1.E6 !-m-3
242        cdnc_pi(i, k) = min(1000.E6, max(20.E6,cdnc_pi(i,k)))
243
244        ! --present-day case
245        rad_chaud(i, k) = 1.1*((pqlwp(i,k)*pplay(i, &
246          k)/(rd*t(i,k)))/(4./3*rpi*1000.*cdnc(i,k)))**(1./3.)
247        rad_chaud(i, k) = max(rad_chaud(i,k)*1.E6, 5.)
248
249        ! --pre-industrial case
250        rad_chaud_pi(i, k) = 1.1*((pqlwp(i,k)*pplay(i, &
251          k)/(rd*t(i,k)))/(4./3.*rpi*1000.*cdnc_pi(i,k)))**(1./3.)
252        rad_chaud_pi(i, k) = max(rad_chaud_pi(i,k)*1.E6, 5.)
253
254        ! --pre-industrial case
255        ! --liquid/ice cloud water paths:
256        IF (pclc(i,k)<=seuil_neb) THEN
257
258          pcldtaupi(i, k) = 0.0
259
260        ELSE
261
262          zflwp_var = 1000.*(1.-zfice(i,k))*pqlwp(i, k)/pclc(i, k)* &
263            rhodz(i, k)
264          zfiwp_var = 1000.*zfice(i, k)*pqlwp(i, k)/pclc(i, k)*rhodz(i, k)
265          tc = t(i, k) - 273.15
266          rei = d_rei_dt*tc + rei_max
267          IF (tc<=-81.4) rei = rei_min
268
269          ! -- cloud optical thickness :
270          ! [for liquid clouds, traditional formula,
271          ! for ice clouds, Ebert & Curry (1992)]
272
273          IF (zfiwp_var==0. .OR. rei<=0.) rei = 1.
274          pcldtaupi(i, k) = 3.0/2.0*zflwp_var/rad_chaud_pi(i, k) + &
275            zfiwp_var*(3.448E-03+2.431/rei)
276
277        END IF
278
279      END DO
280    END DO
281
282  ELSE !--not ok_cdnc
283
284    ! -prescribed cloud droplet radius
285
286    DO k = 1, min(3, klev)
287      DO i = 1, klon
288        rad_chaud(i, k) = rad_chau2
289        rad_chaud_pi(i, k) = rad_chau2
290      END DO
291    END DO
292    DO k = min(3, klev) + 1, klev
293      DO i = 1, klon
294        rad_chaud(i, k) = rad_chau1
295        rad_chaud_pi(i, k) = rad_chau1
296      END DO
297    END DO
298
299  END IF !--ok_cdnc
300
301  ! --computation of cloud optical depth and emissivity
302  ! --in the general case
303
304  DO k = 1, klev
305    DO i = 1, klon
306
307      IF (pclc(i,k)<=seuil_neb) THEN
308
309        ! effective cloud droplet radius (microns) for liquid water clouds:
310        ! For output diagnostics cloud droplet effective radius [um]
311        ! we multiply here with f * xl (fraction of liquid water
312        ! clouds in the grid cell) to avoid problems in the averaging of the
313        ! output.
314        ! In the output of IOIPSL, derive the REAL cloud droplet
315        ! effective radius as re/fl
316
317        fl(i, k) = seuil_neb*(1.-zfice(i,k))
318        re(i, k) = rad_chaud(i, k)*fl(i, k)
319        rel = 0.
320        rei = 0.
321        pclc(i, k) = 0.0
322        pcltau(i, k) = 0.0
323        pclemi(i, k) = 0.0
324
325      ELSE
326
327        ! -- liquid/ice cloud water paths:
328
329        zflwp_var = 1000.*(1.-zfice(i,k))*pqlwp(i, k)/pclc(i, k)*rhodz(i, k)
330        zfiwp_var = 1000.*zfice(i, k)*pqlwp(i, k)/pclc(i, k)*rhodz(i, k)
331
332        ! effective cloud droplet radius (microns) for liquid water clouds:
333        ! For output diagnostics cloud droplet effective radius [um]
334        ! we multiply here with f * xl (fraction of liquid water
335        ! clouds in the grid cell) to avoid problems in the averaging of the
336        ! output.
337        ! In the output of IOIPSL, derive the REAL cloud droplet
338        ! effective radius as re/fl
339
340        fl(i, k) = pclc(i, k)*(1.-zfice(i,k))
341        re(i, k) = rad_chaud(i, k)*fl(i, k)
342
343        rel = rad_chaud(i, k)
344
345        ! for ice clouds: as a function of the ambiant temperature
346        ! [formula used by Iacobellis and Somerville (2000), with an
347        ! asymptotical value of 3.5 microns at T<-81.4 C added to be
348        ! consistent with observations of Heymsfield et al. 1986]:
349        ! 2011/05/24 : rei_min = 3.5 becomes a free PARAMETER as well as
350        ! rei_max=61.29
351
352        tc = t(i, k) - 273.15
353        rei = d_rei_dt*tc + rei_max
354        IF (tc<=-81.4) rei = rei_min
355
356        ! -- cloud optical thickness :
357        ! [for liquid clouds, traditional formula,
358        ! for ice clouds, Ebert & Curry (1992)]
359
360        IF (zflwp_var==0.) rel = 1.
361        IF (zfiwp_var==0. .OR. rei<=0.) rei = 1.
362        pcltau(i, k) = 3.0/2.0*(zflwp_var/rel) + zfiwp_var*(3.448E-03+2.431/ &
363          rei)
364
365        ! -- cloud infrared emissivity:
366        ! [the broadband infrared absorption coefficient is PARAMETERized
367        ! as a function of the effective cld droplet radius]
368        ! Ebert and Curry (1992) formula as used by Kiehl & Zender (1995):
369
370        k_ice = k_ice0 + 1.0/rei
371
372        pclemi(i, k) = 1.0 - exp(-coef_chau*zflwp_var-df*k_ice*zfiwp_var)
373
374      END IF
375
376      reice(i, k) = rei
377
378      xflwp(i) = xflwp(i) + xflwc(i, k)*rhodz(i, k)
379      xfiwp(i) = xfiwp(i) + xfiwc(i, k)*rhodz(i, k)
380
381    END DO
382  END DO
383
384  ! --if cloud droplet radius is fixed, then pcldtaupi=pcltau
385
386  IF (.NOT. ok_cdnc) THEN
387    DO k = 1, klev
388      DO i = 1, klon
389        pcldtaupi(i, k) = pcltau(i, k)
390        reice_pi(i, k) = reice(i, k)
391      END DO
392    END DO
393  END IF
394
395  DO k = 1, klev
396    DO i = 1, klon
397      reliq(i, k) = rad_chaud(i, k)
398      reliq_pi(i, k) = rad_chaud_pi(i, k)
399      reice_pi(i, k) = reice(i, k)
400    END DO
401  END DO
402
403  ! COMPUTE CLOUD LIQUID PATH AND TOTAL CLOUDINESS
404  ! IM cf. CR:test: calcul prenant ou non en compte le recouvrement
405  ! initialisations
406
407  DO i = 1, klon
408    zclear(i) = 1.
409    zcloud(i) = 0.
410    zcloudh(i) = 0.
411    zcloudm(i) = 0.
412    zcloudl(i) = 0.
413    pch(i) = 1.0
414    pcm(i) = 1.0
415    pcl(i) = 1.0
416    pctlwp(i) = 0.0
417  END DO
418
419  ! --calculation of liquid water path
420
421  DO k = klev, 1, -1
422    DO i = 1, klon
423      pctlwp(i) = pctlwp(i) + pqlwp(i, k)*rhodz(i, k)
424    END DO
425  END DO
426
427  ! --calculation of cloud properties with cloud overlap
428
429  IF (novlp==1) THEN
430    DO k = klev, 1, -1
431      DO i = 1, klon
432        zclear(i) = zclear(i)*(1.-max(pclc(i,k),zcloud(i)))/(1.-min(real( &
433          zcloud(i),kind=8),1.-zepsec))
434        pct(i) = 1. - zclear(i)
435        IF (paprs(i,k)<prmhc) THEN
436          pch(i) = pch(i)*(1.-max(pclc(i,k),zcloudh(i)))/(1.-min(real(zcloudh &
437            (i),kind=8),1.-zepsec))
438          zcloudh(i) = pclc(i, k)
439        ELSE IF (paprs(i,k)>=prmhc .AND. paprs(i,k)<prlmc) THEN
440          pcm(i) = pcm(i)*(1.-max(pclc(i,k),zcloudm(i)))/(1.-min(real(zcloudm &
441            (i),kind=8),1.-zepsec))
442          zcloudm(i) = pclc(i, k)
443        ELSE IF (paprs(i,k)>=prlmc) THEN
444          pcl(i) = pcl(i)*(1.-max(pclc(i,k),zcloudl(i)))/(1.-min(real(zcloudl &
445            (i),kind=8),1.-zepsec))
446          zcloudl(i) = pclc(i, k)
447        END IF
448        zcloud(i) = pclc(i, k)
449      END DO
450    END DO
451  ELSE IF (novlp==2) THEN
452    DO k = klev, 1, -1
453      DO i = 1, klon
454        zcloud(i) = max(pclc(i,k), zcloud(i))
455        pct(i) = zcloud(i)
456        IF (paprs(i,k)<prmhc) THEN
457          pch(i) = min(pclc(i,k), pch(i))
458        ELSE IF (paprs(i,k)>=prmhc .AND. paprs(i,k)<prlmc) THEN
459          pcm(i) = min(pclc(i,k), pcm(i))
460        ELSE IF (paprs(i,k)>=prlmc) THEN
461          pcl(i) = min(pclc(i,k), pcl(i))
462        END IF
463      END DO
464    END DO
465  ELSE IF (novlp==3) THEN
466    DO k = klev, 1, -1
467      DO i = 1, klon
468        zclear(i) = zclear(i)*(1.-pclc(i,k))
469        pct(i) = 1 - zclear(i)
470        IF (paprs(i,k)<prmhc) THEN
471          pch(i) = pch(i)*(1.0-pclc(i,k))
472        ELSE IF (paprs(i,k)>=prmhc .AND. paprs(i,k)<prlmc) THEN
473          pcm(i) = pcm(i)*(1.0-pclc(i,k))
474        ELSE IF (paprs(i,k)>=prlmc) THEN
475          pcl(i) = pcl(i)*(1.0-pclc(i,k))
476        END IF
477      END DO
478    END DO
479  END IF
480
481  DO i = 1, klon
482    pch(i) = 1. - pch(i)
483    pcm(i) = 1. - pcm(i)
484    pcl(i) = 1. - pcl(i)
485  END DO
486
487  ! ========================================================
488  ! DIAGNOSTICS CALCULATION FOR CMIP5 PROTOCOL
489  ! ========================================================
490  ! change by Nicolas Yan (LSCE)
491  ! Cloud Droplet Number Concentration (CDNC) : 3D variable
492  ! Fractionnal cover by liquid water cloud (LCC3D) : 3D variable
493  ! Cloud Droplet Number Concentration at top of cloud (CLDNCL) : 2D variable
494  ! Droplet effective radius at top of cloud (REFFCLWTOP) : 2D variable
495  ! Fractionnal cover by liquid water at top of clouds (LCC) : 2D variable
496
497  IF (ok_cdnc) THEN
498
499    DO k = 1, klev
500      DO i = 1, klon
501        phase3d(i, k) = 1 - zfice(i, k)
502        IF (pclc(i,k)<=seuil_neb) THEN
503          lcc3d(i, k) = seuil_neb*phase3d(i, k)
504        ELSE
505          lcc3d(i, k) = pclc(i, k)*phase3d(i, k)
506        END IF
507        scdnc(i, k) = lcc3d(i, k)*cdnc(i, k) ! m-3
508      END DO
509    END DO
510
511    DO i = 1, klon
512      lcc(i) = 0.
513      reffclwtop(i) = 0.
514      cldncl(i) = 0.
515      IF (random .OR. maximum_random) tcc(i) = 1.
516      IF (maximum) tcc(i) = 0.
517    END DO
518
519    DO i = 1, klon
520      DO k = klev - 1, 1, -1 !From TOA down
521
522          ! Test, if the cloud optical depth exceeds the necessary
523          ! threshold:
524
525        IF (pcltau(i,k)>thres_tau .AND. pclc(i,k)>thres_neb) THEN
526
527          IF (maximum) THEN
528            IF (first) THEN
529              WRITE (*, *) 'Hypothese de recouvrement: MAXIMUM'
530              first = .FALSE.
531            END IF
532            flag_max = -1.
533            ftmp(i) = max(tcc(i), pclc(i,k))
534          END IF
535
536          IF (random) THEN
537            IF (first) THEN
538              WRITE (*, *) 'Hypothese de recouvrement: RANDOM'
539              first = .FALSE.
540            END IF
541            flag_max = 1.
542            ftmp(i) = tcc(i)*(1-pclc(i,k))
543          END IF
544
545          IF (maximum_random) THEN
546            IF (first) THEN
547              WRITE (*, *) 'Hypothese de recouvrement: MAXIMUM_ &
548                &                                             &
549                &                                          RANDOM'
550              first = .FALSE.
551            END IF
552            flag_max = 1.
553            ftmp(i) = tcc(i)*(1.-max(pclc(i,k),pclc(i,k+1)))/(1.-min(pclc(i, &
554              k+1),1.-thres_neb))
555          END IF
556          ! Effective radius of cloud droplet at top of cloud (m)
557          reffclwtop(i) = reffclwtop(i) + rad_chaud(i, k)*1.0E-06*phase3d(i, &
558            k)*(tcc(i)-ftmp(i))*flag_max
559          ! CDNC at top of cloud (m-3)
560          cldncl(i) = cldncl(i) + cdnc(i, k)*phase3d(i, k)*(tcc(i)-ftmp(i))* &
561            flag_max
562          ! Liquid Cloud Content at top of cloud
563          lcc(i) = lcc(i) + phase3d(i, k)*(tcc(i)-ftmp(i))*flag_max
564          ! Total Cloud Content at top of cloud
565          tcc(i) = ftmp(i)
566
567        END IF ! is there a visible, not-too-small cloud?
568      END DO ! loop over k
569
570      IF (random .OR. maximum_random) tcc(i) = 1. - tcc(i)
571
572    END DO ! loop over i
573
574    ! ! Convective and Stratiform Cloud Droplet Effective Radius (REFFCLWC
575    ! REFFCLWS)
576    DO i = 1, klon
577      DO k = 1, klev
578        ! Weight to be used for outputs: eau_liquide*couverture nuageuse
579        lcc3dcon(i, k) = rnebcon(i, k)*phase3d(i, k)*clwcon(i, k) ! eau liquide convective
580        lcc3dstra(i, k) = pclc(i, k)*pqlwp(i, k)*phase3d(i, k)
581        lcc3dstra(i, k) = lcc3dstra(i, k) - lcc3dcon(i, k) ! eau liquide stratiforme
582        lcc3dstra(i, k) = max(lcc3dstra(i,k), 0.0)
583        ! Compute cloud droplet radius as above in meter
584        radius = 1.1*((pqlwp(i,k)*pplay(i,k)/(rd*t(i,k)))/(4./3*rpi*1000.* &
585          cdnc(i,k)))**(1./3.)
586        radius = max(radius, 5.E-6)
587        ! Convective Cloud Droplet Effective Radius (REFFCLWC) : variable 3D
588        reffclwc(i, k) = radius
589        reffclwc(i, k) = reffclwc(i, k)*lcc3dcon(i, k)
590        ! Stratiform Cloud Droplet Effective Radius (REFFCLWS) : variable 3D
591        reffclws(i, k) = radius
592        reffclws(i, k) = reffclws(i, k)*lcc3dstra(i, k)
593      END DO !klev
594    END DO !klon
595
596    ! Column Integrated Cloud Droplet Number (CLDNVI) : variable 2D
597
598    DO i = 1, klon
599      cldnvi(i) = 0.
600      lcc_integrat(i) = 0.
601      height(i) = 0.
602      DO k = 1, klev
603        cldnvi(i) = cldnvi(i) + cdnc(i, k)*lcc3d(i, k)*dh(i, k)
604        lcc_integrat(i) = lcc_integrat(i) + lcc3d(i, k)*dh(i, k)
605        height(i) = height(i) + dh(i, k)
606      END DO ! klev
607      lcc_integrat(i) = lcc_integrat(i)/height(i)
608      IF (lcc_integrat(i)<=1.0E-03) THEN
609        cldnvi(i) = cldnvi(i)*lcc(i)/seuil_neb
610      ELSE
611        cldnvi(i) = cldnvi(i)*lcc(i)/lcc_integrat(i)
612      END IF
613    END DO ! klon
614
615    DO i = 1, klon
616      DO k = 1, klev
617        IF (scdnc(i,k)<=0.0) scdnc(i, k) = 0.0
618        IF (reffclws(i,k)<=0.0) reffclws(i, k) = 0.0
619        IF (reffclwc(i,k)<=0.0) reffclwc(i, k) = 0.0
620        IF (lcc3d(i,k)<=0.0) lcc3d(i, k) = 0.0
621        IF (lcc3dcon(i,k)<=0.0) lcc3dcon(i, k) = 0.0
622        IF (lcc3dstra(i,k)<=0.0) lcc3dstra(i, k) = 0.0
623      END DO
624      IF (reffclwtop(i)<=0.0) reffclwtop(i) = 0.0
625      IF (cldncl(i)<=0.0) cldncl(i) = 0.0
626      IF (cldnvi(i)<=0.0) cldnvi(i) = 0.0
627      IF (lcc(i)<=0.0) lcc(i) = 0.0
628    END DO
629
630  END IF !ok_cdnc
631
632  RETURN
633
634END SUBROUTINE newmicro
Note: See TracBrowser for help on using the repository browser.