44 character(len=*),
intent(in) :: filename
45 type(mct_avect) ,
intent(inout) :: av
46 type(mct_gsmap) ,
intent(in) :: gsmap
47 integer(ip_i4_p),
intent(in) :: mpicom
48 character(len=*),
intent(in) :: avfld
49 character(len=*),
intent(in) :: filefld
50 character(len=*),
intent(in),
optional :: fldtype
53 integer(ip_i4_p) :: n,n1,i,j,fk,fk1
54 integer(ip_i4_p) :: nx
55 integer(ip_i4_p) :: ny
56 type(mct_avect) :: av_g
57 integer(ip_i4_p) :: master_task,iam,ierr
58 integer(ip_i4_p) :: ncid,dimid,dimid2(2),varid
59 integer(ip_i4_p) :: dlen
60 integer(ip_i4_p) :: status
62 real(ip_double_p),
allocatable :: array2(:,:)
63 integer(ip_i4_p) ,
allocatable :: array2i(:,:)
64 integer(ip_i4_p) :: ifldtype
66 character(len=*),
parameter :: subname =
'(oasis_io_read_avfld)'
72 IF (mpicom == mpi_comm_null)
return
77 if (len_trim(filename) < 1)
then
83 call mpi_comm_rank(mpicom,iam,ierr)
86 if (
present(fldtype))
then
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'
96 call mct_avect_gather(av,av_g,gsmap,master_task,mpicom)
98 if (iam == master_task)
then
100 inquire(file=trim(filename),exist=exists)
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))
106 write(nulprt,*) subname,estr,
'file missing ',trim(filename)
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)
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))
120 write(nulprt,*) subname,estr,
'variable ndims ne 2 ',trim(filefld),dlen
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))
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
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))
141 n = mct_avect_indexia(av_g,trim(avfld))
146 av_g%iAttr(n,n1) = array2i(i,j)
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))
156 n = mct_avect_indexra(av_g,trim(avfld))
161 av_g%rAttr(n,n1) = array2(i,j)
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))
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)
192 character(len=*),
intent(in) :: rstfile
193 type(mct_avect) ,
intent(in) :: av
194 type(mct_gsmap) ,
intent(in) :: gsmap
195 integer(ip_i4_p),
intent(in) :: mpicom
196 integer(ip_i4_p),
intent(in) :: nx
197 integer(ip_i4_p),
intent(in) :: ny
198 character(len=*),
intent(in),
optional :: nampre
201 integer(ip_i4_p) :: n,n1,i,j,fk,fk1
202 integer(ip_i4_p) :: nxf,nyf
203 type(mct_avect) :: av_g
204 type(mct_string) :: mstring
205 character(ic_med) :: itemc
206 character(ic_med) :: lnampre
207 character(ic_med) :: lstring
208 integer(ip_i4_p) :: master_task,iam,ierr
209 integer(ip_i4_p) :: ncid,dimid,dimid2(2),varid
210 integer(ip_i4_p) :: dlen
211 integer(ip_i4_p) :: status
213 real(ip_double_p),
allocatable :: array2(:,:)
215 character(len=*),
parameter :: subname =
'(oasis_io_write_avfile)'
221 IF (mpicom == mpi_comm_null)
return
227 if (len_trim(rstfile) < 1)
then
233 if (
present(nampre))
then
234 lnampre = trim(nampre)
238 call mpi_comm_rank(mpicom,iam,ierr)
240 call mct_avect_gather(av,av_g,gsmap,master_task,mpicom)
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
249 inquire(file=trim(rstfile),exist=exists)
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)
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))
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)
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))
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))
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))
281 write(nulprt,*) subname,wstr,
'dlen ne nx ',dlen,nx
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))
290 write(nulprt,*) subname,wstr,
'dlen ne ny ',dlen,ny
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))
304 status = nf90_enddef(ncid)
305 IF (status /= nf90_noerr)
WRITE(nulprt,*) subname,
' model :',compid,
' proc :',&
306 mpi_rank_local,
':',trim(nf90_strerror(status))
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))
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'
332 allocate(array2(nxf,nyf))
339 array2(i,j) = av_g%rAttr(n,n1)
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))
348 call mct_avect_clean(av_g)
350 status = nf90_close(ncid)
351 IF (status /= nf90_noerr)
WRITE(nulprt,*) subname,
' model :',compid,
' proc :',&
352 mpi_rank_local,
':',trim(nf90_strerror(status))
368 character(len=*),
intent(in) :: rstfile
369 type(mct_avect) ,
intent(inout) :: av
370 type(mct_gsmap) ,
intent(in) :: gsmap
371 integer(ip_i4_p),
intent(in) :: mpicom
372 logical ,
intent(in) ,
optional :: abort
373 character(len=*),
intent(in) ,
optional :: nampre
374 logical ,
intent(out),
optional :: didread
377 integer(ip_i4_p) :: n,n1,i,j,fk,fk1
378 integer(ip_i4_p) :: nx
379 integer(ip_i4_p) :: ny
380 type(mct_avect) :: av_g
381 type(mct_string) :: mstring
382 character(ic_med) :: itemc
383 character(ic_med) :: lnampre
384 character(ic_med) :: lstring
385 integer(ip_i4_p) :: master_task,iam,ierr
386 integer(ip_i4_p) :: ncid,dimid,dimid2(2),varid
387 integer(ip_i4_p) :: dlen
388 integer(ip_i4_p) :: status
391 real(ip_double_p),
allocatable :: array2(:,:)
393 character(len=*),
parameter :: subname =
'(oasis_io_read_avfile)'
399 IF (mpicom == mpi_comm_null)
return
403 if (
present(didread)) didread = .false.
407 if (len_trim(rstfile) < 1)
then
413 if (
present(abort))
then
418 if (
present(nampre))
then
419 lnampre = trim(nampre)
423 call mpi_comm_rank(mpicom,iam,ierr)
425 call mct_avect_gather(av,av_g,gsmap,master_task,mpicom)
427 if (iam == master_task)
then
429 inquire(file=trim(rstfile),exist=exists)
430 if (.not.exists)
then
432 write(nulprt,*) subname,estr,
'file missing ',trim(rstfile)
435 write(nulprt,*) subname,wstr,
'file missing ',trim(rstfile)
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))
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)
449 status = nf90_inq_varid(ncid,trim(itemc),varid)
451 if (status /= nf90_noerr)
then
453 write(nulprt,*) subname,estr,
'var missing on file = ',trim(itemc),
':',trim(nf90_strerror(status))
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))
465 write(nulprt,*) subname,estr,
'variable ndims ne 2 on file ',trim(itemc),dlen
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))
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
480 allocate(array2(nx,ny))
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))
490 av_g%rAttr(n,n1) = array2(i,j)
493 if (
present(didread)) didread = .true.
499 status = nf90_close(ncid)
500 IF (status /= nf90_noerr)
WRITE(nulprt,*) subname,
' model :',compid,
' proc :',&
501 mpi_rank_local,
':',trim(nf90_strerror(status))
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)
523 character(len=*),
intent(in) :: rstfile
524 integer(ip_i4_p),
intent(in) :: mpicom
525 integer(ip_i4_p),
intent(inout),
optional :: iarray(:)
526 character(len=*),
intent(in),
optional :: ivarname
527 real(ip_double_p),
intent(inout),
optional :: rarray(:)
528 character(len=*),
intent(in),
optional :: rvarname
529 logical ,
intent(in),
optional :: abort
532 integer(ip_i4_p) :: ncnt
533 integer(ip_i4_p) :: master_task,iam,ierr
534 integer(ip_i4_p) :: ncid,dimid,dimid1(1),varid
535 integer(ip_i4_p) :: dlen
536 integer(ip_i4_p) :: status
540 character(len=*),
parameter :: subname =
'(oasis_io_read_array)'
546 if (mpicom == mpi_comm_null)
return
552 if (len_trim(rstfile) < 1)
then
558 if (
present(abort))
then
563 call mpi_comm_rank(mpicom,iam,ierr)
565 if (iam == master_task)
then
567 inquire(file=trim(rstfile),exist=exists)
568 if (.not.exists)
then
570 write(nulprt,*) subname,estr,
'file missing ',trim(rstfile)
573 write(nulprt,*) subname,wstr,
'file missing ',trim(rstfile)
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))
581 if (
present(iarray))
then
582 if (.not.
present(ivarname))
then
583 write(nulprt,*) subname,estr,
'iarray must have ivarname set'
589 status = nf90_inq_varid(ncid,trim(ivarname),varid)
590 if (status /= nf90_noerr)
then
592 write(nulprt,*) subname,estr,
'var missing on file = ',trim(ivarname),
':',trim(nf90_strerror(status))
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))
603 write(nulprt,*) subname,estr,
'variable ndims ne 1 ',trim(ivarname),dlen
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))
610 if (ncnt /= dlen)
then
611 write(nulprt,*) subname,estr,
'iarray ncnt dlen mismatch ',ncnt,dlen
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))
621 if (
present(rarray))
then
622 if (.not.
present(rvarname))
then
623 write(nulprt,*) subname,estr,
'rarray must have rvarname set'
629 status = nf90_inq_varid(ncid,trim(rvarname),varid)
630 if (status /= nf90_noerr)
then
632 write(nulprt,*) subname,estr,
'var missing on file = ',trim(rvarname),
':',trim(nf90_strerror(status))
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))
643 write(nulprt,*) subname,estr,
'variable ndims ne 1 ',trim(rvarname),dlen
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))
650 if (ncnt /= dlen)
then
651 write(nulprt,*) subname,estr,
'rarray ncnt dlen mismatch ',ncnt,dlen
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))
661 status = nf90_close(ncid)
662 IF (status /= nf90_noerr)
WRITE(nulprt,*) subname,
' model :',compid,
' proc :',&
663 mpi_rank_local,
':',trim(nf90_strerror(status))
668 if (
present(iarray))
then
669 call
oasis_mpi_bcast(iarray,mpicom,trim(subname)//
':iarray',master_task)
672 if (
present(rarray))
then
673 call
oasis_mpi_bcast(rarray,mpicom,trim(subname)//
':rarray',master_task)
688 character(len=*),
intent(in) :: rstfile
689 integer(ip_i4_p),
intent(in) :: mpicom
690 integer(ip_i4_p),
intent(in),
optional :: iarray(:)
691 character(len=*),
intent(in),
optional :: ivarname
692 real(ip_double_p),
intent(in),
optional :: rarray(:)
693 character(len=*),
intent(in),
optional :: rvarname
696 integer(ip_i4_p) :: ncnt
697 integer(ip_i4_p) :: master_task,iam,ierr
698 integer(ip_i4_p) :: ncid,dimid,dimid1(1),varid
699 integer(ip_i4_p) :: dlen
700 integer(ip_i4_p) :: status
703 character(len=*),
parameter :: subname =
'(oasis_io_write_array)'
709 IF (mpicom == mpi_comm_null)
return
715 if (len_trim(rstfile) < 1)
then
721 call mpi_comm_rank(mpicom,iam,ierr)
723 if (iam == master_task)
then
725 inquire(file=trim(rstfile),exist=exists)
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)
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))
737 if (
present(iarray))
then
738 if (.not.
present(ivarname))
then
739 write(nulprt,*) subname,estr,
'iarray must have ivarname set'
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))
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
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))
766 if (
present(rarray))
then
767 if (.not.
present(rvarname))
then
768 write(nulprt,*) subname,estr,
'rarray must have rvarname set'
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))
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
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))
795 status = nf90_enddef(ncid)
796 IF (status /= nf90_noerr)
WRITE(nulprt,*) subname,
' model :',compid,
' proc :',&
797 mpi_rank_local,
':',trim(nf90_strerror(status))
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))
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))
817 status = nf90_close(ncid)
818 IF (status /= nf90_noerr)
WRITE(nulprt,*) subname,
' model :',compid,
' proc :',&
819 mpi_rank_local,
':',trim(nf90_strerror(status))
840 type(mct_avect) ,
intent(in) :: av
841 type(mct_gsmap) ,
intent(in) :: gsmap
842 integer(ip_i4_p),
intent(in) :: mpicom
843 integer(ip_i4_p),
intent(in) :: nx
844 integer(ip_i4_p),
intent(in) :: ny
845 integer(ip_i4_p),
intent(in),
optional :: msec
846 character(len=*),
intent(in),
optional :: f_string
847 character(len=*),
intent(in),
optional :: filename
850 integer(ip_i4_p) :: n,n1,i,j,fk,fk1
851 type(mct_avect) :: av_g
852 type(mct_string) :: mstring
853 character(ic_med) :: itemc
854 character(ic_med) :: lfn
855 character(ic_med) :: lstring
856 integer(ip_i4_p) :: master_task,iam,ierr
857 integer(ip_i4_p) :: ncid,dimid,dimid3(3),varid
858 integer(ip_i4_p) :: start3(3),count3(3)
859 integer(ip_i4_p) :: start1(1),count1(1)
860 integer(ip_i4_p) :: lmsec(1)
861 integer(ip_i4_p) :: dlen
862 integer(ip_i4_p) :: status
864 logical :: whead,wdata
865 real(ip_double_p),
allocatable :: array3(:,:,:)
866 real(ip_double_p) :: tbnds(2)
868 character(len=*),
parameter :: subname =
'(oasis_io_write_avfbf)'
874 IF (mpicom == mpi_comm_null)
return
879 if (
present(msec))
then
884 if (
present(f_string))
then
885 lstring = trim(f_string)
889 call mpi_comm_rank(mpicom,iam,ierr)
893 call oasis_ioshr_wopen(lfn,clobber=.true.,cdf64=.true.)
899 elseif (fk == 2)
then
903 WRITE(nulprt,*) subname,estr,
'fk illegal'
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)
911 call oasis_ioshr_write(lfn, gsmap, av,
'pout', &
912 whead=whead,wdata=wdata,nx=nx,ny=ny,nt=1, &
915 if (fk == 1) call oasis_ioshr_enddef(lfn)
918 call oasis_ioshr_close(lfn)
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
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
936 lfn = trim(itemc)//trim(lstring)//
'.nc'
942 array3(i,j,1) = av_g%rAttr(n,n1)
953 inquire(file=trim(lfn),exist=exists)
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))
965 start3(3) = start1(1)
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))
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))
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))
1008 call mct_avect_clean(av_g)
1031 type(mct_avect) ,
intent(inout) :: av
1032 type(mct_gsmap) ,
intent(in) :: gsmap
1033 integer(ip_i4_p),
intent(in) :: mpicom
1034 integer(ip_i4_p),
intent(in),
optional :: msec
1035 character(len=*),
intent(in),
optional :: f_string
1036 character(len=*),
intent(in),
optional :: filename
1039 integer(ip_i4_p) :: n,n1,i,j,fk,fk1
1040 integer(ip_i4_p) :: nx,ny
1041 type(mct_avect) :: av_g
1042 type(mct_string) :: mstring
1043 character(ic_med) :: itemc
1044 character(ic_med) :: lfn
1045 character(ic_med) :: lstring
1046 integer(ip_i4_p) :: master_task,iam,ierr
1047 integer(ip_i4_p) :: ncid,dimid,dimid3(3),varid
1048 integer(ip_i4_p) :: start3(3),count3(3)
1049 integer(ip_i4_p) :: lmsec(1)
1050 integer(ip_i4_p) :: dlen
1051 integer(ip_i4_p) :: status
1053 logical :: whead,wdata
1054 real(ip_double_p),
allocatable :: array3(:,:,:)
1055 integer(ip_i4_p) ,
allocatable :: time(:)
1056 real(ip_double_p) :: tbnds(2)
1058 character(len=*),
parameter :: subname =
'(oasis_io_read_avfbf)'
1064 IF (mpicom == mpi_comm_null)
return
1069 if (
present(msec))
then
1074 if (
present(f_string))
then
1075 lstring = trim(f_string)
1079 call mpi_comm_rank(mpicom,iam,ierr)
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)
1090 lfn = trim(itemc)//trim(lstring)//
'.nc'
1093 inquire(file=trim(lfn),exist=exists)
1094 if (.not.exists)
then
1095 write(nulprt,*) subname,estr,
'file not found ',trim(lfn)
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))
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))
1118 if (time(j) == lmsec(1)) n1 = j
1122 write(nulprt,*) subname,estr,
'time not found on file ',trim(lfn),lmsec
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))
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
1149 allocate(array3(nx,ny,1))
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))
1162 av_g%rAttr(n,n1) = array3(i,j,1)
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)
1192 character(len=*) ,
intent(in) :: filename
1193 character(len=*) ,
intent(in) :: fldname
1194 integer(ip_i4_p) ,
intent(inout),
optional :: ifld2(:,:)
1195 real(ip_realwp_p),
intent(inout),
optional :: fld2(:,:)
1196 real(ip_realwp_p),
intent(inout),
optional :: fld3(:,:,:)
1197 integer(ip_i4_p) ,
intent(inout),
optional :: nx
1198 integer(ip_i4_p) ,
intent(inout),
optional :: ny
1199 integer(ip_i4_p) ,
intent(inout),
optional :: nz
1202 integer(ip_i4_p) :: ncid,varid
1203 integer(ip_i4_p) :: n,ndims,xtype
1204 integer(ip_i4_p),
allocatable :: dimid(:),nd(:)
1205 integer(ip_i4_p) :: status
1206 integer(ip_i4_p) :: ind
1208 character(len=ic_med) :: gridname
1210 character(len=*),
parameter :: subname =
'(oasis_io_read_field_fromroot)'
1221 if (oasis_debug >= 5)
then
1222 write(nulprt,*) subname,
' read ',trim(filename),
' ',trim(fldname)
1225 inquire(file=trim(filename),exist=exists)
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))
1231 write(nulprt,*) subname,estr,
'in filename ',trim(filename)
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)
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))
1245 allocate(dimid(ndims),nd(ndims))
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))
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))
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))
1270 write(nulprt,*) subname,estr,
'mismatch in field and data'
1275 status = nf90_close(ncid)
1276 IF (status /= nf90_noerr)
WRITE(nulprt,*) subname,
' model :',compid,
' proc :',&
1277 mpi_rank_local,
':',trim(nf90_strerror(status))
1279 if (
present(nx))
then
1283 if (
present(ny))
then
1287 if (
present(nz))
then
1291 deallocate(dimid,nd)
1311 character(len=*),
intent(in) :: filename
1312 character(len=*),
intent(in) :: fldname
1313 real(ip_realwp_p),
intent(in) :: fld(:,:)
1314 integer(ip_i4_p),
intent(in) :: nx
1315 integer(ip_i4_p),
intent(in) :: ny
1318 integer(ip_i4_p) :: ncid,dimid,dimid2(2),varid
1319 integer(ip_i4_p) :: status
1320 integer(ip_i4_p) :: ind
1322 character(len=ic_med) :: gridname
1324 character(len=*),
parameter :: subname =
'(oasis_io_write_2dgridfld_fromroot)'
1335 if (oasis_debug >= 5)
then
1336 write(nulprt,*) subname,
' write ',trim(filename),
' ',trim(fldname),nx,ny
1339 ind = index(trim(fldname),
'.')
1341 write(nulprt,*) subname,estr,
'in fldname ',trim(fldname)
1344 gridname = fldname(1:ind-1)
1346 inquire(file=trim(filename),exist=exists)
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)
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))
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))
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))
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))
1387 status = nf90_enddef(ncid)
1388 IF (status /= nf90_noerr)
WRITE(nulprt,*) subname,
' model :',compid,
' proc :',&
1389 mpi_rank_local,
':',trim(nf90_strerror(status))
1392 status = nf90_close(ncid)
1393 IF (status /= nf90_noerr)
WRITE(nulprt,*) subname,
' model :',compid,
' proc :',&
1394 mpi_rank_local,
':',trim(nf90_strerror(status))
1414 character(len=*),
intent(in) :: filename
1415 character(len=*),
intent(in) :: fldname
1416 integer(ip_i4_p),
intent(in) :: fld(:,:)
1417 integer(ip_i4_p),
intent(in) :: nx
1418 integer(ip_i4_p),
intent(in) :: ny
1421 integer(ip_i4_p) :: ncid,dimid,dimid2(2),varid
1422 integer(ip_i4_p) :: status
1423 integer(ip_i4_p) :: ind
1425 character(len=ic_med) :: gridname
1427 character(len=*),
parameter :: subname =
'(oasis_io_write_2dgridint_fromroot)'
1438 if (oasis_debug >= 5)
then
1439 write(nulprt,*) subname,
' write ',trim(filename),
' ',trim(fldname),nx,ny
1442 ind = index(trim(fldname),
'.')
1444 write(nulprt,*) subname,estr,
'in fldname ',trim(fldname)
1447 gridname = fldname(1:ind-1)
1449 inquire(file=trim(filename),exist=exists)
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)
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))
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))
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))
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))
1490 status = nf90_enddef(ncid)
1491 IF (status /= nf90_noerr)
WRITE(nulprt,*) subname,
' model :',compid,
' proc :',&
1492 mpi_rank_local,
':',trim(nf90_strerror(status))
1495 status = nf90_close(ncid)
1496 IF (status /= nf90_noerr)
WRITE(nulprt,*) subname,
' model :',compid,
' proc :',&
1497 mpi_rank_local,
':',trim(nf90_strerror(status))
1517 character(len=*),
intent(in) :: filename
1518 character(len=*),
intent(in) :: fldname
1519 real(ip_realwp_p),
intent(in) :: fld(:,:,:)
1520 integer(ip_i4_p),
intent(in) :: nx
1521 integer(ip_i4_p),
intent(in) :: ny
1522 integer(ip_i4_p),
intent(in) :: nc
1525 integer(ip_i4_p) :: ncid,dimid,dimid3(3),varid
1526 integer(ip_i4_p) :: status
1527 integer(ip_i4_p) :: ind
1529 character(len=ic_med) :: gridname
1531 character(len=*),
parameter :: subname =
'(oasis_io_write_3dgridfld_fromroot)'
1542 if (oasis_debug >= 5)
then
1543 write(nulprt,*) subname,
' write ',trim(filename),
' ',trim(fldname),nx,ny,nc
1546 ind = index(trim(fldname),
'.')
1548 write(nulprt,*) subname,estr,
'in fldname ',trim(fldname)
1551 gridname = fldname(1:ind-1)
1553 inquire(file=trim(filename),exist=exists)
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)
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))
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))
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))
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))
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))
1602 status = nf90_enddef(ncid)
1603 IF (status /= nf90_noerr)
WRITE(nulprt,*) subname,
' model :',compid,
' proc :',&
1604 mpi_rank_local,
':',trim(nf90_strerror(status))
1607 status = nf90_close(ncid)
1608 IF (status /= nf90_noerr)
WRITE(nulprt,*) subname,
' model :',compid,
' proc :',&
1609 mpi_rank_local,
':',trim(nf90_strerror(status))
subroutine, public oasis_debug_exit(string)
Used when a subroutine is exited, write info to log file at some debug level.
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.
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.
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.