source: trunk/src/concoutput_surf_nest.f90 @ 28

Last change on this file since 28 was 20, checked in by igpis, 10 years ago

move version 9.1.8 form branches to trunk. Contributions from HSO, saeck, pesei, NIK, RT, XKF, IP and others

File size: 23.6 KB
Line 
1!**********************************************************************
2! Copyright 1998,1999,2000,2001,2002,2005,2007,2008,2009,2010         *
3! Andreas Stohl, Petra Seibert, A. Frank, Gerhard Wotawa,             *
4! Caroline Forster, Sabine Eckhardt, John Burkhart, Harald Sodemann   *
5!                                                                     *
6! This file is part of FLEXPART.                                      *
7!                                                                     *
8! FLEXPART is free software: you can redistribute it and/or modify    *
9! it under the terms of the GNU General Public License as published by*
10! the Free Software Foundation, either version 3 of the License, or   *
11! (at your option) any later version.                                 *
12!                                                                     *
13! FLEXPART is distributed in the hope that it will be useful,         *
14! but WITHOUT ANY WARRANTY; without even the implied warranty of      *
15! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the       *
16! GNU General Public License for more details.                        *
17!                                                                     *
18! You should have received a copy of the GNU General Public License   *
19! along with FLEXPART.  If not, see <http://www.gnu.org/licenses/>.   *
20!**********************************************************************
21
22subroutine concoutput_surf_nest(itime,outnum)
23  !                        i     i
24  !*****************************************************************************
25  !                                                                            *
26  !     Output of the concentration grid and the receptor concentrations.      *
27  !                                                                            *
28  !     Author: A. Stohl                                                       *
29  !                                                                            *
30  !     24 May 1995                                                            *
31  !                                                                            *
32  !     13 April 1999, Major update: if output size is smaller, dump output    *
33  !                    in sparse matrix format; additional output of           *
34  !                    uncertainty                                             *
35  !                                                                            *
36  !     05 April 2000, Major update: output of age classes; output for backward*
37  !                    runs is time spent in grid cell times total mass of     *
38  !                    species.                                                *
39  !                                                                            *
40  !     17 February 2002, Appropriate dimensions for backward and forward runs *
41  !                       are now specified in file par_mod                    *
42  !                                                                            *
43  !     June 2006, write grid in sparse matrix with a single write command     *
44  !                in order to save disk space                                 *
45  !                                                                            *
46  !     2008 new sparse matrix format                                          *
47  !                                                                            *
48  !*****************************************************************************
49  !                                                                            *
50  ! Variables:                                                                 *
51  ! outnum          number of samples                                          *
52  ! ncells          number of cells with non-zero concentrations               *
53  ! sparse          .true. if in sparse matrix format, else .false.            *
54  ! tot_mu          1 for forward, initial mass mixing ration for backw. runs  *
55  !                                                                            *
56  !*****************************************************************************
57
58  use unc_mod
59  use point_mod
60  use outg_mod
61  use par_mod
62  use com_mod
63
64  implicit none
65
66  real(kind=dp) :: jul
67  integer :: itime,i,ix,jy,kz,ks,kp,l,iix,jjy,kzz,nage,jjjjmmdd,ihmmss
68  integer :: sp_count_i,sp_count_r
69  real :: sp_fact
70  real :: outnum,densityoutrecept(maxreceptor),xl,yl
71
72  !real densityoutgrid(0:numxgrid-1,0:numygrid-1,numzgrid),
73  !    +grid(0:numxgrid-1,0:numygrid-1,numzgrid,maxspec,maxpointspec_act,
74  !    +    maxageclass)
75  !real wetgrid(0:numxgrid-1,0:numygrid-1,maxspec,maxpointspec_act,
76  !    +       maxageclass)
77  !real drygrid(0:numxgrid-1,0:numygrid-1,maxspec,
78  !    +       maxpointspec_act,maxageclass)
79  !real gridsigma(0:numxgrid-1,0:numygrid-1,numzgrid,maxspec,
80  !    +       maxpointspec_act,maxageclass),
81  !    +     drygridsigma(0:numxgrid-1,0:numygrid-1,maxspec,
82  !    +     maxpointspec_act,maxageclass),
83  !    +     wetgridsigma(0:numxgrid-1,0:numygrid-1,maxspec,
84  !    +     maxpointspec_act,maxageclass)
85  !real factor(0:numxgrid-1,0:numygrid-1,numzgrid)
86  !real sparse_dump_r(numxgrid*numygrid*numzgrid)
87  !integer sparse_dump_i(numxgrid*numygrid*numzgrid)
88
89  !real sparse_dump_u(numxgrid*numygrid*numzgrid)
90  real :: auxgrid(nclassunc)
91  real :: halfheight,dz,dz1,dz2,tot_mu(maxspec,maxpointspec_act)
92  real,parameter :: smallnum = tiny(0.0) ! smallest number that can be handled
93  real,parameter :: weightair=28.97
94  logical :: sp_zer
95  character :: adate*8,atime*6
96  character(len=3) :: anspec
97
98
99  ! Determine current calendar date, needed for the file name
100  !**********************************************************
101
102  jul=bdate+real(itime,kind=dp)/86400._dp
103  call caldate(jul,jjjjmmdd,ihmmss)
104  write(adate,'(i8.8)') jjjjmmdd
105  write(atime,'(i6.6)') ihmmss
106
107
108  ! For forward simulations, output fields have dimension MAXSPEC,
109  ! for backward simulations, output fields have dimension MAXPOINT.
110  ! Thus, make loops either about nspec, or about numpoint
111  !*****************************************************************
112
113
114    if (ldirect.eq.1) then
115       do ks=1,nspec
116         do kp=1,maxpointspec_act
117           tot_mu(ks,kp)=1
118         end do
119       end do
120   else
121      do ks=1,nspec
122             do kp=1,maxpointspec_act
123               tot_mu(ks,kp)=xmass(kp,ks)
124             end do
125      end do
126    endif
127
128
129  !*******************************************************************
130  ! Compute air density: sufficiently accurate to take it
131  ! from coarse grid at some time
132  ! Determine center altitude of output layer, and interpolate density
133  ! data to that altitude
134  !*******************************************************************
135
136  do kz=1,numzgrid
137    if (kz.eq.1) then
138      halfheight=outheight(1)/2.
139    else
140      halfheight=(outheight(kz)+outheight(kz-1))/2.
141    endif
142    do kzz=2,nz
143      if ((height(kzz-1).lt.halfheight).and. &
144           (height(kzz).gt.halfheight)) goto 46
145    end do
14646   kzz=max(min(kzz,nz),2)
147    dz1=halfheight-height(kzz-1)
148    dz2=height(kzz)-halfheight
149    dz=dz1+dz2
150    do jy=0,numygridn-1
151      do ix=0,numxgridn-1
152        xl=outlon0n+real(ix)*dxoutn
153        yl=outlat0n+real(jy)*dyoutn
154        xl=(xl-xlon0)/dx
155        yl=(yl-ylat0)/dy
156        iix=max(min(nint(xl),nxmin1),0)
157        jjy=max(min(nint(yl),nymin1),0)
158        densityoutgrid(ix,jy,kz)=(rho(iix,jjy,kzz,2)*dz1+ &
159             rho(iix,jjy,kzz-1,2)*dz2)/dz
160      end do
161    end do
162  end do
163
164    do i=1,numreceptor
165      xl=xreceptor(i)
166      yl=yreceptor(i)
167      iix=max(min(nint(xl),nxmin1),0)
168      jjy=max(min(nint(yl),nymin1),0)
169      densityoutrecept(i)=rho(iix,jjy,1,2)
170    end do
171
172
173  ! Output is different for forward and backward simulations
174    do kz=1,numzgrid
175      do jy=0,numygridn-1
176        do ix=0,numxgridn-1
177          if (ldirect.eq.1) then
178            factor3d(ix,jy,kz)=1.e12/volumen(ix,jy,kz)/outnum
179          else
180            factor3d(ix,jy,kz)=real(abs(loutaver))/outnum
181          endif
182        end do
183      end do
184    end do
185
186  !*********************************************************************
187  ! Determine the standard deviation of the mean concentration or mixing
188  ! ratio (uncertainty of the output) and the dry and wet deposition
189  !*********************************************************************
190
191  do ks=1,nspec
192
193  write(anspec,'(i3.3)') ks
194  if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5)) then
195    if (ldirect.eq.1) then
196      open(unitoutgrid,file=path(2)(1:length(2))//'grid_conc_nest_' &
197           //adate// &
198           atime//'_'//anspec,form='unformatted')
199    else
200      open(unitoutgrid,file=path(2)(1:length(2))//'grid_time_nest_' &
201           //adate// &
202           atime//'_'//anspec,form='unformatted')
203    endif
204     write(unitoutgrid) itime
205   endif
206
207  if ((iout.eq.2).or.(iout.eq.3)) then      ! mixing ratio
208   open(unitoutgridppt,file=path(2)(1:length(2))//'grid_pptv_nest_' &
209        //adate// &
210        atime//'_'//anspec,form='unformatted')
211
212    write(unitoutgridppt) itime
213  endif
214
215  do kp=1,maxpointspec_act
216  do nage=1,nageclass
217
218    do jy=0,numygridn-1
219      do ix=0,numxgridn-1
220
221  ! WET DEPOSITION
222        if ((WETDEP).and.(ldirect.gt.0)) then
223            do l=1,nclassunc
224              auxgrid(l)=wetgriduncn(ix,jy,ks,kp,l,nage)
225            end do
226            call mean(auxgrid,wetgrid(ix,jy), &
227                 wetgridsigma(ix,jy),nclassunc)
228  ! Multiply by number of classes to get total concentration
229            wetgrid(ix,jy)=wetgrid(ix,jy) &
230                 *nclassunc
231  ! Calculate standard deviation of the mean
232            wetgridsigma(ix,jy)= &
233                 wetgridsigma(ix,jy)* &
234                 sqrt(real(nclassunc))
235        endif
236
237  ! DRY DEPOSITION
238        if ((DRYDEP).and.(ldirect.gt.0)) then
239            do l=1,nclassunc
240              auxgrid(l)=drygriduncn(ix,jy,ks,kp,l,nage)
241            end do
242            call mean(auxgrid,drygrid(ix,jy), &
243                 drygridsigma(ix,jy),nclassunc)
244  ! Multiply by number of classes to get total concentration
245            drygrid(ix,jy)=drygrid(ix,jy)* &
246                 nclassunc
247  ! Calculate standard deviation of the mean
248            drygridsigma(ix,jy)= &
249                 drygridsigma(ix,jy)* &
250                 sqrt(real(nclassunc))
251        endif
252
253  ! CONCENTRATION OR MIXING RATIO
254        do kz=1,numzgrid
255            do l=1,nclassunc
256              auxgrid(l)=griduncn(ix,jy,kz,ks,kp,l,nage)
257            end do
258            call mean(auxgrid,grid(ix,jy,kz), &
259                 gridsigma(ix,jy,kz),nclassunc)
260  ! Multiply by number of classes to get total concentration
261            grid(ix,jy,kz)= &
262                 grid(ix,jy,kz)*nclassunc
263  ! Calculate standard deviation of the mean
264            gridsigma(ix,jy,kz)= &
265                 gridsigma(ix,jy,kz)* &
266                 sqrt(real(nclassunc))
267        end do
268      end do
269    end do
270
271
272  !*******************************************************************
273  ! Generate output: may be in concentration (ng/m3) or in mixing
274  ! ratio (ppt) or both
275  ! Output the position and the values alternated multiplied by
276  ! 1 or -1, first line is number of values, number of positions
277  ! For backward simulations, the unit is seconds, stored in grid_time
278  !*******************************************************************
279
280  ! Concentration output
281  !*********************
282  if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5)) then
283
284  ! Wet deposition
285         sp_count_i=0
286         sp_count_r=0
287         sp_fact=-1.
288         sp_zer=.true.
289         if ((ldirect.eq.1).and.(WETDEP)) then
290         do jy=0,numygridn-1
291            do ix=0,numxgridn-1
292  !oncentraion greater zero
293              if (wetgrid(ix,jy).gt.smallnum) then
294                 if (sp_zer.eqv..true.) then ! first non zero value
295                    sp_count_i=sp_count_i+1
296                    sparse_dump_i(sp_count_i)=ix+jy*numxgridn
297                    sp_zer=.false.
298                    sp_fact=sp_fact*(-1.)
299                 endif
300                 sp_count_r=sp_count_r+1
301                 sparse_dump_r(sp_count_r)= &
302                      sp_fact*1.e12*wetgrid(ix,jy)/arean(ix,jy)
303                 sparse_dump_u(sp_count_r)= &
304                      1.e12*wetgridsigma(ix,jy)/area(ix,jy)
305              else ! concentration is zero
306                  sp_zer=.true.
307              endif
308            end do
309         end do
310         else
311            sp_count_i=0
312            sp_count_r=0
313         endif
314         write(unitoutgrid) sp_count_i
315         write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
316         write(unitoutgrid) sp_count_r
317         write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
318         write(unitoutgrid) sp_count_r
319         write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
320
321  ! Dry deposition
322         sp_count_i=0
323         sp_count_r=0
324         sp_fact=-1.
325         sp_zer=.true.
326         if ((ldirect.eq.1).and.(DRYDEP)) then
327          do jy=0,numygridn-1
328            do ix=0,numxgridn-1
329              if (drygrid(ix,jy).gt.smallnum) then
330                 if (sp_zer.eqv..true.) then ! first non zero value
331                    sp_count_i=sp_count_i+1
332                    sparse_dump_i(sp_count_i)=ix+jy*numxgridn
333                    sp_zer=.false.
334                    sp_fact=sp_fact*(-1.)
335                 endif
336                 sp_count_r=sp_count_r+1
337                 sparse_dump_r(sp_count_r)= &
338                      sp_fact* &
339                      1.e12*drygrid(ix,jy)/arean(ix,jy)
340                 sparse_dump_u(sp_count_r)= &
341                      1.e12*drygridsigma(ix,jy)/area(ix,jy)
342              else ! concentration is zero
343                  sp_zer=.true.
344              endif
345            end do
346          end do
347         else
348            sp_count_i=0
349            sp_count_r=0
350         endif
351         write(unitoutgrid) sp_count_i
352         write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
353         write(unitoutgrid) sp_count_r
354         write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
355         write(unitoutgrid) sp_count_r
356         write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
357
358
359
360  ! Concentrations
361
362  ! if surf_only write only 1st layer
363
364         if(surf_only.eq.1) then
365         sp_count_i=0
366         sp_count_r=0
367         sp_fact=-1.
368         sp_zer=.true.
369          do kz=1,1
370            do jy=0,numygridn-1
371              do ix=0,numxgridn-1
372                if (grid(ix,jy,kz).gt.smallnum) then
373                  if (sp_zer.eqv..true.) then ! first non zero value
374                    sp_count_i=sp_count_i+1
375                    sparse_dump_i(sp_count_i)= &
376                         ix+jy*numxgridn+kz*numxgridn*numygridn
377                    sp_zer=.false.
378                    sp_fact=sp_fact*(-1.)
379                   endif
380                   sp_count_r=sp_count_r+1
381                   sparse_dump_r(sp_count_r)= &
382                        sp_fact* &
383                        grid(ix,jy,kz)* &
384                        factor3d(ix,jy,kz)/tot_mu(ks,kp)
385  !                 if ((factor(ix,jy,kz)/tot_mu(ks,kp)).eq.0)
386  !    +              write (*,*) factor(ix,jy,kz),tot_mu(ks,kp),ks,kp
387                   sparse_dump_u(sp_count_r)= &
388                        gridsigma(ix,jy,kz)* &
389                        factor3d(ix,jy,kz)/tot_mu(ks,kp)
390              else ! concentration is zero
391                  sp_zer=.true.
392              endif
393              end do
394            end do
395          end do
396         write(unitoutgrid) sp_count_i
397         write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
398         write(unitoutgrid) sp_count_r
399         write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
400         write(unitoutgrid) sp_count_r
401         write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
402         else
403
404  ! write full vertical resolution
405
406         sp_count_i=0
407         sp_count_r=0
408         sp_fact=-1.
409         sp_zer=.true.
410          do kz=1,numzgrid
411            do jy=0,numygridn-1
412              do ix=0,numxgridn-1
413                if (grid(ix,jy,kz).gt.smallnum) then
414                  if (sp_zer.eqv..true.) then ! first non zero value
415                    sp_count_i=sp_count_i+1
416                    sparse_dump_i(sp_count_i)= &
417                         ix+jy*numxgridn+kz*numxgridn*numygridn
418                    sp_zer=.false.
419                    sp_fact=sp_fact*(-1.)
420                   endif
421                   sp_count_r=sp_count_r+1
422                   sparse_dump_r(sp_count_r)= &
423                        sp_fact* &
424                        grid(ix,jy,kz)* &
425                        factor3d(ix,jy,kz)/tot_mu(ks,kp)
426  !                 if ((factor(ix,jy,kz)/tot_mu(ks,kp)).eq.0)
427  !    +              write (*,*) factor(ix,jy,kz),tot_mu(ks,kp),ks,kp
428                   sparse_dump_u(sp_count_r)= &
429                        gridsigma(ix,jy,kz)* &
430                        factor3d(ix,jy,kz)/tot_mu(ks,kp)
431              else ! concentration is zero
432                  sp_zer=.true.
433              endif
434              end do
435            end do
436          end do
437         write(unitoutgrid) sp_count_i
438         write(unitoutgrid) (sparse_dump_i(i),i=1,sp_count_i)
439         write(unitoutgrid) sp_count_r
440         write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
441         write(unitoutgrid) sp_count_r
442         write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
443         endif ! surf_only
444
445
446    endif !  concentration output
447
448  ! Mixing ratio output
449  !********************
450
451  if ((iout.eq.2).or.(iout.eq.3)) then      ! mixing ratio
452
453  ! Wet deposition
454         sp_count_i=0
455         sp_count_r=0
456         sp_fact=-1.
457         sp_zer=.true.
458         if ((ldirect.eq.1).and.(WETDEP)) then
459          do jy=0,numygridn-1
460            do ix=0,numxgridn-1
461                if (wetgrid(ix,jy).gt.smallnum) then
462                  if (sp_zer.eqv..true.) then ! first non zero value
463                    sp_count_i=sp_count_i+1
464                    sparse_dump_i(sp_count_i)= &
465                         ix+jy*numxgridn
466                    sp_zer=.false.
467                    sp_fact=sp_fact*(-1.)
468                 endif
469                 sp_count_r=sp_count_r+1
470                 sparse_dump_r(sp_count_r)= &
471                      sp_fact* &
472                      1.e12*wetgrid(ix,jy)/arean(ix,jy)
473                 sparse_dump_u(sp_count_r)= &
474                      1.e12*wetgridsigma(ix,jy)/area(ix,jy)
475              else ! concentration is zero
476                  sp_zer=.true.
477              endif
478            end do
479          end do
480         else
481           sp_count_i=0
482           sp_count_r=0
483         endif
484         write(unitoutgridppt) sp_count_i
485         write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
486         write(unitoutgridppt) sp_count_r
487         write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
488         write(unitoutgridppt) sp_count_r
489         write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
490
491
492  ! Dry deposition
493         sp_count_i=0
494         sp_count_r=0
495         sp_fact=-1.
496         sp_zer=.true.
497         if ((ldirect.eq.1).and.(DRYDEP)) then
498          do jy=0,numygridn-1
499            do ix=0,numxgridn-1
500                if (drygrid(ix,jy).gt.smallnum) then
501                  if (sp_zer.eqv..true.) then ! first non zero value
502                    sp_count_i=sp_count_i+1
503                    sparse_dump_i(sp_count_i)= &
504                         ix+jy*numxgridn
505                    sp_zer=.false.
506                    sp_fact=sp_fact*(-1)
507                 endif
508                 sp_count_r=sp_count_r+1
509                 sparse_dump_r(sp_count_r)= &
510                      sp_fact* &
511                      1.e12*drygrid(ix,jy)/arean(ix,jy)
512                 sparse_dump_u(sp_count_r)= &
513                      1.e12*drygridsigma(ix,jy)/area(ix,jy)
514              else ! concentration is zero
515                  sp_zer=.true.
516              endif
517            end do
518          end do
519         else
520           sp_count_i=0
521           sp_count_r=0
522         endif
523         write(unitoutgridppt) sp_count_i
524         write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
525         write(unitoutgridppt) sp_count_r
526         write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
527         write(unitoutgridppt) sp_count_r
528         write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
529
530
531  ! Mixing ratios
532
533    ! if surf_only write only 1st layer
534
535         if(surf_only.eq.1) then
536         sp_count_i=0
537         sp_count_r=0
538         sp_fact=-1.
539         sp_zer=.true.
540          do kz=1,1
541            do jy=0,numygridn-1
542              do ix=0,numxgridn-1
543                if (grid(ix,jy,kz).gt.smallnum) then
544                  if (sp_zer.eqv..true.) then ! first non zero value
545                    sp_count_i=sp_count_i+1
546                    sparse_dump_i(sp_count_i)= &
547                         ix+jy*numxgridn+kz*numxgridn*numygridn
548                    sp_zer=.false.
549                    sp_fact=sp_fact*(-1.)
550                 endif
551                 sp_count_r=sp_count_r+1
552                 sparse_dump_r(sp_count_r)= &
553                      sp_fact* &
554                      1.e12*grid(ix,jy,kz) &
555                      /volumen(ix,jy,kz)/outnum* &
556                      weightair/weightmolar(ks)/densityoutgrid(ix,jy,kz)
557                 sparse_dump_u(sp_count_r)= &
558                      1.e12*gridsigma(ix,jy,kz)/volumen(ix,jy,kz)/ &
559                      outnum*weightair/weightmolar(ks)/ &
560                      densityoutgrid(ix,jy,kz)
561              else ! concentration is zero
562                  sp_zer=.true.
563              endif
564              end do
565            end do
566          end do
567         write(unitoutgridppt) sp_count_i
568         write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
569         write(unitoutgridppt) sp_count_r
570         write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
571         write(unitoutgridppt) sp_count_r
572         write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
573         else
574
575  ! write full vertical resolution
576
577         sp_count_i=0
578         sp_count_r=0
579         sp_fact=-1.
580         sp_zer=.true.
581          do kz=1,numzgrid
582            do jy=0,numygridn-1
583              do ix=0,numxgridn-1
584                if (grid(ix,jy,kz).gt.smallnum) then
585                  if (sp_zer.eqv..true.) then ! first non zero value
586                    sp_count_i=sp_count_i+1
587                    sparse_dump_i(sp_count_i)= &
588                         ix+jy*numxgridn+kz*numxgridn*numygridn
589                    sp_zer=.false.
590                    sp_fact=sp_fact*(-1.)
591                 endif
592                 sp_count_r=sp_count_r+1
593                 sparse_dump_r(sp_count_r)= &
594                      sp_fact* &
595                      1.e12*grid(ix,jy,kz) &
596                      /volumen(ix,jy,kz)/outnum* &
597                      weightair/weightmolar(ks)/densityoutgrid(ix,jy,kz)
598                 sparse_dump_u(sp_count_r)= &
599                      1.e12*gridsigma(ix,jy,kz)/volumen(ix,jy,kz)/ &
600                      outnum*weightair/weightmolar(ks)/ &
601                      densityoutgrid(ix,jy,kz)
602              else ! concentration is zero
603                  sp_zer=.true.
604              endif
605              end do
606            end do
607          end do
608         write(unitoutgridppt) sp_count_i
609         write(unitoutgridppt) (sparse_dump_i(i),i=1,sp_count_i)
610         write(unitoutgridppt) sp_count_r
611         write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
612         write(unitoutgridppt) sp_count_r
613         write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
614         endif ! surf_only
615
616      endif ! output for ppt
617
618  end do
619  end do
620
621    close(unitoutgridppt)
622    close(unitoutgrid)
623
624  end do
625
626
627
628  ! Reinitialization of grid
629  !*************************
630
631  do ks=1,nspec
632  do kp=1,maxpointspec_act
633    do i=1,numreceptor
634      creceptor(i,ks)=0.
635    end do
636    do jy=0,numygridn-1
637      do ix=0,numxgridn-1
638        do l=1,nclassunc
639          do nage=1,nageclass
640            do kz=1,numzgrid
641              griduncn(ix,jy,kz,ks,kp,l,nage)=0.
642            end do
643          end do
644        end do
645      end do
646    end do
647  end do
648  end do
649
650
651end subroutine concoutput_surf_nest
652
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG