source: branches/jerome/src_flexwrf_v3.1/write_ncconc.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: 12.9 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!* Adam Dingwell                                                      *
8!*                                                                    *
9!* This file is part of FLEXPART WRF                                  *
10!                                                                     *
11! FLEXPART is free software: you can redistribute it and/or modify    *
12! it under the terms of the GNU General Public License as published by*
13! the Free Software Foundation, either version 3 of the License, or   *
14! (at your option) any later version.                                 *
15!                                                                     *
16! FLEXPART is distributed in the hope that it will be useful,         *
17! but WITHOUT ANY WARRANTY; without even the implied warranty of      *
18! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the       *
19! GNU General Public License for more details.                        *
20!                                                                     *
21! You should have received a copy of the GNU General Public License   *
22! along with FLEXPART.  If not, see <http://www.gnu.org/licenses/>.   *
23!**********************************************************************
24
25subroutine write_ncconc(itime,outnum,ks,kp,nage,tot_mu_scalar,nesting_level)
26 
27  !*****************************************************************************
28  !                                                                            *
29  !  This routine writes concentration, mixing ratio and deposition fields     *
30  !  to a netcdf file defined by flex_ncheader.                                *
31  !                                                                            *
32  !  flex_ncheader is called from within write_ncconc when it's time for a new *
33  !  output file.                                                              *
34  !                                                                            *
35  !  write_ncconc should be called by concoutput_irreg and concoutput_reg      *
36  !  it is separate from the binary and ascii output routines to avoid mixing  *
37  !  of sparse and full grid approaches.  Netcdf will output the full grid.    *
38  !                                                                            *
39  !      Author: A. Dingwell                                                   *
40  !                                                                            *
41  !      29 May 2013                                                           *
42  !                                                                            *
43  ! Modifications:                                                             *
44  ! June 5 2013: J. Brioude: compression using deflate level, optimization of  *
45  !  the writing procedure. bug fixes for backtrajectory mode                  *
46  !*****************************************************************************
47 
48  use point_mod
49  use outg_mod
50  use par_mod
51  use com_mod
52
53  implicit none
54
55  include 'netcdf.inc'
56
57  real    :: outnum         ! Number of samples for each concentration calculation
58  integer :: itime          ! Current simulation time [s]
59  integer :: ks,kp,nage     ! species, maxpointspec_act and nageclass indices resp.
60  real    :: tot_mu_scalar  ! total mass per source and species (backward)
61                            ! or unity (forward).  Should probably be sent as
62                            ! tot_mu(ks,kp) from concoutput*.f90
63  integer :: nesting_level  ! 0 for main (mother) grid, 1 for nest (child)
64
65  real(kind=dp) :: jul          ! Julian date
66  integer   :: jjjjmmdd,ihmmss  ! date & time as integer
67  character :: adate*8,atime*6  ! date and time strings, used for filename
68
69  integer :: ncid           ! Pointer to netcdf file, depends on nesting level
70  integer :: grid_nx,grid_ny! outgrid dimensions, depend on the nesting level
71  integer :: ncret          ! Netcdf:  return code
72  integer :: ix,jy,kz       ! iterators
73  character :: datestr*15   ! For the Times variable
74  integer :: deflate_level=5 ! compression level
75   
76  if (option_verbose.ge.1) then
77    write(*,*)'write_ncconc: writing netcdf output for: domain,kp,nage =',&
78      nesting_level+1,kp,nage
79  endif
80
81  ! Determine which nest/outfile we are writing to
82  !***********************************************
83  if (nesting_level.eq.0) then
84    ncid    = ncout
85    grid_nx = numxgrid
86    grid_ny = numygrid
87  elseif (nesting_level.eq.1) then
88    ncid    = ncoutn
89    grid_nx = numxgridn
90    grid_ny = numygridn
91  else
92    write(*,*) '***write_ncconc error: nesting level  must be 0 or 1'
93    ! Note for future development: If additional output nests are to be
94    ! supported for netcdf output, modification must be made here as well as in
95    ! the respective nesting_level if-block in write_ncheader
96  endif
97  ! Update/Initialize record index
98  !*******************************
99   if ((ks.eq.1).and.(kp.eq.1).and.(nage.eq.1)) then
100!   print*,'ncirec',ncirec,ncnumrec
101  if (nesting_level.eq.0) then  ! Only update for first domain
102    if (itime.eq.loutstep) then  ! first output
103      ncirec = 1  ! initialize record index
104    elseif (ncirec.eq.ncnumrec) then  ! new file
105!      print*,'file is closing'
106      ncirec = 1  ! reset record index
107      ncret=nf_close(ncid)      ! close the old file
108      call check_ncerror(ncret)
109!      print*,'file is closed'
110    else
111      ncirec=ncirec+1 ! move on to next record
112    endif
113  endif
114!   print*,'ncirec',ncirec,ncnumrec
115  endif
116
117  ! Check if it's time to create a new netcdf file
118  !***********************************************
119  if (ncirec.eq.1) then         ! First output in current file?
120!   write(*,*) 'itime=',itime
121   if ((ks.eq.1).and.(kp.eq.1).and.(nage.eq.1)) then
122!    if (itime.ne.loutstep) then ! Not the first output file?
123!      ncret=nf_close(ncid)      ! close the old file
124!      call check_ncerror(ncret)
125!     print*,'file is closed'
126!    endif
127!   call write_ncheader(itime,nesting_level)  ! Create new file
128    if (option_verbose.ge.1) &
129      write(*,*)'write_ncconc: calling write_ncinfo'
130    call write_ncinfo(itime,nesting_level)  ! Create new file
131    ! Reassign file handle to the newly created file:
132     endif
133    if (nesting_level.eq.0) ncid=ncout
134    if (nesting_level.eq.1) ncid=ncoutn
135  endif
136
137  if (option_verbose.ge.10) &
138    write(*,*) 'ncid,nccovid=',ncid,nccovid
139
140
141  ! Create output for the current record index
142  !*******************************************
143  jul=bdate+real(itime,kind=dp)/86400._dp
144  call caldate(jul,jjjjmmdd,ihmmss)
145  write(adate,'(i8.8)') jjjjmmdd
146  write(atime,'(i6.6)') ihmmss
147 
148   if ((ks.eq.1).and.(kp.eq.1).and.(nage.eq.1)) then
149
150  if (option_verbose.ge.10) write(*,*)'write_ncconc: record index',ncirec
151  write(datestr,'(I8.8,A1,I6.6)') jjjjmmdd,'_',ihmmss
152  ncret = nf_put_vara_text(ncid,ncrecvid,(/1,ncirec/),(/15,1/),datestr)
153  call check_ncerror(ncret)
154   endif
155
156  if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5)) then ! concentration
157    if (option_verbose.ge.1)  &
158      write(*,*)'write_ncconc: concentration output',kp,nage,ncirec,nccovid,ncid
159
160
161    do kz=1,numzgrid
162    do jy=0,grid_ny-1
163    do ix=0,grid_nx-1
164      grid2(ix,jy,kz,kp)= grid(ix,jy,kz)*factor3d(ix,jy,kz)/tot_mu_scalar
165    enddo ! ix=1,grid_nx-1
166    enddo ! jy=1,grid_ny-1
167    enddo ! kz=1,numzgrid
168
169    if (kp.eq.maxpointspec_act) then
170    if (ldirect.eq.-1) then
171      ncret = nf_put_vara_real(ncid,nccovid, &
172        (/1,1,1,1,nage,ncirec/),(/grid_nx,grid_ny,numzgrid,kp,1,1/), &
173        grid2(0:grid_nx-1,0:grid_ny-1,1:numzgrid,1:kp))
174      call check_ncerror(ncret)
175    else
176    if (kp.gt.1) then
177      ncret = nf_put_vara_real(ncid,nccovid, &
178        (/1,1,1,1,ks,nage,ncirec/),(/grid_nx,grid_ny,numzgrid,kp,1,1,1/), &
179        grid2(0:grid_nx-1,0:grid_ny-1,1:numzgrid,1:kp))
180    else
181      ncret = nf_put_vara_real(ncid,nccovid, &
182        (/1,1,1,ks,nage,ncirec/),(/grid_nx,grid_ny,numzgrid,1,1,1/), &
183        grid2(0:grid_nx-1,0:grid_ny-1,1:numzgrid,1))
184    endif
185      call check_ncerror(ncret)
186    endif
187  endif
188!    do kz=1,numzgrid
189!    do jy=0,grid_ny-1
190!    do ix=0,grid_nx-1
191!      ncret = nf_put_vara_real(ncid,nccovid, &
192!        (/ix+1,jy+1,kz,kp,nage,ncirec/),(/1,1,1,1,1,1/), &
193!        grid(ix,jy,kz)*factor3d(ix,jy,kz)/tot_mu_scalar)
194!      call check_ncerror(ncret)
195!    enddo ! ix=1,grid_nx-1
196!    enddo ! jy=1,grid_ny-1
197!    enddo ! kz=1,numzgrid
198  endif ! concentraion
199
200  if ((iout.eq.2).or.(iout.eq.3)) then  ! mixing ratio
201    if (option_verbose.ge.1)write(*,*)'write_ncconc: mixing ratio output'
202    do kz=1,numzgrid
203    do jy=0,grid_ny-1
204    do ix=0,grid_nx-1
205      grid3(ix,jy,kz,kp)= 1.e12*grid(ix,jy,kz)/volume(ix,jy,kz)/outnum*  &
206        weightair/weightmolar(ks)/densityoutgrid(ix,jy,kz)
207    enddo ! ix=1,grid_nx-1
208    enddo ! jy=1,grid_ny-1
209    enddo ! kz=1,numzgrid
210    if (kp.eq.maxpointspec_act) then
211    if (ldirect.eq.-1) then
212      ncret = nf_put_vara_real(ncid,ncravid, &
213        (/1,1,1,1,nage,ncirec/),(/grid_nx,grid_ny,numzgrid,kp,1,1/), &
214        grid3(0:grid_nx-1,0:grid_ny-1,1:numzgrid,1:kp))
215      call check_ncerror(ncret)
216  else
217    if (kp.gt.1) then
218      ncret = nf_put_vara_real(ncid,ncravid, &
219        (/1,1,1,1,ks,nage,ncirec/),(/grid_nx,grid_ny,numzgrid,kp,1,1,1/), &
220        grid3(0:grid_nx-1,0:grid_ny-1,1:numzgrid,1:kp))
221    else
222      ncret = nf_put_vara_real(ncid,ncravid, &
223        (/1,1,1,ks,nage,ncirec/),(/grid_nx,grid_ny,numzgrid,1,1,1/), &
224        grid3(0:grid_nx-1,0:grid_ny-1,1:numzgrid,1))
225    endif
226      call check_ncerror(ncret)
227  endif
228    endif
229
230!    do kz=1,numzgrid
231!    do jy=0,grid_ny-1
232!    do ix=0,grid_nx-1
233!      ncret = nf_put_vara_real(ncid,ncravid, &
234!        (/ix+1,jy+1,kz,kp,nage,ncirec/),(/1,1,1,1,1,1/), &
235!        1.e12*grid(ix,jy,kz)/volume(ix,jy,kz)/outnum*  &
236!        weightair/weightmolar(ks)/densityoutgrid(ix,jy,kz))
237!      call check_ncerror(ncret)
238!    enddo ! ix=1,grid_nx-1
239!    enddo ! jy=1,numygrid-1
240!    enddo ! kz=1,numzgrid
241  endif ! mixing ratio
242
243  if ((ldirect.eq.1).and.(WETDEP)) then ! WETDEP
244    if (option_verbose.ge.1)write(*,*)'write_ncconc: wet deposition output'
245    do jy=0,grid_ny-1
246    do ix=0,grid_nx-1
247    if (nesting_level.eq.0)  wetgrid2(ix,jy,kp)=1.e12*wetgrid(ix,jy)/area(ix,jy)
248    if (nesting_level.eq.1)  wetgrid2(ix,jy,kp)=1.e12*wetgrid(ix,jy)/arean(ix,jy)
249    enddo ! ix=1,grid_nx-1
250    enddo ! jy=1,grid_ny-1
251    if (kp.eq.maxpointspec_act) then
252  if (ldirect.eq.-1) then
253      ncret = nf_put_vara_real(ncid,ncwdvid, &
254        (/1,1,1,nage,ncirec/),(/grid_nx,grid_ny,kp,1,1/), &
255        wetgrid2(0:grid_nx-1,0:grid_ny-1,1:kp))
256      call check_ncerror(ncret)
257  else
258    if (kp.gt.1) then
259      ncret = nf_put_vara_real(ncid,ncwdvid, &
260        (/1,1,1,ks,nage,ncirec/),(/grid_nx,grid_ny,kp,1,1,1/), &
261        wetgrid2(0:grid_nx-1,0:grid_ny-1,1:kp))
262    else
263      ncret = nf_put_vara_real(ncid,ncwdvid, &
264        (/1,1,ks,nage,ncirec/),(/grid_nx,grid_ny,1,1,1/), &
265        wetgrid2(0:grid_nx-1,0:grid_ny-1,1))
266    endif
267      call check_ncerror(ncret)
268    endif
269  endif
270!    do jy=0,grid_ny-1
271!    do ix=0,grid_nx-1
272!      ncret = nf_put_vara_real(ncid,ncwdvid, &
273!        (/ix+1,jy+1,kp,nage,ncirec/),(/1,1,1,1,1/), &
274!        1.e12*wetgrid(ix,jy)/area(ix,jy))
275!      call check_ncerror(ncret)
276!    enddo ! ix=1,grid_nx-1
277!    enddo ! jy=1,numygrid-1
278  endif ! WETDEP
279
280  if ((ldirect.eq.1).and.(DRYDEP)) then ! DRYDEP
281    if (option_verbose.ge.1)write(*,*)'write_ncconc: dry deposition output'
282    do jy=0,grid_ny-1
283    do ix=0,grid_nx-1
284    if (nesting_level.eq.0)  drygrid2(ix,jy,kp)=1.e12*drygrid(ix,jy)/area(ix,jy)
285    if (nesting_level.eq.1)  drygrid2(ix,jy,kp)=1.e12*drygrid(ix,jy)/arean(ix,jy)
286    enddo ! ix=1,grid_nx-1
287    enddo ! jy=1,grid_ny-1
288    if (kp.eq.maxpointspec_act) then
289  if (ldirect.eq.-1) then
290      ncret = nf_put_vara_real(ncid,ncddvid, &
291        (/1,1,1,nage,ncirec/),(/grid_nx,grid_ny,kp,1,1/), &
292        drygrid2(0:grid_nx-1,0:grid_ny-1,1:kp))
293      call check_ncerror(ncret)
294  else
295    if (kp.gt.1) then
296      ncret = nf_put_vara_real(ncid,ncddvid, &
297        (/1,1,1,ks,nage,ncirec/),(/grid_nx,grid_ny,kp,1,1,1/), &
298        drygrid2(0:grid_nx-1,0:grid_ny-1,1:kp))
299    else
300      ncret = nf_put_vara_real(ncid,ncddvid, &
301        (/1,1,ks,nage,ncirec/),(/grid_nx,grid_ny,1,1,1/), &
302        drygrid2(0:grid_nx-1,0:grid_ny-1,1))
303    endif
304      call check_ncerror(ncret)
305    endif
306  endif
307
308!    do jy=0,grid_ny-1
309!    do ix=0,grid_nx-1
310!      ncret = nf_put_vara_real(ncid,ncddvid, &
311!        (/ix+1,jy+1,kp,nage,ncirec/),(/1,1,1,1,1/), &
312!        1.e12*drygrid(ix,jy)/area(ix,jy))
313!      call check_ncerror(ncret)
314!    enddo ! ix=1,grid_nx-1
315!    enddo ! jy=1,numygrid-1
316  endif ! DRYDEP
317
318   ncret=nf_sync(ncid)
319   call check_ncerror(ncret)
320
321end subroutine write_ncconc
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG