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

Annotation of /trunk/Sources/phylmd/Radlwsw/sw2s.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: 14078 byte(s)
Sources inside, compilation outside.
1 guez 81 SUBROUTINE sw2s(knu, paer, flag_aer, tauae, pizae, cgae, paki, palbd, palbp, &
2     pcg, pcld, pclear, pcldsw, pdsig, pomega, poz, prmu, psec, ptau, pud, &
3     pwv, pqs, pfdown, pfup)
4     USE dimens_m
5     USE dimphy
6     USE raddim
7     USE radepsi
8     IMPLICIT NONE
9 guez 24
10 guez 81 ! ------------------------------------------------------------------
11     ! PURPOSE.
12     ! --------
13    
14     ! THIS ROUTINE COMPUTES THE SHORTWAVE RADIATION FLUXES IN THE
15     ! SECOND SPECTRAL INTERVAL FOLLOWING FOUQUART AND BONNEL (1980).
16    
17     ! METHOD.
18     ! -------
19    
20     ! 1. COMPUTES REFLECTIVITY/TRANSMISSIVITY CORRESPONDING TO
21     ! CONTINUUM SCATTERING
22     ! 2. COMPUTES REFLECTIVITY/TRANSMISSIVITY CORRESPONDING FOR
23     ! A GREY MOLECULAR ABSORPTION
24     ! 3. LAPLACE TRANSFORM ON THE PREVIOUS TO GET EFFECTIVE AMOUNTS
25     ! OF ABSORBERS
26     ! 4. APPLY H2O AND U.M.G. TRANSMISSION FUNCTIONS
27     ! 5. MULTIPLY BY OZONE TRANSMISSION FUNCTION
28    
29     ! REFERENCE.
30     ! ----------
31    
32     ! SEE RADIATION'S PART OF THE ECMWF RESEARCH DEPARTMENT
33     ! DOCUMENTATION, AND FOUQUART AND BONNEL (1980)
34    
35     ! AUTHOR.
36     ! -------
37     ! JEAN-JACQUES MORCRETTE *ECMWF*
38    
39     ! MODIFICATIONS.
40     ! --------------
41     ! ORIGINAL : 89-07-14
42     ! 94-11-15 J.-J. MORCRETTE DIRECT/DIFFUSE ALBEDO
43     ! ------------------------------------------------------------------
44     ! * ARGUMENTS:
45    
46     INTEGER knu
47     ! -OB
48     DOUBLE PRECISION flag_aer
49     DOUBLE PRECISION tauae(kdlon, kflev, 2)
50     DOUBLE PRECISION pizae(kdlon, kflev, 2)
51     DOUBLE PRECISION cgae(kdlon, kflev, 2)
52     DOUBLE PRECISION paer(kdlon, kflev, 5)
53     DOUBLE PRECISION paki(kdlon, 2)
54     DOUBLE PRECISION palbd(kdlon, 2)
55     DOUBLE PRECISION palbp(kdlon, 2)
56     DOUBLE PRECISION pcg(kdlon, 2, kflev)
57     DOUBLE PRECISION pcld(kdlon, kflev)
58     DOUBLE PRECISION pcldsw(kdlon, kflev)
59     DOUBLE PRECISION pclear(kdlon)
60     DOUBLE PRECISION pdsig(kdlon, kflev)
61     DOUBLE PRECISION pomega(kdlon, 2, kflev)
62     DOUBLE PRECISION poz(kdlon, kflev)
63     DOUBLE PRECISION pqs(kdlon, kflev)
64     DOUBLE PRECISION prmu(kdlon)
65     DOUBLE PRECISION psec(kdlon)
66     DOUBLE PRECISION ptau(kdlon, 2, kflev)
67     DOUBLE PRECISION pud(kdlon, 5, kflev+1)
68     DOUBLE PRECISION pwv(kdlon, kflev)
69    
70     DOUBLE PRECISION pfdown(kdlon, kflev+1)
71     DOUBLE PRECISION pfup(kdlon, kflev+1)
72    
73     ! * LOCAL VARIABLES:
74    
75     INTEGER iind2(2), iind3(3)
76     DOUBLE PRECISION zcgaz(kdlon, kflev)
77     DOUBLE PRECISION zfd(kdlon, kflev+1)
78     DOUBLE PRECISION zfu(kdlon, kflev+1)
79     DOUBLE PRECISION zg(kdlon)
80     DOUBLE PRECISION zgg(kdlon)
81     DOUBLE PRECISION zpizaz(kdlon, kflev)
82     DOUBLE PRECISION zrayl(kdlon)
83     DOUBLE PRECISION zray1(kdlon, kflev+1)
84     DOUBLE PRECISION zray2(kdlon, kflev+1)
85     DOUBLE PRECISION zref(kdlon)
86     DOUBLE PRECISION zrefz(kdlon, 2, kflev+1)
87     DOUBLE PRECISION zre1(kdlon)
88     DOUBLE PRECISION zre2(kdlon)
89     DOUBLE PRECISION zrj(kdlon, 6, kflev+1)
90     DOUBLE PRECISION zrj0(kdlon, 6, kflev+1)
91     DOUBLE PRECISION zrk(kdlon, 6, kflev+1)
92     DOUBLE PRECISION zrk0(kdlon, 6, kflev+1)
93     DOUBLE PRECISION zrl(kdlon, 8)
94     DOUBLE PRECISION zrmue(kdlon, kflev+1)
95     DOUBLE PRECISION zrmu0(kdlon, kflev+1)
96     DOUBLE PRECISION zrmuz(kdlon)
97     DOUBLE PRECISION zrneb(kdlon)
98     DOUBLE PRECISION zruef(kdlon, 8)
99     DOUBLE PRECISION zr1(kdlon)
100     DOUBLE PRECISION zr2(kdlon, 2)
101     DOUBLE PRECISION zr3(kdlon, 3)
102     DOUBLE PRECISION zr4(kdlon)
103     DOUBLE PRECISION zr21(kdlon)
104     DOUBLE PRECISION zr22(kdlon)
105     DOUBLE PRECISION zs(kdlon)
106     DOUBLE PRECISION ztauaz(kdlon, kflev)
107     DOUBLE PRECISION zto1(kdlon)
108     DOUBLE PRECISION ztr(kdlon, 2, kflev+1)
109     DOUBLE PRECISION ztra1(kdlon, kflev+1)
110     DOUBLE PRECISION ztra2(kdlon, kflev+1)
111     DOUBLE PRECISION ztr1(kdlon)
112     DOUBLE PRECISION ztr2(kdlon)
113     DOUBLE PRECISION zw(kdlon)
114     DOUBLE PRECISION zw1(kdlon)
115     DOUBLE PRECISION zw2(kdlon, 2)
116     DOUBLE PRECISION zw3(kdlon, 3)
117     DOUBLE PRECISION zw4(kdlon)
118     DOUBLE PRECISION zw5(kdlon)
119    
120     INTEGER jl, jk, k, jaj, ikm1, ikl, jn, jabs, jkm1
121     INTEGER jref, jkl, jklp1, jajp, jkki, jkkp4, jn2j, iabs
122     DOUBLE PRECISION zrmum1, zwh2o, zcneb, zaa, zbb, zrki, zre11
123    
124     ! * Prescribed Data:
125    
126     DOUBLE PRECISION rsun(2)
127     SAVE rsun
128     DOUBLE PRECISION rray(2, 6)
129     SAVE rray
130     DATA rsun(1)/0.441676/
131     DATA rsun(2)/0.558324/
132     DATA (rray(1,k), k=1, 6)/.428937E-01, .890743E+00, -.288555E+01, &
133     .522744E+01, -.469173E+01, .161645E+01/
134     DATA (rray(2,k), k=1, 6)/.697200E-02, .173297E-01, -.850903E-01, &
135     .248261E+00, -.302031E+00, .129662E+00/
136    
137     ! ------------------------------------------------------------------
138    
139     ! * 1. SECOND SPECTRAL INTERVAL (0.68-4.00 MICRON)
140     ! -------------------------------------------
141    
142    
143    
144     ! * 1.1 OPTICAL THICKNESS FOR RAYLEIGH SCATTERING
145     ! -----------------------------------------
146    
147    
148     DO jl = 1, kdlon
149     zrmum1 = 1. - prmu(jl)
150     zrayl(jl) = rray(knu, 1) + zrmum1*(rray(knu,2)+zrmum1*(rray(knu, &
151     3)+zrmum1*(rray(knu,4)+zrmum1*(rray(knu,5)+zrmum1*rray(knu,6)))))
152     END DO
153    
154    
155     ! ------------------------------------------------------------------
156    
157     ! * 2. CONTINUUM SCATTERING CALCULATIONS
158     ! ---------------------------------
159    
160    
161     ! * 2.1 CLEAR-SKY FRACTION OF THE COLUMN
162     ! --------------------------------
163    
164    
165     CALL swclr(knu, paer, flag_aer, tauae, pizae, cgae, palbp, pdsig, zrayl, &
166     psec, zcgaz, zpizaz, zray1, zray2, zrefz, zrj0, zrk0, zrmu0, ztauaz, &
167     ztra1, ztra2)
168    
169    
170     ! * 2.2 CLOUDY FRACTION OF THE COLUMN
171     ! -----------------------------
172    
173    
174     CALL swr(knu, palbd, pcg, pcld, pdsig, pomega, zrayl, psec, ptau, zcgaz, &
175     zpizaz, zray1, zray2, zrefz, zrj, zrk, zrmue, ztauaz, ztra1, ztra2)
176    
177    
178     ! ------------------------------------------------------------------
179    
180     ! * 3. SCATTERING CALCULATIONS WITH GREY MOLECULAR ABSORPTION
181     ! ------------------------------------------------------
182    
183    
184     jn = 2
185    
186     DO jabs = 1, 2
187    
188    
189     ! * 3.1 SURFACE CONDITIONS
190     ! ------------------
191    
192    
193     DO jl = 1, kdlon
194     zrefz(jl, 2, 1) = palbd(jl, knu)
195     zrefz(jl, 1, 1) = palbd(jl, knu)
196     END DO
197    
198    
199     ! * 3.2 INTRODUCING CLOUD EFFECTS
200     ! -------------------------
201    
202    
203     DO jk = 2, kflev + 1
204     jkm1 = jk - 1
205     ikl = kflev + 1 - jkm1
206     DO jl = 1, kdlon
207     zrneb(jl) = pcld(jl, jkm1)
208     IF (jabs==1 .AND. zrneb(jl)>2.*zeelog) THEN
209     zwh2o = max(pwv(jl,jkm1), zeelog)
210     zcneb = max(zeelog, min(zrneb(jl),1.-zeelog))
211     zbb = pud(jl, jabs, jkm1)*pqs(jl, jkm1)/zwh2o
212     zaa = max((pud(jl,jabs,jkm1)-zcneb*zbb)/(1.-zcneb), zeelog)
213     ELSE
214     zaa = pud(jl, jabs, jkm1)
215     zbb = zaa
216     END IF
217     zrki = paki(jl, jabs)
218     zs(jl) = exp(-zrki*zaa*1.66)
219     zg(jl) = exp(-zrki*zaa/zrmue(jl,jk))
220     ztr1(jl) = 0.
221     zre1(jl) = 0.
222     ztr2(jl) = 0.
223     zre2(jl) = 0.
224    
225     zw(jl) = pomega(jl, knu, jkm1)
226     zto1(jl) = ptau(jl, knu, jkm1)/zw(jl) + ztauaz(jl, jkm1)/zpizaz(jl, &
227     jkm1) + zbb*zrki
228    
229     zr21(jl) = ptau(jl, knu, jkm1) + ztauaz(jl, jkm1)
230     zr22(jl) = ptau(jl, knu, jkm1)/zr21(jl)
231     zgg(jl) = zr22(jl)*pcg(jl, knu, jkm1) + (1.-zr22(jl))*zcgaz(jl, jkm1)
232     zw(jl) = zr21(jl)/zto1(jl)
233     zref(jl) = zrefz(jl, 1, jkm1)
234     zrmuz(jl) = zrmue(jl, jk)
235     END DO
236    
237     CALL swde(zgg, zref, zrmuz, zto1, zw, zre1, zre2, ztr1, ztr2)
238    
239     DO jl = 1, kdlon
240    
241     zrefz(jl, 2, jk) = (1.-zrneb(jl))*(zray1(jl,jkm1)+zrefz(jl,2,jkm1)* &
242     ztra1(jl,jkm1)*ztra2(jl,jkm1))*zg(jl)*zs(jl) + zrneb(jl)*zre1(jl)
243    
244     ztr(jl, 2, jkm1) = zrneb(jl)*ztr1(jl) + (ztra1(jl,jkm1))*zg(jl)*(1.- &
245     zrneb(jl))
246    
247     zrefz(jl, 1, jk) = (1.-zrneb(jl))*(zray1(jl,jkm1)+zrefz(jl,1,jkm1)* &
248     ztra1(jl,jkm1)*ztra2(jl,jkm1)/(1.-zray2(jl,jkm1)*zrefz(jl,1, &
249     jkm1)))*zg(jl)*zs(jl) + zrneb(jl)*zre2(jl)
250    
251     ztr(jl, 1, jkm1) = zrneb(jl)*ztr2(jl) + (ztra1(jl,jkm1)/(1.-zray2(jl, &
252     jkm1)*zrefz(jl,1,jkm1)))*zg(jl)*(1.-zrneb(jl))
253    
254     END DO
255     END DO
256    
257     ! * 3.3 REFLECT./TRANSMISSIVITY BETWEEN SURFACE AND LEVEL
258     ! -------------------------------------------------
259    
260    
261     DO jref = 1, 2
262    
263     jn = jn + 1
264    
265     DO jl = 1, kdlon
266     zrj(jl, jn, kflev+1) = 1.
267     zrk(jl, jn, kflev+1) = zrefz(jl, jref, kflev+1)
268     END DO
269    
270     DO jk = 1, kflev
271     jkl = kflev + 1 - jk
272     jklp1 = jkl + 1
273     DO jl = 1, kdlon
274     zre11 = zrj(jl, jn, jklp1)*ztr(jl, jref, jkl)
275     zrj(jl, jn, jkl) = zre11
276     zrk(jl, jn, jkl) = zre11*zrefz(jl, jref, jkl)
277     END DO
278     END DO
279     END DO
280     END DO
281    
282    
283     ! ------------------------------------------------------------------
284    
285     ! * 4. INVERT GREY AND CONTINUUM FLUXES
286     ! --------------------------------
287    
288    
289    
290     ! * 4.1 UPWARD (ZRK) AND DOWNWARD (ZRJ) PSEUDO-FLUXES
291     ! ---------------------------------------------
292    
293    
294     DO jk = 1, kflev + 1
295     DO jaj = 1, 5, 2
296     jajp = jaj + 1
297     DO jl = 1, kdlon
298     zrj(jl, jaj, jk) = zrj(jl, jaj, jk) - zrj(jl, jajp, jk)
299     zrk(jl, jaj, jk) = zrk(jl, jaj, jk) - zrk(jl, jajp, jk)
300     zrj(jl, jaj, jk) = max(zrj(jl,jaj,jk), zeelog)
301     zrk(jl, jaj, jk) = max(zrk(jl,jaj,jk), zeelog)
302     END DO
303     END DO
304     END DO
305    
306     DO jk = 1, kflev + 1
307     DO jaj = 2, 6, 2
308     DO jl = 1, kdlon
309     zrj(jl, jaj, jk) = max(zrj(jl,jaj,jk), zeelog)
310     zrk(jl, jaj, jk) = max(zrk(jl,jaj,jk), zeelog)
311     END DO
312     END DO
313     END DO
314    
315     ! * 4.2 EFFECTIVE ABSORBER AMOUNTS BY INVERSE LAPLACE
316     ! ---------------------------------------------
317    
318    
319     DO jk = 1, kflev + 1
320     jkki = 1
321     DO jaj = 1, 2
322     iind2(1) = jaj
323     iind2(2) = jaj
324     DO jn = 1, 2
325     jn2j = jn + 2*jaj
326     jkkp4 = jkki + 4
327    
328     ! * 4.2.1 EFFECTIVE ABSORBER AMOUNTS
329     ! --------------------------
330    
331    
332     DO jl = 1, kdlon
333     zw2(jl, 1) = log(zrj(jl,jn,jk)/zrj(jl,jn2j,jk))/paki(jl, jaj)
334     zw2(jl, 2) = log(zrk(jl,jn,jk)/zrk(jl,jn2j,jk))/paki(jl, jaj)
335     END DO
336    
337     ! * 4.2.2 TRANSMISSION FUNCTION
338     ! ---------------------
339    
340    
341     CALL swtt1(knu, 2, iind2, zw2, zr2)
342    
343     DO jl = 1, kdlon
344     zrl(jl, jkki) = zr2(jl, 1)
345     zruef(jl, jkki) = zw2(jl, 1)
346     zrl(jl, jkkp4) = zr2(jl, 2)
347     zruef(jl, jkkp4) = zw2(jl, 2)
348     END DO
349    
350     jkki = jkki + 1
351     END DO
352     END DO
353    
354     ! * 4.3 UPWARD AND DOWNWARD FLUXES WITH H2O AND UMG ABSORPTION
355     ! ------------------------------------------------------
356    
357    
358     DO jl = 1, kdlon
359     pfdown(jl, jk) = zrj(jl, 1, jk)*zrl(jl, 1)*zrl(jl, 3) + &
360     zrj(jl, 2, jk)*zrl(jl, 2)*zrl(jl, 4)
361     pfup(jl, jk) = zrk(jl, 1, jk)*zrl(jl, 5)*zrl(jl, 7) + &
362     zrk(jl, 2, jk)*zrl(jl, 6)*zrl(jl, 8)
363     END DO
364     END DO
365    
366    
367     ! ------------------------------------------------------------------
368    
369     ! * 5. MOLECULAR ABSORPTION ON CLEAR-SKY FLUXES
370     ! ----------------------------------------
371    
372    
373    
374     ! * 5.1 DOWNWARD FLUXES
375     ! ---------------
376    
377    
378     jaj = 2
379     iind3(1) = 1
380     iind3(2) = 2
381     iind3(3) = 3
382    
383     DO jl = 1, kdlon
384     zw3(jl, 1) = 0.
385     zw3(jl, 2) = 0.
386     zw3(jl, 3) = 0.
387     zw4(jl) = 0.
388     zw5(jl) = 0.
389     zr4(jl) = 1.
390     zfd(jl, kflev+1) = zrj0(jl, jaj, kflev+1)
391     END DO
392     DO jk = 1, kflev
393     ikl = kflev + 1 - jk
394     DO jl = 1, kdlon
395     zw3(jl, 1) = zw3(jl, 1) + pud(jl, 1, ikl)/zrmu0(jl, ikl)
396     zw3(jl, 2) = zw3(jl, 2) + pud(jl, 2, ikl)/zrmu0(jl, ikl)
397     zw3(jl, 3) = zw3(jl, 3) + poz(jl, ikl)/zrmu0(jl, ikl)
398     zw4(jl) = zw4(jl) + pud(jl, 4, ikl)/zrmu0(jl, ikl)
399     zw5(jl) = zw5(jl) + pud(jl, 5, ikl)/zrmu0(jl, ikl)
400     END DO
401    
402     CALL swtt1(knu, 3, iind3, zw3, zr3)
403    
404     DO jl = 1, kdlon
405     ! ZR4(JL) = EXP(-RSWCE*ZW4(JL)-RSWCP*ZW5(JL))
406     zfd(jl, ikl) = zr3(jl, 1)*zr3(jl, 2)*zr3(jl, 3)*zr4(jl)* &
407     zrj0(jl, jaj, ikl)
408     END DO
409     END DO
410    
411    
412     ! * 5.2 UPWARD FLUXES
413     ! -------------
414    
415    
416     DO jl = 1, kdlon
417     zfu(jl, 1) = zfd(jl, 1)*palbp(jl, knu)
418     END DO
419    
420     DO jk = 2, kflev + 1
421     ikm1 = jk - 1
422     DO jl = 1, kdlon
423     zw3(jl, 1) = zw3(jl, 1) + pud(jl, 1, ikm1)*1.66
424     zw3(jl, 2) = zw3(jl, 2) + pud(jl, 2, ikm1)*1.66
425     zw3(jl, 3) = zw3(jl, 3) + poz(jl, ikm1)*1.66
426     zw4(jl) = zw4(jl) + pud(jl, 4, ikm1)*1.66
427     zw5(jl) = zw5(jl) + pud(jl, 5, ikm1)*1.66
428     END DO
429    
430     CALL swtt1(knu, 3, iind3, zw3, zr3)
431    
432     DO jl = 1, kdlon
433     ! ZR4(JL) = EXP(-RSWCE*ZW4(JL)-RSWCP*ZW5(JL))
434     zfu(jl, jk) = zr3(jl, 1)*zr3(jl, 2)*zr3(jl, 3)*zr4(jl)* &
435     zrk0(jl, jaj, jk)
436     END DO
437     END DO
438    
439    
440     ! ------------------------------------------------------------------
441    
442     ! * 6. INTRODUCTION OF OZONE AND H2O CONTINUUM ABSORPTION
443     ! --------------------------------------------------
444    
445     iabs = 3
446    
447     ! * 6.1 DOWNWARD FLUXES
448     ! ---------------
449    
450     DO jl = 1, kdlon
451     zw1(jl) = 0.
452     zw4(jl) = 0.
453     zw5(jl) = 0.
454     zr1(jl) = 0.
455     pfdown(jl, kflev+1) = ((1.-pclear(jl))*pfdown(jl,kflev+1)+pclear(jl)*zfd( &
456     jl,kflev+1))*rsun(knu)
457     END DO
458    
459     DO jk = 1, kflev
460     ikl = kflev + 1 - jk
461     DO jl = 1, kdlon
462     zw1(jl) = zw1(jl) + poz(jl, ikl)/zrmue(jl, ikl)
463     zw4(jl) = zw4(jl) + pud(jl, 4, ikl)/zrmue(jl, ikl)
464     zw5(jl) = zw5(jl) + pud(jl, 5, ikl)/zrmue(jl, ikl)
465     ! ZR4(JL) = EXP(-RSWCE*ZW4(JL)-RSWCP*ZW5(JL))
466     END DO
467    
468     CALL swtt(knu, iabs, zw1, zr1)
469    
470     DO jl = 1, kdlon
471     pfdown(jl, ikl) = ((1.-pclear(jl))*zr1(jl)*zr4(jl)*pfdown(jl,ikl)+ &
472     pclear(jl)*zfd(jl,ikl))*rsun(knu)
473     END DO
474     END DO
475    
476    
477     ! * 6.2 UPWARD FLUXES
478     ! -------------
479    
480     DO jl = 1, kdlon
481     pfup(jl, 1) = ((1.-pclear(jl))*zr1(jl)*zr4(jl)*pfup(jl,1)+pclear(jl)*zfu( &
482     jl,1))*rsun(knu)
483     END DO
484    
485     DO jk = 2, kflev + 1
486     ikm1 = jk - 1
487     DO jl = 1, kdlon
488     zw1(jl) = zw1(jl) + poz(jl, ikm1)*1.66
489     zw4(jl) = zw4(jl) + pud(jl, 4, ikm1)*1.66
490     zw5(jl) = zw5(jl) + pud(jl, 5, ikm1)*1.66
491     ! ZR4(JL) = EXP(-RSWCE*ZW4(JL)-RSWCP*ZW5(JL))
492     END DO
493    
494     CALL swtt(knu, iabs, zw1, zr1)
495    
496     DO jl = 1, kdlon
497     pfup(jl, jk) = ((1.-pclear(jl))*zr1(jl)*zr4(jl)*pfup(jl,jk)+pclear(jl)* &
498     zfu(jl,jk))*rsun(knu)
499     END DO
500     END DO
501    
502     ! ------------------------------------------------------------------
503    
504     RETURN
505     END SUBROUTINE sw2s

  ViewVC Help
Powered by ViewVC 1.1.21