source: flexpart.git/src/writeheader_surf.f90

10.4.1_peseiGFS_025bugfixes+enhancementsdevrelease-10release-10.4.1scaling-bug
Last change on this file was 92fab65, checked in by Ignacio Pisso <ip@…>, 4 years ago

add SPDX-License-Identifier to all .f90 files

  • Property mode set to 100644
File size: 5.5 KB
Line 
1! SPDX-FileCopyrightText: FLEXPART 1998-2019, see flexpart_license.txt
2! SPDX-License-Identifier: GPL-3.0-or-later
3
4subroutine writeheader_surf
5
6  !*****************************************************************************
7  !                                                                            *
8  !  This routine produces a file header containing basic information on the   *
9  !  settings of the FLEXPART run.                                             *
10  !  The header file is essential and must be read in by any postprocessing    *
11  !  program before reading in the output data.                                *
12  !                                                                            *
13  !     Author: A. Stohl                                                       *
14  !                                                                            *
15  !     7 August 2002                                                          *
16  !                                                                            *
17  !*****************************************************************************
18  !                                                                            *
19  !  Modified to remove TRIM around the output of flexversion so that          *
20  !  it will be a constant length (defined in com_mod.f90) in output header    *
21  !                                                                            *
22  !     Don Morton, Boreal Scientific Computing                                *
23  !     07 May 2017                                                            *
24  !                                                                            *
25  !*****************************************************************************
26  !                                                                            *
27  ! Variables:                                                                 *
28  !                                                                            *
29  ! xlon                   longitude                                           *
30  ! xl                     model x coordinate                                  *
31  ! ylat                   latitude                                            *
32  ! yl                     model y coordinate                                  *
33  !                                                                            *
34  !*****************************************************************************
35
36  use point_mod
37  use outg_mod
38  use par_mod
39  use com_mod
40
41  implicit none
42
43  integer :: jjjjmmdd,ihmmss,i,ix,jy,j
44  real :: xp1,yp1,xp2,yp2
45
46
47  !************************
48  ! Open header output file
49  !************************
50
51  open(unitheader,file=path(2)(1:length(2))//'header_grid_time', &
52       form='unformatted',err=998)
53
54
55  ! Write the header information
56  !*****************************
57
58  if (ldirect.eq.1) then
59    write(unitheader) ibdate,ibtime, flexversion
60  else
61    write(unitheader) iedate,ietime, flexversion
62  endif
63
64  ! Write info on output interval, averaging time, sampling time
65  !*************************************************************
66
67  write(unitheader) loutstep,loutaver,loutsample
68
69  ! Write information on output grid setup
70  !***************************************
71
72  write(unitheader) outlon0,outlat0,numxgrid,numygrid, &
73       dxout,dyout
74  write(unitheader) 1,(outheight(1),i=1,1)
75
76  call caldate(bdate,jjjjmmdd,ihmmss)
77  write(unitheader) jjjjmmdd,ihmmss
78
79  ! Write number of species, and name for each species (+extra name for depositions)
80  ! Indicate the dimension of the fields (i.e., 1 for deposition fields, numzgrid for
81  ! concentration fields
82  !*****************************************************************************
83
84  write(unitheader) 3*nspec,maxpointspec_act
85  do i=1,nspec
86    write(unitheader) 1,'WD_'//species(i)(1:7)
87    write(unitheader) 1,'DD_'//species(i)(1:7)
88    write(unitheader) 1,species(i)
89  end do
90
91  ! Write information on release points: total number, then for each point:
92  ! start, end, coordinates, # of particles, name, mass
93  !************************************************************************
94
95  write(unitheader) numpoint
96  do i=1,numpoint
97    write(unitheader) ireleasestart(i),ireleaseend(i),kindz(i)
98    xp1=xpoint1(i)*dx+xlon0
99    yp1=ypoint1(i)*dy+ylat0
100    xp2=xpoint2(i)*dx+xlon0
101    yp2=ypoint2(i)*dy+ylat0
102    write(unitheader) xp1,yp1,xp2,yp2,zpoint1(i),zpoint2(i)
103    write(unitheader) npart(i),1
104    if (numpoint.le.1000) then
105      write(unitheader) compoint(i)
106    else
107      write(unitheader) compoint(1001)
108    endif
109    do j=1,nspec
110      write(unitheader) xmass(i,j)
111      write(unitheader) xmass(i,j)
112      write(unitheader) xmass(i,j)
113    end do
114  end do
115
116  ! Write information on some model switches
117  !*****************************************
118
119  write(unitheader) method,lsubgrid,lconvection, &
120       ind_source,ind_receptor
121
122  ! Write age class information
123  !****************************
124
125  write(unitheader) nageclass,(lage(i),i=1,nageclass)
126
127
128  ! Write topography to output file
129  !********************************
130
131  do ix=0,numxgrid-1
132    write(unitheader) (oroout(ix,jy),jy=0,numygrid-1)
133  end do
134  close(unitheader)
135
136  return
137
138
139998   write(*,*) ' #### FLEXPART MODEL ERROR!   THE FILE         #### '
140  write(*,*) ' #### '//path(2)(1:length(2))//'header'//' #### '
141  write(*,*) ' #### CANNOT BE OPENED. IF A FILE WITH THIS    #### '
142  write(*,*) ' #### NAME ALREADY EXISTS, DELETE IT AND START #### '
143  write(*,*) ' #### THE PROGRAM AGAIN.                       #### '
144  stop
145
146end subroutine writeheader_surf
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG