source: flexpart.git/src/writeheader_txt.f90 @ 4c64400

10.4.1_peseiGFS_025bugfixes+enhancementsdevrelease-10release-10.4.1scaling-bugunivie
Last change on this file since 4c64400 was 8a65cb0, checked in by Espen Sollum ATMOS <espen@…>, 9 years ago

Added code, makefile for dev branch

  • Property mode set to 100644
File size: 7.7 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 writeheader_txt
23
24  !*****************************************************************************
25  !                                                                            *
26  !  This routine produces a file header containing basic information on the   *
27  !  settings of the FLEXPART run.                                             *
28  !  The header file is essential and must be read in by any postprocessing    *
29  !  program before reading in the output data.                                *
30  !                                                                            *
31  !     Author: A. Stohl                                                       *
32  !                                                                            *
33  !     7 August 2002                                                          *
34  !     modified IP 2013 for text output                                       *
35  !*****************************************************************************
36  !                                                                            *
37  ! Variables:                                                                 *
38  !                                                                            *
39  ! xlon                   longitude                                           *
40  ! xl                     model x coordinate                                  *
41  ! ylat                   latitude                                            *
42  ! yl                     model y coordinate                                  *
43  !                                                                            *
44  !*****************************************************************************
45
46  use point_mod
47  use outg_mod
48  use par_mod
49  use com_mod
50
51  implicit none
52
53!  integer :: jjjjmmdd,ihmmss,i,ix,jy,j
54  integer :: jjjjmmdd,ihmmss,i,j
55  real :: xp1,yp1,xp2,yp2
56
57
58  !************************
59  ! Open header output file
60  !************************
61
62  open(unitheader,file=path(2)(1:length(2))//'header_txt', &
63       form='formatted',err=998)
64  open(unitheader_txt,file=path(2)(1:length(2))//'header_txt_releases', &
65       form='formatted',err=998)
66
67
68  ! Write the header information
69  !*****************************
70 
71  write(unitheader,*) '# ibdate,ibtime, iedate, ietime, flexversion'
72  write(unitheader,*) ibdate, ibtime, iedate, ietime, trim(flexversion) !  'FLEXPART V9.0'
73  !if (ldirect.eq.1) then
74  !  write(unitheader,*) ibdate,ibtime,trim(flexversion) !  'FLEXPART V9.0'
75  !else
76  !  write(unitheader,*) iedate,ietime,trim(flexversion) ! 'FLEXPART V9.0'
77  !endif
78
79  ! Write info on output interval, averaging time, sampling time
80  !*************************************************************
81 
82  write(unitheader,*) '# interval, averaging time, sampling time'
83  write(unitheader,*) loutstep,loutaver,loutsample
84
85  ! Write information on output grid setup
86  !***************************************
87 
88  write(unitheader,*) '# information on grid setup    '
89  write(unitheader,*) '#outlon0,outlat0,numxgrid,numygrid,dxout,dyout'
90  write(unitheader,*) outlon0,outlat0,numxgrid,numygrid, &
91       dxout,dyout 
92  write(unitheader,*) '# numzgrid, outheight(.) '
93  write(unitheader,*) numzgrid,(outheight(i),i=1,numzgrid)
94
95  write(unitheader,*) '# jjjjmmdd,ihmmss'
96  call caldate(bdate,jjjjmmdd,ihmmss)
97  write(unitheader,*) jjjjmmdd,ihmmss
98
99  ! Write number of species, and name for each species (+extra name for depositions)
100  ! Indicate the vertical dimension of the fields (i.e., 1 for deposition fields, numzgrid for
101  ! concentration fields
102  !*****************************************************************************
103
104  write(unitheader,*) '# information on species'
105  write(unitheader,*) '# 3*nspec,maxpointspec_act'
106  write(unitheader,*) 3*nspec,maxpointspec_act
107  write(unitheader,*) '# for nspec:'
108  write(unitheader,*) '# 1, WD_ '
109  write(unitheader,*) '# 1, DD_ '
110  write(unitheader,*) '# numzgrid,species'
111  do i=1,nspec
112    write(unitheader,*) 1,'WD_'//species(i)(1:7)
113    write(unitheader,*) 1,'DD_'//species(i)(1:7)
114    write(unitheader,*) numzgrid,species(i)
115  end do
116
117  ! Write information on release points: total number, then for each point:
118  ! start, end, coordinates, # of particles, name, mass
119  !************************************************************************
120
121
122  write(unitheader_txt,*) '# information on release points'
123  write(unitheader_txt,*) '# numpoint'
124  write(unitheader_txt,*) numpoint
125  write(unitheader_txt,*) '# for numpoint:'
126  do i=1,numpoint
127    write(unitheader_txt,*) ireleasestart(i),ireleaseend(i),kindz(i)
128    xp1=xpoint1(i)*dx+xlon0
129    yp1=ypoint1(i)*dy+ylat0
130    xp2=xpoint2(i)*dx+xlon0
131    yp2=ypoint2(i)*dy+ylat0
132    write(unitheader_txt,*) xp1,yp1,xp2,yp2,zpoint1(i),zpoint2(i)
133    write(unitheader_txt,*) npart(i),1
134    if (numpoint.le.1000) then
135      write(unitheader_txt,*) compoint(i)
136    else
137      write(unitheader_txt,*) compoint(1001)
138    endif
139    do j=1,nspec
140      write(unitheader_txt,*) xmass(i,j)
141      write(unitheader_txt,*) xmass(i,j)
142      write(unitheader_txt,*) xmass(i,j)
143    end do
144  end do
145
146  ! Write information on model switches
147  !*****************************************
148
149  write(unitheader,*) '# information on model switches'
150  write(unitheader,*) '# method,lsubgrid,lconvection, ind_source,ind_receptor'
151  write(unitheader,*) method,lsubgrid,lconvection, &
152       ind_source,ind_receptor
153
154  ! Write age class information
155  !****************************
156 
157  write(unitheader,*) '# information on age class     '
158  write(unitheader,*) nageclass,(lage(i),i=1,nageclass)
159
160
161  !Do not write topography to text output file. Keep it on the binary one
162  !********************************
163
164  !do ix=0,numxgrid-1
165  !  write(unitheader,*) (oroout(ix,jy),jy=0,numygrid-1)
166  !end do
167
168
169 
170
171
172  close(unitheader)
173  close(unitheader_txt)
174
175
176!  open(unitheader,file=path(2)(1:length(2))//'header_nml', &
177!        form='formatted',err=998)
178!  write(unitheader,NML=COMMAND)
179!  close(unitheader)
180
181  return
182
183
184998   write(*,*) ' #### FLEXPART MODEL ERROR!   THE FILE         #### '
185  write(*,*) ' #### '//path(2)(1:length(2))//'header_txt'//' #### '
186  write(*,*) ' #### CANNOT BE OPENED. IF A FILE WITH THIS    #### '
187  write(*,*) ' #### NAME ALREADY EXISTS, DELETE IT AND START #### '
188  write(*,*) ' #### THE PROGRAM AGAIN.                       #### '
189  stop
190
191end subroutine writeheader_txt
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG