35 module time_interp_external2_mod
44 use platform_mod,
only : r8_kind, r4_kind
54 use time_interp_mod,
only :
time_interp, time_interp_init
57 use platform_mod,
only: r8_kind, fms_path_len, fms_file_len
58 use horiz_interp_mod,
only :
horiz_interp, horiz_interp_type
59 use horiz_interp_type_mod,
only : conserve
60 use fms2_io_mod,
only : valid_t, fmsnetcdfdomainfile_t,
open_file, get_unlimited_dimension_name, &
61 variable_att_exists, fmsnetcdffile_t, &
62 variable_exists, get_valid, get_variable_num_dimensions,
read_data, &
63 is_valid,
close_file, get_dimension_size, get_variable_dimension_names, &
64 get_variable_size, get_time_calendar, get_variable_missing, get_variable_units
70 #include<file_version.h>
72 integer,
parameter,
public :: NO_REGION=0, inside_region=1, outside_region=2
73 integer,
parameter,
private :: modulo_year= 0001
74 integer,
parameter,
private :: LINEAR_TIME_INTERP = 1
75 integer,
parameter,
public :: SUCCESS = 0, err_field_not_found = 1
76 real(r8_kind),
parameter,
private :: DEFAULT_MISSING_VALUE = -1e20_r8_kind
77 integer,
private :: max_fields = 100, max_files= 40
78 integer,
private :: num_fields = 0, num_files=0
80 integer,
private :: num_io_buffers = 2
81 logical,
private :: module_initialized = .false.
82 logical,
private :: debug_this_module = .false.
87 public get_external_fileobj
97 type(fmsnetcdffile_t),
pointer :: fileobj=>null()
98 character(len=128) :: name, units
99 integer :: siz(4), ndim
100 character(len=32) :: axisname(4)
103 type(
time_type),
dimension(:),
pointer :: start_time =>null(), end_time =>null()
104 type(
time_type),
dimension(:),
pointer :: period =>null()
105 logical :: modulo_time
106 real(r8_kind),
dimension(:,:,:,:),
pointer :: domain_data =>null()
107 logical,
dimension(:,:,:,:),
pointer :: mask =>null()
108 integer,
dimension(:),
pointer :: ibuf =>null()
109 real(r8_kind),
dimension(:,:,:,:),
pointer :: src_data =>null()
110 type(valid_t) :: valid
112 logical :: domain_present
113 real(r8_kind) :: slope, intercept
114 integer :: isc,iec,jsc,jec
115 type(
time_type) :: modulo_time_beg, modulo_time_end
116 logical :: have_modulo_times, correct_leap_year_inconsistency
117 integer :: region_type
118 integer :: is_region, ie_region, js_region, je_region
119 integer :: is_src, ie_src, js_src, je_src
121 integer :: numwindows
122 logical,
dimension(:,:),
pointer :: need_compute=>null()
123 real(r8_kind) :: missing
129 character(len=FMS_FILE_LEN) :: filename =
''
130 type(FmsNetcdfFile_t),
pointer :: fileobj => null()
146 module procedure time_interp_external_0d_r4
147 module procedure time_interp_external_1d_r4
148 module procedure time_interp_external_2d_r4
149 module procedure time_interp_external_3d_r4
150 module procedure time_interp_external_0d_r8
151 module procedure time_interp_external_1d_r8
152 module procedure time_interp_external_2d_r8
153 module procedure time_interp_external_3d_r8
157 module procedure time_interp_external_bridge_0d_r4
158 module procedure time_interp_external_bridge_1d_r4
159 module procedure time_interp_external_bridge_2d_r4
160 module procedure time_interp_external_bridge_3d_r4
161 module procedure time_interp_external_bridge_0d_r8
162 module procedure time_interp_external_bridge_1d_r8
163 module procedure time_interp_external_bridge_2d_r8
164 module procedure time_interp_external_bridge_3d_r8
172 type(ext_fieldtype),
save,
private,
pointer :: loaded_fields(:) => null()
173 type(filetype),
save,
private,
pointer :: opened_files(:) => null()
175 real(r8_kind),
private,
parameter :: time_interp_missing=-1e99_r8_kind
187 integer :: io_status, logunit, ierr
189 namelist /time_interp_external_nml/ num_io_buffers, debug_this_module, &
190 max_fields, max_files
194 if(module_initialized)
return
198 call write_version_number(
"TIME_INTERP_EXTERNAL_MOD", version)
200 read (input_nml_file, time_interp_external_nml, iostat=io_status)
203 write(logunit,time_interp_external_nml)
207 module_initialized = .true.
209 call time_interp_init()
275 verbose,axis_names, axis_sizes,override,correct_leap_year_inconsistency,&
276 permit_calendar_conversion,use_comp_domain,ierr, nwindows, ignore_axis_atts, ongrid )
278 character(len=*),
intent(in) :: file,fieldname
279 logical,
intent(in),
optional :: verbose
280 character(len=*),
intent(in),
optional :: desired_units
281 type(
domain2d),
intent(in),
optional :: domain
282 integer,
intent(inout),
optional :: axis_sizes(4)
283 character(len=*),
intent(out),
optional :: axis_names(4)
284 logical,
intent(in),
optional :: override, correct_leap_year_inconsistency,&
285 permit_calendar_conversion,use_comp_domain
286 integer,
intent(out),
optional :: ierr
287 integer,
intent(in),
optional :: nwindows
288 logical,
optional :: ignore_axis_atts
289 logical,
optional :: ongrid
291 logical :: ongrid_local
295 real(r8_kind) :: slope, intercept
296 integer :: ndim,ntime,i,j
297 integer :: iscomp,iecomp,jscomp,jecomp,isglobal,ieglobal,jsglobal,jeglobal
298 integer :: isdata,iedata,jsdata,jedata, dxsize, dysize,dxsize_max,dysize_max
299 logical :: verb, transpose_xy,use_comp_domain1
300 real(r8_kind),
dimension(:),
allocatable :: tstamp, tstart, tend, tavg
301 character(len=1) :: cart
302 character(len=1),
dimension(4) :: cart_dir
303 character(len=128) :: units, fld_units
304 character(len=128) :: msg, calendar_type, timebeg, timeend
305 character(len=128) :: timename, timeunits
306 character(len=128),
allocatable :: axisname(:)
307 integer,
allocatable :: axislen(:)
308 integer :: siz(4), siz_in(4), gxsize, gysize,gxsize_max, gysize_max
310 integer :: yr, mon, day, hr, minu, sec
311 integer :: len, nfile, nfields_orig, nbuf, nx,ny
312 integer :: numwindows
313 logical :: ignore_axatts
314 logical :: have_modulo_time
315 type(fmsnetcdffile_t),
pointer :: fileobj=>null()
316 integer,
dimension(:),
allocatable :: pes
318 if (.not. module_initialized)
call mpp_error(fatal,
'Must call time_interp_external_init first')
319 if(
present(ierr)) ierr = success
320 ignore_axatts=.false.
321 cart_dir(1)=
'X';cart_dir(2)=
'Y';cart_dir(3)=
'Z';cart_dir(4)=
'T'
322 if(
present(ignore_axis_atts)) ignore_axatts = ignore_axis_atts
323 use_comp_domain1 = .false.
324 if(
PRESENT(use_comp_domain)) use_comp_domain1 = use_comp_domain
326 if (
PRESENT(verbose)) verb=verbose
327 if (debug_this_module) verb = .true.
329 if(
present(nwindows)) numwindows = nwindows
332 if (
PRESENT(desired_units))
then
333 units = desired_units
334 call mpp_error(fatal,
'==> Unit conversion via time_interp_external &
335 &has been temporarily deprecated. Previous versions of&
336 &this module used udunits_mod to perform unit conversion.&
337 & Udunits_mod is in the process of being replaced since &
338 &there were portability issues associated with this code.&
339 & Please remove the desired_units argument from calls to &
344 if(trim(opened_files(i)%filename) == trim(file))
then
350 num_files = num_files + 1
351 if(num_files > max_files)
then
357 call mpp_error(fatal,
"time_interp_external: num_files is greater than max_files, "// &
358 "increase time_interp_external_nml max_files")
360 opened_files(num_files)%filename = trim(file)
361 allocate(opened_files(num_files)%fileobj)
362 fileobj => opened_files(num_files)%fileobj
364 if(.not.
open_file(opened_files(num_files)%fileobj, trim(file),
'read')) &
365 call mpp_error(fatal,
'time_interp_external_mod: Error in opening file '//trim(file))
367 fileobj => opened_files(nfile)%fileobj
370 call get_unlimited_dimension_name(fileobj, timename)
371 call get_dimension_size(fileobj, timename, ntime)
374 write(msg,
'(a15,a,a58)')
'external field ',trim(fieldname),&
375 ' does not have an associated record dimension (REQUIRED) '
380 call get_time_calendar(fileobj, timename, calendar_type)
383 call get_variable_units(fileobj, timename, timeunits)
389 call mpp_get_current_pelist(pes)
390 allocate(tstamp(ntime),tstart(ntime),tend(ntime),tavg(ntime))
393 if (mpp_root_pe() .eq.
mpp_pe())
call read_data(fileobj, timename, tstamp)
394 call mpp_broadcast(tstamp,
size(tstamp), mpp_root_pe(), pelist=pes)
397 transpose_xy = .false.
398 isdata=1; iedata=1; jsdata=1; jedata=1
402 if (
PRESENT(domain))
then
404 nx = iecomp-iscomp+1; ny = jecomp-jscomp+1
405 call mpp_get_data_domain(domain,isdata,iedata,jsdata,jedata,dxsize,dxsize_max,dysize,dysize_max)
406 call mpp_get_global_domain(domain,isglobal,ieglobal,jsglobal,jeglobal,gxsize,gxsize_max,gysize,gysize_max)
407 ongrid_local = .false.
408 if (
present(ongrid)) ongrid_local = ongrid
411 if (ongrid_local)
then
417 elseif(use_comp_domain1)
then
418 call mpp_error(fatal,
"init_external_field:"//&
419 " use_comp_domain=true but domain is not present")
423 nfields_orig = num_fields
425 if (.not. variable_exists(fileobj, fieldname) )
then
426 if (
present(ierr))
then
427 ierr = err_field_not_found
430 call mpp_error(fatal,
'external field "'//trim(fieldname)//
'" not found in file "'//trim(file)//
'"')
437 if(variable_att_exists(fileobj, fieldname,
'time_avg_info'))
then
438 if(variable_exists(fileobj,
'average_T1'))
call read_data(fileobj,
'average_T1', tstart)
439 if(variable_exists(fileobj,
'average_T2'))
call read_data(fileobj,
'average_T2', tend)
440 if(variable_exists(fileobj,
'average_DT'))
call read_data(fileobj,
'average_DT', tavg)
443 num_fields = num_fields + 1
444 if(num_fields > max_fields)
then
450 call mpp_error(fatal,
"time_interp_external: num_fields is greater than max_fields, "// &
451 "increase time_interp_external_nml max_fields")
455 call get_variable_units(fileobj, fieldname, fld_units)
458 loaded_fields(num_fields)%fileobj => fileobj
459 loaded_fields(num_fields)%name = trim(fieldname)
460 loaded_fields(num_fields)%units = trim(fld_units)
461 loaded_fields(num_fields)%isc = 1
462 loaded_fields(num_fields)%iec = 1
463 loaded_fields(num_fields)%jsc = 1
464 loaded_fields(num_fields)%jec = 1
465 loaded_fields(num_fields)%region_type = no_region
466 loaded_fields(num_fields)%is_region = 0
467 loaded_fields(num_fields)%ie_region = -1
468 loaded_fields(num_fields)%js_region = 0
469 loaded_fields(num_fields)%je_region = -1
470 if (
PRESENT(domain))
then
471 loaded_fields(num_fields)%domain_present = .true.
472 loaded_fields(num_fields)%domain = domain
473 loaded_fields(num_fields)%isc=iscomp;loaded_fields(num_fields)%iec = iecomp
474 loaded_fields(num_fields)%jsc=jscomp;loaded_fields(num_fields)%jec = jecomp
476 loaded_fields(num_fields)%domain_present = .false.
479 loaded_fields(num_fields)%valid = get_valid(fileobj, fieldname)
480 ndim = get_variable_num_dimensions(fileobj, fieldname)
482 'invalid array rank <=4d fields supported')
484 loaded_fields(num_fields)%siz = 1
485 loaded_fields(num_fields)%ndim = ndim
486 loaded_fields(num_fields)%tdim = 4
488 loaded_fields(num_fields)%missing = get_variable_missing(fileobj, fieldname)
490 allocate(axisname(ndim), axislen(ndim))
492 call get_variable_dimension_names(fileobj, fieldname, axisname)
493 call get_variable_size(fileobj, fieldname, axislen)
494 do j=1,loaded_fields(num_fields)%ndim
497 if (cart ==
'N' .and. .not. ignore_axatts)
then
498 write(msg,
'(a,"/",a)') trim(file),trim(fieldname)
499 call mpp_error(fatal,
'file/field '//trim(msg)// &
500 ' couldnt recognize axis atts in time_interp_external')
501 else if (cart ==
'N' .and. ignore_axatts)
then
506 if (j.eq.2) transpose_xy = .true.
507 if (.not.
PRESENT(domain) .and. .not.
PRESENT(override))
then
512 loaded_fields(num_fields)%isc=iscomp;loaded_fields(num_fields)%iec=iecomp
513 elseif (
PRESENT(override))
then
515 if (
PRESENT(axis_sizes)) axis_sizes(1) = len
517 loaded_fields(num_fields)%axisname(1) = axisname(j)
518 if(use_comp_domain1)
then
519 loaded_fields(num_fields)%siz(1) = nx
521 loaded_fields(num_fields)%siz(1) = dxsize
523 if (len /= gxsize)
then
524 write(msg,
'(a,"/",a)') trim(file),trim(fieldname)
525 call mpp_error(fatal,
'time_interp_ext, file/field '//trim(msg)//
' x dim doesnt match model')
528 loaded_fields(num_fields)%axisname(2) = axisname(j)
529 if (.not.
PRESENT(domain) .and. .not.
PRESENT(override))
then
534 loaded_fields(num_fields)%jsc=jscomp;loaded_fields(num_fields)%jec=jecomp
535 elseif (
PRESENT(override))
then
537 if (
PRESENT(axis_sizes)) axis_sizes(2) = len
539 if(use_comp_domain1)
then
540 loaded_fields(num_fields)%siz(2) = ny
542 loaded_fields(num_fields)%siz(2) = dysize
544 if (len /= gysize)
then
545 write(msg,
'(a,"/",a)') trim(file),trim(fieldname)
546 call mpp_error(fatal,
'time_interp_ext, file/field '//trim(msg)//
' y dim doesnt match model')
549 loaded_fields(num_fields)%axisname(3) = axisname(j)
550 loaded_fields(num_fields)%siz(3) = len
552 loaded_fields(num_fields)%axisname(4) = axisname(j)
553 loaded_fields(num_fields)%siz(4) = ntime
554 loaded_fields(num_fields)%tdim = j
557 siz = loaded_fields(num_fields)%siz
558 if(
PRESENT(axis_names)) axis_names = loaded_fields(num_fields)%axisname
559 if (
PRESENT(axis_sizes) .and. .not.
PRESENT(override))
then
560 axis_sizes = loaded_fields(num_fields)%siz
563 if (verb)
write(outunit,
'(a,4i6)')
'field x,y,z,t local size= ',siz
564 if (verb)
write(outunit,*)
'field contains data in units = ',trim(loaded_fields(num_fields)%units)
565 if (transpose_xy)
call mpp_error(fatal,
'axis ordering not supported')
566 if (num_io_buffers .le. 1)
call mpp_error(fatal,
'time_interp_ext:num_io_buffers should be at least 2')
567 nbuf = min(num_io_buffers,siz(4))
569 loaded_fields(num_fields)%numwindows = numwindows
570 allocate(loaded_fields(num_fields)%need_compute(nbuf, numwindows))
571 loaded_fields(num_fields)%need_compute = .true.
573 allocate(loaded_fields(num_fields)%domain_data(isdata:iedata,jsdata:jedata,siz(3),nbuf),&
574 loaded_fields(num_fields)%mask(isdata:iedata,jsdata:jedata,siz(3),nbuf) )
575 loaded_fields(num_fields)%mask = .false.
576 loaded_fields(num_fields)%domain_data = 0.0_r8_kind
577 slope=1.0_r8_kind;intercept=0.0_r8_kind
583 loaded_fields(num_fields)%slope = slope
584 loaded_fields(num_fields)%intercept = intercept
585 allocate(loaded_fields(num_fields)%ibuf(nbuf))
586 loaded_fields(num_fields)%ibuf = -1
587 loaded_fields(num_fields)%nbuf = 0
588 if(
PRESENT(override))
then
589 loaded_fields(num_fields)%is_src = 1
590 loaded_fields(num_fields)%ie_src = gxsize
591 loaded_fields(num_fields)%js_src = 1
592 loaded_fields(num_fields)%je_src = gysize
593 allocate(loaded_fields(num_fields)%src_data(gxsize,gysize,siz(3),nbuf))
595 loaded_fields(num_fields)%is_src = isdata
596 loaded_fields(num_fields)%ie_src = iedata
597 loaded_fields(num_fields)%js_src = jsdata
598 loaded_fields(num_fields)%je_src = jedata
599 allocate(loaded_fields(num_fields)%src_data(isdata:iedata,jsdata:jedata,siz(3),nbuf))
602 allocate(loaded_fields(num_fields)%time(ntime))
603 allocate(loaded_fields(num_fields)%period(ntime))
604 allocate(loaded_fields(num_fields)%start_time(ntime))
605 allocate(loaded_fields(num_fields)%end_time(ntime))
608 loaded_fields(num_fields)%time(j) =
get_cal_time(tstamp(j),trim(timeunits),trim(calendar_type), &
609 & permit_calendar_conversion)
610 loaded_fields(num_fields)%start_time(j) =
get_cal_time(tstart(j),trim(timeunits),trim(calendar_type), &
611 & permit_calendar_conversion)
612 loaded_fields(num_fields)%end_time(j) =
get_cal_time( tend(j),trim(timeunits),trim(calendar_type), &
613 & permit_calendar_conversion)
616 if (loaded_fields(num_fields)%modulo_time)
then
617 call set_time_modulo(loaded_fields(num_fields)%Time)
618 call set_time_modulo(loaded_fields(num_fields)%start_time)
619 call set_time_modulo(loaded_fields(num_fields)%end_time)
622 if(
present(correct_leap_year_inconsistency))
then
623 loaded_fields(num_fields)%correct_leap_year_inconsistency = correct_leap_year_inconsistency
625 loaded_fields(num_fields)%correct_leap_year_inconsistency = .false.
628 if(have_modulo_time)
then
630 loaded_fields(num_fields)%modulo_time_beg =
set_time(timebeg)
631 loaded_fields(num_fields)%modulo_time_end =
set_time(timeend)
633 loaded_fields(num_fields)%modulo_time_beg =
set_date(timebeg)
634 loaded_fields(num_fields)%modulo_time_end =
set_date(timeend)
636 loaded_fields(num_fields)%have_modulo_times = .true.
638 loaded_fields(num_fields)%have_modulo_times = .false.
641 call mpp_error(note,
'time_interp_external_mod: file '//trim(file)//
' has only one time level')
644 loaded_fields(num_fields)%period(j) = loaded_fields(num_fields)%end_time(j) &
645 - loaded_fields(num_fields)%start_time(j)
646 if (loaded_fields(num_fields)%period(j) >
set_time(0,0))
then
647 call get_time(loaded_fields(num_fields)%period(j), sec, day)
648 sec = sec/2+mod(day,2)*43200
650 loaded_fields(num_fields)%time(j) = loaded_fields(num_fields)%start_time(j)+&
653 if (j > 1 .and. j < ntime)
then
654 tdiff = loaded_fields(num_fields)%time(j+1) - loaded_fields(num_fields)%time(j-1)
656 sec = sec/2+mod(day,2)*43200
658 loaded_fields(num_fields)%period(j) =
set_time(sec,day)
659 sec = sec/2+mod(day,2)*43200
661 loaded_fields(num_fields)%start_time(j) = loaded_fields(num_fields)%time(j) -
set_time(sec,day)
662 loaded_fields(num_fields)%end_time(j) = loaded_fields(num_fields)%time(j) +
set_time(sec,day)
663 elseif ( j == 1)
then
664 tdiff = loaded_fields(num_fields)%time(2) - loaded_fields(num_fields)%time(1)
666 loaded_fields(num_fields)%period(j) =
set_time(sec,day)
667 sec = sec/2+mod(day,2)*43200
669 loaded_fields(num_fields)%start_time(j) = loaded_fields(num_fields)%time(j) -
set_time(sec,day)
670 loaded_fields(num_fields)%end_time(j) = loaded_fields(num_fields)%time(j) +
set_time(sec,day)
672 tdiff = loaded_fields(num_fields)%time(ntime) - loaded_fields(num_fields)%time(ntime-1)
674 loaded_fields(num_fields)%period(j) =
set_time(sec,day)
675 sec = sec/2+mod(day,2)*43200
677 loaded_fields(num_fields)%start_time(j) = loaded_fields(num_fields)%time(j) -
set_time(sec,day)
678 loaded_fields(num_fields)%end_time(j) = loaded_fields(num_fields)%time(j) +
set_time(sec,day)
685 if (loaded_fields(num_fields)%time(j) >= loaded_fields(num_fields)%time(j+1))
then
686 write(msg,
'(A,i20)')
"times not monotonically increasing. Filename: " &
687 //trim(file)//
" field: "//trim(fieldname)//
" timeslice: ", j
692 loaded_fields(num_fields)%modulo_time =
get_axis_modulo(fileobj, timename)
695 if (loaded_fields(num_fields)%modulo_time)
write(outunit,*)
'data are being treated as modulo in time'
697 write(outunit,*)
'time index, ', j
698 call get_date(loaded_fields(num_fields)%start_time(j),yr,mon,day,hr,minu,sec)
699 write(outunit,
'(a,i4,a,i2,a,i2,1x,i2,a,i2,a,i2)') &
700 'start time: yyyy/mm/dd hh:mm:ss= ',yr,
'/',mon,
'/',day,hr,
':',minu,
':',sec
701 call get_date(loaded_fields(num_fields)%time(j),yr,mon,day,hr,minu,sec)
702 write(outunit,
'(a,i4,a,i2,a,i2,1x,i2,a,i2,a,i2)') &
703 'mid time: yyyy/mm/dd hh:mm:ss= ',yr,
'/',mon,
'/',day,hr,
':',minu,
':',sec
704 call get_date(loaded_fields(num_fields)%end_time(j),yr,mon,day,hr,minu,sec)
705 write(outunit,
'(a,i4,a,i2,a,i2,1x,i2,a,i2,a,i2)') &
706 'end time: yyyy/mm/dd hh:mm:ss= ',yr,
'/',mon,
'/',day,hr,
':',minu,
':',sec
710 deallocate(axisname, axislen)
711 deallocate(tstamp, tstart, tend, tavg)
718 function get_external_fileobj(filename, fileobj)
719 character(len=*),
intent(in) :: filename
720 type(fmsnetcdffile_t),
intent(out) :: fileobj
721 logical :: get_external_fileobj
724 get_external_fileobj = .false.
726 if(trim(opened_files(i)%filename) == trim(filename))
then
727 fileobj = opened_files(i)%fileobj
728 get_external_fileobj = .true.
733 end function get_external_fileobj
736 subroutine set_time_modulo(Time)
738 type(
time_type),
intent(inout),
dimension(:) :: time
741 integer :: yr, mon, dy, hr, minu, sec
743 ntime =
size(time(:))
746 call get_date(time(n), yr, mon, dy, hr, minu, sec)
748 time(n) =
set_date(yr, mon, dy, hr, minu, sec)
752 end subroutine set_time_modulo
757 subroutine load_record(field, rec, interp, is_in, ie_in, js_in, je_in, window_id_in)
759 integer ,
intent(in) :: rec
761 integer,
intent(in),
optional :: is_in, ie_in, js_in, je_in
762 integer,
intent(in),
optional :: window_id_in
766 integer :: isw,iew,jsw,jew
767 integer :: is_region, ie_region, js_region, je_region, i, j
768 integer :: start(4), nread(4)
769 logical :: need_compute
771 real(r8_kind) :: mask_in(size(field%src_data,1),size(field%src_data,2),size(field%src_data,3))
772 real(r8_kind),
allocatable :: mask_out(:,:,:)
773 real(r4_kind),
allocatable :: hi_tmp_data(:,:,:,:)
774 real(r4_kind),
allocatable :: hi_tmp_msk_out(:,:,:)
777 if(
PRESENT(window_id_in) ) window_id = window_id_in
778 need_compute = .true.
785 need_compute = .false.
788 field%nbuf = field%nbuf + 1
789 if(field%nbuf >
size(field%domain_data,4).or.field%nbuf <= 0) field%nbuf = 1
792 field%need_compute(ib,:) = .true.
794 if (debug_this_module)
write(outunit,*)
'reading record without domain for field ',trim(field%name)
796 start(1) = field%is_src; nread(1) = field%ie_src - field%is_src + 1
797 start(2) = field%js_src; nread(2) = field%je_src - field%js_src + 1
798 start(3) = 1; nread(3) =
size(field%src_data,3)
799 start(field%tdim) = rec; nread(field%tdim) = 1
800 call read_data(field%fileobj,field%name,field%src_data(:,:,:,ib:ib),corner=start,edge_lengths=nread)
803 isw=field%isc;iew=field%iec
804 jsw=field%jsc;jew=field%jec
806 if( field%numwindows > 1)
then
807 if( .NOT.
PRESENT(is_in) .OR. .NOT.
PRESENT(ie_in) .OR. .NOT.
PRESENT(js_in) .OR. .NOT.
PRESENT(je_in) )
then
809 &
'time_interp_external(load_record): is_in, ie_in, js_in, je_in must be present when numwindows>1')
811 isw = isw + is_in - 1
812 iew = isw + ie_in - is_in
813 jsw = jsw + js_in - 1
814 jew = jsw + je_in - js_in
819 need_compute = field%need_compute(ib, window_id)
820 if(need_compute)
then
821 if(
PRESENT(interp))
then
822 is_region = field%is_region; ie_region = field%ie_region
823 js_region = field%js_region; je_region = field%je_region
824 mask_in = 0.0_r8_kind
825 where (is_valid(field%src_data(:,:,:,ib), field%valid)) mask_in = 1.0_r8_kind
826 if ( field%region_type .NE. no_region )
then
827 if( any(mask_in == 0.0_r8_kind) )
then
828 call mpp_error(fatal,
"time_interp_external: mask_in should be all 1 when region_type is not NO_REGION")
830 if(interp%interp_method == conserve)
then
832 "conservative interpolation when region_type is not NO_REGION. Make sure interp matches the region")
834 if( field%region_type == outside_region)
then
835 do j = js_region, je_region
836 do i = is_region, ie_region
837 mask_in(i,j,:) = 0.0_r8_kind
841 do j = 1,
size(mask_in,2)
842 do i = 1,
size(mask_in,1)
843 if( j<js_region .OR. j>je_region .OR. i<is_region .OR. i>ie_region ) mask_in(i,j,:) = 0.0_r8_kind
852 if (interp%horizInterpReals4_type%is_allocated)
then
854 allocate(hi_tmp_msk_out(isw:iew,jsw:jew,
SIZE(field%src_data,3)))
855 allocate(hi_tmp_data(lbound(field%domain_data,1):ubound(field%domain_data,1), &
856 lbound(field%domain_data,2):ubound(field%domain_data,2), &
857 lbound(field%domain_data,3):ubound(field%domain_data,3), &
858 lbound(field%domain_data,4):ubound(field%domain_data,4)))
860 hi_tmp_data = real(field%domain_data, r4_kind)
862 if(interp%interp_method == conserve)
then
863 call horiz_interp(interp, real(field%src_data(:,:,:,ib),r4_kind), hi_tmp_data(isw:iew,jsw:jew,:,ib))
864 hi_tmp_msk_out = 1.0_r4_kind
866 call horiz_interp(interp, real(field%src_data(:,:,:,ib),r4_kind), hi_tmp_data(isw:iew,jsw:jew,:,ib), &
867 mask_in=real(mask_in,r4_kind), mask_out=hi_tmp_msk_out)
870 field%domain_data = real(hi_tmp_data, r8_kind)
871 field%mask(isw:iew,jsw:jew,:,ib) = hi_tmp_msk_out(isw:iew,jsw:jew,:) > 0.0_r4_kind
873 if(
allocated(hi_tmp_data))
deallocate(hi_tmp_data)
874 if(
allocated(hi_tmp_msk_out))
deallocate(hi_tmp_msk_out)
876 allocate(mask_out(isw:iew,jsw:jew,
size(field%src_data,3)))
877 if(interp%interp_method == conserve)
then
878 call horiz_interp(interp, field%src_data(:,:,:,ib), field%domain_data(isw:iew,jsw:jew,:,ib))
879 mask_out = 1.0_r8_kind
881 call horiz_interp(interp, field%src_data(:,:,:,ib),field%domain_data(isw:iew,jsw:jew,:,ib), &
882 mask_in=mask_in, mask_out=mask_out)
884 field%mask(isw:iew,jsw:jew,:,ib) = mask_out(isw:iew,jsw:jew,:) > 0.0_r8_kind
889 if ( field%region_type .NE. no_region )
then
890 call mpp_error(fatal,
"time_interp_external: region_type should be NO_REGION when interp is not present")
892 field%domain_data(isw:iew,jsw:jew,:,ib) = field%src_data(isw:iew,jsw:jew,:,ib)
893 field%mask(isw:iew,jsw:jew,:,ib) = is_valid(field%domain_data(isw:iew,jsw:jew,:,ib),field%valid)
896 where(field%mask(isw:iew,jsw:jew,:,ib)) field%domain_data(isw:iew,jsw:jew,:,ib) = &
897 field%domain_data(isw:iew,jsw:jew,:,ib)*field%slope + field%intercept
898 field%need_compute(ib, window_id) = .false.
907 integer ,
intent(in) :: rec
910 integer :: start(4), nread(4)
918 field%nbuf = field%nbuf + 1
919 if(field%nbuf >
size(field%domain_data,4).or.field%nbuf <= 0) field%nbuf = 1
923 if (debug_this_module)
write(outunit,*)
'reading record without domain for field ',trim(field%name)
925 start(3) = 1; nread(3) =
size(field%src_data,3)
926 start(field%tdim) = rec; nread(field%tdim) = 1
927 call read_data(field%fileobj,field%name,field%src_data(:,:,:,ib),corner=start,edge_lengths=nread)
928 if ( field%region_type .NE. no_region )
then
929 call mpp_error(fatal,
"time_interp_external: region_type should be NO_REGION when field is scalar")
931 field%domain_data(1,1,:,ib) = field%src_data(1,1,:,ib)
932 field%mask(1,1,:,ib) = is_valid(field%domain_data(1,1,:,ib),field%valid)
934 where(field%mask(1,1,:,ib)) field%domain_data(1,1,:,ib) = &
935 field%domain_data(1,1,:,ib)*field%slope + field%intercept
942 integer,
intent(in) :: index
943 integer,
intent(in) :: is, ie, js, je
946 if( is == loaded_fields(index)%is_src .AND. ie == loaded_fields(index)%ie_src .AND. &
947 js == loaded_fields(index)%js_src .AND. ie == loaded_fields(index)%je_src )
return
949 if( .NOT.
ASSOCIATED(loaded_fields(index)%src_data) )
call mpp_error(fatal, &
950 "time_interp_external: field(index)%src_data is not associated")
951 nk =
size(loaded_fields(index)%src_data,3)
952 nbuf =
size(loaded_fields(index)%src_data,4)
953 deallocate(loaded_fields(index)%src_data)
954 allocate(loaded_fields(index)%src_data(is:ie,js:je,nk,nbuf))
955 loaded_fields(index)%is_src = is
956 loaded_fields(index)%ie_src = ie
957 loaded_fields(index)%js_src = js
958 loaded_fields(index)%je_src = je
963 subroutine set_override_region(index, region_type, is_region, ie_region, js_region, je_region)
964 integer,
intent(in) :: index, region_type
965 integer,
intent(in) :: is_region, ie_region, js_region, je_region
967 loaded_fields(index)%region_type = region_type
968 loaded_fields(index)%is_region = is_region
969 loaded_fields(index)%ie_region = ie_region
970 loaded_fields(index)%js_region = js_region
971 loaded_fields(index)%je_region = je_region
975 end subroutine set_override_region
979 integer,
intent(in) :: n
984 if (
associated(opened_files))
then
985 if (n <=
size(opened_files))
return
991 if(
Associated(ptr(i)%fileobj))
then
993 deallocate(ptr(i)%fileobj)
997 if (
associated(opened_files))
then
998 ptr(1:
size(opened_files)) = opened_files(:)
999 deallocate(opened_files)
1007 integer,
intent(in) :: n
1012 if (
associated(loaded_fields))
then
1013 if (n <=
size(loaded_fields))
return
1018 ptr(i)%fileobj => null()
1023 ptr(i)%domain = null_domain2d
1024 if (
ASSOCIATED(ptr(i)%time))
DEALLOCATE(ptr(i)%time, stat=ier)
1025 if (
ASSOCIATED(ptr(i)%start_time))
DEALLOCATE(ptr(i)%start_time, stat=ier)
1026 if (
ASSOCIATED(ptr(i)%end_time))
DEALLOCATE(ptr(i)%end_time, stat=ier)
1027 if (
ASSOCIATED(ptr(i)%period))
DEALLOCATE(ptr(i)%period, stat=ier)
1028 ptr(i)%modulo_time=.false.
1029 if (
ASSOCIATED(ptr(i)%domain_data))
DEALLOCATE(ptr(i)%domain_data, stat=ier)
1030 if (
ASSOCIATED(ptr(i)%ibuf))
DEALLOCATE(ptr(i)%ibuf, stat=ier)
1031 if (
ASSOCIATED(ptr(i)%src_data))
DEALLOCATE(ptr(i)%src_data, stat=ier)
1033 ptr(i)%domain_present=.false.
1034 ptr(i)%slope=1.0_r8_kind
1035 ptr(i)%intercept=0.0_r8_kind
1036 ptr(i)%isc=-1;ptr(i)%iec=-1
1037 ptr(i)%jsc=-1;ptr(i)%jec=-1
1039 if (
associated(loaded_fields))
then
1040 ptr(1:
size(loaded_fields)) = loaded_fields(:)
1041 deallocate(loaded_fields)
1051 integer,
dimension(:) :: buf
1061 if (buf(i) == indx)
then
1077 if (index .lt. 1 .or. index .gt. num_fields) &
1078 call mpp_error(fatal,
'invalid index in call to get_external_field_size')
1094 if (index .lt. 1 .or. index .gt. num_fields) &
1095 call mpp_error(fatal,
'invalid index in call to get_external_field_size')
1103 integer ,
intent(in) :: index
1108 if (index < 1.or.index > num_fields) &
1109 call mpp_error(fatal,
'invalid index in call to get_time_axis')
1111 n = min(
size(time),
size(loaded_fields(index)%time))
1113 time(1:n) = loaded_fields(index)%time(1:n)
1126 deallocate(loaded_fields(i)%time,loaded_fields(i)%start_time,loaded_fields(i)%end_time,&
1127 loaded_fields(i)%period,loaded_fields(i)%domain_data,loaded_fields(i)%mask,loaded_fields(i)%ibuf)
1128 if (
ASSOCIATED(loaded_fields(i)%src_data))
deallocate(loaded_fields(i)%src_data)
1129 loaded_fields(i)%domain = null_domain2d
1130 loaded_fields(i)%nbuf = 0
1131 loaded_fields(i)%slope = 0.0_r8_kind
1132 loaded_fields(i)%intercept = 0.0_r8_kind
1138 deallocate(opened_files(i)%fileobj)
1141 deallocate(loaded_fields)
1142 deallocate(opened_files)
1146 module_initialized = .false.
1150 #include "time_interp_external2_r4.fh"
1151 #include "time_interp_external2_r8.fh"
1153 #include "time_interp_external2_bridge_r4.fh"
1154 #include "time_interp_external2_bridge_r8.fh"
1156 end module time_interp_external2_mod
subroutine, public get_axis_cart(fileobj, axisname, cart)
Returns X,Y,Z or T cartesian attribute.
logical function, public get_axis_modulo(fileobj, axisname)
Checks if 'modulo' variable exists for a given axis.
logical function, public get_axis_modulo_times(fileobj, axisname, tbeg, tend)
Close a netcdf or domain file opened with open_file or open_virtual_file.
Opens a given netcdf or domain file.
Read data from a defined field in a file.
integer function, public check_nml_error(IOSTAT, NML_NAME)
Checks the iostat argument that is returned after reading a namelist and determines if the error code...
subroutine, public write_version_number(version, tag, unit)
Prints to the log file (or a specified unit) the version id string and tag name.
Added for mixed precision support. Updates force time_manager math to be done with kind=8 reals _wrap...
Subroutine for performing the horizontal interpolation between two grids.
Holds data pointers and metadata for horizontal interpolations, passed between the horiz_interp modul...
These routines retrieve the axis specifications associated with the compute domains....
These routines retrieve the axis specifications associated with the data domains. The domain is a der...
These routines retrieve the axis specifications associated with the global domains....
The domain2D type contains all the necessary information to define the global, compute and data domai...
integer function stdout()
This function returns the current standard fortran unit numbers for output.
integer function stdlog()
This function returns the current standard fortran unit numbers for log messages. Log messages,...
integer function mpp_npes()
Returns processor count for current pelist.
integer function mpp_pe()
Returns processor ID.
Perform parallel broadcasts.
subroutine, public time_interp_external_exit()
exit time_interp_external_mod. Close all open files and release storage
integer function, dimension(4), public get_external_field_size(index)
Returns size of field after call to init_external_field. Ordering is X/Y/Z/T. This call only makes se...
subroutine, public reset_src_data_region(index, is, ie, js, je)
Reallocates src_data for field from module level loaded_fields array.
subroutine load_record(field, rec, interp, is_in, ie_in, js_in, je_in, window_id_in)
load specified record from file Always loads in r8, casts down for horiz_interp if interp argument is...
subroutine, public time_interp_external_init()
Initialize the time_interp_external_mod module.
real(r8_kind) function, public get_external_field_missing(index)
return missing value
subroutine realloc_files(n)
reallocates array of fields, increasing its size
integer function, public init_external_field(file, fieldname, domain, desired_units, verbose, axis_names, axis_sizes, override, correct_leap_year_inconsistency, permit_calendar_conversion, use_comp_domain, ierr, nwindows, ignore_axis_atts, ongrid)
Initialize an external field. Buffer "num_io_buffers" (default=2) in memory to reduce memory allocati...
integer function, private find_buf_index(indx, buf)
simple linear search for given value in given list TODO should use better search if this list is bigg...
subroutine, public get_time_axis(index, time)
subroutine realloc_fields(n)
reallocates array of fields,increasing its size
subroutine load_record_0d(field, rec)
Given a initialized ext_fieldtype and record number, loads the given index into fieldsrc_data.
Provide data from external file interpolated to current model time. Data may be local to current proc...
Represents external fields.
Holds filename and file object.
Returns a weight and dates or indices for interpolating between two dates. The interface fraction_of_...
subroutine, public get_time(Time, seconds, days, ticks, err_msg)
Returns days and seconds ( < 86400 ) corresponding to a time. err_msg should be checked for any error...
subroutine, public get_date(time, year, month, day, hour, minute, second, tick, err_msg)
Gets the date for different calendar types. Given a time_interval, returns the corresponding date und...
integer function, public days_in_year(Time)
Returns the number of days in the calendar year corresponding to the date represented by time for the...
type(time_type) function, public increment_time(Time, seconds, days, ticks, err_msg, allow_neg_inc)
Increments a time by seconds and days.
integer function, public days_in_month(Time, err_msg)
Given a time, computes the corresponding date given the selected date time mapping algorithm.
integer function, public get_calendar_type()
Returns default calendar type for mapping from time to date.
Given an input date in year, month, days, etc., creates a time_type that represents this time interva...
Given some number of seconds and days, returns the corresponding time_type.
Type to represent amounts of time. Implemented as seconds and days to allow for larger intervals.