source: flexpart.git/src/concoutput_surf_nest.f90 @ 16b61a5

10.4.1_peseiGFS_025bugfixes+enhancementsdevrelease-10release-10.4.1scaling-bugunivie
Last change on this file since 16b61a5 was 6a678e3, checked in by Espen Sollum ATMOS <eso@…>, 8 years ago

Added option to use double precision for calculating the deposition fields

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