source: trunk/src/writeheader_surf.f90 @ 28

Last change on this file since 28 was 20, checked in by igpis, 10 years ago

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

File size: 6.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 writeheader_surf
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  !                                                                            *
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  real :: xp1,yp1,xp2,yp2
55
56
57  !************************
58  ! Open header output file
59  !************************
60
61  open(unitheader,file=path(2)(1:length(2))//'header_grid_time', &
62       form='unformatted',err=998)
63
64
65  ! Write the header information
66  !*****************************
67
68  if (ldirect.eq.1) then
69    write(unitheader) ibdate,ibtime, trim(flexversion)
70  else
71    write(unitheader) iedate,ietime, trim(flexversion)
72  endif
73
74  ! Write info on output interval, averaging time, sampling time
75  !*************************************************************
76
77  write(unitheader) loutstep,loutaver,loutsample
78
79  ! Write information on output grid setup
80  !***************************************
81
82  write(unitheader) outlon0,outlat0,numxgrid,numygrid, &
83       dxout,dyout
84  write(unitheader) 1,(outheight(1),i=1,1)
85
86  call caldate(bdate,jjjjmmdd,ihmmss)
87  write(unitheader) jjjjmmdd,ihmmss
88
89  ! Write number of species, and name for each species (+extra name for depositions)
90  ! Indicate the dimension of the fields (i.e., 1 for deposition fields, numzgrid for
91  ! concentration fields
92  !*****************************************************************************
93
94  write(unitheader) 3*nspec,maxpointspec_act
95  do i=1,nspec
96    write(unitheader) 1,'WD_'//species(i)(1:7)
97    write(unitheader) 1,'DD_'//species(i)(1:7)
98    write(unitheader) 1,species(i)
99  end do
100
101  ! Write information on release points: total number, then for each point:
102  ! start, end, coordinates, # of particles, name, mass
103  !************************************************************************
104
105  write(unitheader) numpoint
106  do i=1,numpoint
107    write(unitheader) ireleasestart(i),ireleaseend(i),kindz(i)
108    xp1=xpoint1(i)*dx+xlon0
109    yp1=ypoint1(i)*dy+ylat0
110    xp2=xpoint2(i)*dx+xlon0
111    yp2=ypoint2(i)*dy+ylat0
112    write(unitheader) xp1,yp1,xp2,yp2,zpoint1(i),zpoint2(i)
113    write(unitheader) npart(i),1
114    if (numpoint.le.1000) then
115      write(unitheader) compoint(i)
116    else
117      write(unitheader) compoint(1001)
118    endif
119    do j=1,nspec
120      write(unitheader) xmass(i,j)
121      write(unitheader) xmass(i,j)
122      write(unitheader) xmass(i,j)
123    end do
124  end do
125
126  ! Write information on some model switches
127  !*****************************************
128
129  write(unitheader) method,lsubgrid,lconvection, &
130       ind_source,ind_receptor
131
132  ! Write age class information
133  !****************************
134
135  write(unitheader) nageclass,(lage(i),i=1,nageclass)
136
137
138  ! Write topography to output file
139  !********************************
140
141  do ix=0,numxgrid-1
142    write(unitheader) (oroout(ix,jy),jy=0,numygrid-1)
143  end do
144  close(unitheader)
145
146  return
147
148
149998   write(*,*) ' #### FLEXPART MODEL ERROR!   THE FILE         #### '
150  write(*,*) ' #### '//path(2)(1:length(2))//'header'//' #### '
151  write(*,*) ' #### CANNOT BE OPENED. IF A FILE WITH THIS    #### '
152  write(*,*) ' #### NAME ALREADY EXISTS, DELETE IT AND START #### '
153  write(*,*) ' #### THE PROGRAM AGAIN.                       #### '
154  stop
155
156end subroutine writeheader_surf
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG