source: LMDZ6/branches/contrails/libf/phylmd/ecrad/ifsrrtm/rrtm_taumol12.F90 @ 5440

Last change on this file since 5440 was 4773, checked in by idelkadi, 12 months ago
  • Update of Ecrad in LMDZ The same organization of the Ecrad offline version is retained in order to facilitate the updating of Ecrad in LMDZ and the comparison between online and offline results. version 1.6.1 of Ecrad (​https://github.com/lguez/ecrad.git)
  • Implementation of the double call of Ecrad in LMDZ


File size: 10.8 KB
Line 
1!----------------------------------------------------------------------------
2SUBROUTINE RRTM_TAUMOL12 (KIDIA,KFDIA,KLEV,P_TAU,&
3 & P_TAUAERL,P_FAC00,P_FAC01,P_FAC10,P_FAC11,P_FORFAC,P_FORFRAC,K_INDFOR,K_JP,K_JT,K_JT1,P_ONEMINUS,&
4 & P_COLH2O,P_COLCO2,K_LAYTROP,P_SELFFAC,P_SELFFRAC,K_INDSELF,PFRAC, & 
5 & PRAT_H2OCO2, PRAT_H2OCO2_1)
6
7!     BAND 12:  1800-2080 cm-1 (low - H2O,CO2; high - nothing)
8
9!     AUTHOR.
10!     -------
11!      JJMorcrette, ECMWF
12
13!     MODIFICATIONS.
14!     --------------
15!      M.Hamrud      01-Oct-2003 CY28 Cleaning
16!      NEC           25-Oct-2007 Optimisations
17!      JJMorcrette 20110613 flexible number of g-points
18!      ABozzo updated to rrtmg v4.85
19! ---------------------------------------------------------------------------
20
21USE PARKIND1  ,ONLY : JPIM     ,JPRB
22USE YOMHOOK   ,ONLY : LHOOK, DR_HOOK, JPHOOK
23
24USE PARRRTM  , ONLY : JPBAND
25USE YOERRTM  , ONLY : JPGPT  ,NG12 ,NGS11
26USE YOERRTWN , ONLY :      NSPA   
27USE YOERRTA12, ONLY : ABSA   ,FRACREFA,SELFREF,FORREF
28USE YOERRTRF, ONLY : CHI_MLS
29
30IMPLICIT NONE
31
32INTEGER(KIND=JPIM),INTENT(IN)    :: KIDIA
33INTEGER(KIND=JPIM),INTENT(IN)    :: KFDIA
34INTEGER(KIND=JPIM),INTENT(IN)    :: KLEV
35REAL(KIND=JPRB)   ,INTENT(OUT)   :: P_TAU(KIDIA:KFDIA,JPGPT,KLEV)
36REAL(KIND=JPRB)   ,INTENT(IN)    :: P_TAUAERL(KIDIA:KFDIA,KLEV,JPBAND)
37REAL(KIND=JPRB)   ,INTENT(IN)    :: P_FAC00(KIDIA:KFDIA,KLEV)
38REAL(KIND=JPRB)   ,INTENT(IN)    :: P_FAC01(KIDIA:KFDIA,KLEV)
39REAL(KIND=JPRB)   ,INTENT(IN)    :: P_FAC10(KIDIA:KFDIA,KLEV)
40REAL(KIND=JPRB)   ,INTENT(IN)    :: P_FAC11(KIDIA:KFDIA,KLEV)
41INTEGER(KIND=JPIM),INTENT(IN)    :: K_JP(KIDIA:KFDIA,KLEV)
42INTEGER(KIND=JPIM),INTENT(IN)    :: K_JT(KIDIA:KFDIA,KLEV)
43INTEGER(KIND=JPIM),INTENT(IN)    :: K_JT1(KIDIA:KFDIA,KLEV)
44REAL(KIND=JPRB)   ,INTENT(IN)    :: P_ONEMINUS
45REAL(KIND=JPRB)   ,INTENT(IN)    :: P_COLH2O(KIDIA:KFDIA,KLEV)
46REAL(KIND=JPRB)   ,INTENT(IN)    :: P_COLCO2(KIDIA:KFDIA,KLEV)
47INTEGER(KIND=JPIM),INTENT(IN)    :: K_LAYTROP(KIDIA:KFDIA)
48REAL(KIND=JPRB)   ,INTENT(IN)    :: P_SELFFAC(KIDIA:KFDIA,KLEV)
49REAL(KIND=JPRB)   ,INTENT(IN)    :: P_SELFFRAC(KIDIA:KFDIA,KLEV)
50INTEGER(KIND=JPIM),INTENT(IN)    :: K_INDSELF(KIDIA:KFDIA,KLEV)
51REAL(KIND=JPRB)   ,INTENT(OUT)   :: PFRAC(KIDIA:KFDIA,JPGPT,KLEV)
52
53REAL(KIND=JPRB)   ,INTENT(IN)   :: PRAT_H2OCO2(KIDIA:KFDIA,KLEV)
54REAL(KIND=JPRB)   ,INTENT(IN)   :: PRAT_H2OCO2_1(KIDIA:KFDIA,KLEV)
55INTEGER(KIND=JPIM),INTENT(IN)   :: K_INDFOR(KIDIA:KFDIA,KLEV)
56REAL(KIND=JPRB)   ,INTENT(IN)   :: P_FORFRAC(KIDIA:KFDIA,KLEV)
57REAL(KIND=JPRB)   ,INTENT(IN)   :: P_FORFAC(KIDIA:KFDIA,KLEV)
58
59
60! ---------------------------------------------------------------------------
61
62REAL(KIND=JPRB) :: Z_SPECCOMB(KLEV),Z_SPECCOMB1(KLEV),Z_SPECCOMB_PLANCK(KLEV)
63INTEGER(KIND=JPIM) :: IND0(KLEV),IND1(KLEV),INDS(KLEV),INDF(KLEV)
64
65INTEGER(KIND=JPIM) :: IG, JS, JLAY,JS1, JPL
66INTEGER(KIND=JPIM) :: JLON
67
68! REAL(KIND=JPRB) :: Z_FAC000, Z_FAC001, Z_FAC010, Z_FAC011, Z_FAC100, Z_FAC101,&
69!  & Z_FAC110, Z_FAC111
70REAL(KIND=JPRB) :: Z_FS, Z_SPECMULT, Z_SPECPARM,  &
71& Z_FS1, Z_SPECMULT1, Z_SPECPARM1, &
72& Z_FPL, Z_SPECMULT_PLANCK, Z_SPECPARM_PLANCK
73
74REAL(KIND=JPRB) ::  Z_FAC000, Z_FAC100, Z_FAC200,&
75 & Z_FAC010, Z_FAC110, Z_FAC210, &
76 & Z_FAC001, Z_FAC101, Z_FAC201, &
77 & Z_FAC011, Z_FAC111, Z_FAC211
78REAL(KIND=JPRB) :: ZP, ZP4, ZFK0, ZFK1, ZFK2
79REAL(KIND=JPRB) :: ZTAUFOR,ZTAUSELF,ZTAU_MAJOR,ZTAU_MAJOR1
80REAL(KIND=JPRB) :: ZREFRAT_PLANCK_A
81
82REAL(KIND=JPHOOK) :: ZHOOK_HANDLE
83
84!  ----------------------------------------------------------
85
86IF (LHOOK) CALL DR_HOOK('RRTM_TAUMOL12',0,ZHOOK_HANDLE)
87
88! Calculate reference ratio to be used in calculation of Planck
89! fraction in lower/upper atmosphere.
90
91! P =   174.164 mb
92      ZREFRAT_PLANCK_A = CHI_MLS(1,10)/CHI_MLS(2,10)
93
94! Compute the optical depth by interpolating in ln(pressure),
95! temperature, and appropriate species.  Below laytrop, the water
96! vapor self-continuum adn foreign continuum is interpolated
97! (in temperature) separately. 
98
99DO JLAY = 1, KLEV
100  DO JLON = KIDIA, KFDIA
101    IF (JLAY <= K_LAYTROP(JLON)) THEN
102      Z_SPECCOMB(JLAY) = P_COLH2O(JLON,JLAY) + PRAT_H2OCO2(JLON,JLAY)*P_COLCO2(JLON,JLAY)
103      Z_SPECPARM = P_COLH2O(JLON,JLAY)/Z_SPECCOMB(JLAY)
104      Z_SPECPARM=MIN(P_ONEMINUS,Z_SPECPARM)
105      Z_SPECMULT = 8._JPRB*(Z_SPECPARM)
106      JS = 1 + INT(Z_SPECMULT)
107      Z_FS = MOD(Z_SPECMULT,1.0_JPRB)
108
109      Z_SPECCOMB1(JLAY) = P_COLH2O(JLON,JLAY) + PRAT_H2OCO2_1(JLON,JLAY)*P_COLCO2(JLON,JLAY)
110      Z_SPECPARM1 = P_COLH2O(JLON,JLAY)/Z_SPECCOMB1(JLAY)
111      IF (Z_SPECPARM1 >= P_ONEMINUS) Z_SPECPARM1 = P_ONEMINUS
112      Z_SPECMULT1 = 8._JPRB*(Z_SPECPARM1)
113      JS1 = 1 + INT(Z_SPECMULT1)
114      Z_FS1 = MOD(Z_SPECMULT1,1.0_JPRB)
115
116      Z_SPECCOMB_PLANCK(JLAY) = P_COLH2O(JLON,JLAY)+ZREFRAT_PLANCK_A*P_COLCO2(JLON,JLAY)
117      Z_SPECPARM_PLANCK = P_COLH2O(JLON,JLAY)/Z_SPECCOMB_PLANCK(JLAY)
118      IF (Z_SPECPARM_PLANCK >= P_ONEMINUS) Z_SPECPARM_PLANCK=P_ONEMINUS
119      Z_SPECMULT_PLANCK = 8._JPRB*Z_SPECPARM_PLANCK
120      JPL= 1 + INT(Z_SPECMULT_PLANCK)
121      Z_FPL = MOD(Z_SPECMULT_PLANCK,1.0_JPRB)
122
123      IND0(JLAY) = ((K_JP(JLON,JLAY)-1)*5+(K_JT(JLON,JLAY)-1))*NSPA(12) + JS
124      IND1(JLAY) = (K_JP(JLON,JLAY)*5+(K_JT1(JLON,JLAY)-1))*NSPA(12) + JS1
125      INDS(JLAY) = K_INDSELF(JLON,JLAY)
126      INDF(JLAY) = K_INDFOR(JLON,JLAY)
127
128  IF (Z_SPECPARM < 0.125_JPRB) THEN
129            ZP = Z_FS - 1
130            ZP4 = ZP**4
131            ZFK0 = ZP4
132            ZFK1 = 1 - ZP - 2.0_JPRB*ZP4
133            ZFK2 = ZP + ZP4
134            Z_FAC000 = ZFK0*P_FAC00(JLON,JLAY)
135            Z_FAC100 = ZFK1*P_FAC00(JLON,JLAY)
136            Z_FAC200 = ZFK2*P_FAC00(JLON,JLAY)
137            Z_FAC010 = ZFK0*P_FAC10(JLON,JLAY)
138            Z_FAC110 = ZFK1*P_FAC10(JLON,JLAY)
139            Z_FAC210 = ZFK2*P_FAC10(JLON,JLAY)
140      ELSEIF (Z_SPECPARM > 0.875_JPRB) THEN
141            ZP = -Z_FS
142            ZP4 = ZP**4
143            ZFK0 = ZP4
144            ZFK1 = 1 - ZP - 2.0_JPRB*ZP4
145            ZFK2 = ZP + ZP4
146            Z_FAC000 = ZFK0*P_FAC00(JLON,JLAY)
147            Z_FAC100 = ZFK1*P_FAC00(JLON,JLAY)
148            Z_FAC200 = ZFK2*P_FAC00(JLON,JLAY)
149            Z_FAC010 = ZFK0*P_FAC10(JLON,JLAY)
150            Z_FAC110 = ZFK1*P_FAC10(JLON,JLAY)
151            Z_FAC210 = ZFK2*P_FAC10(JLON,JLAY)
152      ELSE
153            Z_FAC000 = (1._JPRB - Z_FS) * P_FAC00(JLON,JLAY)
154            Z_FAC010 = (1._JPRB - Z_FS) * P_FAC10(JLON,JLAY)
155            Z_FAC100 = Z_FS * P_FAC00(JLON,JLAY)
156            Z_FAC110 = Z_FS * P_FAC10(JLON,JLAY)
157      ENDIF
158      IF (Z_SPECPARM1 < 0.125_JPRB) THEN
159            ZP = Z_FS1 - 1
160            ZP4 = ZP**4
161            ZFK0 = ZP4
162            ZFK1 = 1 - ZP - 2.0_JPRB*ZP4
163            ZFK2 = ZP + ZP4
164            Z_FAC001 = ZFK0*P_FAC01(JLON,JLAY)
165            Z_FAC101 = ZFK1*P_FAC01(JLON,JLAY)
166            Z_FAC201 = ZFK2*P_FAC01(JLON,JLAY)
167            Z_FAC011 = ZFK0*P_FAC11(JLON,JLAY)
168            Z_FAC111 = ZFK1*P_FAC11(JLON,JLAY)
169            Z_FAC211 = ZFK2*P_FAC11(JLON,JLAY)
170      ELSEIF (Z_SPECPARM1 > 0.875_JPRB) THEN
171            ZP = -Z_FS1
172            ZP4 = ZP**4
173            ZFK0 = ZP4
174            ZFK1 = 1 - ZP - 2.0_JPRB*ZP4
175            ZFK2 = ZP + ZP4
176            Z_FAC001 = ZFK0*P_FAC01(JLON,JLAY)
177            Z_FAC101 = ZFK1*P_FAC01(JLON,JLAY)
178            Z_FAC201 = ZFK2*P_FAC01(JLON,JLAY)
179            Z_FAC011 = ZFK0*P_FAC11(JLON,JLAY)
180            Z_FAC111 = ZFK1*P_FAC11(JLON,JLAY)
181            Z_FAC211 = ZFK2*P_FAC11(JLON,JLAY)
182      ELSE
183            Z_FAC001 = (1._JPRB - Z_FS1) * P_FAC01(JLON,JLAY)
184            Z_FAC011 = (1._JPRB - Z_FS1) * P_FAC11(JLON,JLAY)
185            Z_FAC101 = Z_FS1 * P_FAC01(JLON,JLAY)
186            Z_FAC111 = Z_FS1 * P_FAC11(JLON,JLAY)
187      ENDIF
188
189
190
191!-- DS_000515
192!CDIR UNROLL=NG12
193      DO IG = 1, NG12
194!-- DS_000515
195        ZTAUSELF = P_SELFFAC(JLON,JLAY)* (SELFREF(INDS(JLAY),IG) + P_SELFFRAC(JLON,JLAY) * &
196          &       (SELFREF(INDS(JLAY)+1,IG) - SELFREF(INDS(JLAY),IG)))
197        ZTAUFOR = P_FORFAC(JLON,JLAY) * (FORREF(INDF(JLAY),IG) + P_FORFRAC(JLON,JLAY) * &
198          &       (FORREF(INDF(JLAY)+1,IG) - FORREF(INDF(JLAY),IG)))
199
200            IF (Z_SPECPARM < 0.125_JPRB) THEN
201               ZTAU_MAJOR = Z_SPECCOMB(JLAY) * &
202                 &   (Z_FAC000 * ABSA(IND0(JLAY),IG) + &
203                 &   Z_FAC100 * ABSA(IND0(JLAY)+1,IG) + &
204                 &   Z_FAC200 * ABSA(IND0(JLAY)+2,IG) + &
205                 &   Z_FAC010 * ABSA(IND0(JLAY)+9,IG) + &
206                 &   Z_FAC110 * ABSA(IND0(JLAY)+10,IG) + &
207                 &   Z_FAC210 * ABSA(IND0(JLAY)+11,IG))
208            ELSEIF (Z_SPECPARM > 0.875_JPRB) THEN
209               ZTAU_MAJOR = Z_SPECCOMB(JLAY) * &
210                 &   (Z_FAC200 * ABSA(IND0(JLAY)-1,IG) + &
211                 &   Z_FAC100 * ABSA(IND0(JLAY),IG) + &
212                 &   Z_FAC000 * ABSA(IND0(JLAY)+1,IG) + &
213                 &   Z_FAC210 * ABSA(IND0(JLAY)+8,IG) + &
214                 &   Z_FAC110 * ABSA(IND0(JLAY)+9,IG) + &
215                 &   Z_FAC010 * ABSA(IND0(JLAY)+10,IG))
216            ELSE
217               ZTAU_MAJOR = Z_SPECCOMB(JLAY) * &
218                 &   (Z_FAC000 * ABSA(IND0(JLAY),IG) + &
219                 &   Z_FAC100 * ABSA(IND0(JLAY)+1,IG) + &
220                 &   Z_FAC010 * ABSA(IND0(JLAY)+9,IG) + &
221                 &   Z_FAC110 * ABSA(IND0(JLAY)+10,IG))
222            ENDIF
223
224            IF (Z_SPECPARM1 < 0.125_JPRB) THEN
225               ZTAU_MAJOR1 = Z_SPECCOMB1(JLAY) * &
226                &    (Z_FAC001 * ABSA(IND1(JLAY),IG) + &
227                &    Z_FAC101 * ABSA(IND1(JLAY)+1,IG) + &
228                &    Z_FAC201 * ABSA(IND1(JLAY)+2,IG) + &
229                &    Z_FAC011 * ABSA(IND1(JLAY)+9,IG) + &
230                &    Z_FAC111 * ABSA(IND1(JLAY)+10,IG) + &
231                &    Z_FAC211 * ABSA(IND1(JLAY)+11,IG))
232            ELSEIF (Z_SPECPARM1 > 0.875_JPRB) THEN
233               ZTAU_MAJOR1 = Z_SPECCOMB1(JLAY) * &
234                &    (Z_FAC201 * ABSA(IND1(JLAY)-1,IG) + &
235                &    Z_FAC101 * ABSA(IND1(JLAY),IG) + &
236                &    Z_FAC001 * ABSA(IND1(JLAY)+1,IG) + &
237                &    Z_FAC211 * ABSA(IND1(JLAY)+8,IG) + &
238                &    Z_FAC111 * ABSA(IND1(JLAY)+9,IG) + &
239                &    Z_FAC011 * ABSA(IND1(JLAY)+10,IG))
240            ELSE
241               ZTAU_MAJOR1 = Z_SPECCOMB1(JLAY) * &
242                &    (Z_FAC001 * ABSA(IND1(JLAY),IG) +  &
243                &    Z_FAC101 * ABSA(IND1(JLAY)+1,IG) + &
244                &    Z_FAC011 * ABSA(IND1(JLAY)+9,IG) + &
245                &    Z_FAC111 * ABSA(IND1(JLAY)+10,IG))
246            ENDIF
247
248        P_TAU(JLON,NGS11+IG,JLAY) = ZTAU_MAJOR + ZTAU_MAJOR1 &
249         & + ZTAUSELF + ZTAUFOR &
250         & + P_TAUAERL(JLON,JLAY,12) 
251        PFRAC(JLON,NGS11+IG,JLAY) = FRACREFA(IG,JPL) + Z_FPL *&
252         & (FRACREFA(IG,JPL+1) - FRACREFA(IG,JPL))
253      ENDDO
254    ENDIF
255
256!-- JJM_000517
257    IF (JLAY > K_LAYTROP(JLON)) THEN
258!CDIR UNROLL=NG12
259      DO IG = 1, NG12
260!-- JJM_000517
261        P_TAU(JLON,NGS11+IG,JLAY) = P_TAUAERL(JLON,JLAY,12)
262        PFRAC(JLON,NGS11+IG,JLAY) = 0.0_JPRB
263      ENDDO
264    ENDIF
265  ENDDO
266ENDDO
267
268IF (LHOOK) CALL DR_HOOK('RRTM_TAUMOL12',1,ZHOOK_HANDLE)
269
270END SUBROUTINE RRTM_TAUMOL12
Note: See TracBrowser for help on using the repository browser.