!== dc_date.f90 - 日付および時刻に関する手続きを提供するモジュール
!
! Authors:: Yasuhiro MORIKAWA, Eizi TOYODA
! Version:: $Id: dc_date.f90,v 1.21 2007/06/15 04:59:09 morikawa Exp $
! Tag Name:: $Name: gt4f90io-20070824 $
! Copyright:: Copyright (C) GFD Dennou Club, 2000-2005. All rights reserved.
! License:: See COPYRIGHT[link:../../COPYRIGHT]
module dc_date
!== Overview
!
! 日付および時刻を扱うための手続きを提供します
!
!== List
!
! 以下の手続きは構造型 dc_date_types#DC_DATETIME または
! dc_date_types#DC_DIFFTIME 変数 (日時, 時刻に関する情報を格納)
! を対象とします.
!
! Create :: 初期化
!
! #assignment(=) :: 〃
!
! Eval :: 日時, 時刻情報を個別に取得
!
! toChar :: 日時, 時刻情報を文字型変数へ変換
!
! PutLine :: 日時, 時刻情報の印字を行います.
!
! EvalDay :: 日数 (実数型) に換算して取得
! EvalHour :: 時間 (実数型) に換算して取得
! EvalMin :: 分 (実数型) に換算して取得
! EvalSec :: 秒 (実数型) に換算して取得
! EvalByUnit :: 単位を指定し, 日, 時, 分, 秒のいづれか (実数型)
! に換算して取得
!
! #operator(+) :: 加算 (dc_date_types#DC_DATETIME 型 および
! dc_date_types#DC_DIFFTIME 型 同士)
! #operator(-) :: 減算 (dc_date_types#DC_DATETIME 型 および
! dc_date_types#DC_DIFFTIME 型 同士)
! #operator(*) :: 乗算 (dc_date_types#DC_DIFFTIME 型と数値型)
! #operator(/) :: 除算 (dc_date_types#DC_DIFFTIME 型と数値型)
! mod :: 余り (dc_date_types#DC_DIFFTIME 型同士)
! #operator(==) :: 比較 (dc_date_types#DC_DATETIME 型同士)
! #operator(>) :: 比較 (dc_date_types#DC_DATETIME 型同士)
! #operator(<) :: 比較 (dc_date_types#DC_DATETIME 型同士)
!
! SetZone :: タイムゾーンを変更
!
!
! 以下の手続きは dc_date_types 内部の変数を変更します.
!
! SetCaltype :: 暦法のデフォルトを変更
! SetSecOfDay :: 1 日の秒数のデフォルトを変更
!
! その他の手続き
!
! ValidCaltype :: 暦法が有効なものかをチェック
! ValidZone :: タイムゾーンとして有効化をチェック
! ZoneToDiff :: タイムゾーンを dc_date_types#DC_DIFFTIME 変数へと変換
!
!== Usage
!
!=== 現在時刻の表示
!
! DC_DATETIME 型の変数にサブルーチン Create を用いると, 時刻が設定されます.
! 下記のように特に年月日を指定しないと現在時刻が設定されます.
! 設定された時刻は toChar によって文字型変数へと変換できます.
! サブルーチン Printf に関しては dc_string#Printf を参照ください.
!
! program dc_date_sapmle1
! use dc_string, only: Printf
! use dc_date_types, only: DC_DATETIME
! use dc_date, only: Create, toChar
! implicit none
! type(DC_DATETIME) :: time
!
! call Create(time)
! call Printf(fmt='current date and time is %c', c1=trim(toChar(time)))
! end program dc_date_sapmle1
!
!=== 日時, 時刻情報の加算
!
! DC_DIFFTIME 型の変数は日時差を表現します. 下記の例では,
! 日時差を表現するための変数として *diff*
! を用意し, サブルーチン Create によって 25 日 + 12 時間 + 50 分の日時差
! を設定しています. DC_DATETIME 型の変数 *time_before* と *diff* とを
! #operator(+) によって加算することで *time_before* から
! 25 日 + 12 時間 + 50 分を進めた日時 *time_after* を取得しています.
!
! program dc_date_sapmle2
! use dc_string, only: Printf
! use dc_date_types, only: DC_DATETIME, DC_DIFFTIME
! use dc_date, only: Create, toChar, operator(+)
! implicit none
! type(DC_DATETIME) :: time_before, time_after
! type(DC_DIFFTIME) :: diff
!
! call Create(time_before, year=2006, mon=6, day=10, hour=14, min=15, sec=0.0d0)
! call Create(diff, day=25, hour=12, min=50)
! time_after = time_before + diff
!
! call Printf(fmt='%c + %c = %c', &
! & c1=trim(toChar(time_before)), c2=trim(toChar(diff)), &
! & c3=trim(toChar(time_after)))
! end program dc_date_sapmle2
!
!
!=== 時間積分のループへの応用
!
! 以下は dA/dt = - αA (初期値 1, α=0.0001) を t = 12 (時間)
! まで解くプログラムの例です. 時間積分には前進差分を用いています.
! Δt, データの出力間隔, 計算時間に DC_DIFFTIME を用いることで,
! ループの終了処理や
! データ出力の際の時刻の比較が容易となります.
!
! program dc_date_sapmle3
! use dc_types, only: DP
! use dc_date, only: Create, EvalSec, EvalByUnit, mod, &
! & operator(*), operator(==), operator(>)
! use dc_date_types, only: DC_DIFFTIME
! implicit none
! real(DP) :: func_a = 1.0d0 ! 関数 A の初期値
! real(DP), parameter :: alph = 0.0001d0 ! 係数 α
! character(*), parameter :: out_unit = 'hour' ! 出力される時刻の単位
! type(DC_DIFFTIME):: DelTimef, intervalf, calctimef
! integer :: i
! continue
! call Create(DelTimef, 5.0d0, 'sec') ! Δt = 5.0 (秒)
! call Create(intervalf, 1.0d0, 'min') ! データの出力間隔 = 1.0 (分)
! call Create(calctimef, 12.0d0, 'hour') ! 計算時間 = 12.0 (時間)
!
! open(10, file='dc_date_sample.dat')
! write(10,'(A,A,A)') '# ', out_unit, ' value'
!
! i = 1
! do
! if (DelTimef * i > calctimef) exit ! 計算時間を過ぎたら終了
!
! !---------------------------------------------
! ! A_(n+1) = (1 - αΔt) * A_(n)
! !---------------------------------------------
! func_a = (1.0 - alph * EvalSec(DelTimef)) * func_a
!
! !---------------------------------------------
! ! intervalf (1 分) 毎にデータを出力
! !---------------------------------------------
! if (mod(DelTimef * i, intervalf) == 0) then
! write(10,*) ' ', EvalByUnit(DelTimef * i, out_unit), func_a
! end if
! i = i + 1
! end do
! end program dc_date_sapmle3
!
!
! 日付時刻と時間間隔を区別する。
! 型宣言によって自明に定まるサブルーチンは dc_date_types に置く。
use dc_date_types, only: DC_DATETIME, DC_DIFFTIME
use dc_types, only: DP, STRING
use dc_present, only: present_and_not_empty
implicit none
private
public :: Create, assignment(=), Eval
public :: SetCaltype, SetZone, SetSecOfDay
public :: ValidCaltype, ValidZone, ZoneToDiff
public :: mod, operator(/), operator(-), operator(+), operator(*)
public :: operator(<), operator(>), operator(==)
public :: toChar, EvalDay, EvalHour, EvalMin, EvalSec, EvalByUnit, PutLine
public :: dcdate_normalize, dcdate_parse_unit
interface Create
subroutine DCDateTimeCreate1(time, &
& year, mon, day, hour, min, sec, &
& zone, caltype, day_seconds, err) !:doc-priority 40:
use dc_types, only: DP
use dc_date_types, only: DC_DATETIME
type(DC_DATETIME), intent(out):: time
integer, intent(in), optional:: year, mon, day, hour, min
real(DP),intent(in), optional:: sec, day_seconds
character(*), intent(in), optional :: zone
integer, intent(in), optional:: caltype
logical, intent(out), optional:: err
end subroutine DCDateTimeCreate1
subroutine DCDiffTimeCreate1(diff, &
& year, mon, day, hour, min, sec, day_seconds) !:doc-priority 60:
use dc_types, only: DP
use dc_date_types, only: DC_DIFFTIME
type(DC_DIFFTIME), intent(out) :: diff
integer, intent(in), optional:: year, mon, day, hour, min
real(DP),intent(in), optional:: sec, day_seconds
end subroutine DCDiffTimeCreate1
subroutine DCDiffTimeCreate2(diff, value, unit, err) !:doc-priority 70:
use dc_types, only: DP
use dc_date_types, only: DC_DIFFTIME
type(DC_DIFFTIME), intent(out) :: diff
real(DP), intent(in) :: value
character(*), intent(in) :: unit
logical, intent(out), optional :: err
end subroutine DCDiffTimeCreate2
end interface
interface assignment(=)
subroutine DCDateTimeCreateR(time, sec) !:doc-priority 30:
use dc_date_types, only: DC_DATETIME
type(DC_DATETIME), intent(out):: time
real, intent(in):: sec
end subroutine DCDateTimeCreateR
subroutine DCDateTimeCreateD(time, sec) !:doc-priority 40:
use dc_types, only: DP
use dc_date_types, only: DC_DATETIME
type(DC_DATETIME), intent(out):: time
real(DP), intent(in):: sec
end subroutine DCDateTimeCreateD
subroutine DCDiffTimeCreateR(diff, sec) !:doc-priority 60:
use dc_date_types, only: DC_DIFFTIME
type(DC_DIFFTIME), intent(out):: diff
real, intent(in):: sec
end subroutine DCDiffTimeCreateR
subroutine DCDiffTimeCreateD(diff, sec) !:doc-priority 70:
use dc_types, only: DP
use dc_date_types, only: DC_DIFFTIME
type(DC_DIFFTIME), intent(out):: diff
real(DP), intent(in):: sec
end subroutine DCDiffTimeCreateD
!!$ subroutine DCDateLetFC(diff, string)
!!$ use dc_date_types, only: DC_DIFFTIME
!!$ type(DC_DIFFTIME), intent(out):: diff
!!$ character(len = *), intent(in):: string
!!$ end subroutine DCDateLetFC
!!$
!!$ subroutine DCDateLetTC(time, string)
!!$ use dc_date_types, only: DC_DATETIME
!!$ type(DC_DATETIME), intent(out):: time
!!$ character(len = *), intent(in):: string
!!$ end subroutine DCDateLetTC
end interface
interface SetCaltype
subroutine DCDateTimeSetCaltype(caltype)
integer, intent(in):: caltype
end subroutine DCDateTimeSetCaltype
end interface
interface SetSecOfDay
subroutine DCDateTimeSetSecOfDay(sec)
use dc_types, only: DP
real(DP), intent(in):: sec
end subroutine DCDateTimeSetSecOfDay
end interface
interface ValidCaltype
function DCDateTimeValidCaltype(caltype) result(result)
integer, intent(in):: caltype
logical:: result
end function DCDateTimeValidCaltype
end interface
interface ValidZone
function DCDateTimeValidZone(zone) result(result)
character(*), intent(in):: zone
logical:: result
end function DCDateTimeValidZone
end interface
interface ZoneToDiff
function DCDateTimeZoneToDiff(zone) result(diff)
use dc_date_types, only: DC_DIFFTIME
type(DC_DIFFTIME):: diff
character(*), intent(in):: zone
end function DCDateTimeZoneToDiff
end interface
interface SetZone
subroutine DCDateTimeSetZone(time, zone, err)
use dc_date_types, only: DC_DATETIME
type(DC_DATETIME), intent(inout):: time
character(*), intent(in):: zone
logical, intent(out), optional:: err
end subroutine DCDateTimeSetZone
end interface
interface Eval
subroutine DCDateTimeEval1(time, &
& year, mon, day, hour, min, sec, caltype, zone) !:doc-priority 40:
use dc_types, only: DP
use dc_date_types, only: DC_DATETIME
type(DC_DATETIME), intent(in):: time
integer, intent(out), optional:: year, mon, day, hour, min, caltype
real(DP), intent(out), optional:: sec
character(*), intent(out), optional:: zone
end subroutine DCDateTimeEval1
!!$ subroutine DCDateTimeEval0(time, mon, day, sec)
!!$ use dc_date_types, only: DC_DATETIME
!!$ use dc_types, only: DP
!!$ type(DC_DATETIME), intent(in):: time
!!$ integer, intent(out):: mon, day
!!$ real(DP), intent(out):: sec
!!$ end subroutine DCDateTimeEval0
subroutine DCDiffTimeEval1(diff, year, mon, day, hour, min, sec) !:doc-priority 60:
use dc_types, only: DP
use dc_date_types, only: DC_DIFFTIME
type(DC_DIFFTIME), intent(in):: diff
integer, intent(out), optional:: year, mon, day, hour, min
real(DP), intent(out), optional:: sec
end subroutine DCDiffTimeEval1
end interface
interface EvalDay
function DCDateTimeEvalDay(time) result(result) !:doc-priority 40:
use dc_types, only: DP
use dc_date_types, only: DC_DATETIME
real(DP):: result
type(DC_DATETIME), intent(in):: time
end function DCDateTimeEvalDay
function DCDiffTimeEvalDay(diff) result(result) !:doc-priority 60:
use dc_types, only: DP
use dc_date_types, only: DC_DIFFTIME
real(DP):: result
type(DC_DIFFTIME), intent(in):: diff
end function DCDiffTimeEvalDay
end interface
interface EvalHour
function DCDateTimeEvalHour(time) result(result) !:doc-priority 40:
use dc_types, only: DP
use dc_date_types, only: DC_DATETIME
real(DP):: result
type(DC_DATETIME), intent(in):: time
end function DCDateTimeEvalHour
function DCDifftimeEvalHour(diff) result(result) !:doc-priority 60:
use dc_types, only: DP
use dc_date_types, only: DC_DIFFTIME
real(DP):: result
type(DC_DIFFTIME), intent(in):: diff
end function DCDifftimeEvalHour
end interface
interface EvalMin
function DCDateTimeEvalMin(time) result(result) !:doc-priority 40:
use dc_types, only: DP
use dc_date_types, only: DC_DATETIME
real(DP):: result
type(DC_DATETIME), intent(in):: time
end function DCDateTimeEvalMin
function DCDifftimeEvalMin(diff) result(result) !:doc-priority 60:
use dc_types, only: DP
use dc_date_types, only: DC_DIFFTIME
real(DP):: result
type(DC_DIFFTIME), intent(in):: diff
end function DCDifftimeEvalMin
end interface
interface EvalSec
function DCDateTimeEvalSec(time) result(result) !:doc-priority 40:
use dc_types, only: DP
use dc_date_types, only: DC_DATETIME
real(DP):: result
type(DC_DATETIME), intent(in):: time
end function DCDateTimeEvalSec
function DCDifftimeEvalSec(diff) result(result) !:doc-priority 60:
use dc_types, only: DP
use dc_date_types, only: DC_DIFFTIME
real(DP):: result
type(DC_DIFFTIME), intent(in):: diff
end function DCDifftimeEvalSec
end interface
interface EvalByUnit
function DCDateTimeEvalByUnit(time, unit) result(result)
use dc_types, only: DP, TOKEN
use dc_date_types, only: DC_DATETIME
real(DP):: result
type(DC_DATETIME), intent(in):: time
character(*), intent(in):: unit
end function DCDateTimeEvalByUnit
function DCDiffTimeEvalByUnit(diff, unit) result(result)
use dc_types, only: DP, TOKEN
use dc_date_types, only: DC_DIFFTIME
real(DP):: result
type(DC_DIFFTIME), intent(in):: diff
character(*), intent(in):: unit
end function DCDiffTimeEvalByUnit
end interface
interface toChar
function DCDateTimeToChar(time) result(result) !:doc-priority 40:
use dc_types, only: STRING
use dc_date_types, only: DC_DATETIME
character(STRING) :: result
type(DC_DATETIME), intent(in):: time
end function DCDateTimeToChar
function DCDiffTimeToChar(diff) result(result) !:doc-priority 60:
use dc_types, only: STRING
use dc_date_types, only: DC_DIFFTIME
character(STRING) :: result
type(DC_DIFFTIME), intent(in):: diff
end function DCDiffTimeToChar
end interface
interface PutLine
subroutine DCDateTimePutLine(time, unit)
use dc_date_types, only: DC_DATETIME
type(DC_DATETIME), intent(in) :: time
integer, intent(in), optional :: unit
end subroutine DCDateTimePutLine
subroutine DCDiffTimePutLine(diff, unit)
use dc_date_types, only: DC_DIFFTIME
type(DC_DIFFTIME), intent(in) :: diff
integer, intent(in), optional :: unit
end subroutine DCDiffTimePutLine
end interface
interface operator(+)
module procedure dcdate_add_ft
module procedure dcdate_add_tf
module procedure dcdate_add_ff
end interface
interface operator(-)
module procedure dcdate_sub_tf !:doc-priority 40:
module procedure dcdate_sub_tt
module procedure dcdate_sub_ff
end interface
interface operator(*)
module procedure dcdate_mul_if !:doc-priority 51:
module procedure dcdate_mul_fi !:doc-priority 52:
module procedure dcdate_mul_rf !:doc-priority 61:
module procedure dcdate_mul_fr !:doc-priority 62:
module procedure dcdate_mul_df !:doc-priority 71:
module procedure dcdate_mul_fd !:doc-priority 72:
end interface
interface operator(/)
module procedure dcdate_div_fi
module procedure dcdate_div_fr
module procedure dcdate_div_fd
module procedure dcdate_div_ff
end interface
interface mod
module procedure dcdate_mod_ff
end interface
interface operator(==)
module procedure dcdate_eq_tt !:doc-priority 30:
module procedure dcdate_eq_ff !:doc-priority 40:
module procedure dcdate_eq_if !:doc-priority 51:
module procedure dcdate_eq_fi !:doc-priority 52:
module procedure dcdate_eq_rf !:doc-priority 61:
module procedure dcdate_eq_fr !:doc-priority 62:
module procedure dcdate_eq_df !:doc-priority 71:
module procedure dcdate_eq_fd !:doc-priority 72:
end interface
interface operator(>)
module procedure dcdate_gt_tt !:doc-priority 30:
module procedure dcdate_gt_ff !:doc-priority 40:
end interface
interface operator(<)
module procedure dcdate_lt_tt !:doc-priority 30:
module procedure dcdate_lt_ff !:doc-priority 40:
end interface
contains
subroutine dcdate_normalize(day, sec, day_seconds)
!
!=== 日と秒の正規化
!
! このサブルーチンは内部向けなので dc_date モジュール外では
! 極力使用しないでください.
!
! 日付 *day* と秒数 *sec* の正規化を行います. *sec* が *day_seconds*
! (省略される場合は dc_date_types#day_seconds) を超える場合, *day*
! に繰上げを行います.
! また, *sec* と *day* の符号が逆の場合, 同符号になるよう
! 設定します.
!
use dc_date_types, only: day_seconds_default => day_seconds
implicit none
integer, intent(inout):: day
real(DP), intent(inout):: sec
real(DP), intent(in), optional:: day_seconds
integer:: sgn
real(DP):: day_sec
continue
if (present(day_seconds)) then
day_sec = day_seconds
else
day_sec = day_seconds_default
end if
if (abs(sec) > day_sec) then
day = day + int(sec / day_sec)
sec = modulo(sec, day_sec)
end if
if ((sec > 0.0 .and. day < 0) .or. (sec < 0.0 .and. day > 0)) then
sgn = sign(day, 1)
day = day - sgn
sec = sec + sgn * day_sec
endif
end subroutine dcdate_normalize
character(TOKEN) function dcdate_parse_unit(str) result(unit)
!
! このサブルーチンは内部向けなので dc_date モジュール外では
! 極力使用しないでください.
!
! 引数 *str* に与えられた文字列を解釈し, 日時の単位を
! 返します. それぞれ以下の文字列が日時の単位として解釈されます.
! 大文字と小文字は区別されません.
! 返る文字列は以下の文字型の配列の先頭の文字列です.
! (例: *str* に 'hrs.' が与えられる場合, dc_date_types#UNIT_HOUR
! 配列の先頭の文字列 UNIT_HOUR(1) が返ります.)
!
! 年 :: dc_date_types#UNIT_YEAR
! 月 :: dc_date_types#UNIT_MONTH
! 日 :: dc_date_types#UNIT_DAY
! 時 :: dc_date_types#UNIT_HOUR
! 分 :: dc_date_types#UNIT_MIN
! 秒 :: dc_date_types#UNIT_SEC
!
! これらに該当しない文字列を *str* に与えた場合, 空文字が返ります.
!
use dc_types, only: TOKEN
use dc_date_types, only: UNIT_YEAR, UNIT_MONTH, UNIT_DAY, &
& UNIT_HOUR, UNIT_MIN, UNIT_SEC
use dc_string, only: StriEq
implicit none
character(*), intent(in):: str
integer :: unit_str_size, i
continue
unit = adjustl(str)
unit_str_size = size(UNIT_SEC)
do i = 1, unit_str_size
if (StriEq(trim(unit), trim(UNIT_SEC(i)))) then
unit = UNIT_SEC(1)
return
end if
end do
unit_str_size = size(UNIT_MIN)
do i = 1, unit_str_size
if (StriEq(trim(unit), trim(UNIT_MIN(i)))) then
unit = UNIT_MIN(1)
return
end if
end do
unit_str_size = size(UNIT_HOUR)
do i = 1, unit_str_size
if (StriEq(trim(unit), trim(UNIT_HOUR(i)))) then
unit = UNIT_HOUR(1)
return
end if
end do
unit_str_size = size(UNIT_DAY)
do i = 1, unit_str_size
if (StriEq(trim(unit), trim(UNIT_DAY(i)))) then
unit = UNIT_DAY(1)
return
end if
end do
unit_str_size = size(UNIT_MONTH)
do i = 1, unit_str_size
if (StriEq(trim(unit), trim(UNIT_MONTH(i)))) then
unit = UNIT_MONTH(1)
return
end if
end do
unit_str_size = size(UNIT_YEAR)
do i = 1, unit_str_size
if (StriEq(trim(unit), trim(UNIT_YEAR(i)))) then
unit = UNIT_YEAR(1)
return
end if
end do
unit = ''
end function dcdate_parse_unit
type(DC_DATETIME) function dcdate_add_ft(diff, time) result(result)
!
! 2 つの日時 (DC_DATETIME 型) もしくは
! 日時差 (DC_DIFFTIME 型)の加算を行います.
!
implicit none
type(DC_DIFFTIME), intent(in):: diff
type(DC_DATETIME), intent(in):: time
integer:: time_year, time_mon, time_day, time_caltype
real(DP):: time_sec
character(6):: time_zone
continue
call Eval(time, year = time_year, mon = time_mon, day = time_day, &
& sec = time_sec, caltype = time_caltype, zone = time_zone)
call Create(result, year=time_year, mon = time_mon + diff % mon, &
& day = time_day + diff % day, &
& sec = time_sec + diff % sec, &
& caltype = time_caltype, zone = time_zone)
end function dcdate_add_ft
type(DC_DATETIME) function dcdate_add_tf(time, diff) result(result)
use dc_types, only: DP
implicit none
type(DC_DATETIME), intent(in):: time
type(DC_DIFFTIME), intent(in):: diff
continue
result = dcdate_add_ft(diff, time)
end function dcdate_add_tf
type(DC_DIFFTIME) function dcdate_add_ff(diff1, diff2) result(result)
implicit none
type(DC_DIFFTIME), intent(in):: diff1, diff2
continue
result % mon = diff1 % mon + diff2 % mon
result % day = diff1 % day + diff2 % day
result % sec = diff1 % sec + diff2 % sec
result % day_seconds = diff1 % day_seconds
call dcdate_normalize(result % day, result % sec, result % day_seconds)
end function dcdate_add_ff
type(DC_DATETIME) function dcdate_sub_tf(time, diff) result(result)
!
! 2 つの日時 (DC_DATETIME 型) もしくは
! 日時差 (DC_DIFFTIME 型)の減算を行います.
!
implicit none
type(DC_DATETIME), intent(in):: time
type(DC_DIFFTIME), intent(in):: diff
integer:: time_year, time_mon, time_day, time_caltype
real(DP):: time_sec
character(6):: time_zone
continue
call Eval(time, year = time_year, mon = time_mon, day = time_day, &
& sec = time_sec, caltype = time_caltype, zone = time_zone)
call Create(result, year=time_year, mon = time_mon - diff % mon, &
& day = time_day - diff % day, &
& sec = time_sec - diff % sec, &
& caltype = time_caltype, zone = time_zone)
end function dcdate_sub_tf
type(DC_DIFFTIME) function dcdate_sub_tt(time1, time2) result(result)
implicit none
type(DC_DATETIME), intent(in):: time1, time2
continue
result % day = time1 % day - time2 % day
result % sec = time1 % sec - time2 % sec &
& + EvalSec(ZoneToDiff(time1 % zone) - ZoneToDiff(time2 % zone))
result % day_seconds = time1 % day_seconds
call dcdate_normalize(result % day, result % sec, result % day_seconds)
end function dcdate_sub_tt
type(DC_DIFFTIME) function dcdate_sub_ff(diff1, diff2) result(result)
implicit none
type(DC_DIFFTIME), intent(in):: diff1, diff2
continue
result % mon = diff1 % mon - diff2 % mon
result % day = diff1 % day - diff2 % day
result % sec = diff1 % sec - diff2 % sec
result % day_seconds = diff1 % day_seconds
call dcdate_normalize(result % day, result % sec, result % day_seconds)
end function dcdate_sub_ff
type(DC_DIFFTIME) function dcdate_mul_if(factor, diff) result(result)
!
! 日時差 *diff* と *facter* とを乗算した結果を返します.
!
implicit none
integer, intent(in):: factor
type(DC_DIFFTIME), intent(in):: diff
continue
result % mon = factor * diff % mon
result % day = factor * diff % day
result % sec = factor * diff % sec
result % day_seconds = diff % day_seconds
call dcdate_normalize(result % day, result % sec, result % day_seconds)
end function dcdate_mul_if
type(DC_DIFFTIME) function dcdate_mul_fi(diff, factor) result(result)
implicit none
type(DC_DIFFTIME), intent(in):: diff
integer, intent(in):: factor
continue
result = dcdate_mul_if(factor, diff)
end function dcdate_mul_fi
type(DC_DIFFTIME) function dcdate_mul_rf(factor, diff) result(result)
!
! ※ 注意 : 月差を非整数倍すると近似的結果になるおそれがあります
use dc_types, only: DP
implicit none
real, intent(in):: factor
type(DC_DIFFTIME), intent(in):: diff
continue
result = dcdate_mul_df(real(factor, DP), diff)
end function dcdate_mul_rf
type(DC_DIFFTIME) function dcdate_mul_fr(diff, factor) result(result)
!
! ※ 注意 : 月差を非整数倍すると近似的結果になるおそれがあります
implicit none
type(DC_DIFFTIME), intent(in):: diff
real, intent(in):: factor
continue
result = dcdate_mul_rf(factor, diff)
end function dcdate_mul_fr
type(DC_DIFFTIME) function dcdate_mul_df(factor, diff) result(result)
!
! ※ 注意 : 月差を非整数倍すると近似的結果になるおそれがあります
use dc_types, only: DP
use dc_date_types, only: CYCLIC_MDAYS
implicit none
real(DP), intent(in):: factor
type(DC_DIFFTIME), intent(in):: diff
real(DP):: month, day
continue
month = factor * diff % mon
result % mon = int(month)
day = factor * diff % day + int(CYCLIC_MDAYS * (month - result % mon))
result % day = int(day)
result % sec = &
& factor * diff % sec + (day - result % day) * diff % day_seconds
result % day_seconds = diff % day_seconds
call dcdate_normalize(result % day, result % sec, result % day_seconds)
end function dcdate_mul_df
type(DC_DIFFTIME) function dcdate_mul_fd(diff, factor) result(result)
!
! ※ 注意 : 月差を非整数倍すると近似的結果になるおそれがあります
use dc_types, only: DP
implicit none
type(DC_DIFFTIME), intent(in):: diff
real(DP), intent(in):: factor
continue
result = dcdate_mul_df(factor, diff)
end function dcdate_mul_fd
type(DC_DIFFTIME) function dcdate_div_fi(diff, denominator) result(result)
!
! 日時差 *diff* を *denominator* で除算した結果を返します.
!
! ※ 注意 : 月差を除算すると近似的結果になるおそれがあります
use dc_date_types, only: CYCLIC_MDAYS
implicit none
type(DC_DIFFTIME), intent(in):: diff
integer, intent(in):: denominator
continue
result % mon = diff % mon / denominator
! 月からの近似的繰り下がりは日単位でしか行わない
result % day = diff % day / denominator + &
& int((CYCLIC_MDAYS * mod(diff % mon, denominator)) / &
& denominator)
result % sec = diff % sec / denominator + &
& (diff % day_seconds * mod(diff % day, denominator)) / &
& denominator
end function dcdate_div_fi
type(DC_DIFFTIME) function dcdate_div_fr(diff, denominator) result(result)
!
! ※ 注意 : 月差を除算すると近似的結果になるおそれがあります
implicit none
type(DC_DIFFTIME), intent(in):: diff
real, intent(in):: denominator
continue
result = dcdate_div_fd(diff, real(denominator, DP))
end function dcdate_div_fr
type(DC_DIFFTIME) function dcdate_div_fd(diff, denominator) result(result)
!
! ※ 注意 : 月差を除算すると近似的結果になるおそれがあります
use dc_date_types, only: CYCLIC_MDAYS
implicit none
type(DC_DIFFTIME), intent(in):: diff
real(DP), intent(in):: denominator
real(DP):: month, day
continue
month = diff % mon / denominator
result % mon = int(month)
day = diff % day / denominator + int(CYCLIC_MDAYS * (month - result % mon))
result % day = int(day)
result % sec = &
& diff % sec / denominator + (day - result % day) * diff % day_seconds
result % day_seconds = diff % day_seconds
call dcdate_normalize(result % day, result % sec, result % day_seconds)
end function dcdate_div_fd
real(DP) function dcdate_div_ff(diff1, diff2) result(result)
!
! ※ 注意 : 月差と日時の混在する除算は近似的結果になるおそれがあります
use dc_date_types, only: CYCLIC_MDAYS
implicit none
type(DC_DIFFTIME), intent(in):: diff1, diff2
continue
! ゼロ割対応コードが必要か?
result = &
& (diff1 % day_seconds * (CYCLIC_MDAYS * diff1 % mon + diff1 % day) &
& + diff1 % sec) / &
& (diff2 % day_seconds * (CYCLIC_MDAYS * diff2 % mon + diff2 % day) &
& + diff2 % sec)
end function dcdate_div_ff
type(DC_DIFFTIME) function dcdate_mod_ff(diff1, diff2) result(result)
!
! 引数 diff1 を diff2 で除算した際の余りを返します.
!
! ※ 注意: 月差と日時の混在する除算は近似的結果になるおそれがあります
!
use dc_date_types, only: CYCLIC_MDAYS
implicit none
type(DC_DIFFTIME), intent(in):: diff1, diff2
real(DP):: sec1, sec2
continue
result % day_seconds = diff1 % day_seconds
if (diff1 % day == 0 .and. diff2 % day == 0 .and. &
& diff1 % sec == 0.0 .and. diff2 % sec == 0.0) then
result % mon = mod(diff1 % mon, diff2 % mon)
result % day = 0
result % sec = 0.0
else if (diff1 % sec == 0.0 .and. diff2 % sec == 0.0) then
result % mon = 0
result % day = mod((CYCLIC_MDAYS * diff1 % mon + diff1 % day), &
& (CYCLIC_MDAYS * diff2 % mon + diff2 % day))
result % sec = 0.0
else
sec1 = diff1 % day_seconds * (CYCLIC_MDAYS * diff1 % mon + diff1 % day) &
& + diff1 % sec
sec2 = diff2 % day_seconds * (CYCLIC_MDAYS * diff2 % mon + diff2 % day) &
& + diff2 % sec
result % sec = mod(sec1, sec2)
result % day = 0.0
result % mon = 0.0
call dcdate_normalize(result % day, result % sec, result % day_seconds)
endif
end function dcdate_mod_ff
logical function dcdate_gt_tt(time1, time2) result(result)
!
! 2 つの引数の日時を比較します.
! 1 つ目の引数に格納される日時が 2 つ目の引数に格納される日時
! よりも進んでいる場合, .true. が返ります.
!
implicit none
type(DC_DATETIME), intent(in):: time1, time2
integer:: year1, year2
real(DP):: time1_sec, time2_sec
continue
call Eval(time1, year=year1)
call Eval(time2, year=year2)
if (year1 > year2) then
result = .true.
elseif (year1 < year2) then
result = .false.
else
time1_sec = EvalSec(time1) + EvalSec(ZoneToDiff(time1 % zone))
time2_sec = EvalSec(time2) + EvalSec(ZoneToDiff(time2 % zone))
if (time1_sec > time2_sec) then
result = .true.
else
result = .false.
end if
end if
end function dcdate_gt_tt
logical function dcdate_gt_ff(diff1, diff2) result(result)
!
! 2 つの引数の日時差を比較します.
! 1 つ目の引数に格納される日時差が 2 つ目の引数に格納される日時差
! よりも大きい場合, .true. が返ります.
!
use dc_date_types, only: CYCLIC_MDAYS
implicit none
type(DC_DIFFTIME), intent(in):: diff1, diff2
continue
if (EvalSec(diff1) > EvalSec(diff2)) then
result = .true.
else
result = .false.
end if
end function dcdate_gt_ff
logical function dcdate_lt_tt(time1, time2) result(result)
!
! 2 つの引数の日時を比較します.
! 2 つ目の引数に格納される日時が 1 つ目の引数に格納される日時
! よりも進んでいる場合, .true. が返ります.
!
implicit none
type(DC_DATETIME), intent(in):: time1, time2
continue
result = .not. dcdate_gt_tt(time1, time2)
end function dcdate_lt_tt
logical function dcdate_lt_ff(diff1, diff2) result(result)
!
! 2 つの引数の日時差を比較します.
! 2 つ目の引数に格納される日時差が 1 つ目の引数に格納される日時差
! よりも大きい場合, .true. が返ります.
!
implicit none
type(DC_DIFFTIME), intent(in):: diff1, diff2
continue
result = .not. dcdate_gt_ff(diff1, diff2)
end function dcdate_lt_ff
logical function dcdate_eq_tt(time1, time2) result(result)
!
! 2 つの引数の日時を比較します.
! 1 つ目の引数に格納される日時が 2 つ目の引数に格納される日時
! と同じ場合, .true. が返ります.
!
implicit none
type(DC_DATETIME), intent(in):: time1, time2
integer:: year1, year2
real(DP):: time1_sec, time2_sec
continue
call Eval(time1, year=year1)
call Eval(time2, year=year2)
time1_sec = EvalSec(time1) + EvalSec(ZoneToDiff(time1 % zone))
time2_sec = EvalSec(time2) + EvalSec(ZoneToDiff(time2 % zone))
if (year1 == year2 .and. time1_sec == time2_sec) then
result = .true.
else
result = .false.
end if
end function dcdate_eq_tt
logical function dcdate_eq_ff(diff1, diff2) result(result)
!
! 2 つの引数の日時差を比較します.
! 1 つ目の引数に格納される日時差が 2 つ目の引数に格納される日時差
! と同じ場合, .true. が返ります.
!
implicit none
type(DC_DIFFTIME), intent(in):: diff1, diff2
continue
if (EvalSec(diff1) == EvalSec(diff2)) then
result = .true.
else
result = .false.
end if
end function dcdate_eq_ff
logical function dcdate_eq_if(i, diff) result(result)
!
! 引数 *diff* の日時差が *i* と等しいかどうかを比較します. *diff*
! を秒数に換算した値と *i* とが等しい場合, .true. が返ります.
!
implicit none
type(DC_DIFFTIME), intent(in):: diff
integer, intent(in):: i
continue
result = dcdate_eq_rf(real(i), diff)
end function dcdate_eq_if
logical function dcdate_eq_fi(diff, i) result(result)
implicit none
type(DC_DIFFTIME), intent(in):: diff
integer, intent(in):: i
continue
result = dcdate_eq_if(i, diff)
end function dcdate_eq_fi
logical function dcdate_eq_rf(r, diff) result(result)
!
! 引数 *diff* の日時差が *r* と等しいかどうかを比較します. *diff*
! を秒数に換算した値と *r* とが等しい場合, .true. が返ります.
!
implicit none
type(DC_DIFFTIME), intent(in):: diff
real, intent(in):: r
continue
if (real(EvalSec(diff)) == r) then
result = .true.
else
result = .false.
end if
end function dcdate_eq_rf
logical function dcdate_eq_fr(diff, r) result(result)
implicit none
type(DC_DIFFTIME), intent(in):: diff
real, intent(in):: r
continue
result = dcdate_eq_rf(r, diff)
end function dcdate_eq_fr
logical function dcdate_eq_df(d, diff) result(result)
!
! 引数 *diff* の日時差が *d* と等しいかどうかを比較します. *diff*
! を秒数に換算した値と *d* とが等しい場合, .true. が返ります.
!
use dc_types, only: DP
implicit none
type(DC_DIFFTIME), intent(in):: diff
real(DP), intent(in):: d
continue
if (real(EvalSec(diff)) == d) then
result = .true.
else
result = .false.
end if
end function dcdate_eq_df
logical function dcdate_eq_fd(diff, d) result(result)
use dc_types, only: DP
implicit none
type(DC_DIFFTIME), intent(in):: diff
real(DP), intent(in):: d
continue
result = dcdate_eq_df(d, diff)
end function dcdate_eq_fd
end module dc_date