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 ///////////////////////////////////////////////////////////////////////////
171 #include "AliTimestamp.h"
172 #include "Riostream.h"
174 ClassImp(AliTimestamp) // Class implementation to enable ROOT I/O
176 AliTimestamp::AliTimestamp() : TTimeStamp()
178 // Default constructor
179 // Creation of an AliTimestamp object and initialisation of parameters.
180 // All attributes are initialised to the current date/time as specified
181 // in the docs of TTimeStamp.
186 ///////////////////////////////////////////////////////////////////////////
187 AliTimestamp::AliTimestamp(TTimeStamp& t) : TTimeStamp(t)
189 // Creation of an AliTimestamp object and initialisation of parameters.
190 // All attributes are initialised to the values of the input TTimeStamp.
195 ///////////////////////////////////////////////////////////////////////////
196 AliTimestamp::~AliTimestamp()
198 // Destructor to delete dynamically allocated memory.
200 ///////////////////////////////////////////////////////////////////////////
201 AliTimestamp::AliTimestamp(const AliTimestamp& t) : TTimeStamp(t)
212 ///////////////////////////////////////////////////////////////////////////
213 void AliTimestamp::Date(Int_t mode,Double_t offset)
215 // Print date/time info.
217 // mode = 1 ==> Only the UT yy-mm-dd hh:mm:ss.sss and GMST info is printed
218 // 2 ==> Only the Julian parameter info is printed
219 // 3 ==> Both the UT, GMST and Julian parameter info is printed
220 // -1 ==> Only the UT yy-mm-dd hh:mm:ss.sss and GAST info is printed
221 // -3 ==> Both the UT, GAST and Julian parameter info is printed
223 // offset : Local time offset from UT (and also GMST) in fractional hours.
225 // When an offset value is specified, the corresponding local times
226 // LT and LMST (or LAST) are printed as well.
228 // The default values are mode=3 and offset=0.
230 // Note : In case the (M/T)JD falls outside the TTimeStamp range,
231 // the yy-mm-dd info will be omitted.
233 Int_t mjd,mjsec,mjns,mjps;
234 GetMJD(mjd,mjsec,mjns);
237 TString month[12]={"Jan","Feb","Mar","Apr","May","Jun","Jul","Aug","Sep","Oct","Nov","Dec"};
238 TString day[7]={"Mon","Tue","Wed","Thu","Fri","Sat","Sun"};
240 Int_t hh,mm,ss,ns,ps;
243 if (abs(mode)==1 || abs(mode)==3)
245 if (mjd>=40587 && (mjd<65442 || (mjd==65442 && mjsec<8047)))
247 GetDate(kTRUE,0,&y,&m,&d);
248 wd=GetDayOfWeek(kTRUE,0);
249 cout << " " << day[wd-1].Data() << ", " << setfill('0') << setw(2) << d << " "
250 << setfill(' ') << month[m-1].Data() << " " << y << " ";
256 GetUT(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 << " (UT) ";
262 GetGMST(hh,mm,ss,ns,ps);
267 Convert(gast,hh,mm,ss,ns,ps);
269 cout << setfill('0') << setw(2) << hh << ":"
270 << setw(2) << mm << ":" << setw(2) << ss << "."
271 << setw(9) << ns << setw(3) << ps;
274 cout << " (GMST)" << endl;
278 cout << " (GAST)" << endl;
281 // Local time information
284 // Determine the new date by including the offset
285 AliTimestamp t2(*this);
287 Int_t mjd2,mjsec2,mjns2;
288 t2.GetMJD(mjd2,mjsec2,mjns2);
289 if (mjd2>=40587 && (mjd2<65442 || (mjd2==65442 && mjsec2<8047)))
291 t2.GetDate(kTRUE,0,&y,&m,&d);
292 wd=t2.GetDayOfWeek(kTRUE,0);
293 cout << " " << day[wd-1].Data() << ", " << setfill('0') << setw(2) << d << " "
294 << setfill(' ') << month[m-1].Data() << " " << y << " ";
300 // Determine the local time by including the offset w.r.t. the original timestamp
301 Double_t hlt=GetLT(offset);
305 hlst=GetLMST(offset);
309 hlst=GetLAST(offset);
311 PrintTime(hlt,12); cout << " (LT) "; PrintTime(hlst,12);
314 cout << " (LMST)" << endl;
318 cout << " (LAST)" << endl;
322 if (abs(mode)==2 || abs(mode)==3)
326 Int_t tjd,tjsec,tjns;
327 GetTJD(tjd,tjsec,tjns);
328 cout << " Julian Epoch : " << setprecision(25) << GetJE()
329 << " Besselian Epoch : " << setprecision(25) << GetBE() << endl;
330 cout << " JD : " << jd << " sec : " << jsec << " ns : " << jns << " ps : " << fJps
331 << " Fractional : " << setprecision(25) << GetJD() << endl;
332 cout << " MJD : " << mjd << " sec : " << mjsec << " ns : " << mjns << " ps : " << fJps
333 << " Fractional : " << setprecision(25) << GetMJD() << endl;
334 cout << " TJD : " << tjd << " sec : " << tjsec << " ns : " << tjns << " ps : " << fJps
335 << " Fractional : " << setprecision(25) << GetTJD() << endl;
338 ///////////////////////////////////////////////////////////////////////////
339 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
341 // Provide the (fractional) Julian Date (JD) corresponding to the UT date
342 // and time in the Gregorian calendar as specified by the input arguments.
344 // The input arguments represent the following :
345 // y : year in UT (e.g. 1952, 2003 etc...)
346 // m : month in UT (1=jan 2=feb etc...)
347 // d : day in UT (1-31)
348 // hh : elapsed hours in UT (0-23)
349 // mm : elapsed minutes in UT (0-59)
350 // ss : elapsed seconds in UT (0-59)
351 // ns : remaining fractional elapsed second of UT in nanosecond
353 // This algorithm is valid for all AD dates in the Gregorian calendar
354 // following the recipe of R.W. Sinnott Sky & Telescope 82, (aug. 1991) 183.
355 // See also http://scienceworld.wolfram.com/astronomy/JulianDate.html
357 // In case of invalid input, a value of -1 is returned.
361 // This memberfunction only provides the JD corresponding to the
362 // UT input arguments. It does NOT set the corresponding Julian parameters
363 // for the current AliTimestamp instance.
364 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
365 // To set the Julian parameters for the current AliTimestamp instance,
366 // please use the corresponding SET() memberfunctions of either AliTimestamp
369 if (y<0 || m<1 || m>12 || d<1 || d>31) return -1;
370 if (hh<0 || hh>23 || mm<0 || mm>59 || ss<0 || ss>59 || ns<0 || ns>1e9) return -1;
372 // The UT daytime in fractional hours
373 Double_t ut=double(hh)+double(mm)/60.+(double(ss)+double(ns)*1.e-9)/3600.;
377 JD=367*y-int(7*(y+int((m+9)/12))/4)
378 -int(3*(int((y+(m-9)/7)/100)+1)/4)
379 +int(275*m/9)+d+1721028.5+ut/24.;
383 ///////////////////////////////////////////////////////////////////////////
384 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
386 // Provide the (fractional) Modified Julian Date corresponding to the UT
387 // date and time in the Gregorian calendar as specified by the input arguments.
389 // The input arguments represent the following :
390 // y : year in UT (e.g. 1952, 2003 etc...)
391 // m : month in UT (1=jan 2=feb etc...)
392 // d : day in UT (1-31)
393 // hh : elapsed hours in UT (0-23)
394 // mm : elapsed minutes in UT (0-59)
395 // ss : elapsed seconds in UT (0-59)
396 // ns : remaining fractional elapsed second of UT in nanosecond
398 // This algorithm is valid for all AD dates in the Gregorian calendar
399 // following the recipe of R.W. Sinnott Sky & Telescope 82, (aug. 1991) 183.
400 // See also http://scienceworld.wolfram.com/astronomy/JulianDate.html
402 // In case of invalid input, a value of -1 is returned.
406 // This memberfunction only provides the MJD corresponding to the
407 // UT input arguments. It does NOT set the corresponding Julian parameters
408 // for the current AliTimestamp instance.
409 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
410 // To set the Julian parameters for the current AliTimestamp instance,
411 // please use the corresponding SET() memberfunctions of either AliTimestamp
414 Double_t JD=GetJD(y,m,d,hh,mm,ss,ns);
418 Double_t MJD=JD-2400000.5;
422 ///////////////////////////////////////////////////////////////////////////
423 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
425 // Provide the (fractional) Truncated Julian Date corresponding to the UT
426 // date and time in the Gregorian calendar as specified by the input arguments.
428 // The input arguments represent the following :
429 // y : year in UT (e.g. 1952, 2003 etc...)
430 // m : month in UT (1=jan 2=feb etc...)
431 // d : day in UT (1-31)
432 // hh : elapsed hours in UT (0-23)
433 // mm : elapsed minutes in UT (0-59)
434 // ss : elapsed seconds in UT (0-59)
435 // ns : remaining fractional elapsed second of UT in nanosecond
437 // This algorithm is valid for all AD dates in the Gregorian calendar
438 // following the recipe of R.W. Sinnott Sky & Telescope 82, (aug. 1991) 183.
439 // See also http://scienceworld.wolfram.com/astronomy/JulianDate.html
441 // In case of invalid input, a value of -1 is returned.
445 // This memberfunction only provides the TJD corresponding to the
446 // UT input arguments. It does NOT set the corresponding Julian parameters
447 // for the current AliTimestamp instance.
448 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
449 // To set the Julian parameters for the current AliTimestamp instance,
450 // please use the corresponding SET() memberfunctions of either AliTimestamp
453 Double_t JD=GetJD(y,m,d,hh,mm,ss,ns);
457 Double_t TJD=JD-2440000.5;
461 ///////////////////////////////////////////////////////////////////////////
462 Double_t AliTimestamp::GetJE(Double_t date,TString mode) const
464 // Provide the Julian Epoch (JE) corresponding to the specified date.
465 // The argument "mode" indicates the type of the argument "date".
467 // Available modes are :
468 // mode = "jd" ==> date represents the Julian Date
469 // = "mjd" ==> date represents the Modified Julian Date
470 // = "tjd" ==> date represents the Truncated Julian Date
472 // The default is mode="jd".
474 // In case of invalid input, a value of -99999 is returned.
478 // This memberfunction only provides the JE corresponding to the
479 // input arguments. It does NOT set the corresponding Julian parameters
480 // for the current AliTimestamp instance.
481 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
482 // To set the Julian parameters for the current AliTimestamp instance,
483 // please use the corresponding SET() memberfunctions of either AliTimestamp
486 if ((mode != "jd") && (mode != "mjd") && (mode != "tjd")) return -99999;
489 if (mode=="mjd") jd=date+2400000.5;
490 if (mode=="tjd") jd=date+2440000.5;
492 Double_t je=2000.+(jd-2451545.)/365.25;
496 ///////////////////////////////////////////////////////////////////////////
497 Double_t AliTimestamp::GetBE(Double_t date,TString mode) const
499 // Provide the Besselian Epoch (JE) corresponding to the specified date.
500 // The argument "mode" indicates the type of the argument "date".
502 // Available modes are :
503 // mode = "jd" ==> date represents the Julian Date
504 // = "mjd" ==> date represents the Modified Julian Date
505 // = "tjd" ==> date represents the Truncated Julian Date
507 // The default is mode="jd".
509 // In case of invalid input, a value of -99999 is returned.
513 // This memberfunction only provides the BE corresponding to the
514 // input arguments. It does NOT set the corresponding Julian parameters
515 // for the current AliTimestamp instance.
516 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
517 // To set the Julian parameters for the current AliTimestamp instance,
518 // please use the corresponding SET() memberfunctions of either AliTimestamp
521 if ((mode != "jd") && (mode != "mjd") && (mode != "tjd")) return -99999;
524 if (mode=="mjd") jd=date+2400000.5;
525 if (mode=="tjd") jd=date+2440000.5;
527 Double_t be=1900.+(jd-2415020.31352)/365.242198781;
531 ///////////////////////////////////////////////////////////////////////////
532 void AliTimestamp::Convert(Double_t date,Int_t& days,Int_t& secs,Int_t& ns) const
534 // Convert date as fractional day count into integer days, secs and ns.
536 // Note : Due to computer accuracy the ns value may become inaccurate.
538 // The arguments represent the following :
539 // date : The input date as fractional day count
540 // days : Number of elapsed days
541 // secs : Remaining number of elapsed seconds
542 // ns : Remaining fractional elapsed second in nanoseconds
546 // This memberfunction only converts the input date into the corresponding
547 // integer parameters. It does NOT set the corresponding Julian parameters
548 // for the current AliTimestamp instance.
549 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
550 // To set the Julian parameters for the current AliTimestamp instance,
551 // please use the corresponding SET() memberfunctions of either AliTimestamp
555 date=date-double(days);
556 Int_t daysecs=24*3600;
557 date=date*double(daysecs);
559 date=date-double(secs);
562 ///////////////////////////////////////////////////////////////////////////
563 Double_t AliTimestamp::Convert(Int_t days,Int_t secs,Int_t ns) const
565 // Convert date in integer days, secs and ns into fractional day count.
567 // Note : Due to computer accuracy the ns precision may be lost.
569 // The input arguments represent the following :
570 // days : Number of elapsed days
571 // secs : Remaining number of elapsed seconds
572 // ns : Remaining fractional elapsed second in nanoseconds
576 // This memberfunction only converts the input integer parameters into the
577 // corresponding fractional day count. It does NOT set the corresponding
578 // Julian parameters for the current AliTimestamp instance.
579 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
580 // To set the Julian parameters for the current AliTimestamp instance,
581 // please use the corresponding SET() memberfunctions of either AliTimestamp
584 Double_t frac=double(secs)+double(ns)*1.e-9;
585 Int_t daysecs=24*3600;
586 frac=frac/double(daysecs);
587 Double_t date=double(days)+frac;
590 ///////////////////////////////////////////////////////////////////////////
591 void AliTimestamp::Convert(Double_t h,Int_t& hh,Int_t& mm,Int_t& ss,Int_t& ns,Int_t& ps) const
593 // Convert fractional hour count h into hh:mm:ss:ns:ps.
594 // The sign of the input value will be neglected, so h<0 will result in
595 // the same output values as h>0.
597 // Note : Due to computer accuracy the ps value may become inaccurate.
601 // This memberfunction only converts the input "h" into the corresponding
602 // integer parameters. It does NOT set the corresponding Julian parameters
603 // for the current AliTimestamp instance.
604 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
605 // To set the Julian parameters for the current AliTimestamp instance,
606 // please use the corresponding SET() memberfunctions of either AliTimestamp
626 ///////////////////////////////////////////////////////////////////////////
627 void AliTimestamp::Convert(Double_t h,Int_t& hh,Int_t& mm,Double_t& ss) const
629 // Convert fractional hour count h into hh:mm:ss.s.
630 // The sign of the input value will be neglected, so h<0 will result in
631 // the same output values as h>0.
635 // 1) This memberfunction only converts the input "h" into the corresponding
636 // hh:mm:ss.s values. It does NOT set the corresponding Julian parameters
637 // for the current AliTimestamp instance.
638 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
639 // To set the Julian parameters for the current AliTimestamp instance,
640 // please use the corresponding SET() memberfunctions of either AliTimestamp
642 // 2) This facility can also be used to convert degrees in arcminutes etc...
654 ///////////////////////////////////////////////////////////////////////////
655 Double_t AliTimestamp::Convert(Int_t hh,Int_t mm,Int_t ss,Int_t ns,Int_t ps) const
657 // Convert hh:mm:ss:ns:ps into fractional hour count.
658 // The sign of the input values will be neglected, so the output value
659 // will always correspond to a positive hh:mm:ss:ns:ps specification.
661 // Note : Due to computer accuracy the ps precision may be lost.
665 // This memberfunction only converts the input integer parameters into the
666 // corresponding fractional hour count. It does NOT set the corresponding
667 // Julian parameters for the current AliTimestamp instance.
668 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
669 // To set the Julian parameters for the current AliTimestamp instance,
670 // please use the corresponding SET() memberfunctions of either AliTimestamp
673 // Neglect the sign of the input values
681 h+=double(mm)/60.+(double(ss)+double(ns)*1.e-9+double(ps)*1.e-12)/3600.;
685 ///////////////////////////////////////////////////////////////////////////
686 Double_t AliTimestamp::Convert(Int_t hh,Int_t mm,Double_t ss) const
688 // Convert hh:mm:ss.s into fractional hour count.
689 // The sign of the input values will be neglected, so the output value
690 // will always correspond to a positive hh:mm:ss.s specification.
694 // 1) This memberfunction only converts the input hh:mm:ss.s data into the
695 // corresponding fractional hour count. It does NOT set the corresponding
696 // Julian parameters for the current AliTimestamp instance.
697 // As such the TTimeStamp limitations do NOT apply to this memberfunction.
698 // To set the Julian parameters for the current AliTimestamp instance,
699 // please use the corresponding SET() memberfunctions of either AliTimestamp
701 // 2) This facility can also be used to convert ddd:mm:ss.s into fractional degrees.
703 // Neglect the sign of the input values
709 h+=double(mm)/60.+ss/3600.;
713 ///////////////////////////////////////////////////////////////////////////
714 void AliTimestamp::PrintTime(Double_t h,Int_t ndig) const
716 // Print a fractional hour count in hh:mm:ss.ssss format.
717 // The range of the printed hour value is : -24 < hh < 24.
718 // The argument "ndig" specifies the number of digits for the fractional
719 // seconds (e.g. ndig=6 corresponds to microsecond precision).
720 // No rounding will be performed, so a second count of 3.473 with ndig=1
721 // will appear as 03.4 on the output.
722 // Due to computer accuracy, precision on the picosecond level may get lost.
724 // The default is ndig=1.
726 // Note : The time info is printed without additional spaces or "endline".
727 // This allows the print to be included in various composite output formats.
748 if (h<0) cout << "-";
750 << setw(2) << hh << ":" << setw(2) << mm << ":"
751 << setw(2) << ss << "." << setw(ndig) << sfrac;
753 ///////////////////////////////////////////////////////////////////////////
754 void AliTimestamp::FillJulian()
756 // Calculation and setting of the Julian date/time parameters corresponding
757 // to the current TTimeStamp date/time parameters.
759 UInt_t y,m,d,hh,mm,ss;
761 GetDate(kTRUE,0,&y,&m,&d);
762 GetTime(kTRUE,0,&hh,&mm,&ss);
763 Int_t ns=GetNanoSec();
765 Double_t mjd=GetMJD(y,m,d,hh,mm,ss,ns);
768 fJsec=GetSec()%(24*3600); // Daytime in elapsed seconds
769 fJns=ns; // Remaining fractional elapsed second in nanoseconds
771 // Store the TTimeStamp seconds and nanoseconds values
772 // for which this Julian calculation was performed.
774 fCalcns=GetNanoSec();
776 ///////////////////////////////////////////////////////////////////////////
777 void AliTimestamp::GetMJD(Int_t& mjd,Int_t& sec,Int_t& ns)
779 // Provide the Modified Julian Date (MJD) and time corresponding to the
780 // currently stored AliTimestamp date/time parameters.
782 // The returned arguments represent the following :
783 // mjd : The modified Julian date.
784 // sec : The number of seconds elapsed within the MJD.
785 // ns : The remaining fractional number of seconds (in ns) elapsed within the MJD.
787 if (fCalcs != GetSec() || fCalcns != GetNanoSec()) FillJulian();
793 ///////////////////////////////////////////////////////////////////////////
794 Double_t AliTimestamp::GetMJD()
796 // Provide the (fractional) Modified Julian Date (MJD) corresponding to the
797 // currently stored AliTimestamp date/time parameters.
799 // Due to computer accuracy the ns precision may be lost.
800 // It is advised to use the (mjd,sec,ns) getter instead.
807 Double_t date=Convert(mjd,sec,ns);
811 ///////////////////////////////////////////////////////////////////////////
812 void AliTimestamp::GetTJD(Int_t& tjd,Int_t& sec, Int_t& ns)
814 // Provide the Truncated Julian Date (TJD) and time corresponding to the
815 // currently stored AliTimestamp date/time parameters.
817 // The returned arguments represent the following :
818 // tjd : The modified Julian date.
819 // sec : The number of seconds elapsed within the MJD.
820 // ns : The remaining fractional number of seconds (in ns) elapsed within the MJD.
827 ///////////////////////////////////////////////////////////////////////////
828 Double_t AliTimestamp::GetTJD()
830 // Provide the (fractional) Truncated Julian Date (TJD) corresponding to the
831 // currently stored AliTimestamp date/time parameters.
833 // Due to computer accuracy the ns precision may be lost.
834 // It is advised to use the (mjd,sec,ns) getter instead.
841 Double_t date=Convert(tjd,sec,ns);
845 ///////////////////////////////////////////////////////////////////////////
846 void AliTimestamp::GetJD(Int_t& jd,Int_t& sec, Int_t& ns)
848 // Provide the Julian Date (JD) and time corresponding to the currently
849 // stored AliTimestamp date/time parameters.
851 // The returned arguments represent the following :
852 // jd : The Julian date.
853 // sec : The number of seconds elapsed within the JD.
854 // ns : The remaining fractional number of seconds (in ns) elapsed within the JD.
867 ///////////////////////////////////////////////////////////////////////////
868 Double_t AliTimestamp::GetJD()
870 // Provide the (fractional) Julian Date (JD) corresponding to the currently
871 // stored AliTimestamp date/time parameters.
873 // Due to computer accuracy the ns precision may be lost.
874 // It is advised to use the (jd,sec,ns) getter instead.
881 Double_t date=Convert(jd,sec,ns);
885 ///////////////////////////////////////////////////////////////////////////
886 Double_t AliTimestamp::GetJE()
888 // Provide the Julian Epoch (JE) corresponding to the currently stored
889 // AliTimestamp date/time parameters.
892 Double_t je=GetJE(jd);
895 ///////////////////////////////////////////////////////////////////////////
896 Double_t AliTimestamp::GetBE()
898 // Provide the Besselian Epoch (BE) corresponding to the currently stored
899 // AliTimestamp date/time parameters.
902 Double_t be=GetBE(jd);
905 ///////////////////////////////////////////////////////////////////////////
906 void AliTimestamp::SetMJD(Int_t mjd,Int_t sec,Int_t ns,Int_t ps)
908 // Set the Modified Julian Date (MJD) and time and update the TTimeStamp
909 // parameters accordingly (if possible).
913 // The TTimeStamp EPOCH starts at 01-jan-1970 00:00:00 UTC
914 // which corresponds to the start of MJD=40587.
915 // Using the corresponding MJD of this EPOCH allows construction of
916 // the yy-mm-dd hh:mm:ss:ns TTimeStamp from a given input MJD and time.
917 // Obviously this TTimeStamp implementation would prevent usage of MJD values
918 // smaller than 40587.
919 // Furthermore, due to a limitation on the "seconds since the EPOCH start" count
920 // in TTimeStamp, the latest accessible date/time is 19-jan-2038 02:14:08 UTC.
921 // However, this AliTimestamp facility provides support for the full range
922 // of (M)JD values, but the setting of the corresponding TTimeStamp parameters
923 // is restricted to the values allowed by the TTimeStamp implementation.
924 // For these earlier/later MJD values, the standard TTimeStamp parameters will
925 // be set corresponding to the start of the TTimeStamp EPOCH.
926 // This implies that for these earlier/later MJD values the TTimeStamp parameters
927 // do not match the Julian parameters of AliTimestamp.
929 // The input arguments represent the following :
930 // mjd : The modified Julian date.
931 // sec : The number of seconds elapsed within the MJD.
932 // ns : The remaining fractional number of seconds (in ns) elapsed within the MJD.
933 // ps : The remaining fractional number of nanoseconds (in ps) elapsed within the MJD.
935 // Note : ps=0 is the default value.
937 if (sec<0 || sec>=24*3600 || ns<0 || ns>=1e9 || ps<0 || ps>=1000)
939 cout << " *AliTimestamp::SetMJD* Invalid input."
940 << " sec : " << sec << " ns : " << ns << endl;
949 Int_t epoch=40587; // MJD of the start of the epoch
950 Int_t limit=65442; // MJD of the latest possible TTimeStamp date/time
953 if (mjd<epoch || mjd>limit || (mjd==limit && sec>=8047))
955 Set(0,kFALSE,0,kFALSE);
958 Set(date,time,0,kTRUE,0);
962 // The elapsed time since start of EPOCH
963 Int_t days=mjd-epoch;
964 UInt_t secs=days*24*3600;
966 Set(secs,kFALSE,0,kFALSE);
969 Set(date,time,ns,kTRUE,0);
972 // Denote that the Julian and TTimeStamp parameters are synchronised,
973 // even in the case the MJD falls outside the TTimeStamp validity range.
974 // The latter still allows retrieval of Julian parameters for these
977 fCalcns=GetNanoSec();
979 ///////////////////////////////////////////////////////////////////////////
980 void AliTimestamp::SetMJD(Double_t mjd)
982 // Set the Modified Julian Date (MJD) 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 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 MJD values
992 // smaller than 40587.
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 MJD 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 MJD values the TTimeStamp parameters
1001 // do not match the Julian parameters of AliTimestamp.
1003 // Due to computer accuracy the ns precision may be lost.
1004 // It is advised to use the (mjd,sec,ns) setting instead.
1006 // The input argument represents the following :
1007 // mjd : The modified Julian date as fractional day count.
1012 Convert(mjd,days,secs,ns);
1013 SetMJD(days,secs,ns);
1015 ///////////////////////////////////////////////////////////////////////////
1016 void AliTimestamp::SetJD(Int_t jd,Int_t sec,Int_t ns,Int_t ps)
1018 // Set the Julian Date (JD) and time and update the TTimeStamp
1019 // parameters accordingly (if possible).
1023 // The TTimeStamp EPOCH starts at 01-jan-1970 00:00:00 UTC
1024 // which corresponds to JD=2440587.5 or the start of MJD=40587.
1025 // Using the corresponding MJD of this EPOCH allows construction of
1026 // the yy-mm-dd hh:mm:ss:ns TTimeStamp from a given input MJD and time.
1027 // Obviously this TTimeStamp implementation would prevent usage of values
1028 // smaller than JD=2440587.5.
1029 // Furthermore, due to a limitation on the "seconds since the EPOCH start" count
1030 // in TTimeStamp, the latest accessible date/time is 19-jan-2038 02:14:08 UTC.
1031 // However, this AliTimestamp facility provides support for the full range
1032 // of (M)JD values, but the setting of the corresponding TTimeStamp parameters
1033 // is restricted to the values allowed by the TTimeStamp implementation.
1034 // For these earlier/later JD values, the standard TTimeStamp parameters will
1035 // be set corresponding to the start of the TTimeStamp EPOCH.
1036 // This implies that for these earlier/later (M)JD values the TTimeStamp parameters
1037 // do not match the Julian parameters of AliTimestamp.
1039 // The input arguments represent the following :
1040 // jd : The Julian date.
1041 // sec : The number of seconds elapsed within the JD.
1042 // ns : The remaining fractional number of seconds (in ns) elapsed within the JD.
1043 // ps : The remaining fractional number of nanoseconds (in ps) elapsed within the JD.
1045 // Note : ps=0 is the default value.
1047 Int_t mjd=jd-2400000;
1055 SetMJD(mjd,sec,ns,ps);
1057 ///////////////////////////////////////////////////////////////////////////
1058 void AliTimestamp::SetJD(Double_t jd)
1060 // Set the Julian Date (JD) and time and update the TTimeStamp
1061 // parameters accordingly (if possible).
1065 // The TTimeStamp EPOCH starts at 01-jan-1970 00:00:00 UTC
1066 // which corresponds to JD=2440587.5 or the start of MJD=40587.
1067 // Using the corresponding MJD of this EPOCH allows construction of
1068 // the yy-mm-dd hh:mm:ss:ns TTimeStamp from a given input MJD and time.
1069 // Obviously this TTimeStamp implementation would prevent usage of values
1070 // smaller than JD=2440587.5.
1071 // Furthermore, due to a limitation on the "seconds since the EPOCH start" count
1072 // in TTimeStamp, the latest accessible date/time is 19-jan-2038 02:14:08 UTC.
1073 // However, this AliTimestamp facility provides support for the full range
1074 // of (M)JD values, but the setting of the corresponding TTimeStamp parameters
1075 // is restricted to the values allowed by the TTimeStamp implementation.
1076 // For these earlier/later JD values, the standard TTimeStamp parameters will
1077 // be set corresponding to the start of the TTimeStamp EPOCH.
1078 // This implies that for these earlier/later (M)JD values the TTimeStamp parameters
1079 // do not match the Julian parameters of AliTimestamp.
1081 // Due to computer accuracy the ns precision may be lost.
1082 // It is advised to use the (jd,sec,ns) setting instead.
1084 // The input argument represents the following :
1085 // jd : The Julian date as fractional day count.
1090 Convert(jd,days,secs,ns);
1092 SetJD(days,secs,ns);
1094 ///////////////////////////////////////////////////////////////////////////
1095 void AliTimestamp::SetTJD(Int_t tjd,Int_t sec,Int_t ns,Int_t ps)
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 // The input arguments represent the following :
1119 // tjd : The Truncated Julian date.
1120 // sec : The number of seconds elapsed within the JD.
1121 // ns : The remaining fractional number of seconds (in ns) elapsed within the JD.
1122 // ps : The remaining fractional number of nanoseconds (in ps) elapsed within the JD.
1124 // Note : ps=0 is the default value.
1126 Int_t mjd=tjd+40000;
1128 SetMJD(mjd,sec,ns,ps);
1130 ///////////////////////////////////////////////////////////////////////////
1131 void AliTimestamp::SetTJD(Double_t tjd)
1133 // Set the Truncated Julian Date (TJD) and time and update the TTimeStamp
1134 // parameters accordingly (if possible).
1138 // The TTimeStamp EPOCH starts at 01-jan-1970 00:00:00 UTC
1139 // which corresponds to JD=2440587.5 or the start of TJD=587.
1140 // Using the corresponding MJD of this EPOCH allows construction of
1141 // the yy-mm-dd hh:mm:ss:ns TTimeStamp from a given input MJD and time.
1142 // Obviously this TTimeStamp implementation would prevent usage of values
1143 // smaller than TJD=587.
1144 // Furthermore, due to a limitation on the "seconds since the EPOCH start" count
1145 // in TTimeStamp, the latest accessible date/time is 19-jan-2038 02:14:08 UTC.
1146 // However, this AliTimestamp facility provides support for the full range
1147 // of (T)JD values, but the setting of the corresponding TTimeStamp parameters
1148 // is restricted to the values allowed by the TTimeStamp implementation.
1149 // For these earlier/later JD values, the standard TTimeStamp parameters will
1150 // be set corresponding to the start of the TTimeStamp EPOCH.
1151 // This implies that for these earlier/later (T)JD values the TTimeStamp parameters
1152 // do not match the Julian parameters of AliTimestamp.
1154 // Due to computer accuracy the ns precision may be lost.
1155 // It is advised to use the (jd,sec,ns) setting instead.
1157 // The input argument represents the following :
1158 // tjd : The Truncated Julian date as fractional day count.
1163 Convert(tjd,days,secs,ns);
1165 SetTJD(days,secs,ns);
1167 ///////////////////////////////////////////////////////////////////////////
1168 void AliTimestamp::SetNs(Int_t ns)
1170 // Set the remaining fractional number of seconds in nanosecond precision.
1173 // 1) The allowed range for the argument "ns" is [0,99999999].
1174 // Outside that range no action is performed.
1175 // 2) The ns fraction can also be entered directly via SetMJD() etc...
1176 // 3) For additional accuracy see SetPs().
1178 if (ns>=0 && ns<=99999999) fJns=ns;
1180 ///////////////////////////////////////////////////////////////////////////
1181 Int_t AliTimestamp::GetNs() const
1183 // Provide the remaining fractional number of seconds in nanosecond precision.
1184 // This function allows trigger/timing analysis for (astro)particle physics
1186 // Note : For additional accuracy see also GetPs().
1190 ///////////////////////////////////////////////////////////////////////////
1191 void AliTimestamp::SetPs(Int_t ps)
1193 // Set the remaining fractional number of nanoseconds in picoseconds.
1196 // 1) The allowed range for the argument "ps" is [0,999].
1197 // Outside that range no action is performed.
1198 // 2) The ps fraction can also be entered directly via SetMJD() etc...
1200 if (ps>=0 && ps<=999) fJps=ps;
1202 ///////////////////////////////////////////////////////////////////////////
1203 Int_t AliTimestamp::GetPs() const
1205 // Provide remaining fractional number of nanoseconds in picoseconds.
1206 // This function allows time of flight analysis for particle physics
1211 ///////////////////////////////////////////////////////////////////////////
1212 void AliTimestamp::Add(Int_t d,Int_t s,Int_t ns,Int_t ps)
1214 // Add (or subtract) a certain time difference to the current timestamp.
1215 // Subtraction can be achieved by entering negative values as input arguments.
1217 // The time difference is entered via the following input arguments :
1219 // d : elapsed number of days
1220 // s : (remaining) elapsed number of seconds
1221 // ns : (remaining) elapsed number of nanoseconds
1222 // ps : (remaining) elapsed number of picoseconds
1224 // The specified d, s, ns and ps values will be used in an additive
1225 // way to determine the time difference.
1226 // So, specification of d=1, s=100, ns=0, ps=0 will result in the
1227 // same time difference addition as d=0, s=24*3600+100, ns=0, ps=0.
1228 // However, by making use of the latter the user should take care
1229 // of possible integer overflow problems in the input arguments,
1230 // which obviously will provide incorrect results.
1232 // Note : ps=0 is the default value.
1237 // Use Get functions to ensure updated Julian parameters.
1238 GetMJD(days,secs,nsec);
1254 nsec+=ns%1000000000;
1255 secs+=ns/1000000000;
1261 while (nsec>999999999)
1274 while (secs>=24*3600)
1282 SetMJD(days,secs,nsec,psec);
1284 ///////////////////////////////////////////////////////////////////////////
1285 void AliTimestamp::Add(Double_t hours)
1287 // Add (or subtract) a certain time difference to the current timestamp.
1288 // The time difference is specified as a (fractional) number of hours.
1289 // Subtraction can be achieved by entering a negative value as input argument.
1292 Double_t h=fabs(hours);
1302 if (hours>0) Add(d,s,ns,ps);
1303 if (hours<0) Add(-d,-s,-ns,-ps);
1305 ///////////////////////////////////////////////////////////////////////////
1306 Int_t AliTimestamp::GetDifference(AliTimestamp* t,Int_t& d,Int_t& s,Int_t& ns,Int_t& ps)
1308 // Provide the time difference w.r.t the AliTimestamp specified on the input.
1309 // This memberfunction supports both very small (i.e. time of flight analysis
1310 // for particle physics experiments) and very long (i.e. investigation of
1311 // astrophysical phenomena) timescales.
1313 // The time difference is returned via the following output arguments :
1314 // d : elapsed number of days
1315 // s : remaining elapsed number of seconds
1316 // ns : remaining elapsed number of nanoseconds
1317 // ps : remaining elapsed number of picoseconds
1321 // The calculated time difference is the absolute value of the time interval.
1322 // This implies that the values of d, s, ns and ps are always positive or zero.
1324 // The integer return argument indicates whether the AliTimestamp specified
1325 // on the input argument occurred earlier (-1), simultaneously (0) or later (1).
1329 // Ensure updated Julian parameters for this AliTimestamp instance
1330 if (fCalcs != GetSec() || fCalcns != GetNanoSec()) FillJulian();
1332 // Use Get functions to ensure updated Julian parameters.
1341 if (!d && !s && !ns && !ps) return 0;
1348 if (!sign && s>0) sign=1;
1349 if (!sign && s<0) sign=-1;
1351 if (!sign && ns>0) sign=1;
1352 if (!sign && ns<0) sign=-1;
1354 if (!sign && ps>0) sign=1;
1355 if (!sign && ps<0) sign=-1;
1357 // In case the input stamp was earlier, take the reverse difference
1358 // to simplify the algebra.
1367 // Here we always have a positive time difference
1368 // and can now unambiguously correct for other negative values.
1389 ///////////////////////////////////////////////////////////////////////////
1390 Int_t AliTimestamp::GetDifference(AliTimestamp& t,Int_t& d,Int_t& s,Int_t& ns,Int_t& ps)
1392 // Provide the time difference w.r.t the AliTimestamp specified on the input.
1393 // This memberfunction supports both very small (i.e. time of flight analysis
1394 // for particle physics experiments) and very long (i.e. investigation of
1395 // astrophysical phenomena) timescales.
1397 // The time difference is returned via the following output arguments :
1398 // d : elapsed number of days
1399 // s : remaining elapsed number of seconds
1400 // ns : remaining elapsed number of nanoseconds
1401 // ps : remaining elapsed number of picoseconds
1405 // The calculated time difference is the absolute value of the time interval.
1406 // This implies that the values of d, s, ns and ps are always positive or zero.
1408 // The integer return argument indicates whether the AliTimestamp specified
1409 // on the input argument occurred earlier (-1), simultaneously (0) or later (1).
1411 return GetDifference(&t,d,s,ns,ps);
1413 ///////////////////////////////////////////////////////////////////////////
1414 Double_t AliTimestamp::GetDifference(AliTimestamp* t,TString u,Int_t mode)
1416 // Provide the time difference w.r.t the AliTimestamp specified on the input
1417 // argument in the units as specified by the TString argument.
1418 // A positive return value means that the AliTimestamp specified on the input
1419 // argument occurred later, whereas a negative return value indicates an
1420 // earlier occurence.
1422 // The units may be specified as :
1423 // u = "d" ==> Time difference returned as (fractional) day count
1424 // "s" ==> Time difference returned as (fractional) second count
1425 // "ns" ==> Time difference returned as (fractional) nanosecond count
1426 // "ps" ==> Time difference returned as picosecond count
1428 // It may be clear that for a time difference of several days, the picosecond
1429 // and even the nanosecond accuracy may be lost.
1430 // To cope with this, the "mode" argument has been introduced to allow
1431 // timestamp comparison on only the specified units.
1433 // The following operation modes are supported :
1434 // mode = 1 : Full time difference is returned in specified units
1435 // 2 : Time difference is returned in specified units by
1436 // neglecting the elapsed time for the larger units than the
1438 // 3 : Time difference is returned in specified units by only
1439 // comparing the timestamps on the level of the specified units.
1443 // AliTimestamp t1; // Corresponding to days=3, secs=501, ns=31, ps=7
1444 // AliTimestamp t2; // Corresponding to days=5, secs=535, ns=12, ps=15
1446 // The statement : Double_t val=t1.GetDifference(t2,....)
1447 // would return the following values :
1448 // val=(2*24*3600)+34-(19*1e-9)+(8*1e-12) for u="s" and mode=1
1449 // val=34-(19*1e-9)+(8*1e-12) for u="s" and mode=2
1450 // val=34 for u="s" and mode=3
1451 // val=-19 for u="ns" and mode=3
1453 // The default is mode=1.
1455 if (!t || mode<1 || mode>3) return 0;
1459 // Ensure updated Julian parameters for this AliTimestamp instance
1460 if (fCalcs != GetSec() || fCalcns != GetNanoSec()) FillJulian();
1467 // Use Get functions to ensure updated Julian parameters.
1468 t->GetMJD(dd,ds,dns);
1476 // Time difference for the specified units only
1481 if (u=="ns") dt=dns;
1482 if (u=="ps") dt=dps;
1486 // Suppress elapsed time for the larger units than specified
1503 // Compute the time difference as requested
1504 if (u=="s" || u=="d")
1506 // The time difference in (fractional) seconds
1507 dt=double(dd*24*3600+ds)+(double(dns)*1e-9)+(double(dps)*1e-12);
1508 if (u=="d") dt=dt/double(24*3600);
1510 if (u=="ns") dt=(double(dd*24*3600+ds)*1e9)+double(dns)+(double(dps)*1e-3);
1511 if (u=="ps") dt=(double(dd*24*3600+ds)*1e12)+(double(dns)*1e3)+double(dps);
1515 ///////////////////////////////////////////////////////////////////////////
1516 Double_t AliTimestamp::GetDifference(AliTimestamp& t,TString u,Int_t mode)
1518 // Provide the time difference w.r.t the AliTimestamp specified on the input
1519 // argument in the units as specified by the TString argument.
1520 // A positive return value means that the AliTimestamp specified on the input
1521 // argument occurred later, whereas a negative return value indicates an
1522 // earlier occurence.
1524 // The units may be specified as :
1525 // u = "d" ==> Time difference returned as (fractional) day count
1526 // "s" ==> Time difference returned as (fractional) second count
1527 // "ns" ==> Time difference returned as (fractional) nanosecond count
1528 // "ps" ==> Time difference returned as picosecond count
1530 // It may be clear that for a time difference of several days, the picosecond
1531 // and even the nanosecond accuracy may be lost.
1532 // To cope with this, the "mode" argument has been introduced to allow
1533 // timestamp comparison on only the specified units.
1535 // The following operation modes are supported :
1536 // mode = 1 : Full time difference is returned in specified units
1537 // 2 : Time difference is returned in specified units by
1538 // neglecting the elapsed time for the larger units than the
1540 // 3 : Time difference is returned in specified units by only
1541 // comparing the timestamps on the level of the specified units.
1545 // AliTimestamp t1; // Corresponding to days=3, secs=501, ns=31, ps=7
1546 // AliTimestamp t2; // Corresponding to days=5, secs=535, ns=12, ps=15
1548 // The statement : Double_t val=t1.GetDifference(t2,....)
1549 // would return the following values :
1550 // val=(2*24*3600)+34-(19*1e-9)+(8*1e-12) for u="s" and mode=1
1551 // val=34-(19*1e-9)+(8*1e-12) for u="s" and mode=2
1552 // val=34 for u="s" and mode=3
1553 // val=-19 for u="ns" and mode=3
1555 // The default is mode=1.
1557 return GetDifference(&t,u,mode);
1559 ///////////////////////////////////////////////////////////////////////////
1560 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)
1562 // Set the AliTimestamp parameters corresponding to the UT date and time
1563 // in the Gregorian calendar as specified by the input arguments.
1564 // This facility is exact upto picosecond precision and as such is
1565 // for scientific observations preferable above the corresponding
1566 // Set function(s) of TTimestamp.
1567 // The latter has a random spread in the sub-second part, which
1568 // might be of use in generating distinguishable timestamps while
1569 // still keeping second precision.
1571 // The input arguments represent the following :
1572 // y : year in UT (e.g. 1952, 2003 etc...)
1573 // m : month in UT (1=jan 2=feb etc...)
1574 // d : day in UT (1-31)
1575 // hh : elapsed hours in UT (0-23)
1576 // mm : elapsed minutes in UT (0-59)
1577 // ss : elapsed seconds in UT (0-59)
1578 // ns : remaining fractional elapsed second of UT in nanosecond
1579 // ps : remaining fractional elapsed nanosecond of UT in picosecond
1581 // Note : ns=0 and ps=0 are the default values.
1583 // This facility first determines the elapsed days, seconds etc...
1584 // since the beginning of the specified UT year on basis of the
1585 // input arguments. Subsequently it invokes the SetUT memberfunction
1586 // for the elapsed timespan.
1587 // As such this facility is valid for all AD dates in the Gregorian
1588 // calendar with picosecond precision.
1590 Int_t day=GetDayOfYear(d,m,y);
1591 Int_t secs=hh*3600+mm*60+ss;
1592 SetUT(y,day-1,secs,ns,ps);
1594 ///////////////////////////////////////////////////////////////////////////
1595 void AliTimestamp::SetUT(Int_t y,Int_t d,Int_t s,Int_t ns,Int_t ps)
1597 // Set the AliTimestamp parameters corresponding to the specified elapsed
1598 // timespan since the beginning of the new UT year.
1599 // This facility is exact upto picosecond precision and as such is
1600 // for scientific observations preferable above the corresponding
1601 // Set function(s) of TTimestamp.
1602 // The latter has a random spread in the sub-second part, which
1603 // might be of use in generating distinguishable timestamps while
1604 // still keeping second precision.
1606 // The UT year and elapsed time span is entered via the following input arguments :
1608 // y : year in UT (e.g. 1952, 2003 etc...)
1609 // d : elapsed number of days
1610 // s : (remaining) elapsed number of seconds
1611 // ns : (remaining) elapsed number of nanoseconds
1612 // ps : (remaining) elapsed number of picoseconds
1614 // The specified d, s, ns and ps values will be used in an additive
1615 // way to determine the elapsed timespan.
1616 // So, specification of d=1, s=100, ns=0, ps=0 will result in the
1617 // same elapsed time span as d=0, s=24*3600+100, ns=0, ps=0.
1618 // However, by making use of the latter the user should take care
1619 // of possible integer overflow problems in the input arguments,
1620 // which obviously will provide incorrect results.
1622 // Note : ns=0 and ps=0 are the default values.
1624 // This facility first sets the (M)JD corresponding to the start (01-jan 00:00:00)
1625 // of the specified UT year following the recipe of R.W. Sinnott
1626 // Sky & Telescope 82, (aug. 1991) 183.
1627 // Subsequently the day and (sub)second parts are added to the AliTimestamp.
1628 // As such this facility is valid for all AD dates in the Gregorian calendar.
1630 Double_t jd=GetJD(y,1,1,0,0,0,0);
1634 GetMJD(mjd,sec,nsec);
1638 ///////////////////////////////////////////////////////////////////////////
1639 void AliTimestamp::GetUT(Int_t& hh,Int_t& mm,Int_t& ss,Int_t& ns,Int_t& ps)
1641 // Provide the corrresponding UT as hh:mm:ss:ns:ps.
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 GetMJD(mjd,sec,nsec);
1657 ///////////////////////////////////////////////////////////////////////////
1658 Double_t AliTimestamp::GetUT()
1660 // Provide the corrresponding UT in fractional hours.
1661 // This facility is based on the MJD, so the TTimeStamp limitations
1662 // do not apply here.
1664 Int_t hh,mm,ss,ns,ps;
1666 GetUT(hh,mm,ss,ns,ps);
1668 Double_t ut=Convert(hh,mm,ss,ns,ps);
1672 ///////////////////////////////////////////////////////////////////////////
1673 void AliTimestamp::GetGMST(Int_t& hh,Int_t& mm,Int_t& ss,Int_t& ns,Int_t& ps)
1675 // Provide the corrresponding Greenwich Mean Sideral Time (GMST).
1676 // The algorithm used is the one described at p. 83 of the book
1677 // Astronomy Methods by Hale Bradt.
1678 // This facility is based on the MJD, so the TTimeStamp limitations
1679 // do not apply here.
1681 Int_t mjd,sec,nsec,psec;
1683 // The current UT based timestamp data
1684 GetMJD(mjd,sec,nsec);
1687 // The basis for the daily corrections in units of Julian centuries w.r.t. J2000.
1688 // Note : Epoch J2000 starts at 01-jan-2000 12:00:00 UT.
1689 Double_t tau=(GetJD()-2451545.)/36525.;
1691 // Syncronise sidereal time with current timestamp
1693 sid.SetMJD(mjd,sec,nsec,psec);
1695 // Add offset for GMST start value defined as 06:41:50.54841 at 01-jan 00:00:00 UT
1696 sec=6*3600+41*60+50;
1699 sid.Add(0,sec,nsec,psec);
1701 // Daily correction for precession and polar motion
1702 Double_t addsec=8640184.812866*tau+0.093104*pow(tau,2)-6.2e-6*pow(tau,3);
1704 addsec-=double(sec);
1705 nsec=int(addsec*1.e9);
1706 addsec-=double(nsec)*1.e-9;
1707 psec=int(addsec*1.e12);
1708 sid.Add(0,sec,nsec,psec);
1710 sid.GetMJD(mjd,sec,nsec);
1720 ///////////////////////////////////////////////////////////////////////////
1721 Double_t AliTimestamp::GetGMST()
1723 // Provide the corrresponding Greenwich Mean Sideral Time (GMST)
1724 // in fractional hours.
1725 // This facility is based on the MJD, so the TTimeStamp limitations
1726 // do not apply here.
1728 Int_t hh,mm,ss,ns,ps;
1730 GetGMST(hh,mm,ss,ns,ps);
1732 Double_t gst=Convert(hh,mm,ss,ns,ps);
1736 ///////////////////////////////////////////////////////////////////////////
1737 Double_t AliTimestamp::GetGAST()
1739 // Provide the corrresponding Greenwich Apparent Sideral Time (GAST)
1740 // in fractional hours.
1741 // In case a hh:mm:ss.sss format is needed, please invoke the Convert()
1742 // memberfunction for conversion of the provided fractional hour value.
1744 // The GAST is the GMST corrected for the shift of the vernal equinox
1745 // due to nutation. The right ascension component of the nutation correction
1746 // of the vernal equinox is called the "equation of the equinoxes".
1749 // GAST = GMST + (equation of the equinoxes)
1751 // The equation of the equinoxes is determined via the Almanac() memberfunction.
1753 // Since GMST is based on the MJD, the TTimeStamp limitations do not apply here.
1755 Double_t da=Almanac();
1757 // Convert to fractional hours
1760 Double_t gast=GetGMST()+da;
1773 ///////////////////////////////////////////////////////////////////////////
1774 Double_t AliTimestamp::GetLT(Double_t offset)
1776 // Provide the corresponding local time in fractional hours.
1777 // The "offset" denotes the time difference in (fractional) hours w.r.t. UT.
1778 // A mean solar day lasts 24h (i.e. 86400s).
1780 // In case a hh:mm:ss format is needed, please use the Convert() facility.
1782 // Current UT time in fractional hours
1798 ///////////////////////////////////////////////////////////////////////////
1799 Double_t AliTimestamp::GetLMST(Double_t offset)
1801 // Provide the corresponding Local Mean Sidereal Time (LMST) in fractional hours.
1802 // The "offset" denotes the time difference in (fractional) hours w.r.t. GMST.
1803 // A sidereal day corresponds to 23h 56m 04.09s (i.e. 86164.09s) mean solar time.
1804 // The definition of GMST is such that a sidereal clock corresponds with
1805 // 24 sidereal hours per revolution of the Earth.
1806 // As such, local time offsets w.r.t. UT and GMST can be treated similarly.
1808 // In case a hh:mm:ss format is needed, please use the Convert() facility.
1810 // Current GMST time in fractional hours
1811 Double_t h=GetGMST();
1826 ///////////////////////////////////////////////////////////////////////////
1827 Double_t AliTimestamp::GetLAST(Double_t offset)
1829 // Provide the corresponding Local Apparent Sidereal Time (LAST) in fractional hours.
1830 // The "offset" denotes the time difference in (fractional) hours w.r.t. GAST.
1831 // A sidereal day corresponds to 23h 56m 04.09s (i.e. 86164.09s) mean solar time.
1832 // The definition of GMST and GAST is such that a sidereal clock corresponds with
1833 // 24 sidereal hours per revolution of the Earth.
1834 // As such, local time offsets w.r.t. UT, GMST and GAST can be treated similarly.
1836 // In case a hh:mm:ss.sss format is needed, please use the Convert() facility.
1838 // Current GAST time in fractional hours
1839 Double_t h=GetGAST();
1854 ///////////////////////////////////////////////////////////////////////////
1855 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)
1857 // Set the AliTimestamp parameters corresponding to the LT date and time
1858 // in the Gregorian calendar as specified by the input arguments.
1859 // This facility is exact upto picosecond precision and as such is
1860 // for scientific observations preferable above the corresponding
1861 // Set function(s) of TTimestamp.
1862 // The latter has a random spread in the sub-second part, which
1863 // might be of use in generating distinguishable timestamps while
1864 // still keeping second precision.
1866 // The input arguments represent the following :
1868 // dt : the local time offset in fractional hours w.r.t. UT.
1869 // y : year in LT (e.g. 1952, 2003 etc...)
1870 // m : month in LT (1=jan 2=feb etc...)
1871 // d : day in LT (1-31)
1872 // hh : elapsed hours in LT (0-23)
1873 // mm : elapsed minutes in LT (0-59)
1874 // ss : elapsed seconds in LT (0-59)
1875 // ns : remaining fractional elapsed second of LT in nanosecond
1876 // ps : remaining fractional elapsed nanosecond of LT in picosecond
1878 // Note : ns=0 and ps=0 are the default values.
1880 // This facility first sets the UT as specified by the input arguments
1881 // and then corrects the UT by subtracting the local time offset w.r.t. UT.
1882 // As such this facility is valid for all AD dates in the Gregorian
1883 // calendar with picosecond precision.
1885 SetUT(y,m,d,hh,mm,ss,ns,ps);
1888 ///////////////////////////////////////////////////////////////////////////
1889 void AliTimestamp::SetLT(Double_t dt,Int_t y,Int_t d,Int_t s,Int_t ns,Int_t ps)
1891 // Set the AliTimestamp parameters corresponding to the specified elapsed
1892 // timespan since the beginning of the new LT year.
1893 // This facility is exact upto picosecond precision and as such is
1894 // for scientific observations preferable above the corresponding
1895 // Set function(s) of TTimestamp.
1896 // The latter has a random spread in the sub-second part, which
1897 // might be of use in generating distinguishable timestamps while
1898 // still keeping second precision.
1900 // The LT year and elapsed time span is entered via the following input arguments :
1902 // dt : the local time offset in fractional hours w.r.t. UT.
1903 // y : year in LT (e.g. 1952, 2003 etc...)
1904 // d : elapsed number of days
1905 // s : (remaining) elapsed number of seconds
1906 // ns : (remaining) elapsed number of nanoseconds
1907 // ps : (remaining) elapsed number of picoseconds
1909 // The specified d, s, ns and ps values will be used in an additive
1910 // way to determine the elapsed timespan.
1911 // So, specification of d=1, s=100, ns=0, ps=0 will result in the
1912 // same elapsed time span as d=0, s=24*3600+100, ns=0, ps=0.
1913 // However, by making use of the latter the user should take care
1914 // of possible integer overflow problems in the input arguments,
1915 // which obviously will provide incorrect results.
1917 // Note : ns=0 and ps=0 are the default values.
1919 // This facility first sets the UT as specified by the input arguments
1920 // and then corrects the UT by subtracting the local time offset w.r.t. UT.
1921 // As such this facility is valid for all AD dates in the Gregorian calendar.
1926 ///////////////////////////////////////////////////////////////////////////
1927 Double_t AliTimestamp::GetJD(Double_t e,TString mode) const
1929 // Provide the fractional Julian Date from epoch e.
1930 // The sort of epoch may be specified via the "mode" parameter.
1932 // mode = "J" ==> Julian epoch
1933 // "B" ==> Besselian epoch
1935 // The default value is mode="J".
1939 if (mode=="J" || mode=="j") jd=(e-2000.0)*365.25+2451545.0;
1941 if (mode=="B" || mode=="b") jd=(e-1900.0)*365.242198781+2415020.31352;
1945 ///////////////////////////////////////////////////////////////////////////
1946 Double_t AliTimestamp::GetMJD(Double_t e,TString mode) const
1948 // Provide the fractional Modified Julian Date from epoch e.
1949 // The sort of epoch may be specified via the "mode" parameter.
1951 // mode = "J" ==> Julian epoch
1952 // "B" ==> Besselian epoch
1954 // The default value is mode="J".
1956 Double_t mjd=GetJD(e,mode)-2400000.5;
1960 ///////////////////////////////////////////////////////////////////////////
1961 Double_t AliTimestamp::GetTJD(Double_t e,TString mode) const
1963 // Provide the fractional Truncated Julian Date from epoch e.
1964 // The sort of epoch may be specified via the "mode" parameter.
1966 // mode = "J" ==> Julian epoch
1967 // "B" ==> Besselian epoch
1969 // The default value is mode="J".
1971 Double_t tjd=GetJD(e,mode)-2440000.5;
1975 ///////////////////////////////////////////////////////////////////////////
1976 Double_t AliTimestamp::Almanac(Double_t* dpsi,Double_t* deps,Double_t* eps)
1978 // Determination of some astronomical observables which may be needed
1979 // for further calculations like e.g. precession of coordinates.
1981 // The standard returned value is the "equation of the equinoxes"
1982 // (i.e. the nutational shift of the RA of the vernal equinox) in seconds.
1983 // The memberfunction arguments provide the possibility of retrieving
1984 // optional returned values. The corresponding observables are :
1986 // dpsi : Nutational shift in ecliptic longitude in arcseconds
1987 // deps : Nutational shift in ecliptic obliquity in arcseconds
1988 // eps : Mean obliquity of the ecliptic in arcseconds
1990 // All shifts are determined for the current timestamp with
1991 // J2000.0 (i.e. 01-jan-2000 12:00:00 UT) as the reference epoch.
1993 // Invokation example :
1994 // --------------------
1996 // Double_t da,dpsi,deps,eps;
1997 // da=t.Almanac(&dpsi,&deps,&eps);
1999 // The nutation model used is the new one as documented in :
2000 // "The IAU Resolutions on Astronomical Reference Systems,
2001 // Time Scales and Earth Rotation Models".
2002 // This document is freely available as Circular 179 (2005) of the
2003 // United States Naval Observatory (USNO).
2004 // (See : http://aa.usno.navy.mil/publications/docs).
2006 // The change in ecliptic longitude (dpsi) and ecliptic obliquity (deps)
2007 // are evaluated using the IAU 2000A nutation series expansion
2008 // as provided in the USNO Circular 179.
2009 // The new expression for the equation of the equinoxes is based on a series
2010 // expansion and is the most accurate one known to date.
2011 // The components are documented on p.17 of the USNO Circular 179.
2013 // In the current implementation only the first 28 terms of the nutation series
2014 // are used. This provides an accuracy of about 0.01 arcsec corresponding to 0.001 sec.
2015 // In case a better accuracy is required, the series can be extended.
2016 // The total series expansion consists of 1365 terms.
2018 // Since all calculations are based on the JD, the TTimeStamp limitations
2019 // do not apply here.
2021 Double_t pi=acos(-1.);
2023 Double_t t; // Time difference in fractional Julian centuries w.r.t. the start of J2000.
2024 Double_t epsilon; // Mean obliquity of the ecliptic
2025 Double_t l; // Mean anomaly of the Moon
2026 Double_t lp; // Mean anomaly of the Sun
2027 Double_t f; // Mean argument of latitude of the moon
2028 Double_t d; // Mean elongation of the Moon from the Sun
2029 Double_t om; // Mean longitude of the Moon's mean ascending mode
2031 t=(GetJD()-2451545.0)/36525.;
2033 // Values of epsilon and the fundamental luni-solar arguments in arcseconds
2034 epsilon=84381.406-46.836769*t-0.0001831*pow(t,2)+0.00200340*pow(t,3)
2035 -0.000000576*pow(t,4)-0.0000000434*pow(t,5);
2036 l=485868.249036+1717915923.2178*t+31.8792*pow(t,2)+0.051635*pow(t,3)-0.00024470*pow(t,4);
2037 lp=1287104.79305+129596581.0481*t-0.5532*pow(t,2)+0.000136*pow(t,3)-0.00001149*pow(t,4);
2038 f=335779.526232+1739527262.8478*t-12.7512*pow(t,2)-0.001037*pow(t,3)+0.00000417*pow(t,4);
2039 d=1072260.70369+1602961601.2090*t-6.3706*pow(t,2)+0.006593*pow(t,3)-0.00003169*pow(t,4);
2040 om=450160.398036-6962890.5431*t+7.4722*pow(t,2)+0.007702*pow(t,3)-0.00005939*pow(t,4);
2042 if (eps) *eps=epsilon;
2044 // Convert to radians
2045 epsilon=epsilon*pi/(180.*3600.);
2046 f=f*pi/(180.*3600.);
2047 d=d*pi/(180.*3600.);
2048 l=l*pi/(180.*3600.);
2049 lp=lp*pi/(180.*3600.);
2050 om=om*pi/(180.*3600.);
2052 //The IAU 2000A nutation series expansion.
2053 Double_t phi[28]={om,2.*(f-d+om),2.*(f+om),2.*om,lp,lp+2.*(f-d+om),l,
2054 2.*f+om,l+2.*(f+om),2.*(f-d+om)-lp,2.*(f-d)+om,2.*(f+om)-l,2.*d-l,l+om,
2055 om-l,2.*(f+d+om)-l,l+2.*f+om,2.*(f-l)+om,2.*d,2.*(f+d+om),2.*(f-d+om-lp),
2056 2.*(d-l),2.*(l+d+om),l+2.*(f-d+om),2.*f+om-l,2.*l,2.*f,lp+om};
2057 Double_t s[28]={-17.2064161,-1.3170907,-0.2276413, 0.2074554, 0.1475877,-0.0516821, 0.0711159,
2058 -0.0387298,-0.0301461, 0.0215829, 0.0128227, 0.0123457, 0.0156994, 0.0063110,
2059 -0.0057976,-0.0059641,-0.0051613, 0.0045893, 0.0063384,-0.0038571, 0.0032481,
2060 -0.0047722,-0.0031046, 0.0028593, 0.0020441, 0.0029243, 0.0025887,-0.0014053};
2061 Double_t sd[28]={-0.0174666,-0.0001675,-0.0000234, 0.0000207,-0.0003633, 0.0001226, 0.0000073,
2062 -0.0000367,-0.0000036,-0.0000494, 0.0000137, 0.0000011, 0.0000010, 0.0000063,
2063 -0.0000063,-0.0000011,-0.0000042, 0.0000050, 0.0000011,-0.0000001, 0.0000000,
2064 0.0000000,-0.0000001, 0.0000000, 0.0000021, 0.0000000, 0.0000000,-0.0000025};
2065 Double_t cp[28]={ 0.0033386,-0.0013696, 0.0002796,-0.0000698, 0.0011817,-0.0000524,-0.0000872,
2066 0.0000380, 0.0000816, 0.0000111, 0.0000181, 0.0000019,-0.0000168, 0.0000027,
2067 -0.0000189, 0.0000149, 0.0000129, 0.0000031,-0.0000150, 0.0000158, 0.0000000,
2068 -0.0000018, 0.0000131,-0.0000001, 0.0000010,-0.0000074,-0.0000066, 0.0000079};
2069 Double_t c[28]= { 9.2052331, 0.5730336, 0.0978459,-0.0897492, 0.0073871, 0.0224386,-0.0006750,
2070 0.0200728, 0.0129025,-0.0095929,-0.0068982,-0.0053311,-0.0001235,-0.0033228,
2071 0.0031429, 0.0025543, 0.0026366,-0.0024236,-0.0001220, 0.0016452,-0.0013870,
2072 0.0000477, 0.0013238,-0.0012338,-0.0010758,-0.0000609,-0.0000550, 0.0008551};
2073 Double_t cd[28]={ 0.0009086,-0.0003015,-0.0000485, 0.0000470,-0.0000184,-0.0000677, 0.0000000,
2074 0.0000018,-0.0000063, 0.0000299,-0.0000009, 0.0000032, 0.0000000, 0.0000000,
2075 0.0000000,-0.0000011, 0.0000000,-0.0000010, 0.0000000,-0.0000011, 0.0000000,
2076 0.0000000,-0.0000011, 0.0000010, 0.0000000, 0.0000000, 0.0000000,-0.0000002};
2077 Double_t sp[28]={ 0.0015377,-0.0004587, 0.0001374,-0.0000291,-0.0001924,-0.0000174, 0.0000358,
2078 0.0000318, 0.0000367, 0.0000132, 0.0000039,-0.0000004, 0.0000082,-0.0000009,
2079 -0.0000075, 0.0000066, 0.0000078, 0.0000020, 0.0000029, 0.0000068, 0.0000000,
2080 -0.0000025, 0.0000059,-0.0000003,-0.0000003, 0.0000013, 0.0000011,-0.0000045};
2082 Double_t dp=0,de=0,da=0;
2083 for (Int_t i=0; i<28; i++)
2085 dp+=(s[i]+sd[i]*t)*sin(phi[i])+cp[i]*cos(phi[i]);
2086 de+=(c[i]+cd[i]*t)*cos(phi[i])+sp[i]*sin(phi[i]);
2089 da=dp*cos(epsilon)+0.00264096*sin(om)+0.00006352*sin(2.*om)
2090 +0.00001175*sin(2.*f-2.*d+3.*om)+0.00001121*sin(2.*f-2.*d+om)
2091 -0.00000455*sin(2.*f-2.*d+2.*om)+0.00000202*sin(2.*f+3.*om)+0.00000198*sin(2.*f+om)
2092 -0.00000172*sin(3.*om)-0.00000087*t*sin(om);
2097 // Convert to seconds
2102 ///////////////////////////////////////////////////////////////////////////
2103 void AliTimestamp::SetEpoch(Double_t e,TString mode)
2105 // Set the timestamp parameters according to the epoch as specified by
2106 // the input argument "e".
2107 // Via the input argument "mode" the user can specify the type of epoch
2109 // mode = "B" ==> Besselian epoch
2110 // "J" ==> Julian epoch
2112 Double_t jd=GetJD(e,mode);
2115 ///////////////////////////////////////////////////////////////////////////
2116 Double_t AliTimestamp::GetEpoch(TString mode)
2118 // Provide the corresponding epoch value.
2119 // Via the input argument "mode" the user can specify the type of epoch
2121 // mode = "B" ==> Besselian epoch
2122 // "J" ==> Julian epoch
2125 if (mode=="B" || mode=="b") e=GetBE();
2126 if (mode=="J" || mode=="j") e=GetJE();
2129 ///////////////////////////////////////////////////////////////////////////