source: branches/petra/src/readwind_nests.f90 @ 37

Last change on this file since 37 was 37, checked in by pesei, 9 years ago

Wet dep quick fix and other small changes. Wet depo quick fix not final yet.

File size: 19.2 KB
Line 
1!**********************************************************************
2! Copyright 1998-2015                                                 *
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 readwind_nests(indj,n,uuhn,vvhn,wwhn)
23  !                           i   i  o    o    o
24  !*****************************************************************************
25  !                                                                            *
26  !     This routine reads the wind fields for the nested model domains.       *
27  !     It is similar to subroutine readwind, which reads the mother domain.   *
28  !                                                                            *
29  !     Authors: A. Stohl, G. Wotawa                                           *
30  !                                                                            *
31  !     8 February 1999                                                        *
32  !                                                                            *
33  !     Last update: 17 October 2000, A. Stohl                                 *
34  !                                                                            *
35  !*****************************************************************************
36  !  CHANGES:
37  !
38  !  2/2001, Bernc C. Krueger:                                              *
39  !        Variables tthn and qvhn (on eta coordinates) in common block        *
40  !  11/01/2008, Harald Sodemann, GRIB1/2 input with ECMWF grib_api            *
41  !  03/12/2008, Harald Sodemann, update to f90 with ECMWF grib_api            *
42  !  2/2015, PS: add missing paramter iret in call to grib subr
43  !  3/2015, PS: add verbosity output
44  !  3/2015, PS: add GRIB2 coding for deta/dt dp/deta used by ZAMG
45  !
46  !*****************************************************************************
47
48  use grib_api
49  use par_mod
50  use com_mod
51
52  implicit none
53
54  !HSO  parameters for grib_api
55  integer :: ifile
56  integer :: iret
57  integer :: igrib
58  integer :: gribVer,parCat,parNum,typSurf,valSurf,discipl
59  integer :: parId !!added by mc for making it consistent with new readwind.f90
60  integer :: gotGrid
61  !HSO  end
62
63  real :: uuhn(0:nxmaxn-1,0:nymaxn-1,nuvzmax,maxnests)
64  real :: vvhn(0:nxmaxn-1,0:nymaxn-1,nuvzmax,maxnests)
65  real :: wwhn(0:nxmaxn-1,0:nymaxn-1,nwzmax,maxnests)
66  integer :: indj,i,j,k,n,levdiff2,ifield,iumax,iwmax,l
67
68  ! VARIABLES AND ARRAYS NEEDED FOR GRIB DECODING
69
70  ! dimension of isec2 at least (22+n), where n is the number of parallels or
71  ! meridians in a quasi-regular (reduced) Gaussian or lat/long grid
72
73  ! dimension of zsec2 at least (10+nn), where nn is the number of vertical
74  ! coordinate parameters
75
76  integer :: isec1(56),isec2(22+nxmaxn+nymaxn)
77  real(kind=4) :: zsec4(jpunp)
78  real(kind=4) :: xaux,yaux
79  real(kind=8) :: xauxin,yauxin
80  real,parameter :: eps=1.e-4
81  real :: ewss(0:nxmaxn-1,0:nymaxn-1),nsss(0:nxmaxn-1,0:nymaxn-1)
82  real :: plev1,pmean,tv,fu,hlev1,ff10m,fflev1
83  real :: conversion_factor !added by mc to make it consistent with new gridchek.f90
84
85  logical :: hflswitch,strswitch
86
87  !HSO  grib api error messages
88  character(len=24) :: gribErrorMsg = 'Error reading grib file'
89  character(len=20) :: gribFunction = 'readwind_nests'
90 
91  call verboutput(verbosity,1,'      entered readwind_nests')
92
93  do l=1,numbnests
94    hflswitch=.false.
95    strswitch=.false.
96    levdiff2=nlev_ec-nwz+1
97    iumax=0
98    iwmax=0
99
100    ifile=0
101    igrib=0
102    iret=0
103
104  !
105  ! OPENING OF DATA FILE (GRIB CODE)
106  !
107
1085 continue
109  call verboutput(verbosity,1,'        call grib_open_file '//&
110    path(numpath+2*(l-1)+1)(1:length(numpath+2*(l-1)+1))//trim(wfnamen(l,indj)))
111  call grib_open_file(ifile,path(numpath+2*(l-1)+1) &
112         (1:length(numpath+2*(l-1)+1))//trim(wfnamen(l,indj)),'r')
113  if (iret.ne.GRIB_SUCCESS) then
114    goto 888   ! ERROR DETECTED
115  endif
116  !turn on support for multi fields messages */
117  !call grib_multi_support_on
118
119    gotGrid=0
120    ifield=0
12110   ifield=ifield+1
122  !
123  ! GET NEXT FIELDS
124  !
125  call grib_new_from_file(ifile,igrib,iret)
126  if (iret.eq.GRIB_END_OF_FILE)  then
127    goto 50    ! EOF DETECTED
128  elseif (iret.ne.GRIB_SUCCESS) then
129    goto 888   ! ERROR DETECTED
130  endif
131
132  !first see if we read GRIB1 or GRIB2
133  call grib_get_int(igrib,'editionNumber',gribVer,iret)
134  call grib_check(iret,gribFunction,gribErrorMsg)
135
136  if (gribVer.eq.1) then ! GRIB Edition 1
137
138    !print*,'GRiB Edition 1'
139   
140    ! read the grib2 identifiers
141   
142    call grib_get_int(igrib,'indicatorOfParameter',isec1(6),iret)
143    call grib_check(iret,gribFunction,gribErrorMsg)
144    call grib_get_int(igrib,'level',isec1(8),iret)
145    call grib_check(iret,gribFunction,gribErrorMsg)
146
147    !change code for etadot to code for omega
148    if (isec1(6).eq.77) then
149      isec1(6)=135
150    endif
151
152    conversion_factor=1.
153
154  else
155
156    !print*,'GRiB Edition 2'
157
158    ! read the grib2 identifiers
159
160    call grib_get_int(igrib,'discipline',discipl,iret)
161    call grib_check(iret,gribFunction,gribErrorMsg)
162    call grib_get_int(igrib,'parameterCategory',parCat,iret)
163    call grib_check(iret,gribFunction,gribErrorMsg)
164    call grib_get_int(igrib,'parameterNumber',parNum,iret)
165    call grib_check(iret,gribFunction,gribErrorMsg)
166    call grib_get_int(igrib,'typeOfFirstFixedSurface',typSurf,iret)
167    call grib_check(iret,gribFunction,gribErrorMsg)
168    call grib_get_int(igrib,'level',valSurf,iret)
169    call grib_check(iret,gribFunction,gribErrorMsg)
170    call grib_get_int(igrib,'paramId',parId,iret) !added by mc to make it consistent with new readwind.f90
171    call grib_check(iret,gribFunction,gribErrorMsg) !added by mc to make it consistent with new readwind.f90
172
173    !print*,discipl,parCat,parNum,typSurf,valSurf
174
175    !convert to grib1 identifiers
176    isec1(6)=-1
177    isec1(7)=-1
178    isec1(8)=-1
179    isec1(8)=valSurf     ! level
180     conversion_factor=1.
181    if ((parCat.eq.0).and.(parNum.eq.0).and.(typSurf.eq.105)) then ! T
182      isec1(6)=130         ! indicatorOfParameter
183    elseif ((parCat.eq.2).and.(parNum.eq.2).and.(typSurf.eq.105)) then ! U
184      isec1(6)=131         ! indicatorOfParameter
185    elseif ((parCat.eq.2).and.(parNum.eq.3).and.(typSurf.eq.105)) then ! V
186      isec1(6)=132         ! indicatorOfParameter
187    elseif ((parCat.eq.1).and.(parNum.eq.0).and.(typSurf.eq.105)) then ! Q
188      isec1(6)=133         ! indicatorOfParameter
189    elseif ((parCat.eq.3).and.(parNum.eq.0).and.(typSurf.eq.1)) then !SP
190      isec1(6)=134         ! indicatorOfParameter
191    elseif ((parCat.eq.2).and.(parNum.eq.32)) then ! W, actually detadtdpdeta
192      isec1(6)=135         ! indicatorOfParameter
193    elseif ((parCat.eq.128).and.(parNum.eq.77)) then ! W, actually detadtdpdeta
194      isec1(6)=135         ! indicatorOfParameter !added by mc to make it consistent with new readwind.f90
195    elseif ((parCat.eq.2).and.(parNum.eq.8)) then ! Omega, actually detadtdpdeta
196      isec1(6)=135         ! indicatorOfParameter !added by mc to make it consistent with new readwind.f90
197    elseif ((parCat.eq.3).and.(parNum.eq.0).and.(typSurf.eq.101)) then !SLP
198      isec1(6)=151         ! indicatorOfParameter
199    elseif ((parCat.eq.2).and.(parNum.eq.2).and.(typSurf.eq.103)) then ! 10U
200      isec1(6)=165         ! indicatorOfParameter
201    elseif ((parCat.eq.2).and.(parNum.eq.3).and.(typSurf.eq.103)) then ! 10V
202      isec1(6)=166         ! indicatorOfParameter
203    elseif ((parCat.eq.0).and.(parNum.eq.0).and.(typSurf.eq.103)) then ! 2T
204      isec1(6)=167         ! indicatorOfParameter
205    elseif ((parCat.eq.0).and.(parNum.eq.6).and.(typSurf.eq.103)) then ! 2D
206      isec1(6)=168         ! indicatorOfParameter
207    elseif ((parCat.eq.1).and.(parNum.eq.11).and.(typSurf.eq.1)) then ! SD
208      isec1(6)=141         ! indicatorOfParameter
209      conversion_factor=1000. !added by mc to make it consistent with new readwind.f90
210    elseif ((parCat.eq.6).and.(parNum.eq.1) .or. parId .eq. 164) then ! CC !added by mc to make it consistent with new readwind.f90
211      isec1(6)=164         ! indicatorOfParameter
212    elseif ((parCat.eq.1).and.(parNum.eq.9) .or. parId .eq. 142) then ! LSP !added by mc to make it consistent with new readwind.f90
213      isec1(6)=142         ! indicatorOfParameter
214    elseif ((parCat.eq.1).and.(parNum.eq.10)) then ! CP
215      isec1(6)=143         ! indicatorOfParameter
216      conversion_factor=1000. !added by mc to make it consistent with new readwind.f90
217    elseif ((parCat.eq.0).and.(parNum.eq.11).and.(typSurf.eq.1)) then ! SHF
218      isec1(6)=146         ! indicatorOfParameter
219    elseif ((parCat.eq.4).and.(parNum.eq.9).and.(typSurf.eq.1)) then ! SR
220      isec1(6)=176         ! indicatorOfParameter
221    elseif ((parCat.eq.2).and.(parNum.eq.17) .or. parId .eq. 180) then ! EWSS !added by mc to make it consistent with new readwind.f90
222      isec1(6)=180         ! indicatorOfParameter
223    elseif ((parCat.eq.2).and.(parNum.eq.18) .or. parId .eq. 181) then ! NSSS !added by mc to make it consistent with new readwind.f90
224      isec1(6)=181         ! indicatorOfParameter
225    elseif ((parCat.eq.3).and.(parNum.eq.4)) then ! ORO
226      isec1(6)=129         ! indicatorOfParameter
227     elseif ((parCat.eq.3).and.(parNum.eq.7) .or. parId .eq. 160) then ! SDO !added by mc to make it consistent with new readwind.f90
228      isec1(6)=160         ! indicatorOfParameter
229    elseif ((discipl.eq.2).and.(parCat.eq.0).and.(parNum.eq.0).and. &
230         (typSurf.eq.1)) then ! LSM
231      isec1(6)=172         ! indicatorOfParameter
232    else
233      print*,'***WARNING: undefined GRiB2 message found!',discipl, &
234           parCat,parNum,typSurf
235    endif
236    if(parId .ne. isec1(6) .and. parId .ne. 77) then !added by mc to make it consistent with new readwind.f90
237      write(*,*) 'parId',parId, 'isec1(6)',isec1(6)  !
238  !    stop
239    endif
240
241  endif
242
243  !HSO  get the size and data of the values array
244  if (isec1(6).ne.-1) then
245    call grib_get_real4_array(igrib,'values',zsec4,iret)
246    call grib_check(iret,gribFunction,gribErrorMsg)
247  endif
248
249  !HSO  get the required fields from section 2 in a gribex compatible manner
250  if(ifield.eq.1) then
251  call grib_get_int(igrib,'numberOfPointsAlongAParallel', &
252       isec2(2),iret)
253  call grib_check(iret,gribFunction,gribErrorMsg)
254  call grib_get_int(igrib,'numberOfPointsAlongAMeridian', &
255       isec2(3),iret)
256  call grib_check(iret,gribFunction,gribErrorMsg)
257  call grib_get_int(igrib,'numberOfVerticalCoordinateValues', &
258       isec2(12),iret)
259  call grib_check(iret,gribFunction,gribErrorMsg)
260  ! CHECK GRID SPECIFICATIONS
261  if(isec2(2).ne.nxn(l)) stop &
262    'READWIND: NX NOT CONSISTENT FOR A NESTING LEVEL'
263  if(isec2(3).ne.nyn(l)) stop &
264    'READWIND: NY NOT CONSISTENT FOR A NESTING LEVEL'
265  if(isec2(12)/2-1.ne.nlev_ec) stop &
266    'READWIND: VERTICAL DISCRETIZATION NOT CONSISTENT FOR A NESTING LEVEL'
267  endif ! ifield
268
269  !HSO  get the second part of the grid dimensions only from GRiB1 messages
270 if (isec1(6) .eq. 167 .and. (gotGrid.eq.0)) then ! !added by mc to make it consistent with new readwind.f90
271    call grib_get_real8(igrib,'longitudeOfFirstGridPointInDegrees', &
272         xauxin,iret)
273    call grib_check(iret,gribFunction,gribErrorMsg)
274    call grib_get_real8(igrib,'latitudeOfLastGridPointInDegrees', &
275         yauxin,iret)
276    call grib_check(iret,gribFunction,gribErrorMsg)
277    if (xauxin.gt.180.) xauxin=xauxin-360.0
278    if (xauxin.lt.-180.) xauxin=xauxin+360.0
279
280    xaux=xauxin
281    yaux=yauxin
282    if (abs(xaux-xlon0n(l)).gt.eps) &
283    stop 'READWIND: LOWER LEFT LONGITUDE NOT CONSISTENT FOR A NESTING LEVEL'
284    if (abs(yaux-ylat0n(l)).gt.eps) &
285    stop 'READWIND: LOWER LEFT LATITUDE NOT CONSISTENT FOR A NESTING LEVEL'
286    gotGrid=1
287  endif
288
289    do j=0,nyn(l)-1
290      do i=0,nxn(l)-1
291        k=isec1(8)
292        if(isec1(6).eq.130) tthn(i,j,nlev_ec-k+2,n,l)= &!! TEMPERATURE
293             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
294        if(isec1(6).eq.131) uuhn(i,j,nlev_ec-k+2,l)= &!! U VELOCITY
295             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
296        if(isec1(6).eq.132) vvhn(i,j,nlev_ec-k+2,l)= &!! V VELOCITY
297             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
298        if(isec1(6).eq.133) then                         !! SPEC. HUMIDITY
299          qvhn(i,j,nlev_ec-k+2,n,l)=zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
300          if (qvhn(i,j,nlev_ec-k+2,n,l) .lt. 0.) &
301               qvhn(i,j,nlev_ec-k+2,n,l) = 0.
302  !          this is necessary because the gridded data may contain
303  !          spurious negative values
304        endif
305        if(isec1(6).eq.134) psn(i,j,1,n,l)= &!! SURF. PRESS.
306             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
307        if(isec1(6).eq.135) wwhn(i,j,nlev_ec-k+1,l)= &!! W VELOCITY
308             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
309        if(isec1(6).eq.141) sdn(i,j,1,n,l)= &!! SNOW DEPTH
310             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)/conversion_factor !added by mc to make it consistent with new readwind.f90!
311        if(isec1(6).eq.151) msln(i,j,1,n,l)= &!! SEA LEVEL PRESS.
312             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
313        if(isec1(6).eq.164) tccn(i,j,1,n,l)= &!! CLOUD COVER
314             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
315        if(isec1(6).eq.165) u10n(i,j,1,n,l)= &!! 10 M U VELOCITY
316             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
317        if(isec1(6).eq.166) v10n(i,j,1,n,l)= &!! 10 M V VELOCITY
318             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
319        if(isec1(6).eq.167) tt2n(i,j,1,n,l)= &!! 2 M TEMPERATURE
320             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
321        if(isec1(6).eq.168) td2n(i,j,1,n,l)= &!! 2 M DEW POINT
322             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
323        if(isec1(6).eq.142) then                         !! LARGE SCALE PREC.
324          lsprecn(i,j,1,n,l)=zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
325          if (lsprecn(i,j,1,n,l).lt.0.) lsprecn(i,j,1,n,l)=0.
326        endif
327        if(isec1(6).eq.143) then                         !! CONVECTIVE PREC.
328          convprecn(i,j,1,n,l)=zsec4(nxn(l)*(nyn(l)-j-1)+i+1)/conversion_factor !added by mc to make it consistent with new readwind.f90
329          if (convprecn(i,j,1,n,l).lt.0.) convprecn(i,j,1,n,l)=0.
330        endif
331        if(isec1(6).eq.146) sshfn(i,j,1,n,l)= &!! SENS. HEAT FLUX
332             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
333        if((isec1(6).eq.146).and. &
334             (zsec4(nxn(l)*(nyn(l)-j-1)+i+1).ne.0.)) hflswitch=.true.    ! Heat flux available
335        if(isec1(6).eq.176) then                         !! SOLAR RADIATION
336          ssrn(i,j,1,n,l)=zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
337          if (ssrn(i,j,1,n,l).lt.0.) ssrn(i,j,1,n,l)=0.
338        endif
339        if(isec1(6).eq.180) ewss(i,j)= &!! EW SURFACE STRESS
340             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
341        if(isec1(6).eq.181) nsss(i,j)= &!! NS SURFACE STRESS
342             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
343        if(((isec1(6).eq.180).or.(isec1(6).eq.181)).and. &
344             (zsec4(nxn(l)*(nyn(l)-j-1)+i+1).ne.0.)) strswitch=.true.    ! stress available
345        if(isec1(6).eq.129) oron(i,j,l)= &!! ECMWF OROGRAPHY
346             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)/ga
347        if(isec1(6).eq.160) excessoron(i,j,l)= &!! STANDARD DEVIATION OF OROGRAPHY
348             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
349        if(isec1(6).eq.172) lsmn(i,j,l)= &!! ECMWF LAND SEA MASK
350             zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
351        if(isec1(6).eq.131) iumax=max(iumax,nlev_ec-k+1)
352        if(isec1(6).eq.135) iwmax=max(iwmax,nlev_ec-k+1)
353
354      end do
355    end do
356
357  call grib_release(igrib)
358  goto 10                      !! READ NEXT LEVEL OR PARAMETER
359  !
360  ! CLOSING OF INPUT DATA FILE
361  !
36250   call grib_close_file(ifile)
363
364  !error message if no fields found with correct first longitude in it
365  if (gotGrid.eq.0) then
366    print*,'***ERROR: input file needs to contain GRiB1 formatted'// &
367         'messages'
368    stop
369  endif
370
371  if(levdiff2.eq.0) then
372    iwmax=nlev_ec+1
373    do i=0,nxn(l)-1
374      do j=0,nyn(l)-1
375        wwhn(i,j,nlev_ec+1,l)=0.
376      end do
377    end do
378  endif
379
380  do i=0,nxn(l)-1
381    do j=0,nyn(l)-1
382      surfstrn(i,j,1,n,l)=sqrt(ewss(i,j)**2+nsss(i,j)**2)
383    end do
384  end do
385
386  if ((.not.hflswitch).or.(.not.strswitch)) then
387    write(*,*) 'WARNING: No flux data contained in GRIB file ', &
388         wfnamen(l,indj)
389
390  ! CALCULATE USTAR AND SSHF USING THE PROFILE METHOD
391  ! As ECMWF has increased the model resolution, such that now the first model
392  ! level is at about 10 m (where 10-m wind is given), use the 2nd ECMWF level
393  ! (3rd model level in FLEXPART) for the profile method
394  !***************************************************************************
395
396    do i=0,nxn(l)-1
397      do j=0,nyn(l)-1
398        plev1=akz(3)+bkz(3)*psn(i,j,1,n,l)
399        pmean=0.5*(psn(i,j,1,n,l)+plev1)
400        tv=tthn(i,j,3,n,l)*(1.+0.61*qvhn(i,j,3,n,l))
401        fu=-r_air*tv/ga/pmean
402        hlev1=fu*(plev1-psn(i,j,1,n,l))   ! HEIGTH OF FIRST MODEL LAYER
403        ff10m= sqrt(u10n(i,j,1,n,l)**2+v10n(i,j,1,n,l)**2)
404        fflev1=sqrt(uuhn(i,j,3,l)**2+vvhn(i,j,3,l)**2)
405        call pbl_profile(psn(i,j,1,n,l),td2n(i,j,1,n,l),hlev1, &
406             tt2n(i,j,1,n,l),tthn(i,j,3,n,l),ff10m,fflev1, &
407             surfstrn(i,j,1,n,l),sshfn(i,j,1,n,l))
408        if(sshfn(i,j,1,n,l).gt.200.) sshfn(i,j,1,n,l)=200.
409        if(sshfn(i,j,1,n,l).lt.-400.) sshfn(i,j,1,n,l)=-400.
410      end do
411    end do
412  endif
413
414
415  ! Assign 10 m wind to model level at eta=1.0 to have one additional model
416  ! level at the ground
417  ! Specific humidity is taken the same as at one level above
418  ! Temperature is taken as 2 m temperature
419  !**************************************************************************
420
421    do i=0,nxn(l)-1
422      do j=0,nyn(l)-1
423        uuhn(i,j,1,l)=u10n(i,j,1,n,l)
424        vvhn(i,j,1,l)=v10n(i,j,1,n,l)
425        qvhn(i,j,1,n,l)=qvhn(i,j,2,n,l)
426        tthn(i,j,1,n,l)=tt2n(i,j,1,n,l)
427      end do
428    end do
429
430    if(iumax.ne.nuvz-1) stop &
431         'READWIND: NUVZ NOT CONSISTENT FOR A NESTING LEVEL'
432    if(iwmax.ne.nwz) stop &
433         'READWIND: NWZ NOT CONSISTENT FOR A NESTING LEVEL'
434
435  end do
436
437  call verboutput(verbosity,1,'      leaving readwind_nests')
438
439  return
440
441888   write(*,*) ' #### FLEXPART MODEL ERROR! WINDFIELD         #### '
442  write(*,*) ' #### ',wfnamen(l,indj),' FOR NESTING LEVEL  #### '
443  write(*,*) ' #### ',l,' IS NOT GRIB FORMAT !!!           #### '
444  stop 'Execution terminated'
445
446
447999   write(*,*) ' #### FLEXPART MODEL ERROR! WINDFIELD         #### '
448  write(*,*) ' #### ',wfnamen(l,indj),'                    #### '
449  write(*,*) ' #### CANNOT BE OPENED FOR NESTING LEVEL ',l,'####'
450  call verboutput(verbosity,1,'      leaving readwind')
451
452end subroutine readwind_nests
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG