GammaLib  1.7.0.dev
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
GTime.cpp
Go to the documentation of this file.
1 /***************************************************************************
2  * GTime.cpp - Time class *
3  * ----------------------------------------------------------------------- *
4  * copyright (C) 2010-2019 by Juergen Knoedlseder *
5  * ----------------------------------------------------------------------- *
6  * *
7  * This program is free software: you can redistribute it and/or modify *
8  * it under the terms of the GNU General Public License as published by *
9  * the Free Software Foundation, either version 3 of the License, or *
10  * (at your option) any later version. *
11  * *
12  * This program is distributed in the hope that it will be useful, *
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of *
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
15  * GNU General Public License for more details. *
16  * *
17  * You should have received a copy of the GNU General Public License *
18  * along with this program. If not, see <http://www.gnu.org/licenses/>. *
19  * *
20  ***************************************************************************/
21 /**
22  * @file GTime.cpp
23  * @brief Time class implementation
24  * @author Juergen Knoedlseder
25  */
26 
27 /* __ Includes ___________________________________________________________ */
28 #ifdef HAVE_CONFIG_H
29 #include <config.h>
30 #endif
31 #include <ctime>
32 #include <cstring> // std::memcpy
33 #include <cstdio>
34 #include "GTools.hpp"
35 #include "GMath.hpp"
36 #include "GException.hpp"
37 #include "GTime.hpp"
38 #include "GTimeReference.hpp"
39 
40 /* __ Constants __________________________________________________________ */
41 const double mjd_ref = 55197.000766018518519; //!< MJD of time=0
42 const double jd_ref = mjd_ref + 2400000.5; //!< JD of time=0
43 
44 /* __ Method name definitions ____________________________________________ */
45 #define G_CONSTRUCT "GTime::GTime(double&, std::string&)"
46 #define G_SECS_GET "GTime::secs(std::string&)"
47 #define G_SECS_SET "GTime::secs(double&, std::string&)"
48 #define G_UTC "GTime::utc(std::string&)"
49 #define G_UTC_GET "GTime::utc(int&)"
50 
51 /* __ Macros _____________________________________________________________ */
52 
53 /* __ Coding definitions _________________________________________________ */
54 
55 /* __ Debug definitions __________________________________________________ */
56 
57 
58 /*==========================================================================
59  = =
60  = Constructors/destructors =
61  = =
62  ==========================================================================*/
63 
64 /***********************************************************************//**
65  * @brief Void constructor
66  ***************************************************************************/
68 {
69  // Initialise private members
70  init_members();
71 
72  // Return
73  return;
74 }
75 
76 
77 /***********************************************************************//**
78  * @brief Copy constructor
79  *
80  * @param[in] time Time.
81  ***************************************************************************/
82 GTime::GTime(const GTime& time)
83 {
84  // Initialise private members
85  init_members();
86 
87  // Copy members
88  copy_members(time);
89 
90  // Return
91  return;
92 }
93 
94 
95 /***********************************************************************//**
96  * @brief Time constructor
97  *
98  * @param[in] time Time value in TT (seconds or days).
99  * @param[in] unit Time unit string.
100  *
101  * @exception GException::time_invalid_unit
102  * Invalid time unit specified.
103  *
104  * Constructs a GTime object by setting the time in the native reference
105  * in the TT time system in units of seconds (default) or days.
106  ***************************************************************************/
107 GTime::GTime(const double& time, const std::string& unit)
108 {
109  // Initialise private members
110  init_members();
111 
112  // Set time according to timeunit string
113  std::string timeunit = gammalib::tolower(unit);
114  if (timeunit == "d" || timeunit == "day" || timeunit == "days") {
115  days(time);
116  }
117  else if (timeunit == "s" || timeunit == "sec" || timeunit == "secs") {
118  secs(time);
119  }
120  else {
122  "Valid timeunit values are: \"d\", \"day\", \"days\","
123  " \"s\", \"sec\" or \"secs\"");
124  }
125 
126  // Return
127  return;
128 }
129 
130 
131 /***********************************************************************//**
132  * @brief Time constructor
133  *
134  * @param[in] time Time in given reference system.
135  * @param[in] ref Reference system.
136  *
137  * Constructs a GTime object by setting the time to a value given in a
138  * specific reference system.
139  ***************************************************************************/
140 GTime::GTime(const double& time, const GTimeReference& ref)
141 {
142  // Initialise private members
143  init_members();
144 
145  // Set time
146  set(time, ref);
147 
148  // Return
149  return;
150 }
151 
152 
153 /***********************************************************************//**
154  * @brief Time constructor
155  *
156  * @param[in] time Time string in given reference system.
157  * @param[in] ref Reference system.
158  *
159  * Constructs a GTime object by setting the time to a string value. See the
160  * set(std::string&, GTimeReference&) method for valid time strings.
161  ***************************************************************************/
162 GTime::GTime(const std::string& time, const GTimeReference& ref)
163 {
164  // Initialise private members
165  init_members();
166 
167  // Set time
168  set(time, ref);
169 
170  // Return
171  return;
172 }
173 
174 
175 /***********************************************************************//**
176  * @brief Time constructor
177  *
178  * @param[in] time Time string.
179  *
180  * Constructs a GTime object by setting the time to a string value. See the
181  * set(const std::string& time) method for valid time strings.
182  ***************************************************************************/
183 GTime::GTime(const std::string& time)
184 {
185  // Initialise private members
186  init_members();
187 
188  // Set time
189  set(time);
190 
191  // Return
192  return;
193 }
194 
195 
196 /***********************************************************************//**
197  * @brief Destructor
198  ***************************************************************************/
200 {
201  // Free members
202  free_members();
203 
204  // Return
205  return;
206 }
207 
208 
209 /*==========================================================================
210  = =
211  = Operators =
212  = =
213  ==========================================================================*/
214 
215 /***********************************************************************//**
216  * @brief Assignment operator
217  *
218  * @param[in] time Time.
219  * @return Time.
220  ***************************************************************************/
222 {
223  // Execute only if object is not identical
224  if (this != &time) {
225 
226  // Free members
227  free_members();
228 
229  // Initialise private members for clean destruction
230  init_members();
231 
232  // Copy members
233  copy_members(time);
234 
235  } // endif: object was not identical
236 
237  // Return
238  return *this;
239 }
240 
241 
242 /*==========================================================================
243  = =
244  = Public methods =
245  = =
246  ==========================================================================*/
247 
248 /***********************************************************************//**
249  * @brief Clear time
250  ***************************************************************************/
251 void GTime::clear(void)
252 {
253  // Free members
254  free_members();
255 
256  // Initialise members
257  init_members();
258 
259  // Return
260  return;
261 }
262 
263 
264 /***********************************************************************//**
265  * @brief Clone time
266  *
267  * @return Pointer to deep copy of time.
268  ***************************************************************************/
269 GTime* GTime::clone(void) const
270 {
271  // Clone this image
272  return new GTime(*this);
273 }
274 
275 
276 /***********************************************************************//**
277  * @brief Return time in Julian Days (TT)
278  *
279  * @return Time in Julian Days (TT) (days).
280  *
281  * Returns the time in Julian Days (JD) in the Terrestrial Time (TT) system.
282  ***************************************************************************/
283 double GTime::jd(void) const
284 {
285  // Convert time from MET to JD
286  double jd = m_time * gammalib::sec2day + jd_ref;
287 
288  // Return Julian Days
289  return jd;
290 }
291 
292 
293 /***********************************************************************//**
294  * @brief Return time in Julian Days for time system
295  *
296  * @param[in] timesys Time system (one of "TT", "TAI", "UTC")
297  * @return Time in Julian Days (days).
298  *
299  * Returns the time in Julian Days (JD) for the specified time system.
300  ***************************************************************************/
301 double GTime::jd(const std::string& timesys) const
302 {
303  // Convert time to Julian Days
304  double jd = secs(timesys) * gammalib::sec2day + jd_ref;
305 
306  // Return Julian Days
307  return jd;
308 }
309 
310 
311 /***********************************************************************//**
312  * @brief Return time in Modified Julian Days (TT)
313  *
314  * @return Time in Modified Julian Days (TT) (days).
315  *
316  * Returns the time in Modified Julian Days (MJD) in the Terrestrial Time
317  * (TT) system.
318  ***************************************************************************/
319 double GTime::mjd(void) const
320 {
321  // Convert time to Modified Julian Days
322  double mjd = m_time * gammalib::sec2day + mjd_ref;
323 
324  // Return Modified Julian Days
325  return mjd;
326 }
327 
328 
329 /***********************************************************************//**
330  * @brief Return time in Modified Julian Days for time system
331  *
332  * @param[in] timesys Time system (one of "TT", "TAI", "UTC")
333  * @return Time in Modified Julian Days (days).
334  *
335  * Returns the time in Modified Julian Days (JD) for the specified time
336  * system.
337  ***************************************************************************/
338 double GTime::mjd(const std::string& timesys) const
339 {
340  // Convert time to Modified Julian Days
341  double mjd = secs(timesys) * gammalib::sec2day + mjd_ref;
342 
343  // Return Modified Julian Days
344  return mjd;
345 }
346 
347 
348 /***********************************************************************//**
349  * @brief Return time in seconds in native reference for time system
350  *
351  * @param[in] timesys Time system (one of "TT", "TAI", "UTC")
352  * @return Time (seconds).
353  ***************************************************************************/
354 double GTime::secs(const std::string& timesys) const
355 {
356  // Initialise time
357  double time;
358 
359  // If system is TT then return the stored time ...
360  if (timesys == "TT") {
361  time = m_time;
362  }
363 
364  // ... otherwise if the system if TAI then subtract an offset ...
365  else if (timesys == "TAI") {
366  time = m_time - gammalib::tai2tt;
367  }
368 
369  // ... otherwise if the system is UTC then convert the time from TT to TAI
370  // and subtract the leap seconds. The repeated calling of the leap second
371  // method is a kluge to converge as the argument is given in the UTC system
372  // but upon start we're in the TAI system.
373  else if (timesys == "UTC") {
374  time = m_time - gammalib::tai2tt; // Time in TAI
375  double mjd = time * gammalib::sec2day + mjd_ref;
376  double leaps = leap_seconds(mjd);
377  leaps = leap_seconds(mjd - leaps * gammalib::sec2day);
378  leaps = leap_seconds(mjd - leaps * gammalib::sec2day);
379  time -= leaps;
380  }
381 
382  // ... otherwise throw an exception
383  else {
384  std::string msg = "Unknown time system \""+timesys+"\". Either specify "
385  "\"TT\", \"TAI\" or \"UTC\".";
387  }
388 
389  // Return time
390  return time;
391 }
392 
393 
394 /***********************************************************************//**
395  * @brief Return time in days in native reference (TT)
396  *
397  * @return Time in native reference (days).
398  ***************************************************************************/
399 double GTime::days(void) const
400 {
401  // Return time
402  return (m_time * gammalib::sec2day);
403 }
404 
405 
406 /***********************************************************************//**
407  * @brief Return time in days in native reference for time system
408  *
409  * @param[in] timesys Time system (one of "TT", "TAI", "UTC")
410  * @return Time in native reference (days).
411  ***************************************************************************/
412 double GTime::days(const std::string& timesys) const
413 {
414  // Return time
415  return (secs(timesys) * gammalib::sec2day);
416 }
417 
418 
419 /***********************************************************************//**
420  * @brief Return Julian epoch in native reference (TT)
421  *
422  * @return Julian epoch (years).
423  ***************************************************************************/
424 double GTime::julian_epoch(void) const
425 {
426  // Compute Julian epoch
427  double julian_epoch = 2000.0 + (jd() - 2451545.0) / 365.25;
428 
429  // Return Julian epoch
430  return (julian_epoch);
431 }
432 
433 
434 /***********************************************************************//**
435  * @brief Return Julian epoch in native reference for time system
436  *
437  * @param[in] timesys Time system (one of "TT", "TAI", "UTC")
438  * @return Julian epoch (years).
439  ***************************************************************************/
440 double GTime::julian_epoch(const std::string& timesys) const
441 {
442  // Compute Julian epoch
443  double julian_epoch = 2000.0 + (jd(timesys) - 2451545.0) / 365.25;
444 
445  // Return Julian epoch
446  return (julian_epoch);
447 }
448 
449 
450 /***********************************************************************//**
451  * @brief Return time as string in UTC time system
452  *
453  * @param[in] precision Digits of precision to show in the seconds field
454  * @return Time as string in UTC time system.
455  *
456  * Returns time in the format YYYY-MM-DDThh:mm:ss, where YYYY is a four-digit
457  * year, MM a two-digit month, DD a two-digit day of month, hh two digits of
458  * hour (0 through 23), mm two digits of minutes, and ss two digits of
459  * second (ISO 8601 time standard). In case that @p precision > 0 digits
460  * in the second after the comma will be returned.
461  *
462  * The method is only valid for dates from year 1972 on.
463  ***************************************************************************/
464 std::string GTime::utc(const int& precision) const
465 {
466  // Check argument
467  if (precision < 0) {
468  std::string msg = "Invalid precision \""+gammalib::str(precision)+"\""
469  " encountered. Please specify a precision >= 0.";
471  }
472 
473  // Define number of days per month
474  static int daymonth[12] = {31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31};
475 
476  // Get MJD in TAI time system
477  double mjd = this->mjd() - gammalib::tai2tt * gammalib::sec2day;
478 
479  // Correct for leap seconds (repeat is a kluge to converge as the
480  // argument is given in the UTC system but upon start we're in the TAI
481  // system
482  double leaps = leap_seconds(mjd) * gammalib::sec2day;
483  leaps = leap_seconds(mjd - leaps) * gammalib::sec2day;
484  leaps = leap_seconds(mjd - leaps) * gammalib::sec2day;
485  mjd -= leaps;
486 
487  // Split in day and fraction
488  int day = (int)mjd;
489  double fraction = mjd - (double)day;
490 
491  // Compute time in day. We add a margin of 0.5 to the seconds and
492  // subtract it later to avoid rounding of 59 to 60.
493  double second = fraction * gammalib::sec_in_day + 0.5;
494  int hour = (int)second / 3600;
495  second -= hour * 3600.0;
496  int minute = (int)second / 60;
497  second -= minute * 60.0;
498  if (hour > 23) {
499  hour -= 24;
500  day++;
501  }
502  second -= 0.5;
503  if (second < 0.0) {
504  second = 0.0;
505  }
506 
507  // Compute year and day in the year
508  int year = 1972; // Set year to 1972
509  day -= 41317; // Subtract MJD of 1-1-1972
510  day++; // Day 0 is 1 January
511  int days = days_in_year(year);
512  while (day > days) {
513  day -= days;
514  year++;
515  days = days_in_year(year);
516  }
517 
518  // Adjust number of days per month for leap years
519  if (is_leap_year(year)) {
520  daymonth[1] = 29;
521  }
522  else {
523  daymonth[1] = 28;
524  }
525 
526  // Compute month and day in the month
527  int month = 0;
528  while (month < 12) {
529  if (day <= daymonth[month]) {
530  break;
531  }
532  day -= daymonth[month];
533  month++;
534  }
535  month++;
536 
537  // Calculate the width of the seconds' pattern, plus two for the integer
538  // seconds, plus one if a decimal is needed
539  int sec_width = precision + 2;
540  if (precision > 0) {
541  sec_width += 1;
542  }
543 
544  // Format pattern for variable seconds precision
545  char utc_pattern[50];
546  sprintf(utc_pattern, "%%4.4d-%%2.2d-%%2.2dT%%2.2d:%%2.2d:%%0%d.0%df",
547  sec_width, precision);
548 
549  // Create string
550  char utc[32];
551  sprintf(utc, utc_pattern,
552  year, month, day, hour, minute, second);
553 
554  // Return
555  return (std::string(utc));
556 }
557 
558 
559 /***********************************************************************//**
560  * @brief Return Greenwich mean sidereal time in hours in a day
561  *
562  * @return Greenwich mean sidereal time (hours).
563  *
564  * See http://aa.usno.navy.mil/faq/docs/GAST.php
565  ***************************************************************************/
566 double GTime::gmst(void) const
567 {
568  // Get days since 1 January 2000, 12h in Universal Time
569  double d = days("UTC") + 3652.50076602;
570 
571  // Compute Greenwich mean sidereal time in hours
572  double gmst = 18.697374558 + 24.06570982441908 * d;
573 
574  // Put into [0,24]
575  gmst -= floor(gmst/24.0) * 24.0;
576 
577  // Return Greenwich mean sidereal time in hours
578  return gmst;
579 }
580 
581 
582 /***********************************************************************//**
583  * @brief Return Greenwich apparent sidereal time in hours in a day
584  *
585  * @return Greenwich apparent sidereal time (hours).
586  *
587  * See http://aa.usno.navy.mil/faq/docs/GAST.php
588  ***************************************************************************/
589 double GTime::gast(void) const
590 {
591  // Get days since 1 January 2000, 12h in Universal Time
592  double d = days("UTC") + 3652.50076602;
593 
594  // Compute longitude of the ascending node of the Moon in degrees
595  double Omega = 125.04 - 0.052954 * d;
596 
597  // Compute mean longitude of the Sun in degrees
598  double L = 280.47 + 0.98565 * d;
599 
600  // Compute nutation in longitude in hours
601  double DeltaPsi = -0.000319 * std::sin(Omega * gammalib::deg2rad) -
602  0.000024 * std::sin(2.0 * L * gammalib::deg2rad);
603 
604  // Compute the obliquity in degrees
605  double epsilon = 23.4393 - 0.0000004 * d;
606 
607  // Compute equation of equinoxes
608  double eqeq = DeltaPsi * std::cos(epsilon * gammalib::deg2rad);
609 
610  // Compute Greenwich apparent sidereal time in hours
611  double gast = gmst() + eqeq;
612 
613  // Put into [0,24]
614  gast -= floor(gast/24.0) * 24.0;
615 
616  // Return Greenwich apparent sidereal time in hours
617  return gast;
618 }
619 
620 
621 /***********************************************************************//**
622  * @brief Return local mean sidereal time in hours in a day
623  *
624  * @param[in] geolon Geographic longitude West of Greenwich (degrees).
625  * @return Local mean sidereal time (hours).
626  *
627  * See http://aa.usno.navy.mil/faq/docs/GAST.php
628  ***************************************************************************/
629 double GTime::lmst(const double& geolon) const
630 {
631  // Compute local mean siderial time
632  double lmst = gmst() - geolon/15.0;
633 
634  // Put into [0,24]
635  lmst -= floor(lmst/24.0) * 24.0;
636 
637  // Return local mean sidereal time in hours
638  return lmst;
639 }
640 
641 
642 /***********************************************************************//**
643  * @brief Return local apparent sidereal time in hours in a day
644  *
645  * @param[in] geolon Geographic longitude West of Greenwich (degrees).
646  * @return Local apparent sidereal time (hours).
647  *
648  * See http://aa.usno.navy.mil/faq/docs/GAST.php
649  ***************************************************************************/
650 double GTime::last(const double& geolon) const
651 {
652  // Compute local mean siderial time
653  double last = gast() - geolon/15.0;
654 
655  // Put into [0,24]
656  last -= floor(last/24.0) * 24.0;
657 
658  // Return local mean sidereal time in hours
659  return last;
660 }
661 
662 
663 /***********************************************************************//**
664  * @brief Return time in specified reference
665  *
666  * @return Time in specified reference.
667  *
668  * Convert the time from the native reference system into the specified
669  * reference system.
670  ***************************************************************************/
671 double GTime::convert(const GTimeReference& ref) const
672 {
673  // Retrieve time in native reference (TT in seconds)
674  double time = m_time;
675 
676  // Compute time offset in seconds
677  double offset = (mjd_ref - ref.mjdref()) * gammalib::sec_in_day;
678 
679  // Add time offset in seconds
680  time += offset;
681 
682  // Subtract leap seconds in case that time is requested in UTC time
683  // system
684  if (gammalib::toupper(ref.timesys()) == "UTC") {
685 
686  // Get MJD in TAI time system
687  double mjd = this->mjd() - gammalib::tai2tt * gammalib::sec2day;
688 
689  // Get leap seconds (repeat is a kluge to converge as the
690  // argument is given in the UTC system but upon start we're in the TAI
691  // system
692  double leaps = leap_seconds(mjd);
693  leaps = leap_seconds(mjd - leaps * gammalib::sec2day);
694  leaps = leap_seconds(mjd - leaps * gammalib::sec2day);
695 
696  // Subtract leap seconds and TAI offset
697  time -= (leaps + gammalib::tai2tt);
698 
699  } // endif: UTC time system has been requested
700 
701  // ... otherwise, if time is requested in TAI system then subtract the
702  // TAI offset
703  else if (gammalib::toupper(ref.timesys()) == "TAI") {
704  time -= gammalib::tai2tt;
705  }
706 
707  // Convert to specified time unit
708  double to_unit = ref.unitseconds();
709  if (to_unit != 1.0) {
710  time /= to_unit;
711  }
712 
713  // Return time
714  return time;
715 }
716 
717 
718 /***********************************************************************//**
719  * @brief Set time in Julian Days in native reference (TT)
720  *
721  * @param[in] time Time in Julian Days (TT) (days).
722  ***************************************************************************/
723 void GTime::jd(const double& time)
724 {
725  // Convert time from Julian Days to seconds in native reference
726  m_time = (time - jd_ref) * gammalib::sec_in_day;
727 
728  // Return
729  return;
730 }
731 
732 
733 /***********************************************************************//**
734  * @brief Set time in Julian Days in native reference for time system
735  *
736  * @param[in] time Time in Julian Days (days).
737  * @param[in] timesys Time system (one of "TT", "TAI", "UTC")
738  ***************************************************************************/
739 void GTime::jd(const double& time, const std::string& timesys)
740 {
741  // Convert time from Julian Days to seconds in native reference
742  double seconds = (time - jd_ref) * gammalib::sec_in_day;
743 
744  // Set time according to the specified time system
745  secs(seconds, timesys);
746 
747  // Return
748  return;
749 }
750 
751 
752 /***********************************************************************//**
753  * @brief Set time in Modified Julian Days in native reference (TT)
754  *
755  * @param[in] time Time in Modified Julian Days (TT) (days).
756  ***************************************************************************/
757 void GTime::mjd(const double& time)
758 {
759  // Convert time from Modified Julian Days to native (seconds)
760  m_time = (time - mjd_ref) * gammalib::sec_in_day;
761 
762  // Return
763  return;
764 }
765 
766 
767 /***********************************************************************//**
768  * @brief Set time in Modified Julian Days in native reference for time
769  * system
770  *
771  * @param[in] time Time in Modified Julian Days (days).
772  * @param[in] timesys Time system (one of "TT", "TAI", "UTC")
773  ***************************************************************************/
774 void GTime::mjd(const double& time, const std::string& timesys)
775 {
776  // Convert time from Modified Julian Days to seconds in native reference
777  double seconds = (time - mjd_ref) * gammalib::sec_in_day;
778 
779  // Set time according to the specified time system
780  secs(seconds, timesys);
781 
782  // Return
783  return;
784 }
785 
786 
787 /***********************************************************************//**
788  * @brief Set time in seconds in native reference for time system
789  *
790  * @param[in] seconds Time in native reference (seconds).
791  * @param[in] timesys Time system (one of "TT", "TAI", "UTC")
792  ***************************************************************************/
793 void GTime::secs(const double& seconds, const std::string& timesys)
794 {
795  // If system is TT then simply set time ...
796  if (timesys == "TT") {
797  m_time = seconds;
798  }
799 
800  // ... otherwise if the system if TAI then add an offset ...
801  else if (timesys == "TAI") {
802  m_time = seconds + gammalib::tai2tt;
803  }
804 
805  // ... otherwise if the system is UTC then convert the time from TT to TAI
806  // and add the leap seconds ...
807  else if (timesys == "UTC") {
808  double mjd = seconds * gammalib::sec2day + mjd_ref; // Time in UTC
809  double leaps = leap_seconds(mjd);
810  m_time = seconds + gammalib::tai2tt + leaps;
811  }
812 
813  // ... otherwise throw an exception
814  else {
815  std::string msg = "Unknown time system \""+timesys+"\". Either specify "
816  "\"TT\", \"TAI\" or \"UTC\".";
818  }
819 
820  // Return
821  return;
822 }
823 
824 
825 /***********************************************************************//**
826  * @brief Set time in days in native reference (TT)
827  *
828  * @param[in] days Time (TT) (days).
829  ***************************************************************************/
830 void GTime::days(const double& days)
831 {
832  // Set time
833  m_time = days * gammalib::sec_in_day;
834 
835  // Return
836  return;
837 }
838 
839 
840 /***********************************************************************//**
841  * @brief Set time in days in native reference for time system
842  *
843  * @param[in] days Time (TT) (days).
844  * @param[in] timesys Time system (one of "TT", "TAI", "UTC")
845  ***************************************************************************/
846 void GTime::days(const double& days, const std::string& timesys)
847 {
848  // Convert time from days to seconds in native reference
849  double seconds = days * gammalib::sec_in_day;
850 
851  // Set time according to the specified time system
852  secs(seconds, timesys);
853 
854  // Return
855  return;
856 }
857 
858 
859 /***********************************************************************//**
860  * @brief Set time as string in UTC time system
861  *
862  * @param[in] time Time string (UTC).
863  *
864  * @exception GException::invalid_argument
865  * Invalid time string specified.
866  *
867  * The time has to be given in the format YYYY-MM-DDThh:mm:ss.s, where
868  * YYYY is a four-digit year, MM a two-digit month, DD a two-digit day of
869  * month, hh two digits of hour (0 through 23), mm two digits of minutes,
870  * ss two digits of second and s one or more digits representing a
871  * decimal fraction of a second (ISO 8601 time standard).
872  *
873  * The method is only valid for dates from year 1972 on.
874  ***************************************************************************/
875 void GTime::utc(const std::string& time)
876 {
877  // Define number of days per month
878  static int daymonth[12] = {31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31};
879 
880  // Analyse the UTC string
881  long year = 0;
882  int month = 0;
883  long day = 0;
884  long hour = 0;
885  long minute = 0;
886  double second = 0.0;
887  int n = sscanf(time.c_str(), "%ld-%d-%ldT%ld:%ld:%lg",
888  &year, &month, &day, &hour, &minute, &second);
889 
890  // Adjust number of days per month for leap years
891  if (is_leap_year(year)) {
892  daymonth[1] = 29;
893  }
894  else {
895  daymonth[1] = 28;
896  }
897 
898  // Check time string
899  if (n != 3 && n != 6) {
900  std::string msg = "Invalid time string \""+time+"\" encountered. "
901  "Please specify the time in the format YYYY-MM-DD "
902  "or YYYY-MM-DDThh:mm:ss.s.";
904  }
905  if (year < 1000 || year > 9999) {
906  std::string msg = "Invalid year "+gammalib::str(year)+" specified. "
907  "The year needs to be a four-digit year.";
909  }
910  if (month < 1 || month > 12) {
911  std::string msg = "Invalid month "+gammalib::str(month)+" specified. "
912  "The month needs to be comprised between 01 and 12";
914  }
915  if (day < 1 || day > daymonth[month-1]) {
916  std::string msg = "Invalid day "+gammalib::str(day)+" specified. "
917  "The day for month "+gammalib::str(month)+" needs "
918  "to be comprised between 01 and "+
919  gammalib::str(daymonth[month-1])+".";
921  }
922  if (hour < 0 || hour > 23) {
923  std::string msg = "Invalid hour "+gammalib::str(hour)+" specified. "
924  "The hour needs to be comprised between 00 and 23";
926  }
927  if (minute < 0 || minute > 59) {
928  std::string msg = "Invalid minute "+gammalib::str(minute)+" specified. "
929  "The minute needs to be comprised between 00 and 59";
931  }
932  if (second < 0 || second >= 60.0) {
933  std::string msg = "Invalid second "+gammalib::str(second)+" specified. "
934  "The second needs to be comprised between 00 and <60";
936  }
937 
938  // Compute MJD day
939  month--;
940  for (int i = 0; i < month; ++i) { // Add days passed per month
941  day += daymonth[i];
942  }
943  day += (year - 1972) * 365 - 1; // Add days passed per year
944  day += (year - 1969) / 4; // Add leap days passed (every 4 years)
945  day -= (year - 1901) / 100; // Add leap days passed (not every 100 years)
946  day += (year - 1601) / 400; // Add leap days passed (every 400 years)
947  day += 41317; // Add MJD at 1972
948 
949  // Compute MJD fraction (UTC)
950  double fraction = ((double)hour * 3600.0 + (double)minute * 60.0 + second) *
952 
953  // Compute MJD (UTC)
954  double mjd = (double)day + fraction;
955 
956  // Get conversion from UTC to TAI to TT
957  double correction = leap_seconds(mjd) + gammalib::tai2tt;
958 
959  // Convert MJD from UTC to TT
960  mjd += correction * gammalib::sec2day;
961 
962  // Set time
963  this->mjd(mjd);
964 
965  // Return
966  return;
967 }
968 
969 
970 /***********************************************************************//**
971  * @brief Set time given in specified reference
972  *
973  * @param[in] time Time in given reference system.
974  * @param[in] ref Reference system.
975  *
976  * Set the time to a value given in a specific reference system.
977  ***************************************************************************/
978 void GTime::set(const double& time, const GTimeReference& ref)
979 {
980  // Convert time to specified time unit
981  m_time = time * ref.unitseconds();
982 
983  // Compute time offset in seconds
984  double offset = (mjd_ref - ref.mjdref()) * gammalib::sec_in_day;
985 
986  // Subtract time offset in seconds
987  m_time -= offset;
988 
989  // Add leap seconds in case that time was given in UTC time system
990  if (gammalib::toupper(ref.timesys()) == "UTC") {
991 
992  // Add leap seconds and offset between TAI and TT time system
993  m_time += leap_seconds(this->mjd());
995 
996  } // endif: Time was given in UTC time system
997 
998  // ... otherwise if system is TAI system then add TAI offset
999  else if (gammalib::toupper(ref.timesys()) == "TAI") {
1001  }
1002 
1003  // Return
1004  return;
1005 }
1006 
1007 
1008 /***********************************************************************//**
1009  * @brief Set time from string
1010  *
1011  * @param[in] time Time string.
1012  * @param[in] ref Reference system.
1013  *
1014  * Sets the time from a string for a given reference system. The following
1015  * strings are valid:
1016  *
1017  * "2016-10-05T15:08:56" (UTC string)
1018  * "1800.0" (MET seconds in specified reference system)
1019  * "1800.0 (TT)" (MET seconds in specified reference, TT system)
1020  * "1800.0 (UTC)" (MET seconds in specified reference, UTC time system)
1021  * "1800.0 (TAI)" (MET seconds in specified reference, TAI time system)
1022  * "MJD 54609" (Modified Julian Days, TT system)
1023  * "MJD 54609 (TT)" (Modified Julian Days, TT system)
1024  * "MJD 54609 (UTC)" (Modified Julian Days, UTC system)
1025  * "MJD 54609 (TAI)" (Modified Julian Days, TAI system)
1026  * "JD 54609" (Julian Days, TT system)
1027  * "JD 54609 (TT)" (Julian Days, TT system)
1028  * "JD 54609 (UTC)" (Julian Days, UTC system)
1029  * "JD 54609 (TAI)" (Julian Days, TAI system)
1030  *
1031  * If any other string is encountered, the numerical value is interpreted as
1032  * time is seconds using the TT system.
1033  *
1034  * Note that the TT, UTC or TAI attributes overwrite the values contained in
1035  * the specified reference system. The reference system is only used to
1036  * convert MET times in seconds.
1037  ***************************************************************************/
1038 void GTime::set(const std::string& time, const GTimeReference& ref)
1039 {
1040  // Strip any whitespace from string and convert it to upper case
1041  std::string str = gammalib::toupper(gammalib::strip_whitespace(time));
1042 
1043  // First check if the string is a UTC string
1044  long year = 0;
1045  int month = 0;
1046  long day = 0;
1047  long hour = 0;
1048  long minute = 0;
1049  double second = 0.0;
1050  int n = sscanf(str.c_str(), "%ld-%d-%ldT%ld:%ld:%lg",
1051  &year, &month, &day, &hour, &minute, &second);
1052  if (n == 3 || n == 6) {
1053  utc(str);
1054  }
1055 
1056  // ... otherwise check for Modified Julian Days
1057  else if (str.find("MJD") == 0) {
1058  double timeval = extract_timeval(str);
1059  std::string timesys = extract_timesys(str);
1060  if (timesys.empty()) {
1061  timesys = "TT";
1062  }
1063  mjd(timeval, timesys);
1064  }
1065 
1066  // ... otherwise check for Julian Days
1067  else if (str.find("JD") == 0) {
1068  double timeval = extract_timeval(str);
1069  std::string timesys = extract_timesys(str);
1070  if (timesys.empty()) {
1071  timesys = "TT";
1072  }
1073  jd(timeval, timesys);
1074  }
1075 
1076  // ... otherwise take time as seconds and use the specified reference
1077  // system. If no TT, UTC or TAI arrtibutes are specified use the value
1078  // specified by the reference system.
1079  else {
1080  double timeval = extract_timeval(str);
1081  std::string timesys = extract_timesys(str);
1082  if (timesys.empty()) {
1083  timesys = ref.timesys();
1084  }
1085  GTimeReference timeref(ref.mjdref(), ref.timeunit(), timesys, ref.timeref());
1086  set(timeval, timeref);
1087  }
1088 
1089  // Return
1090  return;
1091 }
1092 
1093 
1094 /***********************************************************************//**
1095  * @brief Set time to current time
1096  *
1097  * Sets time to current time.
1098  ***************************************************************************/
1099 void GTime::now(void)
1100 {
1101  // Allocate variables
1102  struct std::tm timeStruct;
1103  std::time_t now;
1104  char buffer[100];
1105 
1106  // Get time
1107  now = std::time(NULL);
1108  #ifdef HAVE_GMTIME_R
1109  std::gmtime_r(&now, &timeStruct);
1110  #else
1111  std::memcpy(&timeStruct, gmtime(&now), sizeof(struct tm));
1112  #endif
1113 
1114  // Write message type, time and task name to buffer
1115  std::sprintf(buffer, "%04d-%02d-%02dT%02d:%02d:%02d",
1116  timeStruct.tm_year + 1900,
1117  timeStruct.tm_mon + 1,
1118  timeStruct.tm_mday,
1119  timeStruct.tm_hour,
1120  timeStruct.tm_min,
1121  timeStruct.tm_sec);
1122 
1123  // Build string from buffer
1124  std::string date = buffer;
1125 
1126  // Set UTC time
1127  utc(date);
1128 
1129  // Return
1130  return;
1131 }
1132 
1133 
1134 /***********************************************************************//**
1135  * @brief Returns native time reference
1136  *
1137  * @return Native time reference.
1138  *
1139  * Returns the native GammaLib time reference. The GammaLib native time
1140  * reference (i.e. time=0) is defined as January 1, 2010, 00:00:00 (TT).
1141  * The time system is Terrestrial Time (TT). Time is stored in seconds.
1142  ***************************************************************************/
1144 {
1145  // Allocate native time reference
1146  GTimeReference reference(mjd_ref, "s", "TT", "LOCAL");
1147 
1148  // Return reference
1149  return reference;
1150 }
1151 
1152 
1153 /***********************************************************************//**
1154  * @brief Print time
1155  *
1156  * @param[in] chatter Chattiness.
1157  * @return String containing time in seconds in native reference.
1158  *
1159  * Prints time in seconds in the native reference.
1160  ***************************************************************************/
1161 std::string GTime::print(const GChatter& chatter) const
1162 {
1163  // Initialise result string
1164  std::string result;
1165 
1166  // Continue only if chatter is not silent
1167  if (chatter != SILENT) {
1168 
1169  // Append time
1170  result.append(gammalib::str(m_time)+" s (TT)");
1171 
1172  } // endif: chatter was not silent
1173 
1174  // Return
1175  return result;
1176 }
1177 
1178 
1179 /*==========================================================================
1180  = =
1181  = Private methods =
1182  = =
1183  ==========================================================================*/
1184 
1185 /***********************************************************************//**
1186  * @brief Initialise class members
1187  ***************************************************************************/
1189 {
1190  // Initialise members
1191  m_time = 0.0;
1192 
1193  // Return
1194  return;
1195 }
1196 
1197 
1198 /***********************************************************************//**
1199  * @brief Copy class members
1200  *
1201  * @param[in] time Time.
1202  ***************************************************************************/
1203 void GTime::copy_members(const GTime& time)
1204 {
1205  // Copy time
1206  m_time = time.m_time;
1207 
1208  // Return
1209  return;
1210 }
1211 
1212 
1213 /***********************************************************************//**
1214  * @brief Delete class members
1215  ***************************************************************************/
1217 {
1218  // Return
1219  return;
1220 }
1221 
1222 
1223 /***********************************************************************//**
1224  * @brief Returns number of leap seconds for a given MJD
1225  *
1226  * @param[in] mjd Modified Julian Day in UTC time system.
1227  * @return Number of lead seconds.
1228  *
1229  * Return the number of leap seconds for a given MJD specified in the UTC
1230  * time system. This method returns valid number of leap seconds for the
1231  * years 1972-2017.
1232  *
1233  * See http://www.nist.gov/pml/div688/grp50/leapsecond.cfm for a table of
1234  * leap seconds.
1235  ***************************************************************************/
1236 double GTime::leap_seconds(const double& mjd) const
1237 {
1238  // Leap second table from 1972 on
1239  // see http://www.nist.gov/pml/div688/grp50/leapsecond.cfm
1240  // The first entry is 1-Jan-1972
1241  const long leapsmjd[] = {41317, // 1972-01-01
1242  41498, // 1972-06-30
1243  41682, // 1972-12-31
1244  42047, // 1973-12-31
1245  42412, // 1974-12-31
1246  42777, // 1975-12-31
1247  43143, // 1976-12-31
1248  43508, // 1977-12-31
1249  43873, // 1978-12-31
1250  44238, // 1979-12-31
1251  44785, // 1981-06-30
1252  45150, // 1982-06-30
1253  45515, // 1983-06-30
1254  46246, // 1985-06-30
1255  47160, // 1987-12-31
1256  47891, // 1989-12-31
1257  48256, // 1990-12-31
1258  48803, // 1992-06-30
1259  49168, // 1993-06-30
1260  49533, // 1994-06-30
1261  50082, // 1995-12-31
1262  50629, // 1997-06-30
1263  51178, // 1998-12-31
1264  53735, // 2005-12-31
1265  54831, // 2008-12-31
1266  56108, // 2012-06-30
1267  57203, // 2015-06-30
1268  57753}; // 2016-12-31
1269  const double leapsecs[] = {10.0, 11.0, 12.0, 13.0, 14.0,
1270  15.0, 16.0, 17.0, 18.0, 19.0,
1271  20.0, 21.0, 22.0, 23.0, 24.0,
1272  25.0, 26.0, 27.0, 28.0, 29.0,
1273  30.0, 31.0, 32.0, 33.0, 34.0,
1274  35.0, 36.0, 37.0};
1275  const int n_leapsecs = sizeof(leapsmjd)/sizeof(long);
1276 
1277  // Extract MJD day and MJD fraction
1278  long day = (long)mjd;
1279 
1280  // Find the leap second MJD that is equal or larger to the specified
1281  // day
1282  int i = n_leapsecs - 1;
1283  while ((day <= leapsmjd[i]) && i > 0) {
1284  i--;
1285  }
1286 
1287  // Return leap seconds
1288  return (leapsecs[i]);
1289 }
1290 
1291 
1292 /***********************************************************************//**
1293  * @brief Extract time value from time string
1294  *
1295  * @param[in] time Time string.
1296  * @return Time value.
1297  *
1298  * Extracts the time value from a time string. The method strips any prefix
1299  * such as "MJD" or "JD" and any suffix starting with a left parentheses "("
1300  * and converts the remainder into a double precision value.
1301  ***************************************************************************/
1302 double GTime::extract_timeval(const std::string& time) const
1303 {
1304  // Find time
1305  size_t length = time.length();
1306  size_t start = time.find_first_of("0123456789+-.");
1307  size_t stop = time.find("(");
1308 
1309  // If there is a stop position then set length to stop
1310  if (stop != std::string::npos) {
1311  length = stop;
1312  }
1313 
1314  // Reduce length by start
1315  if (start != std::string::npos) {
1316  length -= start;
1317  }
1318 
1319  // Get substring containing the time value
1320  std::string value = time.substr(start, length);
1321 
1322  // Convert value into double precision
1323  double result = gammalib::todouble(value);
1324 
1325  // Return result
1326  return result;
1327 }
1328 
1329 
1330 /***********************************************************************//**
1331  * @brief Extract time system from time string
1332  *
1333  * @param[in] time Time string.
1334  * @return Time system.
1335  *
1336  * Extracts the time system from a time string. Valid time systems are:
1337  *
1338  * "(TT)" (TT system)
1339  * "(UTC)" (UTC system)
1340  * "(TAI)" (TAI system)
1341  *
1342  * If no time system is found a blank string is returned.
1343  ***************************************************************************/
1344 std::string GTime::extract_timesys(const std::string& time) const
1345 {
1346  // Initialise time system with blank string
1347  std::string timesys = "";
1348 
1349  // Strip any whitespace from string and convert it to upper case
1350  std::string str = gammalib::toupper(gammalib::strip_whitespace(time));
1351 
1352  // Check for UTC
1353  if (str.find("(UTC)") != std::string::npos) {
1354  timesys = "UTC";
1355  }
1356  else if (str.find("(TAI)") != std::string::npos) {
1357  timesys = "TAI";
1358  }
1359  else if (str.find("(TT)") != std::string::npos) {
1360  timesys = "TT";
1361  }
1362 
1363  // Return time system
1364  return timesys;
1365 }
std::string extract_timesys(const std::string &time) const
Extract time system from time string.
Definition: GTime.cpp:1344
#define G_SECS_GET
Definition: GTime.cpp:46
void init_members(void)
Initialise class members.
Definition: GTime.cpp:1188
void now(void)
Set time to current time.
Definition: GTime.cpp:1099
GTimeReference reference(void) const
Returns native time reference.
Definition: GTime.cpp:1143
const std::string & timeunit(void) const
Return time unit.
const double mjd_ref
MJD of time=0.
Definition: GTime.cpp:41
double lmst(const double &geolon) const
Return local mean sidereal time in hours in a day.
Definition: GTime.cpp:629
void copy_members(const GTime &time)
Copy class members.
Definition: GTime.cpp:1203
GVector cos(const GVector &vector)
Computes cosine of vector elements.
Definition: GVector.cpp:1100
double unitseconds(void) const
Return the time unit in seconds.
void clear(void)
Clear time.
Definition: GTime.cpp:251
const double sec2day
Definition: GTools.hpp:50
Time class.
Definition: GTime.hpp:54
Gammalib tools definition.
std::string strip_whitespace(const std::string &arg)
Strip leading and trailing whitespace from string.
Definition: GTools.cpp:73
#define G_SECS_SET
Definition: GTime.cpp:47
double extract_timeval(const std::string &time) const
Extract time value from time string.
Definition: GTime.cpp:1302
double jd(void) const
Return time in Julian Days (TT)
Definition: GTime.cpp:283
bool is_leap_year(const int &year) const
Signals if year is a leap year.
Definition: GTime.hpp:179
#define G_UTC_GET
Definition: GTime.cpp:49
const double sec_in_day
Definition: GTools.hpp:49
virtual ~GTime(void)
Destructor.
Definition: GTime.cpp:199
const double deg2rad
Definition: GMath.hpp:43
const double & mjdref(void) const
Return MJD reference (units: days)
const double jd_ref
JD of time=0.
Definition: GTime.cpp:42
void free_members(void)
Delete class members.
Definition: GTime.cpp:1216
GTime(void)
Void constructor.
Definition: GTime.cpp:67
const double tai2tt
Definition: GTools.hpp:51
Time reference class interface definition.
const std::string & timesys(void) const
Return time system.
GChatter
Definition: GTypemaps.hpp:33
std::string print(const GChatter &chatter=NORMAL) const
Print time.
Definition: GTime.cpp:1161
double leap_seconds(const double &mjd) const
Returns number of leap seconds for a given MJD.
Definition: GTime.cpp:1236
const std::string & timeref(void) const
Return time reference.
GTime * clone(void) const
Clone time.
Definition: GTime.cpp:269
double gmst(void) const
Return Greenwich mean sidereal time in hours in a day.
Definition: GTime.cpp:566
const double & secs(void) const
Return time in seconds in native reference (TT)
Definition: GTime.hpp:153
double last(const double &geolon) const
Return local apparent sidereal time in hours in a day.
Definition: GTime.cpp:650
GTime & operator=(const GTime &time)
Assignment operator.
Definition: GTime.cpp:221
double gast(void) const
Return Greenwich apparent sidereal time in hours in a day.
Definition: GTime.cpp:589
GVector sin(const GVector &vector)
Computes sine of vector elements.
Definition: GVector.cpp:1226
Exception handler interface definition.
std::string toupper(const std::string &s)
Convert string to upper case.
Definition: GTools.cpp:820
Implements a time reference.
std::string tolower(const std::string &s)
Convert string to lower case.
Definition: GTools.cpp:834
double mjd(void) const
Return time in Modified Julian Days (TT)
Definition: GTime.cpp:319
#define G_CONSTRUCT
Definition: GTime.cpp:45
int days_in_year(const int &year) const
Returns number of days in year.
Definition: GTime.hpp:194
#define G_UTC
Definition: GTime.cpp:48
double julian_epoch(void) const
Return Julian epoch in native reference (TT)
Definition: GTime.cpp:424
std::string utc(const int &precision=0) const
Return time as string in UTC time system.
Definition: GTime.cpp:464
Time class interface definition.
void set(const double &time, const GTimeReference &ref)
Set time given in specified reference.
Definition: GTime.cpp:978
double convert(const GTimeReference &ref) const
Return time in specified reference.
Definition: GTime.cpp:671
double days(void) const
Return time in days in native reference (TT)
Definition: GTime.cpp:399
double m_time
Time in seconds in native reference (TT)
Definition: GTime.hpp:131
Mathematical function definitions.
double todouble(const std::string &arg)
Convert string into double precision value.
Definition: GTools.cpp:805
std::string str(const unsigned short int &value)
Convert unsigned short integer value into string.
Definition: GTools.cpp:413