source: trunk/MagicSoft/Mars/mastro/MAstroCatalog.cc@ 3568

Last change on this file since 3568 was 3568, checked in by tbretz, 21 years ago
*** empty log message ***
File size: 26.4 KB
Line 
1/* ======================================================================== *\
2!
3! *
4! * This file is part of MARS, the MAGIC Analysis and Reconstruction
5! * Software. It is distributed to you in the hope that it can be a useful
6! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
7! * It is distributed WITHOUT ANY WARRANTY.
8! *
9! * Permission to use, copy, modify and distribute this software and its
10! * documentation for any purpose is hereby granted without fee,
11! * provided that the above copyright notice appear in all copies and
12! * that both that copyright notice and this permission notice appear
13! * in supporting documentation. It is provided "as is" without express
14! * or implied warranty.
15! *
16!
17!
18! Author(s): Thomas Bretz, 03/2004 <mailto:tbretz@astro.uni-wuerzburg.de>
19!
20! Copyright: MAGIC Software Development, 2002-2004
21!
22!
23\* ======================================================================== */
24
25//////////////////////////////////////////////////////////////////////////////
26//
27// MAstroCatalog
28//
29// THIS IMPLEMENTATION IS PRELIMINARY AND WILL BE MERGED WITH
30// SOME PARTS OF THE DRIVE SOFTWARE SOON!
31//
32//////////////////////////////////////////////////////////////////////////////
33#include "MAstroCatalog.h"
34
35#include <errno.h>
36#include <fstream>
37#include <stdlib.h>
38#include <limits.h> // INT_MAX (Suse 7.3/gcc 2.95)
39
40#include <KeySymbols.h>
41
42#include <TPad.h> // TPad::GetMaxPickDistance
43#include <TLine.h>
44#include <TMarker.h>
45#include <TCanvas.h>
46#include <TArrayI.h>
47#include <TGToolTip.h>
48#include <TRotation.h>
49#include <TPaveText.h>
50#include <TStopwatch.h>
51
52#include "MLog.h"
53#include "MLogManip.h"
54
55#include "MTime.h"
56#include "MAstro.h"
57#include "MAstroSky2Local.h"
58#include "MObservatory.h"
59
60ClassImp(MVector3);
61ClassImp(MAstroCatalog);
62
63using namespace std;
64/*
65class MRotation : public TRotation
66{
67public:
68 MRotation(Double_t gmst, const MObservatory &obs) : TRotation(1, 0, 0, 0, -1, 0, 0, 0, 1)
69 {
70 RotateZ(gmst + obs.GetElong());
71 RotateY(obs.GetPhi()-TMath::Pi()/2);
72 RotateZ(TMath::Pi());
73 }
74 MRotation(const MTime &t, const MObservatory &obs) : TRotation(1, 0, 0, 0, -1, 0, 0, 0, 1)
75 {
76 RotateZ(t.GetGmst() + obs.GetElong());
77 RotateY(obs.GetPhi()-TMath::Pi()/2);
78 RotateZ(TMath::Pi());
79 }
80};
81*/
82/*
83MVector3 MVector3::GetZdAz(const MObservatory &obs, Double_t gmst) const
84{
85 if (!fType==kIsRaDec)
86 return MVector3();
87
88 MVector3 v(*this);
89 v *= MAstroSky2Local(gmst, obs);
90
91 return v;
92
93 // ------ Using vectors -------
94 // v(1) = -v(1); // phi --> -phi
95 // v.RotateZ(gmst + obs.GetElong()); // -phi --> alpha-phi
96 // v.RotateY(obs.GetPhi()-TMath::Pi()/2);
97 // v.RotateZ(TMath::Pi());
98
99 // ------ The same using slalib, tested in the drive system -------
100 // const Double_t alpha = slaGmst(mjd) + obs.GetElong();
101 // Double_t el;
102 // slaDe2h(fAlpha-ra, dec, obs.GetPhi(), &az, &el);
103 // zd = TMath::Pi()/2-el;
104}
105
106MVector3 MVector3::GetZdAz(const MTime &time, MObservatory &obs) const
107{
108 return GetZdAz(obs, time.GetGmst());
109}
110
111MVector3 MVector3::GetRaDec(const MObservatory &obs, Double_t gmst) const
112{
113 if (!fType==kIsZdAz)
114 return MVector3();
115
116 MVector3 v(*this);
117 v *= MAstroSky2Local(gmst, obs).Inverse();
118
119 return v;
120
121 // ------ Using vectors -------
122 // v.RotateZ(-TMath::Pi());
123 // v.RotateY(TMath::Pi()/2-obs.GetPhi());
124 // v.RotateZ(-gmst - obs.GetElong()); // alpha-phi --> -phi
125 // v(1) = -v(1); // -phi --> phi
126
127 // ------ The same using slalib, tested in the drive system -------
128 // const Double_t alpha = slaGmst(mjd) + obs.GetElong();
129 // Double_t el;
130 // slaDe2h(fAlpha-ra, dec, obs.GetPhi(), &az, &el);
131 // zd = TMath::Pi()/2-el;
132}
133
134MVector3 MVector3::GetRaDec(const MTime &time, MObservatory &obs) const
135{
136 return GetRaDec(obs, time.GetGmst());
137}
138*/
139// --------------------------------------------------------------------------
140MAstroCatalog::MAstroCatalog() : fLimMag(99), fRadiusFOV(90), fToolTip(0), fObservatory(0), fTime(0)
141{
142 fList.SetOwner();
143 fToolTip = new TGToolTip(0, "", 0);
144}
145
146// --------------------------------------------------------------------------
147MAstroCatalog::~MAstroCatalog()
148{
149 if (fTime)
150 delete fTime;
151 if (fObservatory)
152 delete fObservatory;
153
154 fToolTip->Hide();
155 delete fToolTip;
156
157 DeleteMap();
158
159 // FIXME: There must be an easier way!
160 TIter Next(gROOT->GetListOfCanvases());
161 TCanvas *c;
162 while ((c=(TCanvas*)Next()))
163 c->Disconnect("ProcessedEvent(Int_t,Int_t,Int_t,TObject*)", this,
164 "EventInfo(Int_t,Int_t,Int_t,TObject*)");
165
166}
167
168// --------------------------------------------------------------------------
169TString MAstroCatalog::FindToken(TString &line, Char_t tok)
170{
171 Ssiz_t token = line.First(tok);
172 if (token<0)
173 {
174 const TString copy(line);
175 line = "";
176 return copy;
177 }
178
179 const TString res = line(0, token);
180 line.Remove(0, token+1);
181 return res;
182}
183
184// --------------------------------------------------------------------------
185Int_t MAstroCatalog::atoi(const TSubString &sub)
186{
187 return atoi(TString(sub));
188}
189
190// --------------------------------------------------------------------------
191Float_t MAstroCatalog::atof(const TSubString &sub)
192{
193 return atof(TString(sub));
194}
195
196// --------------------------------------------------------------------------
197Int_t MAstroCatalog::atoi(const TString &s)
198{
199 return std::atoi(s);
200}
201
202// --------------------------------------------------------------------------
203Float_t MAstroCatalog::atof(const TString &s)
204{
205 return std::atof(s);
206}
207
208// --------------------------------------------------------------------------
209Int_t MAstroCatalog::ReadXephem(TString catalog)
210{
211 SetBit(kHasChanged);
212
213 gLog << inf << "Reading Xephem catalog: " << catalog << endl;
214
215 ifstream fin(catalog);
216 if (!fin)
217 {
218 gLog << err << "Cannot open file " << catalog << ": ";
219 gLog << strerror(errno) << endl;
220 return 0;
221 }
222
223 Int_t add =0;
224 Int_t cnt =0;
225 Int_t line=0;
226
227 Double_t maxmag=0;
228
229 while (1)
230 {
231 TString row;
232 row.ReadLine(fin);
233 if (!fin)
234 break;
235
236 line++;
237
238 if (row[0]=='#')
239 continue;
240
241 TString line(row);
242
243 TString name = FindToken(line);
244 TString dummy = FindToken(line);
245 TString r = FindToken(line);
246 TString d = FindToken(line);
247 TString m = FindToken(line);
248 TString epoch = FindToken(line);
249
250 if (name.IsNull() || r.IsNull() || d.IsNull() || m.IsNull() || epoch.IsNull())
251 {
252 gLog << warn << "Invalid Entry Line #" << line << ": " << row << endl;
253 continue;
254 }
255
256 cnt++;
257
258 const Double_t mag = atof(m);
259
260 maxmag = TMath::Max(maxmag, mag);
261
262 if (mag>fLimMag)
263 continue;
264
265 if (epoch!="2000")
266 {
267 gLog << warn << "Epoch != 2000... skipped." << endl;
268 continue;
269 }
270
271 Double_t ra0, dec0;
272 MAstro::Coordinate2Angle(r, ra0);
273 MAstro::Coordinate2Angle(d, dec0);
274
275 ra0 *= TMath::Pi()/12;
276 dec0 *= TMath::Pi()/180;
277
278 MVector3 *star=new MVector3;
279 star->SetRaDec(ra0, dec0, mag);
280 star->SetName(name);
281 if (star->Angle(fRaDec)*TMath::RadToDeg()>fRadiusFOV)
282 {
283 delete star;
284 continue;
285 }
286
287 fList.Add(star);
288 add++;
289 }
290 gLog << inf << "Read " << add << " out of " << cnt << " (Total max mag=" << maxmag << ")" << endl;
291
292 return add;
293}
294
295// --------------------------------------------------------------------------
296Int_t MAstroCatalog::ReadNGC2000(TString catalog)
297{
298 SetBit(kHasChanged);
299
300 gLog << inf << "Reading NGC2000 catalog: " << catalog << endl;
301
302 ifstream fin(catalog);
303 if (!fin)
304 {
305 gLog << err << "Cannot open file " << catalog << ": ";
306 gLog << strerror(errno) << endl;
307 return 0;
308 }
309
310 Int_t add=0;
311 Int_t cnt=0;
312 Int_t n =0;
313
314 Double_t maxmag=0;
315
316 while (1)
317 {
318 TString row;
319 row.ReadLine(fin);
320 if (!fin)
321 break;
322
323 cnt++;
324
325 const Int_t rah = atoi(row(13, 2));
326 const Float_t ram = atof(row(16, 4));
327 const Char_t decs = row(22);
328 const Int_t decd = atoi(row(23, 2));
329 const Int_t decm = atoi(row(26, 2));
330 const TString m = row(43, 4);
331
332 if (m.Strip().IsNull())
333 continue;
334
335 n++;
336
337 const Double_t mag = atof(m);
338
339 maxmag = TMath::Max(maxmag, mag);
340
341 if (mag>fLimMag)
342 continue;
343
344 const Double_t ra = MAstro::Hms2Rad(rah, (int)ram, fmod(ram, 1)*60);
345 const Double_t dec = MAstro::Dms2Rad(decd, decm, 0, decs);
346
347 MVector3 *star=new MVector3;
348 star->SetRaDec(ra, dec, mag);
349 star->SetName(row(0, 8));
350 if (star->Angle(fRaDec)*TMath::RadToDeg()>fRadiusFOV)
351 {
352 delete star;
353 continue;
354 }
355
356 fList.Add(star);
357 add++;
358 }
359
360 gLog << inf << "Read " << add << " out of " << n << " (Total max mag=" << maxmag << ")" << endl;
361
362 return add;
363}
364
365// --------------------------------------------------------------------------
366Int_t MAstroCatalog::ReadBSC(TString catalog)
367{
368 SetBit(kHasChanged);
369
370 gLog << inf << "Reading Bright Star Catalog (BSC5) catalog: " << catalog << endl;
371
372 ifstream fin(catalog);
373 if (!fin)
374 {
375 gLog << err << "Cannot open file " << catalog << ": ";
376 gLog << strerror(errno) << endl;
377 return 0;
378 }
379
380 Int_t add=0;
381 Int_t cnt=0;
382 Int_t n =0;
383
384 Double_t maxmag=0;
385
386 while (1)
387 {
388 TString row;
389 row.ReadLine(fin);
390 if (!fin)
391 break;
392
393 cnt++;
394
395 const Int_t rah = atoi(row(75, 2));
396 const Int_t ram = atoi(row(77, 2));
397 const Float_t ras = atof(row(79, 4));
398 const Char_t decsgn = row(83);
399 const Int_t decd = atoi(row(84, 2));
400 const Int_t decm = atoi(row(86, 2));
401 const Int_t decs = atoi(row(88, 2));
402 const TString m = row(102, 5);
403
404 if (m.Strip().IsNull())
405 continue;
406
407 n++;
408
409 const Double_t mag = atof(m.Data());
410
411 maxmag = TMath::Max(maxmag, mag);
412
413 if (mag>fLimMag)
414 continue;
415
416 const Double_t ra = MAstro::Hms2Rad(rah, ram, ras);
417 const Double_t dec = MAstro::Dms2Rad(decd, decm, decs, decsgn);
418
419 MVector3 *star=new MVector3;
420 star->SetRaDec(ra, dec, mag);
421 star->SetName(row(4,9));
422 if (star->Angle(fRaDec)*TMath::RadToDeg()>fRadiusFOV)
423 {
424 delete star;
425 continue;
426 }
427
428 fList.Add(star);
429 add++;
430 }
431
432 gLog << inf << "Read " << add << " out of " << n << " (Total max mag=" << maxmag << ")" << endl;
433
434 return add;
435}
436
437// --------------------------------------------------------------------------
438void MAstroCatalog::Paint(Option_t *o)
439{
440 SetRangePad();
441
442 if (TestBit(kHasChanged))
443 DrawPrimitives(o);
444}
445
446// --------------------------------------------------------------------------
447void MAstroCatalog::DrawStar(Double_t x, Double_t y, const TVector3 &v, Bool_t transparent, const char *txt)
448{
449 const Double_t ra = v.Phi()*TMath::RadToDeg()/15;
450 const Double_t dec = (TMath::Pi()/2-v.Theta())*TMath::RadToDeg();
451
452 TString str = v.GetName();
453 str += Form(": Ra=%.2fh", ra);
454 str += Form(" Dec=%.1fd", dec);
455 str += Form(" Mag=%.1f", -2.5*log10(v.Mag()));
456 if (txt)
457 str += Form(" (%s)", txt);
458
459 // draw star on the camera display
460 TMarker *tip=new TMarker(x, y, transparent ? kDot : kFullDotLarge);;
461 tip->SetMarkerColor(kBlack);
462 AddMap(tip, new TString(str));
463}
464
465// --------------------------------------------------------------------------
466void MAstroCatalog::Update(Bool_t upd)
467{
468 if (gPad && TestBit(kMustCleanup))
469 {
470 SetBit(kHasChanged);
471 gPad->Modified();
472 if (upd)
473 gPad->Update();
474 }
475}
476
477// --------------------------------------------------------------------------
478//
479// Set the observation time. Necessary to use local coordinate
480// system.
481//
482void MAstroCatalog::SetTime(const MTime &time)
483{
484 if (fTime)
485 delete fTime;
486 fTime=(MTime*)time.Clone();
487}
488
489// --------------------------------------------------------------------------
490//
491// Set the observatory location. Necessary to use local coordinate
492// system.
493//
494void MAstroCatalog::SetObservatory(const MObservatory &obs)
495{
496 if (fObservatory)
497 delete fObservatory;
498 fObservatory=(MObservatory*)obs.Clone();
499}
500
501// --------------------------------------------------------------------------
502//
503// Convert the vector to pad coordinates. After conversion
504// the x- coordinate of the vector must be the x coordinate
505// of the pad - the same for y. If the coordinate is inside
506// the current draw area return kTRUE, otherwise kFALSE.
507// If it is an invalid coordinate return kERROR
508//
509Int_t MAstroCatalog::ConvertToPad(const TVector3 &w0, TVector2 &v) const
510{
511 TVector3 w(w0);
512
513 // Stretch such, that the Z-component is alwas the same. Now
514 // X and Y contains the intersection point between the star-light
515 // and the plain of a virtual plain screen (ccd...)
516 if (TestBit(kPlainScreen))
517 w *= 1./w(2);
518
519 w *= TMath::RadToDeg();
520 v.Set(w(0), w(1));
521
522 if (w(2)<0)
523 return kERROR;
524
525 return w(0)>gPad->GetX1() && w(1)>gPad->GetY1() &&
526 w(0)<gPad->GetX2() && w(1)<gPad->GetY2();
527}
528
529// --------------------------------------------------------------------------
530Int_t MAstroCatalog::Convert(const TRotation &rot, TVector2 &v) const
531{
532 MVector3 w;
533 w.SetMagThetaPhi(1, v.Y(), v.X());
534 w *= rot;
535
536 return ConvertToPad(w, v);
537}
538
539// --------------------------------------------------------------------------
540Bool_t MAstroCatalog::DrawLine(const TVector2 &v, Double_t dx, Double_t dy, const TRotation &rot, Int_t type)
541{
542 const TVector2 add(dx*TMath::DegToRad(), dy*TMath::DegToRad());
543
544 TVector2 v0 = v;
545 TVector2 v1 = v+add;
546
547 const Int_t rc0 = Convert(rot, v0);
548 const Int_t rc1 = Convert(rot, v1);
549
550 // Both are kFALSE or both are kERROR
551 if ((rc0|rc1)==kFALSE || (rc0&rc1)==kERROR)
552 return kFALSE;
553
554 TLine *line = new TLine(v0.X(), v0.Y(), v1.X(), v1.Y());
555 line->SetLineStyle(kDashDotted); //kDashed, kDotted, kDashDotted
556 line->SetLineColor(kWhite+type*2);
557 AddMap(line);
558
559 const TVector2 deg = v*TMath::RadToDeg();
560 TString txt = type==1 ?
561 Form("Ra=%.2fh Dec=%.1fd", fmod(deg.X()/15+48, 24), fmod(90-deg.Y()+270,180)-90) :
562 Form("Zd=%.1fd Az=%.1fd", fmod(deg.Y()+270,180)-90, fmod(deg.X()+720, 360));
563
564 TMarker *tip=new TMarker(v0.X(), v0.Y(), kDot);
565 tip->SetMarkerColor(kWhite+type*2);
566 AddMap(tip, new TString(txt));
567
568 return kTRUE;
569}
570
571// --------------------------------------------------------------------------
572//
573// Use "local" draw option to align the display to the local
574// coordinate system instead of the sky coordinate system.
575//
576void MAstroCatalog::Draw(const TVector2 &v0, const TRotation &rot, TArrayI &dx, TArrayI &dy, Int_t stepx, Int_t stepy, Int_t type)
577{
578 const TVector2 v1 = v0 + TVector2(dx[0]*TMath::DegToRad(), dy[0]*TMath::DegToRad());
579
580 Int_t idx[] = {1, 1, 1, 1};
581
582 Int_t dirs[4][2] = { {0, stepy}, {stepx, 0}, {0, -stepy}, {-stepx, 0} };
583
584 for (int i=0; i<dx.GetSize(); i++)
585 {
586 for (int j=0; j<4; j++)
587 {
588 const Bool_t rcx0 = (dx[i]+720)%360==(dx[0]+dirs[j][0]+720)%360;
589 const Bool_t rcy0 = (dy[i]+360)%180==(dy[0]+dirs[j][1]+360)%180;
590 if (rcx0&&rcy0)
591 idx[j] = 0;
592 }
593 }
594
595 // Enhance size of array by 1, copy current
596 // position as last entry
597 dx.Set(dx.GetSize()+1);
598 dy.Set(dy.GetSize()+1);
599
600 dx[dx.GetSize()-1] = dx[0];
601 dy[dy.GetSize()-1] = dy[0];
602
603 // Store current positon
604 const Int_t d[2] = { dx[0], dy[0] };
605
606 for (int i=0; i<4; i++)
607 if (idx[i])
608 {
609 // Calculate new position
610 dx[0] = d[0]+dirs[i][0];
611 dy[0] = d[1]+dirs[i][1];
612
613 //cout << dx[0] << " " << dy[0] << endl;
614
615 // Draw corresponding line and iterate through grid
616 if (DrawLine(v1, dirs[i][0], dirs[i][1], rot, type))
617 Draw(v0, rot, dx, dy, stepx, stepy, type);
618
619 dx[0]=d[0];
620 dy[0]=d[1];
621 }
622}
623
624// --------------------------------------------------------------------------
625void MAstroCatalog::DrawGrid(const TVector3 &v0, const TRotation &rot, Int_t type)
626{
627 TArrayI dx(1);
628 TArrayI dy(1);
629
630 // align to 1deg boundary
631 TVector2 v(v0.Phi()*TMath::RadToDeg(), v0.Theta()*TMath::RadToDeg());
632 v.Set((Float_t)TMath::Nint(v.X()), (Float_t)TMath::Nint(v.Y()));
633
634 // calculate stepsizes based on visible FOV
635 Int_t stepx = 1;
636
637 if (v.Y()<fRadiusFOV || v.Y()>180-fRadiusFOV)
638 stepx=36;
639 else
640 {
641 // This is a rough estimate how many degrees are visible
642 const Float_t m = log(fRadiusFOV/180.)/log(90./(fRadiusFOV+1)+1);
643 const Float_t t = log(180.)-m*log(fRadiusFOV);
644 const Float_t f = m*log(90-fabs(90-v.Y()))+t;
645 const Int_t nx = (Int_t)(exp(f)+0.5);
646 stepx = nx<4 ? 1 : nx/4;
647 if (stepx>36)
648 stepx=36;
649 }
650
651 const Int_t ny = (Int_t)(fRadiusFOV+1);
652 Int_t stepy = ny<4 ? 1 : ny/4;
653
654 // align stepsizes to be devisor or 180 and 90
655 while (180%stepx)
656 stepx++;
657 while (90%stepy)
658 stepy++;
659
660 // align to step-size boundary (search for the nearest one)
661 Int_t dv = 1;
662 while ((int)(v.X())%stepx)
663 {
664 v.Set(v.X()+dv, v.Y());
665 dv = -TMath::Sign(TMath::Abs(dv)+1, dv);
666 }
667
668 dv = 1;
669 while ((int)(v.Y())%stepy)
670 {
671 v.Set(v.X(), v.Y()+dv);
672 dv = -TMath::Sign(TMath::Abs(dv)+1, dv);
673 }
674
675 // draw...
676 v *= TMath::DegToRad();
677
678 Draw(v, rot, dx, dy, stepx, stepy, type);
679}
680
681// --------------------------------------------------------------------------
682//
683// Get a rotation matrix which aligns the pointing position
684// to the center of the x,y plain
685//
686TRotation MAstroCatalog::AlignCoordinates(const TVector3 &v) const
687{
688 TRotation trans;
689 trans.RotateZ(-v.Phi());
690 trans.RotateY(-v.Theta());
691 trans.RotateZ(-TMath::Pi()/2);
692 return trans;
693}
694
695// --------------------------------------------------------------------------
696//
697// Return the rotation matrix which converts either sky or
698// local coordinates to coordinates which pole is the current
699// pointing direction.
700//
701TRotation MAstroCatalog::GetGrid(Bool_t local)
702{
703 const Bool_t enable = fTime && fObservatory;
704
705 if (!local)
706 {
707 const TRotation trans(AlignCoordinates(fRaDec));
708
709 DrawGrid(fRaDec, trans, 1);
710
711 if (enable)
712 {
713 const MAstroSky2Local rot(*fTime, *fObservatory);
714 DrawGrid(rot*fRaDec, trans*rot.Inverse(), 2);
715 }
716
717 return trans;
718 }
719
720 if (local && enable)
721 {
722 const MAstroSky2Local rot(*fTime, *fObservatory);
723
724 const TRotation trans(AlignCoordinates(rot*fRaDec));
725
726 DrawGrid(fRaDec, trans*rot, 1);
727 DrawGrid(rot*fRaDec, trans, 2);
728
729 return trans*rot;
730 }
731
732 return TRotation();
733}
734
735// --------------------------------------------------------------------------
736TString MAstroCatalog::GetPadTitle() const
737{
738 const Double_t ra = fRaDec.Phi()*TMath::RadToDeg();
739 const Double_t dec = (TMath::Pi()/2-fRaDec.Theta())*TMath::RadToDeg();
740
741 TString txt;
742 txt += Form("\\alpha=%.2fh ", fmod(ra/15+48, 24));
743 txt += Form("\\delta=%.1f\\circ ", fmod(dec+270,180)-90);
744 txt += Form("/ FOV=%.1f\\circ", fRadiusFOV);
745
746 if (!fTime || !fObservatory)
747 return txt;
748
749 const MAstroSky2Local rot(*fTime, *fObservatory);
750 const TVector3 loc = rot*fRaDec;
751
752 const Double_t rho = rot.RotationAngle(fRaDec.Phi(), TMath::Pi()/2-fRaDec.Theta());
753
754 const Double_t zd = TMath::RadToDeg()*loc.Theta();
755 const Double_t az = TMath::RadToDeg()*loc.Phi();
756
757 txt.Prepend("#splitline{");
758
759 txt += Form(" \\theta=%.1fh ", fmod(zd+270,180)-90);
760 txt += Form("\\phi=%.1f\\circ ", fmod(az+720, 360));
761 txt += Form(" / \\rho=%.1f\\circ", rho*TMath::RadToDeg());
762 txt += "}{<";
763 txt += fTime->GetSqlDateTime();
764 txt += ">}";
765 return txt;
766}
767
768// --------------------------------------------------------------------------
769void MAstroCatalog::AddPrimitives(Option_t *o)
770{
771 const Bool_t local = TString(o).Contains("local", TString::kIgnoreCase);
772
773 cout << "Opt: " << o << endl;
774
775 const TRotation rot(GetGrid(local));
776
777 TIter Next(&fList);
778 TVector3 *v=0;
779 while ((v=(TVector3*)Next()))
780 {
781 // FIXME: Check Magnitude!
782 TVector2 s(v->Phi(), v->Theta());
783 if (Convert(rot, s)==kTRUE)
784 DrawStar(s.X(), TMath::Pi()/2-s.Y(), *v, kFALSE);
785 }
786
787 TPaveText *pv = new TPaveText(0.01, 0.90, 0.63, 0.99, "brNDC");
788 pv->AddText(GetPadTitle());
789 AddMap(pv);
790
791 TMarker *mk=new TMarker(0, 0, kMultiply);
792 mk->SetMarkerColor(kBlack);
793 mk->SetMarkerSize(1.5);
794 AddMap(mk);
795}
796
797// --------------------------------------------------------------------------
798void MAstroCatalog::SetRangePad()
799{
800 const Double_t edge = fRadiusFOV/TMath::Sqrt(2.);
801 gPad->Range(-edge, -edge, edge, edge);
802
803 const Float_t w = gPad->GetWw();
804 const Float_t h = gPad->GetWh();
805
806 if (w<h)
807 gPad->Range(-edge, -edge*h/w, edge, edge*h/w);
808 else
809 gPad->Range(-edge*w/h, -edge, edge*w/h, edge);
810}
811
812// --------------------------------------------------------------------------
813void MAstroCatalog::DrawPrimitives(Option_t *o)
814{
815 DeleteMap();
816
817 SetRangePad();
818
819 TStopwatch clk;
820 clk.Start();
821 AddPrimitives(o);
822 clk.Stop();
823 clk.Print();
824
825 // Append to a possible second pad
826 if (!gPad->GetListOfPrimitives()->FindObject(this))
827 AppendPad(o);
828
829 // Append all objects to pad
830 DrawMap();
831
832 ResetBit(kHasChanged);
833}
834
835// --------------------------------------------------------------------------
836void MAstroCatalog::DrawMap()
837{
838 Long_t key, val;
839 TExMapIter map(&fMapG);
840 while (map.Next(key, val))
841 ((TObject*)key)->Draw();
842}
843
844// --------------------------------------------------------------------------
845void MAstroCatalog::Draw(Option_t *o)
846{
847 // Append to first pad
848 AppendPad(o);
849
850 // If contents have not previously changed make sure that
851 // all primitives are recreated.
852 if (!TestBit(kHasChanged))
853 DrawPrimitives(o);
854
855 // Connect all TCanvas::ProcessedEvent to this->EventInfo
856 // This means, that after TCanvas has processed an event
857 // EventInfo of this class is called, see TCanvas::HandleInput
858 gPad->GetCanvas()->Connect("ProcessedEvent(Int_t,Int_t,Int_t,TObject*)",
859 "MAstroCatalog", this,
860 "EventInfo(Int_t,Int_t,Int_t,TObject*)");
861
862 // Do this instead of fListG.Draw, because
863 // TCollection overwrites Draw
864 // Would be nice, but doesn't work because the single
865 // graphical object are not handled by TPad anymore...
866 // fListG.AppendPad();
867}
868
869// --------------------------------------------------------------------------
870//
871// This function was connected to all created canvases. It is used
872// to redirect GetObjectInfo into our own status bar.
873//
874// The 'connection' is done in AddTab
875//
876void MAstroCatalog::EventInfo(Int_t event, Int_t px, Int_t py, TObject *selected)
877{
878 TCanvas *c = (TCanvas*)gTQSender;
879
880 gPad = c ? c->GetSelectedPad() : NULL;
881 if (!gPad)
882 return;
883
884 // Try to find a corresponding object with kCannotPick set and
885 // an available TString (for a tool tip)
886 TString *str=0;
887 if (!selected || selected==this)
888 {
889 Long_t key, val;
890 TExMapIter map(&fMapG);
891 while (map.Next(key, val))
892 {
893 if (!val)
894 continue;
895
896 TObject *o=(TObject*)key;
897 if (o->DistancetoPrimitive(px, py)>TPad::GetMaxPickDistance())
898 continue;
899
900 selected = o;
901 str = (TString*)val;
902 break;
903 }
904 }
905
906 if (!selected)
907 return;
908
909 switch (event)
910 {
911 case kMouseMotion:
912 if (!fToolTip->IsMapped() && str)
913 ShowToolTip(px, py, *str);
914 break;
915
916 case kMouseLeave:
917 if (fToolTip->IsMapped())
918 fToolTip->Hide();
919 break;
920
921 case kKeyPress:
922 ExecuteEvent(kKeyPress, px, py);
923 break;
924 }
925}
926
927// --------------------------------------------------------------------------
928void MAstroCatalog::ExecuteEventKbd(Int_t keycode, Int_t keysym)
929{
930 Double_t dra =0;
931 Double_t ddec=0;
932
933 switch (keysym)
934 {
935 case kKey_Left:
936 dra = -TMath::DegToRad();
937 break;
938 case kKey_Right:
939 dra = +TMath::DegToRad();
940 break;
941 case kKey_Up:
942 ddec = +TMath::DegToRad();
943 break;
944 case kKey_Down:
945 ddec = -TMath::DegToRad();
946 break;
947 case kKey_Plus:
948 SetRadiusFOV(fRadiusFOV+1);
949 break;
950 case kKey_Minus:
951 SetRadiusFOV(fRadiusFOV-1);
952 break;
953
954 default:
955 return;
956 }
957
958 const Double_t r = fRaDec.Phi();
959 const Double_t d = TMath::Pi()/2-fRaDec.Theta();
960
961 SetRaDec(r+dra, d+ddec);
962
963 gPad->Update();
964}
965
966// ------------------------------------------------------------------------
967//
968// Execute a gui event on the camera
969//
970void MAstroCatalog::ExecuteEvent(Int_t event, Int_t mp1, Int_t mp2)
971{
972 if (!TestBit(kGuiActive))
973 return;
974
975 if (event==kKeyPress)
976 ExecuteEventKbd(mp1, mp2);
977}
978
979// --------------------------------------------------------------------------
980Int_t MAstroCatalog::DistancetoPrimitive(Int_t px, Int_t py)
981{
982 Int_t min = INT_MAX;
983
984 Long_t key, val;
985 TExMapIter map(&fMapG);
986 while (map.Next(key, val))
987 {
988 TObject *o=(TObject*)key;
989
990 const Int_t d = o->DistancetoPrimitive(px, py);
991
992 if (d<TPad::GetMaxPickDistance())
993 return 0;
994
995 if (d<min)
996 min=d;
997 }
998
999 return min;
1000}
1001
1002// --------------------------------------------------------------------------
1003void MAstroCatalog::ShowToolTip(Int_t px, Int_t py, const char *txt)
1004{
1005 Int_t x=0;
1006 Int_t y=0;
1007
1008 const Window_t id1 = gVirtualX->GetWindowID(gPad->GetCanvasID());
1009 const Window_t id2 = fToolTip->GetParent()->GetId();
1010
1011 Window_t id3;
1012 gVirtualX->TranslateCoordinates(id1, id2, px, py, x, y, id3);
1013
1014 // Show tool tip
1015 fToolTip->SetText(txt);
1016 fToolTip->Show(x+4, y+4);
1017}
Note: See TracBrowser for help on using the repository browser.