/[lmdze]/trunk/phylmd/Radlwsw/swclr.f90
ViewVC logotype

Contents of /trunk/phylmd/Radlwsw/swclr.f90

Parent Directory Parent Directory | Revision Log Revision Log


Revision 328 - (show annotations)
Thu Jun 13 14:40:06 2019 UTC (4 years, 11 months ago) by guez
File size: 7201 byte(s)
Change all `.f` suffixes to `.f90`. (The opposite was done in revision
82.)  Because of change of philosopy in GNUmakefile: we already had a
rewritten rule for `.f`, so it does not make the makefile longer to
replace it by a rule for `.f90`. And it spares us options of
makedepf90 and of the compiler. Also we prepare the way for a simpler
`CMakeLists.txt`.

1 module swclr_m
2
3 IMPLICIT NONE
4
5 contains
6
7 SUBROUTINE swclr(knu, palbp, pdsig, prayl, psec, ppizaz, pray1, pray2, &
8 prefz, prj, prk, prmu0, ptauaz, ptra1, ptra2)
9
10 ! PURPOSE.
11 ! COMPUTES THE REFLECTIVITY AND TRANSMISSIVITY IN CASE OF
12 ! CLEAR-SKY COLUMN
13
14 ! REFERENCE.
15 ! SEE RADIATION'S PART OF THE ECMWF RESEARCH DEPARTMENT
16 ! DOCUMENTATION, AND FOUQUART AND BONNEL (1980)
17
18 ! AUTHOR.
19 ! JEAN-JACQUES MORCRETTE *ECMWF*
20
21 ! MODIFICATIONS.
22 ! ORIGINAL : 94-11-15
23
24 USE raddim, only: kdlon, kflev
25 USE radepsi, only: zepsec
26 USE radopt, only: novlp
27
28 ! ARGUMENTS:
29
30 INTEGER knu
31 DOUBLE PRECISION palbp(kdlon, 2)
32 DOUBLE PRECISION, intent(in):: pdsig(kdlon, kflev)
33 DOUBLE PRECISION, intent(in):: prayl(kdlon)
34 DOUBLE PRECISION psec(kdlon)
35
36 DOUBLE PRECISION, intent(out):: ppizaz(kdlon, kflev)
37 DOUBLE PRECISION pray1(kdlon, kflev + 1)
38 DOUBLE PRECISION pray2(kdlon, kflev + 1)
39 DOUBLE PRECISION prefz(kdlon, 2, kflev + 1)
40 DOUBLE PRECISION prj(kdlon, 6, kflev + 1)
41 DOUBLE PRECISION prk(kdlon, 6, kflev + 1)
42 DOUBLE PRECISION prmu0(kdlon, kflev + 1)
43 DOUBLE PRECISION, intent(out):: ptauaz(kdlon, kflev)
44 DOUBLE PRECISION ptra1(kdlon, kflev + 1)
45 DOUBLE PRECISION ptra2(kdlon, kflev + 1)
46
47 ! LOCAL VARIABLES:
48 DOUBLE PRECISION zc0i(kdlon, kflev + 1)
49 DOUBLE PRECISION zclear(kdlon)
50 DOUBLE PRECISION zr21(kdlon)
51 DOUBLE PRECISION zss0(kdlon)
52 DOUBLE PRECISION zscat(kdlon)
53 DOUBLE PRECISION ztr(kdlon, 2, kflev + 1)
54 INTEGER jl, jk, ja, jkl, jklp1, jaj, jkm1
55 DOUBLE PRECISION zfacoa, zcorae
56 DOUBLE PRECISION zmue, zgap, zww, zto, zden, zmu1, zden1
57 DOUBLE PRECISION zbmu0, zbmu1, zre11
58 double precision, parameter:: REPSCT = 1d-10
59
60 !------------------------------------------------------------------
61
62 ! 1. OPTICAL PARAMETERS FOR AEROSOLS AND RAYLEIGH
63
64 DO jk = 1, kflev + 1
65 DO ja = 1, 6
66 DO jl = 1, kdlon
67 prj(jl, ja, jk) = 0d0
68 prk(jl, ja, jk) = 0d0
69 END DO
70 END DO
71 END DO
72
73 DO jk = 1, kflev
74 DO jl = 1, kdlon
75 ptauaz(jl, jk) = prayl(jl) * pdsig(jl, jk)
76 ppizaz(jl, jk) = 1d0 - repsct
77 END DO
78 END DO
79
80 ! 2. TOTAL EFFECTIVE CLOUDINESS ABOVE A GIVEN LEVEL
81
82 DO jl = 1, kdlon
83 zc0i(jl, kflev + 1) = 0d0
84 zclear(jl) = 1d0
85 zscat(jl) = 0d0
86 END DO
87
88 jk = 1
89 jkl = kflev + 1 - jk
90 jklp1 = jkl + 1
91 DO jl = 1, kdlon
92 zfacoa = 1d0
93 zcorae = zfacoa * ptauaz(jl, jkl) * psec(jl)
94 zr21(jl) = exp(- zcorae)
95 zss0(jl) = 1d0 - zr21(jl)
96
97 IF (novlp == 1) THEN
98 ! maximum-random
99 zclear(jl) = zclear(jl) * (1d0 - max(zss0(jl), zscat(jl))) / (1d0 &
100 - min(zscat(jl), 1d0 - zepsec))
101 zc0i(jl, jkl) = 1d0 - zclear(jl)
102 zscat(jl) = zss0(jl)
103 ELSE IF (novlp == 2) THEN
104 ! maximum
105 zscat(jl) = max(zss0(jl), zscat(jl))
106 zc0i(jl, jkl) = zscat(jl)
107 ELSE IF (novlp == 3) THEN
108 ! random
109 zclear(jl) = zclear(jl) * (1d0 - zss0(jl))
110 zscat(jl) = 1d0 - zclear(jl)
111 zc0i(jl, jkl) = zscat(jl)
112 END IF
113 END DO
114
115 DO jk = 2, kflev
116 jkl = kflev + 1 - jk
117 jklp1 = jkl + 1
118 DO jl = 1, kdlon
119 zfacoa = 1d0
120 zcorae = zfacoa * ptauaz(jl, jkl) * psec(jl)
121 zr21(jl) = exp(- zcorae)
122 zss0(jl) = 1d0 - zr21(jl)
123
124 IF (novlp == 1) THEN
125 ! maximum-random
126 zclear(jl) = zclear(jl) * (1d0 - max(zss0(jl), zscat(jl))) &
127 / (1d0 - min(zscat(jl), 1d0 - zepsec))
128 zc0i(jl, jkl) = 1d0 - zclear(jl)
129 zscat(jl) = zss0(jl)
130 ELSE IF (novlp == 2) THEN
131 ! maximum
132 zscat(jl) = max(zss0(jl), zscat(jl))
133 zc0i(jl, jkl) = zscat(jl)
134 ELSE IF (novlp == 3) THEN
135 ! random
136 zclear(jl) = zclear(jl) * (1d0 - zss0(jl))
137 zscat(jl) = 1d0 - zclear(jl)
138 zc0i(jl, jkl) = zscat(jl)
139 END IF
140 END DO
141 END DO
142
143 ! 3. REFLECTIVITY/TRANSMISSIVITY FOR PURE SCATTERING
144
145 DO jl = 1, kdlon
146 pray1(jl, kflev + 1) = 0d0
147 pray2(jl, kflev + 1) = 0d0
148 prefz(jl, 2, 1) = palbp(jl, knu)
149 prefz(jl, 1, 1) = palbp(jl, knu)
150 ptra1(jl, kflev + 1) = 1d0
151 ptra2(jl, kflev + 1) = 1d0
152 END DO
153
154 DO jk = 2, kflev + 1
155 jkm1 = jk - 1
156 DO jl = 1, kdlon
157
158 ! 3.1 EQUIVALENT ZENITH ANGLE
159
160 zmue = (1d0 - zc0i(jl, jk)) * psec(jl) + zc0i(jl, jk) * 1.66d0
161 prmu0(jl, jk) = 1d0 / zmue
162
163 ! 3.2 REFLECT./TRANSMISSIVITY DUE TO RAYLEIGH AND AEROSOLS
164
165 zgap = 0d0
166 zbmu0 = 0.5d0 - 0.75d0 * zgap / zmue
167 zww = ppizaz(jl, jkm1)
168 zto = ptauaz(jl, jkm1)
169 zden = 1d0 + (1d0 - zww + zbmu0 * zww) * zto * zmue + (1d0 - zww) &
170 * (1d0 - zww + 2d0 * zbmu0 * zww) * zto * zto * zmue * zmue
171 pray1(jl, jkm1) = zbmu0 * zww * zto * zmue / zden
172 ptra1(jl, jkm1) = 1d0 / zden
173
174 zmu1 = 0.5d0
175 zbmu1 = 0.5d0 - 0.75d0 * zgap * zmu1
176 zden1 = 1d0 + (1d0 - zww + zbmu1 * zww) * zto / zmu1 + (1d0 - zww) &
177 * (1d0 - zww + 2d0 * zbmu1 * zww &
178 ) * zto * zto / zmu1 / zmu1
179 pray2(jl, jkm1) = zbmu1 * zww * zto / zmu1 / zden1
180 ptra2(jl, jkm1) = 1d0 / zden1
181
182 prefz(jl, 1, jk) = (pray1(jl, jkm1) + prefz(jl, 1, jkm1) &
183 * ptra1(jl, jkm1)* ptra2(jl, jkm1) / (1d0 - pray2(jl, jkm1) &
184 * prefz(jl, 1, jkm1)))
185
186 ztr(jl, 1, jkm1) = (ptra1(jl, jkm1) / (1d0 - pray2(jl, jkm1) &
187 * prefz(jl, 1, jkm1)))
188
189 prefz(jl, 2, jk) = (pray1(jl, jkm1) + prefz(jl, 2, jkm1) &
190 * ptra1(jl, jkm1) * ptra2(jl, jkm1))
191
192 ztr(jl, 2, jkm1) = ptra1(jl, jkm1)
193
194 END DO
195 END DO
196 DO jl = 1, kdlon
197 zmue = (1d0 - zc0i(jl, 1)) * psec(jl) + zc0i(jl, 1) * 1.66d0
198 prmu0(jl, 1) = 1d0 / zmue
199 END DO
200
201 ! 3.5 REFLECT./TRANSMISSIVITY BETWEEN SURFACE AND LEVEL
202
203 IF (knu == 1) THEN
204 jaj = 2
205 DO jl = 1, kdlon
206 prj(jl, jaj, kflev + 1) = 1d0
207 prk(jl, jaj, kflev + 1) = prefz(jl, 1, kflev + 1)
208 END DO
209
210 DO jk = 1, kflev
211 jkl = kflev + 1 - jk
212 jklp1 = jkl + 1
213 DO jl = 1, kdlon
214 zre11 = prj(jl, jaj, jklp1) * ztr(jl, 1, jkl)
215 prj(jl, jaj, jkl) = zre11
216 prk(jl, jaj, jkl) = zre11 * prefz(jl, 1, jkl)
217 END DO
218 END DO
219 ELSE
220 DO jaj = 1, 2
221 DO jl = 1, kdlon
222 prj(jl, jaj, kflev + 1) = 1d0
223 prk(jl, jaj, kflev + 1) = prefz(jl, jaj, kflev + 1)
224 END DO
225
226 DO jk = 1, kflev
227 jkl = kflev + 1 - jk
228 jklp1 = jkl + 1
229 DO jl = 1, kdlon
230 zre11 = prj(jl, jaj, jklp1) * ztr(jl, jaj, jkl)
231 prj(jl, jaj, jkl) = zre11
232 prk(jl, jaj, jkl) = zre11 * prefz(jl, jaj, jkl)
233 END DO
234 END DO
235 END DO
236 END IF
237
238 END SUBROUTINE swclr
239
240 end module swclr_m

  ViewVC Help
Powered by ViewVC 1.1.21