SCALE-RM
Data Types | Functions/Subroutines
scale_external_input Module Reference

module EXTERNAL INPUT More...

Functions/Subroutines

subroutine, public extin_setup (model)
 Setup. More...
 
subroutine, public extin_regist (basename, varname, axistype, enable_periodic_year, enable_periodic_month, enable_periodic_day, step_fixed, offset, defval, check_coordinates, step_limit, exist)
 Regist data. More...
 
subroutine extin_update_1d (var, varname, current_time, error)
 Read data. More...
 
subroutine extin_update_2d (var, varname, current_time, error)
 Read data. More...
 
subroutine extin_update_3d (var, varname, current_time, error)
 Read data. More...
 

Detailed Description

module EXTERNAL INPUT

Description
External file input module
Author
Team SCALE
NAMELIST
  • EXTITEM
    nametypedefault valuecomment
    BASENAME character(len=*)
    VARNAME character(len=*) item name
    AXISTYPE character(len=*)
    STEP_LIMIT integer limit number for reading data
    STEP_FIXED integer fixed step position to read
    ENABLE_PERIODIC_YEAR logical treat as yearly periodic data?
    ENABLE_PERIODIC_MONTH logical treat as yearly,monthly periodic data?
    ENABLE_PERIODIC_DAY logical treat as yearly,monthly,daily periodic data?
    OFFSET real(RP)
    DEFVAL real(RP)
    CHECK_COORDINATES logical

History Output
No history output

Function/Subroutine Documentation

◆ extin_setup()

subroutine, public scale_external_input::extin_setup ( character(len=*), intent(in)  model)

Setup.

Definition at line 160 of file scale_external_input.F90.

References scale_const::const_undef, extin_regist(), scale_external_input_rm::extin_rm_get_dims_1d(), scale_external_input_rm::extin_rm_get_dims_2d(), scale_external_input_rm::extin_rm_get_dims_3d(), scale_stdio::io_fid_conf, scale_stdio::io_fid_log, scale_stdio::io_fid_nml, scale_stdio::io_l, scale_stdio::io_nml, and scale_process::prc_mpistop().

Referenced by mod_rm_driver::scalerm().

160  use scale_process, only: &
162  use scale_const, only: &
163  undef => const_undef
164  use scale_external_input_rm, only: &
168  implicit none
169  character(len=*), intent(in) :: model
170 
171  character(len=H_LONG) :: basename
172  character(len=H_SHORT) :: varname
173  character(len=H_SHORT) :: axistype
174  integer :: step_limit ! limit number for reading data
175  integer :: step_fixed ! fixed step position to read
176  logical :: enable_periodic_year ! treat as yearly periodic data?
177  logical :: enable_periodic_month ! treat as yearly,monthly periodic data?
178  logical :: enable_periodic_day ! treat as yearly,monthly,daily periodic data?
179  real(RP) :: offset
180  real(RP) :: defval
181  logical :: check_coordinates
182 
183  namelist /extitem/ &
184  basename, &
185  varname, &
186  axistype, &
187  step_limit, &
188  step_fixed, &
189  enable_periodic_year, &
190  enable_periodic_month, &
191  enable_periodic_day, &
192  offset, &
193  defval, &
194  check_coordinates
195 
196  integer :: count
197  integer :: ierr
198  !---------------------------------------------------------------------------
199 
200  if( io_l ) write(io_fid_log,*)
201  if( io_l ) write(io_fid_log,*) '++++++ Module[EXTIN] / Categ[ATMOS-RM IO] / Origin[SCALElib]'
202 
203  select case ( model )
204  case ( 'RM' )
205  extin_get_dims_1d => extin_rm_get_dims_1d
206  extin_get_dims_2d => extin_rm_get_dims_2d
207  extin_get_dims_3d => extin_rm_get_dims_3d
208 ! case ( 'GM' )
209  case default
210  write(*,*) 'xxx EXTIN is not support for the model: ', trim(model)
211  call prc_mpistop
212  end select
213 
214  ! count external data from namelist
215  rewind(io_fid_conf)
216  do count = 1, extin_item_limit
217 
218  ! set default
219  step_limit = extin_step_limit
220  basename = ''
221  varname = ''
222  axistype = ''
223  step_fixed = -1
224  enable_periodic_year = .false.
225  enable_periodic_month = .false.
226  enable_periodic_day = .false.
227  offset = 0.0_rp
228  defval = undef
229  check_coordinates = .false.
230 
231  ! read namelist
232  read(io_fid_conf,nml=extitem,iostat=ierr)
233 
234  if( ierr < 0 ) then !--- no more items
235  exit
236  elseif( ierr > 0 ) then !--- fatal error
237  write(*,*) 'xxx Not appropriate names in namelist EXTITEM. Check!', count
238  call prc_mpistop
239  endif
240 
241  if( io_nml .AND. io_fid_nml /= io_fid_log ) write(io_fid_nml,nml=extitem)
242 
243  call extin_regist( basename, & ! [IN]
244  varname, & ! [IN]
245  axistype, & ! [IN]
246  enable_periodic_year, & ! [IN]
247  enable_periodic_month, & ! [IN]
248  enable_periodic_day, & ! [IN]
249  step_fixed, & ! [IN]
250  offset, & ! [IN]
251  defval, & ! [IN]
252  check_coordinates, & ! [IN]
253  step_limit ) ! [IN]
254 
255  enddo
256 
257  return
subroutine, public prc_mpistop
Abort MPI.
subroutine, public extin_rm_get_dims_2d(dim1_max, dim1_S, dim1_E, dim2_max, dim2_S, dim2_E, transpose, varname, axistype)
subroutine, public extin_rm_get_dims_3d(dim1_max, dim1_S, dim1_E, dim2_max, dim2_S, dim2_E, dim3_max, dim3_S, dim3_E, transpose, varname, axistype)
module PROCESS
module CONSTANT
Definition: scale_const.F90:14
subroutine, public extin_rm_get_dims_1d(dim1_max, dim1_S, dim1_E, varname, axistype)
get_dims
integer, public io_fid_conf
Config file ID.
Definition: scale_stdio.F90:55
Here is the call graph for this function:
Here is the caller graph for this function:

