source: flexpart.git/src/readoutgrid.f90 @ 3481cc1

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

move license from headers to a different file
  • Property mode set to 100644
File size: 7.3 KB
Line 
1subroutine readoutgrid
2
3  !*****************************************************************************
4  !                                                                            *
5  !     This routine reads the user specifications for the output grid.        *
6  !                                                                            *
7  !     Author: A. Stohl                                                       *
8  !                                                                            *
9  !     4 June 1996                                                            *
10  !     HSO, 1 July 2014
11  !     Added optional namelist input
12  !                                                                            *
13  !*****************************************************************************
14  !                                                                            *
15  ! Variables:                                                                 *
16  ! dxout,dyout          grid distance                                         *
17  ! numxgrid,numygrid,numzgrid    grid dimensions                              *
18  ! outlon0,outlat0      lower left corner of grid                             *
19  ! outheight(maxzgrid)  height levels of output grid [m]                      *
20  !                                                                            *
21  ! Constants:                                                                 *
22  ! unitoutgrid          unit connected to file OUTGRID                        *
23  !                                                                            *
24  !*****************************************************************************
25
26  use outg_mod
27  use par_mod
28  use com_mod
29
30  implicit none
31
32  integer :: i,j,stat
33  real :: outhelp,xr,xr1,yr,yr1
34  real,parameter :: eps=1.e-4
35
36  ! namelist variables
37  integer, parameter :: maxoutlev=500
38  integer :: readerror
39  real,allocatable, dimension (:) :: outheights
40
41  ! declare namelist
42  namelist /outgrid/ &
43    outlon0,outlat0, &
44    numxgrid,numygrid, &
45    dxout,dyout, &
46    outheights
47
48  ! allocate large array for reading input
49  allocate(outheights(maxoutlev),stat=stat)
50  if (stat.ne.0) write(*,*)'ERROR: could not allocate outheights'
51
52  ! helps identifying failed namelist input
53  dxout=-1.0
54  outheights=-1.0
55
56  ! Open the OUTGRID file and read output grid specifications
57  !**********************************************************
58
59  open(unitoutgrid,file=path(1)(1:length(1))//'OUTGRID',status='old',form='formatted',err=999)
60
61  ! try namelist input
62  read(unitoutgrid,outgrid,iostat=readerror)
63  close(unitoutgrid)
64
65  if ((dxout.le.0).or.(readerror.ne.0)) then
66
67    readerror=1
68
69    open(unitoutgrid,file=path(1)(1:length(1))//'OUTGRID',status='old',err=999)
70
71    call skplin(5,unitoutgrid)
72
73    ! 1.  Read horizontal grid specifications
74    !****************************************
75
76    call skplin(3,unitoutgrid)
77    read(unitoutgrid,'(4x,f11.4)') outlon0
78    call skplin(3,unitoutgrid)
79    read(unitoutgrid,'(4x,f11.4)') outlat0
80    call skplin(3,unitoutgrid)
81    read(unitoutgrid,'(4x,i5)') numxgrid
82    call skplin(3,unitoutgrid)
83    read(unitoutgrid,'(4x,i5)') numygrid
84    call skplin(3,unitoutgrid)
85    read(unitoutgrid,'(4x,f12.5)') dxout
86    call skplin(3,unitoutgrid)
87    read(unitoutgrid,'(4x,f12.5)') dyout
88
89  endif
90
91  ! Check validity of output grid (shall be within model domain)
92  !*************************************************************
93
94  xr=outlon0+real(numxgrid)*dxout
95  yr=outlat0+real(numygrid)*dyout
96  xr1=xlon0+real(nxmin1)*dx
97  yr1=ylat0+real(nymin1)*dy
98  if ((outlon0+eps.lt.xlon0).or.(outlat0+eps.lt.ylat0) &
99       .or.(xr.gt.xr1+eps).or.(yr.gt.yr1+eps)) then
100    write(*,*) outlon0,outlat0
101    write(*,*) xr1,yr1,xlon0,ylat0,xr,yr,dxout,dyout
102    write(*,*) ' #### FLEXPART MODEL ERROR! PART OF OUTPUT    ####'
103    write(*,*) ' #### GRID IS OUTSIDE MODEL DOMAIN. CHANGE    ####'
104    write(*,*) ' #### FILE OUTGRID IN DIRECTORY               ####'
105    write(*,'(a)') path(1)(1:length(1))
106    stop
107  endif
108
109  ! 2. Count Vertical levels of output grid
110  !****************************************
111
112  if (readerror.ne.0) then
113    j=0
114100 j=j+1
115    do i=1,3
116      read(unitoutgrid,*,end=99)
117    end do
118    read(unitoutgrid,'(4x,f7.1)',end=99) outhelp
119    if (outhelp.eq.0.) goto 99
120    goto 100
12199  numzgrid=j-1
122  else
123    do i=1,maxoutlev
124      if (outheights(i).lt.0) exit
125    end do
126    numzgrid=i-1
127  end if
128
129  allocate(outheight(numzgrid),stat=stat)
130  if (stat.ne.0) write(*,*)'ERROR: could not allocate outheight'
131  allocate(outheighthalf(numzgrid),stat=stat)
132  if (stat.ne.0) write(*,*)'ERROR: could not allocate outheighthalf'
133
134  ! 2. Vertical levels of output grid
135  !**********************************
136
137  if (readerror.ne.0) then
138
139    rewind(unitoutgrid)
140    call skplin(29,unitoutgrid)
141
142    do j=1,numzgrid
143      do i=1,3
144        read(unitoutgrid,*)
145      end do
146      read(unitoutgrid,'(4x,f7.1)') outhelp
147      outheight(j)=outhelp
148      outheights(j)=outhelp
149    end do
150    close(unitoutgrid)
151
152  else
153
154    do j=1,numzgrid
155      outheight(j)=outheights(j)
156    end do
157
158  endif
159
160  ! write outgrid file in namelist format to output directory if requested
161  if (nmlout.and.lroot) then
162    ! reallocate outheights with actually required dimension for namelist writing
163    deallocate(outheights)
164    allocate(outheights(numzgrid),stat=stat)
165    if (stat.ne.0) write(*,*)'ERROR: could not allocate outheights'
166
167    do j=1,numzgrid
168      outheights(j)=outheight(j)
169    end do
170
171    open(unitoutgrid,file=path(2)(1:length(2))//'OUTGRID.namelist',err=1000)
172    write(unitoutgrid,nml=outgrid)
173    close(unitoutgrid)
174  endif
175
176  ! Check whether vertical levels are specified in ascending order
177  !***************************************************************
178
179  do j=2,numzgrid
180    if (outheight(j).le.outheight(j-1)) then
181    write(*,*) ' #### FLEXPART MODEL ERROR! YOUR SPECIFICATION#### '
182    write(*,*) ' #### OF OUTPUT LEVELS IS CORRUPT AT LEVEL    #### '
183    write(*,*) ' #### ',j,'                              #### '
184    write(*,*) ' #### PLEASE MAKE CHANGES IN FILE OUTGRID.    #### '
185    endif
186  end do
187
188  ! Determine the half levels, i.e. middle levels of the output grid
189  !*****************************************************************
190
191  outheighthalf(1)=outheight(1)/2.
192  do j=2,numzgrid
193    outheighthalf(j)=(outheight(j-1)+outheight(j))/2.
194  end do
195
196  xoutshift=xlon0-outlon0
197  youtshift=ylat0-outlat0
198
199  allocate(oroout(0:numxgrid-1,0:numygrid-1),stat=stat)
200  if (stat.ne.0) write(*,*)'ERROR: could not allocate oroout'
201  allocate(area(0:numxgrid-1,0:numygrid-1),stat=stat)
202  if (stat.ne.0) write(*,*)'ERROR: could not allocate area'
203  allocate(volume(0:numxgrid-1,0:numygrid-1,numzgrid),stat=stat)
204  if (stat.ne.0) write(*,*)'ERROR: could not allocate volume'
205  allocate(areaeast(0:numxgrid-1,0:numygrid-1,numzgrid),stat=stat)
206  if (stat.ne.0) write(*,*)'ERROR: could not allocate areaeast'
207  allocate(areanorth(0:numxgrid-1,0:numygrid-1,numzgrid),stat=stat)
208  if (stat.ne.0) write(*,*)'ERROR: could not allocate areanorth'
209  return
210
211
212999 write(*,*) ' #### FLEXPART MODEL ERROR! FILE "OUTGRID"    #### '
213  write(*,*) ' #### CANNOT BE OPENED IN THE DIRECTORY       #### '
214  write(*,'(a)') path(1)(1:length(1))
215  stop
216
2171000 write(*,*) ' #### FLEXPART MODEL ERROR! FILE "OUTGRID"    #### '
218  write(*,*) ' #### CANNOT BE OPENED IN THE DIRECTORY       #### '
219  write(*,'(a)') path(2)(1:length(2))
220  stop
221
222end subroutine readoutgrid
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG