Class | read_time_series |
In: |
io/read_time_series.f90
|
Note that Japanese and English are described in parallel.
海表面温度, O3 量などの時系列データを NetCDF ファイルから読み込む.
Reading time series data, such as sea surface temperature, O3, and so on, from NetCDF file.
!$ ! GroundFileGet : | 地表面データファイルの入力 |
!$ !—
!$ ! GroundFileOpen : | 地表面データファイルのオープン |
!$ ! GroundFileOutput : | 地表面データファイルへのデータ出力 |
!$ ! GroundFileClose : | 地表面データファイルのクローズ |
++
———— : | ———— |
!$ ! GroundFileGet : | Input ground data file |
!$ !—
!$ ! GroundFileOpen : | Open ground data file |
!$ ! GroundFileOutput : | Data output to ground data file |
!$ ! GroundFileClose : | Close ground data file |
++
!$ ! NAMELIST#ground_file_io_nml
Subroutine : |
subroutine ReadTimeSeriesInit ! ファイル入出力補助 ! File I/O support ! use dc_iounit, only: FileOpen ! NAMELIST ファイル入力に関するユーティリティ ! Utilities for NAMELIST file input ! use namelist_util, only: namelist_filename, NmlutilMsg, NmlutilAryValid ! 宣言文 ; Declaration statements ! logical :: flag_mpi_init integer:: unit_nml ! NAMELIST ファイルオープン用装置番号. ! Unit number for NAMELIST file open integer:: iostat_nml ! NAMELIST 読み込み時の IOSTAT. ! IOSTAT of NAMELIST read ! NAMELIST 変数群 ! NAMELIST group name ! !!$ namelist /read_time_series_nml/ & !!$ & ! ! デフォルト値については初期化手続 "read_time_series#ReadTimeSeriesInit" ! のソースコードを参照のこと. ! ! Refer to source codes in the initialization procedure ! "read_time_series#ReadTimeSeriesInit" for the default values. ! if ( read_time_series_inited ) return ! デフォルト値の設定 ! Default values settings ! ! NAMELIST の読み込み ! NAMELIST is input ! !!$ if ( trim(namelist_filename) /= '' ) then !!$ call FileOpen( unit_nml, & ! (out) !!$ & namelist_filename, mode = 'r' ) ! (in) !!$ !!$ rewind( unit_nml ) !!$ read( unit_nml, & ! (in) !!$ & nml = read_time_series_nml, & ! (out) !!$ & iostat = iostat_nml ) ! (out) !!$ close( unit_nml ) !!$ !!$ call NmlutilMsg( iostat_nml, module_name ) ! (in) !!$ end if ! 印字 ; Print ! call MessageNotify( 'M', module_name, '----- Initialization Messages -----' ) !!$ call MessageNotify( 'M', module_name, 'InFileName = %c', c1 = trim(InFileName) ) call MessageNotify( 'M', module_name, '-- version = %c', c1 = trim(version) ) read_time_series_inited = .true. end subroutine ReadTimeSeriesInit
Subroutine : | |||
FileName : | character(*), intent(in )
| ||
VarName : | character(*), intent(in )
| ||
xy_Var(0:imax-1, 1:jmax) : | real(DP) , intent(inout)
| ||
keyword : | character(*), intent(in ) |
Alias for SetValFromTimeSeriesWrapper2D
Subroutine : | |||
FileName : | character(*), intent(in )
| ||
VarName : | character(*), intent(in )
| ||
xyz_Press(0:imax-1, 1:jmax, 1:kmax) : | real(DP) , intent(in )
| ||
xyz_Var(0:imax-1, 1:jmax, 1:kmax) : | real(DP) , intent(inout)
| ||
keyword : | character(*), intent(in ) |
Alias for SetValFromTimeSeriesWrapper3D
Subroutine : | |||
FileName : | character(*), intent(in )
| ||
VarName : | character(*), intent(in )
| ||
xy_Var(0:imax-1, 1:jmax) : | real(DP) , intent(inout)
| ||
keyword : | character(*), intent(in ) |
subroutine SetValFromTimeSeriesWrapper2D( FileName, VarName, xy_Var, keyword ) ! 宣言文 ; Declaration statements ! implicit none character(*), intent(in ):: FileName ! ファイル名. ! File name character(*), intent(in ):: VarName ! 変数名. ! Variable name real(DP) , intent(inout):: xy_Var(0:imax-1, 1:jmax) ! 地表面温度. ! Surface temperature character(*), intent(in ):: keyword ! 初期化確認 ! Initialization check ! if ( .not. read_time_series_inited ) then call MessageNotify( 'E', module_name, 'This module has not been initialized.' ) end if if ( keyword == 'SST' ) then call SetValuesFromTimeSeries( TSDataInfoSST, FileName, VarName, xy_Var = xy_Var ) else if ( keyword == 'SIC' ) then call SetValuesFromTimeSeries( TSDataInfoSIC, FileName, VarName, xy_Var = xy_Var ) else if ( keyword == 'DOD' ) then call SetValuesFromTimeSeries( TSDataInfoDOD, FileName, VarName, xy_Var = xy_Var ) else stop 'Unsupported keyword' end if end subroutine SetValFromTimeSeriesWrapper2D
Subroutine : | |||
FileName : | character(*), intent(in )
| ||
VarName : | character(*), intent(in )
| ||
xyz_Press(0:imax-1, 1:jmax, 1:kmax) : | real(DP) , intent(in )
| ||
xyz_Var(0:imax-1, 1:jmax, 1:kmax) : | real(DP) , intent(inout)
| ||
keyword : | character(*), intent(in ) |
subroutine SetValFromTimeSeriesWrapper3D( FileName, VarName, xyz_Press, xyz_Var, keyword ) ! 宣言文 ; Declaration statements ! implicit none character(*), intent(in ):: FileName ! ファイル名. ! File name character(*), intent(in ):: VarName ! 変数名. ! Variable name real(DP) , intent(in ):: xyz_Press(0:imax-1, 1:jmax, 1:kmax) ! 地表面温度. ! Surface temperature real(DP) , intent(inout):: xyz_Var(0:imax-1, 1:jmax, 1:kmax) ! 地表面温度. ! Surface temperature character(*), intent(in ):: keyword ! 初期化確認 ! Initialization check ! if ( .not. read_time_series_inited ) then call MessageNotify( 'E', module_name, 'This module has not been initialized.' ) end if if ( keyword == 'O3' ) then call SetValuesFromTimeSeries( TSDataInfoO3, FileName, VarName, xyz_Press = xyz_Press, xyz_Var = xyz_Var, FlagLogVLev = .true., FlagLog = .true., FlagPositiveDef = .true. ) else stop 'Unsupported keyword' end if end subroutine SetValFromTimeSeriesWrapper3D
Subroutine : | |||
TSDataInfo : | type(time_series_data), intent(inout) | ||
FileName : | character(*) , intent(in )
| ||
VarName : | character(*) , intent(in )
| ||
xyz_Press(0:imax-1, 1:jmax, 1:kmax) : | real(DP) , optional, intent(in )
| ||
xy_Var(0:imax-1, 1:jmax) : | real(DP) , optional, intent(inout)
| ||
xyz_Var(0:imax-1, 1:jmax, 1:kmax) : | real(DP) , optional, intent(inout)
| ||
FlagLogVLev : | logical , optional, intent(in ) | ||
FlagLog : | logical , optional, intent(in ) | ||
FlagPositiveDef : | logical , optional, intent(in ) |
地表面の諸々のデータを設定します. xy_SurfTemp 以外は一回目に呼ばれた時のみ設定されます.
Get various data on ground. Arguments excluding "xy_SurfTemp" are configured at first only.
subroutine SetValuesFromTimeSeries( TSDataInfo, FileName, VarName, xyz_Press, xy_Var, xyz_Var, FlagLogVLev, FlagLog, FlagPositiveDef ) ! ! 地表面の諸々のデータを設定します. ! xy_SurfTemp 以外は一回目に呼ばれた時のみ設定されます. ! ! Get various data on ground. ! Arguments excluding "xy_SurfTemp" are configured at first only. ! ! モジュール引用 ; USE statements ! ! gtool4 データ入力 ! Gtool4 data input ! use gtool_history, only: HistoryGet ! 文字列操作 ! Character handling ! use dc_string, only: toChar ! 日付および時刻の取り扱い ! Date and time handler ! use dc_calendar, only: DCCalDateEvalSecOfYear ! 時刻管理 ! Time control ! use timeset, only: TimeN, InitialDate ! 宣言文 ; Declaration statements ! implicit none type(time_series_data), intent(inout):: TSDataInfo character(*) , intent(in ):: FileName ! ファイル名. ! File name character(*) , intent(in ):: VarName ! 変数名. ! Variable name real(DP) , optional, intent(in ):: xyz_Press(0:imax-1, 1:jmax, 1:kmax) ! 3 次元配列. ! 3D array real(DP) , optional, intent(inout):: xy_Var (0:imax-1, 1:jmax) ! 2 次元配列. ! 2D array real(DP) , optional, intent(inout):: xyz_Var(0:imax-1, 1:jmax, 1:kmax) ! 3 次元配列. ! 3D array logical , optional, intent(in ):: FlagLogVLev logical , optional, intent(in ):: FlagLog logical , optional, intent(in ):: FlagPositiveDef ! 作業変数 ! Work variables ! real(DP):: SecOfYear integer :: tindex !!$ real(DP):: xyza_VarInterpolated(0:imax-1, 1:jmax, 1:kmax, 2) integer :: t ! 初期化確認 ! Initialization check ! if ( .not. read_time_series_inited ) then call MessageNotify( 'E', module_name, 'This module has not been initialized.' ) end if if ( ( .not. present( xy_Var ) ) .and. ( .not. present( xyz_Var ) ) ) then call MessageNotify( 'E', module_name, 'xy_Var or xyz_Var have to be given.' ) end if if ( ( present( xyz_Var ) ) .and. ( .not. present( xyz_Press ) ) ) then call MessageNotify( 'E', module_name, 'xyz_Press has to be given, when xyz_Var is given.' ) end if ! 実行文 ; Executable statement ! if ( .not. associated( TSDataInfo % a_time ) ) then if ( present( xy_Var ) ) then call StructureInit( FileName, VarName, 2, TSDataInfo ) else call StructureInit( FileName, VarName, 3, TSDataInfo ) end if end if ! Update values of time paying attention to leap year ! call UpdateTime( TSDataInfo ) if ( TSDataInfo % tmax >= 2 ) then SecOfYear = DCCalDateEvalSecOfYear( TimeN, date = InitialDate ) !!$ write( 6, * ) TSDataInfo%a_time(TSDataInfo%a_tindex(1)), & !!$ & SecOfYear, TSDataInfo%a_time(TSDataInfo%a_tindex(2)) if ( ( SecOfYear < TSDataInfo%a_time(TSDataInfo%a_tindex(1)) ) .or. ( SecOfYear >= TSDataInfo%a_time(TSDataInfo%a_tindex(2)) ) ) then TSDataInfo % a_tindex(1) = 0 TSDataInfo % a_tindex(2) = 1 do t = 1, TSDataInfo%tmax if ( TSDataInfo % a_time(t) <= SecOfYear ) then TSDataInfo % a_tindex(1) = t TSDataInfo % a_tindex(2) = t+1 end if end do !!$ write( 6, * ) '##############################################' !!$ write( 6, * ) '##############################################' !!$ write( 6, * ) '##############################################' !!$ write( 6, * ) '##############################################' !!$ write( 6, * ) '##############################################' !!$ write( 6, * ) 'in if', TSDataInfo % a_tindex(1), TSDataInfo % a_tindex(2) !!$ write( 6, * ) '##############################################' !!$ write( 6, * ) '##############################################' !!$ write( 6, * ) '##############################################' !!$ write( 6, * ) '##############################################' !!$ write( 6, * ) '##############################################' tindex = TSDataInfo % a_tindex(1) if ( tindex == 0 ) then tindex = TSDataInfo % tmax else if ( tindex == TSDataInfo % tmax + 1 ) then tindex = 1 else tindex = tindex end if if ( present( xy_Var ) ) then call HistoryGet( TSDataInfo%FileName, TSDataInfo%VarName, TSDataInfo%xyza_SavedData(:,:,1,1), range = 'time=^'//toChar(tindex), flag_mpi_split = flag_mpi_split ) else call HistoryGet( TSDataInfo%FileName, TSDataInfo%VarName, TSDataInfo%xyza_SavedData(:,:,:,1), range = 'time=^'//toChar(tindex), flag_mpi_split = flag_mpi_split ) end if tindex = TSDataInfo % a_tindex(2) if ( tindex == 0 ) then tindex = TSDataInfo % tmax else if ( tindex == TSDataInfo % tmax + 1 ) then tindex = 1 else tindex = tindex end if if ( present( xy_Var ) ) then call HistoryGet( TSDataInfo%FileName, TSDataInfo%VarName, TSDataInfo%xyza_SavedData(:,:,1,2), range = 'time=^'//toChar(tindex), flag_mpi_split = flag_mpi_split ) else call HistoryGet( TSDataInfo%FileName, TSDataInfo%VarName, TSDataInfo%xyza_SavedData(:,:,:,2), range = 'time=^'//toChar(tindex), flag_mpi_split = flag_mpi_split ) end if end if if ( present( xy_Var ) ) then xy_Var = ( TSDataInfo%xyza_SavedData(:,:,1,2) - TSDataInfo%xyza_SavedData(:,:,1,1) ) / ( TSDataInfo%a_time(TSDataInfo%a_tindex(2)) - TSDataInfo%a_time(TSDataInfo%a_tindex(1)) ) * ( SecOfYear - TSDataInfo%a_time(TSDataInfo%a_tindex(1)) ) + TSDataInfo%xyza_SavedData(:,:,1,1) else TSDataInfo%xyz_VarTimeInterpolated = ( TSDataInfo%xyza_SavedData(:,:,:,2) - TSDataInfo%xyza_SavedData(:,:,:,1) ) / ( TSDataInfo%a_time(TSDataInfo%a_tindex(2)) - TSDataInfo%a_time(TSDataInfo%a_tindex(1)) ) * ( SecOfYear - TSDataInfo%a_time(TSDataInfo%a_tindex(1)) ) + TSDataInfo%xyza_SavedData(:,:,:,1) call VerticalInterpolation( TSDataInfo%NVLevels, TSDataInfo%z_VLevels, TSDataInfo%xyz_VarTimeInterpolated, xyz_Press, xyz_Var, FlagLogVLev, FlagLog, FlagPositiveDef ) end if else if ( present( xy_Var ) ) then xy_Var = TSDataInfo%xyza_SavedData(:,:,1,1) else call VerticalInterpolation( TSDataInfo%NVLevels, TSDataInfo%z_VLevels, TSDataInfo%xyza_SavedData(:,:,:,1), xyz_Press, xyz_Var, FlagLogVLev, FlagLog, FlagPositiveDef ) end if end if end subroutine SetValuesFromTimeSeries
Subroutine : | |||
FileName : | character(*) , intent(in )
| ||
VarName : | character(*) , intent(in )
| ||
NDims : | integer , intent(in )
| ||
TSDataInfo : | type(time_series_data), intent(inout) |
subroutine StructureInit( FileName, VarName, NDims, TSDataInfo ) ! モジュール引用 ; USE statements ! ! gtool4 データ入力 ! Gtool4 data input ! use gtool_history, only: HistoryGet, HistoryGetPointer ! 文字列操作 ! Character handling ! use dc_string, only: toChar ! 日付および時刻の取り扱い ! Date and time handler ! use dc_calendar, only: DCCalInquire, DCCalDateEvalSecOfYear ! 時刻管理 ! Time control ! use timeset, only: TimeN, InitialDate ! NetCDF wrapper ! use netcdf_wrapper, only: NWGetAtt, NWInqDimLen ! 宣言文 ; Declaration statements ! implicit none character(*) , intent(in ):: FileName ! ファイル名. ! File name character(*) , intent(in ):: VarName ! 変数名. ! Variable name integer , intent(in ):: NDims ! 時間以外の変数の次元. ! Number of dimensions except for time dimension type(time_series_data), intent(inout):: TSDataInfo ! 作業変数 ! Work variables ! real(DP) :: SecOfYear !!$ real(DP), pointer:: a_time(:) integer :: tindex integer :: t character(STRING):: VLevelName character(STRING):: attchar integer:: hour_in_day, min_in_hour, day_in_year integer, pointer:: day_in_month_ptr(:) => null() real(DP):: sec_in_min, sec_in_day ! 実行文 ; Executable statement ! ! 初期化確認 ! Initialization check ! if ( .not. read_time_series_inited ) then call MessageNotify( 'E', module_name, 'This module has not been initialized.' ) end if !!$ nullify( a_time ) call DCCalInquire( day_in_month_ptr = day_in_month_ptr , hour_in_day = hour_in_day , min_in_hour = min_in_hour , sec_in_min = sec_in_min ) ! (out) day_in_year = sum( day_in_month_ptr ) deallocate( day_in_month_ptr ) sec_in_day = hour_in_day * min_in_hour * sec_in_min TSDataInfo % FileName = FileName TSDataInfo % VarName = VarName TSDataInfo % NDims = NDims if ( TSDataInfo % NDims == 2 ) then TSDataInfo % NVLevels = 1 else VLevelName = "plev" call HistoryGetPointer( TSDataInfo%FileName, VLevelName, TSDataInfo % z_VLevels, flag_mpi_split = flag_mpi_split ) TSDataInfo % NVLevels = size( TSDataInfo % z_Vlevels ) end if !----- !!$ call HistoryGetPointer( & !!$ & TSDataInfo%FileName, & !!$ & 'time', & !!$ & a_time, & !!$ & flag_mpi_split = flag_mpi_split & !!$ & ) !!$ TSDataInfo % tmax = size( a_time ) !!$ !!$ call NWGetAtt( TSDataInfo%FileName, 'time', 'units', attchar ) !!$ if ( ( attchar(1:3) == 'day' ) .or. & !!$ & ( attchar(1:3) == 'DAY' ) .or. & !!$ & ( attchar(1:3) == 'Day' ) ) then !!$ a_time = a_time * DAY_SECONDS !!$ else !!$ call MessageNotify( 'E', module_name, 'Unit of time, %c, is not supported.', c1 = trim(attchar) ) !!$ end if !----- call NWInqDimLen( TSDataInfo%FileName, 'time', TSDataInfo%tmax ) if ( TSDataInfo % tmax >= 2 ) then SecOfYear = DCCalDateEvalSecOfYear( TimeN, date = InitialDate ) !----- !!$ allocate( TSDataInfo % a_time(0:TSDataInfo%tmax+1) ) !!$ TSDataInfo % a_time(0) = - ( YEAR_DAYS * DAY_SECONDS - a_time(TSDataInfo%tmax) ) !!$ do t = 1, TSDataInfo%tmax !!$ TSDataInfo % a_time(t) = a_time(t) !!$ end do !!$ TSDataInfo % a_time(TSDataInfo%tmax+1) = YEAR_DAYS * DAY_SECONDS + a_time(1) !!$ deallocate( a_time ) !----- allocate( TSDataInfo % a_time(0:TSDataInfo%tmax+1) ) call HistoryGet( TSDataInfo%FileName, 'time', TSDataInfo%a_time(1:TSDataInfo%tmax), flag_mpi_split = flag_mpi_split ) call NWGetAtt( TSDataInfo%FileName, 'time', 'units', attchar ) if ( ( attchar(1:3) == 'day' ) .or. ( attchar(1:3) == 'DAY' ) .or. ( attchar(1:3) == 'Day' ) ) then TSDataInfo%a_time(1:TSDataInfo%tmax) = TSDataInfo%a_time(1:TSDataInfo%tmax) * sec_in_day else call MessageNotify( 'E', module_name, 'Unit of time, %c, is not supported.', c1 = trim(attchar) ) end if TSDataInfo % a_time(0 ) = - ( day_in_year * sec_in_day - TSDataInfo % a_time(TSDataInfo%tmax) ) TSDataInfo % a_time(TSDataInfo%tmax+1) = day_in_year * sec_in_day + TSDataInfo % a_time(1 ) allocate( TSDataInfo % a_tindex(1:2) ) allocate( TSDataInfo % xyza_SavedData(0:imax-1,1:jmax,1:TSDataInfo%NVLevels,1:2) ) TSDataInfo % a_tindex(1) = 0 TSDataInfo % a_tindex(2) = 1 do t = 1, TSDataInfo%tmax if ( TSDataInfo % a_time(t) <= SecOfYear ) then TSDataInfo % a_tindex(1) = t TSDataInfo % a_tindex(2) = t+1 end if end do tindex = TSDataInfo % a_tindex(1) if ( tindex == 0 ) then tindex = TSDataInfo % tmax else if ( tindex == TSDataInfo % tmax + 1 ) then tindex = 1 else tindex = tindex end if if ( TSDataInfo % NDims == 2 ) then call HistoryGet( TSDataInfo%FileName, TSDataInfo%VarName, TSDataInfo%xyza_SavedData(:,:,1,1), range = 'time=^'//toChar(tindex), flag_mpi_split = flag_mpi_split ) else call HistoryGet( TSDataInfo%FileName, TSDataInfo%VarName, TSDataInfo%xyza_SavedData(:,:,:,1), range = 'time=^'//toChar(tindex), flag_mpi_split = flag_mpi_split ) end if tindex = TSDataInfo % a_tindex(2) if ( tindex == 0 ) then tindex = TSDataInfo % tmax else if ( tindex == TSDataInfo % tmax + 1 ) then tindex = 1 else tindex = tindex end if if ( TSDataInfo % NDims == 2 ) then call HistoryGet( TSDataInfo%FileName, TSDataInfo%VarName, TSDataInfo%xyza_SavedData(:,:,1,2), range = 'time=^'//toChar(tindex), flag_mpi_split = flag_mpi_split ) else call HistoryGet( TSDataInfo%FileName, TSDataInfo%VarName, TSDataInfo%xyza_SavedData(:,:,:,2), range = 'time=^'//toChar(tindex), flag_mpi_split = flag_mpi_split ) end if allocate( TSDataInfo%xyz_VarTimeInterpolated(0:imax-1, 1:jmax, 1:TSDataInfo%NVLevels) ) else !----- !!$ allocate( TSDataInfo % a_time(TSDataInfo%tmax) ) !!$ do t = 1, TSDataInfo%tmax !!$ TSDataInfo % a_time(t) = a_time(t) !!$ end do !----- allocate( TSDataInfo % a_time(TSDataInfo%tmax) ) call HistoryGet( TSDataInfo%FileName, 'time', TSDataInfo%a_time(1:TSDataInfo%tmax), flag_mpi_split = flag_mpi_split ) call NWGetAtt( TSDataInfo%FileName, 'time', 'units', attchar ) if ( ( attchar(1:3) == 'day' ) .or. ( attchar(1:3) == 'DAY' ) .or. ( attchar(1:3) == 'Day' ) ) then TSDataInfo%a_time(1:TSDataInfo%tmax) = TSDataInfo%a_time(1:TSDataInfo%tmax) * sec_in_day else call MessageNotify( 'E', module_name, 'Unit of time, %c, is not supported.', c1 = trim(attchar) ) end if allocate( TSDataInfo % a_tindex(1:1) ) allocate( TSDataInfo % xyza_SavedData(0:imax-1,1:jmax,1:TSDataInfo%NVLevels,1:1) ) TSDataInfo % a_tindex(1) = 1 if ( TSDataInfo % NDims == 2 ) then call HistoryGet( TSDataInfo%FileName, TSDataInfo%VarName, TSDataInfo%xyza_SavedData(:,:,1,1), range = 'time=^'//toChar(TSDataInfo%a_tindex(1)), flag_mpi_split = flag_mpi_split ) else call HistoryGet( TSDataInfo%FileName, TSDataInfo%VarName, TSDataInfo%xyza_SavedData(:,:,:,1), range = 'time=^'//toChar(TSDataInfo%a_tindex(1)), flag_mpi_split = flag_mpi_split ) end if end if end subroutine StructureInit
Subroutine : | |
TSDataInfo : | type(time_series_data), intent(inout) |
subroutine UpdateTime( TSDataInfo ) ! モジュール引用 ; USE statements ! ! 暦と日時の取り扱い ! Calendar and Date handler ! use dc_calendar, only: DC_CAL_DATE ! 日時を表現するデータ型. ! Data type for date and time ! 日付および時刻の取り扱い ! Date and time handler ! use dc_calendar, only: DCCalInquire, DCCalDateCreate, DCCalDateInquire, DCCalDateDifference, DCCalDateChkLeapYear ! 時刻管理 ! Time control ! use timeset, only: TimeN, InitialDate ! 宣言文 ; Declaration statements ! implicit none type(time_series_data), intent(inout):: TSDataInfo ! 作業変数 ! Work variables ! integer :: hour_in_day, min_in_hour, day_in_year integer, pointer :: day_in_month_ptr(:) => null() real(DP) :: sec_in_min, sec_in_day integer :: year type(DC_CAL_DATE) :: PreYearDate ! 実行文 ; Executable statement ! ! 初期化確認 ! Initialization check ! if ( .not. read_time_series_inited ) then call MessageNotify( 'E', module_name, 'This module has not been initialized.' ) end if if ( TSDataInfo % tmax >= 2 ) then call DCCalInquire( day_in_month_ptr = day_in_month_ptr , hour_in_day = hour_in_day , min_in_hour = min_in_hour , sec_in_min = sec_in_min ) ! (out) day_in_year = sum( day_in_month_ptr ) deallocate( day_in_month_ptr ) sec_in_day = hour_in_day * min_in_hour * sec_in_min ! Set TSDataInfo % a_time(0) ! call DCCalDateInquire( year = year, elapse_sec = TimeN, date = InitialDate ) call DCCalDateCreate( year = year - 1, month = 1, day = 1, hour = 0, min = 0, sec = 0.0_DP, date = PreYearDate ) if ( DCCalDateChkLeapYear( 0.0_DP, PreYearDate ) ) then TSDataInfo % a_time(0 ) = - ( ( day_in_year + 1 ) * sec_in_day - TSDataInfo % a_time(TSDataInfo%tmax) ) else TSDataInfo % a_time(0 ) = - ( day_in_year * sec_in_day - TSDataInfo % a_time(TSDataInfo%tmax) ) end if ! Set TSDataInfo % a_time(TSDataInfo%tmax+1) ! if ( DCCalDateChkLeapYear( TimeN, InitialDate ) ) then TSDataInfo % a_time(TSDataInfo%tmax+1) = ( day_in_year + 1 ) * sec_in_day + TSDataInfo % a_time(1 ) else TSDataInfo % a_time(TSDataInfo%tmax+1) = day_in_year * sec_in_day + TSDataInfo % a_time(1 ) end if end if end subroutine UpdateTime
Subroutine : | |
NVLevels : | integer , intent(in ) |
z_VLevelsIn(1:NVLevels) : | real(DP), intent(in ) |
xyz_VarIn(0:imax-1, 1:jmax, 1:NVLevels) : | real(DP), intent(in ) |
xyz_VarOut(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(out) |
FlagLogVLev : | logical , intent(in ), optional |
FlagLog : | logical , intent(in ), optional |
FlagPositiveDef : | logical , intent(in ), optional |
subroutine VerticalInterpolation( NVLevels, z_VLevelsIn, xyz_VarIn, xyz_VLevelsOut, xyz_VarOut, FlagLogVLev, FlagLog, FlagPositiveDef ) ! モジュール引用 ; USE statements ! ! 宣言文 ; Declaration statements ! implicit none integer , intent(in ):: NVLevels real(DP), intent(in ):: z_VLevelsIn(1:NVLevels) real(DP), intent(in ):: xyz_VarIn (0:imax-1, 1:jmax, 1:NVLevels) real(DP), intent(in ):: xyz_VlevelsOut(0:imax-1, 1:jmax, 1:kmax) real(DP), intent(out):: xyz_VarOut (0:imax-1, 1:jmax, 1:kmax) logical , intent(in ), optional:: FlagLogVLev logical , intent(in ), optional:: FlagLog logical , intent(in ), optional:: FlagPositiveDef ! 作業変数 ! Work variables ! logical:: FlagLogVLevLV logical:: FlagLogLV logical:: FlagPositiveDefLV integer:: i integer:: j integer:: k integer:: kk integer:: xyz_kk(0:imax-1, 1:jmax, 1:kmax) ! 実行文 ; Executable statement ! ! 初期化確認 ! Initialization check ! if ( .not. read_time_series_inited ) then call MessageNotify( 'E', module_name, 'This module has not been initialized.' ) end if ! Check order of vertical levels ! if ( z_VlevelsIn(1) < z_VlevelsIn(2) ) then call MessageNotify( 'E', module_name, 'The order of vertical levels is inappropriate.' ) end if ! Flags are checked ! if ( present( FlagLogVLev ) ) then FlagLogVLevLV = FlagLogVLev else FlagLogVLevLV = .false. end if if ( present( FlagLog ) ) then FlagLogLV = FlagLog else FlagLogLV = .false. end if if ( present( FlagPositiveDef ) ) then FlagPositiveDefLV = FlagPositiveDef else FlagPositiveDefLV = .false. end if xyz_kk = 1 do kk = 1, NVLevels-1 do k = 1, kmax do j = 1, jmax do i = 0, imax-1 if ( xyz_VLevelsOut(i,j,k) < z_VLevelsIn(kk) ) then xyz_kk(i,j,k) = kk end if end do end do end do end do if ( FlagLogVLevLV ) then if ( FlagLogLV ) then ! log-log do k = 1, kmax do j = 1, jmax do i = 0, imax-1 !!$ if ( xyz_VLevelsOut(i,j,k) > z_VLevelsIn(1) ) then !!$ xyz_VarOut(i,j,k) = xyz_VarIn(i,j,1) !!$ !!$! else if ( xyz_VLevelsOut(i,j,k) < z_VLevelsIn(NVLevels) ) then !!$! call MessageNotify( 'E', module_name, 'Vertical level is out of given range.' ) !!$!! xyz_VarOut(i,j,k) = 0.0d0 !!$ !!$ else !!$ !!$ xyz_VarOut(i,j,k) = & !!$ & ( xyz_VarIn(i,j,xyz_kk(i,j,k)+1) - xyz_VarIn(i,j,xyz_kk(i,j,k)) ) & !!$ & / log( z_VLevelsIn (xyz_kk(i,j,k)+1) / z_VLevelsIn(xyz_kk(i,j,k)) ) & !!$ & * log( xyz_VLevelsOut(i,j,k) / z_VLevelsIn(xyz_kk(i,j,k)) ) & !!$ & + xyz_VarIn(i,j,xyz_kk(i,j,k)) !!$ !!$ end if xyz_VarOut(i,j,k) = log( ( xyz_VarIn(i,j,xyz_kk(i,j,k)+1) + 1.0d-100 ) / ( xyz_VarIn(i,j,xyz_kk(i,j,k) ) + 1.0d-100 ) ) / log( z_VLevelsIn (xyz_kk(i,j,k)+1) / z_VLevelsIn(xyz_kk(i,j,k)) ) * log( xyz_VLevelsOut(i,j,k) / z_VLevelsIn(xyz_kk(i,j,k)) ) + log( xyz_VarIn(i,j,xyz_kk(i,j,k)) + 1.0d-100 ) xyz_VarOut(i,j,k) = exp( xyz_VarOut(i,j,k) ) end do end do end do else ! log-linear do k = 1, kmax do j = 1, jmax do i = 0, imax-1 xyz_VarOut(i,j,k) = ( xyz_VarIn(i,j,xyz_kk(i,j,k)+1) - xyz_VarIn(i,j,xyz_kk(i,j,k) ) ) / log( z_VLevelsIn (xyz_kk(i,j,k)+1) / z_VLevelsIn(xyz_kk(i,j,k)) ) * log( xyz_VLevelsOut(i,j,k) / z_VLevelsIn(xyz_kk(i,j,k)) ) + xyz_VarIn(i,j,xyz_kk(i,j,k)) end do end do end do end if else if ( FlagLogLV ) then ! linear-log do k = 1, kmax do j = 1, jmax do i = 0, imax-1 xyz_VarOut(i,j,k) = log( ( xyz_VarIn(i,j,xyz_kk(i,j,k)+1) + 1.0d-100 ) / ( xyz_VarIn(i,j,xyz_kk(i,j,k) ) + 1.0d-100 ) ) / ( z_VLevelsIn(xyz_kk(i,j,k)+1) - z_VLevelsIn(xyz_kk(i,j,k) ) ) * ( xyz_VLevelsOut(i,j,k) - z_VLevelsIn(xyz_kk(i,j,k)) ) + log( xyz_VarIn(i,j,xyz_kk(i,j,k)) + 1.0d-100 ) xyz_VarOut(i,j,k) = exp( xyz_VarOut(i,j,k) ) end do end do end do else ! linear-linear do k = 1, kmax do j = 1, jmax do i = 0, imax-1 xyz_VarOut(i,j,k) = ( xyz_VarIn(i,j,xyz_kk(i,j,k)+1) - xyz_VarIn(i,j,xyz_kk(i,j,k) ) ) / ( z_VLevelsIn(xyz_kk(i,j,k)+1) - z_VLevelsIn(xyz_kk(i,j,k) ) ) * ( xyz_VLevelsOut(i,j,k) - z_VLevelsIn(xyz_kk(i,j,k)) ) + xyz_VarIn(i,j,xyz_kk(i,j,k)) end do end do end do end if end if if ( FlagPositiveDefLV ) then do k = 1, kmax do j = 1, jmax do i = 0, imax-1 if ( xyz_VarOut(i,j,k) < 0.0d0 ) then xyz_VarOut(i,j,k) = 0.0d0 end if end do end do end do end if !!$ i = 0 !!$ j = jmax / 2 !!$ do k = 1, NVLevels !!$ write( 92, * ) z_VLevelsIn(k), xyz_VarIn(i,j,k) !!$ end do !!$ write( 92, * ) !!$ call flush( 92 ) !!$ do k = 1, kmax !!$ write( 93, * ) xyz_VLevelsOut(i,j,k), xyz_VarOut(i,j,k) !!$ end do !!$ write( 93, * ) !!$ call flush( 93 ) !!$ stop end subroutine VerticalInterpolation
Constant : | |||
module_name = ‘read_time_series_from_file‘ : | character(*), parameter
|
Variable : | |||
read_time_series_inited = .false. : | logical, save
|
Derived Type : | |||
FileName : | character(STRING)
| ||
VarName : | character(STRING)
| ||
NDims : | integer | ||
tmax : | integer | ||
NVLevels : | integer | ||
z_VLevels(:) : | real(DP), pointer | ||
a_time(:) : | real(DP), pointer | ||
a_tindex(:) : | integer , pointer | ||
xyza_SavedData(:,:,:,:) : | real(DP), pointer | ||
xyz_VarTimeInterpolated(:,:,:) : | real(DP), pointer |
Constant : | |||
version = ’$Name: dcpam5-20120922 $’ // ’$Id: read_time_series.f90,v 1.14 2012-07-27 00:52:39 yot Exp $’ : | character(*), parameter
|