◆ extin_regist()

subroutine, public scale_external_input::extin_regist ( character(len=*), intent(in)  basename,
character(len=*), intent(in)  varname,
character(len=*), intent(in)  axistype,
logical, intent(in)  enable_periodic_year,
logical, intent(in)  enable_periodic_month,
logical, intent(in)  enable_periodic_day,
integer, intent(in)  step_fixed,
real(rp), intent(in)  offset,
real(rp), intent(in)  defval,
logical, intent(in), optional  check_coordinates,
integer, intent(in), optional  step_limit,
logical, intent(out), optional  exist 
)

Regist data.

Definition at line 275 of file scale_external_input.F90.

References scale_calendar::calendar_adjust_daysec(), scale_calendar::calendar_cfunits2sec(), scale_calendar::calendar_combine_daysec(), scale_calendar::calendar_date2daysec(), scale_calendar::calendar_daysec2date(), gtool_file_h::file_fread, gtool_file::fileopen(), scale_calendar::i_day, scale_calendar::i_month, scale_calendar::i_year, scale_stdio::io_fid_log, scale_stdio::io_l, scale_process::prc_mpistop(), scale_process::prc_myrank, scale_time::time_nowdaysec, scale_time::time_offset_year, and scale_time::time_startdaysec.

Referenced by scale_atmos_phy_rd_offline::atmos_phy_rd_offline_setup(), extin_setup(), and scale_ocean_phy_file::ocean_phy_file_setup().

