source: trunk/src/readoutgrid.f90 @ 24

Last change on this file since 24 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: 7.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 readoutgrid
23
24  !*****************************************************************************
25  !                                                                            *
26  !     This routine reads the user specifications for the output grid.        *
27  !                                                                            *
28  !     Author: A. Stohl                                                       *
29  !                                                                            *
30  !     4 June 1996                                                            *
31  !                                                                            *
32  !*****************************************************************************
33  !                                                                            *
34  ! Variables:                                                                 *
35  ! dxout,dyout          grid distance                                         *
36  ! numxgrid,numygrid,numzgrid    grid dimensions                              *
37  ! outlon0,outlat0      lower left corner of grid                             *
38  ! outheight(maxzgrid)  height levels of output grid [m]                      *
39  !                                                                            *
40  ! Constants:                                                                 *
41  ! unitoutgrid          unit connected to file OUTGRID                        *
42  !                                                                            *
43  !*****************************************************************************
44
45  use outg_mod
46  use par_mod
47  use com_mod
48
49  implicit none
50
51  integer :: i,j,stat
52  real :: outhelp,xr,xr1,yr,yr1
53  real,parameter :: eps=1.e-4
54
55
56
57  ! Open the OUTGRID file and read output grid specifications
58  !**********************************************************
59
60  open(unitoutgrid,file=path(1)(1:length(1))//'OUTGRID',status='old', &
61       err=999)
62
63
64  call skplin(5,unitoutgrid)
65
66
67  ! 1.  Read horizontal grid specifications
68  !****************************************
69
70  call skplin(3,unitoutgrid)
71  read(unitoutgrid,'(4x,f11.4)') outlon0
72  call skplin(3,unitoutgrid)
73  read(unitoutgrid,'(4x,f11.4)') outlat0
74  call skplin(3,unitoutgrid)
75  read(unitoutgrid,'(4x,i5)') numxgrid
76  call skplin(3,unitoutgrid)
77  read(unitoutgrid,'(4x,i5)') numygrid
78  call skplin(3,unitoutgrid)
79  read(unitoutgrid,'(4x,f12.5)') dxout
80  call skplin(3,unitoutgrid)
81  read(unitoutgrid,'(4x,f12.5)') dyout
82
83
84  ! Check validity of output grid (shall be within model domain)
85  !*************************************************************
86
87  xr=outlon0+real(numxgrid)*dxout
88  yr=outlat0+real(numygrid)*dyout
89  xr1=xlon0+real(nxmin1)*dx
90  yr1=ylat0+real(nymin1)*dy
91  if ((outlon0+eps.lt.xlon0).or.(outlat0+eps.lt.ylat0) &
92       .or.(xr.gt.xr1+eps).or.(yr.gt.yr1+eps)) then
93    write(*,*) 'outlon0,outlat0:'
94    write(*,*) outlon0,outlat0
95    write(*,*) 'xr1,yr1,xlon0,ylat0,xr,yr,dxout,dyout:'
96    write(*,*) xr1,yr1,xlon0,ylat0,xr,yr,dxout,dyout
97    write(*,*) ' #### FLEXPART MODEL ERROR! PART OF OUTPUT    ####'
98    write(*,*) ' #### GRID IS OUTSIDE MODEL DOMAIN. CHANGE    ####'
99    write(*,*) ' #### FILE OUTGRID IN DIRECTORY               ####'
100    write(*,'(a)') path(1)(1:length(1))
101    stop
102  endif
103
104  ! 2. Count Vertical levels of output grid
105  !****************************************
106  j=0
107100   j=j+1
108    do i=1,3
109      read(unitoutgrid,*,end=99)
110    end do
111    read(unitoutgrid,'(4x,f7.1)',end=99) outhelp
112    if (outhelp.eq.0.) goto 99
113    goto 100
11499   numzgrid=j-1
115
116    allocate(outheight(numzgrid) &
117         ,stat=stat)
118    if (stat.ne.0) write(*,*)'ERROR: could not allocate outh'
119    allocate(outheighthalf(numzgrid) &
120         ,stat=stat)
121    if (stat.ne.0) write(*,*)'ERROR: could not allocate outh'
122
123
124  rewind(unitoutgrid)
125  call skplin(29,unitoutgrid)
126
127  ! 2. Vertical levels of output grid
128  !**********************************
129
130  j=0
1311000   j=j+1
132    do i=1,3
133      read(unitoutgrid,*,end=990)
134    end do
135    read(unitoutgrid,'(4x,f7.1)',end=990) outhelp
136    if (outhelp.eq.0.) goto 99
137    outheight(j)=outhelp
138    goto 1000
139990   numzgrid=j-1
140
141
142  ! Check whether vertical levels are specified in ascending order
143  !***************************************************************
144
145  do j=2,numzgrid
146    if (outheight(j).le.outheight(j-1)) then
147    write(*,*) ' #### FLEXPART MODEL ERROR! YOUR SPECIFICATION#### '
148    write(*,*) ' #### OF OUTPUT LEVELS IS CORRUPT AT LEVEL    #### '
149    write(*,*) ' #### ',j,'                              #### '
150    write(*,*) ' #### PLEASE MAKE CHANGES IN FILE OUTGRID.    #### '
151    endif
152  end do
153
154  ! Determine the half levels, i.e. middle levels of the output grid
155  !*****************************************************************
156
157  outheighthalf(1)=outheight(1)/2.
158  do j=2,numzgrid
159    outheighthalf(j)=(outheight(j-1)+outheight(j))/2.
160  end do
161
162
163  xoutshift=xlon0-outlon0
164  youtshift=ylat0-outlat0
165  close(unitoutgrid)
166
167    allocate(oroout(0:numxgrid-1,0:numygrid-1) &
168         ,stat=stat)
169    if (stat.ne.0) write(*,*)'ERROR: could not allocate outh'
170    allocate(area(0:numxgrid-1,0:numygrid-1) &
171         ,stat=stat)
172    if (stat.ne.0) write(*,*)'ERROR: could not allocate outh'
173    allocate(volume(0:numxgrid-1,0:numygrid-1,numzgrid) &
174         ,stat=stat)
175    if (stat.ne.0) write(*,*)'ERROR: could not allocate outh'
176    allocate(areaeast(0:numxgrid-1,0:numygrid-1,numzgrid) &
177         ,stat=stat)
178    if (stat.ne.0) write(*,*)'ERROR: could not allocate outh'
179    allocate(areanorth(0:numxgrid-1,0:numygrid-1,numzgrid) &
180         ,stat=stat)
181    if (stat.ne.0) write(*,*)'ERROR: could not allocate outh'
182  return
183
184
185999   write(*,*) ' #### FLEXPART MODEL ERROR! FILE "OUTGRID"    #### '
186  write(*,*) ' #### CANNOT BE OPENED IN THE DIRECTORY       #### '
187  write(*,*) ' #### xxx/flexpart/options                    #### '
188  stop
189
190end subroutine readoutgrid
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG