FMS 2025.01.02-dev
Flexible Modeling System
All Classes Namespaces Files Functions Variables Modules Pages
time_manager.F90
1!***********************************************************************
2!* GNU Lesser General Public License
3!*
4!* This file is part of the GFDL Flexible Modeling System (FMS).
5!*
6!* FMS is free software: you can redistribute it and/or modify it under
7!* the terms of the GNU Lesser General Public License as published by
8!* the Free Software Foundation, either version 3 of the License, or (at
9!* your option) any later version.
10!*
11!* FMS is distributed in the hope that it will be useful, but WITHOUT
12!* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13!* FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
14!* for more details.
15!*
16!* You should have received a copy of the GNU Lesser General Public
17!* License along with FMS. If not, see <http://www.gnu.org/licenses/>.
18!***********************************************************************
19!> @defgroup time_manager_mod time_manager_mod
20!> @ingroup time_manager
21!> @brief A software package that provides a set of simple interfaces for
22!! modelers to perform computations related to time and dates.
23!!
24!! Optional error flag can be used in calling arguments of public routines.
25!! This allows the using routine to terminate the program. It is likely that more
26!! diagnostic information is available from the user than from time_manager alone.
27!! If the error flag is present then it is the responsibility of the using
28!! routine to test it and add additional information to the error message.
29!!
30!! Calendar specific routines are private.
31!! They are not used, and should not be used, by any using code.
32!!
33!! The module defines a type that can be used to represent discrete
34!! times (accurate to one second) and to map these times into dates
35!! using a variety of calendars. A time is mapped to a date by
36!! representing the time with respect to an arbitrary base date (refer
37!! to <B>NOTES</B> section for the base date setting).
38!!
39!! The time_manager provides a single defined type, time_type, which is
40!! used to store time and date quantities. A time_type is a positive
41!! definite quantity that represents an interval of time. It can be
42!! most easily thought of as representing the number of seconds in some
43!! time interval. A time interval can be mapped to a date under a given
44!! calendar definition by using it to represent the time that has passed
45!! since some base date. A number of interfaces are provided to operate
46!! on time_type variables and their associated calendars. Time intervals
47!! can be as large as n days where n is the largest number represented by
48!! the default integer type on a compiler. This is typically considerably
49!! greater than 10 million years (assuming 32 bit integer representation)
50!! which is likely to be adequate for most applications. The description
51!! of the interfaces is separated into two sections. The first deals with
52!! operations on time intervals while the second deals with operations
53!! that convert time intervals to dates for a given calendar.
54!!
55!! The smallest increment of time is referred to as a tick.
56!! A tick cannot be larger than 1 second, which also is the default.
57!! The number of ticks per second is set via pubic subroutine set_ticks_per_second.
58!! For example, ticks_per_second = 1000 will set the tick to one millisecond.
59
60!> @addtogroup time_manager_mod
61!> @{
62module time_manager_mod
63
64
65use platform_mod, only: r8_kind, r4_kind
66use constants_mod, only: rseconds_per_day=>seconds_per_day
67use fms_mod, only: error_mesg, fatal, warning, write_version_number, stdout
68
69implicit none
70private
71
72! Module defines a single type
73public time_type
74
75! Operators defined on time_type
76public operator(+), operator(-), operator(*), operator(/), &
77 operator(>), operator(>=), operator(==), operator(/=), &
78 operator(<), operator(<=), operator(//), assignment(=)
79
80! Subroutines and functions operating on time_type
83public time_list_error
84
85! List of available calendar types
86public thirty_day_months, julian, gregorian, noleap, no_calendar, invalid_calendar
87
88! Subroutines and functions involving relations between time and calendar
93public set_date
94public get_date
95public increment_date
96public decrement_date
97public days_in_month
98public leap_year
99public length_of_year
100public days_in_year
101public day_of_year
102public month_name
103
105
106! Subroutines for printing version number and time type
108
109! The following exist only for interpolator.F90
110! interpolator.F90 uses them to do a calendar conversion,
111! which is also done by get_cal_time. interpolator.F90
112! should be modified to use get_cal_time instead.
113! After interpolator.F90 is fixed, these can be removed
114! and the corresponding private routines can be renamed.
115! (e.g., rename set_date_julian_private to be just set_date_julian)
116public :: set_date_julian, set_date_no_leap, get_date_julian, get_date_no_leap
117
118public :: date_to_string
119
120!====================================================================
121
122! Global data to define calendar type
123integer, parameter :: THIRTY_DAY_MONTHS = 1, julian = 2, &
124 gregorian = 3, noleap = 4, &
125 no_calendar = 0, invalid_calendar =-1
126integer, private :: calendar_type = no_calendar
127integer, parameter :: max_type = 4
128
129! Define number of days per month
130integer, private :: days_per_month(12) = (/31,28,31,30,31,30,31,31,30,31,30,31/)
131integer, parameter :: seconds_per_day = rseconds_per_day ! This should automatically cast real to integer
132integer, parameter :: days_in_400_year_period = 146097 !< Used only for gregorian
133integer,parameter :: do_floor = 0
134integer,parameter :: do_nearest = 1
135
136!> @}
137
138!> @brief Type to represent amounts of time.
139!> Implemented as seconds and days to allow for larger intervals.
140!> @ingroup time_manager_mod
142 private
143 integer:: seconds
144 integer:: days
145 integer:: ticks
146end type time_type
147
148!> Operator override interface for use with @ref time_type
149!> @ingroup time_manager_mod
150interface operator (+); module procedure time_plus; end interface
151!> Operator override interface for use with @ref time_type
152!> @ingroup time_manager_mod
153interface operator (-); module procedure time_minus; end interface
154!> Operator override interface for use with @ref time_type
155!> @ingroup time_manager_mod
156interface operator (*); module procedure time_scalar_mult
157 module procedure scalar_time_mult; end interface
158!> Operator override interface for use with @ref time_type
159!> @ingroup time_manager_mod
160interface operator (/); module procedure time_scalar_divide
161 module procedure time_divide; end interface
162!> Operator override interface for use with @ref time_type
163!> @ingroup time_manager_mod
164interface operator (>); module procedure time_gt; end interface
165!> Operator override interface for use with @ref time_type
166!> @ingroup time_manager_mod
167interface operator (>=); module procedure time_ge; end interface
168!> Operator override interface for use with @ref time_type
169!> @ingroup time_manager_mod
170interface operator (<); module procedure time_lt; end interface
171!> Operator override interface for use with @ref time_type
172!> @ingroup time_manager_mod
173interface operator (<=); module procedure time_le; end interface
174!> Operator override interface for use with @ref time_type
175!> @ingroup time_manager_mod
176interface operator (==); module procedure time_eq; end interface
177!> Operator override interface for use with @ref time_type
178!> @ingroup time_manager_mod
179interface operator (/=); module procedure time_ne; end interface
180!> Operator override interface for use with @ref time_type
181!> @ingroup time_manager_mod
182interface operator (//); module procedure time_real_divide; end interface
183!> Operator override interface for use with @ref time_type
184!> @ingroup time_manager_mod
185interface assignment(=); module procedure time_assignment; end interface
186
187!======================================================================
188
189!> @brief Given some number of seconds and days, returns the
190!! corresponding time_type.
191!!
192!> Given some number of seconds and days, returns the
193!! corresponding time_type. set_time has two forms;
194!! one accepts integer input, the other a character string with the day and second counts.
195!! For the first form, there are no restrictions on the range of the inputs,
196!! except that the result must be positive time.
197!! e.g. days=-1, seconds=86401 is acceptable.
198!! For the second form, days and seconds must both be positive.
199!!
200!! <br>Example usage:
201!! @code{.F90}
202!! type(time_type) :: time1, time2
203!! time1 = set_time(seconds, days, ticks, err_msg)
204!! time2 = set_time("100 43200", err_msg, allow_rounding)
205!! @endcode
206!> @ingroup time_manager_mod
207interface set_time
208 module procedure set_time_i, set_time_c
209end interface
210
211!> @brief Given an input date in year, month, days, etc., creates a
212!! time_type that represents this time interval from the
213!! internally defined base date.
214!!
215!> Given a date, computes the corresponding time given the selected
216!! date time mapping algorithm. Note that it is possible to specify
217!! any number of illegal dates; these should be checked for and generate
218!! errors as appropriate.
219!!
220!! <br>Example usage:
221!! <br> Integer input
222!! @code{.F90} time = set_date(year, month, day, hours, minute, second, tick, err_msg) @endcode
223!! <br> String input
224!! @code{.F90} time = set_date(time_string, zero_year_warning, err_msg, allow_rounding) @endcode
225!!
226!! @param time_string A character string containing a date formatted
227!! according to CF conventions. e.g. '1980-12-31 23:59:59.9'
228!!
229!! @param zero_year_warning If the year number is zero, it will be silently changed to one,
230!! unless zero_year_warning=.true., in which case a WARNING message
231!! will also be issued.
232!!
233!! @param allow_rounding When .true., any fractions of a second will be rounded off to the nearest
234!! tick. When .false., it is a fatal error if the second fraction cannot be exactly
235!! represented by a number of ticks.
236!!
237!! @param err_msg When present, and when non-blank, a fatal error condition as been detected.
238!! The string itself is an error message.
239!! It is recommended that, when err_msg is present in the call
240!! to this routine, the next line of code should be something
241!! similar to this:
242!! @code{.F90}
243!! if(err_msg /= '') call error_mesg('my_routine','additional info: '//trim(err_msg) ,FATAL)
244!! @endcode
245!!
246!> @ingroup time_manager_mod
247interface set_date
248 module procedure set_date_i, set_date_c
249end interface
250
251!> Wrapper for the real to time interface
252!! Takes seconds as reals to convert to a time_type representation of an interval
253!! r4 versions just casts to r8
255 module procedure real4_to_time_type
256 module procedure real8_to_time_type
257end interface
258
259!> @addtogroup time_manager_mod
260!> @{
261
262!======================================================================
263
264! Include variable "version" to be written to log file.
265#include<file_version.h>
266logical :: module_is_initialized = .false.
267
268!======================================================================
269
270! A tick is the smallest increment of time.
271! That is, smallest increment of time = (1/ticks_per_second) seconds
272
273integer :: ticks_per_second = 1
274
275!======================================================================
276contains
277
278! First define all operations on time intervals independent of calendar
279
280!> Returns a time interval corresponding to this number of days, seconds, and ticks.
281!! days, seconds and ticks may be negative, but resulting time must be positive.
282 function set_time_private(seconds, days, ticks, Time_out, err_msg)
283
284! -- pjp --
285! To understand why inputs may be negative,
286! one needs to understand the intrinsic function "modulo".
287! The expanation below is copied from a web page on fortran 90
288
289! In addition, CEILING, FLOOR and MODULO have been added to Fortran 90.
290! Only the last one is difficult to explain, which is most easily done with the examples from ISO (1991)
291
292! MOD (8,5) gives 3 MODULO (8,5) gives 3
293! MOD (-8,5) gives -3 MODULO (-8,5) gives 2
294! MOD (8,-5) gives 3 MODULO (8,-5) gives -2
295! MOD (-8,-5) gives -3 MODULO (-8,-5) gives -3
296
297! I don't think it is difficult to explain.
298! I think that is it sufficient to say this:
299! "The result of modulo(n,m) has the sign of m"
300! -- pjp --
301
302 logical :: set_time_private
303 integer, intent(in) :: seconds, days, ticks
304 type(time_type), intent(out) :: time_out
305 character(len=*), intent(out) :: err_msg
306 integer :: seconds_new, days_new, ticks_new
307
308 seconds_new = seconds + floor(real(ticks, r8_kind)/real(ticks_per_second, r8_kind))
309 ticks_new = modulo(ticks,ticks_per_second)
310 days_new = days + floor(real(seconds_new, r8_kind)/real(seconds_per_day, r8_kind))
311 seconds_new = modulo(seconds_new,seconds_per_day)
312
313 if ( seconds_new < 0 .or. ticks_new < 0) then
314 call error_mesg('function set_time_i','Bad result for time. Contact those responsible for maintaining time_manager'&
315 & ,fatal)
316 endif
317
318 if(days_new < 0) then
319 write(err_msg,'(a,i6,a,i6,a,i6)') 'time is negative. days=',days_new,' seconds=',seconds_new,' ticks=',ticks_new
320 set_time_private = .false.
321 else
322 time_out%days = days_new
323 time_out%seconds = seconds_new
324 time_out%ticks = ticks_new
325 err_msg = ''
326 set_time_private = .true.
327 endif
328
329 end function set_time_private
330!---------------------------------------------------------------------------
331 !> @brief Returns a time_type set to the given amount of time via integer amounts.
332 function set_time_i(seconds, days, ticks, err_msg)
333 type(time_type) :: set_time_i
334 integer, intent(in) :: seconds !< A number of seconds
335 integer, intent(in), optional :: days !< A number of days
336 integer, intent(in), optional :: ticks !< A number of ticks
337 character(len=*), intent(out), optional :: err_msg !< When present, and when non-blank, a fatal
338 !! error condition as been detected. The string itself is an error message.
339 !! It is recommended that, when err_msg is present in the call
340 !! to this routine, the next line of code should be something
341 !! similar to this:
342 !! if(err_msg /= '') call error_mesg('my_routine','additional info: '//trim(err_msg),FATAL)
343 character(len=128) :: err_msg_local
344 integer :: odays, oticks
345
346 if(.not.module_is_initialized) call time_manager_init
347
348 odays = 0; if(present(days)) odays = days
349 oticks = 0; if(present(ticks)) oticks = ticks
350 if(present(err_msg)) err_msg = ''
351
352 if(.not.set_time_private(seconds, odays, oticks, set_time_i, err_msg_local)) then
353 if(error_handler('function set_time_i', trim(err_msg_local), err_msg)) return
354 endif
355
356 end function set_time_i
357!---------------------------------------------------------------------------
358
359 !> @brief Returns a time_type set to the given amount of time via a string
360 function set_time_c(string, err_msg, allow_rounding)
361
362 type(time_type) :: set_time_c
363 character(len=*), intent(in) :: string !< Contains days and seconds separated by a single blank.
364 !! days must be integer, seconds may be integer or real.
365 !! Examples: '100 43200' '100 43200.50'
366 character(len=*), intent(out), optional :: err_msg !< When present, and when non-blank, a fatal
367 !! error condition as been detected. The string itself is an error message.
368 !! It is recommended that, when err_msg is present in the call
369 !! to this routine, the next line of code should be something
370 !! similar to this:
371 !! if(err_msg /= '') call error_mesg('my_routine','additional info: '//trim(err_msg),FATAL)
372 logical, intent(in), optional :: allow_rounding !< When .true., any fractions of a second will be
373 !! rounded off to the nearest tick. When .false., it is a fatal error
374 !! if the second fraction cannot be exactly represented by a number of ticks.
375
376 character(len=4) :: formt='(i )'
377 integer :: i1, i2, i3, day, second, tick, nsps
378 character(len=32) :: string_sifted_left
379 character(len=128) :: err_msg_local
380 logical :: allow_rounding_local
381
382 if(.not.module_is_initialized) call time_manager_init
383 if(present(err_msg)) err_msg = ''
384 allow_rounding_local=.true.; if(present(allow_rounding)) allow_rounding_local=allow_rounding
385
386 err_msg_local = 'Form of character time stamp is incorrect. The character time stamp is: '//trim(string)
387
388 string_sifted_left = adjustl(string)
389 i1 = index(trim(string_sifted_left),' ')
390 if(i1 == 0) then
391 if(error_handler('function set_time_c', err_msg_local, err_msg)) return
392 endif
393 if(index(string,'-') /= 0 .or. index(string,':') /= 0) then
394 if(error_handler('function set_time_c', err_msg_local, err_msg)) return
395 endif
396
397 i2 = index(trim(string_sifted_left),'.')
398 i3 = len_trim(cut0(string_sifted_left))
399
400 if(i2 /= 0) then ! There is no decimal point
401 ! Check that decimal is on seconds (not days)
402 if(i2 < i1) then
403 if(error_handler('function set_time_c', err_msg_local, err_msg)) return
404 endif
405 endif
406 write(formt(3:3),'(i1)') i1-1
407 read(string_sifted_left(1:i1-1),formt) day
408
409 if(i2 == 0) then ! There is no decimal point
410 write(formt(3:3),'(i1)') i3-i1
411 read(string_sifted_left(i1+1:i3),formt) second
412 tick = 0
413 else ! There is a decimal point
414 ! nsps = spaces occupied by whole number of seconds
415 nsps = i2-i1-1
416 if(nsps == 0) then
417 second = 0
418 else
419 write(formt(3:3),'(i1)') nsps
420 read(string_sifted_left(i1+1:i2-1),formt) second
421 endif
422
423 if(.not.get_tick_from_string(string_sifted_left(i2:i3), err_msg_local, allow_rounding_local, tick)) then
424 if(error_handler('function set_time_c', err_msg_local, err_msg)) return
425 endif
426 ! If tick has been rounded up to ticks_per_second, then bump up second.
427 if(tick == ticks_per_second) then
428 second = second + 1
429 tick = 0
430 endif
431 endif
432
433 if(.not.set_time_private(second, day, tick, set_time_c, err_msg_local)) then
434 if(error_handler('function set_time_c', err_msg_local, err_msg)) return
435 endif
436
437 end function set_time_c
438!---------------------------------------------------------------------------
439
440 function get_tick_from_string(string, err_msg, allow_rounding, tick)
441
442 logical :: get_tick_from_string
443 character(len=*), intent(in) :: string
444 character(len=*), intent(out) :: err_msg
445 logical, intent(in) :: allow_rounding
446 integer, intent(out) :: tick
447
448 character(len=4) :: formt='(i )'
449 integer :: i3, nspf, fraction, magnitude, tpsfrac
450
451 err_msg = ''
452 get_tick_from_string = .true.
453 i3 = len_trim(string)
454 nspf = i3 - 1 ! nspf = spaces occupied by fractional seconds, excluding decimal point
455 if(nspf == 0) then
456 tick = 0 ! Nothing to the right of the decimal point
457 else
458 write(formt(3:3),'(i1)') nspf
459 read(string(2:i3),formt) fraction
460 if(fraction == 0) then
461 tick = 0 ! All zeros to the right of the decimal point
462 else
463 magnitude = 10**nspf
464 tpsfrac = ticks_per_second*fraction
465 if(allow_rounding) then
466 tick = nint((real(tpsfrac, r8_kind)/real(magnitude, r8_kind)))
467 else
468 if(modulo(tpsfrac,magnitude) == 0) then
469 tick = tpsfrac/magnitude
470 else
471 write(err_msg,'(a,i6)') 'Second fraction cannot be exactly represented with ticks. '// &
472 'fraction='//trim(string)//' ticks_per_second=',ticks_per_second
473 get_tick_from_string = .false.
474 endif
475 endif
476 endif
477 endif
478
479 end function get_tick_from_string
480
481!> @brief Returns days and seconds ( < 86400 ) corresponding to a time.
482!! <TT>err_msg</TT> should be checked for any errors.
483!!
484!> @param time A @ref time_type interval to get days and seconds from
485!! @param [out] seconds The number of seconds
486!! @param [out] days The number of seconds
487!! @param [out] ticks The number of ticks
488!! @param [out] err_msg Contains an error message on failure
489!! <br>Example usage:
490!! @code{.F90} get_time(time, seconds, days, ticks, err_msg) @endcode
491subroutine get_time(Time, seconds, days, ticks, err_msg)
492
493! Returns days and seconds ( < 86400 ) corresponding to a time.
494
495type(time_type), intent(in) :: time
496integer, intent(out) :: seconds
497integer, intent(out), optional :: days, ticks
498character(len=*), intent(out), optional :: err_msg
499character(len=128) :: err_msg_local
500
501if(.not.module_is_initialized) call time_manager_init
502if(present(err_msg)) err_msg = ''
503
504seconds = time%seconds
505
506if(present(ticks)) then
507 ticks = time%ticks
508else
509 if(time%ticks /= 0) then
510 err_msg_local = 'subroutine get_time: ticks must be present when time has a second fraction'
511 if(error_handler('subroutine get_time', err_msg_local, err_msg)) return
512 endif
513endif
514
515if (present(days)) then
516 days = time%days
517else
518 if (time%days > (huge(seconds) - seconds)/seconds_per_day) then
519 err_msg_local = 'Integer overflow in seconds. Optional argument days must be present.'
520 if(error_handler('subroutine get_time', err_msg_local, err_msg)) return
521 endif
522 seconds = seconds + time%days * seconds_per_day
523endif
524
525end subroutine get_time
526
527 !> @brief Increments a time by seconds and days.
528 !!
529 !> Given a time and an increment of days and seconds, returns
530 !! a new time_type that represents the given time after the given increment.
531 !! @returns incremented @ref time_type
532 function increment_time(Time, seconds, days, ticks, err_msg, allow_neg_inc)
533
535 type(time_type), intent(in) :: time !< A time interval
536 integer, intent(in) :: seconds !< Increment of seconds
537 integer, intent(in), optional :: days, ticks !< Increment of days and ticks
538 character(len=*), intent(out), optional :: err_msg !< When present and non-blank, a fatal error
539 !! condition has been detected, with the string itself
540 !! as the error message.
541 logical, intent(in), optional :: allow_neg_inc !< When false, negative increments give fatal errors
542 !! Defaults to true.
543 integer :: odays, oticks
544 character(len=128) :: err_msg_local
545 logical :: allow_neg_inc_local
546
547 odays = 0; if(present(days)) odays = days
548 oticks = 0; if(present(ticks)) oticks = ticks
549 allow_neg_inc_local=.true.; if(present(allow_neg_inc)) allow_neg_inc_local=allow_neg_inc
550
551 if(.not.allow_neg_inc_local) then
552 if(seconds < 0 .or. odays < 0 .or. oticks < 0) then
553 write(err_msg_local,10) seconds, odays, oticks
554 10 format('One or more time increments are negative: seconds=',i6,' days=',i6,' ticks=',i6)
555 if(error_handler('function increment_time', err_msg_local, err_msg)) return
556 endif
557 endif
558
559 if(.not.increment_time_private(time, seconds, odays, oticks, increment_time, err_msg_local)) then
560 if(error_handler('function increment_time', err_msg_local, err_msg)) return
561 endif
562
563 end function increment_time
564!--------------------------------------------------------------------------
565
566 !> Increments a time by seconds, days and ticks.
567 function increment_time_private(Time_in, seconds, days, ticks, Time_out, err_msg)
568
569
570 logical :: increment_time_private
571 type(time_type), intent(in) :: time_in
572 integer, intent(in) :: seconds, days, ticks
573 type(time_type), intent(out) :: time_out
574 character(len=*), intent(out) :: err_msg
575
576! Watch for immediate overflow on days or seconds
577 if(days >= huge(days) - time_in%days) then
578 err_msg = 'Integer overflow in days in increment_time'
579 increment_time_private = .false.
580 return
581 endif
582 if(seconds >= huge(seconds) - time_in%seconds) then
583 err_msg = 'Integer overflow in seconds in increment_time'
584 increment_time_private = .false.
585 return
586 endif
587
588 increment_time_private = set_time_private(time_in%seconds+seconds, time_in%days+days, time_in%ticks+ticks, &
589 & time_out, err_msg)
590
591 end function increment_time_private
592
593!--------------------------------------------------------------------------
594
595!> @brief Decrements a time by seconds and days.
596!!
597!> Given a time and a decrement of days and seconds, returns
598!! a time that subtracts this decrement from an input time.
599!> @returns A time that suvtracts this decrement from an input time. A negative result is a fatal error.
600function decrement_time(Time, seconds, days, ticks, err_msg, allow_neg_inc)
601
602! Decrements a time by seconds, days and ticks.
603
605type(time_type), intent(in) :: time !< A time interval
606integer, intent(in) :: seconds !< Decrement of seconds
607integer, intent(in), optional :: days, ticks !< Decrement of days and ticks
608character(len=*), intent(out), optional :: err_msg !< Present and non-blank when a fatal error has
609 !! occured, holds the error message.
610logical, intent(in), optional :: allow_neg_inc !< Throws fatal warning when set to false if
611 !! negative values are used to decrement. Default
612 !! is true.
613
614integer :: odays, oticks
615character(len=128) :: err_msg_local
616logical :: allow_neg_inc_local
617
618odays = 0; if (present(days)) odays = days
619oticks = 0; if (present(ticks)) oticks = ticks
620allow_neg_inc_local=.true.; if(present(allow_neg_inc)) allow_neg_inc_local=allow_neg_inc
621
622if(.not.allow_neg_inc_local) then
623 if(seconds < 0 .or. odays < 0 .or. oticks < 0) then
624 write(err_msg_local,10) seconds,odays,oticks
625 10 format('One or more time increments are negative: seconds=',i6,' days=',i6,' ticks=',i6)
626 if(error_handler('function decrement_time', err_msg_local, err_msg)) return
627 endif
628endif
629
630 if(.not.increment_time_private(time, -seconds, -odays, -oticks, decrement_time, err_msg_local)) then
631 if(error_handler('function decrement_time', err_msg_local, err_msg)) return
632 endif
633
634end function decrement_time
635
636!--------------------------------------------------------------------------
637
638!> @brief Returns true if time1 > time2
639function time_gt(time1, time2)
640
641logical :: time_gt
642type(time_type), intent(in) :: time1 !< time interval to compare
643type(time_type), intent(in) :: time2 !< time interval to compare
644
645time_gt = (time1%days > time2%days)
646if(time1%days == time2%days) then
647 if(time1%seconds == time2%seconds) then
648 time_gt = (time1%ticks > time2%ticks)
649 else
650 time_gt = (time1%seconds > time2%seconds)
651 endif
652endif
653
654end function time_gt
655
656!--------------------------------------------------------------------------
657
658!> Returns true if time1 >= time2.
659function time_ge(time1, time2)
660logical :: time_ge
661type(time_type), intent(in) :: time1 !< time interval to compare
662type(time_type), intent(in) :: time2 !< time interval to compare
663
664time_ge = (time_gt(time1, time2) .or. time_eq(time1, time2))
665
666end function time_ge
667
668!--------------------------------------------------------------------------
669
670!> Returns true if time1 < time2
671function time_lt(time1, time2)
672logical :: time_lt
673type(time_type), intent(in) :: time1 !< time interval to compare
674type(time_type), intent(in) :: time2 !< time interval to compare
675
676time_lt = (time1%days < time2%days)
677if(time1%days == time2%days)then
678 if(time1%seconds == time2%seconds) then
679 time_lt = (time1%ticks < time2%ticks)
680 else
681 time_lt = (time1%seconds < time2%seconds)
682 endif
683endif
684end function time_lt
685
686!--------------------------------------------------------------------------
687
688!> Returns true if time1 <= time2
689function time_le(time1, time2)
690logical :: time_le
691type(time_type), intent(in) :: time1 !< time interval to compare
692type(time_type), intent(in) :: time2 !< time interval to compare
693
694time_le = (time_lt(time1, time2) .or. time_eq(time1, time2))
695
696end function time_le
697
698!--------------------------------------------------------------------------
699
700!> Returns true if time1 == time2
701function time_eq(time1, time2)
702logical :: time_eq
703type(time_type), intent(in) :: time1 !< time interval to compare
704type(time_type), intent(in) :: time2 !< time interval to compare
705
706if(.not.module_is_initialized) call time_manager_init
707
708time_eq = (time1%seconds == time2%seconds .and. time1%days == time2%days &
709 .and. time1%ticks == time2%ticks)
710
711end function time_eq
712
713!--------------------------------------------------------------------------
714
715!> Returns true if time1 /= time2
716function time_ne(time1, time2)
717logical :: time_ne
718type(time_type), intent(in) :: time1 !< time interval to compare
719type(time_type), intent(in) :: time2 !< time interval to compare
720
721time_ne = (.not. time_eq(time1, time2))
722
723end function time_ne
724
725!-------------------------------------------------------------------------
726
727!> Returns sum of two time_types
728function time_plus(time1, time2)
729type(time_type) :: time_plus
730type(time_type), intent(in) :: time1 !< time interval to add
731type(time_type), intent(in) :: time2 !< time interval to add
732
733if(.not.module_is_initialized) call time_manager_init
734
735time_plus = increment_time(time1, time2%seconds, time2%days, time2%ticks)
736
737end function time_plus
738
739!-------------------------------------------------------------------------
740
741!> Returns difference of two time_types. WARNING: a time type is positive
742!! so by definition time1 - time2 is the same as time2 - time1.
743function time_minus(time1, time2)
744type(time_type) :: time_minus
745type(time_type), intent(in) :: time1 !< time interval to subtract
746type(time_type), intent(in) :: time2 !< time interval to subtract
747
748if(.not.module_is_initialized) call time_manager_init
749
750if(time1 > time2) then
751 time_minus = decrement_time(time1, time2%seconds, time2%days, time2%ticks)
752else
753 time_minus = decrement_time(time2, time1%seconds, time1%days, time1%ticks)
754endif
755
756end function time_minus
757
758!--------------------------------------------------------------------------
759
760!> Returns time multiplied by integer factor n
761function time_scalar_mult(time, n)
763type(time_type), intent(in) :: time !< time interval to multply
764integer, intent(in) :: n !< factor to multiply by
765integer :: days, seconds, ticks, num_sec
766real(r8_kind) :: sec_prod, tick_prod
767
768if(.not.module_is_initialized) call time_manager_init
769
770!! Multiplying here in a reasonable fashion to avoid overflow is tricky
771!! Could multiply by some large factor n, and seconds could be up to 86399
772!! Need to avoid overflowing integers and wrapping around to negatives
773!! ticks could be up to ticks_per_second-1
774
775tick_prod = real(time%ticks, r8_kind) * real(n, r8_kind)
776num_sec = int(tick_prod/real(ticks_per_second, r8_kind))
777sec_prod = real(time%seconds, r8_kind) * real(n, r8_kind) + real(num_sec, r8_kind)
778ticks = int(tick_prod) - (num_sec * ticks_per_second)
779
780!! If sec_prod is large compared to precision of double precision, things
781!! can go bad. Need to warn and abort on this.
782!! The same is true of tick_prod but is is more likely to happen to sec_prod,
783!! so let's just test sec_prod. (A test of tick_prod would be necessary only
784!! if ticks_per_second were greater than seconds_per_day)
785if(sec_prod /= 0.0_r8_kind) then
786 if(log10(sec_prod) > precision(sec_prod) - 3) call error_mesg('time_scalar_mult', &
787 'Insufficient precision to handle scalar product in time_scalar_mult; contact developer',fatal)
788end if
789
790days = int(sec_prod / real(seconds_per_day, r8_kind))
791seconds = int(sec_prod - real(days, r8_kind) * real(seconds_per_day, r8_kind))
792
793time_scalar_mult = set_time(seconds, time%days * n + days, ticks)
794
795end function time_scalar_mult
796
797!-------------------------------------------------------------------------
798
799!> Returns time multipled by integer factor n
800function scalar_time_mult(n, time)
801
803type(time_type), intent(in) :: time !< a time interval
804integer, intent(in) :: n !< factor to mulitply by
805
807
808end function scalar_time_mult
809
810!-------------------------------------------------------------------------
811
812!> Returns the largest integer, n, for which time1 >= time2 * n.
813function time_divide(time1, time2)
814
815integer :: time_divide
816type(time_type), intent(in) :: time1 !< a time interval (dividend)
817type(time_type), intent(in) :: time2 !< a time interval (divisor)
818real(r8_kind) :: d1, d2
819
820if(.not.module_is_initialized) call time_manager_init
821
822! Convert time intervals to floating point days; risky for general performance?
823d1 = real(time1%days, r8_kind) * real(seconds_per_day, r8_kind) + real(time1%seconds, r8_kind) &
824 + real(time1%ticks, r8_kind)/real(ticks_per_second, r8_kind)
825d2 = real(time2%days, r8_kind) * real(seconds_per_day, r8_kind) + real(time2%seconds, r8_kind) &
826 + real(time2%ticks,r8_kind)/real(ticks_per_second, r8_kind)
827
828! Get integer quotient of this, check carefully to avoid round-off problems.
829time_divide = int(d1 / d2)
830
831! Verify time_divide*time2 is <= time1 and (time_divide + 1)*time2 is > time1
832if(time_divide * time2 > time1 .or. (time_divide + 1) * time2 <= time1) &
833 call error_mesg('time_divide',' quotient error :: notify developer',fatal)
834
835end function time_divide
836
837!> Returns the double precision quotient of two times
838function time_real_divide(time1, time2)
839
840real(r8_kind) :: time_real_divide
841type(time_type), intent(in) :: time1 !< a time interval (dividend)
842type(time_type), intent(in) :: time2 !< a time interval (divisor)
843real(r8_kind) :: d1, d2
844
845if(.not.module_is_initialized) call time_manager_init
846
847! Convert time intervals to floating point seconds; risky for general performance?
848d1 = real(time1%days, r8_kind) * real(seconds_per_day, r8_kind) + real(time1%seconds, r8_kind) + &
849 real(time1%ticks, r8_kind)/real(ticks_per_second, r8_kind)
850d2 = real(time2%days, r8_kind) * real(seconds_per_day, r8_kind) + real(time2%seconds, r8_kind) + &
851 real(time2%ticks, r8_kind)/real(ticks_per_second, r8_kind)
852
853time_real_divide = d1 / d2
854
855end function time_real_divide
856
857!-------------------------------------------------------------------------
858
859!> Assigns all components of the time_type variable on
860!! RHS to same components of time_type variable on LHS.
861subroutine time_assignment(time1, time2)
862type(time_type), intent(out) :: time1
863type(time_type), intent(in) :: time2
864 time1%seconds = time2%seconds
865 time1%days = time2%days
866 time1%ticks = time2%ticks
867end subroutine time_assignment
868
869!> Converts time to seconds and returns it as a real number
870function time_type_to_real(time)
871
872real(kind=r8_kind) :: time_type_to_real
873type(time_type), intent(in) :: time
874
875if(.not.module_is_initialized) call time_manager_init
876
877time_type_to_real = real(time%days, r8_kind) * 86400.e0_r8_kind + real(time%seconds, r8_kind) + &
878 real(time%ticks, r8_kind)/real(ticks_per_second, r8_kind)
879
880end function time_type_to_real
881
882!> @brief Convert a real number of seconds into a time_type variable.
883!! @return A filled time type variable, and an error message if an
884!! error occurs.
885function real8_to_time_type(x,err_msg) result(t)
886 real(r8_kind),intent(in) :: x !< Number of seconds.
887 character(len=*),intent(out),optional :: err_msg !< Error message.
888 type(time_type) :: t
889 integer :: days
890 integer :: seconds
891 integer :: ticks
892 character(len=128) :: err_msg_local
893 real(r8_kind),parameter :: spd = 86400.0_r8_kind
894 real(r8_kind) :: tps
895 real(r8_kind) :: a
896 tps = real(ticks_per_second, r8_kind)
897 a = x/spd
898 days = safe_rtoi(a,do_floor)
899 a = x - real(days, r8_kind)*spd
900 seconds = safe_rtoi(a,do_floor)
901 a = (a - real(seconds, r8_kind))*tps
902 ticks = safe_rtoi(a,do_nearest)
903 if (.not. set_time_private(seconds,days,ticks,t,err_msg_local)) then
904 if (error_handler('function real8_to_time_type',err_msg_local,err_msg)) then
905 return
906 endif
907 endif
908end function real8_to_time_type
909
910function real4_to_time_type(x, err_msg) result(t)
911 real(r4_kind), intent(in) :: x !< number of seconds
912 character(len=*),intent(out),optional :: err_msg !< Error message.
913 type(time_type) :: t
914 t = real_to_time_type(real(x, r8_kind), err_msg)
915end function
916
917!> @brief Convert a floating point value to an integer value.
918!! @return The integer value, using the input rounding mode.
919function safe_rtoi(rval,mode) result(ival)
920 real(r8_kind),intent(in) :: rval !< A floating point value.
921 integer,intent(in) :: mode !< A rouding mode (either "do_floor" or
922 !! "do_nearest")
923 integer :: ival
924 real(r8_kind) :: big
925 big = real(huge(ival), r8_kind)
926 if (rval .le. big .and. -1.0_r8_kind*rval .ge. -1.0_r8_kind*big) then
927 if (mode .eq. do_floor) then
928 ival = floor(rval)
929 elseif (mode .eq. do_nearest) then
930 ival = nint(rval)
931 else
932 call error_mesg("safe_rtoi","mode must be either do_floor" &
933 //" or do_nearest.",fatal)
934 endif
935 else
936 call error_mesg("safe_rtoi","input value cannot be safely" &
937 //" converted to a 32-bit integer.",fatal)
938 endif
939end function safe_rtoi
940
941!-------------------------------------------------------------------------
942
943!> Returns the largest time, t, for which n * t <= time.
944function time_scalar_divide(time, n)
946type(time_type), intent(in) :: time !< time interval to divide
947integer, intent(in) :: n !< divisor
948real(r8_kind) :: d, div, dseconds_per_day, dticks_per_second
949integer :: days, seconds, ticks
950type(time_type) :: prod1, prod2
951character(len=128) tmp1,tmp2
952logical :: ltmp
953
954! Convert time interval to floating point days; risky for general performance?
955dseconds_per_day = real(seconds_per_day, r8_kind)
956dticks_per_second = real(ticks_per_second, r8_kind)
957d = real(time%days,r8_kind)*dseconds_per_day*dticks_per_second + real(time%seconds, r8_kind)*dticks_per_second + &
958 real(time%ticks, r8_kind)
959div = d/real(n, r8_kind)
960
961days = int(div/(dseconds_per_day*dticks_per_second))
962seconds = int(div/dticks_per_second - real(days, r8_kind)*dseconds_per_day)
963ticks = int(div - (real(days, r8_kind)*dseconds_per_day + real(seconds, r8_kind))*dticks_per_second)
964time_scalar_divide = set_time(seconds, days, ticks)
965
966! Need to make sure that roundoff isn't killing this
967prod1 = n * time_scalar_divide
968prod2 = n * (increment_time(time_scalar_divide, days=0, seconds=0, ticks=1))
969if(prod1 > time .or. prod2 <= time) then
970 call get_time(time, seconds, days, ticks)
971 write(tmp1,20) days,seconds,ticks
972 call get_time(time_scalar_divide, seconds, days, ticks)
973 write(tmp2,30) n,days,seconds,ticks
974 ltmp = error_handler('time_scalar_divide',' quotient error:'//trim(tmp1)//trim(tmp2))
975 20 format('time=',i7,' days, ',i6,' seconds, ',i6,' ticks')
976 30 format(' time divided by',i6,'=',i7,' days, ',i6,' seconds, ',i6,' ticks')
977endif
978
979end function time_scalar_divide
980
981!> Supports a commonly used type of test on times for models. Given the
982!! current time, and a time for an alarm, determines if this is the closest
983!! time to the alarm time given a time step of time_interval. If this
984!! is the closest time (alarm - time <= time_interval/2), the function
985!! returns true and the alarm is incremented by the alarm_interval. Watch
986!! for problems if the new alarm time is less than time + time_interval
987!!
988!! This is a specialized operation that is frequently performed in models.
989!! Given a time, and a time interval, this function is true if this is the
990!! closest time step to the alarm time. The actual computation is:
991!!
992!! if((alarm_time - time) &#60;&#61; (time_interval / 2))
993!!
994!! If the function is true, the alarm time is incremented by the
995!! alarm_interval; WARNING, this is a featured side effect. Otherwise, the
996!! function is false and there are no other effects. CAUTION: if the
997!! alarm_interval is smaller than the time_interval, the alarm may fail to
998!! return true ever again. Watch
999!! for problems if the new alarm time is less than time + time_interval
1000function interval_alarm(time, time_interval, alarm, alarm_interval)
1001logical :: interval_alarm
1002type(time_type), intent(in) :: time !< current time
1003type(time_type), intent(in) :: time_interval !< a time interval
1004type(time_type), intent(in) :: alarm_interval !< a time interval
1005type(time_type), intent(inout) :: alarm !< An alarm time, which is incremented by the alarm_interval if
1006 !! the function is true.
1007
1008if((alarm - time) <= (time_interval / 2)) then
1009 interval_alarm = .true.
1010 alarm = alarm + alarm_interval
1011else
1012 interval_alarm = .false.
1013end if
1014
1015end function interval_alarm
1016
1017!--------------------------------------------------------------------------
1018
1019!> Repeat_alarm supports an alarm that goes off with alarm_frequency and
1020!! lasts for alarm_length. If the nearest occurence of an alarm time
1021!! is less than half an alarm_length from the input time, repeat_alarm
1022!! is true. For instance, if the alarm_frequency is 1 day, and the
1023!! alarm_length is 2 hours, then repeat_alarm is true from time 2300 on
1024!! day n to time 0100 on day n + 1 for all n.
1025function repeat_alarm(time, alarm_frequency, alarm_length)
1026logical :: repeat_alarm
1027type(time_type), intent(in) :: time !< current time
1028type(time_type), intent(in) :: alarm_frequency !< a time interval for time in between alarm activations
1029type(time_type), intent(in) :: alarm_length !< a time interval for amount of time alarm is active for
1030type(time_type) :: prev, next
1031
1032prev = (time / alarm_frequency) * alarm_frequency
1033next = prev + alarm_frequency
1034if(time - prev <= alarm_length / 2 .or. next - time <= alarm_length / 2) then
1035 repeat_alarm = .true.
1036else
1037 repeat_alarm = .false.
1038endif
1039
1040end function repeat_alarm
1041
1042!--------------------------------------------------------------------------
1043
1044!=========================================================================
1045! CALENDAR OPERATIONS BEGIN HERE
1046!=========================================================================
1047
1048! <SUBROUTINE NAME="set_calendar_type">
1049
1050! <OVERVIEW>
1051! Sets the default calendar type for mapping time intervals to dates.
1052! </OVERVIEW>
1053! <DESCRIPTION>
1054! A constant number for setting the calendar type.
1055! </DESCRIPTION>
1056! <TEMPLATE> set_calendar_type(type, err_msg) </TEMPLATE>
1057
1058! <IN NAME="type" TYPE="integer" DIM="(scalar)" DEFAULT="NO_CALENDAR">
1059! A constant number for setting the calendar type.
1060! </IN>
1061! <OUT NAME="err_msg" TYPE="character, optional" DIM="(scalar)">
1062! When present, and when non-blank, a fatal error condition as been detected.
1063! The string itself is an error message.
1064! It is recommended that, when err_msg is present in the call
1065! to this routine, the next line of code should be something
1066! similar to this:
1067! if(err_msg /= '') call error_mesg('my_routine','additional info: '//trim(err_msg),FATAL)
1068! </OUT>
1069
1070!> @brief Sets calendar_type for mapping an interval to a date.
1071!! For the Gregorian calendar, negative years and the proleptic calendar are not used;
1072!! and the discontinuity of days in October 1582 (when the Gregorian calendar was adopted by select groups in Europe)
1073!! is also not taken into account.
1074subroutine set_calendar_type(type, err_msg)
1075
1076! Selects calendar for default mapping from time to date.
1077
1078integer, intent(in) :: type !< constant parameter value (ie. one NO_CALENDAR, )
1079character(len=*), intent(out), optional :: err_msg
1080character(len=256) :: err_msg_local
1081
1082if(.not.module_is_initialized) call time_manager_init()
1083
1084if(present(err_msg)) err_msg = ''
1085
1086if(type < 0 .or. type > max_type) then
1087 err_msg_local = 'Illegal calendar type'
1088 if(error_handler('subroutine set_calendar_type', err_msg_local, err_msg)) return
1089endif
1090
1091if(seconds_per_day /= 86400 .and. type /= no_calendar ) then
1092 err_msg_local = 'Only calendar type NO_CALENDAR is allowed when seconds_per_day is not 86400.'// &
1093 ' You are using '//trim(valid_calendar_types(type))//' and seconds_per_day='
1094 write(err_msg_local(len_trim(err_msg_local)+1:len_trim(err_msg_local)+8),'(i8)') seconds_per_day
1095 if(error_handler('subroutine set_calendar_type', err_msg_local, err_msg)) return
1096endif
1097
1098calendar_type = type
1099
1100end subroutine set_calendar_type
1101
1102!------------------------------------------------------------------------
1103
1104!> Returns default calendar type for mapping from time to date.
1106
1107integer :: get_calendar_type
1108
1109get_calendar_type = calendar_type
1110
1111end function get_calendar_type
1112
1113!------------------------------------------------------------------------
1114
1115!> Sets the number of ticks per second.
1117integer, intent(in) :: tps
1118
1119ticks_per_second = tps
1120
1121end subroutine set_ticks_per_second
1122
1123!------------------------------------------------------------------------
1124
1125!> Returns the number of ticks per second.
1127integer :: get_ticks_per_second
1128
1129get_ticks_per_second = ticks_per_second
1130
1131end function get_ticks_per_second
1132
1133!------------------------------------------------------------------------
1134
1135 !> @brief Gets the date for different calendar types.
1136 !! Given a time_interval, returns the corresponding date under
1137 !! the selected calendar.
1138 !! When err_msg present, and when non-blank, a fatal error condition as been detected.
1139 !! The string itself is an error message.
1140 !! It is recommended that, when err_msg is present in the call
1141 !! to this routine, the next line of code should be something
1142 !! similar to this: <br>
1143 !!
1144 !! if(err_msg /= '') call error_mesg('my_routine','additional info: '//trim(err_msg),FATAL)
1145 subroutine get_date(time, year, month, day, hour, minute, second, tick, err_msg)
1146
1147! Given a time, computes the corresponding date given the selected calendar
1148
1149 type(time_type), intent(in) :: time
1150 integer, intent(out) :: second, minute, hour, day, month, year
1151 integer, intent(out), optional :: tick
1152 character(len=*), intent(out), optional :: err_msg
1153 character(len=128) :: err_msg_local
1154 integer :: tick1
1155
1156 if(.not.module_is_initialized) call time_manager_init
1157 if(present(err_msg)) err_msg = ''
1158
1159 select case(calendar_type)
1160 case(thirty_day_months)
1161 call get_date_thirty(time, year, month, day, hour, minute, second, tick1)
1162 case(gregorian)
1163 call get_date_gregorian(time, year, month, day, hour, minute, second, tick1)
1164 case(julian)
1165 call get_date_julian_private(time, year, month, day, hour, minute, second, tick1)
1166 case(noleap)
1167 call get_date_no_leap_private(time, year, month, day, hour, minute, second, tick1)
1168 case(no_calendar)
1169 err_msg_local = 'Cannot produce a date when the calendar type is NO_CALENDAR'
1170 if(error_handler('subroutine get_date', err_msg_local, err_msg)) return
1171 case default
1172 err_msg_local = 'Invalid calendar type'
1173 if(error_handler('subroutine get_date', err_msg_local, err_msg)) return
1174 end select
1175
1176 if(present(tick)) then
1177 tick = tick1
1178 else
1179 if(tick1 /= 0) then
1180 err_msg_local = 'tick must be present when time has a second fraction'
1181 if(error_handler('subroutine get_date', err_msg_local, err_msg)) return
1182 endif
1183 endif
1184
1185 end subroutine get_date
1186
1187!------------------------------------------------------------------------
1188
1189!> @brief Gets the date on a Gregorian calendar.
1190!! Computes the year, month, day on the fly from the quantity time%days
1191 subroutine get_date_gregorian(time, year, month, day, hour, minute, second, tick)
1192
1193 type(time_type), intent(in) :: time
1194 integer, intent(out) :: year, month, day, hour, minute, second
1195 integer, intent(out) :: tick
1196 integer :: iday, isec
1197
1198 integer :: l !< value of 1 if leap year; value of 0 otherwise
1199 integer :: ncenturies !< number of centuries passed in the current 400 year period
1200 integer :: nlpyrs !< number of leap years in the current century
1201 integer :: yearx, monthx, dayx, idayx !< temporary values for year, month, day
1202 integer :: i !< counter, dummy variable
1203
1204 ! Computes date corresponding to time for gregorian calendar
1205
1206 !Carried over from the old subroutine
1207 if(time%seconds >= 86400) then ! This check appears to be unecessary.
1208 call error_mesg('get_date',.ge.'Time%seconds 86400 in subroutine get_date_gregorian',fatal)
1209 endif
1210
1211 iday = mod(time%days+1,days_in_400_year_period)
1212
1213 yearx = 1
1214 idayx = 0
1215 if( iday.eq.0 ) then !year 400
1216 yearx = 0
1217 idayx = -366
1218 else if( iday.gt.365 ) then
1219 yearx = int(iday/365) - 1 !approximation off by -1 year at most
1220 ncenturies = int(yearx/100)
1221 nlpyrs = int((yearx-ncenturies*100)/4)
1222 idayx = ncenturies*36524 + (yearx-ncenturies*100)*365 + nlpyrs !36524 days in a century
1223 if( ncenturies.eq.4 ) idayx = idayx + 1 !year 400 is a leap year
1224 l = 0 ; if ( leap_year_gregorian_int(yearx+1) ) l = 1
1225 if ( (iday-idayx).gt.365+l ) then !off by -1 year
1226 yearx = yearx + 1
1227 idayx = idayx + 365 + l
1228 end if
1229 yearx = yearx + 1
1230 end if
1231
1232 year = 400*int((time%days+1)/days_in_400_year_period) + yearx
1233
1234 l = 0 ; if( leap_year_gregorian_int(year) ) l = 1
1235 dayx = iday - idayx
1236 if( dayx.le.31 ) then
1237 month = 1
1238 day = dayx
1239 else
1240 monthx = int(dayx/30)
1241 if( l.eq.1 ) then
1242 do i=1, monthx
1243 dayx = dayx - days_per_month(i)
1244 if(i.eq.2) dayx = dayx - l
1245 end do
1246 else
1247 do i=1, monthx
1248 dayx = dayx - days_per_month(i)
1249 end do
1250 end if
1251 month = monthx + 1
1252 day = dayx
1253 if( dayx.le.0 ) then
1254 month = monthx
1255 day = dayx + days_per_month(monthx)
1256 if(monthx.eq.2) day = day + l
1257 end if
1258 end if
1259
1260 hour = time%seconds / 3600
1261 isec = time%seconds - 3600*hour
1262 minute = isec / 60
1263 second = isec - 60*minute
1264 tick = time%ticks
1265
1266 end subroutine get_date_gregorian
1267!------------------------------------------------------------------------
1268
1269 function cut0(string)
1270 character(len=256) :: cut0
1271 character(len=*), intent(in) :: string
1272 integer :: i
1273
1274 cut0 = string
1275
1276 do i=1,len(string)
1277 if(ichar(string(i:i)) == 0 ) then
1278 cut0(i:i) = ' '
1279 endif
1280 enddo
1281
1282 return
1283 end function cut0
1284!------------------------------------------------------------------------
1285
1286!> Base date for Julian calendar is year 1 with all multiples of 4
1287!! years being leap years.
1288 subroutine get_date_julian_private(time, year, month, day, hour, minute, second, tick)
1289
1290
1291 type(time_type), intent(in) :: time
1292 integer, intent(out) :: second, minute, hour, day, month, year
1293 integer, intent(out) :: tick
1294 integer :: m, t, nfour, nex, days_this_month
1295 logical :: leap
1296
1297! find number of four year periods; also get modulo number of days
1298 nfour = time%days / (4 * 365 + 1)
1299 day = modulo(time%days, (4 * 365 + 1))
1300
1301! Find out what year in four year chunk
1302 nex = day / 365
1303 if(nex == 4) then
1304 nex = 3
1305 day = 366
1306 else
1307 day=modulo(day, 365) + 1
1308 endif
1309
1310! Is this a leap year?
1311 leap = (nex == 3)
1312
1313 year = 1 + 4 * nfour + nex
1314
1315! find month and day
1316 do m = 1, 12
1317 month = m
1318 days_this_month = days_per_month(m)
1319 if(leap .and. m == 2) days_this_month = 29
1320 if(day <= days_this_month) exit
1321 day = day - days_this_month
1322 end do
1323
1324! find hour,minute and second
1325 t = time%seconds
1326 hour = t / (60 * 60)
1327 t = t - hour * (60 * 60)
1328 minute = t / 60
1329 second = t - 60 * minute
1330 tick = time%ticks
1331 end subroutine get_date_julian_private
1332
1333!------------------------------------------------------------------------
1334 subroutine get_date_julian(time, year, month, day, hour, minute, second)
1335
1336! No need to include tick in argument list because this routine
1337! exists only for interpolator.F90, which does not need it.
1338
1339 type(time_type), intent(in) :: time
1340 integer, intent(out) :: second, minute, hour, day, month, year
1341 integer :: tick
1342
1343 call get_date_julian_private(time, year, month, day, hour, minute, second, tick)
1344
1345 end subroutine get_date_julian
1346
1347!------------------------------------------------------------------------
1348
1349 !> Computes date corresponding to time interval for 30 day months, 12
1350 !! month years.
1351 subroutine get_date_thirty(time, year, month, day, hour, minute, second, tick)
1352
1353
1354 type(time_type), intent(in) :: time
1355 integer, intent(out) :: second, minute, hour, day, month, year
1356 integer, intent(out) :: tick
1357 integer :: t, dmonth, dyear
1358
1359 t = time%days
1360 dyear = t / (30 * 12)
1361 year = dyear + 1
1362 t = t - dyear * (30 * 12)
1363 dmonth = t / 30
1364 month = 1 + dmonth
1365 day = t -dmonth * 30 + 1
1366
1367 t = time%seconds
1368 hour = t / (60 * 60)
1369 t = t - hour * (60 * 60)
1370 minute = t / 60
1371 second = t - 60 * minute
1372 tick = time%ticks
1373
1374 end subroutine get_date_thirty
1375!------------------------------------------------------------------------
1376
1377 subroutine get_date_no_leap_private(time, year, month, day, hour, minute, second, tick)
1378
1379! Base date for NOLEAP calendar is year 1.
1380
1381 type(time_type), intent(in) :: time
1382 integer, intent(out) :: second, minute, hour, day, month, year
1383 integer, intent(out) :: tick
1384 integer :: m, t
1385
1386! get modulo number of days
1387 year = time%days / 365 + 1
1388 day = modulo(time%days, 365) + 1
1389
1390! find month and day
1391 do m = 1, 12
1392 month = m
1393 if(day <= days_per_month(m)) exit
1394 day = day - days_per_month(m)
1395 end do
1396
1397! find hour,minute and second
1398 t = time%seconds
1399 hour = t / (60 * 60)
1400 t = t - hour * (60 * 60)
1401 minute = t / 60
1402 second = t - 60 * minute
1403 tick = time%ticks
1404
1405 end subroutine get_date_no_leap_private
1406
1407!------------------------------------------------------------------------
1408 subroutine get_date_no_leap(time, year, month, day, hour, minute, second)
1409
1410! No need to include tick in argument list because this routine
1411! exists only for interpolator.F90, which does not need it.
1412
1413 type(time_type), intent(in) :: time
1414 integer, intent(out) :: second, minute, hour, day, month, year
1415 integer :: tick
1416
1417 call get_date_no_leap_private(time, year, month, day, hour, minute, second, tick)
1418
1419 end subroutine get_date_no_leap
1420!------------------------------------------------------------------------
1421
1422 !> @brief Sets days for different calendar types.
1423 !! Given an input date in year, month, days, etc., creates a
1424 !! time_type that represents this time interval from the
1425 !! internally defined base date.
1426 !!
1427 !! @note that it is possible to specify
1428 !! any number of illegal dates; these are checked for and generate
1429 !! errors as appropriate.
1430 function set_date_private(year, month, day, hour, minute, second, tick, Time_out, err_msg)
1431
1432! Given a date, computes the corresponding time given the selected
1433! date time mapping algorithm.
1434 logical :: set_date_private
1435 integer, intent(in) :: year, month, day, hour, minute, second, tick
1436 type(time_type) :: time_out
1437 character(len=*), intent(out) :: err_msg !< error message, if non-empty an error has occured
1438
1439 if(.not.module_is_initialized) call time_manager_init
1440
1441 err_msg = ''
1442
1443 select case(calendar_type)
1444 case(thirty_day_months)
1445 set_date_private = set_date_thirty(year, month, day, hour, minute, second, tick, time_out, err_msg)
1446 case(gregorian)
1447 set_date_private = set_date_gregorian(year, month, day, hour, minute, second, tick, time_out, err_msg)
1448 case(julian)
1449 set_date_private = set_date_julian_private(year, month, day, hour, minute, second, tick, time_out, err_msg)
1450 case(noleap)
1451 set_date_private = set_date_no_leap_private(year, month, day, hour, minute, second, tick, time_out, err_msg)
1452 case (no_calendar)
1453 err_msg = 'Cannot produce a date when calendar type is NO_CALENDAR'
1454 set_date_private = .false.
1455 case default
1456 err_msg = 'Invalid calendar type'
1457 set_date_private = .false.
1458 end select
1459
1460 end function set_date_private
1461
1462!------------------------------------------------------------------------
1463
1464 !> @brief Calls set_date_private to set days for different calendar types.
1465 function set_date_i(year, month, day, hour, minute, second, tick, err_msg)
1466 type(time_type) :: set_date_i
1467 integer, intent(in) :: day, month, year
1468 integer, intent(in), optional :: second, minute, hour, tick
1469 character(len=*), intent(out), optional :: err_msg
1470 integer :: osecond, ominute, ohour, otick
1471 character(len=128) :: err_msg_local
1472
1473 if(.not.module_is_initialized) call time_manager_init
1474 if(present(err_msg)) err_msg = ''
1475
1476! Missing optionals are set to 0
1477 osecond = 0; if(present(second)) osecond = second
1478 ominute = 0; if(present(minute)) ominute = minute
1479 ohour = 0; if(present(hour)) ohour = hour
1480 otick = 0; if(present(tick)) otick = tick
1481
1482 if(.not.set_date_private(year, month, day, ohour, ominute, osecond, otick, set_date_i, err_msg_local)) then
1483 if(error_handler('function set_date_i', err_msg_local, err_msg)) return
1484 end if
1485
1486 end function set_date_i
1487!------------------------------------------------------------------------
1488
1489 !> @brief Calls set_date_private for different calendar types when given a string input.
1490 !! Examples of acceptable forms of string:
1491 !!
1492 !! 1980-01-01 00:00:00
1493 !! 1980-01-01 00:00:00.50
1494 !! 1980-1-1 0:0:0
1495 !! 1980-1-1
1496 !!
1497 !! year number must occupy 4 spaces.
1498 !! months, days, hours, minutes, seconds may occupy 1 or 2 spaces
1499 !! year, month and day must be separated by a '-'
1500 !! hour, minute, second must be separated by a ':'
1501 !! hour, minute, second are optional. If not present then zero is assumed.
1502 !! second may be a real number.
1503 !!
1504 !! zero_year_warning:
1505 !! If the year number is zero, it will be silently changed to one,
1506 !! unless zero_year_warning=.true., in which case a WARNING message
1507 !! will also be issued
1508 function set_date_c(string, zero_year_warning, err_msg, allow_rounding)
1509 type(time_type) :: set_date_c
1510 character(len=*), intent(in) :: string
1511 logical, intent(in), optional :: zero_year_warning
1512 character(len=*), intent(out), optional :: err_msg
1513 logical, intent(in), optional :: allow_rounding
1514 character(len=4) :: formt='(i )'
1515 logical :: correct_form, zero_year_warning_local, allow_rounding_local
1516 integer :: i1, i2, i3, i4, i5, i6, i7
1517 character(len=32) :: string_sifted_left
1518 integer :: year, month, day, hour, minute, second, tick
1519 character(len=128) :: err_msg_local
1520
1521 if(.not.module_is_initialized) call time_manager_init()
1522 if(present(err_msg)) err_msg = ''
1523 if(present(zero_year_warning)) then
1524 zero_year_warning_local = zero_year_warning
1525 else
1526 zero_year_warning_local = .true.
1527 endif
1528 if(present(allow_rounding)) then
1529 allow_rounding_local = allow_rounding
1530 else
1531 allow_rounding_local = .true.
1532 endif
1533
1534 string_sifted_left = adjustl(string)
1535 i1 = index(string_sifted_left,'-')
1536 i2 = index(string_sifted_left,'-',back=.true.)
1537 i3 = index(string_sifted_left,':')
1538 i4 = index(string_sifted_left,':',back=.true.)
1539 i5 = len_trim(cut0(string_sifted_left))
1540 i6 = index(string_sifted_left,'.',back=.true.)
1541 correct_form = (i1 > 1) ! year number must occupy at least 1 space
1542 correct_form = correct_form .and. (i2-i1 == 2 .or. i2-i1 == 3) ! month number must occupy 1 or 2 spaces
1543 if(.not.correct_form) then
1544 err_msg_local = 'Form of character time stamp is incorrect. The character time stamp is: '//trim(string)
1545 if(error_handler('function set_date_c', err_msg_local, err_msg)) return
1546 endif
1547 write(formt(3:3),'(i1)') i1-1
1548 read(string_sifted_left(1:i1-1),formt) year
1549 if(year == 0) then
1550 year = 1
1551 if(zero_year_warning_local) then
1552 call error_mesg('set_date_c','Year zero is invalid. Resetting year to 1', warning)
1553 endif
1554 endif
1555 write(formt(3:3),'(i1)') i2-i1-1
1556 read(string_sifted_left(i1+1:i2-1),formt) month
1557 i7 = min(i2+2,i5)
1558 read(string_sifted_left(i2+1:i7),'(i2)') day
1559
1560 if(i3 == 0) then
1561! There are no minutes or seconds in the string
1562 minute = 0
1563 second = 0
1564 tick = 0
1565 if(i5 <= i2+2) then
1566 ! There is no clocktime in the string at all
1567 hour = 0
1568 else
1569 ! The clocktime includes only hours
1570 read(string_sifted_left(i5-1:i5),'(i2)') hour
1571 endif
1572 else if(i3 == i4) then
1573 ! The string includes hours and minutes, but no seconds
1574 read(string_sifted_left(i3-2:i3-1),'(i2)') hour
1575 write(formt(3:3),'(i1)') i5-i3
1576 read(string_sifted_left(i3+1:i5),formt) minute
1577 second = 0
1578 tick = 0
1579 else
1580 ! The string includes hours, minutes, and seconds
1581 read(string_sifted_left(i3-2:i3-1),'(i2)') hour
1582 write(formt(3:3),'(i1)') i4-i3-1
1583 read(string_sifted_left(i3+1:i4-1),formt) minute
1584 write(formt(3:3),'(i1)') i5-i4
1585 if(i6 == 0) then
1586 ! There are no fractional seconds
1587 read(string_sifted_left(i4+1:i5),formt) second
1588 tick = 0
1589 else
1590 read(string_sifted_left(i4+1:i6-1),formt) second
1591 if(.not.get_tick_from_string(string_sifted_left(i6:i5), err_msg_local, allow_rounding_local, tick)) then
1592 if(error_handler('function set_date_c', err_msg_local, err_msg)) return
1593 endif
1594 ! If tick has been rounded up to ticks_per_second, then bump up second.
1595 if(tick == ticks_per_second) then
1596 second = second + 1
1597 tick = 0
1598 endif
1599 endif
1600 endif
1601
1602 if(.not.set_date_private(year, month, day, hour, minute, second, tick, set_date_c, err_msg_local)) then
1603 if(error_handler('function set_date_c', err_msg_local, err_msg)) return
1604 end if
1605
1606 end function set_date_c
1607
1608!------------------------------------------------------------------------
1609
1610!> @brief Sets Time_out%days on a Gregorian calendar
1611!! Computes the total number of days between 1/1/0001 to the current month/day/year
1612 function set_date_gregorian(year, month, day, hour, minute, second, tick, Time_out, err_msg)
1613 logical :: set_date_gregorian
1614
1615! Computes time corresponding to date for gregorian calendar.
1616
1617 integer, intent(in) :: year, month, day, hour, minute, second, tick
1618 type(time_type), intent(out) :: time_out
1619 character(len=*), intent(out) :: err_msg
1620 integer :: yearx, monthx, dayx, hrx, minx, secx, tickx, ncenturies, nlpyrs, l
1621
1622 if( .not.valid_increments(year,month,day,hour,minute,second,tick,err_msg) ) then
1623 set_date_gregorian = .false.
1624 return
1625 endif
1626
1627 l = 0 ; if( leap_year_gregorian_int(year) ) l = 1
1628
1629 ! Check if date is invalid
1630 if(month.eq.2) then
1631 if(day.gt.days_per_month(month)+l .or. day.lt.1) then
1632 err_msg = 'Invalid_date. Date='//convert_integer_date_to_char(year,month,day,hour,minute,second)
1633 set_date_gregorian = .false.
1634 return
1635 end if
1636 else
1637 if(day.gt.days_per_month(month) .or. day.lt.1) then
1638 err_msg = 'Invalid_date. Date='//convert_integer_date_to_char(year,month,day,hour,minute,second)
1639 set_date_gregorian = .false.
1640 return
1641 end if
1642 end if
1643
1644 time_out%seconds = second + 60*(minute + 60*hour)
1645
1646 yearx = mod(year-1,400)
1647 dayx = 0
1648 if(yearx.gt.0) then
1649 ncenturies = int( yearx/100 )
1650 nlpyrs = int( (yearx-ncenturies*100)/4 )
1651 dayx = ncenturies*36524 + (yearx-ncenturies*100)*365 + nlpyrs ! 36524 days in 100 years, year 100 not
1652 !! leap year
1653 if(ncenturies.eq.4) dayx = dayx + 1 ! year 400 is a leap year
1654 end if
1655
1656 select case( month )
1657 case(1) ; dayx = dayx
1658 case(2) ; dayx = dayx + 31
1659 case(3) ; dayx = dayx + 59 + l
1660 case(4) ; dayx = dayx + 90 + l
1661 case(5) ; dayx = dayx + 120 + l
1662 case(6) ; dayx = dayx + 151 + l
1663 case(7) ; dayx = dayx + 181 + l
1664 case(8) ; dayx = dayx + 212 + l
1665 case(9) ; dayx = dayx + 243 + l
1666 case(10) ; dayx = dayx + 273 + l
1667 case(11) ; dayx = dayx + 304 + l
1668 case(12) ; dayx = dayx + 334 + l
1669 end select
1670
1671 dayx = int((year-1)/400)*days_in_400_year_period + dayx + day - 1
1672 time_out%days = dayx
1673 time_out%ticks = tick
1674
1675 err_msg = ''
1676 set_date_gregorian = .true.
1677
1678 ! check
1679 yearx = year ; monthx = month ; dayx = day
1680 hrx = hour ; minx = minute ; secx = second ; tickx = tick
1681 call get_date_gregorian(time_out, yearx, monthx, dayx, hrx, minx, secx, tickx)
1682 l = 0 ; if( leap_year_gregorian_int(yearx) ) l = 1
1683 if( monthx.lt.1 .or. monthx.gt.12 ) then
1684 err_msg = 'Invalid_date. Date='//convert_integer_date_to_char(yearx,monthx,dayx,hour,minute,second)
1685 set_date_gregorian = .false.
1686 else if( dayx.lt.1 .or. dayx.gt.days_per_month(monthx) ) then
1687 if( monthx.eq.2 .and. dayx.le.days_per_month(monthx)+l ) return
1688 err_msg = 'Invalid_date. Date='//convert_integer_date_to_char(yearx,monthx,dayx,hour,minute,second)
1689 set_date_gregorian = .false.
1690 end if
1691
1692 end function set_date_gregorian
1693
1694!------------------------------------------------------------------------
1695
1696 function set_date_julian_private(year, month, day, hour, minute, second, tick, Time_out, err_msg)
1697 logical :: set_date_julian_private
1698
1699! Returns time corresponding to date for julian calendar.
1700
1701 integer, intent(in) :: year, month, day, hour, minute, second, tick
1702 type(time_type), intent(out) :: time_out
1703 character(len=*), intent(out) :: err_msg
1704 integer :: ndays, m, nleapyr
1705 logical :: leap
1706
1707 if( .not.valid_increments(year,month,day,hour,minute,second,tick,err_msg) ) then
1708 set_date_julian_private = .false.
1709 return
1710 endif
1711
1712 if(month /= 2 .and. day > days_per_month(month)) then
1713 err_msg = 'Invalid date. Date='//convert_integer_date_to_char(year,month,day,hour,minute,second)
1714 set_date_julian_private = .false.
1715 return
1716 endif
1717
1718! Is this a leap year?
1719 leap = (modulo(year,4) == 0)
1720! compute number of complete leap years from year 1
1721 nleapyr = (year - 1) / 4
1722
1723! Finish checking for day specication errors
1724 if(month == 2 .and. (day > 29 .or. ((.not. leap) .and. day > 28))) then
1725 err_msg = 'Invalid date. Date='//convert_integer_date_to_char(year,month,day,hour,minute,second)
1726 set_date_julian_private = .false.
1727 return
1728 endif
1729
1730 ndays = 0
1731 do m = 1, month - 1
1732 ndays = ndays + days_per_month(m)
1733 if(leap .and. m == 2) ndays = ndays + 1
1734 enddo
1735
1736 time_out%seconds = second + 60 * (minute + 60 * hour)
1737 time_out%days = day -1 + ndays + 365*(year - nleapyr - 1) + 366*(nleapyr)
1738 time_out%ticks = tick
1739 err_msg = ''
1740 set_date_julian_private = .true.
1741
1742 end function set_date_julian_private
1743
1744!------------------------------------------------------------------------
1745 function set_date_julian(year, month, day, hour, minute, second)
1746
1747! No need to include tick or err_msg in argument list because this
1748! routine exists only for interpolator.F90, which does not need them.
1749
1750 type(time_type) :: set_date_julian
1751 integer, intent(in) :: year, month, day, hour, minute, second
1752 character(len=128) :: err_msg
1753
1754 if(.not.set_date_julian_private(year, month, day, hour, minute, second, 0, set_date_julian, err_msg)) then
1755 call error_mesg('set_date_julian',trim(err_msg),fatal)
1756 endif
1757
1758 end function set_date_julian
1759!------------------------------------------------------------------------
1760
1761 function set_date_thirty(year, month, day, hour, minute, second, tick, Time_out, err_msg)
1762 logical :: set_date_thirty
1763
1764! Computes time corresponding to date for thirty day months.
1765
1766 integer, intent(in) :: year, month, day, hour, minute, second, tick
1767 type(time_type), intent(out) :: time_out
1768 character(len=*), intent(out) :: err_msg
1769
1770 if( .not.valid_increments(year,month,day,hour,minute,second,tick,err_msg) ) then
1771 set_date_thirty = .false.
1772 return
1773 endif
1774
1775 if(day > 30) then
1776 err_msg = 'Invalid date. Date='//convert_integer_date_to_char(year,month,day,hour,minute,second)
1777 set_date_thirty = .false.
1778 return
1779 endif
1780
1781 time_out%days = (day - 1) + 30 * ((month - 1) + 12 * (year - 1))
1782 time_out%seconds = second + 60 * (minute + 60 * hour)
1783 time_out%ticks = tick
1784 err_msg = ''
1785 set_date_thirty = .true.
1786
1787 end function set_date_thirty
1788
1789!------------------------------------------------------------------------
1790
1791 function set_date_no_leap_private(year, month, day, hour, minute, second, tick, Time_out, err_msg)
1792 logical :: set_date_no_leap_private
1793
1794! Computes time corresponding to date for fixed 365 day year calendar.
1795
1796 integer, intent(in) :: year, month, day, hour, minute, second, tick
1797 type(time_type), intent(out) :: time_out
1798 character(len=*), intent(out) :: err_msg
1799 integer :: ndays, m
1800
1801 if( .not.valid_increments(year,month,day,hour,minute,second,tick,err_msg) ) then
1802 set_date_no_leap_private = .false.
1803 return
1804 endif
1805
1806 if(day > days_per_month(month)) then
1807 err_msg = 'Invalid date. Date='//convert_integer_date_to_char(year,month,day,hour,minute,second)
1808 set_date_no_leap_private = .false.
1809 return
1810 endif
1811
1812 ndays = 0
1813 do m = 1, month - 1
1814 ndays = ndays + days_per_month(m)
1815 enddo
1816
1817! No need for err_msg in call to set_time because previous checks ensure positive value of time.
1818 time_out = set_time(second + 60 * (minute + 60 * hour), day -1 + ndays + 365 * (year - 1), tick)
1819 err_msg = ''
1820 set_date_no_leap_private = .true.
1821
1822 end function set_date_no_leap_private
1823!------------------------------------------------------------------------
1824
1825 function set_date_no_leap(year, month, day, hour, minute, second)
1826
1827! No need to include tick or err_msg in argument list because this
1828! routine exists only for interpolator.F90, which does not need them.
1829
1830 type(time_type) :: set_date_no_leap
1831 integer, intent(in) :: year, month, day, hour, minute, second
1832 character(len=128) :: err_msg
1833
1834 if(.not.set_date_no_leap_private(year, month, day, hour, minute, second, 0, set_date_no_leap, err_msg)) then
1835 call error_mesg('set_date_no_leap',trim(err_msg),fatal)
1836 endif
1837
1838 end function set_date_no_leap
1839
1840!=========================================================================
1841
1842 function valid_increments(year, month, day, hour, minute, second, tick, err_msg)
1843 logical :: valid_increments
1844 integer, intent(in) :: year, month, day, hour, minute, second, tick
1845 character(len=128), intent(out) :: err_msg
1846
1847! Check for invalid values
1848
1849 err_msg = ''
1850 valid_increments = .true.
1851 if(second > 59 .or. second < 0 .or. minute > 59 .or. minute < 0 &
1852 .or. hour > 23 .or. hour < 0 .or. day > 31 .or. day < 1 &
1853 .or. month > 12 .or. month < 1 .or. year < 1) then
1854 err_msg = 'Invalid date. Date='//convert_integer_date_to_char(year,month,day,hour,minute,second)
1855 valid_increments = .false.
1856 return
1857 endif
1858 if(tick < 0 .or. tick >= ticks_per_second) then
1859 write(err_msg,'(a,i6)') 'Invalid number of ticks. tick=',tick
1860 valid_increments = .false.
1861 endif
1862
1863 end function valid_increments
1864
1865!=========================================================================
1866
1867 function convert_integer_date_to_char(year, month, day, hour, minute, second)
1868 character(len=19) :: convert_integer_date_to_char
1869 integer, intent(in) :: year, month, day
1870 integer, intent(in) :: hour, minute, second
1871
1872 write(convert_integer_date_to_char,10) year,month,day,hour,minute,second
1873 10 format(i4.4, '-', i2.2, '-', i2.2, ' ', i2.2, ':', i2.2, ':', i2.2)
1874
1875 end function convert_integer_date_to_char
1876
1877!=========================================================================
1878! END OF set_date BLOCK
1879!=========================================================================
1880
1881! <FUNCTION NAME="increment_date">
1882
1883! <OVERVIEW>
1884! Increments the date represented by a time interval and the
1885! default calendar type by a number of seconds, etc.
1886! </OVERVIEW>
1887! <DESCRIPTION>
1888! Given a time and some date increment, computes a new time. Depending
1889! on the mapping algorithm from date to time, it may be possible to specify
1890! undefined increments (i.e. if one increments by 68 days and 3 months in
1891! a Julian calendar, it matters which order these operations are done and
1892! we don't want to deal with stuff like that, make it an error).
1893! </DESCRIPTION>
1894! <TEMPLATE>
1895! increment_date(time, years, months, days, hours, minutes, seconds, ticks, err_msg)
1896! </TEMPLATE>
1897! <IN NAME="time" TYPE="time_type"> A time interval.</IN>
1898! <IN NAME="years" TYPE="integer">An increment of years.</IN>
1899! <IN NAME="months" TYPE="integer">An increment of months.</IN>
1900! <IN NAME="days" TYPE="integer">An increment of days.</IN>
1901! <IN NAME="hours" TYPE="integer">An increment of hours.</IN>
1902! <IN NAME="minutes" TYPE="integer">An increment of minutes.</IN>
1903! <IN NAME="seconds" TYPE="integer">An increment of seconds.</IN>
1904! <IN NAME="ticks" TYPE="integer">An increment of ticks.</IN>
1905! <OUT NAME="err_msg" TYPE="character, optional" DIM="(scalar)">
1906! When present, and when non-blank, a fatal error condition as been detected.
1907! The string itself is an error message.
1908! It is recommended that, when err_msg is present in the call
1909! to this routine, the next line of code should be something
1910! similar to this:
1911! if(err_msg /= '') call error_mesg('my_routine','additional info: '//trim(err_msg),FATAL)
1912! </OUT>
1913! <OUT NAME="increment_date" TYPE="time_type"> A new time based on the input
1914! time interval and the calendar type.
1915! </OUT>
1916! <IN NAME="allow_neg_inc" TYPE="logical, optional" DIM="(scalar)" DEFAULT=".true.">
1917! When .false., it is a fatal error if any of the input time increments are negative.
1918! This mimics the behavior of lima and earlier revisions.
1919! </IN>
1920! <NOTE>
1921! For all but the thirty_day_months calendar, increments to months
1922! and years must be made separately from other units because of the
1923! non-associative nature of addition.
1924! If the result is a negative time (i.e. date before the base date)
1925! it is considered a fatal error.
1926! </NOTE>
1927
1928 function increment_date(Time, years, months, days, hours, minutes, seconds, ticks, err_msg, allow_neg_inc)
1929
1930! Given a time and some date increment, computes a new time. Depending
1931! on the mapping algorithm from date to time, it may be possible to specify
1932! undefined increments (i.e. if one increments by 68 days and 3 months in
1933! a Julian calendar, it matters which order these operations are done and
1934! we don't want to deal with stuff like that, make it an error).
1935
1936! This routine operates in one of two modes.
1937! 1. days, hours, minutes, seconds, ticks are incremented, years and months must be zero or absent arguments.
1938! 2. years and/or months are incremented, other time increments must be zero or absent arguments.
1939
1940 type(time_type) :: increment_date
1941 type(time_type), intent(in) :: time
1942 integer, intent(in), optional :: years, months, days, hours, minutes, seconds, ticks
1943 character(len=*), intent(out), optional :: err_msg
1944 logical, intent(in), optional :: allow_neg_inc
1945
1946 integer :: oyears, omonths, odays, ohours, ominutes, oseconds, oticks
1947 character(len=128) :: err_msg_local
1948 logical :: allow_neg_inc_local
1949
1950 if(.not.module_is_initialized) call time_manager_init
1951 if(present(err_msg)) err_msg = ''
1952
1953! Missing optionals are set to 0
1954 oseconds = 0; if(present(seconds)) oseconds = seconds
1955 ominutes = 0; if(present(minutes)) ominutes = minutes
1956 ohours = 0; if(present(hours)) ohours = hours
1957 odays = 0; if(present(days)) odays = days
1958 omonths = 0; if(present(months)) omonths = months
1959 oyears = 0; if(present(years)) oyears = years
1960 oticks = 0; if(present(ticks)) oticks = ticks
1961 allow_neg_inc_local=.true.; if(present(allow_neg_inc)) allow_neg_inc_local=allow_neg_inc
1962
1963 if(.not.allow_neg_inc_local) then
1964 if(oyears < 0 .or. omonths < 0 .or. odays < 0 .or. ohours < 0 .or. ominutes < 0 .or. oseconds < 0 .or. &
1965 & oticks < 0) then
1966 write(err_msg_local,10) oyears, omonths, odays, ohours, ominutes, oseconds, oticks
1967 if(error_handler('function increment_time', err_msg_local, err_msg)) return
1968 endif
1969 endif
1970 10 format('One or more time increments are negative: '// &
1971 'years=',i6,' months=',i6,' days=',i6,' hours=',i6,' minutes=',i6,' seconds=',i6,' ticks=',i6)
1972
1973 if(.not.increment_date_private( &
1974 time, oyears, omonths, odays, ohours, ominutes, oseconds, oticks, increment_date, err_msg_local)) then
1975 if(error_handler('function increment_date', err_msg_local, err_msg)) return
1976 endif
1977
1978 end function increment_date
1979
1980!=======================================================================
1981
1982 !> Given a time and some date increment, computes a new time. Depending
1983 !! on the mapping algorithm from date to time, it may be possible to specify
1984 !! undefined increments (i.e. if one increments by 68 days and 3 months in
1985 !! a Julian calendar, it matters which order these operations are done and
1986 !! we don't want to deal with stuff like that, make it an error).
1987 !!
1988 !! This routine operates in one of two modes.
1989 !! 1. days, hours, minutes, seconds, ticks are incremented, years and months must be zero or absent arguments.
1990 !! 2. years and/or months are incremented, other time increments must be zero or absent arguments.
1991 !!
1992 !! Negative increments are always allowed in the private version of this routine.
1993 function increment_date_private(Time, years, months, days, hours, minutes, seconds, ticks, Time_out, err_msg)
1994
1995
1996 logical :: increment_date_private
1997 type(time_type), intent(in) :: time
1998 integer, intent(in) :: years, months, days, hours, minutes, seconds, ticks
1999 type(time_type), intent(out) :: time_out
2000 character(len=*), intent(out) :: err_msg
2001 integer :: cyear , cmonth , cday , chour , cminute , csecond , ctick
2002 logical :: mode_1, mode_2
2003
2004 err_msg = ''
2005 increment_date_private = .true.
2006
2007 mode_1 = days /= 0 .or. hours /= 0 .or. minutes /= 0 .or. seconds /= 0 .or. ticks /= 0
2008 mode_2 = years /= 0 .or. months /= 0
2009
2010 if(.not.mode_1 .and. .not.mode_2) then
2011 ! All time increments are zero
2012 time_out = time
2013 return
2014 endif
2015
2016 if(mode_1 .and. mode_2) then
2017 err_msg = 'years and/or months must not be incremented with other time units'
2018 increment_date_private = .false.
2019 return
2020 endif
2021
2022 if(mode_1) then
2023 csecond = seconds + 60 * (minutes + 60 * hours)
2024 increment_date_private = increment_time_private(time, csecond, days, ticks, time_out, err_msg)
2025 endif
2026
2027 if(mode_2) then
2028 ! Convert Time to a date
2029 select case(calendar_type)
2030 case(thirty_day_months)
2031 call get_date_thirty (time, cyear, cmonth, cday, chour, cminute, csecond, ctick)
2032 case(noleap)
2033 call get_date_no_leap_private (time, cyear, cmonth, cday, chour, cminute, csecond, ctick)
2034 case(julian)
2035 call get_date_julian_private (time, cyear, cmonth, cday, chour, cminute, csecond, ctick)
2036 case(gregorian)
2037 call get_date_gregorian(time, cyear, cmonth, cday, chour, cminute, csecond, ctick)
2038 case(no_calendar)
2039 err_msg = 'Cannot increment a date when the calendar type is NO_CALENDAR'
2040 increment_date_private = .false.
2041 return
2042 case default
2043 err_msg = 'Invalid calendar type'
2044 increment_date_private = .false.
2045 return
2046 end select
2047
2048 ! Add month increment
2049 cmonth = cmonth + months
2050
2051 ! Adjust year and month number when cmonth falls outside the range 1 to 12
2052 cyear = cyear + floor(real(cmonth-1,r8_kind)/12.0_r8_kind)
2053 cmonth = modulo((cmonth-1),12) + 1
2054
2055 ! Add year increment
2056 cyear = cyear + years
2057
2058 ! Convert this back into a time.
2059 select case(calendar_type)
2060 case(thirty_day_months)
2061 increment_date_private = set_date_thirty(cyear, cmonth, cday, chour, cminute, csecond, ctick, time_out, err_msg)
2062 case(noleap)
2063 increment_date_private = set_date_no_leap_private(cyear, cmonth, cday, chour, cminute, csecond, ctick, &
2064 & time_out, err_msg)
2065 case(julian)
2066 increment_date_private = set_date_julian_private(cyear, cmonth, cday, chour, cminute, csecond, ctick, &
2067 & time_out, err_msg)
2068 case(gregorian)
2069 increment_date_private = set_date_gregorian(cyear, cmonth, cday, chour, cminute, csecond, ctick, time_out, &
2070 & err_msg)
2071 end select
2072 endif ! if(mode_2)
2073
2074 end function increment_date_private
2075
2076!=========================================================================
2077
2078 !> Given a time and some date decrement, computes a new time. Depending
2079 !! on the mapping algorithm from date to time, it may be possible to specify
2080 !! undefined decrements (i.e. if one decrements by 68 days and 3 months in
2081 !! a Julian calendar, it matters which order these operations are done and
2082 !! we don't want to deal with stuff like that, make it an error).
2083 !!
2084 !! @note For all but the thirty_day_months calendar, decrements to months
2085 !! and years must be made separately from other units because of the
2086 !! non-associative nature of addition.
2087 !! If the result is a negative time (i.e. date before the base date)
2088 !! it is considered a fatal error.
2089 function decrement_date(Time, years, months, days, hours, minutes, seconds, ticks, err_msg, allow_neg_inc)
2090
2091 type(time_type) :: decrement_date !< Time after the given decrement is applied
2092 type(time_type), intent(in) :: time !< time interval to decrement
2093 integer, intent(in), optional :: seconds, minutes, hours, days, months, years, ticks !< amount of time to decrement by
2094 !! units should not exceed next largest unit (ie. 61 seconds
2095 !! should be 1 min 1 sec )
2096 character(len=*), intent(out), optional :: err_msg
2097 logical, intent(in), optional :: allow_neg_inc
2098
2099 integer :: oseconds, ominutes, ohours, odays, omonths, oyears, oticks
2100 character(len=128) :: err_msg_local
2101 logical :: allow_neg_inc_local
2102
2103 if(present(err_msg)) err_msg = ''
2104
2105 ! Missing optionals are set to 0
2106 oseconds = 0; if(present(seconds)) oseconds = seconds
2107 ominutes = 0; if(present(minutes)) ominutes = minutes
2108 ohours = 0; if(present(hours)) ohours = hours
2109 odays = 0; if(present(days)) odays = days
2110 omonths = 0; if(present(months)) omonths = months
2111 oyears = 0; if(present(years)) oyears = years
2112 oticks = 0; if(present(ticks)) oticks = ticks
2113 allow_neg_inc_local=.true.; if(present(allow_neg_inc)) allow_neg_inc_local=allow_neg_inc
2114
2115 if(.not.allow_neg_inc_local) then
2116 if(oyears < 0 .or. omonths < 0 .or. odays < 0 .or. ohours < 0 .or. ominutes < 0 .or. oseconds < 0 .or. &
2117 & oticks < 0) then
2118 write(err_msg_local,10) oyears, omonths, odays, ohours, ominutes, oseconds, oticks
2119 if(error_handler('function decrement_date', err_msg_local, err_msg)) return
2120 endif
2121 endif
2122 10 format('One or more time increments are negative: '// &
2123 'years=',i6,' months=',i6,' days=',i6,' hours=',i6,' minutes=',i6,' seconds=',i6,' ticks=',i6)
2124
2125 if(.not.increment_date_private( &
2126 time, -oyears, -omonths, -odays, -ohours, -ominutes, -oseconds, -oticks, decrement_date, err_msg_local)) then
2127 if(error_handler('function decrement_date', err_msg_local, err_msg)) return
2128 endif
2129
2130 end function decrement_date
2131
2132!--------------------------------------------------------------------------
2133
2134!> Given a time, computes the corresponding date given the selected
2135!! date time mapping algorithm
2136function days_in_month(Time, err_msg)
2137integer :: days_in_month !< number of days in month given the current selected calendar type
2138type(time_type), intent(in) :: time !< a time interval
2139character(len=*), intent(out), optional :: err_msg
2140
2141if(.not.module_is_initialized) call time_manager_init
2142if(present(err_msg)) err_msg = ''
2143
2144select case(calendar_type)
2145case(thirty_day_months)
2147case(gregorian)
2149case(julian)
2151case(noleap)
2153case(no_calendar)
2154 if(error_handler('function days_in_month', &
2155 'days_in_month makes no sense when the calendar type is NO_CALENDAR', err_msg)) return
2156case default
2157 if(error_handler('function days_in_month', 'Invalid calendar type', err_msg)) return
2158end select
2159end function days_in_month
2160
2161!--------------------------------------------------------------------------
2162
2163!> Returns the number of days in a gregorian month.
2165integer :: days_in_month_gregorian
2166type(time_type), intent(in) :: time
2167integer :: year, month, day, hour, minute, second, ticks
2168
2169call get_date_gregorian(time, year, month, day, hour, minute, second, ticks)
2170days_in_month_gregorian = days_per_month(month)
2171if(leap_year_gregorian_int(year) .and. month == 2) days_in_month_gregorian = 29
2172
2173end function days_in_month_gregorian
2174
2175!--------------------------------------------------------------------------
2176
2177!> Returns the number of days in a julian month.
2179integer :: days_in_month_julian
2180type(time_type), intent(in) :: time
2181integer :: year, month, day, hour, minute, second, ticks
2182
2183call get_date_julian_private(time, year, month, day, hour, minute, second, ticks)
2184days_in_month_julian = days_per_month(month)
2185if(leap_year_julian(time) .and. month == 2) days_in_month_julian = 29
2186
2187end function days_in_month_julian
2188
2189!--------------------------------------------------------------------------
2190
2191!> Returns the number of days in a thirty day month (needed for transparent
2192!! changes to calendar type).
2194integer :: days_in_month_thirty
2195type(time_type), intent(in) :: time
2196
2198
2199end function days_in_month_thirty
2200
2201!--------------------------------------------------------------------------
2202
2203!> Returns the number of days in a 365 day year month.
2205integer :: days_in_month_no_leap
2206type(time_type), intent(in) :: time
2207integer :: year, month, day, hour, minute, second, ticks
2208
2209call get_date_no_leap_private(time, year, month, day, hour, minute, second, ticks)
2210days_in_month_no_leap= days_per_month(month)
2211
2212end function days_in_month_no_leap
2213
2214!> Returns true if the year corresponding to the input time is
2215!! a leap year (for default calendar). Always returns false for THIRTY_DAY_MONTHS and NOLEAP.
2216function leap_year(Time, err_msg)
2217logical :: leap_year
2218type(time_type), intent(in) :: time !< a time interval to check if leap year
2219character(len=*), intent(out), optional :: err_msg
2220
2221if(.not.module_is_initialized) call time_manager_init
2222if(present(err_msg)) err_msg=''
2223
2224select case(calendar_type)
2225case(thirty_day_months)
2227case(gregorian)
2228 leap_year = leap_year_gregorian(time)
2229case(julian)
2231case(noleap)
2233case default
2234 if(error_handler('function leap_year', 'Invalid calendar type in leap_year', err_msg)) return
2235end select
2236end function leap_year
2237
2238!--------------------------------------------------------------------------
2239
2240function leap_year_gregorian(Time)
2241
2242! Is this a leap year for gregorian calendar?
2243
2244logical :: leap_year_gregorian
2245type(time_type), intent(in) :: time
2246integer :: seconds, minutes, hours, day, month, year
2247
2248call get_date(time, year, month, day, hours, minutes, seconds)
2249leap_year_gregorian = leap_year_gregorian_int(year)
2250
2251end function leap_year_gregorian
2252
2253!--------------------------------------------------------------------------
2254
2255function leap_year_gregorian_int(year)
2256logical :: leap_year_gregorian_int
2257integer, intent(in) :: year
2258
2259leap_year_gregorian_int = mod(year,4) == 0
2260leap_year_gregorian_int = leap_year_gregorian_int .and. .not.mod(year,100) == 0
2261leap_year_gregorian_int = leap_year_gregorian_int .or. mod(year,400) == 0
2262
2263end function leap_year_gregorian_int
2264
2265!--------------------------------------------------------------------------
2266
2267!> Returns the number of days in a julian month.
2268function leap_year_julian(Time)
2269logical :: leap_year_julian
2270type(time_type), intent(in) :: time
2271integer :: seconds, minutes, hours, day, month, year
2272
2273call get_date(time, year, month, day, hours, minutes, seconds)
2274leap_year_julian = ((year / 4 * 4) == year)
2275
2276end function leap_year_julian
2277
2278!--------------------------------------------------------------------------
2279
2280!> No leap years in thirty day months, included for transparency.
2281function leap_year_thirty(Time)
2282logical :: leap_year_thirty
2283type(time_type), intent(in) :: time
2284
2285leap_year_thirty = .false.
2286
2287end function leap_year_thirty
2288
2289!--------------------------------------------------------------------------
2290
2291!> Another tough one; no leap year returns false for leap year inquiry.
2292function leap_year_no_leap(Time)
2293logical :: leap_year_no_leap
2294type(time_type), intent(in) :: time
2295
2296leap_year_no_leap = .false.
2297
2298end function leap_year_no_leap
2299
2300!> @brief Returns the mean length of the year in the default calendar setting.
2301!!
2302!> There are no arguments in this function. It returns the mean
2303!! length of the year for the default calendar.
2305! What is the length of the year for the default calendar type
2306
2308
2309if(.not.module_is_initialized) call time_manager_init
2310
2311select case(calendar_type)
2312case(thirty_day_months)
2313 length_of_year = length_of_year_thirty()
2314case(gregorian)
2315 length_of_year = length_of_year_gregorian()
2316case(julian)
2317 length_of_year = length_of_year_julian()
2318case(noleap)
2319 length_of_year = length_of_year_no_leap()
2320case default
2321 call error_mesg('length_of_year','Invalid calendar type in length_of_year',fatal)
2322end select
2323end function length_of_year
2324
2325!--------------------------------------------------------------------------
2326
2327function length_of_year_thirty()
2328
2329type(time_type) :: length_of_year_thirty
2330
2331length_of_year_thirty = set_time(0, 360)
2332
2333end function length_of_year_thirty
2334
2335!---------------------------------------------------------------------------
2336
2337function length_of_year_gregorian()
2338
2339type(time_type) :: length_of_year_gregorian
2340
2341length_of_year_gregorian = set_time(20952, 365) !20952 = 86500 * (days_in_400_yrs/400. - (days_in_400_yrs/400))
2342
2343end function length_of_year_gregorian
2344
2345!--------------------------------------------------------------------------
2346
2347function length_of_year_julian()
2348
2349type(time_type) :: length_of_year_julian
2350
2351length_of_year_julian = set_time(21600, 365) !21600 = (24/4) * 60 * 60
2352
2353end function length_of_year_julian
2354
2355!--------------------------------------------------------------------------
2356
2357function length_of_year_no_leap()
2358
2359type(time_type) :: length_of_year_no_leap
2360
2361length_of_year_no_leap = set_time(0, 365)
2362
2363end function length_of_year_no_leap
2364
2365!--------------------------------------------------------------------------
2366
2367!> Returns number of day in year for given time. Jan 1st is day 1, not zero!
2368function day_of_year(time)
2369 integer :: day_of_year
2370 type(time_type), intent(in) :: time
2371
2372 integer :: second, minute, hour, day, month, year
2373 type(time_type) :: t
2374
2375 call get_date(time,year,month,day,hour,minute,second)
2376 t = time-set_date(year,1,1,0,0,0)
2377 day_of_year = t%days + 1
2378end
2379
2380!> @brief Returns the number of days in the calendar year corresponding to the date represented by
2381!! time for the default calendar.
2382!> @returns The number of days in this year for the default calendar type.
2383function days_in_year(Time)
2384
2385! What is the number of days in this year for the default calendar type
2386
2387integer :: days_in_year
2388type(time_type), intent(in) :: time !< A time interval
2389
2390if(.not.module_is_initialized) call time_manager_init
2391
2392select case(calendar_type)
2393case(thirty_day_months)
2394 days_in_year = days_in_year_thirty(time)
2395case(gregorian)
2396 days_in_year = days_in_year_gregorian(time)
2397case(julian)
2398 days_in_year = days_in_year_julian(time)
2399case(noleap)
2400 days_in_year = days_in_year_no_leap(time)
2401case default
2402 call error_mesg('days_in_year','Invalid calendar type in days_in_year',fatal)
2403end select
2404end function days_in_year
2405
2406!--------------------------------------------------------------------------
2407
2408function days_in_year_thirty(Time)
2409
2410integer :: days_in_year_thirty
2411type(time_type), intent(in) :: time
2412
2413days_in_year_thirty = 360
2414
2415end function days_in_year_thirty
2416
2417!---------------------------------------------------------------------------
2418
2419function days_in_year_gregorian(Time)
2420
2421integer :: days_in_year_gregorian
2422type(time_type), intent(in) :: time
2423
2424if(leap_year_gregorian(time)) then
2425 days_in_year_gregorian = 366
2426else
2427 days_in_year_gregorian = 365
2428endif
2429
2430end function days_in_year_gregorian
2431
2432!--------------------------------------------------------------------------
2433function days_in_year_julian(Time)
2434
2435integer :: days_in_year_julian
2436type(time_type), intent(in) :: time
2437
2438if(leap_year_julian(time)) then
2439 days_in_year_julian = 366
2440else
2441 days_in_year_julian = 365
2442endif
2443
2444end function days_in_year_julian
2445
2446!--------------------------------------------------------------------------
2447
2448function days_in_year_no_leap(Time)
2449
2450integer :: days_in_year_no_leap
2451type(time_type), intent(in) :: time
2452
2453days_in_year_no_leap = 365
2454
2455end function days_in_year_no_leap
2456
2457!--------------------------------------------------------------------------
2458
2459!> @brief Returns a character string containing the name of the
2460!! month corresponding to month number n.
2461!!
2462!> Definition is the same for all calendar types.
2463!! @returns The character string associated with a month. All calendars have 12 months and return
2464!! full month names, not abreviations.
2465function month_name(n)
2466
2467! Returns character string associated with a month, for now, all calendars
2468! have 12 months and will return standard names.
2469
2470character (len=9) :: month_name
2471integer, intent(in) :: n !< Month number
2472character (len = 9), dimension(12) :: months = (/'January ', 'February ', &
2473 'March ', 'April ', 'May ', 'June ', 'July ', &
2474 'August ', 'September', 'October ', 'November ', 'December '/)
2475
2476if(.not.module_is_initialized) call time_manager_init
2477
2478if(n < 1 .or. n > 12) call error_mesg('month_name','Illegal month index',fatal)
2479
2480month_name = months(n)
2481
2482end function month_name
2483
2484!==========================================================================
2485
2486!> The purpose of this routine is to prevent the addition of an excessive amount of code in order to implement
2487!! the error handling scheme involving an optional error flag of type character.
2488!! It allows one line of code to accomplish what would otherwise require 6 lines.
2489!! A value of .true. for this function is a flag to the caller that it should immediately return to it's caller.
2490 function error_handler(routine, err_msg_local, err_msg)
2491
2492
2493 logical :: error_handler
2494 character(len=*), intent(in) :: routine, err_msg_local
2495 character(len=*), intent(out), optional :: err_msg
2496
2497 error_handler = .false.
2498 if(present(err_msg)) then
2499 err_msg = err_msg_local
2500 error_handler = .true.
2501 else
2502 call error_mesg(trim(routine),trim(err_msg_local),fatal)
2503 endif
2504
2505 end function error_handler
2506
2507!------------------------------------------------------------------------
2508
2509!> Initialization routine. Writes the version information to the log file
2510subroutine time_manager_init ( )
2511
2512 if (module_is_initialized) return ! silent return if already called
2513
2514 call write_version_number("TIME_MANAGER_MOD", version)
2515 module_is_initialized = .true.
2516
2517end subroutine time_manager_init
2518
2519!------------------------------------------------------------------------
2520
2521!> @brief Prints the given time_type argument as a time (using days, seconds and ticks)
2522!!
2523!> @note There is no check for PE number.
2524subroutine print_time (Time,str,unit)
2525type(time_type) , intent(in) :: time !< Time that will be printed
2526character (len=*), intent(in), optional :: str !< Character string that precedes the printed time
2527integer , intent(in), optional :: unit !< Unit number for printed output, defaults to stdout
2528integer :: s,d,ticks, ns,nd,nt, unit_in
2529character(len=19) :: fmt
2530
2531! prints the time to standard output (or optional unit) as days and seconds
2532! NOTE: there is no check for PE number
2533
2534 unit_in = stdout()
2535 if (present(unit)) unit_in = unit
2536
2537 call get_time (time,s,d,ticks)
2538
2539! format output
2540! get number of digits for days and seconds strings
2541 nd = int(log10(real(max(1,d))))+1
2542 ns = int(log10(real(max(1,s))))+1
2543 nt = int(log10(real(max(1,ticks))))+1
2544 write (fmt,10) nd, ns, nt
254510 format ('(a,i',i2.2,',a,i',i2.2,',a,i',i2.2,')')
2546
2547 if (present(str)) then
2548 write (unit_in,fmt) trim(str)//' day=', d, ', sec=', s, ', ticks=', ticks
2549 else
2550 write (unit_in,fmt) 'TIME: day=', d, ', sec=', s, ', ticks=', ticks
2551 endif
2552
2553end subroutine print_time
2554
2555!> @brief Prints the time to standard output (or optional unit) as a date.
2556!!
2557!! Prints the given time_type argument as a date (using year, month, day,
2558!! hour, minutes, seconds and ticks). NOTE: there is no check for PE number.
2559subroutine print_date (Time,str,unit)
2560type(time_type) , intent(in) :: time !< Time that will be printed
2561character (len=*), intent(in), optional :: str !< Character string that precedes the printed time
2562integer , intent(in), optional :: unit !< Unit number for printed output, defaults to stdout
2563integer :: y,mo,d,h,m,s, unit_in
2564character(len=9) :: mon
2565
2566! prints the time to standard output (or optional unit) as a date
2567! NOTE: there is no check for PE number
2568
2569 unit_in = stdout()
2570 if (present(unit)) unit_in = unit
2571
2572 call get_date (time,y,mo,d,h,m,s)
2573 mon = month_name(mo)
2574 if (present(str)) then
2575 write (unit_in,10) trim(str)//' ', y,mon(1:3),' ',d,' ',h,':',m,':',s
2576 else
2577 write (unit_in,10) 'DATE: ', y,mon(1:3),' ',d,' ',h,':',m,':',s
2578 endif
257910 format (a,i4,1x,a3,4(a1,i2.2))
2580
2581end subroutine print_date
2582
2583!------------------------------------------------------------------------
2584
2585!> @brief Returns a character string that describes the
2586!! calendar type corresponding to the input integer.
2587!!
2588!> @returns A character string describing the calendar type
2589function valid_calendar_types(ncal, err_msg)
2590integer, intent(in) :: ncal !< Integer corresponding to a valid calendar type
2591character(len=*), intent(out), optional :: err_msg !< Holds an error message when present
2592character(len=24) :: valid_calendar_types
2593character(len=128) :: err_msg_local
2594
2595if(.not.module_is_initialized) call time_manager_init
2596
2597if(present(err_msg)) err_msg = ''
2598
2599if(ncal == no_calendar) then
2600 valid_calendar_types = 'NO_CALENDAR '
2601else if(ncal == thirty_day_months) then
2602 valid_calendar_types = '360_DAY '
2603else if(ncal == julian) then
2604 valid_calendar_types = 'JULIAN '
2605else if(ncal == gregorian) then
2606 valid_calendar_types = 'GREGORIAN '
2607else if(ncal == noleap) then
2608 valid_calendar_types = 'NOLEAP '
2609else
2610 write(err_msg_local,'(a,i4,a)') 'calendar type=',ncal,' is invalid.'
2611 if(error_handler('function valid_calendar_types', err_msg_local, err_msg)) return
2612endif
2613end function valid_calendar_types
2614
2615!------------------------------------------------------------------------
2616
2617!> Get the a character string that represents the time. The format will be
2618!! yyyymmdd.hhmmss
2619function date_to_string(time, err_msg)
2620 type(time_type), intent(in) :: time
2621 character(len=*), intent(out), optional :: err_msg
2622 character(len=128) :: err_msg_local
2623 character(len=15) :: date_to_string
2624 integer :: yr,mon,day,hr,min,sec
2625
2626 if(present(err_msg)) err_msg = ''
2627 call get_date(time,yr,mon,day,hr,min,sec)
2628 if (yr <= 9999) then
2629 write(date_to_string,'(I4.4,I2.2,I2.2,".",I2.2,I2.2,I2.2)') yr, mon, day, hr, min, sec
2630 else
2631 write(err_msg_local, '(a,i4.4,a)') 'year = ', yr, ' should be less than 10000'
2632 if(error_handler('function date_to_string', err_msg_local, err_msg)) return
2633 endif
2634
2635end function date_to_string
2636
2637!> \author Tom Robinson thomas.robinson@noaa.gov
2638!! \brief This routine converts the integer t%days to a string
2639subroutine time_list_error (T,Terr)
2640 type(time_type), intent(in) :: t !< time_type input
2641 character(len=:), allocatable :: terr !< String holding the t%days
2642!> Allocate the string
2643 allocate (character(len=10) :: terr)
2644!> Write the integer to the string
2645 write (terr,'(I0)') t%days
2646end subroutine time_list_error
2647
2648end module time_manager_mod
type(time_type) function real8_to_time_type(x, err_msg)
Convert a real number of seconds into a time_type variable.
subroutine, public print_date(time, str, unit)
Prints the time to standard output (or optional unit) as a date.
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 decrement_time(time, seconds, days, ticks, err_msg, allow_neg_inc)
Decrements a time by seconds and days.
character(len=24) function, public valid_calendar_types(ncal, err_msg)
Returns a character string that describes the calendar type corresponding to the input integer.
type(time_type) function scalar_time_mult(n, time)
Returns time multipled by integer factor n.
integer function, public get_ticks_per_second()
Returns the number of ticks per second.
subroutine, public time_list_error(t, terr)
This routine converts the integer tdays to a string.
logical function increment_date_private(time, years, months, days, hours, minutes, seconds, ticks, time_out, err_msg)
Given a time and some date increment, computes a new time. Depending on the mapping algorithm from da...
integer function, public get_calendar_type()
Returns default calendar type for mapping from time to date.
logical function leap_year_julian(time)
Returns the number of days in a julian month.
type(time_type) function, public increment_time(time, seconds, days, ticks, err_msg, allow_neg_inc)
Increments a time by seconds and days.
logical function increment_time_private(time_in, seconds, days, ticks, time_out, err_msg)
Increments a time by seconds, days and ticks.
subroutine, public set_ticks_per_second(tps)
Sets the number of ticks per second.
logical function set_date_gregorian(year, month, day, hour, minute, second, tick, time_out, err_msg)
Sets Time_outdays on a Gregorian calendar Computes the total number of days between 1/1/0001 to the c...
type(time_type) function, public length_of_year()
Returns the mean length of the year in the default calendar setting.
logical function set_date_private(year, month, day, hour, minute, second, tick, time_out, err_msg)
Sets days for different calendar types. Given an input date in year, month, days, etc....
integer function days_in_month_no_leap(time)
Returns the number of days in a 365 day year month.
subroutine, public set_calendar_type(type, err_msg)
Sets calendar_type for mapping an interval to a date. For the Gregorian calendar, negative years and ...
integer function, public days_in_month(time, err_msg)
Given a time, computes the corresponding date given the selected date time mapping algorithm.
type(time_type) function set_date_i(year, month, day, hour, minute, second, tick, err_msg)
Calls set_date_private to set days for different calendar types.
subroutine get_date_gregorian(time, year, month, day, hour, minute, second, tick)
Gets the date on a Gregorian calendar. Computes the year, month, day on the fly from the quantity tim...
type(time_type) function set_date_c(string, zero_year_warning, err_msg, allow_rounding)
Calls set_date_private for different calendar types when given a string input. Examples of acceptable...
character(len=9) function, public month_name(n)
Returns a character string containing the name of the month corresponding to month number n.
logical function set_time_private(seconds, days, ticks, time_out, err_msg)
Returns a time interval corresponding to this number of days, seconds, and ticks. days,...
logical function leap_year_no_leap(time)
Another tough one; no leap year returns false for leap year inquiry.
integer function days_in_month_gregorian(time)
Returns the number of days in a gregorian month.
type(time_type) function time_scalar_divide(time, n)
Returns the largest time, t, for which n * t <= time.
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...
type(time_type) function set_time_i(seconds, days, ticks, err_msg)
Returns a time_type set to the given amount of time via integer amounts.
integer, parameter days_in_400_year_period
Used only for gregorian.
type(time_type) function time_scalar_mult(time, n)
Returns time multiplied by integer factor n.
integer function days_in_month_julian(time)
Returns the number of days in a julian month.
logical function, public leap_year(time, err_msg)
Returns true if the year corresponding to the input time is a leap year (for default calendar)....
character(len=15) function, public date_to_string(time, err_msg)
Get the a character string that represents the time. The format will be yyyymmdd.hhmmss.
logical function, public repeat_alarm(time, alarm_frequency, alarm_length)
Repeat_alarm supports an alarm that goes off with alarm_frequency and lasts for alarm_length....
integer function time_divide(time1, time2)
Returns the largest integer, n, for which time1 >= time2 * n.
logical function leap_year_thirty(time)
No leap years in thirty day months, included for transparency.
subroutine get_date_julian_private(time, year, month, day, hour, minute, second, tick)
Base date for Julian calendar is year 1 with all multiples of 4 years being leap years.
type(time_type) function set_time_c(string, err_msg, allow_rounding)
Returns a time_type set to the given amount of time via a string.
integer function days_in_month_thirty(time)
Returns the number of days in a thirty day month (needed for transparent changes to calendar type).
integer function, public day_of_year(time)
Returns number of day in year for given time. Jan 1st is day 1, not zero!
real(kind=r8_kind) function, public time_type_to_real(time)
Converts time to seconds and returns it as a real number.
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...
type(time_type) function real4_to_time_type(x, err_msg)
subroutine get_date_thirty(time, year, month, day, hour, minute, second, tick)
Computes date corresponding to time interval for 30 day months, 12 month years.
subroutine, public print_time(time, str, unit)
Prints the given time_type argument as a time (using days, seconds and ticks)
subroutine, public time_manager_init()
Initialization routine. Writes the version information to the log file.
integer function safe_rtoi(rval, mode)
Convert a floating point value to an integer value.
logical function, public interval_alarm(time, time_interval, alarm, alarm_interval)
Supports a commonly used type of test on times for models. Given the current time,...
type(time_type) function, public decrement_date(time, years, months, days, hours, minutes, seconds, ticks, err_msg, allow_neg_inc)
Given a time and some date decrement, computes a new time. Depending on the mapping algorithm from da...
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.
Wrapper for the real to time interface Takes seconds as reals to convert to a time_type representatio...