275  use gtool_file_h, only: &
276  file_fread
277  use gtool_file, only: &
278  fileopen, &
279  filegetalldatainfo, &
280  fileread
281  use scale_process, only: &
282  prc_myrank, &
284  use scale_fileio, only: &
285  fileio_check_coordinates
286  use scale_calendar, only: &
292  i_year, &
293  i_month, &
294  i_day
295  use scale_time, only: &
297  time_nowdaysec, &
299  implicit none
300 
301  character(len=*) , intent(in) :: basename
302  character(len=*) , intent(in) :: varname
303  character(len=*) , intent(in) :: axistype
304  integer , intent(in) :: step_fixed ! fixed step position to read
305  logical , intent(in) :: enable_periodic_year ! treat as yearly periodic data?
306  logical , intent(in) :: enable_periodic_month ! treat as yearly,monthly periodic data?
307  logical , intent(in) :: enable_periodic_day ! treat as yearly,monthly,daily periodic data?
308  real(RP) , intent(in) :: offset
309  real(RP) , intent(in) :: defval
310  logical, optional, intent(in) :: check_coordinates
311  integer, optional, intent(in) :: step_limit ! limit number for reading data
312  logical, optional, intent(out) :: exist
313 
314  integer :: step_nmax
315  character(len=H_LONG) :: description
316  character(len=H_SHORT) :: unit
317  integer :: datatype
318  integer :: dim_rank
319  character(len=H_SHORT) :: dim_name (3)
320  integer :: dim_size (3)
321  real(DP) :: time_start(EXTIN_step_limit)
322  real(DP) :: time_end (EXTIN_step_limit)
323  character(len=H_MID) :: time_units
324 
325  integer :: datadate(6)
326  real(DP) :: datasubsec
327  integer :: dataday
328  real(DP) :: datasec
329  integer :: offset_year
330 
331  integer :: dim1_max, dim1_S, dim1_E
332  integer :: dim2_max, dim2_S, dim2_E
333  integer :: dim3_max, dim3_S, dim3_E
334 
335  integer :: step_limit_
336 
337  integer :: fid
338  integer :: nid, n
339  !---------------------------------------------------------------------------
340 
341  if ( present(step_limit) ) then
342  if ( step_limit > 0 ) then
343  step_limit_ = step_limit
344  else
345  step_limit_ = extin_step_limit
346  end if
347  else
348  step_limit_ = extin_step_limit
349  end if
350 
351  do nid = 1, extin_item_count
352  if ( extin_item(nid)%varname == varname ) then
353  write(*,*) 'xxx Data is already registered! basename,varname = ', trim(basename), ', ', trim(varname)
354  call prc_mpistop
355  end if
356  end do
357 
358  extin_item_count = extin_item_count + 1
359 
360  if ( extin_item_count > extin_item_limit ) then
361  write(*,*) 'xxx Number of EXT data exceedes the limit', extin_item_count, extin_item_limit
362  call prc_mpistop
363  end if
364 
365  call fileopen( fid, basename, file_fread, myrank=prc_myrank )
366 
367  ! read from file
368  call filegetalldatainfo( step_limit_, & ! [IN]
369  extin_dim_limit, & ! [IN]
370  fid, & ! [IN]
371  varname, & ! [IN]
372  step_nmax, & ! [OUT]
373  description, & ! [OUT]
374  unit, & ! [OUT]
375  datatype, & ! [OUT]
376  dim_rank, & ! [OUT]
377  dim_name(:), & ! [OUT]
378  dim_size(:), & ! [OUT]
379  time_start(1:step_limit_), & ! [OUT]
380  time_end(1:step_limit_), & ! [OUT]
381  time_units ) ! [OUT]
382 
383  if ( step_nmax == 0 ) then
384  if ( present(exist) ) then
385  exist = .false.
386  return
387  end if
388  write(*,*) 'xxx Data not found! basename,varname = ', trim(basename), ', ', trim(varname)
389  call prc_mpistop
390  endif
391 
392  if ( present(exist) ) exist = .true.
393 
394  do n = dim_rank+1, 3
395  dim_size(n) = 1
396  end do
397 
398  nid = extin_item_count
399 
400  ! setup item
401  extin_item(nid)%fid = fid
402  extin_item(nid)%varname = varname
403  extin_item(nid)%dim_size(:) = dim_size(:)
404  extin_item(nid)%step_num = step_nmax
405 
406  if ( enable_periodic_day ) then
407  extin_item(nid)%flag_periodic = i_periodic_day
408  elseif( enable_periodic_month ) then
409  extin_item(nid)%flag_periodic = i_periodic_month
410  elseif( enable_periodic_year ) then
411  extin_item(nid)%flag_periodic = i_periodic_year
412  else
413  extin_item(nid)%flag_periodic = 0
414  endif
415 
416  allocate( extin_item(nid)%value(dim_size(1),dim_size(2),dim_size(3),2) )
417  extin_item(nid)%value(:,:,:,:) = defval
418  extin_item(nid)%offset = offset
419 
420  allocate( extin_item(nid)%time(step_nmax) )
421  extin_item(nid)%time(:) = 0.0_dp
422 
423  do n = 1, extin_item(nid)%step_num
424  extin_item(nid)%time(n) = calendar_cfunits2sec( time_end(n), time_units, time_offset_year, time_startdaysec )
425  enddo
426 
427  if ( extin_item(nid)%step_num == 1 ) then
428 
429  extin_item(nid)%fixed_step = .true.
430  extin_item(nid)%data_steppos(i_prev) = 1
431  extin_item(nid)%data_steppos(i_next) = 1
432 
433  else if ( step_fixed > 0 ) then ! fixed time step mode
434 
435  extin_item(nid)%fixed_step = .true.
436  extin_item(nid)%data_steppos(i_prev) = step_fixed
437  extin_item(nid)%data_steppos(i_next) = step_fixed
438 
439  else
440 
441  extin_item(nid)%fixed_step = .false.
442 
443  ! seek start position
444  extin_item(nid)%data_steppos(i_next) = 1
445  do n = 1, extin_item(nid)%step_num
446  if ( extin_item(nid)%time(n) > time_nowdaysec ) exit
447  extin_item(nid)%data_steppos(i_next) = n + 1
448  enddo
449 
450  extin_item(nid)%data_steppos(i_prev) = extin_item(nid)%data_steppos(i_next) - 1
451 
452  if ( extin_item(nid)%flag_periodic > 0 ) then ! periodic time step mode
453 
454  if ( extin_item(nid)%data_steppos(i_next) == 1 ) then ! between first-1 and first
455 
456  ! first-1 = last
457  extin_item(nid)%data_steppos(i_prev) = extin_item(nid)%step_num
458 
459  elseif( extin_item(nid)%data_steppos(i_next) == extin_item(nid)%step_num+1 ) then ! between last and last+1
460 
461  ! last+1 = first
462  extin_item(nid)%data_steppos(i_next) = 1
463 
464  ! update data time in periodic condition
465  do n = 1, extin_item(nid)%step_num
466  dataday = 0
467  datasec = extin_item(nid)%time(n)
468  offset_year = 0
469  call calendar_adjust_daysec( dataday, datasec ) ! [INOUT]
470 
471  call calendar_daysec2date( datadate(:), & ! [OUT]
472  datasubsec, & ! [OUT]
473  dataday, & ! [IN]
474  datasec, & ! [IN]
475  offset_year ) ! [IN]
476 
477  if ( extin_item(nid)%flag_periodic == i_periodic_day ) then
478  datadate(i_day) = datadate(i_day) + 1
479  elseif( extin_item(nid)%flag_periodic == i_periodic_month ) then
480  datadate(i_month) = datadate(i_month) + 1
481  elseif( extin_item(nid)%flag_periodic == i_periodic_year ) then
482  datadate(i_year) = datadate(i_year) + 1
483  endif
484 
485  call calendar_date2daysec( dataday, & ! [OUT]
486  datasec, & ! [OUT]
487  datadate(:), & ! [IN]
488  datasubsec, & ! [IN]
489  offset_year ) ! [IN]
490 
491  extin_item(nid)%time(n) = calendar_combine_daysec( dataday, datasec )
492  enddo
493 
494  if( io_l ) write(io_fid_log,*) '*** data time is updated.'
495  endif
496 
497  else ! normal mode
498 
499  if ( extin_item(nid)%data_steppos(i_next) == 1 &
500  .OR. extin_item(nid)%data_steppos(i_next) == extin_item(nid)%step_num+1 ) then
501  write(*,*) 'xxx Current time is out of period of external data! ', trim(varname)
502  call prc_mpistop
503  endif
504 
505  endif
506 
507  endif
508 
509  !--- read first data
510  if( io_l ) write(io_fid_log,'(1x,A,A15)') '*** Initial read of external data : ', trim(varname)
511 
512  if ( dim_size(1) >= 1 &
513  .AND. dim_size(2) == 1 &
514  .AND. dim_size(3) == 1 ) then ! 1D
515 
516  call extin_get_dims_1d( &
517  dim1_max, dim1_s, dim1_e, & ! (out)
518  varname, axistype ) ! (in)
519 
520  extin_item(nid)%ndim = 1
521  extin_item(nid)%transpose = .false.
522  allocate( extin_item(nid)%dim_start(1) )
523  extin_item(nid)%dim_start(1) = dim1_s
524 
525  if ( dim1_max /= dim_size(1) ) then
526  write(*,*) 'xxx data length does not match! ', trim(axistype), ' item:', trim(varname)
527  write(*,*) 'xxx dim 1 (data,requested) : ', dim_size(1), dim1_max
528  call prc_mpistop
529  endif
530 
531  ! read prev
532  if( io_l ) write(io_fid_log,'(1x,A,A,A,I4,A)') &
533  '*** Read 1D var : ', trim(extin_item(nid)%varname), &
534  ' (step= ', extin_item(nid)%data_steppos(i_prev), ')'
535  call fileread( extin_item(nid)%value(:,1,1,i_prev), &
536  extin_item(nid)%fid, &
537  extin_item(nid)%varname, &
538  extin_item(nid)%data_steppos(i_prev) )
539  ! read next
540  if( io_l ) write(io_fid_log,'(1x,A,A,A,I4,A)') &
541  '*** Read 1D var : ', trim(extin_item(nid)%varname), &
542  ' (step= ', extin_item(nid)%data_steppos(i_next), ')'
543  call fileread( extin_item(nid)%value(:,1,1,i_next), &
544  extin_item(nid)%fid, &
545  extin_item(nid)%varname, &
546  extin_item(nid)%data_steppos(i_next) )
547 
548  elseif ( dim_size(1) >= 1 &
549  .AND. dim_size(2) > 1 &
550  .AND. dim_size(3) == 1 ) then ! 2D
551 
552  call extin_get_dims_2d( &
553  dim1_max, dim1_s, dim1_e, & ! (out)
554  dim2_max, dim2_s, dim2_e, & ! (out)
555  extin_item(nid)%transpose, & ! (out)
556  varname, axistype ) ! (in)
557 
558  extin_item(nid)%ndim = 2
559  allocate( extin_item(nid)%dim_start(2) )
560  extin_item(nid)%dim_start(1) = dim1_s
561  extin_item(nid)%dim_start(2) = dim2_s
562 
563  if ( dim1_max /= dim_size(1) &
564  .OR. dim2_max /= dim_size(2) ) then
565  write(*,*) 'xxx data length does not match! ', trim(axistype), ' item:', trim(varname)
566  write(*,*) 'xxx dim 1 (data,requested) : ', dim_size(1), dim1_max
567  write(*,*) 'xxx dim 2 (data,requested) : ', dim_size(2), dim2_max
568  call prc_mpistop
569  endif
570 
571  ! read prev
572  if( io_l ) write(io_fid_log,'(1x,A,A,A,I4,A)') &
573  '*** Read 2D var : ', trim(extin_item(nid)%varname), &
574  ' (step= ', extin_item(nid)%data_steppos(i_prev), ')'
575  call fileread( extin_item(nid)%value(:,:,1,i_prev), &
576  extin_item(nid)%fid, &
577  extin_item(nid)%varname, &
578  extin_item(nid)%data_steppos(i_prev) )
579  ! read next
580  if( io_l ) write(io_fid_log,'(1x,A,A,A,I4,A)') &
581  '*** Read 2D var : ', trim(extin_item(nid)%varname), &
582  ' (step= ', extin_item(nid)%data_steppos(i_next), ')'
583  call fileread( extin_item(nid)%value(:,:,1,i_next), &
584  extin_item(nid)%fid, &
585  extin_item(nid)%varname, &
586  extin_item(nid)%data_steppos(i_next) )
587 
588  elseif ( dim_size(1) >= 1 &
589  .AND. dim_size(2) > 1 &
590  .AND. dim_size(3) > 1 ) then ! 3D
591 
592  call extin_get_dims_3d( &
593  dim1_max, dim1_s, dim1_e, & ! (out)
594  dim2_max, dim2_s, dim2_e, & ! (out)
595  dim3_max, dim3_s, dim3_e, & ! (out)
596  extin_item(nid)%transpose, & ! (out)
597  varname, axistype ) ! (in)
598 
599  extin_item(nid)%ndim = 3
600  allocate( extin_item(nid)%dim_start(3) )
601  extin_item(nid)%dim_start(1) = dim1_s
602  extin_item(nid)%dim_start(2) = dim2_s
603  extin_item(nid)%dim_start(3) = dim3_s
604 
605  if ( dim1_max /= dim_size(1) &
606  .OR. dim2_max /= dim_size(2) &
607  .OR. dim3_max /= dim_size(3) ) then
608  write(*,*) 'xxx data length does not match! ', trim(axistype), ' item:', trim(varname)
609  write(*,*) 'xxx dim 1 (data,requested) : ', dim_size(1), dim1_max
610  write(*,*) 'xxx dim 2 (data,requested) : ', dim_size(2), dim2_max
611  write(*,*) 'xxx dim 3 (data,requested) : ', dim_size(3), dim3_max
612  call prc_mpistop
613  endif
614 
615  ! read prev
616  if( io_l ) write(io_fid_log,'(1x,A,A,A,I4,A)') &
617  '*** Read 3D var : ', trim(extin_item(nid)%varname), &
618  ' (step= ', extin_item(nid)%data_steppos(i_prev), ')'
619  call fileread( extin_item(nid)%value(:,:,:,i_prev), &
620  extin_item(nid)%fid, &
621  extin_item(nid)%varname, &
622  extin_item(nid)%data_steppos(i_prev) )
623 
624  ! read next
625  if( io_l ) write(io_fid_log,'(1x,A,A,A,I4,A)') &
626  '*** Read 3D var : ', trim(extin_item(nid)%varname), &
627  ' (step= ', extin_item(nid)%data_steppos(i_next), ')'
628  call fileread( extin_item(nid)%value(:,:,:,i_next), &
629  extin_item(nid)%fid, &
630  extin_item(nid)%varname, &
631  extin_item(nid)%data_steppos(i_next) )
632 
633  else
634  write(*,*) 'xxx Unexpected dimsize: ', dim_size(:)
635  call prc_mpistop
636  endif
637 
638  if ( present(check_coordinates) ) then
639  if ( check_coordinates ) &
640  call fileio_check_coordinates( fid, &
641  atmos = extin_item(nid)%ndim==3, &
642  transpose = extin_item(nid)%transpose )
643  end if
644 
645  return
module GTOOL_FILE
Definition: gtool_file.f90:17
subroutine, public prc_mpistop
Abort MPI.
real(dp) function, public calendar_combine_daysec(absday, abssec)
Combine day and second.
real(dp), public time_nowdaysec
second of current time [sec]
Definition: scale_time.F90:69
real(dp), public time_startdaysec
second of start time [sec]
Definition: scale_time.F90:74
module FILE I/O (netcdf)
subroutine, public calendar_adjust_daysec(absday, abssec)
Adjust day and second.
real(dp) function, public calendar_cfunits2sec(cftime, cfunits, offset_year, startdaysec)
Convert time in units of the CF convention to second.
integer, public time_offset_year
time offset [year]
Definition: scale_time.F90:73
module TIME
Definition: scale_time.F90:15
module PROCESS
integer, parameter, public file_fread
integer, public prc_myrank
process num in local communicator
subroutine, public calendar_daysec2date(ymdhms, subsec, absday, abssec, offset_year)
Convert from gregorian date to absolute day/second.
module CALENDAR
module FILE I/O HEADER
subroutine, public fileopen(fid, basename, mode, single, comm, myrank)
Definition: gtool_file.f90:495
subroutine, public calendar_date2daysec(absday, abssec, ymdhms, subsec, offset_year)
Convert from gregorian date to absolute day/second.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ extin_update_1d()

