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 ///////////////////////////////////////////////////////////////////////////
130 #include "AliAstrolab.h"
131 #include "Riostream.h"
133 ClassImp(AliAstrolab) // Class implementation to enable ROOT I/O
135 AliAstrolab::AliAstrolab(const char* name,const char* title) : TTask(name,title),AliTimestamp()
137 // Default constructor
151 ///////////////////////////////////////////////////////////////////////////
152 AliAstrolab::~AliAstrolab()
154 // Destructor to delete all allocated memory.
187 ///////////////////////////////////////////////////////////////////////////
188 AliAstrolab::AliAstrolab(const AliAstrolab& t) : TTask(t),AliTimestamp(t)
195 if (t.fXsig) fXsig=new AliSignal(*(t.fXsig));
199 Int_t size=t.fRefs->GetSize();
200 fRefs=new TObjArray(size);
201 for (Int_t i=0; i<size; i++)
203 AliSignal* sx=(AliSignal*)t.fRefs->At(i);
204 if (sx) fRefs->AddAt(sx->Clone(),i);
216 ///////////////////////////////////////////////////////////////////////////
217 void AliAstrolab::Data(Int_t mode,TString u)
219 // Provide lab information.
221 // "mode" indicates the mode of the timestamp info (see AliTimestamp::Date).
223 // The string argument "u" allows to choose between different angular units
224 // in case e.g. a spherical frame is selected.
225 // u = "rad" : angles provided in radians
226 // "deg" : angles provided in degrees
227 // "dms" : angles provided in ddd:mm:ss.sss
228 // "hms" : angles provided in hh:mm:ss.sss
230 // The defaults are mode=1 and u="deg".
232 const char* name=GetName();
233 const char* title=GetTitle();
234 cout << " *" << ClassName() << "::Data*";
235 if (strlen(name)) cout << " Name : " << GetName();
236 if (strlen(title)) cout << " Title : " << GetTitle();
240 GetLabPosition(l,b,"deg");
241 cout << " Lab position longitude : "; PrintAngle(l,"deg",u,2);
242 cout << " latitude : "; PrintAngle(b,"deg",u,2);
244 cout << " Lab time offset w.r.t. UT : "; PrintTime(fToffset,12); cout << endl;
246 // UT and Local time info
249 ///////////////////////////////////////////////////////////////////////////
250 void AliAstrolab::SetLabPosition(Ali3Vector& p)
252 // Set the lab position in the terrestrial coordinates.
253 // The right handed reference frame is defined such that the North Pole
254 // corresponds to a polar angle theta=0 and the Greenwich meridian corresponds
255 // to an azimuth angle phi=0, with phi increasing eastwards.
257 fLabPos.SetPosition(p);
259 // Determine local time offset in fractional hours w.r.t. UT.
261 p.GetVector(vec,"sph","deg");
265 ///////////////////////////////////////////////////////////////////////////
266 void AliAstrolab::SetLabPosition(Double_t l,Double_t b,TString u)
268 // Set the lab position in the terrestrial longitude (l) and latitude (b).
269 // Positions north of the equator have b>0, whereas b<0 indicates
270 // locations south of the equator.
271 // Positions east of Greenwich have l>0, whereas l<0 indicates
272 // locations west of Greenwich.
274 // The string argument "u" allows to choose between different angular units
275 // u = "rad" : angles provided in radians
276 // "deg" : angles provided in degrees
277 // "dms" : angles provided in dddmmss.sss
278 // "hms" : angles provided in hhmmss.sss
280 // The default is u="deg".
282 Double_t r=1,theta=0,phi=0;
284 l=ConvertAngle(l,u,"deg");
285 b=ConvertAngle(b,u,"deg");
292 Double_t p[3]={r,theta,phi};
293 fLabPos.SetPosition(p,"sph","deg");
295 // Local time offset in fractional hours w.r.t. UT.
298 ///////////////////////////////////////////////////////////////////////////
299 AliPosition AliAstrolab::GetLabPosition() const
301 // Provide the lab position in the terrestrial coordinates.
302 // The right handed reference frame is defined such that the North Pole
303 // corresponds to a polar angle theta=0 and the Greenwich meridian corresponds
304 // to an azimuth angle phi=0, with phi increasing eastwards.
308 ///////////////////////////////////////////////////////////////////////////
309 void AliAstrolab::GetLabPosition(Double_t& l,Double_t& b,TString u) const
311 // Provide the lab position in the terrestrial longitude (l) and latitude (b).
312 // Positions north of the equator have b>0, whereas b<0 indicates
313 // locations south of the equator.
314 // Positions east of Greenwich have l>0, whereas l<0 indicates
315 // locations west of Greenwich.
317 // The string argument "u" allows to choose between different angular units
318 // u = "rad" : angles provided in radians
319 // "deg" : angles provided in degrees
321 // The default is u="deg".
323 Double_t pi=acos(-1.);
326 if (u=="rad") offset=pi/2.;
329 fLabPos.GetPosition(p,"sph",u);
333 ///////////////////////////////////////////////////////////////////////////
334 Double_t AliAstrolab::GetLT()
336 // Provide the Lab's local time in fractional hours.
337 // A mean solar day lasts 24h (i.e. 86400s).
339 // In case a hh:mm:ss format is needed, please use the Convert() facility.
341 Double_t h=GetLT(fToffset);
344 ///////////////////////////////////////////////////////////////////////////
345 Double_t AliAstrolab::GetLMST()
347 // Provide the Lab's Local Mean Sidereal Time (LMST) in fractional hours.
348 // A sidereal day corresponds to 23h 56m 04.09s (i.e. 86164.09s) mean solar time.
349 // The definition of GMST is such that a sidereal clock corresponds with
350 // 24 sidereal hours per revolution of the Earth.
351 // As such, local time offsets w.r.t. UT and GMST can be treated similarly.
353 // In case a hh:mm:ss format is needed, please use the Convert() facility.
355 Double_t h=GetLMST(fToffset);
358 ///////////////////////////////////////////////////////////////////////////
359 Double_t AliAstrolab::GetLAST()
361 // Provide the Lab's Local Apparent Sidereal Time (LAST) in fractional hours.
362 // A sidereal day corresponds to 23h 56m 04.09s (i.e. 86164.09s) mean solar time.
363 // The definition of GMST and GAST is such that a sidereal clock corresponds with
364 // 24 sidereal hours per revolution of the Earth.
365 // As such, local time offsets w.r.t. UT, GMST and GAST can be treated similarly.
367 // In case a hh:mm:ss format is needed, please use the Convert() facility.
369 Double_t h=GetLAST(fToffset);
372 ///////////////////////////////////////////////////////////////////////////
373 void AliAstrolab::PrintAngle(Double_t a,TString in,TString out,Int_t ndig) const
375 // Printing of angles in various formats.
377 // The input argument "a" denotes the angle to be printed.
378 // The string arguments "in" and "out" specify the angular I/O formats.
380 // in = "rad" : input angle provided in radians
381 // "deg" : input angle provided in degrees
382 // "dms" : input angle provided in dddmmss.sss
383 // "hms" : input angle provided in hhmmss.sss
385 // out = "rad" : output angle provided in radians
386 // "deg" : output angle provided in degrees
387 // "dms" : output angle provided in dddmmss.sss
388 // "hms" : output angle provided in hhmmss.sss
390 // The argument "ndig" specifies the number of digits for the fractional
391 // part (e.g. ndig=6 for "dms" corresponds to micro-arcsecond precision).
392 // No rounding will be performed, so an arcsecond count of 3.473 with ndig=1
393 // will appear as 03.4 on the output.
394 // Due to computer accuracy, precision on the pico-arcsecond level may get lost.
396 // The default is ndig=1.
398 // Note : The angle info is printed without additional spaces or "endline".
399 // This allows the print to be included in various composite output formats.
401 Double_t b=ConvertAngle(a,in,out);
403 if (out=="deg" || out=="rad")
405 cout.setf(ios::fixed,ios::floatfield);
406 cout << setprecision(ndig) << b << " " << out.Data();
407 cout.unsetf(ios::fixed);
411 Double_t epsilon=1.e-12; // Accuracy in (arc)seconds
412 Int_t word=0,ddd=0,hh=0,mm=0,ss=0;
424 s=fabs(b)-Double_t(ddd*10000+mm*100+ss);
446 if (b<0) cout << "-";
447 cout << ddd << "d " << mm << "' " << ss << "."
448 << setfill('0') << setw(ndig) << sfrac << "\"";
460 s=fabs(b)-Double_t(hh*10000+mm*100+ss);
482 if (b<0) cout << "-";
483 cout << hh << "h " << mm << "m " << ss << "."
484 << setfill('0') << setw(ndig) << sfrac << "s";
488 ///////////////////////////////////////////////////////////////////////////
489 void AliAstrolab::SetSignal(Ali3Vector* r,TString frame,TString mode,AliTimestamp* ts,Int_t jref,TString name)
491 // Store a signal as specified by the position r and the timestamp ts.
492 // The position is stored in International Celestial Reference System (ICRS) coordinates.
493 // The ICRS is a fixed, time independent frame and as such provides a unique reference
494 // frame without the need of specifying any epoch etc...
495 // The ICRS coordinate definitions match within 20 mas with the mean ones of the J2000.0
496 // equatorial system. Nevertheless, to obtain the highest accuracy, the slight
497 // coordinate correction between J2000 and ICRS is performed here via the
498 // so-called frame bias matrix.
499 // For further details see the U.S. Naval Observatory (USNO) circular 179 (2005),
500 // which is available on http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
502 // The input parameter "frame" allows the user to specify the frame to which
503 // the components of r refer. Available options are :
505 // frame = "equ" ==> Equatorial coordinates with right ascension (a) and declination (d),
506 // where the "sph" components of r correspond to theta=(pi/2)-d and phi=a.
507 // "gal" ==> Galactic coordinates with longitude (l) and latitude (b).
508 // where the "sph" components of r correspond to theta=(pi/2)-b and phi=l.
509 // "ecl" ==> Ecliptic coordinates with longitude (l) and latitude (b),
510 // where the "sph" components of r correspond to theta=(pi/2)-b and phi=l.
511 // "hor" ==> Horizontal coordinates at the AliAstrolab location, where the "sph"
512 // components of r correspond to theta=zenith angle and phi=pi-azimuth.
513 // "icr" ==> ICRS coordinates with longitude (l) and latitude (b),
514 // where the "sph" components of r correspond to theta=(pi/2)-b and phi=l.
515 // "loc" ==> Local coordinates at the AliAstrolab location, where the "sph"
516 // components of r correspond to the usual theta and phi angles.
518 // In case the coordinates are the equatorial right ascension and declination (a,d),
519 // they can represent so-called "mean" and "true" values.
520 // The distinction between these two representations is the following :
522 // mean values : (a,d) are only corrected for precession and not for nutation
523 // true values : (a,d) are corrected for both precession and nutation
525 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
526 // values for the input in case of equatorial (a,d) coordinates.
528 // mode = "M" --> Input coordinates are the mean values
529 // "T" --> Input coordinates are the true values
531 // The input parameter "jref" allows the user to store so-called "reference" signals.
532 // These reference signals may be used to check space-time event coincidences with the
533 // stored measurement (e.g. coincidence of the measurement with transient phenomena).
535 // jref = 0 --> Storage of the measurement
536 // j --> Storage of a reference signal at the j-th position (j=1 is first)
537 // < 0 --> Add a reference signal at the next available position
539 // Via the input argument "name" the user can give the stored signal also a name.
541 // The default values are jref=0 and name="".
543 // Note : In case ts=0 the current timestamp of the lab will be taken.
555 if (frame!="equ" && frame!="gal" && frame!="ecl" && frame!="hor" && frame!="icr" && frame!="loc")
565 if (frame=="equ" && mode!="M" && mode!="m" && mode!="T" && mode!="t")
575 if (!ts) ts=(AliTimestamp*)this;
577 Double_t vec[3]={1,0,0};
578 vec[1]=r->GetX(2,"sph","rad");
579 vec[2]=r->GetX(3,"sph","rad");
581 q.SetVector(vec,"sph","rad");
585 if (!jref) // Storage of the measurement
589 fXsig=new AliSignal();
595 if (name != "") fXsig->SetName(name);
596 fXsig->SetTitle("Event in ICRS coordinates");
597 fXsig->SetTimestamp(*ts);
599 else // Storage of a reference signal
603 fRefs=new TObjArray();
606 // Expand array size if needed
607 if (jref>0 && jref>=fRefs->GetSize()) fRefs->Expand(jref+1);
608 sxref=new AliSignal();
609 if (name != "") sxref->SetName(name);
610 sxref->SetTitle("Reference event in ICRS coordinates");
611 sxref->SetTimestamp(*ts);
616 // Convert to horizontal coordinates
617 q=q.GetUnprimed(&fL);
620 SetSignal(&q,"hor",mode,ts,jref);
626 // Convert to "mean" values at specified epoch
627 if (mode=="T" || mode=="t")
630 q=q.GetUnprimed(&fN);
633 // Convert to "mean" values at J2000
635 q=q.GetUnprimed(&fP);
637 // Convert to ICRS values
638 if (!fBias) SetBmatrix();
639 q=q.GetUnprimed(&fB);
644 // Convert to J2000 equatorial mean coordinates
645 if (fGal != 2) SetGmatrix("J");
646 q=q.GetUnprimed(&fG);
648 // Convert to ICRS values
649 if (!fBias) SetBmatrix();
650 q=q.GetUnprimed(&fB);
655 // Convert to mean equatorial values at specified epoch
657 q=q.GetUnprimed(&fE);
659 // Convert to "mean" values at J2000
661 q=q.GetUnprimed(&fP);
663 // Convert to ICRS values
664 if (!fBias) SetBmatrix();
665 q=q.GetUnprimed(&fB);
670 // Convert to "true" equatorial values at the specified timestamp
672 q=q.GetUnprimed(&fH);
674 // Convert to "mean" values at specified timestamp
676 q=q.GetUnprimed(&fN);
678 // Convert to "mean" values at J2000
680 q=q.GetUnprimed(&fP);
682 // Convert to ICRS values
683 if (!fBias) SetBmatrix();
684 q=q.GetUnprimed(&fB);
687 // Store the signal in ICRS coordinates
688 if (!jref) // Storage of a regular signal
690 fXsig->SetPosition(q);
692 else // Storage of a reference signal
694 sxref->SetPosition(q);
701 fRefs->AddAt(sxref,jref-1);
705 ///////////////////////////////////////////////////////////////////////////
706 void AliAstrolab::SetSignal(Double_t a,Double_t d,TString s,Double_t e,TString mode,Int_t jref,TString name)
708 // Store a signal with right ascension (a) and declination (d) given for epoch e.
709 // The position is stored in International Celestial Reference System (ICRS) coordinates.
710 // The ICRS is a fixed, time independent frame and as such provides a unique reference
711 // frame without the need of specifying any epoch etc...
712 // The ICRS coordinate definitions match within 20 mas the mean ones of the J2000.0
713 // equatorial system. Nevertheless, to obtain the highest accuracy, the slight
714 // coordinate correction between J2000 and ICRS is performed here via the
715 // so-called frame bias matrix.
716 // For further details see the U.S. Naval Observatory (USNO) circular 179 (2005),
717 // which is available on http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
719 // The coordinates (a,d) can represent so-called "mean" and "true" values.
720 // The distinction between these two representations is the following :
722 // mean values : (a,d) are only corrected for precession and not for nutation
723 // true values : (a,d) are corrected for both precession and nutation
725 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
726 // values for the input (a,d) coordinates.
728 // a : Right ascension in hhmmss.sss
729 // d : Declination in dddmmss.sss
730 // s = "B" --> Besselian reference epoch.
731 // "J" --> Julian reference epoch.
732 // e : Reference epoch for the input coordinates (e.g. 1900, 1950, 2000,...)
733 // mode = "M" --> Input coordinates are the mean values
734 // "T" --> Input coordinates are the true values
736 // The input parameter "jref" allows the user to store so-called "reference" signals.
737 // These reference signals may be used to check space-time event coincidences with the
738 // stored measurement (e.g. coincidence of the measurement with transient phenomena).
740 // jref = 0 --> Storage of the measurement
741 // j --> Storage of a reference signal at the j-th position (j=1 is first)
742 // < 0 --> Add a reference signal at the next available position
744 // Via the input argument "name" the user can give the stored signal also a name.
746 // The default values are jref=0 and name="".
748 if (s!="B" && s!="b" && s!="J" && s!="j")
758 if (mode!="M" && mode!="m" && mode!="T" && mode!="t")
768 // Convert coordinates to fractional degrees.
769 a=ConvertAngle(a,"hms","deg");
770 d=ConvertAngle(d,"dms","deg");
777 Double_t vec[3]={1.,90.-d,a};
778 r.SetVector(vec,"sph","deg");
780 SetSignal(&r,"equ",mode,&tx,jref,name);
782 ///////////////////////////////////////////////////////////////////////////
783 void AliAstrolab::SetSignal(Double_t a,Double_t d,TString mode,AliTimestamp* ts,Int_t jref,TString name)
785 // Store a signal with right ascension (a) and declination (d) given for timestamp ts.
786 // The position is stored in International Celestial Reference System (ICRS) coordinates.
787 // The ICRS is a fixed, time independent frame and as such provides a unique reference
788 // frame without the need of specifying any epoch etc...
789 // The ICRS coordinate definitions match within 20 mas the mean ones of the J2000.0
790 // equatorial system. Nevertheless, to obtain the highest accuracy, the slight
791 // coordinate correction between J2000 and ICRS is performed here via the
792 // so-called frame bias matrix.
793 // For further details see the U.S. Naval Observatory (USNO) circular 179 (2005),
794 // which is available on http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
796 // The coordinates (a,d) can represent so-called "mean" and "true" values.
797 // The distinction between these two representations is the following :
799 // mean values : (a,d) are only corrected for precession and not for nutation
800 // true values : (a,d) are corrected for both precession and nutation
802 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
803 // values for the input (a,d) coordinates.
805 // a : Right ascension in hhmmss.sss
806 // d : Declination in dddmmss.sss
807 // mode = "M" --> Input coordinates are the mean values
808 // "T" --> Input coordinates are the true values
810 // The input parameter "jref" allows the user to store so-called "reference" signals.
811 // These reference signals may be used to check space-time event coincidences with the
812 // stored measurement (e.g. coincidence of the measurement with transient phenomena).
814 // jref = 0 --> Storage of the measurement
815 // j --> Storage of a reference signal at the j-th position (j=1 is first)
816 // < 0 --> Add a reference signal at the next available position
818 // Via the input argument "name" the user can give the stored signal also a name.
820 // The default values are jref=0 and name="".
822 // Note : In case ts=0 the current timestamp of the lab will be taken.
824 // Convert coordinates to fractional degrees.
825 a=ConvertAngle(a,"hms","deg");
826 d=ConvertAngle(d,"dms","deg");
829 Double_t vec[3]={1.,90.-d,a};
830 r.SetVector(vec,"sph","deg");
832 SetSignal(&r,"equ",mode,ts,jref,name);
834 ///////////////////////////////////////////////////////////////////////////
835 AliSignal* AliAstrolab::GetSignal(Ali3Vector& r,TString frame,TString mode,AliTimestamp* ts,Int_t jref)
837 // Provide the user specified coordinates of a signal at the specific timestamp ts.
838 // The coordinates are returned via the vector argument "r".
839 // In addition also a pointer to the stored signal object is provided.
840 // In case no stored signal was available or one of the input arguments was
841 // invalid, the returned pointer will be 0.
843 // The input parameter "frame" allows the user to specify the frame to which
844 // the components of r refer. Available options are :
846 // frame = "equ" ==> Equatorial coordinates with right ascension (a) and declination (d),
847 // where the "sph" components of r correspond to theta=(pi/2)-d and phi=a.
848 // "gal" ==> Galactic coordinates with longitude (l) and latitude (b).
849 // where the "sph" components of r correspond to theta=(pi/2)-b and phi=l.
850 // "ecl" ==> Ecliptic coordinates with longitude (l) and latitude (b),
851 // where the "sph" components of r correspond to theta=(pi/2)-b and phi=l.
852 // "hor" ==> Horizontal coordinates at the AliAstrolab location, where the "sph"
853 // components of r correspond to theta=zenith angle and phi=pi-azimuth.
854 // "icr" ==> ICRS coordinates with longitude (l) and latitude (b),
855 // where the "sph" components of r correspond to theta=(pi/2)-b and phi=l.
856 // "loc" ==> Local coordinates at the AliAstrolab location, where the "sph"
857 // components of r correspond to the usual theta and phi angles.
859 // In case the coordinates are the equatorial right ascension and declination (a,d),
860 // they can represent so-called "mean" and "true" values.
861 // The distinction between these two representations is the following :
863 // mean values : (a,d) are only corrected for precession and not for nutation
864 // true values : (a,d) are corrected for both precession and nutation
866 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
867 // values for the input in case of equatorial (a,d) coordinates.
869 // mode = "M" --> Input coordinates are the mean values
870 // "T" --> Input coordinates are the true values
872 // The input parameter "jref" allows the user to access so-called "reference" signals.
873 // These reference signals may be used to check space-time event coincidences with the
874 // stored measurement (e.g. coincidence of the measurement with transient phenomena).
876 // jref = 0 --> Access to the measurement
877 // j --> Access to the reference signal at the j-th position (j=1 is first)
879 // Default value is jref=0.
881 // Note : In case ts=0 the current timestamp of the lab will be taken.
885 if (frame!="equ" && frame!="gal" && frame!="ecl" && frame!="hor" && frame!="icr" && frame!="loc") return 0;
887 if (frame=="equ" && mode!="M" && mode!="m" && mode!="T" && mode!="t") return 0;
889 AliSignal* sx=GetSignal(jref);
893 if (!ts) ts=(AliTimestamp*)this;
896 sx->GetPosition(vec,"sph","rad");
898 q.SetVector(vec,"sph","rad");
906 // Convert from ICRS to equatorial J2000 coordinates
907 if (!fBias) SetBmatrix();
912 // Precess to specified timestamp
914 ts1.SetEpoch(2000,"J");
917 // Nutation correction if requested
918 if (mode=="T" || mode=="t") Nutate(q,ts);
923 // Convert from equatorial J2000 to galactic
924 if (fGal != 2) SetGmatrix("J");
930 // Precess to specified timestamp
932 ts1.SetEpoch(2000,"J");
935 // Convert from equatorial to ecliptic coordinates
942 // Precess to specified timestamp
944 ts1.SetEpoch(2000,"J");
947 // Nutation correction
950 // Convert from equatorial to horizontal coordinates
957 // Get the signal in horizontal coordinates
958 GetSignal(q,"hor",mode,ts);
960 // Convert from horizontal local-frame coordinates
967 ///////////////////////////////////////////////////////////////////////////
968 AliSignal* AliAstrolab::GetSignal(Ali3Vector& r,TString frame,TString mode,AliTimestamp* ts,TString name)
970 // Provide the user specified coordinates of the signal with the specified
971 // name at the specific timestamp ts.
972 // The coordinates are returned via the vector argument "r".
973 // In addition also a pointer to the stored signal object is provided.
974 // In case no such stored signal was available or one of the input arguments was
975 // invalid, the returned pointer will be 0.
977 // The input parameter "frame" allows the user to specify the frame to which
978 // the components of r refer. Available options are :
980 // frame = "equ" ==> Equatorial coordinates with right ascension (a) and declination (d),
981 // where the "sph" components of r correspond to theta=(pi/2)-d and phi=a.
982 // "gal" ==> Galactic coordinates with longitude (l) and latitude (b).
983 // where the "sph" components of r correspond to theta=(pi/2)-b and phi=l.
984 // "ecl" ==> Ecliptic coordinates with longitude (l) and latitude (b),
985 // where the "sph" components of r correspond to theta=(pi/2)-b and phi=l.
986 // "hor" ==> Horizontal coordinates at the AliAstrolab location, where the "sph"
987 // components of r correspond to theta=zenith angle and phi=pi-azimuth.
988 // "icr" ==> ICRS coordinates with longitude (l) and latitude (b),
989 // where the "sph" components of r correspond to theta=(pi/2)-b and phi=l.
990 // "loc" ==> Local coordinates at the AliAstrolab location, where the "sph"
991 // components of r correspond to the usual theta and phi angles.
993 // In case the coordinates are the equatorial right ascension and declination (a,d),
994 // they can represent so-called "mean" and "true" values.
995 // The distinction between these two representations is the following :
997 // mean values : (a,d) are only corrected for precession and not for nutation
998 // true values : (a,d) are corrected for both precession and nutation
1000 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
1001 // values for the input in case of equatorial (a,d) coordinates.
1003 // mode = "M" --> Input coordinates are the mean values
1004 // "T" --> Input coordinates are the true values
1006 // Note : In case ts=0 the current timestamp of the lab will be taken.
1009 Int_t j=GetSignalIndex(name);
1010 if (j>=0) sx=GetSignal(r,frame,mode,ts,j);
1013 ///////////////////////////////////////////////////////////////////////////
1014 AliSignal* AliAstrolab::GetSignal(Double_t& a,Double_t& d,TString mode,AliTimestamp* ts,Int_t jref)
1016 // Provide precession (and nutation) corrected right ascension (a) and
1017 // declination (d) of the stored signal object at the specified timestamp.
1018 // In addition also a pointer to the stored signal object is provided.
1019 // In case no stored signal was available or one of the input arguments was
1020 // invalid, the returned pointer will be 0.
1022 // The coordinates (a,d) can represent so-called "mean" and "true" values.
1023 // The distinction between these two representations is the following :
1025 // mean values : (a,d) are only corrected for precession and not for nutation
1026 // true values : (a,d) are corrected for both precession and nutation
1028 // The input parameter "mode" allows the user to select either
1029 // "mean" or "true" values for (a,d).
1031 // The correction methods used are the new IAU 2000 ones as described in the
1032 // U.S. Naval Observatory (USNO) circular 179 (2005), which is available on
1033 // http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
1035 // a : Right ascension in hhmmss.sss
1036 // d : Declination in dddmmss.sss
1037 // mode = "M" --> Output coordinates are the mean values
1038 // "T" --> Output coordinates are the true values
1039 // ts : Timestamp at which the corrected coordinate values are requested.
1041 // The input parameter "jref" allows the user to access so-called "reference" signals.
1042 // These reference signals may be used to check space-time event coincidences with the
1043 // stored regular signal (e.g. coincidence of the stored signal with transient phenomena).
1045 // jref = 0 --> Access to the measurement
1046 // j --> Access to the reference signal at the j-th position (j=1 is first)
1048 // Default value is jref=0.
1050 // Note : In case ts=0 the current timestamp of the lab will be taken.
1056 AliSignal* sx=GetSignal(r,"equ",mode,ts,jref);
1060 // Retrieve the requested (a,d) values
1062 r.GetVector(vec,"sph","deg");
1083 // Convert coordinates to appropriate format
1084 a=ConvertAngle(a,"deg","hms");
1085 d=ConvertAngle(d,"deg","dms");
1089 ///////////////////////////////////////////////////////////////////////////
1090 AliSignal* AliAstrolab::GetSignal(Double_t& a,Double_t& d,TString mode,AliTimestamp* ts,TString name)
1092 // Provide precession (and nutation) corrected right ascension (a) and
1093 // declination (d) of the stored signal object with the specified name
1094 // at the specific timestamp ts.
1095 // In addition also a pointer to the stored signal object is provided.
1096 // In case no stored signal was available or one of the input arguments was
1097 // invalid, the returned pointer will be 0.
1099 // The coordinates (a,d) can represent so-called "mean" and "true" values.
1100 // The distinction between these two representations is the following :
1102 // mean values : (a,d) are only corrected for precession and not for nutation
1103 // true values : (a,d) are corrected for both precession and nutation
1105 // The input parameter "mode" allows the user to select either
1106 // "mean" or "true" values for (a,d).
1108 // The correction methods used are the new IAU 2000 ones as described in the
1109 // U.S. Naval Observatory (USNO) circular 179 (2005), which is available on
1110 // http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
1112 // a : Right ascension in hhmmss.sss
1113 // d : Declination in dddmmss.sss
1114 // mode = "M" --> Output coordinates are the mean values
1115 // "T" --> Output coordinates are the true values
1116 // ts : Timestamp at which the corrected coordinate values are requested.
1117 // name : Name of the requested signal object
1119 // Note : In case ts=0 the current timestamp of the lab will be taken.
1122 Int_t j=GetSignalIndex(name);
1123 if (j>=0) sx=GetSignal(a,d,mode,ts,j);
1126 ///////////////////////////////////////////////////////////////////////////
1127 AliSignal* AliAstrolab::GetSignal(Double_t& a,Double_t& d,TString s,Double_t e,TString mode,Int_t jref)
1129 // Provide precession (and nutation) corrected right ascension (a) and
1130 // declination (d) of the stored signal object at the specified epoch e.
1131 // In addition also a pointer to the stored signal object is provided.
1132 // In case no stored signal was available or one of the input arguments was
1133 // invalid, the returned pointer will be 0.
1135 // The coordinates (a,d) can represent so-called "mean" and "true" values.
1136 // The distinction between these two representations is the following :
1138 // mean values : (a,d) are only corrected for precession and not for nutation
1139 // true values : (a,d) are corrected for both precession and nutation
1141 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
1142 // values for the input (a,d) coordinates.
1144 // a : Right ascension in hhmmss.sss
1145 // d : Declination in dddmmss.sss
1146 // s = "B" --> Besselian reference epoch.
1147 // "J" --> Julian reference epoch.
1148 // e : Reference epoch for the input coordinates (e.g. 1900, 1950, 2000,...)
1149 // mode = "M" --> Input coordinates are the mean values
1150 // "T" --> Input coordinates are the true values
1152 // The input parameter "jref" allows the user to access so-called "reference" signals.
1153 // These reference signals may be used to check space-time event coincidences with the
1154 // stored measurement (e.g. coincidence of the measurement with transient phenomena).
1156 // jref = 0 --> Access to the measurement
1157 // j --> Access to the reference signal at the j-th position (j=1 is first)
1159 // Default value is jref=0.
1164 if (s!="B" && s!="b" && s!="J" && s!="j") return 0;
1166 if (mode!="M" && mode!="m" && mode!="T" && mode!="t") return 0;
1168 // Convert coordinates to fractional degrees.
1169 a=ConvertAngle(a,"hms","deg");
1170 d=ConvertAngle(d,"dms","deg");
1176 AliSignal* sx=GetSignal(a,d,mode,&tx,jref);
1179 ///////////////////////////////////////////////////////////////////////////
1180 AliSignal* AliAstrolab::GetSignal(Double_t& a,Double_t& d,TString s,Double_t e,TString mode,TString name)
1182 // Provide precession (and nutation) corrected right ascension (a) and
1183 // declination (d) of the stored signal object with the specified name
1184 // at the specific epoch e.
1185 // In addition also a pointer to the stored signal object is provided.
1186 // In case no stored signal was available or one of the input arguments was
1187 // invalid, the returned pointer will be 0.
1189 // The coordinates (a,d) can represent so-called "mean" and "true" values.
1190 // The distinction between these two representations is the following :
1192 // mean values : (a,d) are only corrected for precession and not for nutation
1193 // true values : (a,d) are corrected for both precession and nutation
1195 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
1196 // values for the input (a,d) coordinates.
1198 // a : Right ascension in hhmmss.sss
1199 // d : Declination in dddmmss.sss
1200 // s = "B" --> Besselian reference epoch.
1201 // "J" --> Julian reference epoch.
1202 // e : Reference epoch for the input coordinates (e.g. 1900, 1950, 2000,...)
1203 // mode = "M" --> Input coordinates are the mean values
1204 // "T" --> Input coordinates are the true values
1205 // name : Name of the requested signal object
1208 Int_t j=GetSignalIndex(name);
1209 if (j>=0) sx=GetSignal(a,d,s,e,mode,j);
1212 ///////////////////////////////////////////////////////////////////////////
1213 AliSignal* AliAstrolab::GetSignal(Int_t jref)
1215 // Provide the pointer to a stored signal object.
1217 // The input parameter "jref" allows the user to access so-called "reference" signals.
1218 // These reference signals may be used to check space-time event coincidences with the
1219 // stored measurement (e.g. coincidence of the measurement with transient phenomena).
1221 // jref = 0 --> Access to the measurement
1222 // j --> Access to the reference signal at the j-th position (j=1 is first)
1224 // Default value is jref=0.
1233 if (jref>0 && jref<fRefs->GetSize()) sx=(AliSignal*)fRefs->At(jref-1);
1237 ///////////////////////////////////////////////////////////////////////////
1238 AliSignal* AliAstrolab::GetSignal(TString name)
1240 // Provide the pointer to the stored signal object with the specified name.
1243 Int_t j=GetSignalIndex(name);
1244 if (j>=0) sx=GetSignal(j);
1247 ///////////////////////////////////////////////////////////////////////////
1248 void AliAstrolab::RemoveRefSignal(Int_t j,Int_t compress)
1250 // Remove the reference signal which was stored at the j-th position (j=1 is first).
1251 // Note : j=0 means that all stored reference signals will be removed.
1252 // j<0 allows array compression (see below) without removing any signals.
1254 // The "compress" parameter allows compression of the ref. signal storage array.
1256 // compress = 1 --> Array will be compressed
1257 // 0 --> Array will not be compressed
1259 // Note : Compression of the storage array means that the indices of the
1260 // reference signals in the storage array will change.
1264 // Clearing of the complete storage
1272 // Removing a specific reference signal
1273 if (j>0 && j<fRefs->GetSize())
1275 TObject* obj=fRefs->RemoveAt(j-1);
1276 if (obj) delete obj;
1279 // Compression of the storage array
1280 if (compress) fRefs->Compress();
1282 ///////////////////////////////////////////////////////////////////////////
1283 void AliAstrolab::RemoveRefSignal(TString name,Int_t compress)
1285 // Remove the reference signal with the specified name.
1287 // The "compress" parameter allows compression of the ref. signal storage array.
1289 // compress = 1 --> Array will be compressed
1290 // 0 --> Array will not be compressed
1292 // Note : Compression of the storage array means that the indices of the
1293 // reference signals in the storage array will change.
1295 Int_t j=GetSignalIndex(name);
1296 if (j>0) RemoveRefSignal(j,compress);
1298 ///////////////////////////////////////////////////////////////////////////
1299 Int_t AliAstrolab::GetSignalIndex(TString name)
1301 // Provide storage index of the signal with the specified name.
1302 // In case the name matches with the stored measurement,
1303 // the value 0 is returned.
1304 // In case no signal with the specified name was found, the value -1 is returned.
1310 if (name==fXsig->GetName()) return 0;
1313 if (!fRefs) return -1;
1315 for (Int_t i=0; i<fRefs->GetSize(); i++)
1317 AliSignal* sx=(AliSignal*)fRefs->At(i);
1320 if (name==sx->GetName())
1329 ///////////////////////////////////////////////////////////////////////////
1330 void AliAstrolab::PrintSignal(TString frame,TString mode,AliTimestamp* ts,Int_t ndig,Int_t jref)
1332 // Print data of a stored signal in user specified coordinates at the specific timestamp ts.
1333 // In case no stored signal was available or one of the input arguments was
1334 // invalid, no printout is produced.
1336 // The argument "ndig" specifies the number of digits for the fractional
1337 // part (e.g. ndig=6 for "dms" corresponds to micro-arcsecond precision).
1338 // No rounding will be performed, so an arcsecond count of 3.473 with ndig=1
1339 // will appear as 03.4 on the output.
1340 // Due to computer accuracy, precision on the pico-arcsecond level may get lost.
1342 // Note : The angle info is printed without additional spaces or "endline".
1343 // This allows the print to be included in various composite output formats.
1345 // The input parameter "frame" allows the user to specify the frame to which
1346 // the coordinates refer. Available options are :
1348 // frame = "equ" ==> Equatorial coordinates with right ascension (a) and declination (d).
1350 // "gal" ==> Galactic coordinates with longitude (l) and latitude (b).
1352 // "ecl" ==> Ecliptic coordinates with longitude (l) and latitude (b).
1354 // "hor" ==> Horizontal azimuth and altitude coordinates at the AliAstrolab location.
1356 // "icr" ==> ICRS coordinates with longitude (l) and latitude (b).
1358 // "loc" ==> Local spherical angles theta and phi at the AliAstrolab location.
1360 // In case the coordinates are the equatorial right ascension and declination (a,d),
1361 // they can represent so-called "mean" and "true" values.
1362 // The distinction between these two representations is the following :
1364 // mean values : (a,d) are only corrected for precession and not for nutation
1365 // true values : (a,d) are corrected for both precession and nutation
1367 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
1368 // values for the input in case of equatorial (a,d) coordinates.
1370 // mode = "M" --> Input coordinates are the mean values
1371 // "T" --> Input coordinates are the true values
1373 // The input parameter "mode" also determines which type of time and
1374 // local hour angle will appear in the printout.
1376 // mode = "M" --> Mean Sidereal Time (MST) and Local Mean Hour Angle (LMHA)
1377 // "T" --> Apparent Sidereal Time (AST) and Local Apparent Hour Angle (LAHA)
1379 // The input parameter "jref" allows printing of a so-called "reference" signal.
1380 // These reference signals may serve to check space-time event coincidences with the
1381 // stored measurement (e.g. coincidence of the measurement with transient phenomena).
1383 // jref = 0 --> Printing of the measurement
1384 // j --> Printing of the j-th reference signal
1386 // Default value is jref=0.
1388 // Note : In case ts=0 the current timestamp of the lab will be taken.
1391 AliSignal* sx=GetSignal(r,frame,mode,ts,jref);
1395 // Local Hour Angle of the signal
1396 Double_t lha=GetHourAngle("A",ts,jref);
1397 TString slha="LAHA";
1398 if (mode=="M" || mode=="m")
1400 lha=GetHourAngle("M",ts,jref);
1404 TString name=sx->GetName();
1405 if (name != "") cout << name.Data() << " ";
1410 d=90.-r.GetX(2,"sph","deg");
1411 a=r.GetX(3,"sph","rad");
1412 cout << "Equatorial (" << mode.Data() <<") a : "; PrintAngle(a,"rad","hms",ndig);
1413 cout << " d : "; PrintAngle(d,"deg","dms",ndig);
1414 cout << " " << slha.Data() << " : "; PrintAngle(lha,"deg","hms",ndig);
1421 b=90.-r.GetX(2,"sph","deg");
1422 l=r.GetX(3,"sph","deg");
1423 cout << "Galactic l : "; PrintAngle(l,"deg","deg",ndig);
1424 cout << " b : "; PrintAngle(b,"deg","deg",ndig);
1425 cout << " " << slha.Data() << " : "; PrintAngle(lha,"deg","hms",ndig);
1432 d=90.-r.GetX(2,"sph","deg");
1433 a=r.GetX(3,"sph","rad");
1434 cout << "ICRS l : "; PrintAngle(a,"rad","hms",ndig);
1435 cout << " b : "; PrintAngle(d,"deg","dms",ndig);
1436 cout << " " << slha.Data() << " : "; PrintAngle(lha,"deg","hms",ndig);
1443 d=90.-r.GetX(2,"sph","deg");
1444 a=r.GetX(3,"sph","deg");
1445 cout << "Ecliptic l : "; PrintAngle(a,"deg","deg",ndig);
1446 cout << " b : "; PrintAngle(d,"deg","deg",ndig);
1447 cout << " " << slha.Data() << " : "; PrintAngle(lha,"deg","hms",ndig);
1453 Double_t alt=90.-r.GetX(2,"sph","deg");
1454 Double_t azi=180.-r.GetX(3,"sph","deg");
1463 cout << "Horizontal azi : "; PrintAngle(azi,"deg","deg",ndig);
1464 cout << " alt : "; PrintAngle(alt,"deg","deg",ndig);
1465 cout << " " << slha.Data() << " : "; PrintAngle(lha,"deg","hms",ndig);
1471 Double_t theta=r.GetX(2,"sph","deg");
1472 Double_t phi=r.GetX(3,"sph","deg");
1473 cout << "Local-frame phi : "; PrintAngle(phi,"deg","deg",ndig);
1474 cout << " theta : "; PrintAngle(theta,"deg","deg",ndig);
1475 cout << " " << slha.Data() << " : "; PrintAngle(lha,"deg","hms",ndig);
1479 ///////////////////////////////////////////////////////////////////////////
1480 void AliAstrolab::PrintSignal(TString frame,TString mode,AliTimestamp* ts,Int_t ndig,TString name)
1482 // Print data of the stored signal with the specified name in user specified coordinates
1483 // at the specific timestamp ts.
1484 // In case no such stored signal was available or one of the input arguments was
1485 // invalid, no printout is produced.
1487 // The argument "ndig" specifies the number of digits for the fractional
1488 // part (e.g. ndig=6 for "dms" corresponds to micro-arcsecond precision).
1489 // No rounding will be performed, so an arcsecond count of 3.473 with ndig=1
1490 // will appear as 03.4 on the output.
1491 // Due to computer accuracy, precision on the pico-arcsecond level may get lost.
1493 // Note : The angle info is printed without additional spaces or "endline".
1494 // This allows the print to be included in various composite output formats.
1496 // The input parameter "frame" allows the user to specify the frame to which
1497 // the coordinates refer. Available options are :
1499 // frame = "equ" ==> Equatorial coordinates with right ascension (a) and declination (d).
1501 // "gal" ==> Galactic coordinates with longitude (l) and latitude (b).
1503 // "ecl" ==> Ecliptic coordinates with longitude (l) and latitude (b).
1505 // "hor" ==> Horizontal azimuth and altitude coordinates at the AliAstrolab location.
1507 // "icr" ==> ICRS coordinates with longitude (l) and latitude (b).
1509 // "loc" ==> Local spherical angles theta and phi at the AliAstrolab location.
1511 // In case the coordinates are the equatorial right ascension and declination (a,d),
1512 // they can represent so-called "mean" and "true" values.
1513 // The distinction between these two representations is the following :
1515 // mean values : (a,d) are only corrected for precession and not for nutation
1516 // true values : (a,d) are corrected for both precession and nutation
1518 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
1519 // values for the input in case of equatorial (a,d) coordinates.
1521 // mode = "M" --> Input coordinates are the mean values
1522 // "T" --> Input coordinates are the true values
1524 // The input parameter "mode" also determines which type of time and
1525 // local hour angle will appear in the printout.
1527 // mode = "M" --> Mean Sidereal Time (MST) and Local Mean Hour Angle (LMHA)
1528 // "T" --> Apparent Sidereal Time (AST) and Local Apparent Hour Angle (LAHA)
1530 // Note : In case ts=0 the current timestamp of the lab will be taken.
1532 Int_t j=GetSignalIndex(name);
1533 if (j>=0) PrintSignal(frame,mode,ts,ndig,j);
1535 ///////////////////////////////////////////////////////////////////////////
1536 void AliAstrolab::ListSignals(TString frame,TString mode,Int_t ndig)
1538 // List all stored signals in user specified coordinates at the timestamp
1539 // of the actual recording of the stored measurement under investigation.
1540 // In case no (timestamp of the) actual measurement is available,
1541 // the current timestamp of the lab will be taken.
1542 // In case no stored signal is available or one of the input arguments is
1543 // invalid, no printout is produced.
1545 // The argument "ndig" specifies the number of digits for the fractional
1546 // part (e.g. ndig=6 for "dms" corresponds to micro-arcsecond precision).
1547 // No rounding will be performed, so an arcsecond count of 3.473 with ndig=1
1548 // will appear as 03.4 on the output.
1549 // Due to computer accuracy, precision on the pico-arcsecond level may get lost.
1551 // The default value is ndig=1.
1553 // Note : The angle info is printed without additional spaces or "endline".
1554 // This allows the print to be included in various composite output formats.
1556 // The input parameter "frame" allows the user to specify the frame to which
1557 // the coordinates refer. Available options are :
1559 // frame = "equ" ==> Equatorial coordinates with right ascension (a) and declination (d).
1561 // "gal" ==> Galactic coordinates with longitude (l) and latitude (b).
1563 // "ecl" ==> Ecliptic coordinates with longitude (l) and latitude (b).
1565 // "hor" ==> Horizontal azimuth and altitude coordinates at the AliAstrolab location.
1567 // "icr" ==> ICRS coordinates with longitude (l) and latitude (b).
1569 // "loc" ==> Local spherical angles theta and phi at the AliAstrolab location.
1571 // In case the coordinates are the equatorial right ascension and declination (a,d),
1572 // they can represent so-called "mean" and "true" values.
1573 // The distinction between these two representations is the following :
1575 // mean values : (a,d) are only corrected for precession and not for nutation
1576 // true values : (a,d) are corrected for both precession and nutation
1578 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
1579 // values for the input in case of equatorial (a,d) coordinates.
1581 // mode = "M" --> Input coordinates are the mean values
1582 // "T" --> Input coordinates are the true values
1584 // The input parameter "mode" also determines which type of time and
1585 // local hour angle will appear in the listing.
1587 // mode = "M" --> Mean Sidereal Time (MST) and Local Mean Hour Angle (LMHA)
1588 // "T" --> Apparent Sidereal Time (AST) and Local Apparent Hour Angle (LAHA)
1595 if (mode=="T" || mode=="t") dform=-1;
1599 tx=fXsig->GetTimestamp();
1600 if (!tx) tx=(AliTimestamp*)this;
1601 cout << " *AliAstrolab::ListSignals* List of all stored signals." << endl;
1602 cout << " === The measurement under investigation ===" << endl;
1603 cout << " Timestamp of the actual observation" << endl;
1604 tx->Date(dform,fToffset);
1605 cout << " Location of the actual observation" << endl;
1606 cout << " "; PrintSignal(frame,mode,tx,ndig); cout << endl;
1612 for (Int_t i=1; i<=fRefs->GetSize(); i++)
1614 AliSignal* sx=GetSignal(i);
1619 cout << " *AliAstrolab::ListRefSignals* List of all stored signals." << endl;
1620 tx=(AliTimestamp*)this;
1621 cout << " Current timestamp of the laboratory" << endl;
1622 tx->Date(dform,fToffset);
1627 cout << " === All stored reference signals according to the above timestamp ===" << endl;
1630 cout << " Index : " << i << " "; PrintSignal(frame,mode,tx,ndig,i); cout << endl;
1633 ///////////////////////////////////////////////////////////////////////////
1634 void AliAstrolab::Precess(Ali3Vector& r,AliTimestamp* ts1,AliTimestamp* ts2)
1636 // Correct mean right ascension and declination given for Julian date "jd"
1637 // for the earth's precession corresponding to the specified timestamp.
1638 // The results are the so-called "mean" (i.e. precession corrected) values,
1639 // corresponding to the specified timestamp.
1640 // The method used is the new IAU 2000 one as described in the
1641 // U.S. Naval Observatory (USNO) circular 179 (2005), which is available on
1642 // http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
1643 // Since the standard reference epoch is J2000, this implies that all
1644 // input (a,d) coordinates will be first internally converted to the
1645 // corresponding J2000 values before the precession correction w.r.t. the
1646 // specified lab timestamp will be applied.
1648 // r : Input vector containing the right ascension and declination information
1649 // in the form of standard Ali3Vector coordinates.
1650 // In spherical coordinates the phi corresponds to the right ascension,
1651 // whereas the declination corresponds to (pi/2)-theta.
1652 // jd : Julian date corresponding to the input coordinate values.
1653 // ts : Timestamp corresponding to the requested corrected coordinate values.
1655 // Note : In case ts=0 the current timestamp of the lab will be taken.
1657 // Convert back to J2000 values
1660 r0=r.GetUnprimed(&fP);
1662 // Precess to the specified timestamp
1663 if (!ts2) ts2=(AliTimestamp*)this;
1665 r=r0.GetPrimed(&fP);
1667 ///////////////////////////////////////////////////////////////////////////
1668 void AliAstrolab::Nutate(Ali3Vector& r,AliTimestamp* ts)
1670 // Correct mean right ascension and declination for the earth's nutation
1671 // corresponding to the specified timestamp.
1672 // The results are the so-called "true" (i.e. nutation corrected) values,
1673 // corresponding to the specified timestamp.
1674 // The method used is the new IAU 2000 one as described in the
1675 // U.S. Naval Observatory (USNO) circular 179 (2005), which is available on
1676 // http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
1678 // r : Input vector containing the right ascension and declination information
1679 // in the form of standard Ali3Vector coordinates.
1680 // In spherical coordinates the phi corresponds to the right ascension,
1681 // whereas the declination corresponds to (pi/2)-theta.
1682 // ts : Timestamp for which the corrected coordinate values are requested.
1684 // Note : In case ts=0 the current timestamp of the lab will be taken.
1686 // Nutation correction for the specified timestamp
1687 if (!ts) ts=(AliTimestamp*)this;
1691 ///////////////////////////////////////////////////////////////////////////
1692 void AliAstrolab::SetBmatrix()
1694 // Set the frame bias matrix elements.
1695 // The formulas and numerical constants used are the ones from the
1696 // U.S. Naval Observatory (USNO) circular 179 (2005), which is available on
1697 // http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
1699 Double_t pi=acos(-1.);
1701 // Parameters in mas
1703 Double_t x=-16.6170;
1706 // Convert to radians
1707 a*=pi/(180.*3600.*1000.);
1708 x*=pi/(180.*3600.*1000.);
1709 e*=pi/(180.*3600.*1000.);
1712 mat[0]=1.-0.5*(a*a+x*x);
1716 mat[4]=1.-0.5*(a*a+e*e);
1720 mat[8]=1.-0.5*(e*e+x*x);
1725 ///////////////////////////////////////////////////////////////////////////
1726 void AliAstrolab::SetPmatrix(AliTimestamp* ts)
1728 // Set precession matrix elements for Julian date jd w.r.t. J2000.
1729 // The formulas and numerical constants used are the ones from the
1730 // U.S. Naval Observatory (USNO) circular 179 (2005), which is available on
1731 // http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
1732 // All numerical constants refer to the standard reference epoch J2000.
1734 Double_t mat[9]={0,0,0,0,0,0,0,0,0};
1741 Double_t pi=acos(-1.);
1743 Double_t t=(ts->GetJD()-2451545.0)/36525.; // Julian centuries since J2000.0
1745 // Parameters for the precession matrix in arcseconds
1746 Double_t eps0=84381.406; // Mean ecliptic obliquity at J2000.0
1747 Double_t psi=5038.481507*t-1.0790069*pow(t,2)-0.00114045*pow(t,3)+0.000132851*pow(t,4)
1748 -0.0000000951*pow(t,4);
1749 Double_t om=eps0-0.025754*t+0.0512623*pow(t,2)-0.00772503*pow(t,3)-0.000000467*pow(t,4)
1750 +0.0000003337*pow(t,5);
1751 Double_t chi=10.556403*t-2.3814292*pow(t,2)-0.00121197*pow(t,3)+0.000170663*pow(t,4)
1752 -0.0000000560*pow(t,5);
1754 // Convert to radians
1755 eps0*=pi/(180.*3600.);
1756 psi*=pi/(180.*3600.);
1757 om*=pi/(180.*3600.);
1758 chi*=pi/(180.*3600.);
1760 Double_t s1=sin(eps0);
1761 Double_t s2=sin(-psi);
1762 Double_t s3=sin(-om);
1763 Double_t s4=sin(chi);
1764 Double_t c1=cos(eps0);
1765 Double_t c2=cos(-psi);
1766 Double_t c3=cos(-om);
1767 Double_t c4=cos(chi);
1769 mat[0]=c4*c2-s2*s4*c3;
1770 mat[1]=c4*s2*c1+s4*c3*c2*c1-s1*s4*s3;
1771 mat[2]=c4*s2*s1+s4*c3*c2*s1+c1*s4*s3;
1772 mat[3]=-s4*c2-s2*c4*c3;
1773 mat[4]=-s4*s2*c1+c4*c3*c2*c1-s1*c4*s3;
1774 mat[5]=-s4*s2*s1+c4*c3*c2*s1+c1*c4*s3;
1776 mat[7]=-s3*c2*c1-s1*c3;
1777 mat[8]=-s3*c2*s1+c3*c1;
1781 ///////////////////////////////////////////////////////////////////////////
1782 void AliAstrolab::SetNmatrix(AliTimestamp* ts)
1784 // Set nutation matrix elements for the specified Julian date jd.
1785 // The formulas and numerical constants used are the ones from the
1786 // U.S. Naval Observatory (USNO) circular 179 (2005), which is available on
1787 // http://aa.usno,navy.mil/publications/docs/Circular_179.pdf.
1789 Double_t mat[9]={0,0,0,0,0,0,0,0,0};
1796 Double_t pi=acos(-1.);
1798 Double_t dpsi,deps,eps;
1799 ts->Almanac(&dpsi,&deps,&eps);
1801 // Convert to radians
1802 dpsi*=pi/(180.*3600.);
1803 deps*=pi/(180.*3600.);
1804 eps*=pi/(180.*3600.);
1806 Double_t s1=sin(eps);
1807 Double_t s2=sin(-dpsi);
1808 Double_t s3=sin(-(eps+deps));
1809 Double_t c1=cos(eps);
1810 Double_t c2=cos(-dpsi);
1811 Double_t c3=cos(-(eps+deps));
1817 mat[4]=c3*c2*c1-s1*s3;
1818 mat[5]=c3*c2*s1+c1*s3;
1820 mat[7]=-s3*c2*c1-s1*c3;
1821 mat[8]=-s3*c2*s1+c3*c1;
1825 ///////////////////////////////////////////////////////////////////////////
1826 void AliAstrolab::SetGmatrix(TString mode)
1828 // Set the mean equatorial to galactic coordinate conversion matrix.
1829 // The B1950 parameters were taken from section 22.3 of the book
1830 // "An Introduction to Modern Astrophysics" by Carrol and Ostlie (1996).
1831 // The J2000 parameters are obtained by precession of the B1950 values.
1833 // Via the input argument "mode" the required epoch can be selected
1834 // mode = "B" ==> B1950
1837 Ali3Vector x; // The Galactic x-axis in the equatorial frame
1838 Ali3Vector y; // The Galactic y-axis in the equatorial frame
1839 Ali3Vector z; // The Galactic z-axis in the equatorial frame
1842 Double_t vec[3]={1,0,0};
1844 fGal=1; // Set flag to indicate B1950 matrix values
1846 // B1950 equatorial coordinates of the North Galactic Pole (NGP)
1849 a=ConvertAngle(a,"hms","deg");
1850 d=ConvertAngle(d,"dms","deg");
1853 z.SetVector(vec,"sph","deg");
1855 // B1950 equatorial coordinates of the Galactic l=b=0 point
1858 a=ConvertAngle(a,"hms","deg");
1859 d=ConvertAngle(d,"dms","deg");
1862 x.SetVector(vec,"sph","deg");
1864 // Precess to the corresponding J2000 values if requested
1867 fGal=2; // Set flag to indicate J2000 matrix values
1869 t1.SetEpoch(1950,"B");
1871 t2.SetEpoch(2000,"J");
1876 // The Galactic y-axis is determined for the right handed frame
1879 fG.SetAngles(x.GetX(2,"sph","deg"),x.GetX(3,"sph","deg"),
1880 y.GetX(2,"sph","deg"),y.GetX(3,"sph","deg"),
1881 z.GetX(2,"sph","deg"),z.GetX(3,"sph","deg"));
1883 ///////////////////////////////////////////////////////////////////////////
1884 void AliAstrolab::SetEmatrix(AliTimestamp* ts)
1886 // Set the mean equatorial to ecliptic coordinate conversion matrix
1887 // for the specified timestamp.
1888 // A nice sketch and explanation of the two frames can be found
1889 // in chapter 3 of the book "Astronomy Methods" by Hale Bradt (2004).
1891 Double_t dpsi,deps,eps;
1892 ts->Almanac(&dpsi,&deps,&eps);
1894 // Convert to degrees
1897 // Positions of the ecliptic axes w.r.t. the equatorial ones
1898 // at the moment of the specified timestamp
1899 Double_t theta1=90; // Ecliptic x-axis
1901 Double_t theta2=90.-eps; //Ecliptic y-axis
1903 Double_t theta3=eps; // Ecliptic z-axis
1906 fE.SetAngles(theta1,phi1,theta2,phi2,theta3,phi3);
1908 ///////////////////////////////////////////////////////////////////////////
1909 void AliAstrolab::SetHmatrix(AliTimestamp* ts)
1911 // Set the mean equatorial to horizontal coordinate conversion matrix
1912 // for the specified timestamp.
1913 // A nice sketch and explanation of the two frames can be found
1914 // in chapter 3 of the book "Astronomy Methods" by Hale Bradt (2004).
1916 // Note : In order to simplify the calculations, we use here a
1917 // right-handed horizontal frame.
1919 Ali3Vector x; // The (South pointing) horizontal x-axis in the equatorial frame
1920 Ali3Vector y; // The (East pointing) horizontal y-axis in the equatorial frame
1921 Ali3Vector z; // The (Zenith pointing) horizontal z-axis in the equatorial frame
1924 GetLabPosition(l,b,"deg");
1927 Double_t vec[3]={1,0,0};
1929 // Equatorial coordinates of the horizontal z-axis
1930 // at the moment of the specified timestamp
1931 a=ts->GetLAST(fToffset);
1932 a*=15.; // Convert fractional hours to degrees
1935 z.SetVector(vec,"sph","deg");
1937 // Equatorial coordinates of the horizontal x-axis
1938 // at the moment of the specified timestamp
1941 x.SetVector(vec,"sph","deg");
1943 // The horizontal y-axis is determined for the right handed frame
1946 fH.SetAngles(x.GetX(2,"sph","deg"),x.GetX(3,"sph","deg"),
1947 y.GetX(2,"sph","deg"),y.GetX(3,"sph","deg"),
1948 z.GetX(2,"sph","deg"),z.GetX(3,"sph","deg"));
1950 ///////////////////////////////////////////////////////////////////////////
1951 void AliAstrolab::SetLocalFrame(Double_t t1,Double_t p1,Double_t t2,Double_t p2,Double_t t3,Double_t p3)
1953 // Specification of the orientations of the local-frame axes.
1954 // The input arguments represent the angles (in degrees) of the local-frame axes
1955 // w.r.t. a so called Master Reference Frame (MRF), with the same convention
1956 // as the input arguments of TRotMatrix::SetAngles.
1958 // The right handed Master Reference Frame is defined as follows :
1959 // Z-axis : Points to the local Zenith
1960 // X-axis : Makes an angle of 90 degrees with the Z-axis and points South
1961 // Y-axis : Makes an angle of 90 degrees with the Z-axis and points East
1963 // Once the user has specified the local reference frame, any observed event
1964 // can be related to astronomical space-time locations via the SetSignal
1965 // and GetSignal memberfunctions.
1967 // Set the matrix for the conversion of our reference frame coordinates
1968 // into the local-frame ones.
1970 fL.SetAngles(t1,p1,t2,p2,t3,p3);
1972 ///////////////////////////////////////////////////////////////////////////
1973 Double_t AliAstrolab::ConvertAngle(Double_t a,TString in,TString out) const
1975 // Conversion of various angular formats.
1977 // The input argument "a" denotes the angle to be converted.
1978 // The string arguments "in" and "out" specify the angular I/O formats.
1980 // in = "rad" : input angle provided in radians
1981 // "deg" : input angle provided in degrees
1982 // "dms" : input angle provided in dddmmss.sss
1983 // "hms" : input angle provided in hhmmss.sss
1984 // "hrs" : input angle provided in fractional hours
1986 // out = "rad" : output angle provided in radians
1987 // "deg" : output angle provided in degrees
1988 // "dms" : output angle provided in dddmmss.sss
1989 // "hms" : output angle provided in hhmmss.sss
1990 // "hrs" : output angle provided in fractional hours
1992 if (in==out) return a;
1994 // Convert input to its absolute value in (fractional) degrees.
1995 Double_t pi=acos(-1.);
1996 Double_t epsilon=1.e-12; // Accuracy in (arc)seconds
1997 Int_t word=0,ddd=0,hh=0,mm=0,ss=0;
2002 if (in=="rad") b*=180./pi;
2004 if (in=="hrs") b*=15.;
2013 s=b-Double_t(ddd*10000+mm*100+ss);
2014 b=Double_t(ddd)+Double_t(mm)/60.+(Double_t(ss)+s)/3600.;
2024 s=b-Double_t(hh*10000+mm*100+ss);
2025 b=15.*(Double_t(hh)+Double_t(mm)/60.+(Double_t(ss)+s)/3600.);
2033 if (out=="rad") b*=pi/180.;
2035 if (out=="hrs") b/=15.;
2066 b=Double_t(10000*ddd+100*mm+ss)+s;
2099 b=Double_t(10000*hh+100*mm+ss)+s;
2106 ///////////////////////////////////////////////////////////////////////////
2107 Double_t AliAstrolab::GetHourAngle(TString mode,AliTimestamp* ts,Int_t jref)
2109 // Provide the Local Hour Angle (in fractional degrees) of a stored signal
2110 // object at the specified timestamp.
2112 // The input parameter "mode" allows the user to select either the
2113 // "mean" or "apparent" value for the returned Hour Angle.
2115 // mode = "M" --> Output is the Mean Hour Angle
2116 // "A" --> Output is the Apparent Hour Angle
2117 // ts : Timestamp at which the hour angle is requested.
2119 // The input parameter "jref" allows the user to specify so-called "reference" signals.
2120 // These reference signals may be used to check space-time event coincidences with the
2121 // stored measurement (e.g. coincidence of the measurement with transient phenomena).
2123 // jref = 0 --> Use the stored measurement
2124 // j --> Use the reference signal at the j-th position (j=1 is first)
2126 // Default value is jref=0.
2128 // Note : In case ts=0 the current timestamp of the lab will be taken.
2130 if (!ts) ts=(AliTimestamp*)this;
2132 // Get corrected right ascension and declination for the specified timestamp.
2134 if (mode=="M" || mode=="m") GetSignal(a,d,"M",ts,jref);
2135 if (mode=="A" || mode=="a") GetSignal(a,d,"T",ts,jref);
2137 // Convert coordinates to fractional degrees.
2138 a=ConvertAngle(a,"hms","deg");
2139 d=ConvertAngle(d,"dms","deg");
2141 a/=15.; // Convert a to fractional hours
2143 if (mode=="M" || mode=="m") ha=ts->GetLMST(fToffset)-a;
2144 if (mode=="A" || mode=="a") ha=ts->GetLAST(fToffset)-a;
2145 ha*=15.; // Convert to (fractional) degrees
2149 ///////////////////////////////////////////////////////////////////////////
2150 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)
2152 // Set the AliTimestamp parameters corresponding to the LT date and time
2153 // in the Gregorian calendar as specified by the input arguments.
2155 // The input arguments represent the following :
2156 // y : year in LT (e.g. 1952, 2003 etc...)
2157 // m : month in LT (1=jan 2=feb etc...)
2158 // d : day in LT (1-31)
2159 // hh : elapsed hours in LT (0-23)
2160 // mm : elapsed minutes in LT (0-59)
2161 // ss : elapsed seconds in LT (0-59)
2162 // ns : remaining fractional elapsed second of LT in nanosecond
2163 // ps : remaining fractional elapsed nanosecond of LT in picosecond
2165 // Note : ns=0 and ps=0 are the default values.
2168 SetLT(fToffset,y,m,d,hh,mm,ss,ns,ps);
2170 ///////////////////////////////////////////////////////////////////////////
2171 void AliAstrolab::SetLT(Int_t y,Int_t d,Int_t s,Int_t ns,Int_t ps)
2173 // Set the AliTimestamp parameters corresponding to the specified elapsed
2174 // timespan since the beginning of the new LT year.
2176 // The LT year and elapsed time span is entered via the following input arguments :
2178 // y : year in LT (e.g. 1952, 2003 etc...)
2179 // d : elapsed number of days
2180 // s : (remaining) elapsed number of seconds
2181 // ns : (remaining) elapsed number of nanoseconds
2182 // ps : (remaining) elapsed number of picoseconds
2184 // The specified d, s, ns and ps values will be used in an additive
2185 // way to determine the elapsed timespan.
2186 // So, specification of d=1, s=100, ns=0, ps=0 will result in the
2187 // same elapsed time span as d=0, s=24*3600+100, ns=0, ps=0.
2188 // However, by making use of the latter the user should take care
2189 // of possible integer overflow problems in the input arguments,
2190 // which obviously will provide incorrect results.
2192 // Note : ns=0 and ps=0 are the default values.
2194 SetLT(fToffset,y,d,s,ns,ps);
2196 ///////////////////////////////////////////////////////////////////////////
2197 Double_t AliAstrolab::GetDifference(Int_t j,TString au,Double_t& dt,TString tu,Int_t mode,Int_t* ia,Int_t* it)
2199 // Provide space and time difference between the j-th reference signal
2200 // (j=1 indicates first) and the stored measurement.
2202 // The return value of this memberfunction provides the positional angular
2203 // difference, whereas the output argument "dt" provides the time difference.
2205 // The units of the angular difference can be specified via the the "au"
2206 // input argument, where
2208 // au = "rad" --> Angular difference in (fractional) radians
2209 // "deg" --> Angular difference in (fractional) degrees
2211 // The units of the time difference can be specified via the "tu" and "mode"
2212 // input arguments. For details please refer to AliTimestamp::GetDifference().
2213 // Also here mode=1 is the default value.
2215 // For the time difference the reference signal is used as the standard.
2216 // This means that in case of a positive time difference, the stored
2217 // measurement occurred later than the reference signal.
2219 // In case j=0, the stored measurement will be compared with each
2220 // reference signal and the returned angular and time differences are
2221 // the minimal differences which were encountered.
2222 // In this case the user may obtain the indices of the two stored reference signals
2223 // which had the minimal angular and minimal time difference via the output
2224 // arguments "ia" and "it" as follows :
2226 // ia = Index of the stored reference signal with minimial angular difference
2227 // it = Index of the stored reference signal with minimial time difference
2229 // In case these indices are the same, there obviously was 1 single reference signal
2230 // which showed both the minimal angular and time difference.
2232 // The default values are mode=1, ia=0 and it=0;
2239 if (!fRefs) return dang;
2241 Ali3Vector rx; // Position of the measurement
2242 Ali3Vector r0; // Position of the reference signal
2244 AliSignal* sx=GetSignal(rx,"icr","M",0);
2245 if (!sx) return dang;
2247 AliTimestamp* tx=sx->GetTimestamp();
2248 if (!tx) return dang;
2250 // Space and time difference w.r.t. a specific reference signal
2253 AliSignal* s0=GetSignal(r0,"icr","M",0,j);
2254 if (!s0) return dang;
2255 AliTimestamp* t0=s0->GetTimestamp();
2257 if (!t0) return dang;
2259 dang=r0.GetOpeningAngle(rx,au);
2260 dt=t0->GetDifference(tx,tu,mode);
2264 // Minimal space and time difference encountered over all reference signals
2265 Double_t dangmin=dang;
2267 for (Int_t i=1; i<=fRefs->GetSize(); i++)
2269 AliSignal* s0=GetSignal(r0,"icr","M",0,i);
2271 AliTimestamp* t0=s0->GetTimestamp();
2273 dang=r0.GetOpeningAngle(rx,au);
2274 dt=t0->GetDifference(tx,tu,mode);
2275 if (fabs(dang)<dangmin)
2290 ///////////////////////////////////////////////////////////////////////////
2291 Double_t AliAstrolab::GetDifference(TString name,TString au,Double_t& dt,TString tu,Int_t mode)
2293 // Provide space and time difference between the reference signal
2294 // with the specified name and the stored measurement.
2296 // The return value of this memberfunction provides the positional angular
2297 // difference, whereas the output argument "dt" provides the time difference.
2299 // The units of the angular difference can be specified via the the "au"
2300 // input argument, where
2302 // au = "rad" --> Angular difference in (fractional) radians
2303 // "deg" --> Angular difference in (fractional) degrees
2305 // The units of the time difference can be specified via the "tu" and "mode"
2306 // input arguments. For details please refer to AliTimestamp::GetDifference().
2307 // Also here mode=1 is the default value.
2309 // For the time difference the reference signal is used as the standard.
2310 // This means that in case of a positive time difference, the stored
2311 // measurement occurred later than the reference signal.
2316 Int_t j=GetSignalIndex(name);
2317 if (j>0) dang=GetDifference(j,au,dt,tu,mode);
2320 ///////////////////////////////////////////////////////////////////////////
2321 TArrayI* AliAstrolab::MatchRefSignal(Double_t da,TString au,Double_t dt,TString tu,Int_t mode)
2323 // Provide the storage indices of the reference signals which match in space
2324 // and time with the stored measurement.
2325 // The indices are returned via a pointer to a TArrayI object.
2326 // In case no matches were found, the null pointer is returned.
2327 // A reference signal is regarded as matching with the stored measurement
2328 // if the positional angular difference doesn't exceed "da" and the absolute
2329 // value of the time difference doesn't exceed "dt".
2331 // The units of the angular difference "da" can be specified via the the "au"
2332 // input argument, where
2334 // au = "rad" --> Angular difference in (fractional) radians
2335 // "deg" --> Angular difference in (fractional) degrees
2337 // The units of the time difference "dt" can be specified via the "tu" and "mode"
2338 // input arguments. For details please refer to AliTimestamp::GetDifference().
2339 // Also here mode=1 is the default value.
2341 if (!fXsig || !fRefs) return 0;
2343 Int_t nrefs=fRefs->GetEntries();
2345 if (fIndices) delete fIndices;
2346 fIndices=new TArrayI(nrefs);
2348 Double_t dang,dtime;
2350 for (Int_t i=1; i<=fRefs->GetSize(); i++)
2352 dang=GetDifference(i,au,dtime,tu,mode);
2353 if (fabs(dang)<=da && fabs(dtime)<=dt)
2355 fIndices->AddAt(i,jfill);
2360 fIndices->Set(jfill);
2363 ///////////////////////////////////////////////////////////////////////////
2364 void AliAstrolab::DisplaySignal(TString frame,TString mode,AliTimestamp* ts,Int_t jref,TString proj,Int_t clr)
2366 // Display a stored signal in a user specified coordinate projection
2367 // at the specific timestamp ts.
2369 // In case no stored signal was available or one of the input arguments was
2370 // invalid, no display is produced.
2372 // Note : In case ts=0 the current timestamp of the lab will be taken.
2374 // The input parameter "frame" allows the user to specify the frame to which
2375 // the coordinates refer. Available options are :
2377 // frame = "equ" ==> Equatorial coordinates with right ascension (a) and declination (d).
2379 // "gal" ==> Galactic coordinates with longitude (l) and latitude (b).
2381 // "ecl" ==> Ecliptic coordinates with longitude (l) and latitude (b).
2383 // "hor" ==> Horizontal azimuth and altitude coordinates at the AliAstrolab location.
2385 // "icr" ==> ICRS coordinates with longitude (l) and latitude (b).
2387 // "loc" ==> Local spherical angles theta and phi at the AliAstrolab location.
2389 // In case the coordinates are the equatorial right ascension and declination (a,d),
2390 // they can represent so-called "mean" and "true" values.
2391 // The distinction between these two representations is the following :
2393 // mean values : (a,d) are only corrected for precession and not for nutation
2394 // true values : (a,d) are corrected for both precession and nutation
2396 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
2397 // values for the input in case of equatorial (a,d) coordinates.
2399 // mode = "M" --> Input coordinates are the mean values
2400 // "T" --> Input coordinates are the true values
2402 // The input parameter "jref" allows display of a so-called "reference" signal.
2403 // These reference signals may serve to check space-time event coincidences with the
2404 // stored measurement (e.g. coincidence of the measurement with transient phenomena).
2406 // jref = 0 --> Display of the measurement
2407 // j --> Display of the j-th reference signal
2409 // Default value is jref=0.
2411 // The input parameter "proj" allows the user to specify the desired projection.
2412 // The available projection modes are :
2414 // cyl : Cylindrical projection plotted with colored markers
2415 // cylh : Cylindrical projection plotted in a 2-D histogram
2416 // ham : Hammer equal area projection plotted with colored markers
2417 // hamh : Hammer equal area projection plotted in a 2-D histogram
2418 // ait : Aitoff projection plotted with colored markers
2419 // aith : Aitoff projection plotted in a 2-D histogram
2420 // mer : Mercator projection plotted with colored markers
2421 // merh : Mercator projection plotted in a 2-D histogram
2422 // ang : Straight cos(b) vs. l plot with colored markers
2423 // angh : Straight cos(b) vs. l plot in a 2-D histogram
2425 // Note : The ang(h) plot allows for easy identification of a homogeneous distribution.
2427 // The input argument "clr" allows to clear (1) the display before drawing or not (0).
2429 // The default values are : jref=0, proj="ham" and clr=0.
2431 // This routine is based on the work by Garmt de Vries-Uiterweerd.
2434 AliSignal* sx=GetSignal(r,frame,mode,ts,jref);
2438 // The generic input angles (in rad) for the projections
2442 Double_t pi=acos(-1.);
2444 if (frame=="equ" || frame=="gal" || frame=="icr" || frame=="ecl" || frame=="loc")
2446 theta=(pi/2.)-r.GetX(2,"sph","rad");
2447 phi=r.GetX(3,"sph","rad");
2452 theta=(pi/2.)-r.GetX(2,"sph","rad");
2453 phi=pi-r.GetX(3,"sph","rad");
2456 // Automatic choice of central meridian if not selected by the user
2469 // Obtain the projected (x,y) position
2472 Project(phi,theta,proj,x,y);
2475 if (proj=="hamh" || proj=="aith" || proj=="merh" || proj=="cylh" || proj=="angh") hist=1;
2477 // Update the display for this signal position
2479 // Create a new canvas if needed
2480 if (!fCanvas) fCanvas=new TCanvas("AliAstrolab","Skymap");
2482 // Construct the title string for this map
2483 TString title="Projection : ";
2487 title+=" Center : ";
2489 if (frame=="equ" || frame=="icr")
2491 ang=int(ConvertAngle(fMeridian,"rad","hms"));
2505 ang=int(ConvertAngle(fMeridian,"rad","dms"));
2518 if (!hist) // 2-D Marker display (i.e. not a histogram)
2520 // Remove existing markers, grid and outline from display if needed
2521 if (clr==1 || proj!=fProj)
2534 fMarkers=new TObjArray();
2535 fMarkers->SetOwner();
2536 // Set canvas range, header and axes
2539 fCanvas->Range(-2.2,-6.,2.2,6.);
2540 TText* header=new TText;
2541 fMarkers->Add(header);
2542 header->DrawText(-1.4,5.5,title.Data());
2543 TLine* line=new TLine(-2,0,2,0);
2544 fMarkers->Add(line);
2546 line=new TLine(0,5.5,0,-5.5);
2547 fMarkers->Add(line);
2552 fCanvas->Range(-2.2,-1.2,2.2,1.2);
2553 TText* header=new TText;
2554 fMarkers->Add(header);
2555 header->DrawText(-1.4,1.1,title.Data());
2556 TLine* line=new TLine(-2,0,2,0);
2557 fMarkers->Add(line);
2559 line=new TLine(0,1,0,-1);
2560 fMarkers->Add(line);
2564 if (proj=="ham" || proj=="ait")
2566 // Draw ellips outline with axes
2567 TEllipse* outline=new TEllipse(0,0,2,1);
2568 fMarkers->Add(outline);
2573 if (!jref) // The measurement
2575 marker=new TMarker(x,y,29);
2576 marker->SetMarkerColor(kRed);
2578 else // The reference signal(s)
2580 marker=new TMarker(x,y,20);
2581 marker->SetMarkerColor(kBlue);
2583 marker->SetMarkerSize(1);
2584 fMarkers->Add(marker);
2587 else if (hist==1) // 2-D display via histogram
2589 // Reset the histogram if needed
2590 if (clr==1 || proj!=fProj)
2593 if (!fHist) fHist=new TH2F();
2595 fHist->SetMarkerStyle(20);
2596 fHist->SetMarkerSize(1);
2597 fHist->SetMarkerColor(kBlue);
2598 fHist->SetNameTitle("SkyMap",title.Data());
2601 fHist->SetBins(100,-2.2,2.2,100,-5.5,5.5);
2605 fHist->SetBins(100,-2.2,2.2,100,-1.1,1.1);
2613 ///////////////////////////////////////////////////////////////////////////
2614 void AliAstrolab::DisplaySignal(TString frame,TString mode,AliTimestamp* ts,TString name,TString proj,Int_t clr)
2616 // Display the stored signal with the specified name in a user specified
2617 // coordinate projection at the specific timestamp ts.
2619 // Note : In case ts=0 the current timestamp of the lab will be taken.
2621 // In case no such stored signal was available or one of the input arguments was
2622 // invalid, no display is produced.
2624 // The input parameter "frame" allows the user to specify the frame to which
2625 // the coordinates refer. Available options are :
2627 // frame = "equ" ==> Equatorial coordinates with right ascension (a) and declination (d).
2629 // "gal" ==> Galactic coordinates with longitude (l) and latitude (b).
2631 // "ecl" ==> Ecliptic coordinates with longitude (l) and latitude (b).
2633 // "hor" ==> Horizontal azimuth and altitude coordinates at the AliAstrolab location.
2635 // "icr" ==> ICRS coordinates with longitude (l) and latitude (b).
2637 // "loc" ==> Local spherical angles theta and phi at the AliAstrolab location.
2639 // In case the coordinates are the equatorial right ascension and declination (a,d),
2640 // they can represent so-called "mean" and "true" values.
2641 // The distinction between these two representations is the following :
2643 // mean values : (a,d) are only corrected for precession and not for nutation
2644 // true values : (a,d) are corrected for both precession and nutation
2646 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
2647 // values for the input in case of equatorial (a,d) coordinates.
2649 // mode = "M" --> Input coordinates are the mean values
2650 // "T" --> Input coordinates are the true values
2652 // The input parameter "proj" allows the user to specify the desired projection.
2653 // The available projection modes are :
2655 // cyl : Cylindrical projection plotted with colored markers
2656 // cylh : Cylindrical projection plotted in a 2-D histogram
2657 // ham : Hammer equal area projection plotted with colored markers
2658 // hamh : Hammer equal area projection plotted in a 2-D histogram
2659 // ait : Aitoff projection plotted with colored markers
2660 // aith : Aitoff projection plotted in a 2-D histogram
2661 // mer : Mercator projection plotted with colored markers
2662 // merh : Mercator projection plotted in a 2-D histogram
2663 // ang : Straight cos(b) vs. l plot with colored markers
2664 // angh : Straight cos(b) vs. l plot in a 2-D histogram
2666 // Note : The ang(h) plot allows for easy identification of a homogeneous distribution.
2668 // The input argument "clr" allows to clear (1) the display before drawing or not (0).
2670 // The default values are : proj="ham" and clr=0.
2672 // This routine is based on the work by Garmt de Vries-Uiterweerd.
2674 Int_t j=GetSignalIndex(name);
2675 if (j>=0) DisplaySignal(frame,mode,ts,j,proj,clr);
2677 ///////////////////////////////////////////////////////////////////////////
2678 void AliAstrolab::DisplaySignals(TString frame,TString mode,AliTimestamp* ts,TString proj,Int_t clr)
2680 // Display of all stored signals in user specified coordinate projection
2681 // at the specific timestamp ts.
2682 // In case ts=0, the timestamp of the actual recording of the stored measurement
2683 // under investigation will be used.
2684 // In case no (timestamp of the) actual measurement is available,
2685 // the current timestamp of the lab will be taken.
2686 // In case no stored signal is available or one of the input arguments is
2687 // invalid, no display is produced.
2689 // The input parameter "frame" allows the user to specify the frame to which
2690 // the coordinates refer. Available options are :
2692 // frame = "equ" ==> Equatorial coordinates with right ascension (a) and declination (d).
2694 // "gal" ==> Galactic coordinates with longitude (l) and latitude (b).
2696 // "ecl" ==> Ecliptic coordinates with longitude (l) and latitude (b).
2698 // "hor" ==> Horizontal azimuth and altitude coordinates at the AliAstrolab location.
2700 // "icr" ==> ICRS coordinates with longitude (l) and latitude (b).
2702 // "loc" ==> Local spherical angles theta and phi at the AliAstrolab location.
2704 // In case the coordinates are the equatorial right ascension and declination (a,d),
2705 // they can represent so-called "mean" and "true" values.
2706 // The distinction between these two representations is the following :
2708 // mean values : (a,d) are only corrected for precession and not for nutation
2709 // true values : (a,d) are corrected for both precession and nutation
2711 // The input parameter "mode" allows the user to specifiy either "mean" or "true"
2712 // values for the input in case of equatorial (a,d) coordinates.
2714 // mode = "M" --> Input coordinates are the mean values
2715 // "T" --> Input coordinates are the true values
2717 // The input parameter "proj" allows the user to specify the desired projection.
2718 // The available projection modes are :
2720 // cyl : Cylindrical projection plotted with colored markers
2721 // cylh : Cylindrical projection plotted in a 2-D histogram
2722 // ham : Hammer equal area projection plotted with colored markers
2723 // hamh : Hammer equal area projection plotted in a 2-D histogram
2724 // ait : Aitoff projection plotted with colored markers
2725 // aith : Aitoff projection plotted in a 2-D histogram
2726 // mer : Mercator projection plotted with colored markers
2727 // merh : Mercator projection plotted in a 2-D histogram
2728 // ang : Straight cos(b) vs. l plot with colored markers
2729 // angh : Straight cos(b) vs. l plot in a 2-D histogram
2731 // Note : The ang(h) plot allows for easy identification of a homogeneous distribution.
2733 // The input argument "clr" allows to clear (1) the display before drawing or not (0).
2735 // The default values are : proj="ham" and clr=0.
2737 // This routine is based on the work by Garmt de Vries-Uiterweerd.
2739 AliTimestamp* tx=ts;
2740 if (fXsig && !tx) tx=fXsig->GetTimestamp();
2741 if (!tx) tx=(AliTimestamp*)this;
2743 // Display all stored reference signals
2746 for (Int_t i=1; i<=fRefs->GetSize(); i++)
2748 DisplaySignal(frame,mode,tx,i,proj,clr);
2749 clr=0; // No display clear for subsequent signals
2753 // Display the measurement
2754 DisplaySignal(frame,mode,tx,0,proj,clr);
2756 ///////////////////////////////////////////////////////////////////////////
2757 void AliAstrolab::SetCentralMeridian(Double_t phi,TString u)
2759 // Set the central meridian for the sky display.
2760 // Setting a value smaller than -pi will induce automatic meridian setting
2762 // By default the central meridian is set at -999 in the constructor.
2764 // The string argument "u" allows to choose between different angular units
2765 // u = "rad" : angle provided in radians
2766 // "deg" : angle provided in degrees
2767 // "dms" : angle provided in dddmmss.sss
2768 // "hms" : angle provided in hhmmss.sss
2770 // The default is u="deg".
2772 // This routine is based on the work by Garmt de Vries-Uiterweerd.
2774 fMeridian=ConvertAngle(phi,u,"rad");
2776 ///////////////////////////////////////////////////////////////////////////
2777 void AliAstrolab::Project(Double_t l,Double_t b,TString proj,Double_t& x,Double_t& y)
2779 // Generic interface for projection of a (long,lat) pair onto a (x,y) pair.
2781 // Meaning of the arguments :
2782 // --------------------------
2783 // l : Longitude (e.g. right ascension)
2784 // b : Latitude (e.g. declination)
2785 // proj : Projection mode (e.g. "ham")
2786 // x : The resulting x coordinate for the selected projection
2787 // y : The resulting x coordinate for the selected projection
2789 // The available projection modes are :
2791 // cyl : Cylindrical projection plotted with colored markers
2792 // cylh : Cylindrical projection plotted in a 2-D histogram
2793 // ham : Hammer equal area projection plotted with colored markers
2794 // hamh : Hammer equal area projection plotted in a 2-D histogram
2795 // ait : Aitoff projection plotted with colored markers
2796 // aith : Aitoff projection plotted in a 2-D histogram
2797 // mer : Mercator projection plotted with colored markers
2798 // merh : Mercator projection plotted in a 2-D histogram
2799 // ang : Straight cos(b) vs. l plot with colored markers
2800 // angh : Straight cos(b) vs. l plot in a 2-D histogram
2802 // Note : The ang(h) plot allows for easy identification of a homogeneous distribution.
2804 // This routine is based on the work by Garmt de Vries-Uiterweerd.
2806 Double_t pi=acos(-1.);
2808 // Subtract central meridian from longitude
2811 // Take l between -180 and 180 degrees
2824 // Convert (l,b) to (x,y) with -2 < x <= 2
2825 if (proj=="cyl" || proj=="cylh") ProjectCylindrical(l,b,x,y);
2826 if (proj=="ham" || proj=="hamh") ProjectHammer(l,b,x,y);
2827 if (proj=="ait" || proj=="aith") ProjectAitoff(l,b,x,y);
2828 if (proj=="mer" || proj=="merh") ProjectMercator(l,b,x,y);
2829 if (proj=="ang" || proj=="angh")
2835 ///////////////////////////////////////////////////////////////////////////
2836 void AliAstrolab::ProjectCylindrical(Double_t l,Double_t b,Double_t& x, Double_t& y)
2838 // Cylindrical projection of a (long,lat) coordinate pair
2840 // This routine is based on the work by Garmt de Vries-Uiterweerd.
2842 Double_t pi=acos(-1.);
2846 ///////////////////////////////////////////////////////////////////////////
2847 void AliAstrolab::ProjectHammer(Double_t l,Double_t b,Double_t& x,Double_t& y)
2849 // Hammer-Aitoff projection of a (long,lat) coordinate pair.
2850 // This is an equal-area projection.
2852 // This routine is based on the work by Garmt de Vries-Uiterweerd.
2854 Double_t k=1./sqrt(1.+cos(b)*cos(l/2.));
2855 x=2.*k*cos(b)*sin(l/2.);
2858 ///////////////////////////////////////////////////////////////////////////
2859 void AliAstrolab::ProjectAitoff(Double_t l,Double_t b,Double_t& x,Double_t& y)
2861 // Aitoff projection of a (long,lat) coordinate pair.
2863 // This routine is based on the work by Garmt de Vries-Uiterweerd.
2865 Double_t pi=acos(-1.);
2868 Double_t k=acos(cos(b)*cos(l/2.));
2871 x=4.*k*cos(b)*sin(l/2.)/(pi*sin(k));
2872 y=2.*k*sin(b)/(pi*sin(k));
2875 ///////////////////////////////////////////////////////////////////////////
2876 void AliAstrolab::ProjectMercator(Double_t l,Double_t b,Double_t& x,Double_t& y)
2878 // Mercator projection of a (long,lat) coordinate pair
2880 // This routine is based on the work by Garmt de Vries-Uiterweerd.
2882 Double_t pi=acos(-1.);
2885 if (b >= 89.*pi/180.) b=89.*pi/180.;
2886 if (b <= -89.*pi/180.) y=-89.*pi/180.;
2887 if (fabs(y) < 1) y=log((1.+sin(b))/(1.-sin(b)))/2.;
2889 ///////////////////////////////////////////////////////////////////////////
2890 TObject* AliAstrolab::Clone(const char* name) const
2892 // Make a deep copy of the current object and provide the pointer to the copy.
2893 // This memberfunction enables automatic creation of new objects of the
2894 // correct type depending on the object type, a feature which may be very useful
2895 // for containers when adding objects in case the container owns the objects.
2897 AliAstrolab* lab=new AliAstrolab(*this);
2900 if (strlen(name)) lab->SetName(name);
2904 ///////////////////////////////////////////////////////////////////////////