source: trunk/src/concoutput.f90

Last change on this file was 20, checked in by igpis, 11 years ago

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

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