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

Last change on this file since 2153 was 2109, checked in by jbmadeleine, 10 years ago

nettoyage de icefrac_lsc: les paramètres t_glace_min, t_glace_max et
exposant_glace sont à présents chargés à partir de nuage.h et non
donnés en argument;

cleaning of icefrac_lsc: parameters t_glace_min, t_glace_max and
exposant_glace are loaded from the nuage.h file instead of being
given as arguments;
  • 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.1 KB
Line 
1! $Id: newmicro.F90 2109 2014-08-22 10:22:58Z 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 ! computes ice fraction (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),zfice(:,k))
205 
206
207        ! JBM: icefrac_lsc is now contained icefrac_lsc_mod
208!       zfice(i, k) = icefrac_lsc(t(i,k), t_glace_min, &
209!                                 t_glace_max, exposant_glace)
210      DO i = 1, klon
211        ! -layer calculation
212        rhodz(i, k) = (paprs(i,k)-paprs(i,k+1))/rg ! kg/m2
213        zrho(i, k) = pplay(i, k)/t(i, k)/rd ! kg/m3
214        dh(i, k) = rhodz(i, k)/zrho(i, k) ! m
215        ! -IM Total Liquid/Ice water content
216        xflwc(i, k) = (1.-zfice(i,k))*pqlwp(i, k)
217        xfiwc(i, k) = zfice(i, k)*pqlwp(i, k)
218      END DO
219    END DO
220  ENDIF
221
222  IF (ok_cdnc) THEN
223
224    ! --we compute cloud properties as a function of the aerosol load
225
226    DO k = 1, klev
227      DO i = 1, klon
228
229        ! Formula "D" of Boucher and Lohmann, Tellus, 1995
230        ! Cloud droplet number concentration (CDNC) is restricted
231        ! to be within [20, 1000 cm^3]
232
233        ! --present-day case
234        cdnc(i, k) = 10.**(bl95_b0+bl95_b1*log(max(mass_solu_aero(i,k), &
235          1.E-4))/log(10.))*1.E6 !-m-3
236        cdnc(i, k) = min(1000.E6, max(20.E6,cdnc(i,k)))
237
238        ! --pre-industrial case
239        cdnc_pi(i, k) = 10.**(bl95_b0+bl95_b1*log(max(mass_solu_aero_pi(i,k), &
240          1.E-4))/log(10.))*1.E6 !-m-3
241        cdnc_pi(i, k) = min(1000.E6, max(20.E6,cdnc_pi(i,k)))
242
243        ! --present-day case
244        rad_chaud(i, k) = 1.1*((pqlwp(i,k)*pplay(i, &
245          k)/(rd*t(i,k)))/(4./3*rpi*1000.*cdnc(i,k)))**(1./3.)
246        rad_chaud(i, k) = max(rad_chaud(i,k)*1.E6, 5.)
247
248        ! --pre-industrial case
249        rad_chaud_pi(i, k) = 1.1*((pqlwp(i,k)*pplay(i, &
250          k)/(rd*t(i,k)))/(4./3.*rpi*1000.*cdnc_pi(i,k)))**(1./3.)
251        rad_chaud_pi(i, k) = max(rad_chaud_pi(i,k)*1.E6, 5.)
252
253        ! --pre-industrial case
254        ! --liquid/ice cloud water paths:
255        IF (pclc(i,k)<=seuil_neb) THEN
256
257          pcldtaupi(i, k) = 0.0
258
259        ELSE
260
261          zflwp_var = 1000.*(1.-zfice(i,k))*pqlwp(i, k)/pclc(i, k)* &
262            rhodz(i, k)
263          zfiwp_var = 1000.*zfice(i, k)*pqlwp(i, k)/pclc(i, k)*rhodz(i, k)
264          tc = t(i, k) - 273.15
265          rei = d_rei_dt*tc + rei_max
266          IF (tc<=-81.4) rei = rei_min
267
268          ! -- cloud optical thickness :
269          ! [for liquid clouds, traditional formula,
270          ! for ice clouds, Ebert & Curry (1992)]
271
272          IF (zfiwp_var==0. .OR. rei<=0.) rei = 1.
273          pcldtaupi(i, k) = 3.0/2.0*zflwp_var/rad_chaud_pi(i, k) + &
274            zfiwp_var*(3.448E-03+2.431/rei)
275
276        END IF
277
278      END DO
279    END DO
280
281  ELSE !--not ok_cdnc
282
283    ! -prescribed cloud droplet radius
284
285    DO k = 1, min(3, klev)
286      DO i = 1, klon
287        rad_chaud(i, k) = rad_chau2
288        rad_chaud_pi(i, k) = rad_chau2
289      END DO
290    END DO
291    DO k = min(3, klev) + 1, klev
292      DO i = 1, klon
293        rad_chaud(i, k) = rad_chau1
294        rad_chaud_pi(i, k) = rad_chau1
295      END DO
296    END DO
297
298  END IF !--ok_cdnc
299
300  ! --computation of cloud optical depth and emissivity
301  ! --in the general case
302
303  DO k = 1, klev
304    DO i = 1, klon
305
306      IF (pclc(i,k)<=seuil_neb) THEN
307
308        ! effective cloud droplet radius (microns) for liquid water clouds:
309        ! For output diagnostics cloud droplet effective radius [um]
310        ! we multiply here with f * xl (fraction of liquid water
311        ! clouds in the grid cell) to avoid problems in the averaging of the
312        ! output.
313        ! In the output of IOIPSL, derive the REAL cloud droplet
314        ! effective radius as re/fl
315
316        fl(i, k) = seuil_neb*(1.-zfice(i,k))
317        re(i, k) = rad_chaud(i, k)*fl(i, k)
318        rel = 0.
319        rei = 0.
320        pclc(i, k) = 0.0
321        pcltau(i, k) = 0.0
322        pclemi(i, k) = 0.0
323
324      ELSE
325
326        ! -- liquid/ice cloud water paths:
327
328        zflwp_var = 1000.*(1.-zfice(i,k))*pqlwp(i, k)/pclc(i, k)*rhodz(i, k)
329        zfiwp_var = 1000.*zfice(i, k)*pqlwp(i, k)/pclc(i, k)*rhodz(i, k)
330
331        ! effective cloud droplet radius (microns) for liquid water clouds:
332        ! For output diagnostics cloud droplet effective radius [um]
333        ! we multiply here with f * xl (fraction of liquid water
334        ! clouds in the grid cell) to avoid problems in the averaging of the
335        ! output.
336        ! In the output of IOIPSL, derive the REAL cloud droplet
337        ! effective radius as re/fl
338
339        fl(i, k) = pclc(i, k)*(1.-zfice(i,k))
340        re(i, k) = rad_chaud(i, k)*fl(i, k)
341
342        rel = rad_chaud(i, k)
343
344        ! for ice clouds: as a function of the ambiant temperature
345        ! [formula used by Iacobellis and Somerville (2000), with an
346        ! asymptotical value of 3.5 microns at T<-81.4 C added to be
347        ! consistent with observations of Heymsfield et al. 1986]:
348        ! 2011/05/24 : rei_min = 3.5 becomes a free PARAMETER as well as
349        ! rei_max=61.29
350
351        tc = t(i, k) - 273.15
352        rei = d_rei_dt*tc + rei_max
353        IF (tc<=-81.4) rei = rei_min
354
355        ! -- cloud optical thickness :
356        ! [for liquid clouds, traditional formula,
357        ! for ice clouds, Ebert & Curry (1992)]
358
359        IF (zflwp_var==0.) rel = 1.
360        IF (zfiwp_var==0. .OR. rei<=0.) rei = 1.
361        pcltau(i, k) = 3.0/2.0*(zflwp_var/rel) + zfiwp_var*(3.448E-03+2.431/ &
362          rei)
363
364        ! -- cloud infrared emissivity:
365        ! [the broadband infrared absorption coefficient is PARAMETERized
366        ! as a function of the effective cld droplet radius]
367        ! Ebert and Curry (1992) formula as used by Kiehl & Zender (1995):
368
369        k_ice = k_ice0 + 1.0/rei
370
371        pclemi(i, k) = 1.0 - exp(-coef_chau*zflwp_var-df*k_ice*zfiwp_var)
372
373      END IF
374
375      reice(i, k) = rei
376
377      xflwp(i) = xflwp(i) + xflwc(i, k)*rhodz(i, k)
378      xfiwp(i) = xfiwp(i) + xfiwc(i, k)*rhodz(i, k)
379
380    END DO
381  END DO
382
383  ! --if cloud droplet radius is fixed, then pcldtaupi=pcltau
384
385  IF (.NOT. ok_cdnc) THEN
386    DO k = 1, klev
387      DO i = 1, klon
388        pcldtaupi(i, k) = pcltau(i, k)
389        reice_pi(i, k) = reice(i, k)
390      END DO
391    END DO
392  END IF
393
394  DO k = 1, klev
395    DO i = 1, klon
396      reliq(i, k) = rad_chaud(i, k)
397      reliq_pi(i, k) = rad_chaud_pi(i, k)
398      reice_pi(i, k) = reice(i, k)
399    END DO
400  END DO
401
402  ! COMPUTE CLOUD LIQUID PATH AND TOTAL CLOUDINESS
403  ! IM cf. CR:test: calcul prenant ou non en compte le recouvrement
404  ! initialisations
405
406  DO i = 1, klon
407    zclear(i) = 1.
408    zcloud(i) = 0.
409    zcloudh(i) = 0.
410    zcloudm(i) = 0.
411    zcloudl(i) = 0.
412    pch(i) = 1.0
413    pcm(i) = 1.0
414    pcl(i) = 1.0
415    pctlwp(i) = 0.0
416  END DO
417
418  ! --calculation of liquid water path
419
420  DO k = klev, 1, -1
421    DO i = 1, klon
422      pctlwp(i) = pctlwp(i) + pqlwp(i, k)*rhodz(i, k)
423    END DO
424  END DO
425
426  ! --calculation of cloud properties with cloud overlap
427
428  IF (novlp==1) THEN
429    DO k = klev, 1, -1
430      DO i = 1, klon
431        zclear(i) = zclear(i)*(1.-max(pclc(i,k),zcloud(i)))/(1.-min(real( &
432          zcloud(i),kind=8),1.-zepsec))
433        pct(i) = 1. - zclear(i)
434        IF (paprs(i,k)<prmhc) THEN
435          pch(i) = pch(i)*(1.-max(pclc(i,k),zcloudh(i)))/(1.-min(real(zcloudh &
436            (i),kind=8),1.-zepsec))
437          zcloudh(i) = pclc(i, k)
438        ELSE IF (paprs(i,k)>=prmhc .AND. paprs(i,k)<prlmc) THEN
439          pcm(i) = pcm(i)*(1.-max(pclc(i,k),zcloudm(i)))/(1.-min(real(zcloudm &
440            (i),kind=8),1.-zepsec))
441          zcloudm(i) = pclc(i, k)
442        ELSE IF (paprs(i,k)>=prlmc) THEN
443          pcl(i) = pcl(i)*(1.-max(pclc(i,k),zcloudl(i)))/(1.-min(real(zcloudl &
444            (i),kind=8),1.-zepsec))
445          zcloudl(i) = pclc(i, k)
446        END IF
447        zcloud(i) = pclc(i, k)
448      END DO
449    END DO
450  ELSE IF (novlp==2) THEN
451    DO k = klev, 1, -1
452      DO i = 1, klon
453        zcloud(i) = max(pclc(i,k), zcloud(i))
454        pct(i) = zcloud(i)
455        IF (paprs(i,k)<prmhc) THEN
456          pch(i) = min(pclc(i,k), pch(i))
457        ELSE IF (paprs(i,k)>=prmhc .AND. paprs(i,k)<prlmc) THEN
458          pcm(i) = min(pclc(i,k), pcm(i))
459        ELSE IF (paprs(i,k)>=prlmc) THEN
460          pcl(i) = min(pclc(i,k), pcl(i))
461        END IF
462      END DO
463    END DO
464  ELSE IF (novlp==3) THEN
465    DO k = klev, 1, -1
466      DO i = 1, klon
467        zclear(i) = zclear(i)*(1.-pclc(i,k))
468        pct(i) = 1 - zclear(i)
469        IF (paprs(i,k)<prmhc) THEN
470          pch(i) = pch(i)*(1.0-pclc(i,k))
471        ELSE IF (paprs(i,k)>=prmhc .AND. paprs(i,k)<prlmc) THEN
472          pcm(i) = pcm(i)*(1.0-pclc(i,k))
473        ELSE IF (paprs(i,k)>=prlmc) THEN
474          pcl(i) = pcl(i)*(1.0-pclc(i,k))
475        END IF
476      END DO
477    END DO
478  END IF
479
480  DO i = 1, klon
481    pch(i) = 1. - pch(i)
482    pcm(i) = 1. - pcm(i)
483    pcl(i) = 1. - pcl(i)
484  END DO
485
486  ! ========================================================
487  ! DIAGNOSTICS CALCULATION FOR CMIP5 PROTOCOL
488  ! ========================================================
489  ! change by Nicolas Yan (LSCE)
490  ! Cloud Droplet Number Concentration (CDNC) : 3D variable
491  ! Fractionnal cover by liquid water cloud (LCC3D) : 3D variable
492  ! Cloud Droplet Number Concentration at top of cloud (CLDNCL) : 2D variable
493  ! Droplet effective radius at top of cloud (REFFCLWTOP) : 2D variable
494  ! Fractionnal cover by liquid water at top of clouds (LCC) : 2D variable
495
496  IF (ok_cdnc) THEN
497
498    DO k = 1, klev
499      DO i = 1, klon
500        phase3d(i, k) = 1 - zfice(i, k)
501        IF (pclc(i,k)<=seuil_neb) THEN
502          lcc3d(i, k) = seuil_neb*phase3d(i, k)
503        ELSE
504          lcc3d(i, k) = pclc(i, k)*phase3d(i, k)
505        END IF
506        scdnc(i, k) = lcc3d(i, k)*cdnc(i, k) ! m-3
507      END DO
508    END DO
509
510    DO i = 1, klon
511      lcc(i) = 0.
512      reffclwtop(i) = 0.
513      cldncl(i) = 0.
514      IF (random .OR. maximum_random) tcc(i) = 1.
515      IF (maximum) tcc(i) = 0.
516    END DO
517
518    DO i = 1, klon
519      DO k = klev - 1, 1, -1 !From TOA down
520
521          ! Test, if the cloud optical depth exceeds the necessary
522          ! threshold:
523
524        IF (pcltau(i,k)>thres_tau .AND. pclc(i,k)>thres_neb) THEN
525
526          IF (maximum) THEN
527            IF (first) THEN
528              WRITE (*, *) 'Hypothese de recouvrement: MAXIMUM'
529              first = .FALSE.
530            END IF
531            flag_max = -1.
532            ftmp(i) = max(tcc(i), pclc(i,k))
533          END IF
534
535          IF (random) THEN
536            IF (first) THEN
537              WRITE (*, *) 'Hypothese de recouvrement: RANDOM'
538              first = .FALSE.
539            END IF
540            flag_max = 1.
541            ftmp(i) = tcc(i)*(1-pclc(i,k))
542          END IF
543
544          IF (maximum_random) THEN
545            IF (first) THEN
546              WRITE (*, *) 'Hypothese de recouvrement: MAXIMUM_ &
547                &                                             &
548                &                                          RANDOM'
549              first = .FALSE.
550            END IF
551            flag_max = 1.
552            ftmp(i) = tcc(i)*(1.-max(pclc(i,k),pclc(i,k+1)))/(1.-min(pclc(i, &
553              k+1),1.-thres_neb))
554          END IF
555          ! Effective radius of cloud droplet at top of cloud (m)
556          reffclwtop(i) = reffclwtop(i) + rad_chaud(i, k)*1.0E-06*phase3d(i, &
557            k)*(tcc(i)-ftmp(i))*flag_max
558          ! CDNC at top of cloud (m-3)
559          cldncl(i) = cldncl(i) + cdnc(i, k)*phase3d(i, k)*(tcc(i)-ftmp(i))* &
560            flag_max
561          ! Liquid Cloud Content at top of cloud
562          lcc(i) = lcc(i) + phase3d(i, k)*(tcc(i)-ftmp(i))*flag_max
563          ! Total Cloud Content at top of cloud
564          tcc(i) = ftmp(i)
565
566        END IF ! is there a visible, not-too-small cloud?
567      END DO ! loop over k
568
569      IF (random .OR. maximum_random) tcc(i) = 1. - tcc(i)
570
571    END DO ! loop over i
572
573    ! ! Convective and Stratiform Cloud Droplet Effective Radius (REFFCLWC
574    ! REFFCLWS)
575    DO i = 1, klon
576      DO k = 1, klev
577        ! Weight to be used for outputs: eau_liquide*couverture nuageuse
578        lcc3dcon(i, k) = rnebcon(i, k)*phase3d(i, k)*clwcon(i, k) ! eau liquide convective
579        lcc3dstra(i, k) = pclc(i, k)*pqlwp(i, k)*phase3d(i, k)
580        lcc3dstra(i, k) = lcc3dstra(i, k) - lcc3dcon(i, k) ! eau liquide stratiforme
581        lcc3dstra(i, k) = max(lcc3dstra(i,k), 0.0)
582        ! Compute cloud droplet radius as above in meter
583        radius = 1.1*((pqlwp(i,k)*pplay(i,k)/(rd*t(i,k)))/(4./3*rpi*1000.* &
584          cdnc(i,k)))**(1./3.)
585        radius = max(radius, 5.E-6)
586        ! Convective Cloud Droplet Effective Radius (REFFCLWC) : variable 3D
587        reffclwc(i, k) = radius
588        reffclwc(i, k) = reffclwc(i, k)*lcc3dcon(i, k)
589        ! Stratiform Cloud Droplet Effective Radius (REFFCLWS) : variable 3D
590        reffclws(i, k) = radius
591        reffclws(i, k) = reffclws(i, k)*lcc3dstra(i, k)
592      END DO !klev
593    END DO !klon
594
595    ! Column Integrated Cloud Droplet Number (CLDNVI) : variable 2D
596
597    DO i = 1, klon
598      cldnvi(i) = 0.
599      lcc_integrat(i) = 0.
600      height(i) = 0.
601      DO k = 1, klev
602        cldnvi(i) = cldnvi(i) + cdnc(i, k)*lcc3d(i, k)*dh(i, k)
603        lcc_integrat(i) = lcc_integrat(i) + lcc3d(i, k)*dh(i, k)
604        height(i) = height(i) + dh(i, k)
605      END DO ! klev
606      lcc_integrat(i) = lcc_integrat(i)/height(i)
607      IF (lcc_integrat(i)<=1.0E-03) THEN
608        cldnvi(i) = cldnvi(i)*lcc(i)/seuil_neb
609      ELSE
610        cldnvi(i) = cldnvi(i)*lcc(i)/lcc_integrat(i)
611      END IF
612    END DO ! klon
613
614    DO i = 1, klon
615      DO k = 1, klev
616        IF (scdnc(i,k)<=0.0) scdnc(i, k) = 0.0
617        IF (reffclws(i,k)<=0.0) reffclws(i, k) = 0.0
618        IF (reffclwc(i,k)<=0.0) reffclwc(i, k) = 0.0
619        IF (lcc3d(i,k)<=0.0) lcc3d(i, k) = 0.0
620        IF (lcc3dcon(i,k)<=0.0) lcc3dcon(i, k) = 0.0
621        IF (lcc3dstra(i,k)<=0.0) lcc3dstra(i, k) = 0.0
622      END DO
623      IF (reffclwtop(i)<=0.0) reffclwtop(i) = 0.0
624      IF (cldncl(i)<=0.0) cldncl(i) = 0.0
625      IF (cldnvi(i)<=0.0) cldnvi(i) = 0.0
626      IF (lcc(i)<=0.0) lcc(i) = 0.0
627    END DO
628
629  END IF !ok_cdnc
630
631  RETURN
632
633END SUBROUTINE newmicro
Note: See TracBrowser for help on using the repository browser.