source: LMDZ6/branches/contrails/libf/phylmd/rrtm/rrtm_ecrt_140gp.F90 @ 5452

Last change on this file since 5452 was 5294, checked in by Laurent Fairhead, 2 months ago

Keeping clesphys.h was not the right solution
LF
  • 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:keywords set to Author Date Id Revi
File size: 13.4 KB
Line 
1!
2! $Id: rrtm_ecrt_140gp.F90 5294 2024-10-29 18:35:00Z aborella $
3!
4!****************** SUBROUTINE RRTM_ECRT_140GP **************************
5
6SUBROUTINE RRTM_ECRT_140GP &
7 & ( K_IPLON, klon , klev, kcld,&
8 & paer , paph , pap,&
9 & pts  , pth  , pt,&
10 & P_ZEMIS, P_ZEMIW,&
11 & pq   , pcco2, pozn, pcldf, ptaucld, ptclear,&
12 & P_CLDFRAC,P_TAUCLD,&
13 & PTAU_LW,&
14 & P_COLDRY,P_WKL,P_WX,&
15 & P_TAUAERL,PAVEL,P_TAVEL,PZ,P_TZ,P_TBOUND,K_NLAYERS,P_SEMISS,K_IREFLECT ) 
16
17!     Reformatted for F90 by JJMorcrette, ECMWF, 980714
18
19!     Read in atmospheric profile from ECMWF radiation code, and prepare it
20!     for use in RRTM.  Set other RRTM input parameters.  Values are passed
21!     back through existing RRTM arrays and commons.
22
23!- Modifications
24
25!     2000-05-15 Deborah Salmond  Speed-up
26
27USE PARKIND1  ,ONLY : JPIM     ,JPRB
28USE YOMHOOK   ,ONLY : LHOOK,   DR_HOOK
29
30USE PARRRTM  , ONLY : JPBAND   ,JPXSEC   ,JPLAY   ,&
31 & JPINPX 
32USE YOERAD   , ONLY : NLW      ,NOVLP
33!MPL/IM 20160915 on prend GES de phylmd USE YOERDI   , ONLY :    RCH4     ,RN2O    ,RCFC11  ,RCFC12
34USE YOESW    , ONLY : RAER
35! Temporary fix waiting for cleaner interface (or not)
36USE clesphys_mod_h, ONLY: NSW, rcfc11, rcfc12, rch4, rn2o
37
38!------------------------------Arguments--------------------------------
39
40IMPLICIT NONE
41
42
43INTEGER(KIND=JPIM),INTENT(IN)    :: KLON! Number of atmospheres (longitudes)
44INTEGER(KIND=JPIM),INTENT(IN)    :: KLEV! Number of atmospheric layers
45INTEGER(KIND=JPIM),INTENT(IN)    :: K_IPLON
46INTEGER(KIND=JPIM),INTENT(OUT)   :: KCLD
47REAL(KIND=JPRB)   ,INTENT(IN)    :: PAER(KLON,6,KLEV) ! Aerosol optical thickness
48REAL(KIND=JPRB)   ,INTENT(IN)    :: PAPH(KLON,KLEV+1) ! Interface pressures (Pa)
49REAL(KIND=JPRB)   ,INTENT(IN)    :: PAP(KLON,KLEV) ! Layer pressures (Pa)
50REAL(KIND=JPRB)   ,INTENT(IN)    :: PTS(KLON) ! Surface temperature (K)
51REAL(KIND=JPRB)   ,INTENT(IN)    :: PTH(KLON,KLEV+1) ! Interface temperatures (K)
52REAL(KIND=JPRB)   ,INTENT(IN)    :: PT(KLON,KLEV) ! Layer temperature (K)
53REAL(KIND=JPRB)   ,INTENT(IN)    :: P_ZEMIS(KLON) ! Non-window surface emissivity
54REAL(KIND=JPRB)   ,INTENT(IN)    :: P_ZEMIW(KLON) ! Window surface emissivity
55REAL(KIND=JPRB)   ,INTENT(IN)    :: PQ(KLON,KLEV) ! H2O specific humidity (mmr)
56REAL(KIND=JPRB)   ,INTENT(IN)    :: PCCO2 ! CO2 mass mixing ratio
57REAL(KIND=JPRB)   ,INTENT(IN)    :: POZN(KLON,KLEV) ! O3 mass mixing ratio
58REAL(KIND=JPRB)   ,INTENT(IN)    :: PCLDF(KLON,KLEV) ! Cloud fraction
59REAL(KIND=JPRB)   ,INTENT(IN)    :: PTAUCLD(KLON,KLEV,JPBAND) ! Cloud optical depth
60!--C.Kleinschmitt
61REAL(KIND=JPRB)   ,INTENT(IN)    :: PTAU_LW(KLON,KLEV,NLW) ! LW Optical depth of aerosols 
62!--end
63REAL(KIND=JPRB)   ,INTENT(OUT)   :: PTCLEAR
64REAL(KIND=JPRB)   ,INTENT(OUT)   :: P_CLDFRAC(JPLAY) ! Cloud fraction
65REAL(KIND=JPRB)   ,INTENT(OUT)   :: P_TAUCLD(JPLAY,JPBAND) ! Spectral optical thickness
66REAL(KIND=JPRB)   ,INTENT(OUT)   :: P_COLDRY(JPLAY)
67REAL(KIND=JPRB)   ,INTENT(OUT)   :: P_WKL(JPINPX,JPLAY)
68REAL(KIND=JPRB)   ,INTENT(OUT)   :: P_WX(JPXSEC,JPLAY) ! Amount of trace gases
69REAL(KIND=JPRB)   ,INTENT(OUT)   :: P_TAUAERL(JPLAY,JPBAND)
70REAL(KIND=JPRB)   ,INTENT(OUT)   :: PAVEL(JPLAY)
71REAL(KIND=JPRB)   ,INTENT(OUT)   :: P_TAVEL(JPLAY)
72REAL(KIND=JPRB)   ,INTENT(OUT)   :: PZ(0:JPLAY)
73REAL(KIND=JPRB)   ,INTENT(OUT)   :: P_TZ(0:JPLAY)
74REAL(KIND=JPRB)   ,INTENT(OUT)   :: P_TBOUND
75INTEGER(KIND=JPIM),INTENT(OUT)   :: K_NLAYERS
76REAL(KIND=JPRB)   ,INTENT(OUT)   :: P_SEMISS(JPBAND)
77INTEGER(KIND=JPIM),INTENT(OUT)   :: K_IREFLECT
78!      real rch4                       ! CH4 mass mixing ratio
79!      real rn2o                       ! N2O mass mixing ratio
80!      real rcfc11                     ! CFC11 mass mixing ratio
81!      real rcfc12                     ! CFC12 mass mixing ratio
82!- from AER
83!- from PROFILE             
84!- from SURFACE             
85REAL(KIND=JPRB) :: ztauaer(5)
86REAL(KIND=JPRB) :: zc1j(0:klev)               ! total cloud from top and level k
87REAL(KIND=JPRB) :: Z_AMD                  ! Effective molecular weight of dry air (g/mol)
88REAL(KIND=JPRB) :: Z_AMW                  ! Molecular weight of water vapor (g/mol)
89REAL(KIND=JPRB) :: Z_AMCO2                ! Molecular weight of carbon dioxide (g/mol)
90REAL(KIND=JPRB) :: Z_AMO                  ! Molecular weight of ozone (g/mol)
91REAL(KIND=JPRB) :: Z_AMCH4                ! Molecular weight of methane (g/mol)
92REAL(KIND=JPRB) :: Z_AMN2O                ! Molecular weight of nitrous oxide (g/mol)
93REAL(KIND=JPRB) :: Z_AMC11                ! Molecular weight of CFC11 (g/mol) - CFCL3
94REAL(KIND=JPRB) :: Z_AMC12                ! Molecular weight of CFC12 (g/mol) - CF2CL2
95REAL(KIND=JPRB) :: Z_AVGDRO               ! Avogadro's number (molecules/mole)
96REAL(KIND=JPRB) :: Z_GRAVIT               ! Gravitational acceleration (cm/sec2)
97
98! Atomic weights for conversion from mass to volume mixing ratios; these
99!  are the same values used in ECRT to assure accurate conversion to vmr
100data Z_AMD   /  28.970_JPRB    /
101data Z_AMW   /  18.0154_JPRB   /
102data Z_AMCO2 /  44.011_JPRB    /
103data Z_AMO   /  47.9982_JPRB   /
104data Z_AMCH4 /  16.043_JPRB    /
105data Z_AMN2O /  44.013_JPRB    /
106data Z_AMC11 / 137.3686_JPRB   /
107data Z_AMC12 / 120.9140_JPRB   /
108data Z_AVGDRO/ 6.02214E23_JPRB /
109data Z_GRAVIT/ 9.80665E02_JPRB /
110
111INTEGER(KIND=JPIM) :: IATM, IMOL, IXMAX, J1, J2, JAE, JB, JK, JL, I_L
112INTEGER(KIND=JPIM) :: I_NMOL, I_NXMOL
113
114REAL(KIND=JPRB) :: Z_AMM, ZCLDLY, ZCLEAR, ZCLOUD, ZEPSEC
115REAL(KIND=JPRB) :: ZHOOK_HANDLE
116
117!MPL/IM 20160915 on prend GES de phylmd
118!!#include "clesphys.h"
119! ***
120
121! *** mji
122! Initialize all molecular amounts and aerosol optical depths to zero here,
123! then pass ECRT amounts into RRTM arrays below.
124
125!      DATA ZWKL /MAXPRDW*0.0/
126!      DATA ZWX  /MAXPROD*0.0/
127!      DATA KREFLECT /0/
128
129! Activate cross section molecules:
130!     NXMOL     - number of cross-sections input by user
131!     IXINDX(I) - index of cross-section molecule corresponding to Ith
132!                 cross-section specified by user
133!                 = 0 -- not allowed in RRTM
134!                 = 1 -- CCL4
135!                 = 2 -- CFC11
136!                 = 3 -- CFC12
137!                 = 4 -- CFC22
138!      DATA KXMOL  /2/
139!      DATA KXINDX /0,2,3,0,31*0/
140
141!      IREFLECT=KREFLECT
142!      NXMOL=KXMOL
143
144IF (LHOOK) CALL DR_HOOK('RRTM_ECRT_140GP',0,ZHOOK_HANDLE)
145K_IREFLECT=0
146I_NXMOL=2
147
148DO J1=1,35
149! IXINDX(J1)=0
150  DO J2=1,KLEV
151    P_WKL(J1,J2)=0.0_JPRB
152  ENDDO
153ENDDO
154!IXINDX(2)=2
155!IXINDX(3)=3
156
157!     Set parameters needed for RRTM execution:
158IATM    = 0
159!      IXSECT  = 1
160!      NUMANGS = 0
161!      IOUT    = -1
162IXMAX   = 4
163
164!     Bands 6,7,8 are considered the 'window' and allowed to have a
165!     different surface emissivity (as in ECMWF).  Eli wrote this part....
166P_SEMISS(1)  = P_ZEMIS(K_IPLON)
167P_SEMISS(2)  = P_ZEMIS(K_IPLON)
168P_SEMISS(3)  = P_ZEMIS(K_IPLON)
169P_SEMISS(4)  = P_ZEMIS(K_IPLON)
170P_SEMISS(5)  = P_ZEMIS(K_IPLON)
171P_SEMISS(6)  = P_ZEMIW(K_IPLON)
172P_SEMISS(7)  = P_ZEMIW(K_IPLON)
173P_SEMISS(8)  = P_ZEMIW(K_IPLON)
174P_SEMISS(9)  = P_ZEMIS(K_IPLON)
175P_SEMISS(10) = P_ZEMIS(K_IPLON)
176P_SEMISS(11) = P_ZEMIS(K_IPLON)
177P_SEMISS(12) = P_ZEMIS(K_IPLON)
178P_SEMISS(13) = P_ZEMIS(K_IPLON)
179P_SEMISS(14) = P_ZEMIS(K_IPLON)
180P_SEMISS(15) = P_ZEMIS(K_IPLON)
181P_SEMISS(16) = P_ZEMIS(K_IPLON)
182
183!     Set surface temperature. 
184
185P_TBOUND = pts(K_IPLON)
186
187!     Install ECRT arrays into RRTM arrays for pressure, temperature,
188!     and molecular amounts.  Pressures are converted from Pascals
189!     (ECRT) to mb (RRTM).  H2O, CO2, O3 and trace gas amounts are
190!     converted from mass mixing ratio to volume mixing ratio.  CO2
191!     converted with same dry air and CO2 molecular weights used in
192!     ECRT to assure correct conversion back to the proper CO2 vmr.
193!     The dry air column COLDRY (in molec/cm2) is calculated from
194!     the level pressures PZ (in mb) based on the hydrostatic equation
195!     and includes a correction to account for H2O in the layer.  The
196!     molecular weight of moist air (amm) is calculated for each layer.
197!     Note: RRTM levels count from bottom to top, while the ECRT input
198!     variables count from the top down and must be reversed here.
199
200K_NLAYERS = klev
201I_NMOL = 6
202PZ(0) = paph(K_IPLON,klev+1)/100._JPRB
203P_TZ(0) = pth(K_IPLON,klev+1)
204DO I_L = 1, KLEV
205  PAVEL(I_L) = pap(K_IPLON,KLEV-I_L+1)/100._JPRB
206  P_TAVEL(I_L) = pt(K_IPLON,KLEV-I_L+1)
207  PZ(I_L) = paph(K_IPLON,KLEV-I_L+1)/100._JPRB
208  P_TZ(I_L) = pth(K_IPLON,KLEV-I_L+1)
209  P_WKL(1,I_L) = pq(K_IPLON,KLEV-I_L+1)*Z_AMD/Z_AMW
210  P_WKL(2,I_L) = pcco2*Z_AMD/Z_AMCO2
211  P_WKL(3,I_L) = pozn(K_IPLON,KLEV-I_L+1)*Z_AMD/Z_AMO
212  P_WKL(4,I_L) = rn2o*Z_AMD/Z_AMN2O
213  P_WKL(6,I_L) = rch4*Z_AMD/Z_AMCH4
214  Z_AMM = (1-P_WKL(1,I_L))*Z_AMD + P_WKL(1,I_L)*Z_AMW
215  P_COLDRY(I_L) = (PZ(I_L-1)-PZ(I_L))*1.E3_JPRB*Z_AVGDRO/(Z_GRAVIT*Z_AMM*(1+P_WKL(1,I_L)))
216ENDDO
217
218!- Fill RRTM aerosol arrays with operational ECMWF aerosols,
219!  do the mixing and distribute over the 16 spectral intervals
220
221DO I_L=1,KLEV
222  JK=KLEV-I_L+1
223!       DO JAE=1,5
224  JAE=1
225  ZTAUAER(JAE) =&
226   & RAER(JAE,1)*PAER(K_IPLON,1,JK)+RAER(JAE,2)*PAER(K_IPLON,2,JK)&
227   & +RAER(JAE,3)*PAER(K_IPLON,3,JK)+RAER(JAE,4)*PAER(K_IPLON,4,JK)&
228   & +RAER(JAE,5)*PAER(K_IPLON,5,JK)+RAER(JAE,6)*PAER(K_IPLON,6,JK) 
229  P_TAUAERL(I_L, 1)=ZTAUAER(1)
230  P_TAUAERL(I_L, 2)=ZTAUAER(1)
231  JAE=2
232  ZTAUAER(JAE) =&
233   & RAER(JAE,1)*PAER(K_IPLON,1,JK)+RAER(JAE,2)*PAER(K_IPLON,2,JK)&
234   & +RAER(JAE,3)*PAER(K_IPLON,3,JK)+RAER(JAE,4)*PAER(K_IPLON,4,JK)&
235   & +RAER(JAE,5)*PAER(K_IPLON,5,JK)+RAER(JAE,6)*PAER(K_IPLON,6,JK) 
236  P_TAUAERL(I_L, 3)=ZTAUAER(2)
237  P_TAUAERL(I_L, 4)=ZTAUAER(2)
238  P_TAUAERL(I_L, 5)=ZTAUAER(2)
239  JAE=3
240  ZTAUAER(JAE) =&
241   & RAER(JAE,1)*PAER(K_IPLON,1,JK)+RAER(JAE,2)*PAER(K_IPLON,2,JK)&
242   & +RAER(JAE,3)*PAER(K_IPLON,3,JK)+RAER(JAE,4)*PAER(K_IPLON,4,JK)&
243   & +RAER(JAE,5)*PAER(K_IPLON,5,JK)+RAER(JAE,6)*PAER(K_IPLON,6,JK) 
244  P_TAUAERL(I_L, 6)=ZTAUAER(3)
245  P_TAUAERL(I_L, 8)=ZTAUAER(3)
246  P_TAUAERL(I_L, 9)=ZTAUAER(3)
247  JAE=4
248  ZTAUAER(JAE) =&
249   & RAER(JAE,1)*PAER(K_IPLON,1,JK)+RAER(JAE,2)*PAER(K_IPLON,2,JK)&
250   & +RAER(JAE,3)*PAER(K_IPLON,3,JK)+RAER(JAE,4)*PAER(K_IPLON,4,JK)&
251   & +RAER(JAE,5)*PAER(K_IPLON,5,JK)+RAER(JAE,6)*PAER(K_IPLON,6,JK) 
252  P_TAUAERL(I_L, 7)=ZTAUAER(4)
253  JAE=5
254  ZTAUAER(JAE) =&
255   & RAER(JAE,1)*PAER(K_IPLON,1,JK)+RAER(JAE,2)*PAER(K_IPLON,2,JK)&
256   & +RAER(JAE,3)*PAER(K_IPLON,3,JK)+RAER(JAE,4)*PAER(K_IPLON,4,JK)&
257   & +RAER(JAE,5)*PAER(K_IPLON,5,JK)+RAER(JAE,6)*PAER(K_IPLON,6,JK) 
258!       END DO
259  P_TAUAERL(I_L,10)=ZTAUAER(5)
260  P_TAUAERL(I_L,11)=ZTAUAER(5)
261  P_TAUAERL(I_L,12)=ZTAUAER(5)
262  P_TAUAERL(I_L,13)=ZTAUAER(5)
263  P_TAUAERL(I_L,14)=ZTAUAER(5)
264  P_TAUAERL(I_L,15)=ZTAUAER(5)
265  P_TAUAERL(I_L,16)=ZTAUAER(5)
266ENDDO
267!--Use LW AOD from own Mie calculations (C. Kleinschmitt)
268DO I_L=1,KLEV
269  JK=KLEV-I_L+1
270  DO JAE=1, NLW
271    P_TAUAERL(I_L,JAE) = MAX( PTAU_LW(K_IPLON, JK, JAE), 1e-30 )
272  ENDDO
273ENDDO
274!--end C. Kleinschmitt
275
276DO J2=1,KLEV
277  DO J1=1,JPXSEC
278    P_WX(J1,J2)=0.0_JPRB
279  ENDDO
280ENDDO
281
282DO I_L = 1, KLEV
283!- Set cross section molecule amounts from ECRT; convert to vmr
284  P_WX(2,I_L) = rcfc11*Z_AMD/Z_AMC11
285  P_WX(3,I_L) = rcfc12*Z_AMD/Z_AMC12
286  P_WX(2,I_L) = P_COLDRY(I_L) * P_WX(2,I_L) * 1.E-20_JPRB
287  P_WX(3,I_L) = P_COLDRY(I_L) * P_WX(3,I_L) * 1.E-20_JPRB
288
289!- Here, all molecules in WKL and WX are in volume mixing ratio; convert to
290!  molec/cm2 based on COLDRY for use in RRTM
291
292  DO IMOL = 1, I_NMOL
293    P_WKL(IMOL,I_L) = P_COLDRY(I_L) * P_WKL(IMOL,I_L)
294  ENDDO 
295 
296! DO IX = 1,JPXSEC
297! IF (IXINDX(IX)  /=  0) THEN
298!     WX(IXINDX(IX),L) = COLDRY(L) * WX(IX,L) * 1.E-20_JPRB
299! ENDIF
300! END DO 
301
302ENDDO
303
304!- Approximate treatment for various cloud overlaps
305ZCLEAR=1.0_JPRB
306ZCLOUD=0.0_JPRB
307ZC1J(0)=0.0_JPRB
308ZEPSEC=1.E-03_JPRB
309JL=K_IPLON
310
311!++MODIFCODE
312IF ((NOVLP == 1).OR.(NOVLP ==6).OR.(NOVLP ==8)) THEN
313!--MODIFCODE
314
315  DO JK=1,KLEV
316    IF (pcldf(JL,JK) > ZEPSEC) THEN
317      ZCLDLY=pcldf(JL,JK)
318      ZCLEAR=ZCLEAR &
319       & *(1.0_JPRB-MAX( ZCLDLY , ZCLOUD ))&
320       & /(1.0_JPRB-MIN( ZCLOUD , 1.0_JPRB-ZEPSEC )) 
321      ZCLOUD = ZCLDLY
322      ZC1J(JK)= 1.0_JPRB - ZCLEAR
323    ELSE
324      ZCLDLY=0.0_JPRB
325      ZCLEAR=ZCLEAR &
326       & *(1.0_JPRB-MAX( ZCLDLY , ZCLOUD ))&
327       & /(1.0_JPRB-MIN( ZCLOUD , 1.0_JPRB-ZEPSEC )) 
328      ZCLOUD = ZCLDLY
329      ZC1J(JK)= 1.0_JPRB - ZCLEAR
330    ENDIF
331  ENDDO
332
333!++MODIFCODE
334ELSEIF ((NOVLP == 2).OR.(NOVLP ==7)) THEN
335!--MODIFCODE
336
337  DO JK=1,KLEV
338    IF (pcldf(JL,JK) > ZEPSEC) THEN
339      ZCLDLY=pcldf(JL,JK)
340      ZCLOUD = MAX( ZCLDLY , ZCLOUD )
341      ZC1J(JK) = ZCLOUD
342    ELSE
343      ZCLDLY=0.0_JPRB
344      ZCLOUD = MAX( ZCLDLY , ZCLOUD )
345      ZC1J(JK) = ZCLOUD
346    ENDIF
347  ENDDO
348
349!++MODIFCODE
350ELSEIF ((NOVLP == 3).OR.(NOVLP ==5)) THEN
351!--MODIFCODE
352
353  DO JK=1,KLEV
354    IF (pcldf(JL,JK) > ZEPSEC) THEN
355      ZCLDLY=pcldf(JL,JK)
356      ZCLEAR = ZCLEAR * (1.0_JPRB-ZCLDLY)
357      ZCLOUD = 1.0_JPRB - ZCLEAR
358      ZC1J(JK) = ZCLOUD
359    ELSE
360      ZCLDLY=0.0_JPRB
361      ZCLEAR = ZCLEAR * (1.0_JPRB-ZCLDLY)
362      ZCLOUD = 1.0_JPRB - ZCLEAR
363      ZC1J(JK) = ZCLOUD
364    ENDIF
365  ENDDO
366
367ELSEIF (NOVLP == 4) THEN
368
369ENDIF
370PTCLEAR=1.0_JPRB-ZC1J(KLEV)
371
372! Transfer cloud fraction and cloud optical depth to RRTM arrays;
373! invert array index for pcldf to go from bottom to top for RRTM
374
375!- clear-sky column
376IF (PTCLEAR  >  1.0_JPRB-ZEPSEC) THEN
377  KCLD=0
378  DO I_L = 1, KLEV
379    P_CLDFRAC(I_L) = 0.0_JPRB
380  ENDDO
381  DO JB=1,JPBAND
382    DO I_L=1,KLEV
383      P_TAUCLD(I_L,JB) = 0.0_JPRB
384    ENDDO
385  ENDDO
386
387ELSE
388
389!- cloudy column
390!   The diffusivity factor (Savijarvi, 1997) on the cloud optical
391!   thickness TAUCLD has already been applied in RADLSW
392
393  KCLD=1
394  DO I_L=1,KLEV
395    P_CLDFRAC(I_L) = pcldf(K_IPLON,I_L)
396  ENDDO
397  DO JB=1,JPBAND
398    DO I_L=1,KLEV
399      P_TAUCLD(I_L,JB) = ptaucld(K_IPLON,I_L,JB)
400    ENDDO
401  ENDDO
402
403ENDIF
404
405!     ------------------------------------------------------------------
406 
407IF (LHOOK) CALL DR_HOOK('RRTM_ECRT_140GP',1,ZHOOK_HANDLE)
408END SUBROUTINE RRTM_ECRT_140GP
Note: See TracBrowser for help on using the repository browser.