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

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

Reworked the domain-filling option (MPI). Fixed a slow loop which had errors in loop counter (MPI)

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