#ifndef COORD_H
#define COORD_H

#include <math.h>          // floor
#include <iostream>

#include "MAGIC.h"

/* pi/180:  degrees to radians */
const double kDeg2Rad = 0.017453292519943295769236907684886127134428718885417;

/* 180/pi:  radians to degrees */
//const double kRad2Deg = 57.295779513082320876798154814105170332405472466564;

/* pi/2:  90 degrees in radians */
const double kPiDiv2 = 1.5707963267948966192313216916397514420985846996876;

/* 2pi */
const double k2Pi = 6.2831853071795864769252867665590057683943387987502;

class Deg
{
protected:
    double fDeg;

public:
    Deg(const double d) : fDeg(d) {}
    Deg(const Deg &c) { fDeg = c.fDeg; }

    void Set(double d) { fDeg=d; }

    double operator()() { return fDeg; }

    operator double() const { return fDeg*kDeg2Rad; }
};

class XY
{
    friend istream& operator>>(istream &in, XY &xy);
    friend ostream& operator<<(ostream &in, XY &xy);

protected:
    double fX;
    double fY;

public:
    XY(double x=0, double y=0) : fX(x), fY(y) {}
    XY(const XY &c) { fX = c.fX; fY = c.fY; }

    void Set(double x, double y) { fX=x; fY=y; }

    double X() const { return fX; }
    double Y() const { return fY; }

    void X(double x) { fX=x; }
    void Y(double y) { fY=y; }

    void operator/=(double c) { fX/=c; fY/=c; }
    void operator*=(double c) { fX*=c; fY*=c; }

    XY operator/(double c) const { return XY(fX/c, fY/c); }
    XY operator*(double c) const { return XY(fX*c, fY*c); }
    XY operator*(const XY &c) const { return XY(fX*c.fX, fY*c.fY); }
    XY operator/(const XY &c) const { return XY(fX/c.fX, fY/c.fY); }
    XY operator+(const XY &c) const { return XY(fX+c.fX, fY+c.fY); }
    XY operator-(const XY &c) const { return XY(fX-c.fX, fY-c.fY); }
    XY operator-() const { return XY(-fX, -fY); }

    double Sqr()   const { return fX*fX + fY*fY; }
    double Sqrt()  const { return sqrt(Sqr()); }
    double Ratio() const { return fX/fY; }
    void Round()         { fX=(int)(floor(fX+.5)); fY=(int)(floor(fY+.5)); }
};

inline std::istream& operator>>(std::istream &in,  XY &xy) { in  >> xy.fX; in  >> xy.fY; return in; }
inline std::ostream& operator<<(std::ostream &out, const XY &xy) { out << xy.X() << " " << xy.Y(); return out; }

class AltAz : public XY
{
public:
    AltAz(double alt=0, double az=0) : XY(alt, az) {}

    double Alt() const { return fX; }
    double Az()  const { return fY; }

    void operator*=(double c) { fX*=c; fY*=c; }
    void operator/=(double c) { fX*=c; fY*=c; }

    void Alt(double d) { fX=d; }
    void Az(double d)  { fY=d; }
    void operator*=(const XY &c)    { fX*=c.X(); fY*=c.Y(); }
    void operator/=(const XY &c)    { fX/=c.X(); fY/=c.Y(); }
    void operator-=(const AltAz &c) { fX-=c.fX; fY-=c.fY; }
    void operator+=(const AltAz &c) { fX+=c.fX; fY+=c.fY; }

    AltAz operator/(double c) const { return AltAz(fX/c, fY/c); }
    AltAz operator*(double c) const { return AltAz(fX*c, fY*c); }
    AltAz operator*(const XY &c) const { return AltAz(fX*c.X(), fY*c.Y()); }
    AltAz operator/(const XY &c) const { return AltAz(fX/c.X(), fY/c.Y()); }
    AltAz operator+(const AltAz &c) const { return AltAz(fX+c.fX, fY+c.fY); }
    AltAz operator-(const AltAz &c) const { return AltAz(fX-c.fX, fY-c.fY); }
    AltAz operator-() const { return AltAz(-fX, -fY); }
};

class ZdAz : public XY
{
public:
    ZdAz(double zd=0, double az=0) : XY(zd, az) {}
    ZdAz(const ZdAz &c) : XY(c) {}

    void operator*=(double c) { fX*=c; fY*=c; }
    void operator/=(double c) { fX*=c; fY*=c; }

    double Zd() const { return fX; }
    double Az() const { return fY; }

    void Zd(double d) { fX=d; }
    void Az(double d) { fY=d; }
    void operator*=(const XY &c)   { fX*=c.X(); fY*=c.Y(); }
    void operator/=(const XY &c)   { fX/=c.X(); fY/=c.Y(); }
    void operator-=(const ZdAz &c) { fX-=c.fX; fY-=c.fY; }
    void operator+=(const ZdAz &c) { fX+=c.fX; fY+=c.fY; }

    ZdAz operator/(double c) const { return ZdAz(fX/c, fY/c); }
    ZdAz operator*(double c) const { return ZdAz(fX*c, fY*c); }
    ZdAz operator*(const XY &c) const { return ZdAz(fX*c.X(), fY*c.Y()); }
    ZdAz operator/(const XY &c) const { return ZdAz(fX/c.X(), fY/c.Y()); }
    ZdAz operator+(const ZdAz &c) const { return ZdAz(fX+c.fX, fY+c.fY); }
    ZdAz operator-(const ZdAz &c) const { return ZdAz(fX-c.fX, fY-c.fY); }
    ZdAz operator-() const { return ZdAz(-fX, -fY); }
};

class RaDec : public XY
{
public:
    RaDec(double ra=0, double dec=0) : XY(ra, dec) {}

    double Ra()  const { return fX; }
    double Dec() const { return fY; }

    void operator*=(double c) { fX*=c; fY*=c; }
    void operator/=(double c) { fX*=c; fY*=c; }

    void Ra(double x)  { fX = x; }
    void Dec(double y) { fY = y; }

    RaDec operator/(double c) const { return RaDec(fX/c, fY/c); }
    RaDec operator*(double c) const { return RaDec(fX*c, fY*c); }
    RaDec operator*(const XY &c) const { return RaDec(fX*c.X(), fY*c.Y()); }
    RaDec operator+(const RaDec &c) const { return RaDec(fX+c.fX, fY+c.fY); }
    RaDec operator-(const RaDec &c) const { return RaDec(fX-c.fX, fY-c.fY); }
    RaDec operator-() const { return RaDec(-fX, -fY); }
};

inline double Rad2Deg(double rad)
{
    return kRad2Deg*rad;
}

inline double Deg2Rad(double rad)
{
    return kDeg2Rad*rad;
}
#endif
