49 integer :: gribver,parcat,parnum,typsurf,valsurf,discipl
50 #if defined WITH_CTBTO_PATCHES
55 integer :: i,j,k,l,ifn,ifield,iumax,iwmax,numskip,nlev_ecn
57 real :: akmn(nwzmax),bkmn(nwzmax),akzn(nuvzmax),bkzn(nuvzmax)
58 real(kind=4) :: xaux1,xaux2,yaux1,yaux2
59 real(kind=8) :: xaux1in,xaux2in,yaux1in,yaux2in
69 integer :: isec1(56),isec2(22+nxmaxn+nymaxn)
70 real(kind=4) :: zsec2(60+2*nuvzmax),zsec4(jpunp)
71 #if defined WITH_CTBTO_PATCHES
72 real(kind=4) ::conversion_factor
76 character(len=24) :: griberrormsg =
'Error reading grib file'
77 character(len=20) :: gribfunction =
'gridcheck_nests'
102 5 call grib_open_file(ifile,path(numpath+2*(l-1)+1) &
103 (1:length(numpath+2*(l-1)+1))//trim(wfnamen(l,ifn)),
'r',iret)
104 if (iret.ne.grib_success)
then
117 call grib_new_from_file(ifile,igrib,iret)
118 if (iret.eq.grib_end_of_file)
then
120 elseif (iret.ne.grib_success)
then
124 #if defined WITH_CTBTO_PATCHES
125 conversion_factor=1.0
129 call grib_get_int(igrib,
'editionNumber',gribver,iret)
130 call grib_check(iret,gribfunction,griberrormsg)
132 if (gribver.eq.1)
then
136 call grib_get_int(igrib,
'indicatorOfParameter',isec1(6),iret)
137 call grib_check(iret,gribfunction,griberrormsg)
138 call grib_get_int(igrib,
'level',isec1(8),iret)
139 call grib_check(iret,gribfunction,griberrormsg)
142 #if defined WITH_CTBTO_PATCHES
143 if (isec1(6).eq.77.or.isec1(6).eq.87)
then
145 if (isec1(6).eq.77)
then
155 #if defined WITH_CTBTO_PATCHES
156 call
grib2check(igrib,isec1,conversion_factor)
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)
178 if ((parcat.eq.0).and.(parnum.eq.0).and.(typsurf.eq.105))
then
180 elseif ((parcat.eq.2).and.(parnum.eq.2).and.(typsurf.eq.105))
then
182 elseif ((parcat.eq.2).and.(parnum.eq.3).and.(typsurf.eq.105))
then
184 elseif ((parcat.eq.1).and.(parnum.eq.0).and.(typsurf.eq.105))
then
186 elseif ((parcat.eq.3).and.(parnum.eq.0).and.(typsurf.eq.1))
then
188 elseif ((parcat.eq.2).and.(parnum.eq.32))
then
190 elseif ((parcat.eq.3).and.(parnum.eq.0).and.(typsurf.eq.101))
then
192 elseif ((parcat.eq.2).and.(parnum.eq.2).and.(typsurf.eq.103))
then
194 elseif ((parcat.eq.2).and.(parnum.eq.3).and.(typsurf.eq.103))
then
196 elseif ((parcat.eq.0).and.(parnum.eq.0).and.(typsurf.eq.103))
then
198 elseif ((parcat.eq.0).and.(parnum.eq.6).and.(typsurf.eq.103))
then
200 elseif ((parcat.eq.1).and.(parnum.eq.11).and.(typsurf.eq.1))
then
202 elseif ((parcat.eq.6).and.(parnum.eq.1))
then
204 elseif ((parcat.eq.1).and.(parnum.eq.9))
then
206 elseif ((parcat.eq.1).and.(parnum.eq.10))
then
208 elseif ((parcat.eq.0).and.(parnum.eq.11).and.(typsurf.eq.1))
then
210 elseif ((parcat.eq.4).and.(parnum.eq.9).and.(typsurf.eq.1))
then
212 elseif ((parcat.eq.2).and.(parnum.eq.17))
then
214 elseif ((parcat.eq.2).and.(parnum.eq.18))
then
216 elseif ((parcat.eq.3).and.(parnum.eq.4))
then
218 elseif ((parcat.eq.3).and.(parnum.eq.7))
then
220 elseif ((discipl.eq.2).and.(parcat.eq.0).and.(parnum.eq.0).and. &
224 print*,
'***ERROR: undefined GRiB2 message found!',discipl, &
225 parcat,parnum,typsurf
232 if (isec1(6).ne.-1)
then
233 call grib_get_real4_array(igrib,
'values',zsec4,iret)
234 #if defined WITH_CTBTO_PATCHES
235 zsec4=zsec4/conversion_factor
237 call grib_check(iret,gribfunction,griberrormsg)
241 if (ifield.eq.1)
then
242 call grib_get_int(igrib,
'numberOfPointsAlongAParallel', &
244 call grib_check(iret,gribfunction,griberrormsg)
245 call grib_get_int(igrib,
'numberOfPointsAlongAMeridian', &
247 call grib_check(iret,gribfunction,griberrormsg)
248 call grib_get_int(igrib,
'numberOfVerticalCoordinateValues', &
250 call grib_check(iret,gribfunction,griberrormsg)
252 call grib_get_real4_array(igrib,
'pv',zsec2,iret)
253 call grib_check(iret,gribfunction,griberrormsg)
257 nlev_ecn=isec2(12)/2-1
261 #if defined WITH_CTBTO_PATCHES
262 if (gotgrib.eq.0)
then
264 if ((gribver.eq.1).and.(gotgrib.eq.0))
then
266 call grib_get_real8(igrib,
'longitudeOfFirstGridPointInDegrees', &
268 call grib_check(iret,gribfunction,griberrormsg)
269 call grib_get_real8(igrib,
'longitudeOfLastGridPointInDegrees', &
271 call grib_check(iret,gribfunction,griberrormsg)
272 call grib_get_real8(igrib,
'latitudeOfLastGridPointInDegrees', &
274 call grib_check(iret,gribfunction,griberrormsg)
275 call grib_get_real8(igrib,
'latitudeOfFirstGridPointInDegrees', &
277 call grib_check(iret,gribfunction,griberrormsg)
282 if(xaux1.gt.180) xaux1=xaux1-360.0
283 if(xaux2.gt.180) xaux2=xaux2-360.0
284 if(xaux1.lt.-180) xaux1=xaux1+360.0
285 if(xaux2.lt.-180) xaux2=xaux2+360.0
286 if (xaux2.lt.xaux1) xaux2=xaux2+360.
289 dxn(l)=(xaux2-xaux1)/
real(nxn(l)-1)
290 dyn(l)=(yaux2-yaux1)/
real(nyn(l)-1)
295 if(isec1(6).eq.131) iumax=max(iumax,nlev_ec-k+1)
296 if(isec1(6).eq.135) iwmax=max(iwmax,nlev_ec-k+1)
298 if(isec1(6).eq.129)
then
301 oron(i,j,l)=zsec4(nxn(l)*(nyn(l)-j-1)+i+1)/ga
305 if(isec1(6).eq.172)
then
308 lsmn(i,j,l)=zsec4(nxn(l)*(nyn(l)-j-1)+i+1)/ga
312 if(isec1(6).eq.160)
then
315 excessoron(i,j,l)=zsec4(nxn(l)*(nyn(l)-j-1)+i+1)/ga
320 call grib_release(igrib)
326 30 call grib_close_file(ifile)
329 if (gotgrib.eq.0)
then
330 print*,
'***ERROR: input file needs to contain GRiB1 formatted'// &
337 if(nuvzn.eq.nlev_ec) nwzn=nlev_ecn+1
339 if (nxn(l).gt.nxmaxn)
then
340 write(*,*)
'FLEXPART error: Too many grid points in x direction.'
341 write(*,*)
'Reduce resolution of wind fields (file GRIDSPEC)'
342 write(*,*)
'for nesting level ',l
343 write(*,*)
'Or change parameter settings in file par_mod.'
344 write(*,*) nxn(l),nxmaxn
348 if (nyn(l).gt.nymaxn)
then
349 write(*,*)
'FLEXPART error: Too many grid points in y direction.'
350 write(*,*)
'Reduce resolution of wind fields (file GRIDSPEC)'
351 write(*,*)
'for nesting level ',l
352 write(*,*)
'Or change parameter settings in file par_mod.'
353 write(*,*) nyn(l),nymaxn
357 if ((nuvzn.gt.nuvzmax).or.(nwzn.gt.nwzmax))
then
358 write(*,*)
'FLEXPART error: Nested wind fields have too many'// &
360 write(*,*)
'Problem was encountered for nesting level ',l
368 write(*,
'(a,i2)')
'Nested domain #: ',l
369 write(*,
'(a,f10.2,a1,f10.2,a,f10.2)')
' Longitude range: ', &
370 xlon0n(l),
' to ',xlon0n(l)+(nxn(l)-1)*dxn(l), &
371 ' Grid distance: ',dxn(l)
372 write(*,
'(a,f10.2,a1,f10.2,a,f10.2)')
' Latitude range: ', &
373 ylat0n(l),
' to ',ylat0n(l)+(nyn(l)-1)*dyn(l), &
374 ' Grid distance: ',dyn(l)
390 xaux2=xlon0n(l)+
real(nxn(l)-1)*dxn(l)
392 yaux2=ylat0n(l)+
real(nyn(l)-1)*dyn(l)
394 xln(l)=(xaux1-xlon0)/dx
395 xrn(l)=(xaux2-xlon0)/dx
396 yln(l)=(yaux1-ylat0)/dy
397 yrn(l)=(yaux2-ylat0)/dy
400 if ((xln(l).lt.0.).or.(yln(l).lt.0.).or. &
401 (xrn(l).gt.
real(nxmin1)).or.(yrn(l).gt.
real(nymin1))) then
402 write(*,*)
'Nested domain does not fit into mother domain'
403 write(*,*)
'For global mother domain fields, you can shift'
404 write(*,*)
'shift the mother domain into x-direction'
405 write(*,*)
'by setting nxshift (file par_mod) to a'
406 write(*,*)
'positive value. Execution is terminated.'
414 numskip=nlev_ecn-nuvzn
417 k=nlev_ecn+1+numskip+i
418 akmn(nwzn-i+1)=zsec2(j)
419 bkmn(nwzn-i+1)=zsec2(k)
434 akzn(i+1)=0.5*(akmn(i+1)+akmn(i))
435 bkzn(i+1)=0.5*(bkmn(i+1)+bkmn(i))
445 if ((akzn(i).ne.akz(i)).or.(bkzn(i).ne.bkz(i)))
then
446 write(*,*)
'FLEXPART error: The wind fields of nesting level',l
447 write(*,*)
'are not consistent with the mother domain:'
448 write(*,*)
'Differences in vertical levels detected.'
454 if ((akmn(i).ne.akm(i)).or.(bkmn(i).ne.bkm(i)))
then
455 write(*,*)
'FLEXPART error: The wind fields of nesting level',l
456 write(*,*)
'are not consistent with the mother domain:'
457 write(*,*)
'Differences in vertical levels detected.'
467 write(*,*)
' ###########################################'// &
469 write(*,*)
' FLEXPART SUBROUTINE GRIDCHECK:'
470 write(*,*)
' CAN NOT OPEN INPUT DATA FILE '//wfnamen(l,ifn)
471 write(*,*)
' FOR NESTING LEVEL ',k
472 write(*,*)
' ###########################################'// &
subroutine gridcheck_nests
subroutine grib2check(igrib, isec1, conversion_factor)