Context Navigation

sugas_corrk.F90 @ 1956

Last change on this file since 1956 was 1822, checked in by jvatant, 8 years ago
Preliminary modifs for the optical coupling of haze + Moved inits of setspi/v before init of mufi + Added access to tarcers in optci/v + Some coherence in call to directories JVO
File size: 22.0 KB

Line
1	subroutine sugas_corrk
2
3	!==================================================================
4	!
5	! Purpose
6	! -------
7	! Set up gaseous absorption parameters used by the radiation code.
8	! This subroutine is a replacement for the old 'setrad', which contained
9	! both absorption and scattering data.
10	!
11	! Authors
12	! -------
13	! Adapted and generalised from the NASA Ames code by Robin Wordsworth (2009)
14	! Added double gray case by Jeremy Leconte (2012)
15	! New HITRAN continuum data section by RW (2012)
16	!
17	! Summary
18	! -------
19	!
20	!==================================================================
21
22	use radinc_h
23	use radcommon_h, only : pgasref,pfgasref,pgasmin,pgasmax
24	use radcommon_h, only : tgasref,tgasmin,tgasmax
25	use radcommon_h, only : gasv,gasi,FZEROI,FZEROV,gweight
26	use radcommon_h, only : WNOI,WNOV
27	use datafile_mod, only: datadir, corrkdir, banddir
28	use comcstfi_mod, only: mugaz
29	use gases_h
30	use ioipsl_getin_p_mod, only: getin_p
31	use callkeys_mod, only: graybody,callgasvis, continuum
32	implicit none
33
34	!==================================================================
35
36	logical file_ok
37
38	integer n, nt, np, nh, ng, nw, m, i
39	integer L_NGAUSScheck
40
41	character(len=200) :: file_id
42	character(len=500) :: file_path
43
44	! ALLOCATABLE ARRAYS -- AS 12/2011
45	REAL*8, DIMENSION(:,:,:,:,:), ALLOCATABLE,SAVE :: gasi8, gasv8 !read by master
46
47	real*8 x, xi(4), yi(4), ans, p
48	! For gray case (JL12)
49	real kappa_IR, kappa_VI, IR_VI_wnlimit
50	integer nVI_limit,nIR_limit
51
52	integer ngas, igas, jgas
53
54	double precision testcont ! for continuum absorption initialisation
55
56	integer :: dummy
57
58	!=======================================================================
59	! Load variable species data, exit if we have wrong database
60	file_id='/corrk_data/' // TRIM(corrkdir) // '/Q.dat'
61	file_path=TRIM(datadir)//TRIM(file_id)
62
63	! check that the file exists
64	inquire(FILE=file_path,EXIST=file_ok)
65	if(.not.file_ok) then
66	write(,)'The file ',TRIM(file_path)
67	write(,)'was not found by sugas_corrk.F90, exiting.'
68	write(,)'Check that your path to datagcm:',trim(datadir)
69	write(,)' is correct. You can change it in callphys.def with:'
70	write(,)' datadir = /absolute/path/to/datagcm'
71	write(,)'Also check that the corrkdir you chose in callphys.def exists.'
72	call abort
73	endif
74
75	!$OMP MASTER
76	! check that database matches varactive toggle
77	open(111,file=TRIM(file_path),form='formatted')
78	read(111,*) ngas
79
80	if(ngas.gt.5 .or. ngas.lt.1)then
81	print*,ngas,' species in database [', &
82	corrkdir(1:LEN_TRIM(corrkdir)), &
83	'], radiative code cannot handle this.'
84	call abort
85	endif
86
87	L_REFVAR = 1 ! JVO 2017 : set to 1 to keep the code running until the new variable species treatment
88
89	!=======================================================================
90	! Set the weighting in g-space
91
92	file_id='/corrk_data/' // TRIM(corrkdir) // '/g.dat'
93	file_path=TRIM(datadir)//TRIM(file_id)
94
95	! check that the file exists
96	inquire(FILE=file_path,EXIST=file_ok)
97	if(.not.file_ok) then
98	write(,)'The file ',TRIM(file_path)
99	write(,)'was not found by sugas_corrk.F90, exiting.'
100	write(,)'Check that your path to datagcm:',trim(datadir)
101	write(,)' is correct. You can change it in callphys.def with:'
102	write(,)' datadir = /absolute/path/to/datagcm'
103	write(,)'Also check that the corrkdir you chose in callphys.def exists.'
104	call abort
105	endif
106
107	! check the array size is correct, load the coefficients
108	open(111,file=TRIM(file_path),form='formatted')
109	read(111,*) L_NGAUSScheck
110	if(.not.(L_NGAUSScheck.eq.L_NGAUSS)) then
111	print*,'The size of your radiative transfer g-space array does '
112	print*,'not match the value given in g.dat, exiting.'
113	call abort
114	endif
115	read(111,*) gweight
116	close(111)
117
118	! display the values
119	print*,'Correlated-k g-space grid:'
120	do n=1,L_NGAUSS
121	print*,n,'.',gweight(n)
122	end do
123	print*,''
124
125	!=======================================================================
126	! Set the reference pressure and temperature arrays. These are
127	! the pressures and temperatures at which we have k-coefficients.
128
129	!-----------------------------------------------------------------------
130	! pressure
131
132	file_id='/corrk_data/' // TRIM(corrkdir) // '/p.dat'
133	file_path=TRIM(datadir)//TRIM(file_id)
134
135	! check that the file exists
136	inquire(FILE=file_path,EXIST=file_ok)
137	if(.not.file_ok) then
138	write(,)'The file ',TRIM(file_path)
139	write(,)'was not found by sugas_corrk.F90, exiting.'
140	write(,)'Check that your path to datagcm:',trim(datadir)
141	write(,)' is correct. You can change it in callphys.def with:'
142	write(,)' datadir = /absolute/path/to/datagcm'
143	write(,)'Also check that the corrkdir you chose in callphys.def exists.'
144	call abort
145	endif
146
147	! get array size, load the coefficients
148	open(111,file=TRIM(file_path),form='formatted')
149	read(111,*) L_NPREF
150	IF( .NOT. ALLOCATED( pgasref ) ) ALLOCATE( PGASREF(L_NPREF) )
151	read(111,*) pgasref
152	close(111)
153	L_PINT = (L_NPREF-1)*5+1
154	IF( .NOT. ALLOCATED( pfgasref ) ) ALLOCATE( PFGASREF(L_PINT) )
155
156	! display the values
157	print*,'Correlated-k pressure grid (mBar):'
158	do n=1,L_NPREF
159	print*,n,'. 1 x 10^',pgasref(n),' mBar'
160	end do
161	print*,''
162
163	! save the min / max matrix values
164	pgasmin = 10.0**pgasref(1)
165	pgasmax = 10.0**pgasref(L_NPREF)
166
167	! interpolate to finer grid, adapted to uneven grids
168	do n=1,L_NPREF-1
169	do m=1,5
170	pfgasref((n-1)5+m) = pgasref(n)+(m-1)(pgasref(n+1) - pgasref(n))/5.
171	end do
172	end do
173	pfgasref(L_PINT) = pgasref(L_NPREF)
174
175	!-----------------------------------------------------------------------
176	! temperature
177
178	file_id='/corrk_data/' // TRIM(corrkdir) // '/T.dat'
179	file_path=TRIM(datadir)//TRIM(file_id)
180
181	! check that the file exists
182	inquire(FILE=file_path,EXIST=file_ok)
183	if(.not.file_ok) then
184	write(,)'The file ',TRIM(file_path)
185	write(,)'was not found by sugas_corrk.F90, exiting.'
186	write(,)'Check that your path to datagcm:',trim(datadir)
187	write(,)' is correct. You can change it in callphys.def with:'
188	write(,)' datadir = /absolute/path/to/datagcm'
189	write(,)'Also check that the corrkdir you chose in callphys.def exists.'
190	call abort
191	endif
192
193	! get array size, load the coefficients
194	open(111,file=TRIM(file_path),form='formatted')
195	read(111,*) L_NTREF
196	IF( .NOT. ALLOCATED( tgasref ) ) ALLOCATE( TGASREF(L_NTREF) )
197	read(111,*) tgasref
198	close(111)
199
200	! display the values
201	print*,'Correlated-k temperature grid:'
202	do n=1,L_NTREF
203	print*,n,'.',tgasref(n),' K'
204	end do
205
206	! save the min / max matrix values
207	tgasmin = tgasref(1)
208	tgasmax = tgasref(L_NTREF)
209
210	IF( .NOT. ALLOCATED( gasi8 ) ) ALLOCATE( gasi8(L_NTREF,L_NPREF,L_REFVAR,L_NSPECTI,L_NGAUSS) )
211	IF( .NOT. ALLOCATED( gasv8 ) ) ALLOCATE( gasv8(L_NTREF,L_NPREF,L_REFVAR,L_NSPECTV,L_NGAUSS) )
212	!$OMP END MASTER
213	!$OMP BARRIER
214
215	!-----------------------------------------------------------------------
216	! allocate the multidimensional arrays in radcommon_h
217	IF( .NOT. ALLOCATED( gasi ) ) ALLOCATE( gasi(L_NTREF,L_PINT,L_REFVAR,L_NSPECTI,L_NGAUSS) )
218	IF( .NOT. ALLOCATED( gasv ) ) ALLOCATE( gasv(L_NTREF,L_PINT,L_REFVAR,L_NSPECTV,L_NGAUSS) )
219
220	! display the values
221	print*,''
222	print*,'Correlated-k matrix size:'
223	print*,'[',L_NTREF,',',L_NPREF,',',L_REFVAR,',',L_NGAUSS,']'
224
225	!=======================================================================
226	! Get gaseous k-coefficients and interpolate onto finer pressure grid
227
228
229	! wavelength used to separate IR from VI in graybody. We will need that anyway
230	IR_VI_wnlimit=3000.
231	write(,)"graybody: Visible / Infrared separation set at",10000./IR_VI_wnlimit,"um"
232
233	nVI_limit=0
234	Do nw=1,L_NSPECTV
235	if ((WNOV(nw).gt.IR_VI_wnlimit).and.(L_NSPECTV.gt.1)) then
236	nVI_limit=nw-1
237	exit
238	endif
239	End do
240	nIR_limit=L_NSPECTI
241	Do nw=1,L_NSPECTI
242	if ((WNOI(nw).gt.IR_VI_wnlimit).and.(L_NSPECTI.gt.1)) then
243	nIR_limit=nw-1
244	exit
245	endif
246	End do
247
248	if (graybody) then
249	! constant absorption coefficient in visible
250	write(,)"graybody: constant absorption coefficient in visible:"
251	kappa_VI=-100000.
252	call getin_p("kappa_VI",kappa_VI)
253	write(,)" kappa_VI = ",kappa_VI
254	kappa_VI=kappa_VI1.e4 mugaz * 1.672621e-27 ! conversion from m^2/kg to cm^2/molecule
255
256	! constant absorption coefficient in IR
257	write(,)"graybody: constant absorption coefficient in InfraRed:"
258	kappa_IR=-100000.
259	call getin_p("kappa_IR",kappa_IR)
260	write(,)" kappa_IR = ",kappa_IR
261	kappa_IR=kappa_IR1.e4 mugaz * 1.672621e-27 ! conversion from m^2/kg to cm^2/molecule
262
263	write(,)"graybody: Visible / Infrared separation set at band: IR=",nIR_limit,", VI=",nVI_limit
264
265	Else
266	kappa_VI=1.e-30
267	kappa_IR=1.e-30
268	End if
269
270	!$OMP MASTER
271	! print*,corrkdir(1:4)
272	! VISIBLE
273	if (callgasvis) then
274	if ((corrkdir(1:4).eq.'null'))then !(TRIM(corrkdir).eq.'null_LowTeffStar')) then
275	gasv8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,1:L_NSPECTV,1:L_NGAUSS)=0.0
276	print*,'using no corrk data'
277	print*,'Visible corrk gaseous absorption is set to zero if graybody=F'
278	else
279	file_id='/corrk_data/'//trim(adjustl(banddir))//'/corrk_gcm_VI.dat'
280	file_path=TRIM(datadir)//TRIM(file_id)
281
282	! check that the file exists
283	inquire(FILE=file_path,EXIST=file_ok)
284	if(.not.file_ok) then
285	write(,)'The file ',TRIM(file_path)
286	write(,)'was not found by sugas_corrk.F90.'
287	write(,)'Are you sure you have absorption data for these bands?'
288	call abort
289	endif
290
291	open(111,file=TRIM(file_path),form='formatted')
292	read(111,*) gasv8
293	close(111)
294	end if
295
296	if(nVI_limit.eq.0) then
297	gasv8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,1:L_NSPECTV,1:L_NGAUSS)= &
298	gasv8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,1:L_NSPECTV,1:L_NGAUSS)+kappa_VI
299	else if (nVI_limit.eq.L_NSPECTV) then
300	gasv8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,1:L_NSPECTV,1:L_NGAUSS)= &
301	gasv8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,1:L_NSPECTV,1:L_NGAUSS)+kappa_IR
302	else
303	gasv8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,1:nVI_limit,1:L_NGAUSS)= &
304	gasv8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,1:nVI_limit,1:L_NGAUSS)+kappa_IR
305	gasv8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,nVI_limit+1:L_NSPECTV,1:L_NGAUSS)= &
306	gasv8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,nVI_limit+1:L_NSPECTV,1:L_NGAUSS)+kappa_VI
307	end if
308	else
309	print*,'Visible corrk gaseous absorption is set to zero.'
310	gasv8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,1:L_NSPECTV,1:L_NGAUSS)=0.0
311	endif
312	!$OMP END MASTER
313	!$OMP BARRIER
314
315	! INFRA-RED
316	if ((corrkdir(1:4).eq.'null'))then !.or.(TRIM(corrkdir).eq.'null_LowTeffStar')) then
317	print*,'Infrared corrk gaseous absorption is set to zero if graybody=F'
318	!$OMP MASTER
319	gasi8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,1:L_NSPECTI,1:L_NGAUSS)=0.0
320	!$OMP END MASTER
321	!$OMP BARRIER
322	else
323	file_id='/corrk_data/'//trim(adjustl(banddir))//'/corrk_gcm_IR.dat'
324	file_path=TRIM(datadir)//TRIM(file_id)
325
326	! check that the file exists
327	inquire(FILE=file_path,EXIST=file_ok)
328	if(.not.file_ok) then
329	write(,)'The file ',TRIM(file_path)
330	write(,)'was not found by sugas_corrk.F90.'
331	write(,)'Are you sure you have absorption data for these bands?'
332	call abort
333	endif
334
335	!$OMP MASTER
336	open(111,file=TRIM(file_path),form='formatted')
337	read(111,*) gasi8
338	close(111)
339	!$OMP END MASTER
340	!$OMP BARRIER
341
342	! 'fzero' is a currently unused feature that allows optimisation
343	! of the radiative transfer by neglecting bands where absorption
344	! is close to zero. As it could be useful in the future, this
345	! section of the code has been kept commented and not erased.
346	! RW 7/3/12.
347
348	do nw=1,L_NSPECTI
349	fzeroi(nw) = 0.d0
350	! do nt=1,L_NTREF
351	! do np=1,L_NPREF
352	! do nh=1,L_REFVAR
353	! do ng = 1,L_NGAUSS
354	! if(gasi8(nt,np,nh,nw,ng).lt.1.0e-25)then
355	! fzeroi(nw)=fzeroi(nw)+1.d0
356	! endif
357	! end do
358	! end do
359	! end do
360	! end do
361	! fzeroi(nw)=fzeroi(nw)/dble(L_NTREFL_NPREFL_REFVAR*L_NGAUSS)
362	end do
363
364	do nw=1,L_NSPECTV
365	fzerov(nw) = 0.d0
366	! do nt=1,L_NTREF
367	! do np=1,L_NPREF
368	! do nh=1,L_REFVAR
369	! do ng = 1,L_NGAUSS
370	! if(gasv8(nt,np,nh,nw,ng).lt.1.0e-25)then
371	! fzerov(nw)=fzerov(nw)+1.d0
372	! endif
373	! end do
374	! end do
375	! end do
376	! end do
377	! fzerov(nw)=fzerov(nw)/dble(L_NTREFL_NPREFL_REFVAR*L_NGAUSS)
378	end do
379
380	endif
381
382	!$OMP MASTER
383	if(nIR_limit.eq.0) then
384	gasi8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,1:L_NSPECTI,1:L_NGAUSS)= &
385	gasi8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,1:L_NSPECTI,1:L_NGAUSS)+kappa_VI
386	else if (nIR_limit.eq.L_NSPECTI) then
387	gasi8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,1:L_NSPECTI,1:L_NGAUSS)= &
388	gasi8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,1:L_NSPECTI,1:L_NGAUSS)+kappa_IR
389	else
390	gasi8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,1:nIR_limit,1:L_NGAUSS)= &
391	gasi8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,1:nIR_limit,1:L_NGAUSS)+kappa_IR
392	gasi8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,nIR_limit+1:L_NSPECTI,1:L_NGAUSS)= &
393	gasi8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,nIR_limit+1:L_NSPECTI,1:L_NGAUSS)+kappa_VI
394	end if
395
396
397	! Take log10 of the values - this is what we will interpolate.
398	! Smallest value is 1.0E-200.
399
400	do nt=1,L_NTREF
401	do np=1,L_NPREF
402	do nh=1,L_REFVAR
403	do ng = 1,L_NGAUSS
404
405	do nw=1,L_NSPECTV
406	if(gasv8(nt,np,nh,nw,ng).gt.1.0d-200) then
407	gasv8(nt,np,nh,nw,ng) = log10(gasv8(nt,np,nh,nw,ng))
408	else
409	gasv8(nt,np,nh,nw,ng) = -200.0
410	end if
411	end do
412
413	do nw=1,L_NSPECTI
414	if(gasi8(nt,np,nh,nw,ng).gt.1.0d-200) then
415	gasi8(nt,np,nh,nw,ng) = log10(gasi8(nt,np,nh,nw,ng))
416	else
417	gasi8(nt,np,nh,nw,ng) = -200.0
418	end if
419	end do
420
421	end do
422	end do
423	end do
424	end do
425	!$OMP END MASTER
426	!$OMP BARRIER
427
428	! Interpolate the values: first the longwave
429
430	do nt=1,L_NTREF
431	do nh=1,L_REFVAR
432	do nw=1,L_NSPECTI
433	do ng=1,L_NGAUSS
434
435	! First, the initial interval
436
437	n = 1
438	do m=1,5
439	x = pfgasref(m)
440	xi(1) = pgasref(n)
441	xi(2) = pgasref(n+1)
442	xi(3) = pgasref(n+2)
443	xi(4) = pgasref(n+3)
444	yi(1) = gasi8(nt,n,nh,nw,ng)
445	yi(2) = gasi8(nt,n+1,nh,nw,ng)
446	yi(3) = gasi8(nt,n+2,nh,nw,ng)
447	yi(4) = gasi8(nt,n+3,nh,nw,ng)
448	call lagrange(x,xi,yi,ans)
449	gasi(nt,m,nh,nw,ng) = 10.0**ans
450	end do
451
452	do n=2,L_NPREF-2
453	do m=1,5
454	i = (n-1)*5+m
455	x = pfgasref(i)
456	xi(1) = pgasref(n-1)
457	xi(2) = pgasref(n)
458	xi(3) = pgasref(n+1)
459	xi(4) = pgasref(n+2)
460	yi(1) = gasi8(nt,n-1,nh,nw,ng)
461	yi(2) = gasi8(nt,n,nh,nw,ng)
462	yi(3) = gasi8(nt,n+1,nh,nw,ng)
463	yi(4) = gasi8(nt,n+2,nh,nw,ng)
464	call lagrange(x,xi,yi,ans)
465	gasi(nt,i,nh,nw,ng) = 10.0**ans
466	end do
467	end do
468
469	! Now, get the last interval
470
471	n = L_NPREF-1
472	do m=1,5
473	i = (n-1)*5+m
474	x = pfgasref(i)
475	xi(1) = pgasref(n-2)
476	xi(2) = pgasref(n-1)
477	xi(3) = pgasref(n)
478	xi(4) = pgasref(n+1)
479	yi(1) = gasi8(nt,n-2,nh,nw,ng)
480	yi(2) = gasi8(nt,n-1,nh,nw,ng)
481	yi(3) = gasi8(nt,n,nh,nw,ng)
482	yi(4) = gasi8(nt,n+1,nh,nw,ng)
483	call lagrange(x,xi,yi,ans)
484	gasi(nt,i,nh,nw,ng) = 10.0**ans
485	end do
486
487	! Fill the last pressure point
488
489	gasi(nt,L_PINT,nh,nw,ng) = &
490	10.0**gasi8(nt,L_NPREF,nh,nw,ng)
491
492	end do
493	end do
494	end do
495	end do
496
497	! Interpolate the values: now the shortwave
498
499	do nt=1,L_NTREF
500	do nh=1,L_REFVAR
501	do nw=1,L_NSPECTV
502	do ng=1,L_NGAUSS
503
504	! First, the initial interval
505
506	n = 1
507	do m=1,5
508	x = pfgasref(m)
509	xi(1) = pgasref(n)
510	xi(2) = pgasref(n+1)
511	xi(3) = pgasref(n+2)
512	xi(4) = pgasref(n+3)
513	yi(1) = gasv8(nt,n,nh,nw,ng)
514	yi(2) = gasv8(nt,n+1,nh,nw,ng)
515	yi(3) = gasv8(nt,n+2,nh,nw,ng)
516	yi(4) = gasv8(nt,n+3,nh,nw,ng)
517	call lagrange(x,xi,yi,ans)
518	gasv(nt,m,nh,nw,ng) = 10.0**ans
519	end do
520
521	do n=2,L_NPREF-2
522	do m=1,5
523	i = (n-1)*5+m
524	x = pfgasref(i)
525	xi(1) = pgasref(n-1)
526	xi(2) = pgasref(n)
527	xi(3) = pgasref(n+1)
528	xi(4) = pgasref(n+2)
529	yi(1) = gasv8(nt,n-1,nh,nw,ng)
530	yi(2) = gasv8(nt,n,nh,nw,ng)
531	yi(3) = gasv8(nt,n+1,nh,nw,ng)
532	yi(4) = gasv8(nt,n+2,nh,nw,ng)
533	call lagrange(x,xi,yi,ans)
534	gasv(nt,i,nh,nw,ng) = 10.0**ans
535	end do
536	end do
537
538	! Now, get the last interval
539
540	n = L_NPREF-1
541	do m=1,5
542	i = (n-1)*5+m
543	x = pfgasref(i)
544	xi(1) = pgasref(n-2)
545	xi(2) = pgasref(n-1)
546	xi(3) = pgasref(n)
547	xi(4) = pgasref(n+1)
548	yi(1) = gasv8(nt,n-2,nh,nw,ng)
549	yi(2) = gasv8(nt,n-1,nh,nw,ng)
550	yi(3) = gasv8(nt,n,nh,nw,ng)
551	yi(4) = gasv8(nt,n+1,nh,nw,ng)
552	call lagrange(x,xi,yi,ans)
553	gasv(nt,i,nh,nw,ng) = 10.0**ans
554	end do
555
556	! Fill the last pressure point
557
558	gasv(nt,L_PINT,nh,nw,ng) = &
559	10.0**gasv8(nt,L_NPREF,nh,nw,ng)
560
561	end do
562	end do
563	end do
564	end do
565
566
567	!=======================================================================
568	! Initialise the continuum absorption data
569	if(continuum)then
570	do igas=1,ngasmx
571
572	if (igas .eq. igas_N2) then
573
574	dummy = -9999
575	call interpolateN2N2(100.D+0,250.D+0,17500.D+0,testcont,.true.,dummy)
576
577	elseif (igas .eq. igas_H2) then
578
579	! first do self-induced absorption
580	dummy = -9999
581	call interpolateH2H2(500.D+0,250.D+0,17500.D+0,testcont,.true.,dummy)
582	! then cross-interactions with other gases
583	do jgas=1,ngasmx
584	if (jgas .eq. igas_N2) then
585	dummy = -9999
586	call interpolateN2H2(592.D+0,278.15D+0,200000.D+0,10000.D+0,testcont,.true.,dummy)
587	endif
588	enddo
589
590	elseif (igas .eq. igas_CH4) then
591
592	! first do self-induced absorption
593	dummy = -9999
594	call interpolateCH4CH4(200.D+0,200.D+0,7500.D+0,testcont,.true.,dummy)
595	! then cross-interactions with other gases
596	do jgas=1,ngasmx
597	if (jgas .eq. igas_N2) then
598	dummy = -9999
599	call interpolateN2CH4(200.D+0,250.0D+0,100000.D+0,5000.D+0,testcont,.true.,dummy)
600	endif
601	enddo
602
603	endif
604
605	enddo
606	endif
607
608	print*,'----------------------------------------------------'
609	print*,'And that`s all we have. It`s possible that other'
610	print*,'continuum absorption may be present, but if it is we'
611	print*,'don`t yet have data for it...'
612	print*,''
613
614	! Deallocate local arrays
615	!$OMP BARRIER
616	!$OMP MASTER
617	IF( ALLOCATED( gasi8 ) ) DEALLOCATE( gasi8 )
618	IF( ALLOCATED( gasv8 ) ) DEALLOCATE( gasv8 )
619	IF( ALLOCATED( pgasref ) ) DEALLOCATE( pgasref )
620	!$OMP END MASTER
621	!$OMP BARRIER
622
623	return
624	end subroutine sugas_corrk

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

Download in other formats:

Original Format