source: trunk/MagicSoft/Cosy/catalog/SlaPlanets.cc@ 1757

#include "SlaPlanets.h"

#include <iostream.h>

#include "slalib.h"

ClassImp(SlaPlanets);

SlaPlanets::SlaPlanets() : Slalib(), fDt(slaDt(2000.0)/60./60./24.)
{
}

SlaPlanets::~SlaPlanets()
{
}

void SlaPlanets::SetMjd(const double mjd)
{
    Slalib::SetMjd(mjd);

    //
    // coordinates of earth: heliocentric, equatoril, J2000
    //
    double veb[3]; // barycentric velocity
    double ceb[3]; // barycentric coordinates
    double veh[3]; // heliocentric velocity

    slaEvp(mjd + fDt, 2000.0, veb, ceb, veh, fEarth);
}

#include "SlaStars.h"

void SlaPlanets::UpdatePlanetPos(ePlanets_t planet)
{
    /*
    double cp[6] = { 0, 0, 0, 0, 0, 0 };

    double coord[6];

    switch (planet)
    {
    case kEMoon:
        //
        // -- approximate geocentric pos and vel --
        //
        slaDmoon(GetMjd() + fDt, coord);
        break;

    case kEMercury:
    case kEVenus:
    case kEMars:
    case kEJupiter:
    case kESaturn:
    case kEUranus:
    case kENeptune:
    case kEPluto:
    case kESun:
        {
        */
    //
    // coordinates of planet: heliocentric, equatoril, J2000
    //
    double ra,dec,diam;
    // TDB, planet, elong, phi, *ra, *dec, *diam

    // One can use TT instead of TDB for all plenets
    // (except the moon)

    const double tt = GetMjd() + slaDtt(GetMjd())/60./60./24.;

    // FIXME: Shell be UT1
    const double ut1 = fmod(GetMjd(), 1);

    // tt should be tdb, but tt can be used
    const double tdb = tt + slaRcc(tt, ut1, GetElong(), 8300*cos(GetPhi()), 8300*sin(GetPhi() ))/60./60./24.;
    /*
     cout << "Dt: " << slaDtt(GetMjd()) << endl;
     cout << "Rc: " << slaRcc(tt, ut1, GetElong(), 8300*cos(GetPhi()), 8300*sin(GetPhi() )) << endl;
     */
    slaRdplan(tt, planet, GetElong(), GetPhi(),
              &ra, &dec, &diam);
    /*
     switch (planet)
     {
     case kEMoon: cout << "Moon: "; break;
     case kEMercury: cout << "Merkur: "; break;
     case kEVenus: cout << "Venus: "; break;
     case kEMars: cout << "Mars: "; break;
     case kEJupiter: cout << "Jupiter: "; break;
     case kESaturn: cout << "Saturn: "; break;
     case kEUranus: cout << "Uranus: "; break;
     case kENeptune: cout << "Neptun: "; break;
     case kEPluto: cout << "Pluto: "; break;
     case kESun: cout << "Sonne: "; break;
     }

     cout << (int)(ra*12/TMath::Pi()) << "h"<< (int)(ra*12/TMath::Pi()*60)%60 <<"min ";
     cout << (int)(dec*kRad2Deg) << "deg" << (int)(dec*kRad2Deg*60)%60 << "min" << endl;
     */
    /*
     SlaStars sla;
     sla.SetMjd(GetMjd());
     ZdAz za = sla.CalcZdAz(RaDec(ra,dec));
     cout << za.Zd()*kRad2Deg << " " << fmod(za.Az()*kRad2Deg+720,360.) << endl;
     fZdAz[planet] = za;
     */

    double az, el;
    slaDe2h(GetAlpha()-ra, dec, GetPhi(), &az, &el);

    cout << 90-el*kRad2Deg << " " << fmod(az*kRad2Deg+720,360.) << endl;
    fZdAz[planet].Set(TMath::Pi()/2-el, az);
    /*
     }
     // !! Fall through!
     break;
     case kESun:
     coord[0] = cp[0] - fEarth[0];
     coord[1] = cp[1] - fEarth[1];
     coord[2] = cp[2] - fEarth[2];
     break;
     }
     fZdAz[planet] = XYZ2ZdAz(coord);

     if (planet==kESun)
     {
     cout << fZdAz[kESun].Zd()*kRad2Deg << " " <<  fZdAz[kESun].Az()*kRad2Deg << endl;
     }*/
}