Ignore:
Timestamp:
05/15/04 16:46:27 (21 years ago)
Author:
tbretz
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Cosy/catalog
Files:
3 edited

Legend:

Unmodified
Added
Removed

  • trunk/MagicSoft/Cosy/catalog/SlaStars.cc

    r3897 r4076  
    159159    double az;
    160160    slaAopqk (r, d, (double*)fAoprms,
    161               &az,    // observed azimuth (radians: N=0,E=90)
    162               &zd,    // observed zenith distance (radians) [-pi/2, pi/2]
     161              &az,    // observed azimuth (radians: N=0,E=90) [-pi, pi]
     162              &zd,    // observed zenith distance (radians)   [-pi/2, pi/2]
    163163              &h0,    // observed hour angle (radians)
    164164              &d1,    // observed declination (radians)
    … …  
    175175ZdAz SlaStars::GetApproxVel(const RaDec &radec) const // [rad/rad]
    176176{
     177    // radec        [rad]
     178    // GetApproxVel [rad/rad]
    177179    double az, vaz, aaz;
    178180    double el, vel, ael;

  • trunk/MagicSoft/Cosy/catalog/StarCatalog.cc

    r2407 r4076  
    55
    66#include <TSystem.h>
     7#include <TRotation.h>
    78
    89#include "slalib.h"
    … …  
    1112
    1213#include "MStarList.h"
     14#include "MAstroCatalog.h"
    1315
    1416ClassImp(StarCatalog);
    1517
    16 StarCatalog::StarCatalog(MObservatory::LocationName_t key) : SlaStars(key), fSao(NULL), fSrt(NULL), fEntries(0), fSinAngle(0), fCosAngle(1)
    17 {
    18     // p = pointer to MainFrame (not owner)
    19 
    20     //
    21     // read index file
    22     //
    23     File idx("sao/sao-sort.idx", "r");
    24     if (!idx)
    25         return;
    26 
    27     while (!idx.Eof())
    28     {
    29         idx.Newline();
    30         fEntries++;
    31     }
    32 
    33     idx.Reset();
    34 
    35     fSrt = new sort_t[fEntries];
    36 
    37     for (int i=0; i<fEntries; i++)
    38     {
    39         fSrt[i].ra  = idx.Geti(4);
    40         fSrt[i].dec = idx.Geti(4);
    41         fSrt[i].nr  = idx.Geti(10);
    42         idx.Newline();
    43     }
    44 
    45     //
    46     // open catalog
    47     //
    48     fSao = new SaoFile("sao/sao-sort.cmp");
     18StarCatalog::StarCatalog(MObservatory::LocationName_t key) : SlaStars(key), fAstro(0), /*fSao(NULL), fSrt(NULL), fEntries(0),*/ fSinAngle(0), fCosAngle(1)
     19{
     20    fAstro = new MAstroCatalog;
     21    fAstro->SetObservatory(*this);
     22    fAstro->SetPlainScreen();
    4923}
    5024
    5125StarCatalog::~StarCatalog()
    5226{
    53     if (fSrt)
    54         delete fSrt;
    55     if (fSao)
    56         delete fSao;
     27    delete fAstro;
    5728}
    5829
    5930void StarCatalog::SetPixSize(const double pixsize)
    6031{
    61     fPixSize = D2PI/360.0 * pixsize;
    62 
    63     fWidth  = fPixSize * 768/2;
    64     fHeight = fPixSize * 576/2;
     32    // pixsize [arcsec/pixel]
     33    fPixSize = D2PI/360.0*pixsize/3600; // [rad / (deg*pixel)]
     34    fAstro->SetRadiusFOV(pixsize, 768, 576);
    6535}
    6636
    … …  
    7040}
    7141
     42void StarCatalog::SetLimitMag(const float mag)
     43{
     44    fLimitMag = mag;
     45    fAstro->SetLimMag(mag);
     46}
     47
     48void StarCatalog::SetMjd(double mjd)
     49{
     50    SlaStars::SetMjd(mjd);
     51    fAstro->SetTime(MTime(mjd));
     52}
    7253
    7354void StarCatalog::SetAltAz(const AltAz &altaz)
    7455{
    7556    fAltAz = altaz * D2PI/360.0;
    76 
    77     cout << "Set --> Alt: " << 360.0/D2PI*fAltAz.Alt();
    78     cout << "  Az: " << fAltAz.Az() << endl;
    79 
    8057    fRaDec = CalcRaDec(fAltAz);
    8158
    82     cout << "Ra: " << 360.0/D2PI*fRaDec.Ra();
    83     cout << "  Dec: " << 360.0/D2PI*fRaDec.Dec() << endl;
    84 
    85     CalcAltAzRange();
    86     CalcRaDecRange();
     59    fAstro->SetRaDec(fRaDec.Ra(), fRaDec.Dec());
     60}
     61
     62void StarCatalog::Reload()
     63{
     64    fAstro->SetLimMag(99);
     65    //fAstro->ReadBSC("bsc5.dat");
     66    //fAstro->ReadHeasarcPPM("heasarc_ppm.tdat");
     67    fAstro->ReadCompressed("ppm9.bin");
     68    fAstro->SetLimMag(fLimitMag);
    8769}
    8870
    8971void StarCatalog::SetRaDec(const RaDec &radec)
    9072{
    91     fRaDec = radec;
    92     fRaDec *= D2PI/360.0;
    93 
     73    const RaDec rd = fRaDec*360.0/D2PI;;
     74
     75    const Bool_t same =
     76        rd.Ra() >radec.Ra() -1e-5 && rd.Ra() <radec.Ra() +1e-5 &&
     77        rd.Dec()>radec.Dec()-1e-5 && rd.Dec()<radec.Dec()+1e-5;
     78
     79    fRaDec = radec * D2PI/360.0;
    9480    fAltAz = CalcAltAz(fRaDec);
    9581
    96     //cout << "Alt: " << 360.0/D2PI*fAltAz.Alt() << "  ";
    97     //cout << "Az: "  << 360.0/D2PI*fAltAz.Az()  << endl;
    98 
    99     CalcRaDecRange();
    100     CalcAltAzRange();
    101 }
    102 
    103 void StarCatalog::CalcAltAzRange()
    104 {
    105     byte fAlt0[180];
    106 
    107     for (int h=0; h<180; h++)
    108         fAlt0[h] = kFALSE;
    109 
    110     for (int h=0; h<360; h++)
    111         fAz0[h] = kFALSE;
    112 
    113     double az0, alt0;
    114     double az1, alt1;
    115     //
    116     // scan horizontal border
    117     //
    118     for (int x=-768/2; x<768/2+1; x++)
    119     {
    120         slaDh2e(DPI+x*fPixSize, -fHeight, DPI/2-fAltAz.Alt(), &az0, &alt0);
    121         slaDh2e(DPI+x*fPixSize, +fHeight, DPI/2-fAltAz.Alt(), &az1, &alt1);
    122 
    123         const int z0 = ((int)(360.0/D2PI*(az0+fAltAz.Az()))+360)%360;
    124         const int t0 = (int)(360.0/D2PI*alt0);
    125 
    126         fAz0[z0] = kTRUE;
    127 
    128         if (-89<=t0 && t0<=90)
    129             fAlt0[90-t0] = kTRUE;
    130 
    131         const int z1 = ((int)(360.0/D2PI*(az1+fAltAz.Az()))+360)%360;
    132         const int t1 = (int)(360.0/D2PI*alt1);
    133 
    134         fAz0[z1] = kTRUE;
    135 
    136         if (-89<=t1 && t1<=90)
    137             fAlt0[90-t1] = kTRUE;
    138     }
    139 
    140     //
    141     // scan vertical border
    142     //
    143     for (int y=-576/2; y<576/2+1; y++)
    144     {
    145         slaDh2e(DPI-fWidth, y*fPixSize, DPI/2-fAltAz.Alt(), &az0, &alt0);
    146         slaDh2e(DPI+fWidth, y*fPixSize, DPI/2-fAltAz.Alt(), &az1, &alt1);
    147 
    148         const int z0 = ((int)(360.0/D2PI*(az0+fAltAz.Az()))+360)%360;
    149         const int t0 = (int)(360.0/D2PI*alt0);
    150 
    151         fAz0[z0] = kTRUE;
    152 
    153         if (-89<=t0 && t0<=90)
    154             fAlt0[90-t0] = kTRUE;
    155 
    156         const int z1 = ((int)(360.0/D2PI*(az1+fAltAz.Az()))+360)%360;
    157         const int t1 = (int)(360.0/D2PI*alt1);
    158 
    159         fAz0[z1] = kTRUE;
    160 
    161         if (-89<=t1 && t1<=90)
    162             fAlt0[90-t1] = kTRUE;
    163     }
    164 
    165     //
    166     // count degrees of azimut
    167     //
    168     fAzCnt=0;
    169     for (int x=0; x<360; x++)
    170         if (fAz0[x])
    171             fAzCnt++;
    172 
    173     //cout << "fAzCnt: " << setw(3) << fAzCnt << "  " << flush;
    174 
    175     //
    176     // calculate min and max of altitude
    177     //
    178     fAltMin=0;
    179     fAltMax=0;
    180     for (int y=0; y<180; y++)
    181     {
    182         if (fAlt0[y])
    183             fAltMax = y;
    184 
    185         if (fAlt0[179-y])
    186             fAltMin = 179-y;
    187     }
    188 
    189     fAltMin -= 90;
    190     fAltMax -= 90;
    191 
    192     //
    193     // check whether altaz north- or south-pole is in the visible region
    194     //
    195     byte img[768*576];
    196     if (DrawAltAz(0, img, 90, 0))
    197     {
    198         fAltMax=89;
    199         cout << "Alt Az Pole1 Inside!" << endl;
    200     }
    201     if (DrawAltAz(0, img, -90, 0))
    202     {
    203         fAltMin=-90;
    204         cout << "Alt Az Pole2 Inside!" << endl;
    205     }
    206 
    207     //cout << "fAltMin: " << setw(3) << fAltMin << "  ";
    208     //cout << "fAltMax: " << setw(3) << fAltMax << endl;
    209 }
    210 
    211 void StarCatalog::CalcRaDecRange()
    212 {
    213     //
    214     // calculate range to search in
    215     //
    216     byte fDec[180];
    217 
    218     for (int h=0; h<180; h++)
    219         fDec[h] = kFALSE;
    220 
    221     for (int h=0; h<360; h++)
    222         fRa0[h] = kFALSE;
    223 
    224     double ha0, ha1;
    225     double de0, de1;
    226 
    227     const double phi   = GetPhi();
    228     const double alpha = GetAlpha();
    229     //
    230     // scan horizontal border
    231     //
    232     for (int x=-768/2; x<768/2+1; x++)
    233     {
    234         double dx, dy;
    235         slaDh2e(DPI-x*fPixSize, -fHeight, DPI/2-fAltAz.Alt(), &dx, &dy);
    236         slaDh2e(fAltAz.Az()+dx, -dy, phi, &ha0, &de0);
    237 
    238         slaDh2e(DPI-x*fPixSize, +fHeight, DPI/2-fAltAz.Alt(), &dx, &dy);
    239         slaDh2e(fAltAz.Az()+dx, -dy, phi, &ha1, &de1);
    240 
    241         const int h0 = ((int)(360.0/D2PI*(alpha-ha0))+360)%360;
    242         const int d0 = (int)(360.0/D2PI*de0);
    243 
    244         fRa0[h0] = kTRUE;
    245 
    246         if (-90<=d0 && d0<=89)
    247             fDec[d0+90] = kTRUE;
    248 
    249         const int h1 = ((int)(360.0/D2PI*(alpha-ha1))+360)%360;
    250         const int d1 = (int)(360.0/D2PI*de1);
    251 
    252         fRa0[h1] = kTRUE;
    253 
    254         if (-90<=d1 && d1<=89)
    255             fDec[d1+90] = kTRUE;
    256     }
    257 
    258     //
    259     // scan vertical border
    260     //
    261     for (int y=-576/2; y<576/2+1; y++)
    262     {
    263         double dx, dy;
    264         slaDh2e(DPI-fWidth, -y*fPixSize, DPI/2-fAltAz.Alt(), &dx, &dy);
    265         slaDh2e(fAltAz.Az()+dx, -dy, phi, &ha0, &de0);
    266 
    267         slaDh2e(DPI+fWidth, -y*fPixSize, DPI/2-fAltAz.Alt(), &dx, &dy);
    268         slaDh2e(fAltAz.Az()+dx, -dy, phi, &ha1, &de1);
    269 
    270         const int h0 = ((int)(360.0/D2PI*(alpha-ha0))+360)%360;
    271         const int d0 = (int)(360.0/D2PI*de0);
    272 
    273         fRa0[h0] = kTRUE;
    274 
    275         if (-90<=d0 && d0<=89)
    276             fDec[d0+90] = kTRUE;
    277 
    278         const int h1 = ((int)(360.0/D2PI*(alpha-ha1))+360)%360;
    279         const int d1 = (int)(360.0/D2PI*de1);
    280 
    281         fRa0[h1] = kTRUE;
    282 
    283         if (-90<=d1 && d1<=89)
    284             fDec[d1+90] = kTRUE;
    285     }
    286 
    287     //
    288     // count degrees of right ascension
    289     //
    290     fRaCnt=0;
    291     for (int x=0; x<360; x++)
    292         if (fRa0[x])
    293             fRaCnt++;
    294     //cout << "fRaCnt: " << setw(3) << fRaCnt << "  " << flush;
    295 
    296     //
    297     // calculate min and max of declination
    298     //
    299     for (int y=0; y<180; y++)
    300     {
    301         if (fDec[y])
    302             fDecMax = y;
    303 
    304         if (fDec[179-y])
    305             fDecMin = 179-y;
    306     }
    307 
    308     fDecMin -= 90;
    309     fDecMax -= 90;
    310 
    311     //
    312     // check whether radec north- or south-pole is in the visible region
    313     //
    314     byte img[768*576];
    315     if (DrawRaDec(0, img, 0, 90))
    316     {
    317         fDecMax=89;
    318         cout << "Ra Dec Pole1 Inside!" << endl;
    319     }
    320     if (DrawRaDec(0, img, 0, -90))
    321     {
    322         fDecMin=-90;
    323         cout << "Ra Dec Pole1 Inside!" << endl;
    324     }
    325 
    326     //cout << "fDecMin: " << setw(3) << fDecMin << "  ";
    327     //cout << "fDecMax: " << setw(3) << fDecMax << endl;
    328 }
    329 
    330 void StarCatalog::DrawSCAltAz(byte *img, const int color) const
    331 {
    332     //
    333     // ------------ draw az lines ---------------
    334     //
    335     for (int az=0; az<360; az++)
    336     {
    337         if (!fAz0[az])
    338             continue;
    339 
    340         for (double alt=fAltMin-1; alt<fAltMax+1; alt+=0.006*(fAltMax-fAltMin))
    341         {
    342             if ((alt>88 && az%5) || alt>89.5)
    343                 continue;
    344 
    345             DrawAltAz(color, img, alt, az);
    346         }
    347     }
    348 
    349     //
    350     // ------------ draw alt lines ---------------
    351     //
    352     for (int alt=fAltMin; alt<fAltMax+1; alt++)
    353     {
    354         for (double az=0; az<360; az+=0.004*fAzCnt)
    355         {
    356             if (!fAz0[(int)az] && !fAz0[(int)(az+359)%360] && !fAz0[(int)(az+1)%360])
    357                 continue;
    358 
    359             DrawAltAz(color, img, alt, az);
    360         }
    361     }
    362 }
    363 
    364 void StarCatalog::DrawSCRaDec(byte *img, const int color) const
    365 {
    366     //
    367     // ------------ draw ra lines ---------------
    368     //
    369     for (int ra=0; ra<360; ra++)
    370     {
    371         if (!fRa0[ra])
    372             continue;
    373 
    374         for (double dec=fDecMin-1; dec<fDecMax+1; dec+=0.005*(fDecMax-fDecMin))
    375         {
    376             if ((dec>88 && ra%5) || dec>89.5)
    377                 continue;
    378 
    379             DrawRaDec(color, img, ra, dec, ra==0||ra==90);
    380         }
    381     }
    382 
    383     //
    384     // ------------ draw dec lines ---------------
    385     //
    386     for (int dec=fDecMin; dec<fDecMax+1; dec++)
    387     {
    388         for (double ra=0; ra<360; ra+=0.003*fRaCnt)
    389         {
    390             if (!fRa0[(int)ra])
    391                 continue;
    392 
    393             DrawRaDec(color, img, ra, dec, dec==89);
    394         }
    395     }
     82    fAstro->SetRaDec(fRaDec.Ra(), fRaDec.Dec());
     83    if (!same)
     84        Reload();
    39685}
    39786
    … …  
    40998    memset(cimg, 0, 768*576);
    41099
    411     DrawSCAltAz(cimg, 2<<4);
    412     DrawSCRaDec(cimg, 2);
    413 
    414100    DrawStars(list, cimg);
    415101    DrawCross(img, 768/2, 576/2);
    416 }
    417 
    418 void StarCatalog::GetImg(byte *img, byte *cimg, const double utc,
    419                          const RaDec &radec)
    420 {
    421     MStarList list;
    422     GetStars(list, utc, radec);
    423     GetImg(img, cimg, list);
    424     /*
    425      // memset(img,  0, 768*576);
    426      SetMjd(utc);
    427      //fAlpha = sla.GetAlpha();
    428      SetRaDec(radec);
    429      //CalcImg(cimg);
    430      */
    431 }
    432 
    433 void StarCatalog::GetImg(byte *img, byte *cimg, const double utc,
    434                          const AltAz &altaz)
    435 {
    436     MStarList list;
    437     GetStars(list, utc, altaz);
    438     GetImg(img, cimg, list);
    439     /*
    440      // memset(img,  0, 768*576);
    441 
    442      SetMjd(utc);
    443      //fAlpha = sla.GetAlpha();
    444      SetAltAz(altaz);
    445 
    446      CalcRaDecRange();
    447 
    448      //CalcImg(img);
    449      */
    450 
    451 }
    452 
    453 void StarCatalog::GetStars(MStarList &list, const double utc, const RaDec &radec)
    454 {
    455     SetMjd(utc);
    456     SetRaDec(radec);
    457 
    458     CalcStars(list);
    459 }
    460 
    461 void StarCatalog::GetStars(MStarList &list, const double utc, const AltAz &altaz)
    462 {
    463     SetMjd(utc);
    464     SetAltAz(altaz);
    465 
    466     CalcRaDecRange();
    467     CalcStars(list);
    468102}
    469103
    … …  
    489123}
    490124
    491 Bool_t StarCatalog::DrawAltAz(const int color, byte *img, double alt, double az, int size) const
    492 {
    493     //
    494     // alt/az[deg] -> alt/az[rad]
    495     //
    496     alt *= D2PI/360.0;
    497     az  *= D2PI/360.0;
    498 
    499     //
    500     // alt/az[rad] -> alt/az[pix]
    501     //
    502     double dx, dy;
    503     slaDe2h(az-fAltAz.Az(), -alt, DPI/2-fAltAz.Alt(), &dx, &dy);
    504 
    505     //
    506     // Align alt/az[pix]
    507     //
    508     const int xx = (int)(((dx-DPI)*fCosAngle - dy*fSinAngle + fWidth)/fPixSize);
    509     const int yy = (int)(((dx-DPI)*fSinAngle + dy*fCosAngle + fHeight)/fPixSize);
    510     //const int xx = 767-(int)((fWidth-dx+DPI)/fPixSize);
    511     //const int yy =     (int)((fHeight+dy)/fPixSize);
    512 
    513     //
    514     // Range Check
    515     //
    516     if (!(0<=xx && xx<768 && 0<=yy && yy<576))
    517         return kFALSE;
    518 
    519     //
    520     // Draw
    521     //
    522     DrawCircle(color, img, xx, yy, size);
    523 
    524     return kTRUE;
    525 }
    526 
    527 Bool_t StarCatalog::Draw(const int color, byte *img, const AltAz &altaz)
    528 {
    529     return DrawAltAz(color, img, altaz.Alt(), altaz.Az());
    530 }
    531 
    532 /*
    533 Bool_t StarCatalog::Draw(const int color, byte *img, const SaoFile *sao)
    534 {
    535     if (sao->MagV() > fLimitMag)
    536         return kFALSE;
    537 
    538     //
    539     // ---- mean to observed ---
    540     //
    541     AltAz altaz=CalcAltAz(sao->GetRaDec()) * 360.0/D2PI;
    542 
    543     const int mag = (10 - (sao->MagV()>1 ? (int)sao->MagV() : 1))/2;
    544 
    545     //
    546     // ---- imaging -----
    547     //
    548     return DrawAltAz(color, img, altaz.Alt(), altaz.Az(), mag);
    549 }
    550 */
    551 
    552 Bool_t StarCatalog::DrawRaDec(const int color, byte *img, double ra, double dec, int size) const
    553 {
    554     //
    555     // radec[deg] -> radec[rad]
    556     //
    557     ra  *= D2PI/360.0;
    558     dec *= D2PI/360.0;
    559 
    560     //
    561     // radec[rad] -> hadec[rad]
    562     //
    563     const double ha = GetAlpha()-ra;
    564 
    565     //
    566     // hadec[rad] -> altaz[rad]
    567     //
    568     double alt, az;
    569     slaDe2h(ha, dec, GetPhi(), &az, &alt);
    570 
    571     //
    572     // altaz[rad] -> altaz[deg]
    573     //
    574     alt *= 360.0/D2PI;
    575     az  *= 360.0/D2PI;
    576 
    577     return DrawAltAz(color, img, alt, az, size);
    578 }
    579 
    580 Bool_t StarCatalog::Draw(const int color, byte *img, const RaDec &radec)
    581 {
    582     return DrawRaDec(color, img, radec.Ra(), radec.Dec());
    583 }
    584 /*
    585 void StarCatalog::CalcImg(byte *img)
    586 {
    587 
    588     //
    589     // --------- search for stars in catalog ----------
    590     //
    591     int count   = 0;
    592     int deleted = 0;
    593 
    594     int idx     = 0;
    595 
    596     while (fSrt[idx].dec<fDecMin)
    597         idx++;
    598 
    599     idx--;
    600     while (++idx<fEntries && fSrt[idx].dec<fDecMax+1)
    601     {
    602         const int ra = fSrt[idx].ra;
    603 
    604         if (!fRa0[ra])
    605             continue;
    606 
    607         int nr = fSrt[idx].nr;
    608         do
    609         {
    610             //
    611             // Get entry from catalog
    612             //
    613             fSao->GetEntry(nr++);
    614 
    615             //
    616             // Try to draw star into the image
    617             //  white = 0xff
    618             //
    619             if (!Draw(0x0f, img, fSao))
    620                 deleted++;
    621 
    622             count++;
    623         }
    624         while ((int)(360.0/D2PI*fSao->Ra())==ra);
    625     }
    626 
    627     cout << " " << count << "-" << deleted << "=" << count-deleted << " " << flush;
    628 }
    629 */
     125void StarCatalog::PaintImg(unsigned char *buf, int w, int h)
     126{
     127    fAstro->PaintImg(buf, w, h);
     128}
     129
    630130void StarCatalog::DrawStars(MStarList &list, byte *img)
    631131{
    … …  
    637137        const int mag = (10 - (star->GetMag()>1 ? (int)star->GetMag() : 1))/2;
    638138
    639         Double_t color = 0x0f;
    640 
     139        Double_t color = 0xf0; //0x0f;
    641140        DrawCircle(color, img, (int)star->GetX(), (int)star->GetY(), mag);
    642141    }
    … …  
    645144void StarCatalog::CalcStars(MStarList &list) const
    646145{
    647     //
    648     // --------- search for stars in catalog ----------
    649     //
    650     if (fEntries==0)
    651         return;
    652 
    653     int count   = 0;
    654     int deleted = 0;
    655 
    656     int idx     = 0;
    657 
    658     while (fSrt[idx].dec<fDecMin)
    659         idx++;
    660 
    661     idx--;
    662     while (++idx<fEntries && fSrt[idx].dec<fDecMax+1)
     146    // Align stars into telescope system
     147    // (Move the telescope to pointing position)
     148    TRotation align;
     149    align.RotateZ(-fAltAz.Az());
     150    align.RotateY(-(TMath::Pi()/2-fAltAz.Alt()));
     151    align.RotateZ(TMath::Pi()/2);
     152
     153    // For an apropriate unit conversion to pixels
     154    const Double_t scale = TMath::RadToDeg()*sqrt(768*768 + 576*576)/(fAstro->GetRadiusFOV()*2);
     155
     156    // Get List of stars from catalog
     157    TIter Next(fAstro->GetList());
     158    TVector3 *star=0;
     159
     160    const Double_t limmag = pow(10, -fLimitMag/2.5);
     161
     162    while ((star=(TVector3*)Next()))
    663163    {
    664         const int ra = fSrt[idx].ra;
    665 
    666         if (!fRa0[ra])
     164        // Check for limiting magnitude
     165        const Double_t mag = star->Mag();
     166        if (mag < limmag)
    667167            continue;
    668168
    669         int nr = fSrt[idx].nr;
    670         do
    671         {
    672             //
    673             // Get entry from catalog
    674             //
    675             fSao->GetEntry(nr++);
    676 
    677             if (fSao->MagV() > fLimitMag)
    678                 continue;
    679 
    680             //
    681             // ---- mean to observed ---
    682             //
    683             AltAz altaz=CalcAltAz(fSao->GetRaDec());
    684 
    685             //
    686             // alt/az[rad] -> alt/az[pix]
    687             //
    688             double dx, dy;
    689             slaDe2h(altaz.Az()-fAltAz.Az(), -altaz.Alt(),
    690                     DPI/2-fAltAz.Alt(), &dx, &dy);
    691 
    692             //
    693             // Align and rotate alt/az[pix]
    694             //
    695             float xx = ((dx-DPI)*fCosAngle - dy*fSinAngle + fWidth)/fPixSize;
    696             float yy = ((dx-DPI)*fSinAngle + dy*fCosAngle + fHeight)/fPixSize;
    697 
    698             //
    699             // Range Check, add stars to the list
    700             //
    701             if (!(0<=xx && xx<768 && 0<=yy && yy<576))
    702             {
    703                 deleted++;
    704                 continue;
    705             }
    706 
    707             list.Add(xx, yy, fSao->MagV());
    708             count++;
    709         }
    710         while ((int)(360.0/D2PI*fSao->Ra())==ra);
     169        // Get star position and do an apropiate
     170        // conversion to local coordinates
     171        const RaDec rd(star->Phi(), TMath::Pi()/2-star->Theta());
     172        const ZdAz  za(CalcZdAz(rd));
     173
     174        // Virtually move telescope to pointing position
     175        TVector3 loc;
     176        loc.SetMagThetaPhi(1, za.Zd(), za.Az());
     177        loc *= align;
     178
     179        // Sanity check
     180        if (loc(2)<0)
     181            continue;
     182
     183        // Stretch such, that the Z-component is alwas the same. Now
     184        // X and Y contains the intersection point between the star-light
     185        // and the plain of a virtual plain screen (ccd...)
     186        loc *= 1./loc(2);
     187
     188        // Do an apropriate unit conversion to pixels
     189        loc *= scale;
     190
     191        // if (loc.Mod2()>fRadiusFOV*fRadiusFOV)
     192        //     continue;
     193
     194        // Rotate by the rotation angle of the video camera
     195        Float_t xx = loc.X()*fCosAngle - loc.Y()*fSinAngle;
     196        Float_t yy = loc.X()*fSinAngle + loc.Y()*fCosAngle;
     197
     198        // Store pixel coordinates of star in list
     199        list.Add(xx+768/2, yy+576/2, -2.5*log10(mag));
    711200    }
    712 
    713     cout << "Showing " << count+deleted << "-" << deleted << "=" << count << " stars." << endl;
    714201}
    715202
    716203AltAz StarCatalog::CalcAltAzFromPix(Double_t pixx, Double_t pixy) const
    717204{
    718     pixx *= fPixSize;
    719     pixy *= fPixSize;
    720 
    721     const double dx =  (pixx-fWidth)*fCosAngle + (pixy-fHeight)*fSinAngle;
    722     const double dy = -(pixx-fWidth)*fSinAngle + (pixy-fHeight)*fCosAngle;
     205    double dx =  (pixx-576/2)*fCosAngle + (pixy-768/2)*fSinAngle;
     206    double dy = -(pixx-576/2)*fSinAngle + (pixy-768/2)*fCosAngle;
     207
     208    dx *= fPixSize;
     209    dy *= fPixSize;
    723210
    724211    //const double dx = (pixx-768.0)*fPixSize + fWidth+DPI;

  • trunk/MagicSoft/Cosy/catalog/StarCatalog.h

    r1953 r4076  
    44#ifndef ROOT_TROOT
    55#include <TROOT.h>
     6#endif
     7#ifndef ROOT_GuiTypes
     8#include <GuiTypes.h>
    69#endif
    710#ifndef SAOFILE_H
    … …  
    1720
    1821class MStarList;
     22class MAstroCatalog;
    1923
    2024class StarCatalog : public SlaStars
    2125{
    2226private:
    23     SaoFile *fSao;
    24     sort_t  *fSrt;
    25     int      fEntries;
     27    MAstroCatalog *fAstro;
    2628
    2729    double   fPixSize;  // [rad/pix] size of one pixel
    28     double   fWidth;    // size of display
    29     double   fHeight;   //
    3030    double   fSinAngle;
    3131    double   fCosAngle;
    … …  
    3434
    3535    AltAz    fAltAz;    // [rad]
    36     byte     fAz0[360];
    37     int      fAltMin;
    38     int      fAltMax;
    39     int      fAzCnt;
    40 
    4136    RaDec    fRaDec;    // [rad]
    42     byte     fRa0[360];
    43     int      fRaCnt;
    44     int      fDecMin;
    45     int      fDecMax;
    4637
    4738    static void DrawCross(byte *img, const int x, const int y);
    4839    static void DrawCircle(int color, byte *img, int xx, int yy, int size);
    4940
    50     Bool_t DrawAltAz(const int color, byte *img, double alt, double az,  int size=0) const;
    51     Bool_t DrawRaDec(const int color, byte *img, double ra,  double dec, int size=0) const;
    52 
    53     Bool_t Draw(const int color, byte *img, const AltAz &altaz);
    54     Bool_t Draw(const int color, byte *img, const RaDec &radec);
    55 
    56     //Bool_t Draw(const int color, byte *img, const SaoFile *sao);
    57     //void   CalcImg(byte *);
    58 
    59     void   CalcStars(MStarList &list) const;
    60 
    61     static void DrawStars(MStarList &list, byte *img);
    62 
    6341    void   SetRaDec(const RaDec &radec);
    6442    void   SetAltAz(const AltAz &altaz);
    65     void   DrawSCAltAz(byte *img, const int color) const;
    66     void   DrawSCRaDec(byte *img, const int color) const;
    67  
    68     void   CalcRaDecRange();
    69     void   CalcAltAzRange();
    70 
    71 //    RaDec  AltAz2RaDec(const AltAz &altaz) const;
    72 //    AltAz  RaDec2AltAz(const RaDec &radec, const RaDec &rdpm) const;
    7343
    7444public:
    … …  
    7747
    7848    void GetImg(byte *img, byte *cimg, MStarList &list) const;
    79     void GetImg(byte *img, byte *cimg, const double utc, const RaDec &radec);
    80     void GetImg(byte *img, byte *cimg, const double utc, const AltAz &altaz);
    81 
    82     void GetStars(MStarList &list, const double utc, const RaDec &radec);
    83     void GetStars(MStarList &list, const double utc, const AltAz &altaz);
     49    void PaintImg(unsigned char *buf, int w, int h);
    8450
    8551    const AltAz GetAltAz() const { return fAltAz*kRad2Deg; }
    … …  
    8854
    8955    void  SetPixSize(const double pixsize);
    90     void  SetLimitMag(const float mag) { fLimitMag = mag; }
     56    void  SetLimitMag(const float mag);
    9157    void  SetRotationAngle(const float angle) { fSinAngle = sin(angle/kRad2Deg); fCosAngle = cos(angle/kRad2Deg); }
     58    void  Reload();
    9259
    9360    double GetPixSize() const;
    … …  
    9562    AltAz CalcAltAzFromPix(Double_t pixx, Double_t pixy) const;
    9663
     64    virtual void SetMjd(double mjd);
     65
     66    void SetPointing(double mjd, const RaDec &radec)
     67    {
     68        SetMjd(mjd); SetRaDec(radec);
     69    }
     70
     71    void   CalcStars(MStarList &list) const;
     72    static void DrawStars(MStarList &list, byte *img);
     73
    9774    ClassDef(StarCatalog, 0)
    9875};
Note: See TracChangeset for help on using the changeset viewer.