1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
18 ///////////////////////////////////////////////////////////////////////////
20 // Handling of timestamps for (astro)particle physics reserach.
22 // This class is derived from TTimeStamp and provides additional
23 // facilities (e.g. Julian date) which are commonly used in the
24 // field of (astro)particle physics.
26 // The Julian Date (JD) indicates the number of days since noon (UT) on
27 // 01 jan -4712 (i.e. noon 01 jan 4713 BC), being day 0 of the Julian calendar.
29 // The Modified Julian Date (MJD) indicates the number of days since midnight
30 // (UT) on 17-nov-1858, which corresponds to 2400000.5 days after day 0 of the
33 // The Truncated Julian Date (TJD) corresponds to 2440000.5 days after day 0
34 // of the Julian calendar and consequently TJD=MJD-40000.
35 // This TJD date indication was used by the Vela and Batse missions in
36 // view of Gamma Ray Burst investigations.
38 // The Julian Epoch (JE) indicates the fractional elapsed Julian year count
39 // since the start of the Gregorian year count.
40 // A Julian year is defined to be 365.25 days and starts at 01-jan 12:00:00 UT.
41 // As such, the integer part of JE corresponds to the usual Gregorian year count,
42 // apart from 01-jan before 12:00:00 UT.
43 // So, 01-jan-1965 12:00:00 UT corresponds to JE=1965.0
45 // The Besselian Epoch (BE) indicates the fractional elapsed Besselian year count
46 // since the start of the Gregorian year count.
47 // A Besselian (or tropical) year is defined to be 365.242198781 days.
48 // The date 31-dec-1949 22:09:46.862 UT corresponds to BE=1950.0
50 // The Besselian and Julian epochs are used in astronomical catalogs
51 // to denote values of time varying observables like e.g. right ascension.
53 // Because of the fact that the Julian date indicators are all w.r.t. UT
54 // they provide an absolute timescale irrespective of timezone or daylight
57 // In view of astronomical observations and positioning it is convenient
58 // to have also a UT equivalent related to stellar meridian transitions.
59 // This is achieved by the Greenwich Sidereal Time (GST).
60 // The GST is defined as the right ascension of the objects passing
61 // the Greenwich meridian at 00:00:00 UT.
62 // Due to the rotation of the Earth around the Sun, a sidereal day
63 // lasts 86164.09 seconds (23h 56m 04.09s) compared to the mean solar
64 // day of 86400 seconds (24h).
65 // Furthermore, precession of the earth's spin axis results in the fact
66 // that the zero point of right ascension (vernal equinox) gradually
67 // moves along the celestial equator.
68 // In addition, tidal friction and ocean and atmospheric effects will
69 // induce seasonal variations in the earth's spin rate and polar motion
70 // of the earth's spin axis.
71 // Taking the above effects into account leads to what is called
72 // the Greenwich Mean Sidereal Time (GMST).
73 // In case also the nutation of the earth's spin axis is taken into
74 // account we speak of the Greenwich Apparent Sidereal Time (GAST).
76 // This AliTimestamp facility allows for picosecond precision, in view
77 // of time of flight analyses for particle physics experiments.
78 // For normal date/time indication the standard nanosecond precision
79 // will in general be sufficient.
80 // Note that when the fractional JD, MJD and TJD counts are used instead
81 // of the integer (days,sec,ns) specification, the nanosecond precision
82 // may be lost due to computer accuracy w.r.t. floating point operations.
84 // The TTimeStamp EPOCH starts at 01-jan-1970 00:00:00 UTC
85 // which corresponds to JD=2440587.5 or the start of MJD=40587 or TJD=587.
86 // Using the corresponding MJD of this EPOCH allows construction of
87 // the yy-mm-dd hh:mm:ss:ns TTimeStamp from a given input (M/T)JD and time.
88 // Obviously this TTimeStamp implementation would prevent usage of values
89 // smaller than JD=2440587.5 or MJD=40587 or TJD=587.
90 // Furthermore, due to a limitation on the "seconds since the EPOCH start" count
91 // in TTimeStamp, the latest accessible date/time is 19-jan-2038 02:14:08 UTC.
92 // However, this AliTimestamp facility provides support for the full range
93 // of (M/T)JD values, but the setting of the corresponding TTimeStamp parameters
94 // is restricted to the values allowed by the TTimeStamp implementation.
95 // For these earlier/later (M/T)JD values, the standard TTimeStamp parameters will
96 // be set corresponding to the start of the TTimeStamp EPOCH.
97 // This implies that for these earlier/later (M/T)JD values the TTimeStamp parameters
98 // do not match the Julian parameters of AliTimestamp.
99 // As such the standard TTimeStamp parameters do not appear on the print output
100 // when invoking the Date() memberfunction for these earlier/later (M/T)JD values.
105 // Note : All TTimeStamp functionality is available as well.
111 // // Retrieve Julian Date
112 // Int_t jd,jsec,jns;
113 // t.GetJD(jd,jsec,jns);
115 // // Retrieve fractional Truncated Julian Date
116 // Double_t tjd=t.GetTJD();
118 // // Retrieve fractional Julian Epoch
119 // Double_t je=t.GetJE();
121 // // Set to a specific Modified Julian Date
124 // Int_t mjns=185643;
125 // t.SetMJD(mjd,mjsec,mjns);
129 // // Time intervals for e.g. trigger or TOF analysis
131 // AliTrack* tx=evt.GetTrack(5);
132 // AliTimestamp* timex=tx->GetTimestamp();
133 // Double_t dt=evt.GetDifference(timex,"ps");
134 // AliTimestamp trig((AliTimestamp)evt);
135 // trig.Add(0,0,2,173);
136 // AliSignal* sx=evt.GetHit(23);
137 // AliTimestamp* timex=sx->GetTimestamp();
138 // Double_t dt=trig.GetDifference(timex,"ps");
140 // trig.GetDifference(timex,d,s,ns,ps);
142 // // Some practical conversion facilities
143 // // Note : They don't influence the actual date/time settings
144 // // and as such can also be invoked as AliTimestamp::Convert(...) etc...
152 // Double_t jdate=t.GetJD(y,m,d,hh,mm,ss,ns);
154 // Int_t days,secs,nsecs;
155 // Double_t date=421.1949327;
156 // t.Convert(date,days,secs,nsecs);
161 // date=t.Convert(days,secs,nsecs);
163 // Double_t mjdate=40563.823744;
164 // Double_t epoch=t.GetJE(mjdate,"mjd");
166 //--- Author: Nick van Eijndhoven 28-jan-2005 Utrecht University.
167 //- Modified: NvE $Date$ Utrecht University.
168 ///////////////////////////////////////////////////////////////////////////
170 #include "AliTimestamp.h"
171 #include "Riostream.h"
173 ClassImp(AliTimestamp) // Class implementation to enable ROOT I/O
175 AliTimestamp::AliTimestamp() : TTimeStamp()
177 // Default constructor
178 // Creation of an AliTimestamp object and initialisation of parameters.
179 // All attributes are initialised to the current date/time as specified
180 // in the docs of TTimeStamp.
185 ///////////////////////////////////////////////////////////////////////////
186 AliTimestamp::AliTimestamp(TTimeStamp& t) : TTimeStamp(t)
188 // Creation of an AliTimestamp object and initialisation of parameters.
189 // All attributes are initialised to the values of the input TTimeStamp.
194 ///////////////////////////////////////////////////////////////////////////
195 AliTimestamp::~AliTimestamp()
197 // Destructor to delete dynamically allocated memory.
199 ///////////////////////////////////////////////////////////////////////////
200 AliTimestamp::AliTimestamp(const AliTimestamp& t) : TTimeStamp(t)
211 ///////////////////////////////////////////////////////////////////////////
212 void AliTimestamp::Date(Int_t mode,Double_t offset)
214 // Print date/time info.
216 // mode = 1 ==> Only the UT yy-mm-dd hh:mm:ss.sss and GMST info is printed
217 // 2 ==> Only the Julian parameter info is printed
218 // 3 ==> Both the UT, GMST and Julian parameter info is printed
220 // offset : Local time offset from UT (and also GMST) in fractional hours.
222 // When an offset value is specified, the corresponding local times
223 // LT and LMST are printed as well.
225 // The default values are mode=3 and offset=0.
227 // Note : In case the (M/T)JD falls outside the TTimeStamp range,
228 // the yy-mm-dd info will be omitted.
230 Int_t mjd,mjsec,mjns,mjps;
231 GetMJD(mjd,mjsec,mjns);
234 TString month[12]={"Jan","Feb","Mar","Apr","May","Jun","Jul","Aug","Sep","Oct","Nov","Dec"};
235 TString day[7]={"Mon","Tue","Wed","Thu","Fri","Sat","Sun"};
237 Int_t hh,mm,ss,ns,ps;
239 if (mode==1 || mode==3)
241 if (mjd>=40587 && (mjd<65442 || (mjd==65442 && mjsec<8047)))
243 GetDate(kTRUE,0,&y,&m,&d);
244 wd=GetDayOfWeek(kTRUE,0);
245 cout << " " << day[wd-1].Data() << ", " << setfill('0') << setw(2) << d << " "
246 << setfill(' ') << month[m-1].Data() << " " << y << " ";
252 GetUT(hh,mm,ss,ns,ps);
253 cout << setfill('0') << setw(2) << hh << ":"
254 << setw(2) << mm << ":" << setw(2) << ss << "."
255 << setw(9) << ns << setw(3) << ps << " (UT) ";
256 GetGMST(hh,mm,ss,ns,ps);
257 cout << setfill('0') << setw(2) << hh << ":"
258 << setw(2) << mm << ":" << setw(2) << ss << "."
259 << setw(9) << ns << setw(3) << ps << " (GMST)" << endl;
261 // Local time information
264 // Determine the new date by including the offset
265 AliTimestamp t2(*this);
267 Int_t mjd2,mjsec2,mjns2;
268 t2.GetMJD(mjd2,mjsec2,mjns2);
269 if (mjd2>=40587 && (mjd2<65442 || (mjd2==65442 && mjsec2<8047)))
271 t2.GetDate(kTRUE,0,&y,&m,&d);
272 wd=t2.GetDayOfWeek(kTRUE,0);
273 cout << " " << day[wd-1].Data() << ", " << setfill('0') << setw(2) << d << " "
274 << setfill(' ') << month[m-1].Data() << " " << y << " ";
280 // Determine the local time by including the offset w.r.t. the original timestamp
281 Double_t hlt=GetLT(offset);
282 Double_t hlst=GetLMST(offset);
283 PrintTime(hlt,12); cout << " (LT) "; PrintTime(hlst,12); cout << " (LMST)" << endl;
286 if (mode==2 || mode==3)
290 Int_t tjd,tjsec,tjns;
291 GetTJD(tjd,tjsec,tjns);
292 cout << " Julian Epoch : " << setprecision(25) << GetJE()
293 << " Besselian Epoch : " << setprecision(25) << GetBE() << endl;
294 cout << " JD : " << jd << " sec : " << jsec << " ns : " << jns << " ps : " << fJps
295 << " Fractional : " << setprecision(25) << GetJD() << endl;
296 cout << " MJD : " << mjd << " sec : " << mjsec << " ns : " << mjns << " ps : " << fJps
297 << " Fractional : " << setprecision(25) << GetMJD() << endl;
298 cout << " TJD : " << tjd << " sec : " << tjsec << " ns : " << tjns << " ps : " << fJps
299 << " Fractional : " << setprecision(25) << GetTJD() << endl;
302 ///////////////////////////////////////////////////////////////////////////
303 Double_t AliTimestamp::GetJD(Int_t y,Int_t m,Int_t d,Int_t hh,Int_t mm,Int_t ss,Int_t ns) const
305 // Provide the (fractional) Julian Date (JD) corresponding to the UT date
306 // and time in the Gregorian calendar as specified by the input arguments.
308 // The input arguments represent the following :
309 // y : year in UT (e.g. 1952, 2003 etc...)
310 // m : month in UT (1=jan 2=feb etc...)
311 // d : day in UT (1-31)
312 // hh : elapsed hours in UT (0-23)
313 // mm : elapsed minutes in UT (0-59)
314 // ss : elapsed seconds in UT (0-59)
315 // ns : remaining fractional elapsed second of UT in nanosecond
317 // This algorithm is valid for all AD dates in the Gregorian calendar
318 // following the recipe of R.W. Sinnott Sky & Telescope 82, (aug. 1991) 183.
319 // See also http://scienceworld.wolfram.com/astronomy/JulianDate.html
321 // In case of invalid input, a value of -1 is returned.
325 // This memberfunction only provides the JD corresponding to the
326 // UT input arguments. It does NOT set the corresponding Julian parameters
327 // for the current AliTimestamp instance.
328 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
329 // To set the Julian parameters for the current AliTimestamp instance,
330 // please use the corresponding SET() memberfunctions of either AliTimestamp
333 if (y<0 || m<1 || m>12 || d<1 || d>31) return -1;
334 if (hh<0 || hh>23 || mm<0 || mm>59 || ss<0 || ss>59 || ns<0 || ns>1e9) return -1;
336 // The UT daytime in fractional hours
337 Double_t ut=double(hh)+double(mm)/60.+(double(ss)+double(ns)*1.e-9)/3600.;
341 JD=367*y-int(7*(y+int((m+9)/12))/4)
342 -int(3*(int((y+(m-9)/7)/100)+1)/4)
343 +int(275*m/9)+d+1721028.5+ut/24.;
347 ///////////////////////////////////////////////////////////////////////////
348 Double_t AliTimestamp::GetMJD(Int_t y,Int_t m,Int_t d,Int_t hh,Int_t mm,Int_t ss,Int_t ns) const
350 // Provide the (fractional) Modified Julian Date corresponding to the UT
351 // date and time in the Gregorian calendar as specified by the input arguments.
353 // The input arguments represent the following :
354 // y : year in UT (e.g. 1952, 2003 etc...)
355 // m : month in UT (1=jan 2=feb etc...)
356 // d : day in UT (1-31)
357 // hh : elapsed hours in UT (0-23)
358 // mm : elapsed minutes in UT (0-59)
359 // ss : elapsed seconds in UT (0-59)
360 // ns : remaining fractional elapsed second of UT in nanosecond
362 // This algorithm is valid for all AD dates in the Gregorian calendar
363 // following the recipe of R.W. Sinnott Sky & Telescope 82, (aug. 1991) 183.
364 // See also http://scienceworld.wolfram.com/astronomy/JulianDate.html
366 // In case of invalid input, a value of -1 is returned.
370 // This memberfunction only provides the MJD corresponding to the
371 // UT input arguments. It does NOT set the corresponding Julian parameters
372 // for the current AliTimestamp instance.
373 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
374 // To set the Julian parameters for the current AliTimestamp instance,
375 // please use the corresponding SET() memberfunctions of either AliTimestamp
378 Double_t JD=GetJD(y,m,d,hh,mm,ss,ns);
382 Double_t MJD=JD-2400000.5;
386 ///////////////////////////////////////////////////////////////////////////
387 Double_t AliTimestamp::GetTJD(Int_t y,Int_t m,Int_t d,Int_t hh,Int_t mm,Int_t ss,Int_t ns) const
389 // Provide the (fractional) Truncated Julian Date corresponding to the UT
390 // date and time in the Gregorian calendar as specified by the input arguments.
392 // The input arguments represent the following :
393 // y : year in UT (e.g. 1952, 2003 etc...)
394 // m : month in UT (1=jan 2=feb etc...)
395 // d : day in UT (1-31)
396 // hh : elapsed hours in UT (0-23)
397 // mm : elapsed minutes in UT (0-59)
398 // ss : elapsed seconds in UT (0-59)
399 // ns : remaining fractional elapsed second of UT in nanosecond
401 // This algorithm is valid for all AD dates in the Gregorian calendar
402 // following the recipe of R.W. Sinnott Sky & Telescope 82, (aug. 1991) 183.
403 // See also http://scienceworld.wolfram.com/astronomy/JulianDate.html
405 // In case of invalid input, a value of -1 is returned.
409 // This memberfunction only provides the TJD corresponding to the
410 // UT input arguments. It does NOT set the corresponding Julian parameters
411 // for the current AliTimestamp instance.
412 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
413 // To set the Julian parameters for the current AliTimestamp instance,
414 // please use the corresponding SET() memberfunctions of either AliTimestamp
417 Double_t JD=GetJD(y,m,d,hh,mm,ss,ns);
421 Double_t TJD=JD-2440000.5;
425 ///////////////////////////////////////////////////////////////////////////
426 Double_t AliTimestamp::GetJE(Double_t date,TString mode) const
428 // Provide the Julian Epoch (JE) corresponding to the specified date.
429 // The argument "mode" indicates the type of the argument "date".
431 // Available modes are :
432 // mode = "jd" ==> date represents the Julian Date
433 // = "mjd" ==> date represents the Modified Julian Date
434 // = "tjd" ==> date represents the Truncated Julian Date
436 // The default is mode="jd".
438 // In case of invalid input, a value of -99999 is returned.
442 // This memberfunction only provides the JE corresponding to the
443 // input arguments. It does NOT set the corresponding Julian parameters
444 // for the current AliTimestamp instance.
445 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
446 // To set the Julian parameters for the current AliTimestamp instance,
447 // please use the corresponding SET() memberfunctions of either AliTimestamp
450 if ((mode != "jd") && (mode != "mjd") && (mode != "tjd")) return -99999;
453 if (mode=="mjd") jd=date+2400000.5;
454 if (mode=="tjd") jd=date+2440000.5;
456 Double_t je=2000.+(jd-2451545.)/365.25;
460 ///////////////////////////////////////////////////////////////////////////
461 Double_t AliTimestamp::GetBE(Double_t date,TString mode) const
463 // Provide the Besselian Epoch (JE) corresponding to the specified date.
464 // The argument "mode" indicates the type of the argument "date".
466 // Available modes are :
467 // mode = "jd" ==> date represents the Julian Date
468 // = "mjd" ==> date represents the Modified Julian Date
469 // = "tjd" ==> date represents the Truncated Julian Date
471 // The default is mode="jd".
473 // In case of invalid input, a value of -99999 is returned.
477 // This memberfunction only provides the BE corresponding to the
478 // input arguments. It does NOT set the corresponding Julian parameters
479 // for the current AliTimestamp instance.
480 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
481 // To set the Julian parameters for the current AliTimestamp instance,
482 // please use the corresponding SET() memberfunctions of either AliTimestamp
485 if ((mode != "jd") && (mode != "mjd") && (mode != "tjd")) return -99999;
488 if (mode=="mjd") jd=date+2400000.5;
489 if (mode=="tjd") jd=date+2440000.5;
491 Double_t be=1900.+(jd-2415020.31352)/365.242198781;
495 ///////////////////////////////////////////////////////////////////////////
496 void AliTimestamp::Convert(Double_t date,Int_t& days,Int_t& secs,Int_t& ns) const
498 // Convert date as fractional day count into integer days, secs and ns.
500 // Note : Due to computer accuracy the ns value may become inaccurate.
502 // The arguments represent the following :
503 // date : The input date as fractional day count
504 // days : Number of elapsed days
505 // secs : Remaining number of elapsed seconds
506 // ns : Remaining fractional elapsed second in nanoseconds
510 // This memberfunction only converts the input date into the corresponding
511 // integer parameters. It does NOT set the corresponding Julian parameters
512 // for the current AliTimestamp instance.
513 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
514 // To set the Julian parameters for the current AliTimestamp instance,
515 // please use the corresponding SET() memberfunctions of either AliTimestamp
519 date=date-double(days);
520 Int_t daysecs=24*3600;
521 date=date*double(daysecs);
523 date=date-double(secs);
526 ///////////////////////////////////////////////////////////////////////////
527 Double_t AliTimestamp::Convert(Int_t days,Int_t secs,Int_t ns) const
529 // Convert date in integer days, secs and ns into fractional day count.
531 // Note : Due to computer accuracy the ns precision may be lost.
533 // The input arguments represent the following :
534 // days : Number of elapsed days
535 // secs : Remaining number of elapsed seconds
536 // ns : Remaining fractional elapsed second in nanoseconds
540 // This memberfunction only converts the input integer parameters into the
541 // corresponding fractional day count. It does NOT set the corresponding
542 // Julian parameters for the current AliTimestamp instance.
543 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
544 // To set the Julian parameters for the current AliTimestamp instance,
545 // please use the corresponding SET() memberfunctions of either AliTimestamp
548 Double_t frac=double(secs)+double(ns)*1.e-9;
549 Int_t daysecs=24*3600;
550 frac=frac/double(daysecs);
551 Double_t date=double(days)+frac;
554 ///////////////////////////////////////////////////////////////////////////
555 void AliTimestamp::Convert(Double_t h,Int_t& hh,Int_t& mm,Int_t& ss,Int_t& ns,Int_t& ps) const
557 // Convert fractional hour count h into hh:mm:ss:ns:ps.
558 // The sign of the input value will be neglected, so h<0 will result in
559 // the same output values as h>0.
561 // Note : Due to computer accuracy the ps value may become inaccurate.
565 // This memberfunction only converts the input "h" into the corresponding
566 // integer parameters. It does NOT set the corresponding Julian parameters
567 // for the current AliTimestamp instance.
568 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
569 // To set the Julian parameters for the current AliTimestamp instance,
570 // please use the corresponding SET() memberfunctions of either AliTimestamp
590 ///////////////////////////////////////////////////////////////////////////
591 void AliTimestamp::Convert(Double_t h,Int_t& hh,Int_t& mm,Double_t& ss) const
593 // Convert fractional hour count h into hh:mm:ss.s.
594 // The sign of the input value will be neglected, so h<0 will result in
595 // the same output values as h>0.
599 // 1) This memberfunction only converts the input "h" into the corresponding
600 // hh:mm:ss.s values. It does NOT set the corresponding Julian parameters
601 // for the current AliTimestamp instance.
602 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
603 // To set the Julian parameters for the current AliTimestamp instance,
604 // please use the corresponding SET() memberfunctions of either AliTimestamp
606 // 2) This facility can also be used to convert degrees in arcminutes etc...
618 ///////////////////////////////////////////////////////////////////////////
619 Double_t AliTimestamp::Convert(Int_t hh,Int_t mm,Int_t ss,Int_t ns,Int_t ps) const
621 // Convert hh:mm:ss:ns:ps into fractional hour count.
622 // The sign of the input values will be neglected, so the output value
623 // will always correspond to a positive hh:mm:ss:ns:ps specification.
625 // Note : Due to computer accuracy the ps precision may be lost.
629 // This memberfunction only converts the input integer parameters into the
630 // corresponding fractional hour count. It does NOT set the corresponding
631 // Julian parameters for the current AliTimestamp instance.
632 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
633 // To set the Julian parameters for the current AliTimestamp instance,
634 // please use the corresponding SET() memberfunctions of either AliTimestamp
637 // Neglect the sign of the input values
645 h+=double(mm)/60.+(double(ss)+double(ns)*1.e-9+double(ps)*1.e-12)/3600.;
649 ///////////////////////////////////////////////////////////////////////////
650 Double_t AliTimestamp::Convert(Int_t hh,Int_t mm,Double_t ss) const
652 // Convert hh:mm:ss.s into fractional hour count.
653 // The sign of the input values will be neglected, so the output value
654 // will always correspond to a positive hh:mm:ss.s specification.
658 // 1) This memberfunction only converts the input hh:mm:ss.s data into the
659 // corresponding fractional hour count. It does NOT set the corresponding
660 // Julian parameters for the current AliTimestamp instance.
661 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
662 // To set the Julian parameters for the current AliTimestamp instance,
663 // please use the corresponding SET() memberfunctions of either AliTimestamp
665 // 2) This facility can also be used to convert ddd:mm:ss.s into fractional degrees.
667 // Neglect the sign of the input values
673 h+=double(mm)/60.+ss/3600.;
677 ///////////////////////////////////////////////////////////////////////////
678 void AliTimestamp::PrintTime(Double_t h,Int_t ndig) const
680 // Print a fractional hour count in hh:mm:ss.ssss format.
681 // The range of the printed hour value is : -24 < hh < 24.
682 // The argument "ndig" specifies the number of digits for the fractional
683 // seconds (e.g. ndig=6 corresponds to microsecond precision).
684 // No rounding will be performed, so a second count of 3.473 with ndig=1
685 // will appear as 03.4 on the output.
686 // Due to computer accuracy, precision on the picosecond level may get lost.
688 // The default is ndig=1.
690 // Note : The time info is printed without additional spaces or "endline".
691 // This allows the print to be included in various composite output formats.
712 if (h<0) cout << "-";
714 << setw(2) << hh << ":" << setw(2) << mm << ":"
715 << setw(2) << ss << "." << setw(ndig) << sfrac;
717 ///////////////////////////////////////////////////////////////////////////
718 void AliTimestamp::FillJulian()
720 // Calculation and setting of the Julian date/time parameters corresponding
721 // to the current TTimeStamp date/time parameters.
723 UInt_t y,m,d,hh,mm,ss;
725 GetDate(kTRUE,0,&y,&m,&d);
726 GetTime(kTRUE,0,&hh,&mm,&ss);
727 Int_t ns=GetNanoSec();
729 Double_t mjd=GetMJD(y,m,d,hh,mm,ss,ns);
732 fJsec=GetSec()%(24*3600); // Daytime in elapsed seconds
733 fJns=ns; // Remaining fractional elapsed second in nanoseconds
735 // Store the TTimeStamp seconds and nanoseconds values
736 // for which this Julian calculation was performed.
738 fCalcns=GetNanoSec();
740 ///////////////////////////////////////////////////////////////////////////
741 void AliTimestamp::GetMJD(Int_t& mjd,Int_t& sec,Int_t& ns)
743 // Provide the Modified Julian Date (MJD) and time corresponding to the
744 // currently stored AliTimestamp date/time parameters.
746 // The returned arguments represent the following :
747 // mjd : The modified Julian date.
748 // sec : The number of seconds elapsed within the MJD.
749 // ns : The remaining fractional number of seconds (in ns) elapsed within the MJD.
751 if (fCalcs != GetSec() || fCalcns != GetNanoSec()) FillJulian();
757 ///////////////////////////////////////////////////////////////////////////
758 Double_t AliTimestamp::GetMJD()
760 // Provide the (fractional) Modified Julian Date (MJD) corresponding to the
761 // currently stored AliTimestamp date/time parameters.
763 // Due to computer accuracy the ns precision may be lost.
764 // It is advised to use the (mjd,sec,ns) getter instead.
771 Double_t date=Convert(mjd,sec,ns);
775 ///////////////////////////////////////////////////////////////////////////
776 void AliTimestamp::GetTJD(Int_t& tjd,Int_t& sec, Int_t& ns)
778 // Provide the Truncated Julian Date (TJD) and time corresponding to the
779 // currently stored AliTimestamp date/time parameters.
781 // The returned arguments represent the following :
782 // tjd : The modified Julian date.
783 // sec : The number of seconds elapsed within the MJD.
784 // ns : The remaining fractional number of seconds (in ns) elapsed within the MJD.
791 ///////////////////////////////////////////////////////////////////////////
792 Double_t AliTimestamp::GetTJD()
794 // Provide the (fractional) Truncated Julian Date (TJD) corresponding to the
795 // currently stored AliTimestamp date/time parameters.
797 // Due to computer accuracy the ns precision may be lost.
798 // It is advised to use the (mjd,sec,ns) getter instead.
805 Double_t date=Convert(tjd,sec,ns);
809 ///////////////////////////////////////////////////////////////////////////
810 void AliTimestamp::GetJD(Int_t& jd,Int_t& sec, Int_t& ns)
812 // Provide the Julian Date (JD) and time corresponding to the currently
813 // stored AliTimestamp date/time parameters.
815 // The returned arguments represent the following :
816 // jd : The Julian date.
817 // sec : The number of seconds elapsed within the JD.
818 // ns : The remaining fractional number of seconds (in ns) elapsed within the JD.
831 ///////////////////////////////////////////////////////////////////////////
832 Double_t AliTimestamp::GetJD()
834 // Provide the (fractional) Julian Date (JD) corresponding to the currently
835 // stored AliTimestamp date/time parameters.
837 // Due to computer accuracy the ns precision may be lost.
838 // It is advised to use the (jd,sec,ns) getter instead.
845 Double_t date=Convert(jd,sec,ns);
849 ///////////////////////////////////////////////////////////////////////////
850 Double_t AliTimestamp::GetJE()
852 // Provide the Julian Epoch (JE) corresponding to the currently stored
853 // AliTimestamp date/time parameters.
856 Double_t je=GetJE(jd);
859 ///////////////////////////////////////////////////////////////////////////
860 Double_t AliTimestamp::GetBE()
862 // Provide the Besselian Epoch (BE) corresponding to the currently stored
863 // AliTimestamp date/time parameters.
866 Double_t be=GetBE(jd);
869 ///////////////////////////////////////////////////////////////////////////
870 void AliTimestamp::SetMJD(Int_t mjd,Int_t sec,Int_t ns,Int_t ps)
872 // Set the Modified Julian Date (MJD) and time and update the TTimeStamp
873 // parameters accordingly (if possible).
877 // The TTimeStamp EPOCH starts at 01-jan-1970 00:00:00 UTC
878 // which corresponds to the start of MJD=40587.
879 // Using the corresponding MJD of this EPOCH allows construction of
880 // the yy-mm-dd hh:mm:ss:ns TTimeStamp from a given input MJD and time.
881 // Obviously this TTimeStamp implementation would prevent usage of MJD values
882 // smaller than 40587.
883 // Furthermore, due to a limitation on the "seconds since the EPOCH start" count
884 // in TTimeStamp, the latest accessible date/time is 19-jan-2038 02:14:08 UTC.
885 // However, this AliTimestamp facility provides support for the full range
886 // of (M)JD values, but the setting of the corresponding TTimeStamp parameters
887 // is restricted to the values allowed by the TTimeStamp implementation.
888 // For these earlier/later MJD values, the standard TTimeStamp parameters will
889 // be set corresponding to the start of the TTimeStamp EPOCH.
890 // This implies that for these earlier/later MJD values the TTimeStamp parameters
891 // do not match the Julian parameters of AliTimestamp.
893 // The input arguments represent the following :
894 // mjd : The modified Julian date.
895 // sec : The number of seconds elapsed within the MJD.
896 // ns : The remaining fractional number of seconds (in ns) elapsed within the MJD.
897 // ps : The remaining fractional number of nanoseconds (in ps) elapsed within the MJD.
899 // Note : ps=0 is the default value.
901 if (sec<0 || sec>=24*3600 || ns<0 || ns>=1e9 || ps<0 || ps>=1000)
903 cout << " *AliTimestamp::SetMJD* Invalid input."
904 << " sec : " << sec << " ns : " << ns << endl;
913 Int_t epoch=40587; // MJD of the start of the epoch
914 Int_t limit=65442; // MJD of the latest possible TTimeStamp date/time
917 if (mjd<epoch || (mjd>=limit && sec>=8047))
919 Set(0,kFALSE,0,kFALSE);
922 Set(date,time,0,kTRUE,0);
926 // The elapsed time since start of EPOCH
927 Int_t days=mjd-epoch;
928 UInt_t secs=days*24*3600;
930 Set(secs,kFALSE,0,kFALSE);
933 Set(date,time,ns,kTRUE,0);
936 // Denote that the Julian and TTimeStamp parameters are synchronised,
937 // even in the case the MJD falls outside the TTimeStamp validity range.
938 // The latter still allows retrieval of Julian parameters for these
941 fCalcns=GetNanoSec();
943 ///////////////////////////////////////////////////////////////////////////
944 void AliTimestamp::SetMJD(Double_t mjd)
946 // Set the Modified Julian Date (MJD) and time and update the TTimeStamp
947 // parameters accordingly (if possible).
951 // The TTimeStamp EPOCH starts at 01-jan-1970 00:00:00 UTC
952 // which corresponds to the start of MJD=40587.
953 // Using the corresponding MJD of this EPOCH allows construction of
954 // the yy-mm-dd hh:mm:ss:ns TTimeStamp from a given input MJD and time.
955 // Obviously this TTimeStamp implementation would prevent usage of MJD values
956 // smaller than 40587.
957 // Furthermore, due to a limitation on the "seconds since the EPOCH start" count
958 // in TTimeStamp, the latest accessible date/time is 19-jan-2038 02:14:08 UTC.
959 // However, this AliTimestamp facility provides support for the full range
960 // of (M)JD values, but the setting of the corresponding TTimeStamp parameters
961 // is restricted to the values allowed by the TTimeStamp implementation.
962 // For these earlier/later MJD values, the standard TTimeStamp parameters will
963 // be set corresponding to the start of the TTimeStamp EPOCH.
964 // This implies that for these earlier/later MJD values the TTimeStamp parameters
965 // do not match the Julian parameters of AliTimestamp.
967 // Due to computer accuracy the ns precision may be lost.
968 // It is advised to use the (mjd,sec,ns) setting instead.
970 // The input argument represents the following :
971 // mjd : The modified Julian date as fractional day count.
976 Convert(mjd,days,secs,ns);
977 SetMJD(days,secs,ns);
979 ///////////////////////////////////////////////////////////////////////////
980 void AliTimestamp::SetJD(Int_t jd,Int_t sec,Int_t ns,Int_t ps)
982 // Set the Julian Date (JD) and time and update the TTimeStamp
983 // parameters accordingly (if possible).
987 // The TTimeStamp EPOCH starts at 01-jan-1970 00:00:00 UTC
988 // which corresponds to JD=2440587.5 or the start of MJD=40587.
989 // Using the corresponding MJD of this EPOCH allows construction of
990 // the yy-mm-dd hh:mm:ss:ns TTimeStamp from a given input MJD and time.
991 // Obviously this TTimeStamp implementation would prevent usage of values
992 // smaller than JD=2440587.5.
993 // Furthermore, due to a limitation on the "seconds since the EPOCH start" count
994 // in TTimeStamp, the latest accessible date/time is 19-jan-2038 02:14:08 UTC.
995 // However, this AliTimestamp facility provides support for the full range
996 // of (M)JD values, but the setting of the corresponding TTimeStamp parameters
997 // is restricted to the values allowed by the TTimeStamp implementation.
998 // For these earlier/later JD values, the standard TTimeStamp parameters will
999 // be set corresponding to the start of the TTimeStamp EPOCH.
1000 // This implies that for these earlier/later (M)JD values the TTimeStamp parameters
1001 // do not match the Julian parameters of AliTimestamp.
1003 // The input arguments represent the following :
1004 // jd : The Julian date.
1005 // sec : The number of seconds elapsed within the JD.
1006 // ns : The remaining fractional number of seconds (in ns) elapsed within the JD.
1007 // ps : The remaining fractional number of nanoseconds (in ps) elapsed within the JD.
1009 // Note : ps=0 is the default value.
1011 Int_t mjd=jd-2400000;
1019 SetMJD(mjd,sec,ns,ps);
1021 ///////////////////////////////////////////////////////////////////////////
1022 void AliTimestamp::SetJD(Double_t jd)
1024 // Set the Julian Date (JD) and time and update the TTimeStamp
1025 // parameters accordingly (if possible).
1029 // The TTimeStamp EPOCH starts at 01-jan-1970 00:00:00 UTC
1030 // which corresponds to JD=2440587.5 or the start of MJD=40587.
1031 // Using the corresponding MJD of this EPOCH allows construction of
1032 // the yy-mm-dd hh:mm:ss:ns TTimeStamp from a given input MJD and time.
1033 // Obviously this TTimeStamp implementation would prevent usage of values
1034 // smaller than JD=2440587.5.
1035 // Furthermore, due to a limitation on the "seconds since the EPOCH start" count
1036 // in TTimeStamp, the latest accessible date/time is 19-jan-2038 02:14:08 UTC.
1037 // However, this AliTimestamp facility provides support for the full range
1038 // of (M)JD values, but the setting of the corresponding TTimeStamp parameters
1039 // is restricted to the values allowed by the TTimeStamp implementation.
1040 // For these earlier/later JD values, the standard TTimeStamp parameters will
1041 // be set corresponding to the start of the TTimeStamp EPOCH.
1042 // This implies that for these earlier/later (M)JD values the TTimeStamp parameters
1043 // do not match the Julian parameters of AliTimestamp.
1045 // Due to computer accuracy the ns precision may be lost.
1046 // It is advised to use the (jd,sec,ns) setting instead.
1048 // The input argument represents the following :
1049 // jd : The Julian date as fractional day count.
1054 Convert(jd,days,secs,ns);
1056 SetJD(days,secs,ns);
1058 ///////////////////////////////////////////////////////////////////////////
1059 void AliTimestamp::SetTJD(Int_t tjd,Int_t sec,Int_t ns,Int_t ps)
1061 // Set the Truncated Julian Date (TJD) and time and update the TTimeStamp
1062 // parameters accordingly (if possible).
1066 // The TTimeStamp EPOCH starts at 01-jan-1970 00:00:00 UTC
1067 // which corresponds to JD=2440587.5 or the start of TJD=587.
1068 // Using the corresponding MJD of this EPOCH allows construction of
1069 // the yy-mm-dd hh:mm:ss:ns TTimeStamp from a given input MJD and time.
1070 // Obviously this TTimeStamp implementation would prevent usage of values
1071 // smaller than TJD=587.
1072 // Furthermore, due to a limitation on the "seconds since the EPOCH start" count
1073 // in TTimeStamp, the latest accessible date/time is 19-jan-2038 02:14:08 UTC.
1074 // However, this AliTimestamp facility provides support for the full range
1075 // of (T)JD values, but the setting of the corresponding TTimeStamp parameters
1076 // is restricted to the values allowed by the TTimeStamp implementation.
1077 // For these earlier/later JD values, the standard TTimeStamp parameters will
1078 // be set corresponding to the start of the TTimeStamp EPOCH.
1079 // This implies that for these earlier/later (T)JD values the TTimeStamp parameters
1080 // do not match the Julian parameters of AliTimestamp.
1082 // The input arguments represent the following :
1083 // tjd : The Truncated Julian date.
1084 // sec : The number of seconds elapsed within the JD.
1085 // ns : The remaining fractional number of seconds (in ns) elapsed within the JD.
1086 // ps : The remaining fractional number of nanoseconds (in ps) elapsed within the JD.
1088 // Note : ps=0 is the default value.
1090 Int_t mjd=tjd+40000;
1092 SetMJD(mjd,sec,ns,ps);
1094 ///////////////////////////////////////////////////////////////////////////
1095 void AliTimestamp::SetTJD(Double_t tjd)
1097 // Set the Truncated Julian Date (TJD) and time and update the TTimeStamp
1098 // parameters accordingly (if possible).
1102 // The TTimeStamp EPOCH starts at 01-jan-1970 00:00:00 UTC
1103 // which corresponds to JD=2440587.5 or the start of TJD=587.
1104 // Using the corresponding MJD of this EPOCH allows construction of
1105 // the yy-mm-dd hh:mm:ss:ns TTimeStamp from a given input MJD and time.
1106 // Obviously this TTimeStamp implementation would prevent usage of values
1107 // smaller than TJD=587.
1108 // Furthermore, due to a limitation on the "seconds since the EPOCH start" count
1109 // in TTimeStamp, the latest accessible date/time is 19-jan-2038 02:14:08 UTC.
1110 // However, this AliTimestamp facility provides support for the full range
1111 // of (T)JD values, but the setting of the corresponding TTimeStamp parameters
1112 // is restricted to the values allowed by the TTimeStamp implementation.
1113 // For these earlier/later JD values, the standard TTimeStamp parameters will
1114 // be set corresponding to the start of the TTimeStamp EPOCH.
1115 // This implies that for these earlier/later (T)JD values the TTimeStamp parameters
1116 // do not match the Julian parameters of AliTimestamp.
1118 // Due to computer accuracy the ns precision may be lost.
1119 // It is advised to use the (jd,sec,ns) setting instead.
1121 // The input argument represents the following :
1122 // tjd : The Truncated Julian date as fractional day count.
1127 Convert(tjd,days,secs,ns);
1129 SetTJD(days,secs,ns);
1131 ///////////////////////////////////////////////////////////////////////////
1132 void AliTimestamp::SetNs(Int_t ns)
1134 // Set the remaining fractional number of seconds in nanosecond precision.
1137 // 1) The allowed range for the argument "ns" is [0,99999999].
1138 // Outside that range no action is performed.
1139 // 2) The ns fraction can also be entered directly via SetMJD() etc...
1140 // 3) For additional accuracy see SetPs().
1142 if (ns>=0 && ns<=99999999) fJns=ns;
1144 ///////////////////////////////////////////////////////////////////////////
1145 Int_t AliTimestamp::GetNs() const
1147 // Provide the remaining fractional number of seconds in nanosecond precision.
1148 // This function allows trigger/timing analysis for (astro)particle physics
1150 // Note : For additional accuracy see also GetPs().
1154 ///////////////////////////////////////////////////////////////////////////
1155 void AliTimestamp::SetPs(Int_t ps)
1157 // Set the remaining fractional number of nanoseconds in picoseconds.
1160 // 1) The allowed range for the argument "ps" is [0,999].
1161 // Outside that range no action is performed.
1162 // 2) The ps fraction can also be entered directly via SetMJD() etc...
1164 if (ps>=0 && ps<=999) fJps=ps;
1166 ///////////////////////////////////////////////////////////////////////////
1167 Int_t AliTimestamp::GetPs() const
1169 // Provide remaining fractional number of nanoseconds in picoseconds.
1170 // This function allows time of flight analysis for particle physics
1175 ///////////////////////////////////////////////////////////////////////////
1176 void AliTimestamp::Add(Int_t d,Int_t s,Int_t ns,Int_t ps)
1178 // Add (or subtract) a certain time difference to the current timestamp.
1179 // Subtraction can be achieved by entering negative values as input arguments.
1181 // The time difference is entered via the following input arguments :
1183 // d : elapsed number of days
1184 // s : (remaining) elapsed number of seconds
1185 // ns : (remaining) elapsed number of nanoseconds
1186 // ps : (remaining) elapsed number of picoseconds
1188 // The specified d, s, ns and ps values will be used in an additive
1189 // way to determine the time difference.
1190 // So, specification of d=1, s=100, ns=0, ps=0 will result in the
1191 // same time difference addition as d=0, s=24*3600+100, ns=0, ps=0.
1192 // However, by making use of the latter the user should take care
1193 // of possible integer overflow problems in the input arguments,
1194 // which obviously will provide incorrect results.
1196 // Note : ps=0 is the default value.
1201 // Use Get functions to ensure updated Julian parameters.
1202 GetMJD(days,secs,nsec);
1218 nsec+=ns%1000000000;
1219 secs+=ns/1000000000;
1225 while (nsec>999999999)
1238 while (secs>=24*3600)
1246 SetMJD(days,secs,nsec,psec);
1248 ///////////////////////////////////////////////////////////////////////////
1249 void AliTimestamp::Add(Double_t hours)
1251 // Add (or subtract) a certain time difference to the current timestamp.
1252 // The time difference is specified as a (fractional) number of hours.
1253 // Subtraction can be achieved by entering a negative value as input argument.
1256 Double_t h=fabs(hours);
1266 if (hours>0) Add(d,s,ns,ps);
1267 if (hours<0) Add(-d,-s,-ns,-ps);
1269 ///////////////////////////////////////////////////////////////////////////
1270 Int_t AliTimestamp::GetDifference(AliTimestamp* t,Int_t& d,Int_t& s,Int_t& ns,Int_t& ps)
1272 // Provide the time difference w.r.t the AliTimestamp specified on the input.
1273 // This memberfunction supports both very small (i.e. time of flight analysis
1274 // for particle physics experiments) and very long (i.e. investigation of
1275 // astrophysical phenomena) timescales.
1277 // The time difference is returned via the following output arguments :
1278 // d : elapsed number of days
1279 // s : remaining elapsed number of seconds
1280 // ns : remaining elapsed number of nanoseconds
1281 // ps : remaining elapsed number of picoseconds
1285 // The calculated time difference is the absolute value of the time interval.
1286 // This implies that the values of d, s, ns and ps are always positive or zero.
1288 // The integer return argument indicates whether the AliTimestamp specified
1289 // on the input argument occurred earlier (-1), simultaneously (0) or later (1).
1293 // Ensure updated Julian parameters for this AliTimestamp instance
1294 if (fCalcs != GetSec() || fCalcns != GetNanoSec()) FillJulian();
1296 // Use Get functions to ensure updated Julian parameters.
1305 if (!d && !s && !ns && !ps) return 0;
1312 if (!sign && s>0) sign=1;
1313 if (!sign && s<0) sign=-1;
1315 if (!sign && ns>0) sign=1;
1316 if (!sign && ns<0) sign=-1;
1318 if (!sign && ps>0) sign=1;
1319 if (!sign && ps<0) sign=-1;
1321 // In case the input stamp was earlier, take the reverse difference
1322 // to simplify the algebra.
1331 // Here we always have a positive time difference
1332 // and can now unambiguously correct for other negative values.
1353 ///////////////////////////////////////////////////////////////////////////
1354 Int_t AliTimestamp::GetDifference(AliTimestamp& t,Int_t& d,Int_t& s,Int_t& ns,Int_t& ps)
1356 // Provide the time difference w.r.t the AliTimestamp specified on the input.
1357 // This memberfunction supports both very small (i.e. time of flight analysis
1358 // for particle physics experiments) and very long (i.e. investigation of
1359 // astrophysical phenomena) timescales.
1361 // The time difference is returned via the following output arguments :
1362 // d : elapsed number of days
1363 // s : remaining elapsed number of seconds
1364 // ns : remaining elapsed number of nanoseconds
1365 // ps : remaining elapsed number of picoseconds
1369 // The calculated time difference is the absolute value of the time interval.
1370 // This implies that the values of d, s, ns and ps are always positive or zero.
1372 // The integer return argument indicates whether the AliTimestamp specified
1373 // on the input argument occurred earlier (-1), simultaneously (0) or later (1).
1375 return GetDifference(&t,d,s,ns,ps);
1377 ///////////////////////////////////////////////////////////////////////////
1378 Double_t AliTimestamp::GetDifference(AliTimestamp* t,TString u,Int_t mode)
1380 // Provide the time difference w.r.t the AliTimestamp specified on the input
1381 // argument in the units as specified by the TString argument.
1382 // A positive return value means that the AliTimestamp specified on the input
1383 // argument occurred later, whereas a negative return value indicates an
1384 // earlier occurence.
1386 // The units may be specified as :
1387 // u = "d" ==> Time difference returned as (fractional) day count
1388 // "s" ==> Time difference returned as (fractional) second count
1389 // "ns" ==> Time difference returned as (fractional) nanosecond count
1390 // "ps" ==> Time difference returned as picosecond count
1392 // It may be clear that for a time difference of several days, the picosecond
1393 // and even the nanosecond accuracy may be lost.
1394 // To cope with this, the "mode" argument has been introduced to allow
1395 // timestamp comparison on only the specified units.
1397 // The following operation modes are supported :
1398 // mode = 1 : Full time difference is returned in specified units
1399 // 2 : Time difference is returned in specified units by
1400 // neglecting the elapsed time for the larger units than the
1402 // 3 : Time difference is returned in specified units by only
1403 // comparing the timestamps on the level of the specified units.
1407 // AliTimestamp t1; // Corresponding to days=3, secs=501, ns=31, ps=7
1408 // AliTimestamp t2; // Corresponding to days=5, secs=535, ns=12, ps=15
1410 // The statement : Double_t val=t1.GetDifference(t2,....)
1411 // would return the following values :
1412 // val=(2*24*3600)+34-(19*1e-9)+(8*1e-12) for u="s" and mode=1
1413 // val=34-(19*1e-9)+(8*1e-12) for u="s" and mode=2
1414 // val=34 for u="s" and mode=3
1415 // val=-19 for u="ns" and mode=3
1417 // The default is mode=1.
1419 if (!t || mode<1 || mode>3) return 0;
1423 // Ensure updated Julian parameters for this AliTimestamp instance
1424 if (fCalcs != GetSec() || fCalcns != GetNanoSec()) FillJulian();
1431 // Use Get functions to ensure updated Julian parameters.
1432 t->GetMJD(dd,ds,dns);
1440 // Time difference for the specified units only
1445 if (u=="ns") dt=dns;
1446 if (u=="ps") dt=dps;
1450 // Suppress elapsed time for the larger units than specified
1467 // Compute the time difference as requested
1468 if (u=="s" || u=="d")
1470 // The time difference in (fractional) seconds
1471 dt=double(dd*24*3600+ds)+(double(dns)*1e-9)+(double(dps)*1e-12);
1472 if (u=="d") dt=dt/double(24*3600);
1474 if (u=="ns") dt=(double(dd*24*3600+ds)*1e9)+double(dns)+(double(dps)*1e-3);
1475 if (u=="ps") dt=(double(dd*24*3600+ds)*1e12)+(double(dns)*1e3)+double(dps);
1479 ///////////////////////////////////////////////////////////////////////////
1480 Double_t AliTimestamp::GetDifference(AliTimestamp& t,TString u,Int_t mode)
1482 // Provide the time difference w.r.t the AliTimestamp specified on the input
1483 // argument in the units as specified by the TString argument.
1484 // A positive return value means that the AliTimestamp specified on the input
1485 // argument occurred later, whereas a negative return value indicates an
1486 // earlier occurence.
1488 // The units may be specified as :
1489 // u = "d" ==> Time difference returned as (fractional) day count
1490 // "s" ==> Time difference returned as (fractional) second count
1491 // "ns" ==> Time difference returned as (fractional) nanosecond count
1492 // "ps" ==> Time difference returned as picosecond count
1494 // It may be clear that for a time difference of several days, the picosecond
1495 // and even the nanosecond accuracy may be lost.
1496 // To cope with this, the "mode" argument has been introduced to allow
1497 // timestamp comparison on only the specified units.
1499 // The following operation modes are supported :
1500 // mode = 1 : Full time difference is returned in specified units
1501 // 2 : Time difference is returned in specified units by
1502 // neglecting the elapsed time for the larger units than the
1504 // 3 : Time difference is returned in specified units by only
1505 // comparing the timestamps on the level of the specified units.
1509 // AliTimestamp t1; // Corresponding to days=3, secs=501, ns=31, ps=7
1510 // AliTimestamp t2; // Corresponding to days=5, secs=535, ns=12, ps=15
1512 // The statement : Double_t val=t1.GetDifference(t2,....)
1513 // would return the following values :
1514 // val=(2*24*3600)+34-(19*1e-9)+(8*1e-12) for u="s" and mode=1
1515 // val=34-(19*1e-9)+(8*1e-12) for u="s" and mode=2
1516 // val=34 for u="s" and mode=3
1517 // val=-19 for u="ns" and mode=3
1519 // The default is mode=1.
1521 return GetDifference(&t,u,mode);
1523 ///////////////////////////////////////////////////////////////////////////
1524 void AliTimestamp::SetUT(Int_t y,Int_t m,Int_t d,Int_t hh,Int_t mm,Int_t ss,Int_t ns,Int_t ps)
1526 // Set the AliTimestamp parameters corresponding to the UT date and time
1527 // in the Gregorian calendar as specified by the input arguments.
1528 // This facility is exact upto picosecond precision and as such is
1529 // for scientific observations preferable above the corresponding
1530 // Set function(s) of TTimestamp.
1531 // The latter has a random spread in the sub-second part, which
1532 // might be of use in generating distinguishable timestamps while
1533 // still keeping second precision.
1535 // The input arguments represent the following :
1536 // y : year in UT (e.g. 1952, 2003 etc...)
1537 // m : month in UT (1=jan 2=feb etc...)
1538 // d : day in UT (1-31)
1539 // hh : elapsed hours in UT (0-23)
1540 // mm : elapsed minutes in UT (0-59)
1541 // ss : elapsed seconds in UT (0-59)
1542 // ns : remaining fractional elapsed second of UT in nanosecond
1543 // ps : remaining fractional elapsed nanosecond of UT in picosecond
1545 // Note : ns=0 and ps=0 are the default values.
1547 // This facility first determines the elapsed days, seconds etc...
1548 // since the beginning of the specified UT year on basis of the
1549 // input arguments. Subsequently it invokes the SetUT memberfunction
1550 // for the elapsed timespan.
1551 // As such this facility is valid for all AD dates in the Gregorian
1552 // calendar with picosecond precision.
1554 Int_t day=GetDayOfYear(d,m,y);
1555 Int_t secs=hh*3600+mm*60+ss;
1556 SetUT(y,day-1,secs,ns,ps);
1558 ///////////////////////////////////////////////////////////////////////////
1559 void AliTimestamp::SetUT(Int_t y,Int_t d,Int_t s,Int_t ns,Int_t ps)
1561 // Set the AliTimestamp parameters corresponding to the specified elapsed
1562 // timespan since the beginning of the new UT year.
1563 // This facility is exact upto picosecond precision and as such is
1564 // for scientific observations preferable above the corresponding
1565 // Set function(s) of TTimestamp.
1566 // The latter has a random spread in the sub-second part, which
1567 // might be of use in generating distinguishable timestamps while
1568 // still keeping second precision.
1570 // The UT year and elapsed time span is entered via the following input arguments :
1572 // y : year in UT (e.g. 1952, 2003 etc...)
1573 // d : elapsed number of days
1574 // s : (remaining) elapsed number of seconds
1575 // ns : (remaining) elapsed number of nanoseconds
1576 // ps : (remaining) elapsed number of picoseconds
1578 // The specified d, s, ns and ps values will be used in an additive
1579 // way to determine the elapsed timespan.
1580 // So, specification of d=1, s=100, ns=0, ps=0 will result in the
1581 // same elapsed time span as d=0, s=24*3600+100, ns=0, ps=0.
1582 // However, by making use of the latter the user should take care
1583 // of possible integer overflow problems in the input arguments,
1584 // which obviously will provide incorrect results.
1586 // Note : ns=0 and ps=0 are the default values.
1588 // This facility first sets the (M)JD corresponding to the start (01-jan 00:00:00)
1589 // of the specified UT year following the recipe of R.W. Sinnott
1590 // Sky & Telescope 82, (aug. 1991) 183.
1591 // Subsequently the day and (sub)second parts are added to the AliTimestamp.
1592 // As such this facility is valid for all AD dates in the Gregorian calendar.
1594 Double_t jd=GetJD(y,1,1,0,0,0,0);
1598 GetMJD(mjd,sec,nsec);
1602 ///////////////////////////////////////////////////////////////////////////
1603 void AliTimestamp::GetUT(Int_t& hh,Int_t& mm,Int_t& ss,Int_t& ns,Int_t& ps)
1605 // Provide the corrresponding UT as hh:mm:ss:ns:ps.
1606 // This facility is based on the MJD, so the TTimeStamp limitations
1607 // do not apply here.
1609 Int_t mjd,sec,nsec,psec;
1611 GetMJD(mjd,sec,nsec);
1621 ///////////////////////////////////////////////////////////////////////////
1622 Double_t AliTimestamp::GetUT()
1624 // Provide the corrresponding UT in fractional hours.
1625 // This facility is based on the MJD, so the TTimeStamp limitations
1626 // do not apply here.
1628 Int_t hh,mm,ss,ns,ps;
1630 GetUT(hh,mm,ss,ns,ps);
1632 Double_t ut=Convert(hh,mm,ss,ns,ps);
1636 ///////////////////////////////////////////////////////////////////////////
1637 void AliTimestamp::GetGMST(Int_t& hh,Int_t& mm,Int_t& ss,Int_t& ns,Int_t& ps)
1639 // Provide the corrresponding Greenwich Mean Sideral Time (GMST).
1640 // The algorithm used is the one described at p. 83 of the book
1641 // Astronomy Methods by Hale Bradt.
1642 // This facility is based on the MJD, so the TTimeStamp limitations
1643 // do not apply here.
1645 Int_t mjd,sec,nsec,psec;
1647 // The current UT based timestamp data
1648 GetMJD(mjd,sec,nsec);
1651 // The basis for the daily corrections in units of Julian centuries w.r.t. J2000.
1652 // Note : Epoch J2000 starts at 01-jan-2000 12:00:00 UT.
1653 Double_t tau=(GetJD()-2451545.)/36525.;
1655 // Syncronise sidereal time with current timestamp
1657 sid.SetMJD(mjd,sec,nsec,psec);
1659 // Add offset for GMST start value defined as 06:41:50.54841 at 01-jan 00:00:00 UT
1660 sec=6*3600+41*60+50;
1663 sid.Add(0,sec,nsec,psec);
1665 // Daily correction for precession and polar motion
1666 Double_t addsec=8640184.812866*tau+0.093104*pow(tau,2)-6.2e-6*pow(tau,3);
1668 addsec-=double(sec);
1669 nsec=int(addsec*1.e9);
1670 addsec-=double(nsec)*1.e-9;
1671 psec=int(addsec*1.e12);
1672 sid.Add(0,sec,nsec,psec);
1674 sid.GetMJD(mjd,sec,nsec);
1684 ///////////////////////////////////////////////////////////////////////////
1685 Double_t AliTimestamp::GetGMST()
1687 // Provide the corrresponding Greenwich Mean Sideral Time (GMST)
1688 // in fractional hours.
1689 // This facility is based on the MJD, so the TTimeStamp limitations
1690 // do not apply here.
1692 Int_t hh,mm,ss,ns,ps;
1694 GetGMST(hh,mm,ss,ns,ps);
1696 Double_t gst=Convert(hh,mm,ss,ns,ps);
1700 ///////////////////////////////////////////////////////////////////////////
1701 Double_t AliTimestamp::GetGAST()
1703 // Provide the corrresponding Greenwich Apparent Sideral Time (GAST)
1704 // in fractional hours.
1705 // In case a hh:mm:ss.sss format is needed, please invoke the Convert()
1706 // memberfunction for conversion of the provided fractional hour value.
1708 // The GAST is the GMST corrected for the shift of the vernal equinox
1709 // due to nutation. The right ascension component of the nutation correction
1710 // of the vernal equinox is called the "equation of the equinoxes".
1713 // GAST = GMST + (equation of the equinoxes)
1715 // With the right ascension and declination values to be (0,0) for the
1716 // vernal equinox, we can workout the shift in right ascension of the vernal
1717 // equinox due to nutation.
1718 // The nutation model used is the new one as documented in :
1719 // "The IAU Resolutions on Astronomical Reference Systems,
1720 // Time Scales and Earth Rotation Models".
1721 // This document is freely available as Circular 179 (2005) of the
1722 // United States Naval Observatory (USNO).
1723 // (See : http://aa.usno.navy.mil/publications/docs).
1724 // The new expression for the equation of the equinoxes is based on a series
1725 // expansion and is the most accurate one known to date.
1726 // The components are documented on p.17 of the USNO Circular 179.
1728 // Since GMST is based on the MJD, the TTimeStamp limitations
1729 // do not apply here.
1731 Double_t pi=acos(-1.);
1733 Double_t gmst=GetGMST();
1735 Double_t days; // Time difference in fractional Julian days w.r.t. the start of J2000.
1736 Double_t t; // Time difference in fractional Julian centuries w.r.t. the start of J2000.
1737 Double_t epsilon; // Mean obliquity of the ecliptic
1738 Double_t omega; // Mean longitude of the Moon's mean ascending mode
1739 Double_t f; // Mean argument of latitude of the moon
1740 Double_t d; // Mean elongation of the Moon from the Sun
1742 days=GetJD()-2451545.0;
1745 // Values in arcseconds
1746 epsilon=84381.406-46.836769*t-0.0001831*pow(t,2)+0.00200340*pow(t,3)
1747 -0.000000576*pow(t,4)-0.0000000434*pow(t,5);
1748 omega=450160.398036-6962890.5431*t+7.4722*pow(t,2)+0.007702*pow(t,3)-0.00005939*pow(t,4);
1749 f=335779.526232+1739527262.8478*t-12.7512*pow(t,2)-0.001037*pow(t,3)+0.00000417*pow(t,4);
1750 d=1072260.70369+1602961601.2090*t-6.3706*pow(t,2)+0.006593*pow(t,3)-0.00003169*pow(t,4);
1752 // Convert to radians
1753 epsilon=epsilon*pi/(180.*3600.);
1754 omega=omega*pi/(180.*3600.);
1755 f=f*pi/(180.*3600.);
1756 d=d*pi/(180.*3600.);
1758 Double_t dpsi; // Nutation in longitude in arcseconds
1759 Double_t beta=(125.04-0.052954*days)*pi/180.;
1760 Double_t gamma=(200.94+1.97130*days)*pi/180.;
1761 dpsi=-17.226*sin(beta)-1.296*sin(gamma);
1763 Double_t da; // Right ascension shift of the vernal equinox in arcseconds
1764 da=dpsi*cos(epsilon)+0.00264096*sin(omega)+0.00006352*sin(2.*omega)
1765 +0.00001175*sin(2.*f-2.*d+3.*omega)+0.00001121*sin(2.*f-2.*d+omega)
1766 -0.00000455*sin(2.*f-2.*d+2.*omega)+0.00000202*sin(2.*f+3.*omega)+0.00000198*sin(2.*f+omega)
1767 -0.00000172*sin(3.*omega)-0.00000087*t*sin(omega);
1769 // Convert to fractional hours
1772 Double_t gast=gmst+da;
1785 ///////////////////////////////////////////////////////////////////////////
1786 Double_t AliTimestamp::GetLT(Double_t offset)
1788 // Provide the corresponding local time in fractional hours.
1789 // The "offset" denotes the time difference in (fractional) hours w.r.t. UT.
1790 // A mean solar day lasts 24h (i.e. 86400s).
1792 // In case a hh:mm:ss format is needed, please use the Convert() facility.
1794 // Current UT time in fractional hours
1810 ///////////////////////////////////////////////////////////////////////////
1811 Double_t AliTimestamp::GetLMST(Double_t offset)
1813 // Provide the corresponding Local Mean Sidereal Time (LMST) in fractional hours.
1814 // The "offset" denotes the time difference in (fractional) hours w.r.t. GMST.
1815 // A sidereal day corresponds to 23h 56m 04.09s (i.e. 86164.09s) mean solar time.
1816 // The definition of GMST is such that a sidereal clock corresponds with
1817 // 24 sidereal hours per revolution of the Earth.
1818 // As such, local time offsets w.r.t. UT and GMST can be treated similarly.
1820 // In case a hh:mm:ss format is needed, please use the Convert() facility.
1822 // Current GMST time in fractional hours
1823 Double_t h=GetGMST();
1838 ///////////////////////////////////////////////////////////////////////////
1839 Double_t AliTimestamp::GetLAST(Double_t offset)
1841 // Provide the corresponding Local Apparent Sidereal Time (LAST) in fractional hours.
1842 // The "offset" denotes the time difference in (fractional) hours w.r.t. GAST.
1843 // A sidereal day corresponds to 23h 56m 04.09s (i.e. 86164.09s) mean solar time.
1844 // The definition of GMST and GAST is such that a sidereal clock corresponds with
1845 // 24 sidereal hours per revolution of the Earth.
1846 // As such, local time offsets w.r.t. UT, GMST and GAST can be treated similarly.
1848 // In case a hh:mm:ss.sss format is needed, please use the Convert() facility.
1850 // Current GAST time in fractional hours
1851 Double_t h=GetGAST();
1866 ///////////////////////////////////////////////////////////////////////////
1867 void AliTimestamp::SetLT(Double_t dt,Int_t y,Int_t m,Int_t d,Int_t hh,Int_t mm,Int_t ss,Int_t ns,Int_t ps)
1869 // Set the AliTimestamp parameters corresponding to the LT date and time
1870 // in the Gregorian calendar as specified by the input arguments.
1871 // This facility is exact upto picosecond precision and as such is
1872 // for scientific observations preferable above the corresponding
1873 // Set function(s) of TTimestamp.
1874 // The latter has a random spread in the sub-second part, which
1875 // might be of use in generating distinguishable timestamps while
1876 // still keeping second precision.
1878 // The input arguments represent the following :
1880 // dt : the local time offset in fractional hours w.r.t. UT.
1881 // y : year in LT (e.g. 1952, 2003 etc...)
1882 // m : month in LT (1=jan 2=feb etc...)
1883 // d : day in LT (1-31)
1884 // hh : elapsed hours in LT (0-23)
1885 // mm : elapsed minutes in LT (0-59)
1886 // ss : elapsed seconds in LT (0-59)
1887 // ns : remaining fractional elapsed second of LT in nanosecond
1888 // ps : remaining fractional elapsed nanosecond of LT in picosecond
1890 // Note : ns=0 and ps=0 are the default values.
1892 // This facility first sets the UT as specified by the input arguments
1893 // and then corrects the UT by subtracting the local time offset w.r.t. UT.
1894 // As such this facility is valid for all AD dates in the Gregorian
1895 // calendar with picosecond precision.
1897 SetUT(y,m,d,hh,mm,ss,ns,ps);
1900 ///////////////////////////////////////////////////////////////////////////
1901 void AliTimestamp::SetLT(Double_t dt,Int_t y,Int_t d,Int_t s,Int_t ns,Int_t ps)
1903 // Set the AliTimestamp parameters corresponding to the specified elapsed
1904 // timespan since the beginning of the new LT year.
1905 // This facility is exact upto picosecond precision and as such is
1906 // for scientific observations preferable above the corresponding
1907 // Set function(s) of TTimestamp.
1908 // The latter has a random spread in the sub-second part, which
1909 // might be of use in generating distinguishable timestamps while
1910 // still keeping second precision.
1912 // The LT year and elapsed time span is entered via the following input arguments :
1914 // dt : the local time offset in fractional hours w.r.t. UT.
1915 // y : year in LT (e.g. 1952, 2003 etc...)
1916 // d : elapsed number of days
1917 // s : (remaining) elapsed number of seconds
1918 // ns : (remaining) elapsed number of nanoseconds
1919 // ps : (remaining) elapsed number of picoseconds
1921 // The specified d, s, ns and ps values will be used in an additive
1922 // way to determine the elapsed timespan.
1923 // So, specification of d=1, s=100, ns=0, ps=0 will result in the
1924 // same elapsed time span as d=0, s=24*3600+100, ns=0, ps=0.
1925 // However, by making use of the latter the user should take care
1926 // of possible integer overflow problems in the input arguments,
1927 // which obviously will provide incorrect results.
1929 // Note : ns=0 and ps=0 are the default values.
1931 // This facility first sets the UT as specified by the input arguments
1932 // and then corrects the UT by subtracting the local time offset w.r.t. UT.
1933 // As such this facility is valid for all AD dates in the Gregorian calendar.
1938 ///////////////////////////////////////////////////////////////////////////
1939 Double_t AliTimestamp::GetJD(Double_t e,TString mode) const
1941 // Provide the fractional Julian Date from epoch e.
1942 // The sort of epoch may be specified via the "mode" parameter.
1944 // mode = "J" ==> Julian epoch
1945 // "B" ==> Besselian epoch
1947 // The default value is mode="J".
1951 if (mode=="J" || mode=="j") jd=(e-2000.0)*365.25+2451545.0;
1953 if (mode=="B" || mode=="b") jd=(e-1900.0)*365.242198781+2415020.31352;
1957 ///////////////////////////////////////////////////////////////////////////
1958 Double_t AliTimestamp::GetMJD(Double_t e,TString mode) const
1960 // Provide the fractional Modified Julian Date from epoch e.
1961 // The sort of epoch may be specified via the "mode" parameter.
1963 // mode = "J" ==> Julian epoch
1964 // "B" ==> Besselian epoch
1966 // The default value is mode="J".
1968 Double_t mjd=GetJD(e,mode)-2400000.5;
1972 ///////////////////////////////////////////////////////////////////////////
1973 Double_t AliTimestamp::GetTJD(Double_t e,TString mode) const
1975 // Provide the fractional Truncated Julian Date from epoch e.
1976 // The sort of epoch may be specified via the "mode" parameter.
1978 // mode = "J" ==> Julian epoch
1979 // "B" ==> Besselian epoch
1981 // The default value is mode="J".
1983 Double_t tjd=GetJD(e,mode)-2440000.5;
1987 ///////////////////////////////////////////////////////////////////////////