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

Last change on this file since 3709 was 3707, checked in by tbretz, 21 years ago
*** empty log message ***
File size: 31.7 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// To display a starfield you must have a supported catalog, then do:
33//
34// MTime time;
35// // Time for which to get the picture
36// time.Set(2004, 2, 28, 20, 14, 7);
37// // Current observatory
38// MObservatory magic1;
39// // Right Ascension [h] and declination [deg] of source
40// // Currently 'perfect' pointing is assumed
41// const Double_t ra = MAstro::Hms2Rad(5, 34, 31.9);
42// const Double_t dec = MAstro::Dms2Rad(22, 0, 52.0);
43// MAstroCatalog stars;
44// // Magnitude up to which the stars are loaded from the catalog
45// stars.SetLimMag(6);
46// // Radius of FOV around the source position to load the stars
47// stars.SetRadiusFOV(3);
48// // Source position
49// stars.SetRaDec(ra, dec);
50// // Catalog to load (here: Bright Star Catalog V5)
51// stars.ReadBSC("bsc5.dat");
52// // Obersavatory and time to also get local coordinate information
53// stars.SetObservatory(magic1);
54// stars.SetTime(time);
55// // Enable interactive GUI
56// stars.SetGuiActive();
57// //Clone the catalog due to the validity range of the instance
58// TObject *o = stars.Clone();
59// o->SetBit(kCanDelete);
60// o->Draw();
61//
62// If no time and/or Obervatory location is given no local coordinate
63// information is displayed.
64//
65// The conversion from sky coordinates to local coordinates is done using
66// MAstroSky2Local which does a simple rotation of the coordinate system.
67// This is inaccurate in the order of 30arcsec due to ignorance of all
68// astrometrical corrections (nutation, precission, abberation, ...)
69//
70// GUI:
71// * If the gui is interactive you can use the cursor keys to change
72// the position you are looking at and with plus/minus you
73// can (un)zoom the FOV (Field Of View)
74//
75// * The displayed values mean the following:
76// + alpha: Right Ascension
77// + delta: Declination
78// + theta: zenith distance / zenith angle
79// + phi: azimuth angle
80// + rho: angle of rotation sky-coordinate system vs local-
81// coordinate system
82// + time of display
83//
84// * Move the mouse on top of the grid points or the stars to get
85// more setailed information.
86//
87// ToDo:
88// - replace MVetcor3 by a more convinient class. Maybe use TExMap, too.
89//
90//////////////////////////////////////////////////////////////////////////////
91#include "MAstroCatalog.h"
92
93#include <errno.h>
94#include <fstream>
95#include <stdlib.h>
96#include <limits.h> // INT_MAX (Suse 7.3/gcc 2.95)
97
98#include <KeySymbols.h>
99
100#include <TPad.h> // TPad::GetMaxPickDistance
101#include <TLine.h>
102#include <TMarker.h>
103#include <TCanvas.h>
104#include <TArrayI.h>
105#include <TGToolTip.h>
106#include <TRotation.h>
107#include <TPaveText.h>
108
109#include "MLog.h"
110#include "MLogManip.h"
111
112#include "MTime.h"
113#include "MAstro.h"
114#include "MAstroSky2Local.h"
115#include "MObservatory.h"
116
117#undef DEBUG
118//#define DEBUG
119
120#ifdef DEBUG
121#include <TStopwatch.h>
122#endif
123
124ClassImp(MVector3);
125ClassImp(MAstroCatalog);
126
127using namespace std;
128
129// --------------------------------------------------------------------------
130//
131// Default Constructor. Set Default values:
132// fLimMag = 99
133// fRadiusFOV = 90
134//
135MAstroCatalog::MAstroCatalog() : fLimMag(99), fRadiusFOV(90), fToolTip(0), fObservatory(0), fTime(0)
136{
137 fList.SetOwner();
138 fToolTip = new TGToolTip(0, "", 0);
139}
140
141// --------------------------------------------------------------------------
142//
143// Destructor. Delete fTime, fObservatory. Disconnect signal. delete tooltip.
144// Delete Map with gui primitives
145//
146MAstroCatalog::~MAstroCatalog()
147{
148 // First disconnect the EventInfo...
149 // FIXME: There must be an easier way!
150 TIter Next(gROOT->GetListOfCanvases());
151 TCanvas *c;
152 while ((c=(TCanvas*)Next()))
153 c->Disconnect("ProcessedEvent(Int_t,Int_t,Int_t,TObject*)", this,
154 "EventInfo(Int_t,Int_t,Int_t,TObject*)");
155
156 // Now delete the data members
157 if (fTime)
158 delete fTime;
159 if (fObservatory)
160 delete fObservatory;
161
162 fToolTip->Hide();
163 delete fToolTip;
164
165 DeleteMap();
166}
167
168// --------------------------------------------------------------------------
169//
170// Snippet to for reading ctalog files.
171//
172TString MAstroCatalog::FindToken(TString &line, Char_t tok)
173{
174 Ssiz_t token = line.First(tok);
175 if (token<0)
176 {
177 const TString copy(line);
178 line = "";
179 return copy;
180 }
181
182 const TString res = line(0, token);
183 line.Remove(0, token+1);
184 return res;
185}
186
187// --------------------------------------------------------------------------
188//
189// return int correspoding to TSubString
190//
191Int_t MAstroCatalog::atoi(const TSubString &sub)
192{
193 return atoi(TString(sub));
194}
195
196// --------------------------------------------------------------------------
197//
198// return float correspoding to TSubString
199//
200Float_t MAstroCatalog::atof(const TSubString &sub)
201{
202 return atof(TString(sub));
203}
204
205// --------------------------------------------------------------------------
206//
207// return int correspoding to TString
208//
209Int_t MAstroCatalog::atoi(const TString &s)
210{
211 return std::atoi(s);
212}
213
214// --------------------------------------------------------------------------
215//
216// return float correspoding to TString
217//
218Float_t MAstroCatalog::atof(const TString &s)
219{
220 return std::atof(s);
221}
222
223// --------------------------------------------------------------------------
224//
225// Read from a xephem catalog, set bit kHasChahanged.
226// Already read data is not deleted. To delete the stored data call
227// Delete().
228//
229Int_t MAstroCatalog::ReadXephem(TString catalog)
230{
231 SetBit(kHasChanged);
232
233 gLog << inf << "Reading Xephem catalog: " << catalog << endl;
234
235 ifstream fin(catalog);
236 if (!fin)
237 {
238 gLog << err << "Cannot open file " << catalog << ": ";
239 gLog << strerror(errno) << endl;
240 return 0;
241 }
242
243 Int_t add =0;
244 Int_t cnt =0;
245 Int_t line=0;
246
247 Double_t maxmag=0;
248
249 while (1)
250 {
251 TString row;
252 row.ReadLine(fin);
253 if (!fin)
254 break;
255
256 line++;
257
258 if (row[0]=='#')
259 continue;
260
261 TString line(row);
262
263 TString name = FindToken(line);
264 TString dummy = FindToken(line);
265 TString r = FindToken(line);
266 TString d = FindToken(line);
267 TString m = FindToken(line);
268 TString epoch = FindToken(line);
269
270 if (name.IsNull() || r.IsNull() || d.IsNull() || m.IsNull() || epoch.IsNull())
271 {
272 gLog << warn << "Invalid Entry Line #" << line << ": " << row << endl;
273 continue;
274 }
275
276 cnt++;
277
278 const Double_t mag = atof(m);
279
280 maxmag = TMath::Max(maxmag, mag);
281
282 if (mag>fLimMag)
283 continue;
284
285 if (epoch!="2000")
286 {
287 gLog << warn << "Epoch != 2000... skipped." << endl;
288 continue;
289 }
290
291 Double_t ra0, dec0;
292 MAstro::Coordinate2Angle(r, ra0);
293 MAstro::Coordinate2Angle(d, dec0);
294
295 ra0 *= TMath::Pi()/12;
296 dec0 *= TMath::Pi()/180;
297
298 MVector3 *star=new MVector3;
299 star->SetRaDec(ra0, dec0, mag);
300 star->SetName(name);
301 if (star->Angle(fRaDec)*TMath::RadToDeg()>fRadiusFOV)
302 {
303 delete star;
304 continue;
305 }
306
307 fList.Add(star);
308 add++;
309 }
310 gLog << inf << "Read " << add << " out of " << cnt << " (Total max mag=" << maxmag << ")" << endl;
311
312 return add;
313}
314
315// --------------------------------------------------------------------------
316//
317// Read from a NGC2000 catalog. set bit kHasChanged
318// Already read data is not deleted. To delete the stored data call
319// Delete().
320//
321Int_t MAstroCatalog::ReadNGC2000(TString catalog)
322{
323 SetBit(kHasChanged);
324
325 gLog << inf << "Reading NGC2000 catalog: " << catalog << endl;
326
327 ifstream fin(catalog);
328 if (!fin)
329 {
330 gLog << err << "Cannot open file " << catalog << ": ";
331 gLog << strerror(errno) << endl;
332 return 0;
333 }
334
335 Int_t add=0;
336 Int_t cnt=0;
337 Int_t n =0;
338
339 Double_t maxmag=0;
340
341 while (1)
342 {
343 TString row;
344 row.ReadLine(fin);
345 if (!fin)
346 break;
347
348 cnt++;
349
350 const Int_t rah = atoi(row(13, 2));
351 const Float_t ram = atof(row(16, 4));
352 const Char_t decs = row(22);
353 const Int_t decd = atoi(row(23, 2));
354 const Int_t decm = atoi(row(26, 2));
355 const TString m = row(43, 4);
356
357 if (m.Strip().IsNull())
358 continue;
359
360 n++;
361
362 const Double_t mag = atof(m);
363
364 maxmag = TMath::Max(maxmag, mag);
365
366 if (mag>fLimMag)
367 continue;
368
369 const Double_t ra = MAstro::Hms2Rad(rah, (int)ram, fmod(ram, 1)*60);
370 const Double_t dec = MAstro::Dms2Rad(decd, decm, 0, decs);
371
372 MVector3 *star=new MVector3;
373 star->SetRaDec(ra, dec, mag);
374 star->SetName(row(0, 8));
375 if (star->Angle(fRaDec)*TMath::RadToDeg()>fRadiusFOV)
376 {
377 delete star;
378 continue;
379 }
380
381 fList.Add(star);
382 add++;
383 }
384
385 gLog << inf << "Read " << add << " out of " << n << " (Total max mag=" << maxmag << ")" << endl;
386
387 return add;
388}
389
390// --------------------------------------------------------------------------
391//
392// Read from a Bright Star Catalog catalog. set bit kHasChanged
393// Already read data is not deleted. To delete the stored data call
394// Delete().
395//
396Int_t MAstroCatalog::ReadBSC(TString catalog)
397{
398 SetBit(kHasChanged);
399
400 gLog << inf << "Reading Bright Star Catalog (BSC5) catalog: " << catalog << endl;
401
402 ifstream fin(catalog);
403 if (!fin)
404 {
405 gLog << err << "Cannot open file " << catalog << ": ";
406 gLog << strerror(errno) << endl;
407 return 0;
408 }
409
410 Int_t add=0;
411 Int_t cnt=0;
412 Int_t n =0;
413
414 Double_t maxmag=0;
415
416 while (1)
417 {
418 TString row;
419 row.ReadLine(fin);
420 if (!fin)
421 break;
422
423 cnt++;
424
425 const Int_t rah = atoi(row(75, 2));
426 const Int_t ram = atoi(row(77, 2));
427 const Float_t ras = atof(row(79, 4));
428 const Char_t decsgn = row(83);
429 const Int_t decd = atoi(row(84, 2));
430 const Int_t decm = atoi(row(86, 2));
431 const Int_t decs = atoi(row(88, 2));
432 const TString m = row(102, 5);
433
434 if (m.Strip().IsNull())
435 continue;
436
437 n++;
438
439 const Double_t mag = atof(m.Data());
440
441 maxmag = TMath::Max(maxmag, mag);
442
443 if (mag>fLimMag)
444 continue;
445
446 const Double_t ra = MAstro::Hms2Rad(rah, ram, ras);
447 const Double_t dec = MAstro::Dms2Rad(decd, decm, decs, decsgn);
448
449 MVector3 *star=new MVector3;
450 star->SetRaDec(ra, dec, mag);
451 star->SetName(row(4,9));
452 if (star->Angle(fRaDec)*TMath::RadToDeg()>fRadiusFOV)
453 {
454 delete star;
455 continue;
456 }
457
458 fList.Add(star);
459 add++;
460 }
461
462 gLog << inf << "Read " << add << " out of " << n << " (Total max mag=" << maxmag << ")" << endl;
463
464 return add;
465}
466
467// --------------------------------------------------------------------------
468//
469// Set Range of pad. If something has changed create and draw new primitives.
470// Paint all gui primitives.
471//
472void MAstroCatalog::Paint(Option_t *o)
473{
474 SetRangePad(o);
475
476 if (TestBit(kHasChanged))
477 DrawPrimitives(o);
478
479 PaintMap();
480}
481
482// --------------------------------------------------------------------------
483//
484// Draw a black marker at the position of the star. Create a corresponding
485// tooltip with the coordinates.
486// x, y: Pad Coordinates to draw star
487// v: Sky position (Ra/Dec) of the star
488// transparent: Draw marker or tooltip only
489// txt: additional tooltip text
490//
491void MAstroCatalog::DrawStar(Double_t x, Double_t y, const TVector3 &v, Bool_t transparent, const char *txt)
492{
493 const Double_t ra = v.Phi()*TMath::RadToDeg()/15;
494 const Double_t dec = (TMath::Pi()/2-v.Theta())*TMath::RadToDeg();
495
496 TString str = v.GetName();
497 str += Form(": Ra=%.2fh", ra);
498 str += Form(" Dec=%.1fd", dec);
499 str += Form(" Mag=%.1f", -2.5*log10(v.Mag()));
500 if (txt)
501 str += Form(" (%s)", txt);
502
503 // draw star on the camera display
504 TMarker *tip=new TMarker(x, y, transparent ? kDot : kFullDotMedium);;
505 tip->SetMarkerColor(kBlack);
506 AddMap(tip, new TString(str));
507}
508
509// --------------------------------------------------------------------------
510//
511// Set pad as modified.
512//
513void MAstroCatalog::Update(Bool_t upd)
514{
515 if (gPad && TestBit(kMustCleanup))
516 {
517 SetBit(kHasChanged);
518 gPad->Modified();
519 if (upd)
520 gPad->Update();
521 }
522}
523
524// --------------------------------------------------------------------------
525//
526// Set the observation time. Necessary to use local coordinate
527// system. The MTime object is cloned.
528//
529void MAstroCatalog::SetTime(const MTime &time)
530{
531 if (fTime)
532 delete fTime;
533 fTime=(MTime*)time.Clone();
534}
535
536// --------------------------------------------------------------------------
537//
538// Set the observatory location. Necessary to use local coordinate
539// system. The MObservatory object is cloned.
540//
541void MAstroCatalog::SetObservatory(const MObservatory &obs)
542{
543 if (fObservatory)
544 delete fObservatory;
545 fObservatory=(MObservatory*)obs.Clone();
546}
547
548// --------------------------------------------------------------------------
549//
550// Convert the vector to pad coordinates. After conversion
551// the x- coordinate of the vector must be the x coordinate
552// of the pad - the same for y. If the coordinate is inside
553// the current draw area return kTRUE, otherwise kFALSE.
554// If it is an invalid coordinate return kERROR
555//
556Int_t MAstroCatalog::ConvertToPad(const TVector3 &w0, TVector2 &v) const
557{
558 TVector3 w(w0);
559
560 // Stretch such, that the Z-component is alwas the same. Now
561 // X and Y contains the intersection point between the star-light
562 // and the plain of a virtual plain screen (ccd...)
563 if (TestBit(kPlainScreen))
564 w *= 1./w(2);
565
566 w *= TMath::RadToDeg(); // FIXME: *conversion factor?
567 v.Set(TestBit(kMirrorX) ? -w(0) : w(0),
568 TestBit(kMirrorY) ? -w(1) : w(1));
569
570 if (w(2)<0)
571 return kERROR;
572
573 return v.X()>gPad->GetX1() && v.Y()>gPad->GetY1() &&
574 v.X()<gPad->GetX2() && v.Y()<gPad->GetY2();
575}
576
577// --------------------------------------------------------------------------
578//
579// Convert theta/phi coordinates of v by TRotation into new coordinate
580// system and convert the coordinated to pad by ConvertToPad.
581// The result is retunred in v.
582//
583Int_t MAstroCatalog::Convert(const TRotation &rot, TVector2 &v) const
584{
585 MVector3 w;
586 w.SetMagThetaPhi(1, v.Y(), v.X());
587 w *= rot;
588
589 return ConvertToPad(w, v);
590}
591
592// --------------------------------------------------------------------------
593//
594// Draw a line from v to v+(dx,dy) using Convert/ConvertToPad to get the
595// corresponding pad coordinates.
596//
597Bool_t MAstroCatalog::DrawLine(const TVector2 &v, Double_t dx, Double_t dy, const TRotation &rot, Int_t type)
598{
599 const TVector2 add(dx*TMath::DegToRad(), dy*TMath::DegToRad());
600
601 TVector2 v0 = v;
602 TVector2 v1 = v+add;
603
604 const Int_t rc0 = Convert(rot, v0);
605 const Int_t rc1 = Convert(rot, v1);
606
607 // Both are kFALSE or both are kERROR
608 if ((rc0|rc1)==kFALSE || (rc0&rc1)==kERROR)
609 return kFALSE;
610
611 TLine *line = new TLine(v0.X(), v0.Y(), v1.X(), v1.Y());
612 line->SetLineStyle(kDashDotted); //kDashed, kDotted, kDashDotted
613 line->SetLineColor(kWhite+type*2);
614 AddMap(line);
615
616 const TVector2 deg = v*TMath::RadToDeg();
617 TString txt = type==1 ?
618 Form("Ra=%.2fh Dec=%.1fd", fmod(deg.X()/15+48, 24), fmod(90-deg.Y()+270,180)-90) :
619 Form("Zd=%.1fd Az=%.1fd", fmod(deg.Y()+270,180)-90, fmod(deg.X()+720, 360));
620
621 TMarker *tip=new TMarker(v0.X(), v0.Y(), kDot);
622 tip->SetMarkerColor(kWhite+type*2);
623 AddMap(tip, new TString(txt));
624
625 return kTRUE;
626}
627
628// --------------------------------------------------------------------------
629//
630// Use "local" draw option to align the display to the local
631// coordinate system instead of the sky coordinate system.
632// dx, dy are arrays storing recuresively all touched points
633// stepx, stepy are the step-size of the current grid.
634//
635void MAstroCatalog::Draw(const TVector2 &v0, const TRotation &rot, TArrayI &dx, TArrayI &dy, Int_t stepx, Int_t stepy, Int_t type)
636{
637 // Calculate the end point
638 const TVector2 v1 = v0 + TVector2(dx[0]*TMath::DegToRad(), dy[0]*TMath::DegToRad());
639
640 // Check whether the point has already been touched.
641 Int_t idx[] = {1, 1, 1, 1};
642
643 Int_t dirs[4][2] = { {0, stepy}, {stepx, 0}, {0, -stepy}, {-stepx, 0} };
644
645 // Check for ambiguities.
646 for (int i=0; i<dx.GetSize(); i++)
647 {
648 for (int j=0; j<4; j++)
649 {
650 const Bool_t rcx0 = (dx[i]+720)%360==(dx[0]+dirs[j][0]+720)%360;
651 const Bool_t rcy0 = (dy[i]+360)%180==(dy[0]+dirs[j][1]+360)%180;
652 if (rcx0&&rcy0)
653 idx[j] = 0;
654 }
655 }
656
657 // Enhance size of array by 1, copy current
658 // position as last entry
659 dx.Set(dx.GetSize()+1);
660 dy.Set(dy.GetSize()+1);
661
662 dx[dx.GetSize()-1] = dx[0];
663 dy[dy.GetSize()-1] = dy[0];
664
665 // Store current positon
666 const Int_t d[2] = { dx[0], dy[0] };
667
668 for (int i=0; i<4; i++)
669 if (idx[i])
670 {
671 // Calculate new position
672 dx[0] = d[0]+dirs[i][0];
673 dy[0] = d[1]+dirs[i][1];
674
675 // Draw corresponding line and iterate through grid
676 if (DrawLine(v1, dirs[i][0], dirs[i][1], rot, type))
677 Draw(v0, rot, dx, dy, stepx, stepy, type);
678
679 dx[0]=d[0];
680 dy[0]=d[1];
681 }
682}
683
684// --------------------------------------------------------------------------
685//
686// Draw a grid recursively around the point v0 (either Ra/Dec or Zd/Az)
687// The points in the grid are converted by a TRotation and CovertToPad
688// to pad coordinates. The type arguemnts is neccessary to create the
689// correct tooltip (Ra/Dec, Zd/Az) at the grid-points.
690// From the pointing position the step-size of teh gris is caluclated.
691//
692void MAstroCatalog::DrawGrid(const TVector3 &v0, const TRotation &rot, Int_t type)
693{
694 TArrayI dx(1);
695 TArrayI dy(1);
696
697 // align to 1deg boundary
698 TVector2 v(v0.Phi()*TMath::RadToDeg(), v0.Theta()*TMath::RadToDeg());
699 v.Set((Float_t)TMath::Nint(v.X()), (Float_t)TMath::Nint(v.Y()));
700
701 // calculate stepsizes based on visible FOV
702 Int_t stepx = 1;
703
704 if (v.Y()<fRadiusFOV || v.Y()>180-fRadiusFOV)
705 stepx=36;
706 else
707 {
708 // This is a rough estimate how many degrees are visible
709 const Float_t m = log(fRadiusFOV/180.)/log(90./(fRadiusFOV+1)+1);
710 const Float_t t = log(180.)-m*log(fRadiusFOV);
711 const Float_t f = m*log(90-fabs(90-v.Y()))+t;
712 const Int_t nx = (Int_t)(exp(f)+0.5);
713 stepx = nx<4 ? 1 : nx/4;
714 if (stepx>36)
715 stepx=36;
716 }
717
718 const Int_t ny = (Int_t)(fRadiusFOV+1);
719 Int_t stepy = ny<4 ? 1 : ny/4;
720
721 // align stepsizes to be devisor or 180 and 90
722 while (180%stepx)
723 stepx++;
724 while (90%stepy)
725 stepy++;
726
727 // align to step-size boundary (search for the nearest one)
728 Int_t dv = 1;
729 while ((int)(v.X())%stepx)
730 {
731 v.Set(v.X()+dv, v.Y());
732 dv = -TMath::Sign(TMath::Abs(dv)+1, dv);
733 }
734
735 dv = 1;
736 while ((int)(v.Y())%stepy)
737 {
738 v.Set(v.X(), v.Y()+dv);
739 dv = -TMath::Sign(TMath::Abs(dv)+1, dv);
740 }
741
742 // draw...
743 v *= TMath::DegToRad();
744
745 Draw(v, rot, dx, dy, stepx, stepy, type);
746}
747
748// --------------------------------------------------------------------------
749//
750// Get a rotation matrix which aligns the pointing position
751// to the center of the x,y plain
752//
753TRotation MAstroCatalog::AlignCoordinates(const TVector3 &v) const
754{
755 TRotation trans;
756 trans.RotateZ(-v.Phi());
757 trans.RotateY(-v.Theta());
758 trans.RotateZ(-TMath::Pi()/2);
759 return trans;
760}
761
762// --------------------------------------------------------------------------
763//
764// Return the rotation matrix which converts either sky or
765// local coordinates to coordinates which pole is the current
766// pointing direction.
767//
768TRotation MAstroCatalog::GetGrid(Bool_t local)
769{
770 const Bool_t enable = fTime && fObservatory;
771
772 // If sky coordinate view is requested get rotation matrix and
773 // draw corresponding sky-grid and if possible local grid
774 if (!local)
775 {
776 const TRotation trans(AlignCoordinates(fRaDec));
777
778 DrawGrid(fRaDec, trans, 1);
779
780 if (enable)
781 {
782 const MAstroSky2Local rot(*fTime, *fObservatory);
783 DrawGrid(rot*fRaDec, trans*rot.Inverse(), 2);
784 }
785
786 // Return the correct rotation matrix
787 return trans;
788 }
789
790 // If local coordinate view is requested get rotation matrix and
791 // draw corresponding sky-grid and if possible local grid
792 if (local && enable)
793 {
794 const MAstroSky2Local rot(*fTime, *fObservatory);
795
796 const TRotation trans(AlignCoordinates(rot*fRaDec));
797
798 DrawGrid(fRaDec, trans*rot, 1);
799 DrawGrid(rot*fRaDec, trans, 2);
800
801 // Return the correct rotation matrix
802 return trans*rot;
803 }
804
805 return TRotation();
806}
807
808// --------------------------------------------------------------------------
809//
810// Create the title for the pad.
811//
812TString MAstroCatalog::GetPadTitle() const
813{
814 const Double_t ra = fRaDec.Phi()*TMath::RadToDeg();
815 const Double_t dec = (TMath::Pi()/2-fRaDec.Theta())*TMath::RadToDeg();
816
817 TString txt;
818 txt += Form("\\alpha=%.2fh ", fmod(ra/15+48, 24));
819 txt += Form("\\delta=%.1f\\circ ", fmod(dec+270,180)-90);
820 txt += Form("/ FOV=%.1f\\circ", fRadiusFOV);
821
822 if (!fTime || !fObservatory)
823 return txt;
824
825 const MAstroSky2Local rot(*fTime, *fObservatory);
826 const TVector3 loc = rot*fRaDec;
827
828 const Double_t rho = rot.RotationAngle(fRaDec.Phi(), TMath::Pi()/2-fRaDec.Theta());
829
830 const Double_t zd = TMath::RadToDeg()*loc.Theta();
831 const Double_t az = TMath::RadToDeg()*loc.Phi();
832
833 txt.Prepend("#splitline{");
834 txt += Form(" \\theta=%.1f\\circ ", fmod(zd+270,180)-90);
835 txt += Form("\\phi=%.1f\\circ ", fmod(az+720, 360));
836 txt += Form(" / \\rho=%.1f\\circ", rho*TMath::RadToDeg());
837 txt += "}{<";
838 txt += fTime->GetSqlDateTime();
839 txt += ">}";
840 return txt;
841}
842
843// --------------------------------------------------------------------------
844//
845// To overlay the catalog make sure, that in any case you are using
846// the 'same' option.
847//
848// If you want to overlay this on top of any picture which is created
849// by derotation of the camera plain you have to use the 'mirror' option
850// the compensate the mirroring of the image in the camera plain.
851//
852// If you have already compensated this by x=-x and y=-y when creating
853// the histogram you can simply overlay the catalog.
854//
855// To overlay the catalog on a 2D histogram the histogram must have
856// units of degrees (which are plain, like you directly convert the
857// camera units by multiplication to degrees)
858//
859// To be 100% exact you must use the option 'plain' which assumes a plain
860// screen. This is not necessary for the MAGIC-camera because the
861// difference between both is less than 1e-3.
862//
863// You should always be aware of the fact, that the shown stars and the
864// displayed grid is the ideal case, like a reflection on a virtual
865// perfectly aligned central mirror. In reality the star-positions are
866// smeared to the edge of the camera the more the distance to the center
867// is, such that the center of gravity of the light distribution might
868// be more far away from the center than the display shows.
869//
870//
871void MAstroCatalog::AddPrimitives(TString o)
872{
873 const Bool_t same = o.Contains("same", TString::kIgnoreCase);
874 const Bool_t local = o.Contains("local", TString::kIgnoreCase);
875 const Bool_t mirx = o.Contains("mirrorx", TString::kIgnoreCase);
876 const Bool_t miry = o.Contains("mirrory", TString::kIgnoreCase);
877 const Bool_t mirror = o.Contains("mirror", TString::kIgnoreCase) && !mirx && !miry;
878
879 // X is vice versa, because ra is defined anti-clockwise
880 mirx || mirror ? ResetBit(kMirrorX) : SetBit(kMirrorX);
881 miry || mirror ? SetBit(kMirrorY) : ResetBit(kMirrorY);
882
883 const TRotation rot(GetGrid(local));
884
885 TIter Next(&fList);
886 TVector3 *v=0;
887 while ((v=(TVector3*)Next()))
888 {
889 // FIXME: Check Magnitude!
890 TVector2 s(v->Phi(), v->Theta());
891 if (Convert(rot, s)==kTRUE)
892 DrawStar(s.X(), s.Y(), *v, kFALSE);
893 }
894
895 if (!same)
896 {
897 TPaveText *pv = new TPaveText(0.01, 0.90, 0.63, 0.99, "brNDC");
898 pv->AddText(GetPadTitle());
899 AddMap(pv);
900 }
901
902 TMarker *mk=new TMarker(0, 0, kMultiply);
903 mk->SetMarkerColor(kBlack);
904 mk->SetMarkerSize(1.5);
905 AddMap(mk);
906}
907
908// --------------------------------------------------------------------------
909//
910// Do nothing if 'same' option given.
911// Otherwise set pad-range such that x- and y- coordinates have the same
912// step-size
913//
914void MAstroCatalog::SetRangePad(Option_t *o)
915{
916 if (TString(o).Contains("same", TString::kIgnoreCase))
917 return;
918
919 const Double_t edge = fRadiusFOV/TMath::Sqrt(2.);
920 gPad->Range(-edge, -edge, edge, edge);
921
922 const Float_t w = gPad->GetWw();
923 const Float_t h = gPad->GetWh();
924
925 if (w<h)
926 gPad->Range(-edge, -edge*h/w, edge, edge*h/w);
927 else
928 gPad->Range(-edge*w/h, -edge, edge*w/h, edge);
929}
930
931// --------------------------------------------------------------------------
932//
933// First delete all gui elements.
934// Set the correct range of the pad.
935// Create all gui primitives
936// If "this" is not in the current pad add it to the current pad.
937// Reset vit kHasChanged
938//
939void MAstroCatalog::DrawPrimitives(Option_t *o)
940{
941 DeleteMap();
942
943 SetRangePad(o);
944
945#ifdef DEBUG
946 TStopwatch clk;
947 clk.Start();
948#endif DEBUG
949 AddPrimitives(o);
950#ifdef DEBUG
951 clk.Stop();
952 clk.Print();
953#endif DEBUG
954
955 // Append to a possible second pad
956 if (!gPad->GetListOfPrimitives()->FindObject(this))
957 AppendPad(o);
958
959 ResetBit(kHasChanged);
960}
961
962// --------------------------------------------------------------------------
963//
964// Call Paint() of all gui elements
965//
966void MAstroCatalog::PaintMap()
967{
968 Long_t key, val;
969 TExMapIter map(&fMapG);
970 while (map.Next(key, val))
971 ((TObject*)key)->Paint();
972}
973
974// --------------------------------------------------------------------------
975//
976// Append "this" to current pad
977// set bit kHasChanged to recreate all gui elements
978// Connect signal
979//
980void MAstroCatalog::Draw(Option_t *o)
981{
982 // Append to first pad
983 AppendPad(o);
984
985 // If contents have not previously changed make sure that
986 // all primitives are recreated.
987 SetBit(kHasChanged);
988
989 // Connect all TCanvas::ProcessedEvent to this->EventInfo
990 // This means, that after TCanvas has processed an event
991 // EventInfo of this class is called, see TCanvas::HandleInput
992 gPad->GetCanvas()->Connect("ProcessedEvent(Int_t,Int_t,Int_t,TObject*)",
993 "MAstroCatalog", this,
994 "EventInfo(Int_t,Int_t,Int_t,TObject*)");
995}
996
997// --------------------------------------------------------------------------
998//
999// This function was connected to all created canvases. It is used
1000// to redirect GetObjectInfo into our own status bar.
1001//
1002// The 'connection' is done in AddTab
1003//
1004void MAstroCatalog::EventInfo(Int_t event, Int_t px, Int_t py, TObject *selected)
1005{
1006 TCanvas *c = (TCanvas*)gTQSender;
1007
1008 gPad = c ? c->GetSelectedPad() : NULL;
1009 if (!gPad)
1010 return;
1011
1012 // Try to find a corresponding object with kCannotPick set and
1013 // an available TString (for a tool tip)
1014 TString *str=0;
1015 if (!selected || selected==this)
1016 {
1017 Long_t key, val;
1018 TExMapIter map(&fMapG);
1019 while (map.Next(key, val))
1020 {
1021 if (!val)
1022 continue;
1023
1024 TObject *o=(TObject*)key;
1025 if (o->DistancetoPrimitive(px, py)>TPad::GetMaxPickDistance())
1026 continue;
1027
1028 selected = o;
1029 str = (TString*)val;
1030 break;
1031 }
1032 }
1033
1034 if (!selected)
1035 return;
1036
1037 // Handle some gui events
1038 switch (event)
1039 {
1040 case kMouseMotion:
1041 if (!fToolTip->IsMapped() && str)
1042 ShowToolTip(px, py, *str);
1043 break;
1044
1045 case kMouseLeave:
1046 if (fToolTip->IsMapped())
1047 fToolTip->Hide();
1048 break;
1049
1050 case kKeyPress:
1051 ExecuteEvent(kKeyPress, px, py);
1052 break;
1053 }
1054}
1055
1056// --------------------------------------------------------------------------
1057//
1058// Handle keyboard events.
1059//
1060void MAstroCatalog::ExecuteEventKbd(Int_t keycode, Int_t keysym)
1061{
1062 Double_t dra =0;
1063 Double_t ddec=0;
1064
1065 switch (keysym)
1066 {
1067 case kKey_Left:
1068 dra = -TMath::DegToRad();
1069 break;
1070 case kKey_Right:
1071 dra = +TMath::DegToRad();
1072 break;
1073 case kKey_Up:
1074 ddec = +TMath::DegToRad();
1075 break;
1076 case kKey_Down:
1077 ddec = -TMath::DegToRad();
1078 break;
1079 case kKey_Plus:
1080 SetRadiusFOV(fRadiusFOV+1);
1081 break;
1082 case kKey_Minus:
1083 SetRadiusFOV(fRadiusFOV-1);
1084 break;
1085
1086 default:
1087 return;
1088 }
1089
1090 const Double_t r = fRaDec.Phi();
1091 const Double_t d = TMath::Pi()/2-fRaDec.Theta();
1092
1093 SetRaDec(r+dra, d+ddec);
1094
1095 gPad->Update();
1096}
1097
1098// ------------------------------------------------------------------------
1099//
1100// Execute a gui event on the camera
1101//
1102void MAstroCatalog::ExecuteEvent(Int_t event, Int_t mp1, Int_t mp2)
1103{
1104 if (!TestBit(kGuiActive))
1105 return;
1106
1107 if (event==kKeyPress)
1108 ExecuteEventKbd(mp1, mp2);
1109}
1110
1111// --------------------------------------------------------------------------
1112//
1113// Calculate distance to primitive by checking all gui elements
1114//
1115Int_t MAstroCatalog::DistancetoPrimitive(Int_t px, Int_t py)
1116{
1117 Int_t min = INT_MAX;
1118
1119 Long_t key, val;
1120 TExMapIter map(&fMapG);
1121 while (map.Next(key, val))
1122 {
1123 TObject *o=(TObject*)key;
1124
1125 const Int_t d = o->DistancetoPrimitive(px, py);
1126
1127 if (d<TPad::GetMaxPickDistance())
1128 return 0;
1129
1130 if (d<min)
1131 min=d;
1132 }
1133
1134 return min;
1135}
1136
1137// --------------------------------------------------------------------------
1138void MAstroCatalog::ShowToolTip(Int_t px, Int_t py, const char *txt)
1139{
1140 Int_t x=0;
1141 Int_t y=0;
1142
1143 const Window_t id1 = gVirtualX->GetWindowID(gPad->GetCanvasID());
1144 const Window_t id2 = fToolTip->GetParent()->GetId();
1145
1146 Window_t id3;
1147 gVirtualX->TranslateCoordinates(id1, id2, px, py, x, y, id3);
1148
1149 // Show tool tip
1150 fToolTip->SetText(txt);
1151 fToolTip->Show(x+4, y+4);
1152}
1153
1154/*
1155void MAstroCatalog::RecursiveRemove(TObject *obj)
1156{
1157 ULong_t hash;
1158 Long_t key, val;
1159
1160 TExMapIter map(&fMapG);
1161 while (map.Next(hash, key, val))
1162 {
1163 if (key != (Long_t)obj)
1164 continue;
1165
1166 fMapG.Remove(hash, key);
1167 delete (TObject*)(key);
1168 if (val)
1169 delete (TString*)(val);
1170 break;
1171 }
1172}
1173*/
Note: See TracBrowser for help on using the repository browser.