52#if defined PIO_LIB && defined DISTRIBUTE
61 & ncvname, ncvarid, Nx, Ny, &
62 & Xmin, Xmax, X, Ymin, Ymax, Y, &
70 integer,
intent(in) :: ng, model, ncid, ncvarid
71 integer,
intent(in) :: Nx, Ny
73 character (len=*),
intent(in) :: ncname
74 character (len=*),
intent(in) :: ncvname
76 logical,
intent(out) :: rectangular
78 real(r8),
intent(out) :: Xmin, Xmax, Ymin, Ymax
80 real(r8),
intent(out) :: X(Nx,Ny)
81 real(r8),
intent(out) :: Y(Nx,Ny)
87 integer :: i, ic, j, jc
88 integer :: alen, blank1, blank2, my_varid
90 real(r8),
dimension(Nx) :: Xwrk
91 real(r8),
dimension(Ny) :: Ywrk
93 character (len=20) :: Xname, Yname
95 character (len=*),
parameter :: MyFile = &
96 & __FILE__//
", get_varcoords_nf90"
108 & myvarname = ncvname)
117 IF (trim(
var_aname(i)).eq.
'coordinates')
THEN
120 blank2=index(
var_achar(i)(blank1+1:alen),
' ')
122 IF (blank2.gt.0)
THEN
123 yname=
var_achar(i)(blank1+1:blank1+blank2-1)
131 IF (.not.foundit)
THEN
132 IF (
master)
WRITE (
stdout,10) trim(ncvname), trim(ncname)
141 & myvarname = xname, &
170 IF (x(1,j).eq.x(1,1))
THEN
180 & myvarname = yname, &
209 IF (y(i,1).eq.y(1,1))
THEN
217 IF (((ic.ne.0).and.(ic.eq.nx)).and. &
218 & ((jc.ne.0).and.(jc.eq.ny)))
THEN
230 xmin=min(xmin,x(i,j))
231 xmax=max(xmax,x(i,j))
232 ymin=min(ymin,y(i,j))
233 ymax=max(ymax,y(i,j))
237 10
FORMAT (/,
' GET_VARCOORDS_NF90 - Cannot find "coordinates" ', &
238 &
'attribute for variable:',2x,a,/,22x,
'in file:',2x,a,/, &
239 & /,22x,
'This attribute is needed to interpolate input data', &
240 & /,22x,
'to model grid. Following CF compliance, we need:',/, &
241 & /,22x,
'float my_var(time, lat, lon) ;', &
242 & /,22x,
' my_var:long_name = "my variable long name" ;', &
243 & /,22x,
' my_var:units = "my variable units" ;', &
244 & /,22x,
' my_var:coordinates = "lon lat" ;', &
245 & /,22x,
' my_var:time = "my_var_time" ;',/)
250#if defined PIO_LIB && defined DISTRIBUTE
254 & ncvname, pioVar, Nx, Ny, &
255 & Xmin, Xmax, X, Ymin, Ymax, Y, &
263 logical,
intent(out) :: rectangular
265 integer,
intent(in) :: ng, model
266 integer,
intent(in) :: Nx, Ny
268 real(r8),
intent(out) :: Xmin, Xmax, Ymin, Ymax
270 real(r8),
intent(out) :: X(Nx,Ny)
271 real(r8),
intent(out) :: Y(Nx,Ny)
273 character (len=*),
intent(in) :: ncname
274 character (len=*),
intent(in) :: ncvname
276 TYPE (File_desc_t),
intent(inout) :: pioFile
277 TYPE (My_VarDesc),
intent(inout) :: pioVar
283 integer :: i, ic, j, jc
284 integer :: alen, blank1, blank2
286 real(r8),
dimension(Nx) :: Xwrk
287 real(r8),
dimension(Ny) :: Ywrk
289 character (len=20) :: Xname, Yname
291 character (len=*),
parameter :: MyFile = &
292 & __FILE__//
", get_varcoords_pio"
294 TYPE (Var_desc_t) :: my_pioVar
305 & piofile = piofile, &
306 & myvarname = ncvname)
315 IF (trim(var_aname(i)).eq.
'coordinates')
THEN
316 alen=len_trim(var_achar(i))
317 blank1=index(var_achar(i)(1:alen),
' ')
318 blank2=index(var_achar(i)(blank1+1:alen),
' ')
319 xname=var_achar(i)(1:blank1-1)
320 IF (blank2.gt.0)
THEN
321 yname=var_achar(i)(blank1+1:blank1+blank2-1)
323 yname=var_achar(i)(blank1+1:alen)
329 IF (.not.foundit)
THEN
330 IF (
master)
WRITE (
stdout,10) trim(ncvname), trim(ncname)
338 & piofile = piofile, &
339 & myvarname = xname, &
340 & piovar = my_piovar)
343 IF (n_vdim.eq.1)
THEN
346 & piofile = piofile, &
361 & piofile = piofile, &
368 IF (x(1,j).eq.x(1,1))
THEN
377 & piofile = piofile, &
378 & myvarname = yname, &
379 & piovar = my_piovar)
382 IF (n_vdim.eq.1)
THEN
385 & piofile = piofile, &
400 & piofile = piofile, &
407 IF (y(i,1).eq.y(1,1))
THEN
415 IF (((ic.ne.0).and.(ic.eq.nx)).and. &
416 & ((jc.ne.0).and.(jc.eq.ny)))
THEN
428 xmin=min(xmin,x(i,j))
429 xmax=max(xmax,x(i,j))
430 ymin=min(ymin,y(i,j))
431 ymax=max(ymax,y(i,j))
435 10
FORMAT (/,
' GET_VARCOORDS_PIO - Cannot find "coordinates" ', &
436 &
'attribute for variable:',2x,a,/,21x,
'in file:',2x,a,/, &
437 & /,21x,
'This attribute is needed to interpolate input data', &
438 & /,21x,
'to model grid. Following CF compliance, we need:',/, &
439 & /,21x,
'float my_var(time, lat, lon) ;', &
440 & /,21x,
' my_var:long_name = "my variable long name" ;', &
441 & /,21x,
' my_var:units = "my variable units" ;', &
442 & /,21x,
' my_var:coordinates = "lon lat" ;', &
443 & /,21x,
' my_var:time = "my_var_time" ;',/)
subroutine get_varcoords_pio(ng, model, ncname, piofile, ncvname, piovar, nx, ny, xmin, xmax, x, ymin, ymax, y, rectangular)
subroutine get_varcoords_nf90(ng, model, ncname, ncid, ncvname, ncvarid, nx, ny, xmin, xmax, x, ymin, ymax, y, rectangular)
character(len=256) sourcefile
character(len=1024), dimension(nvara) var_achar
character(len=100), dimension(nvara) var_aname
subroutine, public netcdf_inq_var(ng, model, ncname, ncid, myvarname, searchvar, varid, nvardim, nvaratt)
subroutine, public pio_netcdf_inq_var(ng, model, ncname, piofile, myvarname, searchvar, piovar, nvardim, nvaratt)
real(dp), parameter spval
logical function, public founderror(flag, noerr, line, routine)