Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

concoutput_reg.f90 @ 16

Last change on this file since 16 was 16, checked in by jebri, 11 years ago
sources for flexwrf v3.1
File size: 28.0 KB

Line
1	!***********************************************************************
2	!* Copyright 2012,2013 *
3	!* Jerome Brioude, Delia Arnold, Andreas Stohl, Wayne Angevine, *
4	!* John Burkhart, Massimo Cassiani, Adam Dingwell, Richard C Easter, Sabine Eckhardt,*
5	!* Stephanie Evan, Jerome D Fast, Don Morton, Ignacio Pisso, *
6	!* Petra Seibert, Gerard Wotawa, Caroline Forster, Harald Sodemann, *
7	!* *
8	!* This file is part of FLEXPART WRF *
9	!* *
10	!* FLEXPART is free software: you can redistribute it and/or modify *
11	!* it under the terms of the GNU General Public License as published by*
12	!* the Free Software Foundation, either version 3 of the License, or *
13	!* (at your option) any later version. *
14	!* *
15	!* FLEXPART is distributed in the hope that it will be useful, *
16	!* but WITHOUT ANY WARRANTY; without even the implied warranty of *
17	!* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
18	!* GNU General Public License for more details. *
19	!* *
20	!* You should have received a copy of the GNU General Public License *
21	!* along with FLEXPART. If not, see <http://www.gnu.org/licenses/>. *
22	!***********************************************************************
23	subroutine concoutput_reg(itime,outnum,gridtotalunc,wetgridtotalunc, &
24	drygridtotalunc)
25	! i i o o
26	! o
27	!*******************************************************************************
28	! *
29	! Note: This is the FLEXPART_WRF version of subroutine concoutput *
30	! *
31	! Output of the concentration grid and the receptor concentrations. *
32	! *
33	! Author: A. Stohl *
34	! *
35	! 24 May 1995 *
36	! *
37	! 13 April 1999, Major update: if output size is smaller, dump output *
38	! in sparse matrix format; additional output of uncertainty *
39	! *
40	! 05 April 2000, Major update: output of age classes; output for backward *
41	! runs is time spent in grid cell times total mass of *
42	! species. *
43	! *
44	! 17 February 2002, Appropriate dimensions for backward and forward runs *
45	! are now specified in file includepar *
46	! *
47	! Dec 2005, J. Fast - Output files can be either binary or ascii. *
48	! Sparse output option is turned off. *
49	! Dec 2005, R. Easter - changed names of "lon0" & "lat0" variables *
50	! 2012, J. Brioude- modify output format to flexpart 8, latlon regular output
51	! *
52	!*******************************************************************************
53	! *
54	! Variables: *
55	! outnum number of samples *
56	! ncells number of cells with non-zero concentrations *
57	! sparse .true. if in sparse matrix format, else .false. *
58	! nspeciesdim either nspec (forward runs), or numpoint (backward runs) *
59	! tot_mu 1 for forward, initial mass mixing ration for backw. runs *
60	! maxpointspec maxspec for forward runs, maxpoint for backward runs *
61	! *
62	!*******************************************************************************
63
64	! include 'includepar'
65	! include 'includecom'
66	!
67	! double precision jul
68	! integer itime,i,ix,jy,kz,k,l,iix,jjy,kzz,nage,jjjjmmdd,ihmmss
69	! integer ncells(maxpointspec,maxageclass)
70	! integer ncellsd(maxpointspec,maxageclass)
71	! integer ncellsw(maxpointspec,maxageclass),nspeciesdim
72	! real outnum,weightair,densityoutrecept(maxreceptor),xl,yl
73	! real densityoutgrid(0:maxxgrid-1,0:maxygrid-1,maxzgrid),
74	! +grid(0:maxxgrid-1,0:maxygrid-1,maxzgrid,maxpointspec,maxageclass)
75	! real wetgrid(0:maxxgrid-1,0:maxygrid-1,maxpointspec,maxageclass)
76	! real drygrid(0:maxxgrid-1,0:maxygrid-1,maxpointspec,maxageclass)
77	! real gridsigma(0:maxxgrid-1,0:maxygrid-1,maxzgrid,maxpointspec,
78	! +maxageclass),
79	! +drygridsigma(0:maxxgrid-1,0:maxygrid-1,maxpointspec,maxageclass),
80	! +wetgridsigma(0:maxxgrid-1,0:maxygrid-1,maxpointspec,maxageclass)
81	! real auxgrid(nclassunc),gridtotal,gridsigmatotal,gridtotalunc
82	! real wetgridtotal,wetgridsigmatotal,wetgridtotalunc
83	! real drygridtotal,drygridsigmatotal,drygridtotalunc
84	! real factor(0:maxxgrid-1,0:maxygrid-1,maxzgrid)
85	! real halfheight,dz,dz1,dz2,tot_mu(maxpointspec)
86	! real xnelat,xnelon
87	! real xsw,xne,ysw,yne,tmpx,tmpy,tmplon,tmplat
88	! parameter(weightair=28.97)
89	! logical sparse(maxpointspec,maxageclass)
90	! logical sparsed(maxpointspec,maxageclass)
91	! logical sparsew(maxpointspec,maxageclass)
92	! character adate8,atime6
93	use unc_mod
94	use point_mod
95	use outg_mod
96	use par_mod
97	use com_mod
98
99	implicit none
100
101	real(kind=dp) :: jul
102	integer :: itime,i,ix,jy,kz,ks,kp,l,iix,jjy,kzz,nage,jjjjmmdd,ihmmss
103	integer :: sp_count_i,sp_count_r
104	real :: sp_fact
105	real :: outnum,densityoutrecept(maxreceptor),xl,yl,xl2,yl2
106	real :: auxgrid(nclassunc),gridtotal,gridsigmatotal,gridtotalunc
107	real :: wetgridtotal,wetgridsigmatotal,wetgridtotalunc
108	real :: drygridtotal,drygridsigmatotal,drygridtotalunc
109	real :: halfheight,dz,dz1,dz2,tot_mu(maxspec,maxpointspec_act)
110	real :: xsw,xne,ysw,yne,tmpx,tmpy,tmplon,tmplat
111	real,parameter :: smallnum = tiny(0.0) ! smallest number that can be handled
112	! real,parameter :: weightair=28.97 !AD: moved this to par_mod.f90
113	logical :: sp_zer
114	character :: adate8,atime6
115	character(len=3) :: anspec
116
117
118	! Determine current calendar date, needed for the file name
119	!**********************************************************
120
121	jul=bdate+real(itime,kind=dp)/86400._dp
122
123	call caldate(jul,jjjjmmdd,ihmmss)
124	write(adate,'(i8.8)') jjjjmmdd
125	write(atime,'(i6.6)') ihmmss
126	write(unitdates,'(a)') adate//atime
127
128
129	! For forward simulations, output fields have dimension MAXSPEC,
130	! for backward simulations, output fields have dimension MAXPOINT.
131	! Thus, make loops either about nspec, or about numpoint
132	!*****************************************************************
133
134	if (ldirect.eq.1) then
135	do ks=1,nspec
136	do kp=1,maxpointspec_act
137	tot_mu(ks,kp)=1
138	end do
139	end do
140	else
141	do ks=1,nspec
142	do kp=1,maxpointspec_act
143	tot_mu(ks,kp)=xmass(kp,ks)
144	end do
145	end do
146	endif
147
148	!*******************************************************************
149	! Compute air density: sufficiently accurate to take it
150	! from coarse grid at some time
151	! Determine center altitude of output layer, and interpolate density
152	! data to that altitude
153	!*******************************************************************
154
155	do kz=1,numzgrid
156	if (kz.eq.1) then
157	halfheight=outheight(1)/2.
158	else
159	halfheight=(outheight(kz)+outheight(kz-1))/2.
160	endif
161	do kzz=2,nz
162	if ((height(kzz-1).lt.halfheight).and. &
163	(height(kzz).gt.halfheight)) goto 46
164	end do
165	46 kzz=max(min(kzz,nz),2)
166	dz1=halfheight-height(kzz-1)
167	dz2=height(kzz)-halfheight
168	dz=dz1+dz2
169	do jy=0,numygrid-1
170	do ix=0,numxgrid-1
171	! xl=out_xm0+float(ix)*dxout
172	! yl=out_ym0+float(jy)*dyout
173	! xl=(xl-xmet0)/dx
174	! yl=(yl-ymet0)/dx
175	xl2=outlon0+float(ix)*dxoutl !long
176	yl2=outlat0+float(jy)*dyoutl !lat
177	call ll_to_xymeter_wrf(xl2,yl2,xl,yl) !xl is coord
178	xl=(xl-xmet0)/dx
179	yl=(yl-ymet0)/dy
180
181	iix=max(min(nint(xl),nxmin1),0)
182	jjy=max(min(nint(yl),nymin1),0)
183	densityoutgrid(ix,jy,kz)=(rho(iix,jjy,kzz,2)*dz1+ &
184	rho(iix,jjy,kzz-1,2)*dz2)/dz
185	end do
186	end do
187	end do
188
189	do i=1,numreceptor
190	xl=xreceptor(i)
191	yl=yreceptor(i)
192	iix=max(min(nint(xl),nxmin1),0)
193	jjy=max(min(nint(yl),nymin1),0)
194	densityoutrecept(i)=rho(iix,jjy,1,2)
195	end do
196
197	! Output is different for forward and backward simulations
198	do kz=1,numzgrid
199	do jy=0,numygrid-1
200	do ix=0,numxgrid-1
201	if (ldirect.eq.1) then
202	factor3d(ix,jy,kz)=1.e12/volume(ix,jy,kz)/outnum
203	else
204	factor3d(ix,jy,kz)=real(abs(loutaver))/outnum
205	endif
206	end do
207	end do
208	end do
209
210	!*********************************************************************
211	! Determine the standard deviation of the mean concentration or mixing
212	! ratio (uncertainty of the output) and the dry and wet deposition
213	!*********************************************************************
214
215	gridtotal=0.
216	gridsigmatotal=0.
217	gridtotalunc=0.
218	wetgridtotal=0.
219	wetgridsigmatotal=0.
220	wetgridtotalunc=0.
221	drygridtotal=0.
222	drygridsigmatotal=0.
223	drygridtotalunc=0.
224
225	!*******************************************************************
226	! Generate output: may be in concentration (ng/m3) or in mixing
227	! ratio (ppt) or both
228	! Output either in full grid dump or sparse matrix format
229	! For backward simulations, the unit is seconds, stored in grid_conc
230	!*******************************************************************
231
232	! Concentration output
233	!*********************
234	!
235	! open(23,file=path(1)(1:length(1))//'latlon.txt' &
236	! ,form='formatted')
237	! open(24,file=path(1)(1:length(1))//'latlon_corner.txt' &
238	! ,form='formatted')
239	!
240	!! xnelat=outgrid_nelat
241	!! xnelon=outgrid_nelon
242	! call ll_to_xymeter_wrf(outgrid_swlon,outgrid_swlat,xsw,ysw)
243	! call ll_to_xymeter_wrf(outgrid_nelon,outgrid_nelat,xne,yne)
244	! do jy=1,numygrid
245	! do ix=1,numxgrid
246	!! tmpx=out_xm0+(ix-1)*dxout
247	!! tmpy=out_ym0+(jy-1)*dyout
248	!! tmpx=out_xm0+(float(ix)-0.5)*dxout
249	!! tmpy=out_ym0+(float(jy)-0.5)*dyout
250	! tmpx=xsw+(xne-xsw)*float(ix-1)/float(numxgrid-1)
251	! tmpy=ysw+(yne-ysw)*float(jy-1)/float(numygrid-1)
252	!! print*,'jb','tmpx','tmpy',dxout,dyout,ix,jy
253	! call xymeter_to_ll_wrf(tmpx,tmpy,tmplon,tmplat)
254	! xl2=outlon0+(float(ix)-0.5)*dxoutl !long
255	! yl2=outlat0+(float(jy)-0.5)*dyoutl !lat
256	!
257	!!jb if(iouttype.eq.0) write(unitoutgrid) tmplon,tmplat
258	!! if(iouttype.eq.1) write(unitoutgrid,*) tmplon,tmplat
259	!! write(23,*) tmplon,tmplat
260	! write(23,*) xl2,yl2
261	!! tmpx=out_xm0+(ix-1-0.5)*dxout
262	!! tmpy=out_ym0+(jy-1-0.5)*dyout
263	!! tmpx=out_xm0+(float(ix)-1.)*dxout
264	!! tmpy=out_ym0+(float(jy)-1.)*dyout
265	! xl2=outlon0+float(ix-1)*dxoutl !long
266	! yl2=outlat0+float(jy-1)*dyoutl !lat
267	! write(24,*) xl2,yl2
268	!! tmpx=xsw+(xne-xsw)*float(ix-1)/float(numxgrid-1)
269	!! tmpy=ysw+(yne-ysw)*float(jy-1)/float(numygrid-1)
270	!! print*,'jb2','tmpx','tmpy',dxout,dyout,ix,jy
271	!! call xymeter_to_ll_wrf(tmpx,tmpy,tmplon,tmplat)
272	!! call xymeter_to_ll_wrf_out(tmpx,tmpy,tmplon,tmplat)
273	!! write(24,*) tmplon,tmplat
274	! enddo
275	! enddo
276	! close(23)
277	! close(24)
278
279	do ks=1,nspec
280	write(anspec,'(i3.3)') ks
281	if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5)) then
282	if (ldirect.eq.1) then
283	if (iouttype.eq.0) &
284	open(unitoutgrid,file=path(1)(1:length(1))//'grid_conc_'//adate// &
285	atime//'_'//anspec,form='unformatted')
286	if (iouttype.eq.1) &
287	open(unitoutgrid,file=path(1)(1:length(1))//'grid_conc_'//adate// &
288	atime//'_'//anspec,form='formatted')
289	else
290	if (iouttype.eq.0) &
291	open(unitoutgrid,file=path(1)(1:length(1))//'grid_time_'//adate// &
292	atime//'_'//anspec,form='unformatted')
293	if (iouttype.eq.1) &
294	open(unitoutgrid,file=path(1)(1:length(1))//'grid_time_'//adate// &
295	atime//'_'//anspec,form='formatted')
296	endif
297	write(unitoutgrid) itime
298	endif
299
300	if ((iout.eq.2).or.(iout.eq.3)) then ! mixing ratio
301	if (iouttype.eq.0) &
302	open(unitoutgridppt,file=path(1)(1:length(1))//'grid_pptv_'//adate// &
303	atime//'_'//anspec,form='unformatted')
304	if (iouttype.eq.1) &
305	open(unitoutgridppt,file=path(1)(1:length(1))//'grid_pptv_'//adate// &
306	atime//'_'//anspec,form='formatted')
307
308	write(unitoutgridppt) itime
309	endif
310
311
312	do kp=1,maxpointspec_act
313	do nage=1,nageclass
314
315	do jy=0,numygrid-1
316	do ix=0,numxgrid-1
317
318	! WET DEPOSITION
319	if ((WETDEP).and.(ldirect.gt.0)) then
320	do l=1,nclassunc
321	auxgrid(l)=wetgridunc(ix,jy,ks,kp,l,nage)
322	end do
323	call mean(auxgrid,wetgrid(ix,jy), &
324	wetgridsigma(ix,jy),nclassunc)
325	! Multiply by number of classes to get total concentration
326	wetgrid(ix,jy)=wetgrid(ix,jy) &
327	*nclassunc
328	wetgridtotal=wetgridtotal+wetgrid(ix,jy)
329	! Calculate standard deviation of the mean
330	wetgridsigma(ix,jy)= &
331	wetgridsigma(ix,jy)* &
332	sqrt(real(nclassunc))
333	wetgridsigmatotal=wetgridsigmatotal+ &
334	wetgridsigma(ix,jy)
335	endif
336
337	! DRY DEPOSITION
338	if ((DRYDEP).and.(ldirect.gt.0)) then
339	do l=1,nclassunc
340	auxgrid(l)=drygridunc(ix,jy,ks,kp,l,nage)
341	end do
342	call mean(auxgrid,drygrid(ix,jy), &
343	drygridsigma(ix,jy),nclassunc)
344	! Multiply by number of classes to get total concentration
345	drygrid(ix,jy)=drygrid(ix,jy)* &
346	nclassunc
347	drygridtotal=drygridtotal+drygrid(ix,jy)
348	! Calculate standard deviation of the mean
349	drygridsigma(ix,jy)= &
350	drygridsigma(ix,jy)* &
351	sqrt(real(nclassunc))
352	125 drygridsigmatotal=drygridsigmatotal+ &
353	drygridsigma(ix,jy)
354	endif
355	! CONCENTRATION OR MIXING RATIO
356	do kz=1,numzgrid
357	do l=1,nclassunc
358	auxgrid(l)=gridunc(ix,jy,kz,ks,kp,l,nage)
359	end do
360	call mean(auxgrid,grid(ix,jy,kz), &
361	gridsigma(ix,jy,kz),nclassunc)
362	! Multiply by number of classes to get total concentration
363	grid(ix,jy,kz)= &
364	grid(ix,jy,kz)*nclassunc
365	gridtotal=gridtotal+grid(ix,jy,kz)
366	! Calculate standard deviation of the mean
367	gridsigma(ix,jy,kz)= &
368	gridsigma(ix,jy,kz)* &
369	sqrt(real(nclassunc))
370	gridsigmatotal=gridsigmatotal+ &
371	gridsigma(ix,jy,kz)
372	end do
373	end do
374	end do
375
376	!*******************************************************************
377	! Generate output: may be in concentration (ng/m3) or in mixing
378	! ratio (ppt) or both
379	! Output the position and the values alternated multiplied by
380	! 1 or -1, first line is number of values, number of positions
381	! For backward simulations, the unit is seconds, stored in grid_time
382	!*******************************************************************
383
384	if (iouttype.eq.2) then ! netcdf output
385	if (option_verbose.ge.1) then
386	write(,) 'concoutput_irreg: Calling write_ncconc for main outgrid'
387	endif
388	call write_ncconc(itime,outnum,ks,kp,nage,tot_mu(ks,kp),0) ! 0= nest level
389	else ! binary or ascii output
390
391	! Concentration output
392	!*********************
393	if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5)) then
394
395	! Wet deposition
396	sp_count_i=0
397	sp_count_r=0
398	sp_fact=-1.
399	sp_zer=.true.
400	if ((ldirect.eq.1).and.(WETDEP)) then
401	do jy=0,numygrid-1
402	do ix=0,numxgrid-1
403	! concentraion greater zero
404	if (wetgrid(ix,jy).gt.smallnum) then
405	if (sp_zer.eqv..true.) then ! first non zero value
406	sp_count_i=sp_count_i+1
407	sparse_dump_i(sp_count_i)=ix+jy*numxgrid
408	sp_zer=.false.
409	sp_fact=sp_fact*(-1.)
410	endif
411	sp_count_r=sp_count_r+1
412	sparse_dump_r(sp_count_r)= &
413	sp_fact1.e12wetgrid(ix,jy)/area(ix,jy)
414	! sparse_dump_u(sp_count_r)=
415	!+ 1.e12*wetgridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
416	else ! concentration is zero
417	sp_zer=.true.
418	endif
419	end do
420	end do
421	else
422	sp_count_i=0
423	sp_count_r=0
424	endif
425	if (iouttype.eq.0) then
426	write(unitoutgrid) sp_count_i
427	write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
428	write(unitoutgrid) sp_count_r
429	write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
430	endif
431	if (iouttype.eq.1) then
432	write(unitoutgrid,*) sp_count_i
433	write(unitoutgrid,*) (sparse_dump_i(i),i=1,sp_count_i)
434	write(unitoutgrid,*) sp_count_r
435	write(unitoutgrid,*) (sparse_dump_r(i),i=1,sp_count_r)
436	endif
437	! write(unitoutgrid) sp_count_u
438	! write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
439
440	! Dry deposition
441	sp_count_i=0
442	sp_count_r=0
443	sp_fact=-1.
444	sp_zer=.true.
445	if ((ldirect.eq.1).and.(DRYDEP)) then
446	do jy=0,numygrid-1
447	do ix=0,numxgrid-1
448	if (drygrid(ix,jy).gt.smallnum) then
449	if (sp_zer.eqv..true.) then ! first non zero value
450	sp_count_i=sp_count_i+1
451	sparse_dump_i(sp_count_i)=ix+jy*numxgrid
452	sp_zer=.false.
453	sp_fact=sp_fact*(-1.)
454	endif
455	sp_count_r=sp_count_r+1
456	sparse_dump_r(sp_count_r)= &
457	sp_fact* &
458	1.e12*drygrid(ix,jy)/area(ix,jy)
459	! sparse_dump_u(sp_count_r)=
460	!+ 1.e12*drygridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
461	else ! concentration is zero
462	sp_zer=.true.
463	endif
464	end do
465	end do
466	else
467	sp_count_i=0
468	sp_count_r=0
469	endif
470	if (iouttype.eq.0) then
471	write(unitoutgrid) sp_count_i
472	write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
473	write(unitoutgrid) sp_count_r
474	write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
475	endif
476	if (iouttype.eq.1) then
477	write(unitoutgrid,*) sp_count_i
478	write(unitoutgrid,*) (sparse_dump_i(i),i=1,sp_count_i)
479	write(unitoutgrid,*) sp_count_r
480	write(unitoutgrid,*) (sparse_dump_r(i),i=1,sp_count_r)
481	endif
482	! write(,) sp_count_u
483	! write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
484
485	! Concentrations
486	sp_count_i=0
487	sp_count_r=0
488	sp_fact=-1.
489	sp_zer=.true.
490	do kz=1,numzgrid
491	do jy=0,numygrid-1
492	do ix=0,numxgrid-1
493	if (grid(ix,jy,kz).gt.smallnum) then
494	if (sp_zer.eqv..true.) then ! first non zero value
495	sp_count_i=sp_count_i+1
496	sparse_dump_i(sp_count_i)= &
497	ix+jynumxgrid+kznumxgrid*numygrid
498	sp_zer=.false.
499	sp_fact=sp_fact*(-1.)
500	endif
501	sp_count_r=sp_count_r+1
502	sparse_dump_r(sp_count_r)= &
503	sp_fact* &
504	grid(ix,jy,kz)* &
505	factor3d(ix,jy,kz)/tot_mu(ks,kp)
506	! if ((factor(ix,jy,kz)/tot_mu(ks,kp)).eq.0)
507	! + write (,) factor(ix,jy,kz),tot_mu(ks,kp),ks,kp
508	! sparse_dump_u(sp_count_r)=
509	!+ ,gridsigma(ix,jy,kz,ks,kp,nage)*
510	!+ factor(ix,jy,kz)/tot_mu(ks,kp)
511	else ! concentration is zero
512	sp_zer=.true.
513	endif
514	end do
515	end do
516	end do
517	if (iouttype.eq.0) then
518	write(unitoutgrid) sp_count_i
519	write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
520	write(unitoutgrid) sp_count_r
521	write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
522	endif
523	if (iouttype.eq.1) then
524	write(unitoutgrid,*) sp_count_i
525	write(unitoutgrid,*) (sparse_dump_i(i),i=1,sp_count_i)
526	write(unitoutgrid,*) sp_count_r
527	write(unitoutgrid,*) (sparse_dump_r(i),i=1,sp_count_r)
528	endif
529	! write(unitoutgrid) sp_count_u
530	! write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
531
532
533
534	endif ! concentration output
535	! Mixing ratio output
536	!********************
537
538	if ((iout.eq.2).or.(iout.eq.3)) then ! mixing ratio
539
540	! Wet deposition
541	sp_count_i=0
542	sp_count_r=0
543	sp_fact=-1.
544	sp_zer=.true.
545	if ((ldirect.eq.1).and.(WETDEP)) then
546	do jy=0,numygrid-1
547	do ix=0,numxgrid-1
548	if (wetgrid(ix,jy).gt.smallnum) then
549	if (sp_zer.eqv..true.) then ! first non zero value
550	sp_count_i=sp_count_i+1
551	sparse_dump_i(sp_count_i)= &
552	ix+jy*numxgrid
553	sp_zer=.false.
554	sp_fact=sp_fact*(-1.)
555	endif
556	sp_count_r=sp_count_r+1
557	sparse_dump_r(sp_count_r)= &
558	sp_fact* &
559	1.e12*wetgrid(ix,jy)/area(ix,jy)
560	! sparse_dump_u(sp_count_r)=
561	! + ,1.e12*wetgridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
562	else ! concentration is zero
563	sp_zer=.true.
564	endif
565	end do
566	end do
567	else
568	sp_count_i=0
569	sp_count_r=0
570	endif
571	if (iouttype.eq.0) then
572	write(unitoutgridppt) sp_count_i
573	write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
574	write(unitoutgridppt) sp_count_r
575	write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
576	endif
577	if (iouttype.eq.1) then
578	write(unitoutgridppt,*) sp_count_i
579	write(unitoutgridppt,*) (sparse_dump_i(i),i=1,sp_count_i)
580	write(unitoutgridppt,*) sp_count_r
581	write(unitoutgridppt,*) (sparse_dump_r(i),i=1,sp_count_r)
582	endif
583	! write(unitoutgridppt) sp_count_u
584	! write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
585
586	! Dry deposition
587	sp_count_i=0
588	sp_count_r=0
589	sp_fact=-1.
590	sp_zer=.true.
591	if ((ldirect.eq.1).and.(DRYDEP)) then
592	do jy=0,numygrid-1
593	do ix=0,numxgrid-1
594	if (drygrid(ix,jy).gt.smallnum) then
595	if (sp_zer.eqv..true.) then ! first non zero value
596	sp_count_i=sp_count_i+1
597	sparse_dump_i(sp_count_i)= &
598	ix+jy*numxgrid
599	sp_zer=.false.
600	sp_fact=sp_fact*(-1)
601	endif
602	sp_count_r=sp_count_r+1
603	sparse_dump_r(sp_count_r)= &
604	sp_fact* &
605	1.e12*drygrid(ix,jy)/area(ix,jy)
606	! sparse_dump_u(sp_count_r)=
607	! + ,1.e12*drygridsigma(ix,jy,ks,kp,nage)/area(ix,jy)
608	else ! concentration is zero
609	sp_zer=.true.
610	endif
611	end do
612	end do
613	else
614	sp_count_i=0
615	sp_count_r=0
616	endif
617	if (iouttype.eq.0) then
618	write(unitoutgridppt) sp_count_i
619	write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
620	write(unitoutgridppt) sp_count_r
621	write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
622	endif
623	if (iouttype.eq.1) then
624	write(unitoutgridppt,*) sp_count_i
625	write(unitoutgridppt,*) (sparse_dump_i(i),i=1,sp_count_i)
626	write(unitoutgridppt,*) sp_count_r
627	write(unitoutgridppt,*) (sparse_dump_r(i),i=1,sp_count_r)
628	endif
629	! write(unitoutgridppt) sp_count_u
630	! write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
631
632
633	! Mixing ratios
634	sp_count_i=0
635	sp_count_r=0
636	sp_fact=-1.
637	sp_zer=.true.
638	do kz=1,numzgrid
639	do jy=0,numygrid-1
640	do ix=0,numxgrid-1
641	if (grid(ix,jy,kz).gt.smallnum) then
642	if (sp_zer.eqv..true.) then ! first non zero value
643	sp_count_i=sp_count_i+1
644	sparse_dump_i(sp_count_i)= &
645	ix+jynumxgrid+kznumxgrid*numygrid
646	sp_zer=.false.
647	sp_fact=sp_fact*(-1.)
648	endif
649	sp_count_r=sp_count_r+1
650	sparse_dump_r(sp_count_r)= &
651	sp_fact* &
652	1.e12*grid(ix,jy,kz) &
653	/volume(ix,jy,kz)/outnum* &
654	weightair/weightmolar(ks)/densityoutgrid(ix,jy,kz)
655	! sparse_dump_u(sp_count_r)=
656	!+ ,1.e12*gridsigma(ix,jy,kz,ks,kp,nage)/volume(ix,jy,kz)/
657	!+ outnum*weightair/weightmolar(ks)/
658	!+ densityoutgrid(ix,jy,kz)
659	else ! concentration is zero
660	sp_zer=.true.
661	endif
662	end do
663	end do
664	end do
665	if (iouttype.eq.0) then
666	write(unitoutgridppt) sp_count_i
667	write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
668	write(unitoutgridppt) sp_count_r
669	write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
670	endif
671	if (iouttype.eq.1) then
672	write(unitoutgridppt,*) sp_count_i
673	write(unitoutgridppt,*) (sparse_dump_i(i),i=1,sp_count_i)
674	write(unitoutgridppt,*) sp_count_r
675	write(unitoutgridppt,*) (sparse_dump_r(i),i=1,sp_count_r)
676	endif
677	! write(unitoutgridppt) sp_count_u
678	! write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
679
680	endif ! output for ppt
681
682	endif ! iouttype.eq.2
683
684	end do
685	end do
686
687	if((iouttype.eq.0).or.(iouttype.eq.1)) then ! binary or ascii output
688	close(unitoutgridppt)
689	close(unitoutgrid)
690	endif
691
692	end do
693
694	if (gridtotal.gt.0.) gridtotalunc=gridsigmatotal/gridtotal
695	if (wetgridtotal.gt.0.) wetgridtotalunc=wetgridsigmatotal/ &
696	wetgridtotal
697	if (drygridtotal.gt.0.) drygridtotalunc=drygridsigmatotal/ &
698	drygridtotal
699
700	! Dump of receptor concentrations
701
702	if (numreceptor.gt.0 .and. (iout.eq.2 .or. iout.eq.3) ) then
703	write(unitoutreceptppt) itime
704	do ks=1,nspec
705	write(unitoutreceptppt) (1.e12creceptor(i,ks)/outnum &
706	weightair/weightmolar(ks)/densityoutrecept(i),i=1,numreceptor)
707	end do
708	endif
709
710	! Dump of receptor concentrations
711
712	if (numreceptor.gt.0) then
713	write(unitoutrecept) itime
714	do ks=1,nspec
715	write(unitoutrecept) (1.e12*creceptor(i,ks)/outnum, &
716	i=1,numreceptor)
717	end do
718	endif
719
720	do ks=1,nspec
721	do kp=1,maxpointspec_act
722	do i=1,numreceptor
723	creceptor(i,ks)=0.
724	end do
725	do jy=0,numygrid-1
726	do ix=0,numxgrid-1
727	do l=1,nclassunc
728	do nage=1,nageclass
729	do kz=1,numzgrid
730	gridunc(ix,jy,kz,ks,kp,l,nage)=0.
731	end do
732	end do
733	end do
734	end do
735	end do
736	end do
737	end do
738
739
740	end subroutine concoutput_reg
741

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

Download in other formats:

Original Format