Context Navigation

← Previous Revision
Next Revision →
Blame
Revision Log

optcv.F90

Last change on this file was 3959, checked in by debatzbr, 7 weeks ago
Pluto PCM: Take clouds into account in radiative transfer. BBT
Property svn:executable set to ``*
File size: 19.0 KB

Line
1	MODULE optcv_mod
2
3	IMPLICIT NONE
4
5	CONTAINS
6
7	SUBROUTINE OPTCV(DTAUV,TAUV,TAUCUMV,PLEV, &
8	QXVAER,QSVAER,GVAER,WBARV,COSBV, &
9	TAURAY,TAUAERO,TMID,PMID,TAUGSURF,QVAR,MUVAR,&
10	DTAUV_AER,WBARV_AER)
11
12	use radinc_h, only: L_NLAYRAD, L_NLEVRAD, L_LEVELS, L_NSPECTV, L_NGAUSS, L_REFVAR, NAERKIND
13	use radcommon_h, only: gasv, tlimit, wrefVAR, Cmk, tgasref, pfgasref,wnov,scalep,indv,glat_ig
14	use gases_h, only: gfrac, ngasmx, igas_H2, igas_H2O, igas_He, igas_N2, &
15	igas_CH4, igas_N2
16	use comcstfi_mod, only: g, r, mugaz
17	use callkeys_mod, only: kastprof,continuum,graybody,callgasvis,callmufi,callmuclouds, &
18	Fabs_aers_VI,Fabs_aerf_VI,Fabs_cldd_VI
19	use recombin_corrk_mod, only: corrk_recombin, gasv_recomb
20	use tpindex_mod, only: tpindex
21
22	implicit none
23
24	!==================================================================
25	!
26	! Purpose
27	! -------
28	! Calculates shortwave optical constants at each level.
29	!
30	! Authors
31	! -------
32	! Adapted from the NASA Ames code by R. Wordsworth (2009)
33	!
34	!==================================================================
35	!
36	! THIS SUBROUTINE SETS THE OPTICAL CONSTANTS IN THE VISUAL
37	! IT CALCULATES FOR EACH LAYER, FOR EACH SPECTRAL INTERVAL IN THE VISUAL
38	! LAYER: WBAR, DTAU, COSBAR
39	! LEVEL: TAU
40	!
41	! TAUV(L,NW,NG) is the cumulative optical depth at the top of radiation code
42	! layer L. NW is spectral wavelength interval, ng the Gauss point index.
43	!
44	! TLEV(L) - Temperature at the layer boundary
45	! PLEV(L) - Pressure at the layer boundary (i.e. level)
46	! GASV(NT,NPS,NW,NG) - Visible k-coefficients
47	!
48	!-------------------------------------------------------------------
49
50
51	real*8,intent(out) :: DTAUV(L_NLAYRAD,L_NSPECTV,L_NGAUSS)
52	real*8 DTAUKV(L_LEVELS,L_NSPECTV,L_NGAUSS)
53	real*8,intent(out) :: TAUV(L_NLEVRAD,L_NSPECTV,L_NGAUSS)
54	real*8,intent(out) :: TAUCUMV(L_LEVELS,L_NSPECTV,L_NGAUSS)
55	real*8,intent(in) :: PLEV(L_LEVELS)
56	real*8,intent(in) :: TMID(L_LEVELS), PMID(L_LEVELS)
57	real*8,intent(out) :: COSBV(L_NLAYRAD,L_NSPECTV,L_NGAUSS)
58	real*8,intent(out) :: WBARV(L_NLAYRAD,L_NSPECTV,L_NGAUSS)
59
60	! for aerosols
61	real*8,intent(in) :: QXVAER(L_LEVELS,L_NSPECTV,NAERKIND)
62	real*8,intent(in) :: QSVAER(L_LEVELS,L_NSPECTV,NAERKIND)
63	real*8,intent(in) :: GVAER(L_LEVELS,L_NSPECTV,NAERKIND)
64	real*8,intent(in) :: TAUAERO(L_LEVELS,NAERKIND)
65	real*8,intent(out) :: DTAUV_AER(L_NLAYRAD,L_NSPECTV,NAERKIND)
66	real*8,intent(out) :: WBARV_AER(L_NLAYRAD,L_NSPECTV,NAERKIND)
67
68	! local arrays (saved for convenience as need be allocated)
69	real*8,save,allocatable :: TAUAEROLK(:,:,:)
70	real*8,save,allocatable :: TAEROS(:,:,:)
71	!$OMP THREADPRIVATE(TAUAEROLK,TAEROS)
72
73	integer L, NW, NG, K, LK, IAER
74	integer MT(L_LEVELS), MP(L_LEVELS), NP(L_LEVELS)
75	real*8 ANS, TAUGAS
76	real*8,intent(in) :: TAURAY(L_NSPECTV) ! Rayleigh scattering
77	real*8 TRAY(L_LEVELS,L_NSPECTV) ! Rayleigh scattering
78	real*8 DPR(L_LEVELS), U(L_LEVELS)
79	real*8 LCOEF(4), LKCOEF(L_LEVELS,4)
80
81	real*8,intent(out) :: taugsurf(L_NSPECTV,L_NGAUSS-1)
82	real*8 DCONT,DAERO
83	real*8 DRAYAER
84	double precision wn_cont, p_cont, p_air, T_cont, dtemp, dtempc
85	double precision p_cross
86
87	! variable species mixing ratio variables
88	real*8,intent(in) :: QVAR(L_LEVELS)
89	real*8,intent(in) :: MUVAR(L_LEVELS)
90	real*8 :: WRATIO(L_LEVELS)
91	real*8 KCOEF(4)
92	integer NVAR(L_LEVELS)
93
94	! temporary variables to reduce memory access time to gasv
95	real*8 tmpk(2,2)
96	real*8 tmpkvar(2,2,2)
97
98	! temporary variables for multiple aerosol calculation
99	real*8 atemp(L_NLAYRAD,L_NSPECTV)
100	real*8 btemp(L_NLAYRAD,L_NSPECTV)
101	real*8 ctemp(L_NLAYRAD,L_NSPECTV)
102
103	! variables for k in units m^-1
104	real*8 dz(L_LEVELS)
105
106	! Variables for aerosol absorption
107	real*8 Fabs_aer(NAERKIND)
108	real*8 wbarv_prime
109
110	integer igas, jgas
111
112	integer interm
113
114	logical :: firstcall=.true.
115	!$OMP THREADPRIVATE(firstcall)
116
117	if (firstcall) then
118	! allocate local arrays of size "naerkind" (which are also
119	! "saved" so that this is done only once in for all even if
120	! we don't need to store the value from a time step to the next)
121	allocate(TAUAEROLK(L_LEVELS,L_NSPECTV,NAERKIND))
122	allocate(TAEROS(L_LEVELS,L_NSPECTV,NAERKIND))
123	firstcall=.false.
124	endif ! of if (firstcall)
125
126	!! AS: to save time in computing continuum (see bilinearbig)
127	IF (.not.ALLOCATED(indv)) THEN
128	ALLOCATE(indv(L_NSPECTV,ngasmx,ngasmx))
129	indv = -9999 ! this initial value means "to be calculated"
130	ENDIF
131
132	!=======================================================================
133	! Determine the total gas opacity throughout the column, for each
134	! spectral interval, NW, and each Gauss point, NG.
135	! Calculate the continuum opacities, i.e., those that do not depend on
136	! NG, the Gauss index.
137
138	dtauv(:,:,:) = 0.0
139	tauv(:,:,:) = 0.0
140	taucumv(:,:,:) = 0.0
141
142	taugsurf(:,:) = 0.0
143	dpr(:) = 0.0 ! pressure difference between levels
144	lkcoef(:,:) = 0.0
145	DTAUKV(:,:,:) = 0.0
146	dtauv_aer(:,:,:) = 0.0
147	wbarv_aer(:,:,:) = 1.0
148
149	if (callmufi) then
150	Fabs_aer(1) = Fabs_aers_VI
151	Fabs_aer(2) = Fabs_aerf_VI
152	if (callmuclouds) then
153	Fabs_aer(3) = Fabs_cldd_VI
154	endif
155	else
156	Fabs_aer(:) = 1.0
157	endif
158
159	do K=2,L_LEVELS
160	DPR(k) = PLEV(K)-PLEV(K-1)
161
162	! if we have continuum opacities, we need dz
163	if(kastprof)then
164	dz(k) = dpr(k)(1000.0d08.3145d0/muvar(k))TMID(K)/(gPMID(K))
165	U(k) = CmkDPR(k)mugaz/muvar(k)
166	else
167	dz(k) = dpr(k)RTMID(K)/(glat_igPMID(K))mugaz/muvar(k)
168	U(k) = CmkDPR(k)mugaz/muvar(k) ! only Cmk line in optci.F
169	!JL13 the mugaz/muvar factor takes into account water meanmolecular weight if water is present
170	endif
171
172	call tpindex(PMID(K),TMID(K),QVAR(K),pfgasref,tgasref,WREFVAR, &
173	LCOEF,MT(K),MP(K),NVAR(K),WRATIO(K))
174
175	do LK=1,4
176	LKCOEF(K,LK) = LCOEF(LK)
177	end do
178	end do ! levels
179
180	! Spectral dependance of aerosol absorption
181	!JL18 It seems to be good to have aerosols in the first "radiative layer" of the gcm in the IR
182	! but visible does not handle very well diffusion in first layer.
183	! The tauaero and tauray are thus set to 0 (a small value for rayleigh because the code crashes otherwise)
184	! in the 4 first semilayers in optcv, but not optci.
185	! This solves random variations of the sw heating at the model top.
186	do iaer=1,naerkind
187	do NW=1,L_NSPECTV
188	TAEROS(1:4,NW,IAER)=0.d0
189	do K=5,L_LEVELS
190	TAEROS(K,NW,IAER) = TAUAERO(K,IAER) * QXVAER(K,NW,IAER)
191	! Aerosol/Cloud absorption correction --> impact on opacity.
192	if (callmufi .and. (TAEROS(K,NW,IAER).gt.0.d0)) then
193	TAEROS(K,NW,IAER) = TAEROS(K,NW,IAER) * &
194	((QSVAER(K,NW,IAER)/QXVAER(K,NW,IAER)) + Fabs_aer(IAER)*(1.-(QSVAER(K,NW,IAER)/QXVAER(K,NW,IAER))))
195	endif
196	end do ! L_LEVELS
197	end do ! L_NSPECTV
198	end do ! naerkind
199
200	! Rayleigh scattering
201	do NW=1,L_NSPECTV
202	TRAY(1:4,NW) = 1d-30
203	do K=5,L_LEVELS
204	TRAY(K,NW) = TAURAY(NW) * DPR(K)
205	end do ! levels
206	end do
207
208	! we ignore K=1...
209	do K=2,L_LEVELS
210
211	do NW=1,L_NSPECTV
212
213	DRAYAER = TRAY(K,NW)
214	! DRAYAER is Tau RAYleigh scattering, plus AERosol opacity
215	do iaer=1,naerkind
216	DRAYAER = DRAYAER + TAEROS(K,NW,IAER)
217	end do
218
219	DCONT = 0.0 ! continuum absorption
220
221	if(continuum.and.(.not.graybody).and.callgasvis)then
222	! include continua if necessary
223	wn_cont = dble(wnov(nw))
224	T_cont = dble(TMID(k))
225	do igas=1,ngasmx
226
227	if(gfrac(igas).eq.-1)then ! variable
228	p_cont = dble(PMID(k)scalepQVAR(k)) ! qvar = mol/mol
229	else
230	p_cont = dble(PMID(k)scalepgfrac(igas)*(1.-QVAR(k)))
231	endif
232
233	dtemp=0.0
234	if(igas.eq.igas_N2)then
235
236	interm = indv(nw,igas,igas)
237	! call interpolateN2N2(wn_cont,T_cont,p_cont,dtemp,.false.,interm)
238	indv(nw,igas,igas) = interm
239	! only goes to 500 cm^-1, so unless we're around a cold brown dwarf, this is irrelevant in the visible
240
241	! elseif(igas.eq.igas_H2)then !AF24: removed
242
243	elseif(igas.eq.igas_CH4)then
244
245	! first do self-induced absorption
246	interm = indv(nw,igas,igas)
247	call interpolateCH4CH4(wn_cont,T_cont,p_cont,dtemp,.false.,interm)
248	indv(nw,igas,igas) = interm
249
250	! then cross-interactions with other gases !AF24: removed
251
252	! elseif(igas.eq.igas_H2O.and.T_cont.gt.100.0)then !AF24: removed
253
254	endif
255	DCONT = DCONT + dtemp
256
257	enddo
258
259	DCONT = DCONT*dz(k)
260
261	endif
262
263	do ng=1,L_NGAUSS-1
264
265	! Now compute TAUGAS
266
267	! Interpolate between water mixing ratios
268	! WRATIO = 0.0 if the requested water amount is equal to, or outside the
269	! the water data range
270
271	if(L_REFVAR.eq.1)then ! added by RW for special no variable case
272
273	! JVO 2017 : added tmpk because the repeated calls to gasi/v increased dramatically
274	! the execution time of optci/v -> ~ factor 2 on the whole radiative
275	! transfer on the tested simulations !
276
277	IF (corrk_recombin) THEN ! Added by JVO
278	tmpk = GASV_RECOMB(MT(K):MT(K)+1,MP(K):MP(K)+1,1,NW,NG) ! contains the mix of recombined species
279	ELSE
280	tmpk = GASV(MT(K):MT(K)+1,MP(K):MP(K)+1,1,NW,NG)
281	ENDIF
282
283	KCOEF(1) = tmpk(1,1) ! KCOEF(1) = GASV(MT(K),MP(K),1,NW,NG)
284	KCOEF(2) = tmpk(1,2) ! KCOEF(2) = GASV(MT(K),MP(K)+1,1,NW,NG)
285	KCOEF(3) = tmpk(2,2) ! KCOEF(3) = GASV(MT(K)+1,MP(K)+1,1,NW,NG)
286	KCOEF(4) = tmpk(2,1) ! KCOEF(4) = GASV(MT(K)+1,MP(K),1,NW,NG)
287
288	else
289
290	IF (corrk_recombin) THEN
291	tmpkvar = GASV_RECOMB(MT(K):MT(K)+1,MP(K):MP(K)+1,NVAR(K):NVAR(K)+1,NW,NG)
292	ELSE
293	tmpkvar = GASV(MT(K):MT(K)+1,MP(K):MP(K)+1,NVAR(K):NVAR(K)+1,NW,NG)
294	ENDIF
295
296	KCOEF(1) = tmpkvar(1,1,1) + WRATIO(K) * &
297	( tmpkvar(1,1,2)-tmpkvar(1,1,1) )
298
299	KCOEF(2) = tmpkvar(1,2,1) + WRATIO(K) * &
300	( tmpkvar(1,2,2)-tmpkvar(1,2,1) )
301
302	KCOEF(3) = tmpkvar(2,2,1) + WRATIO(K) * &
303	( tmpkvar(2,2,2)-tmpkvar(2,2,1) )
304
305	KCOEF(4) = tmpkvar(2,1,1) + WRATIO(K) * &
306	( tmpkvar(2,1,2)-tmpkvar(2,1,1) )
307
308
309	endif
310
311	! Interpolate the gaseous k-coefficients to the requested T,P values
312
313	ANS = LKCOEF(K,1)KCOEF(1) + LKCOEF(K,2)KCOEF(2) + &
314	LKCOEF(K,3)KCOEF(3) + LKCOEF(K,4)KCOEF(4)
315
316	TAUGAS = U(k)*ANS
317
318	TAUGSURF(NW,NG) = TAUGSURF(NW,NG) + TAUGAS + DCONT
319	DTAUKV(K,nw,ng) = TAUGAS &
320	+ DRAYAER & ! DRAYAER includes all scattering contributions
321	+ DCONT ! For parameterized continuum aborption
322
323	end do
324
325	! Now fill in the "clear" part of the spectrum (NG = L_NGAUSS),
326	! which holds continuum opacity only
327
328	NG = L_NGAUSS
329	DTAUKV(K,nw,ng) = DRAYAER + DCONT ! Scattering + parameterized continuum absorption
330
331	end do
332	end do
333
334
335	!=======================================================================
336	! Now the full treatment for the layers, where besides the opacity
337	! we need to calculate the scattering albedo and asymmetry factors
338
339	!JL18 It seems to be good to have aerosols in the first "radiative layer" of the gcm in the IR
340	! but not in the visible
341	! The tauaero is thus set to 0 in the 4 first semilayers in optcv, but not optci.
342	! This solves random variations of the sw heating at the model top.
343	do iaer=1,naerkind
344	DO NW=1,L_NSPECTV
345	TAUAEROLK(1:4,NW,IAER)=0.d0
346	DO K=5,L_LEVELS
347	TAUAEROLK(K,NW,IAER) = TAUAERO(K,IAER) * QSVAER(K,NW,IAER) ! effect of scattering albedo
348	ENDDO
349	ENDDO
350	end do
351
352	DO NW=1,L_NSPECTV
353	DO L=1,L_NLAYRAD-1
354	K = 2*L+1
355	atemp(L,NW) = SUM(GVAER(K,NW,1:naerkind) * TAUAEROLK(K,NW,1:naerkind)) + &
356	SUM(GVAER(K+1,NW,1:naerkind) * TAUAEROLK(K+1,NW,1:naerkind))
357	! Aerosol absorption correction --> impact on single scattering albedo.
358	if (callmufi) then
359	! Spherical aerosols
360	if ((TAEROS(K,NW,1).gt.0.d0) .and. TAEROS(K+1,NW,1).gt.0.d0) then
361	btemp(L,NW) = TAUAEROLK(K,NW,1) / &
362	(QSVAER(K,NW,1)/QXVAER(K,NW,1) + &
363	Fabs_aer(1)*(1.-QSVAER(K,NW,1)/QXVAER(K,NW,1))) + &
364	TAUAEROLK(K+1,NW,1) / &
365	(QSVAER(K+1,NW,1)/QXVAER(K+1,NW,1) + &
366	Fabs_aer(1)*(1.-QSVAER(K+1,NW,1)/QXVAER(K+1,NW,1)))
367	else
368	btemp(L,NW) = TAUAEROLK(K,NW,1) + TAUAEROLK(K+1,NW,1)
369	endif
370	! Fractal aerosols
371	if ((TAEROS(K,NW,2).gt.0.d0) .and. TAEROS(K+1,NW,2).gt.0.d0) then
372	btemp(L,NW) = btemp(L,NW) + &
373	(TAUAEROLK(K,NW,2) / &
374	(QSVAER(K,NW,2)/QXVAER(K,NW,2) + &
375	Fabs_aer(2)*(1.-QSVAER(K,NW,2)/QXVAER(K,NW,2)))) + &
376	(TAUAEROLK(K+1,NW,2) / &
377	(QSVAER(K+1,NW,2)/QXVAER(K+1,NW,2) + &
378	Fabs_aer(2)*(1.-QSVAER(K+1,NW,2)/QXVAER(K+1,NW,2))))
379	else
380	btemp(L,NW) = btemp(L,NW) + TAUAEROLK(K,NW,2) + TAUAEROLK(K+1,NW,2)
381	endif
382
383	! Cloud absorption correction --> impact on single scattering albedo.
384	if (callmuclouds) then
385	if ((TAEROS(K,NW,3).gt.0.d0) .and. TAEROS(K+1,NW,3).gt.0.d0) then
386	btemp(L,NW) = btemp(L,NW) + &
387	(TAUAEROLK(K,NW,3) / &
388	(QSVAER(K,NW,3)/QXVAER(K,NW,3) + &
389	Fabs_aer(3)*(1.-QSVAER(K,NW,3)/QXVAER(K,NW,3)))) + &
390	(TAUAEROLK(K+1,NW,3) / &
391	(QSVAER(K+1,NW,3)/QXVAER(K+1,NW,3) + &
392	Fabs_aer(3)*(1.-QSVAER(K+1,NW,3)/QXVAER(K+1,NW,3))))
393	else
394	btemp(L,NW) = btemp(L,NW) + TAUAEROLK(K,NW,3) + TAUAEROLK(K+1,NW,3)
395	endif
396	endif ! end callmuclouds
397	else
398	btemp(L,NW) = SUM(TAUAEROLK(K,NW,1:naerkind)) + SUM(TAUAEROLK(K+1,NW,1:naerkind))
399	endif ! callmufi
400	ctemp(L,NW) = btemp(L,NW) + 0.9999*(TRAY(K,NW) + TRAY(K+1,NW)) ! JVO 2017 : does this 0.999 is really meaningful ?
401
402	btemp(L,NW) = btemp(L,NW) + TRAY(K,NW) + TRAY(K+1,NW)
403	COSBV(L,NW,1:L_NGAUSS) = atemp(L,NW)/btemp(L,NW)
404	END DO ! L vertical loop
405
406	! Last level
407	!-----------
408	L = L_NLAYRAD
409	K = 2*L+1
410	atemp(L,NW) = SUM(GVAER(K,NW,1:naerkind) * TAUAEROLK(K,NW,1:naerkind))
411	! Aerosol absorption correction --> impact on single scattering albedo.
412	if (callmufi) then
413	! Spherical aerosols
414	if (TAEROS(K,NW,1).gt.0.d0) then
415	btemp(L,NW) = TAUAEROLK(K,NW,1) / &
416	(QSVAER(K,NW,1)/QXVAER(K,NW,1) + &
417	Fabs_aer(1)*(1.-QSVAER(K,NW,1)/QXVAER(K,NW,1)))
418	else
419	btemp(L,NW) = TAUAEROLK(K,NW,1)
420	endif
421	! Fractal aerosols
422	if (TAEROS(K,NW,2).gt.0.d0) then
423	btemp(L,NW) = btemp(L,NW) + &
424	(TAUAEROLK(K,NW,2) / &
425	(QSVAER(K,NW,2)/QXVAER(K,NW,2) + &
426	Fabs_aer(2)*(1.-QSVAER(K,NW,2)/QXVAER(K,NW,2))))
427	else
428	btemp(L,NW) = btemp(L,NW) + TAUAEROLK(K,NW,2)
429	endif
430	! Cloud absorption correction --> impact on single scattering albedo.
431	if (callmuclouds) then
432	if (TAEROS(K,NW,3).gt.0.d0) then
433	btemp(L,NW) = btemp(L,NW) + &
434	(TAUAEROLK(K,NW,3) / &
435	(QSVAER(K,NW,3)/QXVAER(K,NW,3) + &
436	Fabs_aer(3)*(1.-QSVAER(K,NW,3)/QXVAER(K,NW,3))))
437	else
438	btemp(L,NW) = btemp(L,NW) + TAUAEROLK(K,NW,3)
439	endif
440	endif ! end callmuclouds
441	else
442	btemp(L,NW) = SUM(TAUAEROLK(K,NW,1:naerkind))
443	endif ! callmufi
444	ctemp(L,NW) = btemp(L,NW) + 0.9999*TRAY(K,NW) ! JVO 2017 : does this 0.999 is really meaningful ?
445
446	btemp(L,NW) = btemp(L,NW) + TRAY(K,NW)
447	COSBV(L,NW,1:L_NGAUSS) = atemp(L,NW)/btemp(L,NW)
448	END DO ! NW spectral loop
449
450	DO NG=1,L_NGAUSS
451	DO NW=1,L_NSPECTV
452	DO L=1,L_NLAYRAD-1
453	K = 2*L+1
454	DTAUV(L,nw,ng) = DTAUKV(K,NW,NG) + DTAUKV(K+1,NW,NG)
455	WBARV(L,nw,ng) = ctemp(L,NW) / DTAUV(L,nw,ng)
456	END DO ! L vertical loop
457
458	! Last level
459	!-----------
460	L = L_NLAYRAD
461	K = 2*L+1
462
463	DTAUV(L,nw,ng) = DTAUKV(K,NW,NG)
464	WBARV(L,NW,NG) = ctemp(L,NW) / DTAUV(L,NW,NG)
465	END DO ! NW spectral loop
466	END DO ! NG Gauss loop
467
468	! Aerosols extinction optical depths
469	DO iaer = 1, naerkind
470	DO nw = 1, L_NSPECTV
471	DO L = 1, L_NLAYRAD-1
472	K = 2*L+1
473	DTAUV_AER(L,nw,iaer) = TAEROS(K,nw,iaer) + TAEROS(K+1,nw,iaer)
474
475	IF (QXVAER(K,nw,iaer) > 0.0D0) THEN
476	wbarv_prime = (QSVAER(K,nw,iaer)/QXVAER(K,nw,iaer)) / &
477	(QSVAER(K,nw,iaer)/QXVAER(K,nw,iaer) + Fabs_aer(iaer)*(1.-QSVAER(K,nw,iaer)/QXVAER(K,nw,iaer)))
478	ELSE
479	wbarv_prime = 1.0
480	ENDIF
481	WBARV_AER(L,nw,iaer) = wbarv_prime * TAEROS(K,nw,iaer)
482	IF (QXVAER(K+1,nw,iaer) > 0.0D0) THEN
483	wbarv_prime = (QSVAER(K+1,nw,iaer)/QXVAER(K+1,nw,iaer)) / &
484	(QSVAER(K+1,nw,iaer)/QXVAER(K+1,nw,iaer) + Fabs_aer(iaer)*(1.-QSVAER(K+1,nw,iaer)/QXVAER(K+1,nw,iaer)))
485	ELSE
486	wbarv_prime = 1.0
487	ENDIF
488	WBARV_AER(L,nw,iaer) = WBARV_AER(L,nw,iaer) + (wbarv_prime * TAEROS(K+1,nw,iaer))
489	IF (DTAUV_AER(L,nw,iaer) > 0.0D0) THEN
490	WBARV_AER(L,nw,iaer) = WBARV_AER(L,nw,iaer) / DTAUV_AER(L,nw,iaer)
491	ELSE
492	WBARV_AER(L,nw,iaer) = 1.0
493	ENDIF
494	END DO ! L vertical loop
495
496	! Last level
497	!-----------
498	L = L_NLAYRAD
499	K = 2*L+1
500	DTAUV_AER(L,nw,iaer) = TAEROS(K,nw,iaer)
501	IF (QXVAER(K,nw,iaer) > 0.0D0) THEN
502	WBARV_AER(L,nw,iaer) = (QSVAER(K,nw,iaer)/QXVAER(K,nw,iaer)) / &
503	(QSVAER(K,nw,iaer)/QXVAER(K,nw,iaer) + Fabs_aer(iaer)*(1.-QSVAER(K,nw,iaer)/QXVAER(K,nw,iaer)))
504	ELSE
505	WBARV_AER(L,nw,iaer) = 1.0
506	ENDIF
507	END DO ! nw spectral loop
508	END DO ! iaer Gauss loop
509
510	! Total extinction optical depths
511	DO NG=1,L_NGAUSS ! full gauss loop
512	DO NW=1,L_NSPECTV
513	TAUV(1,NW,NG)=0.0D0
514	TAUCUMV(1,NW,NG)=0.0D0
515
516	DO K=2,L_LEVELS
517	TAUCUMV(K,NW,NG)=TAUCUMV(K-1,NW,NG)+DTAUKV(K,NW,NG)
518	END DO
519
520	DO L=1,L_NLAYRAD
521	TAUV(L,NW,NG)=TAUCUMV(2*L,NW,NG)
522	END DO
523	TAUV(L,NW,NG)=TAUCUMV(2*L_NLAYRAD+1,NW,NG)
524	END DO
525	END DO ! end full gauss loop
526
527	end subroutine optcv
528
529	END MODULE optcv_mod

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format