subroutine scale_external_input::extin_update_1d ( real(rp), dimension(:), intent(out)  var,
character(len=*), intent(in)  varname,
real(dp), intent(in)  current_time,
logical, intent(out)  error 
)

Read data.

Definition at line 655 of file scale_external_input.F90.

References scale_stdio::io_fid_log, scale_stdio::io_l, and scale_process::prc_mpistop().

655  use gtool_file, only: &
656  fileread
657  use scale_process, only: &
659  implicit none
660 
661  real(RP), intent(out) :: var(:) ! variable
662  character(len=*), intent(in) :: varname ! item name
663  real(DP), intent(in) :: current_time ! current time
664  logical, intent(out) :: error ! error code
665 
666  integer :: nid
667  real(RP) :: weight
668  logical :: do_readfile
669 
670  integer :: n
671  integer :: n1
672  integer :: nn1
673  !---------------------------------------------------------------------------
674 
675  error = .true.
676 
677  ! searching the data ID
678  nid = -1
679  do n = 1, extin_item_count
680  if( varname == extin_item(n)%varname ) nid = n
681  enddo
682 
683  if ( nid == 0 ) then
684  if( io_l ) write(io_fid_log,*) 'xxx Variable was not registered: ', trim(varname)
685  return
686  end if
687 
688  if ( extin_item(nid)%ndim /= 1 ) then
689  write(*,*) 'xxx Data is not 1D var: ', trim(extin_item(nid)%varname)
690  call prc_mpistop
691  end if
692 
693  call extin_time_advance( nid, & ! [IN]
694  current_time, & ! [IN]
695  weight, & ! [OUT]
696  do_readfile ) ! [OUT]
697 
698  if ( do_readfile ) then
699 
700  if( io_l ) write(io_fid_log,'(1x,A,A,A,I4,A)') &
701  '*** Read 1D var : ', trim(extin_item(nid)%varname), &
702  ' (step= ', extin_item(nid)%data_steppos(i_next), ')'
703 
704  ! next -> prev
705  extin_item(nid)%value(:,:,:,i_prev) = extin_item(nid)%value(:,:,:,i_next)
706 
707  ! read next
708  call fileread( extin_item(nid)%value(:,1,1,i_next), & ! [OUT]
709  extin_item(nid)%fid, & ! [IN]
710  extin_item(nid)%varname, & ! [IN]
711  extin_item(nid)%data_steppos(i_next) )
712  endif
713 
714  ! store data with weight
715  do n1 = 1, extin_item(nid)%dim_size(1)
716  nn1 = n1 + extin_item(nid)%dim_start(1) - 1
717 
718  var(nn1) = ( 1.0_rp-weight ) * extin_item(nid)%value(n1,1,1,i_prev) &
719  + ( weight ) * extin_item(nid)%value(n1,1,1,i_next)
720  enddo
721 
722  error = .false.
723 
724  return
module GTOOL_FILE
Definition: gtool_file.f90:17
subroutine, public prc_mpistop
Abort MPI.
module PROCESS
Here is the call graph for this function:

