source: LMDZ6/branches/IPSLCM6.0.15/libf/phylmd/newmicro.F90 @ 3248

Last change on this file since 3248 was 3248, checked in by Laurent Fairhead, 6 years ago

Inclusion of r3245 from trunk

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