Oasis3-MCT
 All Classes Files Functions Variables Macros Pages
mod_oasis_io.F90
Go to the documentation of this file.
1 
2 !> Provides reusable IO routines for OASIS
3 
5 
9  USE mod_oasis_mpi
10  USE mod_oasis_sys
11  USE mod_oasis_ioshr
12  USE mct_mod
13  USE netcdf
14 
15  implicit none
16 
17  private
18 
19  !--- interfaces ---
20  public :: oasis_io_read_avfld
21  public :: oasis_io_read_avfile
22  public :: oasis_io_write_avfile
23  public :: oasis_io_read_array
24  public :: oasis_io_write_array
25  public :: oasis_io_read_avfbf
26  public :: oasis_io_write_avfbf
31 
32 !===========================================================================
33 CONTAINS
34 !===========================================================================
35 
36 !===============================================================================
37 
38 !> Reads single field from a file into an attribute Vector
39 
40 subroutine oasis_io_read_avfld(filename,av,gsmap,mpicom,avfld,filefld,fldtype)
41 
42  implicit none
43 
44  character(len=*), intent(in) :: filename !< filename
45  type(mct_avect) , intent(inout) :: av !< avect
46  type(mct_gsmap) , intent(in) :: gsmap !< gsmap
47  integer(ip_i4_p), intent(in) :: mpicom !< mpicom
48  character(len=*), intent(in) :: avfld !< av field name
49  character(len=*), intent(in) :: filefld !< file field name
50  character(len=*), intent(in),optional :: fldtype !< kind
51 
52  !--- local ---
53  integer(ip_i4_p) :: n,n1,i,j,fk,fk1 ! index
54  integer(ip_i4_p) :: nx ! 2d global size nx
55  integer(ip_i4_p) :: ny ! 2d global size ny
56  type(mct_avect) :: av_g ! avect global data
57  integer(ip_i4_p) :: master_task,iam,ierr ! mpi info
58  integer(ip_i4_p) :: ncid,dimid,dimid2(2),varid ! netcdf info
59  integer(ip_i4_p) :: dlen ! dimension length
60  integer(ip_i4_p) :: status ! error code
61  logical :: exists ! file existance
62  real(ip_double_p),allocatable :: array2(:,:)
63  integer(ip_i4_p) ,allocatable :: array2i(:,:)
64  integer(ip_i4_p) :: ifldtype ! field type int (1) or real (2)
65 
66  character(len=*),parameter :: subname = '(oasis_io_read_avfld)'
67 
68 !-------------------------------------------------------------------------------
69 !
70 !-------------------------------------------------------------------------------
71 
72  IF (mpicom == mpi_comm_null) return
73 
74  ! empty filename, just return
75 
76  call oasis_debug_enter(subname)
77  if (len_trim(filename) < 1) then
78  call oasis_debug_exit(subname)
79  return
80  endif
81 
82  master_task = 0
83  call mpi_comm_rank(mpicom,iam,ierr)
84 
85  ifldtype = 2 ! real default
86  if (present(fldtype)) then
87  ifldtype = 0
88  if (trim(fldtype) == 'int') ifldtype = 1
89  if (trim(fldtype) == 'real') ifldtype = 2
90  if (ifldtype == 0) then
91  WRITE(nulprt,*) subname,estr,'in fldtype argument'
92  call oasis_abort()
93  endif
94  endif
95 
96  call mct_avect_gather(av,av_g,gsmap,master_task,mpicom)
97 
98  if (iam == master_task) then
99 
100  inquire(file=trim(filename),exist=exists)
101  if (exists) then
102  status = nf90_open(trim(filename),nf90_nowrite,ncid)
103  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
104  mpi_rank_local,':',trim(nf90_strerror(status))
105  else
106  write(nulprt,*) subname,estr,'file missing ',trim(filename)
107  call oasis_abort()
108  endif
109 
110  status = nf90_inq_varid(ncid,trim(filefld),varid)
111  if (status /= nf90_noerr) then
112  write(nulprt,*) subname,':',trim(nf90_strerror(status))
113  WRITE(nulprt,*) subname,estr,'filefld variable not found '//trim(filefld)
114  call oasis_abort()
115  endif
116  status = nf90_inquire_variable(ncid,varid,ndims=dlen,dimids=dimid2)
117  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
118  mpi_rank_local,':',trim(nf90_strerror(status))
119  if (dlen /= 2) then
120  write(nulprt,*) subname,estr,'variable ndims ne 2 ',trim(filefld),dlen
121  call oasis_abort()
122  endif
123  status = nf90_inquire_dimension(ncid,dimid2(1),len=nx)
124  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
125  mpi_rank_local,':',trim(nf90_strerror(status))
126  status = nf90_inquire_dimension(ncid,dimid2(2),len=ny)
127  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
128  mpi_rank_local,':',trim(nf90_strerror(status))
129 
130  if (size(av_g%rAttr,dim=2) /= nx*ny) then
131  write(nulprt,*) subname,estr,'av gsize nx ny mismatch ',size(av_g%rAttr,dim=2),nx,ny
132  call oasis_abort()
133  endif
134 
135  if (ifldtype == 1) then
136  allocate(array2i(nx,ny))
137  status = nf90_get_var(ncid,varid,array2i)
138  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
139  mpi_rank_local,':',trim(nf90_strerror(status))
140 
141  n = mct_avect_indexia(av_g,trim(avfld))
142  n1 = 0
143  do j = 1,ny
144  do i = 1,nx
145  n1 = n1 + 1
146  av_g%iAttr(n,n1) = array2i(i,j)
147  enddo
148  enddo
149  deallocate(array2i)
150  else
151  allocate(array2(nx,ny))
152  status = nf90_get_var(ncid,varid,array2)
153  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
154  mpi_rank_local,':',trim(nf90_strerror(status))
155 
156  n = mct_avect_indexra(av_g,trim(avfld))
157  n1 = 0
158  do j = 1,ny
159  do i = 1,nx
160  n1 = n1 + 1
161  av_g%rAttr(n,n1) = array2(i,j)
162  enddo
163  enddo
164  deallocate(array2)
165  endif
166 
167  status = nf90_close(ncid)
168  IF (status /= nf90_noerr) THEN
169  WRITE(nulprt,*) subname,' model :',compid,' proc :',mpi_rank_local,':',&
170  trim(nf90_strerror(status))
171  ENDIF
172 
173  endif
174 
175  call mct_avect_scatter(av_g,av,gsmap,master_task,mpicom)
176  if (iam == master_task) then
177  call mct_avect_clean(av_g)
178  endif
179 
180  call oasis_debug_exit(subname)
181 
182 end subroutine oasis_io_read_avfld
183 
184 !===============================================================================
185 
186 !> Writes all fields from an attribute vector to a file
187 
188 subroutine oasis_io_write_avfile(rstfile,av,gsmap,mpicom,nx,ny,nampre)
189 
190  implicit none
191 
192  character(len=*), intent(in) :: rstfile !< filename
193  type(mct_avect) , intent(in) :: av !< avect
194  type(mct_gsmap) , intent(in) :: gsmap !< gsmap
195  integer(ip_i4_p), intent(in) :: mpicom !< mpicom
196  integer(ip_i4_p), intent(in) :: nx !< 2d global nx size
197  integer(ip_i4_p), intent(in) :: ny !< 2d global ny size
198  character(len=*), intent(in),optional :: nampre !< field name prepend string on file
199 
200  !--- local ---
201  integer(ip_i4_p) :: n,n1,i,j,fk,fk1 ! index
202  integer(ip_i4_p) :: nxf,nyf ! field size on file
203  type(mct_avect) :: av_g ! avect global data
204  type(mct_string) :: mstring ! mct char type
205  character(ic_med) :: itemc ! string converted to char
206  character(ic_med) :: lnampre ! local nampre
207  character(ic_med) :: lstring ! local filename
208  integer(ip_i4_p) :: master_task,iam,ierr ! mpi info
209  integer(ip_i4_p) :: ncid,dimid,dimid2(2),varid ! netcdf info
210  integer(ip_i4_p) :: dlen ! dimension length
211  integer(ip_i4_p) :: status ! error code
212  logical :: exists ! file existance
213  real(ip_double_p),allocatable :: array2(:,:)
214 
215  character(len=*),parameter :: subname = '(oasis_io_write_avfile)'
216 
217 !-------------------------------------------------------------------------------
218 !
219 !-------------------------------------------------------------------------------
220 
221  IF (mpicom == mpi_comm_null) return
222 
223  call oasis_debug_enter(subname)
224 
225  ! empty filename, just return
226 
227  if (len_trim(rstfile) < 1) then
228  call oasis_debug_exit(subname)
229  return
230  endif
231 
232  lnampre = ""
233  if (present(nampre)) then
234  lnampre = trim(nampre)
235  endif
236 
237  master_task = 0
238  call mpi_comm_rank(mpicom,iam,ierr)
239 
240  call mct_avect_gather(av,av_g,gsmap,master_task,mpicom)
241 
242  if (iam == master_task) then
243  if (size(av_g%rAttr,dim=2) /= nx*ny) then
244  write(nulprt,*) subname,estr,'av gsize nx ny mismatch ',&
245  size(av_g%rAttr,dim=2),nx,ny
246  call oasis_abort()
247  endif
248 
249  inquire(file=trim(rstfile),exist=exists)
250  if (exists) then
251  status = nf90_open(trim(rstfile),nf90_write,ncid)
252  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
253  mpi_rank_local,':',trim(nf90_strerror(status))
254  status = nf90_redef(ncid)
255  else
256  status = nf90_create(trim(rstfile),nf90_clobber,ncid)
257  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
258  mpi_rank_local,':',trim(nf90_strerror(status))
259  endif
260 
261  do n = 1,mct_avect_nrattr(av_g)
262  call mct_avect_getrlist(mstring,n,av_g)
263  itemc = mct_string_tochar(mstring)
264  itemc = trim(lnampre)//trim(itemc)
265  call mct_string_clean(mstring)
266 
267  status = nf90_inq_dimid(ncid,trim(itemc)//'_nx',dimid2(1))
268  if (status /= nf90_noerr) then
269  status = nf90_def_dim(ncid,trim(itemc)//'_nx',nx,dimid2(1))
270  endif
271 
272  status = nf90_inq_dimid(ncid,trim(itemc)//'_ny',dimid2(2))
273  if (status /= nf90_noerr) then
274  status = nf90_def_dim(ncid,trim(itemc)//'_ny',ny,dimid2(2))
275  endif
276 
277  status = nf90_inquire_dimension(ncid,dimid2(1),len=dlen)
278  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
279  mpi_rank_local,':',trim(nf90_strerror(status))
280  if (dlen /= nx) then
281  write(nulprt,*) subname,wstr,'dlen ne nx ',dlen,nx
282  CALL oasis_flush(nulprt)
283 ! call oasis_abort()
284  endif
285 
286  status = nf90_inquire_dimension(ncid,dimid2(2),len=dlen)
287  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
288  mpi_rank_local,':',trim(nf90_strerror(status))
289  if (dlen /= ny) then
290  write(nulprt,*) subname,wstr,'dlen ne ny ',dlen,ny
291  CALL oasis_flush(nulprt)
292 ! call oasis_abort()
293  endif
294 
295  status = nf90_inq_varid(ncid,trim(itemc),varid)
296  if (status /= nf90_noerr) then
297  status = nf90_def_var(ncid,trim(itemc),nf90_double,dimid2,varid)
298  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
299  mpi_rank_local,':',trim(nf90_strerror(status))
300  endif
301 
302  enddo
303 
304  status = nf90_enddef(ncid)
305  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
306  mpi_rank_local,':',trim(nf90_strerror(status))
307 
308  nxf = 0
309  nyf = 0
310  do n = 1,mct_avect_nrattr(av_g)
311  call mct_avect_getrlist(mstring,n,av_g)
312  itemc = mct_string_tochar(mstring)
313  itemc = trim(lnampre)//trim(itemc)
314  call mct_string_clean(mstring)
315  status = nf90_inq_varid(ncid,trim(itemc),varid)
316  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
317  mpi_rank_local,':',trim(nf90_strerror(status))
318  if (n == 1) then
319  status = nf90_inquire_variable(ncid,varid,ndims=dlen,dimids=dimid2)
320  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
321  mpi_rank_local,':',trim(nf90_strerror(status))
322  status = nf90_inquire_dimension(ncid,dimid2(1),len=nxf)
323  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
324  mpi_rank_local,':',trim(nf90_strerror(status))
325  status = nf90_inquire_dimension(ncid,dimid2(2),len=nyf)
326  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
327  mpi_rank_local,':',trim(nf90_strerror(status))
328  if (dlen /= 2 .or. nx*ny /= nxf*nyf) then
329  WRITE(nulprt,*) subname,estr,'ndims and size does not match on file'
330  call oasis_abort()
331  endif
332  allocate(array2(nxf,nyf))
333  endif
334 
335  n1 = 0
336  do j = 1,nyf
337  do i = 1,nxf
338  n1 = n1 + 1
339  array2(i,j) = av_g%rAttr(n,n1)
340  enddo
341  enddo
342 
343  status = nf90_put_var(ncid,varid,array2)
344  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
345  mpi_rank_local,':',trim(nf90_strerror(status))
346  enddo
347  deallocate(array2)
348  call mct_avect_clean(av_g)
349 
350  status = nf90_close(ncid)
351  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
352  mpi_rank_local,':',trim(nf90_strerror(status))
353 
354  endif
355 
356  call oasis_debug_exit(subname)
357 
358 end subroutine oasis_io_write_avfile
359 
360 !===============================================================================
361 
362 !> Reads all fields for an attribute vector from a file
363 
364 subroutine oasis_io_read_avfile(rstfile,av,gsmap,mpicom,abort,nampre,didread)
365 
366  implicit none
367 
368  character(len=*), intent(in) :: rstfile !< filename
369  type(mct_avect) , intent(inout) :: av !< avect
370  type(mct_gsmap) , intent(in) :: gsmap !< gsmap
371  integer(ip_i4_p), intent(in) :: mpicom !< mpicom
372  logical , intent(in) ,optional :: abort !< abort on error flag, default is true
373  character(len=*), intent(in) ,optional :: nampre !< name prepend string for fields on file
374  logical , intent(out),optional :: didread !< flag indicating that read was successful
375 
376  !--- local ---
377  integer(ip_i4_p) :: n,n1,i,j,fk,fk1 ! index
378  integer(ip_i4_p) :: nx ! 2d global size nx
379  integer(ip_i4_p) :: ny ! 2d global size ny
380  type(mct_avect) :: av_g ! avect global data
381  type(mct_string) :: mstring ! mct char type
382  character(ic_med) :: itemc ! string converted to char
383  character(ic_med) :: lnampre ! local nampre
384  character(ic_med) :: lstring ! local filename
385  integer(ip_i4_p) :: master_task,iam,ierr ! mpi info
386  integer(ip_i4_p) :: ncid,dimid,dimid2(2),varid ! netcdf info
387  integer(ip_i4_p) :: dlen ! dimension length
388  integer(ip_i4_p) :: status ! error code
389  logical :: exists ! file existance
390  logical :: labort ! local abort flag
391  real(ip_double_p),allocatable :: array2(:,:)
392 
393  character(len=*),parameter :: subname = '(oasis_io_read_avfile)'
394 
395 !-------------------------------------------------------------------------------
396 !
397 !-------------------------------------------------------------------------------
398 
399  IF (mpicom == mpi_comm_null) return
400 
401  call oasis_debug_enter(subname)
402 
403  if (present(didread)) didread = .false.
404 
405  ! empty filename, just return
406 
407  if (len_trim(rstfile) < 1) then
408  call oasis_debug_exit(subname)
409  return
410  endif
411 
412  labort = .true.
413  if (present(abort)) then
414  labort = abort
415  endif
416 
417  lnampre = ""
418  if (present(nampre)) then
419  lnampre = trim(nampre)
420  endif
421 
422  master_task = 0
423  call mpi_comm_rank(mpicom,iam,ierr)
424 
425  call mct_avect_gather(av,av_g,gsmap,master_task,mpicom)
426 
427  if (iam == master_task) then
428 
429  inquire(file=trim(rstfile),exist=exists)
430  if (.not.exists) then
431  IF (labort) THEN
432  write(nulprt,*) subname,estr,'file missing ',trim(rstfile)
433  CALL oasis_abort()
434  ELSE
435  write(nulprt,*) subname,wstr,'file missing ',trim(rstfile)
436  CALL oasis_flush(nulprt)
437  ENDIF
438  else
439  status = nf90_open(trim(rstfile),nf90_nowrite,ncid)
440  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
441  mpi_rank_local,':',trim(nf90_strerror(status))
442 
443  do n = 1,mct_avect_nrattr(av_g)
444  call mct_avect_getrlist(mstring,n,av_g)
445  itemc = mct_string_tochar(mstring)
446  itemc = trim(lnampre)//trim(itemc)
447  call mct_string_clean(mstring)
448 
449  status = nf90_inq_varid(ncid,trim(itemc),varid)
450 
451  if (status /= nf90_noerr) then
452  IF (labort) THEN
453  write(nulprt,*) subname,estr,'var missing on file = ',trim(itemc),':',trim(nf90_strerror(status))
454  CALL oasis_abort()
455 ! ELSE
456 ! write(nulprt,*) subname,wstr,'var missing on file = ',trim(itemc),':',trim(nf90_strerror(status))
457 ! CALL oasis_flush(nulprt)
458  ENDIF
459 
460  else
461  status = nf90_inquire_variable(ncid,varid,ndims=dlen,dimids=dimid2)
462  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
463  mpi_rank_local,':',trim(nf90_strerror(status))
464  if (dlen /= 2) then
465  write(nulprt,*) subname,estr,'variable ndims ne 2 on file ',trim(itemc),dlen
466  call oasis_abort()
467  endif
468  status = nf90_inquire_dimension(ncid,dimid2(1),len=nx)
469  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
470  mpi_rank_local,':',trim(nf90_strerror(status))
471  status = nf90_inquire_dimension(ncid,dimid2(2),len=ny)
472  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
473  mpi_rank_local,':',trim(nf90_strerror(status))
474 
475  if (size(av_g%rAttr,dim=2) /= nx*ny) then
476  write(nulprt,*) subname,estr,'av gsize nx ny mismatch = ',size(av_g%rAttr,dim=2),nx,ny
477  call oasis_abort()
478  endif
479 
480  allocate(array2(nx,ny))
481 
482  status = nf90_get_var(ncid,varid,array2)
483  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
484  mpi_rank_local,':',trim(nf90_strerror(status))
485 
486  n1 = 0
487  do j = 1,ny
488  do i = 1,nx
489  n1 = n1 + 1
490  av_g%rAttr(n,n1) = array2(i,j)
491  enddo
492  enddo
493  if (present(didread)) didread = .true.
494 
495  deallocate(array2)
496  endif ! varid valid
497  enddo
498 
499  status = nf90_close(ncid)
500  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
501  mpi_rank_local,':',trim(nf90_strerror(status))
502 
503  endif ! file exists
504  endif
505 
506  call mct_avect_scatter(av_g,av,gsmap,master_task,mpicom)
507  if (iam == master_task) then
508  call mct_avect_clean(av_g)
509  endif
510 
511  call oasis_debug_exit(subname)
512 
513 end subroutine oasis_io_read_avfile
514 
515 !===============================================================================
516 
517 !> Reads an integer or real field from a file into an array
518 
519 subroutine oasis_io_read_array(rstfile,mpicom,iarray,ivarname,rarray,rvarname,abort)
520 
521  implicit none
522 
523  character(len=*), intent(in) :: rstfile !< filename
524  integer(ip_i4_p), intent(in) :: mpicom !< mpicom
525  integer(ip_i4_p), intent(inout),optional :: iarray(:) !< integer data on root
526  character(len=*), intent(in),optional :: ivarname !< integer variable name on file
527  real(ip_double_p),intent(inout),optional :: rarray(:) !< real data on root
528  character(len=*), intent(in),optional :: rvarname !< real variable name on file
529  logical , intent(in),optional :: abort !< abort on error flag, default is true
530 
531  !--- local ---
532  integer(ip_i4_p) :: ncnt
533  integer(ip_i4_p) :: master_task,iam,ierr ! mpi info
534  integer(ip_i4_p) :: ncid,dimid,dimid1(1),varid ! netcdf info
535  integer(ip_i4_p) :: dlen ! dimension length
536  integer(ip_i4_p) :: status ! error code
537  logical :: exists ! file existance
538  logical :: labort ! local abort flag
539 
540  character(len=*),parameter :: subname = '(oasis_io_read_array)'
541 
542 !-------------------------------------------------------------------------------
543 !
544 !-------------------------------------------------------------------------------
545 
546  if (mpicom == mpi_comm_null) return
547 
548  call oasis_debug_enter(subname)
549 
550  ! empty filename, just return
551 
552  if (len_trim(rstfile) < 1) then
553  call oasis_debug_exit(subname)
554  return
555  endif
556 
557  labort = .true.
558  if (present(abort)) then
559  labort = abort
560  endif
561 
562  master_task = 0
563  call mpi_comm_rank(mpicom,iam,ierr)
564 
565  if (iam == master_task) then
566 
567  inquire(file=trim(rstfile),exist=exists)
568  if (.not.exists) then
569  IF (labort) THEN
570  write(nulprt,*) subname,estr,'file missing ',trim(rstfile)
571  CALL oasis_abort()
572  ELSE
573  write(nulprt,*) subname,wstr,'file missing ',trim(rstfile)
574  CALL oasis_flush(nulprt)
575  ENDIF
576  else
577  status = nf90_open(trim(rstfile),nf90_nowrite,ncid)
578  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
579  mpi_rank_local,':',trim(nf90_strerror(status))
580 
581  if (present(iarray)) then
582  if (.not. present(ivarname)) then
583  write(nulprt,*) subname,estr,'iarray must have ivarname set'
584  call oasis_abort()
585  endif
586 
587  ncnt = size(iarray)
588 
589  status = nf90_inq_varid(ncid,trim(ivarname),varid)
590  if (status /= nf90_noerr) then
591  IF (labort) THEN
592  write(nulprt,*) subname,estr,'var missing on file = ',trim(ivarname),':',trim(nf90_strerror(status))
593  CALL oasis_abort()
594 ! ELSE
595 ! write(nulprt,*) subname,wstr,'var missing on file = ',trim(ivarname),':',trim(nf90_strerror(status))
596 ! CALL oasis_flush(nulprt)
597  ENDIF
598  else
599  status = nf90_inquire_variable(ncid,varid,ndims=dlen,dimids=dimid1)
600  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
601  mpi_rank_local,':',trim(nf90_strerror(status))
602  if (dlen /= 1) then
603  write(nulprt,*) subname,estr,'variable ndims ne 1 ',trim(ivarname),dlen
604  call oasis_abort()
605  endif
606  status = nf90_inquire_dimension(ncid,dimid1(1),len=dlen)
607  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
608  mpi_rank_local,':',trim(nf90_strerror(status))
609 
610  if (ncnt /= dlen) then
611  write(nulprt,*) subname,estr,'iarray ncnt dlen mismatch ',ncnt,dlen
612  call oasis_abort()
613  endif
614 
615  status = nf90_get_var(ncid,varid,iarray)
616  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
617  mpi_rank_local,':',trim(nf90_strerror(status))
618  endif
619  endif
620 
621  if (present(rarray)) then
622  if (.not. present(rvarname)) then
623  write(nulprt,*) subname,estr,'rarray must have rvarname set'
624  call oasis_abort()
625  endif
626 
627  ncnt = size(rarray)
628 
629  status = nf90_inq_varid(ncid,trim(rvarname),varid)
630  if (status /= nf90_noerr) then
631  IF (labort) THEN
632  write(nulprt,*) subname,estr,'var missing on file = ',trim(rvarname),':',trim(nf90_strerror(status))
633  CALL oasis_abort()
634 ! ELSE
635 ! write(nulprt,*) subname,wstr,'var missing on file = ',trim(rvarname),':',trim(nf90_strerror(status))
636 ! CALL oasis_flush(nulprt)
637  ENDIF
638  else
639  status = nf90_inquire_variable(ncid,varid,ndims=dlen,dimids=dimid1)
640  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
641  mpi_rank_local,':',trim(nf90_strerror(status))
642  if (dlen /= 1) then
643  write(nulprt,*) subname,estr,'variable ndims ne 1 ',trim(rvarname),dlen
644  call oasis_abort()
645  endif
646  status = nf90_inquire_dimension(ncid,dimid1(1),len=dlen)
647  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
648  mpi_rank_local,':',trim(nf90_strerror(status))
649 
650  if (ncnt /= dlen) then
651  write(nulprt,*) subname,estr,'rarray ncnt dlen mismatch ',ncnt,dlen
652  call oasis_abort()
653  endif
654 
655  status = nf90_get_var(ncid,varid,rarray)
656  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
657  mpi_rank_local,':',trim(nf90_strerror(status))
658  endif
659  endif
660 
661  status = nf90_close(ncid)
662  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
663  mpi_rank_local,':',trim(nf90_strerror(status))
664 
665  endif
666  endif
667 
668  if (present(iarray)) then
669  call oasis_mpi_bcast(iarray,mpicom,trim(subname)//':iarray',master_task)
670  endif
671 
672  if (present(rarray)) then
673  call oasis_mpi_bcast(rarray,mpicom,trim(subname)//':rarray',master_task)
674  endif
675 
676  call oasis_debug_exit(subname)
677 
678 end subroutine oasis_io_read_array
679 
680 !===============================================================================
681 
682 !> Writes a real or integer array to a file
683 
684 subroutine oasis_io_write_array(rstfile,mpicom,iarray,ivarname,rarray,rvarname)
685 
686  implicit none
687 
688  character(len=*), intent(in) :: rstfile !< filename
689  integer(ip_i4_p), intent(in) :: mpicom !< mpicom
690  integer(ip_i4_p), intent(in),optional :: iarray(:) !< integer data on root
691  character(len=*), intent(in),optional :: ivarname !< integer variable name on file
692  real(ip_double_p),intent(in),optional :: rarray(:) !< real data on root
693  character(len=*), intent(in),optional :: rvarname !< real variable name on file
694 
695  !--- local ---
696  integer(ip_i4_p) :: ncnt
697  integer(ip_i4_p) :: master_task,iam,ierr ! mpi info
698  integer(ip_i4_p) :: ncid,dimid,dimid1(1),varid ! netcdf info
699  integer(ip_i4_p) :: dlen ! dimension length
700  integer(ip_i4_p) :: status ! error code
701  logical :: exists ! file existance
702 
703  character(len=*),parameter :: subname = '(oasis_io_write_array)'
704 
705 !-------------------------------------------------------------------------------
706 !
707 !-------------------------------------------------------------------------------
708 
709  IF (mpicom == mpi_comm_null) return
710 
711  call oasis_debug_enter(subname)
712 
713  ! empty filename, just return
714 
715  if (len_trim(rstfile) < 1) then
716  call oasis_debug_exit(subname)
717  return
718  endif
719 
720  master_task = 0
721  call mpi_comm_rank(mpicom,iam,ierr)
722 
723  if (iam == master_task) then
724 
725  inquire(file=trim(rstfile),exist=exists)
726  if (exists) then
727  status = nf90_open(trim(rstfile),nf90_write,ncid)
728  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
729  mpi_rank_local,':',trim(nf90_strerror(status))
730  status = nf90_redef(ncid)
731  else
732  status = nf90_create(trim(rstfile),nf90_clobber,ncid)
733  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
734  mpi_rank_local,':',trim(nf90_strerror(status))
735  endif
736 
737  if (present(iarray)) then
738  if (.not. present(ivarname)) then
739  write(nulprt,*) subname,estr,'iarray must have ivarname set'
740  call oasis_abort()
741  endif
742 
743  ncnt = size(iarray)
744 
745  status = nf90_inq_dimid(ncid,trim(ivarname)//'_ncnt',dimid1(1))
746  if (status /= nf90_noerr) then
747  status = nf90_def_dim(ncid,trim(ivarname)//'_ncnt',ncnt,dimid1(1))
748  endif
749 
750  status = nf90_inquire_dimension(ncid,dimid1(1),len=dlen)
751  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
752  mpi_rank_local,':',trim(nf90_strerror(status))
753  if (dlen /= ncnt) then
754  write(nulprt,*) subname,estr,'iarray dlen ne ncnt ',dlen,ncnt
755  call oasis_abort()
756  endif
757 
758  status = nf90_inq_varid(ncid,trim(ivarname),varid)
759  if (status /= nf90_noerr) then
760  status = nf90_def_var(ncid,trim(ivarname),nf90_int,dimid1,varid)
761  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
762  mpi_rank_local,':',trim(nf90_strerror(status))
763  endif
764  endif
765 
766  if (present(rarray)) then
767  if (.not. present(rvarname)) then
768  write(nulprt,*) subname,estr,'rarray must have rvarname set'
769  call oasis_abort()
770  endif
771 
772  ncnt = size(rarray)
773 
774  status = nf90_inq_dimid(ncid,trim(rvarname)//'_ncnt',dimid1(1))
775  if (status /= nf90_noerr) then
776  status = nf90_def_dim(ncid,trim(rvarname)//'_ncnt',ncnt,dimid1(1))
777  endif
778 
779  status = nf90_inquire_dimension(ncid,dimid1(1),len=dlen)
780  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
781  mpi_rank_local,':',trim(nf90_strerror(status))
782  if (dlen /= ncnt) then
783  write(nulprt,*) subname,estr,'rarray dlen ne ncnt ',dlen,ncnt
784  call oasis_abort()
785  endif
786 
787  status = nf90_inq_varid(ncid,trim(rvarname),varid)
788  if (status /= nf90_noerr) then
789  status = nf90_def_var(ncid,trim(rvarname),nf90_double,dimid1,varid)
790  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
791  mpi_rank_local,':',trim(nf90_strerror(status))
792  endif
793  endif
794 
795  status = nf90_enddef(ncid)
796  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
797  mpi_rank_local,':',trim(nf90_strerror(status))
798 
799  if (present(iarray)) then
800  status = nf90_inq_varid(ncid,trim(ivarname),varid)
801  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
802  mpi_rank_local,':',trim(nf90_strerror(status))
803  status = nf90_put_var(ncid,varid,iarray)
804  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
805  mpi_rank_local,':',trim(nf90_strerror(status))
806  endif
807 
808  if (present(rarray)) then
809  status = nf90_inq_varid(ncid,trim(rvarname),varid)
810  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
811  mpi_rank_local,':',trim(nf90_strerror(status))
812  status = nf90_put_var(ncid,varid,rarray)
813  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
814  mpi_rank_local,':',trim(nf90_strerror(status))
815  endif
816 
817  status = nf90_close(ncid)
818  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
819  mpi_rank_local,':',trim(nf90_strerror(status))
820 
821  endif
822 
823  call oasis_debug_exit(subname)
824 
825 end subroutine oasis_io_write_array
826 
827 !===============================================================================
828 
829 !> Write each field in an attribute vector to an individual files
830 
831 subroutine oasis_io_write_avfbf(av,gsmap,mpicom,nx,ny,msec,f_string,filename)
832 
833  ! ---------------------------------------
834  ! This works only for a single av to a file
835  ! Optionally can specify time info, and filename info
836  ! ---------------------------------------
837 
838  implicit none
839 
840  type(mct_avect) , intent(in) :: av !< avect
841  type(mct_gsmap) , intent(in) :: gsmap !< gsmap
842  integer(ip_i4_p), intent(in) :: mpicom !< mpicom
843  integer(ip_i4_p), intent(in) :: nx !< 2d global nx size
844  integer(ip_i4_p), intent(in) :: ny !< 2d global ny size
845  integer(ip_i4_p), intent(in),optional :: msec !< optional time info in seconds
846  character(len=*), intent(in),optional :: f_string !< optional f_string to append to filename
847  character(len=*), intent(in),optional :: filename !< optional output filename
848 
849  !--- local ---
850  integer(ip_i4_p) :: n,n1,i,j,fk,fk1 ! index
851  type(mct_avect) :: av_g ! avect global data
852  type(mct_string) :: mstring ! mct char type
853  character(ic_med) :: itemc ! f_string converted to char
854  character(ic_med) :: lfn ! local filename
855  character(ic_med) :: lstring ! local filename
856  integer(ip_i4_p) :: master_task,iam,ierr ! mpi info
857  integer(ip_i4_p) :: ncid,dimid,dimid3(3),varid ! netcdf info
858  integer(ip_i4_p) :: start3(3),count3(3) ! netcdf info
859  integer(ip_i4_p) :: start1(1),count1(1) ! netcdf info
860  integer(ip_i4_p) :: lmsec(1) ! local msec value
861  integer(ip_i4_p) :: dlen ! dimension length
862  integer(ip_i4_p) :: status ! error code
863  logical :: exists ! file existance
864  logical :: whead,wdata ! for writing restart/history cdf files
865  real(ip_double_p),allocatable :: array3(:,:,:)
866  real(ip_double_p) :: tbnds(2)
867 
868  character(len=*),parameter :: subname = '(oasis_io_write_avfbf)'
869 
870 !-------------------------------------------------------------------------------
871 !
872 !-------------------------------------------------------------------------------
873 
874  IF (mpicom == mpi_comm_null) return
875 
876  call oasis_debug_enter(subname)
877 
878  lmsec = 0
879  if (present(msec)) then
880  lmsec = msec
881  endif
882 
883  lstring = " "
884  if (present(f_string)) then
885  lstring = trim(f_string)
886  endif
887 
888  master_task = 0
889  call mpi_comm_rank(mpicom,iam,ierr)
890 
891 #if (PIO_DEFINED)
892 ! tcraig, not working as of Oct 2011
893  call oasis_ioshr_wopen(lfn,clobber=.true.,cdf64=.true.)
894 
895  do fk = fk1,2
896  if (fk == 1) then
897  whead = .true.
898  wdata = .false.
899  elseif (fk == 2) then
900  whead = .false.
901  wdata = .true.
902  else
903  WRITE(nulprt,*) subname,estr,'fk illegal'
904  call oasis_abort()
905  end if
906 
907  call oasis_ioshr_write(lfn,&
908  time_units='seconds',time_cal='seconds',time_val=real(lmsec,ip_double_p),&
909  nt=1,whead=whead,wdata=wdata)
910 
911  call oasis_ioshr_write(lfn, gsmap, av, 'pout', &
912  whead=whead,wdata=wdata,nx=nx,ny=ny,nt=1, &
913  use_float=.false.)
914 
915  if (fk == 1) call oasis_ioshr_enddef(lfn)
916  enddo
917 
918  call oasis_ioshr_close(lfn)
919 #else
920 
921  call mct_avect_gather(av,av_g,gsmap,master_task,mpicom)
922  if (iam == master_task) then
923  if (size(av_g%rAttr,dim=2) /= nx*ny) then
924  write(nulprt,*) subname,estr,'av gsize nx ny mismatch ',size(av_g%rAttr,dim=2),nx,ny
925  call oasis_abort()
926  endif
927 
928  allocate(array3(nx,ny,1))
929  do n = 1,mct_avect_nrattr(av_g)
930  call mct_avect_getrlist(mstring,n,av_g)
931  itemc = mct_string_tochar(mstring)
932  call mct_string_clean(mstring)
933  if (present(filename)) then
934  lfn = trim(filename)
935  else
936  lfn = trim(itemc)//trim(lstring)//'.nc'
937  endif
938  n1 = 0
939  do j = 1,ny
940  do i = 1,nx
941  n1 = n1 + 1
942  array3(i,j,1) = av_g%rAttr(n,n1)
943  enddo
944  enddo
945 
946  start1 = 1
947  count1 = 1
948  start3 = 1
949  count3(1) = nx
950  count3(2) = ny
951  count3(3) = 1
952 
953  inquire(file=trim(lfn),exist=exists)
954  if (exists) then
955  status = nf90_open(lfn,nf90_write,ncid)
956  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
957  mpi_rank_local,':',trim(nf90_strerror(status))
958  status = nf90_inq_dimid(ncid,'time',dimid)
959  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
960  mpi_rank_local,':',trim(nf90_strerror(status))
961  status = nf90_inquire_dimension(ncid,dimid,len=dlen)
962  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
963  mpi_rank_local,':',trim(nf90_strerror(status))
964  start1(1) = dlen + 1
965  start3(3) = start1(1)
966  else
967  status = nf90_create(lfn,nf90_clobber,ncid)
968  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
969  mpi_rank_local,':',trim(nf90_strerror(status))
970  status = nf90_def_dim(ncid,'nx',nx,dimid3(1))
971  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
972  mpi_rank_local,':',trim(nf90_strerror(status))
973  status = nf90_def_dim(ncid,'ny',ny,dimid3(2))
974  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
975  mpi_rank_local,':',trim(nf90_strerror(status))
976  status = nf90_def_dim(ncid,'time',nf90_unlimited,dimid)
977  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
978  mpi_rank_local,':',trim(nf90_strerror(status))
979  dimid3(3) = dimid
980  status = nf90_def_var(ncid,'time',nf90_int,dimid,varid)
981  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
982  mpi_rank_local,':',trim(nf90_strerror(status))
983  status = nf90_def_var(ncid,trim(itemc),nf90_double,dimid3,varid)
984  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
985  mpi_rank_local,':',trim(nf90_strerror(status))
986  status = nf90_enddef(ncid)
987  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
988  mpi_rank_local,':',trim(nf90_strerror(status))
989  endif
990 
991  status = nf90_inq_varid(ncid,'time',varid)
992  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
993  mpi_rank_local,':',trim(nf90_strerror(status))
994  status = nf90_put_var(ncid,varid,lmsec,start1,count1)
995  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
996  mpi_rank_local,':',trim(nf90_strerror(status))
997  status = nf90_inq_varid(ncid,trim(itemc),varid)
998  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
999  mpi_rank_local,':',trim(nf90_strerror(status))
1000  status = nf90_put_var(ncid,varid,array3,start3,count3)
1001  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1002  mpi_rank_local,':',trim(nf90_strerror(status))
1003  status = nf90_close(ncid)
1004  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1005  mpi_rank_local,':',trim(nf90_strerror(status))
1006  enddo
1007  deallocate(array3)
1008  call mct_avect_clean(av_g)
1009  endif
1010 
1011 
1012 #endif
1013 
1014  call oasis_debug_exit(subname)
1015 
1016 end subroutine oasis_io_write_avfbf
1017 
1018 !===============================================================================
1019 
1020 !> Read each field in an attribute vector from individual files
1021 
1022 subroutine oasis_io_read_avfbf(av,gsmap,mpicom,msec,f_string,filename)
1023 
1024  ! ---------------------------------------
1025  ! This works only for a single av from a file
1026  ! Optionally can specify time info, and filename info
1027  ! ---------------------------------------
1028 
1029  implicit none
1030 
1031  type(mct_avect) , intent(inout) :: av !< avect
1032  type(mct_gsmap) , intent(in) :: gsmap !< gsmap
1033  integer(ip_i4_p), intent(in) :: mpicom !< mpicom
1034  integer(ip_i4_p), intent(in),optional :: msec !< optional time info in seconds
1035  character(len=*), intent(in),optional :: f_string !< optional f_string to append to filename
1036  character(len=*), intent(in),optional :: filename !< optional input filename
1037 
1038  !--- local ---
1039  integer(ip_i4_p) :: n,n1,i,j,fk,fk1 ! index
1040  integer(ip_i4_p) :: nx,ny ! grid size from file
1041  type(mct_avect) :: av_g ! avect global data
1042  type(mct_string) :: mstring ! mct char type
1043  character(ic_med) :: itemc ! f_string converted to char
1044  character(ic_med) :: lfn ! local filename
1045  character(ic_med) :: lstring ! local filename
1046  integer(ip_i4_p) :: master_task,iam,ierr ! mpi info
1047  integer(ip_i4_p) :: ncid,dimid,dimid3(3),varid ! netcdf info
1048  integer(ip_i4_p) :: start3(3),count3(3) ! netcdf info
1049  integer(ip_i4_p) :: lmsec(1) ! local msec value
1050  integer(ip_i4_p) :: dlen ! dimension length
1051  integer(ip_i4_p) :: status ! error code
1052  logical :: exists ! file existance
1053  logical :: whead,wdata ! for writing restart/history cdf files
1054  real(ip_double_p),allocatable :: array3(:,:,:)
1055  integer(ip_i4_p) ,allocatable :: time(:)
1056  real(ip_double_p) :: tbnds(2)
1057 
1058  character(len=*),parameter :: subname = '(oasis_io_read_avfbf)'
1059 
1060 !-------------------------------------------------------------------------------
1061 !
1062 !-------------------------------------------------------------------------------
1063 
1064  IF (mpicom == mpi_comm_null) return
1065 
1066  call oasis_debug_enter(subname)
1067 
1068  lmsec = 0
1069  if (present(msec)) then
1070  lmsec = msec
1071  endif
1072 
1073  lstring = " "
1074  if (present(f_string)) then
1075  lstring = trim(f_string)
1076  endif
1077 
1078  master_task = 0
1079  call mpi_comm_rank(mpicom,iam,ierr)
1080 
1081  call mct_avect_gather(av,av_g,gsmap,master_task,mpicom)
1082  if (iam == master_task) then
1083  do n = 1,mct_avect_nrattr(av_g)
1084  call mct_avect_getrlist(mstring,n,av_g)
1085  itemc = mct_string_tochar(mstring)
1086  call mct_string_clean(mstring)
1087  if (present(filename)) then
1088  lfn = trim(filename)
1089  else
1090  lfn = trim(itemc)//trim(lstring)//'.nc'
1091  endif
1092 
1093  inquire(file=trim(lfn),exist=exists)
1094  if (.not.exists) then
1095  write(nulprt,*) subname,estr,'file not found ',trim(lfn)
1096  call oasis_abort()
1097  endif
1098 
1099  status = nf90_open(lfn,nf90_nowrite,ncid)
1100  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1101  mpi_rank_local,':',trim(nf90_strerror(status))
1102 
1103  status = nf90_inq_dimid(ncid,'time',dimid)
1104  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1105  mpi_rank_local,':',trim(nf90_strerror(status))
1106  status = nf90_inquire_dimension(ncid,dimid,len=dlen)
1107  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1108  mpi_rank_local,':',trim(nf90_strerror(status))
1109  allocate(time(dlen))
1110  status = nf90_inq_varid(ncid,'time',varid)
1111  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1112  mpi_rank_local,':',trim(nf90_strerror(status))
1113  status = nf90_get_var(ncid,varid,time)
1114  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1115  mpi_rank_local,':',trim(nf90_strerror(status))
1116  n1 = 0
1117  do j = 1,dlen
1118  if (time(j) == lmsec(1)) n1 = j
1119  enddo
1120  deallocate(time)
1121  if (n1 < 1) then
1122  write(nulprt,*) subname,estr,'time not found on file ',trim(lfn),lmsec
1123  call oasis_abort()
1124  endif
1125 
1126  status = nf90_inq_varid(ncid,trim(itemc),varid)
1127  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1128  mpi_rank_local,':',trim(nf90_strerror(status))
1129  status = nf90_inquire_variable(ncid,varid,dimids=dimid3)
1130  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1131  mpi_rank_local,':',trim(nf90_strerror(status))
1132  status = nf90_inquire_dimension(ncid,dimid3(1),len=nx)
1133  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1134  mpi_rank_local,':',trim(nf90_strerror(status))
1135  status = nf90_inquire_dimension(ncid,dimid3(2),len=ny)
1136  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1137  mpi_rank_local,':',trim(nf90_strerror(status))
1138 
1139  if (size(av_g%rAttr,dim=2) /= nx*ny) then
1140  write(nulprt,*) subname,estr,'av gsize nx ny mismatch ',size(av_g%rAttr,dim=2),nx,ny
1141  call oasis_abort()
1142  endif
1143 
1144  start3 = 1
1145  count3(1) = nx
1146  count3(2) = ny
1147  count3(3) = 1
1148  start3(3) = n1
1149  allocate(array3(nx,ny,1))
1150 
1151  status = nf90_get_var(ncid,varid,array3,start3,count3)
1152  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1153  mpi_rank_local,':',trim(nf90_strerror(status))
1154  status = nf90_close(ncid)
1155  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1156  mpi_rank_local,':',trim(nf90_strerror(status))
1157 
1158  n1 = 0
1159  do j = 1,ny
1160  do i = 1,nx
1161  n1 = n1 + 1
1162  av_g%rAttr(n,n1) = array3(i,j,1)
1163  enddo
1164  enddo
1165 
1166  deallocate(array3)
1167 
1168  enddo
1169  endif
1170 
1171  call mct_avect_scatter(av_g,av,gsmap,master_task,mpicom)
1172  if (iam == master_task) then
1173  call mct_avect_clean(av_g)
1174  endif
1175 
1176  call oasis_debug_exit(subname)
1177 
1178 end subroutine oasis_io_read_avfbf
1179 
1180 !===============================================================================
1181 
1182 !> Read a field on the root task from a file into an array.
1183 
1184 subroutine oasis_io_read_field_fromroot(filename,fldname,ifld2,fld2,fld3,nx,ny,nz)
1185 
1186  ! ---------------------------------------
1187  ! Supports 2d integer or real and 3d real arrays.
1188  ! ---------------------------------------
1189 
1190  implicit none
1191 
1192  character(len=*) , intent(in) :: filename !< filename
1193  character(len=*) , intent(in) :: fldname !< field name
1194  integer(ip_i4_p) , intent(inout),optional :: ifld2(:,:) !< 2d integer array
1195  real(ip_realwp_p), intent(inout),optional :: fld2(:,:) !< 2d real array
1196  real(ip_realwp_p), intent(inout),optional :: fld3(:,:,:) !< 3d real array
1197  integer(ip_i4_p) , intent(inout),optional :: nx !< global nx size
1198  integer(ip_i4_p) , intent(inout),optional :: ny !< global ny size
1199  integer(ip_i4_p) , intent(inout),optional :: nz !< global nz size
1200 
1201  !--- local ---
1202  integer(ip_i4_p) :: ncid,varid ! cdf info
1203  integer(ip_i4_p) :: n,ndims,xtype
1204  integer(ip_i4_p),allocatable :: dimid(:),nd(:)
1205  integer(ip_i4_p) :: status ! error code
1206  integer(ip_i4_p) :: ind ! string index
1207  logical :: exists ! file existance
1208  character(len=ic_med) :: gridname ! grid name derived from fldname
1209 
1210  character(len=*),parameter :: subname = '(oasis_io_read_field_fromroot)'
1211 
1212 !-------------------------------------------------------------------------------
1213 !
1214 !-------------------------------------------------------------------------------
1215 
1216  call oasis_debug_enter(subname)
1217 
1218 ! expects to run only on 1 pe.
1219 ! if (iam == master_task) then
1220 
1221  if (oasis_debug >= 5) then
1222  write(nulprt,*) subname,' read ',trim(filename),' ',trim(fldname)
1223  endif
1224 
1225  inquire(file=trim(filename),exist=exists)
1226  if (exists) then
1227  status = nf90_open(filename,nf90_nowrite,ncid)
1228  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1229  mpi_rank_local,':',trim(nf90_strerror(status))
1230  else
1231  write(nulprt,*) subname,estr,'in filename ',trim(filename)
1232  call oasis_abort()
1233  endif
1234 
1235  status = nf90_inq_varid(ncid,trim(fldname),varid)
1236  if (status /= nf90_noerr) then
1237  write(nulprt,*) subname,estr,'in variable name ',trim(fldname)
1238  call oasis_abort()
1239  endif
1240 
1241  status = nf90_inquire_variable(ncid,varid,ndims=ndims,xtype=xtype)
1242  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1243  mpi_rank_local,':',trim(nf90_strerror(status))
1244 
1245  allocate(dimid(ndims),nd(ndims))
1246 
1247  status = nf90_inquire_variable(ncid,varid,dimids=dimid)
1248  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1249  mpi_rank_local,':',trim(nf90_strerror(status))
1250  do n = 1,ndims
1251  status = nf90_inquire_dimension(ncid,dimid(n),len=nd(n))
1252  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1253  mpi_rank_local,':',trim(nf90_strerror(status))
1254  enddo
1255 
1256  if (present(ifld2) .or. present(fld2) .or. present(fld3)) then
1257  if (xtype == nf90_int .and. ndims == 2 .and. present(ifld2)) then
1258  status = nf90_get_var(ncid,varid,ifld2)
1259  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1260  mpi_rank_local,':',trim(nf90_strerror(status))
1261  elseif (xtype /= nf90_int .and. ndims == 2 .and. present(fld2)) then
1262  status = nf90_get_var(ncid,varid,fld2)
1263  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1264  mpi_rank_local,':',trim(nf90_strerror(status))
1265  elseif (xtype /= nf90_int .and. ndims == 3 .and. present(fld3)) then
1266  status = nf90_get_var(ncid,varid,fld3)
1267  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1268  mpi_rank_local,':',trim(nf90_strerror(status))
1269  else
1270  write(nulprt,*) subname,estr,'mismatch in field and data'
1271  call oasis_abort()
1272  endif
1273  endif
1274 
1275  status = nf90_close(ncid)
1276  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1277  mpi_rank_local,':',trim(nf90_strerror(status))
1278 
1279  if (present(nx)) then
1280  nx = nd(1)
1281  endif
1282 
1283  if (present(ny)) then
1284  ny = nd(2)
1285  endif
1286 
1287  if (present(nz)) then
1288  nz = nd(3)
1289  endif
1290 
1291  deallocate(dimid,nd)
1292 
1293 ! endif
1294 
1295  call oasis_debug_exit(subname)
1296 
1297 end subroutine oasis_io_read_field_fromroot
1298 
1299 !===============================================================================
1300 
1301 !> Write a real array named field from the root task to a file.
1302 
1303 subroutine oasis_io_write_2dgridfld_fromroot(filename,fldname,fld,nx,ny)
1304 
1305  ! ---------------------------------------
1306  ! Designed to work with oasis3 write_grid .
1307  ! ---------------------------------------
1308 
1309  implicit none
1310 
1311  character(len=*), intent(in) :: filename !< file name
1312  character(len=*), intent(in) :: fldname !< field name
1313  real(ip_realwp_p),intent(in) :: fld(:,:) !< 2d real field
1314  integer(ip_i4_p), intent(in) :: nx !< 2d global nx size
1315  integer(ip_i4_p), intent(in) :: ny !< 2d global ny size
1316 
1317  !--- local ---
1318  integer(ip_i4_p) :: ncid,dimid,dimid2(2),varid ! cdf info
1319  integer(ip_i4_p) :: status ! error code
1320  integer(ip_i4_p) :: ind ! string index
1321  logical :: exists ! file existance
1322  character(len=ic_med) :: gridname ! grid name derived from fldname
1323 
1324  character(len=*),parameter :: subname = '(oasis_io_write_2dgridfld_fromroot)'
1325 
1326 !-------------------------------------------------------------------------------
1327 !
1328 !-------------------------------------------------------------------------------
1329 
1330  call oasis_debug_enter(subname)
1331 
1332 ! expects to run only on 1 pe.
1333 ! if (iam == master_task) then
1334 
1335  if (oasis_debug >= 5) then
1336  write(nulprt,*) subname,' write ',trim(filename),' ',trim(fldname),nx,ny
1337  endif
1338 
1339  ind = index(trim(fldname),'.')
1340  if (ind < 2) then
1341  write(nulprt,*) subname,estr,'in fldname ',trim(fldname)
1342  call oasis_abort()
1343  endif
1344  gridname = fldname(1:ind-1)
1345 
1346  inquire(file=trim(filename),exist=exists)
1347  if (exists) then
1348  status = nf90_open(filename,nf90_write,ncid)
1349  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1350  mpi_rank_local,':',trim(nf90_strerror(status))
1351  status = nf90_redef(ncid)
1352  else
1353  status = nf90_create(filename,nf90_clobber,ncid)
1354  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1355  mpi_rank_local,':',trim(nf90_strerror(status))
1356  endif
1357 
1358  ! define x dimension if it doesn't exist
1359  status = nf90_inq_dimid(ncid,'x_'//trim(gridname),dimid2(1))
1360  if (status /= nf90_noerr) then
1361  status = nf90_def_dim(ncid,'x_'//trim(gridname),nx,dimid2(1))
1362  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1363  mpi_rank_local,':',trim(nf90_strerror(status))
1364  endif
1365 
1366  ! define y dimension if it doesn't exist
1367  status = nf90_inq_dimid(ncid,'y_'//trim(gridname),dimid2(2))
1368  if (status /= nf90_noerr) then
1369  status = nf90_def_dim(ncid,'y_'//trim(gridname),ny,dimid2(2))
1370  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1371  mpi_rank_local,':',trim(nf90_strerror(status))
1372  endif
1373 
1374  ! define var if it doesn't exist
1375  status = nf90_inq_varid(ncid,trim(fldname),varid)
1376  if (status /= nf90_noerr) then
1377  status = nf90_def_var(ncid,trim(fldname),nf90_double,dimid2,varid)
1378  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1379  mpi_rank_local,':',trim(nf90_strerror(status))
1380  status = nf90_enddef(ncid)
1381  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1382  mpi_rank_local,':',trim(nf90_strerror(status))
1383  status = nf90_put_var(ncid,varid,fld)
1384  if (status /= nf90_noerr) write(nulprt,*) subname,' model :',compid,' proc :',&
1385  mpi_rank_local,':',trim(nf90_strerror(status))
1386  else
1387  status = nf90_enddef(ncid)
1388  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1389  mpi_rank_local,':',trim(nf90_strerror(status))
1390  endif
1391 
1392  status = nf90_close(ncid)
1393  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1394  mpi_rank_local,':',trim(nf90_strerror(status))
1395 
1396 ! endif
1397 
1398  call oasis_debug_exit(subname)
1399 
1400 end subroutine oasis_io_write_2dgridfld_fromroot
1401 
1402 !===============================================================================
1403 
1404 !> Write an integer array named field from the root task to a file.
1405 
1406 subroutine oasis_io_write_2dgridint_fromroot(filename,fldname,fld,nx,ny)
1407 
1408  ! ---------------------------------------
1409  ! Designed to work with oasis3 write_grid .
1410  ! ---------------------------------------
1411 
1412  implicit none
1413 
1414  character(len=*), intent(in) :: filename !< file name
1415  character(len=*), intent(in) :: fldname !< field name
1416  integer(ip_i4_p), intent(in) :: fld(:,:) !< integer field
1417  integer(ip_i4_p), intent(in) :: nx !< 2d global nx size
1418  integer(ip_i4_p), intent(in) :: ny !< 2d global ny size
1419 
1420  !--- local ---
1421  integer(ip_i4_p) :: ncid,dimid,dimid2(2),varid ! cdf info
1422  integer(ip_i4_p) :: status ! error code
1423  integer(ip_i4_p) :: ind ! string index
1424  logical :: exists ! file existance
1425  character(len=ic_med) :: gridname ! grid name derived from fldname
1426 
1427  character(len=*),parameter :: subname = '(oasis_io_write_2dgridint_fromroot)'
1428 
1429 !-------------------------------------------------------------------------------
1430 !
1431 !-------------------------------------------------------------------------------
1432 
1433  call oasis_debug_enter(subname)
1434 
1435 ! expects to run only on 1 pe.
1436 ! if (iam == master_task) then
1437 
1438  if (oasis_debug >= 5) then
1439  write(nulprt,*) subname,' write ',trim(filename),' ',trim(fldname),nx,ny
1440  endif
1441 
1442  ind = index(trim(fldname),'.')
1443  if (ind < 2) then
1444  write(nulprt,*) subname,estr,'in fldname ',trim(fldname)
1445  call oasis_abort()
1446  endif
1447  gridname = fldname(1:ind-1)
1448 
1449  inquire(file=trim(filename),exist=exists)
1450  if (exists) then
1451  status = nf90_open(filename,nf90_write,ncid)
1452  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1453  mpi_rank_local,':',trim(nf90_strerror(status))
1454  status = nf90_redef(ncid)
1455  else
1456  status = nf90_create(filename,nf90_clobber,ncid)
1457  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1458  mpi_rank_local,':',trim(nf90_strerror(status))
1459  endif
1460 
1461  ! define x dimension if it doesn't exist
1462  status = nf90_inq_dimid(ncid,'x_'//trim(gridname),dimid2(1))
1463  if (status /= nf90_noerr) then
1464  status = nf90_def_dim(ncid,'x_'//trim(gridname),nx,dimid2(1))
1465  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1466  mpi_rank_local,':',trim(nf90_strerror(status))
1467  endif
1468 
1469  ! define y dimension if it doesn't exist
1470  status = nf90_inq_dimid(ncid,'y_'//trim(gridname),dimid2(2))
1471  if (status /= nf90_noerr) then
1472  status = nf90_def_dim(ncid,'y_'//trim(gridname),ny,dimid2(2))
1473  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1474  mpi_rank_local,':',trim(nf90_strerror(status))
1475  endif
1476 
1477  ! define var if it doesn't exist
1478  status = nf90_inq_varid(ncid,trim(fldname),varid)
1479  if (status /= nf90_noerr) then
1480  status = nf90_def_var(ncid,trim(fldname),nf90_int,dimid2,varid)
1481  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1482  mpi_rank_local,':',trim(nf90_strerror(status))
1483  status = nf90_enddef(ncid)
1484  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1485  mpi_rank_local,':',trim(nf90_strerror(status))
1486  status = nf90_put_var(ncid,varid,fld)
1487  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1488  mpi_rank_local,':',trim(nf90_strerror(status))
1489  else
1490  status = nf90_enddef(ncid)
1491  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1492  mpi_rank_local,':',trim(nf90_strerror(status))
1493  endif
1494 
1495  status = nf90_close(ncid)
1496  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1497  mpi_rank_local,':',trim(nf90_strerror(status))
1498 
1499 ! endif
1500 
1501  call oasis_debug_exit(subname)
1502 
1503 end subroutine oasis_io_write_2dgridint_fromroot
1504 
1505 !===============================================================================
1506 
1507 !> Write a 3d real array named field from the root task to a file.
1508 
1509 subroutine oasis_io_write_3dgridfld_fromroot(filename,fldname,fld,nx,ny,nc)
1510 
1511  ! ---------------------------------------
1512  ! Designed to work with oasis3 write_grid.
1513  ! ---------------------------------------
1514 
1515  implicit none
1516 
1517  character(len=*), intent(in) :: filename !< file name
1518  character(len=*), intent(in) :: fldname !< field name
1519  real(ip_realwp_p), intent(in) :: fld(:,:,:)!< 3d real array
1520  integer(ip_i4_p), intent(in) :: nx !< 3d global nx size
1521  integer(ip_i4_p), intent(in) :: ny !< 3d global ny size
1522  integer(ip_i4_p), intent(in) :: nc !< 3d global nz size or nc size for corners
1523 
1524  !--- local ---
1525  integer(ip_i4_p) :: ncid,dimid,dimid3(3),varid ! cdf info
1526  integer(ip_i4_p) :: status ! error code
1527  integer(ip_i4_p) :: ind ! string index
1528  logical :: exists ! file existance
1529  character(len=ic_med) :: gridname ! grid name derived from fldname
1530 
1531  character(len=*),parameter :: subname = '(oasis_io_write_3dgridfld_fromroot)'
1532 
1533 !-------------------------------------------------------------------------------
1534 !
1535 !-------------------------------------------------------------------------------
1536 
1537  call oasis_debug_enter(subname)
1538 
1539 ! expects to run only on 1 pe.
1540 ! if (iam == master_task) then
1541 
1542  if (oasis_debug >= 5) then
1543  write(nulprt,*) subname,' write ',trim(filename),' ',trim(fldname),nx,ny,nc
1544  endif
1545 
1546  ind = index(trim(fldname),'.')
1547  if (ind < 2) then
1548  write(nulprt,*) subname,estr,'in fldname ',trim(fldname)
1549  call oasis_abort()
1550  endif
1551  gridname = fldname(1:ind-1)
1552 
1553  inquire(file=trim(filename),exist=exists)
1554  if (exists) then
1555  status = nf90_open(filename,nf90_write,ncid)
1556  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1557  mpi_rank_local,':',trim(nf90_strerror(status))
1558  status = nf90_redef(ncid)
1559  else
1560  status = nf90_create(filename,nf90_clobber,ncid)
1561  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1562  mpi_rank_local,':',trim(nf90_strerror(status))
1563  endif
1564 
1565  ! define x dimension if it doesn't exist
1566  status = nf90_inq_dimid(ncid,'x_'//trim(gridname),dimid3(1))
1567  if (status /= nf90_noerr) then
1568  status = nf90_def_dim(ncid,'x_'//trim(gridname),nx,dimid3(1))
1569  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1570  mpi_rank_local,':',trim(nf90_strerror(status))
1571  endif
1572 
1573  ! define y dimension if it doesn't exist
1574  status = nf90_inq_dimid(ncid,'y_'//trim(gridname),dimid3(2))
1575  if (status /= nf90_noerr) then
1576  status = nf90_def_dim(ncid,'y_'//trim(gridname),ny,dimid3(2))
1577  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1578  mpi_rank_local,':',trim(nf90_strerror(status))
1579  endif
1580 
1581  ! define crn dimension if it doesn't exist
1582  status = nf90_inq_dimid(ncid,'crn_'//trim(gridname),dimid3(3))
1583  if (status /= nf90_noerr) then
1584  status = nf90_def_dim(ncid,'crn_'//trim(gridname),nc,dimid3(3))
1585  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1586  mpi_rank_local,':',trim(nf90_strerror(status))
1587  endif
1588 
1589  ! define var if it doesn't exist
1590  status = nf90_inq_varid(ncid,trim(fldname),varid)
1591  if (status /= nf90_noerr) then
1592  status = nf90_def_var(ncid,trim(fldname),nf90_double,dimid3,varid)
1593  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1594  mpi_rank_local,':',trim(nf90_strerror(status))
1595  status = nf90_enddef(ncid)
1596  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1597  mpi_rank_local,':',trim(nf90_strerror(status))
1598  status = nf90_put_var(ncid,varid,fld)
1599  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1600  mpi_rank_local,':',trim(nf90_strerror(status))
1601  else
1602  status = nf90_enddef(ncid)
1603  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1604  mpi_rank_local,':',trim(nf90_strerror(status))
1605  endif
1606 
1607  status = nf90_close(ncid)
1608  IF (status /= nf90_noerr) WRITE(nulprt,*) subname,' model :',compid,' proc :',&
1609  mpi_rank_local,':',trim(nf90_strerror(status))
1610 
1611 ! endif
1612 
1613  call oasis_debug_exit(subname)
1614 
1615 end subroutine oasis_io_write_3dgridfld_fromroot
1616 
1617 !-------------------------------------------------------------------
1618 
1619 END MODULE mod_oasis_io
subroutine, public oasis_debug_exit(string)
Used when a subroutine is exited, write info to log file at some debug level.
System type methods.
subroutine, public oasis_io_write_array(rstfile, mpicom, iarray, ivarname, rarray, rvarname)
Writes a real or integer array to a file.
Provides a generic and simpler interface into MPI calls for OASIS.
Generic overloaded interface into MPI broadcast.
subroutine, public oasis_io_read_avfile(rstfile, av, gsmap, mpicom, abort, nampre, didread)
Reads all fields for an attribute vector from a file.
Provides a common location for several OASIS variables.
subroutine, public oasis_io_write_avfbf(av, gsmap, mpicom, nx, ny, msec, f_string, filename)
Write each field in an attribute vector to an individual files.
subroutine, public oasis_io_write_3dgridfld_fromroot(filename, fldname, fld, nx, ny, nc)
Write a 3d real array named field from the root task to a file.
subroutine, public oasis_debug_enter(string)
Used when a subroutine is entered, write info to log file at some debug level.
Defines kinds for OASIS.
subroutine, public oasis_flush(nu)
Flushes output to file.
subroutine, public oasis_io_read_avfld(filename, av, gsmap, mpicom, avfld, filefld, fldtype)
Reads single field from a file into an attribute Vector.
IO interfaces based on pio (not supported yet)
subroutine, public oasis_io_write_avfile(rstfile, av, gsmap, mpicom, nx, ny, nampre)
Writes all fields from an attribute vector to a file.
subroutine, public oasis_abort(id_compid, cd_routine, cd_message)
OASIS abort method, publically available to users.
subroutine, public oasis_io_write_2dgridfld_fromroot(filename, fldname, fld, nx, ny)
Write a real array named field from the root task to a file.
subroutine, public oasis_io_write_2dgridint_fromroot(filename, fldname, fld, nx, ny)
Write an integer array named field from the root task to a file.
subroutine, public oasis_io_read_field_fromroot(filename, fldname, ifld2, fld2, fld3, nx, ny, nz)
Read a field on the root task from a file into an array.
subroutine, public oasis_io_read_array(rstfile, mpicom, iarray, ivarname, rarray, rvarname, abort)
Reads an integer or real field from a file into an array.
Provides reusable IO routines for OASIS.
Definition: mod_oasis_io.F90:4
Defines parameters for OASIS.
subroutine, public oasis_io_read_avfbf(av, gsmap, mpicom, msec, f_string, filename)
Read each field in an attribute vector from individual files.