◆ extin_update_2d()

subroutine scale_external_input::extin_update_2d ( real(rp), dimension(:,:), intent(out)  var,
character(len=*), intent(in)  varname,
real(dp), intent(in)  current_time,
logical, intent(out)  error 
)

Read data.

Definition at line 734 of file scale_external_input.F90.

References scale_stdio::io_fid_log, scale_stdio::io_l, and scale_process::prc_mpistop().

734  use gtool_file, only: &
735  fileread
736  use scale_process, only: &
738  implicit none
739 
740  real(RP), intent(out) :: var(:,:) ! variable
741  character(len=*), intent(in) :: varname ! item name
742  real(DP), intent(in) :: current_time ! current time
743  logical, intent(out) :: error ! error code
744 
745  integer :: nid
746  real(RP) :: weight
747  logical :: do_readfile
748 
749  integer :: n
750  integer :: n1, n2
751  integer :: nn1, nn2
752  !---------------------------------------------------------------------------
753 
754  error = .true.
755 
756  ! searching the data ID
757  nid = -1
758  do n = 1, extin_item_count
759  if( varname == extin_item(n)%varname ) nid = n
760  enddo
761 
762  if ( nid == 0 ) then
763  if( io_l ) write(io_fid_log,*) 'xxx variable was not registered: ', trim(varname)
764  return
765  end if
766 
767  call extin_time_advance( nid, & ! [IN]
768  current_time, & ! [IN]
769  weight, & ! [OUT]
770  do_readfile ) ! [OUT]
771 
772  if ( do_readfile ) then
773 
774  if( io_l ) write(io_fid_log,'(1x,A,A,A,I4,A)') &
775  '*** Read 2D var : ', trim(extin_item(nid)%varname), &
776  ' (step= ', extin_item(nid)%data_steppos(i_next), ')'
777 
778  ! next -> prev
779  extin_item(nid)%value(:,:,:,i_prev) = extin_item(nid)%value(:,:,:,i_next)
780 
781  ! read next
782  call fileread( extin_item(nid)%value(:,:,1,i_next), & ! [OUT]
783  extin_item(nid)%fid, & ! [IN]
784  extin_item(nid)%varname, & ! [IN]
785  extin_item(nid)%data_steppos(i_next) ) ! [IN]
786  endif
787 
788  if ( extin_item(nid)%transpose ) then
789  ! store data with weight (x,z)->(z,x)
790  do n1 = 1, extin_item(nid)%dim_size(1)
791  nn1 = n1 + extin_item(nid)%dim_start(1) - 1
792  do n2 = 1, extin_item(nid)%dim_size(2)
793 
794  nn2 = n2 + extin_item(nid)%dim_start(2) - 1
795  var(nn2,nn1) = ( 1.0_rp-weight ) * extin_item(nid)%value(n1,n2,1,i_prev) &
796  + ( weight ) * extin_item(nid)%value(n1,n2,1,i_next)
797  enddo
798  enddo
799  else
800  ! store data with weight
801  do n2 = 1, extin_item(nid)%dim_size(2)
802  nn2 = n2 + extin_item(nid)%dim_start(2) - 1
803  do n1 = 1, extin_item(nid)%dim_size(1)
804  nn1 = n1 + extin_item(nid)%dim_start(1) - 1
805 
806  var(nn1,nn2) = ( 1.0_rp-weight ) * extin_item(nid)%value(n1,n2,1,i_prev) &
807  + ( weight ) * extin_item(nid)%value(n1,n2,1,i_next)
808  enddo
809  enddo
810  end if
811 
812  error = .false.
813 
814  return
module GTOOL_FILE
Definition: gtool_file.f90:17
subroutine, public prc_mpistop
Abort MPI.
module PROCESS
Here is the call graph for this function:

