Ignore:
Timestamp:
03/22/04 09:25:49 (21 years ago)
Author:
tbretz
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars/mastro
Files:
9 edited

Legend:

Unmodified
Added
Removed

  • trunk/MagicSoft/Mars/mastro/MAstro.cc

    r3537 r3568  
    310310    return MTime(ut1).GetGmst();
    311311}
     312
     313// --------------------------------------------------------------------------
     314//
     315// RotationAngle
     316//
     317// calculates the angle for the rotation of the sky coordinate system
     318// with respect to the local coordinate system. This is identical
     319// to the rotation angle of the sky image in the camera.
     320//
     321//  sinl  [rad]: sine of observers latitude
     322//  cosl  [rad]: cosine of observers latitude
     323//  theta [rad]: polar angle/zenith distance
     324//  phi   [rad]: rotation angle/azimuth
     325//
     326// Return sin/cos component of angle
     327//
     328// The convention is such, that the rotation angle is -pi/pi if
     329// right ascension and local rotation angle are counted in the
     330// same direction, 0 if counted in the opposite direction.
     331//
     332// (In other words: The rotation angle is 0 when the source culminates)
     333//
     334// Using vectors it can be done like:
     335//    TVector3 v, p;
     336//    v.SetMagThetaPhi(1, theta, phi);
     337//    p.SetMagThetaPhi(1, TMath::Pi()/2-latitude, 0);
     338//    v = v.Cross(l));
     339//    v.RotateZ(-phi);
     340//    v.Rotate(-theta)
     341//    rho = TMath::ATan2(v(2), v(1));
     342//
     343// For more information see TDAS 00-11, eqs. (18) and (20)
     344//
     345void MAstro::RotationAngle(Double_t sinl, Double_t cosl, Double_t theta, Double_t phi, Double_t &sin, Double_t &cos)
     346{
     347    const Double_t sint = TMath::Sin(theta);
     348    const Double_t cost = TMath::Cos(theta);
     349
     350    const Double_t snlt = sinl*sint;
     351    const Double_t cslt = cosl*cost;
     352
     353    const Double_t sinp = TMath::Sin(phi);
     354    const Double_t cosp = TMath::Cos(phi);
     355
     356    const Double_t v1 = sint*sinp;
     357    const Double_t v2 = cosl - snlt*cosp;
     358
     359    const Double_t denom = TMath::Sqrt(v1*v1 + v2*v2);
     360
     361    cos =  -cosl*sinp      / denom;
     362    sin = (snlt-cslt*cosp) / denom;
     363}
     364
     365// --------------------------------------------------------------------------
     366//
     367// RotationAngle
     368//
     369// calculates the angle for the rotation of the sky coordinate system
     370// with respect to the local coordinate system. This is identical
     371// to the rotation angle of the sky image in the camera.
     372//
     373//  sinl  [rad]: sine of observers latitude
     374//  cosl  [rad]: cosine of observers latitude
     375//  theta [rad]: polar angle/zenith distance
     376//  phi   [rad]: rotation angle/azimuth
     377//
     378// Return angle [rad] in the range -pi, pi
     379//
     380// The convention is such, that the rotation angle is -pi/pi if
     381// right ascension and local rotation angle are counted in the
     382// same direction, 0 if counted in the opposite direction.
     383//
     384// (In other words: The rotation angle is 0 when the source culminates)
     385//
     386// Using vectors it can be done like:
     387//    TVector3 v, p;
     388//    v.SetMagThetaPhi(1, theta, phi);
     389//    p.SetMagThetaPhi(1, TMath::Pi()/2-latitude, 0);
     390//    v = v.Cross(l));
     391//    v.RotateZ(-phi);
     392//    v.Rotate(-theta)
     393//    rho = TMath::ATan2(v(2), v(1));
     394//
     395// For more information see TDAS 00-11, eqs. (18) and (20)
     396//
     397Double_t MAstro::RotationAngle(Double_t sinl, Double_t cosl, Double_t theta, Double_t phi)
     398{
     399    const Double_t snlt = sinl*TMath::Sin(theta);
     400    const Double_t cslt = cosl*TMath::Cos(theta);
     401
     402    const Double_t sinp = TMath::Sin(phi);
     403    const Double_t cosp = TMath::Cos(phi);
     404
     405    return TMath::ATan2(-cosl*sinp, snlt-cslt*cosp);
     406}

  • trunk/MagicSoft/Mars/mastro/MAstro.h

    r3371 r3568  
    5050    static Double_t UT2GMST(Double_t ut1);
    5151
     52    // Rotation angle between local and sky coordinate system
     53    static void     RotationAngle(Double_t sinl, Double_t cosl, Double_t theta, Double_t phi, Double_t &sin, Double_t &cos);
     54    static Double_t RotationAngle(Double_t sinl, Double_t cosl, Double_t theta, Double_t phi);
     55
    5256    ClassDef(MAstro, 0)
    5357};

  • trunk/MagicSoft/Mars/mastro/MAstroCamera.cc

    r3537 r3568  
    5757using namespace std;
    5858
     59// --------------------------------------------------------------------------
    5960MAstroCamera::MAstroCamera() : fGeom(0), fMirrors(0)
    6061{
    … …  
    6364}
    6465
     66// --------------------------------------------------------------------------
    6567MAstroCamera::~MAstroCamera()
    6668{
    … …  
    7375}
    7476
    75 void MAstroCamera::SetMirrors(TClonesArray *arr)
    76 {
    77     if (!arr || arr->GetClass()!=MGeomMirror::Class())
     77// --------------------------------------------------------------------------
     78void MAstroCamera::SetMirrors(TClonesArray &arr)
     79{
     80    if (arr.GetClass()!=MGeomMirror::Class())
    7881        return;
    7982
    80     const Int_t n = arr->GetSize();
     83    const Int_t n = arr.GetSize();
    8184
    8285    if (!fMirrors)
    … …  
    8689
    8790    for (int i=0; i<n; i++)
    88         memcpy((*fMirrors)[i], (*arr)[i], sizeof(MGeomMirror));
    89 
    90 }
    91 
     91        memcpy((*fMirrors)[i], arr[i], sizeof(MGeomMirror));
     92
     93}
     94
     95// --------------------------------------------------------------------------
    9296void MAstroCamera::SetGeom(const MGeomCam &cam)
    9397{
    … …  
    98102}
    99103
    100 Int_t MAstroCamera::ConvertToPad(const TVector3 &w, TVector2 &v)
    101 {
    102     const TVector3 spot = fMirror0->GetReflection(w, fGeom->GetCameraDist())*1000;
    103 
     104// --------------------------------------------------------------------------
     105Int_t MAstroCamera::ConvertToPad(const TVector3 &w, TVector2 &v) const
     106{
    104107    /*
    105108     --- Use this to plot the 'mean grid' instead of the grid of a
    … …  
    113116        */
    114117
    115 
     118    const TVector3 spot = fMirror0->GetReflection(w, fGeom->GetCameraDist())*1000;
    116119    v.Set(spot(0), spot(1));
    117120
    … …  
    120123}
    121124
    122 void MAstroCamera::DrawNet(const TRotation &trans)
    123 {
    124     const TRotation rot(MAstroSky2Local(*fTime, *fObservatory));
    125 
    126     TVector2 radec(fRaDec.Phi(), fRaDec.Theta());
    127     MAstroCatalog::DrawNet(radec, trans*rot, 2);
    128 
    129     const TVector3 zdaz0 = MAstroSky2Local(*fTime, *fObservatory)*fRaDec;
    130     TVector2 zdaz(zdaz0.Phi(), zdaz0.Theta());
    131     MAstroCatalog::DrawNet(zdaz, trans, 1);
    132 }
    133 
     125// --------------------------------------------------------------------------
    134126TObject *FindObjectInPad(const char *name, TVirtualPad *pad)
    135127{
    … …  
    155147}
    156148
     149// --------------------------------------------------------------------------
    157150void MAstroCamera::AddPrimitives(Option_t *o)
    158151{
    … …  
    171164    const Bool_t hasdot  = opt.Contains(".", TString::kIgnoreCase);
    172165    const Bool_t usecam  = opt.Contains("c", TString::kIgnoreCase);
    173 
    174     MAstroSky2Local rot(*fTime, *fObservatory);
    175 
    176     const Float_t rho = rot.RotationAngle(fRaDec.Phi(), TMath::Pi()/2-fRaDec.Theta());
    177 
    178     TString str = fTime->GetSqlDateTime();
    179     str += Form("  (\\alpha=%.1fh \\delta=%.1f\\circ)  \\rho=%.1f\\circ",
    180                 fRaDec.Phi()/TMath::Pi()*12, 90-fRaDec.Theta()*TMath::RadToDeg(),
    181                 rho *TMath::RadToDeg());
    182166
    183167    // Get camera
    … …  
    198182    }
    199183
    200     camera->SetTitle(str);
     184    camera->SetTitle(GetPadTitle());
    201185
    202186    gPad->cd(1);
    … …  
    221205    }
    222206
    223     TVector3 zdaz0 = fRaDec;
    224     zdaz0 *= rot;
    225 
    226     cout << zdaz0.Phi()*TMath::RadToDeg() << " " << zdaz0.Theta()*TMath::RadToDeg() << endl;
    227 
    228     TVector3 test = zdaz0;
    229     test *= rot.Inverse();
    230 
    231     cout << test.Phi()*TMath::RadToDeg()/15 << " " << test.Theta()*TMath::RadToDeg() << endl;
    232 
    233     TRotation rot2;
    234     rot2.RotateZ(-zdaz0.Phi());
    235     rot2.RotateY(-zdaz0.Theta());
    236     rot2.RotateZ(-TMath::Pi()/2); // align coordinate system
    237 
    238     DrawNet(rot2);
     207    const TRotation rot(GetGrid(kTRUE));
    239208
    240209    MVector3 *radec;
    … …  
    247216        TVector3 star(*radec);
    248217
    249         // Calculate local coordinates
     218        // Rotate Star into telescope system
    250219        star *= rot;
    251         // Rotate Star into telescope system
    252         star *= rot2;
    253220
    254221        TVector3 mean;
    … …  
    270237            {
    271238                TMarker *m=new TMarker(spot(0), spot(1), 1);
    272                 m->SetBit(kCannotPick);
    273                 m->SetBit(kCanDelete);
    274239                m->SetMarkerColor(kMagenta);
    275240                m->SetMarkerStyle(kDot);
    … …  
    303268        case kKey_Plus:
    304269            fTime->SetMjd(fTime->GetMjd()+0.25/24);
    305             SetBit(kHasChanged);
    306             gPad->Modified();
    307             gPad->Update();
     270            Update(kTRUE);
    308271            return;
    309272
    310273        case kKey_Minus:
    311274            fTime->SetMjd(fTime->GetMjd()-0.25/24);
    312             SetBit(kHasChanged);
    313             gPad->Modified();
    314             gPad->Update();
     275            Update(kTRUE);
    315276            return;
    316277        }

  • trunk/MagicSoft/Mars/mastro/MAstroCamera.h

    r3537 r3568  
    2121    MGeomMirror  *fMirror0;     //!
    2222
    23     Int_t  ConvertToPad(const TVector3 &w, TVector2 &v);
     23    Int_t  ConvertToPad(const TVector3 &w, TVector2 &v) const;
    2424    void   AddPrimitives(Option_t *o);
    2525    void   SetRangePad() { }
    26     void   DrawNet(const TRotation &rot);
    2726    void   ExecuteEvent(Int_t event, Int_t mp1, Int_t mp2);
    2827
    … …  
    3130    ~MAstroCamera();
    3231
    33     void SetMirrors(TClonesArray *arr);
     32    void SetMirrors(TClonesArray &arr);
    3433    void SetGeom(const MGeomCam &cam);
    3534

  • trunk/MagicSoft/Mars/mastro/MAstroCatalog.cc

    r3537 r3568  
    3333#include "MAstroCatalog.h"
    3434
     35#include <errno.h>
    3536#include <fstream>
    3637#include <stdlib.h>
     38#include <limits.h> // INT_MAX (Suse 7.3/gcc 2.95)
    3739
    3840#include <KeySymbols.h>
    … …  
    4547#include <TGToolTip.h>
    4648#include <TRotation.h>
     49#include <TPaveText.h>
    4750#include <TStopwatch.h>
    4851
    … …  
    5962
    6063using namespace std;
    61 
     64/*
    6265class MRotation : public TRotation
    6366{
    … …  
    7679    }
    7780};
    78 
     81*/
    7982/*
    8083MVector3 MVector3::GetZdAz(const MObservatory &obs, Double_t gmst) const
    … …  
    134137}
    135138*/
    136 MAstroCatalog::MAstroCatalog() : fLimMag(99), fRadiusFOV(99), fToolTip(0), fObservatory(0), fTime(0)
     139// --------------------------------------------------------------------------
     140MAstroCatalog::MAstroCatalog() : fLimMag(99), fRadiusFOV(90), fToolTip(0), fObservatory(0), fTime(0)
    137141{
    138142    fList.SetOwner();
    … …  
    140144}
    141145
     146// --------------------------------------------------------------------------
    142147MAstroCatalog::~MAstroCatalog()
    143148{
    … …  
    161166}
    162167
     168// --------------------------------------------------------------------------
    163169TString MAstroCatalog::FindToken(TString &line, Char_t tok)
    164170{
    … …  
    176182}
    177183
     184// --------------------------------------------------------------------------
    178185Int_t MAstroCatalog::atoi(const TSubString &sub)
    179186{
    … …  
    181188}
    182189
     190// --------------------------------------------------------------------------
    183191Float_t MAstroCatalog::atof(const TSubString &sub)
    184192{
    … …  
    186194}
    187195
     196// --------------------------------------------------------------------------
    188197Int_t MAstroCatalog::atoi(const TString &s)
    189198{
    … …  
    191200}
    192201
     202// --------------------------------------------------------------------------
    193203Float_t MAstroCatalog::atof(const TString &s)
    194204{
    … …  
    196206}
    197207
     208// --------------------------------------------------------------------------
    198209Int_t MAstroCatalog::ReadXephem(TString catalog)
    199210{
    … …  
    203214
    204215    ifstream fin(catalog);
     216    if (!fin)
     217    {
     218        gLog << err << "Cannot open file " << catalog << ": ";
     219        gLog << strerror(errno) << endl;
     220        return 0;
     221    }
    205222
    206223    Int_t add =0;
    … …  
    276293}
    277294
     295// --------------------------------------------------------------------------
    278296Int_t MAstroCatalog::ReadNGC2000(TString catalog)
    279297{
    … …  
    283301
    284302    ifstream fin(catalog);
     303    if (!fin)
     304    {
     305        gLog << err << "Cannot open file " << catalog << ": ";
     306        gLog << strerror(errno) << endl;
     307        return 0;
     308    }
    285309
    286310    Int_t add=0;
    … …  
    339363}
    340364
     365// --------------------------------------------------------------------------
    341366Int_t MAstroCatalog::ReadBSC(TString catalog)
    342367{
    … …  
    346371
    347372    ifstream fin(catalog);
     373    if (!fin)
     374    {
     375        gLog << err << "Cannot open file " << catalog << ": ";
     376        gLog << strerror(errno) << endl;
     377        return 0;
     378    }
    348379
    349380    Int_t add=0;
    … …  
    404435}
    405436
     437// --------------------------------------------------------------------------
    406438void MAstroCatalog::Paint(Option_t *o)
    407439{
    408 //    if (!gPad->IsBatch())
    409 //        gVirtualX->ClearWindow();
     440    SetRangePad();
    410441
    411442    if (TestBit(kHasChanged))
    … …  
    413444}
    414445
     446// --------------------------------------------------------------------------
    415447void MAstroCatalog::DrawStar(Double_t x, Double_t y, const TVector3 &v, Bool_t transparent, const char *txt)
    416448{
    … …  
    428460    TMarker *tip=new TMarker(x, y, transparent ? kDot : kFullDotLarge);;
    429461    tip->SetMarkerColor(kBlack);
    430     tip->SetBit(kCanDelete);
    431     tip->SetBit(kCannotPick);
    432462    AddMap(tip, new TString(str));
    433463}
    434464
    435 void MAstroCatalog::Update()
     465// --------------------------------------------------------------------------
     466void MAstroCatalog::Update(Bool_t upd)
    436467{
    437468    if (gPad && TestBit(kMustCleanup))
    … …  
    439470        SetBit(kHasChanged);
    440471        gPad->Modified();
    441     }
    442 }
    443 
     472        if (upd)
     473            gPad->Update();
     474    }
     475}
     476
     477// --------------------------------------------------------------------------
     478//
     479// Set the observation time. Necessary to use local coordinate
     480// system.
     481//
    444482void MAstroCatalog::SetTime(const MTime &time)
    445483{
    … …  
    449487}
    450488
     489// --------------------------------------------------------------------------
     490//
     491// Set the observatory location. Necessary to use local coordinate
     492// system.
     493//
    451494void MAstroCatalog::SetObservatory(const MObservatory &obs)
    452495{
    … …  
    456499}
    457500
    458 Int_t MAstroCatalog::ConvertToPad(const TVector3 &w0, TVector2 &v)
     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
    459510{
    460511    TVector3 w(w0);
    461512
    462     // Stretch such, that the X-component is alwas the same. Now
    463     // Y and Z contains the crossing point between the star-light
    464     // and the plain of a virtual screen (ccd...)
     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...)
    465516    if (TestBit(kPlainScreen))
    466         w *= TMath::RadToDeg()/w.X();
    467     else
    468         w.SetMag(TMath::RadToDeg());
    469 
    470     v.Set(w(1), w(2));
    471 
    472     if (w(0)<0)
     517        w *= 1./w(2);
     518
     519    w *= TMath::RadToDeg();
     520    v.Set(w(0), w(1));
     521
     522    if (w(2)<0)
    473523        return kERROR;
    474524
    475     return w(1)>gPad->GetX1() && w(2)>gPad->GetY1() &&
    476            w(1)<gPad->GetX2() && w(2)<gPad->GetY2();
    477 }
    478 
    479 Int_t MAstroCatalog::Convert(const TRotation &rot, TVector2 &v)
     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
    480531{
    481532    MVector3 w;
    … …  
    486537}
    487538
     539// --------------------------------------------------------------------------
    488540Bool_t MAstroCatalog::DrawLine(const TVector2 &v, Double_t dx, Double_t dy, const TRotation &rot, Int_t type)
    489541{
    … …  
    501553
    502554    TLine *line = new TLine(v0.X(), v0.Y(), v1.X(), v1.Y());
    503     line->SetBit(kCanDelete);
    504555    line->SetLineStyle(kDashDotted); //kDashed, kDotted, kDashDotted
    505556    line->SetLineColor(kWhite+type*2);
    506     line->SetBit(kCannotPick);
    507557    AddMap(line);
    508558
    509559    const TVector2 deg = v*TMath::RadToDeg();
    510560    TString txt = type==1 ?
    511         Form("Ra=%.1fh  Dec=%.1fd", fmod(deg.X()/15+48, 24),  fmod(90-deg.Y()+270,180)-90) :
     561        Form("Ra=%.2fh  Dec=%.1fd", fmod(deg.X()/15+48, 24),  fmod(90-deg.Y()+270,180)-90) :
    512562        Form("Zd=%.1fd  Az=%.1fd",  fmod(deg.Y()+270,180)-90, fmod(deg.X()+720, 360));
    513563
    514564    TMarker *tip=new TMarker(v0.X(), v0.Y(), kDot);
    515     tip->SetBit(kCanDelete);
    516     tip->SetBit(kCannotPick);
    517565    tip->SetMarkerColor(kWhite+type*2);
    518566    AddMap(tip, new TString(txt));
    … …  
    521569}
    522570
    523 
     571// --------------------------------------------------------------------------
     572//
     573// Use "local" draw option to align the display to the local
     574// coordinate system instead of the sky coordinate system.
     575//
    524576void MAstroCatalog::Draw(const TVector2 &v0, const TRotation &rot, TArrayI &dx, TArrayI &dy, Int_t stepx, Int_t stepy, Int_t type)
    525577{
    … …  
    570622}
    571623
    572 void MAstroCatalog::DrawNet(const TVector2 &v0, const TRotation &rot, Int_t type)
     624// --------------------------------------------------------------------------
     625void MAstroCatalog::DrawGrid(const TVector3 &v0, const TRotation &rot, Int_t type)
    573626{
    574627    TArrayI dx(1);
    … …  
    576629
    577630    // align to 1deg boundary
    578     TVector2 v = v0*TMath::RadToDeg();
     631    TVector2 v(v0.Phi()*TMath::RadToDeg(), v0.Theta()*TMath::RadToDeg());
    579632    v.Set((Float_t)TMath::Nint(v.X()), (Float_t)TMath::Nint(v.Y()));
    580633
    581634    // calculate stepsizes based on visible FOV
    582     Int_t stepx=1;
    583 
    584     if (fabs(90-v.Y())>90-fRadiusFOV || fabs(90-v.Y())<fRadiusFOV)
    585         stepx = 180/10;
     635    Int_t stepx = 1;
     636
     637    if (v.Y()<fRadiusFOV || v.Y()>180-fRadiusFOV)
     638        stepx=36;
    586639    else
    587640    {
    588641        // This is a rough estimate how many degrees are visible
    589         const Float_t m = log(fRadiusFOV/180.)/log(90./fRadiusFOV-1);
     642        const Float_t m = log(fRadiusFOV/180.)/log(90./(fRadiusFOV+1)+1);
    590643        const Float_t t = log(180.)-m*log(fRadiusFOV);
    591         const Int_t n = (Int_t)(exp(m*log(90-fabs(90-v.Y()))+t)+0.5);
    592         stepx = n<6 ? 1 : n/6;
    593     }
    594 
    595     const Int_t n = (Int_t)(fRadiusFOV+1);
    596     Int_t stepy = n<4 ? 1 : n/4;
     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;
    597653
    598654    // align stepsizes to be devisor or 180 and 90
    … …  
    602658        stepy++;
    603659
    604     // align to step-size boundary
     660    // align to step-size boundary (search for the nearest one)
     661    Int_t dv = 1;
    605662    while ((int)(v.X())%stepx)
    606         v.Set(v.X()+1, v.Y());
    607 
     663    {
     664        v.Set(v.X()+dv, v.Y());
     665        dv = -TMath::Sign(TMath::Abs(dv)+1, dv);
     666    }
     667
     668    dv = 1;
    608669    while ((int)(v.Y())%stepy)
    609         v.Set(v.X(), v.Y()+1);
     670    {
     671        v.Set(v.X(), v.Y()+dv);
     672        dv = -TMath::Sign(TMath::Abs(dv)+1, dv);
     673    }
    610674
    611675    // draw...
    … …  
    615679}
    616680
     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// --------------------------------------------------------------------------
    617769void MAstroCatalog::AddPrimitives(Option_t *o)
    618770{
    619     // Precalc Sin/Cos...
    620     TRotation trans;
    621     trans.RotateZ(-fRaDec.Phi());
    622     trans.Rotate(TMath::Pi()/2-fRaDec.Theta(), TVector3(0, 1, 0));
    623 
    624     if (fTime && fObservatory)
    625     {
    626         const TRotation rot(MAstroSky2Local(*fTime, *fObservatory));
    627         const TVector3 zdaz0 = rot*fRaDec;
    628         const TVector2 zdaz(zdaz0.Phi(), zdaz0.Theta());
    629         DrawNet(zdaz, trans*rot.Inverse(), 2);
    630     }
    631 
    632     const TVector2 radec(fRaDec.Phi(), fRaDec.Theta());
    633     DrawNet(radec, trans, 1);
     771    const Bool_t local = TString(o).Contains("local", TString::kIgnoreCase);
     772
     773    cout << "Opt: " << o << endl;
     774
     775    const TRotation rot(GetGrid(local));
    634776
    635777    TIter Next(&fList);
    … …  
    639781        // FIXME: Check Magnitude!
    640782        TVector2 s(v->Phi(), v->Theta());
    641         if (Convert(trans, s)==kTRUE)
     783        if (Convert(rot, s)==kTRUE)
    642784            DrawStar(s.X(), TMath::Pi()/2-s.Y(), *v, kFALSE);
    643785    }
    644 }
    645 
     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// --------------------------------------------------------------------------
    646798void MAstroCatalog::SetRangePad()
    647799{
    … …  
    658810}
    659811
     812// --------------------------------------------------------------------------
    660813void MAstroCatalog::DrawPrimitives(Option_t *o)
    661814{
    … …  
    680833}
    681834
     835// --------------------------------------------------------------------------
    682836void MAstroCatalog::DrawMap()
    683837{
    … …  
    688842}
    689843
     844// --------------------------------------------------------------------------
    690845void MAstroCatalog::Draw(Option_t *o)
    691846{
    … …  
    770925}
    771926
     927// --------------------------------------------------------------------------
    772928void MAstroCatalog::ExecuteEventKbd(Int_t keycode, Int_t keysym)
    773929{
    … …  
    821977}
    822978
     979// --------------------------------------------------------------------------
    823980Int_t MAstroCatalog::DistancetoPrimitive(Int_t px, Int_t py)
    824981{
    … …  
    8431000}
    8441001
     1002// --------------------------------------------------------------------------
    8451003void MAstroCatalog::ShowToolTip(Int_t px, Int_t py, const char *txt)
    8461004{
    … …  
    8581016    fToolTip->Show(x+4, y+4);
    8591017}
    860 

  • trunk/MagicSoft/Mars/mastro/MAstroCatalog.h

    r3537 r3568  
    7474{
    7575private:
    76     Double_t fLimMag;    // [1]   Limiting Magnitude
    77     Double_t fRadiusFOV; // [deg] Radius of Field of View
    78 
    79     TGToolTip *fToolTip; //!
     76    Double_t   fLimMag;    // [1]   Limiting Magnitude
     77    Double_t   fRadiusFOV; // [deg] Radius of Field of View
     78
     79    TExMap     fMapG;      //! A map with all gui primitives and tooltips
     80    TGToolTip *fToolTip;   //! The tooltip currently displayed
    8081
    8182    void ShowToolTip(Int_t px, Int_t py, const char *txt);
    … …  
    9293//#endif
    9394
    94 protected:
    95     enum {
    96         kHasChanged  = BIT(15),
    97         kGuiActive   = BIT(16),
    98         kPlainScreen = BIT(17)
    99     };
    100 
    101     TExMap   fMapG;
    102     TList    fList;      // List of stars loaded
    103     MVector3 fRaDec;     // pointing position
    104 
    105     MObservatory *fObservatory; // Possible obervatora location
    106     MTime        *fTime;        // Possible observation time
    107 
    108     virtual Int_t  ConvertToPad(const TVector3 &w, TVector2 &v);
    109     virtual Int_t  Convert(const TRotation &rot, TVector2 &v);
    110     virtual Bool_t DrawLine(const TVector2 &v0, Double_t dx, Double_t dy, const TRotation &rot, Int_t type);
    111     virtual void   AddPrimitives(Option_t *o);
    112     virtual void   DrawPrimitives(Option_t *o);
    113     virtual void   SetRangePad();
    114     void  Draw(const TVector2 &v0, const TRotation &rot, TArrayI &dx, TArrayI &dy, Int_t stepx, Int_t stepy, Int_t type);
    115     void  DrawNet(const TVector2 &v0, const TRotation &rot, Int_t type);
    116     void  DrawStar(Double_t x, Double_t y, const TVector3 &v, Bool_t t, const char *txt=0);
    117     void  Paint(Option_t *o="");
    118     Int_t DistancetoPrimitive(Int_t px, Int_t py);
    119     void  Update();
    120 
    121     void ExecuteEventKbd(Int_t keycode, Int_t keysym);
    122     void ExecuteEvent(Int_t event, Int_t mp1, Int_t mp2);
    123 
    124     void DrawMap();
    125     void AddMap(void *k, void *v=0) { fMapG.Add(fMapG.GetSize(), (Long_t)k, (Long_t)v); }
    126     void DeleteMap()
     95    virtual Int_t ConvertToPad(const TVector3 &w, TVector2 &v) const;
     96    virtual void  AddPrimitives(Option_t *o);
     97    virtual void  SetRangePad();
     98
     99    Int_t     Convert(const TRotation &rot, TVector2 &v) const;
     100    void      Draw(const TVector2 &v0, const TRotation &rot, TArrayI &dx, TArrayI &dy, Int_t stepx, Int_t stepy, Int_t type);
     101    void      DrawPrimitives(Option_t *o);
     102    Bool_t    DrawLine(const TVector2 &v0, Double_t dx, Double_t dy, const TRotation &rot, Int_t type);
     103    void      DrawGrid(const TVector3 &v0, const TRotation &rot, Int_t type);
     104    TRotation AlignCoordinates(const TVector3 &v) const;
     105    void      Paint(Option_t *o="");
     106    Int_t     DistancetoPrimitive(Int_t px, Int_t py);
     107    void      DrawMap();
     108    void      DeleteMap()
    127109    {
    128110        Long_t key, val;
    … …  
    148130    }
    149131
     132protected:
     133    enum {
     134        kHasChanged  = BIT(15),
     135        kGuiActive   = BIT(16),
     136        kPlainScreen = BIT(17)
     137    };
     138
     139    TList    fList;      // List of stars loaded
     140    MVector3 fRaDec;     // pointing position
     141
     142    MObservatory *fObservatory; // Possible obervatora location
     143    MTime        *fTime;        // Possible observation time
     144
     145    virtual TString GetPadTitle() const;
     146    TRotation GetGrid(Bool_t local);
     147    void      DrawStar(Double_t x, Double_t y, const TVector3 &v, Bool_t t, const char *txt=0);
     148    void      Update(Bool_t upd=kFALSE);
     149
     150    void ExecuteEventKbd(Int_t keycode, Int_t keysym);
     151    void ExecuteEvent(Int_t event, Int_t mp1, Int_t mp2);
     152
     153    void AddMap(TObject *k, void *v=0)
     154    {
     155        k->SetBit(kCanDelete);
     156        k->SetBit(kCannotPick);
     157        fMapG.Add(fMapG.GetSize(), (Long_t)k, (Long_t)v);
     158    }
     159
    150160public:
    151161    MAstroCatalog();
    … …  
    161171        if (deg>max)
    162172            deg=max;
    163         if (deg<=0)
     173        if (deg<1)
    164174            deg=1;
    165175
    … …  
    189199
    190200    void Draw(Option_t *o="");
     201    void SetDrawOption(Option_t *option="")
     202    {
     203        TObject::SetDrawOption(option);
     204        Update(kTRUE);
     205    } //*MENU*
    191206
    192207    virtual void EventInfo(Int_t event, Int_t px, Int_t py, TObject *selected=0);

  • trunk/MagicSoft/Mars/mastro/MAstroSky2Local.cc

    r3537 r3568  
    6262#include "MAstroSky2Local.h"
    6363
     64#include "MAstro.h"
    6465#include "MTime.h"
    6566#include "MObservatory.h"
    … …  
    128129// seen with an Alt/Az telescope.
    129130//
     131// For more information see MAstro::RotationAngle
     132//
    130133Double_t MAstroSky2Local::RotationAngle(Double_t ra, Double_t dec) const
    131134{
    132     TVector3 loc;
    133     loc.SetMagThetaPhi(1, TMath::Pi()/2-dec, ra);
    134     loc *= *this;
     135    TVector3 v;
     136    v.SetMagThetaPhi(1, TMath::Pi()/2-dec, ra);
     137    v *= *this;
    135138
    136     TRotation rot;
    137     rot.RotateZ(-loc.Phi());
    138     rot.RotateY(-loc.Theta());
    139 
    140     TVector3 v(1, 0, 0);
    141     v *= *this;
    142     v *= rot;
    143 
    144     return TMath::ATan2(v(1), v(0));
     139    return MAstro::RotationAngle(ZZ(), XZ(), v.Theta(), v.Phi());
    145140}

  • trunk/MagicSoft/Mars/mastro/MObservatory.cc

    r3498 r3568  
    101101{
    102102    *fLog << all;
    103     *fLog << fObservatoryName << endl;
    104     *fLog << "Latitude " << (fLatitude > 0 ? (fLatitude*kRad2Deg) : -(fLatitude*kRad2Deg)) << " deg " << (fLatitude > 0 ? "W" : "E") << endl;
    105     *fLog << "Longitude " << (fLongitude > 0 ? (fLongitude*kRad2Deg) : -(fLongitude*kRad2Deg)) <<" deg " << (fLongitude < 0 ? "N" : "S") << endl;
    106     *fLog << "Height " << fHeight << "m" << endl;
     103    *fLog << underline << fObservatoryName << ":" << endl;
     104    *fLog << "Latitude:  " << TMath::Abs(fLatitude*kRad2Deg)  << " deg " << (fLatitude > 0 ? "W" : "E") << endl;
     105    *fLog << "Longitude: " << TMath::Abs(fLongitude*kRad2Deg) << " deg " << (fLongitude < 0 ? "N" : "S") << endl;
     106    *fLog << "Height:    " << fHeight << "m" << endl;
    107107}
    108108
    109109// --------------------------------------------------------------------------
    110110//
    111 // RotationAngle
    112 //
    113 // calculates the angle for the rotation of the sky image in the camera;
    114 // this angle is a function of the local coordinates
     111// Get the corresponding rotation angle of the sky coordinate system
     112// seen with an Alt/Az telescope.
    115113//
    116 //  theta [rad]: polar angle/zenith distance
    117 //  phi   [rad]: rotation angle/azimuth
    118 //
    119 // Return sin/cos component of angle
    120 //
    121 // calculate rotation angle alpha of sky image in camera
    122 // (see TDAS 00-11, eqs. (18) and (20))
     114// For more information see MAstro::RotationAngle
    123115//
    124116void MObservatory::RotationAngle(Double_t theta, Double_t phi, Double_t &sin, Double_t &cos) const
    125117{
    126     const Double_t sint = TMath::Sin(theta);
    127     const Double_t cost = TMath::Cos(theta);
    128 
    129     const Double_t sinl = fSinLatitude*sint;
    130     const Double_t cosl = fCosLatitude*cost;
    131 
    132     const Double_t sinp = TMath::Sin(phi);
    133     const Double_t cosp = TMath::Cos(phi);
    134 
    135     const Double_t v1 = sint*sinp;
    136     const Double_t v2 = cosl - sinl*cosp;
    137 
    138     const Double_t denom = TMath::Sqrt(v1*v1 + v2*v2);
    139 
    140     sin =  fCosLatitude*sinp / denom;
    141     cos = (sinl - cosl*cosp) / denom;
     118    MAstro::RotationAngle(fSinLatitude, fCosLatitude, theta, phi, sin, cos);
    142119}
    143120
    144121// --------------------------------------------------------------------------
    145122//
    146 // RotationAngle
    147 //
    148 // calculates the angle for the rotation of the sky image in the camera;
    149 // this angle is a function of the local coordinates
     123// Get the corresponding rotation angle of the sky coordinate system
     124// seen with an Alt/Az telescope.
    150125//
    151 //  theta [rad]: polar angle/zenith distance
    152 //  phi   [rad]: rotation angle/azimuth
    153 //
    154 // Return RotationAngle in rad
    155 //
    156 // calculate rotation angle alpha of sky image in camera
    157 // (see TDAS 00-11, eqs. (18) and (20))
     126// For more information see MAstro::RotationAngle
    158127//
    159128Double_t MObservatory::RotationAngle(Double_t theta, Double_t phi) const
    160129{
    161     const Double_t sint = TMath::Sin(theta);
    162     const Double_t cost = TMath::Cos(theta);
    163 
    164     const Double_t sinp = TMath::Sin(phi);
    165     const Double_t cosp = TMath::Cos(phi);
    166 
    167     const Double_t v1 = sint*sinp;
    168     const Double_t v2 = fCosLatitude*cost - fSinLatitude*sint*cosp;
    169 
    170     const Double_t denom = TMath::Sqrt(v1*v1 + v2*v2);
    171 
    172     return TMath::ASin((fCosLatitude*sinp) / denom);
     130    return MAstro::RotationAngle(fSinLatitude, fCosLatitude, theta, phi);
    173131}
    174 
    175 // --------------------------------------------------------------------------
    176 //
    177 // RotationAngle
    178 //
    179 // calculates the angle for the rotation of the sky image in the camera;
    180 // this angle is a function of the sky coordinates, the observatory
    181 // location and the time
    182 //
    183 //  ra  [rad]: Right ascension
    184 //  dec [rad]: Declination
    185 //
    186 // Return RotationAngle in rad
    187 //
    188 Double_t MObservatory::RotationAngle(Double_t ra, Double_t dec, const MTime &t) const
    189 {
    190     const Double_t alpha = t.GetGmst() + GetElong();
    191 
    192     TVector3 v;
    193     v.SetMagThetaPhi(1, TMath::Pi()/2-dec, alpha-ra);
    194     v.RotateY(GetPhi()-TMath::Pi()/2);
    195 
    196     return RotationAngle(v.Theta(), v.Phi());
    197 }

  • trunk/MagicSoft/Mars/mastro/MObservatory.h

    r3497 r3568  
    5858    void RotationAngle(Double_t theta, Double_t phi, Double_t &sin, Double_t &cos) const;
    5959    Double_t RotationAngle(Double_t theta, Double_t phi) const;
    60     Double_t RotationAngle(Double_t ra, Double_t dec, const MTime &t) const;
    6160
    6261    LocationName_t GetObservatoryKey() const { return fObservatoryKey; }
Note: See TracChangeset for help on using the changeset viewer.