source: palm/trunk/SOURCE/date_and_time_mod.f90 @ 3458

Last change on this file since 3458 was 3458, checked in by kanani, 3 years ago

Reintegrated fixes/changes from branch chemistry

  • Property svn:keywords set to Id
File size: 20.7 KB
Line 
1!> @file date_and_time_mod.f90
2!------------------------------------------------------------------------------!
3! This file is part of the PALM model system.
4!
5! PALM is free software: you can redistribute it and/or modify it under the
6! terms of the GNU General Public License as published by the Free Software
7! Foundation, either version 3 of the License, or (at your option) any later
8! version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
17! Copyright 1997-2018 Leibniz Universitaet Hannover
18!------------------------------------------------------------------------------!
19!
20! Current revisions:
21! ------------------
22!
23!
24! Former revisions:
25! -----------------
26! $Id: date_and_time_mod.f90 3458 2018-10-30 14:51:23Z kanani $
27! from chemistry branch r3443, banzhafs:
28! Added initial hour_of_day, hour_of_year, day_of_year and month_of_year to
29! init_date_and_time
30!
31! 3298 2018-10-02 12:21:11Z kanani
32! - Minor formatting (kanani)
33! - Added Routines for DEFAULT mode of chemistry emissions (Russo)
34! - Added routine for reading-in real dates: format has to be in DDMMYYYY and
35!   passed in the namelist parameter date_init (Russo)
36! - Added calculation of several time indices useful in several routines
37!   of the model (Russo)
38!
39! 2718 2018-01-02 08:49:38Z maronga
40! Corrected "Former revisions" section
41!
42! 2701 2017-12-15 15:40:50Z suehring
43! Changes from last commit documented
44!
45! 2698 2017-12-14 18:46:24Z suehring
46! Bugfix in definition of d_seconds_year.
47!
48! 2696 2017-12-14 17:12:51Z kanani
49! Change in file header (GPL part)
50!
51! 2544 2017-10-13 18:09:32Z maronga
52! Initial revision
53!
54!
55
56!
57! Description:
58! ------------
59!> This routine calculates all needed information on date and time used by
60!> other modules
61!> @todo Further testing and revision of routines for updating indices of
62!>       emissions in the default mode
63!> @todo Add routine for recognizing leap years
64!> @todo Add recognition of exact days of week (Monday, Tuesday, etc.)
65!> @todo Reconsider whether to remove day_of_year_init from the namelist: we
66!>       already implemented changes for calculating it from date_init in
67!>       calc_date_and_time
68!> @todo time_utc during spin-up 
69!------------------------------------------------------------------------------!
70 MODULE date_and_time_mod
71
72    USE control_parameters,                                                   &
73        ONLY:  coupling_start_time, days_since_reference_point,               &
74               message_string, simulated_time, time_since_reference_point
75
76    USE kinds
77
78
79    IMPLICIT NONE
80
81    PRIVATE
82
83!-- Variables Declaration
84
85    INTEGER(iwp)        ::  day_of_year      = 0   !< day of the year (DOY)
86    INTEGER(iwp)        ::  day_of_year_init = 0   !< DOY at model start (default: 0)
87
88    ! --- Most of these indices are updated by the routine calc_date_and_time according to the current date and time of the simulation
89    INTEGER(iwp)        ::  hour_of_year = 1                        !< hour of the current year (1:8760(8784))
90    INTEGER(iwp)        ::  hour_of_day=1                           !< hour of the current day (1:24)
91    INTEGER(iwp)        ::  day_of_month=0                          !< day of the current month (1:31)
92    INTEGER(iwp)        ::  month_of_year=0                         !< month of the current year (1:12)
93    INTEGER(iwp)        ::  current_year=0                          !< current year
94    INTEGER(iwp)        ::  hour_call_emis=0                        !< index used to call the emissions just once every hour
95
96    INTEGER(iwp)        ::  ihour           
97    INTEGER(iwp)        ::  index_mm                                !< index months of the default emission mode
98    INTEGER(iwp)        ::  index_dd                                !< index days of the default emission mode
99    INTEGER(iwp)        ::  index_hh                                !< index hours of the emission mode
100    INTEGER(iwp)        ::  hours_since_reference_point             !< hours of current simulation
101
102    REAL(wp)            ::  time_utc                     !< current model time in UTC
103    REAL(wp)            ::  time_utc_emis                !< current emission module time in UTC
104    REAL(wp)            ::  time_utc_init = 43200.0_wp   !< UTC time at model start
105    REAL(wp)            ::  time_update                  !< used to calculate actual second of the simulation
106 
107    REAL(wp), PARAMETER ::  d_hours_day    = 1.0_wp / 24.0_wp       !< inverse of hours per day (1/24)
108    REAL(wp), PARAMETER ::  d_seconds_hour = 1.0_wp / 3600.0_wp     !< inverse of seconds per hour (1/3600)
109    REAL(wp), PARAMETER ::  d_seconds_year = 1.0_wp / 31536000.0_wp !< inverse of the seconds per year (1/(365*86400))
110   
111    CHARACTER(len=8)    ::  date_init = "21062017"                  !< Starting date of simulation: We selected this because it was a monday
112 
113    !> --- Parameters
114    INTEGER, PARAMETER, DIMENSION(12) :: days = (/31,28,31,30,31,30,31,31,30,31,30,31/) ! total number of days for each month (no leap year)
115
116    SAVE
117
118!-- INTERFACES PART
119    !-- Read initial day and time of simulation
120    INTERFACE init_date_and_time
121       MODULE PROCEDURE init_date_and_time
122    END INTERFACE init_date_and_time
123
124    !-- Get hour index in the DEAFULT case of chemistry emissions :
125    INTERFACE time_default_indices
126       MODULE PROCEDURE time_mdh_indices
127       MODULE PROCEDURE time_hour_indices
128    END INTERFACE time_default_indices
129
130    !-- Get hour index in the PRE-PROCESSED case of chemistry emissions :
131    INTERFACE time_preprocessed_indices
132       MODULE PROCEDURE time_preprocessed_indices
133    END INTERFACE time_preprocessed_indices
134
135
136    !-- Calculate current date and time
137    INTERFACE calc_date_and_time
138       MODULE PROCEDURE calc_date_and_time
139    END INTERFACE
140
141
142    !-- Public Interfaces
143    PUBLIC calc_date_and_time, time_default_indices, init_date_and_time, time_preprocessed_indices
144
145    !-- Public Variables
146    PUBLIC date_init, d_hours_day, d_seconds_hour, d_seconds_year,               &
147           day_of_year, day_of_year_init, time_utc, time_utc_init, day_of_month, &
148           month_of_year, index_mm, index_dd, index_hh, hour_of_day, hour_of_year, &
149           hour_call_emis
150
151 CONTAINS
152
153
154 !------------------------------------------------------------------------------!
155 !> Reads starting date from namelist
156 !------------------------------------------------------------------------------!
157 
158    SUBROUTINE init_date_and_time
159
160       IMPLICIT NONE
161
162       !--    Variables Definition
163       INTEGER ::  i_mon       !< Index for going through the different months
164
165       IF  (day_of_year_init == 0) THEN
166          ! Day of the month at starting time
167          READ(UNIT=date_init(1:2),fmt=*)day_of_month
168
169          ! Month of the year at starting time
170          READ(UNIT=date_init(3:4),fmt=*)month_of_year
171
172          ! Year at starting time
173          READ(UNIT=date_init(5:8),fmt=*)current_year
174       
175       ENDIF
176
177
178       !-- Calculate initial hour of the day: the first hour of the day is from 00:00:00 to 00:59:59.
179
180       hour_of_day = INT( FLOOR( time_utc_init/3600.0_wp ) ) + 1
181
182       !-- Calculate initial day day_of_year_init in case date_init is given or day_of_year_init is given
183       IF ( day_of_year_init == 0 ) THEN
184
185          !> Condition for printing an error when date_init is not provided when day_of_year_init is not given in the namelist or when the format of the date is not the one required by PALM.
186          IF ( day_of_month .GT. 0 .AND. day_of_month .LE. 31 .AND. month_of_year .GT. 0 .AND. month_of_year .LE. 12) THEN
187       
188             IF ( month_of_year == 1 ) THEN  !!month of year is read in input
189
190                day_of_year_init = day_of_month
191
192             ELSE
193
194                day_of_year_init= SUM(days( 1:(month_of_year-1) )) + day_of_month  !day_of_month is read in input in this case
195
196             ENDIF
197!kanani: Revise, we cannot force users to provide date_init, maybe set a default value?
198!           ELSE
199!
200!              message_string = 'date_init not provided in the namelist or'             //          &
201!                               ' given in the wrong format: MUST BE DDMMYYYY'                       
202!              CALL message( 'calc_date_and_time', 'DT0100', 2, 2, 0, 6, 0 )
203     
204          ENDIF
205
206       ENDIF
207
208       
209       !-- Initial day of the year
210       day_of_year = day_of_year_init
211
212       !-- Initial hour of the year
213       hour_of_year = ( (day_of_year-1) * 24 ) + hour_of_day
214
215       !--Initial day of the month and month of the year
216       !> --------------------------------------------------------------------------------
217       !> The first case is when date_init is not provided: we only know day_of_year_init     
218       IF ( month_of_year == 0 .AND. day_of_month == 0) THEN
219
220         
221          IF ( day_of_year .LE. 31 ) THEN
222
223             month_of_year=1
224             day_of_month=day_of_year
225
226          ELSE
227
228             DO i_mon=2,12   !january is considered in the first case
229                IF ( day_of_year .LE. SUM(days(1:i_mon)) .AND. day_of_year .GT. SUM(days(1:(i_mon-1))) ) THEN
230           
231                   month_of_year=i_mon
232
233                   day_of_month=INT(MOD(day_of_year, SUM(days(1:(i_mon-1))))) 
234
235                   GOTO 38
236
237                ENDIF
238
239             38 ENDDO
240          ENDIF
241       !> --------------------------------------------------------------------------------
242       !> in the second condition both day of month and month_of_year are either given in input (passed to date_init) or we are in some day successive to the initial one, so that day_of_month has already be computed in previous step
243       !>TBD: something to calculate the current year is missing
244       ELSEIF ( day_of_month .GT. 0 .AND. day_of_month .LE. 31 .AND. month_of_year .GT. 0 .AND. month_of_year .LE. 12) THEN
245 
246          !> calculate month_of_year. TBD: test the condition when day_of_year==31
247 
248          IF (day_of_year==1) THEN  !> this allows to turn from december to January when passing from a year to another
249 
250             month_of_year = 1
251       
252          ELSE IF (day_of_year .GT. 1 .AND. day_of_year .GT. SUM(days(1:month_of_year))) THEN
253
254             month_of_year = month_of_year + 1
255
256          ENDIF
257
258          !> calculate day_of_month
259          IF ( month_of_year == 1 ) THEN
260           
261            day_of_month=day_of_year
262
263          ELSE
264
265            day_of_month=INT(MOD(day_of_year, SUM(days(1:(month_of_year-1)))))
266
267          ENDIF
268
269
270       ELSE
271
272          !> Condition when date_init is provided but it is given in the wrong format
273          message_string = 'date_init not provided in the namelist or'            //          &
274                              ' given in the wrong format: MUST BE DDMMYYYY'                 
275          CALL message( 'init_date_and_time', 'DT0102', 2, 2, 0, 6, 0 ) 
276
277       ENDIF
278
279
280    END SUBROUTINE init_date_and_time
281
282!------------------------------------------------------------------------------!
283! Description:
284! ------------
285!> Calculate current date and time of the simulation
286!------------------------------------------------------------------------------!
287 
288    SUBROUTINE calc_date_and_time
289
290       IMPLICIT NONE
291
292!--    Variables Definition
293       INTEGER                          :: i_mon       !< Index for going through the different months
294
295       !> Update simulation time in seconds
296       time_update = simulated_time-coupling_start_time
297
298!--    Calculate current day of the simulated time
299       days_since_reference_point=INT(FLOOR( (time_utc_init + time_update) &
300                               / 86400.0_wp ) )
301
302!--    Calculate actual UTC time                       
303       time_utc = MOD((time_utc_init + time_since_reference_point), 86400.0_wp)
304       
305!sB    PRILIMINARY workaround for time_utc changes due to changes in time_since_reference_point in
306!sB    radiation_model_mod during runtime:
307       time_utc_emis = MOD((time_utc_init + time_update), 86400.0_wp)     
308
309!--    Calculate initial day of the year: it is calculated only once. In fact, day_of_year_init is initialized to 0 and then a positive value is passed. This condition is also called only when day_of_year_init is not given in the namelist.
310
311       IF ( day_of_year_init == 0 ) THEN
312
313          !> Condition for printing an error when date_init is not provided when day_of_year_init is not given in the namelist or when the format of the date is not the one required by PALM.
314          IF ( day_of_month .GT. 0 .AND. day_of_month .LE. 31 .AND. month_of_year .GT. 0 .AND. month_of_year .LE. 12) THEN
315       
316             IF ( month_of_year == 1 ) THEN  !!month of year is read in input
317
318                day_of_year_init = day_of_month
319
320             ELSE
321
322                day_of_year_init= SUM(days( 1:(month_of_year-1) )) + day_of_month  !day_of_month is read in input in this case
323
324             ENDIF
325!kanani: Revise, we cannot force users to provide date_init, maybe set a default value?
326!           ELSE
327!
328!              message_string = 'date_init not provided in the namelist or'             //          &
329!                               ' given in the wrong format: MUST BE DDMMYYYY'                       
330!              CALL message( 'calc_date_and_time', 'DT0100', 2, 2, 0, 6, 0 )
331     
332          ENDIF
333
334       ENDIF
335
336      !-- Calculate actual hour of the day: the first hour of the day is from 00:00:00 to 00:59:59.
337
338      hour_of_day = INT( FLOOR( time_utc_emis/3600.0_wp ) ) + 1
339
340!--    Calculate current day of the year !TBD: considetr leap years
341       IF ( (day_of_year_init + days_since_reference_point)  .GT. 365 ) THEN
342
343          day_of_year=INT(MOD((day_of_year_init + days_since_reference_point), 365.0_wp))
344
345       ELSE
346         
347          day_of_year = day_of_year_init + days_since_reference_point
348
349       ENDIF
350
351!
352!--    Calculate current hour of the year
353       hour_of_year = ( (day_of_year-1) * 24 ) + hour_of_day  !> actual hour of the year
354       
355
356!
357!--    UPDATE actual day of the month and month of the year
358       !> --------------------------------------------------------------------------------
359       !> The first case is when date_init is not provided: we only know day_of_year_init     
360       IF ( month_of_year == 0 .AND. day_of_month == 0) THEN
361
362          !> The first case is when date_init is not provided: we only know day_of_year_init
363          !DO i_mon=1,12
364             !IF (day_of_year .LE. SUM(days(1:i_mon))) THEN
365          IF ( day_of_year .LE. 31 ) THEN
366
367             month_of_year=1
368             day_of_month=day_of_year
369
370          ELSE
371
372             DO i_mon=2,12   !january is considered in the first case
373                IF ( day_of_year .LE. SUM(days(1:i_mon)) .AND. day_of_year .GT. SUM(days(1:(i_mon-1))) ) THEN
374           
375                   month_of_year=i_mon
376
377                   day_of_month=INT(MOD(day_of_year, SUM(days(1:(i_mon-1))))) 
378
379                   GOTO 38
380
381                ENDIF
382
383             38 ENDDO
384          ENDIF
385       !> --------------------------------------------------------------------------------
386       !> in the second condition both day of month and month_of_year are either given in input (passed to date_init) or we are in some day successive to the initial one, so that day_of_month has already be computed in previous step
387       !>TBD: something to calculate the current year is missing
388       ELSEIF ( day_of_month .GT. 0 .AND. day_of_month .LE. 31 .AND. month_of_year .GT. 0 .AND. month_of_year .LE. 12) THEN
389 
390          !> calculate month_of_year. TBD: test the condition when day_of_year==31
391 
392          IF (day_of_year==1) THEN  !> this allows to turn from december to January when passing from a year to another
393 
394             month_of_year = 1
395       
396          ELSE IF (day_of_year .GT. 1 .AND. day_of_year .GT. SUM(days(1:month_of_year))) THEN
397
398             month_of_year = month_of_year + 1
399
400          ENDIF
401
402          !> calculate day_of_month
403          IF ( month_of_year == 1 ) THEN
404           
405            day_of_month=day_of_year
406
407          ELSE
408
409            day_of_month=INT(MOD(day_of_year, SUM(days(1:(month_of_year-1)))))
410
411          ENDIF
412
413
414
415
416      ELSE
417
418          !> Condition when date_init is provided but it is given in the wrong format
419          message_string = 'date_init not provided in the namelist or'            //          &
420                              ' given in the wrong format: MUST BE DDMMYYYY'                 
421          CALL message( 'calc_date_and_time', 'DT0101', 2, 2, 0, 6, 0 ) 
422
423      ENDIF       
424
425    END SUBROUTINE calc_date_and_time
426
427
428!------------------------------------------------------------------------------!
429! Description:
430! ------------
431!> This routine determines the time factor index in the PRE-PROCESSED emissions mode.
432!------------------------------------------------------------------------------!
433
434 SUBROUTINE time_preprocessed_indices(index_hh)
435
436    USE indices
437
438    IMPLICIT NONE
439
440!
441!-- In/output
442    INTEGER, INTENT(INOUT) ::  index_hh    !> Index Hour
443!
444!-- Additional Variables for calculateing indices
445!-- Constants
446    INTEGER, PARAMETER ::  nhour = 24
447
448    IF (days_since_reference_point == 0) THEN
449
450       index_hh=hour_of_day
451
452    ELSE
453
454       index_hh=(days_since_reference_point*nhour)+(hour_of_day)
455
456    ENDIF
457
458
459 END SUBROUTINE time_preprocessed_indices
460
461
462!------------------------------------------------------------------------------!
463! Description:
464! ------------
465!> This routine determines the time factor index in the mdh case of the DEFAULT chemistry emissions mode.
466!------------------------------------------------------------------------------!
467
468 SUBROUTINE time_mdh_indices(daytype_mdh,mo, dd, hh, index_mm, index_dd, index_hh)
469
470    USE indices
471
472    IMPLICIT NONE
473
474    !> IN/OUTPUT
475    INTEGER, INTENT(INOUT) :: mo          !> Month of year
476    INTEGER, INTENT(INOUT) :: dd          !> Day of month
477    INTEGER, INTENT(INOUT) :: hh          !> Hour of day
478    INTEGER, INTENT(INOUT) :: index_mm    !> Index Month
479    INTEGER, INTENT(INOUT) :: index_dd    !> Index Day
480    INTEGER, INTENT(INOUT) :: index_hh    !> Index Hour
481
482    CHARACTER(len=80), INTENT(INOUT) :: daytype_mdh !> type of the day in mdh mode: one of 1-WORKDAY
483                                                    !                                      2-WEEKEND
484                                                    !                                      3-HOLIDAY
485
486    REAL(wp)                         :: frac_day=0 
487
488    !> ------------------------------------------------------------------------
489
490    INTEGER                          :: weekday
491
492    !> CONSTANTS
493    INTEGER, PARAMETER               :: nmonth = 12
494    INTEGER, PARAMETER               :: nday = 7
495    INTEGER, PARAMETER               :: nhour = 24
496
497    frac_day= (dd-1)/nday    !> indicates the week of the month, supposing the month starts on monday
498
499   ! 1:7      1:31  7                   (0:30)/7
500    weekday = dd-( nday * (INT( CEILING( frac_day ) ) ) ) ! for now we let the year start on Monday.
501
502    !TBD: set weekday correct based on date
503    index_mm = mo
504    index_dd = nmonth + weekday  !> Index of the days in the mdh mode (13:20)
505
506    SELECT CASE(TRIM(daytype_mdh))
507
508                   CASE ("workday")
509               
510                   index_hh = nmonth+ nday + hh
511
512                   CASE ("weekend")
513               
514                   index_hh = nmonth+ nday + nhour + hh
515
516                   CASE ("holiday")
517               
518                   index_hh = nmonth+ nday + 2*nhour + hh
519    END SELECT
520
521
522 END SUBROUTINE time_mdh_indices
523
524!------------------------------------------------------------------------------!
525! Description:
526! ------------
527!> This routine determines the time factor index in the HOURLY case of the DEFAULT emissions mode.
528!------------------------------------------------------------------------------!
529
530 SUBROUTINE time_hour_indices(mo,dd,hh,index_hh)
531
532    USE indices
533
534    IMPLICIT NONE
535
536    !> IN/OUTPUT
537    INTEGER, INTENT(INOUT)              :: mo          !> Month
538    INTEGER, INTENT(INOUT)              :: hh          !> Hour
539    INTEGER, INTENT(INOUT)              :: dd          !> Day
540    INTEGER, INTENT(INOUT)              :: index_hh    !> Index Hour
541 
542    !> Additional Variables for calculateing indices
543    INTEGER                             :: index_mm    !> Index Month
544    INTEGER                             :: index_dd    !> Index Day
545    INTEGER                             :: i_mon       !> Index for going through the different months
546    INTEGER                             :: sum_dd      !> Sum days
547
548    !> CONSTANTS
549    INTEGER, PARAMETER                  :: nmonth=12
550    INTEGER, PARAMETER                  :: nday = 7
551    INTEGER, PARAMETER                  :: nhour = 24
552    INTEGER, PARAMETER, DIMENSION(12)   :: days = (/31,28,31,30,31,30,31,31,30,31,30,31/) ! no leap year
553   
554
555    index_mm = mo-1
556    index_dd = dd-1
557    sum_dd=0
558
559    IF (mo == 1) THEN
560
561       index_hh=(index_dd*nhour)+hh
562
563    ELSE
564
565        DO i_mon=1,index_mm
566       
567            sum_dd=sum_dd+days(i_mon)
568
569        ENDDO
570     
571    index_hh=(sum_dd*nhour)+(index_dd*nhour)+(hh)
572
573    ENDIF
574 
575
576 END SUBROUTINE time_hour_indices
577
578
579END MODULE date_and_time_mod
Note: See TracBrowser for help on using the repository browser.