/[lmdze]/trunk/Sources/phylmd/Radlwsw/swr.f
ViewVC logotype

Annotation of /trunk/Sources/phylmd/Radlwsw/swr.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 134 - (hide annotations)
Wed Apr 29 15:47:56 2015 UTC (9 years, 1 month ago) by guez
File size: 9586 byte(s)
Sources inside, compilation outside.
1 guez 81 SUBROUTINE swr(knu, palbd, pcg, pcld, pdsig, pomega, prayl, psec, ptau, &
2     pcgaz, ppizaz, pray1, pray2, prefz, prj, prk, prmue, ptauaz, ptra1, &
3     ptra2)
4     USE dimens_m
5     USE dimphy
6     USE raddim
7     USE radepsi
8     USE radopt
9     IMPLICIT NONE
10    
11     ! ------------------------------------------------------------------
12     ! PURPOSE.
13     ! --------
14     ! COMPUTES THE REFLECTIVITY AND TRANSMISSIVITY IN CASE OF
15     ! CONTINUUM SCATTERING
16    
17     ! METHOD.
18     ! -------
19    
20     ! 1. COMPUTES CONTINUUM FLUXES CORRESPONDING TO AEROSOL
21     ! OR/AND RAYLEIGH SCATTERING (NO MOLECULAR GAS ABSORPTION)
22    
23     ! REFERENCE.
24     ! ----------
25    
26     ! SEE RADIATION'S PART OF THE ECMWF RESEARCH DEPARTMENT
27     ! DOCUMENTATION, AND FOUQUART AND BONNEL (1980)
28    
29     ! AUTHOR.
30     ! -------
31     ! JEAN-JACQUES MORCRETTE *ECMWF*
32    
33     ! MODIFICATIONS.
34     ! --------------
35     ! ORIGINAL : 89-07-14
36     ! ------------------------------------------------------------------
37     ! * ARGUMENTS:
38    
39     INTEGER knu
40     DOUBLE PRECISION palbd(kdlon, 2)
41     DOUBLE PRECISION pcg(kdlon, 2, kflev)
42     DOUBLE PRECISION pcld(kdlon, kflev)
43     DOUBLE PRECISION pdsig(kdlon, kflev)
44     DOUBLE PRECISION pomega(kdlon, 2, kflev)
45     DOUBLE PRECISION prayl(kdlon)
46     DOUBLE PRECISION psec(kdlon)
47     DOUBLE PRECISION ptau(kdlon, 2, kflev)
48    
49     DOUBLE PRECISION pray1(kdlon, kflev+1)
50     DOUBLE PRECISION pray2(kdlon, kflev+1)
51     DOUBLE PRECISION prefz(kdlon, 2, kflev+1)
52     DOUBLE PRECISION prj(kdlon, 6, kflev+1)
53     DOUBLE PRECISION prk(kdlon, 6, kflev+1)
54     DOUBLE PRECISION prmue(kdlon, kflev+1)
55     DOUBLE PRECISION pcgaz(kdlon, kflev)
56     DOUBLE PRECISION ppizaz(kdlon, kflev)
57     DOUBLE PRECISION ptauaz(kdlon, kflev)
58     DOUBLE PRECISION ptra1(kdlon, kflev+1)
59     DOUBLE PRECISION ptra2(kdlon, kflev+1)
60    
61     ! * LOCAL VARIABLES:
62    
63     DOUBLE PRECISION zc1i(kdlon, kflev+1)
64     DOUBLE PRECISION zcleq(kdlon, kflev)
65     DOUBLE PRECISION zclear(kdlon)
66     DOUBLE PRECISION zcloud(kdlon)
67     DOUBLE PRECISION zgg(kdlon)
68     DOUBLE PRECISION zref(kdlon)
69     DOUBLE PRECISION zre1(kdlon)
70     DOUBLE PRECISION zre2(kdlon)
71     DOUBLE PRECISION zrmuz(kdlon)
72     DOUBLE PRECISION zrneb(kdlon)
73     DOUBLE PRECISION zr21(kdlon)
74     DOUBLE PRECISION zr22(kdlon)
75     DOUBLE PRECISION zr23(kdlon)
76     DOUBLE PRECISION zss1(kdlon)
77     DOUBLE PRECISION zto1(kdlon)
78     DOUBLE PRECISION ztr(kdlon, 2, kflev+1)
79     DOUBLE PRECISION ztr1(kdlon)
80     DOUBLE PRECISION ztr2(kdlon)
81     DOUBLE PRECISION zw(kdlon)
82    
83     INTEGER jk, jl, ja, jkl, jklp1, jkm1, jaj
84     DOUBLE PRECISION zfacoa, zfacoc, zcorae, zcorcd
85     DOUBLE PRECISION zmue, zgap, zww, zto, zden, zden1
86     DOUBLE PRECISION zmu1, zre11, zbmu0, zbmu1
87    
88     ! ------------------------------------------------------------------
89    
90     ! * 1. INITIALIZATION
91     ! --------------
92    
93    
94     DO jk = 1, kflev + 1
95     DO ja = 1, 6
96     DO jl = 1, kdlon
97     prj(jl, ja, jk) = 0.
98     prk(jl, ja, jk) = 0.
99     END DO
100     END DO
101     END DO
102    
103    
104     ! ------------------------------------------------------------------
105    
106     ! * 2. TOTAL EFFECTIVE CLOUDINESS ABOVE A GIVEN LEVEL
107     ! ----------------------------------------------
108    
109    
110     DO jl = 1, kdlon
111     zr23(jl) = 0.
112     zc1i(jl, kflev+1) = 0.
113     zclear(jl) = 1.
114     zcloud(jl) = 0.
115     END DO
116    
117     jk = 1
118     jkl = kflev + 1 - jk
119     jklp1 = jkl + 1
120     DO jl = 1, kdlon
121     zfacoa = 1. - ppizaz(jl, jkl)*pcgaz(jl, jkl)*pcgaz(jl, jkl)
122     zfacoc = 1. - pomega(jl, knu, jkl)*pcg(jl, knu, jkl)*pcg(jl, knu, jkl)
123     zcorae = zfacoa*ptauaz(jl, jkl)*psec(jl)
124     zcorcd = zfacoc*ptau(jl, knu, jkl)*psec(jl)
125     zr21(jl) = exp(-zcorae)
126     zr22(jl) = exp(-zcorcd)
127     zss1(jl) = pcld(jl, jkl)*(1.0-zr21(jl)*zr22(jl)) + &
128     (1.0-pcld(jl,jkl))*(1.0-zr21(jl))
129     zcleq(jl, jkl) = zss1(jl)
130    
131     IF (novlp==1) THEN
132     ! * maximum-random
133     zclear(jl) = zclear(jl)*(1.0-max(zss1(jl),zcloud(jl)))/ &
134     (1.0-min(zcloud(jl),1.-zepsec))
135     zc1i(jl, jkl) = 1.0 - zclear(jl)
136     zcloud(jl) = zss1(jl)
137     ELSE IF (novlp==2) THEN
138     ! * maximum
139     zcloud(jl) = max(zss1(jl), zcloud(jl))
140     zc1i(jl, jkl) = zcloud(jl)
141     ELSE IF (novlp==3) THEN
142     ! * random
143     zclear(jl) = zclear(jl)*(1.0-zss1(jl))
144     zcloud(jl) = 1.0 - zclear(jl)
145     zc1i(jl, jkl) = zcloud(jl)
146     END IF
147     END DO
148    
149     DO jk = 2, kflev
150     jkl = kflev + 1 - jk
151     jklp1 = jkl + 1
152     DO jl = 1, kdlon
153     zfacoa = 1. - ppizaz(jl, jkl)*pcgaz(jl, jkl)*pcgaz(jl, jkl)
154     zfacoc = 1. - pomega(jl, knu, jkl)*pcg(jl, knu, jkl)*pcg(jl, knu, jkl)
155     zcorae = zfacoa*ptauaz(jl, jkl)*psec(jl)
156     zcorcd = zfacoc*ptau(jl, knu, jkl)*psec(jl)
157     zr21(jl) = exp(-zcorae)
158     zr22(jl) = exp(-zcorcd)
159     zss1(jl) = pcld(jl, jkl)*(1.0-zr21(jl)*zr22(jl)) + &
160     (1.0-pcld(jl,jkl))*(1.0-zr21(jl))
161     zcleq(jl, jkl) = zss1(jl)
162    
163     IF (novlp==1) THEN
164     ! * maximum-random
165     zclear(jl) = zclear(jl)*(1.0-max(zss1(jl),zcloud(jl)))/ &
166     (1.0-min(zcloud(jl),1.-zepsec))
167     zc1i(jl, jkl) = 1.0 - zclear(jl)
168     zcloud(jl) = zss1(jl)
169     ELSE IF (novlp==2) THEN
170     ! * maximum
171     zcloud(jl) = max(zss1(jl), zcloud(jl))
172     zc1i(jl, jkl) = zcloud(jl)
173     ELSE IF (novlp==3) THEN
174     ! * random
175     zclear(jl) = zclear(jl)*(1.0-zss1(jl))
176     zcloud(jl) = 1.0 - zclear(jl)
177     zc1i(jl, jkl) = zcloud(jl)
178 guez 24 END IF
179 guez 81 END DO
180     END DO
181    
182     ! ------------------------------------------------------------------
183    
184     ! * 3. REFLECTIVITY/TRANSMISSIVITY FOR PURE SCATTERING
185     ! -----------------------------------------------
186    
187    
188     DO jl = 1, kdlon
189     pray1(jl, kflev+1) = 0.
190     pray2(jl, kflev+1) = 0.
191     prefz(jl, 2, 1) = palbd(jl, knu)
192     prefz(jl, 1, 1) = palbd(jl, knu)
193     ptra1(jl, kflev+1) = 1.
194     ptra2(jl, kflev+1) = 1.
195     END DO
196    
197     DO jk = 2, kflev + 1
198     jkm1 = jk - 1
199     DO jl = 1, kdlon
200     zrneb(jl) = pcld(jl, jkm1)
201     zre1(jl) = 0.
202     ztr1(jl) = 0.
203     zre2(jl) = 0.
204     ztr2(jl) = 0.
205    
206    
207     ! ------------------------------------------------------------------
208    
209     ! * 3.1 EQUIVALENT ZENITH ANGLE
210     ! -----------------------
211    
212    
213     zmue = (1.-zc1i(jl,jk))*psec(jl) + zc1i(jl, jk)*1.66
214     prmue(jl, jk) = 1./zmue
215    
216    
217     ! ------------------------------------------------------------------
218    
219     ! * 3.2 REFLECT./TRANSMISSIVITY DUE TO RAYLEIGH AND AEROSOLS
220     ! ----------------------------------------------------
221    
222    
223     zgap = pcgaz(jl, jkm1)
224     zbmu0 = 0.5 - 0.75*zgap/zmue
225     zww = ppizaz(jl, jkm1)
226     zto = ptauaz(jl, jkm1)
227     zden = 1. + (1.-zww+zbmu0*zww)*zto*zmue + (1-zww)*(1.-zww+2.*zbmu0*zww) &
228     *zto*zto*zmue*zmue
229     pray1(jl, jkm1) = zbmu0*zww*zto*zmue/zden
230     ptra1(jl, jkm1) = 1./zden
231     ! PRINT *,' LOOP 342 ** 3 ** JL=',JL,PRAY1(JL,JKM1),PTRA1(JL,JKM1)
232    
233     zmu1 = 0.5
234     zbmu1 = 0.5 - 0.75*zgap*zmu1
235     zden1 = 1. + (1.-zww+zbmu1*zww)*zto/zmu1 + (1-zww)*(1.-zww+2.*zbmu1*zww &
236     )*zto*zto/zmu1/zmu1
237     pray2(jl, jkm1) = zbmu1*zww*zto/zmu1/zden1
238     ptra2(jl, jkm1) = 1./zden1
239    
240    
241     ! ------------------------------------------------------------------
242    
243     ! * 3.3 EFFECT OF CLOUD LAYER
244     ! ---------------------
245    
246    
247     zw(jl) = pomega(jl, knu, jkm1)
248     zto1(jl) = ptau(jl, knu, jkm1)/zw(jl) + ptauaz(jl, jkm1)/ppizaz(jl, &
249     jkm1)
250     zr21(jl) = ptau(jl, knu, jkm1) + ptauaz(jl, jkm1)
251     zr22(jl) = ptau(jl, knu, jkm1)/zr21(jl)
252     zgg(jl) = zr22(jl)*pcg(jl, knu, jkm1) + (1.-zr22(jl))*pcgaz(jl, jkm1)
253     ! Modif PhD - JJM 19/03/96 pour erreurs arrondis
254     ! machine
255     ! PHD PROTECTION ZW(JL) = ZR21(JL) / ZTO1(JL)
256     IF (zw(jl)==1. .AND. ppizaz(jl,jkm1)==1.) THEN
257     zw(jl) = 1.
258 guez 24 ELSE
259 guez 81 zw(jl) = zr21(jl)/zto1(jl)
260 guez 24 END IF
261 guez 81 zref(jl) = prefz(jl, 1, jkm1)
262     zrmuz(jl) = prmue(jl, jk)
263     END DO
264    
265     CALL swde(zgg, zref, zrmuz, zto1, zw, zre1, zre2, ztr1, ztr2)
266    
267     DO jl = 1, kdlon
268    
269     prefz(jl, 1, jk) = (1.-zrneb(jl))*(pray1(jl,jkm1)+prefz(jl,1,jkm1)* &
270     ptra1(jl,jkm1)*ptra2(jl,jkm1)/(1.-pray2(jl,jkm1)*prefz(jl,1, &
271     jkm1))) + zrneb(jl)*zre2(jl)
272    
273     ztr(jl, 1, jkm1) = zrneb(jl)*ztr2(jl) + (ptra1(jl,jkm1)/(1.-pray2(jl, &
274     jkm1)*prefz(jl,1,jkm1)))*(1.-zrneb(jl))
275    
276     prefz(jl, 2, jk) = (1.-zrneb(jl))*(pray1(jl,jkm1)+prefz(jl,2,jkm1)* &
277     ptra1(jl,jkm1)*ptra2(jl,jkm1)) + zrneb(jl)*zre1(jl)
278    
279     ztr(jl, 2, jkm1) = zrneb(jl)*ztr1(jl) + ptra1(jl, jkm1)*(1.-zrneb(jl))
280    
281     END DO
282     END DO
283     DO jl = 1, kdlon
284     zmue = (1.-zc1i(jl,1))*psec(jl) + zc1i(jl, 1)*1.66
285     prmue(jl, 1) = 1./zmue
286     END DO
287    
288    
289     ! ------------------------------------------------------------------
290    
291     ! * 3.5 REFLECT./TRANSMISSIVITY BETWEEN SURFACE AND LEVEL
292     ! -------------------------------------------------
293    
294    
295     IF (knu==1) THEN
296     jaj = 2
297     DO jl = 1, kdlon
298     prj(jl, jaj, kflev+1) = 1.
299     prk(jl, jaj, kflev+1) = prefz(jl, 1, kflev+1)
300     END DO
301    
302     DO jk = 1, kflev
303     jkl = kflev + 1 - jk
304     jklp1 = jkl + 1
305     DO jl = 1, kdlon
306     zre11 = prj(jl, jaj, jklp1)*ztr(jl, 1, jkl)
307     prj(jl, jaj, jkl) = zre11
308     prk(jl, jaj, jkl) = zre11*prefz(jl, 1, jkl)
309     END DO
310     END DO
311    
312     ELSE
313    
314     DO jaj = 1, 2
315     DO jl = 1, kdlon
316     prj(jl, jaj, kflev+1) = 1.
317     prk(jl, jaj, kflev+1) = prefz(jl, jaj, kflev+1)
318     END DO
319    
320     DO jk = 1, kflev
321     jkl = kflev + 1 - jk
322     jklp1 = jkl + 1
323     DO jl = 1, kdlon
324     zre11 = prj(jl, jaj, jklp1)*ztr(jl, jaj, jkl)
325     prj(jl, jaj, jkl) = zre11
326     prk(jl, jaj, jkl) = zre11*prefz(jl, jaj, jkl)
327     END DO
328     END DO
329     END DO
330    
331     END IF
332    
333     ! ------------------------------------------------------------------
334    
335     RETURN
336     END SUBROUTINE swr

  ViewVC Help
Powered by ViewVC 1.1.21