CTBTO FLEXPART WO4 (2015-10-15)
 All Classes Files Functions Variables
readwind_nests.f90
Go to the documentation of this file.
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 
22 subroutine 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, Bernd C. Krueger, Feb. 2001: *
37  ! Variables tthn and qvhn (on eta coordinates) in common block *
38  ! CHANGE: 11/01/2008, Harald Sodemann, GRIB1/2 input with ECMWF grib_api *
39  ! CHANGE: 03/12/2008, Harald Sodemann, update to f90 with ECMWF grib_api *
40  !*****************************************************************************
41  ! Changes Arnold, D. and Morton, D. (2015): *
42  ! -- description of local and common variables *
43  !*****************************************************************************
44 
45  use grib_api
46  use par_mod
47  use com_mod
48 
49  implicit none
50 
51  !***********************************************************************
52  ! Subroutine Parameters: *
53  ! input *
54  ! indj indicates number of the wind field to be read in *
55  ! n temporal index for meteorological fields (1 to 3)*
56  ! uuhn,vvhn, wwhn wind components for the input nest in ECMWF *
57  ! model levels *
58  integer :: indj,n
59  real :: uuhn(0:nxmaxn-1,0:nymaxn-1,nuvzmax,maxnests)
60  real :: vvhn(0:nxmaxn-1,0:nymaxn-1,nuvzmax,maxnests)
61  real :: wwhn(0:nxmaxn-1,0:nymaxn-1,nwzmax,maxnests)
62  !***********************************************************************
63 
64  !***********************************************************************
65  ! Local variables *
66  ! *
67  ! HSO variables for grib_api: *
68  ! ifile grib file to be opened and read in *
69  ! iret is a return code for successful or not open *
70  ! igrib grib edition number (whether GRIB1 or GRIB2) *
71  ! gribVer where info on igrib is kept, 1 for GRIB1 and *
72  ! 2 for GRIB2 *
73  ! parCat parameter category ( = number) , how FLEXPART *
74  ! identifies fields *
75  ! parNum parameter number by product discipline and *
76  ! parameter category *
77  ! typSurf type of first fixed surface *
78  ! valSurf level *
79  ! discipl discipline of processed data contained in a *
80  ! GRIB message *
81  integer :: ifile
82  integer :: iret
83  integer :: igrib
84  integer :: gribver,parcat,parnum,typsurf,valsurf,discipl
85  integer :: gotgrid
86  !
87  ! Variables and arrays for grib decoding: *
88  ! dimension of isec2 at least (22+n), where n is the number of *
89  ! parallels or meridians in a quasi-regular (reduced) Gaussian or *
90  ! lat/long grid *
91  ! dimension of zsec2 at least (10+nn), where nn is the number of *
92  ! vertical coordinate parameters *
93  ! *
94  ! isec1 grib definition section (version, center, ... *
95  ! isec2 grid description section *
96  ! zsec4 the binary data section *
97  ! xaux
98  ! yaux
99  ! xauxin
100  ! yauxin
101  ! xaux0 auxiliary variable for xlon0 *
102  ! yaux0 auxiliary variable for xlat0 *
103  ! nsss NS surface stress *
104  ! ewss EW surface stress *
105  ! plev1 pressure of the first model layer *
106  ! pmean mean pressure between ?? *
107  ! tv virtual temperature *
108  ! fu for conversion from pressure to meters *
109  ! hlev1 height of the first model layer *
110  ! ff10m wind speed at 10 m *
111  ! fflev1 wind speed at the first model layere *
112  ! hflswitch logical variable to check existence of flux data *
113  ! strswitch logical variable to check existence of stress *
114  ! data *
115  integer :: isec1(56),isec2(22+nxmaxn+nymaxn)
116  real(kind=4) :: zsec4(jpunp)
117  real(kind=8) :: xauxin,yauxin
118  real(kind=4) :: xaux,yaux,xaux0,yaux0
119  real :: ewss(0:nxmaxn-1,0:nymaxn-1),nsss(0:nxmaxn-1,0:nymaxn-1)
120  real :: plev1,pmean,tv,fu,hlev1,ff10m,fflev1
121  logical :: hflswitch,strswitch
122  !
123  ! Other variables:
124  ! i,j,k loop control indices in each direction *
125  ! levdiff2 number of model levels - number of model levels *
126  ! for the staggered vertical wind +1 *
127  ! ifield index to control field read in *
128  ! iumax
129  ! iwmax
130  integer ::i,j,k,levdiff2,ifield,iumax,iwmax,l
131  !***********************************************************************
132 
133  !***********************************************************************
134  ! Local constants *
135  !HSO grib api error messages:
136  character(len=24) :: griberrormsg = 'Error reading grib file'
137  character(len=20) :: gribfunction = 'readwind_nests'
138  !***********************************************************************
139  !***********************************************************************
140  ! Global variables *
141  ! from par_mod and com_mod *
142  ! wfname File name of data to be read in *
143  ! maxnests maximum number of nests *
144  ! nxn,nyn,nuvz,nwz expected field dimensions *
145  ! nxmaxn,nymaxn maximum dimension of nested wind fields in *
146  ! x and y direction, respectively *
147  ! nuvzmax,nwzmax maximum dimension of (u,v) and (w) wind fields in z
148  ! direction (for fields on eta levels) *
149  ! nlev_ec number of vertical levels ecmwf model *
150  ! u10n,v10n wind fields at 10 m *
151  ! tt2n,td2n temperature and dew point temperature at 2 m *
152  ! tthn,qvhn temperature and specific humidity on original *
153  ! eta levels *
154  ! psn surface pressure *
155  ! sdn snow depth (but not used afterwards!) *
156  ! msln mean sea level pressure *
157  ! ttcn total cloud cover *
158  ! lsprecn large scale precipitation *
159  ! cprecn convective precipitation *
160  ! sshfn sensible heat flux *
161  ! ssrn solar radiation *
162  ! surfstrn surface stress *
163  ! oron orography *
164  ! excessoron standard deviation of orography *
165  ! lsmn land sea mask *
166  ! r_air individual gas constant for dry air [J/kg/K] *
167  ! ga gravity acceleration of earth [m/s**2] *
168  !***********************************************************************
169 
170 !-----------------------------------------------------------------------------
171 
172  do l=1,numbnests
173  hflswitch=.false.
174  strswitch=.false.
175  levdiff2=nlev_ec-nwz+1
176  iumax=0
177  iwmax=0
178 
179  ifile=0
180  igrib=0
181  iret=0
182 
183  !
184  ! OPENING OF DATA FILE (GRIB CODE)
185  !
186 
187 5 call grib_open_file(ifile,path(numpath+2*(l-1)+1) &
188  (1:length(numpath+2*(l-1)+1))//trim(wfnamen(l,indj)),'r')
189  if (iret.ne.grib_success) then
190  goto 888 ! ERROR DETECTED
191  endif
192  !turn on support for multi fields messages */
193  !call grib_multi_support_on
194 
195  gotgrid=0
196  ifield=0
197 10 ifield=ifield+1
198  !
199  ! GET NEXT FIELDS
200  !
201  call grib_new_from_file(ifile,igrib,iret)
202  if (iret.eq.grib_end_of_file) then
203  goto 50 ! EOF DETECTED
204  elseif (iret.ne.grib_success) then
205  goto 888 ! ERROR DETECTED
206  endif
207 
208  !first see if we read GRIB1 or GRIB2
209  call grib_get_int(igrib,'editionNumber',gribver,iret)
210  call grib_check(iret,gribfunction,griberrormsg)
211 
212  if (gribver.eq.1) then ! GRIB Edition 1
213 
214  !print*,'GRiB Edition 1'
215  !read the grib2 identifiers
216  call grib_get_int(igrib,'indicatorOfParameter',isec1(6),iret)
217  call grib_check(iret,gribfunction,griberrormsg)
218  call grib_get_int(igrib,'level',isec1(8),iret)
219  call grib_check(iret,gribfunction,griberrormsg)
220 
221  !change code for etadot to code for omega
222  if (isec1(6).eq.77) then
223  isec1(6)=135
224  endif
225 
226  else
227 
228  !print*,'GRiB Edition 2'
229  !read the grib2 identifiers
230  call grib_get_int(igrib,'discipline',discipl,iret)
231  call grib_check(iret,gribfunction,griberrormsg)
232  call grib_get_int(igrib,'parameterCategory',parcat,iret)
233  call grib_check(iret,gribfunction,griberrormsg)
234  call grib_get_int(igrib,'parameterNumber',parnum,iret)
235  call grib_check(iret,gribfunction,griberrormsg)
236  call grib_get_int(igrib,'typeOfFirstFixedSurface',typsurf,iret)
237  call grib_check(iret,gribfunction,griberrormsg)
238  call grib_get_int(igrib,'level',valsurf,iret)
239  call grib_check(iret,gribfunction,griberrormsg)
240 
241  !print*,discipl,parCat,parNum,typSurf,valSurf
242 
243  !convert to grib1 identifiers
244  isec1(6)=-1
245  isec1(7)=-1
246  isec1(8)=-1
247  isec1(8)=valsurf ! level
248  if ((parcat.eq.0).and.(parnum.eq.0).and.(typsurf.eq.105)) then ! T
249  isec1(6)=130 ! indicatorOfParameter
250  elseif ((parcat.eq.2).and.(parnum.eq.2).and.(typsurf.eq.105)) then ! U
251  isec1(6)=131 ! indicatorOfParameter
252  elseif ((parcat.eq.2).and.(parnum.eq.3).and.(typsurf.eq.105)) then ! V
253  isec1(6)=132 ! indicatorOfParameter
254  elseif ((parcat.eq.1).and.(parnum.eq.0).and.(typsurf.eq.105)) then ! Q
255  isec1(6)=133 ! indicatorOfParameter
256  elseif ((parcat.eq.3).and.(parnum.eq.0).and.(typsurf.eq.1)) then !SP
257  isec1(6)=134 ! indicatorOfParameter
258  elseif ((parcat.eq.2).and.(parnum.eq.32)) then ! W, actually eta dot
259  isec1(6)=135 ! indicatorOfParameter
260  elseif ((parcat.eq.3).and.(parnum.eq.0).and.(typsurf.eq.101)) then !SLP
261  isec1(6)=151 ! indicatorOfParameter
262  elseif ((parcat.eq.2).and.(parnum.eq.2).and.(typsurf.eq.103)) then ! 10U
263  isec1(6)=165 ! indicatorOfParameter
264  elseif ((parcat.eq.2).and.(parnum.eq.3).and.(typsurf.eq.103)) then ! 10V
265  isec1(6)=166 ! indicatorOfParameter
266  elseif ((parcat.eq.0).and.(parnum.eq.0).and.(typsurf.eq.103)) then ! 2T
267  isec1(6)=167 ! indicatorOfParameter
268  elseif ((parcat.eq.0).and.(parnum.eq.6).and.(typsurf.eq.103)) then ! 2D
269  isec1(6)=168 ! indicatorOfParameter
270  elseif ((parcat.eq.1).and.(parnum.eq.11).and.(typsurf.eq.1)) then ! SD
271  isec1(6)=141 ! indicatorOfParameter
272  elseif ((parcat.eq.6).and.(parnum.eq.1)) then ! CC
273  isec1(6)=164 ! indicatorOfParameter
274  elseif ((parcat.eq.1).and.(parnum.eq.9)) then ! LSP
275  isec1(6)=142 ! indicatorOfParameter
276  elseif ((parcat.eq.1).and.(parnum.eq.10)) then ! CP
277  isec1(6)=143 ! indicatorOfParameter
278  elseif ((parcat.eq.0).and.(parnum.eq.11).and.(typsurf.eq.1)) then ! SHF
279  isec1(6)=146 ! indicatorOfParameter
280  elseif ((parcat.eq.4).and.(parnum.eq.9).and.(typsurf.eq.1)) then ! SR
281  isec1(6)=176 ! indicatorOfParameter
282  elseif ((parcat.eq.2).and.(parnum.eq.17)) then ! EWSS
283  isec1(6)=180 ! indicatorOfParameter
284  elseif ((parcat.eq.2).and.(parnum.eq.18)) then ! NSSS
285  isec1(6)=181 ! indicatorOfParameter
286  elseif ((parcat.eq.3).and.(parnum.eq.4)) then ! ORO
287  isec1(6)=129 ! indicatorOfParameter
288  elseif ((parcat.eq.3).and.(parnum.eq.7)) then ! SDO
289  isec1(6)=160 ! indicatorOfParameter
290  elseif ((discipl.eq.2).and.(parcat.eq.0).and.(parnum.eq.0).and. &
291  (typsurf.eq.1)) then ! LSM
292  isec1(6)=172 ! indicatorOfParameter
293  else
294  print*,'***ERROR: undefined GRiB2 message found!',discipl, &
295  parcat,parnum,typsurf
296  endif
297 
298  endif
299 
300  !HSO get the size and data of the values array
301  if (isec1(6).ne.-1) then
302  call grib_get_real4_array(igrib,'values',zsec4,iret)
303  call grib_check(iret,gribfunction,griberrormsg)
304  endif
305 
306  !HSO get the required fields from section 2 in a gribex compatible manner
307  if(ifield.eq.1) then
308  call grib_get_int(igrib,'numberOfPointsAlongAParallel', &
309  isec2(2),iret)
310  call grib_check(iret,gribfunction,griberrormsg)
311  call grib_get_int(igrib,'numberOfPointsAlongAMeridian', &
312  isec2(3),iret)
313  call grib_check(iret,gribfunction,griberrormsg)
314  call grib_get_int(igrib,'numberOfVerticalCoordinateValues', &
315  isec2(12))
316  call grib_check(iret,gribfunction,griberrormsg)
317  ! CHECK GRID SPECIFICATIONS
318  if(isec2(2).ne.nxn(l)) stop &
319  'READWIND: NX NOT CONSISTENT FOR A NESTING LEVEL'
320  if(isec2(3).ne.nyn(l)) stop &
321  'READWIND: NY NOT CONSISTENT FOR A NESTING LEVEL'
322  if(isec2(12)/2-1.ne.nlev_ec) stop 'READWIND: VERTICAL DISCRET&
323  &IZATION NOT CONSISTENT FOR A NESTING LEVEL'
324  endif ! ifield
325 
326  !HSO get the second part of the grid dimensions only from GRiB1 messages
327  if ((gribver.eq.1).and.(gotgrid.eq.0)) then
328  call grib_get_real8(igrib,'longitudeOfFirstGridPointInDegrees', &
329  xauxin,iret)
330  call grib_check(iret,gribfunction,griberrormsg)
331  call grib_get_real8(igrib,'latitudeOfLastGridPointInDegrees', &
332  yauxin,iret)
333  call grib_check(iret,gribfunction,griberrormsg)
334  xaux=xauxin
335  yaux=yauxin
336  xaux0=xlon0n(l)
337  yaux0=ylat0n(l)
338  if(xaux.lt.0.) xaux=xaux+360.
339  if(yaux.lt.0.) yaux=yaux+360.
340  if(xaux0.lt.0.) xaux0=xaux0+360.
341  if(yaux0.lt.0.) yaux0=yaux0+360.
342  if(xaux.ne.xaux0) &
343  stop 'READWIND: LOWER LEFT LONGITUDE NOT CONSISTENT FOR A NES&
344  &TING LEVEL'
345  if(yaux.ne.yaux0) &
346  stop 'READWIND: LOWER LEFT LATITUDE NOT CONSISTENT FOR A NEST&
347  &ING LEVEL'
348  gotgrid=1
349  endif
350 
351  do j=0,nyn(l)-1
352  do i=0,nxn(l)-1
353  k=isec1(8)
354  if(isec1(6).eq.130) tthn(i,j,nlev_ec-k+2,n,l)= &!! TEMPERATURE
355  zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
356  if(isec1(6).eq.131) uuhn(i,j,nlev_ec-k+2,l)= &!! U VELOCITY
357  zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
358  if(isec1(6).eq.132) vvhn(i,j,nlev_ec-k+2,l)= &!! V VELOCITY
359  zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
360  if(isec1(6).eq.133) then !! SPEC. HUMIDITY
361  qvhn(i,j,nlev_ec-k+2,n,l)=zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
362  if (qvhn(i,j,nlev_ec-k+2,n,l) .lt. 0.) &
363  qvhn(i,j,nlev_ec-k+2,n,l) = 0.
364  ! this is necessary because the gridded data may contain
365  ! spurious negative values
366  endif
367  if(isec1(6).eq.134) psn(i,j,1,n,l)= &!! SURF. PRESS.
368  zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
369  if(isec1(6).eq.135) wwhn(i,j,nlev_ec-k+1,l)= &!! W VELOCITY
370  zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
371  if(isec1(6).eq.141) sdn(i,j,1,n,l)= &!! SNOW DEPTH
372  zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
373  if(isec1(6).eq.151) msln(i,j,1,n,l)= &!! SEA LEVEL PRESS.
374  zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
375  if(isec1(6).eq.164) tccn(i,j,1,n,l)= &!! CLOUD COVER
376  zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
377  if(isec1(6).eq.165) u10n(i,j,1,n,l)= &!! 10 M U VELOCITY
378  zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
379  if(isec1(6).eq.166) v10n(i,j,1,n,l)= &!! 10 M V VELOCITY
380  zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
381  if(isec1(6).eq.167) tt2n(i,j,1,n,l)= &!! 2 M TEMPERATURE
382  zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
383  if(isec1(6).eq.168) td2n(i,j,1,n,l)= &!! 2 M DEW POINT
384  zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
385  if(isec1(6).eq.142) then !! LARGE SCALE PREC.
386  lsprecn(i,j,1,n,l)=zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
387  if (lsprecn(i,j,1,n,l).lt.0.) lsprecn(i,j,1,n,l)=0.
388  endif
389  if(isec1(6).eq.143) then !! CONVECTIVE PREC.
390  convprecn(i,j,1,n,l)=zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
391  if (convprecn(i,j,1,n,l).lt.0.) convprecn(i,j,1,n,l)=0.
392  endif
393  if(isec1(6).eq.146) sshfn(i,j,1,n,l)= &!! SENS. HEAT FLUX
394  zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
395  if((isec1(6).eq.146).and. &
396  (zsec4(nxn(l)*(nyn(l)-j-1)+i+1).ne.0.)) hflswitch=.true. ! Heat flux available
397  if(isec1(6).eq.176) then !! SOLAR RADIATION
398  ssrn(i,j,1,n,l)=zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
399  if (ssrn(i,j,1,n,l).lt.0.) ssrn(i,j,1,n,l)=0.
400  endif
401  if(isec1(6).eq.180) ewss(i,j)= &!! EW SURFACE STRESS
402  zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
403  if(isec1(6).eq.181) nsss(i,j)= &!! NS SURFACE STRESS
404  zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
405  if(((isec1(6).eq.180).or.(isec1(6).eq.181)).and. &
406  (zsec4(nxn(l)*(nyn(l)-j-1)+i+1).ne.0.)) strswitch=.true. ! stress available
407  if(isec1(6).eq.129) oron(i,j,l)= &!! ECMWF OROGRAPHY
408  zsec4(nxn(l)*(nyn(l)-j-1)+i+1)/ga
409  if(isec1(6).eq.160) excessoron(i,j,l)= &!! STANDARD DEVIATION OF OROGRAPHY
410  zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
411  if(isec1(6).eq.172) lsmn(i,j,l)= &!! ECMWF LAND SEA MASK
412  zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
413  if(isec1(6).eq.131) iumax=max(iumax,nlev_ec-k+1)
414  if(isec1(6).eq.135) iwmax=max(iwmax,nlev_ec-k+1)
415 
416  end do
417  end do
418 
419  call grib_release(igrib)
420  goto 10 !! READ NEXT LEVEL OR PARAMETER
421  !
422  ! CLOSING OF INPUT DATA FILE
423  !
424 50 call grib_close_file(ifile)
425 
426  !error message if no fields found with correct first longitude in it
427  if (gotgrid.eq.0) then
428  print*,'***ERROR: input file needs to contain GRiB1 formatted'// &
429  'messages'
430  stop
431  endif
432 
433  if(levdiff2.eq.0) then
434  iwmax=nlev_ec+1
435  do i=0,nxn(l)-1
436  do j=0,nyn(l)-1
437  wwhn(i,j,nlev_ec+1,l)=0.
438  end do
439  end do
440  endif
441 
442  do i=0,nxn(l)-1
443  do j=0,nyn(l)-1
444  surfstrn(i,j,1,n,l)=sqrt(ewss(i,j)**2+nsss(i,j)**2)
445  end do
446  end do
447 
448  if ((.not.hflswitch).or.(.not.strswitch)) then
449  write(*,*) 'WARNING: No flux data contained in GRIB file ', &
450  wfnamen(l,indj)
451 
452  ! CALCULATE USTAR AND SSHF USING THE PROFILE METHOD
453  ! As ECMWF has increased the model resolution, such that now the first model
454  ! level is at about 10 m (where 10-m wind is given), use the 2nd ECMWF level
455  ! (3rd model level in FLEXPART) for the profile method
456  !***************************************************************************
457 
458  do i=0,nxn(l)-1
459  do j=0,nyn(l)-1
460  plev1=akz(3)+bkz(3)*psn(i,j,1,n,l)
461  pmean=0.5*(psn(i,j,1,n,l)+plev1)
462  tv=tthn(i,j,3,n,l)*(1.+0.61*qvhn(i,j,3,n,l))
463  fu=-r_air*tv/ga/pmean
464  hlev1=fu*(plev1-psn(i,j,1,n,l)) ! HEIGTH OF FIRST MODEL LAYER
465  ff10m= sqrt(u10n(i,j,1,n,l)**2+v10n(i,j,1,n,l)**2)
466  fflev1=sqrt(uuhn(i,j,3,l)**2+vvhn(i,j,3,l)**2)
467  call pbl_profile(psn(i,j,1,n,l),td2n(i,j,1,n,l),hlev1, &
468  tt2n(i,j,1,n,l),tthn(i,j,3,n,l),ff10m,fflev1, &
469  surfstrn(i,j,1,n,l),sshfn(i,j,1,n,l))
470  if(sshfn(i,j,1,n,l).gt.200.) sshfn(i,j,1,n,l)=200.
471  if(sshfn(i,j,1,n,l).lt.-400.) sshfn(i,j,1,n,l)=-400.
472  end do
473  end do
474  endif
475 
476 
477  ! Assign 10 m wind to model level at eta=1.0 to have one additional model
478  ! level at the ground
479  ! Specific humidity is taken the same as at one level above
480  ! Temperature is taken as 2 m temperature
481  !**************************************************************************
482 
483  do i=0,nxn(l)-1
484  do j=0,nyn(l)-1
485  uuhn(i,j,1,l)=u10n(i,j,1,n,l)
486  vvhn(i,j,1,l)=v10n(i,j,1,n,l)
487  qvhn(i,j,1,n,l)=qvhn(i,j,2,n,l)
488  tthn(i,j,1,n,l)=tt2n(i,j,1,n,l)
489  end do
490  end do
491 
492  if(iumax.ne.nuvz-1) stop &
493  'READWIND: NUVZ NOT CONSISTENT FOR A NESTING LEVEL'
494  if(iwmax.ne.nwz) stop &
495  'READWIND: NWZ NOT CONSISTENT FOR A NESTING LEVEL'
496 
497  end do
498 
499  return
500 888 write(*,*) ' #### FLEXPART MODEL ERROR! WINDFIELD #### '
501  write(*,*) ' #### ',wfnamen(l,indj),' FOR NESTING LEVEL #### '
502  write(*,*) ' #### ',l,' IS NOT GRIB FORMAT !!! #### '
503  stop 'Execution terminated'
504 
505 
506 999 write(*,*) ' #### FLEXPART MODEL ERROR! WINDFIELD #### '
507  write(*,*) ' #### ',wfnamen(l,indj),' #### '
508  write(*,*) ' #### CANNOT BE OPENED FOR NESTING LEVEL ',l,'####'
509 
510 end subroutine readwind_nests
par_mod
Definition: par_mod.F90:34
pbl_profile
subroutine pbl_profile(ps, td2m, zml1, t2m, tml1, u10m, uml1, stress, hf)
Definition: pbl_profile.f90:22
readwind_nests
subroutine readwind_nests(indj, n, uuhn, vvhn, wwhn)
Definition: readwind_nests.f90:22
com_mod
Definition: com_mod.f90:13
Generated on Thu Oct 15 2015 06:20:58 for CTBTO FLEXPART WO4 (2015-10-15) by   doxygen 1.8.6