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 // Virtual lab to correlate measurements with astrophysical phenomena.
22 // This class is derived from TTask, but the only reason for this
23 // is to enable this class to serve as a base class for other TTask
24 // derived classes (e.g. AliEventSelector) without the need for
25 // multiple virtual inheritance.
26 // So, this AliAstrolab class itself does not provide any TTask
27 // related functionality.
29 // The lab can be given a terrestrial location via the usual longitude
30 // and latitude specifications.
31 // Since this class is derived from AliTimestamp, a lab can also be
32 // given a specific timestamp. Together with the terrestrial location
33 // this provides access to local (sidereal) times etc...
34 // In addition to the usual astronomical reference frames, a local
35 // lab reference frame can also be specified. Together with the lab's
36 // timestamp this uniquely defines all the coordinate transformations
37 // between the various reference frames.
38 // These lab characteristics allow a space-time correlation of lab
39 // observations with external (astrophysical) phenomena.
41 // Observations are entered as generic signals containing a position,
42 // reference frame specification and a timestamp.
43 // These observations can then be analysed in various reference frames
44 // via the available GET functions.
46 // Various external (astrophysical) phenomena may be entered as
47 // so-called reference signals.
48 // This class provides facilities (e.g. MatchRefSignal) to check
49 // correlations of the stored measurement with these reference signals.
50 // Also graphical facilities (e.g. DisplaySignals) are available to
51 // provide skymaps in various projections.
55 // gSystem->Load("ralice");
57 // AliAstrolab lab("IceCube","The South Pole Neutrino Observatory");
58 // lab.SetLabPosition(0,-90,"deg"); // South Pole
60 // lab.SetLocalFrame(90,180,90,270,0,0); // Local frame has X-axis to the North
62 // lab.Data(1,"dms"); // Print laboratory parameters
64 // // Enter some observed event to be investigated
66 // ts.SetUT(1989,7,30,8,14,23,738504,0);
67 // Float_t vec[3]={1,23.8,118.65};
69 // r.SetVector(vec,"sph","deg");
70 // lab.SetSignal(&r,"loc","M",&ts,0,"Event10372");
72 // // Enter some reference signals
73 // Float_t alpha=194818.0;
74 // Float_t delta=84400.;
75 // lab.SetSignal(alpha,delta,"B",1950,"M",-1,"Altair");
78 // lab.SetSignal(alpha,delta,"B",1950,"M",-1,"NGP");
81 // lab.SetSignal(alpha,delta,"J",2000,"M",-1,"Sirius");
84 // lab.SetSignal(alpha,delta,"J",2000,"M",-1,"Polaris");
87 // lab.SetSignal(alpha,delta,"J",2000,"M",-1,"Aldebaran");
89 // Float_t b=-35.8903;
90 // Float_t pos[3]={1,90.-b,l};
91 // r.SetVector(pos,"sph","deg");
92 // lab.SetUT(1989,7,30,8,14,16,0,0);
93 // lab.SetSignal(&r,"gal","M",0,-1,"GRB890730");
95 // // List all stored objects
96 // lab.ListSignals("equ","M",5);
98 // // Determine minimal space and time differences with reference signals
101 // da=lab.GetDifference(0,"deg",dt,"s",1,&ia,&it);
102 // cout << " Minimal differences damin (deg) : " << da << " dtmin (s) : " << dt
103 // << " damin-index : " << ia << " dtmin-index : " << it << endl;
104 // cout << " damin for "; lab->PrintSignal("equ","T",&ts,5,ia); cout << endl;
105 // cout << " dtmin for "; lab->PrintSignal("equ","T",&ts,5,it); cout << endl;
107 // // Search for space and time match with the reference signals
110 // TArrayI* arr=lab.MatchRefSignal(da,"deg",dt,"s");
114 // for (Int_t i=0; i<arr->GetSize(); i++)
117 // cout << " Match found for index : " << index << endl;
118 // cout << " Corresponding ref. object "; lab->PrintSignal("equ","T",&ts,5,index); cout << endl;
122 // // Display all stored objects in a skymap (Hammer projection)
123 // lab.DisplaySignals("equ","M",&ts,"ham",1);
125 //--- Author: Nick van Eijndhoven 15-mar-2007 Utrecht University
126 //- Modified: NvE $Date$ Utrecht University
127 ///////////////////////////////////////////////////////////////////////////
129 #include "AliAstrolab.h"
130 #include "Riostream.h"
132 ClassImp(AliAstrolab) // Class implementation to enable ROOT I/O
134 AliAstrolab::AliAstrolab(const char* name,const char* title) : TTask(name,title),AliTimestamp()
136 // Default constructor
150 ///////////////////////////////////////////////////////////////////////////
151 AliAstrolab::~AliAstrolab()
153 // Destructor to delete all allocated memory.
186 ///////////////////////////////////////////////////////////////////////////
187 AliAstrolab::AliAstrolab(const AliAstrolab& t) : TTask(t),AliTimestamp(t)
194 if (t.fXsig) fXsig=new AliSignal(*(t.fXsig));
198 Int_t size=t.fRefs->GetSize();
199 fRefs=new TObjArray(size);
200 for (Int_t i=0; i<size; i++)
202 AliSignal* sx=(AliSignal*)t.fRefs->At(i);
203 if (sx) fRefs->AddAt(sx->Clone(),i);
215 ///////////////////////////////////////////////////////////////////////////
216 void AliAstrolab::Data(Int_t mode,TString u)
218 // Provide lab information.
220 // "mode" indicates the mode of the timestamp info (see AliTimestamp::Date).
222 // The string argument "u" allows to choose between different angular units
223 // in case e.g. a spherical frame is selected.
224 // u = "rad" : angles provided in radians
225 // "deg" : angles provided in degrees
226 // "dms" : angles provided in ddd:mm:ss.sss
227 // "hms" : angles provided in hh:mm:ss.sss
229 // The defaults are mode=1 and u="deg".
231 const char* name=GetName();
232 const char* title=GetTitle();
233 cout << " *" << ClassName() << "::Data*";
234 if (strlen(name)) cout << " Name : " << GetName();
235 if (strlen(title)) cout << " Title : " << GetTitle();
239 GetLabPosition(l,b,"deg");
240 cout << " Lab position longitude : "; PrintAngle(l,"deg",u,2);
241 cout << " latitude : "; PrintAngle(b,"deg",u,2);
243 cout << " Lab time offset w.r.t. UT : "; PrintTime(fToffset,12); cout << endl;
245 // UT and Local time info
248 ///////////////////////////////////////////////////////////////////////////
249 void AliAstrolab::SetLabPosition(Ali3Vector& p)
251 // Set the lab position in the terrestrial coordinates.
252 // The right handed reference frame is defined such that the North Pole
253 // corresponds to a polar angle theta=0 and the Greenwich meridian corresponds
254 // to an azimuth angle phi=0, with phi increasing eastwards.
256 fLabPos.SetPosition(p);
258 // Determine local time offset in fractional hours w.r.t. UT.
260 p.GetVector(vec,"sph","deg");
264 ///////////////////////////////////////////////////////////////////////////
265 void AliAstrolab::SetLabPosition(Double_t l,Double_t b,TString u)
267 // Set the lab position in the terrestrial longitude (l) and latitude (b).
268 // Positions north of the equator have b>0, whereas b<0 indicates
269 // locations south of the equator.
270 // Positions east of Greenwich have l>0, whereas l<0 indicates
271 // locations west of Greenwich.
273 // The string argument "u" allows to choose between different angular units
274 // u = "rad" : angles provided in radians
275 // "deg" : angles provided in degrees
276 // "dms" : angles provided in dddmmss.sss
277 // "hms" : angles provided in hhmmss.sss
279 // The default is u="deg".
281 Double_t r=1,theta=0,phi=0;
283 l=ConvertAngle(l,u,"deg");
284 b=ConvertAngle(b,u,"deg");
291 Double_t p[3]={r,theta,phi};
292 fLabPos.SetPosition(p,"sph","deg");
294 // Local time offset in fractional hours w.r.t. UT.
297 ///////////////////////////////////////////////////////////////////////////
298 AliPosition AliAstrolab::GetLabPosition() const
300 // Provide the lab position in the terrestrial coordinates.
301 // The right handed reference frame is defined such that the North Pole
302 // corresponds to a polar angle theta=0 and the Greenwich meridian corresponds
303 // to an azimuth angle phi=0, with phi increasing eastwards.
307 ///////////////////////////////////////////////////////////////////////////
308 void AliAstrolab::GetLabPosition(Double_t& l,Double_t& b,TString u) const
310 // Provide the lab position in the terrestrial longitude (l) and latitude (b).
311 // Positions north of the equator have b>0, whereas b<0 indicates
312 // locations south of the equator.
313 // Positions east of Greenwich have l>0, whereas l<0 indicates
314 // locations west of Greenwich.
316 // The string argument "u" allows to choose between different angular units
317 // u = "rad" : angles provided in radians
318 // "deg" : angles provided in degrees
320 // The default is u="deg".
322 Double_t pi=acos(-1.);
325 if (u=="rad") offset=pi/2.;
328 fLabPos.GetPosition(p,"sph",u);
332 ///////////////////////////////////////////////////////////////////////////
333 Double_t AliAstrolab::GetLT()
335 // Provide the Lab's local time in fractional hours.
336 // A mean solar day lasts 24h (i.e. 86400s).
338 // In case a hh:mm:ss format is needed, please use the Convert() facility.
340 Double_t h=GetLT(fToffset);
343 ///////////////////////////////////////////////////////////////////////////
344 Double_t AliAstrolab::GetLMST()
346 // Provide the Lab's Local Mean Sidereal Time (LMST) in fractional hours.
347 // A sidereal day corresponds to 23h 56m 04.09s (i.e. 86164.09s) mean solar time.
348 // The definition of GMST is such that a sidereal clock corresponds with
349 // 24 sidereal hours per revolution of the Earth.
350 // As such, local time offsets w.r.t. UT and GMST can be treated similarly.
352 // In case a hh:mm:ss format is needed, please use the Convert() facility.
354 Double_t h=GetLMST(fToffset);
357 ///////////////////////////////////////////////////////////////////////////
358 Double_t AliAstrolab::GetLAST()
360 // Provide the Lab's Local Apparent Sidereal Time (LAST) in fractional hours.
361 // A sidereal day corresponds to 23h 56m 04.09s (i.e. 86164.09s) mean solar time.
362 // The definition of GMST and GAST is such that a sidereal clock corresponds with
363 // 24 sidereal hours per revolution of the Earth.
364 // As such, local time offsets w.r.t. UT, GMST and GAST can be treated similarly.
366 // In case a hh:mm:ss format is needed, please use the Convert() facility.
368 Double_t h=GetLAST(fToffset);
371 ///////////////////////////////////////////////////////////////////////////
372 void AliAstrolab::PrintAngle(Double_t a,TString in,TString out,Int_t ndig) const
374 // Printing of angles in various formats.
376 // The input argument "a" denotes the angle to be printed.
377 // The string arguments "in" and "out" specify the angular I/O formats.
379 // in = "rad" : input angle provided in radians
380 // "deg" : input angle provided in degrees
381 // "dms" : input angle provided in dddmmss.sss
382 // "hms" : input angle provided in hhmmss.sss
384 // out = "rad" : output angle provided in radians
385 // "deg" : output angle provided in degrees
386 // "dms" : output angle provided in dddmmss.sss
387 // "hms" : output angle provided in hhmmss.sss
389 // The argument "ndig" specifies the number of digits for the fractional
390 // part (e.g. ndig=6 for "dms" corresponds to micro-arcsecond precision).
391 // No rounding will be performed, so an arcsecond count of 3.473 with ndig=1
392 // will appear as 03.4 on the output.
393 // Due to computer accuracy, precision on the pico-arcsecond level may get lost.
395 // The default is ndig=1.
397 // Note : The angle info is printed without additional spaces or "endline".
398 // This allows the print to be included in various composite output formats.
400 Double_t b=ConvertAngle(a,in,out);
402 if (out=="deg" || out=="rad")
404 cout.setf(ios::fixed,ios::floatfield);
405 cout << setprecision(ndig) << b << " " << out.Data();
406 cout.unsetf(ios::fixed);
410 Double_t epsilon=1.e-12; // Accuracy in (arc)seconds
411 Int_t word=0,ddd=0,hh=0,mm=0,ss=0;
423 s=fabs(b)-Double_t(ddd*10000+mm*100+ss);
445 if (b<0) cout << "-";
446 cout << ddd << "d " << mm << "' " << ss << "."
447 << setfill('0') << setw(ndig) << sfrac << "\"";
459 s=fabs(b)-Double_t(hh*10000+mm*100+ss);
481 if (b<0) cout << "-";
482 cout << hh << "h " << mm << "m " << ss << "."
483 << setfill('0') << setw(ndig) << sfrac << "s";
487 ///////////////////////////////////////////////////////////////////////////
488 void AliAstrolab::SetSignal(Ali3Vector* r,TString frame,TString mode,AliTimestamp* ts,Int_t jref,TString name)
490 // Store a signal as specified by the position r and the timestamp ts.
491 // The position is stored in International Celestial Reference System (ICRS) coordinates.
492 // The ICRS is a fixed, time independent frame and as such provides a unique reference
493 // frame without the need of specifying any epoch etc...
494 // The ICRS coordinate definitions match within 20 mas with the mean ones of the J2000.0
495 // equatorial system. Nevertheless, to obtain the highest accuracy, the slight
496 // coordinate correction between J2000 and ICRS is performed here via the
497 // so-called frame bias matrix.
498 // For further details see the U.S. Naval Observatory (USNO) circular 179 (2005),
499 // which is available on http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
501 // The input parameter "frame" allows the user to specify the frame to which
502 // the components of r refer. Available options are :
504 // frame = "equ" ==> Equatorial coordinates with right ascension (a) and declination (d),
505 // where the "sph" components of r correspond to theta=(pi/2)-d and phi=a.
506 // "gal" ==> Galactic coordinates with longitude (l) and latitude (b).
507 // where the "sph" components of r correspond to theta=(pi/2)-b and phi=l.
508 // "ecl" ==> Ecliptic coordinates with longitude (l) and latitude (b),
509 // where the "sph" components of r correspond to theta=(pi/2)-b and phi=l.
510 // "hor" ==> Horizontal coordinates at the AliAstrolab location, where the "sph"
511 // components of r correspond to theta=zenith angle and phi=pi-azimuth.
512 // "icr" ==> ICRS coordinates with longitude (l) and latitude (b),
513 // where the "sph" components of r correspond to theta=(pi/2)-b and phi=l.
514 // "loc" ==> Local coordinates at the AliAstrolab location, where the "sph"
515 // components of r correspond to the usual theta and phi angles.
517 // In case the coordinates are the equatorial right ascension and declination (a,d),
518 // they can represent so-called "mean" and "true" values.
519 // The distinction between these two representations is the following :
521 // mean values : (a,d) are only corrected for precession and not for nutation
522 // true values : (a,d) are corrected for both precession and nutation
524 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
525 // values for the input in case of equatorial (a,d) coordinates.
527 // mode = "M" --> Input coordinates are the mean values
528 // "T" --> Input coordinates are the true values
530 // The input parameter "jref" allows the user to store so-called "reference" signals.
531 // These reference signals may be used to check space-time event coincidences with the
532 // stored measurement (e.g. coincidence of the measurement with transient phenomena).
534 // jref = 0 --> Storage of the measurement
535 // j --> Storage of a reference signal at the j-th position (j=1 is first)
536 // < 0 --> Add a reference signal at the next available position
538 // Via the input argument "name" the user can give the stored signal also a name.
540 // The default values are jref=0 and name="".
542 // Note : In case ts=0 the current timestamp of the lab will be taken.
554 if (frame!="equ" && frame!="gal" && frame!="ecl" && frame!="hor" && frame!="icr" && frame!="loc")
564 if (frame=="equ" && mode!="M" && mode!="m" && mode!="T" && mode!="t")
574 if (!ts) ts=(AliTimestamp*)this;
576 Double_t vec[3]={1,0,0};
577 vec[1]=r->GetX(2,"sph","rad");
578 vec[2]=r->GetX(3,"sph","rad");
580 q.SetVector(vec,"sph","rad");
584 if (!jref) // Storage of the measurement
588 fXsig=new AliSignal();
594 if (name != "") fXsig->SetName(name);
595 fXsig->SetTitle("Event in ICRS coordinates");
596 fXsig->SetTimestamp(*ts);
598 else // Storage of a reference signal
602 fRefs=new TObjArray();
605 // Expand array size if needed
606 if (jref>0 && jref>=fRefs->GetSize()) fRefs->Expand(jref+1);
607 sxref=new AliSignal();
608 if (name != "") sxref->SetName(name);
609 sxref->SetTitle("Reference event in ICRS coordinates");
610 sxref->SetTimestamp(*ts);
615 // Convert to horizontal coordinates
616 q=q.GetUnprimed(&fL);
619 SetSignal(&q,"hor",mode,ts,jref);
625 // Convert to "mean" values at specified epoch
626 if (mode=="T" || mode=="t")
629 q=q.GetUnprimed(&fN);
632 // Convert to "mean" values at J2000
634 q=q.GetUnprimed(&fP);
636 // Convert to ICRS values
637 if (!fBias) SetBmatrix();
638 q=q.GetUnprimed(&fB);
643 // Convert to J2000 equatorial mean coordinates
644 if (fGal != 2) SetGmatrix("J");
645 q=q.GetUnprimed(&fG);
647 // Convert to ICRS values
648 if (!fBias) SetBmatrix();
649 q=q.GetUnprimed(&fB);
654 // Convert to mean equatorial values at specified epoch
656 q=q.GetUnprimed(&fE);
658 // Convert to "mean" values at J2000
660 q=q.GetUnprimed(&fP);
662 // Convert to ICRS values
663 if (!fBias) SetBmatrix();
664 q=q.GetUnprimed(&fB);
669 // Convert to "true" equatorial values at the specified timestamp
671 q=q.GetUnprimed(&fH);
673 // Convert to "mean" values at specified timestamp
675 q=q.GetUnprimed(&fN);
677 // Convert to "mean" values at J2000
679 q=q.GetUnprimed(&fP);
681 // Convert to ICRS values
682 if (!fBias) SetBmatrix();
683 q=q.GetUnprimed(&fB);
686 // Store the signal in ICRS coordinates
687 if (!jref) // Storage of a regular signal
689 fXsig->SetPosition(q);
691 else // Storage of a reference signal
693 sxref->SetPosition(q);
700 fRefs->AddAt(sxref,jref-1);
704 ///////////////////////////////////////////////////////////////////////////
705 void AliAstrolab::SetSignal(Double_t a,Double_t d,TString s,Double_t e,TString mode,Int_t jref,TString name)
707 // Store a signal with right ascension (a) and declination (d) given for epoch e.
708 // The position is stored in International Celestial Reference System (ICRS) coordinates.
709 // The ICRS is a fixed, time independent frame and as such provides a unique reference
710 // frame without the need of specifying any epoch etc...
711 // The ICRS coordinate definitions match within 20 mas the mean ones of the J2000.0
712 // equatorial system. Nevertheless, to obtain the highest accuracy, the slight
713 // coordinate correction between J2000 and ICRS is performed here via the
714 // so-called frame bias matrix.
715 // For further details see the U.S. Naval Observatory (USNO) circular 179 (2005),
716 // which is available on http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
718 // The coordinates (a,d) can represent so-called "mean" and "true" values.
719 // The distinction between these two representations is the following :
721 // mean values : (a,d) are only corrected for precession and not for nutation
722 // true values : (a,d) are corrected for both precession and nutation
724 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
725 // values for the input (a,d) coordinates.
727 // a : Right ascension in hhmmss.sss
728 // d : Declination in dddmmss.sss
729 // s = "B" --> Besselian reference epoch.
730 // "J" --> Julian reference epoch.
731 // e : Reference epoch for the input coordinates (e.g. 1900, 1950, 2000,...)
732 // mode = "M" --> Input coordinates are the mean values
733 // "T" --> Input coordinates are the true values
735 // The input parameter "jref" allows the user to store so-called "reference" signals.
736 // These reference signals may be used to check space-time event coincidences with the
737 // stored measurement (e.g. coincidence of the measurement with transient phenomena).
739 // jref = 0 --> Storage of the measurement
740 // j --> Storage of a reference signal at the j-th position (j=1 is first)
741 // < 0 --> Add a reference signal at the next available position
743 // Via the input argument "name" the user can give the stored signal also a name.
745 // The default values are jref=0 and name="".
747 if (s!="B" && s!="b" && s!="J" && s!="j")
757 if (mode!="M" && mode!="m" && mode!="T" && mode!="t")
767 // Convert coordinates to fractional degrees.
768 a=ConvertAngle(a,"hms","deg");
769 d=ConvertAngle(d,"dms","deg");
776 Double_t vec[3]={1.,90.-d,a};
777 r.SetVector(vec,"sph","deg");
779 SetSignal(&r,"equ",mode,&tx,jref,name);
781 ///////////////////////////////////////////////////////////////////////////
782 void AliAstrolab::SetSignal(Double_t a,Double_t d,TString mode,AliTimestamp* ts,Int_t jref,TString name)
784 // Store a signal with right ascension (a) and declination (d) given for timestamp ts.
785 // The position is stored in International Celestial Reference System (ICRS) coordinates.
786 // The ICRS is a fixed, time independent frame and as such provides a unique reference
787 // frame without the need of specifying any epoch etc...
788 // The ICRS coordinate definitions match within 20 mas the mean ones of the J2000.0
789 // equatorial system. Nevertheless, to obtain the highest accuracy, the slight
790 // coordinate correction between J2000 and ICRS is performed here via the
791 // so-called frame bias matrix.
792 // For further details see the U.S. Naval Observatory (USNO) circular 179 (2005),
793 // which is available on http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
795 // The coordinates (a,d) can represent so-called "mean" and "true" values.
796 // The distinction between these two representations is the following :
798 // mean values : (a,d) are only corrected for precession and not for nutation
799 // true values : (a,d) are corrected for both precession and nutation
801 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
802 // values for the input (a,d) coordinates.
804 // a : Right ascension in hhmmss.sss
805 // d : Declination in dddmmss.sss
806 // mode = "M" --> Input coordinates are the mean values
807 // "T" --> Input coordinates are the true values
809 // The input parameter "jref" allows the user to store so-called "reference" signals.
810 // These reference signals may be used to check space-time event coincidences with the
811 // stored measurement (e.g. coincidence of the measurement with transient phenomena).
813 // jref = 0 --> Storage of the measurement
814 // j --> Storage of a reference signal at the j-th position (j=1 is first)
815 // < 0 --> Add a reference signal at the next available position
817 // Via the input argument "name" the user can give the stored signal also a name.
819 // The default values are jref=0 and name="".
821 // Note : In case ts=0 the current timestamp of the lab will be taken.
823 // Convert coordinates to fractional degrees.
824 a=ConvertAngle(a,"hms","deg");
825 d=ConvertAngle(d,"dms","deg");
828 Double_t vec[3]={1.,90.-d,a};
829 r.SetVector(vec,"sph","deg");
831 SetSignal(&r,"equ",mode,ts,jref,name);
833 ///////////////////////////////////////////////////////////////////////////
834 AliSignal* AliAstrolab::GetSignal(Ali3Vector& r,TString frame,TString mode,AliTimestamp* ts,Int_t jref)
836 // Provide the user specified coordinates of a signal at the specific timestamp ts.
837 // The coordinates are returned via the vector argument "r".
838 // In addition also a pointer to the stored signal object is provided.
839 // In case no stored signal was available or one of the input arguments was
840 // invalid, the returned pointer will be 0.
842 // The input parameter "frame" allows the user to specify the frame to which
843 // the components of r refer. Available options are :
845 // frame = "equ" ==> Equatorial coordinates with right ascension (a) and declination (d),
846 // where the "sph" components of r correspond to theta=(pi/2)-d and phi=a.
847 // "gal" ==> Galactic coordinates with longitude (l) and latitude (b).
848 // where the "sph" components of r correspond to theta=(pi/2)-b and phi=l.
849 // "ecl" ==> Ecliptic coordinates with longitude (l) and latitude (b),
850 // where the "sph" components of r correspond to theta=(pi/2)-b and phi=l.
851 // "hor" ==> Horizontal coordinates at the AliAstrolab location, where the "sph"
852 // components of r correspond to theta=zenith angle and phi=pi-azimuth.
853 // "icr" ==> ICRS coordinates with longitude (l) and latitude (b),
854 // where the "sph" components of r correspond to theta=(pi/2)-b and phi=l.
855 // "loc" ==> Local coordinates at the AliAstrolab location, where the "sph"
856 // components of r correspond to the usual theta and phi angles.
858 // In case the coordinates are the equatorial right ascension and declination (a,d),
859 // they can represent so-called "mean" and "true" values.
860 // The distinction between these two representations is the following :
862 // mean values : (a,d) are only corrected for precession and not for nutation
863 // true values : (a,d) are corrected for both precession and nutation
865 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
866 // values for the input in case of equatorial (a,d) coordinates.
868 // mode = "M" --> Input coordinates are the mean values
869 // "T" --> Input coordinates are the true values
871 // The input parameter "jref" allows the user to access so-called "reference" signals.
872 // These reference signals may be used to check space-time event coincidences with the
873 // stored measurement (e.g. coincidence of the measurement with transient phenomena).
875 // jref = 0 --> Access to the measurement
876 // j --> Access to the reference signal at the j-th position (j=1 is first)
878 // Default value is jref=0.
880 // Note : In case ts=0 the current timestamp of the lab will be taken.
884 if (frame!="equ" && frame!="gal" && frame!="ecl" && frame!="hor" && frame!="icr" && frame!="loc") return 0;
886 if (frame=="equ" && mode!="M" && mode!="m" && mode!="T" && mode!="t") return 0;
888 AliSignal* sx=GetSignal(jref);
892 if (!ts) ts=(AliTimestamp*)this;
895 sx->GetPosition(vec,"sph","rad");
897 q.SetVector(vec,"sph","rad");
905 // Convert from ICRS to equatorial J2000 coordinates
906 if (!fBias) SetBmatrix();
911 // Precess to specified timestamp
913 ts1.SetEpoch(2000,"J");
916 // Nutation correction if requested
917 if (mode=="T" || mode=="t") Nutate(q,ts);
922 // Convert from equatorial J2000 to galactic
923 if (fGal != 2) SetGmatrix("J");
929 // Precess to specified timestamp
931 ts1.SetEpoch(2000,"J");
934 // Convert from equatorial to ecliptic coordinates
941 // Precess to specified timestamp
943 ts1.SetEpoch(2000,"J");
946 // Nutation correction
949 // Convert from equatorial to horizontal coordinates
956 // Get the signal in horizontal coordinates
957 GetSignal(q,"hor",mode,ts);
959 // Convert from horizontal local-frame coordinates
966 ///////////////////////////////////////////////////////////////////////////
967 AliSignal* AliAstrolab::GetSignal(Ali3Vector& r,TString frame,TString mode,AliTimestamp* ts,TString name)
969 // Provide the user specified coordinates of the signal with the specified
970 // name at the specific timestamp ts.
971 // The coordinates are returned via the vector argument "r".
972 // In addition also a pointer to the stored signal object is provided.
973 // In case no such stored signal was available or one of the input arguments was
974 // invalid, the returned pointer will be 0.
976 // The input parameter "frame" allows the user to specify the frame to which
977 // the components of r refer. Available options are :
979 // frame = "equ" ==> Equatorial coordinates with right ascension (a) and declination (d),
980 // where the "sph" components of r correspond to theta=(pi/2)-d and phi=a.
981 // "gal" ==> Galactic coordinates with longitude (l) and latitude (b).
982 // where the "sph" components of r correspond to theta=(pi/2)-b and phi=l.
983 // "ecl" ==> Ecliptic coordinates with longitude (l) and latitude (b),
984 // where the "sph" components of r correspond to theta=(pi/2)-b and phi=l.
985 // "hor" ==> Horizontal coordinates at the AliAstrolab location, where the "sph"
986 // components of r correspond to theta=zenith angle and phi=pi-azimuth.
987 // "icr" ==> ICRS coordinates with longitude (l) and latitude (b),
988 // where the "sph" components of r correspond to theta=(pi/2)-b and phi=l.
989 // "loc" ==> Local coordinates at the AliAstrolab location, where the "sph"
990 // components of r correspond to the usual theta and phi angles.
992 // In case the coordinates are the equatorial right ascension and declination (a,d),
993 // they can represent so-called "mean" and "true" values.
994 // The distinction between these two representations is the following :
996 // mean values : (a,d) are only corrected for precession and not for nutation
997 // true values : (a,d) are corrected for both precession and nutation
999 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
1000 // values for the input in case of equatorial (a,d) coordinates.
1002 // mode = "M" --> Input coordinates are the mean values
1003 // "T" --> Input coordinates are the true values
1005 // Note : In case ts=0 the current timestamp of the lab will be taken.
1008 Int_t j=GetSignalIndex(name);
1009 if (j>=0) sx=GetSignal(r,frame,mode,ts,j);
1012 ///////////////////////////////////////////////////////////////////////////
1013 AliSignal* AliAstrolab::GetSignal(Double_t& a,Double_t& d,TString mode,AliTimestamp* ts,Int_t jref)
1015 // Provide precession (and nutation) corrected right ascension (a) and
1016 // declination (d) of the stored signal object at the specified timestamp.
1017 // In addition also a pointer to the stored signal object is provided.
1018 // In case no stored signal was available or one of the input arguments was
1019 // invalid, the returned pointer will be 0.
1021 // The coordinates (a,d) can represent so-called "mean" and "true" values.
1022 // The distinction between these two representations is the following :
1024 // mean values : (a,d) are only corrected for precession and not for nutation
1025 // true values : (a,d) are corrected for both precession and nutation
1027 // The input parameter "mode" allows the user to select either
1028 // "mean" or "true" values for (a,d).
1030 // The correction methods used are the new IAU 2000 ones as described in the
1031 // U.S. Naval Observatory (USNO) circular 179 (2005), which is available on
1032 // http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
1034 // a : Right ascension in hhmmss.sss
1035 // d : Declination in dddmmss.sss
1036 // mode = "M" --> Output coordinates are the mean values
1037 // "T" --> Output coordinates are the true values
1038 // ts : Timestamp at which the corrected coordinate values are requested.
1040 // The input parameter "jref" allows the user to access so-called "reference" signals.
1041 // These reference signals may be used to check space-time event coincidences with the
1042 // stored regular signal (e.g. coincidence of the stored signal with transient phenomena).
1044 // jref = 0 --> Access to the measurement
1045 // j --> Access to the reference signal at the j-th position (j=1 is first)
1047 // Default value is jref=0.
1049 // Note : In case ts=0 the current timestamp of the lab will be taken.
1055 AliSignal* sx=GetSignal(r,"equ",mode,ts,jref);
1059 // Retrieve the requested (a,d) values
1061 r.GetVector(vec,"sph","deg");
1082 // Convert coordinates to appropriate format
1083 a=ConvertAngle(a,"deg","hms");
1084 d=ConvertAngle(d,"deg","dms");
1088 ///////////////////////////////////////////////////////////////////////////
1089 AliSignal* AliAstrolab::GetSignal(Double_t& a,Double_t& d,TString mode,AliTimestamp* ts,TString name)
1091 // Provide precession (and nutation) corrected right ascension (a) and
1092 // declination (d) of the stored signal object with the specified name
1093 // at the specific timestamp ts.
1094 // In addition also a pointer to the stored signal object is provided.
1095 // In case no stored signal was available or one of the input arguments was
1096 // invalid, the returned pointer will be 0.
1098 // The coordinates (a,d) can represent so-called "mean" and "true" values.
1099 // The distinction between these two representations is the following :
1101 // mean values : (a,d) are only corrected for precession and not for nutation
1102 // true values : (a,d) are corrected for both precession and nutation
1104 // The input parameter "mode" allows the user to select either
1105 // "mean" or "true" values for (a,d).
1107 // The correction methods used are the new IAU 2000 ones as described in the
1108 // U.S. Naval Observatory (USNO) circular 179 (2005), which is available on
1109 // http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
1111 // a : Right ascension in hhmmss.sss
1112 // d : Declination in dddmmss.sss
1113 // mode = "M" --> Output coordinates are the mean values
1114 // "T" --> Output coordinates are the true values
1115 // ts : Timestamp at which the corrected coordinate values are requested.
1116 // name : Name of the requested signal object
1118 // Note : In case ts=0 the current timestamp of the lab will be taken.
1121 Int_t j=GetSignalIndex(name);
1122 if (j>=0) sx=GetSignal(a,d,mode,ts,j);
1125 ///////////////////////////////////////////////////////////////////////////
1126 AliSignal* AliAstrolab::GetSignal(Double_t& a,Double_t& d,TString s,Double_t e,TString mode,Int_t jref)
1128 // Provide precession (and nutation) corrected right ascension (a) and
1129 // declination (d) of the stored signal object at the specified epoch e.
1130 // In addition also a pointer to the stored signal object is provided.
1131 // In case no stored signal was available or one of the input arguments was
1132 // invalid, the returned pointer will be 0.
1134 // The coordinates (a,d) can represent so-called "mean" and "true" values.
1135 // The distinction between these two representations is the following :
1137 // mean values : (a,d) are only corrected for precession and not for nutation
1138 // true values : (a,d) are corrected for both precession and nutation
1140 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
1141 // values for the input (a,d) coordinates.
1143 // a : Right ascension in hhmmss.sss
1144 // d : Declination in dddmmss.sss
1145 // s = "B" --> Besselian reference epoch.
1146 // "J" --> Julian reference epoch.
1147 // e : Reference epoch for the input coordinates (e.g. 1900, 1950, 2000,...)
1148 // mode = "M" --> Input coordinates are the mean values
1149 // "T" --> Input coordinates are the true values
1151 // The input parameter "jref" allows the user to access so-called "reference" signals.
1152 // These reference signals may be used to check space-time event coincidences with the
1153 // stored measurement (e.g. coincidence of the measurement with transient phenomena).
1155 // jref = 0 --> Access to the measurement
1156 // j --> Access to the reference signal at the j-th position (j=1 is first)
1158 // Default value is jref=0.
1163 if (s!="B" && s!="b" && s!="J" && s!="j") return 0;
1165 if (mode!="M" && mode!="m" && mode!="T" && mode!="t") return 0;
1167 // Convert coordinates to fractional degrees.
1168 a=ConvertAngle(a,"hms","deg");
1169 d=ConvertAngle(d,"dms","deg");
1175 AliSignal* sx=GetSignal(a,d,mode,&tx,jref);
1178 ///////////////////////////////////////////////////////////////////////////
1179 AliSignal* AliAstrolab::GetSignal(Double_t& a,Double_t& d,TString s,Double_t e,TString mode,TString name)
1181 // Provide precession (and nutation) corrected right ascension (a) and
1182 // declination (d) of the stored signal object with the specified name
1183 // at the specific epoch e.
1184 // In addition also a pointer to the stored signal object is provided.
1185 // In case no stored signal was available or one of the input arguments was
1186 // invalid, the returned pointer will be 0.
1188 // The coordinates (a,d) can represent so-called "mean" and "true" values.
1189 // The distinction between these two representations is the following :
1191 // mean values : (a,d) are only corrected for precession and not for nutation
1192 // true values : (a,d) are corrected for both precession and nutation
1194 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
1195 // values for the input (a,d) coordinates.
1197 // a : Right ascension in hhmmss.sss
1198 // d : Declination in dddmmss.sss
1199 // s = "B" --> Besselian reference epoch.
1200 // "J" --> Julian reference epoch.
1201 // e : Reference epoch for the input coordinates (e.g. 1900, 1950, 2000,...)
1202 // mode = "M" --> Input coordinates are the mean values
1203 // "T" --> Input coordinates are the true values
1204 // name : Name of the requested signal object
1207 Int_t j=GetSignalIndex(name);
1208 if (j>=0) sx=GetSignal(a,d,s,e,mode,j);
1211 ///////////////////////////////////////////////////////////////////////////
1212 AliSignal* AliAstrolab::GetSignal(Int_t jref)
1214 // Provide the pointer to a stored signal object.
1216 // The input parameter "jref" allows the user to access so-called "reference" signals.
1217 // These reference signals may be used to check space-time event coincidences with the
1218 // stored measurement (e.g. coincidence of the measurement with transient phenomena).
1220 // jref = 0 --> Access to the measurement
1221 // j --> Access to the reference signal at the j-th position (j=1 is first)
1223 // Default value is jref=0.
1232 if (jref>0 && jref<fRefs->GetSize()) sx=(AliSignal*)fRefs->At(jref-1);
1236 ///////////////////////////////////////////////////////////////////////////
1237 AliSignal* AliAstrolab::GetSignal(TString name)
1239 // Provide the pointer to the stored signal object with the specified name.
1242 Int_t j=GetSignalIndex(name);
1243 if (j>=0) sx=GetSignal(j);
1246 ///////////////////////////////////////////////////////////////////////////
1247 void AliAstrolab::RemoveRefSignal(Int_t j,Int_t compress)
1249 // Remove the reference signal which was stored at the j-th position (j=1 is first).
1250 // Note : j=0 means that all stored reference signals will be removed.
1251 // j<0 allows array compression (see below) without removing any signals.
1253 // The "compress" parameter allows compression of the ref. signal storage array.
1255 // compress = 1 --> Array will be compressed
1256 // 0 --> Array will not be compressed
1258 // Note : Compression of the storage array means that the indices of the
1259 // reference signals in the storage array will change.
1263 // Clearing of the complete storage
1271 // Removing a specific reference signal
1272 if (j>0 && j<fRefs->GetSize())
1274 TObject* obj=fRefs->RemoveAt(j-1);
1275 if (obj) delete obj;
1278 // Compression of the storage array
1279 if (compress) fRefs->Compress();
1281 ///////////////////////////////////////////////////////////////////////////
1282 void AliAstrolab::RemoveRefSignal(TString name,Int_t compress)
1284 // Remove the reference signal with the specified name.
1286 // The "compress" parameter allows compression of the ref. signal storage array.
1288 // compress = 1 --> Array will be compressed
1289 // 0 --> Array will not be compressed
1291 // Note : Compression of the storage array means that the indices of the
1292 // reference signals in the storage array will change.
1294 Int_t j=GetSignalIndex(name);
1295 if (j>0) RemoveRefSignal(j,compress);
1297 ///////////////////////////////////////////////////////////////////////////
1298 Int_t AliAstrolab::GetSignalIndex(TString name)
1300 // Provide storage index of the signal with the specified name.
1301 // In case the name matches with the stored measurement,
1302 // the value 0 is returned.
1303 // In case no signal with the specified name was found, the value -1 is returned.
1309 if (name==fXsig->GetName()) return 0;
1312 if (!fRefs) return -1;
1314 for (Int_t i=0; i<fRefs->GetSize(); i++)
1316 AliSignal* sx=(AliSignal*)fRefs->At(i);
1319 if (name==sx->GetName())
1328 ///////////////////////////////////////////////////////////////////////////
1329 void AliAstrolab::PrintSignal(TString frame,TString mode,AliTimestamp* ts,Int_t ndig,Int_t jref)
1331 // Print data of a stored signal in user specified coordinates at the specific timestamp ts.
1332 // In case no stored signal was available or one of the input arguments was
1333 // invalid, no printout is produced.
1335 // The argument "ndig" specifies the number of digits for the fractional
1336 // part (e.g. ndig=6 for "dms" corresponds to micro-arcsecond precision).
1337 // No rounding will be performed, so an arcsecond count of 3.473 with ndig=1
1338 // will appear as 03.4 on the output.
1339 // Due to computer accuracy, precision on the pico-arcsecond level may get lost.
1341 // Note : The angle info is printed without additional spaces or "endline".
1342 // This allows the print to be included in various composite output formats.
1344 // The input parameter "frame" allows the user to specify the frame to which
1345 // the coordinates refer. Available options are :
1347 // frame = "equ" ==> Equatorial coordinates with right ascension (a) and declination (d).
1349 // "gal" ==> Galactic coordinates with longitude (l) and latitude (b).
1351 // "ecl" ==> Ecliptic coordinates with longitude (l) and latitude (b).
1353 // "hor" ==> Horizontal azimuth and altitude coordinates at the AliAstrolab location.
1355 // "icr" ==> ICRS coordinates with longitude (l) and latitude (b).
1357 // "loc" ==> Local spherical angles theta and phi at the AliAstrolab location.
1359 // In case the coordinates are the equatorial right ascension and declination (a,d),
1360 // they can represent so-called "mean" and "true" values.
1361 // The distinction between these two representations is the following :
1363 // mean values : (a,d) are only corrected for precession and not for nutation
1364 // true values : (a,d) are corrected for both precession and nutation
1366 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
1367 // values for the input in case of equatorial (a,d) coordinates.
1369 // mode = "M" --> Input coordinates are the mean values
1370 // "T" --> Input coordinates are the true values
1372 // The input parameter "mode" also determines which type of time and
1373 // local hour angle will appear in the printout.
1375 // mode = "M" --> Mean Sidereal Time (MST) and Local Mean Hour Angle (LMHA)
1376 // "T" --> Apparent Sidereal Time (AST) and Local Apparent Hour Angle (LAHA)
1378 // The input parameter "jref" allows printing of a so-called "reference" signal.
1379 // These reference signals may serve to check space-time event coincidences with the
1380 // stored measurement (e.g. coincidence of the measurement with transient phenomena).
1382 // jref = 0 --> Printing of the measurement
1383 // j --> Printing of the j-th reference signal
1385 // Default value is jref=0.
1387 // Note : In case ts=0 the current timestamp of the lab will be taken.
1390 AliSignal* sx=GetSignal(r,frame,mode,ts,jref);
1394 // Local Hour Angle of the signal
1395 Double_t lha=GetHourAngle("A",ts,jref);
1396 TString slha="LAHA";
1397 if (mode=="M" || mode=="m")
1399 lha=GetHourAngle("M",ts,jref);
1403 TString name=sx->GetName();
1404 if (name != "") cout << name.Data() << " ";
1409 d=90.-r.GetX(2,"sph","deg");
1410 a=r.GetX(3,"sph","rad");
1411 cout << "Equatorial (" << mode.Data() <<") a : "; PrintAngle(a,"rad","hms",ndig);
1412 cout << " d : "; PrintAngle(d,"deg","dms",ndig);
1413 cout << " " << slha.Data() << " : "; PrintAngle(lha,"deg","hms",ndig);
1420 b=90.-r.GetX(2,"sph","deg");
1421 l=r.GetX(3,"sph","deg");
1422 cout << "Galactic l : "; PrintAngle(l,"deg","deg",ndig);
1423 cout << " b : "; PrintAngle(b,"deg","deg",ndig);
1424 cout << " " << slha.Data() << " : "; PrintAngle(lha,"deg","hms",ndig);
1431 d=90.-r.GetX(2,"sph","deg");
1432 a=r.GetX(3,"sph","rad");
1433 cout << "ICRS l : "; PrintAngle(a,"rad","hms",ndig);
1434 cout << " b : "; PrintAngle(d,"deg","dms",ndig);
1435 cout << " " << slha.Data() << " : "; PrintAngle(lha,"deg","hms",ndig);
1442 d=90.-r.GetX(2,"sph","deg");
1443 a=r.GetX(3,"sph","deg");
1444 cout << "Ecliptic l : "; PrintAngle(a,"deg","deg",ndig);
1445 cout << " b : "; PrintAngle(d,"deg","deg",ndig);
1446 cout << " " << slha.Data() << " : "; PrintAngle(lha,"deg","hms",ndig);
1452 Double_t alt=90.-r.GetX(2,"sph","deg");
1453 Double_t azi=180.-r.GetX(3,"sph","deg");
1462 cout << "Horizontal azi : "; PrintAngle(azi,"deg","deg",ndig);
1463 cout << " alt : "; PrintAngle(alt,"deg","deg",ndig);
1464 cout << " " << slha.Data() << " : "; PrintAngle(lha,"deg","hms",ndig);
1470 Double_t theta=r.GetX(2,"sph","deg");
1471 Double_t phi=r.GetX(3,"sph","deg");
1472 cout << "Local-frame phi : "; PrintAngle(phi,"deg","deg",ndig);
1473 cout << " theta : "; PrintAngle(theta,"deg","deg",ndig);
1474 cout << " " << slha.Data() << " : "; PrintAngle(lha,"deg","hms",ndig);
1478 ///////////////////////////////////////////////////////////////////////////
1479 void AliAstrolab::PrintSignal(TString frame,TString mode,AliTimestamp* ts,Int_t ndig,TString name)
1481 // Print data of the stored signal with the specified name in user specified coordinates
1482 // at the specific timestamp ts.
1483 // In case no such stored signal was available or one of the input arguments was
1484 // invalid, no printout is produced.
1486 // The argument "ndig" specifies the number of digits for the fractional
1487 // part (e.g. ndig=6 for "dms" corresponds to micro-arcsecond precision).
1488 // No rounding will be performed, so an arcsecond count of 3.473 with ndig=1
1489 // will appear as 03.4 on the output.
1490 // Due to computer accuracy, precision on the pico-arcsecond level may get lost.
1492 // Note : The angle info is printed without additional spaces or "endline".
1493 // This allows the print to be included in various composite output formats.
1495 // The input parameter "frame" allows the user to specify the frame to which
1496 // the coordinates refer. Available options are :
1498 // frame = "equ" ==> Equatorial coordinates with right ascension (a) and declination (d).
1500 // "gal" ==> Galactic coordinates with longitude (l) and latitude (b).
1502 // "ecl" ==> Ecliptic coordinates with longitude (l) and latitude (b).
1504 // "hor" ==> Horizontal azimuth and altitude coordinates at the AliAstrolab location.
1506 // "icr" ==> ICRS coordinates with longitude (l) and latitude (b).
1508 // "loc" ==> Local spherical angles theta and phi at the AliAstrolab location.
1510 // In case the coordinates are the equatorial right ascension and declination (a,d),
1511 // they can represent so-called "mean" and "true" values.
1512 // The distinction between these two representations is the following :
1514 // mean values : (a,d) are only corrected for precession and not for nutation
1515 // true values : (a,d) are corrected for both precession and nutation
1517 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
1518 // values for the input in case of equatorial (a,d) coordinates.
1520 // mode = "M" --> Input coordinates are the mean values
1521 // "T" --> Input coordinates are the true values
1523 // The input parameter "mode" also determines which type of time and
1524 // local hour angle will appear in the printout.
1526 // mode = "M" --> Mean Sidereal Time (MST) and Local Mean Hour Angle (LMHA)
1527 // "T" --> Apparent Sidereal Time (AST) and Local Apparent Hour Angle (LAHA)
1529 // Note : In case ts=0 the current timestamp of the lab will be taken.
1531 Int_t j=GetSignalIndex(name);
1532 if (j>=0) PrintSignal(frame,mode,ts,ndig,j);
1534 ///////////////////////////////////////////////////////////////////////////
1535 void AliAstrolab::ListSignals(TString frame,TString mode,Int_t ndig)
1537 // List all stored signals in user specified coordinates at the timestamp
1538 // of the actual recording of the stored measurement under investigation.
1539 // In case no (timestamp of the) actual measurement is available,
1540 // the current timestamp of the lab will be taken.
1541 // In case no stored signal is available or one of the input arguments is
1542 // invalid, no printout is produced.
1544 // The argument "ndig" specifies the number of digits for the fractional
1545 // part (e.g. ndig=6 for "dms" corresponds to micro-arcsecond precision).
1546 // No rounding will be performed, so an arcsecond count of 3.473 with ndig=1
1547 // will appear as 03.4 on the output.
1548 // Due to computer accuracy, precision on the pico-arcsecond level may get lost.
1550 // The default value is ndig=1.
1552 // Note : The angle info is printed without additional spaces or "endline".
1553 // This allows the print to be included in various composite output formats.
1555 // The input parameter "frame" allows the user to specify the frame to which
1556 // the coordinates refer. Available options are :
1558 // frame = "equ" ==> Equatorial coordinates with right ascension (a) and declination (d).
1560 // "gal" ==> Galactic coordinates with longitude (l) and latitude (b).
1562 // "ecl" ==> Ecliptic coordinates with longitude (l) and latitude (b).
1564 // "hor" ==> Horizontal azimuth and altitude coordinates at the AliAstrolab location.
1566 // "icr" ==> ICRS coordinates with longitude (l) and latitude (b).
1568 // "loc" ==> Local spherical angles theta and phi at the AliAstrolab location.
1570 // In case the coordinates are the equatorial right ascension and declination (a,d),
1571 // they can represent so-called "mean" and "true" values.
1572 // The distinction between these two representations is the following :
1574 // mean values : (a,d) are only corrected for precession and not for nutation
1575 // true values : (a,d) are corrected for both precession and nutation
1577 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
1578 // values for the input in case of equatorial (a,d) coordinates.
1580 // mode = "M" --> Input coordinates are the mean values
1581 // "T" --> Input coordinates are the true values
1583 // The input parameter "mode" also determines which type of time and
1584 // local hour angle will appear in the listing.
1586 // mode = "M" --> Mean Sidereal Time (MST) and Local Mean Hour Angle (LMHA)
1587 // "T" --> Apparent Sidereal Time (AST) and Local Apparent Hour Angle (LAHA)
1594 if (mode=="T" || mode=="t") dform=-1;
1598 tx=fXsig->GetTimestamp();
1599 if (!tx) tx=(AliTimestamp*)this;
1600 cout << " *AliAstrolab::ListSignals* List of all stored signals." << endl;
1601 cout << " === The measurement under investigation ===" << endl;
1602 cout << " Timestamp of the actual observation" << endl;
1603 tx->Date(dform,fToffset);
1604 cout << " Location of the actual observation" << endl;
1605 cout << " "; PrintSignal(frame,mode,tx,ndig); cout << endl;
1611 for (Int_t i=1; i<=fRefs->GetSize(); i++)
1613 AliSignal* sx=GetSignal(i);
1618 cout << " *AliAstrolab::ListRefSignals* List of all stored signals." << endl;
1619 tx=(AliTimestamp*)this;
1620 cout << " Current timestamp of the laboratory" << endl;
1621 tx->Date(dform,fToffset);
1626 cout << " === All stored reference signals according to the above timestamp ===" << endl;
1629 cout << " Index : " << i << " "; PrintSignal(frame,mode,tx,ndig,i); cout << endl;
1632 ///////////////////////////////////////////////////////////////////////////
1633 void AliAstrolab::Precess(Ali3Vector& r,AliTimestamp* ts1,AliTimestamp* ts2)
1635 // Correct mean right ascension and declination given for Julian date "jd"
1636 // for the earth's precession corresponding to the specified timestamp.
1637 // The results are the so-called "mean" (i.e. precession corrected) values,
1638 // corresponding to the specified timestamp.
1639 // The method used is the new IAU 2000 one as described in the
1640 // U.S. Naval Observatory (USNO) circular 179 (2005), which is available on
1641 // http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
1642 // Since the standard reference epoch is J2000, this implies that all
1643 // input (a,d) coordinates will be first internally converted to the
1644 // corresponding J2000 values before the precession correction w.r.t. the
1645 // specified lab timestamp will be applied.
1647 // r : Input vector containing the right ascension and declination information
1648 // in the form of standard Ali3Vector coordinates.
1649 // In spherical coordinates the phi corresponds to the right ascension,
1650 // whereas the declination corresponds to (pi/2)-theta.
1651 // jd : Julian date corresponding to the input coordinate values.
1652 // ts : Timestamp corresponding to the requested corrected coordinate values.
1654 // Note : In case ts=0 the current timestamp of the lab will be taken.
1656 // Convert back to J2000 values
1659 r0=r.GetUnprimed(&fP);
1661 // Precess to the specified timestamp
1662 if (!ts2) ts2=(AliTimestamp*)this;
1664 r=r0.GetPrimed(&fP);
1666 ///////////////////////////////////////////////////////////////////////////
1667 void AliAstrolab::Nutate(Ali3Vector& r,AliTimestamp* ts)
1669 // Correct mean right ascension and declination for the earth's nutation
1670 // corresponding to the specified timestamp.
1671 // The results are the so-called "true" (i.e. nutation corrected) values,
1672 // corresponding to the specified timestamp.
1673 // The method used is the new IAU 2000 one as described in the
1674 // U.S. Naval Observatory (USNO) circular 179 (2005), which is available on
1675 // http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
1677 // r : Input vector containing the right ascension and declination information
1678 // in the form of standard Ali3Vector coordinates.
1679 // In spherical coordinates the phi corresponds to the right ascension,
1680 // whereas the declination corresponds to (pi/2)-theta.
1681 // ts : Timestamp for which the corrected coordinate values are requested.
1683 // Note : In case ts=0 the current timestamp of the lab will be taken.
1685 // Nutation correction for the specified timestamp
1686 if (!ts) ts=(AliTimestamp*)this;
1690 ///////////////////////////////////////////////////////////////////////////
1691 void AliAstrolab::SetBmatrix()
1693 // Set the frame bias matrix elements.
1694 // The formulas and numerical constants used are the ones from the
1695 // U.S. Naval Observatory (USNO) circular 179 (2005), which is available on
1696 // http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
1698 Double_t pi=acos(-1.);
1700 // Parameters in mas
1702 Double_t x=-16.6170;
1705 // Convert to radians
1706 a*=pi/(180.*3600.*1000.);
1707 x*=pi/(180.*3600.*1000.);
1708 e*=pi/(180.*3600.*1000.);
1711 mat[0]=1.-0.5*(a*a+x*x);
1715 mat[4]=1.-0.5*(a*a+e*e);
1719 mat[8]=1.-0.5*(e*e+x*x);
1724 ///////////////////////////////////////////////////////////////////////////
1725 void AliAstrolab::SetPmatrix(AliTimestamp* ts)
1727 // Set precession matrix elements for Julian date jd w.r.t. J2000.
1728 // The formulas and numerical constants used are the ones from the
1729 // U.S. Naval Observatory (USNO) circular 179 (2005), which is available on
1730 // http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
1731 // All numerical constants refer to the standard reference epoch J2000.
1733 Double_t mat[9]={0,0,0,0,0,0,0,0,0};
1740 Double_t pi=acos(-1.);
1742 Double_t t=(ts->GetJD()-2451545.0)/36525.; // Julian centuries since J2000.0
1744 // Parameters for the precession matrix in arcseconds
1745 Double_t eps0=84381.406; // Mean ecliptic obliquity at J2000.0
1746 Double_t psi=5038.481507*t-1.0790069*pow(t,2)-0.00114045*pow(t,3)+0.000132851*pow(t,4)
1747 -0.0000000951*pow(t,4);
1748 Double_t om=eps0-0.025754*t+0.0512623*pow(t,2)-0.00772503*pow(t,3)-0.000000467*pow(t,4)
1749 +0.0000003337*pow(t,5);
1750 Double_t chi=10.556403*t-2.3814292*pow(t,2)-0.00121197*pow(t,3)+0.000170663*pow(t,4)
1751 -0.0000000560*pow(t,5);
1753 // Convert to radians
1754 eps0*=pi/(180.*3600.);
1755 psi*=pi/(180.*3600.);
1756 om*=pi/(180.*3600.);
1757 chi*=pi/(180.*3600.);
1759 Double_t s1=sin(eps0);
1760 Double_t s2=sin(-psi);
1761 Double_t s3=sin(-om);
1762 Double_t s4=sin(chi);
1763 Double_t c1=cos(eps0);
1764 Double_t c2=cos(-psi);
1765 Double_t c3=cos(-om);
1766 Double_t c4=cos(chi);
1768 mat[0]=c4*c2-s2*s4*c3;
1769 mat[1]=c4*s2*c1+s4*c3*c2*c1-s1*s4*s3;
1770 mat[2]=c4*s2*s1+s4*c3*c2*s1+c1*s4*s3;
1771 mat[3]=-s4*c2-s2*c4*c3;
1772 mat[4]=-s4*s2*c1+c4*c3*c2*c1-s1*c4*s3;
1773 mat[5]=-s4*s2*s1+c4*c3*c2*s1+c1*c4*s3;
1775 mat[7]=-s3*c2*c1-s1*c3;
1776 mat[8]=-s3*c2*s1+c3*c1;
1780 ///////////////////////////////////////////////////////////////////////////
1781 void AliAstrolab::SetNmatrix(AliTimestamp* ts)
1783 // Set nutation matrix elements for the specified Julian date jd.
1784 // The formulas and numerical constants used are the ones from the
1785 // U.S. Naval Observatory (USNO) circular 179 (2005), which is available on
1786 // http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
1788 Double_t mat[9]={0,0,0,0,0,0,0,0,0};
1795 Double_t pi=acos(-1.);
1797 Double_t dpsi,deps,eps;
1798 ts->Almanac(&dpsi,&deps,&eps);
1800 // Convert to radians
1801 dpsi*=pi/(180.*3600.);
1802 deps*=pi/(180.*3600.);
1803 eps*=pi/(180.*3600.);
1805 Double_t s1=sin(eps);
1806 Double_t s2=sin(-dpsi);
1807 Double_t s3=sin(-(eps+deps));
1808 Double_t c1=cos(eps);
1809 Double_t c2=cos(-dpsi);
1810 Double_t c3=cos(-(eps+deps));
1816 mat[4]=c3*c2*c1-s1*s3;
1817 mat[5]=c3*c2*s1+c1*s3;
1819 mat[7]=-s3*c2*c1-s1*c3;
1820 mat[8]=-s3*c2*s1+c3*c1;
1824 ///////////////////////////////////////////////////////////////////////////
1825 void AliAstrolab::SetGmatrix(TString mode)
1827 // Set the mean equatorial to galactic coordinate conversion matrix.
1828 // The B1950 parameters were taken from section 22.3 of the book
1829 // "An Introduction to Modern Astrophysics" by Carrol and Ostlie (1996).
1830 // The J2000 parameters are obtained by precession of the B1950 values.
1832 // Via the input argument "mode" the required epoch can be selected
1833 // mode = "B" ==> B1950
1836 Ali3Vector x; // The Galactic x-axis in the equatorial frame
1837 Ali3Vector y; // The Galactic y-axis in the equatorial frame
1838 Ali3Vector z; // The Galactic z-axis in the equatorial frame
1841 Double_t vec[3]={1,0,0};
1843 fGal=1; // Set flag to indicate B1950 matrix values
1845 // B1950 equatorial coordinates of the North Galactic Pole (NGP)
1848 a=ConvertAngle(a,"hms","deg");
1849 d=ConvertAngle(d,"dms","deg");
1852 z.SetVector(vec,"sph","deg");
1854 // B1950 equatorial coordinates of the Galactic l=b=0 point
1857 a=ConvertAngle(a,"hms","deg");
1858 d=ConvertAngle(d,"dms","deg");
1861 x.SetVector(vec,"sph","deg");
1863 // Precess to the corresponding J2000 values if requested
1866 fGal=2; // Set flag to indicate J2000 matrix values
1868 t1.SetEpoch(1950,"B");
1870 t2.SetEpoch(2000,"J");
1875 // The Galactic y-axis is determined for the right handed frame
1878 fG.SetAngles(x.GetX(2,"sph","deg"),x.GetX(3,"sph","deg"),
1879 y.GetX(2,"sph","deg"),y.GetX(3,"sph","deg"),
1880 z.GetX(2,"sph","deg"),z.GetX(3,"sph","deg"));
1882 ///////////////////////////////////////////////////////////////////////////
1883 void AliAstrolab::SetEmatrix(AliTimestamp* ts)
1885 // Set the mean equatorial to ecliptic coordinate conversion matrix
1886 // for the specified timestamp.
1887 // A nice sketch and explanation of the two frames can be found
1888 // in chapter 3 of the book "Astronomy Methods" by Hale Bradt (2004).
1890 Double_t dpsi,deps,eps;
1891 ts->Almanac(&dpsi,&deps,&eps);
1893 // Convert to degrees
1896 // Positions of the ecliptic axes w.r.t. the equatorial ones
1897 // at the moment of the specified timestamp
1898 Double_t theta1=90; // Ecliptic x-axis
1900 Double_t theta2=90.-eps; //Ecliptic y-axis
1902 Double_t theta3=eps; // Ecliptic z-axis
1905 fE.SetAngles(theta1,phi1,theta2,phi2,theta3,phi3);
1907 ///////////////////////////////////////////////////////////////////////////
1908 void AliAstrolab::SetHmatrix(AliTimestamp* ts)
1910 // Set the mean equatorial to horizontal coordinate conversion matrix
1911 // for the specified timestamp.
1912 // A nice sketch and explanation of the two frames can be found
1913 // in chapter 3 of the book "Astronomy Methods" by Hale Bradt (2004).
1915 // Note : In order to simplify the calculations, we use here a
1916 // right-handed horizontal frame.
1918 Ali3Vector x; // The (South pointing) horizontal x-axis in the equatorial frame
1919 Ali3Vector y; // The (East pointing) horizontal y-axis in the equatorial frame
1920 Ali3Vector z; // The (Zenith pointing) horizontal z-axis in the equatorial frame
1923 GetLabPosition(l,b,"deg");
1926 Double_t vec[3]={1,0,0};
1928 // Equatorial coordinates of the horizontal z-axis
1929 // at the moment of the specified timestamp
1930 a=ts->GetLAST(fToffset);
1931 a*=15.; // Convert fractional hours to degrees
1934 z.SetVector(vec,"sph","deg");
1936 // Equatorial coordinates of the horizontal x-axis
1937 // at the moment of the specified timestamp
1940 x.SetVector(vec,"sph","deg");
1942 // The horizontal y-axis is determined for the right handed frame
1945 fH.SetAngles(x.GetX(2,"sph","deg"),x.GetX(3,"sph","deg"),
1946 y.GetX(2,"sph","deg"),y.GetX(3,"sph","deg"),
1947 z.GetX(2,"sph","deg"),z.GetX(3,"sph","deg"));
1949 ///////////////////////////////////////////////////////////////////////////
1950 void AliAstrolab::SetLocalFrame(Double_t t1,Double_t p1,Double_t t2,Double_t p2,Double_t t3,Double_t p3)
1952 // Specification of the orientations of the local-frame axes.
1953 // The input arguments represent the angles (in degrees) of the local-frame axes
1954 // w.r.t. a so called Master Reference Frame (MRF), with the same convention
1955 // as the input arguments of TRotMatrix::SetAngles.
1957 // The right handed Master Reference Frame is defined as follows :
1958 // Z-axis : Points to the local Zenith
1959 // X-axis : Makes an angle of 90 degrees with the Z-axis and points South
1960 // Y-axis : Makes an angle of 90 degrees with the Z-axis and points East
1962 // Once the user has specified the local reference frame, any observed event
1963 // can be related to astronomical space-time locations via the SetSignal
1964 // and GetSignal memberfunctions.
1966 // Set the matrix for the conversion of our reference frame coordinates
1967 // into the local-frame ones.
1969 fL.SetAngles(t1,p1,t2,p2,t3,p3);
1971 ///////////////////////////////////////////////////////////////////////////
1972 Double_t AliAstrolab::ConvertAngle(Double_t a,TString in,TString out) const
1974 // Conversion of various angular formats.
1976 // The input argument "a" denotes the angle to be converted.
1977 // The string arguments "in" and "out" specify the angular I/O formats.
1979 // in = "rad" : input angle provided in radians
1980 // "deg" : input angle provided in degrees
1981 // "dms" : input angle provided in dddmmss.sss
1982 // "hms" : input angle provided in hhmmss.sss
1983 // "hrs" : input angle provided in fractional hours
1985 // out = "rad" : output angle provided in radians
1986 // "deg" : output angle provided in degrees
1987 // "dms" : output angle provided in dddmmss.sss
1988 // "hms" : output angle provided in hhmmss.sss
1989 // "hrs" : output angle provided in fractional hours
1991 if (in==out) return a;
1993 // Convert input to its absolute value in (fractional) degrees.
1994 Double_t pi=acos(-1.);
1995 Double_t epsilon=1.e-12; // Accuracy in (arc)seconds
1996 Int_t word=0,ddd=0,hh=0,mm=0,ss=0;
2001 if (in=="rad") b*=180./pi;
2003 if (in=="hrs") b*=15.;
2012 s=b-Double_t(ddd*10000+mm*100+ss);
2013 b=Double_t(ddd)+Double_t(mm)/60.+(Double_t(ss)+s)/3600.;
2023 s=b-Double_t(hh*10000+mm*100+ss);
2024 b=15.*(Double_t(hh)+Double_t(mm)/60.+(Double_t(ss)+s)/3600.);
2032 if (out=="rad") b*=pi/180.;
2034 if (out=="hrs") b/=15.;
2065 b=Double_t(10000*ddd+100*mm+ss)+s;
2098 b=Double_t(10000*hh+100*mm+ss)+s;
2105 ///////////////////////////////////////////////////////////////////////////
2106 Double_t AliAstrolab::GetHourAngle(TString mode,AliTimestamp* ts,Int_t jref)
2108 // Provide the Local Hour Angle (in fractional degrees) of a stored signal
2109 // object at the specified timestamp.
2111 // The input parameter "mode" allows the user to select either the
2112 // "mean" or "apparent" value for the returned Hour Angle.
2114 // mode = "M" --> Output is the Mean Hour Angle
2115 // "A" --> Output is the Apparent Hour Angle
2116 // ts : Timestamp at which the hour angle is requested.
2118 // The input parameter "jref" allows the user to specify so-called "reference" signals.
2119 // These reference signals may be used to check space-time event coincidences with the
2120 // stored measurement (e.g. coincidence of the measurement with transient phenomena).
2122 // jref = 0 --> Use the stored measurement
2123 // j --> Use the reference signal at the j-th position (j=1 is first)
2125 // Default value is jref=0.
2127 // Note : In case ts=0 the current timestamp of the lab will be taken.
2129 if (!ts) ts=(AliTimestamp*)this;
2131 // Get corrected right ascension and declination for the specified timestamp.
2133 if (mode=="M" || mode=="m") GetSignal(a,d,"M",ts,jref);
2134 if (mode=="A" || mode=="a") GetSignal(a,d,"T",ts,jref);
2136 // Convert coordinates to fractional degrees.
2137 a=ConvertAngle(a,"hms","deg");
2138 d=ConvertAngle(d,"dms","deg");
2140 a/=15.; // Convert a to fractional hours
2142 if (mode=="M" || mode=="m") ha=ts->GetLMST(fToffset)-a;
2143 if (mode=="A" || mode=="a") ha=ts->GetLAST(fToffset)-a;
2144 ha*=15.; // Convert to (fractional) degrees
2148 ///////////////////////////////////////////////////////////////////////////
2149 void AliAstrolab::SetLT(Int_t y,Int_t m,Int_t d,Int_t hh,Int_t mm,Int_t ss,Int_t ns,Int_t ps)
2151 // Set the AliTimestamp parameters corresponding to the LT date and time
2152 // in the Gregorian calendar as specified by the input arguments.
2154 // The input arguments represent the following :
2155 // y : year in LT (e.g. 1952, 2003 etc...)
2156 // m : month in LT (1=jan 2=feb etc...)
2157 // d : day in LT (1-31)
2158 // hh : elapsed hours in LT (0-23)
2159 // mm : elapsed minutes in LT (0-59)
2160 // ss : elapsed seconds in LT (0-59)
2161 // ns : remaining fractional elapsed second of LT in nanosecond
2162 // ps : remaining fractional elapsed nanosecond of LT in picosecond
2164 // Note : ns=0 and ps=0 are the default values.
2167 SetLT(fToffset,y,m,d,hh,mm,ss,ns,ps);
2169 ///////////////////////////////////////////////////////////////////////////
2170 void AliAstrolab::SetLT(Int_t y,Int_t d,Int_t s,Int_t ns,Int_t ps)
2172 // Set the AliTimestamp parameters corresponding to the specified elapsed
2173 // timespan since the beginning of the new LT year.
2175 // The LT year and elapsed time span is entered via the following input arguments :
2177 // y : year in LT (e.g. 1952, 2003 etc...)
2178 // d : elapsed number of days
2179 // s : (remaining) elapsed number of seconds
2180 // ns : (remaining) elapsed number of nanoseconds
2181 // ps : (remaining) elapsed number of picoseconds
2183 // The specified d, s, ns and ps values will be used in an additive
2184 // way to determine the elapsed timespan.
2185 // So, specification of d=1, s=100, ns=0, ps=0 will result in the
2186 // same elapsed time span as d=0, s=24*3600+100, ns=0, ps=0.
2187 // However, by making use of the latter the user should take care
2188 // of possible integer overflow problems in the input arguments,
2189 // which obviously will provide incorrect results.
2191 // Note : ns=0 and ps=0 are the default values.
2193 SetLT(fToffset,y,d,s,ns,ps);
2195 ///////////////////////////////////////////////////////////////////////////
2196 Double_t AliAstrolab::GetDifference(Int_t j,TString au,Double_t& dt,TString tu,Int_t mode,Int_t* ia,Int_t* it)
2198 // Provide space and time difference between the j-th reference signal
2199 // (j=1 indicates first) and the stored measurement.
2201 // The return value of this memberfunction provides the positional angular
2202 // difference, whereas the output argument "dt" provides the time difference.
2204 // The units of the angular difference can be specified via the the "au"
2205 // input argument, where
2207 // au = "rad" --> Angular difference in (fractional) radians
2208 // "deg" --> Angular difference in (fractional) degrees
2210 // The units of the time difference can be specified via the "tu" and "mode"
2211 // input arguments. For details please refer to AliTimestamp::GetDifference().
2212 // Also here mode=1 is the default value.
2214 // For the time difference the reference signal is used as the standard.
2215 // This means that in case of a positive time difference, the stored
2216 // measurement occurred later than the reference signal.
2218 // In case j=0, the stored measurement will be compared with each
2219 // reference signal and the returned angular and time differences are
2220 // the minimal differences which were encountered.
2221 // In this case the user may obtain the indices of the two stored reference signals
2222 // which had the minimal angular and minimal time difference via the output
2223 // arguments "ia" and "it" as follows :
2225 // ia = Index of the stored reference signal with minimial angular difference
2226 // it = Index of the stored reference signal with minimial time difference
2228 // In case these indices are the same, there obviously was 1 single reference signal
2229 // which showed both the minimal angular and time difference.
2231 // The default values are mode=1, ia=0 and it=0;
2238 if (!fRefs) return dang;
2240 Ali3Vector rx; // Position of the measurement
2241 Ali3Vector r0; // Position of the reference signal
2243 AliSignal* sx=GetSignal(rx,"icr","M",0);
2244 if (!sx) return dang;
2246 AliTimestamp* tx=sx->GetTimestamp();
2247 if (!tx) return dang;
2249 // Space and time difference w.r.t. a specific reference signal
2252 AliSignal* s0=GetSignal(r0,"icr","M",0,j);
2253 if (!s0) return dang;
2254 AliTimestamp* t0=s0->GetTimestamp();
2256 if (!t0) return dang;
2258 dang=r0.GetOpeningAngle(rx,au);
2259 dt=t0->GetDifference(tx,tu,mode);
2263 // Minimal space and time difference encountered over all reference signals
2264 Double_t dangmin=dang;
2266 for (Int_t i=1; i<=fRefs->GetSize(); i++)
2268 AliSignal* s0=GetSignal(r0,"icr","M",0,i);
2270 AliTimestamp* t0=s0->GetTimestamp();
2272 dang=r0.GetOpeningAngle(rx,au);
2273 dt=t0->GetDifference(tx,tu,mode);
2274 if (fabs(dang)<dangmin)
2289 ///////////////////////////////////////////////////////////////////////////
2290 Double_t AliAstrolab::GetDifference(TString name,TString au,Double_t& dt,TString tu,Int_t mode)
2292 // Provide space and time difference between the reference signal
2293 // with the specified name and the stored measurement.
2295 // The return value of this memberfunction provides the positional angular
2296 // difference, whereas the output argument "dt" provides the time difference.
2298 // The units of the angular difference can be specified via the the "au"
2299 // input argument, where
2301 // au = "rad" --> Angular difference in (fractional) radians
2302 // "deg" --> Angular difference in (fractional) degrees
2304 // The units of the time difference can be specified via the "tu" and "mode"
2305 // input arguments. For details please refer to AliTimestamp::GetDifference().
2306 // Also here mode=1 is the default value.
2308 // For the time difference the reference signal is used as the standard.
2309 // This means that in case of a positive time difference, the stored
2310 // measurement occurred later than the reference signal.
2315 Int_t j=GetSignalIndex(name);
2316 if (j>0) dang=GetDifference(j,au,dt,tu,mode);
2319 ///////////////////////////////////////////////////////////////////////////
2320 TArrayI* AliAstrolab::MatchRefSignal(Double_t da,TString au,Double_t dt,TString tu,Int_t mode)
2322 // Provide the storage indices of the reference signals which match in space
2323 // and time with the stored measurement.
2324 // The indices are returned via a pointer to a TArrayI object.
2325 // In case no matches were found, the null pointer is returned.
2326 // A reference signal is regarded as matching with the stored measurement
2327 // if the positional angular difference doesn't exceed "da" and the absolute
2328 // value of the time difference doesn't exceed "dt".
2330 // The units of the angular difference "da" can be specified via the the "au"
2331 // input argument, where
2333 // au = "rad" --> Angular difference in (fractional) radians
2334 // "deg" --> Angular difference in (fractional) degrees
2336 // The units of the time difference "dt" can be specified via the "tu" and "mode"
2337 // input arguments. For details please refer to AliTimestamp::GetDifference().
2338 // Also here mode=1 is the default value.
2340 if (!fXsig || !fRefs) return 0;
2342 Int_t nrefs=fRefs->GetEntries();
2344 if (fIndices) delete fIndices;
2345 fIndices=new TArrayI(nrefs);
2347 Double_t dang,dtime;
2349 for (Int_t i=1; i<=fRefs->GetSize(); i++)
2351 dang=GetDifference(i,au,dtime,tu,mode);
2352 if (fabs(dang)<=da && fabs(dtime)<=dt)
2354 fIndices->AddAt(i,jfill);
2359 fIndices->Set(jfill);
2362 ///////////////////////////////////////////////////////////////////////////
2363 void AliAstrolab::DisplaySignal(TString frame,TString mode,AliTimestamp* ts,Int_t jref,TString proj,Int_t clr)
2365 // Display a stored signal in a user specified coordinate projection
2366 // at the specific timestamp ts.
2368 // In case no stored signal was available or one of the input arguments was
2369 // invalid, no display is produced.
2371 // Note : In case ts=0 the current timestamp of the lab will be taken.
2373 // The input parameter "frame" allows the user to specify the frame to which
2374 // the coordinates refer. Available options are :
2376 // frame = "equ" ==> Equatorial coordinates with right ascension (a) and declination (d).
2378 // "gal" ==> Galactic coordinates with longitude (l) and latitude (b).
2380 // "ecl" ==> Ecliptic coordinates with longitude (l) and latitude (b).
2382 // "hor" ==> Horizontal azimuth and altitude coordinates at the AliAstrolab location.
2384 // "icr" ==> ICRS coordinates with longitude (l) and latitude (b).
2386 // "loc" ==> Local spherical angles theta and phi at the AliAstrolab location.
2388 // In case the coordinates are the equatorial right ascension and declination (a,d),
2389 // they can represent so-called "mean" and "true" values.
2390 // The distinction between these two representations is the following :
2392 // mean values : (a,d) are only corrected for precession and not for nutation
2393 // true values : (a,d) are corrected for both precession and nutation
2395 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
2396 // values for the input in case of equatorial (a,d) coordinates.
2398 // mode = "M" --> Input coordinates are the mean values
2399 // "T" --> Input coordinates are the true values
2401 // The input parameter "jref" allows display of a so-called "reference" signal.
2402 // These reference signals may serve to check space-time event coincidences with the
2403 // stored measurement (e.g. coincidence of the measurement with transient phenomena).
2405 // jref = 0 --> Display of the measurement
2406 // j --> Display of the j-th reference signal
2408 // Default value is jref=0.
2410 // The input parameter "proj" allows the user to specify the desired projection.
2411 // The available projection modes are :
2413 // cyl : Cylindrical projection plotted with colored markers
2414 // cylh : Cylindrical projection plotted in a 2-D histogram
2415 // ham : Hammer equal area projection plotted with colored markers
2416 // hamh : Hammer equal area projection plotted in a 2-D histogram
2417 // ait : Aitoff projection plotted with colored markers
2418 // aith : Aitoff projection plotted in a 2-D histogram
2419 // mer : Mercator projection plotted with colored markers
2420 // merh : Mercator projection plotted in a 2-D histogram
2421 // ang : Straight cos(b) vs. l plot with colored markers
2422 // angh : Straight cos(b) vs. l plot in a 2-D histogram
2424 // Note : The ang(h) plot allows for easy identification of a homogeneous distribution.
2426 // The input argument "clr" allows to clear (1) the display before drawing or not (0).
2428 // The default values are : jref=0, proj="ham" and clr=0.
2430 // This routine is based on the work by Garmt de Vries-Uiterweerd.
2433 AliSignal* sx=GetSignal(r,frame,mode,ts,jref);
2437 // The generic input angles (in rad) for the projections
2441 Double_t pi=acos(-1.);
2443 if (frame=="equ" || frame=="gal" || frame=="icr" || frame=="ecl" || frame=="loc")
2445 theta=(pi/2.)-r.GetX(2,"sph","rad");
2446 phi=r.GetX(3,"sph","rad");
2451 theta=(pi/2.)-r.GetX(2,"sph","rad");
2452 phi=pi-r.GetX(3,"sph","rad");
2455 // Automatic choice of central meridian if not selected by the user
2468 // Obtain the projected (x,y) position
2471 Project(phi,theta,proj,x,y);
2474 if (proj=="hamh" || proj=="aith" || proj=="merh" || proj=="cylh" || proj=="angh") hist=1;
2476 // Update the display for this signal position
2478 // Create a new canvas if needed
2479 if (!fCanvas) fCanvas=new TCanvas("AliAstrolab","Skymap");
2481 // Construct the title string for this map
2482 TString title="Projection : ";
2486 title+=" Center : ";
2488 if (frame=="equ" || frame=="icr")
2490 ang=int(ConvertAngle(fMeridian,"rad","hms"));
2504 ang=int(ConvertAngle(fMeridian,"rad","dms"));
2517 if (!hist) // 2-D Marker display (i.e. not a histogram)
2519 // Remove existing markers, grid and outline from display if needed
2520 if (clr==1 || proj!=fProj)
2533 fMarkers=new TObjArray();
2534 fMarkers->SetOwner();
2535 // Set canvas range, header and axes
2538 fCanvas->Range(-2.2,-6.,2.2,6.);
2539 TText* header=new TText;
2540 fMarkers->Add(header);
2541 header->DrawText(-1.4,5.5,title.Data());
2542 TLine* line=new TLine(-2,0,2,0);
2543 fMarkers->Add(line);
2545 line=new TLine(0,5.5,0,-5.5);
2546 fMarkers->Add(line);
2551 fCanvas->Range(-2.2,-1.2,2.2,1.2);
2552 TText* header=new TText;
2553 fMarkers->Add(header);
2554 header->DrawText(-1.4,1.1,title.Data());
2555 TLine* line=new TLine(-2,0,2,0);
2556 fMarkers->Add(line);
2558 line=new TLine(0,1,0,-1);
2559 fMarkers->Add(line);
2563 if (proj=="ham" || proj=="ait")
2565 // Draw ellips outline with axes
2566 TEllipse* outline=new TEllipse(0,0,2,1);
2567 fMarkers->Add(outline);
2572 if (!jref) // The measurement
2574 marker=new TMarker(x,y,29);
2575 marker->SetMarkerColor(kRed);
2577 else // The reference signal(s)
2579 marker=new TMarker(x,y,20);
2580 marker->SetMarkerColor(kBlue);
2582 marker->SetMarkerSize(1);
2583 fMarkers->Add(marker);
2586 else if (hist==1) // 2-D display via histogram
2588 // Reset the histogram if needed
2589 if (clr==1 || proj!=fProj)
2592 if (!fHist) fHist=new TH2F();
2594 fHist->SetMarkerStyle(20);
2595 fHist->SetMarkerSize(1);
2596 fHist->SetMarkerColor(kBlue);
2597 fHist->SetNameTitle("SkyMap",title.Data());
2600 fHist->SetBins(100,-2.2,2.2,100,-5.5,5.5);
2604 fHist->SetBins(100,-2.2,2.2,100,-1.1,1.1);
2612 ///////////////////////////////////////////////////////////////////////////
2613 void AliAstrolab::DisplaySignal(TString frame,TString mode,AliTimestamp* ts,TString name,TString proj,Int_t clr)
2615 // Display the stored signal with the specified name in a user specified
2616 // coordinate projection at the specific timestamp ts.
2618 // Note : In case ts=0 the current timestamp of the lab will be taken.
2620 // In case no such stored signal was available or one of the input arguments was
2621 // invalid, no display is produced.
2623 // The input parameter "frame" allows the user to specify the frame to which
2624 // the coordinates refer. Available options are :
2626 // frame = "equ" ==> Equatorial coordinates with right ascension (a) and declination (d).
2628 // "gal" ==> Galactic coordinates with longitude (l) and latitude (b).
2630 // "ecl" ==> Ecliptic coordinates with longitude (l) and latitude (b).
2632 // "hor" ==> Horizontal azimuth and altitude coordinates at the AliAstrolab location.
2634 // "icr" ==> ICRS coordinates with longitude (l) and latitude (b).
2636 // "loc" ==> Local spherical angles theta and phi at the AliAstrolab location.
2638 // In case the coordinates are the equatorial right ascension and declination (a,d),
2639 // they can represent so-called "mean" and "true" values.
2640 // The distinction between these two representations is the following :
2642 // mean values : (a,d) are only corrected for precession and not for nutation
2643 // true values : (a,d) are corrected for both precession and nutation
2645 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
2646 // values for the input in case of equatorial (a,d) coordinates.
2648 // mode = "M" --> Input coordinates are the mean values
2649 // "T" --> Input coordinates are the true values
2651 // The input parameter "proj" allows the user to specify the desired projection.
2652 // The available projection modes are :
2654 // cyl : Cylindrical projection plotted with colored markers
2655 // cylh : Cylindrical projection plotted in a 2-D histogram
2656 // ham : Hammer equal area projection plotted with colored markers
2657 // hamh : Hammer equal area projection plotted in a 2-D histogram
2658 // ait : Aitoff projection plotted with colored markers
2659 // aith : Aitoff projection plotted in a 2-D histogram
2660 // mer : Mercator projection plotted with colored markers
2661 // merh : Mercator projection plotted in a 2-D histogram
2662 // ang : Straight cos(b) vs. l plot with colored markers
2663 // angh : Straight cos(b) vs. l plot in a 2-D histogram
2665 // Note : The ang(h) plot allows for easy identification of a homogeneous distribution.
2667 // The input argument "clr" allows to clear (1) the display before drawing or not (0).
2669 // The default values are : proj="ham" and clr=0.
2671 // This routine is based on the work by Garmt de Vries-Uiterweerd.
2673 Int_t j=GetSignalIndex(name);
2674 if (j>=0) DisplaySignal(frame,mode,ts,j,proj,clr);
2676 ///////////////////////////////////////////////////////////////////////////
2677 void AliAstrolab::DisplaySignals(TString frame,TString mode,AliTimestamp* ts,TString proj,Int_t clr)
2679 // Display of all stored signals in user specified coordinate projection
2680 // at the specific timestamp ts.
2681 // In case ts=0, the timestamp of the actual recording of the stored measurement
2682 // under investigation will be used.
2683 // In case no (timestamp of the) actual measurement is available,
2684 // the current timestamp of the lab will be taken.
2685 // In case no stored signal is available or one of the input arguments is
2686 // invalid, no display is produced.
2688 // The input parameter "frame" allows the user to specify the frame to which
2689 // the coordinates refer. Available options are :
2691 // frame = "equ" ==> Equatorial coordinates with right ascension (a) and declination (d).
2693 // "gal" ==> Galactic coordinates with longitude (l) and latitude (b).
2695 // "ecl" ==> Ecliptic coordinates with longitude (l) and latitude (b).
2697 // "hor" ==> Horizontal azimuth and altitude coordinates at the AliAstrolab location.
2699 // "icr" ==> ICRS coordinates with longitude (l) and latitude (b).
2701 // "loc" ==> Local spherical angles theta and phi at the AliAstrolab location.
2703 // In case the coordinates are the equatorial right ascension and declination (a,d),
2704 // they can represent so-called "mean" and "true" values.
2705 // The distinction between these two representations is the following :
2707 // mean values : (a,d) are only corrected for precession and not for nutation
2708 // true values : (a,d) are corrected for both precession and nutation
2710 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
2711 // values for the input in case of equatorial (a,d) coordinates.
2713 // mode = "M" --> Input coordinates are the mean values
2714 // "T" --> Input coordinates are the true values
2716 // The input parameter "proj" allows the user to specify the desired projection.
2717 // The available projection modes are :
2719 // cyl : Cylindrical projection plotted with colored markers
2720 // cylh : Cylindrical projection plotted in a 2-D histogram
2721 // ham : Hammer equal area projection plotted with colored markers
2722 // hamh : Hammer equal area projection plotted in a 2-D histogram
2723 // ait : Aitoff projection plotted with colored markers
2724 // aith : Aitoff projection plotted in a 2-D histogram
2725 // mer : Mercator projection plotted with colored markers
2726 // merh : Mercator projection plotted in a 2-D histogram
2727 // ang : Straight cos(b) vs. l plot with colored markers
2728 // angh : Straight cos(b) vs. l plot in a 2-D histogram
2730 // Note : The ang(h) plot allows for easy identification of a homogeneous distribution.
2732 // The input argument "clr" allows to clear (1) the display before drawing or not (0).
2734 // The default values are : proj="ham" and clr=0.
2736 // This routine is based on the work by Garmt de Vries-Uiterweerd.
2738 AliTimestamp* tx=ts;
2739 if (fXsig && !tx) tx=fXsig->GetTimestamp();
2740 if (!tx) tx=(AliTimestamp*)this;
2742 // Display all stored reference signals
2745 for (Int_t i=1; i<=fRefs->GetSize(); i++)
2747 DisplaySignal(frame,mode,tx,i,proj,clr);
2748 clr=0; // No display clear for subsequent signals
2752 // Display the measurement
2753 DisplaySignal(frame,mode,tx,0,proj,clr);
2755 ///////////////////////////////////////////////////////////////////////////
2756 void AliAstrolab::SetCentralMeridian(Double_t phi,TString u)
2758 // Set the central meridian for the sky display.
2759 // Setting a value smaller than -pi will induce automatic meridian setting
2761 // By default the central meridian is set at -999 in the constructor.
2763 // The string argument "u" allows to choose between different angular units
2764 // u = "rad" : angle provided in radians
2765 // "deg" : angle provided in degrees
2766 // "dms" : angle provided in dddmmss.sss
2767 // "hms" : angle provided in hhmmss.sss
2769 // The default is u="deg".
2771 // This routine is based on the work by Garmt de Vries-Uiterweerd.
2773 fMeridian=ConvertAngle(phi,u,"rad");
2775 ///////////////////////////////////////////////////////////////////////////
2776 void AliAstrolab::Project(Double_t l,Double_t b,TString proj,Double_t& x,Double_t& y)
2778 // Generic interface for projection of a (long,lat) pair onto a (x,y) pair.
2780 // Meaning of the arguments :
2781 // --------------------------
2782 // l : Longitude (e.g. right ascension)
2783 // b : Latitude (e.g. declination)
2784 // proj : Projection mode (e.g. "ham")
2785 // x : The resulting x coordinate for the selected projection
2786 // y : The resulting x coordinate for the selected projection
2788 // The available projection modes are :
2790 // cyl : Cylindrical projection plotted with colored markers
2791 // cylh : Cylindrical projection plotted in a 2-D histogram
2792 // ham : Hammer equal area projection plotted with colored markers
2793 // hamh : Hammer equal area projection plotted in a 2-D histogram
2794 // ait : Aitoff projection plotted with colored markers
2795 // aith : Aitoff projection plotted in a 2-D histogram
2796 // mer : Mercator projection plotted with colored markers
2797 // merh : Mercator projection plotted in a 2-D histogram
2798 // ang : Straight cos(b) vs. l plot with colored markers
2799 // angh : Straight cos(b) vs. l plot in a 2-D histogram
2801 // Note : The ang(h) plot allows for easy identification of a homogeneous distribution.
2803 // This routine is based on the work by Garmt de Vries-Uiterweerd.
2805 Double_t pi=acos(-1.);
2807 // Subtract central meridian from longitude
2810 // Take l between -180 and 180 degrees
2823 // Convert (l,b) to (x,y) with -2 < x <= 2
2824 if (proj=="cyl" || proj=="cylh") ProjectCylindrical(l,b,x,y);
2825 if (proj=="ham" || proj=="hamh") ProjectHammer(l,b,x,y);
2826 if (proj=="ait" || proj=="aith") ProjectAitoff(l,b,x,y);
2827 if (proj=="mer" || proj=="merh") ProjectMercator(l,b,x,y);
2828 if (proj=="ang" || proj=="angh")
2834 ///////////////////////////////////////////////////////////////////////////
2835 void AliAstrolab::ProjectCylindrical(Double_t l,Double_t b,Double_t& x, Double_t& y)
2837 // Cylindrical projection of a (long,lat) coordinate pair
2839 // This routine is based on the work by Garmt de Vries-Uiterweerd.
2841 Double_t pi=acos(-1.);
2845 ///////////////////////////////////////////////////////////////////////////
2846 void AliAstrolab::ProjectHammer(Double_t l,Double_t b,Double_t& x,Double_t& y)
2848 // Hammer-Aitoff projection of a (long,lat) coordinate pair.
2849 // This is an equal-area projection.
2851 // This routine is based on the work by Garmt de Vries-Uiterweerd.
2853 Double_t k=1./sqrt(1.+cos(b)*cos(l/2.));
2854 x=2.*k*cos(b)*sin(l/2.);
2857 ///////////////////////////////////////////////////////////////////////////
2858 void AliAstrolab::ProjectAitoff(Double_t l,Double_t b,Double_t& x,Double_t& y)
2860 // Aitoff projection of a (long,lat) coordinate pair.
2862 // This routine is based on the work by Garmt de Vries-Uiterweerd.
2864 Double_t pi=acos(-1.);
2867 Double_t k=acos(cos(b)*cos(l/2.));
2870 x=4.*k*cos(b)*sin(l/2.)/(pi*sin(k));
2871 y=2.*k*sin(b)/(pi*sin(k));
2874 ///////////////////////////////////////////////////////////////////////////
2875 void AliAstrolab::ProjectMercator(Double_t l,Double_t b,Double_t& x,Double_t& y)
2877 // Mercator projection of a (long,lat) coordinate pair
2879 // This routine is based on the work by Garmt de Vries-Uiterweerd.
2881 Double_t pi=acos(-1.);
2884 if (b >= 89.*pi/180.) b=89.*pi/180.;
2885 if (b <= -89.*pi/180.) y=-89.*pi/180.;
2886 if (fabs(y) < 1) y=log((1.+sin(b))/(1.-sin(b)))/2.;
2888 ///////////////////////////////////////////////////////////////////////////
2889 TObject* AliAstrolab::Clone(const char* name) const
2891 // Make a deep copy of the current object and provide the pointer to the copy.
2892 // This memberfunction enables automatic creation of new objects of the
2893 // correct type depending on the object type, a feature which may be very useful
2894 // for containers when adding objects in case the container owns the objects.
2896 AliAstrolab* lab=new AliAstrolab(*this);
2899 if (strlen(name)) lab->SetName(name);
2903 ///////////////////////////////////////////////////////////////////////////