◆ extin_update_3d()

subroutine scale_external_input::extin_update_3d ( real(rp), dimension(:,:,:), intent(out)  var,
character(len=*), intent(in)  varname,
real(dp), intent(in)  current_time,
logical, intent(out)  error 
)

Read data.

Definition at line 824 of file scale_external_input.F90.

References scale_calendar::calendar_adjust_daysec(), scale_calendar::calendar_combine_daysec(), scale_calendar::calendar_date2daysec(), scale_calendar::calendar_daysec2date(), scale_calendar::i_day, scale_calendar::i_month, scale_calendar::i_year, scale_stdio::io_fid_log, scale_stdio::io_l, and scale_process::prc_mpistop().

824  use gtool_file, only: &
825  fileread
826  use scale_process, only: &
828  implicit none
829 
830  real(RP), intent(out) :: var(:,:,:) ! variable
831  character(len=*), intent(in) :: varname ! item name
832  real(DP), intent(in) :: current_time ! current time
833  logical, intent(out) :: error ! error code
834 
835  integer :: nid
836  real(RP) :: weight
837  logical :: do_readfile
838 
839  integer :: n
840  integer :: n1, n2, n3
841  integer :: nn1, nn2, nn3
842  !---------------------------------------------------------------------------
843 
844  error = .true.
845 
846  ! searching the data ID
847  nid = -1
848  do n = 1, extin_item_count
849  if( varname == extin_item(n)%varname ) nid = n
850  enddo
851 
852  if ( nid == 0 ) then
853  if( io_l ) write(io_fid_log,*) 'xxx variable was not registered: ', trim(varname)
854  return
855  end if
856 
857  call extin_time_advance( nid, & ! [IN]
858  current_time, & ! [IN]
859  weight, & ! [OUT]
860  do_readfile ) ! [OUT]
861 
862  if ( do_readfile ) then
863 
864  if( io_l ) write(io_fid_log,'(1x,A,A,A,I4,A)') &
865  '*** Read 3D var : ', trim(extin_item(nid)%varname), &
866  ' (step= ', extin_item(nid)%data_steppos(i_next), ')'
867 
868  ! next -> prev
869  extin_item(nid)%value(:,:,:,i_prev) = extin_item(nid)%value(:,:,:,i_next)
870 
871  ! read next
872  call fileread( extin_item(nid)%value(:,:,:,i_next), & ! [OUT]
873  extin_item(nid)%fid, & ! [IN]
874  extin_item(nid)%varname, & ! [IN]
875  extin_item(nid)%data_steppos(i_next) ) ! [IN]
876  endif
877 
878  if ( extin_item(nid)%transpose ) then
879  ! store data with weight (x,y,z)->(z,x,y)
880  do n2 = 1, extin_item(nid)%dim_size(2)
881  nn2 = n2 + extin_item(nid)%dim_start(2) - 1
882  do n1 = 1, extin_item(nid)%dim_size(1)
883  nn1 = n1 + extin_item(nid)%dim_start(1) - 1
884  do n3 = 1, extin_item(nid)%dim_size(3)
885  nn3 = n3 + extin_item(nid)%dim_start(3) - 1
886 
887  var(nn3,nn1,nn2) = ( 1.0_rp-weight ) * extin_item(nid)%value(n1,n2,n3,i_prev) &
888  + ( weight ) * extin_item(nid)%value(n1,n2,n3,i_next)
889  enddo
890  enddo
891  enddo
892  else
893  ! store data with weight (z,x,y)->(z,x,y)
894  do n3 = 1, extin_item(nid)%dim_size(3)
895  nn3 = n3 + extin_item(nid)%dim_start(3) - 1
896  do n2 = 1, extin_item(nid)%dim_size(2)
897  nn2 = n2 + extin_item(nid)%dim_start(2) - 1
898  do n1 = 1, extin_item(nid)%dim_size(1)
899  nn1 = n1 + extin_item(nid)%dim_start(1) - 1
900 
901  var(nn1,nn2,nn3) = ( 1.0_rp-weight ) * extin_item(nid)%value(n1,n2,n3,i_prev) &
902  + ( weight ) * extin_item(nid)%value(n1,n2,n3,i_next)
903  enddo
904  enddo
905  enddo
906  end if
907 
908  error = .false.
909 
910  return
module GTOOL_FILE
Definition: gtool_file.f90:17
subroutine, public prc_mpistop
Abort MPI.
module PROCESS
Here is the call graph for this function: