source: branches/jerome/src_flexwrf_v3.1/concoutput_nest_reg.f90 @ 16

Last change on this file since 16 was 16, checked in by jebri, 10 years ago

sources for flexwrf v3.1

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