Context Navigation

optci.F90 @ 136

Last change on this file since 136 was 135, checked in by aslmd, 14 years ago
CHANGEMENT ARBORESCENCE ETAPE 2 -- NON COMPLET
File size: 9.2 KB

Line
1	subroutine optci(PLEV,TLEV,DTAUI,TAUCUMI, &
2	QXIAER,QSIAER,GIAER,COSBI,WBARI,TAUAERO, &
3	TMID,PMID,TAUGSURF,QVAR)
4
5
6	use radinc_h
7	use radcommon_h, only: gasi, tlimit, wrefVAR, Cmk,tgasref,pfgasref,wnoi
8	implicit none
9
10	!==================================================================
11	!
12	! Purpose
13	! -------
14	! Calculates longwave optical constants at each level. For each
15	! layer and spectral interval in the IR it calculates WBAR, DTAU
16	! and COSBAR. For each level it calculates TAU.
17	!
18	! TAUI(L,LW) is the cumulative optical depth at level L (or alternatively
19	! at the bottom of layer L), LW is the spectral wavelength interval.
20	!
21	! TLEV(L) - Temperature at the layer boundary (i.e., level)
22	! PLEV(L) - Pressure at the layer boundary (i.e., level)
23	!
24	! Authors
25	! -------
26	! Adapted from the NASA Ames code by R. Wordsworth (2009)
27	!
28	!==================================================================
29
30
31	#include "comcstfi.h"
32	#include "callkeys.h"
33
34
35	real*8 DTAUI(L_NLAYRAD,L_NSPECTI,L_NGAUSS)
36	real*8 DTAUKI(L_LEVELS+1,L_NSPECTI,L_NGAUSS)
37	real*8 TAUI(L_NLEVRAD,L_NSPECTI,L_NGAUSS)
38	real*8 TAUCUMI(L_LEVELS,L_NSPECTI,L_NGAUSS)
39	real*8 PLEV(L_LEVELS)
40	real*8 TLEV(L_LEVELS)
41	real*8 TMID(L_LEVELS), PMID(L_LEVELS)
42	real*8 COSBI(L_NLAYRAD,L_NSPECTI,L_NGAUSS)
43	real*8 WBARI(L_NLAYRAD,L_NSPECTI,L_NGAUSS)
44
45	! For aerosols
46	real*8 QXIAER(L_LEVELS+1,L_NSPECTI,NAERKIND)
47	real*8 QSIAER(L_LEVELS+1,L_NSPECTI,NAERKIND)
48	real*8 GIAER(L_LEVELS+1,L_NSPECTI,NAERKIND)
49	real*8 TAUAERO(L_LEVELS+1,NAERKIND)
50	real*8 TAUAEROLK(L_LEVELS+1,L_NSPECTI,NAERKIND)
51	real*8 TAEROS(L_LEVELS,L_NSPECTI,NAERKIND)
52
53	integer L, NW, NG, K, LK, IAER
54	integer MT(L_LEVELS), MP(L_LEVELS), NP(L_LEVELS)
55	real*8 ANS, TAUGAS
56	real*8 DPR(L_LEVELS), U(L_LEVELS)
57	real*8 LCOEF(4), LKCOEF(L_LEVELS,4)
58
59	real*8 taugsurf(L_NSPECTI,L_NGAUSS-1)
60	real*8 dco2
61
62	! Water mixing ratio variables
63	real*8 QVAR(L_LEVELS), WRATIO(L_LEVELS)
64	real*8 KCOEF(4)
65	integer NVAR(L_LEVELS)
66
67	! temporary variables for multiple aerosol calculation
68	real*8 atemp, btemp
69
70	! variables for k in units m^-1
71	real*8 rho, dz
72
73	!=======================================================================
74	! Determine the total gas opacity throughout the column, for each
75	! spectral interval, NW, and each Gauss point, NG.
76
77	DO NG=1,L_NGAUSS-1
78	do NW=1,L_NSPECTI
79	TAUGSURF(NW,NG) = 0.0D0
80	end do
81	end do
82
83	do K=2,L_LEVELS
84	DPR(k) = PLEV(K)-PLEV(K-1)
85
86	rho = PLEV(K)/(R*TMID(K))
87	dz = -DPR(k)/(g*rho)
88	!print*,'rho=',rho
89	!print*,'dz=',dz
90
91	U(k) = Cmk*DPR(k) ! only Cmk line in optci.F
92	! soon to be replaced by m^-1 !!!
93
94	call tpindex(PMID(K),TMID(K),QVAR(K),pfgasref,tgasref,WREFVAR, &
95	LCOEF,MT(K),MP(K),NVAR(K),WRATIO(K))
96
97	do LK=1,4
98	LKCOEF(K,LK) = LCOEF(LK)
99	end do
100
101	DO NW=1,L_NSPECTI
102	do iaer=1,naerkind
103	TAEROS(K,NW,IAER) = TAUAERO(K,IAER) * QXIAER(K,NW,IAER)
104	end do
105	END DO
106	end do ! levels
107
108
109	do K=2,L_LEVELS
110	do nw=1,L_NSPECTI
111
112	DCO2 = 0.0 ! continuum absorption (no longer used)
113
114	!if(varactive)then
115	! call water_cont(PMID(K),TMID(K),WRATIO(K),WNOI(nw),DCO2)
116	!endif
117
118	do ng=1,L_NGAUSS-1
119
120
121	! Now compute TAUGAS
122
123	! Interpolate between water mixing ratios
124	! WRATIO = 0.0 if the requested water amount is equal to, or outside the
125	! the water data range
126
127
128	if (L_REFVAR.eq.1)then ! added by RW for special no variable case
129	KCOEF(1) = GASI(MT(K),MP(K),1,NW,NG)
130	KCOEF(2) = GASI(MT(K),MP(K)+1,1,NW,NG)
131	KCOEF(3) = GASI(MT(K)+1,MP(K)+1,1,NW,NG)
132	KCOEF(4) = GASI(MT(K)+1,MP(K),1,NW,NG)
133	else
134
135	KCOEF(1) = GASI(MT(K),MP(K),NVAR(K),NW,NG) + WRATIO(K)* &
136	(GASI(MT(K),MP(K),NVAR(K)+1,NW,NG) - &
137	GASI(MT(K),MP(K),NVAR(K),NW,NG))
138
139	KCOEF(2) = GASI(MT(K),MP(K)+1,NVAR(K),NW,NG) + WRATIO(K)* &
140	(GASI(MT(K),MP(K)+1,NVAR(K)+1,NW,NG) - &
141	GASI(MT(K),MP(K)+1,NVAR(K),NW,NG))
142
143	KCOEF(3) = GASI(MT(K)+1,MP(K)+1,NVAR(K),NW,NG) + WRATIO(K)* &
144	(GASI(MT(K)+1,MP(K)+1,NVAR(K)+1,NW,NG) - &
145	GASI(MT(K)+1,MP(K)+1,NVAR(K),NW,NG))
146
147	KCOEF(4) = GASI(MT(K)+1,MP(K),NVAR(K),NW,NG) + WRATIO(K)* &
148	(GASI(MT(K)+1,MP(K),NVAR(K)+1,NW,NG) - &
149	GASI(MT(K)+1,MP(K),NVAR(K),NW,NG))
150	endif
151
152	! Interpolate the gaseous k-coefficients to the requested T,P values
153
154	ANS = LKCOEF(K,1)KCOEF(1) + LKCOEF(K,2)KCOEF(2) + &
155	LKCOEF(K,3)KCOEF(3) + LKCOEF(K,4)KCOEF(4)
156
157	TAUGAS = U(k)*ANS
158
159	TAUGSURF(NW,NG) = TAUGSURF(NW,NG) + TAUGAS
160
161	DTAUKI(K,nw,ng) = TAUGAS
162	do iaer=1,naerkind
163	DTAUKI(K,nw,ng) = DTAUKI(K,nw,ng) + TAEROS(K,NW,IAER) &
164	+ DCO2 ! For Kasting CIA
165	end do
166
167	end do
168
169	! Now fill in the "clear" part of the spectrum (NG = L_NGAUSS),
170	! which holds continuum opacity only
171
172	NG = L_NGAUSS
173	DTAUKI(K,nw,ng) = 0.0
174	do iaer=1,naerkind
175	DTAUKI(K,nw,ng) = DTAUKI(K,nw,ng) + TAEROS(K,NW,IAER) &
176	+ DCO2 ! For parameterized continuum absorption
177	end do ! a bug was found here!!
178
179	end do
180	end do
181
182
183	!=======================================================================
184	! Now the full treatment for the layers, where besides the opacity
185	! we need to calculate the scattering albedo and asymmetry factors
186
187	DO NW=1,L_NSPECTI
188	DO K=2,L_LEVELS+1
189	do iaer=1,naerkind
190	TAUAEROLK(K,NW,IAER) = TAUAERO(K,IAER)*QSIAER(K,NW,IAER)
191	end do
192	ENDDO
193	ENDDO
194
195	DO NW=1,L_NSPECTI
196	NG = L_NGAUSS
197	DO L=1,L_NLAYRAD
198
199	K = 2*L+1
200	DTAUI(L,nw,ng) = DTAUKI(K,NW,NG) + DTAUKI(K+1,NW,NG)! + 1.e-50
201
202	atemp = 0.
203	btemp = 0.
204	if(DTAUI(L,NW,NG) .GT. 1.0E-9) then
205	do iaer=1,naerkind
206	atemp = atemp + &
207	GIAER(K,NW,IAER) * TAUAEROLK(K,NW,IAER) + &
208	GIAER(K+1,NW,IAER) * TAUAEROLK(K+1,NW,IAER)
209	btemp = btemp + TAUAEROLK(K,NW,IAER) + TAUAEROLK(K+1,NW,IAER)
210	! * + 1.e-10
211	end do
212	WBARI(L,nw,ng) = btemp / DTAUI(L,NW,NG)
213	else
214	WBARI(L,nw,ng) = 0.0D0
215	DTAUI(L,NW,NG) = 1.0E-9
216	endif
217
218	if(btemp .GT. 0.0) then
219	cosbi(L,NW,NG) = atemp/btemp
220	else
221	cosbi(L,NW,NG) = 0.0D0
222	end if
223
224	END DO ! L vertical loop
225
226	! Now the other Gauss points, if needed.
227
228	DO NG=1,L_NGAUSS-1
229	IF(TAUGSURF(NW,NG) .gt. TLIMIT) THEN
230
231	DO L=1,L_NLAYRAD
232	K = 2*L+1
233	DTAUI(L,nw,ng) = DTAUKI(K,NW,NG)+DTAUKI(K+1,NW,NG)! + 1.e-50
234
235	btemp = 0.
236	if(DTAUI(L,NW,NG) .GT. 1.0E-9) then
237
238	do iaer=1,naerkind
239	btemp = btemp + TAUAEROLK(K,NW,IAER) + TAUAEROLK(K+1,NW,IAER)
240	end do
241	WBARI(L,nw,ng) = btemp / DTAUI(L,NW,NG)
242
243	else
244	WBARI(L,nw,ng) = 0.0D0
245	DTAUI(L,NW,NG) = 1.0E-9
246	endif
247
248	cosbi(L,NW,NG) = cosbi(L,NW,L_NGAUSS)
249	END DO ! L vertical loop
250	END IF
251
252	END DO ! NG Gauss loop
253	END DO ! NW spectral loop
254
255	! Total extinction optical depths
256
257	DO NW=1,L_NSPECTI
258	DO NG=1,L_NGAUSS ! full gauss loop
259	TAUI(1,NW,NG)=0.0D0
260	DO L=1,L_NLAYRAD
261	TAUI(L+1,NW,NG)=TAUI(L,NW,NG)+DTAUI(L,NW,NG)
262	END DO
263
264	TAUCUMI(1,NW,NG)=0.0D0
265	DO K=2,L_LEVELS
266	TAUCUMI(K,NW,NG)=TAUCUMI(K-1,NW,NG)+DTAUKI(K,NW,NG)
267	END DO
268	END DO ! end full gauss loop
269	END DO
270
271
272	return
273	end subroutine optci
274
275
276	!-------------------------------------------------------------------------
277	subroutine water_cont(p,T,wratio,nu,DCO2)
278	! Calculates the continuum opacity for H2O
279	! must check somewhere that it actually _is_ H2O we are using...
280
281	implicit none
282
283	real p, T, wratio, nu, DCO2
284	real h1, h2
285
286	h1 = exp(1800.*(1./T - 0.0034))
287	h2 = 1.25e-22 + 1.67e-19exp(-2.62e-13nu)
288
289	! DCO2 = h1h2p(pwratio)*2/(RT)
290
291	DCO2=0.0 ! to be implemented...
292
293	return
294
295	end subroutine water_cont

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

Download in other formats:

Original Format