source: trunk/src/readreceptors.f90 @ 29

Last change on this file since 29 was 27, checked in by hasod, 10 years ago
  • Implemented optional namelist input for COMMAND, RELEASES, SPECIES, AGECLASSES,OUTGRID,OUTGRID_NEST,RECEPTORS
  • Implemented com_mod switch nmlout to write input files as namelist to the output directory (.true. by default)
  • Proposed updated startup and runtime output (may change back to previous info if desired)
File size: 7.0 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 readreceptors
23
24  !*****************************************************************************
25  !                                                                            *
26  !     This routine reads the user specifications for the receptor points.    *
27  !                                                                            *
28  !     Author: A. Stohl                                                       *
29  !     1 August 1996                                                          *
30  !     HSO, 14 August 2013
31  !     Added optional namelist input
32  !                                                                            *
33  !*****************************************************************************
34  !                                                                            *
35  ! Variables:                                                                 *
36  ! receptorarea(maxreceptor)  area of dx*dy at location of receptor           *
37  ! receptorname(maxreceptor)  names of receptors                              *
38  ! xreceptor,yreceptor  coordinates of receptor points                        *
39  !                                                                            *
40  ! Constants:                                                                 *
41  ! unitreceptor         unit connected to file RECEPTORS                      *
42  !                                                                            *
43  !*****************************************************************************
44
45  use par_mod
46  use com_mod
47
48  implicit none
49
50  integer :: j
51  real :: x,y,xm,ym
52  character(len=16) :: receptor
53
54  integer :: readerror
55  real :: lon,lat   ! for namelist input, lon/lat are used instead of x,y
56  integer,parameter :: unitreceptorout=2
57
58  ! declare namelist
59  namelist /receptors/ &
60    receptor, lon, lat
61
62  lon=-999.9
63
64  ! For backward runs, do not allow receptor output. Thus, set number of receptors to zero
65  !*****************************************************************************
66
67  if (ldirect.lt.0) then
68    numreceptor=0
69    return
70  endif
71
72
73  ! Open the RECEPTORS file and read output grid specifications
74  !************************************************************
75
76  open(unitreceptor,file=path(1)(1:length(1))//'RECEPTORS',form='formatted',status='old',err=999)
77
78  ! try namelist input
79  read(unitreceptor,receptors,iostat=readerror)
80
81  ! prepare namelist output if requested
82  if (nmlout.eqv..true.) then
83    open(unitreceptorout,file=path(2)(1:length(2))//'RECEPTORS.namelist',access='append',status='new',err=1000)
84  endif
85
86  if ((lon.lt.-900).or.(readerror.ne.0)) then
87
88    close(unitreceptor)
89    open(unitreceptor,file=path(1)(1:length(1))//'RECEPTORS',status='old',err=999)
90    call skplin(5,unitreceptor)
91
92    ! Read the names and coordinates of the receptors
93    !************************************************
94
95    j=0
96100 j=j+1
97    read(unitreceptor,*,end=99)
98    read(unitreceptor,*,end=99)
99    read(unitreceptor,*,end=99)
100    read(unitreceptor,'(4x,a16)',end=99) receptor
101    call skplin(3,unitreceptor)
102    read(unitreceptor,'(4x,f11.4)',end=99) x
103    call skplin(3,unitreceptor)
104    read(unitreceptor,'(4x,f11.4)',end=99) y
105    if ((x.eq.0.).and.(y.eq.0.).and. &
106         (receptor.eq.'                ')) then
107      j=j-1
108      goto 100
109    endif
110    if (j.gt.maxreceptor) then
111      write(*,*) ' #### FLEXPART MODEL ERROR! TOO MANY RECEPTOR #### '
112      write(*,*) ' #### POINTS ARE GIVEN.                       #### '
113      write(*,*) ' #### MAXIMUM NUMBER IS ',maxreceptor,'       #### '
114      write(*,*) ' #### PLEASE MAKE CHANGES IN FILE RECEPTORS   #### '
115    endif
116    receptorname(j)=receptor
117    xreceptor(j)=(x-xlon0)/dx       ! transform to grid coordinates
118    yreceptor(j)=(y-ylat0)/dy
119    xm=r_earth*cos(y*pi/180.)*dx/180.*pi
120    ym=r_earth*dy/180.*pi
121    receptorarea(j)=xm*ym
122
123    ! write receptors file in namelist format to output directory if requested
124    if (nmlout.eqv..true.) then
125      lon=x
126      lat=y
127      write(unitreceptorout,nml=receptors)
128    endif
129
130    goto 100
131
13299  numreceptor=j-1
133
134  else ! continue with namelist input
135
136    j=0
137    do while (readerror.eq.0)
138      j=j+1
139      lon=-999.9
140      read(unitreceptor,receptors,iostat=readerror)
141      if ((lon.lt.-900).or.(readerror.ne.0)) then
142        readerror=1
143      else
144        if (j.gt.maxreceptor) then
145          write(*,*) ' #### FLEXPART MODEL ERROR! TOO MANY RECEPTOR #### '
146          write(*,*) ' #### POINTS ARE GIVEN.                       #### '
147          write(*,*) ' #### MAXIMUM NUMBER IS ',maxreceptor,'       #### '
148          write(*,*) ' #### PLEASE MAKE CHANGES IN FILE RECEPTORS   #### '
149        endif
150        receptorname(j)=receptor
151        xreceptor(j)=(lon-xlon0)/dx       ! transform to grid coordinates
152        yreceptor(j)=(lat-ylat0)/dy
153        xm=r_earth*cos(lat*pi/180.)*dx/180.*pi
154        ym=r_earth*dy/180.*pi
155        receptorarea(j)=xm*ym
156      endif
157
158      ! write receptors file in namelist format to output directory if requested
159      if (nmlout.eqv..true.) then
160        write(unitreceptorout,nml=receptors)
161      endif
162
163    end do
164    numreceptor=j-1
165    close(unitreceptor)
166
167  endif
168
169  if (nmlout.eqv..true.) then
170    close(unitreceptorout)
171  endif
172
173  return
174
175
176999 write(*,*) ' #### FLEXPART MODEL ERROR! FILE "RECEPTORS"  #### '
177  write(*,*) ' #### CANNOT BE OPENED IN THE DIRECTORY       #### '
178  write(*,'(a)') path(1)(1:length(1))
179  stop
180
1811000 write(*,*) ' #### FLEXPART MODEL ERROR! FILE "RECEPTORS"    #### '
182  write(*,*) ' #### CANNOT BE OPENED IN THE DIRECTORY       #### '
183  write(*,'(a)') path(2)(1:length(2))
184  stop
185
186end subroutine readreceptors
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG